Journal cover Journal topic
Biogeosciences An interactive open-access journal of the European Geosciences Union
Journal topic
Biogeosciences, 15, 5315-5327, 2018
https://doi.org/10.5194/bg-15-5315-2018
Biogeosciences, 15, 5315-5327, 2018
https://doi.org/10.5194/bg-15-5315-2018

Research article 03 Sep 2018

Research article | 03 Sep 2018

# Drivers of future seasonal cycle changes in oceanic pCO2

Future pCO2 seasonality
M. Angeles Gallego1, Axel Timmermann2,3,4, Tobias Friedrich2, and Richard E. Zeebe1 M. Angeles Gallego et al.
• 1Department of Oceanography, School of Ocean and Earth Sciences and Technology, University of Hawaii at Manoa, Honolulu, Hawaii, USA
• 2International Pacific Research Center, School of Ocean and Earth Sciences and Technology, University of Hawaii at Manoa, Honolulu, Hawaii, USA
• 3Center for Climate Physics, Institute for Basic Science (IBS), Busan, South Korea
• 4Pusan National University, Busan, South Korea
Abstract

Recent observation-based results show that the seasonal amplitude of surface ocean partial pressure of CO2 (pCO2) has been increasing on average at a rate of 2–3 µatm per decade . Future increases in pCO2 seasonality are expected, as marine CO2 concentration ([CO2]) will increase in response to increasing anthropogenic carbon emissions . Here we use seven different global coupled atmosphere–ocean–carbon cycle–ecosystem model simulations conducted as part of the Coupled Model Intercomparison Project Phase 5 (CMIP5) to study future projections of the pCO2 annual cycle amplitude and to elucidate the causes of its amplification. We find that for the RCP8.5 emission scenario the seasonal amplitude (climatological maximum minus minimum) of upper ocean pCO2 will increase by a factor of 1.5 to 3 over the next 60–80 years. To understand the drivers and mechanisms that control the pCO2 seasonal amplification we develop a complete analytical Taylor expansion of pCO2 seasonality in terms of its four drivers: dissolved inorganic carbon (DIC), total alkalinity (TA), temperature (T), and salinity (S). Using this linear approximation we show that the DIC and T terms are the dominant contributors to the total change in pCO2 seasonality. To first order, their future intensification can be traced back to a doubling of the annual mean pCO2, which enhances DIC and alters the ocean carbonate chemistry. Regional differences in the projected seasonal cycle amplitude are generated by spatially varying sensitivity terms. The subtropical and equatorial regions (40 S–40 N) will experience a ≈30–80 µatm increase in seasonal cycle amplitude almost exclusively due to a larger background CO2 concentration that amplifies the T seasonal effect on solubility. This mechanism is further reinforced by an overall increase in the seasonal cycle of T as a result of stronger ocean stratification and a projected shoaling of mean mixed layer depths. The Southern Ocean will experience a seasonal cycle amplification of ≈90–120 µatm in response to the mean pCO2-driven change in the mean DIC contribution and to a lesser extent to the T contribution. However, a decrease in the DIC seasonal cycle amplitude somewhat counteracts this regional amplification mechanism.

1 Introduction

Owing to its large chemical capacity to resist changes in [CO2] (referred to as buffering capacity), the ocean has absorbed nearly half of the anthropogenic CO2 produced by fossil fuel burning and cement production since the industrial revolution . While the ocean's absorption of CO2 lowers the atmospheric concentration, it also increases the ocean's [CO2] and in turn lowers its buffering capacity. This leads to a reduction in the oceanic uptake of CO2 and an intensification of the pCO2 seasonal cycle (from now on referred to as δpCO2) . In a recent key observational study by , it was demonstrated that the δpCO2 amplitude increased at a rate of ≈2–3 µatm per decade from 1982 to 2015.

The pCO2 already experiences large seasonal fluctuations, which in some regions can reach up to 60 % above and below the annual mean pCO2 . An intensification of the δpCO2 amplitude could produce seasonal hypercapnia conditions , which together with increased [H+] seasonality and aragonite undersaturation events could expose marine life to harmful seawater conditions earlier than expected if considering only annual mean values. Moreover, a projected amplification of δpCO2 might increase the net CO2 uptake in some regions, such as the Southern Ocean, thereby further accelerating the decrease in the buffering capacity in that region .

The pCO2 seasonal amplitude is controlled mainly by seasonal changes in temperature (T) and biological activity together with upwelling changes that alter DIC concentrations. Usually, DIC and T changes work in opposite directions . In subtropical regions higher pCO2 values occur in summer when solubility decreases. In subpolar regions, pCO2 increases in winter when waters upwell that are rich in DIC and when respiration of organic matter takes place. Decreased subpolar pCO2 occurs in summer when primary productivity is higher and upwelling diminishes. Therefore, we find close relationships of δpCO2 with the ocean's [CO2], which controls chemical reactions, and with mean pCO2, which moderates exchange with the atmosphere. Both factors are related by the solubility constant that depends on temperature and salinity.

Furthermore, the regional differences in the influence of temperature and biology on δpCO2 are modulated by the ocean's buffering capacity. This is due to the ability of CO2 to react with seawater to form bicarbonate $\left[{\mathrm{HCO}}_{\mathrm{3}}^{-}\right]$ and carbonate $\left[{\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}\right]$, leaving only a small portion of dissolved carbon dioxide in the form of aqueous CO2 ([CO2(aq)]). [CO2(aq)] together with carbonic acid ([H2CO3]) are defined as [CO2]. Therefore, it is useful to define the total amount of carbon as DIC, which is the sum of the three carbon species ($\left[{\mathrm{HCO}}_{\mathrm{3}}^{-}\right]$, $\left[{\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}\right]$, and [CO2]). At current chemical conditions, most of the DIC is in the form of ${\mathrm{HCO}}_{\mathrm{3}}^{-}$, and therefore the buffering capacity is largely controlled by the ${\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}$ capable of transforming CO2 into bicarbonate through the reaction CO2(aq) +${\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}$+H2O= 2${\mathrm{HCO}}_{\mathrm{3}}^{-}$ . The larger the buffering capacity, the larger the pCO2 ability to resist changes in DIC. To quantify this capacity, we can introduce the sensitivity factor γDIC, which is inversely related to the buffering capacity, defined as ${\mathit{\gamma }}_{\mathrm{DIC}}=\partial \mathrm{ln}\left(p{\mathrm{CO}}_{\mathrm{2}}\right)/\partial \mathrm{DIC}$, . Other sensitivity factors are related to the total alkalinity (γTA), salinity (γS), and temperature (γT) changes and are defined in a similar way as $\partial \mathrm{ln}\left(p{\mathrm{CO}}_{\mathrm{2}}\right)/\partial \mathrm{TA}$, $\partial \mathrm{ln}\left(p{\mathrm{CO}}_{\mathrm{2}}\right)/\partial S$, and $\partial \mathrm{ln}\left(p{\mathrm{CO}}_{\mathrm{2}}\right)/\partial T$, respectively. It is important to note that the pCO2 is highly sensitive to temperature due to two factors: first through solubility changes that account for 2  3 of the present-day temperature impact, and second through the dissociation constants that control the carbon system reactions .

