The Madden–Julian oscillation (MJO) skeleton model offers a theoretical prediction of the MJO’s structure. Here, a method is described for identifying this structure in observational data. The method utilizes projections onto equatorial wave structures, and a main question is: Can this method isolate the MJO without using temporal filtering or empirical orthogonal functions? For the data projection, a wide range of data is incorporated: multiple variables (wind, geopotential height, water vapor, and, as a proxy for convective activity, outgoing longwave radiation), multiple pressure levels (850 and 200 hPa), and multiple latitudes (both equatorial and off-equatorial). Such a data variety is combined using a systematic method, and it allows for a distinction between the Kelvin and Rossby components of the MJO’s structure. Results are illustrated for some well-known cases, and statistical measures are presented to quantify the variability of the MJO skeleton signal, MJOS(x, t), and its amplitude, MJOSA(t). The robustness of the methods is demonstrated through a suite of sensitivity studies, including tests with two projection methods. When the projection is based on the skeleton model’s energy, as opposed to the standard L2 energy, water vapor is seen to be of primary importance. Finally, a simple interpretation is given for the MJO skeleton structure: it is related to the wave response to a moving heat source. From either perspective, the methods here identify signals that project onto coupled convection–circulation patterns, and the results suggest that a large portion of the MJO’s structure is consistent with such a coupled pattern.
To identify the Madden–Julian oscillation (MJO), a large number of methods have been used. In their pioneering studies, Madden and Julian (1971, 1972) detected the oscillation in zonal wind, pressure, and temperature data. Subsequently, many MJO indices have been defined using different fields of interest: winds only (e.g., Chen and Del Genio 2009), cloudiness or precipitation only (e.g., Kiladis et al. 2005), or combinations of winds and cloudiness (e.g., Wheeler and Hendon 2004).
In a sense, each MJO index offers a different definition of what the MJO is. A commonality among these MJO indices is that they are all empirical definitions of the MJO. In contrast, by using a model of the MJO, one can offer a theoretical definition of the MJO.
In recent theoretical work, Majda and Stechmann (2009, 2011) proposed the MJO skeleton model, which predicts the fundamental features of the MJO on intraseasonal and planetary scales: slow eastward phase speed of roughly 5 m s−1, peculiar dispersion relation with ω(k) ≈ constant, and horizontal quadrupole vortex structure. In subsequent work, a stochastic version of the skeleton model has also been developed (Thual et al. 2014a,b).
One of the main goals of the present paper is to investigate the MJO skeleton model’s predicted structure in observational data. To what extent is the structure of an individual MJO event consistent with the MJO skeleton structure?
To this end, the methods here will account for multiple variables and multiple basis functions. For example, data from four variables will be used: wind, geopotential height, water vapor, and, as a proxy for convective activity, outgoing longwave radiation (OLR). To combine such a large variety of data, a systematic method is presented, based on equatorial wave theory and the MJO skeleton model. As another example, since the MJO skeleton structure includes both Kelvin and Rossby components, both equatorial and off-equatorial data are needed here to distinguish the MJO’s off-equatorial Rossby gyres. Again, equatorial wave theory provides a systematic method for accounting for equatorial and off-equatorial structures.
Another main question here is: Can the MJO be isolated in observational data without the use of temporal filtering or empirical orthogonal functions (EOFs)? These two techniques are commonly used to identify the MJO, whereas the present paper does not use these methods.
To identify the MJO, many indices have been defined using different fields of interest: winds only (e.g., Knutson and Weickmann 1987; Chen and Del Genio 2009), cloudiness or precipitation only (e.g., Rui and Wang 1990; Kiladis and Weickmann 1992; Hendon and Salby 1994; Wheeler et al. 2000; Kiladis et al. 2005; Matthews 2008; Kiladis et al. 2014), or combinations of winds and cloudiness (e.g., Weare 2003; Wheeler and Hendon 2004; Ventrice et al. 2013). More extensive summaries and comparisons were given recently by Straub (2013) and Kiladis et al. (2014). The use of multiple variables, such as winds and OLR together, is in line with the conceptual viewpoint of the MJO as a coupled convection–circulation phenomenon. Here, our methods build further on this idea of including more variables into the definition of the MJO’s structure.
Equatorial wave theory distinguishes several different wave types, partly through their different meridional structures. The natural meridional basis functions are the parabolic cylinder functions, which have been useful in several observational analyses of convectively coupled equatorial waves (Yang et al. 2003, 2007; Gehne and Kleeman 2012). Here, in addition, the different wave types—Kelvin, Rossby, etc.—are also distinguished by their different eigenvectors. The eigenvectors used here are from the long-wave version of the theory of Matsuno (1966), which is the equatorial analog of the normal mode expansion of the primitive equations on the sphere (Kasahara 1976; Kasahara and Puri 1981).
In the analysis below, an MJO signal MJOS(x, t) is defined, and it could potentially be used as an MJO index. The definition of MJOS(x, t) is a projection onto a theoretical MJO structure, and it is shown that the structure has a more general interpretation beyond the MJO skeleton model. In particular, it will be shown that the MJO skeleton structure is related to the wave response to a moving heat source. As such, MJOS(x, t) provides insight into the moving-heat-source perspective of the MJO as well.
The rest of the paper is organized as follows. The MJO skeleton model is reviewed in section 2, and the data are described in section 3. The methods are described and illustrated in section 4, including the identification of Kelvin, Rossby, and MJO structures. The results are described in more detail in section 5, including case studies, statistics, and sensitivity studies. In section 6, the MJO skeleton structure is related to the wave response to a moving heat source. Finally, conclusions are summarized in section 7.
2. The MJO skeleton model
Since the data analysis requires observational surrogates for the MJO skeleton variables, we first present an overview of the model itself before proceeding to the details of the data and methods.
In its three-dimensional form, the MJO skeleton model of Majda and Stechmann (2009) is
Here, u, υ, and w are the zonal, meridional, and vertical velocity anomalies, respectively; p, θ, and q are the pressure, potential temperature, and water vapor anomalies, respectively; and a is the amplitude of convective activity, which is a nonnegative quantity. The convective heating and drying are taken to be proportional to the envelope of convective activity: . All other sources of cooling and moistening are represented by sθ and sq, respectively. The equations have been nondimensionalized using the reference scales shown here in Table 1 (Stechmann et al. 2008).
Notice that this model contains a minimal number of parameters, summarized in Table 1: , the (nondimensional) mean background vertical moisture gradient, and Γ = 1, or Γ ≈ 0.47 day−1 (g kg−1) in dimensional units. [See Stechmann et al. (2013) for a theoretical estimate of the value of parameter Γ.] These will be the standard parameter values used throughout the paper unless otherwise stated. The source terms sθ and sq must also be specified (see below). The parameter is used to give the units of a heating rate while keeping a nondimensional, and it is actually irrelevant to the dynamics (as can be seen by rescaling a).
To obtain the simplest model for the MJO, truncated vertical and meridional structures are used, as described next.
To arrive at a two-dimensional form of (1), a vertical truncation is used: only the first baroclinic mode is retained so that , etc., with a slight abuse of notation in the reuse of symbol u. The resulting equations resemble a time-dependent version of a Matsuno–Gill model (Matsuno 1966; Gill 1980), but without damping, plus equations for q and a:
Here, the moisture q(x, y, t) represents the lower-tropospheric water vapor, which has been shown to lead the MJO’s heating anomaly in observations (Myers and Waliser 2003; Kikuchi and Takayabu 2004; Kiladis et al. 2005; Tian et al. 2006).
Finally, to arrive at a one-dimensional form of (2), a meridional truncation is used. It is assumed that a, the envelope of convective activity, has a simple equatorial meridional structure proportional to exp(−y2/2). For the long-wave-scaled equations, such a meridional heating structure is known to excite only Kelvin waves and the first symmetric equatorial Rossby waves (Matsuno 1966; Gill 1980; Majda 2003), and the resulting meridionally truncated equations can be written as
The system (3) is written in terms of the anomaly A from a radiative–convective equilibrium state, , uniform in space and time. The total, positive amplitude of convective activity is then , and the quantities K and R are the amplitudes of Kelvin and equatorial Rossby wave structures, respectively.
To identify surrogates of the variables in (1)–(3), two data sources are used here. First, as a surrogate for convective activity a, the National Oceanic and Atmospheric Administration (NOAA) interpolated outgoing longwave radiation (OLR) is used (Liebmann and Smith 1996). While many possible surrogates exist for convective activity a, OLR is used here as a simple choice to initially investigate. Second, National Centers for Environmental Prediction–National Center for Atmospheric Research (NCEP–NCAR) reanalysis data are used for all other variables: zonal wind, geopotential height, and specific humidity (Kalnay et al. 1996). Despite the fact that reanalysis data are a combined product of observational data and model dynamics, we sometimes refer to it as “observational” for simplicity. Both datasets have a horizontal spatial resolution of 2.5° × 2.5° and a daily temporal resolution from 1 January 1979 to 31 December 2011. For comparison with the MJO skeleton signal, data of the Real-time Multivariate MJO (RMM) index of Wheeler and Hendon (2004) and the OLR MJO index (OMI) of Kiladis et al. (2014) will also be used here.
4. Methods: Defining the MJO skeleton variables and eigenmodes
The end goal of this section is to define the MJO skeleton structure and variability in terms of observational data. To this end, the first two steps follow the path from (1) to (2) to (3): a vertical mode truncation to move from 3D (x, y, z) to 2D (x, y) (section 4a), and then a meridional mode truncation to move to 1D (x) (section 4b). A zonal Fourier mode truncation is also used to isolate the planetary spatial scales (section 4c). Furthermore, observational definitions must be chosen for the lower-tropospheric water vapor Q and the amplitude of convective activity A (section 4d). Finally, the linear eigenvectors of the MJO skeleton model are used to provide the ingredients for defining the MJO skeleton signal MJOS(x, t) (section 4e).
In all figures presented here, a smoothed seasonal cycle is removed via the annual mean and the first three harmonics, and a 120-day running time average is subtracted on each day at each spatial point to remove low-frequency variability, as is recommended by the Climate Variability and Predictability (CLIVAR) Madden–Julian Oscillation Working Group (2009; see also Gottschalck et al. 2010).
a. Vertical mode definition
To define the first baroclinic mode variables u(x, y, t) and θ(x, y, t) from (2), a simple method is used here. The method utilizes data from only two pressure levels, 850 and 200 hPa, in order to capture the essence of comprehensive treatments of vertical basis functions (e.g., Fulton and Schubert 1985; Majda 2003). While the model (1)–(2) is formulated in terms of the height coordinate z, the reanalysis data are given in terms of a vertical pressure coordinate p. To connect the two formulations in a simple way, rather than reformulating the model (1)–(2) in terms of a vertical pressure coordinate, an association is made between pressure levels and vertical levels:
where, in nondimensional height units, z = 0 corresponds to a level near the bottom of the free troposphere (and top of the boundary layer) and z = π corresponds to a level near the tropopause.
Given velocity data u at the two levels in (4), it is natural to associate the data with contributions from two vertical modes. With a barotropic contribution uBT and a first baroclinic mode contribution uBC, the velocity is written as
where 1 and are crude representations of the barotropic and baroclinic mode vertical structures, respectively, and the is a normalization constant. At the specific levels of 850 and 200 hPa, this is
By subtracting these two equations, the value of uBC can then be approximated as
This definition of uBC is used to represent the variable u(x, y, t) from (2), where the subscript BC is dropped to ease notation.
For θ, the derivation follows similarly. With a first baroclinic mode contribution only, θ is written as
Instead of obtaining this quantity from temperature data in the middle troposphere (z ≈ π/2), it is instead obtained from the geopotential height Z in the lower and upper troposphere. Using
This provides an expression for the potential temperature anomaly θ from the model, in terms of the geopotential height anomaly Z from the reanalysis data. For nondimensionalization, the reference θ and Z scales are and c2/g ≈ 265 m, respectively. [Alternative parameter choices were explored, but no compelling reason was found for departing from the standard values used in the previous studies of Stechmann et al. (2008) and Majda and Stechmann (2009).] This definition of θBC is used to represent the variable θ(x, y, t) from (2), where the subscript BC is dropped to ease notation.
b. Meridional modes and K and R
Given the first baroclinic mode variables u(x, y, t) and θ(x, y, t) from (2), the next step is to define the Kelvin and equatorial Rossby wave variables, K(x, t) and R(x, t), from (3). This is accomplished using meridional basis functions, ϕm(y), which are the parabolic cylinder functions,
where Hm(y) are the Hermite polynomials:
The functions ϕm(y) form an orthonormal basis, and the variables u(x, y, t) and θ(x, y, t) can then be expanded as
where the quantities um(x, t) are obtained using the projection
These functions are plotted in Fig. 1, which illustrates the significant off-equatorial contribution of ϕ2(y) in the 15°–30° latitude band. Utilizing these basis functions, the definitions of K and R are then
where um and θm are the meridional projections from (15) and (16). See appendix A for more details of the derivation. The important features of the K and R structures are illustrated in Fig. 2. While the K structure is equatorial with an exp(−y2/2) decay away from the equator, the R structure includes off-equatorial gyres that arise from the ϕ2(y) terms in (20).
While (13)–(20) describe the theory on an equatorial β plane with −∞ < y < +∞, the data analysis must use latitude on a sphere. To relate latitude and meridional length, a simple correspondence is used here: 1° latitude ≈ 110 km. Other basis functions, such as spherical harmonics, are applicable for the entire sphere and include variations at all latitudes. Here, the parabolic cylinder functions are most applicable since the focus is on the equatorial regions only.
The observed variability of K(x, t) is shown in Fig. 3a. Also shown are u0(x, t) and θ0(x, t), which are combined to form K(x, t), as in (19). The time period is 1 July 1992–1 July 1993, including the period of the Tropical Ocean and Global Atmosphere Coupled Ocean–Atmosphere Response Experiment (TOGA COARE). The variability of K(x, t) includes any baroclinic variability whose meridional structure projects onto the structure of the Kelvin wave. For example, this should include dry Kelvin waves, convectively coupled Kelvin waves, and the MJO. In Fig. 3a, eastward propagation dominates the variability of K(x, t), and a variety of propagation speeds are visible. Relatively fast propagation can be seen in the eastern Pacific sector (180°–90°W) during many time periods (e.g., 15 December 1992–1 February 1993), reminiscent of that described by Milliff and Madden (1996), that is, of far-field dry Kelvin modes often observed during or after MJO or other large-scale convective events over the warm pool. Also, relatively slow propagation can be seen in the western Pacific sector (120°E–180°) during January 1993, which coincides with the MJO events during TOGA COARE (Chen et al. 1996; Yanai et al. 2000). The signature of these features can also be seen in the variability of u0(x, t) and θ0(x, t) in Figs. 3b,c. Furthermore, in u0(x, t) and θ0(x, t), westward propagation can be seen (e.g., 1 November 1992, 60°E), but it is filtered out by the definition of . In the variability of K(x, t), it is difficult to discern coherent westward-propagating signals. In short, the features in Fig. 3 lend confidence to the data analysis method; they suggest that the simple methods proposed here can distinguish Kelvin wave structures by using a particular combination of winds and geopotential height.
Similarly, the observed variability of R(x, t) is shown in Fig. 4a. Also shown are u2(x, t) and θ2(x, t), which are combined with u0(x, t) and θ0(x, t) to form R(x, t) as in (20). The variability of R(x, t) includes any baroclinic variability whose meridional structure projects onto the structure of the equatorial Rossby wave. Dry or convectively coupled equatorial Rossby waves could potentially have a nonzero projection onto this structure, and these are the source of the westward-propagating signals in R(x, t) in Fig. 4a. Also, the MJO’s structure should have a nonzero projection onto the Rossby wave structure, and it is likely the source of the eastward-propagating signals in R(x, t) in Fig. 4a (e.g., from 60°E to 120°W during January). What is not present in Fig. 4a is any rapid eastward propagation, which suggests that the simple methods proposed here can distinguish equatorial Rossby wave structures from the Kelvin wave structures.
c. Zonal wavenumber selection
To isolate the planetary zonal scales, a Fourier truncation is also applied to retain only zonal wavenumbers k = −3, −2, −1, +1, +2, and +3, which will be referred to loosely as k = 1–3. The effect of this truncation can be seen by comparing the original K(x, t) in Fig. 3a with the filtered K(x, t) in Fig. 5a. Many smaller-scale features have been filtered out, but the main features are present in both plots, which suggests the main features of K(x, t) reside on planetary scales.
Similarly, the original and filtered R(x, t) are shown in Figs. 4a and 5b, respectively. In this case, the filtered and unfiltered plots appear quite different, which suggests that many of the main features of R(x, t) do not reside in planetary wavenumbers k = 1–3. After the filtering, the strongest remaining signal is eastward propagating during the December–February (DJF) season, because of the Rossby gyre structure of the MJO events during TOGA COARE.
By including multiple zonal wavenumbers, a variety of zonal MJO structures can potentially be identified. For example, an MJO structure localized near the western Pacific warm pool—and not extending around the globe as a k = 1 sinusoid—could arise from a combination of multiple zonal wavenumbers. Also, it is possible that the MJO has some nonlinearity to its structure, as suggested in the nonlinear MJO skeleton results of Majda and Stechmann (2011) and Chen and Stechmann (2014, manuscript submitted to Comm. Math. Sci.). In that work, the convectively active region appears to be narrower than the convectively suppressed region. To represent such features with sinusoids would require multiple zonal wavenumbers.
In short, by applying a zonal spatial filter without temporal filtering, the signatures of a variety of Kelvin and Rossby structures can be distinguished. This zonal spatial filter is used in all plots that follow unless indicated otherwise.
d. Moisture and convection: Q and A
The variable Q is meant to represent the lower-tropospheric water vapor anomaly near 850 hPa. In several studies, it has been shown to lead the MJO’s convective anomaly (Myers and Waliser 2003; Kikuchi and Takayabu 2004; Kiladis et al. 2005; Tian et al. 2006). Accordingly, in this first study, we define
Note that only the projection of q onto ϕ0(y) is used here, as described in appendix A [although to ease notation this is not indicated symbolically in (21)]. We have also explored some alternatives to the definition (21) (e.g., alternative weights or alternative pressure levels) but did not notice any significant differences in the results.
The observed variability of Q(x, t) is illustrated in Fig. 5c. Standing oscillations appear more prominently than eastward- or westward-propagating signals, in contrast to the K and R signals. This is possibly because Q variability is affected by such a wide range of dynamical and thermodynamical influences. Also notice that Q has been nondimensionalized by the natural reference scale , where is the reference potential temperature scale. In nondimensional units, the amplitude of Q is comparable to that of K and R, which suggests that the standard reference scales in Table 1 are a reasonable choice.
The variable A is meant to represent the planetary-scale envelope of convective activity. One could imagine many possible surrogates for A. Here, OLR is used as a surrogate for A:
where A and OLR are both anomalies and are projected onto ϕ0(y). Past estimates of this type of relationship have suggested HOLR values in the range of 0.08–0.25 K day−1 (W m−2)−1 (Christy 1991; Yanai and Tomita 1998). In a recent study by Stechmann and Ogrosky (2015), a Walker circulation pattern is driven by the implied diabatic heating, and it is used to obtain an estimate of 0.06 K day−1 (W m−2)−1, which will be adopted here as the standard value.
The observed variability of A(x, t) is shown in Fig. 5d. The strongest signals appear during TOGA COARE, during the period from 15 December 1992 to 15 February 1993. These strong signals are eastward propagating and associated with the MJO events that occurred during this period (Chen et al. 1996; Yanai et al. 2000; Straub 2013). Other weaker signals are also present and are mostly eastward propagating, but some standing oscillations are present (e.g., from 15 October 1992 to 5 November 1992), and some westward-propagating signals are present (e.g., from 15 August to 1 October 1992). In nondimensional units, the amplitude of A is comparable to that of K, R, and Q, which suggests that the standard value of HOLR = 0.06 K day−1 (W m−2)−1 (Stechmann and Ogrosky 2015) is a reasonable choice for the relationship in (22) on the scales of interest here.
e. Eigenvectors of the MJO skeleton model
From the set of variables (K, R, Q, and A) defined above, a particular combination can be used to construct the MJO skeleton structure. Figure 6 shows the contributions of K, R, Q, and A to the MJO structure. The relative contributions are different for each wavenumber k, and the contributions of K and R decrease relative to the contributions of Q and A as wavenumber increases. Also noteworthy is that the R component is substantial for each wavenumber, which suggests the importance of the MJO’s off-equatorial features in defining its structure.
In addition to the amplitudes shown in Fig. 6, the phase relationships among the variables (K, R, Q, and A) are also needed to complete the definition of the MJO skeleton structure. The phase relationships are shown in terms of velocity, pressure, etc. in Fig. 7. Enhanced convective activity appears at the center of zonal convergence and a quadrupole vortex structure. The relative amplitudes of K and R can also be discerned from Fig. 7a, since the pressure and zonal velocity arise from a linear combination of the K and R structures shown in Fig. 2.
The MJO structure in Figs. 6 and 7 has been chosen to match the MJO eigenvector of the linearized MJO skeleton model (Majda and Stechmann 2009). This is a theoretical prediction of the MJO’s structure, and it provides a method for combining multiple fields—winds, geopotential height, moisture, and convective activity—in a systematic way. Such an approach is in line with the viewpoint of the MJO as a coupled convection–circulation phenomenon.
The MJO skeleton eigenvector will be denoted by , and the following description is a brief summary of its derivation and properties (Majda and Stechmann 2009). The starting point is the linearized system (3). The state vector is U = (K, R, Q, A)T, where WT denotes the transpose of W. Using a plane-wave ansatz, , leads to an eigenvalue problem, where one of the four eigenmodes has properties that resemble the MJO. The model predicts an eigenvalue ωMJO(k) that is roughly constant as a function of wavenumber k, consistent with the observed MJO power spectrum (Kiladis et al. 2009). Several values of ωMJO(k) are listed in Table 2. The model also predicts the eigenvector , which was illustrated in Figs. 6 and 7 and is listed for reference in Table 2. It is a four-component column vector for each zonal wavenumber k, and each component describes the magnitude and phase of the contributions of , and . The eigenvector is complex valued and satisfies the condition , where f* denotes the complex conjugate of f, which ensures that the inverse Fourier transform is real valued. Furthermore, the eigenvector has been normalized to have unit magnitude, and an (arbitrary) overall phase factor has been chosen so that the MJO eigenmode is in phase with convective activity A.
The eigenvector can now be used as an identifier of the MJO’s structure in observational data. To do this, the observed state vector is projected onto the MJO eigenmode:
where W† denotes the conjugate transpose of W, and where the time dependence has been suppressed to ease notation. The inverse Fourier transform then leads to the real-valued scalar quantity
which we refer to as the MJO skeleton signal.
where the inner product matrix is defined in terms of the MJO skeleton model’s conserved energy. See appendix B for a description of this orthogonal projection, which is more in line with MJO skeleton theory. In what follows, the standard L2 projection in (23) will be used, due to the simplicity of its formulation.
As a case study for examining MJOS(x, t), the TOGA COARE case is used here. In the period from November 1992 to February 1993, the MJO activity during TOGA COARE is perhaps the most studied of any on record (Lin and Johnson 1996; Wheeler and Kiladis 1999; Yanai et al. 2000; Straub 2013). Based on OLR, two MJO-like convective events can be identified (Yanai et al. 2000): one in December and one in January. These two convective events are associated with rather different features, such as differing circulation patterns (Lin and Johnson 1996; Straub 2013), which makes this an interesting period for examining the different features (e.g., convection and circulation) of large-scale convective activity.
A first example of MJOS(x, t) variability is shown in Fig. 8. The time period 1 July 1992–1 July 1993 includes TOGA COARE and is identical to that used in Figs. 3–5. Eastward-propagating signals are significantly more common than westward-propagating signals, and a strong-amplitude signal appears in the DJF season and reflects the TOGA COARE MJO events. By comparing this signal with the K, R, Q, and A variability from Fig. 5, it is visually clear that the MJOS variability represents not an individual component K, R, Q, or A but a composite of the K, R, Q, and A components. This suggests that MJOS(x, t) represents the coupled convection–circulation variability of the MJO.
To succinctly quantify the strength of MJOS(x, t) variability, its amplitude MJOSA(t) can be computed as a root-mean square of the zonal variations:
where Pe is the circumference (period) of the Earth at the equator. A plot of MJOSA(t) for the time period 1 July 1992–1 July 1993 is shown in Fig. 9. The period from late December to early February is identified in MJOSA(t) as a time of intense MJOS(x, t) amplitude relative to climatological values. In other words, given the definition of MJOS(x, t), this is a period when the coupled convection–circulation structures, or the projections onto such structures, are particularly strong.
In this section, the variability of MJOS(x, t) is explored further with case studies, statistical measures, and sensitivity studies. Among other questions, we build on the TOGA COARE case of section 4 by asking the following question: In other cases of known strong MJO events, how much of the MJO structure projects onto the theoretical, coupled convection–circulation structure?
a. Case studies
Three additional case studies are now considered to illustrate the range of variability of MJOS(x, t). One case is the well-known period of strong MJO activity from 1 July 1987 to 1 July 1988 (Hendon and Liebmann 1994; Wheeler and Hendon 2004). This is the middle case shown in Fig. 10. The other two cases in Fig. 10 are periods when the presence and/or initiation of MJO activity is somewhat ambiguous (Straub 2013). On the left is the period from 1 July 1986 to 1 July 1987, and on the right is the period from 1 July 1997 to 1 July 1998, which includes the time period of May–June 1998 when the South China Sea Monsoon Experiment (SCSMEX) was being conducted. It should also be noted that a portion of 1986–88 was a strong El Niño period, and hence convection was shifted eastward compared to its climatological mean. Three corresponding plots of amplitude MJOSA(t) are shown in Fig. 11.
In the 1987–88 case in Fig. 10b, several periods show strong MJO variability. The MJOS(x, t) signal is strongest in November–December and February–March, with a moderate signal in the interim during January. An interesting feature is the circumnavigating signal through February–April. The Western Hemisphere component of the signal could be a strong instance of the fast “dry” Kelvin wave signals described by Milliff and Madden (1996), although in this case one can discern concomitant signals in all of the four variables K, R, Q, and A (not shown), which suggests this is an instance of coupled convection–circulation variability.
The 1986–87 case has been ambiguous in earlier work. Wavenumber–frequency-filtered OLR indicates a strong active MJO event throughout January and February, whereas the RMM index indicates a strong MJO event from late February to early March (Straub 2013). Here, MJOS(x, t) does not indicate a particularly strong signal for either of those two periods. A short burst appears at the beginning of March, although there is no particularly strong MJO event indicated by MJOS(x, t) during January–March. It seems that the K, R, Q, and A signals are not sufficiently interrelated to indicate a coupled convection–circulation structure at these times.
The final case of 1997–98 has also been ambiguous in earlier work with other MJO indices (Straub 2013). Complicating the picture is a synoptic-scale convectively coupled Kelvin wave that appears in mid-May (Straub et al. 2006; Straub 2013); its rapidly propagating signature is evident in MJOS(x, t) in Fig. 10c near the date line, from 150°E to 150°W. However, in late April and early May, a planetary-scale signal is also present in MJOS(x, t). Its strongest anomalies are negative from the date line to 60°W and positive from 0° to 90°E, and its slower propagation speed is consistent with the MJO. From this information in the MJOS(x, t) signal, it is possible that an MJO is present in early May and then transforms into the convectively coupled Kelvin wave in mid-May, or it is possible that the MJO signal is present throughout May but is contorted by the concomitant signal of the convectively coupled Kelvin wave. The presence of an MJO is also suggested by a strong Rossby signal R(x, t) that is present during May and propagates slowly eastward (not shown).
An interesting feature in Fig. 10 is the mixture of propagating and standing signals. For example, two periods of standing signals are December 1986 in Fig. 10a and October 1997 in Fig. 10c. This variety suggests that MJOS(x, t) includes a range of intraseasonal variability types.
Finally, it should be noted that MJOS(x, t) variability can be present in any season, as indicated by the case studies in Figs. 8–11. The mean MJOSA(t) amplitude has a slight seasonal cycle, as described below in section 5b, where the boreal summer variability is on average slightly weaker than boreal winter variability. With these items in mind, relatively strong activity can be seen outside boreal winter, for example, during late May and early June 1987 (Fig. 11a), as well as during May 1998 (Fig. 11c). It would be interesting to do a more thorough examination of the seasonal changes in MJOS(x, t). Also, it may be advantageous to explore definitions of MJOS(x, t) that account for asymmetry about the equator, which could potentially be accomplished from a moving-heat-source point of view (Biello and Majda 2005) or using MJO skeleton theory itself (Thual et al. 2014b). With the current techniques, the initial results shown here suggest that the structures of intraseasonal variability—in both boreal summer and winter—have discernable projections onto theoretical, coupled convection–circulation structures that are symmetric about the equator.
To quantify the variability of MJOS(x, t), we present several statistics that summarize both its spatial and temporal characteristics.
The wavenumber–frequency power spectrum of MJOS(x, t) is shown in Fig. 12. The largest peak is associated with eastward-propagating signals of zonal wavenumber k = 1 and oscillation periods in the range of 30–90 days. This is consistent with the general view of MJO variability and with the dominance of eastward-propagating intraseasonal signals in the case studies of Figs. 8 and 10. At the same time, a small amount of westward power is also seen in Fig. 12, which is consistent with the occurrence of standing oscillations in Figs. 10a,c. Another feature in Fig. 12 is a weak tail of higher-frequency eastward power, extending along a line from roughly k = 1 and ω = 0.04 cpd to roughly k = 3 and ω = 0.1 cpd, associated with propagation speeds in a range of roughly 12–23 m s−1. This tail is likely associated with the fast signals that appear mostly over the eastern Pacific in Figs. 8 and 10, which are most likely associated with either convectively coupled Kelvin waves or circumnavigating MJO signals (Milliff and Madden 1996). In short, Fig. 12 indicates the statistical dominance of intraseasonal, eastward-propagating, wavenumber-1 signals in MJOS(x, t), although a range of other signals also contributes.
To illustrate the regional variability of MJOS(x, t), its variance is shown for each longitude x in Fig. 13. While nonzero variability exists at all latitudes, it is greatest from 45°E to 180°, in the Indian and western Pacific Ocean regions, as is commonly identified as the regions of greatest MJO variability (Kiladis et al. 2005). In these regions, the variance exhibits a strong seasonal cycle, with the greatest values occurring during DJF, the boreal winter season, and also during March–May (MAM; not shown). In the Western Hemisphere, on the other hand, very little change is seen from season to season.
Finally, several statistics of the MJOSA(t) time series are shown in Figs. 14–16: histogram, autocorrelation function, and seasonal cycle, respectively. The histogram appears to decay roughly exponentially for large MJOSA values [as also indicated by a log–linear plot (not shown)]. The autocorrelation function has a somewhat rapid decay time scale of roughly 12 days. For longer lags, the autocorrelation function oscillates annually, and this behavior persists out to lags of at least 16 yr (not shown). Since MJOSA(t) contains no phase information, there is no oscillation of correlation for intraseasonal lags. Finally, in Fig. 16 the seasonal cycle is displayed in both the mean of MJOSA(t) (middle curve) and the standard deviation, which is indicated by the separation between the curves. July–September is the period with the smallest mean and standard deviation, whereas February–April is the period with the largest mean and standard deviation.
In summary, these statistics are consistent with the generally accepted behavior of MJO activity, which suggests that MJOS(x, t) identifies MJO variability and its regional and seasonal variations.
c. Comparison with RMM and OMI
The time series of MJOSA(t) is highly correlated with other MJO indices, despite the different methods used in each case. For example, for the 1992–93 case in Fig. 9, the pattern correlation with the RMM index of Wheeler and Hendon (2004) is 0.96, and the pattern correlation with the OMI of Kiladis et al. (2014) is 0.89. For comparison, the pattern correlation of RMM and OMI is 0.91. These results are somewhat similar for the three cases in Fig. 11 as well. The slightly higher correlation between RMM and MJOSA is perhaps due to one of the commonalities in their methods that differs from OMI; for instance, both RMM and MJOSA are real-time indices whereas OMI employs some temporal filtering, and both RMM and MJOSA are based on multiple variables whereas OMI uses OLR only. In any case, all three indices are highly correlated. At the same time, it should be noted that a high statistical correlation can occur even when two indices differ in their detailed characterization of some individual events (Straub 2013; Kiladis et al. 2014).
The statistics in Figs. 14–16 can also be computed for the RMM index and OMI (not shown). Overall, the broad features are similar in all cases, despite the different methods used for each index. One slight difference appears in the seasonal cycles, where MJOSA and RMM reach a maximum in February–April whereas OMI reaches a maximum in December–February. If MJOSA is instead computed using an orthogonal projection as described in appendix B, it also reaches a maximum in December–February (not shown).
d. Parameter sensitivity studies
Several sensitivity studies demonstrate that the methods of this paper are robust. The parameter has been tested at values of 0.8 and 0.95, and the parameters Γ = 1 and have been tested at values of 0.5 and 2.0, and the pattern correlation of the modified results with the standard results in Fig. 8 are greater than 0.99 (not shown). Similar changes in the parameter HOLR lead to pattern correlations of 0.97, as reported in Table 3 along with tests of even larger modifications of HOLR. Finally, to test the influence of different spatial scales, the selection of zonal wavenumbers was changed to k = 1–2 or k = 1–4, which leads to pattern correlations of 0.9. These tests lend confidence in the robustness of the present methods and results.
e. Subsets of the variables K, R, Q, and A
As another type of sensitivity study, we now evaluate the contributions of each variable (K, R, Q, and A) in shaping the combined signal of MJOS(x, t). This is analogous to previous evaluations of the contributions of U850, U200, and OLR to the RMM index (Wheeler and Hendon 2004; Straub 2013). The method used here is to compute the projection (23) with only a subset of the four components (K, R, Q, and A). The resulting MJOS(x, t) is then compared with the standard MJOS(x, t) by computing the pattern correlation. All possible subsets of variables are considered: dropping one variable, dropping two variables, or keeping only a single variable. The results are shown in Table 4.
If only one variable is dropped, does the signal MJOS(x, t) change much? If R or A is dropped, slight changes are seen; whereas if K or Q is dropped, almost no change is seen. This suggests, preliminarily, that K and Q offer little additional information beyond that of the other variables (although subtleties can arise with such statements; see paragraphs below). From a physical point of view, this perhaps indicates a high degree of coupling among the convective activity and the K and R circulations.
At the other extreme, one can compare MJOS(x, t) with a single variable. When compared with Q(x, t), the pattern correlation is very low—only 0.18. This and other entries in Table 4 suggest that the Q formulation in (21) may be inadequate. However, when an orthogonal projection method is used (see appendix B), Q is seen to play a leading role. This suggests that the standard L2 projection in (23) may be inadequate. The orthogonal projection of (25) is more in line with MJO skeleton theory.
When MJOS(x, t) is compared with one of K, R, or A, the pattern correlation is in the range 0.74–0.76. This is noteworthy for two reasons. First, the narrow range of 0.74–0.76 indicates that each of K, R, and A is equally correlated with MJOS(x, t), which, in some sense, indicates that no single variable is “more important” to the definition of MJOS(x, t) than the other variables. Second, pattern correlations of 0.74–0.76 are quite far from a 1.00 perfect correlation. In other words, no single variable is perfectly correlated with MJOS(x, t), which indicates that MJOS(x, t) is truly a multivariate signal.
6. Relation to the wave response to a moving heat source
As an alternative method for obtaining a theoretical MJO structure, one can find the wave response to a prescribed, moving heat source (Chao 1987; Majda and Biello 2004; Biello and Majda 2005). In fact, the MJO skeleton eigenvector can be interpreted from this perspective if some additional assumptions are made.
To summarize this alternative perspective, consider a moving heat source:
The Kelvin wave response can be computed using (3a):
where it is understood that the real part of these complex quantities is desired. Similarly, the Rossby wave response is
and the lower-tropospheric moisture response is
The solutions (29)–(31) offer a simple perspective on the MJO skeleton eigenvector : its components of , , and are precisely the response (29)–(31) to a prescribed moving heat source (27). This can be seen by verifying that the values in Table 2 satisfy the relations (29)–(31). While this analogy offers a simple interpretation of , it is not an equivalence of the two perspectives, since the moving heat source requires a prescribed ωMJO(k), whereas the MJO skeleton model predicts ωMJO(k).
To identify the MJO skeleton in observational data, a method has been presented using theoretical wave structures. Instead of defining an empirical MJO structure (i.e., instead of using EOFs), a theoretical MJO structure was defined using equatorial wave theory. The theoretical definition incorporates a wide variety of data: multiple variables (wind, geopotential height, water vapor, and OLR), multiple pressure levels (850 and 200 hPa), and multiple latitudes (both equatorial and off equatorial). Equatorial wave theory provides a systematic way to combine this variety of data, and it allows, among other things, a distinction between the Kelvin and Rossby components of the MJO.
Using the MJOS(x, t) signal, one can explore the following question: To what extent is the structure of an MJO event consistent with the theoretical structure of a coupled convection–circulation anomaly? In several case studies, considerable amounts of consistency were seen. Furthermore, an alternative simple interpretation was also provided for the MJO skeleton structure: it is related to the wave response to a moving heat source. Consequently, the methods and results here are also applicable to a more general class of theories.
Several statistical measures were used to illustrate the spatial and temporal variability of the MJOS(x, t) signal. For example, it was seen that intraseasonal time scales dominate the MJOS(x, t) signal, without the use of an explicit temporal filter. Instead, through the combination of a spatial filter and a projection onto equatorial wave structures, the data are automatically, implicitly filtered onto intraseasonal time scales. As a result, the methods here can be used to filter data in real time.
A suite of studies was used to demonstrate that MJOS(x, t) is truly a multivariate signal. For example, in studies with subsets of the four variables (K, R, Q, and A), it was shown that MJOS(x, t) is only moderately correlated with any single variable (K, R, Q, or A), and each of the variables K, R, and A is equally correlated with MJOS(x, t). The role of moisture Q was seen to be key, provided an orthogonal projection is used ( appendix B), based on the MJO skeleton model’s conserved energy.
While the methods here were illustrated using observational and reanalysis data, they could also be applied to GCM data. For example, the methods could be used to evaluate the abilities of GCMs to simulate the MJO and its proper structure (Lin et al. 2006; Kim et al. 2009; Hung et al. 2013). Since the MJO structures here are defined theoretically, they provide an independent standard for this application.
Finally, we note that the MJOS(x, t) signal could potentially be used as an MJO “index.” Furthermore, the methods here could be combined with statistical methods such as EOFs or nonlinear Laplacian spectral analysis (Giannakis and Majda 2012), or the methods here could be modified to include the effects of a climatological mean state with zonal variations or meridional variations (Ogrosky and Stechmann 2014, manuscript submitted to Quart. J. Roy. Meteor. Soc.; Thual et al. 2014b). The authors are currently pursuing some of these directions and will present the results elsewhere in the near future.
The authors thank G. Kiladis for helpful discussion on the datasets, and they thank G. Kiladis and an anonymous reviewer for helpful comments on the manuscript. The research of SNS is partially supported by ONR MURI Grant N00014-12-1-0912 and by ONR Young Investigator Award N00014-12-1-0744. The research of AJM is partially supported by ONR MURI Grant N00014-12-1-0912.
Derivation of K, R, Q, and A Variables and Evolution Equations
In this appendix, some selected aspects are presented for the derivations of (3) from (2), including the derivations of the variables K and R. A main part of this derivation is the progression from primitive variables (u, θ) to characteristic variables (r, l) to equatorial wave variables (K, R1, R2, …). For more information, see Gill (1980, 1982), Majda (2003), Biello and Majda (2006), and Majda and Stechmann (2009).
a. Characteristic variables l and r
The equatorial long-wave equations can be written naturally in terms of characteristic variables r and l:
where L− and L+ are the raising and lowering operators, respectively:
which operate on parabolic cylinder functions as
Note that the meridional velocity equation has become (A6), which describes the meridional geostrophic balance in terms of the characteristic variables.
b. Defining K, R1, R2, R3, …
Next, we pass from characteristic variables (r, l) to the equatorial wave variables (K, R1, R2, R3, …). This change of variables is facilitated by an expansion in meridional basis functions ϕm(y), similar to (13)–(16):
Replacing the symbol r0 with K leads to (3a).
The second term in the expansion (A9) corresponds to the mixed Rossby–gravity (MRG) wave. With the equatorial long-wave scaling used here, the remnant of the MRG wave is not a dynamical equation but a diagnostic relation: , which follows from (A4). Since the antisymmetric convective activity a1 is not included in the simplest skeleton model with only K and R, the MRG remnant plays no role here.
Finally, in the infinite sum in (A9), each m corresponds to the mth Rossby wave. The sum of all types of Rossby waves can be written as an expansion in parabolic cylinder functions as
These equations can be written in terms of meridional modes as
To get an equation in terms of Rm alone, υm can be eliminated to yield
c. Recovering u and θ from K, R1, R2, R3, …
The variables u, υ, and θ are recovered in a two-step process. First, from K and Rm, m = 1, 2, 3, …, one can recover the characteristic variables l and r using meridional geostrophic balance (A6) and the definition from (A11). Second, from r and l, one can recover u and θ using (A1). With the meridional mode truncations of this paper described above, the resulting formulas are (Majda 2003; Biello and Majda 2006)
where υ is recovered from (A15) and where the standard meridional mode truncations of this paper were applied as described above.
d. Moisture dynamics
The moisture q(x, y, t) can be expanded as , where the sum is truncated here at q0(x, t), which is relabeled as Q(x, t) in (3).
The dynamics of Q(x, t) in (3c) can be derived from the dynamics of q(x, y, t) in (2d) in the following way. Each term in (2d) can be expanded in parabolic cylinder functions and then projected onto ϕ0 in order to find the evolution equation for q0(x, t). First, ux can be expanded in parabolic cylinder functions using (A17). Next, to handle the υy term, the operator ∂y can be written in terms of the ladder operators as , and υm can be written in terms of Rm and the convective activity using (A15) or (A18). Finally, a projection onto ϕ0 leads to (3d).
e. Convective activity dynamics
The convective activity a(x, y, t) can be expanded as , where the sum is truncated here at a0(x, t), which is relabeled as A(x, t) in (3). The dynamics of A(x, t) in (3d) can be derived from the dynamics of a(x, y, t) in (2e) by expanding q and a in the basis and projecting (2e) onto ϕ0. Assuming that q and a are truncated at the q0 and a0 terms, (3d) is obtained. Note that the projection procedure leads to the integral , which takes the value , which for notational simplicity has been absorbed into a new definition of Γ in (3d) compared to the Γ in (2e).
Orthogonal Projection Method Based on MJO Skeleton Energy
Here, we present a second projection method that can be used instead of the projection in (23). One of the main differences is that (23) is not an orthogonal projection, whereas the second method here uses an orthogonal projection. This is achieved by replacing the standard L2 inner product of (23) with an inner product that is based on the MJO skeleton model’s conserved energy. As described below, by using the MJO skeleton energy, two consequences are (i) water vapor plays a leading role, as opposed to its lesser role in section 5, and (ii) fast waves do not project as strongly onto the MJOS(x, t) signal, as compared to their presence in section 5.
The energy equation of the nonlinear skeleton model (2) has the form
which holds when the time-independent sources of heat and moisture are equal: s = sθ = sq. Note that this energy is a convex function of u, θ, q, and a. Upon using a meridional truncation to move from (2) to (3), and linearizing about a state of radiative–convective equilibrium, the conserved energy of the linearized system is
where in radiative–convective equilibrium.
The energy in (B2) can be used to define a natural inner product for the skeleton model. To do this elegantly, one can put the MJO skeleton model into symmetric hyperbolic form, following Friedrichs and Lax (1971). With this approach, the key quantity is the Hessian matrix of the conserved convex energy:
where the Ui are the elements of the state vector U = (K, R, Q, A)T. Explicitly, the matrix is computed from (B2):
For the parameter values used throughout this paper (see Table 1), this matrix takes the numerical values
where the units are nondimensional. Since E is convex, is a positive definite matrix. Furthermore, defines an inner product, and the linear eigenvectors are orthogonal with respect to this inner product. For example, if eMJO(k) and eK(k) are the MJO and dry Kelvin eigenvectors, respectively, for the same wavenumber k, then . Moreover, each eigenvector can be normalized in a natural way via this inner product using, for example, .
As an example of MJOS(x, t) defined through the orthogonal projection (B6), the TOGA COARE case of 1992–93 is shown in Fig. B1. This figure should be compared with Fig. 8, which was computed with the standard L2 projection in (23). The MJO events during January are prominent for both projection methods. One difference is that the rapidly propagating eastward signals in Fig. 8 (e.g., in the region 90°–120°W, from late December to early February) are absent from Fig. B1. The orthogonal projection method isolates the MJO eigenmode only, whereas the standard L2 projection can contain contributions from other modes, such as the fast dry Kelvin mode.
Another significant difference is in the comparison of subsets of variables (as in section 5e and Table 4). With the orthogonal projection method, Q alone produces a pattern correlation of 0.79, and the pair of Q and A together produce a pattern correlation of 0.96. Hence, the pair of Q and A appears to be the most important subset in defining the MJO skeleton signal. This is quite different from the results in section 5e, where Q was seen to play a minor role. This suggests that the Q definition in (21) is appropriately chosen (in contrast to the suggestions of section 5e).
In terms of the time series of the amplitude MJOSA(t) (not shown), the results with the orthogonal inner product are highly correlated with the L2 inner product results. For example, for the TOGA COARE case, the two MJOSA(t) time series have a pattern correlation of 0.95. For the three cases in Fig. 11, the pattern correlations are in the range of 0.95–0.97. Nevertheless, the spatial variability of MJOS(x, t) shows noticeable differences in the appearance of fast wave signals, as mentioned above, and a high statistical correlation can occur even when two indices differ in their detailed characterization of some individual MJO events (Straub 2013; Kiladis et al. 2014).
The two projection methods in (23) and (25) each have their advantages and disadvantages. The L2 projection was given more emphasis here due to its simplicity and due to its potential applicability for the moving-heat-source perspective on the MJO (section 6). It is perhaps the simplest method to use in more general settings, for example, with a zonally varying base state (Ogrosky and Stechmann 2014, manuscript submitted to Quart. J. Roy. Meteor. Soc.) or a meridionally varying base state (Thual et al. 2014b). On the other hand, the orthogonal projection is desirable because it is more in line with the MJO skeleton theory, it provides the mathematical advantages of orthogonality, and it better filters out contributions from other modes of variability.