<figure><img src="image_1.png" /></figure>

HAL
open science

# Response of the Quasi-Biennial Oscillation to a warming climate in global climate models

Jadwiga Richter, Neal Butchart, Yoshio Kawatani, Andrew Bushell, Laura

Holt, Federico Serva, James Anstey, Isla Simpson, Scott Osprey, Kevin

Hamilton, et al.

► To cite this version:

Jadwiga Richter, Neal Butchart, Yoshio Kawatani, Andrew Bushell, Laura Holt, et al. Response of the Quasi-Biennial Oscillation to a warming climate in global climate models. Quarterly Journal of the Royal Meteorological Society, 2020, 10.1002/qj.3749. hal-03049568

HAL Id: hal-03049568

https://hal.science/hal-03049568v1

Submitted on 9 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

ORIGINAL ARTICLE

Journal Section

Response of the quasi-biennial oscillation to a
warming climate in global climate models

Jadwiga H. Richter¹

Neal Butchart²

Yoshio

Kawatani³

Andrew C. Bushell⁴

Laura Holt⁵

Federico Serva⁶

James Anstey⁷

Isla R. Simpson¹

Scott Osprey⁸

Kevin Hamilton⁹, ¹⁰

Peter

Braesicke¹¹

Chiara Cagnazzo⁶

Chih-Chieh Chen¹

Rolando R. Garcia¹

Lesley J. Gray⁸

Tobias

Kerzenmacher¹¹

Francois Lott¹²

Charles

McLandress⁸, ¹³

Hiroaki Naoe¹⁴

John Scinocca⁷

Timothy N. Stockdale¹⁵

Stefan Versick¹¹

Shingo

Watanabe³

Kohei Yoshida¹⁴

Seiji Yukimoto¹⁴

Author
of
Manuscript

This is the author manuscript accepted for publication and has undergone full peer review and has not been through the copyediting, typesetting, pagination and proofreading processes that may lead to differences between this version and the Version of Record. Please cite this as doi: 10.1002/qj.3749
### Correspondence

Jadwiga H. Richter, Climate and Global Dynamics Laboratory, National Center for Atmospheric Research, Boulder, CO, 80307, USA

Email: jrichter@ucar.edu

### Funding information

We compare the response of the quasi-biennial oscillation (QBO) to a warming climate in eleven atmosphere general circulation models that performed time-slice simulations for present-day, doubled, and quadrupled CO₂ climates. No consistency was found among the models for the QBO period response, with the period decreasing by eight months in some models and lengthening by up to thirteen months in others in the doubled CO₂ simulations. In the quadruped CO₂ simulations a reduction in QBO period of 14 months was found in some models, whereas in several others the tropical oscillation no longer resembled the present day QBO, although could still be identified in the deseasonalized zonal mean zonal wind timeseries. In contrast, all the models projected a decrease in the QBO amplitude in a warmer climate with the largest relative decrease near 60 hPa. In simulations with doubled and quadrupled CO₂ the multi-model mean QBO amplitudes decreased by 36% and 51%, respectively. Across the models the differences in the QBO period response were most strongly related to how the gravity wave momentum flux entering the stratosphere and tropical vertical residual velocity responded to the increases in CO₂ amounts. Likewise it was found that the robust decrease in QBO amplitudes was correlated across the models to changes in vertical residual velocity, parameterized gravity wave momentum fluxes, and to some degree the resolved

Author

<footnote n="1">¹Climate and Global Dynamics Laboratory, National Center for Atmospheric Research, Boulder, CO, 80307, USA</footnote>

<footnote n="2">²Met Office Hadley Centre, Exeter, EX1 3PB, UK</footnote>

<footnote n="3">³Japan Agency for Marine-Earth Science and Technology, Yokohama, Japan</footnote>

<footnote n="4">⁴Met Office, Exeter, UK</footnote>

<footnote n="5">⁵North West Research Associates (NWRA), Boulder, USA</footnote>

<footnote n="6">⁶Institute of Marine Sciences, National Research Council (ISMAR-CNR), Rome, Italy</footnote>

<footnote n="7">⁷Canadian Centre for Climate Modelling and Analysis (CCCma), Victoria, Canada</footnote>

<footnote n="8">⁸National Centre for Atmospheric Science (NCAS), University of Oxford, Oxford, UK</footnote>

<footnote n="9">⁹Department of Atmospheric Sciences, University of Hawaii, Honolulu, USA</footnote>

<footnote n="10">¹⁰International Pacific Research Center (IPRC), Honolulu, USA</footnote>

<footnote n="11">¹¹Karlsruher Institut für Technologie (KIT), Karlsruhe, Germany</footnote>

<footnote n="12">¹²Laboratoire de Météorologie Dynamique (LMD), Paris, France</footnote>

<footnote n="13">¹³University of Toronto, Toronto, Canada</footnote>

<footnote n="14">¹⁴Meteorological Research Institute (MRI), Tsukuba, Japan</footnote>

<footnote n="15">¹⁵European Centre for Medium-Range Weather Forecasts (ECMWF)</footnote>
upward wave flux. We argue that uncertainty in the representation of the parameterized gravity waves is the most likely cause of the spread among the eleven models in the QBO's response to climate change.

KEYWORDS

QBO, QBOi, Climate Change, GCMs, Stratosphere, Gravity Waves

# 1 | INTRODUCTION

The quasi-biennial oscillation (QBO) is the most prominent feature of the circulation of the equatorial stratosphere. Alternating descending easterly and westerly shear zones have been observed consistently near the equator between 100 and 10 hPa since 1953 with an average period of 28 months (Naujokat, 1986; Baldwin et al., 2001). The QBO has been shown to affect the strength of the stratospheric polar vortex (Holton, 1980), extratropical surface variability (Anstey and Shepherd, 2014), tropical convection (Collimore et al., 2003), the Madden-Julian Oscillation (MJO) (Son et al., 2017; Yoo and Son, 2016), Pacific tropical cyclones (Ho et al., 2009) and Atlantic tropical cyclones over a portion of the observational record (Camargo and Sobel, 2010). The QBO also modulates the tropical tropopause (Reid and Gage, 1985), stratospheric water vapor, ozone, and methane (Randel et al., 1998). The QBO is understood to be maintained largely by momentum deposition from upward propagating large scale Kelvin and mixed-Rossby gravity waves, gravity waves with horizontal scales of tens to thousands of kilometers, and vertical advection. Variability in these three forcing terms contributes to the variable period of the QBO, which ranges from 22 to 34 months in the observational record (Baldwin et al., 2001). Although the basic mechanism of QBO formation is well understood, many general circulation models (GCMs) have difficulty reproducing this oscillation. For most GCMs, the ability to reproduce the QBO is strongly dependent on adequate vertical resolution, especially in the upper troposphere and lower stratosphere (Giorgetta et al., 2006; Richter et al., 2014; Geller et al., 2016a), generation of Kelvin and mixed-Rossby gravity waves, and the parameters of the gravity wave (GW) scheme, which are poorly constrained by observations. The partitioning of QBO forcing between large scale and smaller-scale convectively generated waves is not precisely known and can only be indirectly inferred from observations and models. Kelvin and mixed-Rossby gravity waves in GCMs are generated
largely by parameterized convection and vary greatly among the different models (Lin et al., 2006; Lott et al., 2014). GW parameterizations are typically tuned to obtain the correct present-day period and amplitude of the oscillation, and hence together with the forcing provided by resolved waves are a very large source of uncertainty in modeling of the QBO.

Due to the numerous impacts of the QBO on the tropospheric and the stratospheric circulation, it is important to understand how the QBO will change in a future climate. This will depend on how vertical advection, convectively coupled waves, and convectively generated gravity waves change and alter the wave-mean flow interactions. There is observational evidence that changes in planetary-scale Rossby wave propagation can lead to surprising effects on the QBO, like the unexpected breakdown of the QBO cycle in 2016 (Osprey et al., 2016; Coy et al., 2017). Several studies have addressed the question of how the QBO is likely to change in a future climate, however the findings have been inconclusive, primarily due to the large uncertainty associated with parameterized gravity waves. The model study of the QBO changes in a warm climate by Giorgetta and Doege (2005) showed that the QBO period decreases in doubled CO₂ simulations from 29 to 26, 22, and 17 months in experiments in which gravity wave source remained the same, and increased by 10%, and then 20% respectively. The QBO period reduction in the warming climate simulation resulted from both the prescribed increase of wave sources and a simulated decrease in the tropical upwelling though the latter is not in agreement with most model projections of increased tropical upwelling in a warmer climate (Butchart et al., 2006). Kawatani et al. (2011) used a model without parameterized non-orographic gravity waves and found that the effect of enhanced mean tropical upwelling in a warming climate overwhelms the counteracting influence from strengthened wave fluxes. Consequently, the amplitude of the QBO becomes smaller, especially in the lower stratosphere and the period becomes longer. Watanabe and Kawatani (2012) also projected that the QBO will lengthen and the amplitude become smaller in a warming climate albeit in a model in which sources of parameterized gravity waves remained fixed and the change in the QBO resulted from changes to the vertical residual velocity. Kawatani and Hamilton (2013) analyzed four Coupled Model Intercomparison Project phase 5 (CMIP5) models that could simulate a spontaneous QBO. They showed that a long term reduction in the QBO amplitude in the lower stratosphere is robust in all four models and they also found this clearly in observations at 70 hPa from 1953 to 2012. On the other hand, the four models produced different projections of the QBO period changes, while the 60-year observational record showed no significant trends in QBO period. Schirber et al. (2015) found that the response of the amplitude and period of the QBO to a warming climate was dependent on the configuration of gravity wave parameterizations in the ECHAM model.

The lack of agreement among these earlier studies on how the QBO may change in the future and the need to under-
stand how assumptions in gravity wave parameterizations influence these findings helped motivate the Stratosphere-troposphere Processes And their Role in Climate (SPARC) Quasi-Biennial Oscillation initiative (QBOi) designed to improve the fidelity of tropical stratosphere variability in present day GCMs (Butchart et al., 2018). For phase 1 of the QBOi, five experiments were designed in combination to: evaluate the ability of GCMs to represent the QBO in present day climate (QBOi Experiments 1 and 2), examine the response of the QBO under climate-change forcings in various GCMs (QBOi Experiments 2, 3, and 4), and evaluate the ability of GCMs to predict the QBO when initialized with observations (QBOi Experiments 5 and 5a). Seventeen GCMs participated in phase 1 of the QBOi and carried out some or all of the experiments (Butchart et al., 2018). In a companion paper, Bushell et al. (2019) describe the characteristics of the QBO in the present day climate among these models. Holt et al. (2019) evaluates the resolved tropical waves and wave forcing of the QBO in these models and QBO teleconnections will be discussed in Anstey et al. (2019). In this paper we discuss the response of the QBO in a doubled and quadrupled CO₂ climate among eleven GCMs that participated in Experiments 2, 3, and 4. We examine here the robustness of changes in the QBO and changes in QBO forcings among these models, as well as discuss the factors that contribute to the largest uncertainty in the response of the QBO. The details of changes in the momentum budget in the three experiments will be reported in a follow-up publication. QBO teleconnections in future climate, although crucial to understand, are beyond the scope of this paper. The paper is organized as follows: the models, experimental set-up, and changes to the mean state of the stratosphere are described in section 2. Section 3 describes changes to the QBO. Section 4 discusses changes in QBO forcing terms, and section 5 presents the summary and conclusions.

## 2 | MODELS, EXPERIMENTAL SETUP, AND MEAN CLIMATE

### 2.1 | Models

We utilize here output from eleven atmosphere only GCMs that participated in the QBOi and performed Experiments 2, 3 and 4. The participating models, GW parametrization specifications, and the length of simulations in each experiment are listed in Table 1. For references to model descriptions, comparison of vertical/horizontal resolutions, time steps, convection parameterizations, and numerical advection schemes among the models, the reader is referred to Butchart et al. (2018). Note that AGCM3-CMAM did not perform Exp 4 because it was discovered that the radiative transfer parameterization in the non-local thermodynamic equilibrium region may be erroneous when CO₂ levels are larger
106

107

than twice the present-day concentration. Although this error might only affect temperatures at very high altitudes, it was decided to not submit model runs in which the radiative physics contained a known error.

108

109

110

111

112

113

114

115

116

117

118

119

120

The majority of the QBOi models have fixed gravity wave sources, meaning that the spectrum of GWs launched at the source level does not change as the climate changes and is fixed (as are all other model characteristics) between Experiments 2, 3, and 4. The only exception to that are models using the Hines (1997a,b) scheme (ECHAM5sh, EMAC, MIROC-ESM, and MRI-ESM2), for which the GW wavenumber spectrum, and in particular the so-called "cutoff wavenumber" depends on the buoyancy frequency at the source level, which varies between the experiments, resulting in a change in total momentum flux at source level. Four of the eleven models (60LCAM5, CESM1(WACCM5-110L), LMDz6, and UMGA7gws) have gravity wave source parameterizations that link the properties of non-orographic GWs to the properties of their sources in the troposphere (described in more detail in Section 4.2). UMGA7 and UMGA7gws differ only in the representation of parameterized gravity waves, with UMGA7 including a fixed and UMGA7gws including a variable GW source spectrum. One of the models, MIROC-AGCM-LL, is able to reproduce the QBO without parameterizing non-orographic GWs. The cumulus parameterization in MIROC-AGCM-LL is based on the method of Arakawa and Schubert (1974), with a relative humidity limit method. Both this method and high vertical resolution enable the MIROC-AGCM-LL to simulate the QBO without a non-orographic GW parameterization.

121

## 2.2 | Experimental Setup

122

The simulation setup for all models followed the QBOi protocol (Butchart et al., 2018):

123

124

125

126

127

128

* Experiment 2 (hereafter 'Exp 2') is a present-day time-slice experiment in which a repeated annual cycle of sea surface temperatures (SSTs) and sea-ice is employed, and 2002 external forcings.

* Experiment 3 (hereafter 'Exp 3') is the same as Exp 2 but with a 2 K temperature perturbation added uniformly to the SSTs and double the CO₂ concentration of Exp 2.

* Experiment 4 (hereafter 'Exp 4') is the same as Exp 2 but with a 4 K temperature perturbation added uniformly to the SSTs and quadruple the CO₂ concentration of Exp 2.

129

130

