TROPHY: A Topologically Robust Physics-Informed Tracking Framework for Tropical Cyclones

TC tracking has been studied over the past decades for weather forecasting and climate analysis to issue early warnings, assess impacted areas, and provide risk assessments for public and critical infrastructures. Most of the traditional TC tracking algorithms require multiple dynamic and thermodynamic variables at different altitudes to detect and track a TC. These algorithms determine a detected feature as a TC candidate by tuning parameter thresholds; however, the choices of thresholds are generally subjective [14]. An example of a traditional TC tracking algorithms is TRACK [20], which uses relative vorticity at 850, 700, 600, 500, and 250 hPa as the key variable. TRACK uses certain criteria (still based on vorticity) during post-tracking to isolate the warm core of TCs. The criteria must be jointly attained for at least one day. Another example of a traditional TC tracker is TempestExtremes [55, 44], which is used in this study for comparison with TROPHY in Section 5. TempestExtremes uses sea level pressure (SLP) as its key feature-tracking variable. TC candidates are initially identified by minima and a closed contour of the SLP field. Next, a geopotential height difference between 250 to 500 hPa is used to refine the candidate definitions. Instead of the living time, TempestExtremes requires that a tracked storm travels at least $8^{\circ}$ between 10 N to 40 N latitudes. Both TRACK and TempestExtremes algorithms and many others (e.g., [3]) require many input variables that are often not readily available from raw model output and need additional calculations, which can involve handling a big amount of data. Also, these algorithms are not exactly tracking the eyes of TCs. They typically identify minimum SLP and maximum vorticity, where the tracked path of a TC is not exactly along the eye of the TC, and sometimes the path can even be biased toward the eyewall (a ring of tall thunderstorms that produce heavy rains and usually the strongest winds).

Recently, several contour-based and topology-based feature tracking approaches have been proposed [51, 12, 38, 31, 11, 35, 6] that have potential use in TC tracking. Most approaches identify features of interest with regions enclosed by streamlines or level sets (i.e., isosurfaces or contours). Wischgol et al. [51] discussed closed streamlines in 2D fields, which can represent the eyes of the TCs and indicate locations and sizes of the eyes in TC tracking. Correspondences between features in consecutive time steps can be identified via spatial overlap [38], critical point tracking [35], tree structure [6], or combinatorial feature flow fields [31]. Critical point tracking is most relevant to our framework, which is reviewed next.

Critical point tracking establishes the correspondences between critical points in successive time steps. It is a key tool from vector field topology and plays an important role in understanding the behavior of time-varying vector fields. Algorithms for critical point tracking may be classified as proximity-, integral-, and interpolation-based methods. Proximity-based methods (e.g., [18, 19]) find correspondences of critical points based on distance proximity in the domain. Integral-based approaches represent the tracking of critical points as streamlines of a higher dimensional field, called the feature flow field (FFF) [40, 49, 31], and compute feature tracks based on tangent curves in FFF. Interpolation-based methods take into account the time as an additional dimension in addition to the space domain [42, 43, 16, 17].

The robustness of critical points has been introduced to quantify the structural stability of critical points [34] and has been used in vector field simplification [36, 37], feature extraction [46], and visualization [48]. Skraba and Wang [35] showed the potential use of robustness in feature tracking, that is, finding correspondences between critical points based on their closeness in stability, measured by robustness, instead of just distance proximity within the domain. Building on the theoretical basis established by [35], Yan et al. [53] proposed a multilevel robustness (MRT) framework to realize critical point tracking in practice for large-scale scientific simulations, see Section 2.2 for details.

A critical point $x\in{\mathbb{X}}$ in $f$ is where the vector vanishes, that is, $|f(x)|=0$. A critical point $x$ can be classified w.r.t. its degree ${\mathrm{deg}}(x)$, defined as the number of field rotations while traveling along a closed curve counterclockwise surrounding $x$ (enclosing no other critical point). A source/sink/center has degree $+1$, whereas a saddle point has degree $-1$. Critical points are important features in studying flow behavior in many applications; see Fig. 1 as an example.

In most cases, the center of a TC can be detected as a center in a vector field when it is intensified into a strong hurricane, with very low wind speed in the eye and extremely high wind speed along the eyewall. During the dissipating phase of a TC, such as at its landfall, the center of the TC can be detected as either a source or a sink. If a center transforms to a source, then it indicates a divergence in meteorology, which means the weather can be clear and calm. If a center transforms to a sink, then it indicates a convergence, which is associated with clouds and precipitation. An example of this phenomenon is Hurricane Florence in 2018: a clear sink forms in the 2D wind field during its landfall, bringing 1-in-500-year expected flooding due to heavy precipitation.

