On the role of oceanic entrainment temperature ( T e ) in decadal changes of El Niño / Southern Oscillation

The role of decadal changes in ocean thermal structure in modulating El Ni ño/Southern Oscillation (ENSO) properties was examined using a hybrid coupled model (HCM), consisting of a statistical atmospheric model and an oceanic general circulation model (OGCM) with an explicitly embedded empirical parameterization for the temperature of subsurface water entrained into the mixed layer (Te), which was constructed via an EOF analysis of model-based historical data. Using the empirical Te models constructed from two subperiods, 1963–1979 ( T 63−79 e ) and 1980–1996 ( T 80−96 e ), the coupled system exhibits striking different properties of interannual variability, including oscillation periods and the propagation characteristic of sea surface temperature anomalies (SSTAs) along the equator. In the T 63−79 e run, the model features a 2–3 yr oscillation and a westward propagation of SSTAs along the equator, while in theT 80−96 e run, it is characterized by a 4–5 yr oscillation and an eastward propagation. Furthermore, a Lag Covariance Analysis (LCOA) was utilized to illustrate the leading physical processes responsible for decadal change in SST. It is shown that the change in the structure of Te acts to modulate the relative strength of the zonal advective and thermocline feedbacks in the coupled system, leading to changes in ENSO properties. Two additional sensitive experiments were conducted to further illustrate the respective roles of the changes in ocean mean states and in Te in modulating ENSO behaviors. These decadal changes in the simulated ENSO properties are consistent with the observed shift occurred in the late 1970s and a previous simulation performed with an intermediate coupled model (ICM) described in Zhang and Busalacchi (2005), indicating a dominant role Te plays in decadal ENSO changes. Correspondence to: J. Zhu (jszhu@mail.iap.ac.cn)

) and 1980-1996 (T 80−96   e   ), the coupled system exhibits striking different properties of interannual variability, including oscillation periods and the propagation characteristic of sea surface temperature anomalies (SSTAs) along the equator.In the T 63−79 e run, the model features a 2-3 yr oscillation and a westward propagation of SSTAs along the equator, while in the T 80−96 e run, it is characterized by a 4-5 yr oscillation and an eastward propagation.Furthermore, a Lag Covariance Analysis (LCOA) was utilized to illustrate the leading physical processes responsible for decadal change in SST.It is shown that the change in the structure of T e acts to modulate the relative strength of the zonal advective and thermocline feedbacks in the coupled system, leading to changes in ENSO properties.Two additional sensitive experiments were conducted to further illustrate the respective roles of the changes in ocean mean states and in T e in modulating ENSO behaviors.These decadal changes in the simulated ENSO properties are consistent with the observed shift occurred in the late 1970s and a previous simulation performed with an intermediate coupled model (ICM) described in Zhang and Busalacchi (2005), indicating a dominant role T e plays in decadal ENSO changes.

Introduction
A striking decadal shift of ENSO behaviors was observed in the tropical Pacific in the late 1970s, including its oscillation period and the way its SST anomalies propagate on the equator.For example, its dominant period increased from 2-3 years during 1960s/70s to 4-5 years during 1980s/90s (e.g., An and Wang, 2000).In addition, before the late 1970s, the ENSO-related SSTAs were initiated in the South American coast and then propagated westward along the equator (Rasmusson and Carpenter, 1982), while beginning with El Niño in 1982, SSTAs propagated eastward from the central basin or developed concurrently in the central and eastern Pacific (Wallace et al., 1998).
Many efforts have been devoted to understanding the decadal changes in ENSO; some specific mechanisms have been identified, including stochastic atmospheric wind forcing (Kirtman and Schopf, 1998;Flügel and Chang, 1999), the influences of extratropical processes in the atmosphere and/or in the ocean (Gu and Philander, 1997;Zhang et al., 1998;Wang and Liu, 2000), the changes in mean state of the tropical climate system on which El Niño evolves (Fedorov and Philander, 2000;Wang and An, 2001) and the nonlinear processes in the climate system (Timmermann and Jin, 2002;Duan and Mu, 2006).Zhang and Busalacchi (2005) identified T e (the temperature of subsurface water entrained into the mixed layer) as an important factor in decadal ENSO variability as observed in the late 1970s.Using an intermediate coupled model (ICM), they demonstrated that the change in subsurface ocean thermal structure was able to modulate ENSO properties.Zhang In this paper we first examine the role of T e in decadal ENSO changes using a different hybrid coupled model (HCM OGCM ), consisting of a statistical atmospheric model and an OGCM with an empirical parameterization for T e embedded.Then, we analyze the possible mechanisms of how decadal T e changes modulate ENSO properties.Two additional sensitivity experiments are further conducted to examine the relative roles of the changes in the ocean background state and T e in decadal ENSO variability.
The paper is arranged as follows.Section 2 briefly describes a hybrid coupled model (HCM) and experiment designs.The results from modeling experiments with the empirical T e models constructed from two different periods are shown in Sect.3. Section 4 examines the mechanisms by which decadal changes in ENSO properties can be modulated by those in T e .Section 5 presents results from two sensitivity experiments, illustrating the relative roles of the changes in the ocean background state vs. in T e in modulating ENSO behaviors.The paper is concluded in Sect.6.