No changes to the sea-ice cover were made in Exp 3 and Exp 4 despite the higher SSTs. Full details of the QBOi Experiments can be found in Butchart et al. (2018).
<table><caption>**TABLE 1** Models participating in QBOi Experiments 2, 3 and 4 and specifications of their non-orographic gravity wave parameterizations (NOGW Param), GW source parameterization (NOGW Source), GW Source Level (GW SL), and number of simulation years in Exp 2, 3, and 4. Some models performed several ensembles indicated by the number 'x' in the Yrs columns. AGCM3-CMAM did not perform Exp 4. Non-orographic GW parameterizations are abbreviated as follows: Hines (1997a,b) [Hi]; Warner and McIntyre (1999) [WM]; Lindzen (1981) [Li]; Lott et al. (2012) [Lo]).</caption><thead><tr><th rowspan="2">Model</th><th>NOGW</th><th>NOGW</th><th>GW</th><th colspan="3">Yrs</th></tr><tr><th>Param</th><th>Source</th><th>SL</th><th>Exp 2</th><th>Exp 3</th><th>Exp 4</th></tr></thead><tbody><tr><td>60LCAM5</td><td>[Li]</td><td>Richter et al. (2010)</td><td>Variable</td><td>3x38</td><td>38</td><td>38</td></tr><tr><td>AGCM3-CMAM</td><td>[WM]</td><td>Fixed</td><td>90 hPa</td><td>3x31</td><td>3x31</td><td>0</td></tr><tr><td>CESM1(WACCM5-110L)</td><td>[Li]</td><td>Richter et al. (2010)</td><td>Variable</td><td>33, 2x35</td><td>31</td><td>31</td></tr><tr><td>ECHAM5sh</td><td>[Hi]</td><td>Fixed</td><td>600 hPa</td><td>30</td><td>30</td><td>30</td></tr><tr><td>EMAC</td><td>[Hi]</td><td>Fixed</td><td>643 hPa</td><td>106</td><td>106</td><td>106</td></tr><tr><td>LMDz6</td><td>[Lo]</td><td>Lott and Guez (2013)</td><td>500 hPa</td><td>70</td><td>70</td><td>70</td></tr><tr><td>MIROC-AGCM-LL</td><td>None</td><td>N/A</td><td>N/A</td><td>3x30</td><td>3x30</td><td>3x30</td></tr><tr><td>MIROC-ESM</td><td>[Hi]</td><td>Fixed</td><td>657 hPa</td><td>100</td><td>3x30</td><td>3x30</td></tr><tr><td>MRI-ESM2</td><td>[Hi]</td><td>Fixed</td><td>Surface</td><td>30</td><td>30</td><td>30</td></tr><tr><td>UMGA7</td><td>[WM]</td><td>Fixed</td><td>4 km</td><td>105</td><td>105</td><td>105</td></tr><tr><td>UMGA7gws</td><td>[WM]</td><td>Bushell et al. (2015)</td><td>4 km</td><td>105</td><td>105</td><td>105</td></tr></tbody></table>

131

## 2.3 | Mean Climate

132

Figure 1 shows the multi-model mean change in the annual- and zonal-mean zonal wind and temperature as a function of latitude and height for Exp 3 and 4 relative to Exp 2. All models are given equal weight in the presented multi-model means. Consistent with the basic radiative effects of CO₂ and increased SSTs, the temperature increases in the troposphere and decreases in the stratosphere in Exp 3 and Exp 4 as compared to Exp 2, both for the multi-model mean (Figure 1a, b) and also for each of the models (not shown). The multi-model response of the zonal mean zonal wind to increased CO₂ concentration and SSTs shows poleward shifting and upward expansion of the mid-latitude jets in the troposphere and strengthening of the extra-tropical jets throughout the stratosphere, especially in the Southern Hemisphere (Figure 1c, d), broadly similar to previous studies (Barnes and Polvani, 2013; Shepherd and McLandress, 2011; Simpson et al., 2014).

133

134

135

136

137

138

139

140

141

In the tropical troposphere, temperature and static stability increase, on average, in Exp 3 and 4 as compared to Exp 2 (Figure 2), while in the lower stratosphere (90 to 40 hPa) temperature changes are small and not significant compared to the inter model spread in either Exp 3 or 4. All the models agree on the warming in the troposphere and cooling of the

142

143

<figure><img src="image_4.png" /></figure>
stratosphere, but the exact altitude where they transition between this warming and cooling varies across models (not shown). Although there is general agreement in how the tropospheric temperature structure will change in a warming climate between the QBOi models, Figure 2b showing the inter-model standard deviation of tropical temperature shows that the spread between models increases vastly between 300 and 100 hPa. In the mean, the decrease in static stability in the equatorial lower stratosphere in Exp 3 and 4 relative to Exp 2, shown in Figure 2c, is robust across models, however as will be discussed later, the increased spread in upper tropospheric temperatures (and hence static stability) in the tropics can be a factor causing increased spread in QBO characteristics, as reduced static stability may affect the wave forcing of the QBO.

Precipitation is a key variable when considering climate change and the QBO as precipitation patterns and variability influence the convectively coupled Kelvin and mixed-Rossby gravity waves. Spatial patterns of precipitation vary a lot between the QBOi models and hence produce very different spectra of large scale waves (discussed in a companion paper, Holt et al. (2019)). Figure 3a shows the zonal mean precipitation rate for Exp 2 for the QBOi models, illustrating substantial differences between models in the present day climate, with peak precipitation latitude not even being in the same hemisphere in all models. However, changes in precipitation from Exp 2 to Exp 3 and to Exp 4 are rather small, 2.4 and 4.1% in a multi-model mean respectively, and the inter-model spread stays very similar: the across model standard deviation of precipitation rate is: 0.30, 0.32, and 0.28 mm day⁻¹ for Exp 2, 3, and 4 respectively.

3 | QBO

Figure 4 shows a 20-year time series of near-equatorial zonal mean winds for Exp 2, 3, and 4. For current SSTs and amounts of CO₂ (Exp 2) the left hand panels confirm that all the models capture the characteristic features of the observed QBO, including alternating layers of descending easterlies and westerlies in the stratosphere. There is still a spread among models in terms of QBO period and amplitude as described in Bushell et al. (2019). Figure 5 shows the results of the same experiments but with the mean seasonal cycle subtracted from the zonal mean winds as was done in the QBO review paper of Baldwin et al. (2001). After removing the seasonal cycle, the vertical structures of the QBOs in Exp 2 in the various models more closely resemble each other, and regular downward phase progression all the way to the tropopause in each model is quite apparent.

As SSTs increase and the amount of CO₂ is doubled (Exp 3, middle column Figure 4), then quadrupled (Exp 4,
right-most column Figure 4), the characteristic QBO signal weakens in all models in the raw time series, and in some cases is no longer discernible for Exp 4 (e.g. ECHAM5sh). However, if the seasonal cycle is subtracted from the zonal mean winds (Figure 5), then a signal of descending layers of eastward and westward wind shear can be identified in all models for the three experiments as confirmed by the lag-correlations (see below for details) shown in Figure 6 for Exp 4. This raises a question of how the QBO should be defined: as the occurrence of downward propagating alternating westerly and easterly shear zones or as downward propagating alternating westerly and easterly shear zones after the seasonal cycle has been removed? In observations and in the present day climate simulations of the QBO, either definition will result in a QBO (with quite similar characteristics), and some past studies of QBO teleconnections have defined QBO easterly and westerly phases based on raw and some on deseasonalized time series of the zonal averaged zonal mean wind (e.g.: Anstey and Shepherd (2014) vs Garfinkel and Hartmann (2010); Yoo and Son (2016)). Here, we consider both definitions of the QBO to be inclusive in our analysis of models that show a QBO (downward propagating easterly and westerly shear zones) only in the deseasonalized time series. The presence of a QBO in the deseasonalized, and not in the raw time series, suggests that QBO generation maybe different in models than in observations, and different from our current understanding of QBO, which relies on absorption of waves, Kelvin, mixed-Rossby GWs and gravity waves, with positive and negative phase speeds (e.g: Lindzen and Holton (1968), Ern et al. (2014)). In simulations in which the sign of the raw mean wind does not change but there is an oscillation in the deseasonalized time series, the downward propagation of QBO phases may no longer be caused by wave absorption with appreciable contributions from both eastward and westward waves, but still can be caused by a wave spectrum with a range of phase speeds not centered near zero phase speed.

A robust and systematic decrease in the QBO amplitude with increasing SSTs and CO₂ concentration can be seen in Figures 4 and 5. The sign of change in the QBO period, on the other hand is model dependent, with some models projecting a decrease, some projecting a lengthening, and others predicting little change. These differences in behaviour are consistent across Exp 3 and Exp 4 and do not appear to be connected to any particular aspect or combination of aspects of model formulation and in particular any of those cataloged in Table 1. However, we note that given the small size of the multi-model ensemble (11 models) as compared to the number of differences in model formulation, it can be difficult to robustly identify any aspect of model formulation as being responsible for the different responses in QBO period.

The qualitative conclusions obtained from Figures 3 and 4 are further quantified below using metrics established in Bushell et al. (2019) for QBOi Experiment 1 (present-day forcings with varying SSTs), augmented by additional

<figure><img src="image_5.png" /></figure>
diagnostics that help elucidate QBO characteristics in simulations for which the oscillation is more erratic than the present-day QBO. The need for additional diagnostics and metrics for Exp 3 and Exp 4 simulations suggests that the QBOs in models are overly tuned for the present day climate.

## 3.1 | QBO Periods:

We derive QBO periods using three methods, defined as follows:

* 'QBO Transition Times (TT)': following Bushell et al. (2019), zero wind transitions from westward to eastward winds are identified at a reference level of 10 hPa and used to define QBO cycles. Wind transitions are identified in an equatorial averaged zonal mean zonal wind time series between 5°S to 5°N at 10 hPa which was smoothed using a 5-month centered running mean. A QBO period and amplitude for each QBO cycle is identified from the zonal mean wind time series between individual transition points. Using this method, the standard deviation of the periods and amplitudes can then be derived from these distributions in addition to the mean.

* 'Fast Fourier Transform (FFT)’: analysis was performed on a 5°S to 5°N averaged zonal mean zonal wind deaseasonalized time series which was padded at the end with zeros to 10 times the length of original time series to increase the spectral resolution. FFT analysis was performed at levels between 50 and 10 hPa. In Exp 2, the FFT spectra typically have one dominant peak, and the Fourier amplitude is highest at exactly the same period at all levels between 50 and 10 hPa. In Exp 3 and Exp 4 the FFT spectra in many models no longer have a well-defined peak and often have several peaks of similar amplitude. In addition, in some simulations the dominant QBO period varies with altitude (not shown). A range of QBO periods is provided for simulations in which the dominant FFT period between 50 and 10 hPa varies by more than a month.

* 'Lag correlations (LAG)’: The correlation of the time series of equatorial zonal mean wind at 10 hPa with the corresponding lagged time series (i.e. at $t \pm lag$) at each level confirms the presence of an oscillation for all models, with increasing phase lag on descending through the stratosphere (see Figure 6 for Exp 4). In most models a well defined mean period can be deduced from the correlations [e.g. CESM1(WACCM5-110L)] while in other models (e.g. ECHAM5sh and UMGA7gws) the point of maximum correlation is followed, or preceded, by more than one local maximum anti-correlation. The larger of these local maxima, as indicated by the black vertical lines in the figure is used to define the mean period used in the rest of the paper. The absence of an unambiguous period in the
lag-correlations is consistent with the less clear QBO signals in Figures 4 and 5 and the absence of a well defined peak in the power spectra.

The findings from the above diagnostics are summarized in Table 2. Downward propagating easterly and westerly shear zones in the deseasonalized zonal mean zonal wind time series can be identified in all the simulations, however in several of the models in Exp 3 and Exp 4 the downward propagation is less regular, and often interrupted. The changed nature of tropical variability and lack of consistent downward propagation in Exp 3 and Exp 4 is reflected by a discrepancy in the mean period derived by the three methods and/or a range of periods between 50 and 10 hPa. These include ECHAM5sh in Exp 3, and ECHAM5sh, MIROC-ESM, and UMGA7gws in Exp 4. Table 2 (and Table 3) include ERA-Interim reanalysis (ERAI) (Dee et al., 2011) derived values for years 1979 - 2008 for reference. In Exp 2 most models simulate a lower range of QBO periods than in ERAI possibly due to the non-annually varying SSTs and fixed external forcings Bushell et al. (2019).

Changes in QBO periods derived using the TT method are visualized in Figure 7. Using this period-derivation method allows for visualizing the range of QBO periods in addition to the mean and standard deviation. Three models (60LCAM5, CESM1(WACCM-110L), and MRI-ESM2) predict much shorter QBO periods with increasing SSTs/CO₂ in Exp 3 and Exp 4 with predicted mean periods between 15 and 18 months in Exp 4. ECHAM5sh, EMAC, UMGA7, UMGA7gws predict no significant change in the mean QBO period in Exp 3, whereas AGCM3-CMAM, LMDz6, MIROC-AGCM-LL, and MIROC-ESM project lengthening of the QBO period. In Exp 4, LMDz6, MIROC-AGCM-LL, and UMGA7gws predict statistically significant lengthening of the QBO ranging from 3 to 10 months relative to present day. ECHAM5sh, EMAC, MIROC-ESM, and UMGA7 predict no significant change in QBO period in Exp 4. In a multi-model mean (Table 2), the QBO period does not change between the experiments, however as seen from Figure 7 that is not representative of most individual models. Figure 7 and standard deviations in Table 2 for individual models and multi-model mean show that in general the range and variability of simulated QBO periods increases in the warming climate, especially in Exp 4, suggesting that the oscillation is more sporadic in a warmer climate. The variability in Exp 4 remains comparable to or slightly smaller than the variability in Exp 2 for CESM1(WACCM-110L) and MRI-ESM2, which are models that have gone to substantially shorter periods in a warming climate. The more variable nature of the near-equatorial oscillation is also revealed by the FFT analysis, which has shown that spectra tend to have multiple peaks, especially in Exp 4, whereas most models had very well defined spectral peaks in Exp 2 (not shown).
<table><caption>**TABLE 2** Summary of estimates of QBO periods in Exp 2, 3, 4 using three different methods: TT, FFT, and LAG. For the TT method the standard deviation follows the ± sign. For simulations in which the dominant period varies by more than 1 month between 50 and 10 hPa a range of dominant periods is specified. Simulations for which the TT and LAG methods of period estimation do not agree to within 4 months are marked with a *. Period values in Exp 3 and Exp 4 derived using the TT method that are statistically different from Exp 2 at the 95% level as defined by the Student t-test for unequal sized distributions are bolded. ** indicates that for ECHAM5sh transition times were calculated at 15 hPa, as the oscillation was poorly defined at 10 hPa and a clear period was difficult to obtain. The multi-model mean is listed in the bottom row of the table. For periods derived using the FFT method in Exp 3 and Exp 4 the multi-model mean is not defined (ND) due to the range of dominant periods for some of the models. Values for ERAI for years 1979 - 2008 are shown for reference in the bottom row.</caption><thead><tr><th rowspan="2">Model/Reanalysis</th><th colspan="3">Exp 2</th><th colspan="3">Exp 3</th><th colspan="3">Exp 4</th></tr><tr><th>TT</th><th>FFT</th><th>LAG</th><th>TT</th><th>FFT</th><th>LAG</th><th>TT</th><th>FFT</th><th>LAG</th></tr></thead><tbody><tr><td>60LCAM5</td><td>24.4±1.7</td><td>24.3</td><td>27</td><td><strong>18.5±5.2</strong></td><td>20.3</td><td>18</td><td><strong>14.2±3.6</strong></td><td>14.8</td><td>16</td></tr><tr><td>AGCM3-CMAM</td><td>27.0±1.8</td><td>27.1</td><td>27</td><td><strong>41.2±8.0</strong></td><td>40.0-53.0</td><td>42</td><td>N/A</td><td>N/A</td><td>N/A</td></tr><tr><td>CESM1(WACCM5-110L)</td><td>27.6±2.5</td><td>27.5</td><td>28</td><td><strong>19.7±1.3</strong></td><td>19.8</td><td>20</td><td><strong>14.0±2.0</strong></td><td>14.5</td><td>14</td></tr><tr><td>ECHAM5sh</td><td>26.6±1.8</td><td>26.9</td><td>26</td><td>27.0*±5.7</td><td>23.4-39.1</td><td>32*</td><td>32.2±13.6</td><td>25.5-44</td><td>32</td></tr><tr><td>EMAC</td><td>26.6±2.0</td><td>26.7</td><td>26</td><td>27.7±3.5</td><td>27.1</td><td>28</td><td>27.6±4.0</td><td>27.8</td><td>28</td></tr><tr><td>LMDz6</td><td>27.4±2.7</td><td>27.3</td><td>26</td><td><strong>33.0±4.8</strong></td><td>34.7</td><td>32</td><td><strong>35.9±17.7</strong></td><td>32.1</td><td>36</td></tr><tr><td>MIROC-AGCM-LL</td><td>19.3±1.5</td><td>19.3</td><td>19</td><td><strong>20.7±2.0</strong></td><td>20.8</td><td>20</td><td><strong>21.9±4.5</strong></td><td>23.3</td><td>21</td></tr><tr><td>MIROC-ESM</td><td>24.2±2.6</td><td>24.1</td><td>24</td><td><strong>26.6±3.0</strong></td><td>25.4</td><td>27</td><td>22.6±8.8</td><td>26.1-34.6</td><td>32*</td></tr><tr><td>MRI-ESM2</td><td>22.5±3.6</td><td>22.5</td><td>22</td><td><strong>15.6±4.1</strong></td><td>19.7</td><td>15</td><td><strong>14.0±3.9</strong></td><td>17.5</td><td>13</td></tr><tr><td>UMGA7</td><td>26.0±2.6</td><td>26.9</td><td>26</td><td>26.0±4.6</td><td>26.6</td><td>26</td><td>28.7±8.6</td><td>30.3</td><td>29</td></tr><tr><td>UMGA7gws</td><td>24.0±2.1</td><td>23.9</td><td>24</td><td>24.2±5.2</td><td>25.3</td><td>28</td><td>22.1±9.5</td><td>32.7-94.0</td><td>34*</td></tr><tr><td>Multi-Model Mean</td><td><strong>25.0±2.5</strong></td><td><strong>25.1</strong></td><td><strong>25</strong></td><td><strong>25.5±4.3</strong></td><td><strong>ND</strong></td><td><strong>26</strong></td><td><strong>23.3±7.6</strong></td><td><strong>ND</strong></td><td><strong>25</strong></td></tr><tr><td>ERAI</td><td><strong>27.8±3.6</strong></td><td>28.1</td><td>28</td><td>N/A</td><td>N/A</td><td>N/A</td><td>N/A</td><td>N/A</td><td>N/A</td></tr></tbody></table>

