Weak ENSO Asymmetry Due to Weak Nonlinear Air-Sea Interaction in CMIP5 Climate Models

The El Niño-Southern Oscillation (ENSO) is the most prominent interannual variability in the climate system. Its climatic impacts are global, so it is very import to understand the simulation and projection of ENSO variability in climate models (Cane, 1983; Philander, 1983; Cai et al., 2014).

ENSO is not a linear system (Jin et al., 1994). The cold SST anomaly center during La Niña usually emerges farther west than the warm SST anomaly center during El Niño (An, 2009). El Niño is usually stronger in amplitude and lasts for a shorter duration than La Niña. This has been referred to as the ENSO asymmetry (Burgers and Stephenson, 1999). The reasons for the asymmetry between the two phases remain elusive, but many studies suggest that it is a consequence of nonlinearity in the dynamics of the ocean and atmosphere (Jin et al., 2003; An and Jin, 2004).

State-of-the-art climate models suffer from large errors in simulating the tropical Pacific mean state and ENSO variability. The pronounced errors in the mean state, which have existed for almost two decades, are the far west cold tongue and the double ITCZ problem (Mechoso et al., 1995; Li and Xie, 2012; Wang et al., 2014). These large errors in the simulation of the Pacific cold tongue and mean ITCZ limit the skill of coupled models in simulating ENSO and its teleconnection (Latif et al., 2001; Guilyardi, 2006; Ham and Kug, 2012). Recent studies have reported some biases in ENSO performance in CMIP5 (Coupled Model Intercomparison Project, Phase 5) climate models. (Zhang and Jin, 2012) noted narrow biases in ENSO width. (Kim et al., 2014) revealed weak ocean dynamic responses to atmospheric changes in the CMIP5 historical experiment, characterized by weak zonal advective feedback, thermocline feedback and thermodynamics damping. (Zhang and Sun, 2014) found a weak ENSO asymmetry in the CMIP5 "piControl" simulation, and their AMIP runs showed that the weak ENSO asymmetry originates at least in part from atmospheric processes. Note, however, that they only examined 14 CMIP5 models and only a single run was analyzed in each model.

ENSO asymmetry could manifest itself as a rectified change in the background state on the decadal time scale (An et al., 2005). The asymmetry between El Niño and La Niña can lead to a collective effect on the tropical Pacific, with a spatial pattern resembling the anomaly in El Niño (Burgers and Stephenson, 1999; Rodgers et al., 2004; Sun et al., 2009; Sun et al., 2013). (Sun et al., 2013) demonstrated a weak ENSO asymmetry and weak residual time mean effect of ENSO in 19 no-flux-adjustment CMIP3 models. Therefore, the bias in ENSO asymmetry may potentially contribute to the tropical Pacific mean state bias, although the residual effect may be a small part as only a rectification. Given the collective time-mean effect of ENSO, understanding the causes and consequences of ENSO asymmetry may be critical to understanding the decadal variability in the tropics and extratropics. Also, to fully capture the role of ENSO in the changing climate, climate models need to simulate the ENSO asymmetry well.

Previous attempts have been made to evaluate the ENSO asymmetry in climate models (An et al., 2005; Leloup et al., 2008, Zhang et al., 2009; Sun et al., 2013, Zhang and Sun, 2014). However, these studies tended to examine only a single run of the models concerned, or the number of models examined was limited. For example, (An et al., 2005) employed 10 models. (Leloup et al., 2008) noted the weak modeled ENSO asymmetry in the 20th century simulation in 23 CMIP3 models; but they only employed the first run of each model. (Zhang et al., 2009) noted an underestimation of ENSO asymmetry in the NCAR CCSM models and explored the causes. (Sun et al., 2013) made use of the CMIP3 archive and found that an underestimation of the asymmetry is a prevalent problem in the 20 ^th century simulation. However, the cause of the ENSO asymmetry was not well understood yet in these studies, and none of them has examined the simulation of the ENSO asymmetry in the 20th century in CMIP5 climate models. Accordingly, the present study evaluated the simulation of ENSO asymmetry in the 20th century simulation (the "historical" experiment) in CMIP5 climate models that were also included in IPCC AR5. Not only did this study examine a large set of CMIP5 models, but it also analyzed the available ensemble runs of individual models. Furthermore, the potential causes of the ENSO asymmetry biases associated with air-sea interaction biases over the tropical Pacific were investigated. This study may help us understanding the relationship between the mean state simulation and ENSO asymmetry simulation in climate models.