Model descriptions
As illustrated in Fig. 1, we develop a HCM consisting of a statistical atmospheric wind stress anomaly model, an OGCM, and an embedded SST anomaly model with an empirical parameterization for T e .Informations between these submodels are exchanged once a day (Zhu et al., 2009).
The atmospheric model is a statistical one constructed from a SVD of the covariance matrix calculated from time series of monthly mean SST and wind stress (τ ) anomaly fields, as in Zhang and Zebiak (2004) and Zhang and Busalacchi (2005).The OGCM used here was firstly developed by Zhang and Endoh (1992), and has been used for ENSO prediction at the Institute of Atmospheric Physics, the Chinese Academy of Sciences (IAP/CAS) (Zhou and Zeng, 2001).The dynamics of the model are governed by primitive equations under hydrostatics and the Boussinesq approximation in σ -coordinates with a free surface.The model domain is confined within the tropical Pacific region (30 • S-30 • N, 121 • E-69 • W) with realistic land-sea boundaries and flat bottom.There are 14 vertical levels with resolution of 20 m in the upper 60 m and of 30 m between 60 m and 240 m depth.The model horizontal resolution was increased from 1 Zhang and Endoh (1992) to 0.5 • × 0.5 • (Fu et al.,  2005; hereafter IAP-TPOGCM), allowing for a proper depictions of equatorial waves.Detailed descriptions of the OGCM can be found in Zeng et al. (1991), Zhang and Endoh (1992) and Fu et al. (2005).
A SSTA model is embedded into the IAP-TPOGCM.Mathematically, the governing equation determining the evolution of interannual SST variability in the surface mixed layer can be written as (Zebiak and Cane, 1987;Zhang et al., 2005) where, T and T e are anomalies of SST and the temperature of subsurface water entrained into the mixed layer, respectively; T and Te are the prescribed seasonally varying mean SST and T e , which are obtained from observations and the OGCM run; ū and v are the prescribed seasonally varying mean zonal and meridional currents in the mixed layer, and w is the prescribed seasonally varying mean entrainment velocity at the base of mixed layer, which are all obtained from the OGCM run forced by prescribed winds from ERA40 (Uppala et al., 2005) anomaly fields; H is the depth of the mixed layer; 5 × 10 3 m 2 s −1 meridionally and 2.5 × 10 4 m 2 s −1 zonally) and κ v (10 −3 m 2 s −1 ) are the horizontal and vertical diffusion coefficients, respectivly; α is the thermal damping coefficient; and other variables are conventional.
In the embedded SSTA model, an empirical relation between T e and sea level (SL) anomalies has been constructed to parameterize T e (Zhang et al., 2005), which acts to improve ENSO simulations in the IAP-TPOGCM (Zhu et al., 2006).The empirical parameterization scheme is developed in two steps.First, an inverse modeling method is adopted to estimate T e anomalies which are estimated from the SSTA Eq. ( 1) with the observed SST fields and the OGCM modeled currents.Secondly, an EOF technique is adopted to construct a relationship between the inversed T e anomalies and the modeled SL anomalies in the OGCM, from which an empirical model for T e can be derived and used to calculate T e anomalies in the coupled ocean-atmosphere runs.Detailed descriptions about the T e model can be found in Zhang et al. (2005) and Zhu et al. (2006).
In addition, the coupled behaviors in our HCM depend on the so-called relative coupling coefficient (α τ ), i.e., the wind stress anomalies calculated from the statistical atmospheric model are multiplied by a scalar parameter before being used to drive ocean models (Barnett et al., 1993).Different values have been tested to get a sustainable oscillation for the HCM.It can be shown that in the HCM we use, the value of α τ does not change the coupled oscillation period and space structure significantly, which is similar to an ICM simulation (Zhang and Busalacchi, 2005).In all of the experiments shown below, the value of α τ has been chosen for the HCM to produce ENSO oscillation with a reasonable amplitude.