252

## 3.2 | QBO Amplitudes

253

QBO amplitudes are derived using two methods:

254

* 'TT': following from the QBO transition times method for calculating QBO periods described in the previous section. Amplitudes associated with each phase at 10 hPa were estimated by taking the maximum or minimum value within the identified time period between zero transitions in the original timeseries with higher frequency variability reduced through a 5-month centered binomial smoothing following Bushell et al. (2019).

255

256

257

* 'Dunkerton and Delisi (DD)’: Dunkerton and Delisi (1985) observed that, when most of the variability in the monthly mean equatorial winds results from the QBO, the mean amplitude of the oscillation can be approximated by $\sqrt{2\sigma}$, where $\sigma$ is the standard deviation of the deseasonalized time series of the monthly mean eastward winds. This

258

259

260
261

262

263

264

allows the mean amplitudes, though not the distribution of amplitudes, to be easily calculated as functions of latitude and height and these agree well with the mean amplitudes obtained from the TT method for simulations with well defined QBOs. Furthermore, even if clearly defined QBO cycles are not present, the diagnostic can still be usefully interpreted as a general quantitative measure of the amount of variability present.

265

266

267

268

269

Amplitudes derived using the above methods at 10 hPa are summarized in Table 3. We chose not to use the FFT method to derive QBO amplitudes due to the change in the nature of the QBO towards more than one dominant spectral component in the warming climate. In present day, amplitudes derived from the peak in the Fourier spectrum are in good agreement with the TT and DD method, however due to the multiple peaks in the spectra in Exp 4, FFT based amplitudes largely underestimate the amplitudes derived from the other methods.

270

271

272

273

274

275

276

277

278

For all Exp 2 and nearly all Exp 3 and 4 models the highest estimate of the mean amplitude is obtained from the DD method, which often exceeds the mean amplitude derived from the TT method by more than one standard deviation of the TT derived amplitudes. Partially this is due to Bushell et al. (2019)'s use of binomial smoothing (see above) which, as well as damping high frequency variations, will damp the amplitude of the QBO—binomial smoothing is used here merely for consistency with Bushell et al. (2019). Models and experiments for which the DD method does not give the highest estimate of the mean amplitude are ECHAM5sh (Exp 3 and 4) and LMDz6 (Exp 4). In these simulations there is significant ambiguity when identifying individual QBO cycles (cf. Figure 4) and it is likely that the QBO phase transition point algorithm combines some cycles with the consequence of excluding some of the smaller local easterly and westerly maxima, hence giving rise to a spuriously high estimate for the mean amplitude.

279

280

281

282

283

284

285

The bottom panel in Figure 7 shows changes in the distribution of QBO amplitudes at 10 hPa calculated using the TT method. Although there is a wide spread in QBO amplitudes among the models for Exp 2, most models show a decrease in mean amplitudes for Exp 3 and Exp 4 in comparison to Exp 2, with the exception of MIROC-AGCM-LL with a slightly larger mean amplitude in Exp 3 than in Exp 2. The multi-model mean (see Table 3) also shows a consistent decrease of QBO amplitude in a warming climate using both the TT and DD methods. The range of amplitudes in Exp 3 and Exp 4 increases in most models, as also demonstrated by the increase in the standard deviation of the amplitudes of individual models and multi-model mean.

286

287

288

Figure 8 shows the separation of QBO amplitude at 10 hPa into the westerly and easterly phases based on the TT method. On average there is little change in the amplitude of the westerly phase with a decrease in some models and an increase in others. In contrast the amplitude of the easterly phase decreases in all models and, on average, decreases
<table><caption>**TABLE 3** Summary of estimates of QBO amplitudes in Exp 2, 3, 4 using two different methods: TT and DD. For the TT method standard deviation is denoted after the ± sign. Ratio of westerly to easterly amplitudes calculated from the TT method is listed in the 'W/E' column. Amplitude values for individual models in Exp 3 and Exp 4 derived using the TT method that are statistically different from Exp 2 at the 95 % level as defined by the Student t-test for unequal sized distributions are in bold. Values for ERAI for years 1979 - 2008 are shown for reference in the bottom row.</caption><thead><tr><th rowspan="2">Model/Reanalysis</th><th colspan="3">Exp 2</th><th colspan="3">Exp 3</th><th colspan="3">Exp 4</th></tr><tr><th>TT</th><th>DD</th><th>W/E</th><th>TT</th><th>DD</th><th>W/E</th><th>TT</th><th>DD</th><th>W/E</th></tr></thead><tbody><tr><td>60LCAM5</td><td>23.0±2.9</td><td>22.4</td><td>0.68</td><td>21.4±4.8</td><td>22.5</td><td>0.71</td><td><strong>18.7±3.9</strong></td><td>19.1</td><td>0.69</td></tr><tr><td>AGCM3-CMAM</td><td>22.4±0.9</td><td>24.5</td><td>0.36</td><td><strong>20.7±3.3</strong></td><td>22.2</td><td>0.41</td><td>N/A</td><td>N/A</td><td>N/A</td></tr><tr><td>CESM1(WACCM5-110L)</td><td>26.0±1.3</td><td>29.4</td><td>0.43</td><td><strong>24.8±1.4</strong></td><td>27.8</td><td>0.47</td><td><strong>20.8±2.3</strong></td><td>24.4</td><td>0.61</td></tr><tr><td>ECHAM5sh</td><td>27.8±2.6</td><td>30.6</td><td>0.89</td><td><strong>16.5±10.9</strong></td><td>23.4</td><td>1.07</td><td><strong>21.1±3.3</strong></td><td>14.4</td><td>1.49</td></tr><tr><td>EMAC</td><td>25.7±1.4</td><td>28.0</td><td>0.65</td><td><strong>24.9±1.4</strong></td><td>26.6</td><td>0.69</td><td><strong>22.3±1.7</strong></td><td>23.0</td><td>0.78</td></tr><tr><td>LMDz6</td><td>22.3±1.4</td><td>23.5</td><td>0.65</td><td><strong>19.9±1.6</strong></td><td>19.8</td><td>0.72</td><td><strong>16.2±3.1</strong></td><td>15.4</td><td>1.02</td></tr><tr><td>MIROC-AGCM-LL</td><td>19.7±0.8</td><td>21.2</td><td>0.61</td><td>20.0±0.9</td><td>20.8</td><td>0.61</td><td><strong>17.1±2.6</strong></td><td>18.0</td><td>0.73</td></tr><tr><td>MIROC-ESM</td><td>22.8±1.4</td><td>24.4</td><td>0.47</td><td><strong>21.1±1.9</strong></td><td>21.7</td><td>0.56</td><td><strong>16.5±4.1</strong></td><td>18.0</td><td>0.77</td></tr><tr><td>MRI-ESM2</td><td>16.6±3.7</td><td>18.8</td><td>0.59</td><td><strong>13.3±2.4</strong></td><td>14.1</td><td>0.71</td><td><strong>13.4±2.7</strong></td><td>11.1</td><td>0.96</td></tr><tr><td>UMGA7</td><td>28.7±0.9</td><td>32.1</td><td>0.6</td><td><strong>24.6±2.6</strong></td><td>26.6</td><td>0.70</td><td><strong>19.7±4.1</strong></td><td>21.0</td><td>0.75</td></tr><tr><td>UMGA7gws</td><td>23.8±1.2</td><td>25.8</td><td>0.56</td><td><strong>19.6±3.7</strong></td><td>20.0</td><td>0.66</td><td><strong>13.1±3.8</strong></td><td>10.5</td><td>0.73</td></tr><tr><td><strong>Multi-model Mean</strong></td><td><strong>23.5±1.7</strong></td><td><strong>25.6</strong></td><td><strong>0.6</strong></td><td><strong>20.6±3.2</strong></td><td><strong>22.5</strong></td><td><strong>0.7</strong></td><td><strong>17.9±3.1</strong></td><td><strong>17.6</strong></td><td><strong>0.85</strong></td></tr><tr><td><strong>ERAI</strong></td><td><strong>23.9±3.0</strong></td><td><strong>25.5</strong></td><td><strong>0.39</strong></td><td>N/A</td><td>N/A</td><td></td><td>N/A</td><td>N/A</td><td>N/A</td></tr></tbody></table>

from 30 $ms^{-1}$ in Exp 2 to 20 $ms^{-1}$ in Exp 4. Hence, there is a decrease in the asymmetry between the amplitude of easterly and westerly phases in the multi-model mean response to a warming climate (Figure 8c).

In most models, in Exp 2 the QBO amplitude peaks at the Equator near 10 hPa (Figure 9), then decreases in both amplitude and width on descending into the lower stratosphere (Bushell et al., 2019). In Exp 3 and Exp 4 there is a robust decrease in the amplitude relative to Exp 2 in all models at all levels above 80 hPa. At each level the decrease is more or less uniform across the full width of the QBO but with a hint of a slightly stronger decrease in the Southern Hemisphere (Figure 9). At 10 hPa the amplitude decrease is strongest away from the local maximum in QBO amplitude at the equator, indicating a narrowing of the QBO. The percentage decrease in QBO amplitude in the multi-model mean declines with height between ~60 hPa and 10 hPa: from 36% in Exp 3 and 51% in Exp 4 at 60 hPa to 12% and 31%, respectively at 10 hPa. Since the amplitude response for Exp 4 is not double that for Exp 3, this implies a non-linearity to the response, though additional simulations would be needed to confirm this.
## 4 | CHANGES IN QBO FORCING TERMS

In the present climate there is a balance between momentum tendencies from vertical and meridional advection and resolved and parameterized gravity wave drag that produce an oscillation in the tropics with a mean period of 28 months, and an amplitude of ~24 m s⁻¹ at the peak of the oscillation at 15 hPa. The tendencies driving the QBO can be described using the Transformed Eulerian Mean (TEM) zonal wind equation (Andrews et al., 1987):

$$ \bar{u}_t = -\bar{v}^* [(a \cos\phi)^{-1} (\bar{u} \cos\phi)_\phi - f] - \bar{w}^* \bar{u}_z + \bar{X} + (\rho_0 a \cos\phi)^{-1} \nabla \cdot \mathbf{F} \quad (1) $$

where $\bar{u}$ is the zonal-mean zonal wind, $a$ is the radius of the earth, $\phi$ is the latitude, $\bar{v}^*$ and $\bar{w}^*$ are the residual meridional and vertical velocities. The first term on the right-hand side of (1) is the meridional advection and the Coriolis torque ($f$ is the Coriolis parameter), the second term is the vertical advection, $\bar{X}$ is the gravity wave drag and the last term is the Eliassen-Palm (EP) flux divergence representing the resolved wave drag. Definition of the EP flux vector can be found in Andrews et al. (1987, pp. 127-130) and precise details of how these diagnostics were calculated for the QBOi archive are given by Butchart et al. (2018). Meridional advection and the Coriolis torque contribution to the driving the QBO in the equatorial region (5°S to 5°N, where the horizontal gradient in the zonal mean wind is weak) were found to be small compared to other terms in the budget and is not examined further here.

In a warming climate, each of the tendencies in (1) are likely to respond differently and, as the tendency terms themselves are dependent on the mean wind which they alter, the prediction of how the QBO responds to a warming climate is not straightforward. In the subsequent four subsections we consider how vertical advection, EP flux divergence, and momentum flux from parameterized GWs change in Exp 3 and Exp 4, and then analyze how these changes correlate with the projected changes in QBO period and amplitude.

### 4.1 | Vertical Advection

Vertical advection, proportional to the vertical wind shear and the vertical residual velocity or upwelling, acts to slow down the downward propagation of QBO easterly shear zones (e.g. Saravanan, 1990; Dunkerton, 1991). An increase (decrease) in upwelling would hence act to lengthen (shorten) the QBO cycle, all else being equal. There is observational
evidence that the Brewer-Dobson circulation, and hence tropical upwelling have been increasing since 1980 (Fu et al., 2015). Nearly all modeling studies also agree that a warming climate would lead to increased tropical upwelling (e.g. Butchart et al., 2006; Lin and Waugh, 2013), though there is some spread in the projected magnitude of the increase (e.g. Austin et al., 2003). This speeding-up of the Brewer Dobson circulation comes from increased gravity wave drag resulting from changes in the subtropical jets (Butchart et al., 2010) and increased extra-tropical EP flux divergence. The increased EP flux divergence results from changes in the thermal structure of the atmosphere in a warming climate and also changes in tropospheric sources of upward propagating waves (e.g. Mclandress and Shepherd, 2009; Winter and Bourqui, 2010).