While the mechanisms controlling the seasonal cycle of pCO2 at present day are well documented, the future evolution of these drivers has not been fully elucidated. Current literature suggests that seasonal amplification is a consequence of an increase in the T and DIC contributions to δpCO2 and an increased sensitivity of the ocean to these variables .

The aim of our paper is to provide an in-depth analysis of the mechanisms controlling the future strength of δpCO2 and its regional differences using seven CMIP5 global Earth system models. Our analysis focuses on 21st century evolution using the Representative Concentration Pathway 8.5 (RCP8.5) scenario. We give a comprehensive analysis of the projected evolution of the DIC, TA, T, and S contributions to pCO2 seasonality. To achieve this goal, we derive explicit analytical expressions for pCO2 sensitivities in terms of γDIC, γTA, γT, and γS, thereby extending previous work done by .

2 Methodology

## 2.1 CMIP5 models

For our analysis, pCO2, DIC, TA, T, and S monthly mean output variables covering the period from 2006–2100 were obtained from future climate change simulations conducted with seven fully coupled Earth system models that participated in the Coupled Model Intercomparison Project, Phase 5 (CMIP5) (Taylor et al., 2012). The following models were selected based on data availability: CanESM2, CESM1-BGC, GFDL-ESM2M, MPI-ESM-LR, MPI-ESM-MR, HadGEM2-ES, and HadGEM2-CC (see the supplementary material of ). For the purpose of this paper, we used the Representative Concentration Pathway 8.5 (RCP8.5) future climate change simulations (IPCC2013). The ocean's surface data sets were regridded onto a $\mathrm{1}{}^{\circ }×\mathrm{1}{}^{\circ }$ grid using climate data operators (CDOs). The Arctic Ocean and the region poleward of 70 S are removed from the analyses because observational data for model validation are scarce.

## 2.2 Analysis of δpCO2

To elucidate the underlying dynamical, thermodynamical, biological, and chemical processes controlling δpCO2 we calculated a first-order Taylor series expansion of δpCO2 in terms of its four drivers, DIC, TA, T, and S. While T and S are controlled only by physics, DIC and TA are controlled by physical, chemical, and biological processes. Throughout this paper we use salinity-normalized DIC and TA using a mean salinity of 35 psu. This effectively removes the concentration–dilution freshwater effect, following the procedure of . The salinity normalized variables are referred to as DICs and TAs, corresponding to DIC $\cdot {S}_{\mathrm{0}}/S$ and TA $\cdot {S}_{\mathrm{0}}/S$, respectively. The freshwater effect on DIC and TA is now included in the S term, renamed as Sfw. For the Taylor series expansion, each variable (X= DIC, TA, T, and S) is decomposed into $X=\stackrel{\mathrm{‾}}{X}+\mathit{\delta }X$. The term $\stackrel{\mathrm{‾}}{X}$ represents the 21-year mean and δX denotes the seasonal cycle (calculated as the monthly mean deviation from the 21-year average). The Taylor expansion is then computed for an initial (2006–2026) and final (2080–2100) period. We use multi-decade means and eventually multi-model ensemble means to remove the effects of interannual variability. The full first-order series expansion is given by

$\begin{array}{ll}\mathit{\delta }p{\mathrm{CO}}_{\mathrm{2}}\approx & \frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial \mathrm{DIC}}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}\mathit{\delta }{\mathrm{DIC}}_{\mathrm{s}}+\frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial \mathrm{TA}}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}\mathit{\delta }{\mathrm{TA}}_{\mathrm{s}}\\ \text{(1)}& & +\frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial T}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}\mathit{\delta }T+\frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial S}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}\mathit{\delta }{S}_{\mathrm{fw}}.\end{array}$

Each term of the right-hand side of Eq. (1) represents the contribution from one of the four drivers of δpCO2. The analytical expressions for the derivatives (without the salinity normalization) are given by the following.

$\begin{array}{ll}\text{(2)}& & \frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial \mathrm{TA}}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}=\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot \frac{-{\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}}{\stackrel{\mathrm{‾}}{\mathrm{DIC}}\cdot \mathrm{\Theta }-{\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}^{\mathrm{2}}}& \frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial \mathrm{DIC}}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}=\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot \frac{\mathrm{\Theta }}{\stackrel{\mathrm{‾}}{\mathrm{DIC}}\cdot \mathrm{\Theta }-{\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}^{\mathrm{2}}}\\ & \frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial T}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}=\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot \frac{\mathrm{1}}{\stackrel{\mathrm{‾}}{\mathrm{DIC}}\cdot \mathrm{\Theta }-{\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}^{\mathrm{2}}}\\ & \left[{\stackrel{\mathrm{‾}}{\mathrm{TA}}}_{\mathrm{c}}\cdot \left(\frac{\partial {\mathrm{Alk}}_{\mathrm{c}}}{\partial T}+\frac{\partial \left[\mathrm{B}\left(\mathrm{OH}{\right)}_{\mathrm{4}}^{-}\right]}{\partial T}+\frac{\partial \left[{\mathrm{OH}}^{-}\right]}{\partial T}\right)\right\\ \text{(3)}& & -\mathrm{\Theta }\cdot \frac{\partial \left(\mathrm{DIC}-\left[{\mathrm{CO}}_{\mathrm{2}}\right]\right)}{\partial T}]-\frac{\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot }{\stackrel{\mathrm{‾}}{{K}_{\mathrm{0}}}\left(T,S\right)}\frac{\partial {K}_{\mathrm{0}}\left(T,S\right)}{\partial T}& \frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial S}{\mathrm{|}}_{\stackrel{\stackrel{\mathrm{‾}}{\mathrm{TA}},\stackrel{\mathrm{‾}}{\mathrm{DIC}}}{\stackrel{\mathrm{‾}}{T},\stackrel{\mathrm{‾}}{S}}}=\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot \frac{\mathrm{1}}{\stackrel{\mathrm{‾}}{\mathrm{DIC}}\cdot \mathrm{\Theta }-{\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}^{\mathrm{2}}}\\ & \left[{\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}\cdot \left(\frac{\partial {\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}}{\partial S}+\frac{\partial \left[\mathrm{B}\left(\mathrm{OH}{\right)}_{\mathrm{4}}^{-}\right]}{\partial S}+\frac{\partial \left[{\mathrm{OH}}^{-}\right]}{\partial S}\right)\right\\ \text{(4)}& & -\mathrm{\Theta }\cdot \frac{\partial \left(\mathrm{DIC}-\left[{\mathrm{CO}}_{\mathrm{2}}\right]\right)}{\partial S}]-\frac{\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot }{\stackrel{\mathrm{‾}}{{K}_{\mathrm{0}}}\left(T,S\right)}\frac{\partial {K}_{\mathrm{0}}\left(T,S\right)}{\partial S}\end{array}$

