Large-scale controls of propagation of the Madden-Julian Oscillation

With widespread influence on global climate and weather extremes, the Madden-Julian Oscillation (MJO) plays a crucial role in subseasonal prediction. Our latest global climate models (GCMs), however, have great difficulty in realistically simulating the MJO. This model inability is largely due to problems in representation of MJO’s cumulus organization. This study, based on a series of idealized aqua-planet model experiments using an atmospheric-only GCM, clearly demonstrates that MJO propagation is strongly modulated by the large-scale background state in which the lower-tropospheric mean moisture gradient and zonal winds are critical. Therefore, when tuning climate models to achieve improved MJO simulations, particular attention needs to be placed on the model large-scale mean state that is also significantly affected by cumulus parameterizations. This study indicates that model biases in representing MJO propagation may be related to the widely reported double-ITCZ (intertropical convergence zone) problem in climate models.


INTRODUCTION
First detected in the 1970s and named after its two discoverers, the now well-known Madden-Julian Oscillation (MJO) 1 is characterized by slow eastward propagating large-scale convective fluctuations along the equator with characteristic periods of 30-60 days, and has been recognized to have tremendous influence on global weather extremes 2 . Due to the MJO's unique role in bridging weather and climate, predictability of the MJO on the intraseasonal time scale enables useful prediction of extreme weather activity beyond the deterministic forecast limit of about 1 week 3 .
Modeling the MJO, however, remains a grand challenge for the climate research community [4][5][6] . Until most recently, the observed eastward propagation of MJO convection and associated circulations could only be simulated by a limited number of global climate models (GCMs). For example, as shown in Fig. 1, the ECHAM atmosphere-only GCM (AGCM) produces weak propagation of MJO convection over the Indian Ocean when forced by observed climatological sea surface temperatures (SSTs). This is in contrast to the observed systematic eastward propagation of MJO convection from the Indian Ocean to western Pacific. Model deficiencies in representing the MJO are generally ascribed to deficiencies in depicting cumulus processes 7 , although inclusion of interactive ocean feedbacks can often improve MJO simulations 8 . A common practice to improve MJO simulations in GCMs is to inhibit triggering of deep convection in cumulus schemes, for example, through enhanced cumulus entrainment rates or increased rain re-evaporation 7 . However, the improved MJO representation achieved by tuning using such methods often occurs at the cost of a degraded model mean state and other climate phenomena 9 .
Meanwhile, a lack of consensus exists on the fundamental physics of the MJO 6,10 . One traditional body of thought considers the MJO to be a couplet of the equatorial Kelvin and Rossby waves. The east-west asymmetry in wave response to MJO convective heating, i.e., the Kelvin wave to the east of MJO convection and Rossby wave to the west leads to eastward propagation of the MJO 11,12 . Under this framework, MJO phase speed is determined by relative intensity of equatorial low-level easterly anomalies that are largely associated with the Kelvin wave and the westerly anomalies that are associated with the Rossby wave 12,13 . A strong Kelvin wave component promotes strong boundary layer convergence and the formation of shallow cumuli to the east of MJO convection, facilitating moisture buildup through the frictional CISK (conditional instability of the second kind) mechanism that supports eastward propagation 11,14 . The inability of many GCMs to simulate MJO eastward propagation under this view has been ascribed to Kelvin wave circulations to the east of MJO convection that are too weak 15 .
In another view where the MJO is considered to be a moisture mode such that MJO convective activity is largely regulated by moisture perturbations [16][17][18] , recent observational and modeling studies illustrate that the propagation of MJO convection is closely associated with lower-tropospheric moistening or drying processes dominated by horizontal moisture advection, with horizontal advection of background moisture by the MJO circulation being the leading term [19][20][21][22] . These findings suggest a crucial role for the large-scale environment for regulating the MJO propagation. This notion has been supported by recent multi-model analyses, which suggest that model performance in representing MJO propagation is closely linked to skill in representing the mean moisture pattern 20,[23][24][25] . Further, fluctuations in MJO characteristics in current climate (e.g., seasonal and interannual variations), and future climate projections, have been linked to changes in the corresponding background state 19,21,26,27 . Therefore, in addition to properly representing model convective organization as mentioned above, a realistic model basic state is necessary to achieve skillful representation of the MJO in GCMs.
In this study, a series of idealized aqua-planet experiments using the ECHAM AGCM that was used for simulations in Fig. 1b (see "Methods") are used to demonstrate how relatively small changes in SST patterns can lead to dramatic changes in model MJO propagation characteristics, lending further support to the crucial role of the large-scale background state for regulating MJO propagation.