Figure 10 shows the change in the multi-model mean as well as the individual model estimates of residual vertical velocity, $\bar{w}^*$, in Exp 3 and Exp 4. On average, there is an overall increase in $\bar{w}^*$ between 100 and 5 hPa and there is an upward shift in the climatological vertical structure of $\bar{w}^*$, with the minimum in the multi-model mean shifting from 60 hPa in Exp 2 to 50 hPa in Exp 4. The upward shift in $\bar{w}^*$ is consistent with the upward shift of the climatological cold point tropopause (Figure 2b) and upward expansion of the subtrobial jet leading to an upward shift of the shallow branch of the Brewer-Dobson circulation (Garcia and Randel, 2008). In Exp 3 and 4 $\bar{w}^*$ increases in all models at most altitudes in the 100 to 20 hPa range, with a few models not showing much change above 50 hPa (Figures 10b, c). Averaged between 100 and 50 hPa, the near-equatorial residual vertical velocity increases relative to Exp 2 by 10-38% for Exp 3, and by 15-75% for Exp 4. The largest increases occur roughly near the altitude where $\bar{w}^*$ is a minimum between 80 and 50 hPa, and are most pronounced in AGCM3-CMAM, ECHAM5sh, and EMAC in Exp 3 and ECHAM5sh and EMAC in Exp 4.

## 4.2 | Gravity Waves

Parameterized gravity wave drag provides at least 50% of the forcing of the QBO in the QBOi models, except for in MIROC-AGCM-LL, which doesn't employ a GW parameterization (Bushell et al., 2019). In the real atmosphere, GW momentum flux entering the stratosphere can change as a result of changes in tropospheric gravity wave sources, and changes in their propagation due to different background wind and temperature profiles. As the majority of QBOi models uses fixed GW sources (Table 1), the spectrum of waves launched at the source level will be exactly the same in Exp 2, 3, and 4, except for a change in the total momentum flux in models utilizing the Hines (1997a,b) scheme (ECHAM5sh, EMAC, MIROC-ESM, and MRI-ESM2) due to the change in buoyancy frequency at the launch level in Exp 3 and 4 (cf. Figure 2). The launch levels of GWs vary from the surface (MRI-ESM2) to 90 hPa (AGCM3-CMAM), hence,
apart from AGCM3-CMAM, the parameterized GW spectrum reaching the stratosphere is affected by propagation through the troposphere and consequentially is likely to be different in Exp 3 and Exp 4. For MRI-ESM2, in particular, the source spectrum is affected by a 8(18)% increase in Exp 3 (Exp 4) in buoyancy frequency at the surface. ECHAM5sh, EMAC, and MIROC-ESM, also using the Hines (1997a,b) scheme, launch GWs near 600 hPa, at a level at which buoyancy frequency changes between Exp 3 and Exp 4 and Exp 2 are small, hence the GW source spectra are very similar in all experiments. Vertical profiles of eastward and westward GW momentum flux averaged between 10°S and 10°N are shown in Figure 11a for the models with fixed GW sources for Exp 2, 3 and 4. For most models with fixed GW sources the GW momentum flux entering the stratosphere changes by less than 5% for Exp 3, with the exception of AGCM3-CMAM which, shows a 10% decrease of westward momentum flux at 85 hPa, and MRI-ESM2 which shows an 8.8% decrease of eastward momentum flux at 100 hPa. In Exp 4, most momentum flux changes are less than 7.5% relative to Exp 2, except for MRI-ESM2 which shows a 12% increase in eastward momentum flux and 14% increase in westward momentum flux at 100 hPa. MRI-ESM2 has the lowest launching level out of all the models and uses the Hines (1997a,b) schemes, hence the GW spectrum is affected by the increase in *N* in the lowermost troposphere.

The small change in the tropospheric GW momentum fluxes in Exp 3 and Exp 4 for models with fixed GW sources is due, on average, to the absence of any significant change in zonal winds through which the GW spectrum propagates before entering the stratosphere. For the zonal mean, Figure 12a shows that in the tropical troposphere the multi-model mean zonal wind changes by less than 2 m s⁻¹ in Exp 3 and 4 relative to Exp 2 and by much less than the differences between the individual models (indicated by the grey lines in the figure showing ± one inter-model standard deviation). However the zonally averaged zonal wind may not be entirely representative of the local wind profile the GWs propagate through as there are quite significant longitudinal variations in the troposphere (Kawatani et al., 2010). Nonetheless, when averaged, for example, over the western and eastern hemispheres separately (Figures 12b and 12c) there is again little change in tropospheric profiles between Exp 3 and 4 and Exp 2 for the multi-model mean and hence it is reasonable to conclude more generally that in these models there will be little change in GW propagation through the tropical troposphere. Figure 12 shows larger changes in the zonal mean wind, up to 5 m s⁻¹, between Exp 2 and Exp 4 between 120 and 80 hPa. These changes may affect the filtering of GWs in the lowermost part of the QBO, especially for AGCM3-CMAM, which due to its high GW launch level of 90 hPa is very sensitive to the winds in this layer.