Here, $\mathrm{\Theta }=\left[{\mathrm{HCO}}_{\mathrm{3}}^{-}\right]+\mathrm{4}\left[{\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}\right]+\frac{\left[\mathrm{B}\left(\mathrm{OH}{\right)}_{\mathrm{4}}^{-}\right]\left[{\mathrm{H}}^{+}\right]}{\left({k}_{\mathrm{b}}+\left[{\mathrm{H}}^{+}\right]\right)}+\left[{\mathrm{H}}^{+}\right]+\left[{\mathrm{OH}}^{-}\right]$ and ${\stackrel{\mathrm{‾}}{\mathrm{Alk}}}_{\mathrm{c}}=\left[{\mathrm{HCO}}_{\mathrm{3}}^{-}\right]+\mathrm{2}\left[{\mathrm{CO}}_{\mathrm{3}}^{\mathrm{2}-}\right]$. The explicit T and S partial derivatives are given in the Supplement (Text S1). The first two derivatives coincide with the results of and , with the exception of the sign of [OH] in the term S. To verify this approach we compared the sum of the Taylor expansion terms with the full simulated range of δpCO2 from the model's output. The Taylor expansion reproduces the full seasonal cycle amplitude of the original climate model simulations well (Fig. S1 in the Supplement). The analytical expressions for temperature and salinity presented here are – to our knowledge – the first ones of their kind. Previously the calculation of these terms was based on the approximation given by or on numerical calculations.

To gain more insight into the processes causing the amplification of δpCO2 we introduce a method based on a second Taylor series expansion described below. Equation (1) can be rewritten using the expressions for the sensitivities γ determined by the relation $\frac{\mathrm{1}}{\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}}\frac{\partial p{\mathrm{CO}}_{\mathrm{2}}}{\partial \mathrm{X}}={\mathit{\gamma }}_{\mathrm{X}}$. These sensitivities have been historically used to represent the percentage of change in pCO2 per unit of DIC, TA, T, or S. With this notation, Eq. (1) can be expressed in the following way.

$\begin{array}{ll}\mathit{\delta }p{\mathrm{CO}}_{\mathrm{2}}\approx & \phantom{\rule{0.125em}{0ex}}\phantom{\rule{0.125em}{0ex}}\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot \left({\mathit{\gamma }}_{\mathrm{DIC}}\cdot \mathit{\delta }{\mathrm{DIC}}_{\mathrm{s}}+{\mathit{\gamma }}_{\mathrm{TA}}\cdot \mathit{\delta }{\mathrm{TA}}_{\mathrm{s}}+{\mathit{\gamma }}_{\mathrm{T}}\right\\ \text{(5)}& & \cdot \phantom{\rule{0.125em}{0ex}}\mathit{\delta }T+{\mathit{\gamma }}_{{S}_{\mathrm{fw}}}\cdot \mathit{\delta }{S}_{\mathrm{fw}})\end{array}$

Each term in Eq. (5) consists of three parts: $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$, the sensitivity γX, and the corresponding seasonal cycle δX. To understand which component is the main driver for δpCO2 changes, we perform a second Taylor expansion of the end of the century's δpCO2 around the initial state of the system in 2006–2026.

To maximize mathematical clarity we will introduce some definitions: first, we introduce the symbol Δ to indicate the difference between the period 2080–2100 and 2006–2026. Therefore, the total future change in δpCO2 is now referred to as ΔδpCO2. In the same manner, the total changes in sensitivities and seasonal cycles are written as $\mathrm{\Delta }{\mathit{\gamma }}_{\mathrm{DIC}{}_{\mathrm{s}},\mathrm{\Delta }{\mathit{\gamma }}_{\mathrm{TA}{}_{\mathrm{s}}},\mathrm{\Delta }{\mathit{\gamma }}_{\mathrm{T}},\mathrm{\Delta }{\mathit{\gamma }}_{\mathrm{S}{}_{\mathrm{fw}}}}$ and $\mathrm{\Delta }\mathit{\delta }{\mathrm{DIC}}_{\mathrm{s}},\mathrm{\Delta }\mathit{\delta }{\mathrm{TA}}_{\mathrm{s}},\mathrm{\Delta }\mathit{\delta }T,\mathrm{\Delta }\mathit{\delta }{S}_{\mathrm{fw}}$, respectively. Finally, we introduce the vector X formed by the four variables DICs, TAs, T, and Sfw as $\mathit{\left\{}{X}_{\mathrm{0}},{X}_{\mathrm{1}},{X}_{\mathrm{2}},{X}_{\mathrm{3}}\mathit{\right\}}=\mathit{\left\{}{\mathrm{DIC}}_{\mathrm{s}},{\mathrm{TA}}_{\mathrm{s}},T,S\mathit{\right\}}$. With this notation, we can write an expansion of Eq. (5) of the final state of the system by 2080–2100 named Xf around the initial state $\mathbit{X}{}^{\mathrm{i}}=\mathit{\left\{}{{\mathrm{DIC}}_{\mathrm{s}}}^{i},{{\mathrm{TA}}_{\mathrm{s}}}^{i},{T}^{i},{S}_{\mathrm{fw}}^{i}\mathit{\right\}}$ by 2006–2026 as

$\begin{array}{ll}& \mathrm{\Delta }\mathit{\delta }p{\mathrm{CO}}_{\mathrm{2}}=\phantom{\rule{0.125em}{0ex}}\phantom{\rule{0.125em}{0ex}}\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\sum _{k=\mathrm{0}}^{\mathrm{3}}{{\mathit{\gamma }}_{{X}_{k}}}^{i}\cdot \mathit{\delta }{{X}_{k}}^{i}\\ & +{\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}}^{i}\sum _{k=\mathrm{0}}^{\mathrm{3}}\mathrm{\Delta }{\mathit{\gamma }}_{{X}_{k}}\cdot \mathit{\delta }{{X}_{k}}^{i}\\ & +{\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}}^{i}\sum _{k=\mathrm{0}}^{\mathrm{3}}{{\mathit{\gamma }}_{{X}_{k}}}^{i}\cdot \mathrm{\Delta }\mathit{\delta }{X}_{k}\\ & +\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\sum _{k=\mathrm{0}}^{\mathrm{3}}\mathrm{\Delta }{\mathit{\gamma }}_{{X}_{k}}\cdot \mathit{\delta }{{X}_{k}}^{i}\text{(second-order\hspace{0.17em}\hspace{0.17em}terms)}\\ & +\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\sum _{k=\mathrm{0}}^{\mathrm{3}}{{\mathit{\gamma }}_{{X}_{k}}}^{i}\cdot \mathrm{\Delta }\mathit{\delta }{X}_{k}\\ \text{(6)}& & +{\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}}^{i}\sum _{k=\mathrm{0}}^{\mathrm{3}}\mathrm{\Delta }{\mathit{\gamma }}_{{X}_{k}}\cdot \mathrm{\Delta }\mathit{\delta }{X}_{k},\end{array}$