To construct an augmented merge tree from a 2D vector field $f$, first we define a scalar field $f_{0}:{\mathbb{X}}\to{\mathbb{R}}$ by assigning the vector magnitude to each point $x\in{\mathbb{X}}$, that is, $f_{0}(x)=||f(x)||_{2}$. In this paper, $f_{0}$ can be expressed as wind speed. Second, instead of using local minima of $f_{0}$ as leaves of the merge tree, the leaves of our augmented merge tree consist of ${\mathbb{X}}_{0}$, which is precisely the set of critical points of $f$. The tracking of the merging behavior of components is the same as classic merge tree construction. Third, once the merge tree is constructed, it can be further augmented with the degrees of critical points (on leaves) and the degrees of components (on internal nodes). The degree of a component is defined as the sum of degrees of critical points the component contains. See Appendix A for an example.

The robustness of a critical point can be calculated as the function value of its lowest zero-degree ancestor in the augmented merge tree [48]. See Appendix B for an example.

In practice, vector fields generated from large-scale ocean and atmospheric datasets contain features at different scales. The drawback of classic robustness comes from building a single merge tree with critical points in the entire domain, which suffer from undesirable boundary effects [53]. To mitigate such drawbacks, Yan et al.[53] introduced a notion of multilevel robustness (reviewed in Section 2.2.1). This notion captures the multiscale nature of the data and mitigates the boundary effects suffered by classic robustness computation. It also shows initial promise in critical point tracking in practice. We review how the notion of multilevel robustness can improve the feature-tracking results in Section 2.2.2, and we give the pipeline to implement the multilevel robustness framework in Section 2.2.3. For simplicity and comparative purposes, we refer to this original multilevel robustness-based tracking [53] as the MRT framework in the remainder of this paper.

We now introduce a physics-informed feature selection strategy during preprocessing, making TROPHY much more efficient than the MRT framework in studying large-scale datasets. The MRT framework computes multilevel robustness for all detected critical points. For the example in Fig. 3 (A), 102,720 critical points and 22,743 tracks are detected with FTK. Suppose we set the number of levels $N=50$, the MRT framework needs to conduct $102720\times 50$ times classic robustness computation. TROPHY instead focuses on tracking cyclone-liked features, which selects a subset of critical points/tracks for in-processing. Considering the physical properties of real-world cyclones, we incorporate two criteria in physics-informed feature selection before multilevel robustness calculations.

First, since the duration of a tropical storm/cyclone is usually larger than 1 day, we require the selected tracks to contain at least one 1-day segment that consists of $+1$ degree critical points only.

Second, since a critical point representing the center of a cyclone usually has a high stability measure across time before it hits the land and dissipates, we require the segment detected from the previous step to have a high average robustness value. Based on domain knowledge, a tropical storm must have maximum sustained winds of at least 17.5 m/s. We set the threshold to 1.75 for track filtering.

With these two requirements, the numbers of critical points and tracks decrease to 9,784 and 86, respectively; compare Fig. 3 (A) and (B). Thus, TROPHY needs to compute multilevel robustness for only $9.5\%$ of critical points compared with the original MRT framework.

We next propose an adaptive level strategy for the multilevel robustness calculation. Such a strategy also considers physical properties of real-world cyclones and leads to more reasonable minimum multilevel robustness values for TROPHY than does the MRT framework. As discussed in Section 2.2.1, the MRT framework approximates the exact multilevel robustness with a set of evenly spaced radii. This strategy may lead to two consequences that are counterintuitive in real-world cyclone analysis. First, a critical point can be canceled with other critical points that are far away in the known data domain due to boundary effects; see partners for $x_{1}$ and $x_{2}$ in Fig. 2 (A) and (B). Although such cancellations are mathematically justifiable, it is almost impossible to happen in real-world scenarios since no perturbation could happen across the whole Atlantic Ocean. Second, the MRT framework may not find the true minimum multilevel robustness value due to sampling. It may also waste computational resources when an enlarged neighborhood does not include potential cancellation partners. For example, increasing the radius of the neighborhood from 60 to 80 for $x_{3}$ in Fig. 2 (C) does not lead to new cancellation candidates.