Figure 11b shows GW momentum fluxes for models with variable GW sources. The change in 100 hPa GW momentum flux in Exp 3 and Exp 4 relative to Exp 2 is the largest in 60LCAM5 and CESM1(WACCM-110L). In 60LCAM5 and CESM1(WACCM-110L), eastward and westward 100 hPa GW momentum flux averaged between 10°S and 10°N is
20 - 25% larger in Exp 3 than in Exp 2. In Exp 4 these fluxes increase by 48 - 58% (40 - 58%) for eastward (westward) propagating GWs. The changes are very similar when averaged over 5°S to 5°N (not shown). Both 60LCAM5 and CESM1(WACCM-110L) use the Beres (2004) parameterization to calculate GW momentum flux from convectively generated GWs, though with some tuning differences between the two models. In the Beres (2004) parameterization GW momentum flux is proportional to the square of the convective heating rate, which is derived from the convection parameterization (Richter et al., 2010). Although convective heating is not part of the QBOi output protocol, it is very likely that convective heating increases when convective precipitation increases (this was verified in CESM1(WACCM-110L). In Exp 3 (Exp 4) the precipitation rate averaged over 10°S and 10°N in 60LCAM5 and CESM1(WACCM-110L) increases by 3.5% (5%), which is consistent with the increase in GW momentum flux. Averaged between 10°N and 10°S UMGA7gws shows a change of less than 0.5% in tropical eastward momentum flux in Exp 3 and 4 and a decrease of 4% and 9% in the westward 100 hPa GW momentum flux in Exp 3 and 4 respectively. In UMGA7gws, eastward and westward waves are generated equally and launched at 4 km. However, after propagating through tropospheric winds they are no longer symmetrical at 100 hPa. The amplitude of GW momentum flux in UMGA7gws is proportional to the square root of precipitation (Bushell et al., 2015). Averaged over the entire 10°S and 10°N domain, monthly mean precipitation increases by ~ 3% in both Exp 3 and Exp 4 in UMA7gws, however there are regions of both increased and decreased precipitation. The GW momentum flux, dependent on the square root of precipitation, increases in the region of enhanced precipitation over the ITCZ, but decreases in the Indian Ocean and east of Brazil (not shown). The net effect is a slight decrease in GW momentum flux averaged over the 10°S and 10°N domain in Exp 4, which is probably a consequence of the threshold dependence of GW momentum flux on precipitation used by the UMGA7gws source parameterization. LMDz6 also shows a decrease in 100 hPa eastward and westward GW momentum fluxes averaged between 10°S and 10°N of 5% and 10% in Exp 3 and 10% and 18% in Exp 4, respectively. In LMDz6, GW momentum flux is proportional to the square of convective heating (Lott and Guez, 2013), which is estimated from the precipitation (and not taken directly from the convection parameterization). LMDz6 shows a 2% and 5% increase in mean annual precipitation between 10°S and 10°N in Exp 3 and Exp 4, however instantaneous data is also not available, which would have allowed for a full understanding of the reason for the decrease of GW momentum fluxes while precipitation increases. In summary, models that directly link the GW sources (60LCAM5 and CESM1(WACCM-110L) to (parameterized) convective heating show an increase in the GW momentum flux entering the stratosphere, while the two models (LMDZ6 and UMGA7gws) that use precipitation as a proxy for convective heating show a decrease and this uncertainty significantly affects how the QBO responds to warming climate as demonstrated below.
## 4.3 | Resolved waves

407

408

409

410

411

412

413

414

415

416

417

418

419

420

421

422

423

424

425

426

427

428

429

430

431

Resolved waves are likely to change in Exp 3 and Exp 4 due to changes in tropospheric convection and latent heating, as well as changes in the zonal mean wind in the troposphere, although the latter was shown to be rather small (Figure 12). Tropical precipitation (zonal average between 10°S and 10°N) increases by 1 to 7% in the QBOi models between Exp 3 and Exp 2, and by 2 to 13% between Exp 4 and Exp 2 (not shown), suggesting corresponding changes in convection and latent heating and therefore the generation of large scale waves. Figure 13 shows frequency/zonal wavenumber spectra of the vertical components of the EP flux calculated from 6-hourly output averaged between 85 and 70 hPa and between 10°S and 10°N for Exp 2, and the ratio of Exp 3 and Exp 4 to Exp 2. The spectra were computed using the fast Fourier transform using a 72-day window and 12-day overlap, exactly the same methodology as in Holt et al. (2019). In Figure 13 westward propagating (easterly) waves are represented by negative flux and eastward propagating (westerly) waves are represented by positive flux. Some positive flux appears at negative zonal wavenumbers and some negative flux appears at positive zonal wavenumbers as phase speeds are defined relative to the ground and not the background wind. Spectra in Figure 13 are averaged over 10 years of simulation, except for MIROC-AGCM-LL, UMGA7, and UMGA7gws and for which only 4 years of data were available. For the models for which 10 years of data was available, 4-year averages were also made to see if a shorter averaging period changed the spectra, and only very small differences were found. In addition, total easterly and westerly vertical components of EP flux (i.e. the vertical EP flux due to westward and eastward propagating waves, respectively) calculated from 6-hourly data were compared to the provided monthly means and were found to be the same for the sampled periods, and interannual variations were found to be small due to the annually repeating external forcings. Hence it is tacitly assumed that the limited number of years availability of the daily data does not affect any of the conclusions presented in this study. In most models the upward flux due to westward propagating waves with low phase speeds, between 1 and 10 m s⁻¹, increases by a factor of 1.5 to 2 in Exp 3, and by a factor of 2 to 16 in Exp 4. A similar increase in the upward flux due to slow westward propagating waves was found by Kawatani et al. (2011). The slow westward propagating waves contribute directly to driving of the westward phase of the QBO; however, their contribution is smaller than that from the parameterized GWs in most of the QBOi models (Holt et al., 2019; Bushell et al., 2019).

432

433

434

Figure 14 shows that across all the models there is a robust increase in the vertical component of the EP flux for both eastward and westward propagating waves in Exp 3 and Exp 4. For the eastward propagating (westerly) waves the increases range from 6% to 22.5% in Exp 3 and from 20% to 75% in Exp 4, with the multi-model mean increasing by 12%
and 29% in Exp 3 and Exp 4, respectively. For the westward propagating (easterly) propagating waves the increases are even larger: from 18% to 62% in Exp 3 and from 52% to 126% in Exp 4, with the multi-model mean now increasing by 38% in Exp 3 and 91% in Exp 4 (Figure 14). Consequently with the warming climate there is a reduction in the ratio of the westerly to easterly vertical component of EP flux in all models except 60LCAM5 and LMDz6. Figure 8 showed that it is the easterly QBO phase at 10 hPa that primarily decreased in amplitude in Exp 3 and Exp 4, however Figure 14 showed that the easterly component of vertical EP flux increases substantially and more than the westerly component of the vertical EP flux. Hence the connection between changes in QBO amplitude and changes in the vertical component of EP flux is not clear.

## 4.4 | Correlations

In this section we investigate the correlations between the changes in the QBO forcings noted above, and the responses of the QBO amplitude and period to a warming climate. Correlation does not imply causality, however we find these helpful in identifying the possible causes of changes in the QBO characteristics. The top two panels of Figure 15 show that the percentage change in QBO period is positively correlated with the percentage change in the vertical residual velocity, $\bar{w}^*$, in the lower stratosphere averaged between 100 to 50 hPa. The linear correlation coefficient is 0.49 in Exp 3 and in Exp 4 based on the QBO period change derived from the lag correlation method (see Section 3.1). Similar values are found using the TT method and $\bar{w}^*$ at any single level between 100 and 50 hPa. Although not significant these correlations suggest a lengthening of the QBO periods as upwelling increases, as might be expected, but do not fully explain the QBO period response. Instead Figure 15c-f shows that changes in the GW momentum flux entering the stratosphere have a stronger influence on (or at least are better correlated with) the QBO period with correlation coefficients for the eastward (westerly, TauW) and westward (easterly, TauE) GW momentum fluxes, respectively, of -0.67 and -0.84 in Exp 3 and -0.86 and -0.90 in Exp 4. All these correlations are significant at the 95% confidence level apart from the -0.67 in Figure 15c. Hence, larger GW momentum flux entering the stratosphere is most likely the key quantity leading to a shorter QBO period among the QBOi models. Previous studies, for example by Geller et al. (2016a), which varied GW momentum flux (keeping resolved waves and vertical advection constant) showed that increase in GW momentum flux entering the stratosphere leads to a shorter QBO period. In our study, vertical advection and resolved wave forcing also vary, but it appears that the change in GW momentum flux is the dominant driver of QBO period changes. In all the panels in the top three rows of Figure 15 there is an evident clustering of models with large changes
462

463

464

465

466

467

468

469

470

471

472

473

474

475

476

477

478

479

480

481

482

483

484

485

486

487

488

489

490

in GW momentum flux (60LCAM5, CESM1(WACCM-110L), and MRI-ESM2) and those that have small changes in GW momentum flux near 100 hPa. Changes in the easterly and westerly vertical component of the EP flux (bottom two rows of Figure 15) are not as well correlated to changes in QBO period with correlation coefficients of less than 0.2 for westerly component of EP flux, and correlations of -0.52 (-0.30) for the easterly component for Exp 3 (Exp 4) which are found not to be significant. Hence, we conclude that resolved waves do not play a large role in driving the QBO period changes in the QBOi models.

The large changes in GW momentum flux at 100 hPa for some models, especially National Center for Atmospheric Research (NCAR)’s 60LCAM5 and CESM1(WACCM-110L), raise the question of how robust are our conclusions regarding the QBO period being mainly dependent on changes in 100 hPa GW momentum flux. We have repeated the calculation of correlations shown in Figure 15 without the NCAR models and found that some of the significant correlations actually increase. The correlations between percentage QBO period change and: TauW for Exp 4 changes from -0.86 to -0.92, TauE for Exp 3 changes from -0.84 to -0.89, and the correlation between TauE for Exp 4 remains at -0.90. The amplitudes of correlations with other variables change slightly as well, but this does not change whether correlations are significant or not. This, along with seeing no change in the significance of the correlation when one model is randomly removed from the calculation gives us confidence in our conclusion that among the QBOi models GW momentum flux entering the stratosphere is the primary driver of QBO period changes in the warming climate.

Figure 9 showed that the QBO amplitude not only changes with warming climate, but that the relative change varies with altitude and is largest between 60 and 50 hPa. Relative changes in the vertical component of residual velocity also vary with height and are strongest between 80 and 30 hPa (Figure 10b and 10c). Hence, in general, the correlations between the percentage change in the QBO forcings and percentage change in the QBO amplitude are altitude dependent as can be seen in Table 4 for all levels between 85 and 10 hPa. Levels at which the correlations or anti-correlations are strongest are illustrated in Figure 16. In Exp 3 the QBO amplitude change is correlated with $\bar{w}^*$ above 15 hPa, where the correlation coefficient is largest, and is anti-correlated below with the strongest anti-correlation at 70 hPa, where the increase in $\bar{w}^*$ is maximum (cf. Figure 10b). In Exp 4, consistent with the upward shift to the vertical profile of $\bar{w}^*$ (Figure 10a) the strongest anti-correlation with the amplitude response occurs higher, between 40 and 20 hPa, and is now statistically significant. These anti-correlations in the region of large $\bar{w}^*$ strongly suggest that the increase in tropical upwelling is a significant contributor to the amplitude decrease among the QBOi models.

A strong and significant correlation between the change in GW momentum flux entering the stratosphere and QBO amplitude (Figure 16c-16f), suggests that the QBO amplitude response is also influenced by changes to the
parameterized GWs, especially above 50 hPa. Changes in GW momentum flux could be affecting the amount of flux at high gravity wave phase speeds. In particular, while changes in $\bar{w}^*$ are most likely the main cause of QBO amplitude changes in the lower stratosphere, the changes in GW momentum flux entering the stratosphere appear to be a more significant contributor to the changes in QBO amplitudes at higher levels. This seems consistent with the results of Saravanan (1990) who used a generalization of the original Holton and Lindzen (1972) idealized model of the QBO. In Saravanan's study the vertical advection behaved as a damping, which competed with the tendency of the wave forcing to produce an oscillation. Near the lower boundary this damping effect dominated over the wave forcing, but became weaker at higher altitudes, as found here. However, based on the sign of the correlations in Table 4, neither the increase in $\bar{w}^*$ nor the increase in GW momentum flux entering the stratosphere can simply explain the decrease in the amplitude at 10 hPa in the QBOi models (Figure 8a).

Correlations with westerly component of the vertical EP flux, $F_z W$, (Figure 16g, h) suggest that the QBO amplitude decreases as $F_z W$ increases. The increase in the westerly component of EP flux primarily occurs at very low phase speeds (less than $5 \text{ m s}^{-1}$ — see right hand column of Figure 13). It is difficult to interpret exactly how these changes influence QBO amplitude, as the relative roles of resolved and parameterized waves in driving the QBO differ among models (Bushell et al., 2019). Correlations with $F_z W$ and $F_z E$ were also calculated separately for the eastward and westward phases but coefficients were even lower and a robust signal could not be identified. We conclude that resolved waves may play a role in the change of QBO amplitude, however each model's differing momentum budget would need to be considered, which will be reported in follow-up studies. In short, the dependence of QBO amplitude changes on changes in tropical upwelling is a robust feature across the QBOi models especially at lower altitudes, parameterized gravity waves appear to play a role in determining the QBO amplitude above 50 hPa, and the role of resolved waves is inconclusive.

In a manner similar to our tests of the robustness of correlations with QBO period change, we repeated the calculation of the correlations shown in Figure 16 without the NCAR models (60LCAM5 and CESM1(WACCM-110L)) and also by randomly removing models from the sample. In this case, some correlations increase and some decrease, but overall the conclusions remain the same. Removal of NCAR models from the percentage change of QBO amplitude and percentage change in $\bar{w}^*$ correlations raises the correlation coefficient from -0.68 to -0.77. On the other hand, the correlations of QBO amplitude change and GW momentum flux at 100 hPa are reduced (in Exp 4 for TauW from 0.83 to 0.68, and for TauW from 0.87 to 0.68). Correlations without the NCAR models between QBO amplitude and FzE and FzW increase for Exp 4 but not for Exp 3, making the interpretation difficult as was already noted above for resolved
<figure><img src="image_6.png" /><figcaption>**TABLE 4** Linear correlation coefficient at different pressure levels between % change in QBO amplitude and % change in the vertical residual velocity at that level, w*, easterly GW momentum flux at 100 hPa, TauE, and westerly gravity wave momentum flux at 100 hPa, TauW, easterly vertical component of EP flux, EpzE, and westerly vertical component of EP flux, EPzW, in Exp 3 and Exp 4 relative to Exp 2. Correlations significant at the 95 % level are bolded. Colours indicate the level of correlation.</figcaption></figure>

520 wave forcing changes.

5 | DISCUSSION AND CONCLUSIONS

521

522 We have presented here potential changes to the QBO in a warming climate based on eleven general circulation models that participated in QBOi Experiments 2, 3, and 4. QBOi Exp 2 is a present day time-slice experiment with annually repeating sea surface temperatures and present day CO₂. In Exp 3 (Exp 4), SSTs were increased by 2 K (4 K) and CO₂ was doubled (quadrupled). The changes in external forcings were substantial, especially in Exp 4 compared to present day, largely to elucidate whether convergence in the response among the different models would occur under such scenarios. All of the models that participated in the experiments were able to simulate realistic characteristics of the QBO for present day conditions including external forcings, as summarized by Bushell et al. (2019). No changes to

523

524

525

527

528
gravity wave parameterization (or other physical parameterizations, or any other model aspect) were made in the models between Exp 2, 3, and 4. The mean climate response of the QBOi models to a warming climate is consistent with previous similar studies, showing an increase in temperature throughout the troposphere, with largest warming in the tropical upper troposphere, and cooling in the stratosphere. The spread between the individual models in tropical upper tropospheric temperature increases substantially in Exp 3 and Exp 4, as compared to Exp 2, which likely contributes to the large spread in QBO responses.

We found a lack of consistency in the response of the QBO period to a warming climate among the QBOi models. Three of the models found a substantial decrease in the QBO period, by 5 to 9 months in Exp 3, and by 7 to 14 months in Exp 4. Three of the models found a lengthening of the QBO period by 2 to 13 months in Exp 3, and by 2 to 10 months in Exp 4, while other models exhibited little change (less then 2 months). In a multi-model mean, there is little change in the QBO period in a warming climate. In several models in Exp 4 the coherent structure of the QBO characterized by alternating downward propagating easterly and westerly shear zones was lost; however, an oscillation in deseasonalized zonal mean wind was still detected. Vertical residual velocity in the tropical lower stratosphere was found to increase in all QBOi models in Exp 3 and Exp 4, and a positive (although statistically insignificant) correlation between its increase and increase in QBO period was found. We found that the change in QBO period had very large and significant correlations with the momentum flux of gravity waves entering the stratosphere, implying that parameterized gravity waves are the most likely primary driver of QBO period changes. A relationship between an increase in GW momentum flux and shortened QBO period (keeping all other factors the same) in the present day climate has been shown before (e.g.: Geller et al. (2016a,b)). Our results, based on experiments in which both vertical advection and resolved wave forcing are also varying suggest that many of the present models capable of producing QBOs are very dependent on GW parameterizations for generating present day QBOs, which may or may not be realistic. Two of the models that predict substantial QBO period shortening in Exp 3 and Exp 4, 60LCAM5 and CESM1(WACCM5-110L), include a parameterization of gravity wave sources based on Beres et al. (2004), which predicts that in a warming climate eastward (westward) momentum flux at 100 hPa would increase by ~25% (20%) in Exp 3 and ~53% (44%) in Exp 4. MRI-ESM2, the third model predicting shortening of the QBO period showed a ~8% (20%) increase in source level gravity wave flux in Exp 3 (Exp 4) due to GW launching at the surface and increased buoyancy frequency at that level. In the other models in which the period either remained unchanged or lengthened the changes in the GW momentum fluxes are 100 hPa were much smaller. We found no relation between changes in QBO period and change in the vertical component of the EP flux, however the driving from resolved waves is smaller than that from parameterized waves in
the QBOi models, except for in MIROC-AGCM-LL (Bushell et al., 2019; Holt et al., 2019).

Unlike changes in the QBO period, the Exp 3 and Exp 4 simulations showed a robust response across the QBOi models in terms of QBO amplitude changes. All of the models project a decrease of the QBO amplitude in a warming climate. The projected QBO amplitude decrease varies with height and is strongest near 60 hPa, reaching a decrease, on average, of 36% in Exp 3 and 51% in Exp 4 at that level. The QBO amplitude reduction at 10 hPa occurred primarily in the easterly phase in all QBOi models. The QBO amplitude decrease in the lower tropical stratosphere is correlated with an increase in vertical residual velocity in that region. The residual vertical velocity increases by 50% in Exp 3 and by 100% in Exp 4 in a multi-model mean between 60 and 80 hPa. The decrease of QBO amplitude with warming climate, and the relation to changes in tropical upwelling, are consistent with findings of Kawatani and Hamilton (2013) who found observational evidence in QBO amplitude weakening. From about 60 hPa upward, the QBO amplitude change is positively correlated with the change in gravity wave momentum flux entering the stratosphere, such that smaller amplitude reductions occur in models that have increased GW flux. These relations are consistent with the findings of Saravanan (1990), who used a generalization of the original Holton and Lindzen (1972) idealized QBO model to show that the QBO amplitude at its lowest altitudes is damped by upwelling at those altitudes. However, it is also possible that the reasons for the amplitude decrease differ between the models. In this study we also found a negative correlation between the westerly component of the vertical EP flux and change in QBO amplitude in Exp 3; however, due to differing relative roles of resolved waves and parametrized waves in driving the QBO among models these changes seemed less robust and more difficult to interpret.

The largest uncertainty in the response of the QBO in a warming climate comes from the representation of parameterized gravity waves in climate models. This was clear from this study utilizing eleven different GCMs, and is consistent with the study of Schirber et al. (2015) who conducted sensitivity tests using a single GCM with four different gravity wave parameterization configurations. Most GCMs at this time still use GW parameterization with gravity wave sources that are fixed and don't change with changes in tropospheric climate. This assumption is likely to be erroneous, and hence four of the models that participated in this study have employed a source-based gravity wave parameterization. However, as was shown in Figure 11b, different GW source parameterizations predict different changes to GW momentum flux spectra entering the stratosphere, creating more uncertainty rather than convergence with regard to how gravity waves are likely to change in the future climate. The uncertainty in parameterized GWs is rooted in the large uncertainties in modeling unresolved convection. A better understanding of how gravity wave sources will change in the future climate is crucial to improving projections of how the QBO will change in the future.
Convection parameterizations are also a source of uncertainty with regard to resolved waves. MIROC-AGCM-LL, the only model in this study without a gravity wave parameterization, predicts a QBO amplitude decrease in a warming climate similarly to all the other QBOi models in Exp 4 (but predicts no change for Exp 3), and predicts a slight QBO period increase in a warming climate, similarly to some models with fixed and similarly to some models with variable gravity wave sources. That response is also uncertain as the resolved wave spectra, primarily Kelvin and mixed-Rossby GWs in MIROC-AGCM-LL differ from observations in several aspects (Holt et al., 2019).

Our findings suggest that the present day models are tuned for the present climate and reliant on GW parameterizations for the forcing of the QBO, and hence may not be representing the future state of the QBO correctly. As described by Holt et al. (2019), most GCMs underestimate Kelvin and mixed-Rossby gravity waves in present climate and there is large uncertainty with regard to how resolved waves will change in a warming climate, which largely stems from the uncertainties in and deficiencies of convective parameterizations. Our study has shown that the westerly and especially the easterly vertical component of the EP flux are likely to increase in a warmer climate, however these estimates are dependent on changes in tropical precipitation in GCMs which carry biases even in the present day climate. Although the contribution of Kelvin and mixed-Rossby gravity waves to the QBO is uncertain and varies greatly between the QBOi models (Holt et al., 2019), changes in these waves in future climate do not seem to affect the QBO as much as the parameterized gravity waves. There are very weak correlations between changes in resolved waves and the QBO period and inconclusive correlations with the QBO amplitude. However, this might be again a signature of the strong dependence of the QBOs in GCMs on gravity wave parameterizations. Hence, although this study has made significant progress in quantifying the changes to the QBO in a warming climate by utilizing multiple models, and is the first to do so, much uncertainty needs to be resolved before confidence in the response of the QBO changes is reached.

## ACKNOWLEDGEMENTS

We acknowledge the scientific guidance of the WCRP for helping motivate this work, coordinated under the framework of the Stratosphere-troposphere Processes and their Role in Climate (SPARC) QBO initiative (QBOi) led by JA, NB, KH and SO. The QBOi data archive was kindly hosted by the Centre for Environmental Data Analysis (CEDA). This material is based upon work supported by the National Center for Atmospheric Research, which is a major facility sponsored by the National Science Foundation under Cooperative Agreement No. 1852977. The CESM project is supported primarily by the National Science Foundation. Portions of this study were supported by the Regional and
614 Global Model Analysis (RGMA) component of the Earth and Environmental System Modeling Program of the U.S.
615 Department of Energy's Office of Biological & Environmental Research (BER) via National Science Foundation IA
616 1947282. 60LCAM5 and CESM1(WACCM-110L) simulations were carried out on the Yellowstone high-performance
617 computing platform (CISL, 2012). NB was supported by the Met Office Hadley Centre Programme funded by BEIS
618 and Defra, and by the European Commission's Seventh Framework Programme StratoClim project 226520. YK was
619 supported by Japan Society for Promotion of Science (JPSP) KAKENHI Grant Numbers JP15KK0178, JP17K18816 and
620 JP18H01286. YK and KH were supported by the Japan Agency for Marine-Earth Science and Technology (JAMSTEC)
621 through its sponsorship of research at the International Pacific Research Center. SW and YK were partly supported by a
622 Japan Science and Technology Agency (JST) as part of the Belmont Forum GOTHAM project, and by the "Integrated
623 Research Program for Advancing Climate Models (TOUGOU program)" from the Ministry of Education, Culture, Sports,
624 Science and Technology (MEXT), Japan. LH was supported by NASA's Atmospheric Composition Modeling and Analysis
625 Program (ACMAP) Grant #80NSSC18K0069, and NASA's Modeling, Analysis and Prediction (MAP) program Grant
626 #80NSSC17K0169. FL, SW, YK, LG and SO acknowledge funding for the Belmont Forum JPI-Climate GOTHAM
627 project from: Agence Nationale Recherche (ANR-15-JCLI-0004-01), the Japan Science and Technology Agency and
628 the Natural Environment Research Council (NE/P006779/1 and NE/M005828/1). LG and SO are supported by the
629 UK Natural Environment Research Council (NERC) through the National Centre for Atmospheric Science (NCAS). The
630 Earth Simulator was used for MIROC-ESM and MIROC-AGCM-LL simulations. For the EMAC simulations, PB, TK, and
631 SV acknowledge support by the state of Baden-Württemberg through bwHPC. CC and FS have been supported by
632 the Copernicus Climate Change Service, funded by the EU and implemented by ECMWF. HN was supported by an
633 Environment Research and Technology Development Fund (S-12-1(3)) from the Ministry of the Environment of Japan
634 and by Japan Society for Promotion of Science (JPSP) KAKENHI Grant Numbers JP18K03748. KY was supported by
635 KAKENHI Grant Number JP19K03977, JP17H01159 and JP15H05816. We thank Elisa Manzini for helpful comments
636 in the early stages of the manuscript, and Marv Geller and anonymous reviewer for comments that improved the final
637 version.

638 REFERENCES

639 Andrews, D. G., Holton, J. R. and Leovy, C. B. (1987) *Middle Atmosphere Dynamics*. Academic Press.

640 Anstey, J. A. and Shepherd, T. G. (2014) Review Article High-latitude influence of the quasi-biennial oscillation. *Quart. J. Roy.
641 Meteor. Soc., **140**, 1–21.

642 Anstey, J. A., Simpson, I. R., Richter, J. H., Gray, L. J., Naoe, H., Taguchi, M., Bushell, A., Butchart, N., Kawatani, Y., Osprey, S., Richter, J., Serva, F., Braesicke, P., Cagnazzo, C., Chen, C.-C., Chun, H.-Y., Garcia, R., Gray, L. J., Hamilton, K., Kerzenmacher, T., Kim, Y.-H., Lott, F., Mclandress, C., Naoe, H., Scinocca, J., Stockdale, T. N., Watanabe, S., Yoshida, K. and Yukimoto, S. (2019) QBO teleconnections in a multi-model ensemble of QBO-resolving models. *Quart. J. Roy. Meteor. Soc.*, In Preparation.

646 Arakawa, A. and Schubert, W. H. (1974) Interactions of cumulus cloud ensemble with the large-scale environment. *J. Atmos. Sci.*, **31**, 674–701.

648 Austin, J., Shindell, D., Beagley, S. R., Dameris, M., Manzini, E., Nagashima, T., Newman, P. A., Pitari, G., Rozanov, E., Schnadt, C. and Shepherd, T. G. (2003) Uncertainties and assessments of chemistry-climate models of the stratosphere. *Atmospheric Chemistry and Physics*, **3**, 1–27.

651 Baldwin, M. P., Gray, L., Dunkerton, T. J., Hamilton, K., Haynes, P. H., Randel, W. J., Holton, J. R., Alexander, M., Hirota, I., Horinouchi, T., Jones, D. B. A., Kinnersley, J. S., Marquardt, C., Sato, K. and Takahashi, M. (2001) The Quasi-Biennial Oscillation. *Review of Geophysics*, **39**, 179–229.

654 Barnes, E. A. and Polvani, L. (2013) Response of the midlatitude jets, and of their variability, to increased greenhouse gases in the CMIP5 models. *Journal of Climate*, **26**, 7117–7135.

656 Beres, J. H. (2004) Gravity Wave Generation by a Three-Dimensional Thermal Forcing. *J. Atmos. Sci.*, **61**, 1805–1815.

657 Beres, J. H., Alexander, M. J. and Holton, J. R. (2004) A method of specifying the gravity wave spectrum above convection based on latent heating properties and background wind. *J. Atmos. Sci.*, **61**, 324–337.

659 Bushell, A., Anstey, J., Butchart, N., Kawatani, Y., Osprey, S., Richter, J., Serva, F., Braesicke, P., Cagnazzo, C., Chen, C.-C., Chun, H.-Y., Garcia, R., Gray, L. J., Hamilton, K., Kerzenmacher, T., Kim, Y.-H., Lott, F., Mclandress, C., Naoe, H., Scinocca, J., Stockdale, T. N., Watanabe, S., Yoshida, K. and Yukimoto, S. (2019) Evaluation of the Quasi-Biennial Oscillation in global climate models for the SPARC QBO-initiative. *Quarterly Journal of the Royal Meteorological Society*, In revision.

663 Bushell, A. C., Butchart, N., Derbyshire, S. H., Jackson, D. R., Shutts, G. J., Vosper, S. B. and Webster, S. (2015) Parameterized Gravity Wave Momentum Fluxes from Sources Related to Convection and Large-Scale Precipitation Processes in a Global Atmosphere Model. *Journal of the Atmospheric Sciences*, **72**, 4349–4371.

666 Butchart, N., Anstey, J., Hamilton, K., Osprey, S., Mclandress, C., Bushell, A., Kawatani, Y., Kim, Y.-H., Lott, F., Scinocca, J., Stockdale, T., Andrews, M., Bellprat, O., Braesicke, P., Cagnazzo, C., Chen, C.-C., Chun, H.-Y., Dobrynin, M., Garcia, R., Garcia-Serrano, J., Gray, L., Holt, L., Kerzenmacher, T., Naoe, H., Pohlmann, H., Richter, J., Scaife, A., Schenzinger, V., Serva, F., Versick, S., Watanabe, S., Yoshida, K. and Yukimoto, S. (2018) Overview of experiment design and comparison of models participating in phase 1 of the SPARC Quasi-Biennial Oscillation initiative (QBOi). *Geoscientific Model Development*, **11**, 1009–1032.
672 Butchart, N., Cionni, I., Eyring, V., Shepherd, T. G., Waugh, D. W., Akiyoshi, H., Austin, J., Brühl, C., Chipperfield, M. P., Cordero, E., Dameris, M., Deckert, R., Dhomse, S., Frith, S. M., Garcia, R. R., Gettelman, A., Giorgetta, M. A., Kinnison, D. E., Li, F., Mancini, E., Mclandress, C., Pawson, S., Pitari, G., Plummer, D. A., Rozanov, E., Sassi, F., Scinocca, J. F., Shibata, K., Steil, B. and Tian, W. (2010) Chemistry-climate model simulations of twenty-first century stratospheric climate and circulation changes. *Journal of Climate*, **23**, 5349–5374.

673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
Butchart, N., Scaife, A. A., Bourqui, M., Grandpré, J., Hare, S. H., Kettleborough, J., Langematz, U., Manzini, E., Sassi, F., Shibata, K., Shindell, D. and Sigmond, M. (2006) Simulations of anthropogenic change in the strength of the Brewer-Dobson circulation. *Climate Dynamics*, **27**, 727–741.

Camargo, S. J. and Sobel, A. H. (2010) Revisiting the influence of the quasi-biennial oscillation on tropical cyclone activity. *J. Climate.*, **23**, 5810–5825.

CISL (2012) Yellowstone: IBM iDataPlex /FDR-IB. Computational and Information Systems Laboratory. URL: http://n2t.net/ark:/85065/d7wd3xhc.

Collimore, C. C., Martin, D. W., Hitchman, M. H., Huesmann, A. and Waliser, D. E. (2003) On The Relationship between the QBO and Tropical Deep Convection. *Journal of Climate*, **16**, 2552–2568.

Coy, L., Newman, P. A., Pawson, S. and Lait, L. R. (2017) Dynamics of the Disrupted 2015/16 Quasi-Biennial Oscillation. *Journal of Climate*, **30**, 5661–5674. URL: https://doi.org/10.1175/JCLI-D-16-0663.1.

Dee, D. P., Uppala, S. M., Simmons, A. J., Berrisford, P., Poli, P., Kobayashi, S., Andrae, U., Balmaseda, M. A., Balsamo, G., Bauer, P., Bechtold, P., Beljaars, A. C. M., van de Berg, L., Bidlot, J., Bormann, N., Delsol, C., Dragani, R., Fuentes, M., Geer, A. J., Haimberger, L., Healy, S. B., Hersbach, H., V., H. E., Isaksen, L., Kallberg, P. W., Kohler, M., Matricardi, M., McNally, A. P., Monge-Sanz, B. M., Morcrette, J.-J., Park, B. K., Peubey, C., P., d. R., Tavolato, C., Thepaut, J. N. and Vitart, F. (2011) The ERA-Interim reanalysis: Configuration and performance of the data assimilation system. *Quarterly Journal of the Royal Meteorological Society*, **656**, 553–597.

Dunkerton, T. J. (1991) Nonlinear propagation of zonal winds in an atmosphere with Newtonian cooling and equatorial wavdriving. *J. Atmos. Sci.*, **48**, 236–263.

Dunkerton, T. J. and Delisi, D. P. (1985) Climatology of the Equatorial Lower Stratosphere. *Journal of the Atmospheric Sciences*, **42**, 376–396.

Ern, M., Ploeger, F., Preusse, P., Gille, J. C., Gray, L. J., Kalisch, S., Mlynczak, M. G., Russell III, J. M. and Riese, M. (2014) Interaction of gravity waves with the qbo: A satellite perspective. *Journal of Geophysical Research: Atmospheres*, **119**, 2329–2355.

Fu, Q., Solomon, S. and Hartmann, D. L. (2015) Observational evidence of strengthening of the Brewer-Dobson circulation since 1980. *J. Geophys. Res.*, **120**, 10214–10228.
702 Garcia, R. R. and Randel, W. J. (2008) Acceleration of the Brewer-Dobson Circulation due to Increases in Greenhouse Gases. *Journal of the Atmospheric Sciences*, **65**, 2731-2739.

703 Garcia, R. R. and Randel, W. J. (2008) Acceleration of the Brewer-Dobson Circulation due to Increases in Greenhouse Gases. *Journal of the Atmospheric Sciences*, **65**, 2731-2739.

704 Garfinkel, C. I. and Hartmann, D. L. (2010) Influence of the quasi-biennial oscillation on the North Pacific and El Niño teleconnections. *J. Geophys. Res.*, **115**, 1-12.

705 Geller, M. A., Zhou, T., Shindell, D., Ruedy, R., Aleinov, I., Nazarenko, L., Tausnev, N. L., Kelley, M., Cheng, Y., Field, R. D. and Faluvegi, G. (2016a) Modeling the QBO-Improvements resulting from higher- model vertical resolution. *Journal of Advances in Modeling Earth Systems*, **8**, 1092-1105.

708 Geller, M. A., Zhou, T., Shindell, D., Ruedy, R., Aleinov, I., Nazarenko, L., Tausnev, N. L., Kelley, M., Cheng, Y., Field, R. D. and Faluvegi, G. (2016a) Modeling the QBO-Improvements resulting from higher- model vertical resolution. *Journal of Advances in Modeling Earth Systems*, **8**, 1092-1105.

709 Geller, M. A., Zhou, T. and Yuan, W. (2016b) The QBO, gravity waves forced by tropical convection, and ENSO. *Journal of Geophysical Research: Atmospheres*, **121**, 8886-8895. URL: https://agupubs.onlinelibrary.wiley.com/doi/abs/10.1002/2015JD024125.

710

711

712 Giorgetta, M. A. and Doege, M. C. (2005) Sensitivity of the quasi-biennial oscillation to CO2 doubling. *Geophys. Res. Lett.*, **32**, 2-5.

713

714 Giorgetta, M. A., Manzini, E., Roeckner, E., Esch, M. and Bengtsson, L. (2006) Climatology and Forcing of the Quasi-Biennial Oscillation in the MACHAM5 Model. *Journal of Climate*, **19**, 3882-3901.

715

716 Hines, C. O. (1997a) Doppler-spread parametization of gravity-wave momentum deposition in the middle atmosphere. Part 1: Basic formulation. *Journal of Atmospheric and Solar-Terrestrial Physics*, **59**, 371-386.

717

718 — (1997b) Doppler-spread parametization of gravity-wave momentum deposition in the middle atmosphere. Part 2: Broad and quasi monochromatic spectra, and implementation. *Journal of Atmospheric and Solar-Terrestrial Physics*, **59**, 387-400.

719

720 Ho, C.-H., Kim, H.-S., Jeong, J.-H. and Son, S.-W. (2009) Influence of stratospheric quasi-biennial oscillation on tropical cyclone tracks in the western north pacific. *Geophysical Research Letters*, **36**.

721

722 Holt, L., Lott, F., Garcia, R., Kiladis, G. N., Anstey, J., Braesicke, P., Bushell, A., Butchart, N., Cagnazzo, C., Chen, C.-C., Chun, H.-Y., Hamilton, K., Kawatani, Y., Kerzenmacher, T., Kim, Y.-H., McLandress, C., Naoe, H., Osprey, S., Richter, J., Scinocca, J., Serva, F., Versick, S., Watanabe, S., Yoshida, K. and Yukimoto, S. (2019) An evaluation of tropical waves and wave forcing of the QBO in the QBOi models. *Quart. J. Roy. Meteor. Soc.*, Submitted.

724

725

726 Holton, J. R. (1980) Wave propagation and transport in the middle atmosphere. *Phil. Trans. R. Soc. Lond.*, **A 296**, 73-85.

727 Holton, J. R. and Lindzen, R. S. (1972) An updated theory for the quasi-biennial cycle of the tropical stratosphere. *J. Atmos. Sci.*, **29**, 1076-1080.

728

729 Kawatani, Y. and Hamilton, K. (2013) Weakened stratospheric quasibiennial oscillation driven by increased tropical mean upwelling. *Nature*, **497**, 478-481.

730
731 Kawatani, Y., Hamilton, K. and Watanabe, S. (2011) The Quasi-Biennial Oscillation in a Double CO2 Climate. *Journal of the Atmospheric Sciences*, **68**, 265–283.

733 Kawatani, Y., Watanabe, S., Sato, K., Dunkerton, T. J., Miyahara, S. and Takahashi, M. (2010) The Roles of Equatorial Trapped Waves and Internal Inertia-Gravity Waves in Driving the Quasi-Biennial Oscillation. Part II: Three-Dimensional Distribution of Wave Forcing. *Journal of the Atmospheric Sciences*, **67**, 981–997.

736 Lin, J.-L., Kiladis, G. N., Mapes, B. E., M, W. K., R., S. K., W., L., C., W. M., D., S. S., Del Genio, A., Donner, L. J., Emori, S., J, G. and Hourdin, P. J. (2006) Tropical Intraseasonal Variability in 14 IPCC AR4 Climate Models. Part I: Convective Signals. *Journal of Climate*, **19**, 2665–2690.

739 Lin, P. and Waugh, D. W. (2013) Changes in various branches of the Brewer-Dobson circulation from an ensemble of chemistry climate models. *J. Geophys. Res.*, **117**, D20119, doi:10.1029/2012JD017905.

741 Lindzen, R. S. (1981) Turbulence and Stress Owing to Gravity Wave and Tidal Breakdown. *J. Geophys. Res.*, **86**, 9707–9714.

742 Lindzen, R. S. and Holton, J. R. (1968) A theory for the quasi-biennial oscillation. *J. Atmos. Sci.*, **25**, 1095–1107.

743 Lott, F., Denvil, S., Butchart, N., Cagnazzo, C., Giorgetta, M. A., Hardiman, S. C., Manzini, E., Krismer, T., Duvel, J. P., Maury, P., Scinocca, J. F., Watanabe, S. and Yukimoto, S. (2014) Kelvin and rossby-gravity wave packets in the lower stratosphere of some high-top CMIP5 models. *Journal of Geophysical Research*, **119**, 2156–2173.

746 Lott, F. and Guez, L. (2013) A stochastic parameterization of the gravity waves due to convection and its impact on the equatorial stratosphere. *Journal of Geophysical Research Atmospheres*, **118**, 8897–8909.

748 Lott, F., Guez, L. and Maury, P. (2012) A stochastic parameterization of non-orographic gravity waves, Formalism and impact on the equatorial stratosphere. *Geophys. Res. Lett.*, **39**, doi:10.1029/2012GL051001.

750 Mclandress, C. and Shepherd, T. G. (2009) Impact of climate change on stratospheric sudden warmings as simulated by the Canadian Middle Atmosphere model. *Journal of Climate*, **22**, 5449–5463.

752 Naujokat, B. (1986) An update to the observed quasi-biennial oscillation of stratospheric winds over the tropics. *J. Atmos. Sci.*, **43**, 1873–1877.

754 Osprey, S. M., Butchart, N., Knight, J. R., Scaife, A. A., Hamilton, K., Anstey, J. A., Schenzinger, V. and Zhang, C. (2016) An unexpected disruption of the atmospheric quasi-biennial oscillation. *Science*, **353**, 1424–1427.

756 Randel, W. J., Wu, F., Russell, J. M., Roche, A. and Waters, J. W. (1998) Seasonal cycles and QBO variations in stratospheric CH4 and H2O observed in UARS HALOE data. *Journal of the Atmospheric Sciences*, **55**, 163–185.

758 Reid, G. C. and Gage, K. S. (1985) Interannual variations in the height of the tropical tropopause. *Journal of Geophysical Research: Atmospheres*, **90**, 5629–5635. URL: https://agupubs.onlinelibrary.wiley.com/doi/abs/10.1029/JD090iD03p05629.
760 Richter, J. H., Bacmeister, J. T. and Solomon, A. (2014) Effects of vertical resolution and non-orographic gravity wave drag on the simulated climate in the Community Atmosphere Model, Version 5. *J. Adv. Model. Earth Systems*, **6**, 357–383.

762 Richter, J. H., Sassi, F. and Garcia, R. R. (2010) Towards a physically based gravity wave source parameterization in a general circulation model. *J. Atmos. Sci.*, **67**, 136–156.

764 Saravanan, R. (1990) A multiwave model of the quasi-biennial oscillation. *J. Atmos. Sci.*, **47**, 2465–2474.

765 Schirber, S., Manzini, E., Krismer, T. and Giorgetta, M. (2015) The quasi-biennial oscillation in a warmer climate: sensitivity to different gravity wave parameterizations. *Climate Dynamics*, **45**, 825–836.

767 Shepherd, T. G. and McLandress, C. (2011) A Robust Mechanism for Strengthening of the Brewer-Dobson Circulation in Response to Climate Change: Critical-Layer Control of Subtropical Wave Breaking. *Journal of the Atmospheric Sciences*, **68**, 784–797. URL: http://journals.ametsoc.org/doi/abs/10.1175/2010JAS3608.1.

770 Simpson, I. R., Shaw, T. A. and Seager, R. (2014) A Diagnosis of the Seasonally and Longitudinally Varying Midlatitude Circulation Response to Global Warming*. *Journal of the Atmospheric Sciences*, **71**, 2489–2515. URL: http://journals.ametsoc.org/doi/abs/10.1175/JAS-D-13-0325.1.

773 Son, S. W., Lim, Y., Yoo, C., Hendon, H. H. and Kim, J. (2017) Stratospheric control of the [M]adden-[J]ulian oscillation. *Journal of Climate*, **30**, 1909–1922.

775 Warner, C. D. and McIntyre, M. E. (1999) Toward an ultra simple spectral gravity wave parameterization for general circulation models. *Earth Planets Space*, **51**, 475–484.

777 Watanabe, S. and Kawatani, Y. (2012) Sensitivity of the QBO to mean tropical upwelling under a changing climate simulated with an Earth System Model. *Journal of the Meteorological Society of Japan, Series II*, **90A**, 351–360.

779 Winter, B. and Bourqui, M. (2010) Wave forcing in the stratosphere under doubled-CO2 conditions in a 100-year coupled chemistry-climate model study. *Journal of Geophysical Research Atmospheres*, **115**, 1–17.

781 Yoo, C. and Son, S. W. (2016) Modulation of the boreal wintertime Madden-Julian oscillation by the stratospheric quasi-biennial oscillation. *Geophysical Research Letters*, **43**, 1392–1398.

<figure><img src="image_7.png" /></figure>
<figure><img src="image_8.png" /><figcaption>**FIGURE 1** Multi-model mean of the change in zonal mean temperature in K (top panels) and zonal mean zonal wind in m s⁻¹ (bottom panels) between Exp 3 (left panels) and Exp 4 (right panels) relative to Exp 2. Shading indicates temperature (zonal wind) change and solid contours indicate the multi-model mean for Exp2 temperature (zonal mean wind) in top (bottom) panels.</figcaption></figure>
Author Manuscript

<figure><img src="image_9.png" /><figcaption>**FIGURE 2** Zonal mean temperature (left), across-model standard deviation of the zonal mean temperature (middle), and buoyancy frequency (right) averaged between 10°S and 10°S for Exp 2 (solid), Exp 3 (dotted), and Exp 4 (dashed). Grey shading around each line in panels (a) and (c) depicts +/- two standard error (2 times the multi-model standard deviation divided by $\sqrt{n}$, where n is the number of models).</figcaption></figure>
.pdf

<figure><img src="image_10.png" /><figcaption>**FIGURE 3** Zonal mean precipitation rate in units of mm day⁻¹ for Exp 2 (left panel) and for Exp 2 (solid), Exp 3 (dotted), and Exp 4 (dashed) (right panel). Colored lines depict different models.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>FIGURE 4 Zonal mean zonal wind averaged between 5°S and 5°N as a function of pressure and time for Exp 2 (left column), Exp 3(middle column) and Exp 4 (right column). First 20 years of the first ensemble for each model is shown. Model names are noted in the panel titles. Bottom panel in left column shows ERAI reanalysis for years 1979 - 1998.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>**FIGURE 5** Same as Figure 4 but after removing the mean seasonal cycle.</figcaption></figure>

Author Manuscript
<figure><img src="image_13.png" /><figcaption>**FIGURE 6** Correlation of monthly mean zonal mean eastward wind, $\bar{u}$, for Exp 4 at a lag as a function of height, $\bar{u}(z, t \pm \text{lag})$ with $\bar{u}$ at 10 hPa, $\bar{u}(z = 10\text{hPa}, t)$ for the first ensemble member of each model. The black vertical lines are at $\text{lag} = 0$ and at the local minimum correlation (maximum anti-correlation) at 10 hPa either side of the $\text{lag} = 0$ line. The large numbers in black denote the number of months between these two lines, interpreted as the period of the oscillating signal seen in the figure.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>**FIGURE 7** Distribution of QBO periods (top) and amplitudes (bottom) in Exp 2 (grey), Exp 3 (blue), and Exp 4 (red) derived using the TT method at 10 hPa. The distribution median is depicted by a horizontal line in each box, box edges mark the lower quartile and upper quartiles, box whiskers mark the minimum and maximum values. Black dots represent mean values. Models are ordered according to Exp 2 mean period (top panel) amplitude (bottom panel). Horizontal grey lines indicate present day ERAI reanalysis values.</figcaption></figure>

Author Manuscript
<figure><img src="image_15.png" /><figcaption>**FIGURE 8** QBO amplitude changes at 10 hPa separated by phase derived using the TT method for Exp 2 (gray), Exp 3 (blue) and Exp 4 (red).</figcaption></figure>
Author Manuscript

<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /></figure>

<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /></figure>

<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /></figure>

<figure><img src="image_24.png" /><figcaption>**FIGURE 9** Multi-model mean of QBO amplitude calculated using the DD method in m s⁻¹ as a function of pressure (top row) and as a function of latitude at 10 hPa (middle row) and 50 hPa (bottom row). In the left column, thin colored lines depict individual model QBO amplitudes for Exp 2 and multi-model means are shown for Exp 2 (thick black solid line), Exp 3 (dotted line), and Exp 4 (dashed line). Middle and right columns depict ratios of QBO amplitudes for Exp 3 to Exp 2 and Exp 4 to Exp 2 respectively, with colored lines depicting individual models and the dotted (dashed) line depicting the multi-model mean of ratios for Exp 3 to Exp 2 (Exp 4 to Exp 2).</figcaption></figure>
<figure><img src="image_25.png" /><figcaption>**FIGURE 10** Multi-model mean of zonally averaged residual vertical velocity between 5°S to 5°N for Exp 2, 3 and 4 (left panel), and its change from Exp 3 to Exp 2 (center panel) and change from Exp 4 to Exp 2 (right panel). In the left panel solid line represents Exp 2, dotted Exp 3, and dashed Exp 4. Grey (peach) shading around the solid (dashed) line depicts +/- two standard error. In the two rightmost panels, thin coloured lines represent individual models, thick black line is the multi-model mean.</figcaption></figure>
<figure><img src="image_26.png" /><figcaption>**FIGURE 11** Vertical profiles of GW momentum flux averaged between 10°S and 10°N for models with fixed GW source parameterizations (left) and interactive GW sources (right). Solid lines represent Exp 2, dotted lines represent Exp 3 and dashed lines are for Exp 4. Eastward (westward) momentum flux is depicted by positive (negative) values respectively. Note that momentum flux was divided by a factor of 10 for 60LCAM5. GW momentum fluxes were not available for EMAC, and MIROC-AGCM-LL does not use a GW parameterization, hence only nine models are shown.</figcaption></figure>

<figure><img src="image_27.png" /><figcaption>**FIGURE 12** Multi-model mean of zonal mean wind averaged between 5°S and 5°N and all longitudes (left), 0 - 180° (middle), and 180 - 360° (right) (Black lines) and multi-model mean +/- one inter-model standard deviation (grey lines). Solid lines are for Exp 2, dotted for Exp 3, and dashed for Exp 4.</figcaption></figure>

<figure><img src="image_28.png" /></figure>
<figure><img src="image_29.png" /><figcaption>**FIGURE 13** Vertical component of EP flux averaged between 85 and 70 hPa and 10°S and 10°N as a function of zonal wavenumber and frequency for Exp 2 (left column). Middle (right) column shows ratio of vertical EP flux for Exp 3 to Exp 2 (Exp 4 to Exp 2). Positive (negative) wavenumbers represent eastward (westward) propagating waves.</figcaption></figure>
<figure><img src="image_30.png" /><figcaption>**FIGURE 14** Vertical component of EP flux averaged between 85 and 70 hPa and over (a) westerly (negative) wave phase speeds, (b) easterly (positive) wave phase speeds. Ratio of westerly to easterly vertical component of EP flux is shown in panel (c). Values are shown for Exp 2 (gray), Exp 3 (blue) and Exp 4 (red) for individual models as well as the multi-model mean.</figcaption></figure>
<figure><img src="image_31.png" /><figcaption>**FIGURE 15** Scatter plot of relationship between changes in QBO period from Exp 3 to Exp 2 (left panels) and Exp 4 to Exp 2 (right panels) vs % change in $\bar{w}^*$ averaged between 50 and 100 hPa (1st row); westerly (eastward) momentum flux at 100 hPa, TauW (2nd row); easterly (westward) momentum flux at 100 hPa, TauE (3rd row); westerly vertical component of EP flux averaged between 85 and 70 hPa, FzW (4th row); and easterly vertical component of EP Flux averaged between 85 and 70 hPa, FzE (5th row). Circles denote QBO period changes using the lag correlation method, and the grey line shows the correlation line using those values. Triangles depict QBO period changes calculated using the TT method and are only shown when the triangles don't overlap with the circles. Linear correlation coefficient, R, is depicted in the top right corner of the plot, with the first value showing correlations with QBO period changes derived using the lag correlation method and correlations with QBO period changes derived using the TT method are shown in parenthesis. Red color indicates that the correlation is significant according to a two-sided student t-test (p-value < 0.05). Model labels are in the top left panel. In the second and third rows EMAC and MIROC-AGCM-LL are missing as TauW/TauE were not available for EMAC and MIROC-AGCM-LL does not use a GW parameterization.</figcaption></figure>
Author Manuscript

<figure><img src="image_32.png" /><figcaption>**FIGURE 16** Scatter plot of relationship between changes in QBO amplitude derived using the DD method from Exp 3 to Exp 2 (left panels) and Exp 4 to Exp 2 (right panels) vs: % change in $w^*$ at the same level as QBO amplitude depicted in the panel title (top row); westerly (eastward) momentum flux at the same level as QBO amplitude, TauW, (2nd row); easterly (westward) momentum flux at the same level as QBO amplitude, TauW, (3rd row); westerly vertical component of EP flux averaged between 85 and 70 hPa (4th row); and easterly vertical component of EP flux averaged between 85 and 70 hPa (5th row). Correlations with amplitude are shown for the level at which maximum correlation exists. Model labels are shown in the top left panel. Linear correlation coefficient, R, is depicted in the top right corner of the plot. Red color indicates that the correlation is significant according to a two-sided student t-test (p-value < 0.05). In the 2nd and 3rd rows EMAC and MIROC-AGCM-LL are missing as TauW/TauE were not available for EMAC and MIROC-AGCM-LL does not use a GW parameterization.</figcaption></figure>
<caption>a)</caption>

<figure><img src="image_33.png" /><figcaption>b)</figcaption></figure>