Experiment design
To examine the sensitivity of coupled behavior to subsurface temperature variability, two T e models are constructed separately from two different subperiods 1963-1979 and 1980-1996 (hereafter T 63−79 e and T 80−96 e , respectively).These two T e models are then utilized to parameterize the T e fields from a SL anomaly in the OGCM (Fig. 1), which is used for the SSTA calculations in the embedded coupled system, with the other model settings being exactly the same.The time cutoff chosen for the two subperiods (1963-1979 and 1980-1996) used to construct the two empirical T e models corresponds to pre-and post-decadal shift that occurred in the late 1970s.
Figure 2 shows the spatial pattern of the first EOF mode of T e derived from the periods 1963-1979 and 1980-1996, re-spectively.As discussed in Zhang and Busalacchi (2005) and Zhang and DeWitt (2006), T e exhibits large decadal changes before and after the late 1970s.In particular, during the former period, T e anomalies are relatively weak both in the eastern equatorial Pacific and in the off-equatorial tropical North Pacific along 10 • N between 140 • E and the date line.After the late 1970s, the T e anomalies are significantly enhanced in both regions.It is expected that there can be large differences in the amplitude of T e anomalies calculated from a given SL anomaly using T e models constructed from these two subperiods.run, the relative coupling coefficient,α τ , is set to 1.15 for the HCM to sustain an interannual variability.Figure 3 shows the simulated anomalies of SST, zonal wind stress and sea level along the equator.The coupled model captures interannual oscillations with about 3-year period.Another striking feature is that the simulated interan-nual variability is characterized by a pronounced westward propagation in SST and wind variations on the equator.That is, SST anomalies appear first in the eastern equatorial Pacific, and then migrate westward along the equator, which is accompanied by wind anomalies which also show coherent westward propagation.When SST anomalies arrive in the central basin, large wind anomalies are induced in the central region west of SST anomalies.The commonly used Niño index more clearly exhibits a phase lag in SST variations, with SST anomalies in the Niño 1+2 regions leading those in Niño 3,4 region.Evidently, the simulated ENSO properties (e.g., coherent westward propagation of SST anomalies and about 3-year oscillation periods) in the T 63−79 e run are consistent with the development and evolution of El Niño events observed before the late 1970s (Rassmusson and Carpenter, 1982).Figure 4 illustrates the results for the T 80−96 e run.Here, the relative coupling coefficient,α τ , is set to 1.2 for the HCM to sustain an interannual oscillation.It is seen that dramatic changes emerge in ENSO properties simulated from the coupled system, including the dominant period and the way El Niño events evolve.The coupled system is now characterized by interannual oscillations with a 4-to 5-year period, with dramatic changes in the phase propagation of coupled atmosphere-ocean anomalies.In particular, SST anomalies first appear in the central basin near the date line and then propagate eastward.These SST anomalies are coupled with wind anomalies which also progress eastward on the equator.These ENSO properties simulated in the T 80−96 e run (e.g., the dominant period and the eastward propagation) are consistent with those observed after the late 1970s.It is evident that the HCM with the two T e models constructed before and after the late 1970s qualitatively capture the changes in ENSO properties as seen in nature.

