{'fraction_non_alphanumeric': 0.047029198837230834, 'fraction_numerical': 0.034920105908274546, 'mean_word_length': 4.2744140625, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
An Eddy-Zonal Flow Feedback Model for Propagating Annular Modes

The variability of the zonal-mean large-scale extratropical circulation is often studied using individual modes obtained from empirical orthogonal function (EOF) analyses. The prevailing reduced-order model of the leading EOF (EOF1) of zonal-mean zonal wind, called the annular mode, consists of an eddy-mean flow interaction mechanism that results in a positive feedback of EOF1 onto itself. However, a few studies have pointed out that under some circumstances in observations and GCMs, strong couplings exist between EOF1 and EOF2 at some lag times, resulting in decaying-oscillatory, or propagating, annular modes. Here, we introduce a reduced-order model for coupled EOF1 and EOF2 that accounts for potential cross-EOF eddy-zonal flow feedbacks. Using the analytical solution of this model, we derive conditions for the existence of the propagating regime based on the feedback strengths. Using this model, and idealized GCMs and stochastic prototypes, we show that cross-EOF feedbacks play an important role in controlling the persistence of the annular modes by setting the frequency of the oscillation. We find that stronger cross-EOF feedbacks lead to less persistent annular modes. Applying the coupled-EOF model to the Southern Hemisphere reanalysis data shows the existence of strong cross-EOF feedbacks. The results highlight the importance of considering the coupling of EOFs and cross-EOF feedbacks to fully understand the natural and forced variability of the zonal-mean large-scale circulation.


Introduction
At the intraseasonal to interannual time scales, the variability of the large-scale atmospheric circulation in the mid-latitudes of both hemispheres is dominated by the "annular modes", which are usually defined based on empirical orthogonal function (EOF) analysis of zonalmean meteorological fields (e.g., Kidson 1988;Wallace 1998, 2000;Lorenz andHartmann 2001, 2003;Thompson and Woodworth 2014;Thompson and Li 2015). The barotropic annular modes are often derived as the first (i.e., leading) EOF (EOF1) of zonal-mean zonal wind, which exhibits a dipolar meridional structure and describes a north-south meandering of the eddy-driven jet. Note that in this paper, the focus is on the barotropic annular modes, hereafter simply called annular modes (see Thompson and Woodworth 2014;Barnes 2014, andLi 2015 for discussions about the "baroclinic annular modes"). The second EOF of zonal-mean zonal wind (EOF2) has a tripolar meridional structure centered on the jet, describing a strengthening and weakening of the eddy-driven jet (i.e., jet pulsation). By construction, EOF1 and EOF2 (and any two EOFs) are orthogonal and their associated time series (i.e., principal components, PCs), sometimes called zonal index, are independent at zero time lag.
The persistence of the annular mode (EOF1) and its underlying dynamics have been the subject of extensive research and debate in the past three decades (e.g., Robinson 1991;Branstator 1995;Feldstein and Lee 1998;Robinson 2000;Lorenz andHartmann 2001, 2003;Gerber and Vallis 2007;Gerber et al. 2008b;Chen and Plumb 2009;Simpson et al. 2013;Zurita-Gotor 2014;Nie et al. 2014;Byrne et al. 2016;Ma et al. 2017;Hassanzadeh and Kuang 2019). Many of the aforementioned studies have pointed to a positive eddy-zonal flow feedback mechanism as the source of the persistence: The zonal wind and temperature anomalies associated with the annular mode (EOF1) modify the generation and/or propagation of the synoptic eddies at the quasi-steady limit (greater than 7 days) in such a way that the resulting eddy fluxes reinforce the annular mode (see Hassanzadeh and Kuang (2019) and the discussion and references therein). Most notably, Lorenz and Hartmann (2001) developed a linear eddy-zonal flow feedback model (LH01 model hereafter) for the annular modes by regressing the anomalous eddy momentum flux divergence onto the zonal index of EOF1 (z 1 ) and interpreting correlations between z 1 (t) and regressed momentum flux divergence (m 1 (t)) at long lags (greater than 7 days) as evidence for eddy-zonal flow feedbacks, i.e., feedbacks of EOF1 onto itself. Lorenz and Hartmann (2001) developed a similar model, separately, for EOF2 and found, respectively, positive and weak eddy-zonal flow feedbacks for EOF1 and EOF2, respectively, consistent with the longer persistence of EOF1 compared to EOF2. Such single-EOF eddy-zonal flow feedback models have been used in most of the subsequent studies of the annular modes (e.g., Lorenz and Hartmann 2003;Simpson et al. 2013;Lorenz 2014;Robert et al. 2017;Ma et al. 2017;Boljka et al. 2018;Hassanzadeh and Kuang 2019;Lindgren et al. 2020).
While EOF1 and EOF2 are independent at zero lag, a few previous studies have pointed out that these two EOFs can be correlated at long lags (e.g., greater than 10 days), and that in fact the combination of these two leading EOFs represents coherent meridional propagations of the zonal-mean flow anomalies. Such propagating regimes have been observed in both hemispheres in reanalysis data (e.g., Feldstein 1998;Feldstein and Lee 1998;Sheshadri and Plumb 2017). Anomalous poleward propagation of zonal wind typically emerges in low latitudes and mainly migrate poleward over a few months, although non-propagating regimes can also appear in some instances (see Fig. 1 of Sheshadri and Plumb 2017 and Fig. 6 in this paper). Similar behaviors have also been reported by in general circulation models (GCMs) (e.g., James and Dodd 1996;Son and Lee 2006;Son et al. 2008;Sheshadri and Plumb 2017). Son and Lee (2006) found that the leading mode of variability in an idealized dry GCM can be either the propagating or non-propagating regime depending on the details of thermal forcing imposed in the model. They also found that unlike the non-propagating regimes, the z 1 and z 2 of the propagating regimes are strongly correlated at long lags, peaking around 50 days (see their Fig. 3; also Figs. 4b of the present paper). Furthermore, Son and Lee (2006) reported that non-propagating regimes are often characterized by a single time-mean jet with a dominant EOF1 (in terms of the explained variance) while the propagating regimes are characterized by a double time-mean jet in the midlatitudes with the variance associated with EOF2 being at least half of the variance of EOF1. Furthermore, Son et al. (2008) found that the e-folding decorrelation time scale of z 1 in the propagating regime to be much shorter than that of the non-propagating regime. The long e-folding decorrelation time scales for the annular modes in the nonpropagating regime were attributed to an unrealistically strong positive EOF1-onto-EOF1 feedback, while the reason behind the reduction in the persistence of the annular modes in the propagating regime remained unclear.
More recently, Sheshadri and Plumb (2017) presented further evidence for the existence of propagating and nonpropagating regimes and strong lagged correlations between z 1 and z 2 in reanalysis data of the Southern Hemisphere (SH) and in idealized GCMs. Moreover, they elegantly showed, using a principal oscillation patterns (POP) analysis (Hasselmann 1988;Penland 1989), that EOF1 and EOF2 are in fact manifestations of a single, decaying-oscillatory coupled mode of the dynamical system. Specifically, they found that EOF1 and EOF2 are, respectively, the real and imaginary parts of a single POP mode, which describes the dominant aspects of the spatiotemporal evolution of zonal wind anomalies. Sheshadri and Plumb (2017) also showed that in the propagating regime, the auto-correlation functions of z 1 and z 2 decay non-exponentially.
Given the above discussion, a single-EOF model is not enough to describe a propagating regime because the EOF1 and EOF2 in this regime are strongly correlated at long lags and that the auto-correlation functions of the associated PCs do not decay exponentially (but rather show some oscillatory behaviors too). From the perspective of eddy-zonal flow feedbacks, one may wonder whether there are cross-EOF feedbacks in addition to the previously studied EOF1 (EOF2) eddy-zonal flow feedback onto EOF1 (EOF2) in the propagating regime. In cross-EOF feedbacks, EOF1 (EOF2) changes the eddy forcing of EOF2 (EOF1) in the quasi-steady limit. Therefore, there is a need to extend the single-EOF model of LH01 and build a model that includes, at a minimum, both leading EOFs and accounts for their cross feedbacks. The objective of the current study is to develop such a model and to use it to estimate effects of the cross-EOF feedbacks on the variability of propagating annular modes.
The paper is structured as follows: Section 2 compares the characteristics of z 1 , z 2 , m 1 , and m 2 for the nonpropagating and propagating annular modes in reanalysis and idealized GCMs. In Section 3 , we develop a linear eddy-zonal flow feedback model that accounts for cross-EOF feedbacks, validate the model using synthetic data from a stochastic prototype, discuss the key properties of the analytical solution of this model, and apply this model to data from reanalysis and an idealized GCM. The paper ends concluding remarks in Section 4.

Propagating annular modes in an idealized GCM and reanalysis
In this section, we will examine the basic characteristics and statistics of propagating annular modes in an idealized GCM (the dry dynamical core) and reanalysis. We focus on the southern annular mode, which makes it easier to compare the results of the reanalysis and the idealized aquaplanet GCM simulations. We will start with the idealized GCM to demonstrate the characteristics of the propagating and non-propagating annular modes.
a. An idealized GCM: The dry dynamical core We use the Geophysical Fluid Dynamics Laboratory (GFDL) dry dynamical core GCM. The GCM is run with a flat, uniform lower boundary (i.e., aquaplanet) with T63 spectral resolution and 40 evenly spaced sigma levels in FIG. 1. One-point lag-correlation maps of the vertically averaged zonal-mean zonal wind anomalies u , reconstructed from projections onto the two leading EOFs of u for the (a) non-propagating regime and (b) propagating regime in two setups of an idealized GCM. The base latitude is at 30 • S and the contour interval is 0.1. Regions enclosed by contour lines denote values significant at the 95% level. the vertical for 50000-day integrations after spinup. The physics of the model is based on Held and Suarez (1994), an idealized configuration for generating a realistic global circulation with minimal parameterization (Held 2005;Jeevanjee et al. 2017). All diabatic processes are represented by Newtonian relaxation of the temperature field toward a prescribed equilibrium profile, and Rayleigh friction is included in the lower atmosphere to mimic the interactions with the boundary layer.
The non-propagating and propagating regimes are produced in two slightly different setups of this model. For the setup with non-propagating regime, we use the standard configuration of Held and Suarez (1994), which employs an analytical profile approximating a troposphere in unstable radiative-convective equilibrium and an isothermal stratosphere for Newtonian relaxation. For the setup with propagating regime, we follow an approach similar to the one used by Sheshadri and Plumb (2017). In this setup, for the equilibrium temperature profile in the troposphere and stratosphere, we use the perpetual-solstice version of the equilibrium temperature specifications used in Lubis et al. (2018a), calculated from a rapid radiative transfer model (RRTM), with winter conditions in the SH. As will be seen later, these choices result in a large-scale circulation with reasonable annular mode time scales in the SH.
In Fig. 1, we show, following Son and Lee (2006), the one-point lag-correlation maps for the vertically averaged zonal-mean zonal wind anomalies u reconstructed from projections onto the two leading EOFs of u for the two setups (hereafter, angle brackets and overbars denote the vertical and zonal averages, respectively). The anomalies are defined as the deviations from the time mean. The non-propagating and propagating regimes are clearly seen in Figs. 1a and 1b, respectively. In the latter, the propagating anomalies emerge in low latitudes and propagate generally poleward over the course of 3-4 months. In contrast, the non-propagating regime is characterized by persistence zonal flow anomalies in the mid-latitude (Fig. 1a).
To understand the relationship between zonal-mean zonal wind and eddy forcing in the non-propagating and propagating annular modes, the vertically averaged zonalmean zonal wind anomalies ( u ) and vertically averaged zonal-mean eddy momentum flux convergence anomalies ( F ) are projected onto the leading EOFs of u following Lorenz and Hartmann (2001). The time series of zonal index (z) and eddy forcing (m) associated with EOF1 and EOF2 are formulated as: where z 1,2 (m 1,2 ) denotes the component of the field u ( F ) that projects onto the latitudinal structure of the two leading EOFs. u (t) and F (t) are u (φ ,t) and F (φ ,t) with their latitude dimension vectorized, W is a diagonal matrix whose elements are the cos(φ ) weighting used when defining the EOF structure e, and φ is latitude (Simpson et al. 2013;Ma et al. 2017). Here, the vertically averaged zonal-mean eddy momentum flux convergence F is calculated in the spherical coordinate as: where u and v are deviations of zonal wind and meridional wind from their respective zonal means, and a is Earth's radius. Figure 2 shows lagged-correlation analysis between z and m in the GCM setup with non-propagating regime. The auto-correlation of z 1 , as discussed in past studies FIG. 2. Lagged-correlation analysis of the GCM setup with nonpropagating regime. (a) Auto-correlation of z 1 (blue) and z 2 (red), (b) cross-correlation z 1 z 2 , (c) cross-correlation m 1 z 1 , (d) cross-correlation m 2 z 2 , (e) cross-correlation m 1 z 2 , and (f) cross-correlation m 2 z 1 . The two leading EOFs contribute 60.2% and 19.2%, respectively, to the total variance. The e-folding decorrelation time scales of z 1 and z 2 are 64.5 and 4.8 days, respectively. Grey shading represents 5% significance level according to the test of Bartlett (Appendix A).
(e.g., Chen and Plumb 2009;Ma et al. 2017), has a noticeable shoulder at around 5-day lags and shows an unrealistically persistence annular mode, well separated from the faster decaying z 2 , which is consistent with the considerable difference in the contribution of the two EOFs to the total variance (60.2% versus 19.2%). The e-folding decorrelation time scales of z 1 and z 2 are 64.5 and 4.8 days, respectively. The strong, positive cross-correlations of m 1 z 1 and insignificant cross-correlations of m 2 z 2 at large positive lags suggest the existence of a positive eddy-zonal flow feedback for EOF1 (from EOF1) but not for EOF2 (see Son et al. (2008) and Ma et al. (2017)). Figure 2b shows that the z 1 z 2 cross-correlations are weak at positive and negative lags, which consistently with the one-point lag-correlation map of Fig. 1a and Fig. 3 (shown later), are indicative of a non-propagating regime, as reported previously for a similar setup (Son and Lee 2006;Son et al. 2008). The m 1 z 2 and m 2 z 1 cross-correlations are small and often insignificant, suggesting the absence of the cross-EOF feedbacks in the non-propagating regime (Figs. 2ef). All together, the above analysis shows that for the non-propagating regime, single-EOF reduced-order models such as LH01 are sufficient. The weak cross-correlations between z 1 and z 2 in the GCM with non-propagating regime (Fig. 2b) can be also seen by regressing the zonal-mean zonal wind anomalies on the zonal index at 0-and 20-day time lag. Figures 3a and 3b show the wind anomalies regressed on z 1 and z 2 at lag 0, yielding approximately the EOF1 and EOF2 patterns, respectively. Twenty days after z 1 leads zonal wind anomalies, the anomalies do not drift poleward or decay, but rather persist (Fig. 3d). In contrast, 20 days after z 2 leads zonal wind anomalies, the anomalies decay and disappear (Fig. 3c). These observations are consistent with the long and short persistence of z 1 and z 2 , respectively, consistent with the weak cross-correlations of z 1 and z 2 at positive or negative lags, and as become clear below, consistent with the non-propagating nature of this setup. Figure 4 shows lagged-correlation analysis between z and m in the GCM setup with propagating regime. The auto-correlation of z 1 , its persistence compared to that of z 2 , and the explained variance by the two EOFs (40.4% versus 32.5%) are much more similar to what is observed in the SH (shown later in Fig. 7). The e-folding decorrelation time scales of z 1 and z 2 are 14.1 and 9.2 days, respectively. Figure 4b shows that z 1 and z 2 are strongly correlated at long lags peaking at around ±20 days. This behavior along with the one-point lag-correlation map of Fig. 1b and regression map of wind anomalies (Fig. 5, shown later) suggests the existence of a propagating regime, as noted by few previous studies (e.g., Son and Lee 2006;Sheshadri and Plumb 2017). It should be noted that Son and Lee (2006) have proposed a rule of thumb based on the ratio of the explained variance of EOF2 to EOF1: A non-propagating (propagating) regimes exists if the ratio FIG. 4. Lagged-correlation analysis of the GCM setup with propagating regime. (a) Auto-correlation of z 1 (blue) and z 2 (red), (b) crosscorrelation z 1 z 2 , (c) cross-correlation m 1 z 1 , (d) cross-correlation m 2 z 2 , (e) cross-correlation m 1 z 2 , and (f) cross-correlation m 2 z 1 . The two leading EOFs contribute 40.4% and 32.5%, respectively, to the total variance. The e-folding decorrelation time scales of z 1 and z 2 are 14.1 and 9.2 days, respectively. Grey shading represents 5% significance level according to the test of Bartlett (Appendix A). is smaller (larger) than 0.5. The regime of our two setups are consistent with this rule of thumb as the ratios are ∼ 0.3 and ∼ 0.8 in our non-propagating and propagating regimes.
Furthermore, Fig. 4c shows that the m 1 z 1 crosscorrelations are positive at long positive lags (5-20 days) and then negative but small. Fig. 4d indicates small and negative cross-correlations between z 2 and m 2 at the times scale of longer than 20 days (Fig. 4c). Overall, the shape of the m 1 z 1 and m 2 z 2 cross-correlation functions are similar between the non-propagating and propagating regimes, although the m 1 z 1 cross-correlations are larger and more persistent in the non-propagating regime. In contrast, the m 1 z 2 and m 2 z 1 cross-correlations are substantially different between the two regimes (Figs. 4e-f). There are statistically significant and large positive m 1 z 2 crosscorrelations at large positive lags (> 5 days) and statistically significant and large negative m 2 z 1 cross-correlations at positive lags up to 30 days. Note that as emphasized in the figures, positive lags here mean that z 1 (z 2 ) is leading m 2 (m 1 ). Therefore, these cross-correlations, as discussed later, indicate the existence of cross-EOF feedbacks in the propagating regime.  Figure 5 shows anomalous zonal-mean zonal wind regressed on z 1 and z 2 at 0-and 20-day time lag in the GCM setup with propagating regime. Figures 5a and 5b show the wind anomalies regressed on z 1 and z 2 at lag 0, again yielding approximately the EOF1 and EOF2 patterns, respectively. As shown in Fig. 5c, 20 days after z 2 leads zonal wind anomalies, the anomalies have drifted poleward and project strongly onto the structure of wind anomalies associated with EOF1 (Figs. 5a,c, pattern correlation = 0.93). This is consistent with positive correlation of z 1 z 2 at lag +20 days when z 1 leads z 2 (Fig. 4b). Likewise, twenty days after z 1 leads zonal wind anomalies, the anomalies (of Fig. 5a) have drifted poleward and project strongly onto the structure of anomalies associated with EOF2, but with an opposite sign (Figs. 5b,d, pattern correlation = -0.85). This is consistent with negative correlation of z 1 z 2 when z 2 leads z 1 by 20 days (Fig. 4b).
Overall, these results suggest the existence of cross-EOF feedbacks in the propagating annular mode. In Section 3, we will developed a model to quantify these four feedbacks and understand the effects of their magnitude and signs on the variability (e.g., persistence) of z 1 and z 2 . But first, we will examine the variability and characteristics of z and m in reanalysis. In particular, we will see that the z and m cross-correlations in the GCM's propagating regime well resemble those in the SH reanalysis data.

b. Reanalysis
We use the 1979-2013 data from the European Centre for Medium-Range Weather Forecasts (ECMWF) interim reanalysis (ERA-Interim; Dee et al. 2011). Zonal and meridional wind components (u, v) are 6 hourly, on 1.5 • FIG. 6. One-point lag-correlation maps of the vertically averaged zonal-mean zonal wind anomalies from year-round ERA-Interim data integrated across the depth of the troposphere (1000-100 hPa) ( u ) in the Southern Hemisphere. (a) Total anomaly fields and (b) reconstructed from projections onto the two leading EOFs of u . The base latitude is at 30 • S and the contour interval is 0.1. Regions enclosed by contour lines denote values significant at the 95% level according to the t-test. latitude × 1.5 • longitude grid, and on 21 vertical levels between 1000 and 100 hPa. Anomalies used for computing correlations and EOF analyses are defined as the deviations from the climatological seasonal cycle. The mean seasonal cycle is defined as the annual average and the first four Fourier harmonics of the 35-yr daily climatology. Figure 6 shows a one-point lag-correlation map of vertically averaged zonal-mean zonal wind u in the SH, where the base latitude is 30 • S. Comparing this figure with Fig. 1, it can be seen that there is an indication of poleward-propagating anomalies in SH, which appear in low latitudes and migrate poleward over the course of 2-3 months (Fig. 6a). However, the poleward-propagating signals are not as clearly as those observed in the GCM setup with the propagating regime (Fig. 1b, or Fig. 2 of Son and Lee 2006). This is consistent with previous studies (e.g. Feldstein 1998;Feldstein and Lee 1998;Sheshadri and Plumb 2017), showing that both propagating and non-propagating anomalies exist in all seasons in the SH, which somehow obscure the propagating signals. Reconstructions based on the projections onto the two lead- FIG. 7. Lagged-correlation analysis for the Southern Hemisphere, calculated from year-round ERA-Interim data. (a) Auto-correlations of z 1 (blue) and z 2 (red), (b) cross-correlation z 1 z 2 , (c) cross-correlation m 1 z 1 , (d) cross-correlation m 2 z 2 , (e) cross-correlation m 1 z 2 , and (f) cross-correlation m 2 z 1 at different lags. The two leading EOFs contribute to 45.1% and 23.2% of the total variance, respectively. The efolding decorrelation time scales of z 1 and z 2 are 10.3 and 8.1 days, respectively. Grey shading represents 5% significance level according to the test of Bartlett (Appendix A).
ing EOFs of zonal-mean zonal wind further show that most of the mid-latitude SH wind variability can be explained by the two leading EOF modes (Fig. 6b). The ratio of the fractional variance of EOF2 (23.2%) to that of EOF1 (45.1%) is 0.51, which is right at the boundary from the rule of thumb. Overall, as already pointed out by Sheshadri and Plumb (2017), a propagating annular mode exists in the SH and is largely explained by the two leading EOF modes. Figure 7a shows the auto-correlations of z 1 and z 2 . Consistent with Lorenz and Hartmann (2001), the estimated decorrelation time scales of these two PCs are 10.3 and 8.1 days, respectively. Figure 7b depicts the cross-correlation z 1 z 2 , showing statistically significant and relatively strong correlations that peak around ±10 days. As discussed in earlier studies, such lagged correlations are a signature of the propagating annular modes (Feldstein and Lee 1998;Son and Lee 2006;Son et al. 2008;Sheshadri and Plumb 2017), implying that the period of the poleward propagation is about 20-30 days in the SH (Fig. 7b), consistent with Sheshadri and Plumb (2017) and with Fig. 6.
To understand the effects of z 1 and z 2 on m 1 and m 2 , we also examine the cross-correlations between z and m at different lags (Figs. 7c-f). The shape and the magnitude of the m 1 z 1 and m 2 z 2 cross-correlations (Figs. 7c-d) are similar to those originally shown by Lorenz and Hartmann (2001) (see their Figs. 5 and 13a) and later by many others using different reanalysis products and time periods. As discussed in Lorenz and Hartmann (2001), the statistically significant positive m 1 z 1 cross-correlations at long positive lags (∼ 8 − 20 days) and the insignificant m 2 z 2 cross-correlations for time scales longer than ∼ 5 days are indicative that a positive eddy-zonal flow feedback exists only for EOF1, but not for EOF2 (also see Byrne et al. (2016) and Ma et al. (2017)). We emphasize that this positive feedback is from EOF1 onto itself.
To see if there are cross-EOF feedbacks, in Figs. 7ef we plot the m 1 z 2 and m 2 z 1 cross-correlations at different lags. The m 1 z 2 cross-correlations show statistical significant positive correlations at large positive lags, signifying that a cross-EOF feedback, i.e, z 2 modifying m 1 , is present. Note that the magnitude of the m 1 z 2 crosscorrelations at positive lags is overall larger than those of m 1 z 1 (Fig. 7c). There are also statistically significant but negative m 2 z 1 correlations at large positive lags, again suggesting the existence of a cross-EOF feedback, i.e, z 1 modifying m 2 . These results indicate that in the presence of propagating regime in the SH, there are indeed cross-EOF feedbacks; however, these feedbacks were always ignored in the previous studies and reduced-order models of the SH extratropical large-scale circulation.

Eddy-zonal flow feedbacks in the propagating annular modes: Model and quantification
In this section, an eddy-zonal flow feedback model that accounts for the coupling of the leading two EOFs and their feedbacks, including the cross-EOF feedbacks will be introduced. Then this model will be validated using synthetic data from a simple stochastic prototype, and from its analytical solution, we will derive conditions for the existence of the propagating regime. Finally, we will use this model to estimate the feedback strengths of the propagating annular modes in data from the reanalysis (SH) and the idealized GCM.
a. Developing an eddy-zonal flow feedback model for propagating annular modes With the same notations as in Lorenz and Hartmann (2001), the time series of zonal indices (z 1 and z 2 ) and eddy forcing (m 1 and m 2 ) associated with the first two leading EOFs are calculated by projecting the vertically averaged zonal-mean zonal wind u and eddy momentum flux convergence F anomalies onto the patterns of the first and second EOFs of u (see Eqs. (1)-(2)). Equations for the tendency of z 1 and z 2 can be then formulated as: where t is time and the last term on the right-hand side in each equation represents damping (mainly due to surface friction) with time scale τ. As discussed in Lorenz and Hartmann (2001), Eqs. (4)-(5) can be interpreted as the zonally and vertically averaged zonal momentum equation: projected into EOF1 and EOF2, respectively. In the above equation, D includes the effects of surface drag and is modeled as Rayleigh drag in Eqs. (4)-(5). Assuming a linear representation for the feedback of an EOF onto itself, Lorenz and Hartmann (2001) and later studies wrote m 1 (t) =m 1 (t)+b 1 z 1 (t) and m 2 (t) =m 2 (t)+ b 2 z 2 (t), where b 1 and b 2 are the feedback strengths (with b j > 0 implying a positive feedback that prolongs the persistence of z j ).m is the random, zonal flow-independent component of the eddy forcing that drives the highfrequency variability of z (Lorenz and Hartmann 2001;Ma et al. 2017).
Here, to account for the cross-EOF feedbacks, i.e., the effect of z 2 on m 1 and z 1 on m 2 , we extend the LH01 model and write With j, k = 1, 2, b jk is the strength of the linearized feedback of z k onto z j through modifying m j in the quasisteady limit; thus the cross-EOF feedbacks are represented by the terms involving b 12 and b 21 . To find the values of b jk , we can use the lagged-regression method of Simpson et al. (2013), which assumes that reg l (m j , z j ) = sum(m j (t + l)z j (t)) ≈ 0 at large positive lags, l. By lagregressing each term in Eqs. (7) onto z 1 and then onto z 2 , we find and similarly, from Eq. (8) we find where we assumed reg l (m j , z k ) ≈ 0 for j, k = 1, 2. TABLE 1. Prescribed and estimated feedback strengths (in day −1 ) in synthetic data for the case without cross EOF-feedbacks. The imposed damping rates of friction are τ 1 =τ 2 = 8 days. The values of b and τ are motivated by the observed ones, see Table 4.
Prescribed 0.040 0.000 0.000 0.000 Estimated (Eqs. (9)-(10)) 0.042 0.001 -0.0006 0.0005 Note that if one attempts to find b 11 using a single-EOF approach such as LH01, then, from Eq. (7), one would be implicitly assuming that reg l (m 1 + b 12 z 2 , z 1 ) = reg l (m 1 , z 1 ) + b 12 reg l (z 2 , z 1 ) ≈ b 12 reg l (z 2 , z 1 ) is zero. However, as shown earlier, in the propagating regime, the z 1 z 2 cross-correlations can be large at long lags, and as discussed below, the range of time lags needed to be used in Eqs. (9) 2017)), i.e., reg l (m j , z k ) ≈ 0 at "largeenough" positive lags. Note that one cannot use a lag that is too long because then even reg l (z j , z j ) would be small and inaccurate. To find the appropriate lag to use, one must look for non-zero m j z k cross-correlations at positive lags beyond an eddy lifetime. In this study, the strengths of the individual feedbacks are averaged over positive lags of 8 to 20 days for both GCM and reanalysis (e.g., Simpson et al. 2013;Burrows et al. 2016). We choose this range in order to avoid the high-frequency variability at short lags (indicated by impulsive and oscillatory characters of thẽ m auto-correlation) and strong damping at the very long lags.
In the following section, we will present a proof of concept for this eddy-zonal flow feedback model using synthetic data obtained from a simple stochastic prototype and show that using Eqs. (9)-(10), the prescribed feedbacks can be accurately backed out.  −1 ) in synthetic data for the case with cross EOF-feedbacks. The imposed damping rates of friction are τ 1 =τ 2 = 8 days. The values of b and τ are motivated by the observed ones, see Table 4.
Prescribed 0.040 0.060 -0.025 0.000 Estimated (Eqs. (9)-(10)) 0.043 0.067 -0.026 -0.002 b. Validation using synthetic data from a simple stochastic prototype We begin by constructing a simple stochastic system to produce synthetic time series z and m in the presence or absence of cross-EOF feedbacks. The equations of this system are the same as Eqs. (4)-(5) and (7)-(8). Following Simpson et al. (2013), we generate a synthetic time series of the random component of the eddy forcing m 1,2 using a second-order autoregressive (AR2) noise process: where t denotes time (in days) and ε is white noise distributed uniformly between -1 and +1. Synthetic time series of z j and m j are produced by numerically integrating Eqs. (4)-(5), (7)-(8), and (11)-(12) forward in time with two different sets of prescribed b jk . In the first set, there is no cross-EOF feedback, i.e., b 12 = b 21 = 0 ( Table 1). In the second set, b 11 and b 22 = 0 are the same as the first set, but here there is cross-EOF feedback, i.e., b 12 and b 21 = 0 ( Table 2). For both sets, we assumed τ 1 = τ 2 = 8 days. The values of b and τ are reasonably chosen based on the observed values in the SH (see Table  4).
Spectral analysis of z 1,2 and m 1,2 shows that the synthetic data indeed have characteristics similar to those of the observed SH. For example, for the case with crossfeedbacks (Fig. 8), we find that consistent with observations (see Fig. 4 of Lorenz and Hartmann (2001) or Fig. 3 of Ma et al. (2017)), the time scales of z 1 and z 2 are much longer (i.e., slower variability) than m 1 and m 2 , and the power spectra of z can be interpreted, to the first order, as reddening of the power spectra of eddy forcing m Lorenz and Hartmann (2001); Ma et al. (2017). The power spectra of eddy forcings m 1 and m 2 have in general a broad maximum centered at the low and synoptic frequency, consistent with observations. Given that the characteristics of the synthetic data mimic the key characteristics of the observed annular modes, we use this idealized framework to FIG. 8. Spectra of z 1,2 and m 1,2 from the synthetic data with cross-EOF feedbacks. Black lines show the power spectra of (a) z 1 , (b) z 2 , (c) m 1 , and (d) m 2 . The red-noise spectra are indicated by the smooth solid red curves, and the smooth dashed blue lines are the 5% and 95% a priori confidence limits.
validate the lagged-correlation approach of Eqs. (9)-(10) for quantifying eddy-zonal flow feedbacks. Figure 9 shows the lagged-correlation analysis of the synthetic data without cross-EOF feedbacks. It is clearly seen that the only noticeable cross-correlations are that of m 1 z 1 , and there are no (statistically significant) crosscorrelations between z 1 z 2 , m 1 z 2 and m 2 z 1 at any lag, consistent with a non-propagating regime and the absence of cross-EOF feedbacks (Fig. 2). Using Eqs. (9)-(10) and lag l=8-20 days, we can closely estimate the prescribed feedback parameters, i.e., b 11 = 0.04 day −1 and b 22 = b 12 = b 21 = 0 (see Table 1). Figure 10 shows the lagged-correlation analysis of the synthetic data with cross-EOF feedbacks. First, we see that there are statistically significant and often large crosscorrelations in z 1 z 2 , m 1 z 1 , m 1 z 2 , and m 2 z 1 , with the shape of the cross-correlation distributions not that different from that of the SH reanalysis and the idealized GCM setup with propagating regime (Figs. 4 and 7). The positive m 1 z 1 and near zero m 2 z 2 cross-correlations at large positive lags signify a positive z 1 -onto-z 1 feedback through m 1 , but no z 2 -onto-z 2 feedback through m 2 , consistent with the prescribed positive value of b 11 and b 22 = 0. In addition, Figs. 10e-f also show that there are statistically significant and large correlations in m 1 z 2 and m 2 z 1 at positive lags, consistent with the introduction of cross-EOF feedbacks by setting b 12 = 0.06 day −1 and b 21 = −0.025 day −1 . The positive m 1 z 2 cross-correlations are positive lags are higher than those of m 1 z 1 (note that FIG. 9. Lagged-correlation analysis of synthetic data without cross-EOF feedback. (a) Auto-correlation of z 1 (blue) and z 2 (red), (b) crosscorrelation z 1 z 2 , (c) cross-correlation m 1 z 1 , (d) cross-correlation m 2 z 2 , (e) cross-correlation m 1 z 2 , and (f) cross-correlation m 2 z 1 . The e-folding decorrelation time scales of z 1 and z 2 are 18.6 and 9.2 days, respectively. Grey shading represents 5% significance level according to the test of Bartlett (Appendix A). b 12 /b 11 ≈ 1.5), and the sign of m 2 z 1 cross-correlations is opposite to the sign of m 1 z 2 cross-correlations (note that b 12 b 21 < 0). Using Eqs. (9)-(10) and lag l=8-20 days, we can again closely estimate the prescribed feedback parameters, including the strength of the cross-EOF feedbacks (see Table 2).
The above analyses validate the approach using Eqs. (9)-(10) for quantifying the feedback strengths b jk in data from both propagating and non-propagating regimes. Furthermore, a closer examination of z 1 and z 2 autocorrelations in Figs. 9a and 10a show that both z 1 and z 2 in the case without cross-EOF feedbacks are more persistence than those in the case with cross-EOF feedbacks; e.g., the e-forcing deccorelation time scale of z 1 is 18.6 days in Fig. 9a while it is 13.9 days in Fig. 10a. This observation might be counter-intuitive because both cases have the same b 11 > 0 while the case with cross-EOFs feedback has b 12 > 0, which might seem like another positive feedback that should further prolong the persistence of z 1 . Finally, we notice that b 12 b 21 < 0 in Table 2 and in the SH reanalysis and idealized GCM setup with the propagating regime (Tables 4 and 5). Synthetic data generated with the same parameters as in Table 2 but with the sign of b 21 flipped results in cross-correlation distributions FIG. 10. Lagged-correlation analysis of synthetic data with cross-EOF feedback. (a) Auto-correlation of z 1 (blue) and z 2 (red), (b) crosscorrelation z 1 z 2 , (c) cross-correlation m 1 z 1 , (d) cross-correlation m 2 z 2 , (e) cross-correlation m 1 z 2 , and (f) cross-correlation m 2 z 1 . The e-folding decorrelation time scales of z 1 and z 2 are 13.9 and 6.5 days, respectively. The regions outside the gray shading indicate 95% significance level according to the test of Bartlett (Appendix A). that are vastly different from those of Fig. 10 and what is seen in the SH reanalysis and idealized GCM. Inspired by these observations, next we examine the analytical solution of the deterministic version of Eqs. (4)-(5 and (7)-(8) to better understand the impacts of the strength and sign of b jk on the variability and in particular the persistence of z 1 and z 2 .
c. Analytical solution of the two-EOF eddy-zonal flow feedback model We focus on the deterministic (i.e.,m j = 0) version of Eqs. (4)-(5) and (7)-(8), which can be re-written as the following system of ordinary differential equations (ODEs): where z = z 1 z 2 , and A = b 11 − 1 The solution to this system is where V and Λ are the eigenvector and eigenvalue matrices of A: To find the eigenvalues λ , we set the determinant of A equal to zero and solve the resulting quadratic equation to obtain: which, in the limit of τ 1 ≈ τ 2 (reasonable given their estimated values in Tables 4 and 5), simplifies to: (18) The solution (Eq. (15)) can be re-written as where c 1 and c 2 depend on the initial condition. This system has a decaying-oscillatory solution, i.e., is in the propagating regime, if and only if the eigenvalues (18) have non-zero imaginary parts, which requires, as a necessary and sufficient condition: Equation (20) also implies that a necessary condition for the existence of propagating regimes is Thus, non-zero cross-EOF feedbacks of opposite signs are essential components of the propagating regime dynamics. The propagating regimes in the stochastic prototype (Table 2), SH reanalysis (Table 4), and idealized GCM (Table 5) satisfy the conditions of Eqs. (20)-(21), while the non-propagating regimes (Tables 1 and 5) do not.
In the non-propagating regime, λ 1,2 = −σ 1,2 < 0 and v 1,2 are real and in this regime, z 1,2 just decay exponentially according to In the propagating regime, λ 1,2 = −σ ± iω and v 1,2 are complex where In this regime, z 1,2 decay and oscillate according to Realizing that in this case v 11 = v 12 are real, and v 21 = v * 22 and c 1 = c * 2 = c, where * means complex conjugate, we can re-write the above equations as σt) . (27) These equations show that z 1 and z 2 have the same decay rate (σ ) but different oscillatory components with frequency ω. These results are consistent with the POP analysis of Sheshadri and Plumb (2017) who showed that EOF1 and EOF2 are, respectively, the real and imaginary parts of a single decaying-oscillatory POP mode (see their Section 4b). As a results, the two modes have the same decay rate and frequency, but have different auto-correlation function decay rates and have strong lag cross-correlations because the oscillations are out of phase. A key contribution of our work is to find the decay rate σ and frequency ω as a function of b jk and τ j (Eqs. (23)- (24)).
To understand the effects of the feedback strength b jk on the persistence of z j , we compute the analytical solutions for 5 systems that have the same b 11 > 0 and b 22 = 0 (Table 3): In EXP1, there is no cross-EOF feedback (b 12 = b 21 = 0), while in EXP2-EXP5, b 12 > 0 and b 21 < 0 and they have been doubled from experiment to experiment. Figure 11 shows the auto-correlation coefficients of z 1 and their e-folding decorrelation time scales for EXP1-EXP5. EXP1, corresponding to non-propagating regimes, has the slowest-decaying auto-correlation function, i.e., longest e-folding decorrelation time scale (Figs. 11a,b). EXP2-EXP5, which all satisfy condition (Eq. 20), have fasterdecaying auto-correlation functions, i.e., shorter e-folding decorrelation time scale, consistent with our earlier results in idealized GCM and stochastic prototype (Figs. 4 and 10). As discussed above, in the propagating regime, the eigenvectors and the corresponding eigenvalue are complex and thus, z 1,2 do not decay just exponentially, but rather show some oscillatory characteristics too (Fig. 11a, Eqs. (26)- (27)). Since the frequency of these oscillations ω (Eq. (24)) increases as the cross-EOF feedback strengths increase, shorter time scales in z 1 are expected in the experiment with stronger b 12 b 21 (Fig. 11b).
The dependence of the e-folding decorrelation time scales of z 1 and z 2 on the feedback strengths, and in particular the cross-EOF feedback strengths, is further evaluated in Fig. 12. In Fig. 12a, it is clearly seen that the impact of increasing b 11 > 0 in the propagating regime (filled symbols) is to increase the persistence, i.e., decorrelation time scale, of z 1 (Fig. 12a), consistent with increasing the positive eddy-zonal flow feedback (z 1 -ontoz 1 through m 1 ). However, when the feedback is further  (20) for the existence of a decaying-oscillatory solution is not satisfied anymore, and consistent with this, we see that the system undergoes a transition to the non-propagating regime. Further increasing b 11 leads to substantially more persistent z 1 and less persistence z 2 . Note that in non-propagating regimes when b 12 b 21 = 0, the decay of z 2 depends on b 11 too (see Eq. (18)). Figure 12b shows that in the propagating regime, unlike increasing b 11 > 0, increasing b 12 > 0 leads to reduction in the persistence of z 1 . This is the counter-intuitive behavior we had observed earlier in the stochastic prototype (Section 3b). Now we understand that this is because increasing b 12 increases the frequency of the oscillation ω in the system, resulting in reduction in the the decorrelation time scale of z 1 (and z 2 ); also see Fig. 11. Such impact can even be more pronounced when both cross-EOF feedbacks b 12 and b 21 are increased (Fig. 12c), leading to shorter decorrelation time scales. Because a positive b 12 decreases the persistence of z 1 , we do not refer to is as a "positive feedback". To understand this behavior, we have to keep in mind that in the eddy forcing of z 1 (z 2 ), i.e., m 1 in Eq. (7) (m 2 in Eq. (8)), b 12 > 0 (b 21 < 0) is the coefficient of z 2 (z 1 ). When z 2 leads z 1 , they are negatively correlated (Figs. 4b,7b, and 10b), thus z 2 multiplied by b 12 > 0 reduces m 1 that is forcing z 1 , decreasing the persistence of z 1 . Similarly, when z 1 leads z 2 , they are positively correlated, thus z 1 multiplied by b 21 < 0 reduces m 2 and thus the persistence of z 2 .
Finally, for the sake of completeness, we also examine the effect of increasing b 11 in the absence of cross-EOF feedback (Fig. 12d). As expected increasing b 11 leads to increasing the persistence of z 1 and has no impact on the persistence of z 2 as now z 1 and z 2 are completely decoupled.