RESULTS
Distinct MJO propagation in responding to SST patterns Forced by various zonally uniform SST distributions with their meridional profiles transitioning from QOBS to FLAT (see Fig. 2a and Methods), a series of idealized experiments based on the ECHAM AGCM are conducted under an aqua-planet configuration. These SST profiles feature a gradual reduction of the SST gradient near the equator following a previous study 28 (see "Methods"). Consistent with previous work 29 , these small changes in SST gradients near the equator produce dramatic changes in model mean states, including precipitation, moisture, and zonal winds ( Fig. 2b-f). In the QOBS experiment, the mean precipitation pattern is characterized by a single intertropical convergence zone (ITCZ; Fig. 2b) and with associated moisture maximum near the equator at all heights in the troposphere (Fig. 2c). Meanwhile, mean zonal winds in QOBS exhibits two off-equatorial easterly maxima below 600 hPa with strong extratropical westerlies in the mid-upper troposphere in both hemispheres (Fig. 2d). In contrast, in response to a weaker equatorial meridional SST gradient in FLAT, the mean precipitation exhibits a double-ITCZ pattern with the maximum rainfall situated near 13°in both hemispheres (Fig. 2b). Correspondingly, two off-equatorial maxima in mean moisture occur in FLAT (Fig. 2e), in contrast to the equatorial moisture peak in QOBS. Additionally, weak equatorial mean westerlies in the mid-upper troposphere between 10°S and 10°N in QOBS are replaced by equatorial easterlies throughout the entire troposphere in FLAT. Note that for a given SST profile, the corresponding model mean state can differ from model to model due to their different sensitivities of model convective mixing and cloud-radiative feedbacks to the SST distribution [29][30][31][32] .
Accompanying the changes in model mean states associated with different SST profiles, the intraseasonal variability of tropical convection exhibits distinct propagation behavior in QOBS and  FLAT (Fig. 3). In QOBS, tropical convection is characterized by frequent circumnavigating eastward propagation with a phase speed of about 6 deg day −1 , close to that of the observed MJO. A smaller spatial scale of intraseasonal convection than the observed is noted in QOBS, possibly due to the lack of zonal SST gradients in the aqua-planet simulations associated with the warm pools and cold tongues as in the reality (cf, Figs. 1b and 3c). In contrast, westward propagation of intraseasonal convection with a weaker amplitude occurs in FLAT with a propagation speed of about 4 deg day −1 . The transition in propagation direction of intraseasonal convection from eastward to westward is clearly evident in model experiments forced by various intermediate SST profiles that lie between QOBS and FLAT ( Supplementary Fig. 1). Since model physics remains exactly the same among these experiments, this behavior clearly indicates that the propagation of tropical intraseasonal variability is strongly modulated by the large-scale mean state. The pronounced MJO-like eastward propagation of intraseasonal convection simulated in QOBS suggests that the model physics in the ECHAM AGCM is capable of representing the convective organization of the MJO. However, biases in the simulated MJO eastward propagation seen in Fig. 1b in ECHAM when forced by the observed SST could also be attributable to model biases in the mean state. Indeed, the lowertropospheric winter mean moisture pattern in the ECHAM AGCM illustrates large biases compared to its observational counterpart ( Supplementary Fig. 2), particularly with a significantly underestimated zonal moisture gradient over the Indian Ocean. While zonally uniform SST profiles are specified in experiments in this study, i.e., only the role of meridional moisture gradient on MJO propagation are represented, the importance of the zonal mean moisture gradient to observed MJO eastward propagation has been previously suggested 20,25 .
Physical mechanisms underlying distinct MJO propagation Consistent with MJO moisture mode theory, the propagation of intraseasonal convection anomalies in QOBS and FLAT are closely associated with propagation of their corresponding lowertropospheric moisture perturbations ( Supplementary Fig. 3). Diagnosis of the moisture budget associated with the intraseasonal variability in these two experiments can thus provide important insights into key processes regulating propagation of intraseasonal convection (see "Methods"). While moisture tendencies due to the column process (sum of vertical moisture advection minus apparent moisture sink; see "Methods") do significantly project onto the total moisture tendencies in both QOBS and FLAT ( Supplementary Fig. 6a), its induced moistening is largely collocated with positive moisture anomalies associated with enhanced intraseasonal convection (Supplementary Figs. 4,5), thus largely contributing to maintenance of the intraseasonal variability. The modest eastward shift of moisture tendencies produced by the column process relative to maximum moisture anomalies (i.e., convection center) in both QOBS and FLAT suggest that the eastward propagation of intraseasonal moisture anomalies and convection in QOBS and westward propagation in FLAT are primarily due to the lower-tropospheric moistening and drying by horizontal moisture advection (Fig. 4a, b, Supplementary Figs.  4-6). Specifically, a key process for the eastward propagation of intraseasonal convection in QOBS is the moistening (drying) to the east (west) of convection (Fig. 4a) through the meridional advection of background moisture by intraseasonal wind anomalies (Fig. 4c), a process has been previously illustrated to be crucial for observed eastward MJO propagation 20,22 . Given the poleward (equatorward) intraseasonal wind anomalies to the east (west) of convection (Fig. 4c), the presence of a moisture maximum near the equator is key for the eastward propagation in QOBS. Anomalous horizontal advection moistening by high-frequency eddies also contributes to eastward propagation of the intraseasonal variability in QOBS ( Supplementary Fig. 6b), in agreement with previous findings for the observed MJO 33 .
In FLAT, as the mean moisture maxima shift to 13°N and 13°S associated with the double-ITCZ in the mean rainfall pattern, advection across this mean moisture pattern by intraseasonal winds leads to very weak drying near the convection center and slightly to its east, with moistening to the west of convection (Fig. 4d). Strong moistening to the west of intraseasonal convection (Fig. 4b), thus the westward propagation of convection in this experiment, is promoted most strongly by advection of intraseasonal moisture anomalies by the mean easterly flow ( Fig. 4f; Supplementary Fig. 6b). By contrast, the moisture tendency due to this process is rather weak in QOBS (Fig. 4e) because of the weak mean easterlies in the lower-troposphere (cf. Figs. 2d, e and 4e, f). The moisture budget analyses thus suggest that the contrasting propagation characteristics of intraseasonal convection in QOBS and FLAT result from their distinct large-scale environments, including differences in mean moisture gradients and equatorial zonal winds, corresponding to the different SST profiles specified in these experiments.
Note that while the wind-induced surface heat exchanges (WISHE) mechanism has been found to be critical for the MJO-like intraseasonal variability in several recent aqua-planet modeling studies 34,35 , this process does not seem essential for the intraseasonal variability simulated in ECHAM AGCM in this study. The amplitude of anomalous surface heat fluxes associated with intraseasonal variability in both QOBS and FLAT are rather weak compared to the vertically integrated moisture tendencies by horizontal advection and column processes ( Supplementary Fig.  5). The role of surface fluxes for growth and propagation of the intraseasonal variability is also not consistent between the QOBS and FLAT (Supplementary Fig. 5). As previously discussed, in both QOBS and FLAT the column process in which the radiative effect is embedded is critical for growth of intraseasonal convection, while horizontal moisture advection is critical for regulating propagation of the intraseasonal convection.
The crucial role of the mean lower-tropospheric moisture gradient for governing propagation of intraseasonal convection is further demonstrated in various ECHAM experiments forced by intermediate SST profiles between QOBS and FLAT. Figure 5a presents the evolution of the mean lower-tropospheric meridional moisture and SST gradients in these experiments. Meridional moisture and SST gradients are defined by their differences between 10°N/S and the equator. The relationship between the propagation of intraseasonal convection, denoted by the ratio of eastward to westward rainfall spectral power on the time and space scales of 10-90 days and wavenumber 1-5 (hereafter E/W ratio 36 ; see "Methods"), and the mean meridional moisture gradient is further illustrated in Fig. 5b. An E/W ratio greater (less) than one corresponds to an eastward (westward) propagation of intraseasonal convection, and is largely consistent with longitude-time rainfall evolution diagrams ( Supplementary Fig. 1). Figure 5b clearly demonstrates that a transition in propagation of intraseasonal convection from eastward in QOBS to westward in FLAT is closely associated with a reversal in the sign of the mean meridional lowertropospheric moisture gradient, and is thus consistent with previous studies on the critical role of the mean moisture pattern in regulating propagation of intraseasonal convection. Note that the transition from positive to negative meridional moisture gradient is not smooth, and occurs rapidly when the meridional SST gradient is reduced to about 0.5 K, i.e., in the experiment e5, when a rapid transition from the eastward to westward propagation of the intraseasonal convection is also found.
Also noteworthy is that MJO propagation in these model simulations does not show a strong relationship with the relative intensity of equatorial low-level easterly anomalies that are largely associated with the Kelvin wave and westerly anomalies associated with the Rossby wave as suggested by the coupled Rossby-Kelvin wave framework for the MJO 12,13 . For example, corresponding to the strong eastward propagating intraseasonal convection in QOBS, equatorial westerly anomalies at 850 hPa in the Rossby-wave component to the west of convection are stronger than easterly anomalies to the east of convection (Fig. 4a). Similarly, associated with westward propagation of intraseasonal convection in FLAT, the 850 hPa anomalous easterlies to the east of intraseasonal convection X. Jiang et al. exhibit stronger amplitude than westerlies with the Rossbywave component (Fig. 4b). There is no significant correlation between the E/W ratio and relative intensity of equatorial lowlevel easterly and westerly anomalies to the east and west of intraseasonal convection, respectively, across these experiments (figure not shown). Also note that a westward-tilted anomalous moisture structure with altitude is evident in both QOBS and FLAT (Fig. 4a, b), associated with both eastward and westward-propagating intraseasonal variability. Therefore, a westward tilting structure, previously considered an important signal of moisture preconditioning responsible for the eastward propagation of the MJO, is not sufficient to drive eastward propagation of intraseasonal convection.

