Testing chemistry-climate models ’ regulation of tropical lower-stratospheric water vapor

Climate models predict that tropical lower stratospheric humidity will increase as the climate warms, with important implications for the chemistry and climate of the atmosphere. We analyze the trend in 21-century simulations from 12 stateof-the-art chemistry-climate models (CCMs) using a linear regression model to determine the factors driving the trends. Within CCMs, the long-term trend in humidity is primarily driven by warming of the troposphere. This is partially offset in most CCMs by an increase in the strength of the Brewer-Dobson circulation, which tends to cool the tropical tropopause layer (TTL). We 5 also apply the regression model to individual decades from the 21 century CCM runs and compare them to observations. Many of the CCMs, but not all, compare well with observations, lending credibility to their predictions. One notable deficiency in most CCMs is that they underestimate the impact of the quasi-biennial oscillation on lower stratospheric humidity. Our analysis provides a new and potentially superior way to evaluate model trends in lower stratospheric humidity.


Introduction
Variations of stratospheric water vapor can impact both the climate and chemistry of the atmosphere.Because of this, understanding the processes that control the humidity of air entering the tropical lower stratosphere (hereafter [H 2 O] entry ) has been a high priority of the scientific community since Brewer (1949) first described the stratospheric circulation.
It is now well established that the fundamental control over [H 2 O] entry comes from the cold temperatures found in the tropical tropopause layer (TTL) (Fueglistaler et al., 2009b) and that variability in these temperatures translates into variability in [H 2 O] entry .The most well-known example of this is the so-called "tape recorder," in which the seasonal cycle in TTL temperatures is imprinted on tropical stratospheric water vapor (Mote et al., 1996).
On interannual time scales, variability in [H 2 O] entry originates from processes such as the Brewer-Dobson Circlation (BDC) and the quasi-biennial oscillation (QBO).More recently, Dessler et al. (2013Dessler et al. ( , 2014) ) has suggested that the temperature of the 1 Atmos.Chem.Phys. Discuss., doi:10.5194/acp-2016-964, 2016 Manuscript under review for journal Atmos.Chem.Phys.Published: 8 November 2016 c Author(s) 2016.CC-BY 3.0 License.troposphere also exerts an influence on [H 2 O] entry .This implies the existence of a stratospheric water vapor feedback, whereby a warming climate would increase stratospheric water vapor, leading to further warming.
Putting these factors together, Dessler et al. (2013Dessler et al. ( , 2014) ) demonstrated that observed [H 2 O] entry could be accurately reproduced with a simple linear model: Where ∆T is the temperature of the troposphere, BDC is the strength of the Brewer-Dobson circulation, QBO represents the phase of the QBO, and epsilon is the residual.As expected, they found that a stronger BDC, which tends to cool the TTL, reduces [H 2 O] entry ; this is consistent with previous analyses (Brewer, 1949;Randel et al., 2006;Castanheira et al., 2012;Fueglistaler et al., 2014;Gilford et al., 2016).They also found that the QBO introduces significant variability with a time scale of a few years, also consistent with previous work (O'Sullivan and Dunkerton, 1997;Randel et al., 1998;Dunkerton, 2001;Fueglistaler and Haynes, 2005;Choiu et al., 2006;Liang et al., 2011;Castanheira et al., 2012;Khosrawi et al., 2013;Kawatani et al., 2014;Tao et al., 2015).Finally, a warmer troposphere tends to increase [H 2 O] entry , although whether this is through influence on TTL temperatures or some other mechanism, such as convective ice lofting, is not clear.
Virtually all climate models predict that [H 2 O] entry will increase as the climate warms.Dessler et al. (2013)  observations.The purpose of this paper is to use this technique to examine a set of CCMs, with the goal of providing insight into the realism of the models.
We use simulations from the REF-B2 scenario of CCMVal-2.In this scenario, greenhouse gas concentrations during the 21 st century come from the A1B scenario, which lies in the middle of the SRES scenarios (IPCC, 2001).Ozone-depleting substance come from the halogen emission scenario A1 described by (WMO, 2007).CCMVal-2 specifics can be found inSPARC (2010) and Morgenstern et al. (2010).We use the refC2 scenario of the CCMI-1.In this scenario, greenhouse gas concentrations come from the RCP 6.0 scenario (Meinshausen et al., 2011) and ozone-depleting substances come from the halogen emission scenario A1 described described by (WMO, 2014).CCMI-1 model specifics can be found in Morgenstern et al. (2016).In order to maintain a consistent reference period between models, our analysis covers 2000-2097, which we will hereafter refer to as "the 21 st century ".
For each model, we fit CCM [H 2 O] entry using the multivariate linear regression (MLR) model described above.We use tropical average 80-hPa water vapor volume mixing ratio anomaly as a proxy for [H 2 O] entry (all tropical averages in this paper are averages over 30°N-30°S; anomalies are calculated by subtracting off the mean annual cycle from the time series).
For our BDC index, we use 80-hPa diabatic heating rate anomalies (see Fueglistaler et al. (2009a) for details).The tropospheric temperature index is the 500-hPa tropical average temperature anomaly, and for the few CCMI-1 simulations that only produce variables on hybrid pressure levels (CMAM, CCSRNIES-MIROC3.2, and MRI-ESM1r1), we choose a hybrid pressure level close to the 500-hPa pressure surface (See Table 1).All of these choices are similar to those used by Dessler et al. (2013Dessler et al. ( , 2014)).
For the QBO index, we take the standardized anomaly of equatorial 50-hPa zonal winds.By examining 21 st century 50 hPa zonal winds (shown in supplement figures), we find that only 5 of the 12 models simulate a QBO (table 1).As a result, we do not expect the QBO to significantly impact [H 2 O] entry in all of the models.
The MLR returns the coefficients for each regression coefficient in Equation 1, along with an uncertainty for each coefficient.
Unless otherwise noted, we use 95%-confidence intervals in this paper.Autocorrelation of the time series is accounted for in the uncertainties following Santer et al. (2000).

Century Analysis
We first analyze the long-term trend in trends in the time series, we will refer to this as the "trended analysis".

Detrended 21 st Century
One concern with the trended analysis is that the [H 2 O] entry time series, the BDC, and ∆T indices are all dominated by long-term trends.In such a case, an MLR may produce a high adjusted R 2 even if there is no actual relation between the variables.To eliminate the influence of long-term trends on adjusted R 2 , we detrend each variable using a Fourier Transform filter (Donnelly, 2006) to remove long-term variability (> 10 years).We then use the MLR on the detrended [H 2 O] entry and the detrended indices.Detrending by removing the long-term linear trend yields similar results.
Figure 1 shows the adjusted R 2 for the detrended calculation.For most of the models, the adjusted R 2 for the detrended MLR is only slightly smaller than that for the trended one.This confirms that the long-term trends in the data tend to inflate the adjusted R 2 , at least a bit, and also that the models' interannual variability and long-term trends are well represented by the same linear model (Equation 1).Large differences do exist for some CCMs.For instance, the CCSRNIES trended century MLR captures approximately 90% of the variance in [H 2 O] entry , while the detrended century MLR only explains about 40% of interannual variance; the CNRM-CM5-3, NIWA-UKCA, and WACCM show something similar.

Physical Process Effects
The coefficients from the trended and detrended calculations are listed in Tables 2 and 3 respectively.The product of the regression coefficient and its index quantifies that process' impact on Figure 3 shows that [H 2 O] entry increases as ∆T increases in all models and that the ∆T regression coefficients are similar for both trended and detrended MLRs.On average, [H 2 O] entry increases by about 0.3±0.1 ppmv K −1 , with individual models yielding values ranging from about 0.1 to 0.6 ppmv K −1 .This confirms that the stratospheric water vapor feedback identified by Dessler et al. (2013) occurs in all CCMs, although the exact magnitude varies.
This figure also shows that the BDC coefficient is generally negative, meaning that a strengthening BDC reduces This is consistent with previous research, which showed that a stronger BDC reduces TTL temperatures and lower stratospheric water vapor (Randel et al., 2006;Gilford et al., 2016).The trended and detrended BDC coefficients are similar in sign and magnitude.Two models (CNRM-CM5-3 and NIWA-UKCA) yield positive BDC coefficients, indicating potential problems with these models.
Figure 3 shows that all QBO regression coefficients are small, generally within ± 0.04 ppmv, with even the sign of the effect in doubt.Interestingly, one of the CCMs not simulating a QBO, CMAM-CCMI, produces the largest QBO regression coefficients of 0.082 ±0.04 and 0.077 ±0.04 ppmv for the trended and detrended calculations, respectively.Among CCMs that do simulate a QBO, the ensemble average QBO regression coefficient does not differ much from the same quantity (approximately 0 ppmv) for the other models.We will discuss this further in the next section.

Decadal Analysis
Ideally, we would compare the results of the last section to observations.Unfortunately, we don't have 100 years of observations to test the models against.Instead, we will compare regressions of 10-year segments from the CCMs to regressions of 10-years of observations.Specifically, we split 21 st century of each CCM run into 10 decades (2000-2010, 2010-2020, 2020-2030, 2040-2050, etc.) and fit each individual decade using the regression model (Equation 1).The regression calculation used on each 10-year segment is identical to the century analysis, except monthly averaged anomalies are used instead of annual mean anomalies.Figure 4 shows the median ± one standard deviation of the ten decadal adjusted R 2 values generated by each CCM.The ensemble average is approximately 0.6±0.25, with some spread among the models.Also plotted are the adjusted R The MLS data covers the time period 2004-2014.
Many of the models have a range of adjusted R 2 values that overlaps with the observational regression However, not all do: the CCSRNIES, CNRM-CM5-3, and NIWA-UKCA have median decadal adjusted R 2 values below 0.4, well below the observational values.It's worth nothing that these models also had issues in the century regressions.The WACCM and LMDZrepro models also have median adjusted R 2 values below the observations.
Figure 5 shows the median and one standard deviation of each coefficient (values are listed in table 4), along with the coefficients from the regression of the MLS data (taken from Dessler et al. (2014)).We find that the CCMs agree unanimously that increases in ∆T are associated with increased [H 2 O] entry .Overall, though, the CCM ensemble tends to underestimate the observational estimate, although most fall within the observation's 95% confidence intervals.The only models that don't fall within both observational ranges are CCSRNIES, CMAM-CCMI, and CNRM-CM5-3.
In addition, the spread between the different decades for a single model tends to be small, with most CCM decadal ∆T regression coefficient distributions confined to a narrow range of ±0.1 ppmv K −1 around the model's median.This gives us some confidence that the comparison between the CCMs and one decade of observations is meaningful.
Figure 5 shows that there exists a high degree of variability in the CCMs' decadal BDC regression coefficients, with a CCM ensemble average value of about -4±2 ppmv (K/Day) −1 , but with individual CCM values ranging between approximately -12 and +5 ppmv (K/Day) −1 .On all timescales, we expect a strengthening BDC should cool the TTL and reduce [H 2 O] entry , so the coefficient should be negative.We see that the median is indeed negative for all CCMs except for the CNRM-CM5-3 and NIWA-UKCA, both of which yield a positive median BDC coefficient (these models also generated positive BDC coefficients for the century analysis).
Comparing to observations, we find that the model ensemble does well.This nonetheless hides a significant spread among As expected, figure 5 shows that, for CCMs not simulating a QBO, the ensemble average decadal QBO coefficient is approximately 0 ppmv.But even those that do simulate a QBO, as seen in the century analysis, see little impact on [H 2 O] entry from it, with an ensemble average of approximately 0.03±0.04ppmv.This is significantly smaller than the response to the QBO in the observations, and this appears to be a clear deficiency in the model ensemble.
Only CCSRNIES-MIROC3.2 and CMAM-CCMI decadal regressions produce QBO coefficients approaching those from both observational regressions.Again, CMAM-CCMI does not simulate a QBO, and it is not clear to us why the model does so well in this aspect of our analysis.
Previous studies found that the QBO significantly influences TTL temperatures and subsequently [H 2 O] entry (Zhou et al., 2001;Geller et al., 2002;Liang et al., 2011), so the lack of response in the model ensemble seems to be a problem for the models.Previous studies have investigated this issue finding that a higher vertical resolution within the stratosphere can help resolve the QBO's impact on the lower stratosphere (Rind et al., 2014;Anstey et al., 2016;Geller et al., 2016).Clearly, this needs to be investigated further.

Century and Decadal Regression Coefficient Comparison
One interesting question is whether the regression coefficients from the decadal analyses are related to regression coefficients from century regressions.To answer this, Figure 6 shows the coefficients from the trended century regressions of each CCM plotted against the median of the decadal regressions from the same CCM.Also shown is a linear least-squares fit to the points.
As in the last section, uncertainties in the observational coefficients are bound by 95% confidence intervals calculated by Dessler et al. (2014).Uncertainty in the slope, intercept, and century regression predictions are constrained by 95% confidence intervals determined using each least-squares fit.
For the ∆T coefficient, the best fit line is: Thus, the ∆T coefficients from the trended MLRs are slightly larger than those from the decadal MLRs.Using values of β(∆T, decade) from decadal observations and this fit, we can predict β(∆T, century).From equation 2, the observed β(∆T, decade) correspond to β(∆T, century) of 0.50 ±0.06 and 0.55 ±0.08 ppmv K −1 for MERRA and ERAI indices, respectively.
For the BDC coefficient, the best fit line is: The BDC coefficients from the trended MLRs are also slightly larger than those from the decadal MLRs.

Conclusions
Climate models predict that tropical lower stratospheric humidity ([H 2 O] entry ) will increase as the climate warms, with important implications for the chemistry and climate of the atmosphere.We described in this paper a new way to evaluate the realism of these model trends.Our method is based on regressing CCM [H 2 O] entry time series against three processes (tropospheric temperature (∆T ), the strength of the Brewer-Dobson circulation (BDC), and the phase of the QBO) that have been shown to be important to [H 2 O] entry .We do this on two separate time-scales: 1) the 21 st century, and 2) on decadal timescales.
We find that long-term increase in [H 2 O] entry in the CCMs is primarily driven by warming of the troposphere.This is partially offset in most CCMs by an increase in the strength of the Brewer-Dobson circulation, which tends to cool the tropical tropopause layer (TTL).The models show little impact from the QBO, in disagreement with observations; this appears to be a deficiency in the models.
The coefficients from regressions of individual decades in the CCMs can be compared to coefficients from regressions of observations covering a decade.Overall, the CCM ensemble seems to reproduce observations well, except for the fact that the CCMs simulate little response of [H 2 O] entry to the QBO, in disagreement with the observations.In addition, the good agreement on average hides some spread among the models, particularly in the response to the BDC.
Our approach provides more insight into model processes than simply comparing [H 2 O] entry or TTL temperatures.Our overall conclusions are encouraging -the models appear do a reasonable job simulating variability in [H 2 O] entry .However, some models have clear problems, e.g., the models that predict [H 2 O] entry will increase with a strengthening BDC.In addition, nearly the entire ensemble does not reproduce the observed variations of [H 2 O] entry with the phase of the QBO.This analysis should help the modeling groups refine their models' simulations of the 21 st century.