To mitigate the drawbacks of the MRT framework, we introduce an adaptive-level strategy. First, we set a physics-informed neighborhood size, where real-world perturbation could happen. We set this radius to be 10 degrees since hurricanes are among the most destructive real-world perturbations to wind field and are typically about 4.7 degrees wide. The TROPHY considers all possible cancellation partners within this neighborhood using varying radii. Fig. 5 (A) illustrates our adaptive level strategy in calculating the multilevel robustness. For a critical point $x$, we consider all critical points within its neighborhood at radius 10. Suppose there are $N$ ($N\geq 10$) critical points in this neighborhood and their Euclidean distances to $x$ are $\{a_{0},\dots,a_{N-1}\}$. We compute the classic robustness of $x$ within neighborhood defined by these radii, giving rise to its multilevel robustness. If $N<10$, we select an additional $10-N$ closest critical points outside of the selected neighborhood for more candidates of cancellation partners; see Fig. 5 (B). However, in most cases, the cancellation partners for the true minimum multilevel robustness are located in our physics-informed neighborhood, for example, $7.85$ degree for $x_{1}$, $4.15$ degree for $x_{2}$, and $6.05$ degree for $x_{3}$ in Fig. 2.

We now present feature selection aided by the minimum multilevel robustness and physical properties of TCs. We inherit one stability filter from the MRT framework. This filter considers the minimum multilevel robustness ${\mathrm{minR}}_{x}$ and its temporal stability in terms of lifespan. Let $l$ denote a logistic transformation of ${\mathrm{minR}}_{x}$, which maps the ${\mathrm{minR}}_{x}\in[0,\infty]$ to $l({\mathrm{minR}}_{x})\in[0,1]$. This normalization is defined as

where $k$ is the logistic growth rate; see [53, Sect. 5.1] for a detailed discussion on the benefit of this normalization and parameter selection for $k$. Now let $\gamma$ denote a track, $|\gamma|$ is its total length. The stability of a track $\gamma$ is defined as

where $T$ is the temporal span of the input dataset and $t_{\gamma}$ is the temporal span of $\gamma$. The first term in Eq. (1) captures the average pointwise stability (in a logistic scale), whereas the second term encodes the lifespan of the track. By definition, $b(\gamma)\in[0,1]$.

Our second feature selection strategy is referred to as maximum wind speed (MWS) filter. Since the storm systems are usually categorized by the Saffir–Simpson Hurricane Wind Scale, which considers only a hurricane’s maximum sustained wind speed, we use this filter to control the categories of storm that users want to visualize. Formally, for a track $\gamma$, its maximum wind speed is

where $\omega(x)$ is the MWS within the two-degree region of the center $x$.

Our third feature selection strategy is referred to as smoothness filter. It is used to filter out tracks that are too tortuous to be considered as TC tracks. Roughly speaking, we define the smoothness of a track by the average distance between the normalized track and its smooth univariate spline. Formally, for a track $\gamma$, its smoothness is

where $J$ is a normalization term mapping $x$ to $[0,1]\times[0,1]$ and $U$ maps $J(x)$ to the point on the smooth univariate spline of the normalized $\gamma$. TROPHY uses the SciPy Python library to calculate the (degree 3) smooth univariate splines for detected tracks. Fig. 6 shows three curves from our climate dataset and their smooth univariate splines. The left track has the lowest average distance and hence the highest smoothness value, whereas the right track has the lowest smoothness value and usually cannot be considered as a TC track.

In Fig. 7, we apply the above three filters to the ${\mathsf{ERA5\,2010}}$ dataset (see Sections 4 and C for details). After integrating multilevel robustness with FTK tracking results, we obtain tracks in Fig. 7 (A), which still suffer from visual clutter. By applying the stability filter with a threshold of $0.029$ for $b(\gamma)$, we obtain tracks in Fig. 7 (B). If we keep enlarging the stability threshold, track $a$ will be filtered out before track $b$, even if track $a$ represents the Category 1 hurricane Otto, and track $c$ cannot be found in the National Hurricane Center’s Tropical Cyclone Reports. Therefore, we use our MSF filter $w(\gamma)$ based on the MSF along the track and set the threshold to be $13.5$. The remaining tracks are shown in Fig. 7 (D). Again, if we keep increasing the threshold for $w(\gamma)$, track $c$ will be filtered out before track $d$, whereas track $c$ represents the Category 1 hurricane Lisa. Therefore, we apply our third smoothness filter $s(\gamma)$ in Fig. 7 (E) and (F). Now, TROPHY highlights 8 tracks representing hurricanes/tropical storms that can all be found in the National Hurricane Center’s Tropical Cyclone Reports.

We review several metrics used for evaluating TCs in climate data. See [56] for detailed descriptions.

