A Dynamical Systems Perspective to Characterize the El Niño Diversity in Spatiotemporal Patterns

El Niño is the long-lasting anomalous warming of sea surface temperature (SST) and surface air temperature (SAT) over the tropical Pacific. Each El Niño event has its unique impact on the overlaying atmosphere, where the warming exhibits diversity in spatiotemporal patterns. It still remains an open question for discriminating the El Niño diversity, since the single area-averaging SST index often fails to distinguish the impact of the event diversity, which is partially due to the nonlinear and non-uniform variations of the warming patterns. Here, we introduced the Dynamical Systems metrics (DSMs) to measure instantaneous dimensions and persistence of the SAT warming patterns over the tropical Pacific. Our results show that different SAT warming patterns can be discriminated by their corresponding values of dimension and persistence, then the central Pacific and eastern Pacific El Niño events can be discriminated by DSM. Particularly, through the analyses of El Niño events, we can interpret the physical meaning of DSM parameters applied to the space-time SAT field: an instantaneous dimension reflects whether the sub-regions of the SAT field are consistently varying and to what degree the spatial pattern of anomalies is homogeneous, while the instantaneous persistence indicates how long an anomalous SAT pattern can be maintained. This work analyzes the spatiotemporal variability of El Niño from a dynamical system perspective, and DSM may also serve as a useful tool to study extreme events related to SST anomalies.