d. Quantifying eddy-zonal flow feedbacks in reanalysis and idealized GCM
The results of Sections 3b and 3c show the importance of carefully quantifying and interpreting the eddy-zonal flow feedbacks, including the cross-EOF feedbacks, to understand the variability of the zonal-mean flow.  Table 3.  (Table 3). The filled indicates that the parameters satisfy the condition for propagating regimes, i.e., existence of decaying-oscillatory solutions (Eq. (20)). Table 4 presents the feedback strengths obtained from applying (9)-(10) with l = 8 − 20 days to the year-round SH reanalysis data. We find b 11 = 0.038 day −1 , a positive feedback from z 1 onto z 1 , consistent with the findings of Lorenz and Hartmann (2001) in their pioneering work. This estimate of b 11 is slightly higher than what we find using the single-EOF approach (b 11 = 0.035 day −1 ), which is the same as what Lorenz and Hartmann (2001) found using their spectral cross-correlation method. We also find non-zero cross-EOF feedbacks: b 12 = 0.059 day −1 and b 21 = −0.020 day −1 . We also estimate b 22 = 0.017 day −1 that is slightly higher from what the single-EOF approach yields (Table 4). The estimated feedback strengths and friction rates (τ) in Table 4 satisfy the condition for propagating regime (Eq. 20). It should be noted that we also extended our approach to include the leading 3 EOFs and quantified the 9 feedback strengths; however, we found the effects of EOF3 on EOF1 and TABLE 4. Feedback strengths (in day −1 ) estimated for year-round ERA-Interim reanalysis. The damping rates of friction are estimated as τ 1 = 8.3 days and τ 2 = 8.4 days following the methodology in Appendix A of Lorenz and Hartmann (2001).
Eqs. (9)-(10) 0.038 0.059 -0.020 0.017 LH01 0.035 --0.002 EOF2 negligible, which suggests that a two-EOF model (9)-(10) is enough to describe the current SH large-scale circulation (not shown). Table 5 presents the feedback strengths obtained from applying (9)-(10) with l = 8 − 20 days to the two setups of the idealized GCM. In the non-propagating regime, we find b 11 = 0.133 day −1 , and small b 22 and negligible b 12 and b 21 , indicating the absence of cross-EOF feedbacks, consistent with insignificant m 1 z 2 and m 2 z 1 crosscorrelations (Figs. 2e-f). The values of b jk do not satisfy the condition for propagating regime, which is consistent with weak cross-correlation between z 1 and z 2 at long lags (Fig. 2b). These results suggest that a strong z 1 -onto-z 1 feedback dominates the dynamics of the annular mode in this setup (the standard Held-Suarez configuration), which leads to an unrealistically persistent annular mode, similar to what is seen in Fig. 12d, and consistent with the findings of previous studies (Son and Lee 2006;Son et al. 2008;Ma et al. 2017). Using the linear response function (LRF) of this setup Kuang 2016b, 2019) showed that this eddy-zonal flow feedback is due to enhanced lowlevel baroclinicity (as proposed by Robinson (2000) and Lorenz and Hartmann (2001)) and estimated, from a budget analysis, that the positive feedback is increasing the persistence of the annular mode by a factor of two.
In the propagating regime, we find b 11 = 0.101 day −1 , which is slightly lower than b 11 of the non-propagating regime. However, in the propagating regime, we also find strong cross-EOF feedbacks b 12 = 0.075 day −1 , b 21 = −0.043 day −1 as well as b 22 = 0.023 day −1 . These feedback strengths satisfy the condition for propagating regime, consistent with strong cross-correlation between z 1 and z 2 at long lags (Fig. 4b). Comparing the two rows of Table 5 and Figs. 2a and 4a with Table 4 and Fig. 7a suggests that while it is true that the b 11 of the the idealized GCM's non-propagating regime is larger than that of the SH reanalysis (by a factor of 3.5), the unrealistic persistence of z 1 in this setup (time scale ≈ 65 days) compared to that of the reanalysis (time scale ≈ 10 days; compare Figs. 2a and 7a) could be, at least partially, due to the absence of cross-EOF feedbacks (thus oscillations), which as we showed earlier in Section 3c, reduce the persistence of the annular modes. The GCM setup with propagating regime has b 11 that is around 2.7 times larger than that of TABLE 5. Feedback strengths (in day −1 ) estimated for the idealized GCM setups with non-propagating and propagating regimes. The estimated damping rates of friction are τ 1 =7.4 days and τ 2 =7.6 days for the GCM setup with non-propagating regime, and τ 1 =7.1 days and τ 2 =7.4 days for the GCM setup with propagating regime (estimated using the methodology in Appendix A of Lorenz and Hartmann (2001) These findings show the importance of quantifying and examining cross-EOF feedbacks to fully understand the dynamics and variability of the annular modes and to better evaluate how well the GCMs simulate the extratopical large-scale circulation.