The models and observational data used in the study are described in section 2. The ENSO asymmetry, air-sea interaction over the tropical Pacific, and the possible causes of the biases are presented in sections 3 and 4. The paper concludes with a summary and further discussion in section 5.

The observational SST data used in this study were from the Hadley Centre Sea Ice and Sea Surface Temperature dataset, version 1 (Rayner et al., 2003). The SST field was built from in-situ and satellite observations, given on a 1^°× 1^° grid and available from 1871 to 2007. The surface wind stress data were obtained from the Simple Ocean Data Assimilation dataset (1900-2008) (Carton and Giese, 2008). The precipitation data were from the Global Precipitation Climatology Project (1979-2013) (Adler et al., 2003). The climate model data examined were the outputs of the `historical' experiment in CMIP5 (45 climate models, with 184 model runs in total). The last 50 years of the 20th century (1950-99) were the focus in this study, since the observational data are more reliable during this period. Noting that the precipitation data focused on the period 1979-99 only, since satellite observational data is only available from 1979. All the model data and observational datasets were monthly mean values and were interpolated to a uniform 1^°× 1^° horizontal grid for analysis. The models used in the study are listed in Fig. 1. Nine climate models were removed from the analysis because they simulated too weak an ENSO or an incorrect ENSO pattern based on the first EOF mode of the tropical Pacific SST (Fig. S1, available only online).

Figure 1.  (a) Frequency distribution (vertical axis) of monthly mean Niño3 SSTA (horizontal axis) in 184 runs of 45 CMIP5 models and the observation (HadISST). (b) PDF (blue curve) for the skewness of Niño3 SSTA in the 184 runs. The red vertical line indicates the observed value. The vertical blue line is the MME value. The short colored bars on the horizontal axis mark the skewness values of Niño3 SSTA in each individual run. (c) Average Niño3 skewness of all runs in individual models. The nine CMIP5 models without model numbers in front of the model names were removed from the analysis based on the first EOF mode of the tropical Pacific SST. The period focused upon is 1950-99.

In observation, El Niño is usually stronger in amplitude and lasts shorter in duration than La Niña; whereas, in climate models, ENSO events are almost symmetric in amplitude and duration. This can be readily seen from the histogram of the Niño3 SST anomaly (Niño3 SSTA) distribution over the last 50 years in the 20th century (Fig. 1a). The Niño3 SSTA varies from -2^°C to 3.5^°C and the cold anomalies occur more or last longer than the warm anomalies in the observation. In contrast, the Niño3 SSTA varies from -5^°C to 5^°C among models, and the distribution of cold and warm anomalies is almost symmetric. The skewness of Niño3 SSTA is typically used as the way to measure ENSO asymmetry (Burgers and Stephenson, 1999; Sun and Zhang, 2006; An et al., 2009). The multi-model ensemble runs (with 184 runs in total) allowed us to construct a probability density function (PDF) for the Niño3 skewness to depict the climate models' performances in ENSO asymmetry over the last 50 years of the 20th century (Fig. 1b). The results indicate that the CMIP5 climate models still generally underestimate the ENSO asymmetry. Most of the modes' skewness are located near zero, and the average skewness of the 184 model runs is 0.17——much weaker than the observed value (0.88). Close to 86% of the model runs (159 runs in number) have skewness between -0.75 and 0.5, while the remaining 25 runs (around 14% of the total number of runs) have skewness larger than 0.5. Furthermore, only seven runs simulate the Niño3 skewness larger than observed (Fig. 1b). The average of the model runs' Niño3 skewness in each model shows that only six models (MIROC5, M23, GFDL-ESM2M, MCESM1-CAM5, CESM-CAM5-1-FV2, CMCC-CM) can simulate the skewness larger than 0.5 (Fig. 1c).

Figure 2.  PDF analysis for the skewness of (a, d) Niño4 index, (b, e) Niño3.4 index and (c, f) Niño3 index in individual runs of (a-c) CMIP3 and (d-f) CMIP5 climate models. The time period is from 1950 to 1999. For the CMIP3 color schemes, refer to Sun et al. (2013, Table 1) and for the CMIP5 color schemes, refer to Fig. 1.

Figure 3.  Scatter plots between the standard deviation (STD) and skewness of Niño3 SSTA from 1950 to 1999. Each color represents a model, and different markers in the same color represent different runs of the same model. The black pentagram is the observational value and the hexagram is the averaged value of all model runs. The iap_fgoals1_0_g model (dark blue) in CMIP3 and the FIO-ESM model (light yellow) in CMIP5 were removed for the calculation of the correlation coefficient, as their STDs were large and their skewness negative. For the CMIP3 color schemes, refer to Sun et al. (2013, Table 1), and for the CMIP5 color schemes, refer to Fig. 1.

Considering there might be systematic pattern shifting bias in space in climate models, the PDF for the skewness of the Niño4 and Niño3.4 indices were also constructed. The weak ENSO asymmetry result is consistent with the Niño3 skewness analysis (Fig. 2). In the observation, the skewness of Niño3 (or Niño3.4) is positive and the skewness of Niño4 is negative. In the models, however, the average skewness of all the model runs is almost zero for all three indices. This weak modeled ENSO asymmetry bias exists from CMIP3 to CMIP5. Nevertheless, CMIP5 shows a slight improvement in the ENSO asymmetry, with the skewness of the three indices being closer to the observation.

Stronger ENSO activity (measured by the standard deviation of Niño3 SSTA) usually corresponds to a stronger ENSO asymmetry in climate models. A high correlation between the ENSO skewness and amplitude, with the correlation coefficient being 0.5 in CMIP3 and 0.53 in CMIP5, was found (Fig. 3). It is more evident in CMIP5 that strong ENSO activities are usually accompanied by strong ENSO asymmetry. CMIP5 shows a reduced diversity of ENSO activity compared with CMIP3. Weak ENSO asymmetry is mainly positively related to El Niño events among the CMIP5 models. The correlation coefficient between the skewness and El Niño events is 0.55 among all the model runs (Fig. 4a). An even more distinct relationship is found between the skewness and extreme El Niño events, with the correlation coefficient reaching 0.76 among all the model runs (Fig. 4b). On the contrary, there is no significant correlation between the skewness and La Niña events (especially the extreme ones) (Fig. 4). The Niño3 SSTAs are significantly underestimated in the El Niño events, especially in the extreme El Niño events. In contrast, the SSTAs in the La Niña phase are closer to the observation than in the El Niño phase (Fig. 4).

Figure 4.  (a) Scatterplots between the skewness of Niño3 SSTA and the average warm (positive) and cold (negative) Niño3 SSTA. The right (left) $r$ is the correlation coefficient between the skewness and the mean warm (cold) Niño3 SSTA. (b) As in (a) but for the relationships between the skewness and the extreme warm (maximum) and cold (minimum) Niño3 SSTA. In total, 94 runs were examined after performing quality control by removing nine climate models.

Figure 5.  The two atmospheric responses to SST changes over the central Pacific region (5$^\circ$S-5$^\circ$N, 180$^\circ$-140$^\circ$W) and their relationships. The variability of the surface zonal wind stress (units: N m$^-2$) (a) and the precipitation (units: Kg m$^-3$ s$^-1$) (b) responses to the variability of Niño3 SSTAs. (c) The relationship between the variability of the two atmospheric responses. (d-f) As in (a-c) but for the two atmospheric responses to SST changes and their relationships in the warm and cold phase, respectively. Niño3 SST was used to reveal the potential relationship between central Pacific air-sea coupling and ENSO. STD is used for the variability of anomalies. $R$ is the linear correlation coefficient among models. The period analyzed is from 1950 to 1999 for wind stress and 1979 to 1999 for precipitation. The period 1979 to 1999 for their correlation was chosen for consistency. The zonal wind stress and precipitation have been amplified by 10$^2$ N m$^-2$ and 10$^4$ kg m$^-3$ s$^-1$, respectively. The model names are indicated by the numbers in Fig. 1.

Figure 6.  The tropical Pacific mean state of the (a) observation (OBS), (b) CMIP5 MME, and (c) their differences (MME minus OBS). The SST is shown as black contours. For (a, b) the vectors are the surface wind stress (units: N m$^-2$), the shading indicates the precipitation (units: kg m$^-3$ s$^-1$), and the contour interval is 2$^\circ$C. For (c) the shading represents the SST biases, the vectors are the wind stress differences, the green (yellow) contours (interval: 0.1 kg m$^-3$ s$^-1$) denote the wetter (drier) biases of models, and the thickened red (black) contours indicate the 28$^\circ$C SST in the observation (MME). The composite anomalies of El Niño (d), La Niña (g) and the ENSO residual effect [the sum of El Niño and La Niña (i)] in the observation. The shading indicates the SSTA, vectors are the wind stress, and green or yellow contours are the precipitation anomalies (interval: 0.1 kg m$^-3$ s$^-1$). (e, h, k) As in (d, g, j) but for the MME of CMIP5. (f, i, l) The differences between the MME and observation. Shading indicates the SST anomalies, vectors are the surface wind stress anomalies, green (yellow) contours represent the wetter (drier) precipitation anomalies (interval: 0.15 kg m$^-3$ s$^-1$), and the thick contours represent the 28$^\circ$C SST. The magenta contour is for El Niño; the cyan contour is for La Niña; and the red (black) contour is for the MME (OBS). The precipitation value has been amplified by 10$^4$ kg m$^-3$ s$^-1$.

Figure 7.  The 28.5$^\circ$C SST during the (a) composite El Niño, (b) composite La Niña, (c) extreme warmest state and (d) extreme coldest state in CMIP5 models and the observation. The black (red) contours are the observation (MME). Colored contours are the nmulti-run-ensemble-means in individual models.

Figure 8.  The relationship between the air-sea interaction (measured as the air-sea coupling coefficient, alpha; units: 10$^-2$ N m$^-2$ k$^-1$) and the longitude of the equatorially averaged 27.5$^\circ$C SST during the (a) El Niño and (b) La Niña events in the observation and CMIP5 models. The relationship between alpha and the SST difference between the averaged North Pacific (0$^\circ$-20$^\circ$N, 130$^\circ$E-150$^\circ$W) and South Pacific (20$^\circ$S-0$^\circ$, 130$^\circ$E-150$^\circ$W) during (c) El Niño and (d) La Niña. $R$ is the linear correlation coefficient. Note that only the models with alpha larger than $0.2\times 10^-2$ N m$^-2$ k$^-1$ are taken into account in (a-d). (e, f) As in (a, c) but for the relationship between the ENSO asymmetry (measured as the skewness of Niño3 SSTA) and the mean state SST in the zonal and meridional direction distribution, respectively. Note that M30, M23, M12 and M19 were not included in the calculation of the linear correlation coefficients $r$ in (e). All models were considered for the $r$ in (f). The relationship between air-sea interaction over the central Pacific and ENSO asymmetry during (g) El Niño and (h) La Niña events. (i) The relationship between the differences of the air-sea interaction in the two phases of ENSO (El Niño minus La Niña) and the ENSO asymmetry. M30 was excluded for the calculation of the correlation coefficient $r$ in (g-i). For the model names of the numbers, refer to Fig. 1.

Figure 9.  The mean state of the (a) large skewness group (LSG) models, (b) small skewness group (SSG) models, and (c) their differences (LSG minus SSG). The composite El Niño in the (d) LSG, (e) SSG, and (f) their differences. The composite La Niña in the (g) LSG, (h) SSG, and (i) their differences. The residual effect (sum of the two phases) in the (j) LSG, (k) SSG, and (l) their differences. The SSTA is shaded and the wind stress anomaly is represented by vectors (units: N m$^-2$). The contours are the 28$^\circ$C SST. The magenta and cyan contours are the 28$^\circ$C SST in El Niño and La Niña, respectively. Black (red) contours are the mean 28$^\circ$C SST for the LSG (SSG) models. The LSG (SSG) consists of the 10 largest (smallest) skewness models.

Figure 10.  Propagation features of the (a) westerly wind anomaly (shaded), (a, b) precipitation anomaly (colored contours; interval: 2 kg m$^-3$ s$^-1$), and (b) SST anomaly (shaded) in the observation. The black contour is the 27.5$^\circ$C SST. Panels (c, d) and (e, f) are the same as (a, b) but for two typical models in the LSG and SSG. The contour interval for the precipitation anomaly is 2 kg m$^-3$ s$^-1$ in (c, d) and 1 kg m$^-3$ s$^-1$ in (e, f). GFDL-ESM2M is the large positive skewness model with a clear eastward-propagating SSTA feature, while Inmcm4 is the small negative-skewness model with a westward-propagating SST feature through the last 50 years of the 20th century. The zonal wind stress and precipitation have been amplified by 10$^2$ N m$^-2$ and 10$^4$ kg m$^-3$ s$^-1$, respectively.

In the tropical Pacific, the air-sea interaction is largely controlled by positive Bjerknes feedback (Bjerknes, 1969). The amplitude of ENSO events in climate models is usually determined by the sensitivity of the surface zonal wind stress and precipitation to SST changes. The linear intermodel correlation coefficients for the two atmospheric responses can reach 0.94 and 0.8, respectively (Figs. 5a and b). Models with stronger wind-SST feedback (model numbers 23, 13, 14, 15, 16, 37 and 6) all simulate stronger precipitation-SST feedback, except for M15 (CMCC-CM). A high positive correlation coefficient (0.87) was found between the precipitation response and surface zonal wind stress response from 1979 to 1999 among the models (Fig. 5c). A positive westerly (negative easterly) surface zonal wind stress response to warm (cold) SST changes is shown in Fig. 5d, with the relationship coefficient being 0.93 and 0.9, respectively. Nonlinear air-sea coupling was found in the seven models with the strongest westerly zonal wind stress response to the warmest SST changes. These seven models also simulate the strongest easterly zonal wind stress response to the coldest SST changes (Fig. 5d). Precipitation-SST feedback is similar to the zonal wind stress-SST feedback, but shows more diversity among the models, with relationship coefficients of 0.82 and 0.76 in the warm and cold SST phase, respectively (Fig. 5e). Consistent positive correlation was found between the two atmospheric responses among the CMIP5 models (Fig. 5f). In general, the modeled air-sea interaction changes in a more nonlinear way during the warm SST phase, compared with the more linear way during the cold SST phase. The precipitation response is nonlinear in the observation (Graham and Barnett, 1987), and the model parameter scheme may be based on that. The nonlinear air-sea interaction was only pronounced in the models with warmer SST changes that were larger than in the observation.

It must to be mentioned that models with the ENSO asymmetry stronger (weaker) than observation appear more nonlinear (linear) air-sea coupling feature. All of the seven strongest air-sea coupling models simulate a comparably larger ENSO asymmetry than the other models, expect M37 (NorESM1-ME) and M6 (BNU-ESM). This is because M37 and M6 simulate an almost equal or stronger La Niña-related response than El Niño-related response; hence, the weak ENSO asymmetry in the two models. The diversity and uncertainty in the precipitation-SST feedback could further amplify the zonal wind stress-SST feedback biases, since the latent heat released by the precipitation response could intensify the surface wind response (Zhang and Sun, 2014). Hence, the precipitation simulation skill in models could be amplified by the surface zonal wind stress, resulting in the SST feedback biases, which in turn force the atmospheric response biases. Thus, a vicious cycle is perpetuated that is hard to break.

Nevertheless, the weak ENSO asymmetry or weak nonlinear air-sea interaction in models is associated with the mean state biases over the tropical Pacific. The excessive cold tongue and double ITCZ problem was the pronounced biases that last for almost two decades in coupled models (Li and Xie, 2014) (Fig. 6c). The multimodel ensemble mean (MME) of CMIP5 climate models shows a very weak composite ENSO-related response (both El Niño and La Niña) and a weak ENSO residual effect compared with the observation (Figs. 6d-i). ENSO (both El Niño and La Niña) is confined to a narrow equatorial band in the MME (Figs. 6e and h). Due to the far westward extension of the excessive cold tongue in the mean state (Fig. 6b), the equatorial deep convection (represented by the precipitation) center shrinks too far west to the western Pacific warm pool region and has difficulties in developing in the central equatorial Pacific and, further, moving to the eastern equatorial Pacific during El Niño (Fig. 6e). Thus, the westerly wind anomaly is confined to the western Pacific, because of the cold tongue biases and the precipitation-related latent heat release biases. The results also provide some clues regarding the weak Pacific-North America (PNA) teleconnection caused by the weak deep convection during El Niño, as indicated by the cold biases along the coast of North America (Fig. 6f). The weak PNA teleconnection may cause the warmer SST biases in the middle latitudes of the North Pacific, which could weaken the colder biases induced by the cold North Atlantic biases through the Atlantic Multidecadal Oscillation (Wang et al., 2014). The weak ENSO asymmetry and the weak collective residual effect of ENSO could also enhance the cold tongue biases and warm pool biases in models (Fig. 6l).

The La Niña-related response pattern is opposite to the El Niño-related response in the models, while it is a more distinctive westward-shifted response pattern of SSTAs and deep convection in the observation (Fig. 6). The composites of SST during El Niño and La Niña in individual models are shown in Fig. 7 to cover the spread of model simulation. The cold tongue problem is more severe during El Niño than La Niña, especially in extreme cases (Figs. 7a and b). The excessive cold tongue problem leaves models little possibility to develop extreme La Niña and makes it difficult to develop extreme El Niño. The cold tongue almost reaches the coastline of the western Pacific in extreme La Niña cases and cannot shrink to the eastern Pacific during extreme El Niño events in the models (Figs. 7c and d). The cold tongue bias is amplified from El Niño to extreme El Niño events in the CMIP5 models (Fig. 7).

Overall, the CMIP5 climate models underestimate the air-sea coupling dynamics in the tropical Pacific during ENSO (both El Niño and La Niña) events. The air-sea coupling coefficient, estimated as the central Pacific positive-westerly (negative-easterly) surface zonal wind stress response to warm (cold) SST changes in the tropical Pacific, is underestimated in the CMIP5 models (Fig. S2). The SST changes in the Niño3 region are used instead of those in the central equatorial Pacific region to observe the direct ENSO-related air-sea interaction dynamics.

The air-sea coupling dynamics are associated with the mean-state SST spatial distribution. There is a positive relationship between the zonal and meridional asymmetry structure of the mean-state tropical Pacific SST. The excessive westward extension of the cold tongue tends to be accompanied by weaker air-sea interaction during ENSO in the models. The correlation coefficient between the longitude of the mean (27.5^°C) SST on the equator with the air-sea coupling intensity (alpha) is 0.44 and 0.52 during El Niño and La Niña, respectively (Figs. 8a and b). Note that we only counted alpha values larger than 0.2× 10^-2 N m^-2 k^-1. The CMIP5 climate models exhibit a westward extension of the mean (27.5^°C) SST compared with the observation, except M19 (NorESM1-ME), M18 (FGOALS-g2) and M15 (CMCC-CM). It is also evident that less colder the northwestern Pacific than the southwestern Pacific facilitates the stronger air-sea interaction during ENSO events, especially El Niño events. The positive linear correlation coefficient between the meridional asymmetry of the SST (north minus south) over the western and central Pacific in the mean state and alpha can reach 0.52 during El Niño, implying that the meridional location of the deep convection could also affect the air-sea interaction dynamics (Fig. 8c).

ENSO asymmetry is also positively associated with the zonal and meridional SST distribution structure in the mean state. Also, the skewness being positively related to the meridional asymmetry structure of SST is even more consistent in the models, with a correlation coefficient of 0.44 (Fig. 8f). The warmer the Northern Hemisphere tropical western and central Pacific is compared with the Southern Hemisphere, the stronger the ENSO asymmetry found in the models. In short, the weak ENSO asymmetry and weak nonlinear air-sea interaction over the tropical Pacific are all related to the mean-state SST biases in the models.

The meridional shift of the zonal wind stress or SST changes due to the seasonal cycle and phase locking of ENSO make the climatological and anomalous westerly wind shift to the south of the equator in late of the year (December-January). This sets up the discharge condition, which is not conducive to the further development of El Niño events, especially eastern Pacific (EP) El Niño events (Harrison and Vecchi, 1999; McGregor et al., 2012, 2013; Kug et al., 2003, 2009). The more asymmetrical the SST between the equatorial Northern and Southern Hemisphere Pacific, the larger the amplitude of El Niño events, especially EP El Niño events, and the larger the accompanying ENSO asymmetry among models.

The relationship between the ENSO asymmetry (Niño3 skewness) and ENSO-related air-sea interaction (alpha) is shown in Figs. 8g-i. The correlation coefficient is 0.47 and 0.21 during El Niño and La Niña, respectively (Figs. 8 g and h). Only the models with alpha larger than 0.2× 10^-2 N m^-2 k^-1 were included, and M30 was also excluded. A positive and slightly nonlinear relationship, with a linear correlation coefficient of 0.39, was found between the ENSO asymmetry and the differences in air-sea interaction between the two phases (El Niño minus La Niña; Fig. 8i). This confirms that the ENSO asymmetry bias is associated with the ENSO-related air-sea interaction biases. The results were similar when all of the model runs were taken into account for the correlation scatterplots in Fig. 8 (not shown).

Two groups of models based on the ENSO asymmetry were divided and their mean state differences and ENSO-related air-sea interaction differences were investigated. The strong (weak) ENSO asymmetry group consisted of the 10 largest (smallest) Niño3 skewness models. A weak ENSO-related response in terms of the SST anomaly, surface wind stress anomaly and precipitation anomaly was found in both the large skewness group (LSG) and the small skewness group (SSG) models, compared with the observation (Figs. 9d, e, g, and h). Furthermore, the LSG showed a stronger El Niño-related (La Niña-related) response than the SSG, with the El Niño-related response being even stronger (Figs. 9f and i). The residual effect in the LSG was an El Niño-like pattern, while in the SSG a La Niña-like pattern was exhibited (Figs. 9j and k). An prominent eastward expansion (shrinking) of the warm pool (cold tongue) from La Niña to El Niño was found in the LSG models, while the eastward shifting of the system was less distinct than the southward shifting in the SSG models (Fig. 9). These southward expanded warm pool biases may have been due to the positive SST-low cloud feedback in the Southern Hemisphere, as previously found by (Li and Xie, 2012), and the southward shifting of the deep convection may not have facilitated the development of stronger El Niño. The mean SST difference between the LSG and SSG also showed the warmer north off-equatorial Pacific SST differences. The collective effect of the ENSO differences between the two groups also showed the strong El Niño-like pattern, consistent with the aforementioned results.

However, it is reasonable that models can show uncertainty and diversity in the relationship between the mean state and the ENSO-related air-sea interaction or ENSO asymmetry. The mean state biases show diversity in CMIP5 models (Ham and Kug, 2014). The stratification of the equatorial Pacific can also show diversity, and the subsurface response was not examined in this study. The subsurface response is also critical to the surface response during ENSO. The overall linear ENSO stability formulated as the Bjerknes Stability (BJ) index was evaluated by (Kim et al., 2014) to quantify the air-sea feedbacks associated with the ENSO variability across the historical simulations in CMIP3 and CMIP5 models. They found weak zonal advective feedback, weak thermocline feedback, and a weak thermodynamic damping term over the last 50 years of the 20th century in CMIP5 models.

Theoretically, a linear system cannot induce asymmetry, even if external stochastic forcing is involved (An and Jin, 2004). Nonlinear air-sea coupling dynamics result in nonlinear dynamical heating (NDH), which is the omitted term in the BJ index that can contribute to the asymmetry of ENSO (An and Jin, 2004). The NDH in the heat budget of the upper ocean can be obtained from the following SST equation: \begin{eqnarray} \label{eq1} \!\dfrac{\partial T'}{\partial t}\!&=&\!-\!\left(\!{u}'\dfrac{\partial\overline{T}}{\partial x}\!+\! {v}'\dfrac{\partial\overline{T}}{\partial y}\!+\!{w}'\dfrac{\partial\overline{T}}{\partial z}\!+\! \overline{{u}}\dfrac{\partial T'}{\partial x}\!+\!\overline{{v}}\dfrac{\partial T'}{\partial y}\!+\! \overline{{w}}\dfrac{\partial T'}{\partial z}\right)\!-\nonumber\\[1mm] &&\left({u}'\dfrac{\partial T'}{\partial x}+{v}'\dfrac{\partial T'}{\partial y}+{w}'\dfrac{\partial T'}{\partial z}\right)+R' ,(1) \end{eqnarray} where T, u, v and w are the SST and the zonal, meridional and vertical velocities, respectively. The overbars and primes denote the climatological means and anomalies, respectively. Surface heat flux and subgrid-scale contributions (e.g., small oceanic diffusion, heat flux due to tropical instability waves) could attribute to the residual term R'. The second set of brackets in Eq. (2) indicates the NDH (An and Jin, 2004).

(Jin et al., 2003) found that the NDH can strengthen El Niño events and weaken subsequent La Niña events, thus leading to strong ENSO asymmetry during strong ENSO events (1982/83 and 1997/98). (An and Jin, 2004) noticed the eastward (westward) propagation of the westerly zonal wind anomaly provided a favorable (unfavorable) phase relationship between temperature and currents that resulted in positive (negative) nonlinear dynamical warming during both El Niño and La Niña, potentially contributing to ENSO asymmetry. In the observation, the ENSO asymmetry shows interdecadal changes synchronized with the interdecadal changes in ENSO characteristics. Pre-1976 (1950-76 in this study), ENSO occurred once every 3-4 years with comparably small amplitude. The SST and westerly wind anomalies propagated westward during that time, leading to the comparably small NDH, consistent with the comparably small skewness of ENSO during that period. Post-1976 (1977-99 in this study), ENSO occurred less frequently than pre-1976, with the period increased to 4-6 years but with larger amplitude. The SST and westerly wind anomalies propagated eastward or remained stationary during that time, leading to the large NDH, which was also consistent with the large skewness of ENSO during this period (Wang and An, 2002; An and Jin, 2004; An, 2009).

CMIP5 models cannot reproduce the observed interdecadal changes of ENSO asymmetry. However, the models with clear eastward propagation of their SST and westerly wind stress anomalies favor stronger ENSO and ENSO asymmetry (Fig. S3). Meanwhile, the CMIP5 models show problems in their precipitation-SST feedback, and the wind stress anomalies are synchronized to the precipitation anomalies well in both the models and observation. The cold tongue biases might suppress the development of deep convection in the central and eastern Pacific. Hence, the underestimation of negative (positive) NDH in the westward (eastward) propagation of SST and wind stress anomalies in climate models, in theory. Therefore, the CMIP5 climate models produce weak ENSO asymmetry over the whole period of 1950-99. Note that the precipitation data were only available from 1979 in the observation. The period 1979-99 is a strong ENSO and ENSO asymmetry period, as well as a period of eastward-propagating or stationary SST and wind stress anomalies (Fig. 10a). The LSG features a strong eastward-propagating SSTA and westerly wind stress anomaly, while there are very weak and less outstanding features in the SSG.

Two representative models were selected from the two groups, considering the warm or cold phase may cancel the other out due to the phase inconsistency in the models. The eastward-propagating SSTA model (GFDL-ESM2M) shows a strong westerly wind burst, strong air-sea coupling, and clear eastward propagation. In contrast, the westward-propagating SSTA model (Inmcm4) shows a weak westerly wind burst and weak air-sea coupling (Fig. 10). The deep convection develops easily in the central Pacific and further propagates to the eastern Pacific in GFDL-ESM2M, but this is not the case in Inmcm4. The former model has the warmer equatorial Pacific mean SST, and the latter has the colder. In addition, the large-skewness model (GFDL-ESM2M) still does not show the decadal variability of ENSO asymmetry found in the observation, implying many problems remain to be addressed regarding the simulation of ENSO asymmetry in climate models.

Evaluations have revealed that the ENSO statistics in CMIP5 climate models are slightly improved compared with CMIP3. The weak ENSO asymmetry in CMIP5 is closer to the observation. Also, the fact that strong ENSO activities tend to correspond to strong ENSO asymmetry is more evident in CMIP5. The diversity of ENSO amplitude is less and the ENSO asymmetry is mainly associated with El Niño events, especially those extreme ones, in CMIP5 models.

CMIP5 climate models show a weaker positive nonlinear atmospheric response to ENSO-related SST changes (wind-SST and precipitation-SST feedback) than the observation. The excessive cold tongue biases in the mean state limit the skill of the air-sea coupling and thus lead to the too weak nonlinearity of air-sea interaction and too weak ENSO asymmetry.

The excessive cold tongue, which extends too far west in CMIP5 models, pushes the deep convection to the western tropical Pacific warm pool region. The cold tongue bias in the mean state is unfavorable for the development of deep convection over the central equatorial Pacific. The deep convection has difficulty in further moving to the eastern equatorial Pacific during El Niño events, especially extreme events, which confines the westerly wind anomaly to the western Pacific. Hence, strong EP El Niño events do not develop as easily as central Pacific El Niño events in the CMIP5 climate models (Kim and Yu, 2012). This weakens the EP El Niño events, especially extreme El Niño events, and thus leads to weakened ENSO asymmetry in climate models.

The weak modeled ENSO asymmetry and weak modeled nonlinearity in air-sea interactions over the tropical Pacific are both associated with the SST distribution bias in the mean state. The results show that a warmer Northern Hemisphere tropical Pacific favors stronger air-sea interaction during El Niño, and thus stronger ENSO asymmetry.

This study has revealed that the cold tongue biases in the mean state play a role in weak positive nonlinear air-sea interaction, weak EP El Niño events, and weak ENSO asymmetry in CMIP5 models. It should be emphasized that the weak ENSO asymmetry and weak EP El Niño may further amplify the bias in the mean state and thus enhance the bias in the ENSO asymmetry. The cold tongue bias is even more prominent during extreme El Niño events, compared with "normal" El Niño events.

CMIP5 climate models also show some diversity and uncertainty in the relationship between the ENSO asymmetry (ENSO-related air-sea interaction) and the mean state biases. CMIP5 models cannot reproduce the decadal variability over the tropical Pacific and thus the decadal variability of ENSO asymmetry. Models tend to simulate one kind of propagation feature instead of the variety of features observed, implying that many problems remain regarding the simulation of ENSO asymmetry in climate models. To improve the simulation of ENSO statistics and ENSO asymmetry, we suggest that an accurate mean state structure, especially a realistic cold tongue and deep convection over the tropical Pacific, should be the priority for state-of-the-art climate models.