DISCUSSION
In this study, how the propagation of tropical intraseasonal variability relates to the large-scale environment is demonstrated based on a series of idealized model experiments using the ECHAM AGCM. Forced by various zonally uniform SST profiles with relatively small differences in their meridional gradients near the equator, pronounced differences are found in both the mean state and propagation characteristics of the intraseasonal convection. When a sharp SST peak is imposed over the equator, a single-ITCZ type mean rainfall pattern prevails in the model along with pronounced eastward propagation of intraseasonal convection. As the equatorial SST gradient flattens, the mean rainfall gradually transitions to a double-ITCZ pattern, while eastward propagation of intraseasonal convection weakens and eventually transitions to westward propagation. Intraseasonal moisture perturbations coincide with convection perturbations and have the same propagation characteristics, in support of moisture mode theory for the MJO. The eastward propagation of intraseasonal convection associated with the single-ITCZ mean rainfall pattern is found to be promoted by meridional advection of mean lowertropospheric moisture by the intraseasonal winds, a mechanism also found critical for the observed MJO propagation. This process becomes weaker and even reverses its sign under a double-ITCZ environment; as a result, westward-propagating intraseasonal convection is favored in these experiments with flattened SST due to advection of intraseasonal moisture perturbation by the mean equatorial easterlies. These results thus suggest a possible link between the great challenges GCMs have in representing the MJO and their tendency to produce double-ITCZ biases 37 .
Under an idealized aqua-planet setting, this study clearly illustrates that propagation of intraseasonal convective variability is strongly modulated by the large-scale background state, in which the lower-tropospheric moisture gradient and mean zonal winds play a crucial role. In order to improve the representation of the MJO in GCMs, previous community efforts largely emphasized the importance of model physics, particularly those depicted in cumulus parameterizations. Results in this study, however, underscore the importance of maintaining a realistic model mean state when tuning these cumulus and other parameters in order to achieve faithful model representation of the MJO.
Note that idealized zonally symmetric SST profiles are specified in the above model experiments, and therefore the influence of mean moisture on MJO propagation is only felt through its meridional distribution. In reality, zonal gradients of lowertropospheric moisture also play a crucial role in promoting MJO eastward propagation, particularly over the Indian Ocean 20 . Sensitivity of intraseasonal variability to the zonal SST gradient, for example, forced by a warm pool-type SST pattern has also been previously reported in aqua-planet model simulations 38,39 .
In a future study, it would be interesting to conduct a multimodel comparison to further understand how model intraseasonal variability responses to the exactly same SST profiles in different GCMs under an aqua-planet configuration, inspired by previous comparisons 29,32,39 . Also note that the experiments conducted here were AGCM experiments with specified SST forcing. Inclusion of atmosphere-ocean interactions can improve representation of the MJO 8 , although the underlying physics still are a matter of debate. For example, the MJO simulation is significantly improved when ECHAM AGCM is coupled to an ocean model 4,40 .

