CO₂ surface variability: from the stratosphere or not?

Fluctuations in atmospheric CO₂ can be measured with great precision and are used to identify human-driven sources as well as natural cycles of ocean and land carbon. One source of variability is the stratosphere, where the influx of aged CO₂-depleted air can produce fluctuations at the surface. This process has been speculated to be a potential source of interannual variability (IAV) in CO₂ that might obscure the quantification of other sources of IAV. Given the recent success in demonstrating that the stratospheric influx of N₂O- and chlorofluorocarbon-depleted air is a dominant source of their surface IAV in the Southern Hemisphere, I apply the same model and measurement analysis here to CO₂. Using chemistry-transport modeling or scaling of the observed N₂O variability, I find that the stratosphere-driven surface variability in CO₂ is at most 10 % of the observed IAV and is not an important source. Diagnosing the amplitude of the CO₂ annual cycle and its increase from 1985 to 2021 through the annual variance gives rates similar to traditional methods in the Northern Hemisphere (BRW, MLO) but can identify the emergence of small trends (0.08 ppm per decade) in the Southern Hemisphere (SMO, CGO).

1 Introduction

The surface abundance of CO₂, also called the Keeling curve (Fig. 1a), is used as the prime example of the human-driven increases in greenhouse gases. It is also used to demonstrate control of CO₂ by the land biosphere and the oceans through its annual cycles and interannual variations (Le Quéré et al., 2016, 2018). The inverse modeling of surface sources based on these CO₂ observations is used to infer regional sources of fossil fuel emissions as well as year-to-year changes in primary productivity of the biosphere or oceanic degassing (e.g., Gurney et al., 2002; Baker et al., 2006; Engelen et al., 2006; Nassar et al., 2011; Peylin et al., 2013; Frankenberg et al., 2016; Pandey et al., 2016; Nakazawa, 2020). There is concern that atmospheric variations in CO₂, and hence the net sources derived from them, may be affected by interannual variations (IAVs) in tropospheric mixing or stratosphere–troposphere exchange (STE) (Gaubert et al., 2019), but there are no definitive studies. For example, Nazakawa's (2020) review of greenhouse gas studies mentions the stratosphere only in connection with CH₄ and N₂O, not with CO₂.

Figure 1 (a) NOAA surface CO₂ monthly data (ppm, mole fraction) from Dlugokencky et al. (2021b). The five sites are as follows: BRW – Barrow, AK (now Utqiaġvik), at 71^∘ N, 156^∘ W; MLO – Mauna Loa, HI, at 20^∘ N, 156^∘ W; SMO – Tutuila, Am. Samoa, at 14^∘ S, 171^∘ W; CGO – Cape Grim, Tasmania, Australia, at 41^∘ S, 144^∘ E; SPO – South Pole at 90^∘ S. Monthly average in situ observations are used with gaps filled by flask data at the same site. Only one point is interpolated (SMO, May 2015). CGO is flask only. (b) The mean annual cycle in surface CO₂ (ppm) at four sites (BRW not shown) using 1985 through 2021 calculated from the residuals after the polynomial fit was removed. The eCO2 model results are shown as dotted lines for SMO and SPO with the same color coding; eCO2 has similar phasing at SPO and is double-peaked at SMO, but the amplitudes (∼0.18 ppm) are much smaller than observed. The sCO2 amplitude is even smaller (∼0.06 ppm) and not shown. (c) Observed CO₂ variability (ppm) derived from the monthly averages of two SH extratropical stations (SPO and CGO) for the period 1985–2021. The poly-fit (solid red curve) shows the residuals after removal of a polynomial fit (average of third- and fourth-order in time). A 12-month running mean (thin dashed blue curve) is derived from the poly-fit residuals and removes the annual cycle. The interannual variability (IAV, solid black curve) is derived from bandpass filtering described in the text. The bandpass limits [0.20 0.80] are set to truncate periods longer than 5 years and shorter than 1.25 years. (d) Surface CO₂ IAV (ppm) for SH extratropics. The SPO + CGO IAV (black solid line) is compared with the SMO IAV (thin dotted black line). Other IAVs shown are (1) model-calculated sCO2 (dashed red line) from the stratosphere-driven influx of aged, low-CO₂ air and (2) N₂O observed IAV (dashed blue line) scaled to match flux of low-CO₂ air. The timing of the QBO phase change in equatorial zonal wind at 40 hPa from negative (easterlies) to positive (westerlies) is denoted with thick vertical gray bars. The timing of moderate to extreme El Niños (red stars) and La Niñas (blue stars) is also shown. (e) CO₂ annual amplitude (ppm, max–min) derived from the variance across 12-monthly values. Each monthly point (centered on the beginning of each month) is the standard deviation of the surrounding ±6-monthly means, scaled by 2×2.5 to give the max–min amplitude as if it were a sine curve. The line fits for BRW and MLO are shown. The slope and standard error (SE) in units of parts per million (ppm) per decade are given in the legend. The SE is calculated conservatively based on the number of years rather than the number of months.