Data availability
This data can be obtained through the British Atmospheric Data Center (BADC) archive.
Author contributions.KS and AD performed this analysis and wrote most of this manuscript.The other authors contributed information pertaining to their individual models and helped revise this paper.
by the NASA Center for Climate Simulation (NCCS).OM acknowledges funding by the New Zealand Royal Society Marsden Fund (grant no.12-NIW-006).OM and GZ wish to acknowledge the contribution of NeSI high-performance computing facilities to the results of this  The units of ∆T , BDC, and QBO are ppmv K −1 , ppmv (K/Day) −1 , and ppmv.The uncertainty is the 95% confidence interval.
Atmos.Chem.Phys.Discuss., doi:10.5194/acp-2016-964,2016   Manuscript under review for journal Atmos.Chem.Phys.Published: 8 November 2016 c Author(s) 2016.CC-BY 3.0 License.FollowingDessler et al. (2014), decadal regression terms are lagged in order to maximize MLR fit: we lag ∆T by 3 months, the BDC by 1 month, and the QBO by 3 months.These lags reflect the time between changes in each index and the impact on [H 2 O] entry .

Figure 1 .Figure 2 .
Figure 1.Bars corresponds to trended (light grey) and detrended (dark grey) adjusted R 2 values for annual-averaged data.The light grey circle represents the CCM ensemble mean trended adjusted R 2 value, while the dark grey circle represents to the CCM ensemble mean detrended adjusted R 2 value.Error bars on ensemble means corresponds to the ± one standard deviation of the CCM ensemble.