Model and experimental design
The model used in this study is the ECHAM AGCM version 4.6 (ref. 41 ) at a resolution of T30 and 19 vertical levels. The convection scheme used in this GCM is the mass flux scheme for penetrative, shallow, and midlevel convection of Tiedtke modified by Nordeng such that the cloud-base mass flux is linked to convective instability for the penetrative convection 42,43 . This AGCM was first integrated for 20 years with lower-boundary forcing from climatological monthly mean SST, and rainfall evolution associated with the MJO and the mean lower-tropospheric moisture pattern are shown in Fig. 1b and Supplementary Fig. 2b, respectively. A series of idealized aqua-planet simulations were then conducted with the ECHAM AGCM to exclude complex influences from boundary conditions, including land-sea contrasts and orography. In these aqua-planet experiments, the diurnal cycle of solar radiation is turned off and set to be perpetual condition of March. Sea ice is set to zero and aerosols are switched off. To further simplify the processes that regulate propagation of tropical convection, zonally uniform SST patterns are used to force the AGCM with various latitudinal profiles. These SST profiles are defined following a previous study 28 , and are specified as follows: where φ is latitude, and a and b are the parameters controlling the meridional SST gradient near the equator (b = 1 − a). Model experiments suggest that propagation of intraseasonal convection exhibits greatest sensitivity when the parameter a varies from 0.5, which is referred as the QOBS SST since it is closest to the observed annual and zonal mean SST, to 0, referred as the FLAT SST. SST profiles for QOBS and FLAT are shown in Fig. 2a. Additional experiments are conducted where the parameter a is incrementally changed by 0.05 between 0.5 and 0 to examine how propagation of intraseasonal convection gradually changes with the SST profiles (Fig. 5).