where the first, second, and third terms represent the contributions to ΔδpCO2 due to changes in the mean pCO2 ($\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$), the pCO2 sensitivities ($\mathrm{\Delta }{\mathit{\gamma }}_{X{}_{k}}$), and the seasonal cycles (ΔδXk), respectively; the fourth to sixth rows are the second-order terms. This method is similar to the one used by .

Figure 1RCP8.5 ensemble mean pCO2 seasonal cycle amplitude. Amplitude is calculated as climatology maximum minus minimum for (a) initial (2006–2026) and (b) final (2080–2100) periods. Initial and final climatologies were calculated as the monthly deviation from the respective 21-year period mean. (c, d) The ratio and difference between the δpCO2 amplitudes for 2080–2100 and 2006–2026, respectively. (e)(h) The zonal means of (a)(d), respectively, with the individual models shown as colored lines and the ensemble mean overlaid in black. Gray shading represents 1 standard deviation across the models.

3 Results and discussion

## 3.1 δpCO2 amplification

Figure 1a shows the ensemble mean δpCO2 amplitude (calculated as climatological maximum minus minimum) for the initial period 2006–2026. The values range from ≈98µatm for the high latitudes (40–70 S, 40–60 N) to ≈60µatm between 40 S and 40 N. The ensemble mean initial seasonal amplitude range is in good agreement with observational estimates calculated for the reference year 2005 and for the 1982–2015 period . The agreement between models and observations is remarkably good in the equatorial regions, but the initial amplitude is slightly overestimated in the middle and high latitudes (see Fig. S3 in the Supplement). The higher amplitude in models than observations is expected, as the initial period 2006–2026 already experienced an amplification compared to previous years. Moreover, found that the ocean's pCO2 historical trend is larger in models than observations when it is estimated in large-scale areas of the ocean. However, they found that model pCO2 trends agree with observations when the trends are subsampled to the locations where the observations were taken, and therefore they do a good job of reproducing well-known time series. Moreover, differences are expected as suggested that CMIP5 models perform well in reproducing the seasonal cycle timing, but still show considerable errors in reproducing the seasonal amplitude of pCO2 due to differences in the mechanisms represented in each model, especially in subpolar biomes.

Figure 2RCP8.5 ensemble mean $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ by (a) 2006–2026 and (b) 2080–2100. (c) Difference between 2080–2100 and 2006–2026. The North Atlantic and subpolar gyres show the largest difference between initial and final periods. The scale is different in each plot to enhance regional features.

By 2080–2100 the annual cycle amplitude attains values of ≈197 and 101 µatm in the high and middle to low latitudes, respectively (Fig. 1b). These seasonal variations correspond to 20 % and 18 % of annual $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ for the initial and final periods, respectively. Figure 1c shows that the global ocean δpCO2 will intensify by a factor of 1.5 to 3 for the 2080–2100 period relative to the 2006–2026 reference period. Figure 1d shows the difference in amplitude (ΔδpCO2); this pattern differs from the ratio because the ratio overestimates the amplification in areas where the initial amplitude is lower than ≈10µatm. used observations and a neural-network clustering algorithm to project that by the year 2100, the δpCO2 amplitude in some regions could be up to 10 times larger than it was in the year 2000. Our mean amplification factor estimation agrees with the mean threefold amplification found for most of the ocean by . However, the high values in this previous study cannot be reproduced here, mainly because we consider 21-year average ratios instead of single-year ratios, which are strongly affected by interannual variability. Using observations, found an increase of 2.2 µatm per decade, which is smaller than our findings of a total 42 µatm increase by the end of the century between 40 S and 40 N and a global mean change of 81 µatm on the high latitudes. This difference is again possibly due to the higher mean pCO2 values in models than observations.

Figure 3RCP8.5 ensemble mean seasonal cycle (δpCO2) and its Taylor decomposition. Colored lines indicate the contributions of DICs (blue), TAs (pink), T (red), and Sfw (green) to δpCO2 reconstructed from its Taylor decomposition (Eq. 1; dashed black). δpCO2 calculated from monthly pCO2 (solid yellow) is shown for comparison with the Taylor expansion. Column (a) shows the period 2006–2026 and column (b) shows the period 2080–2100. Each row represents the global zonal average for a different latitudinal band. Temperature dominates all latitudes except the Southern Ocean. In the 40–60 N band, T contribution is largely compensated for by DIC. The TAs and Sfw effects are rather small in all latitudes.

The global ocean mean amplification factor of δpCO2 roughly coincides with a doubling of $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ (Fig. 2). The direct relationship between these two is explained in Sect. 3.5. Figure 1e–h show the zonal mean of Fig. 1a–d; in general, towards the end of the century the pCO2 amplifies more in high latitudes, but so does the standard deviation uncertainty among models. This regional pattern agrees with the observation-based findings of , which show that high latitudes have already experienced a larger amplification than middle to low latitudes from 1982 to 2015. Furthermore, the same pattern is projected by CMIP5 models for the seasonal amplification of [H+] by the end of the century . This is expected from the near-linear relation between pCO2 and [H+]. These regional differences in amplification for pCO2 can be explained in terms of the relative magnitudes and the phases between the DIC, TA, T, and S contributions, which are explained in subsequent sections.

Figure 4RCP8.5 ensemble mean pCO2 sensitivities for DICs (row a), TAs (row b), T (row c), and Sfw (row d). Row (b) shows the negative γTA. The first and second columns show the sensitivities by 2006–2026 and 2080–2100, respectively. The third column shows the difference between 2080–2100 and 2006–2026 sensitivities. High latitudes show the largest difference between initial and final periods. While DICs and TAs sensitivities increase, the T and Sfw sensitivities decreases, except in the Western Pacific Warm Pool, where ${\mathit{\gamma }}_{S{}_{\mathrm{fw}}}$ increases.

## 3.2 Present and future drivers of δpCO2

To understand the driving factors of δpCO2 and its spatiotemporal differences, we split δpCO2 into the four different contributions from DICs, TAs, T, and Sfw for the initial and final periods, following Eq. (1). The results are shown in Fig. 3. For most of the ocean, the ensemble mean estimated contributions from DICs and T to the present-day δpCO2 are in good agreement with the data-based estimates of and , particularly in the equatorial regions (see Fig. S3 in the Supplement). However, our T and DIC contributions are slightly larger in middle and high latitudes for the same reasons the pCO2 seasonal amplitude is overestimated (see Sect. 3.1). Also, differences arise between our DICs contribution and the observation-based so-called “nonthermal” contribution because the nonthermal contribution also includes the total alkalinity and salinity effects. Nonetheless, between 40 S and 40 N our ensemble mean shows that δpCO2 is dominated by changes in temperature that control CO2 solubility, which decreases in summer, enhancing pCO2; this is in agreement with observations. The Southern Ocean is controlled by DIC, which responds to changes in upwelling and phytoplankton blooms. Both mechanisms act together to decrease (increase) DIC in summer (winter) .