INTRODUCTION
El Niño is a typical extreme event causing a strong influence on the worldwide weather and climate [1][2][3]. Each El Niño event can last for several months, in which sea surface temperature (SST) over the tropical Pacific gets anomalously warm compared to that in the normal periods. The spatial patterns of SST anomaly (SSTA) are found to have the wide diversity in different types of El Niño events [4][5][6], i.e., the central Pacific (CP) El Niño which has strong warm SSTAs only located in the central Pacific and the eastern Pacific (EP) El Niño whose warm SSTAs distribute over both central and eastern Pacific. The difference in the SSTA spatial patterns can initiate the anomalous atmospheric convections with different strengths and locations [7,8], further causing different global climate effects via atmospheric waves [3,9]. Therefore, the El Niño diversity regarding the spatiotemporal patterns and the impact on atmospheric states has been an important topic in the climate study (e.g., [2,[10][11][12][13]).
El Niño is usually quantified by the area-averaged SSTAs over four given regions ( Figure 1A), i.e., Niño1+2, Niño3, Niño4, and Niño3.4 indices [4,14]. These four Niño indices essentially denote the temporal variations of the sea surface temperature over the four given spatial regions, respectively ( Figure 1A). However, these SST-based Niño indices are obscure in terms of dynamical states of atmospheric fields and then they cannot always accurately infer the El Niño impact on the atmosphere [3,8,13,15]. On the other hand, the dynamical evolution features of a meteorological field have been encoded into both of its spatial and temporal structures, whereas the area-averaged procedure misses this spatial information of meteorological fields [13,16]. Thus, the dynamical evolution features of different El Niño events cannot be adequately reflected by Niño indices [8,12,13,15]. Accordingly, further focusing on atmospheric variables (e.g., surface air temperature) and spatiotemporal variability of the meteorological fields can be favorable to an in-depth understanding of the El Niño diversity and its impact on the dynamical states of the atmosphere.
In this study, we will utilize the Dynamical Systems metrics (DSMs, [17]) to analyze the spatiotemporal variability of the surface air temperature (SAT) over the tropical Pacific during El Niño episodes. DSM was proposed to quantify the instantaneous dynamical states of the climate subsystem over a given region, by incorporating both the spatial and temporal structures of the interested meteorological field. Different from the empirical orthogonal function (EOF), which relies on the linear covariance of the data, the DSM algorithm regards spatiotemporal data as the dynamical trajectory in the state space, and thus both linear and nonlinear statistical features can be characterized by DSM [13,17]. In recent studies, DSM has been applied to various issues in the climate study (e.g., [13,[18][19][20]), where the instantaneous dynamical properties of the interested meteorological field can be effectively measured, and DSM is considered to be a promising tool to analyze the extreme events in the ocean and atmosphere.  Figure 1A). The zoomed sphere (black) shows the idealized schematic of the instantaneous persistence and dimension around ζ t0 .
Frontiers in Physics | www.frontiersin.org June 2022 | Volume 10 | Article 919951 Accordingly, two motivations of this study are to identify the dynamical properties of the SAT anomalies within different El Niño states and to interpret the physical meaning of DSM through analysis of the SAT over the tropical Pacific. On the one hand, the DSM applied to analyze SAT over tropical Pacific may address the aforesaid deficiency of Niño indices in the El Niño study. On the other hand, the interpretation of DSM could contribute further to a better understanding of this newly developed approach and provide a technical reference for future related climate studies.
The rest of this article is arranged as follow: Section 2 will describe the used data and DSM. Section 3 is divided into three parts: Section 3.1 will compare the Niño indices with the instantaneous dynamical parameters quantified by DSM, thus understanding how the dynamical states of SATA are influenced by the extreme events. Based on DSM, Section 3.2 will focus on the SAT anomaly (SATA) over the tropical Pacific during El Niño episodes, with regard to the diversity in spatial and temporal patterns of the SATA. Particularly, through analysis of the spatiotemporal patterns of SATA, the physical meaning of the DSM parameters in terms of the meteorological field will be demonstrated. Then, Section 3.3 will utilize the DSM parameters to compare different El Niño events. Further discussions and conclusions will be made in Section 4.

DATA AND DYNAMICAL SYSTEMS METRICS
This study uses the daily mean 2-m temperature data (i.e., SAT) from the European Centre for Medium-Range Weather Forecasts' ERA-Interim reanalysis [21], with a temporal range from 1979-01-01 to 2019-08-31. To focus on the anomalous variations of SAT, here, we subtract the long-term climatological average for each calendar day over the whole span from the 2-m temperature data. Additionally, the daily time series of Niño1+2, Niño3, Niño4, and Niño3.4 indices are obtained from the KNMI Climate Explorer website: http://climexp.knmi.nl. The spatial domains of Niño1+2, Niño3, Niño4, and Niño3.4 regions are, respectively, illustrated in Figure 1A. To cover these regions and measure the SATA dynamical states over the tropical Pacific, the spatial domain for the DSM computation is set at 160°E-80°W and 10°N-10°S ( Figure 1A). Noting that DSM has been demonstrated to be not sensitive when the spatial resolution is set as 0.75°, 1.5°, and 2.5° [19], here, we select a horizontal resolution of 0.75°for the SATA data.
In the methodology of DSM [17], the evolution of a given meteorological field is acting as a dynamical system, in which the temporal evolutions of states of this dynamical system are described by an observed state-space trajectory X(t) (i.e., a matrix embedded with the spatiotemporal data). In this study, X(t) is a succession of daily SATA latitude-longitude maps within the selected spatial domain (see black box lines in Figures 1A), and a certain point ζ t0 on this trajectory is a SATA latitude-longitude map for a specific day ( Figure 1B).
Aiming at ζ t 0 , we can infer the instantaneous dynamical property of the system at this specific day. The beginning step is to calculate the distance function g(X(t), ζ t0 ) −log[dist(X(t), ζ t0 )], in which "dist" represents the Euclidean distance between ζ t 0 and any other point on the trajectory X(t). g(X(t), ζ t 0 ) is a single-column matrix whose collected element takes large values when the point on the trajectory X(t) is close to ζ t 0 . Then, we can utilize this matrix to sort out the analogs of ζ t 0 on the state-space trajectory: by setting the top 98% quantile of g(X(t), ζ t 0 ) as the threshold s X , we retain the exceedances of u(t, ζ t0 ) g(X(t), ζ t0 )-s X but leave out the part at which u(t, ζ t0 ) < 0. It is noted that DSM is not sensitive to the choice of the quantile ranging from 90% to 99%. Thus, u(t, ζ t0 ) collects the analogs of ζ t0 , and its cumulative probability distribution P[u(ζ t0 )] satisfies the generalized Pareto distribution function as [17]: As reported in the previous studies (19,22,23), the two parameters in Eq. 1 can be used to indicate the instantaneous dynamical properties of the system. 0 < θ(ζ t0 ) < 1 is interpreted as the inverse of instantaneous persistence: θ(ζ t0 ) taking lager values means that residence time of the system around ζ t0 is shorter, while its low values denote the long residence time. Another derived parameter d(ζ t0 ) 1/σ(ζ t0 ) is interpreted as the instantaneous dimension, denoting the degree of freedom of the system around ζ t 0 . In terms of the SATA spatiotemporal field, the instantaneous persistence and dimension for each day are estimated, and then two time series θ SATA (t) and d SATA (t) are obtained. The value of θ is estimated by using the Suveges likelihood estimator [24], and the programing code of The aforementioned interpretations for the DSM parameters are in the sense of dynamical systems theory [18,23]; then, the recent studies have demonstrated that the values of instantaneous persistence and dimension are closely related to the intrinsic predictability of a certain climate subsystem ( [19,20]). In this study, we will use θ SATA (t) and d SATA (t) to investigate the underlying variations of the SATA dynamical states during El Niño episodes, and in turn the physical meaning of θ SATA (t) and d SATA (t) in sense of the structural spatiotemporal patterns will be discussed.

The Extreme Surface Air Temperature Anomaly States Measured by Dynamical Systems Metrics
We first analyze how the extreme SSTA states over the tropical Pacific influence the dynamical properties of SATA. The extreme SSTA states are usually inferred by referring to Niño indices. Though the four Niño indices are highly correlated to one another [5,14], extreme events from four different Niño indices are not completely consistent. For example, when the Niño1+2 index exceeds its 90% quantile, all of Niño3, Niño3.4, and Niño4 indices do not reach their own 90% quantiles ( Figures  2A,B). This also means that extreme SSTA events can emerge from different sub-regions of the tropical Pacific [13]. Then, we observe how the dynamical properties of SATA vary with these extreme SSTA states defined by different Niño indices ( Figures  2C,D). Niño1+2 and Niño4 indices denote the SSTA states in the eastern and western tropical Pacific, respectively ( Figure 1A), and the comparison in Figure 2 shows that when either Niño1+2 index or Niño4 index reaches an extreme state, both θ SATA and d SATA will all take their low values. This suggests that the extreme SSTA states over the tropical Pacific can increase the SATA persistence but decrease the SATA dimension.
Taking a segment of the SATA temporal evolution as an example (Figure 3), the covariation between DSM parameters and SATA spatiotemporal pattern can be observed. During the first 200 days of this segment, the zonal distribution of SATA is inhomogeneous and the SATA signals over different sub-regions are weak, at which θ SATA and d SATA remain at high values. Since the 240th day of this segment, the central and eastern sub-regions of equatorial Pacific get anomalously warm, and at this time θ SATA and d SATA start to decrease. During the days from 330 to 480 of this segment, the central and eastern sub-regions of the equatorial Pacific exhibit a consistent warm pattern, at which θ SATA and d SATA remain at low values. This means that the persistence and dimension of SATA are closely related to the spatiotemporal pattern of SATA, and the low values of θ SATA and d SATA may indicate the extreme states of SATA over the tropical Pacific. These results are reasonable. During the extreme SSTA events, the intensified air-sea interactions can cause the extreme SATA states [12,13], and the positive air-sea-coupled feedback can increase the persistence of the atmosphere over the tropical Pacific [13,25,26].

Instantaneous Persistence and Dimension Within Different El Niño States
To further understand the association between the SATA dynamical properties and the different El Niño states, we inspect how θ SATA and d SATA vary with the Niño3.4 index ( Figure 4A). During 1981-2019, every time the Niño3.4 index exceeds its threshold (±0.5°C), the values of θ SATA and d SATA will  decrease. This is consistent with the result shown in Figure 3. However, the degree of changes in θ SATA and d SATA is not identical during different El Niño events. Here, we mark five well-recognized CP and EP El Niño events [5,12,13] in Figure 4A, and it is obvious that the values of θ SATA and d SATA in EP events are smaller than those in CP events. We then pick out all El Niño states by referring that the Niño3.4 index exceeds 0.5°C for at least five consecutive months and analyze these El Niño states in the scatter plot of θ SATA versus d SATA ( Figure 4B). The CP and EP El Niño states, respectively, cluster into two different groups in this scatter plot, and this means that the EP El Niño event can cause the stronger impact on the SATA dynamical properties than the CP.
Focusing on the SATA dynamical properties during El Niño episodes ( Figure 4B), θ SATA and d SATA are found to vary widely, which corresponds to the variety of El Niño states and spatial warming patterns. We can further investigate these samples toward two informative aspects: the first one is to observe how spatiotemporal patterns of the SATA vary with the different El Niño states. The second one is an in-depth understanding in techniques, i.e., what the varied values of DSM parameters mean for the meteorological spatiotemporal field. Here, we divide these θ SATA values into three groups: the smallest 1/3 (named θ Q1 group), the middle 1/3 (the θ Q2 group), and the largest 1/3 (the θ Q3 group), as the histogram of Figure 4B shows. Similarly, we also divide the samples of d SATA into three groups (i.e., the d Q1 , d Q2 , and d Q3 groups) according to the three quantiles of the d SATA values. In the following, we will analyze the different groups of θ SATA and d SATA , respectively.
To show the spatiotemporal patterns of the SATA field corresponding to different values of θ SATA , here, we create the Hovmöller diagrams of SATA by, respectively, compositing the three groups (i.e., the θ Q1 , θ Q2 , and θ Q3 groups), as Figures 5A-C show. As the values of θ SATA increase, the magnitudes of SATA along the equator decrease, and meanwhile the warming signal in the SATA field gets less persistent. We also utilize time-lagged correlation to analyze the persistence of SATA within these three groups. Taking the samples in the θ Q1 group as an example, we calculate spatial correlations between the SATA field at a certain day and its time-lagged ones (see Text. Supplementary Material S1) and then obtain the average result of the θ Q1 group, as shown in Figure 5D. Comparing the correlation curves of the three groups, it is found that as θ SATA increases, the 1-day-lagged correlation coefficient decreases, and the correlation curve decays faster with time lags. The 90% confidence intervals of these three curves (Supplementary Figure S1) are clearly distinguishable from one another (i.e., confidence intervals are not overlapped with one another), which suggests that the SATA persistence magnitudes within these three groups are significantly different. Here, the results demonstrate that θ SATA can seamlessly reflect the varying temporal persistence of the SATA field. Toward the analysis on d SATA , in Figure 6A, we composite the spatial patterns of SATA corresponding to three groups, respectively (i.e., the d Q1 , d Q2 , and d Q3 groups). For the SATA spatial pattern of the d Q1 group, both central and eastern regions of tropical Pacific get significantly warm, while only the central tropical Pacific gets warm for the d Q3 group. This means that the low dimension corresponds to the more homogeneous spatial pattern of the meteorological field. Furthermore, we would investigate whether d SATA is also related to the synchrony of different grids in the SATA field, i.e., whether the temporal evolutions of the different grids are synchronized. Centered on a grid at the central Pacific (150°W, 0°), we calculate the crosscorrelations of SATAs between other grids and this centered grid, in which results are derived for the three groups, respectively ( Figure 6B). This result shows that the synchrony of this SATA field would decrease with the increased d SATA . We also checked the cross-correlations by selecting different centered grids and found that the conclusion was not sensitive to the choice of the centered grid (Figure not shown here). These results demonstrate that d SATA can reflect whether the sub-regions of a meteorological field are consistently varying and to what degree the spatial pattern is homogeneous.
The aforementioned results interpret the physical meaning of θ SATA and d SATA in the meteorological field, i.e., they characterize the variability of spatiotemporal patterns of SATA over the tropical Pacific. Moreover, the results shown in Figures 5, 6 demonstrate the impact of different El Niño states on the SATA spatiotemporal patterns. To observe the corresponding El Niño states, here, we further analyze the SST anomalies downloaded from NOAA High-resolution Blended Analysis of Daily SST dataset [27] and composite the SSTA spatial patterns of the d Q1 , d Q2 , and d Q3 groups, respectively (Supplementary Figure  S2). Corresponding to the d Q1 group, the SSTs are found to be anomalously warm over the central and eastern regions of the tropical Pacific (Supplementary Figure S2A), while the d Q3 group corresponds to the SSTA pattern at which only the central tropical Pacific gets warm (Supplementary Figure  S2C). Since that θ SATA and d SATA are highly correlated over the tropical Pacific ( Figure 4B), the composite SSTA patterns for the θ Q1 , θ Q2 , and θ Q3 groups are similar to those for the three groups of d SATA (Figure not shown here). This reconfirms the conclusion obtained from Figure 4, i.e., the EP and CP types of El Niño can cause different impact on the persistence and dimension of the SATA field over the tropical Pacific.

Discriminating Different El Niño Events
The aforementioned section demonstrates that θ SATA and d SATA can discriminate the different El Niño impacts on the SATA spatiotemporal patterns, and the impact is found to be related to the EP and CP types of El Niño states. In turn, we can employ θ SATA and d SATA to further compare different El Niño events. In Figure 7A, we present the mean values of θ SATA and d SATA within ten El Niño events during 1981-2019. In the scatter plot of θ SATA versus d SATA , the aforementioned typical EP and CP El Niño events are in two distinct divisions, among which the 2015/16 event is between these two divisions. In terms of magnitudes reflected by the Niño3.4 index (e.g., Figure 3A), the 2015/16 event is a super El Niño, and the corresponding SSTA pattern is close to an EP El Niño. However, recent studies found that the air-sea interactions and climate impact of the 2015/16 event are weaker than those of the typical EP El Niño [12,13,28,29]. Here, θ SATA and d SATA confirm that the 2015/16 event has caused the weaker impact on the SATA spatiotemporal pattern than the EP El Niño but stronger than the CP ones. This result further confirms that the 2015/16 event may be a mixed-type El Niño. The previous studies utilized the scatter plot of Niño3 and Niño4 indices to classify the types of El Niño [4,5,30]. For comparison, in Figure 7B, we present the mean values of Niño3 and Niño4 indices within the ten El Niño events. The diagonal line in this scatter plot can divide the typical EP and CP El Niño events into two divisions, whereas the events in 1986/87, 1991/92, 2006/07, and 2015/16 are close to the diagonal line. A similar phenomenon can be observed in the scatter plot of θ SATA versus d SATA (Figure 7A), in which these four events are also between the divisions of EP and CP events. In the two kinds of scatter plots, the 2002/03 event is in the division of CP events, which is consistent with recent studies where this event is found to be closer to the CP type [12,31]. These results show that θ SATA and d SATA can indeed discriminate different El Niño events.

DISCUSSIONS AND CONCLUSION
In this study, we employ the DSM parameters (i.e., θ SATA and d SATA ) to quantify the El Niño impact on the SATA spatiotemporal patterns over the tropical Pacific. The results show that the extreme SSTA events can increase/decrease the persistence/dimension of SATA spatiotemporal patterns, and the impact of EP and CP El Niño events is different. Furthermore, we find that θ SATA and d SATA can adequately discriminate different El Niño events (i.e., the EP, CP, and mixed-type events). This reconfirms that the El Niño impact on the atmosphere is closely related to the properties of events [12,13,31].
Additionally, in Section 3.2, we systematically discuss the physical meaning of θ SATA and d SATA when they are applied to analyzing the meteorological field. θ SATA denotes the persistence of the SATA patterns, while d SATA can reflect whether the sub-regions of the SATA field are consistently varying and to what degree the spatial pattern is homogeneous. This result can be useful for the future climate study when applying the DSM parameters.
In addition to El Niño, the extreme SSTA events over the tropical Pacific can be La Niña events, the coastal warming events, and the Pacific meridional mode events, etc. [4,13,[32][33][34], but they are not analyzed in this study. From the results of DSM, it is shown that La Niña and the coastal warming events have also reduced the values of θ SATA and d SATA ( Figure 4A). Particularly, noting the dots marked by the black box lines in Figures 2A, 4A, these dots are explained by the 1983 coastal warming event [13,35].
Supplementary Figure  S3 illustrates the spatiotemporal evolution of this event, compared with the time series of θ SATA and d SATA . During this event, most parts of the equatorial Pacific are not getting anomalously warm, but only the coastal areas are extremely warm, which is not the canonical El Niño. However, the time series of θ SATA and d SATA show that the impact on SATA of this coastal warming event is comparable to that of the preceding El Niño event (Supplementary Figures S3B, S3C). This is in line with the previous studies about the extreme impact of coastal warming events on the weather and climate [13,35]. These results further confirm that the DSM parameters can faithfully quantify the impact of these extreme events, and the wide application of this method and the relevant extreme events deserve an in-depth study in the future. The current article only focuses on the DSMbased analysis of the El Niño events during 1981-2019; the more detailed studies for different kinds of extreme SSTA events, including their potential changes in the future climate change scenarios, will be addressed in our following work.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material; further inquiries can be directed to the corresponding author.