The possibility of a true STE-driven IAV CO₂ signal, raised by Gaubert et al. (2019), has not been seriously investigated. For the most part, when studies investigate the stratospheric influence on CO₂ source inversions, they are not concerned about STE fluxes directly but other factors that degrade the results: e.g., gradients across the tropopause, the effective tropospheric air mass diluting surface emissions, or the inclusion of CO₂-depleted stratospheric air in column CO₂ calculations (Nassar et al., 2011; Deng et al., 2015; Frankenberg et al., 2016; Pandey et al., 2016). For example, Le Quéré et al. (2018) are concerned with how emissions will mix throughout the troposphere and the stratosphere, but not how stratospheric air will come back down to the surface. Only studies of the CO₂ triple-oxygen isotope signature (Δ¹⁷O) are concerned with accurate STE fluxes, recognizing its importance in the seasonal isotopic signals (Liang et al., 2017; Koren et al., 2019; Laskar et al., 2019).

Both models and observations have shown that the stratospheric quasi-biennial oscillation (QBO) modulates the STE and drives much of the IAV observed in surface N₂O through the stratospheric influx of N₂O-depleted air (Hamilton and Fan, 2000; Nevison et al., 2004, 2011; Ray et al., 2020; Ruiz et al., 2021; Ruiz and Prather, 2022). Here, I use the N₂O studies of Ruiz et al. (2021) with parallel model simulations of CO₂ to place constraints on the CO₂ IAV caused by atmospheric circulation, finding that this effect is a clear but minor perturbation in driving the observed IAV of CO₂.

2 Methods and analysis

I investigate the CO₂ IAV and its causes using surface CO₂ observations from 1985 through 2020, surface N₂O observations from 1997 through 2020, and tracer simulations from the UC Irvine chemistry-transport model (CTM) simulations for the historical period 1990–2017.