For filters introduced in Section 3.3, we recommend a default value as a threshold for each filter function based on the TC tracks provided by IBTrACS. First, we calculate the spatial pattern of cumulative track density using TC tracks detected by TROPHY with thresholds of stability $b(\gamma)$ varying from 0.02 to 0.035, MWS $w(\gamma)$ varying from 10 to 14, and smoothness $s(\gamma)$ varying from 0.95 to 0.97. We also calculate the spatial pattern of cumulative track density from IBTrACS as a baseline. Then, we evaluate the similarity between the density patterns from TROPHY and IBTrACS using Pearson correlation, marked as $r_{xy,track}$, which is the most exhaustive measure of TC activity. We provide a detailed example for $r_{xy,track}$ computation in Section 5.5. We suggest the thresholds for $b(\gamma)$, $w(\gamma)$ , and $s(\gamma)$ with the values when $r_{xy,track}$ reaches its maximum, that is, 0.029 for stability $b(\gamma)$, 13.5 for MWS $w(\gamma)$, and 0.967 for smoothness $s(\gamma)$. TROPHY also allows users to fine-tune these thresholds independently for individual TC tracks of interest. We report performances of TROPHY with both default thresholds and human-in-the-loop thresholds in Section 5. We demonstrate that, although the tracking results using default threshold are encouraging, results with the human-in-the-loop option can further improve the performance of TROPHY.

We apply TROPHY to the ${\mathsf{ERA5\,Wind}}$ dataset on a cluster with 664 nodes (128GB DDR4 and 36 cores per node). We utilize the method of Tricoche et al. [41] to calculate degrees of critical points for the vector field in each time step and parallelize the computation of multilevel robustness with Eden [32], which can schedule and manage a number of tasks on a high-performance computing cluster.

Even though we have obtained the tracking results using TROPHY for all 30 years of data, we demonstrate TC tracking results for 1981, 1990, 2004, and 2009 in Fig. 8 to highlight the tracking differences between various datasets and algorithms. TempestExtremes tracks TC eyes as positions with the local minimum SLP, whereas the TC eyes from TROPHY are located where wind speeds are zeros. Therefore, TC tracks from TempestExtremes and TROPHY do not overlap exactly. Both methods consider the area within a two-degree great circle of a detected TC eye to find a local maximum and use it to measure storm intensity [55, 44]. Tracks in Figs. 8, 9 and 10 are colored by storm intensity w.r.t. this two-degree maximum wind speed, referred to as wind speed for simplicity. We discuss some preliminary findings based on Fig. 8 in the following paragraphs

First, since TempestExtremes requires both dynamic and thermodynamic variables to meet its specified criteria, it usually cannot detect the beginnings and endings of TC where wind speeds are low. TROPHY can detect such tails as long as the centers of flows can be identified as critical points; see tracks $a$, $b$, and $c$ in Fig. 8.

Second, TROPHY can detect some discrete storm events that are not captured by TempestExtremes, such as tracks $d$ and $e$ in Fig. 8. Track $d$ can be found in IBTrACS, whereas track $e$ cannot. To further investigate these discrepancies, including tracks that cannot be captured by TempestExtremes/IBTrACS but are found by TROPHY, and tracks that are in IBTrACS but are undetected by TROPHY, we conducted case studies (reported in Sections 5.2 and 5.3).

Third, using fine-tuned filter thresholds for TROPHY, we can get TC tracking results that are more similar to IBTrACS. For example, track $f$ is missed by TROPHY with a default threshold for smoothness, since $f$ is severely bent and has a relatively low smoothness value. When, however, we decrease the threshold for $s(\gamma)$ from 0.967 to 0.95, track $f$ is shown in the TC tracking result as indicated by the red arrow in Fig. 8. Similarly, track $e$, which is not recorded in IBTrACS, can be filtered out if we increase the threshold of stability $b(\gamma)$ from 0.029 to 0.03. These sensitivities to parameters in the tracking algorithm are, in general, expected [15]. This human-in-the-loop option provides users the opportunity of fine-tuning, especially for short-term forecasting or weather-scale studies, as accuracy is more important at this scale. At the climate scale and for climate change impacts, users may use the same thresholds for all TC cases in both historic and future periods, to avoid the impacts of parameter uncertainties.

We now investigate a Category 2 hurricane track named Hurricane Bonnie. As illustrated in Fig. 9, (A) shows the TC tracks from TempestExtremes (colored by wind speed) and IBTrACS (in black), whereas (B) shows the TC tracks detected by TROPHY. We observe a much longer tail in (B) with low wind speed compared with (A). (C) gives the zoomed-in views of tails of TC tracks detected by TempestExtremes (1st row) and TROPHY (2nd row). Vector fields from different time steps are shown as background. We also highlight the TC eyes detected by each method in corresponding time steps by arrows on (C). We see TC eyes detected by TROPHY are exactly located where the wind vanishes, whereas the TC eyes from TempestExtremes are slightly shifted to the outside of the eyes. The last two columns of (C) show that TROPHY can track the TC eyes during dissipation because the wind flow is still spinning even if the wind speed is low. TROPHY tracks such flow behaviors as sources until the sources disappear.