Figure 5RCP8.5 ensemble zonal mean seasonal cycles of (a) δDICs and (b) δT for different latitudinal bands. Blue lines represent the 2006–2026 period, depicted for comparison with the 2080–2100 period shown by red lines. Different panels represent different latitudinal sections. Black arrows point out that while the T seasonal cycle is projected to increase in most of the ocean, global DICs is projected to decrease. The shading represents 1 standard deviation across the models. It is important to note that the scale is different for some of the latitudinal bands.

The models show that the δpCO2 in the 40 N to 60 N band is controlled by T, which disagrees with the abovementioned observations that show a non-temperature dominance in this band. The difference between models and observations arises from two regions: the North Atlantic basin and the northwestern Pacific, specifically near the Oyashio current and the outflows from the Sea of Okhotsk (see Fig. S3 in the Supplement). Most models show a T dominance in the North Atlantic basin; only CESM1-BGC and GFDL-ESM2M show a DIC dominance (see Fig. S4 in the Supplement). The North Atlantic is one of the major sinks of anthropogenic CO2; however, some models fail to estimate its uptake capacity . found that models with an efficient carbon sequestration present a DIC-dominated pCO2 seasonal cycle in the North Atlantic, but models with low anthropogenic uptake show a T dominance in this region. In the northwestern Pacific, found that coarse models are not able to capture the intricate oceanographic features of this area, and therefore the pCO2 seasonality is not well captured.

Figure 6Contribution of seasonalities, sensitivities, and mean pCO2 changes to ΔδpCO2. (a) Time series for the terms of Eq. (6) for different latitudinal bands. The Δ symbol represents the total century change calculated as the 2080–2100 value minus the 2006–2026 value. The total change in seasonal pCO2 (ΔδpCO2) is depicted as a dashed black. This change is decomposed into changes in seasonalities (ΔδX, purple), sensitivities (ΔγX, green), mean pCO2 ($\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$, red), and second-order terms (blue) summed over the four variables that control pCO2 (DIC, TA, T and S). For comparison with the expansion, ΔδpCO2 is calculated from model output (yellow). Column (b) shows the total change in DIC (dashed red) and T (dashed blue) contributions. Also shown are two components of the total change in these contributions: the $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ effect on the DIC (solid orange) and T (solid blue) contributions and the ΔδDIC (yellow) and ΔδT (light blue) effects. In column (a), the δpCO2 change follows the $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ effect. Column (b) shows that, actually, the leading cause of amplification is the $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ effect on the T contribution. It is important note the different scale between column (a) and (b). Also, the scale was reduced in the 15 S–15 N region to highlight its features.

Towards the end of the century (Fig. 3, right column), the amplification of δpCO2 is caused by an increase in the DICs and T contributions and to a lesser extent due to TAs and Sfw. Only in the high latitudes does the TAs contribution reinforce the DICs effect. The δDICs and δT relative phase and magnitude play an important role in causing regional differences of future δpCO2. For example, between 40 and 60, we find a lower amplification factor than at 30–40 in both hemispheres (Fig. 1c), contrary to what we expected from the general observed larger amplification at higher latitudes. In this band of lower amplification, the warm water from subtropical regions meets the nutrient-rich water from the subpolar regions, but the DICs and T effects are almost 6 months out of phase, and therefore their cancellation is larger than in the 30–40 latitude band where, for example, in the North Atlantic, there is 9-month phase difference between the two contributions. A clear illustration of this phase effect is found in the Supplement (Fig. S5).

In the Southern Ocean there is a shift in the maximum δpCO2 occurring from August–September to March–April (Fig. 3, last row). This shift is generated because the T contribution gains importance over DICs due to a reduction of δDICs magnitude at the same time that δT increases (Fig. 5). In the equatorial Pacific region (Fig. 5), T dominates over DICs but both contributions are small due to their low seasonality. Therefore, this region will experience a low amplification in δpCO2. In this region some models underestimate the pCO2 trend , and therefore the seasonal amplification might be underestimated too. In the following sections we conduct further analysis by decomposing each contribution as the result of three factors: the mean pCO2 ($\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$), the regional pCO2 sensitivities (γDIC, γTA, γT, and ${\mathit{\gamma }}_{{S}_{\mathrm{fw}}}$), and the seasonal cycles (δDICs, δTAs, δT, and δSfw) as determined in Eq. (5).

## 3.3 Future pCO2 sensitivities

The γDIC and γTA are projected to increase by the end of the century due to a lower ocean buffering capacity produced by increasing temperature and larger background concentrations of DIC . This agrees with our results shown in Fig. 4, which shows that all regions will experience an increase in γDIC and γTA. Lower buffer factors (higher sensitivities factors) are found in regions where DIC and TA have similar values, and they will decrease (increase) as the DIC  TA ratio in the oceans increases . The alkalinity sensitivity is negative, as pCO2 decreases with increasing alkalinity, but we show here the negative γTA for better comparison. γTA will increase (with negative values) more than the DIC sensitivity. However, seasonal changes in open-ocean TAs are small, and therefore the total contribution of alkalinity in our analysis is negligible compared to the DICs and T contributions. ${\mathit{\gamma }}_{{\mathrm{S}}_{\mathrm{fw}}}$ decreases everywhere except in the Western Pacific Warm Pool. In this region ${\mathit{\gamma }}_{{\mathrm{S}}_{\mathrm{fw}}}$ increases, probably due to future changes in precipitation that enhance the freshwater effect. In Fig. 4, the sensitivities (γ) are expressed as a percentage change in pCO2 per unit in DIC, TA, T, and S, respectively. This follows the approach of ; however, in their paper the authors compute the Revelle factor, which is related to γDIC as $R=\mathrm{DIC}\cdot {\mathit{\gamma }}_{\mathrm{DIC}}$. To illustrate the meaning of the sensitivities, we will focus on the subtropical North Pacific in the 15–40 N latitudinal band. In this region γDIC indicates an average 0.6 % change in pCO2 per unit of DIC in 2006–2026. Therefore, for a δDICs seasonal cycle amplitude of 40 µmol kg−1 and $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\approx \mathrm{400}$µatm, the total δpCO2 amplitude equals 96 µatm. Following the same reasoning, by 2080–2100, γDIC increases to 0.7 % and δDICs decreases to 30 µmol kg−1; therefore, for a $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ equal to 800 µatm, the δpCO2 amplitude due to δDIC amounts to 168 µatm.