<figure><img src="image_34.png" /><figcaption>Exp 4 - Exp 2</figcaption></figure>

Temperature Difference (K)

<caption>c)</caption>

<figure><img src="image_35.png" /><figcaption>d)</figcaption></figure>

<figure><img src="image_36.png" /><figcaption>Exp 4 - Exp 2</figcaption></figure>

Author Manuscript
<figure><img src="image_37.png" /></figure>
<figure><img src="image_38.png" /></figure>

<figure><img src="image_39.png" /></figure>
<figure><img src="image_40.png" /></figure>
<figure><img src="image_41.png" /></figure>

<figure><img src="image_42.png" /></figure>
(a) Westerly

<figure><img src="image_43.png" /></figure>

(b) Easterly

<figure><img src="image_44.png" /></figure>

(c) Westerly/Easterly

<figure><img src="image_45.png" /></figure>

60LCAM5
AGCM3-CMAM

CESM1(WACCM-110L)
ECHAM5sh

EMAC

LMDz6

MIROC-AGCM-LL
MIROC-ESM

MRI-ESM2

UMGA7

UMGA7gws

Multi-model mean
<figure><img src="image_46.png" /></figure>

<figure><img src="image_47.png" /></figure>

<figure><img src="image_48.png" /></figure>

<figure><img src="image_49.png" /></figure>

