Responses of tropical Pacific air–sea CO2 flux (fCO2) to El Niño–Southern Oscillation (ENSO) events in 14 models from phase 5 of the Coupled Model Intercomparison Project (CMIP5) are examined. The contributions of sea surface temperature (SST), dissolved inorganic carbon in surface seawater (DIC), and total alkalinity of surface seawater (TALK) to interannual variation of ln(pCO2sea) (instead of partial pressure of CO2 in surface seawater pCO2sea) are quantified based on standardized empirical orthogonal function (EOF) results. Results show that six of the models have poor responses because they fail to reproduce observed interannual variation of pCO2sea in the central-eastern tropical Pacific. These six models underestimate the contribution of DIC interannual variation to interannual variation of pCO2sea in the central-eastern tropical Pacific due to a weak relation between interannual variation of upwelling and ENSO events or a weak relation (including no relation) between interannual variation of upwelling and that of DIC. Furthermore, some models have biases in interannual variation of DIC, in terms of both location and period, that are associated with interannual variation of modeled precipitation. It is also found that two models produce unreasonable interannual variation of bioproductivity, which enlarges interannual variation of DIC in the central-eastern tropical Pacific; this may partly explain why the influence of upwelling on interannual variation of DIC is weak in these models, even when the relationship between interannual variation of DIC and ENSO index is reasonable.
The increase of CO2 (an important greenhouse gas) in the atmosphere in recent decades has received considerable attention for its impact on climate change. The annual-mean CO2 uptake in the 1990s by the oceans was about 1.9–2.2 PgC yr−1 (1 Pg = 1015 g, and PgC is petagrams of carbon; Sabine et al. 2004; Takahashi et al. 2009), whereas the tropical Pacific as the biggest natural ocean source of CO2 to the atmosphere discharges about 0.44–0.48 PgC yr−1 (McGillis et al. 2004). These estimates are supported by different sources of data, such as the shipboard data of fugacity of CO2 and sea surface temperature (SST) from November 1981 to June 2004 (Feely et al. 2006), the observational partial pressure of CO2 in the global ocean during 1998–2011 (Landschützer et al. 2014), the observational partial pressure of CO2 in seawater pCO2sea for the reference year of 2000 (Takahashi et al. 2009), an offline model driven by reanalysis ocean products for 1961–2005 (Valsala et al. 2014), and inversion products, observational pCO2sea data, and forward models’ results of 1990–2009 (Ishii et al. 2014). This CO2 outgassing amount is significant compared with the annual-mean CO2 uptake by the global ocean.
Interannual variation of air–sea CO2 flux (fCO2) in the tropical Pacific is important for interannual variation of global fCO2 (Resplandy et al. 2015; Landschützer et al. 2016). Le Quéré et al. (2000) indicated that the variability of fCO2 in the tropical Pacific accounts for 70% of the variability in the global ocean. Large impacts were also noted by other studies (McKinley et al. 2004a,b). The interannual variation of global fCO2 is dominated by that in the tropical Pacific for two reasons. One reason is that El Niño–Southern Oscillation (ENSO) events occur in the tropical Pacific, which are the strongest interannual events on the global scale (Jones et al. 2001; Cai et al. 2014). The other is that the considerable amount of fCO2 in the tropical Pacific itself should not be neglected when considering the global flux (Jones et al. 2001). Therefore, the study of interannual variation of fCO2 in the tropical Pacific is crucial for accurate estimation of ocean’s CO2 uptake capability and for reliable future climate projection.
The ENSO cycle dominates the interannual variation of the tropical Pacific fCO2, which has been studied for several decades. According to the theories based on observations and model results, during El Niño events, the following changes occur between 5°S and 5°N in the Pacific: 1) reduced upwelling of dissolved inorganic carbon (DIC)-rich deep water to the surface, which in turn affects pCO2sea (Keeling and Revelle 1985; Chavez et al. 1999; Le Borgne et al. 2002; Obata and Kitamura 2003; Patra et al. 2005); 2) reduced upwelling of nutrient-rich deep water to the surface, which constrains bioactivities (Barber and Chavez 1983; Feely et al. 1997, 1999, 2002; Chavez et al. 1999; Le Borgne et al. 2002); 3) reduced CO2 solubility due to warm SST anomalies in the central-eastern tropical Pacific (Jones et al. 2001; McKinley et al. 2004a,b; Park et al. 2010; Ishii et al. 2014); and 4) reduced wind speed over the central-eastern tropical Pacific (McKinley et al. 2004b).
Researchers endeavored to quantify the change of fCO2 during El Niño events using observational data and obtained some robust results. For example, the results of Feely et al. (1994, 2004) and Doney et al. (2009a,b) showed that variability in upwelling during ENSO events dominates the tropical Pacific fCO2 due to its regulation of DIC, which accounts for up to 70% of the fCO2 variation in the tropical Pacific (Le Borgne et al. 2002), and the other 30% is from the variation of wind speed, biological processes, and so on. Research using ocean models also reached similar conclusions (Le Quéré et al. 2000; Obata and Kitamura 2003). Although state-of-the-art ocean models are capable of simulating important spatial characters of fCO2 in the tropical Pacific and the dominating role of ENSO events in fCO2 variation (Le Quéré et al. 2000; Jones et al. 2001; Obata and Kitamura 2003; McKinley et al. 2004a,b; Wetzel et al. 2005; Li and Xu 2013), the performances of some coupled models are less than satisfactory (Dong et al. 2016) because coupled models are not constrained as the ocean models, and combination of drivers differs from the observed counterpart (Keller et al. 2012). Investigation of ocean carbon can also help explain ocean dynamics (Keller et al. 2015).
The latest generation of Earth system models (ESMs) is an important tool for studying the linkage between ocean carbon cycle and climate change. Model output used in this study is from phase 5 of the Coupled Model Intercomparison Project (CMIP5), which contributed to the Fifth Assessment Report (AR5) of the Intergovernmental Panel on Climate Change (IPCC 2013). Ocean biogeochemistry components are considered in these ESMs for the first time in the CMIP history (Taylor et al. 2012), with different complexities ranging from a simple nutrient-forcing model to the multinutrient NPZD model, where NPZD stands for nutrient, phytoplankton, zooplankton, and detritus. In our previous studies, we found that most of the CMIP5 models simulated a better annual mean of fCO2 than the interannual variation of the flux, especially the interannual variation in the tropical Pacific (Dong et al. 2016). Here, we investigate the causes for failures in representing the strong interannual variation of fCO2 in the tropical Pacific by some CMIP5 models.
We analyze the responses of fCO2 and its associated variables [e.g., SST, DIC, and total alkalinity of surface seawater (TALK)] to ENSO events in 14 CMIP5 models by using their historical experiment output in the tropical Pacific. Among all the available CMIP5 ESMs, we choose only the models that provide enough output of variables closely associated with the fCO2 in the tropical Pacific. The CMIP5 historical experiments aimed at modeling the historical climate under the influences of natural and anthropogenic forcing derived from observations (Taylor et al. 2012) and were forced by the historical atmospheric CO2 concentration. We first access the overall responses of fCO2 in the tropical Pacific to ENSO events; then, we analyze the underlying causes of the responses. The paper is organized as follows. We describe the data, models, and methodology in section 2. We show the responses of the tropical Pacific fCO2 to ENSO events and the associated variables to the interannual variation of tropical Pacific fCO2 in section 3. We analyze the role of upwelling in interannual variation of DIC in surface seawater in section 4. We discuss likely causes for the shortcoming of modeled interannual variation of DIC in section 5, and provide a summary in section 6.
2. Data, models, and methodology
a. Data and models
The observational SST dataset from the Hadley Centre Sea Ice and Sea Surface Temperature (HadISST; Rayner et al. 2003; http://www.metoffice.gov.uk/hadobs/hadisst/data/download.html) is used to validate model capability of simulating ENSO events in terms of spatial pattern of EOF1 (the first mode of empirical orthogonal function) of SST and periods of ENSO during 1960–2000. The observational fCO2 from Rödenbeck et al. (2013, hereafter R2013; http://www.bgc-jena.mpg.de/CarboScope/?ID=oc) based on Surface Ocean CO2 Atlas (SOCAT), gridded pCO2sea and fCO2 data from Landschützer et al. (2015, hereafter L2015; http://cdiac.ornl.gov/ftp/oceans/SPCO2_1982_2011_ETH_SOM_FFN) based on a statistical model, and wind speed from Co-ordinated Ocean–Ice Reference Experiments, version 2.0 (CORE2.0; Large and Yeager 2004), are used to validate responses of simulated fields (e.g., fCO2, pCO2sea, and surface wind speed) to ENSO events.
The output of the historical runs of the twentieth century (historical; first ensemble member, r1i1p1) is used here, which uses changing conditions consistent with observations (Taylor et al. 2012). The models used in this study are listed in Table 1, together with their names, host institutions, names of ocean model component, and ocean model grids and the corresponding references. The variables include 10-m wind speed (WS10), near-surface relative humidity (RH), precipitation (PR), sea level pressure (SLP), sea surface salinity (SSS), SST, potential temperature of the seawater (PT), fCO2 (positive into the ocean), partial pressure difference of CO2 between air and sea ΔpCO2, pCO2sea, DIC, TALK, net primary productivity of the carbon by all kinds of phytoplankton (INTPP), chlorophyll concentration in surface seawater (CHL), silicate in surface seawater (SiO32−), zonal velocity of the seawater (UO), meridional velocity of the seawater (VO), and upward water mass transport (WMT). We employ the entrainment velocity at the bottom of the mixed layer to indicate the role of upwelling in transporting DIC-rich deep water to the surface in the model. The entrainment velocity was calculated by using the equation from Stevenson and Niiler (1983), based on CMIP5 model output of PT, UO, VO, and vertical velocity calculated by WMT. Because negative entrainment plays no role in the water exchange with the mixed layer, only the positive values are kept (see section 4 for more details). For the models whose output of ΔpCO2 or pCO2sea is not available (HadGEM2-CC, HadGEM2-ES, CNRM-CM5, IPSL-CM5A-LR, IPSL-CM5A-MR, and IPSL-CM5B-LR), we obtain ΔpCO2 or pCO2sea according to the method given in the appendix.
The period from 1960 to 2000 is chosen for analyzing the responses of fCO2 to ENSO events for two reasons. One is that most CMIP5 model output of the historical experiment are from the 1860s to 2000s, whereas many observational results of ocean carbon including the carbon budget of Le Quéré et al. (2015) started from the 1960s. The other is that we try to remove the long-term changes that play large roles in the fCO2 values (Séférian et al. 2013) as much as possible, for projecting the influence of ENSO. Some ESMs conducted more than one ensemble of historical experiment, which started from different integration times of piControl experiments. To reflect the common characters of these ensembles, we chose the latter part of r1i1p1 (1960–2000) to avoid the impact of initial fields. The output of the 14 ESMs in CMIP5 (Table 1) are obtained through data portals of the Earth System Grid Federation (http://pcmdi9.llnl.gov/). Note that the fCO2 from CanESM2 was expressed as CO2 instead of carbon, so the fCO2 of this model was multiplied by 12/44 to convert it from CO2 to carbon units before being used in our study. The results of model fCO2 suggest that the units of ΔpCO2 in IPSL-CM5A-LR, IPSL-CM5A-MR, and IPSL-CM5B-LR be expressed as micro–standard atmosphere (μatm) instead of pascal, so the units of ΔpCO2 in these three models were converted to pascal before being used in our study.
For discussion on the roles of associated variables in fCO2 interannual variability, standardized regression coefficients (RCs) between associated variables (see sections 3 and 4 for more details) are examined. To examine the correlations between fCO2 and associated drivers in the CMIP5 models, correlation analysis (including time-lagged correlation) between the first principal component (PC1) of fCO2 and PC1 of the drivers (DIC, SST, and INTPP) are carried out. In the tropical Pacific, the role of decadal variation in the annual variation of fCO2 is insignificant (Resplandy et al. 2015); therefore, after removing the trend and seasonal cycle, we do not perform further filtering. For further analysis of the roles of these variables in pCO2sea in the tropical Pacific, a partial differential equation of ln(pCO2sea) is used to quantify the contributions of the associated variables to the interannual variation of pCO2sea based on standardized EOF decomposition; in standardized EOF decomposition, the first mode of the EOF (EOF1) decomposition is the result of regressing the original anomalies against the standardized PC corresponding to EOF1, which is handy when we try to identify causes for poor responses of fCO2 to ENSO events in some models. The EOF decomposition in this study is performed for the tropical Pacific Ocean (15°S–15°N), since the EOF1 of global fCO2 in these models can be found in the work of Dong et al. (2016).
Although the periods of interannual SST fluctuation in the CMIP5 models are close to the observed (Dong et al. 2016), the areas of strong interannual SST fluctuation are slightly shifted toward the west in the models with respect to the observation (Bellenger et al. 2014; also, see Fig. S1 in the supplemental material). To include the areas where model ENSO events occur, we treat the ENSO index as the area-mean SST anomalies in the region of 5°S–5°N, 90°–170°W (covering both the Niño-3 and Niño-4 areas). Because of different resolutions in models and observational data, the CMIP5 model output and the observational data are all interpolated to the same grid of 1° × 1°. The continent in R2013 was not masked clearly with distinct boundary between land and ocean, so we added the mask of L2015 to the output of R2013. The monthly data used in our analysis are detrended and deseasonalized first to give anomalies, unless stated otherwise.
3. Responses of tropical Pacific fCO2 to ENSO events
The global distributions of fCO2 are reasonable in most of the CMIP5 models, and the main differences between model results and observation-based results exist in regions with strong vertical movement (Dong et al. 2016). Unfortunately, the difference of fCO2 between models and observations in the tropical Pacific is significant (Figs. S2a,b in the supplemental material), and the CO2 outgassing in some models like HadGEM2-CC and HadGEM-ES is weak. Unlike strong interannual SST fluctuation expressed as one standard deviation (1 SD) in the tropical Pacific (Fig. S1) in all of the 14 models, only CanESM2, CESM1(BGC), GFDL-ESM2G, GFDL-ESM2M, and MRI-ESM1 show strong interannual fluctuation of fCO2 as in R2013 and L2015 (Fig. S2a) and other observational studies (Feely et al. 1997, 2002; Le Borgne et al. 2002; Doney et al. 2009a,b; Valsala et al. 2014), and interannual variation of fCO2 in the other models is weaker than that in the observations, especially in the central tropical Pacific. The area-averaged results of fCO2 in the region 15°S–15°N, 109°E–77°W from 1960 to 2000 show this feature clearly (Fig. S3a in the supplemental material). The fluctuation of fCO2 is largely affected by pCO2sea, and among four factors (SST, SSS, DIC, and TALK) of pCO2sea, the fluctuation of SST in all 14 models is the most significant, and those of TALK and SSS are the weakest (Figs. S3b–d), which is consistent with previous findings by Doney et al. (2009b). For the fluctuation of DIC, besides of the four models mentioned above, obvious fluctuation is also found in the three models of IPSL, in the two models of MPI, and in NorESM1-ME (Fig. S3e). The relationship of interannual variation of fCO2 and that of each of these factors is discussed in the following content.
Figure 1 shows the RCs between fCO2 anomalies and the ENSO index using the observations (1982–2000) and output from 14 models (1960–2000) in the tropical Pacific (15°S–15°N). Stippled areas indicate that RC is significant at the 99% confidence level of the Student’s t test. Results based on R2013 and L2015 show that the positive RCs between fCO2 anomalies and the ENSO index occur in the central Pacific (R2013) or central-eastern Pacific (L2015). Different methods and time periods adopted in these two datasets induce some differences in the value and position of RCs, and the positive RCs are more extensive and significant in L2015.
Of the 14 models, 8 models [CanESM2, CESM1(BGC), GFDL-ESM2G, GFDL-ESM2M, IPSL-CM5A-LR, IPSL-CM5A-MR, IPSL-CM5B-LR, and MRI-ESM1] are capable of producing positive RCs between fCO2 anomalies and the ENSO index in the central-eastern Pacific, as in the observational results and single model results (Feely et al. 1997, 1999; Rayner et al. 1999; Jones et al. 2001; McKinley et al. 2004a,b; Takahashi et al. 2009). The remaining six models fail in this regard. For example, in HadGEM2-CC, HadGEM2-ES, MPI-ESM-LR, MPI-ESM-MR, and NorESM1-ME, the RCs between fCO2 anomalies and the ENSO index are negative in the eastern tropical Pacific and positive in the western tropical Pacific. In CNRM-CM5, the response of fCO2 to ENSO events is very weak according to the RC.
During ENSO events, the variation of fCO2 obviously lags behind that of SST or DIC because of the long time scale of gas exchange (Feely et al. 2002; Li and Xu 2013). With the spatial patterns similar to that of SST, the fCO2 in eight models [CanESM2, CESM1(BGC), GFDL-ESM2G, GFDL-ESM2M, IPSL-CM5A-LR, IPSL-CM5A-MR, IPSL-CM5B-LR, and MRI-ESM1] represents the positive correlations with SST in phase [at 0 lag; in terms of correlation coefficients (CCs) between two PC1s; Fig. 2], with the coefficients ranging from 0.62 to 0.74, which is statistically significant at the confidence level of 99% (492 months of data for each time series). In these eight models, the CCs between PC1 of fCO2 and PC1 of SST (PC1SST) in CESM1(BGC), GFDL-ESM2M, IPSL-CM5A-LR, and IPSL-CM5B-LR reach the largest values at about 5 months, which agrees with the result of Li and Xu (2013), whereas in the other four models (CanESM2, GFDL-ESM2G, IPSL-CM5B-LR, and MRI-ESM1) the CCs decrease with the lag time.
The CCs between PC1INTPP (the PC1 corresponding to EOF1 of INTPP) and PC1SST in GFDL-ESM2G and MRI-ESM1 are opposite those in other models (Fig. S4 in the supplemental material). However, the underlying causes for similar lag CCs in the two models are different: in MRI-ESM1, the CCs between PC1INTPP and PC1SST decrease rapidly as lag time increases (Fig. S4), which may cause the CCs between PC1TALK and PC1SST and between PC1DIC and PC1SST to decrease rapidly as lag time increases (Figs. S5 and S6 in the supplemental material), and further cause its failure in representing the obvious lag CCs between PC1 of fCO2 and PC1SST. In GFDL-ESM2G, the CCs between PC1INTPP and PC1SST reach the largest lag positive value at about 6–8 months (Fig. S4), which makes the correlations among DIC, TALK, and SST more complicated; that is, the CCs between PC1DIC and PC1SST reach the largest value at 2-month lag, whereas the CCs between PC1TALK and PC1SST decrease as time lag increases. The other two models (CanESM2 and IPSL-CM5B-LR) do not show CCs (between PC1 of the fCO2 and PC1SST) increasing as time lag increases, and we could not identify the reason.
a. Responses of air–sea gas transfer coefficient (wind speed) and ΔpCO2 (pCO2sea) to ENSO events
The fCO2 in each of these CMIP5 models is calculated based on the partial pressure difference of CO2 between air and sea, using the following equation:
In the above equation, K is the air–sea gas transfer coefficient and is the product of k and α, where k is the CO2 gas transfer velocity, a product of wind speed and Schmidt number (Sc) (Wanninkhof 1992), and α is CO2 solubility in seawater and is a function of SST and SSS (Weiss 1974). The term γice is the fraction of sea ice. Equation (1) shows that variables influencing the fCO2 in the tropical Pacific, which is ice free, can be divided into two kinds of K and ΔpCO2 (=pCO2air − pCO2sea). Both K and ΔpCO2 can affect the magnitude of fCO2, while ΔpCO2 also determines flux direction.
The term K is mainly associated with wind speed over the sea surface (Doney et al. 2009b; Turk et al. 2010), although controversies remain regarding the relationship between K and wind speed. Figure S7 in the supplemental material shows that the wind speed anomalies over the sea surface display negative RCs against ENSO index in the central tropical Pacific in the observation and 14 models (both during 1960–2000). Except for CESM1(BGC) and NorESM1-ME, in which WS10 output is not available, the other 12 models are capable of producing the weakening of wind speed over the central tropical Pacific during El Niño events. Note that the RCs of wind speed anomalies against ENSO events are much weaker in MPI-ESM-LR than in the other 11 models.
As mentioned above, the historical experiment in the CMIP5 uses time-varying forcing consistent with observations (including the natural and anthropogenic atmospheric CO2 concentrations) (Taylor et al. 2012). The biases in α may influence the flux results, while the temperature dependences of k and α approximately cancel out (Doney et al. 2009b). The influence of SSS on α is relatively small. Therefore, the influences of biases in α on the flux and ΔpCO2 are not significant, and model biases associated with model ΔpCO2 mainly stem from pCO2sea. Figure S8 in the supplemental material shows that in the models that are unable to reproduce the observed responses of fCO2, there are also poor responses of pCO2sea to ENSO events. Observational results show negative correlations between pCO2sea anomalies and ENSO index in the central-eastern tropical Pacific; this characteristic can be seen in all models except for CNRM-CM5, HadGEM2-CC, HadGEM2-ES, MPI-ESM-LR, and MPI-ESM-MR (Fig. S8). The RCs of pCO2sea anomalies against the ENSO index (Fig. S8) and those of fCO2 anomalies against the ENSO index (Fig. 1) are homologous in terms of spatial pattern and signs in both observations and the models, which shows the dominant role of pCO2sea in the interannual variation of fCO2 and the main reason for model deficiency in describing the relation between interannual variation of fCO2 and ENSO events. Previous studies based on observations and ocean models showed that interannual variation of pCO2sea (hence fCO2) is associated with variations of SST, upwelling, and bioactivity during ENSO events in the tropical Pacific (Murray et al. 1994; Feely et al. 1997; McKinley et al. 2004a,b). Among these variations, the variation in upwelling of DIC-rich deep water associated with ENSO events plays a dominant role in the variation of pCO2sea (hence fCO2) in the tropical Pacific (Feely et al. 1999). We want to know how much these variables contribute to pCO2sea in the CMIP5 models and will try to identify the causes for the poor responses of pCO2sea (hence fCO2) to ENSO events in some of these models.
b. Influences of DIC, SST, and TALK on ln(pCO2sea) variation
The variation of pCO2sea and the associated variables affecting such variation could be expressed in terms of Eq. (2) below (Takahashi et al. 1993). In Eqs. (2)–(4), variation , where includes both trend and seasonal cycle. The influence of SSS anomalies is considered to be small (Doney et al. 2009b); therefore, the last term on the right-hand side of Eq. (2) is ignored in this study:
The partial derivatives of pCO2sea with respect to these variables are not constant; they vary with the environment in the ocean. To accurately quantify the contribution of each variable, both sides of Eq. (2) are divided by pCO2sea for constant values of the partial derivatives of ln(pCO2sea) on these variables.
Before the calculation for determining contributions of these variables (SST, DIC, and TALK), the variables from the 14 models are analyzed using standardized EOF analysis, namely, EOF1 is calculated using the regression of the original anomalies against the standardized PC1 to ensure the comparability of EOF1s and PC1s in different models and different variables. The fraction of variance explained by the first three modes of EOF for SST, DIC, TALK, and ln(pCO2sea) are listed in Table 2, which shows that the explained variance fraction of EOF1 is much larger than those of the other two modes for the same variable. In most models, PC1s of these variables in the same model have similar or the same periods, which are about 3–5 years (Table 3). The CCs between PC1s of SST and DIC are 0.53–0.94, and those between PC1s of DIC and TALK are up to 0.81–0.97 (Table 4). All of these CCs are significant at the 99% confidence level of the Student’s t test. If the CCs are high enough, the anomalies of these variables can be changed to their EOF1:
To assure the interannual variations of variables on the right-hand side of Eq. (4) are in close phase so that the sum of them is meaningful, only when the CCs of PC1s of these two variables are over 0.70 and the periods of their PC1s are the same or very close (difference between any two periods is less one year; Table 3) are these two terms of these variables added. Figure 3 (top) lists nine models, with the influences of three variables added together and shown in the fourth column. The influences of the residual term and SSS on ln(pCO2sea) are obtained by the subtraction of EOF1 of ln(pCO2sea) from the accumulated influence of SST, DIC, and TALK on ln(pCO2sea), which is shown in the last column of Fig. 3, top. Figure 3 (bottom) shows three models (HadGEM2-CC, HadGEM2-ES, and IPSL-CM5B-LR) in which the CCs of PC1SST and PC1DIC or the CCs of PC1SST and PC1TALK are lower than 0.70 (Table 4), and only the influences of DIC and TALK are added. As shown in Figs. 3, the prominent contributions of SST to ln(pCO2sea) variation in these models occur in the central-eastern tropical Pacific with positive contributions during ENSO events. During El Niño events, the maximum contribution (judged by its absolute value) of DIC is negative and occurs in the western tropical Pacific in 11 models, where there is also the maximum contribution of TALK with positive values. The locations of maximum contribution correspond to the locations of maximum variation of variables, because the partial derivatives of ln(pCO2sea) for these variables are constant for the whole ocean basin. As illustrated by Fig. 3, top, the sum of three variables’ contributions is similar to the [EOF1 of ln(pCO2sea), unless stated otherwise], which also has a period similar to that of each of the three variables. Except for MPI-ESM-LR and MPI-ESM-MR, the contribution of DIC is larger than those of SST and TALK in the central-eastern Pacific and change of ln(pCO2sea) variation in seven models. Compared with the other models, the contributions of DIC in HadGEM2-CC and HadGEM2-ES are weaker in the whole tropical Pacific, and their remarkable contributions mainly exist south of the equator (Fig. 3, bottom). The spatial correlations of EOF1DIC and , of EOF1TALK and , and of EOF1SST and in HadGEM2-CC are −0.10, 0.00, and 0.85, respectively; in HadGEM2-CC, they are 0.04, 0.13, and 0.81, respectively (Fig. 3, bottom). The weak contributions of DIC in HadGEM2-CC, HadGEM2-ES, MPI-ESM-LR, and MPI-ESM-MR induce opposite characteristics of pCO2sea variation in the central-eastern tropical Pacific.
4. Role of upwelling in transporting DIC-rich deep water
For the models capable of reproducing the observed response of pCO2sea (hence fCO2) to ENSO events, the interannual variation of pCO2sea is dominated by DIC anomalies [e.g., CanESM2, CESM1(BGC), GFDL-ESM2G, GFDL-ESM2M, and MRI-ESM1]. According to McKinley et al. (2004a,b), the interannual variation of DIC in the tropical Pacific is dominated by interannual changes in upwelling during ENSO events. Accordingly, upwelling is examined using entrainment velocity at the bottom of surface mixed layer. Entrainment velocity is calculated using Eq. (5) below, as in Stevenson and Niiler (1983). Because the negative entrainment velocity plays no role in the materials exchange with the mixed layer, only positive values are kept:
where h is mixed layer depth (MLD). In this study, we use the potential temperature change from the surface temperature by 0.5°C to define the MLD. In Eq. (5), v is horizontal velocity, W is vertical velocity, and ∇ ≡ (∂/∂x, ∂/∂y) is the horizontal gradient operator; x and y are the zonal and meridional coordinates, respectively, and t is time. Here, subscript −h indicates the velocity at depth h.
The models chosen in this study all show upwelling, represented by positive entrainment velocity (PEV), in the central-eastern equatorial Pacific (Fig. S9 in the supplemental material). The areas with strong interannual variation of upwelling are not the same in these models. Except for the three IPSL models, MRI-ESM1, and NorESM1-ME, the strongest interannual variation happens in the warm pool (Fig. S10 in the supplemental material), especially in MPI-ESM-LR. The interannual variation of PEV in HadGEM2-CC and HadGEM2-ES is the weakest.
Figure 4 shows the RCs of PEV anomalies against the ENSO index during 1961–99. The stippled area means the RC is significant at the 99% confidence level of the Student’s t test. Except for the model of GFDL-ESM2G, whose PEV output is not available, the other models all display negative RCs of PEV anomalies against ENSO index along the equator, as well as positive RCs in the downwelling regions along both sides off the equator. The significant RC areas are mainly located in the central-western equatorial Pacific. In HadGEM2-CC, HadGEM2-ES, and MRI-ESM1, the RCs of PEV anomalies against ENSO index are weaker than those in the other models.
Although the nonlimiting nutrient of phosphate in surface seawater (PO43−) in the tropical Pacific mainly comes from the transport of local upwelling (Aumont et al. 1999), in some models PO43− is not the limiting nutrient in the tropics (Laufkötter et al. 2015). Schneider et al. (2008) pointed out that silicate (SiO32−) is a good proxy of ENSO to investigate responses of INTPP to ENSO. Hence, the physical role of upwelling in transporting the materials associated with biogeochemical processes can be examined through the relationship between SiO32− anomalies and upwelling anomalies (represented by PEV anomalies). The RCs of SiO32− anomalies against PEV anomalies (Fig. 5) show that these models [CESM1(BGC), GFDL-ESM2M, IPSL-CM5A-LR, IPSL-CM5A-MR, MPI-ESM-LR, and MPI-ESM-MR] are capable of producing the physical role of upwelling in transporting except for CNRM-CM5 and NorESM1-ME. In CNRM-CM5, the biogeochemical carbon process is not fully coupled with physical processes (Fig. 5), although this model is capable of producing the observed variation of PEV (hence upwelling) during ENSO events (Fig. 4).
In CNRM-CM5, a weak regression relationship also appears in the RC of DIC anomalies against PEV anomalies (Fig. 6). In NorESM1-ME, although the RCs of SiO32− anomalies against PEV anomalies are reasonable and significant at the 99% confidence level (of the Student’s t test), the weak regression relationship of DIC anomalies against PEV anomalies implies that the interannual variation of DIC is also affected by other processes besides upwelling. In HadGEM2-CC and HadGEM2-ES, the poor response of pCO2sea to ENSO events is more likely from physical fields; for example, the weak interannual variation of PEV leads to weak interannual variation of DIC. Figures 5 and 6 show that, in MPI-ESM-LR, the interannual variation of PEV and its effect on the interannual variations of biogeochemical materials mainly occur in the western equatorial Pacific. Other models that reproduce observed responses of pCO2sea to ENSO events are also capable of reproducing the relationship between PEV anomalies and DIC anomalies (including both strength and location).
Analyses above illustrate that the influence of upwelling on the interannual variation of DIC is reasonably represented in most of these models. However, in some models, physical transports of biogeochemical variables by upwelling are reasonable, but the interannual variation of DIC is not (e.g., NorESM1-ME), or the RC of upwelling against the ENSO index is relatively weak, but the interannual variation of DIC is reasonable (e.g., MRI-ESM1). We suspect that other processes may control the interannual variation of DIC in some models. If the dominant role of upwelling in interannual variation of DIC is relatively weak in a model, the roles of other variables (e.g., SST, bioactivity, and freshwater) may become more important, which may counteract or enlarge the interannual variability of pCO2sea (hence fCO2), and cause the model to fail in reproducing the observed response of pCO2sea (hence fCO2) to ENSO events.
a. Biogeochemical process
Bioactivity can influence DIC and TALK concentrations in seawater and further impact pCO2sea. As indicated by the observations and model results, there is interannual variation of bioproduction in the western and eastern tropical Pacific, but due to different causes. In the warm pool, the variation of nutrients induced by the change of halocline can affect bioproduction (Feely et al. 2002). The interannual variation of bioproduction in the eastern tropical Pacific is mainly induced by interannual change of upwelling (Feely et al. 2002; Doney et al. 2009b). Different from the previous results showing that the interannual variation of bioproduction in the warm pool is weaker than that in the central-eastern Pacific, some ESMs show the opposite results; for example, GFDL-ESM2G, GFDL-ESM2M, and MRI-ESM1 show too strong variation of INTPP in the warm pool (Fig. S11 in the supplemental material). The relationship between bioproductivity and ENSO events is examined using the RCs of INTPP anomalies against the ENSO index (Fig. 7). Except for CNRM-CM5, GFDL-ESM2G, and MRI-ESM1, the other models are capable of producing the negative relationship between INTPP anomalies and ENSO index in the central-eastern tropical Pacific. The uncorrelated relationship between INTPP anomalies and ENSO index in CNRM-CM5 is expected because the model has an uncorrelated relationship between nutrients and upwelling. The cause for the poor variation of INTPP in GFDL-ESM2G is also investigated. (We cannot get the CHL data of MRI-ESM1 to investigate further.) The CC between CHL anomalies and ENSO index in GFDL-ESM2G is reasonable, with positive correlation in the central-eastern tropical Pacific (Fig. S12 in the supplemental material). However, the CC between CHL and INTPP anomalies is negative, which is unrealistic; the exact cause for this bias is not clear and needs to be investigated in the future. The unrealistic variation of INTPP in GFDL-ESM2G and MRI-ESM1 enlarges the interannual variation of DIC, which may explain why the interannual variation of DIC is reasonable in MRI-ESM1, even though the relativity between upwelling and ENSO index is relatively weak. Studies by Laufkötter et al. (2015) pointed out that responses of INTPP to climate change in different models are different, which is consistent with our results that responses of INTPP to ENSO events show big differences, and that the influences of physical processes on INTPP variation may be important causes for these differences.
The ocean model components in these 14 ESMs all have free surface; therefore, the freshwater flux is likely to affect the concentrations of DIC and TALK at sea surface (Doney et al. 2009b; Turk et al. 2010; Walker Brown et al. 2015). The quantitative analysis in section 3 shows that the main contributions of DIC and TALK to the interannual variation of pCO2sea are located in the warm pool (Fig. 3) in the CMIP5 models, which is similar to the modeled interannual variation of precipitation in terms of EOF1 (Fig. S13 in the supplemental material). The modeled CCs between PC1DIC and PC1PR and between PC1TALK and PC1PR are significant, with values from −0.50 to −0.90 (Table 3). The interannual variation of normalized DIC [NDIC; salinity-normalized at SSS0 = 35; that is, NDIC = (DIC × SSS0/SSS; Ishii et al. 2014] and normalized TALK (NTALK; normalized by SSS, as for NDIC) in most models tends to be located in the central-eastern tropical Pacific in terms of EOF1 (Fig. 8 and Fig. S14 in the supplemental material), which implies that the role of precipitation should not be neglected in the interannual variations of DIC and TALK in the tropical Pacific of these models. Compared to the observed precipitation, most ESMs overestimate the variation of precipitation in the warm pool (Fig. S13), which becomes a source of bias for the interannual variation of DIC.
Turk et al. (2010) pointed out that precipitation could have a significant effect on air–sea CO2 exchange directly through increased transfer velocity and indirectly via carbon–chemical processes, especially over regions of low wind and high rain (such as the western tropical Pacific). In this study, we mainly consider the effect of precipitation on carbon–chemical processes because we find that poor responses of tropical Pacific fCO2 to ENSO events in some CMIP5 models are likely due to poorly simulated interannual variation of pCO2sea instead of poorly simulated interannual variation of wind speed.
CMIP5 models are the latest coupled ESMs and are important tools for studying the ocean carbon cycle (Nevison et al. 2016), which has significant influence on climate change (Friedlingstein et al. 2014). The surface ground temperature would increase 3.2°–5.4°C under the representative concentration pathway 8.5 (RCP8.5) scenario relative to that of 1850–1900 (Le Quéré et al. 2015). We utilized metrics of air–sea CO2 flux (fCO2), partial pressure difference of CO2 between air and sea ΔpCO2, partial pressure of CO2 in the seawater pCO2sea, dissolved inorganic carbon in the surface seawater (DIC), total alkalinity of the surface seawater (TALK), net primary productivity of the carbon by all kinds of phytoplankton (INTPP), chlorophyll concentration in the surface seawater (CHL), silicate in the surface seawater (SiO32−), sea surface temperature (SST), positive entrainment velocity (PEV), sea surface salinity (SSS), 10-m wind speed (WS10), and precipitation (PR) from 14 available CMIP5 models for two objectives: 1) to shed light on the roles of physical and biogeochemical processes in the variation of fCO2 in the tropical Pacific of these models by studying the responses of its associated variables to ENSO events and 2) to identify causes for poor responses of tropical Pacific fCO2 to ENSO events in some CMIP5 models.
Of the 14 models, 6 models (CNRM-CM5, HadGEM2-CC, HadGEM2-ES, MPI-ESM-LR, MPI-ESM-MR, and NorESM1-ME) fail to reproduce reduced CO2 outgassing during El Niño events because they cannot reproduce observed variation of pCO2sea in the central-eastern tropical Pacific (the opposite occurs during La Niña events), even though the 14 models can reproduce the weakening of WS10 over the central tropical Pacific during El Niño events. A close look at surface variables indicate that most CMIP5 ESMs included in this study can reasonably represent the influences of important physical and biological processes on the tropical Pacific pCO2sea during ENSO events (Table 5). Except for CNRM-CM5, the models are all capable of illustrating the role of upwelling in transporting DIC-rich and SiO32−-rich water to the surface ocean in the central-eastern tropical Pacific. However, the responses of associated physical and biological variables to ENSO events vary in terms of strength and amplitude, which makes most of the models fail to represent the lagged correlations between fCO2 variation and SST variation during ENSO events shown by observational and other ocean model results. Analysis based on the linearly expanded function of pCO2sea anomalies was adopted to quantify the contributions of SST, DIC, and TALK to ln(pCO2sea) variation. The results showed that those models that cannot reproduce observed interannual variation of pCO2sea underestimate DIC variation in the central-eastern tropical Pacific (e.g., HadGEM2-CC, HadGEM2-ES, MPI-ESM-MR, and MPI-ESM-LR), which may be due to weaker role of upwelling. In addition, the interannual variation of PR increases the interannual variation of DIC in the western tropical Pacific in most models, and induces model biases of interannual variation of DIC in location and period in some models.
Two models (GFDL-ESM2G and MRI-ESM1) show a good relation between interannual variation of DIC and ENSO index, even when there is a weak influence of upwelling on interannual variation of DIC. The analysis of modeled bioproductivity indicates that one cause is that these models produce the wrong interannual variation of INTPP, which enlarges, instead of reduces, interannual variation of DIC in the central-eastern Pacific.
This work is supported jointly by the National Grand Fundamental Research Program of China (2014CB441302) and the National Natural Science Foundation of China (Grants 41530426 and 41105087). We acknowledge the World Climate Research Programme’s Working Group on Coupled Modelling, which is responsible for CMIP, and we thank the climate modeling groups listed in Table 1 for producing and making available their model output. We thank the CanESM2 team members for the notification and confirmation of correcting carbon expression. We thank Prof. Yongfu Xu, Prof. Yongqiang Yu, and Dr. Jianqiong Zhan for their constructive comments.
Method for Calculating ΔpCO2 or pCO2sea
The values of ΔpCO2 in the two models of HadGEM2-CC and HadGEM2-ES and pCO2sea in the four models of CNRM-CM5, IPSL-CM5A-LR, IPSL-CM5A-MR, and IPSL-CM5B-LR were obtained from the improved Magnus formulas (Dickson et al. 2007), according to near-surface RH and SLP in the models as well as the contemporary observational atmospheric CO2 concentration pCO2air provided by CMIP5 (Taylor et al. 2012). The formulas are as follows:
pCO2sea: partial pressure of CO2 in surface seawater (Pa);
ΔpCO2sea: partial pressure difference of CO2 between air and sea (Pa);
SLP: sea level pressure (Pa);
RH: near-surface relative humidity (%);
pCO2air: observed atmospheric CO2 concentration (10−6);
ea: actual vapor pressure (Pa); and
e0: saturation vapor pressure (Pa).
Supplemental information related to this paper is available at the Journals Online website: https://doi.org/10.1175/JCLI-D-16-0543.s1.