The mechanisms for decadal SST variations
To understand the decadal changes in ENSO from the two runs, the tendency terms for SST evolution are analyzed (Eq. 1 shown above).As demonstrated by An and Wang (2000) and Capotondi et al. (2006), if a tendency term on the right side of Eq. ( 1) is in phase with a SST anomaly, this term favours ENSO growth; if, on the other hand, it is in quadrature with a SST anomaly with leading by a quarter cycle, this tendency term tends to promote a phase transition.As previous studies suggested (Jin and An, 1999;An and Wang, 2000;Wang and McPhaden, 2000;Kang et al., 2001;Capotondi et al., 2006), among various budget terms of the SST tendency equation, two terms are essential for ENSO dynamics: (1) the zonal advection of surface temperature by anomalous current (the second term on the right in Eq. ( 1); referred to as the zonal advective feedback), and (2) the vertical advection of subsurface temperature anomaly by mean upwelling (the sixth term on the right side in Eq. 1; referred to as the thermocline feedback).Furthermore, An (2005) explored the relative roles of these two feedback terms in determining timescale of the tropical coupled ocean-atmosphere system using an eigen analysis.It is shown that the zonal advective feedback favors a high frequency coupled mode, while the thermocline feedback favors a low frequency coupled mode.
In this section, we will analyze the two terms using a LCOA method (Lau et al., 1992;Kang et al., 2001) to identify principle processes that cause differences in the simulated ENSO properties with different T e models.In this LCOA method, the evolutionary feature of a quantity (S) associated with ENSO is described by a lag covariance between the Niño3 SST anomalies and S with a time lag.To measure an actual magnitude of the quantity (S) related to a reference time series, the lag covariance is divided by the standard deviation of Niño3 SST anomalies.The LCOA can, therefore, be calculated by the following equation (Kang et al., 2001): where NT indicates Niño3 SST anomalies; σ (NT) is the standard deviation of NT.The zero lag in the LCOA analysis is equivalent to the mature phase of ENSO.
Figures 5 and 6 show the LCOAs between the Niño3 SST anomalies and those two terms from the T 63−79 e and T 80−96 e runs.In the T 63−79 e run (Fig. 5), the vertical advection term (Fig. 5b) is more or less in phase with the SST anomaly, whereas the zonal advective term (Fig. 5a) is in a quardrature leading with SST.Following An and Wang (2000) and Capotondi et al. (2006), it is thus evident that the thermocline feedback term (Fig. 5b) mostly favors the growth of ENSO, while the zonal advective feedback (Fig. 5a) mostly favors the transition of ENSO.These results in the T 63−79 e run (Fig. 5) are similar to those in An and Wang (2000).
The results for the T 80−96 e run (Fig. 6) differ greatly from those for the T 63−79 e run (Fig. 5).The thermocline feedback term (Fig. 6b) is generally positive during and before the mature phase of an El Niño over the central and eastern Pacific, but becomes negative afterwards.This implies that the thermocline feedback favors both the growth and transition of ENSO.These results are basically similar to those in Kang et al. (2001) who analyzed the NCEP ocean reanalysis data for 16 years (from 1980 to 1995).Nevertheless, there are also some differences between our analysis and that in Kang et al. (2001).For example, in Kang et al. (2001) the zonal advective feedback term shows small positive amplitude before the mature phase of an El Niño, which may favor the transition of ENSO.In our analysis, the contribution of the zonal advective term is small and tends to be negative (Fig. 6a).Thus the dominant thermocline feedback will contribute to a lower frequency oscillation as demonstrated in An (2005).
Comparing the T 63−79 e run and the T 80−96 e run (Figs.5-6), there are significant differences in the roles the zonal advection and vertical advection terms play in the transition period of ENSO cycles.In the T 63−79 e run, the zonal advective feedback term plays a "transition" role (Fig. 5a), while in the T 80−96 e run, the thermocline feedback plays the role (Fig. 6b).Furthermore, the durations when positive values persist before the peak of an El Niño are different.In the T 63−79 e run (Fig. 5a), the "transition" periods are seen intensely within ∼9 months before the mature phase of an El Niño, while in the T 80−96 e run (Fig. 6b), the "transition" periods are seen to extend broadly to at least 18 months before the mature phase of an El Niño.This can explain why short (about 2-3 years) and long (4-5 years) oscillation periods are seen in the two HCM runs which are correspondent to those observed before and after the late-1970s shift.