Temperature sensitivity has been experimentally determined by , who found a value of 0.0423, meaning that pCO2 changes by about 4 % for every C. This value agrees with our global mean ensemble estimate of 0.0428. However, our analytical expression of γT shows that this value varies regionally and, for reasons unknown to us, it might decrease in the future to a global mean value of 0.0415 (Fig. 4c, third column). The T sensitivity is larger in colder regions and lower in the warmer tropics; however, colder regions will experience a larger reduction in γT, which locally prevents a larger amplification of the T contribution to δpCO2. In the next section we show that the T seasonality is projected to increase in high latitudes, strengthening the T contribution.

## 3.4 Future δDICs, δTAs, δT, and δSfw

Towards the end of the century, the global mean amplitude of δDICs is projected to decrease by ≈26 %–28 % in the high latitudes (Fig.5a), according to all the CMIP5 Earth system model simulations used here. In the middle to low latitudinal band there is no agreement between models; while some show an increase, others project a decrease in amplitude. As suggested by , the larger decrease in the Southern Ocean may be the result of changes in the shallow overturning circulation that prevent CO2 accumulation in this region. This reduction may be counteracted by the predicted increase in productivity owing to a suppression of light and temperature limitations .

According to the CMIP5 models, most of the ocean is projected to experience a slight increase in δT, as shown in Fig. 5b. All models show a slight increase in δT; only one model showed a slightly decrease in the southern region, and two models showed a decrease in the equatorial region during October to December. It is important to note that Fig. 5 shows the seasonal values with the mean T removed. Therefore, when considering the positive T trends, the absolute summer values show an increase and the absolute winter values a decrease. This agrees with the results of , who showed that models project a seasonal intensification of T with larger warm extremes and reduced cold extremes. The authors attributed the T seasonality intensification to an increased oceanic stratification and an overall shoaling of the mixed layer depth, which confines seasonal changes in a reduced volume of water, producing larger changes at the surface. They also showed that the intensification trends are stronger in summer than winter, as the mixed layer depth is shallower in summer. Moreover, ice-covered regions will experience the largest increase in T seasonality due to the loss of sea ice because the ice melting and freezing moderates the surface water temperature seasonality .

The TA seasonality is also projected to increase in the high latitudes according to all models, except CESM1-BGC, which shows a decrease. For δS (see Fig. S6 in the Supplement) there is no agreement among the different CMIP5 models, except in the Southern Ocean where all the models show a slight decrease. demonstrated that the seasonality of the drivers is important to determine future changes in [H+] seasonality. In the same fashion, our results show that the four δpCO2 drivers present changes in seasonality, and in particular δDICs and δT changes are important to explain future projections of the δpCO2 amplitude. The increase in δT enhances the δpCO2 amplification, and the reduction of δDICs in the Southern Ocean locally prevents a larger amplification.

## 3.5 Regional dominant factors

To identify the main cause of the δpCO2 amplification we use the Taylor series expansion method. With this method we consider the system's final state (δpCO2 by 2080–2100) as a perturbation of the initial state (δpCO2 by 2006–2026), as shown in Eq. (6). The expansion is done in three groups of variables: the seasonal cycles of DICs, TAs, T, and S (δX), the sensitivities of pCO2 to the same four variables (γx), and the mean pCO2 ($\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$). Therefore, each term of the expansion represents how much of the total δpCO2 change (indicated by ΔδpCO2 and calculated as the 2080–2100 value minus the 2006–2026 value) is due to the change in each of these factors. We also add the second-order terms that come from their combination. The results are shown in Fig. 6a and they indicate that the leading cause of the δpCO2 amplification is the change in $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ ($\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$), which confirms previous findings by .

It is important to note that our linear Taylor expansion approach neglects one aspect of the highly nonlinear carbonate chemistry of the ocean: it assumes $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ and the sensitivities as independent variables and therefore does not include the positive feedback between larger $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ and increasing γDIC (decreasing buffering capacity). Hence in the following, we use changes in $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ and changes in seawater carbonate chemistry synonymously, overall resulting in an enhanced response of δpCO2 to seasonal changes in DIC, TA, T, and S.

Considering regional differences, we note that amplification increases as we move poleward in spite of decreasing $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ (see Figs. 1 and 2). This characteristic geographical pattern of stronger high-latitude amplification is the result of larger present-day sensitivities (${\mathit{\gamma }}_{\mathrm{DIC}{}_{\mathrm{s}}}$, γT) and seasonal amplitudes (δDICs, δT) in the high latitudes that amplify the effect of $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ even when its value is small compared to other regions (see Eq. 6, first row term). Some exceptions can be found south of Greenland and near the subtropical gyres, where $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ reaches higher values and therefore also presents large amplification. We also found spatial differences on smaller scales; for example, the western equatorial Pacific presents lower initial δpCO2 and amplification than the eastern equatorial Pacific (see Fig. 1). This is because the eastern side of the basin has larger DICs and T contributions than the western side (see Fig. S2) as a consequence of the upwelling of cold, CO2-rich waters in the east, which lower the buffering capacity and induce larger δpCO2 amplitude due to the seasonal effects of productivity and solubility .

To further disentangle which of the two main drivers (DICs or T) is most affected by $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$, we decomposed the DICs and T contributions into their sensitivity, seasonal cycle, and $\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ components. Figure 6b shows the total DIC and T components together with the $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ and seasonal cycle effects on them. The effects from the sensitivities are not depicted, as they only play a minor role. Only the ΔγDIC term gains importance in the Southern Ocean (not shown). In most of the ocean, the $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ effect on T contribution is the leading cause of amplification. This effect is the result of seasonal solubility changes acting over a larger [CO2] . In the northern high latitudes, an increase in δT reinforces the amplification. In general, the ΔδT contribution gains importance as we move poleward in both hemispheres and therefore the second-order terms originating from $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}\cdot \mathrm{\Delta }\mathit{\delta }T$ also reinforce the amplification. Interestingly, in the high latitudes, the amplification through second-order terms is as important as the change in the seasonality of the drivers.

The Southern Ocean is an exception to the T dominance; in this region the $\mathrm{\Delta }\stackrel{\mathrm{‾}}{p{\mathrm{CO}}_{\mathrm{2}}}$ effect on the DICs contribution dominates, and the regional amplification is reinforced by low values of the mean buffering capacity (high ${\mathit{\gamma }}_{\mathrm{DIC}{}_{\mathrm{s}}}$). This result agrees with the findings of . In this area the amplification is counteracted by a reduction in δDICs.

4 Conclusions

In this study, we used output from seven CMIP5 global models, subjected to the RCP8.5 radiative forcing scenarios, to provide a comprehensive analysis of the characteristics and drivers of the intensification of the seasonal cycle of pCO2 between present (2006–2026) and future (2080–2100) conditions. By 2080–2100 the δpCO2 will be 1.5–3 times larger compared to 2006–2026. The projected amplification by the Earth system models and the possible causes of it are consistent with observation-based amplification for the period from 1982 to 2015 . However, the models slightly overestimate the present-day amplification, probably due to the larger pCO2 trends in models than observations .