We next investigate a Category 1 hurricane named Hurricane Lisa, which is captured by TempestExtremes but is partially missed by TROPHY. Fig. 10 (A–C) show TCs of 1998, which are recorded in IBTrACS and detected by TempestExtremes and TROPHY with the ${\mathsf{ERA5\,1998}}$ dataset, respectively. We focus on the analysis of Hurricane Lisa, marked as $\gamma_{1}$, $\gamma_{2}$, and $\gamma_{3}$ in different methods. We note that $\gamma_{3}$, detected by TROPHY, has a clearly shorter track compared with the other two datasets because the eyes of Hurricane Lisa were detected as centers only up to 10/08/1998 00:00 AM, after which Hurricane Lisa lost its distinct eye and could not be detected by TROPHY; see Fig. 10 (M), as well as a zoomed-in view in Fig. 10 (d).

Such a process could be related to TC’s extratropical transitions [25], which may occur as Hurricane Lisa moves over cooler water and into areas of stronger wind shear at higher latitudes. During this transition, the cyclone loses its symmetric and distinct eye and, thus, TROPHY terminates its tracking.

In fact, Hurricane Lisa produced distinct eyes again at time step 10/09/1998 at 12:00 PM and 10/10/1998 at 00:00 AM, so TROPHY was able to start the track again. However, because these tracks were short-lived, they were filtered out during the initial feature selection.

This section and the following present the evaluations of TROPHY using the metrics described in Section 4.1. For annual climatology, cumulative statistics are calculated over the entire data period (30 years), which are then normalized to a per-year basis. Table 1 shows annually averaged statistics for TCs detected by TROPHY and TempestExtremes.

According to annual TC count ($\overline{count}$), IBTrACS contains approximately 11 TCs per year within the studied region, whereas TROPHY with default setting produces the least number of TCs. However, when considering annual storm lifetime (i.e., $\overline{tcd}$), TROPHY detects longer storms than TempestExtremes and is closer to IBTrACS on average. This result is consistent with what we observed in our experiments. TROPHY can detect TCs not only during their movements but also, at least partially, during their formation and dissipation periods, even if their wind speed is low; see Fig. 9 for an example. Also, since TempestExtremes requires dynamic and thermodynamic variables to meet specified criteria, TC tracks may break into pieces when some parts of the track do not meet the requirement; see the track indicated by a purple circle from Fig. 8 for an example. Overall, TROPHY produces TC tracks for a longer time than does TempestExtremes, whereas TempestExtremes detects more TC tracks than does TROPHY.

In addition, although TROPHY’s TC tracks have closer hurricane genesis $\overline{lmi}$ w.r.t. the observation, both TROPHY and TempestExtremes show a poleward bias compared with IBTrACS’s genesis. One of the reasons for this bias could be that the reanalysis data including ERA5 are not able to simulate storm structures near the equator well, according to Knaff et al. [25], indicating that the input data to any tracking algorithms is the most important factor when studying TCs.

For spatial pattern climatology, density maps are generated by aggregating occurrences into $4^{\circ}$ by $4^{\circ}$ bins. An example of the spatial track density patterns is shown in Fig. 11. The spatical correlation between IBTrACS and TROPHY (default) is caculated using patterns from Fig. 11 (A) and (C). To quantify the similarity between IBTrACS and TROPHY (or TempeExtremes), we calculate Pearson correlation coefficient $r_{xy}$ for total occurrence ($r_{xy,track}$), genesis occurrence ($r_{xy,gen}$), maximum wind speed, and minimum SLP. Table 2 shows the pattern correlation of TROPHY/TempestExtremes with observations.

Both TROPHY and TempestExtremes can produce reasonable distributions of storm occurrence when compared with observations ($>0.89$). The spatial patterns of genesis ($r_{xy,gen}$) from TROPHY and TempestExtremes show lower correlations ($<0.64$) with observations. This is not surprising because the initialization and development of a clear eye depend on many other factors such as their 3D evolutions, which are beyond the capability of what our 2D vector fields can represent. In fact, predicting the genesis of TC is one of the scientific challenges in TC research fields, indicating that the currently available data cannot capture TC genesis and that new frameworks of such calculations are needed [54]. Both TROPHY and TempestExtremes are able to produce storm strength similar to that of observations ($\geq 0.89$) in terms of maximum $u_{10}$ and minimum sea level pressure.