Interannual variability of air–sea CO2 exchange is an important metric that represents the interaction between the carbon cycle and climate change. Although previous studies report a large bias in the CO2 flux interannual variability in many Earth system models (ESMs), the reason for this bias remains unclear. In this study, the performance of ESMs in phase 5 of the Coupled Model Intercomparison Project (CMIP5) is assessed in the context of the variability of air–sea CO2 flux over the tropical Pacific related to El Niño–Southern Oscillation (ENSO) using an emission-driven historical experiment. Using empirical orthogonal function (EOF) analysis, the first principal component of air–sea CO2 flux shows a significant relationship with the Niño-3.4 index in both the observation-based product and models. In the observation-based product, the spatial pattern of EOF1 shows negative anomalies in the central Pacific, which is, however, in contrast to those in several ESMs, and even opposite in sign to those in HadGEM2-ES and MPI-ESM-LR. The unrealistic response of the air–sea CO2 flux to ENSO mainly originates from the biases in the anomalous surface-water CO2 partial pressure (). A linear Taylor expansion by decomposing the anomalous into contributions from salinity, sea surface temperature, dissolved inorganic carbon (DIC), and alkalinity is applied to diagnose the biases. The results show that decreased during El Niño results from reduced upwelling of high-concentration DIC from deeper layers that overwhelms the increasing caused by warmer sea surface temperature. Overly weak reduction of vertical motion during El Niño and weak vertical gradients of climatological DIC concentration are the main reasons for biases in the negative surface DIC anomalies and eventually the anomalies. This study highlights the importance of both physical ocean responses to El Niño and climatological distributions of carbon-related tracers in the simulation of the interannual variability of air–sea CO2 fluxes.
The concentration of atmospheric CO2 increased by about 45% from 277 ppm in 1750 (Joos and Spahni 2008) to 402.8 ± 0.1 ppm in 2016 (Dlugokencky and Tans 2018). Meanwhile, the ocean has absorbed 160 ± 20 Pg C (1 Pg = 1015 g), which is about one-third of the anthropogenic CO2 emissions (Le Quéré et al. 2016), from the atmosphere via air–sea CO2 exchanges. Based on about 3 million measurements of surface-water CO2 partial pressure () obtained from 1970 to 2007, Takahashi et al. (2009) estimated that the total ocean uptake CO2 flux was 2.0 ± 1.0 Pg C yr−1 in the reference year 2000. However, the tropical Pacific Ocean (14°S–14°N) releases 0.48 Pg C yr−1 to the atmosphere, which accounts for 70% of the carbon released by the global tropical ocean (Le Quéré et al. 2000). The interannual variability of CO2 outgassing over the tropical Pacific Ocean is affected by El Niño–Southern Oscillation (ENSO), which exhibits large variability (Feely et al. 1999, 2002, 2006; Ishii et al. 2004; Wanninkhof et al. 2013). The peak-to-peak amplitude of CO2 flux during the cold and warm phases of ENSO is 0.40 ± 0.09 Pg C yr−1 according to the ocean general circulation models embedded with biogeochemistry models (OBGCMs), and 0.27 ± 0.07 Pg C yr−1 according to the diagnostic models (Ishii et al. 2014).
The mechanism of interannual variability of the air–sea CO2 flux over the tropical Pacific Ocean associated with ENSO has been well established through previous research on both observations and models. In the El Niño phase, the warmer sea surface temperature (SST) helps increase the air–sea CO2 flux, whereas the reduction in dissolved inorganic carbon (DIC) due to weakened upwelling acts to suppress the flux (Feely et al. 2006; Wang et al. 2006, 2015; Long et al. 2013). Among these two competing processes, the reduction in DIC is a dominant factor in controlling the interannual variability of air–sea CO2 flux related to ENSO (McKinley et al. 2004; Li and Xu 2013; Jin et al. 2017).
Earth system models (ESMs) are the latest generation of the state-of-the-art climate models, in which marine and terrestrial biogeochemical processes are added into coupled atmosphere–ocean general circulation models (AOGCMs). Compared with offline OBGCMs, which are driven by prescribed atmospheric data, ESMs characterize the interaction between the carbon cycle and the physical climate system. Owing to these features, the ESMs can simulate the physical, chemical, and biological processes on Earth, and hence ESMs are powerful tools for understanding and projecting changes in the global carbon cycle (Wang et al. 2016; Li et al. 2016; Li and Ilyina 2018; Kwiatkowski and Orr 2018). Because of the nonlinear carbon cycle and climate feedbacks in ESMs, the complex coupling between biogeochemical processes and the physical processes may also cause some unknown biases in the air–sea CO2 flux in ESMs (Boer and Arora 2009; Zhou et al. 2014; Schwinger et al. 2014). Therefore, before using ESMs to perform climate projections, including the study of carbon cycle and climate feedbacks, their accuracy in reproducing the air–sea CO2 fluxes for the present climate should be assessed.
Evaluation of the ocean carbon cycle in 18 ESMs from phase 5 of the Coupled Model Intercomparison Project (CMIP5) showed that models agree on the sign and magnitude of the CO2 flux, but they show weaker year-to-year variability than the observations (Anav et al. 2013). Dong et al. (2016) examined the interannual variability of global air–sea CO2 flux in 18 CMIP5 models and found that 12 models fail to represent the observed ENSO-related pattern over the tropical Pacific owing to the stronger interannual variability in the Southern Ocean and inconsistent period of air–sea CO2 flux with the period range of ENSO events. However, it is unclear whether the mechanism of ENSO-related air–sea CO2 flux interannual variations can be reproduced by CMIP5 ESMs. As a key component of the global carbon cycle, the air–sea CO2 flux over the tropical Pacific Ocean, especially the prominent relation with ENSO, is necessary for system evaluation. In this study, the physical and biogeochemical processes of the interannual variations in the air–sea CO2 flux over the tropical Pacific associated to ENSO were evaluated quantitatively by decomposing the anomalous into contributions from salinity, SST, DIC, and alkalinity (Alk). The objectives of this study are as follows: 1) to evaluate the performances of CMIP5 ESMs in the simulation of the air–sea CO2 flux over the tropical Pacific Ocean; 2) to address whether the coupled ESMs can reproduce the interannual variability of the air–sea CO2 flux related to ENSO; and 3) to understand why some ESMs fail to reproduce the relation between CO2 flux and ENSO.
2. Model, data, and analysis method
a. CMIP5 models
Fourteen ESMs from the CMIP5 were used based on the availability of monthly outputs of the emission-driven historical experiment (esmhistorical). Different from traditional concentration-driven historical experiment, forced by the historical atmospheric CO2 concentration, esmhistorical is forced by spatially distributed CO2 emissions reconstructed using estimated fossil fuel consumption from 1850 to 2005, which is more like the processes of global carbon cycle in nature. The atmospheric CO2 concentration freely varies with human emission and is modulated by carbon exchange between air, ocean, and the terrestrial biosphere (Taylor et al. 2012). Table 1 summarizes the model names, ocean carbon cycle components, horizontal and vertical resolutions of ocean models, and corresponding references. Only one ensemble member for each ESM is used. Seven (bold in Table 1) out of these 14 models are used to analyze sources of model bias depending on whether the data are available.
Observational and reanalysis datasets include 1) a well-structured climatological mean distribution of the air–sea CO2 fluxes with a resolution of 4° × 5° (Takahashi et al. 2009); 2) monthly mean air–sea CO2 fluxes and derived from neural network technique based on observation for the period from 1982 to 2011 and with a resolution of 1° × 1° (Landschützer et al. 2015); 3) surface wind fields from the European Centre for Medium-Range Weather Forecasts (ERA-Interim) (Dee et al. 2011); and 4) the SST dataset from National Oceanic and Atmospheric Administration Extended Reconstructed SST version 3b (ERSSTv3b) (Smith et al. 2008).
Although the observation-based products are used as “observational metrics” to gauge the performance of the models, these types of products rely heavily on data interpolation; for example, the product of Takahashi et al. (2009) employed an advection-based algorithm, whereas Landschützer et al. (2015) used a neural network interpolation.
For evaluating the climatology, since the observation-based products employ the reference year of 2000, the modeled results are averaged from 1990 to 2005 to match with the observation-based product. In the study of interannual variability, we focus on the boreal winter [December–February (DJF)] from 1982 to 2005, which is available in both the observation-based and model results. To focus on interannual variability, we filtered out variations longer than eight years from the original datasets using a Lanczos filter.
The upward ocean mass transport (kg s−1), a standard output of CMIP5 models, is used to analyze physical processes related to model biases. Its unit is different from that of the upward velocity (m s−1) in CMIP3. To make the upward ocean mass transport with different spatial resolutions in CMIP5 ESMs comparable, we transformed the upward ocean mass transport (kg s−1) to upward velocity (m s−1) by dividing the ocean gridcell area (m2) and seawater density (1025 kg m−3) for each model.
c. Analysis method
The air–sea CO2 fluxes (fgCO2) are calculated as follows:
where k is gas transfer coefficient (cm h−1). As in Landschützer et al. (2014, 2015), we now use ERA-Interim 10-m wind speeds U10 and the quadratic gas transfer formulation of Wanninkhof (1992) rescaled to a global mean gas transfer velocity of 16 cm h−1 [to match a recent estimate by Wanninkhof et al. (2013)]. The method used to calculate the coefficient for the gas transfer velocity is described by Wanninkhof et al. (2013) as . Therefore, the optimal coefficient for the gas transfer velocity parameterization is 0.36 cm h−1 (m s−1)−2 and
The Schmidt (Sct) number is calculated as follows:
where α is the CO2 solubility (mol m−3 atm−1; 1 atm = 1013.25 hPa) influenced by SST T (°C) and salinity S (psu) (Weiss 1974):
where A1, A2, A3, B1, B2, and B3 are constants.
The term is the difference in CO2 partial pressure between the surface seawater () and air (). The interannual variability of is controlled primarily by because of the relatively small interannual variations of .
The terms on the right side represent the effect of variability in SST, salinity-normalized DIC (nDIC), salinity-normalized alkalinity (nAlk), and sea surface salinity (SSS). The partial derivatives are estimated from the variation of due to the change of individual forcing variables, while keeping all other forcing variables constant. These partial derivatives are derived from the isolated code for calculating in an ocean carbon cycle model (Ilyina et al. 2013). We perturb a variable 50 times equidistantly around its climatology, keeping the others constant, and then use the least squares fitting to estimate the corresponding partial derivative. The perturbed ranges of variables are listed in Table 2. The climatological winter mean (DJF) during 1982–2005 is used.
Since every grid point of the spatial pattern is independent, when testing statistical significance on the spatial pattern correlation, we evaluate the effective sample sizes Neff for correlations of X and Y as
where N represents data length, and rX(τ) is an autocorrelation of time series X with a lag of τ. In this study, the estimate of the effective sample size is 18.
a. Annual mean
The annual-mean air–sea CO2 fluxes over the tropical Pacific Ocean in the observation-based products and multimodel ensemble (MME) are shown in Fig. 1. Positive values indicate that CO2 fluxes go upward from ocean to the atmosphere. The modeled air–sea CO2 fluxes are averaged in the period 1990–2005, matching with the observation-based product in reference year 2000 in Takahashi et al. (2009). In the observation-based product, the spatial pattern of air–sea CO2 fluxes shows a large contrast between the domains of the eastern and western tropical Pacific Ocean (Fig. 1a). Evidently, the eastern Pacific is characterized by larger efflux than that in the western Pacific. This feature is reproduced reasonably by MME (Fig. 1b). The spatial correlation coefficient between observation-based product and MME is 0.84, which is statistically significant at the 1% level.
In the observation-based product, the climatological annual-mean air–sea CO2 flux over the tropical Pacific (18°S–18°N, 120°E–80°W) is 0.47 Pg C yr−1 (Fig. 1c, brown bar), which mainly originates from the eastern tropical Pacific (0.43 Pg C yr−1). The annual CO2 flux (0.35 Pg C yr−1) simulated by the MME is about 25% less than the observation-based products, owing to less efflux over the eastern tropical Pacific in the MME (0.32 Pg C yr−1). However, the intermodel spread is large. Among the CMIP5 models, some models (BCC-CSM1.1, BCC-CSM1.1-m, BNU-ESM, CanESM2, IPSL-CM5A-LR, MPI-ESM-LR, and NorESM1-ME) reasonably reproduce the observed magnitude, while others cannot.
The Taylor diagram is presented in Fig. 2 to quantitatively assess the individual performance of CMIP5 models in the simulation of annual-mean air–sea CO2 flux. Better model performance is indicated by a higher pattern correlation coefficient and a normalized standard deviation nearer to 1.0. For the tropical Pacific (Fig. 2, black dot), the pattern correlation coefficients of most models with the observation-based product are statistically significant at the 5% level (correlation coefficient > 0.5), except for MIROC-ESM (0.15) and HadGEM2-ES (0.31), indicating that the overall characteristics of the spatial distribution of climatological air–sea CO2 flux over the tropical Pacific Ocean are reasonably reproduced by most of the ESMs; the details, however, are not well captured. Although the ratios of the standard deviations between model and observation-based product show large spread (from 0.55 for HadGEM2-ES to 2.64 for MRI-ESM), they are distributed symmetrically around 1.0. The pattern correlation coefficient (0.84) and the ratio of standard deviation (1.02) in MME are better than the individual ESM, which implies that the ratio of noise-to-signal in MME is lower than that for the individual ESM. For the western tropical Pacific (Fig. 2, red dot), the spatial correlation coefficients exceed 0.65 except for MIROC-ESM (0.30). However, the ratios of the standard deviations are generally higher than 1.0, indicating a stronger spatial variation than the observation-based product. The pattern correlation coefficients in the 14 ESMs of the eastern tropical Pacific (Fig. 2, blue dot) are less than those of the western tropical Pacific, and the MME yields better performances in terms of both the spatial distribution (0.75) and standard deviation (0.97) over the eastern tropical Pacific.
b. Interannual variability
Previous modeling studies have indicated that the interannual variability of the global air–sea CO2 flux is primarily determined by the interannual variability over the tropical Pacific (account for approximately 70% of the global variance; Le Quéré et al. 2000; McKinley et al. 2004; Obata and Kitamura 2003). However, an offline biogeochemical model indicates that the contribution of the tropical Pacific variability is approximately 40% of the global variance (Valsala et al. 2014). Here the averaged contribution from the 14 ESMs is 25%, ranging from 11% in MPI-ESM-LR to 46% in CESM1(BGC), showing large intermodel spread (figure not shown).
The monthly interannual variability intensity of the air–sea CO2 flux over the tropical Pacific Ocean simulated by CMIP5 MME is shown in Fig. 3. The interannual standard deviation for the period 1982–2005 is calculated as the anomalies relative to the mean annual cycle for each month of the year. Two peaks were observed (one each in the winter and summer); they are associated with the peak phase and developing phase of ENSO, respectively. Since the highest interannual variability intensity of the air–sea CO2 flux occurs in ENSO mature winter, we focus on the first leading mode of the air–sea CO2 flux in DJF.
The EOF analysis was performed to examine whether the ESMs can reproduce the leading mode of interannual variability of the air–sea CO2 flux in winter over the tropical Pacific Ocean (20°S–20°N, 120°E–80°W). In the observation-based product, the first leading mode (EOF1) exhibits negative anomalies over the central tropical Pacific, and this explains 49.8% of the variance (Fig. 4a). The first principal component (PC1) correlates with the Niño-3.4 index (the area-averaged SSTA over 5°S–5°N, 120°–170°W) with a high correlation coefficient of 0.92 (Fig. 5a), indicating the dominant influence of ENSO. In the CMIP5 ESMs, the correlation coefficients between the PC1 and Niño-3.4 index calculated based on the model SST range from 0.59 (IPSL-CM5A-LR) to 0.96 (CanESM2) (Fig. 5), and all values are statistically significant at the 1% level. The correlations between PC1 and the Niño-3.4 index are intended to be positive to ensure that the corresponding EOF1 patterns are comparable across models. Hence, the high correlations between PC1 and the Niño-3.4 index calculated based on the model SST show that ENSO also plays dominant role in determining the air–sea CO2 flux variability in the ESMs.
Biases in the simulated atmospheric CO2 concentration in CMIP5 ESMs were reported in several recent studies (Hoffman et al. 2014; Friedlingstein et al. 2014). To exclude the possible influence of biases in atmospheric CO2 on the interannual variability of air–sea CO2 flux, we repeated the EOF analysis on concentration-driven historical simulations in the corresponding 13 models; FIO-ESM was excluded as it did not have any air–sea CO2 flux output (see the online supplemental material, Figs. S1 and S2 therein). The EOF results of emission-driven and concentration-driven simulations differ only slightly. Hence, the influence of biases in simulated atmospheric CO2 concentration in ESMs can be neglected in this study.
The observed negative anomalies in the central Pacific are well captured by several ESMs [CESM1(BGC), FIO-ESM, GFDL-ESM2M, and MRI-ESM], but not so by some other ESMs, and some models (BNU-ESM, HadGEM2-ES, and MPI-ESM-LR) even exhibit positive anomalies in the spatial pattern of EOF1 (Fig. 4). The key question is what causes the biases in response of the air–sea CO2 flux to El Niño in some ESMs. Therefore, based on the availability of the model output fields related to the ocean carbon cycle, seven models are selected to examine the related physical processes and potential causes of biases. The high-skill models [CESM1(BGC) and GFDL-ESM2M], low-skill models (HadGEM2-ES and MPI-ESM-LR), and medium-skill models (CanESM2, GFDL-ESM2G, and NorESM1-ME) are distinguished by the pattern correlation coefficients shown in Fig. 4 between each model and observation-based product (Table 3).
Based on Eq. (1), the variations in the air–sea CO2 flux can be attributed to the gas transfer coefficient k, CO2 solubility α, and . Regression of these three factors with respect to Niño-3.4 index is shown in Fig. 6. In the observations-based result, the negative gas transfer coefficient anomalies locate at the central tropical Pacific (Fig. 6a2), associated with the westerly wind anomalies at the equator induced by El Niño warming. The observed gas transfer coefficient pattern is reproduced reasonably by six [CanESM2, CESM1(BGC), GFDL-ESM2M, GFDL-ESM2G, HadGEM2-ES, and NorESM1-ME] of seven ESMs (Figs. 6b2–h2), whereas weak positive anomalies over the central tropical Pacific are seen in MPI-ESM-LR (Fig. 6g2). For CO2 solubility, the observed negative anomalies (Fig. 6a3) associated with the warm SST anomalies during El Niño [Eq. (4)] are captured well by all the seven ESMs (Fig. 6). The westward shift of the anomalous solubility center in CanESM2, GFDL-ESM2G, MPI-ESM-LR, and NorESM1-ME (Figs. 6b3,e3,g3,h3) is expected to result from biases in ENSO-related anomalous SST and SSS based on Eq. (4). The largest contrast among the seven ESMs, however, is in the anomalies (Figs. 6a4–h4). The spatial patterns of the response to ENSO in the seven ESMs are similar to those of the air–sea CO2 flux, with negative anomalies in CanESM2, CESM1(BGC), and GFDL-ESM2M (Figs. 6c4 and 6d4) and weak anomalies in GFDL-ESM2G and NorESM1-ME (Figs. 6e4 and 6h4) but positive anomalies in HadGEM2-ES and MPI-ESM-LR (Figs. 6f4 and 6g4). Therefore, the biases in the response to ENSO are the dominant factors for the biases of the air–sea CO2 flux variations, especially the reverse sign biases in HadGEM2-ES and MPI-ESM-LR (Figs. 6f1 and 6g1).
To further understand the causes of biases of the response to ENSO in the several ESMs, a linear decomposition [Eq. (5)] was applied to diagnose contributions of SSS, SST, nDIC, and nAlk based on carbonate chemistry. All models indicate that the apparent anomalies result from the compensation between the positive effect of SST and negative effect of DIC; the effects of SSS and Alk are negligible (Fig. 7). For a quantitative comparison, the area-averaged regression coefficients over the region with the maximum variations (10°S–10°N, 180°–70°W) are shown in Fig. 8. In CESM1(BGC) and GFDL-ESM2M, the DIC effect dominates the anomalies; this result is consistent with previous studies (Long et al. 2013). The negative contributions to from DIC are 2.1 and 1.8 times the positive contributions from SST in CESM1(BGC) and GFDL-ESM2M, respectively (Fig. 8). In contrast, the negative effects of DIC in HadGEM2-ES and MPI-ESM-LR are weaker than the positive effect of SST (Fig. 7f3), with 55% and 96% ratios of the contributions of the DIC and SST (Fig. 8); this is the main cause of the biases of (Figs. 6f4 and 6g4). It indicates that the suppressed upwelling of DIC during El Niño in HadGEM2-ES and MPI-ESM-LR is not as evident as that in CESM1(BGC) and GFDL-ESM2M. In NorESM1-ME, although the negative contributions to from DIC are 2.0 times the positive contributions from SST (Fig. 8), the positive contributions from Alk counteract with the negative effect of DIC (Fig. 7h4) and cause the weak anomalies of (Fig. 6h4). The contributions of DIC and SST are equivalent in GFDL-ESM2G (Fig. 8), leading to the weak anomalies of (Fig. 6e4). The biases of the in CanESM2 (Fig. 6b4) result from both the underestimated DIC anomalies in the western Pacific and overestimated contributions of DIC and alkalinity in the eastern Pacific (Figs. 7b3 and 7b4).
Based on previous studies, the suppressed upwelling in the equatorial Pacific during El Niño is an important process that reduces upward transport of high-concentration DIC from deeper layers (McKinley et al. 2004; Long et al. 2013; Jin et al. 2017). This process can be expressed as , where w′ is the anomalous upward motion of seawater related to ENSO and is the vertical gradient of climatological DIC. To understand why negative contribution to from DIC is too weak in HadGEM2-ES and MPI-ESM-LR, we compare the simulated w′ and across different models. Vertical profiles of anomalous equatorial DIC, upward motion and seawater temperature in the upper ocean (above 250 m) over the tropical Pacific associated with Niño-3.4 index are shown in Fig. 9, as well as climatological upward motion, DIC, and its vertical gradient. The specific values of w′, , and their product averaged in the tropical Pacific are shown in Fig. 10. A weakened upward motion (w′ < 0; positive in the upward direction) induced by a deeper thermocline (Fig. 9c3; e.g., Christian et al. 2001) is observed at depths above 250 m in CESM1(BGC) and GFDL-ESM2M (Figs. 9b2 and 9c2). In the same region, the climatological DIC distributions show large vertical gradients () in both models with larger concentrations in the deeper layers (Figs. 9b4 and 9c4). Therefore, the upward transport of high-concentration DIC water to surface is evidently suppressed (), which is expected during an El Niño year [Figs. 9b1 and 9c1, as previously reported by Feely et al. (2006) and Landschützer et al. (2014)]. This process of anomalous upward motion suppressing the transport of high-concentration DIC water to surface simulated by CESM1(BGC) and GFDL-ESM2M is also shown in Fig. 10. In contrast, the negative anomalies of upwelling are weak in HadGEM2-ES (Fig. 9d2) and shift westward in MPI-ESM-LR (Fig. 9e2). Correspondingly, the positive anomalies of seawater temperature over the eastern and central equatorial Pacific in the upper ocean are weak in HadGEM2-ES (Fig. 9d3) and extend westward in MPI-ESM-LR (Fig. 9e3), indicating the weak deepening of the thermocline in HadGEM2-ES and MPI-ESM-LR. In addition, the vertical gradients of climatological DIC are too weak in the region where reduced upward motion occurs (Figs. 9d4 and 9e4). The weak anomalous upwelling and the vertical gradients of climatological DIC simulated by HadGEM2-ES and MPI-ESM-LR are also shown in Fig. 10. Therefore, the concentration of the surface DIC in HadGEM2-ES and MPI-ESM-LR does not decrease as strongly as that in CESM1(BGC) and GFDL-ESM2M. Finally, the negative effect of decreasing DIC on cannot overwhelm the positive effect of increasing SST (Fig. 8). These weak negative contributions arising from the insufficient reduction of vertical DIC transport to are responsible for the positive response of air–sea CO2 flux to El Niño in HadGEM2-ES and MPI-ESM-LR. In the medium-skill models (CanESM2 and NorESM-ME), although the vertical gradients of climatological DIC are as large as CESM1(BGC) and GFDL-ESM2M (Figs. 9a4,b4,c4,f4), the concentration of the surface IDC in CanESM2 and NorESM-ME cannot decrease as strongly as that in CESM1(BGC) and GFDL-ESM2M owing to the relatively weak negative anomalies of upwelling (Figs. 9a2,b2,c2,f2; see also Fig. 10).
4. Conclusions and discussion
The annual-mean climatology and interannual variability of the air–sea CO2 fluxes over the tropical Pacific Ocean in the 14 CMIP5 ESMs were evaluated in this study. We found that more than half of the models cannot reproduce well the spatial pattern of the tropical Pacific CO2 flux response to ENSO, despite the reasonable performances in the simulation of the mean states. The underlying reason was analyzed by comparing the models with different skills.
The spatial pattern of the observed annual-mean air–sea CO2 fluxes in the tropical Pacific can be reasonably reproduced in CMIP5 MME, despite large spreads in the magnitude and spatial gradient across the ESMs. The models also consistently show the controlling influence of ENSO on the interannual variability of the tropical Pacific air–sea CO2 fluxes in DJF. However, the dominant interannual modes (i.e., the spatial pattern) exhibit even stronger outgassing during an El Niño year in several models, such as HadGEM2-ES and MPI-ESM-LR, which is in contrast with the result obtained using observation-based data.
Further analysis show that the positive anomalies of air–sea CO2 fluxes during the warm events of ENSO result from the positive response of surface-water partial pressure (). During El Niño, increasing SST (apt to increase ) and decreasing surface DIC (apt to decrease ) in the central-eastern Pacific play competing roles in determining the , in which the DIC effect should be larger. By calculating the contributions of SST and surface DIC variation to anomalies, we found that the weaker DIC effect on negative pCO2 anomalies associated with ENSO in HadGEM2-ES and MPI-ESM-LR is the main cause of the air–sea CO2 flux biases. For HadGEM2-ES, the main cause is the weak vertical gradient of the climatological DIC, weak negative anomalies of upwelling, and weak deepening of thermocline. For MPI-ESM-LR, the westward shift of the negative anomalies of upwelling limits the DIC vertical transport, leading to weak surface DIC anomalies. Therefore, the negative effect of decreasing DIC on cannot overwhelm the positive effect of increasing SST, leading to positive and outgassing in HadGEM2-ES and MPI-ESM-LR.
In terms of the seasonal cycle, the equatorial Pacific is mainly controlled by the SST evolution rather than by the surface DIC, which shows less fluctuations due to the subtle balance between the compensating individual components of physical transport, biological drawdown, and air–sea CO2 flux (Feely et al. 2002; Jiang and Chai 2006). This subtle balance may be perturbed by tropical oscillations, such as ENSO. This study supports the importance of DIC transport perturbation in response of to ENSO, for which both the physical upwelling response and climatological DIC profile are crucial factors. Model developers should take care to reduce the biases in these two factors in future by observation-based calibration because biological processes can adjust and feedback dramatically to redistribution of carbon and nutrients; this is necessary to ensure good predictions of air–sea CO2 flux and biological activities (Chavez et al. 1999). MPI-ESM-LR and NorESM1-ME, which employ a similar ocean biogeochemical module (Assmann et al. 2010; Ilyina et al. 2013), show similar vertical profiles of the DIC response to ENSO, indicating that the processes and parameters included in the CMIP5 ocean carbon cycle modules (Table 1) partly contribute to DIC and air–sea CO2 flux biases in response to ENSO.
In the observation-based product, the spatial pattern of the first leading EOF of the air–sea CO2 flux clearly exhibits negative anomalies over the central tropical Pacific (Fig. 4a), which appears like the pattern of El Niño Modoki (Ashok et al. 2007). However, we speculate that the physical mechanism of air–sea CO2 flux anomalies over the central tropical Pacific is different from that of El Niño Modoki. The air–sea CO2 flux anomalies over the central tropical Pacific result from the gas transfer coefficient anomalies (Fig. 6a2) and anomalies (Fig. 6a4), which are regressed onto conventional ENSO index. The air–sea CO2 flux anomalies related to ENSO diversity should be another topic and deserve to be studied further.
We wish to thank the modeling groups joined in the CMIP5 and the Program for Climate Model Diagnosis and Intercomparison (PCMDI), which provide the model data. This research was supported by the National Natural Science Foundation of China (Grants 41330423 and 41605057) and National Key Research and Development Program of China (2018YFF0300105).
Supplemental information related to this paper is available at the Journals Online website: https://doi.org/10.1175/JCLI-D-18-0131.s1.