The models confirm the well-established mechanisms controlling present-day δpCO2 . DICs and T contributions are the main counteracting terms dominating the seasonal evolution of δpCO2. Furthermore, the models show that under future conditions the controlling mechanisms remain unchanged. This result confirms the findings of that identified the same regional controlling mechanism for the past 30 years. The relative role of the DIC and T terms is regionally dependent. High latitudes and upwelling regions, such as the California current system and the coast of Chile, are dominated by DICs and the temperate low latitudes are driven by T. Only in the North Atlantic and northwestern Pacific do the models show a dominance of thermal effects over nonthermal effects, which is in disagreement with observations. This further illustrates the urgent need for models to accurately represent regional oceanographic features to accurately reproduce the δpCO2 characteristics.

In agreement with , the model projections towards the end of this century also demonstrate that the global amplification of δpCO2 is due to the overall long-term increase in anthropogenic CO2. A higher oceanic background CO2 concentration enhances the effect of T-driven solubility changes on δpCO2 and alters the seawater carbonate chemistry, also enhancing the DIC seasonality effect. The spatial differences of δpCO2 amplification, however, are determined by the regional sensitivities and seasonality of pCO2 drivers. For example, polar regions show larger sensitivity to DIC and T and larger seasonal cycles of DIC and T. Therefore, these areas present a strong enhancement of δpCO2 in spite of smaller changes in mean pCO2.

Moreover, the pCO2 seasonal cycle amplitude depends on the relative magnitude and phase of the contributions. The models ensemble mean reproduces the highly effective compensation of DICs and T contributions when they are 6 months out of phase, confirming previous studies . The compensation of DIC and T prevents a larger amplification of δpCO2, even when both contributions are largely amplified.

The amplification of the TA and S contributions has a small impact on δpCO2 in most regions, except in the high latitudes at which the TA contribution complements the DIC one, enhancing the nonthermal effect in this region.

The use of Earth system models allowed us to state the importance of including future changes in driver seasonalities for future δpCO2 projections. The T seasonality is projected to increase in most of the ocean basins, thereby reinforcing the δpCO2 amplification. The δT increase is consistent with an increase in stratification that will confine the seasonal changes in net heat fluxes to a shallower mixed layer . The DICs seasonality decreases in some cold areas and its reduction prevents a larger amplification. For the sensitivities, while γDIC increases, γT decreases. The latter phenomenon needs further study.

The increasing amplitude of δpCO2 might have implications for the net air–sea flux of CO2, in particular in regions where there is an imbalance between winter and summer values . Examples of such behavior can be found in the Southern Ocean (between 50 and 60 S) and in the latitude band from 2–40 in both hemispheres . Moreover, seasonal events of high pCO2 could have an impact on acidification, aragonite undersaturation events , and hypercapnia conditions . Therefore, understanding the drivers of future δpCO2 may help to better assess the response of marine ecosystems to future changes in carbonate chemistry. Finally, our complete analytical expansion of δpCO2 in terms of its four variables provides a practical tool to accurately and quickly diagnose temperature and salinity sensitivities from observational or modeling data sets.

Data availability
Data availability.

The CMIP5 data used in the analysis were obtained from https://esgf-node.llnl.gov/projects/esgf-llnl/ (last access: June 2017; Taylor et al., 2012). The pCO2 data product is available at https://www.nodc.noaa.gov/archive/arc0105/0160558/3.3/data/0-data/ (last access: May 2018); pCO2 data estimates from were obtained at http://cdiac.ess-dive.lbl.gov/ftp/oceans/NDP_094/ (last access: May 2018).

Supplement
Supplement.

Competing interests
Competing interests.

The authors declare that they have no conflict of interest.

Acknowledgements
Acknowledgements.

This work was supported by the National Science Foundation under grant no. 1314209 and by the Institute for Basic Science (IBS), South Korea, under IBS-R028-D1. We thank Peter Landschützer for kindly providing his data sets of thermal and nonthermal components of the pCO2 seasonal cycle.

Edited by: Katja Fennel
Reviewed by: two anonymous referees

References

Alexander, M. A., Scott, J. D., Friedland, K. D., Mills, K. E., Nye, J., Pershing, A. J., and Thomas, A. C.: Projected sea surface temperatures over the 21st century: Changes in the mean, variability and extremes for large marine ecosystem regions of Northern Oceans, Elem. Sci. Anth., 6, 9, https://doi.org/10.1525/elementa.191, 2018. a, b

Bopp, L., Resplandy, L., Orr, J. C., Doney, S. C., Dunne, J. P., Gehlen, M., Halloran, P., Heinze, C., Ilyina, T., Séférian, R., Tjiputra, J., and Vichi, M.: Multiple stressors of ocean ecosystems in the 21st century: projections with CMIP5 models, Biogeosciences, 10, 6225–6245, https://doi.org/10.5194/bg-10-6225-2013, 2013. a

Carton, J. A., Ding, Y., and Arrigo, K. R.: The seasonal cycle of the Arctic Ocean under climate change, Geophys. Res. Lett., 42, 7681–7686, 2015. a

Egleston, E. S., Sabine, C. L., and Morel, F. M. M.: Revelle revisited: Buffer factors that quantify the response of ocean chemistry to changes in DIC and alkalinity, Global Biogeochem. Cy., 24, GB1002, https://doi.org/10.1029/2008GB003407, 2010. a, b, c, d, e

Fassbender, A. J., Sabine, C. L., and Palevsky, H. I.: Nonuniform ocean acidification and attenuation of the ocean carbon sink, Geophys. Res. Lett., 44, 8404–8413, 2017. a, b

Fay, A. R. and McKinley, G. A.: Correlations of surface ocean pCO2 to satellite chlorophyll on monthly to interannual timescales, Global Biogeochem. Cy., 31, 436–455, 2017. a, b

Gorgues, T., Aumont, O., and Rodgers, K. B.: A mechanistic account of increasing seasonal variations in the rate of ocean uptake of anthropogenic carbon, Biogeosciences, 7, 2581–2589, https://doi.org/10.5194/bg-7-2581-2010, 2010. a, b

Goris, N., Tjiputra, J., Olsen, A., Schwinger, J., Lauvset, S. K., and Jeansson, E.: Constraining projection-based estimates of the future North Atlantic carbon uptake, J. Climate, 31, 3959–3978, 2018. a, b

Hagens, M. and Middelburg, J. J.: Attributing seasonal pH variability in surface ocean waters to governing factors, Geophys. Res. Lett., 43, 12528–12537, 2016. a, b

Hauck, J. and Völker, C.: Rising atmospheric CO2 leads to large impact of biology on Southern Ocean CO2 uptake via changes of the Revelle factor, Geophys. Res. Lett., 42, 1459–1464, 2015. a, b

Hauri, C., Friedrich, T., and Timmermann, A.: Abrupt onset and prolongation of aragonite undersaturation events in the Southern Ocean, Nat. Clim. Change, 6, 172–176, 2015. a, b