To study the circulation-driven IAV of CO₂, including STE, I focus on the Southern Hemisphere (SH) because fluctuations in the large biosphere-driven seasonality in the Northern Hemisphere (NH) (Fig. 1a and b) will obscure any stratosphere-driven IAV there. The UCI CTM uses ECMWF integrated forecast fields at 1.1^∘ horizontal resolution and has proven quite successful in simulating the historical IAV of surface N₂O, ozone columns, and the Antarctic ozone hole (Ruiz et al., 2021; Ruiz and Prather, 2022; Tang et al., 2021). For CO₂, I develop two model scenarios to highlight the impacts of atmospheric transport. First, I define a surface-emissions-driven eCO2 scenario, in which the total atmosphere increases at a constant rate of 2 ppm yr⁻¹ (parts per million, dry-air mole fraction), driving a flux of about 2 Pg C yr⁻¹ into the SH. This eCO2 scenario is a simple experiment with area-uniform (20–60^∘ N) and time-constant emissions to test how atmospheric circulation driving a NH–SH gradient might affect the seasonal and interannual variability of SH surface CO₂. It is obviously not realistic, lacking the large biospheric and oceanic seasonality. A second stratospheric-driven sCO2 scenario is forced with a net stratospheric flux of CO₂-depleted air being transported into the troposphere and down to the surface. This STE flux is calculated as the equivalent of the aging of stratospheric CO₂ relative to the troposphere (2 ppm yr⁻¹), yielding an apparent negative CO₂ flux of about 0.4 Pg C yr⁻¹ into each hemisphere. This forcing flux is placed in the uppermost model layer (∼80 km altitude) and transported to the surface. In both of these cases, CO₂ changes linearly with a known trend, and I subtract that trend to get the modeled anomalies. The eCO2 scenario effectively forces the stratosphere with a negative flux of 2 ppm yr⁻¹, but most of the SH signal comes from the much larger interhemispheric flux. A third independent method for deriving CO₂ IAV uses the observed SH surface N₂O signal, driven by stratospheric photochemical loss of 13 Tg N (as N₂O) yr⁻¹, as a measure of STE influence. In this case I scale the results to CO₂ using the ratio of the STE fluxes, i.e., 0.15 ppm CO₂ per ppb N₂O. Ruiz et al. (2021, Figs. 3 and S3) show that the tropospheric QBO patterns for N₂O and CFCl₃ are nearly identical despite the different vertical locations and QBO patterns in their stratospheric loss. A species STE flux pattern (i) scales with the total flux out of the stratosphere and (ii) is determined by the dynamics of the lowermost midlatitude stratosphere (Ruiz and Prather, 2022, Fig. 1).

The monthly CO₂ surface observations are gathered from NOAA ESRL (Dlugokencky et al., 2021b). I use five sites: BRW – Barrow, AK (now Utqiaġvik), at 71^∘ N, 156^∘ W; MLO – Mauna Loa, HI, at 20^∘ N, 156^∘ W; SMO – Tutuila, Am. Samoa, at 14^∘ S, 171^∘ W; CGO – Cape Grim, Tasmania, Australia, at 41^∘ S, 144^∘ E; and SPO – South Pole at 90^∘ S. Monthly average in situ observations are used, and gaps are filled by flask data at the same site. CGO is flask only. I have a continuous monthly record from 1985 through 2020 (Fig. 1a). I convert these to a stationary series of residuals by fitting polynomials, assuming the months are equally spaced. The second-, third-, fourth-, and fifth-order polynomials produce almost identical results for each site (not shown), and the average third- and fourth-order fits are subtracted to calculate the residuals. The CO₂ residuals for the average of SPO+CGO are shown in Fig. 1c as the red line, which shows a clear annual cycle plus equally large variability on decadal scales. The annual cycle of CO₂ and its rate of change is a critical metric used to evaluate the carbon cycle in Earth system models (Graven et al., 2013; Zhao and Zeng, 2014; Wenzel et al., 2016). Here, I calculate the cycle simply by averaging each calendar month of the year, with results shown in Fig. 1b. The annual amplitudes (max–min) are 16.4, 6.4, 0.92, 1.02, and 1.14 ppm for BRW, MLO, SMO, CGO, and SPO, respectively. These are consistent with those previous studies, although SH cycles remain understudied and not well evaluated. In some months SMO at 14^∘ S can be north of the South Pacific Convergence Zone and thus influenced by NH air, explaining its non-sinusoidal cycle when compared with SPO and CGO. Also shown is the annual cycle for the modeled eCO2 scenario (∼0.18 ppm, dotted lines for SMO and SPO). That for the sCO2 scenario is even smaller (∼0.06 ppm) and is not shown. It is interesting that the SH annual cycles in eCO2 are similar in shape to those observed, even catching the double peak at SMO, but the magnitude is much smaller. There is no evidence in our direct modeling or analysis that stratosphere–troposphere exchange, which does drive an annual cycle in N₂O, can produce a detectable annual cycle in CO₂ above the large observed cycle.