<figure><img src="image_50.png" /></figure>

<figure><img src="image_51.png" /></figure>

<figure><img src="image_52.png" /></figure>

<figure><img src="image_53.png" /></figure>

<figure><img src="image_54.png" /></figure>
<figure><img src="image_55.png" /></figure>
<figure><img src="image_56.png" /><figcaption>b) Variable GW Sources</figcaption></figure>

<figure><img src="image_57.png" /></figure>
<figure><img src="image_58.png" /><figcaption>(a) All Lons</figcaption></figure>

<caption>(b) 0 - 180°</caption>

<caption>(c) 180 - 360°</caption>
<figure><img src="image_59.png" /></figure>
<figure><img src="image_60.png" /></figure>
<figure><img src="image_61.png" /></figure>

<figure><img src="image_62.png" /></figure>

<figure><img src="image_63.png" /></figure>

<figure><img src="image_64.png" /></figure>

<figure><img src="image_65.png" /></figure>

<figure><img src="image_66.png" /></figure>

<figure><img src="image_67.png" /></figure>

<figure><img src="image_68.png" /></figure>

<figure><img src="image_69.png" /></figure>

<figure><img src="image_70.png" /></figure>
<figure><img src="image_71.png" /></figure>

<figure><img src="image_72.png" /></figure>