Concluding remarks
The low-frequency variability of the extra-tropical large-scale circulation is often studied using a reducedorder model of the leading EOF of zonal-mean zonal wind. The key component of this model (LH01) is an internal-to-troposphere eddy-zonal flow interaction mechanism which leads to a positive feedback of EOF1 onto itself, thus increasing the persistence of the annular mode (Lorenz and Hartmann 2001). However, several studies have showed that under some circumstances, strong couplings exist between EOF1 and EOF2 at some lag times, resulting in decaying-oscillatory, or propagating, annular modes (e.g. Son and Lee 2006;Son et al. 2008;Sheshadri and Plumb 2017). In the current study, following the methodology of Lorenz and Hartmann (2001) and using data from the SH reanalysis and two setups of an idealized GCM that produce circulations with a dominant non-propagating or propagating regime, we first show strong cross-correlations between EOF1 (EOF2) and the eddy forcing of EOF2 (EOF1) at long lags, suggesting that cross-EOF feedbacks might exist in the propagating regimes. These findings together demonstrate that there is a need to extend the single-EOF model of LH01 and build a model that includes, at a minimum, both leading EOFs and accounts for their cross feedbacks.
With similar assumptions and simplifications used in Lorenz and Hartmann (2001), we have developed a two-EOF model for propagating annular modes (consisting of a system of two coupled ODEs, Eqs. (4)-(5) with (7)-(8)) that can account for the cross-EOF feedbacks. In this model, the strength of the feedback of kth EOF onto jth EOF is b jk ( j, k = 1, 2). Using the analytical solution of this model, we derive conditions for the existence of the propagating regime based on the feedback strengths. It is shown that the propagating regime, which requires a decaying-oscillatory solution of the coupled ODEs, can exist only if the cross-EOF feedbacks have opposite signs (b 12 b 21 < 0), and if and only if the following criterion is satisfied: (b 11 − b 22 ) 2 < −4b 12 b 21 . These criteria show that non-zero cross-EOF feedbacks are essential components of the propagating regime dynamics.
Using this model and the idealized GCM and a stochastic prototype, we further show that cross-EOF feedbacks play an important role in controlling the persistence of the propagating annular modes (i.e., the e-folding decorrelation time scale of the zonal index, z j ) by setting the frequency of the oscillation ω (Eq. (24)). Therefore, in this regime, the persistence of the annular mode (EOF1) does not only depend on the feedback of EOF1 onto itself, but also depends on the cross-EOF feedbacks. We find that as a result of the oscillation, the stronger the cross-EOF feedbacks, the less persistent the annular mode.
Applying the coupled-EOF model to the reanalysis data shows the existence of strong cross-EOF feedbacks in the current SH extratropical large-scale circulation. Annular modes have been found to be too persistent compared to observations in GCMs including IPCC AR4 and CMIP5 models (Gerber and Vallis 2007;Gerber et al. 2008a;Bracegirdle et al. 2020). This long persistence has been often attributed to a too strong positive EOF1-onto-EOF1 feedback in the GCMs. The dynamics and strength of this feedback depends on factors such as the mean flow and surface friction (Robinson 2000;Lorenz and Hartmann 2001;Chen and Plumb 2009;Hassanzadeh and Kuang 2019). External (to troposphere) influence, e.g., from the stratospheric polar vortex, has been also suggested to affect the persistence of the annular modes (Byrne et al. 2016;Saggioro and Shepherd 2019). Our results show that the cross-EOF feedbacks play an important role in the dynamics of the annular modes, and in particular, that their absence or weak amplitudes can increase the persistence, offering another explanation for the too-persistent annular modes in GCMs.
Overall, our findings demonstrate that to fully understand the dynamics of the large-scale extratropical circulation and the reason(s) behind the too-persistent annular modes in GCMs, the coupling of the leading EOFs and the cross-EOF feedbacks should be examine using models such as the one introduced in this study.
An important next step is to investigate the underlying dynamics of the cross-EOF feedbacks. So far we have pointed out that cross-EOF feedbacks are essential components of the propagating annular modes; however, the propagation itself is likely essential for the existence of cross-EOF feedbacks. In fact, our preliminary result shows that the cross-EOF feedbacks result from the outof-phase oscillations of EOF1 (north-south jet displacement) and EOF2 (jet pulsation) leading to an orchestrated combination of equatorward propagation of wave activity (a baroclinic process) and nonlinear wave breaking (a barotropic process), which altogether act to reduce the total eddy forcings (not shown). In ongoing work, we aim to explain and quantify the propagating annular modes dynamics using the LRF framework of Hassanzadeh and Kuang (2016a,b) and finite-amplitude waveactivity framework (Nakamura and Zhu 2010;Lubis et al. 2018a,b) that have been proven useful in understanding the dynamics of the non-propagating annular modes (Nie et al. 2014;Ma et al. 2017;Hassanzadeh and Kuang 2019).