Observation data
The European Centre for Medium-Range Weather Forecasts (ECMWF) reanalysis ERA-Interim 44 with a resolution of 1.5°longitude by 1.5°latitude was used for the observed lower-tropospheric moisture pattern shown in Supplementary Fig. 2. Rainfall observations from the Tropical Rainfall Measuring Mission (TRMM, version 3B42) 45 are also used for this study. TRMM3B42 is a global precipitation product based on multi-satellite and rain gauge analyses. It provides long-term precipitation estimates gridded at a three-hourly temporal resolution and 0.25°spatial resolution in a global belt extending from 50°S to 50°N from 1997 to 2015. Phase speed and structure of the intraseasonal variability The evolution of intraseasonal convection in simulations are derived by lag regression of 10-90-day filtered rainfall anomalies against a reference time series of its averaged value over an Indian Ocean box (75-85°E; 5°S-5°N). The propagation speed of intraseasonal convection can be estimated by the slope of features in an equatorially averaged longitude versus time diagram of rainfall anomalies based on lag-regression rainfall evolution (see Fig. 3c, d and Supplementary Fig. 3 for examples). One disadvantage of this approach lies in the uncertainty in deriving propagation speed during the transition of intraseasonal convection from eastward to westward propagation for certain SST patterns (e.g., e5 and e7 in Supplementary Fig. 1). An alternative approach for calculating propagation speed of intraseasonal convection is to use the ratio of eastward to westward rainfall spectral power on 10-90-day timescales and wavenumbers 1-5 spatial scales following previous work 36 , referred as the E/W ratio. An E/W ratio greater (less) than one corresponds to eastward (westward) propagation of intraseasonal convection, which is largely consistent with longitude-time rainfall evolution diagrams with a stronger E/W ratio corresponding to faster eastward propagation speed ( Supplementary Fig.  1). Three-dimensional structure of the intraseasonal variability in model simulations can be derived by lag regression of 10-90-day filtered intraseasonal anomalies of various variables onto the same reference time series of 10-90 day filtered rainfall averaged over an Indian Ocean box. Note that in model simulations with zonally symmetric SST, results will not change significantly if a different averaging box other than 75-85°E; 5°S-5°N is selected.