The roles of oceanic background state vs. T e : two additional sensitivity experiments
In the above experiments, the empirical T 63−79 e and T 80−96 e models are derived from two different periods, which are then utilized to explicitly incorporate a decadal change in subsurface thermal structure into the coupled models.In deriving the two T e models, an "inverse" modeling approach is taken to get T e fields from the SST anomaly equation, using observed (SST fields and their tendency) and modeled data (mean and anomaly currents).Note that the ocean mean states ( ū, v, w, T ,. . . . . . ) used for the inverse modeling are the long-term mean climatology obtained from the OGCM simulation during the entire period 1963-1996 (O  63−96 ).As some recent studies indicate that the decadal shift in ENSO can be attributed to changes in background ocean state (e.g., Fedorov and Philander, 2000;Wang and An, 2001), it is interesting to see if the decadal changes in ENSO seen in our experiments can also be attributed to those in the background ocean fields, instead in T e ., are performed, in which a T e model constructed from the period 1963-1979 (1980-1996) is used with an ocean mean state ( ū, v, w, T ,. . . . . . ) which is obtained from the OGCM simulation for the periods 1980-1996 (1963-1979) model run, respectively.
In the O 80−96 T 63−79 e run, the relative coupling coefficient, α τ , is set to 1.05 for the HCM to sustain an interannual oscillation.Figure 7 shows the simulated anomalies of SST, zonal wind stress, and sea level along the equator from this run.Similar to the T 63−79 e run (Fig. 3), the coupled variability is characterized by interannual oscillations with a 2-3 year period.Moreover, the simulated SST variations are characterized by a westward propagation along the equator (Fig. 7a).In contrast, simulations from the O 63−79 T 80−96 e run in which the relative coupling coefficient (α τ ) is set to 1.03 (Fig. 8) are very similar to those in the T 80−96 e run (Fig. 4).The coupled variability is characterized by interannual oscillations with a 4-to 5-year period; the coherent atmospheric and oceanic anomalies originate from the central basin and then propagate eastward along the equator (Fig. 8).These two sensitivity experiments suggest that decadal ENSO variability in the embedded HCM is mostly attributable to the alteration in T e structure, while the role of the changes in the ocean mean states is less dominant.