IPCC: IPCC Climate Change 2013: The Physical Science Basis, Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, 2013. a

Kwiatkowski, L. and Orr, J.: Diverging seasonal extremes for ocean acidification during the twenty-first century, Nat. Clim. Change, 8, 141–145, 2018. a, b, c

Landschützer, P., Gruber, N., Bakker, D. C. E., and Schuster, U.: Recent variability of the global ocean carbon sink, Global Biogeochem. Cy., 28, 927–949, 2014. a

Landschützer, P., Gruber, N., and Bakker, D.: An updated observation-based global monthly gridded sea surface pCO2 and air-sea CO2 flux product from 1982 through 2015 and its monthly climatology (NCEI Accession 0160558), Version 2.2. NOAA National Centers for Environmental Information, Dataset, https://doi.org/10.7289/V5Z899N6, 2017. a, b, c

Landschützer, P., Gruber, N., Bakker, D. C. E., Stemmler, I., and Six, K. D.: Strengthening seasonal marine CO2 variations due to increasing atmospheric CO2, Nat. Clim. Change, 8, 146–150, 2018. a, b, c, d, e, f, g, h, i, j, k, l

Lovenduski, N. S., Gruber, N., Doney, S. C., and Lima, I. D.: Enhanced CO2 outgassing in the Southern Ocean from a positive phase of the Southern Annular Mode, Global Biogeochem. Cy., 21, GB2026, https://doi.org/10.1029/2006GB002900, 2007. a

Mckinley, G. A., Takahashi, T., Buitenhuis, E., Chai, F., Christian, J. R., Doney, S. C., Jiang, M., Lindsay, K., Moore, J. K., Quere, C. L., Lima, I., Murtugudde, R., Shi, L., and Wetzel, P.: North Pacific carbon cycle response to climate variability on seasonal to decadal timescales, J. Geophys. Res., 111, C07S06, https://doi.org/10.1029/2005JC003173, 2006. a

McNeil, B. I. and Sasse, T. P.: Future ocean hypercapnia driven by anthropogenic amplification of the natural CO2 cycle, Nature, 529, 383–386, 2016. a, b, c, d, e, f

Pilcher, D. J., Brody, S. R., Johnson, L., and Bronselaer, B.: Assessing the abilities of CMIP5 models to represent the seasonal cycle of surface ocean pCO2, J. Geophys. Res. Oceans, 120, 4625–4637, 2015. a

Sabine, C. L., Feely, R. A., Gruber, N., Key, R. M., Lee, K., Bullister, J. L., Wanninkhof, R., Wong, C., Wallace, D. W. R., Tilbrook, B., Millero, F. J., Peng, T. H., Kozyr, A., Ono, T., and Rios, A. F.: The oceanic sink for anthropogenic CO2, Science, 305, 367–371, 2004. a

Sarmiento, J. L. and Gruber, N.: Ocean Biogeochemical Dynamics, Princeton University Press, Princeton, New Jersey, USA, 2006. a, b, c, d

Sasse, T. P., McNeil, B. I., Matear, R. J., and Lenton, A.: Quantifying the influence of CO2 seasonality on future aragonite undersaturation onset, Biogeosciences, 12, 6017–6031, https://doi.org/10.5194/bg-12-6017-2015, 2015. a, b

Shaw, E. C., McNeil, B. I., Tilbrook, B., Matear, R., and Bates, M. L.: Anthropogenic changes to seawater buffer capacity combined with natural reef metabolism induce extreme future coral reef CO2 conditions, Glob. Change Biol., 19, 1632–1641, 2013. a

Steinacher, M., Joos, F., Frölicher, T. L., Bopp, L., Cadule, P., Cocco, V., Doney, S. C., Gehlen, M., Lindsay, K., Moore, J. K., Schneider, B., and Segschneider, J.: Projected 21st century decrease in marine productivity: a multi-model analysis, Biogeosciences, 7, 979–1005, https://doi.org/10.5194/bg-7-979-2010, 2010. a

Takahashi, T., Olafsson, J., Goddard, J. G., Chipman, D. W., and Sutherland, S. C.: Seasonal variation of CO2 and nutrients in the high-latitude surface oceans: A comparative study, Global Biogeochem. Cy., 7, 843–878, 1993. a, b, c

Takahashi, T., Sutherland, S. C., Sweeney, C., Poisson, A., Metzl, N., Tilbrook, B., Bates, N., Wanninkhof, R., Feely, R. A., Sabine, C. L., Olafsson, J., and Nojiri, Y.: Global sea-air CO2 flux based on climatological surface ocean pCO2, and seasonal biological and temperature effects, Deep-Sea Res. Pt. II, 49, 1601–1623, 2002. a, b, c, d

Takahashi, T., Sutherland, S. C., Chipman, D. W., Goddard, J. G., Cheng, H., Newberger, T., Sweeney, C., and Munro, D. R.: Climatological distributions of pH, pCO2, total CO2, alkalinity, and CaCO3 saturation in the global surface ocean, and temporal changes at selected locations, Mar. Chem., 164, 95–125, 2014a. a

Takahashi, T., Sutherland, S. C., Chipman, D. W., Goddard, J. G., Newberger, T., and Sweeney, C.: Climatological Distributions of pH, pCO2, Total CO2, Alkalinity, and CaCO3 Saturation in the Global Surface Ocean, ORNL/CDIAC-160, NDP-094, Carbon Dioxide Information Analysis Center, https://doi.org/10.3334/CDIAC/OTG.NDP094, 2014b. a, b, c

Taylor, K. E., Stouffer, R. J., and Meehl, G. A.: An Overview of CMIP5 and the experiment design, B. Am. Meteorol. Soc., 93, 485–498, https://doi.org/10.1175/BAMS-D-11-00094.1, 2012.

Tjiputra, J. F., Olsen, A. R. E., Bopp, L., Lenton, A., Pfeil, B., Roy, T., Segschneider, J., Totterdell, I. A. N., and Heinze, C.: Long-term surface pCO2 trends from observations and models, Tellus B, 66, 23083, https://doi.org/10.3402/tellusb.v66.23083, 2014. a, b, c

Valsala, V. K., Roxy, M. K., Ashok, K., and Murtugudde, R.: Spatiotemporal characteristics of seasonal to multidecadal variability of pCO2 and air-sea CO2 fluxes in the equatorial Pacific Ocean, J. Geophys. Res.-Oceans, 119, 8987–9012, 2014. a

Völker, C., Wallace, D. W. R., and Wolf-Gladrow, D. A.: On the role of heat fluxes in the uptake of anthropogenic carbon in the North Atlantic, Global Biogeochem. Cy., 16, 1138, https://doi.org/10.1029/2002GB001897, 2002. a

Zeebe, R. E. and Wolf-Gladrow, D.: CO2 in Seawater: Equilibrium, Kinetics, Isotopes, Elsevier Science, Amsterdam, the Netherlands, and Philadelphia, PA, USA, 2001. a