Figure 3 .Figure 4 .Figure 5 .Figure 6 .
Figure 3. Circles represents the detrended (light grey) and trended (dark grey) coefficients for each model, and error bars correspond to 95 th percentile confidence interval bounding each regression coefficient.An asterisk indicates models simulating a QBO.An asterisk on the ensemble mean corresponds to the average QBO coefficient for only models simulating a QBO, while the ensemble mean with no asterisk corresponds to the average of all model coefficients.The ensemble mean coefficients are also represented by a circle, with associated error bars correspond to ±one standard deviation of the ensemble set of coefficients.The units of β∆t,βBDC , and βQBO are ppmv/K, ppmv/(K/Day), and ppmv, respectively.
[H 2 O] entry over the 21 st century.To do this, we calculate annual average values of [H 2 O] entry .As an example, MRI [H 2 O] entry increases by about 1.2 ppmv during the 21 st century (Figure 2).The regression shows that this is the result of a large increase in [H 2 O] (Dessler et al., 2014)ses ( 1.5 ppmv) that is offset by a strengthening BDC, which reduces [H 2 O] entry by approximately 0.3 ppmv; this is consistent with the results of Dessler et al. (2013).The regression finds virtually no change in [H 2 O] entry in response to the QBO, which does not comport with analyses of observations, which suggests that the QBO causes short-term variations in [H 2 O] entry of 0.3 ppmv(Dessler et al., 2014) (Rienecker et al., 2011)f the tropical average Aura Microwave Limb Sounder (MLS) 82-hPa water vapor mixing ratio observations fromDessler et al. (2014).One regression uses Modern-Era Retrospective Analysis for Research and Applications reanalysis (MERRA) data(Rienecker et al., 2011)and the other uses European Centre for Medium-Range Weather Forecasts interim re- (Dee et al., 2011)e et al., 2011)for the ∆T and BDC indices; the QBO index is standardized anomaly of monthly and zonally averaged equatorial 50-hPa winds obtained from the NOAA Climate Prediction Center (http://www.cpc.ncep.noaa.gov/data/indices).

Table 1 .
research.NZ's national facilities are provided by the NZ eScience Infrastructure and funded jointly by NeSI's collaborator institutions and through the Ministry of Business, Innovation & Employment's Research Infrastructure programme (https://www.nesi.org.nz).HA acknowledges the Environment Research and Technology Development Fund, Ministry of Environment, Japan (2-1303) and NEC-SX9/A(ECO) 10 computers at CGER, NIES.The LMDZ-REPRO contribution was supported by the European Project StratoClim (7th framework programme, Grant agreement 603557) and the Grant'SOLSPEC' from the Centre d'Etude Spatiale (CNES).CCMs used in this analysis.The resolution is listed as (lat x lon x number of pressure levels).31 vertical levels indicates CCM data is given on isobaric levels, while CCMs simulating data on >31 levels are given on sigma (hybrid-pressure) levels

Table 2 .
Coefficients from regressions of trended [H2O]entry time series

Table 4 .
Median coefficients from the decadal regressions of [H2O]entry monthly anomalies