Conclusion and discussion
In this paper, a hybrid coupled model (HCM OGCM ), consisting of a statistical atmospheric model and an OGCM with an empirical parameterization for T e , is used to examine the effect of interdecadal changes in the structure of subsurface temperature on ENSO properties in the tropical Pacific.The decadal changes in upper-ocean temperature observed in the late 1970s are incorporated into the coupled model via an empirical T e parameterization which is constructed from two subperiods corresponding to the pre-and post-climate shift.It turns out that the modeled ENSO characteristics obtained using the two sets of T e specifications are qualitatively consistent with those observed in nature.In particular, using the T e model constructed from the period 1963-1979, the system is characterized by a 2-3 yr oscillation and westward propagation of SST anomalies along the equator; on the other hand, the system features a 4-5 yr period oscillation and an eastward phase propagation along the equator when using the T e model constructed from the 1980-1996 period.These results suggest that the interdecadal changes in T e may play an important role in decadal ENSO changes in nature.
The LCOA is further applied to analyze leading advective processes in the SST tendency equation associated with the ENSO variability.The analysis indicates the way the zonal advection and thermocline feedbacks contribute to ENSO evolution is different when using the two T e models constructed during the periods 1963-1979 and 1980-1996.That is, when using the T 63−79 e model (before the late-1970s), the thermocline feedback term mostly favors the growth of ENSO, while the zonal advective feedback term mostly favors the transition of ENSO.When using the T 80−96 e model (after the late-1970s), the thermocline feedback term favors both the growth and transition of ENSO, while the zonal advective feedback term shows a small negative contribution to the transition of ENSO, thus contributing to lower frequency oscillation.Furthermore, the "transition" terms (the zonal advective feedback term in the T 63−79 e run and the thermocline feedback term in the T 80−96 e run) also behave differently.In the T 63−96 e run, the "transition" is seen to be intense within very short time, while in the T 80−79 e run, the "transition" periods are seen to extend broadly to at least 18 months before the mature phase of an El Niño.These differences may explain why higher (lower) frequency oscillations are observed before (after) the late-1970s shift.
In order to confirm the importance of T e in decadal ENSO variability, two additional sensitivity experiments are conducted to isolate the relative effects of the changes in the background ocean state vs. in T e ; these two factors corresponding to the periods 1963-1979/1980-1996  respectively, taken into account in the coupled simulations.It is seen that the effect of the changes in the ocean background is not essential to decadal changes in ENSO.Nevertheless, this result may partially come from the fact that the HCM used here has some potential weakness.For example, in our HCM OGCM , the mean state of the atmospheric model is prescribed to be seasonally varying, as well as some mean fields in SST equation.The effects of ENSO-related interannual coupling on the climatological fields are excluded in the HCM-based simulations.Thus the role that the changes in the background ocean state play in the detected shift in ENSO properties may be underestimated.
It should be pointed out that, as demonstrated by some previous studies, the decadal changes in the surface winds and the associated ocean surface layer dynamics are uncertain in the equatorial Pacific when using current observational stress products (Wittenberg, 2004).Therefore, further experiments taking hybrid coupled modeling approach are needed to examine the sensitivities to different wind stress products, including the Florida State University (FSU) subjective analysis, the NCEP/NCAR reanalysis and others.
It is also desirable to apply the findings of this HCMbased modeling study to CGCM simulations in order to understand ENSO behaviors.For example, many CGCMs simulate ENSO cycles with a 2-3 year oscillation period, which is shorter than that observed (e.g., referred to as the SSTmode; Guilyardi, 2006).Some studies attributed it to the existence of biases in surface wind stress simulated in AGCMs (An and Wang, 2000;Capotondi et al., 2006).From the present studies, we can see that the biases in OGCM can also contribute to ENSO properties.In particular, vertical mixing schemes currently used in OGCMs likely underestimate the connection between surface and sub-surface layers, which likely lead to an underestimation of the thermocline feedback in coupled simulations, and may result in a relatively short ENSO period.
In addition, some questions remain unknown.For example, what specific processes and mechanisms are involved with the changes in the entrainment temperature that can contribute to a selection of dominant feedbacks (i.e., zonal advection feedback vs. thermocline feedback)?More studies are needed to address them.

Fig. 1 AFig. 1 .
Fig.1 A schematic diagram illustrating a hybrid coupled model consisting of an OGCM, an embedded SST anomaly model with an empirical T e parameterization, and a statistical atmospheric wind stress anomaly model.Informations between these submodels are exchanged once a day.

Fig. 2 Fig. 2 .
Fig.2 Spatial patterns of the first EOF mode for T e derived during the periods (a)

Fig. 3 .
Fig.3 Longitude-time sections of simulated anomalies along the equator for (a) SST,
examine the relative roles of the changes in the ocean mean states and in T e in modulating ENSO.The constructions of the two experiments, O 80−96 T 63−79 e and O 63−79 T 80−96 e , attempt to take into account the changes in the mean conditions in two steps as follows.Firstly, the ocean mean climatological state ( ū, v, w, T ,. . . . . . ) is obtained differently from the periods 1963-1979 (O 63−79 ) and 1980-1996 (O 80−96 ), which are referred to as two different background ocean states.Secondly, when deriving T e fields from the "inverse" modeling, the ocean mean states ( ū, v, w, T ,. . . . . . ) used in the SST anomaly Eq. (1) are specified differently for the two different subperiods for 1963-1979 (O 63−79 ) and 1980-1996 (O 80−96 ), respectively.To examine the relative roles of the changes in the ocean mean states vs. in T e , two "CROSS" experiments, i.e.O 80−96 T 63−79 e and O 63−79 T 80−96 e

Fig. 7 .
Fig.7The same as in Fig.3, but for the

Fig. 8 .
Fig.8The same as in Fig.3, but for the