<figure><img src="image_73.png" /></figure>

<figure><img src="image_74.png" /></figure>

<figure><img src="image_75.png" /></figure>

<figure><img src="image_76.png" /></figure>

<figure><img src="image_77.png" /></figure>

<figure><img src="image_78.png" /></figure>

<figure><img src="image_79.png" /></figure>

<figure><img src="image_80.png" /></figure>
<figure><img src="image_81.png" /><figcaption>QJ 3749_Figure 4.png</figcaption></figure>

This article is protected by copyright. All rights reserved.
<figure><img src="image_82.png" /><figcaption>QJ 3749 Figure 5.png</figcaption></figure>

This article is protected by copyright. All rights reserved.
Response of the quasi-biennial oscillation to a
warming climate in global climate models

Jadwiga H. Richter*, Neal Butchart, Yoshio Kawatani, Andrew C. Bushell, Laura Holt, Federico Serva, James Anstey, Isla R. Simpson, Scott Osprey, Kevin Hamilton, Peter Braesicke, Chiara Cagnazzo, Chih-Chieh Chen, Rolando R. Garcia, Lesley J. Gray, Tobias Kerzenmacher, Francois Lott, Charles McLandress, Hiroaki Naoe, John Scinocca, Timothy N. Stockdale, Stefan Versick, Shingo Watanabe, Kohei Yoshida, Seiji Yukimoto

<figure><img src="image_83.png" /></figure>

The response of the quasi-biennial oscillation (QBOi) to a warming climate was examined in eleven general circulation models. No consistency was found among the models for the QBO period response. In contrast, all the models projected a decrease in the QBO amplitude in a warmer climate.
Exp 3

Level

<table><thead><tr><th></th><th>0.16</th><th>0.64</th><th>0.63</th><th>-0.48</th><th>0.01</th></tr></thead><tbody><tr><th>10 hPa</th><td></td><td></td><td></td><td></td><td></td></tr><tr><th>15 hPa</th><td><strong>0.62</strong></td><td><strong>0.74</strong></td><td>0.63</td><td>-0.34</td><td>0.12</td></tr><tr><th>20 hPa</th><td>-0.25</td><td><strong>0.71</strong></td><td>0.62</td><td>-0.31</td><td>0.14</td></tr><tr><th>30 hPa</th><td>-0.31</td><td><strong>0.78</strong></td><td><strong>0.77</strong></td><td>-0.33</td><td>0.27</td></tr><tr><th>40 hPa</th><td>-0.06</td><td><strong>0.69</strong></td><td><strong>0.82</strong></td><td>-0.45</td><td>0.36</td></tr><tr><th>50 hPa</th><td>-0.10</td><td>0.67</td><td><strong>0.82</strong></td><td>-0.46</td><td>0.36</td></tr><tr><th>60 hPa</th><td>-0.34</td><td>0.62</td><td><strong>0.71</strong></td><td>-0.41</td><td>0.17</td></tr><tr><th>70 hPa</th><td>-0.52</td><td>0.25</td><td>0.23</td><td>-0.45</td><td>-0.25</td></tr><tr><th>85 hPa</th><td>-0.46</td><td>-0.08</td><td>0.11</td><td>-0.25</td><td>0.25</td></tr></tbody></table>

Exp 4

<table><thead><tr><th></th><th>0.40</th><th>0.70</th><th>0.64</th><th>-0.03</th><th>-0.11</th></tr></thead><tbody><tr><th>10 hPa</th><td></td><td></td><td></td><td></td><td></td></tr><tr><th>15 hPa</th><td>-0.07</td><td>0.69</td><td>0.60</td><td>0.18</td><td>-0.19</td></tr><tr><th>20 hPa</th><td></td><td></td><td></td><td></td><td></td></tr><tr><th>30 hPa</th><td><strong>-0.68</strong></td><td><strong>0.80</strong></td><td><strong>0.76</strong></td><td>0.18</td><td>-0.04</td></tr><tr><th>40 hPa</th><td><strong>-0.66</strong></td><td><strong>0.83</strong></td><td><strong>0.87</strong></td><td>0.05</td><td>0.12</td></tr><tr><th>50 hPa</th><td><strong>-0.53</strong></td><td><strong>0.80</strong></td><td><strong>0.86</strong></td><td>-0.05</td><td>0.09</td></tr><tr><th>60 hPa</th><td><strong>-0.28</strong></td><td>0.70</td><td><strong>0.72</strong></td><td>0.17</td><td>-0.11</td></tr><tr><th>70 hPa</th><td><strong>-0.35</strong></td><td>0.43</td><td>0.38</td><td>0.41</td><td>-0.53</td></tr><tr><th>85 hPa</th><td><strong>-0.11</strong></td><td>-0.23</td><td>-0.31</td><td>0.02</td><td>-0.54</td></tr><tr><th>90 hPa</th><td><strong>0.12</strong></td><td>-0.63</td><td>-0.42</td><td><strong>-0.75</strong></td><td>0.20</td></tr></tbody></table>

<figure><img src="image_84.png" /></figure>

<figure><img src="image_85.png" /></figure>
# The SPARC Quasi-Biennial Oscillation initiative

James A. Anstey<fnref n="1" /> | Neal Butchart<fnref n="2" /> | Kevin Hamilton<fnref n="3" />
| Scott M. Osprey<fnref n="4" />

Correspondence

Scott Osprey, Atmospheric, Oceanic and Planetary Physics, Clarendon Laboratory, Parks Road, Oxford, OX1 3PU

Email: scott.osprey@physics.ox.uk

KEYWORDS

WCRP-SPARC, stratosphere, *quasi-biennial oscillation*, gravity
waves, climate modelling

1

2

3

Following its establishment in 1992 the Stratosphere-troposphere Processes And their Role in Climate (SPARC) core project of the World Climate Research Programme (WCRP) has included several initiatives designed to evaluate and potentially improve the representation of the stratosphere in general circulation models (GCMs). With the benefit of contributions from these initiatives, rapid advances have been made over the last 25 years in our capability to model the influences of, and changes in, stratospheric processes. Notably, developments in chemistry-climate models (CCMs) first supported by the SPARC CCM-evaluation activity (CCMval) and more recently by the SPARC CCM-Initiative (CCMI)

<footnote n="1">¹Canadian Centre for Climate Modelling and Analysis (CCCma), Victoria, Canada</footnote>

<footnote n="2">²Met Office Hadley Centre (MOHC), Exeter, UK</footnote>

<footnote n="3">³International Pacific Research Center (IPRC), University of Hawaii, Honolulu, USA</footnote>

<footnote n="4">⁴Atmospheric, Oceanic and Planetary Physics (AOPP) & National Centre for Atmospheric Science (NCAS), University of Oxford, Oxford, UK</footnote>
10

11

12

13

14

15

16

have enabled models to play a key role in strengthening the scientific foundations of World Meteorological Organisation (WMO) / United Nations Environment Programme (UNEP) ozone assessments. Nonetheless, shortcomings in the representation of the stratosphere in models remain and likely affect the representation of interannual variability as well as the two-way coupling between the stratosphere and troposphere. The Quasi-Biennial Oscillation initiative (QBOi) was conceived as a grassroots initiative in 2012 to focus on the fidelity of tropical variability in model stratospheres and the possible influence of this variability on other parts of the atmosphere. QBOi was adopted as an approved SPARC activity in 2015.

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

An early SPARC initiative to assess the dynamics underpinning the performance of middle atmosphere climate models was the GCM-Reality Intercomparison Project for SPARC (GRIPS), established in 1996. By the early 2000s GRIPS analysis had shown models were typically deficient in wave drag, leading to climatological biases in mean temperatures and winds throughout the extratropical middle atmosphere. Subsequently, many modelling groups alleviated these circulation deficiencies by implementing non-orographic gravity wave parametrizations to supply the missing wave drag and drive the corresponding adiabatic heating. GRIPS also identified clear issues in the simulation of the tropical middle atmosphere, notably the lack of a QBO in virtually all models at the time. Increasingly, however, QBO-like oscillations of descending eastward and westward winds appeared in model tropical stratospheres as non-orographic gravity wave parametrizations were introduced and provided that suitable parameter settings were used. In fact, due to a lack of strong observational constraints, these wave drag parametrizations can be tuned to simulate oscillations having similar amplitude and period to the observed QBO. Alongside the development of non-orographic gravity wave parametrizations in the 1990s resolved wave forcing in the tropical lower stratosphere was found to improve with increased vertical resolution and the first realistic simulation of the QBO was actually by a GCM without parameterised non-orographic wave drag. Nowadays it is accepted that the QBO is forced by a spectrum of large- and small-scale tropical waves, and that both high vertical resolution and parameterised non-orographic gravity wave drag are important ingredients for simulating a QBO.

33

34

35

36

37

38

Over the last twenty years this rather simple recipe of improving vertical resolution and including parameterised non-orographic gravity wave drag has led to a significant increase in the number of stratosphere-resolving GCMs simulating realistic QBOs. In addition, a few models (including one participating in QBOi) have sufficient resolved wave forcing to simulate a QBO without the use of parameterised drag. While this progress has been very encouraging it also raises important questions about whether these models produce realistic QBOs for realistic reasons, and how sensitive the simulated QBOs are to changes in model characteristics or external forcings (such as increased greenhouse gas
concentrations). Hence there is a need for in-depth analysis and inter-comparison of simulated QBOs in contemporary models using a more extensive range of QBO diagnostics, and this provided the initial motivation for QBOi. Longer term objectives of QBOi are to reduce biases in the simulated circulation of the tropical stratosphere with a focus on those aspects of the QBO identified as being important for teleconnections with extratropical and tropospheric climate.

From its inception QBOi has been community-driven and at the start-up workshop in Victoria, Canada in 2015 the participants, including many of the world's leading QBO experts, prioritised five key areas to focus on during phase one of the initiative:

* Quantitative assessment of simulated QBOs using an extensive suite of metrics

* Response of simulated QBOs to climate change

* Predictability of the QBO on seasonal to annual timescales

* Simulation of vertically propagating equatorial waves contributing to driving the QBOs

* QBO teleconnections

Collaborative research in these priority areas was facilitated by the community agreeing on the design of a suite of coordinated GCM experiments. Details of this suite of experiments are presented in a separate paper published in *Geoscientific Model Development* (Butchart et al., 2018). Seventeen models completed some or all of these experiments, and output was uploaded to the QBOi data archive hosted by the UK's Centre for Environmental Data Analysis (CEDA). Since many of these models are closely related to the coupled-GCMs participating in phase 6 of the Coupled Model Intercomparison Project (CMIP6), analysis of the QBOi database provides an opportunity for interpreting and understanding the behaviour of the CMIP6 models. However, unlike the CMIP6 endorsed MIPs, QBOi is a more focused, flexible and independent community-driven activity, providing an alternative approach for researchers with perhaps limited resources to contribute to coordinated experiments helping to address policy-relevant scientific questions and issues.

The Special Section that follows is guest edited by QBO guru Mark Baldwin. It contains five core papers that provide an overview of the main research findings using the QBOi multi-model ensemble in each of the priority areas identified above, as well as an additional paper analysing the stratospheric semi-annual oscillation (SAO) in the QBOi models. In-depth follow-up analysis, and other studies using the QBOi data archive, will appear in a corresponding Special Collection of the *Quarterly Journal of the Royal Meteorological Society*. The first core paper (Bushell et al.) evaluates
key metrics of the QBO under present-day forcing conditions, providing an overview of the multi-model ensemble performance. The relative roles of resolved and parameterised waves in driving the QBO are also evaluated, but a detailed breakdown of the momentum budget is reserved for a follow-up study. The second paper (Richter *et al.*.) examines the QBO response to idealised climate change scenarios and the correlation to changes in the QBO drivers. Again a detailed analysis of the momentum budget in a changing climate is postponed to a later publication. The QBO's predictability is examined in the third paper of the collection (Stockdale *et al.*). Analysis of equatorial Kelvin and Rossby-gravity waves in the models and their associated forcing of the mean flow is the main focus of the fourth paper (Holt *et al.*). Conclusions relating to the QBOs impact on other parts of the atmosphere through teleconnections are elaborated in the fifth core paper (Anstey *et al.*). Finally, Smith *et al.* examine how the stratospheric SAO, an equatorial oscillation located above the QBO in the region of the stratopause, is represented by the QBOi models.

The QBOi Phase One Coordinators, James Anstey, Neal Butchart, Kevin Hamilton and Scott Osprey, would like to express their appreciation to the six lead authors for agreeing to organise and lead the write-up of the analysis of the Phase One simulations. We gratefully acknowledge the Natural Environment Research Council (NERC), UK, for the original QBOnet grant which funded the first two workshops (March 2015 in Victoria and September 2016 in Oxford). The generosity of the Canadian Centre for Climate Modelling and Analysis at Environment and Climate Change Canada, the sub-department of Atmospheric, Oceanic and Planetary Physics in the Department of Physics, University of Oxford, and the University of Kyoto, for hosting QBOi workshops in Victoria (2015), Oxford (2016) and, Kyoto (2017) – the Kyoto workshop being held jointly with the SPARC Stratospheric And Tropospheric Influences On Tropical Convective Systems (SATIO-TCS) and Fine Scale Atmospheric Processes and Structures (FISAPS) initiatives – was very much appreciated. Special thanks go to WCRP-SPARC for providing support for early career scientists to attend the workshops, and finally to the UK Centre for Environmental Data Analysis (CEDA) for hosting the QBOi Phase One data archive and providing a shared online collaborative work space.

## REFERENCES

Butchart, N., Anstey, J. A., Hamilton, K., Osprey, S., McLandress, C., Bushell, A. C., Kawatani, Y., Kim, Y.-H., Lott, F., Scinocca, J., Stockdale, T. N., Andrews, M., Bellprat, O., Braesicke, P., Cagnazzo, C., Chen, C.-C., Chun, H.-Y., Dobrynin, M., Garcia, R. R., Garcia-Serrano, J., Gray, L. J., Holt, L., Kerzenmacher, T., Naoe, H., Pohlmann, H., Richter, J. H., Scaife, A. A., Schenzinger, V., Serva, F., Versick, S., Watanabe, S., Yoshida, K. and Yukimoto, S. (2018) Overview of experiment design and comparison of models participating in phase 1 of the SPARC Quasi-Biennial Oscillation initiative (QBOi). Geoscientific Model Development
94

opment, 11, 1009–1032. URL: https://www.geosci-model-dev-discuss.net/gmd-2017-187/https://www.geosci-model-

95

dev.net/11/1009/2018/.

<figure><img src="image_86.png" /></figure>