The observed QBO signal in surface N₂O (Ruiz et al., 2021, Fig. 3) has the largest amplitude in the SH extratropics, becoming weaker in the tropics and NH. Because this signal is nearly uniform across the SH extratropics in both observations and models, I combine the SPO and CGO CO₂ data and focus our efforts on that time series. The challenge is to extract the CO₂ IAV signal in the 2–5-year period range. A simple 12-month running mean is great for removing the annual cycle but leaves the large-amplitude decadal periods (blue dashed line in Fig. 1c). I select bandpass filtering, while recognizing that this method can produce spurious results, especially at the edges. After several false starts and with help from the reviewers, I chose the MATLAB bandpass filter. This function is well documented (https://www.mathworks.com/help/signal/ref/bandpass.html, last access: 8 March 2022), and the bandpass is defined by the lower and upper cut-off frequencies. After experimentation with the CO₂ signal to reduce edge effects, I chose the following settings: bandpass frequency (yr⁻¹) range [0.20 0.80]; bandpass applied forward and backward is averaged; ImpulseResponse=iir; Steepness=0.85. A wider filter, e.g., [0.16 0.95], produced similar results for the middle years but large swings for the beginning and end years. IAV signals derived from frequency filtering for the last 2–3 years of the record are not robust. The resulting bandpass IAV (thick black line in Fig. 1c) clearly shows the patterns seen in the 12-month running mean. The SMO IAV is calculated in the same way and plotted alongside the SPO+CGO IAV in Fig. 1d. The IAVs for N₂O observations as well as the modeled eCO2 and sCO2 scenarios use the same processing.

For monthly N₂O surface observations, also from NOAA ESRL (Dlugokencky et al., 2021a), I focus on SH extratropics using SPO and CGO, plus three other sites: S30 – western Pacific cruise at 30^∘ S, 168^∘ E; USH – Tierra del Fuego, Ushuaia, Argentina, at 55^∘ S, 68^∘ W; and PSA – Palmer Station, Antarctica, at 65^∘ S, 64^∘ W. All five sites have nearly identical N₂O records, and I average them to get our SH IAV signal with the same processing as for CO₂. The QBO circulation is known to reach throughout the stratosphere and into the troposphere (Tung and Yang, 1994; Hamilton and Fan, 2000), and multi-model studies have attributed the surface N₂O IAV to the STE flux (Ruiz et al., 2021). I can thus scale the surface N₂O IAV with the ratio of STE fluxes (CO₂:N₂O) to give an observational estimate of the STE-driven CO₂ IAV in the SH extratropics (dashed blue line in Fig. 1d). The IAV in SH (40–90^∘ S) surface CO₂ calculated from the sCO2 model scenario is also shown (dashed red line in Fig. 1d). The modeled sCO2 and observed N₂O-scaled IAVs are not always in phase, but they are in strong agreement in terms of amplitude: the STE IAV in CO₂ is a small fraction of the observed IAV. In addition, the modeled eCO2 IAV shows that tropospheric circulation changes produce small IAV.

I compare CO₂ with well-known interannual cycles in the Earth system in Fig. 1d by plotting (i) the QBO phase change (from easterly to westerly zonal equatorial wind at 40 hPa, see Newman, 2021) as thick gray vertical bars and (ii) the times of moderate to extreme El Niños (red stars) and La Niñas (blue stars) (Trenberth, 2021). From this analysis, I expect minimal contribution of the QBO–driven circulation to the CO₂ IAV and find no obvious connection between the two in this figure. For the El Niño–Southern Oscillation (ENSO), this simple comparison is inadequate. At best it shows that some of the larger positive SH IAVs align with El Niños, whereas I know that ENSO affects ocean upwelling and continental rainfall, and the CO₂ anomalies correlate very well with tropical ocean temperatures (Wang et al., 2021; Keeling and Graven, 2021).

3 Conclusions, speculations, and digressions

I have shown that the STE fluxes of old stratospheric air with “depleted” CO₂ have little influence on the IAV or annual cycles of CO₂ at the surface. The IAV observed for SH stations has a standard deviation of 0.22 to 0.28 ppm, while that for sCO2 is at most 0.02 ppm in both hemispheres, that for eCO2 is less than 0.02 ppm for SPO to SMO, and that for scaled N₂O is 0.03 ppm for SPO+CGO average. The standard deviation for NH CO₂ IAV is larger at 0.4 to 0.6 ppm and thus even less influenced by stratospheric air. Thus, the speculations of Gaubert et al. (2019) regarding atmospheric transport can be dismissed.

The latitudinal pattern of N₂O IAV provides evidence for causes: e.g., the STE-driven signal weakens in the tropics and changes phase in the NH, and QBO composites show a clear separation of hemispheric sources (Ruiz et al., 2021). The latitudinal pattern of CO₂ IAV may similarly provide information on its cause. Comparing tropics to extratropics in the SH (SPO+CGO vs. SMO, solid and dotted black lines in Fig. 1d), I find remarkably similar patterns after 1990, with similar amplitudes and some phase shifts of at most 1 year. If I add the NH tropics MLO IAV (not shown), the pattern and amplitude are similar. When the sites are in sync, one can only presume that the CO₂ perturbation is tied to changes in the growth and decay of tropical biomass transported equally to both hemispheres (Keeling and Graven, 2021). The challenge lies in the phasing and which region leads or lags in change. Unfortunately, the bandpass IAVs in this analysis do not seem able to accurately determine the phase at a level up to 1 year.

The Samoan site SMO provides a valuable but very challenging record of CO₂ and other trace gases having dominant NH emissions, such as chlorofluorocarbons (Cunnold et al., 1994) and N₂O (Nevison et al., 2007). Sometimes SMO is synchronous with the SH extratropics (CGO and SPO, which are almost always synchronous with each other) and at other times it links with MLO and the NH. Thus, to use SMO CO₂ as a metric for carbon cycle models, one must recognize that SMO is not simply representative of the SH tropics. When evaluating carbon cycle models, one should test tracer transport using the IAV for SMO versus SPO+CGO. As shown in Fig. 1d, there are clear times when SMO is distinct from CGO+SPO (e.g., 1994, 1999, 2008 2012, 2015). At these times the SMO IAV matches that of MLO (not shown). These interannual shifts provide an excellent test for CO₂ historical simulations using weather forecasting systems (e.g., McNorton et al., 2020) and realistic sources and sinks (e.g., Piao et al., 2018; Wang et al., 2020).

The rate of increase of the amplitude of the annual cycle of CO₂ is a key measure of changes in the biospheric and oceanic carbon cycles. The amplitude can be measured from the variance across 12 months. If the cycle is sinusoidal, then the max–min amplitude is equal to twice the square root of 2 times the standard deviation as plotted in Fig. 1e. With this approximation, I calculate a mean amplitude of 17.0, 6.4, 1.3, 1.2, and 1.3 ppm for BRW, MLO, SMO, CGO, and SPO, respectively. These results are almost identical to those from the composited annual cycles (Fig. 1b) but disagree at SMO as might be expected because of its double-peaked cycle. A linear fit to the standard deviations gives trends for the period 1985–2020 of 1.06±0.15, 0.142±0.075, 0.082±0.047, 0.079±0.051, and 0.031±0.050 ppm per decade for BRW, MLO, SMO, CGO, and SPO, respectively. The standard errors quoted here come from a standard linear fit of the monthly values shown in Fig. 1e but are calculated more conservatively using 35 years as the degrees of freedom instead of 420 months. Our results for BRW and MLO agree well with other more extensive data analyses (Graven et al., 2013; Zhao and Zeng, 2014; Wenzel et al., 2016; Piao et al., 2018; Wang et al., 2020) but are able to identify emergent trends in the SH, which is not often used for model evaluation. A more serious uncertainty analysis focusing on the SH sources and sinks, the annual and IAV cycles, and their trends would help solidify our knowledge of the carbon cycle.