Moisture tendency analysis for MJO propagation mechanisms
The local time rate of change of specific humidity can be written as where the terms on the right-hand side are the horizontal and vertical moisture advection, and the apparent moisture sink Q 2 as defined in a previous study 46 , which represents the combined effects of evaporation, condensation, sublimation, and deposition within the column and the flux of moisture by unresolved eddies 47 . To diagnose model processes, Q 2 is directly archived from the ECHAM AGCM at daily intervals, which technically include total moisture tendencies from modules of vertical diffusion, subgrid cumulus, and large-scale condensation. Since there is a near cancellation between the vertical moisture advection and Q 2 , the second and third terms on the right-hand side are combined and referred to as the column process following previous work 48 . As previously described, longitude versus pressure profiles of moisture fields and each moisture tendency term in Eq. (2) associated with intraseasonal variability in model simulations can be derived by lag-0 regression onto the Indian Ocean intraseasonal rainfall index (see Fig. 4a, b for an example of anomalous moisture q′, and horizontal moisture advection Àðṽ Á ∇qÞ 0 profiles). The importance of each moisture tendency term to the total moisture tendency can be estimated by the pattern projection of these tendency terms onto the total moisture tendency on a vertical pressurelongitude space of 30-150°E; 900-300 hPa as shown in Supplementary  Fig. 6a. In order to identify the detailed processes associated with the total horizontal moisture advection, a decomposition of this term is further performed by separating daily horizontal wind and specific humidity into three different timescales following previous studies 20 , i.e., low-frequency (period >90 days, with the mean seasonal cycle included), intraseasonal (MJO; 10-90 days), and high-frequency (<10 days) timescales. Similarly, contribution from each of these decomposed horizontal moisture advection terms to the total intraseasonal moisture tendency in each model simulation can also be derived by pattern projection in the pressure-longitude plane, with results shown in Fig. Supplementary 6b.

DATA AVAILABILITY
All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. The TRMM3B42 rainfall data was downloaded from https://disc.gsfc.nasa.gov/datasets/TRMM_3B42_7/summary (https://doi.org/ 10.5067/TRMM/TMPA/3H/7). The ERA-Interim reanalysis data can be downloaded from the website: http://apps.ecmwf.int/datasets/. Output from ECHAM4 model experiments analyzed in this study is available at the website: https://ucla.box.com/v/ echam4-aquaplanet-mjo.