Variability and trends in dynamical forcing of tropical lower stratospheric temperatures

The contribution of dynamical forcing to variations and trends in tropical lower stratospheric 70 hPa temperature for the period 1980–2011 is estimated based on ERA-Interim and Modern-Era Retrospective Analysis for Research and Applications (MERRA) reanalysis data. The dynamical forcing is estimated from the tropical mean residual upwelling calculated with the momentum balance equation, and with a simple proxy based on eddy heat fluxes averaged between 25 and 75 in both hemispheres. The thermodynamic energy equation with Newtonian cooling is used to relate the dynamical forcing to temperature. The deseasonalised, monthly mean time series of all four calculations are highly correlated (∼ 0.85) with temperature for the period 1995–2011 when variations in radiatively active tracers are small. All four calculations provide additional support to previously noted prominent aspects of the temperature evolution 1980–2011: an anomalously strong dynamical cooling (∼−1 to −2 K) following the Pinatubo eruption that partially offsets the warming from enhanced aerosol, and a few years of enhanced dynamical cooling (∼−0.4 K) after October 2000 that contributes to the prominent drop in water entering the stratosphere at that time. The time series of dynamically forced temperature calculated with the same method are more highly correlated and have more similar trends than those from the same reanalysis but with different methods. For 1980–2011 (without volcanic periods), the eddy heat flux calculations give a dynamical cooling of ∼−0.1 to ∼−0.25 Kdecade (magnitude sensitive to latitude belt considered and reanalysis), largely due to increasing high latitude eddy heat flux trends in September and December–January. The eddy heat flux trends also explain the seasonality of temperature trends very well, with maximum cooling in January–February. Trends derived from momentum balance calculations show near-zero annual mean dynamical cooling, with weaker seasonal trends especially in December–January. These contradictory results arising from uncertainties in data and methods are discussed and put in context to previous analyses.


Introduction
Temperature trends in the lower stratosphere over the last few decades are substantially larger than predicted by general circulation models (GCMs) forced with increasing greenhouse gases only (Ramaswamy et al., 2006).Consequently, trends in stratospheric dynamics and other radiatively active trace constituents (such as ozone) may play a major role for temperature trends in this layer, and decomposition into dynamically and radiatively forced trends has attracted much interest (see, e.g., recent analyses of Fu et al., 2010;Ueyama and Wallace, 2010;Polvani and Solomon, 2012;Bohlinger et al., 2014).
The dynamically forced diabatic upwelling in the tropical stratosphere lowers temperatures below radiative equilibrium, and an increase/decrease in the dynamically forced upwelling induces a decrease/increase in tropical lower stratospheric temperatures.In addition to this dynamical control, changes in radiatively active trace constituents, such as ozone, water vapour and aerosol, change the local bal-Published by Copernicus Publications on behalf of the European Geosciences Union.

S. Fueglistaler et al.: Tropical stratosphere trends
ance of emitted and absorbed radiation which then induces changes in temperature.While the interaction between dynamics, composition and radiation is in principle relatively well understood (Yulaeva et al., 1994;Avallone and Prather, 1996;Randel et al., 2006;Fueglistaler, 2012;Abalos et al., 2013), systematic errors in the time series of any of the relevant quantities (temperature, radiatively active trace constituents and atmospheric dynamics) have prevented an unambiguous interpretation of the evolution of the lower stratosphere over the past 3 decades.
Here we focus on quantification of the dynamical forcing of tropical lower stratospheric temperatures at the 70 hPa level over the period 1980-2011 based on meteorological reanalysis data from ERA-Interim (Dee et al., 2011) and Modern-Era Retrospective Analysis for Research and Applications (MERRA) (Rienecker et al., 2011).We use two methods to estimate the dynamically forced temperature variations, namely the momentum balance equation of Randel et al. (2002) that estimates the tropical mean residual upwelling, and the eddy heat flux proxy introduced by Fueglistaler (2012).The residual, i.e. temperature minus estimated dynamically forced temperature, represents the influence from changes in radiatively active trace constituents at the level of interest and changes in up-and downwelling radiative fluxes, and from (i) errors in the reanalysis temperature, (ii) errors in the method to estimate the dynamically forced temperature and (iii) errors in the reanalysis temperature and wind data used to calculate the dynamically forced temperature.
Of particular interest here are the following questions concerning the evolution of the tropical lower stratosphere over the last 3 decades.(i) Is there an increase in tropical upwelling (leading to dynamical cooling) in the observations as reported in climate models (e.g.Butchart et al., 2010)? (ii) To what extent are the observed temperature trends due to changes in dynamical forcing (as suggested for example by Fu et al., 2010), and to what extent due to changes in radiative active trace constituents (as suggested, e.g., by Polvani and Solomon, 2012)?(iii) Is the pronounced seasonality in tropical lower stratospheric temperature trends (e.g.Fu et al., 2010;Free, 2011;Polvani and Solomon, 2012) the result of seasonal changes in the dynamical forcing, or in radiatively active trace constituents?(iv) Are the new calculations presented here supporting previously published results concerning changes in the strength of the residual circulation following the eruption of Pinatubo in 1991 noted by Fueglistaler (2012), and the increase in the strength of the residual circulation after the year 2000 noted by Randel et al. (2006), which would explain the sudden and prolonged drop of water entering the stratosphere in that period?(v) How consistent are results concerning the dynamical forcing between ERA-Interim and MERRA, and between the two methods to estimate the dynamical forcing of temperatures?An important result of our analysis is that the four calculations give very similar time series, but that their trends over the period 1980-2011 differ significantly.While we cannot unambiguously determine which trends are more accurate with the information available here, we document and discuss which parts of the calculations lead to differing trends.
The paper is structured as follows.Section 2 presents the data and methods used to estimate the variations in dynamically forced temperature.The results are shown in Sect.3, with the evolution of the time series discussed in Sect.3.1, and trends discussed in Sect.3.2.Sections 3.3 and 3.4 document and discuss where the two reanalyses, and the two methods, differ.Finally, Sect. 4 presents the conclusions and outlook.

Framework
Consider the thermodynamic energy equation in the quasigeostrophic, transformed Eulerian mean set of equations, and assume that the diabatic heating term can be approximated by Newtonian cooling, where τ is the radiative relaxation timescale, T E is the radiative equilibrium temperature, N 2 is the buoyancy frequency, H is the scale height, R is the gas constant, w * is the diabatic residual velocity, T is temperature, t is time and all quantities are zonal means (tropical means in this work) (see, e.g., Eq. 1 of Yulaeva et al., 1994).
In this framework, we may force temperature T with w * and T E , where variations in T E are a consequence of variations in radiatively active tracers.Note that changes in w * may also affect tracers (here particularly ozone; see, e.g., Fueglistaler et al., 2011) and hence T E .It is well known that Newtonian cooling is a relatively poor approximation for the lower stratosphere (see, e.g., the analysis by Hitchcock et al., 2010), but in practice there exist reasonably linear relations between local heating rate, temperature and radiatively active tracer amounts for smooth perturbations with vertical scales similar to that of the lower stratosphere (Fueglistaler et al., 2011(Fueglistaler et al., , 2014)).Hence, for the conceptual interpretation of the connection between dynamics, tracers and temperatures, Eq. ( 1) is sufficient.
Radiative transfer calculations (not shown) show that the first-order radiatively active trace constituents that need to be considered for the lower stratosphere are aerosol and ozone; variations with observed magnitudes of water vapour, carbon dioxide and other long-lived greenhouse gases typically affect temperature at about less than 1 order of magnitude.
Under the assumption that temperature variations are only dynamically forced (i.e.w * in Eq. 1), the equilibrium temperature is constant and the solution for the difference between temperature and equilibrium temperature, where f (t) is the forcing, * denotes the convolution and the unit pulse response h(t) = 0 for t > 0 and h(t) = e t/τ for t ≤ 0.
We use two methods to estimate the forcing function f (t).First, we use the momentum balance calculation of Randel et al. (2002) to estimate the tropical mean residual vertical velocity w * (Sect.2.3).Second, we use the eddy heat flux proxy of Fueglistaler (2012) (Sect.2.4).For both calculation, the timescale τ is empirically determined by optimising the correlation of deseasonalised, monthly mean temperature and T dyn .The highest correlations (determined for the period 1995-2011 when radiative perturbations are small) are found for τ ∼ 70 days for the momentum balance equation, and τ ∼ 80 days for the eddy heat flux calculation (see Sect. 2.4).We evaluate the convolution (Eq.2) over a length of three times the e-folding timescale, and correspondingly the first year of data (1979) are omitted, and results are shown for 1980-2011.With values around 80 days, the empirically determined timescale τ is longer than the τ ≈ 40 days expected from idealised radiative transfer calculations (e.g.Fueglistaler et al., 2014).The longer-than-expected timescale may reflect variations in T E (t) that are correlated, but phase lagged to the dynamical forcing f (t).Hence, T dyn is to be understood as being the dynamically forced temperature including variations in T E that are correlated to the dynamical forcing.
The residual, thus represents changes in temperature due to (i) radiative changes in T E (t) from tracer changes not correlated with dynamics as represented by T dyn (t), (ii) errors in the reanalyses' temperatures and winds and (iii) errors from dynamical forcing terms not included in the calculation.
With the focus of this analysis on variability and trends over the period 1980-2011, the data are analysed in terms of monthly means.Time series shown are deseasonalised by subtracting the climatological mean of each month of the year determined from detrended data.Since volcanic eruptions lead to major disturbances in the radiative conditions (i.e.T E is not constant), we chose 1995-2011, when stratospheric aerosol variations are small, as reference period for the mean annual cycle.Linear trends are reported for each month of the year, and for annual means of the data.Trends are shown for the full period 1980-2011, and for the period 1980-2011 with the El Chichón and Pinatubo periods excluded (labelled no-volc).For consistency between trends with monthly and annual data, 3 full calendar years are excluded (i.e. no-volc excludes the years 1982, 1983, 1984 and 1991, 1992, 1993).

Data
Temperature and winds are taken from ERA-Interim (Dee et al., 2011) and MERRA (Rienecker et al., 2011), with ERA-Interim data on the original model levels and MERRA on the (standard) pressure levels of the archived data.Results are shown for the 70 hPa level for MERRA, and for the closest model level of ERA-Interim, which is at 67 hPa.We have also calculated the eddy heat flux proxy using NCEP (National Centers for Environmental Prediction) (Kalnay et al., 1996) and DOE (Department of Energy) NCEP-2 data; the results from these calculations are fairly similar to those from ERA-Interim and MERRA, and will be mentioned only in specific places to allow comparison with previously published studies that use NCEP data.All terms are calculated from 6-hourly (00:00/06:00/12:00/18:00 UTC) analysis data, at a horizontal resolution of 1 • /1 • for ERA-Interim and 1.25 • /1.25 • for MERRA (and 2.5 • /2.5 • for NCEP).
The reanalysis temperatures are complemented with measurements from the microwave sounding unit (MSU) lower stratospheric channel (temperature lower stratosphere, TLS, also known as MSU channel 4) as provided by remote sensing systems (RSS) (Mears and Wentz, 2009).The TLS channel's weighting function is centred near 70 hPa, but the deep vertical weighting function is problematic for interpretation of TLS data in the tropics (Fueglistaler, 2012).The data are shown here for context, and because it provides a tropical mean without sampling bias, and has a strong impact on lower stratospheric temperatures in reanalyses.Randel et al. (2002) show that the time-varying zonal mean tropical average upwelling can be calculated from the vertical integral of the Eliassen-Palm flux divergence and acceleration at the bounding latitudes (see Eq. 11 in Randel et al., 2002).Unlike Randel et al. (2002), we calculate the divergence from the complete Eliassen-Palm fluxes.Problems in this calculation may arise from wave activity not resolved in the reanalyses and, probably more importantly, from spurious drifts in the assimilated observational data (analysed meteorological data are not required to be energetically self-consistent; see, e.g., the discussion in Fueglistaler et al., 2009a).The bounding latitudes used here are 30 • S and 30 • N. The calculated estimate of the tropical mean residual velocity w * mb (subscript mb refers to momentum balance) is converted to temperature using the Newtonian cooling approximation (see Eqs. 1 and 2); with the solution

Estimate of temperature variations based on momentum balance
where The correlation with temperature maximises (not shown) for a timescale of τ ∼ 70 days.Interestingly, with τ = 70 days the coefficient K determined with an ordinary leastwww.atmos-chem-phys.net/14/13439/2014/squares linear regression (or total least-squares fit; implications of this calibration on trends are discussed in Sect.3.2) of T mb against actual temperature is very similar to the value obtained by inserting typical values for N 2 and H for the lower stratosphere (N 2 = 6 × 10 −4 s −2 , H = 7000 m).This suggests that the magnitude of the residual vertical velocity calculated with the momentum balance equation is actually quite accurate (but recall the previous discussion that the empirically determined τ may account for a modification of temperature variations due to correlated trace gas changes).

2.4
The eddy heat flux proxy for tropical temperatures Fueglistaler (2012) shows that tropical lower stratospheric temperature variations around 70 hPa are well correlated with eddy heat fluxes on the same level, and averaged over a wide latitudinal belt in both hemispheres.The approach is motivated by the success of a similar approach to predict temperatures in the stratospheric polar vortices by Newman et al. (2001), and a number of related quantities have been used recently to study stratospheric temperatures.However, the fact that the simple calculation of Fueglistaler (2012) explains tropical temperatures very well (with correlation coefficients of ∼ 0.85 for monthly mean data) is an empirical result, indicating that the many neglected aspects (for example where exactly the waves dissipate, and the role of the momentum fluxes) are apparently of secondary importance for variations of monthly mean temperature.Here we discuss and analyse the sensitivities and limitations of this simple proxy in more detail.We note that, while we use an e-folding of heat fluxes rather than period averages as the studies listed above, our results concerning trends also apply to these studies as the e-folding timescale (or averaging period) is much shorter than the 1980-2011 period for trends.
Similarly to the calculation using the momentum balance equation, we consider the thermodynamic energy equation with Newtonian cooling and simply assume that the forcing from the vertical residual velocity is proportional to the average eddy heat flux.Further, as before, the dynamical forcing f (t) is related to tropical average temperatures using a linear regression to calculate the proportionality coefficient α between the temperature and the e-folded average eddy heat flux (based on Eq. 1 and 2): where T p is the proxy time series, v T is the eddy heat flux and the operator [•] 75 25 is the zonal mean averaged between 25 • and 75 • latitude in both hemispheres.As before, * de- notes the convolution, and h(t) = e t/τ for t ≤ 0 and h(t) = 0 for t > 0. Note that α < 0 (larger eddy heat fluxes force lower temperatures in the tropics).The free parameters are the relaxation timescale τ , and the latitude band over which the eddy heat flux is averaged.The values taken here for these parameters are based on the following analysis.
Figure 1a and b show the correlation coefficients between tropical average temperature and T p (t) at the same pressure level, for ERA-Interim (panel a) and MERRA (panel b).The correlations are plotted as function of end latitude of the averaging domain of the eddy heat fluxes, with integration starting at 25 • latitude (diamonds) and integration starting at 80 • latitude (triangles).The figure shows that for 70 hPa, the proxy's skill is primarily due to the heat fluxes between 40 • and 80 • latitude.For this study we use the average between 25 • and 75 • latitude, but correlations for any proxy that includes the range between about 40 • and 70 • latitude are similar.Trends, however, are sensitive to choice of latitudes (further discussion in Sect.3.4).
In addition to the 70 hPa level, the figure also shows results for the 80 hPa and 96 hPa levels for ERA-Interim, and the 100 hPa level for MERRA.While the correlations are very similar for the 80 hPa level, for 100 hPa the correlations between proxy and temperature are, as expected, weaker (the 100 hPa level in the tropics is in the middle of the TTL, and tropical processes play an important role; see discussion in Fueglistaler et al., 2009b).
Model intercomparison projects (e.g.Eyring et al., 2010) often archive (only) the 100 hPa eddy heat fluxes due to their role for polar vortex temperatures (Newman et al., 2001).It is therefore of interest to what extent the 100 hPa eddy heat flux may be used to explain temperatures at 70 hPa.Figure 1a and  b show that for ERA-Interim and MERRA, the 100 hPa eddy heat flux is, in terms of correlation, as good as the 70 hPa eddy heat flux for temperatures at 70 hPa.
Figure 1c shows the correlation between the deseasonalised eddy heat flux proxy and temperature for the period 1995-2011 as function of the relaxation timescale τ .For ERA-Interim and MERRA, correlations between tropical average temperature and the proxy maximise for relaxation timescales around 80 days (i.e. about 15 % longer than for the momentum balance equation) with a correlation coefficients ∼ 0.85 for the period 1995-2011.For NCEP, the corresponding maximum correlation coefficient is 0.75 for a relaxation timescale of ∼ 60 days.

Removal of QBO variability in tropical temperatures
Tropical stratospheric temperatures are strongly affected by the quasi-biennial oscillation (QBO; see, e.g., Baldwin et al., 2001), whose forcing is not well resolved in reanalysis data.The acceleration of the flow in the inner tropics produces a see-saw pattern in temperature between inner and outer tropics that distorts the relation to meridional eddy heat fluxes.The latitudinal scale of the QBO (order 1000 km) is substantially shorter than the scale of the mid-and high latitude eddy forcing, and this scale separation may be used to reduce the variance in tropical temperatures due to the QBO.
Figure 2 shows the correlation and slope of the linear regression of deseasonalised zonal mean zonal wind shear at the equator with the zonal mean deseasonalised temperature for each latitude at 70 hPa.Figure 2b shows the slope of the linear regression integrated from the equator to the poles as function of the end latitude (in both hemispheres) of the integral.The figure shows that there exists no latitude belt (symmetric about the equator) where the effect of the QBO on temperature exactly cancels.For this study we have cho- sen the latitude belt 35 • S-35 • N as a compromise between a minimisation of QBO-related variance and focus on tropical temperatures.The choice of 35 • S-35 • N for temperature increases the correlation with the eddy heat flux proxy slightly, but does not affect trends (i.e.temperature trends differ by no more than a few percent between 30 • S-30 • N and 35 • S-35 • N).

Results
In the following, Sect.3.1 discusses the deseasonalised time series, specifically also the Pinatubo period and the pre-and post-2000 periods.Section 3.2 discusses the linear trends for the period 1980-2011 in annual mean data, and the seasonality of trends.Sections 3.3 and 3.4 discuss differences between ERA-Interim and MERRA, and between the results using the two different methods.

Time series
Figure 3a shows the deseasonalised global average, the tropical average and the combined extratropical average of the RSS MSU lower stratospheric channel TLS.The strong anticorrelation between tropical and combined extratropical av-   c).Correlations between the estimates of dynamically forced temperature are given in Table 1, between temperature and dynamical estimates in Table 2, and between the residuals in Table 3.All panels show deseasonalised data, with the mean annual cycle evaluated over the period 1995-2011.
erage temperatures is readily visible, though the fact that they nearly perfectly cancel in the global average in response to variations in the strength of the residual circulation is a coincidence (Fueglistaler et al., 2011).Also, the apparent stepwise decrease in global average temperatures following the eruptions of El Chichón and Pinatubo noted by Pawson et al. (1998) is well visible.Figure 3b and c show the tropical average temperatures near 70 hPa, and the dynamically forced temperature variations predicted by the momentum balance equation and the eddy heat flux proxy for ERA-Interim (panel b) and MERRA (panel c). Figure 3d shows the residual (temperature minus dynamically forced temperature).
Figure 4 shows the time series of annual means of the dynamically forced temperature variations for both methods, and for both ERA-Interim and MERRA.The time series of annual mean data are used to calculate trends, and some aspects of the time series are more clearly seen in the annual mean data.

Correlations
The correlations summarised in Tables 1, 2 and 3 confirm the visual impression from Fig. 3b and c of a remarkably good agreement between the two methods to estimate the dynamically forced temperature variations, and between the estimates of dynamically forced temperature variations and actual temperatures (with the expected exception during the volcanic periods).
Table 1 shows that the correlations between the four calculations for dynamically forced temperature variations have correlations larger than 0.7, with higher correlations for the (shorter) period 1995-2011.The highest correlations are observed between the estimates that use the same method (but different reanalysis), with correlations between the two eddy heat flux proxy calculations being as high as 0.96 (for 1995-2011) and 0.97 (for 1980-2011).This implies that the dynamical information in ERA-Interim and MERRA is very similar.The correlations between the eddy heat flux proxy and the momentum balance equation based on the same reanalysis are still high (1995-2011: 0.83 for ERA-Interim and 0.81 for MERRA), but the fact that they are smaller than those between calculations using the same method indicates that there is a small but noticeable difference in the quantity represented by the two methods.
Table 2 shows the correlations between the calculations of dynamically forced temperature variations and the actual temperatures for ERA-Interim and MERRA.The table shows that for both the 1995-2011 and 1980-2011 (without volcanoes) periods both calculations of dynamically induced tem-perature variations agree well with the actual temperatures (with correlations as high as 0.87).For the period 1995-2011, the two methods perform equally well, while for the full period 1980-2011 (excluding volcanoes) the eddy heat flux calculations are better correlated (correlations of 0.83 for both ERA-Interim and MERRA).This is both due to different trends in the two methods (discussed below), and due to a generally lower agreement between temperature and the estimated dynamical temperature in the 1980s particularly for the momentum balance equation (for the period 1980-1990 (excluding volcanoes), the correlation for the momentum balance calculation for MERRA is only 0.51).With the information available here, we cannot address the question as to what extent the lower correlation in the 1980s is due to larger errors in temperature, larger errors in the reanalyses' dynamical information or larger variations in temperature due to radiatively active trace constituents.
Finally, Table 3 shows the correlations between the residuals for the period 1995-2011 when variations in T E due to changes in radiatively active trace constituents are expected to be small, and the residuals may largely reflect errors in the estimates of the dynamical forcing and errors in actual temperature.The residuals are correlated the highest between estimates using the same method but different reanalyses (residuals based on eddy heat flux proxy: 0.90; residuals based on momentum balance: 0.67).This result may be expected from the correlations between the dynamical temperature estimates (Table 1), but note that the residuals have additional information, namely the ERA-Interim and MERRA actual temperatures.However, the ERA-Interim and MERRA temperatures are very similar for 1995-2011 and hence introduce only a small amount of additional variance.Conversely, the residuals using different methods show lower correlations (for MERRA, the correlation of the residuals from the momentum balance calculation and the eddy heat flux proxy is as low as 0.36), which supports the previous notion that the main difference between the four calculations is that two methods represent two slightly different quantities.

Volcanic periods
Figure 3d shows that the residual is -as expected -largest during the El Chichón and Pinatubo volcanic periods.The four estimates of the dynamically forced temperature during the Pinatubo period agree well, and all show the anomalously strong forcing noted by Fueglistaler (2012) in ERA-Interim eddy heat flux proxy calculations.Hence, the calculations show additional strong evidence that this response is real and not due to problems in either the eddy heat flux proxy or ERA-Interim data.A similar response is not observed following the eruption of El Chichón, which suggests that there may not be a canonical dynamical response to volcanic eruptions.Further analysis of the dynamical situation following the Pinatubo eruption is beyond the scope of this paper and  will be presented elsewhere.Finally, we note that the uneven distribution of volcanic eruptions substantially distorts linear trend calculations for temperature, dynamically forced temperature and residual.We recommend that stratospheric trends reported in the literature exclude the volcanic periods in order to allow better comparisons between results using slightly different reference periods.

The year 2000 period
Water entering the stratosphere in the tropics shows a pronounced drop around the year 2000 that lasts for several years and is arguably the most prominent change in stratospheric composition that has been observed over the last 3 decades.Fueglistaler (2012), based on eddy heat flux proxy calculations using ERA-Interim, locates the onset of the dry period to October 2000 as a result of a dynamical cooling induced by anomalous southern hemispheric dynamics.The eddy heat flux proxy based on MERRA, as well as the momentum balance calculations of both reanalyses (see Fig. 3b,  c), confirms this result.Randel et al. (2006) attribute the prolonged dryness over several years to an increase in the strength of the residual circulation of the post-2000 period relative to the pre-2000 period.Table 4 summarises our results for the dynamically induced temperature change over that period for the 70 hPa level.The table shows that (also visually clearly seen in Fig. 4) results depend on the length of comparison period chosen, and that the eddy heat flux calculations give a clearer separation for the pre-/post-2000 years.The smallest temperature difference is observed in the ERA-Interim momentum balance calculation, while the other three calculations give a period-averaged decrease in dynamically forced temperature at 70 hPa of about −0.4 K for the period 2001-2005 vs. 1995-2005.Hence, the calculations shown here support the results of Randel et al. (2006).We note that the eddy heat flux results presented here for the 70 hPa level show a much clearer change around the year 2000 than those presented in Fueglistaler (2012) that were based on 96 hPa eddy heat fluxes, and implications for water entering the stratosphere (which is controlled by conditions around 90-80 hPa rather than 70 hPa) remain to be analysed in more detail.Further, we emphasise that the strengthening of the residual circulation is only true for a few years around the year 2000 -for the entire 1990s versus the 2000s the residual circulation has been weakening (i.e.dynamical warming of the tropics) due to the anomalous dynamics following the Pinatubo eruption.

Conclusions concerning the methods
To summarise this section, we note that the dynamical information in ERA-Interim and MERRA used for the two meth-ods to calculate the dynamical temperature variations is very similar, and results using the same method are more similar than using the same reanalysis but different methods.The implications for the residual are that the residual temperature time series represents not only variations in T E not correlated with dynamical forcing but also errors from incorrect representation of the dynamical forcing in one or both methods.Large events, however, such as the Pinatubo period and the pre-/post-2000 periods, are similar for both methods and both reanalyses.

Trends
The evolution of temperature and of the dynamically forced temperatures shown in Fig. 3 shows a large amount of variability even when excluding the volcanic periods, and the time series may be too short to reliably detect a trend.However, linear trends provide an objective metric to compare the time series quantitatively, and for temperature linear trends are routinely reported in the literature.We therefore calculate and discuss the linear trends in our estimates of dynamically forced temperature, but emphasise that we regard these trends simply as the first moments of the time series.
The trends in dynamically forced temperature -and hence also in the residual -depend on the trend in the dynamical forcing f (t) (Eq.2).Further, the magnitude of the trends depends also on the empirically determined proportionality coefficients K (Eq.4) and α (Eq.5).The 1σ uncertainty in the linear regressions to determine K and α is about 10 % in all calculations.The coefficients determined with a total least-squares fit are about 20 % larger than when determined with an ordinary least-squares regression.Both of these uncertainties are smaller than the differences between the four calculations, and we therefore focus the discussion primarily on the latter.

Overview of trends
Figure 4 shows the trends calculated from the annual mean data, and Fig. 5 summarises the trends (annual means and seasonal cycle) for all calculations for the full period 1980-2011 (Fig. 5, panels a-d), and for 1980-2011 without volcanic periods (Fig. 5, panels e-h).For visual clarity, error bars are not shown.Table 5 summarises the numerical values of the trends based on annual mean data including their statistical 1σ uncertainties.
Figure 5 shows that, as argued above, the volcanic periods have a substantial impact on trends.Therefore, we focus below on trends based on data excluding volcanic periods, but the following aspect is noteworthy with respect to previously published temperature trends that did not exclude the volcanic periods: the largest impact of volcanic periods is a change in the annual mean trend in the residual that is dominated by the impact of the volcanic periods on trends in actual temperature.Conversely, the seasonality of trends is less affected, and the anomalous dynamically forced cooling following Pinatubo does not dominate linear trends since the event is roughly in the middle of the period 1980-2011.

Trends in annual means
Table 5 summarises the linear trends based on annual mean data for actual temperatures in ERA-Interim and MERRA.The table and Fig. 5 show that the largest trend differences between ERA-Interim and MERRA are about Temp.−0.35 ± 0.10 −0.18 ± 0.11 T(v T ) −0.24 ± 0.08 −0.12 ± 0.09 T(m-bal) 0.03 ± 0.07 0.00 ± 0.10 0.15 K decade −1 in actual temperature.Trends in dynamically forced temperature are similar between ERA-Interim and MERRA when using the same method.However, the difference in trends between estimates using the momentum balance calculation and the eddy heat flux calculation is about 0.2 K decade −1 .The eddy heat flux calculation gives a statistically significant cooling of the tropical lower stratosphere when excluding volcanic periods in both ERA-Interim and MERRA.Conversely, the momentum balance calculations give no statistically significant trend, or even a slight warming.Similar results are obtained when using a total least-squares fit for α and K, except that in this case the magnitude of the trends is about 20 % larger.Finally, we note that, to the extent that the two calculations are comparable, the cooling implied by the homogenised ozone data of Hassler et al. (2013) as shown by Solomon et al. (2012Solomon et al. ( ) (order −3 K for period 1995Solomon et al. ( -1997Solomon et al. ( vs. 1979Solomon et al. ( -1981) ) is larger than the largest (no-volc) trend in our residual, or the difference between their two periods qualitatively estimated from our Fig.3d.

Seasonality of trends
Temperature trends in the tropical lower stratosphere show a prominent seasonality (seen in Fig. 5 ) that has has attracted much attention (e.g.Fu et al., 2010;Free, 2011;Polvani and Solomon, 2012).Fu et al. (2010) conclude, based on an analysis using eddy heat flux calculations, that the seasonality is largely driven by a seasonality in the dynamical forcing of temperatures.Conversely, Polvani and Solomon (2012) show that in a GCM with prescribed sea surface temperatures, greenhouse gases and ozone, the model's seasonal cycle of temperature trends is similar to observations only when the model is forced with observed ozone changes (which have a pronounced seasonality).The results shown in Fig. 5 show that not only do the annual mean trends differ particularly between the two methods, but also the seasonality.Figure 6a shows the trends in the dynamically forced temperature of the four calculations in one panel (data are identical to that shown in Fig. 5 e-h).It is evident that the largest differences occur during boreal winter, which is also seen in the seasonality in the residual trend (Fig. 5 e-h, shaded in blue/red).Figure 6b shows the fractional trends in the dynamical forcing (f (t), Eq. 2) of the four calculations.The figure shows that from March to October the four estimates roughly agree, but that the magnitude of the trend in the momentum balance calculation for November, December and January is much smaller than in the eddy heat flux calculations.The phase shifting due to the radiative damping then leads to the marked differences in dynamically forced temperatures in January and February.Similarly, the differences in dynamical forcing during December and January largely account for the annual mean trend differences between the two methods.None of the four calculations, however, show the seasonal structure of upwelling trends (indicating an increase in upwelling during boreal summer) in the GCM calculations of Polvani and Solomon (2012) (compare with their Fig.7b).

Eddy heat flux divergence between 100 and 70 hPa
The trends in dynamically forced temperature based on the eddy heat flux proxy are similar in ERA-Interim and MERRA (Table 5), but the difference is large enough to warrant further analysis.Here we focus on the evolution of the time series and identify when the time series diverge.In Sect.3.4 we also show at which latitudes eddy heat flux trends differ most.
Figure 7a shows the difference ERA-Interim minus MERRA eddy heat flux proxy at 70 hPa and at 100 hPa.The figure shows that from about the year 2000 onwards, ERA-Interim drifts (cooling) relative to MERRA. Figure 7b shows the ratio of the proxy calculated from 100 hPa eddy heat fluxes to that calculated from 70 hPa eddy heat fluxes.The figure shows that for MERRA this ratio remains approximately constant, while for ERA-Interim this ratio increases after about the year 2000.That is, in ERA-Interim the eddy heat flux at 70 hPa slightly increases, and at 100 hPa slightly decreases.Hence, in the ERA-Interim world a change in the eddy dissipation between 100 hPa and 70 hPa occurs in the early 2000s, which is not seen in the MERRA data (and is also not seen in NCEP data; result not shown).It is unclear at this point what exactly causes these drifts between the reanalyses.We consider a problem in ERA-Interim data more likely than what would be, arguably, an interesting signal.A hypothesis is that the onset of assimilation of GPS temperature data in ERA-Interim (from the early 2000s onwards CHAMP (CHAllenging Mini-satellite Payload), and after 2006 COSMIC (Constellation Observing System for Meteorology, Ionosphere and Climate); Poli et al., 2010;Dee et al., 2011) may induce a drift in ERA-Interim that is not present in MERRA (which does not assimilate these data Rienecker et al., 2011).

Differences between momentum balance and eddy heat flux proxy
Figures 4, 5 and Table 5 show a cooling due to an increase in dynamical forcing according to the eddy heat flux proxy both in ERA-Interim and MERRA, with very similar annual mean trends, and seasonality in trends.These results are similar to those presented by Fu et al. (2010), whose analysis is also based on eddy heat fluxes (albeit they used NCEP and ERA-40 data, and their latitude range and pressure levels differ from ours).Hence, eddy heat flux trends seem fairly robust across different reanalyses, with an overall increase in eddy heat fluxes over the last 3 decades (qualitatively consistent with the notion of an overall strengthening of the stratospheric residual circulation), and the seasonality of temperature trends being also a consequence of a shift in the time of year of wave activity (note the decreases in November and March bracketing the increases in December and January).
Conversely, the results using the momentum balance calculation indicate basically zero trend in upwelling (and temperature) between 30 • S and 30 • N, leaving a much larger trend residual that contributes also significantly to the seasonality in tropical temperature trends (i.e. according to the momentum balance calculations, the cooling seen in January/February is largely due to radiatively active trace constituents not correlated with dynamical forcing).
Both the eddy heat flux proxy and the momentum balance calculations may suffer from difficult-to-quantify problems in the reanalyses, and with the available data it is not possible to conclusively resolve which of the two calculations give more reliable trends.In the following, we analyse which aspects of the data used for the eddy heat flux calculations lead to the differences.Identification of the cause of the differences may motivate future work with targeted primitive equation model calculations to address the question as to whether the eddy heat flux proxy may be biased for trends, or whether the momentum balance calculations give erroneous trends due to drifts in the reanalysis input data.
Figure 8 shows the trends in monthly mean eddy heat fluxes seasonally resolved (panels a and b for ERA-Interim and MERRA, respectively), and for the annual mean data (panel c).The figure shows that the cooling trend in the eddy heat flux proxy arises from the high latitudes, while over the mid-latitudes the eddy heat fluxes may even decrease (particularly in MERRA).That is, in ERA-Interim and MERRA, the eddy heat fluxes shift polewards and overall intensify.Panels a and b show that the seasonality seen in Fig. 5 is a consequence of a temporal shift in eddy heat fluxes in the hemispheric winter months, in agreement with Fu et al. (2010).In the Northern Hemisphere, the eddy heat fluxes increase in December and January (giving maximum cooling in the tropics in January/February due to the phase lag from radiative temperature damping), and in the Southern Hemisphere the polewards eddy heat fluxes increase in September and, to a lesser extent, in October and December, while they decrease in November.
The latitudinal structure in eddy heat flux trends raises the question as to how robust the trend in dynamically forced temperatures is with respect to the specific latitude band over which the eddy heat fluxes are averaged.We find (not shown) that the trends in average eddy heat fluxes are reasonably robust against further widening of the latitude belt (i.e.extending further equatorwards than 25 • latitude and/or extending further polewards than 75 • latitude).However, trends are sensitive to narrowing the latitude belt (as may be anticipated from Fig. 8c).
Figure 9 shows the trends in dynamically forced temperatures based on the eddy heat flux averaged from 75 • latitude equatorwards, as function of the equatorward bounding latitude.Note that for this calculation, the forcing function (the average eddy heat flux) is dependent on not only the end latitude but also the proportionality coefficient α (Eq.5).
Figure 9 shows that for all latitude bands ending at high latitudes (here 75 • latitude), the eddy heat flux proxy calculations gives a dynamical cooling of the tropical lower stratosphere.The differences between ERA-Interim and MERRA dynamical temperature trends are due to the differences in mid-and low latitude eddy heat fluxes seen in Fig. 8c.The figure also shows that the trend uncertainty arising from the 1σ uncertainty in α is small.When using a total leastsquares fit to determine α, the latitudinal structure and difference between ERA-Interim and MERRA remain very similar, but the amplitudes of the trends are about 20 % larger (not shown).

Conclusions
We have estimated the variations in dynamical forcing of tropical lower stratospheric temperatures at 70 hPa over the period 1980-2011 due to variations in tropical upwelling.The dynamical forcing is calculated with the momentum balance equation for tropical average upwelling following Randel et al. (2002), and a simple proxy using the average eddy heat flux between 25 • and 75 • in both hemispheres following Fueglistaler (2012).All four calculations yield very high correlations (of about 0.85) with actual temperature for the period 1995-2011 when temperature variations due to variations in radiatively active trace constituents not correlated to the dynamical forcing are small.For the period 1980-1990 (excluding volcanic periods), correlations with the eddy heat flux proxy are slightly lower (between 0.77 and 0.80).For the momentum balance calculation, the correlations drop to 0.67 (ERA-Interim) and 0.51 (MERRA) for reasons not clear at present.
All four calculations confirm the previously noted anomalously strong dynamical cooling following the Pinatubo eruption (Fueglistaler, 2012), and the role of the stratospheric residual circulation for the marked drop of water entering the stratosphere at the end of the year 2000, initiated by anomalous southern hemispheric dynamics in October 2000 (Fueglistaler, 2012) and subsequent enhanced upwelling over a few years relative to the years before the year 2000 (Randel et al., 2006).
In all aspects considered, we find that the calculations using the same method, but different reanalysis (i.e.ERA-Interim and MERRA), agree better than the calculations using the same reanalysis but different methods.Hence, we conclude that the dynamical information in the two reanalyses is very similar, but that the two methods evaluate slightly different quantities.Differences between the results of the two methods are most pronounced in trends over the period 1980-2011.
The eddy heat flux proxy indicates a statistically significant cooling of about −0.1 K to −0.25 K for that period (excluding volcanic periods), with annual mean trends, and the seasonal cycle thereof, being strongly influenced by an increase in very high latitude eddy heat fluxes during December and January.The residual trends (due to trends in radiatively active trace species not correlated to dynamics, and/or errors in data and method) are relatively small for these calculations.Our results based on the eddy heat flux proxy are similar to those presented by Fu et al. (2010), which are also based on eddy heat flux calculations.
Conversely, the momentum balance calculations show no statistically significant trend over the period 1980-2011 (excluding volcanoes), mainly because they do not show the strengthening of upwelling in December and January seen in the eddy heat flux calculations.Consequently, with the momentum balance calculation the residual trends are roughly twice as large as those based on the eddy heat flux calculation, and the seasonality of actual temperature trends during December-April are mainly a consequence of the residual, which then has to be explained by trends in radiatively active trace constituents not correlated with dynamics.
The differences in trends of dynamically forced temperature between the two methods, and the discussed sensitivity of the magnitude of trends in the eddy heat flux proxy calculation to the latitude belt used, requires further attention.The results from the eddy heat flux proxy are highly plausible as they explain the observed temperature trends very well, but the tight relation between average eddy heat fluxes and temperature is an empirical result.The physically correct representation of tropical upwelling calculated with the momentum balance equation, in turn, gives a residual with a much larger seasonality (i.e.trends due to radiatively active tracers not correlated with dynamics would have to account for the cooling during boreal winter, while for boreal fall they would induce barely any temperature trend).
We have identified the very high latitude eddy heat flux trends during December and January to be the main source of the disagreement between the two calculations.Primitive equation model calculations with perturbations of high latitude eddy heat fluxes may be used to resolve the question as to whether the proxy calculations overestimate the impact of high latitude eddy heat flux trends on tropical temperatures.
Given the currently unresolved differences between methods, the role of dynamical forcing (and consequently also the role of trends in radiatively active trace constituents) for tropical lower stratospheric temperature trends remains unclear at present.Conversely, the excellent correlation of all four estimates of dynamically forced temperature with actual temperature variations is very encouraging and suggests that these calculations can be used to study the mechanisms of interaction between dynamics, radiatively active trace constituents and temperature in observation-based data; moreover, for specific periods of interest, the temperature evolution can be attributed to dynamics and radiatively active trace constituents with reasonable accuracy.
For example, Fu et al. (2010) use 3 month mean eddy heat fluxes averaged between 50 hPa and 10 hPa, Lin et al. (2009) use 3 month means of eddy heat flux between 45 and 90 • S at 150 hPa, Ueyama et al. (2013) use an index based on the 100 hPa eddy heat flux, and Bohlinger et al. (2014) use the 100 hPa eddy heat flux integrated from 45 to 70 • N, and averaged over 45 days.

31Figure 1 .
Figure 1.Correlation between deseasonalised tropical average temperature and dynamical temperature based on integrated eddy heat flux, convolved with e-folding timescale of 85 days, for pressure levels near 100 hPa, 80 hPa (where available) and 70 hPa (colour coded as labelled in panels).Integrations are started at 25 • (diamonds) and 80 • (triangles), correlations are plotted as function of the end latitude of the integration.Also shown (green) are the correlations between the e-folded eddy heat flux near 100 hPa with the temperatures near 70 hPa.Results are shown for (a) ERA-Interim, and (b) MERRA.(c) Correlation as function of relaxation timescale τ for the eddy heat flux averaged over the latitude range 25-75 • for ERA-Interim, MERRA and NCEP DOE.All calculations refer to the period 1995-2011.

Figure 2 .
Figure 2. (a) Correlation coefficient and slope of linear regression of zonal mean, deseasonalised temperature at 67 hPa and given latitude with deseasonalised equatorial zonal mean zonal wind shear (dU/dln(p)) at same level, calculated from ERA-Interim over the period 1994-2012.(b) Slope of linear regression shown in (a) integrated from equator, note a small residual remaining even when averaging globally (see text for discussion).Also shown is 35 • , taken in this study as a boundary for tropical average temperature (see text).

Figure 3 .
Figure 3. (a) RSS microwave sounding unit (MSU) lower stratospheric temperatures (MSU channel 4, or TLS), global (black), tropical (magenta) and combined extratropical (purple).(b, c) Tropical (35 • S-35 • N) temperature at 70 hPa, and estimates of dynamical temperature T p (based on eddy heat flux calculation; see text) and T mb (based on momentum balance calculation; see text); for ERA-Interim (panel b), and MERRA (panel c).(d) The residuals (tropical average temperature minus dynamical estimates) for ERA-Interim, and MERRA (colour coded as in panels b,c).Correlations between the estimates of dynamically forced temperature are given in Table1, between temperature and dynamical estimates in Table2, and between the residuals in Table3.All panels show deseasonalised data, with the mean annual cycle evaluated over the period 1995-2011.

Figure 4 .
Figure 4. Time series of annual mean dynamically forced temperature variations (T dyn ) from ERA-Interim and MERRA (colour coded, conventions as in Fig. 3) with results of linear regression excluding the volcanic periods (grey shaded).(a) T dyn based on the momentum balance calculation (i.e.T mb ); (b) T dyn based on the eddy heat flux calculation (i.e.T p ).

Figure 5 .
Figure 5. Overview of trends in temperature and in estimates of dynamically forced temperature based on the momentum balance calculations (T mb , panels a, b, e, f) and eddy heat flux calculations (T p , panels c, d, g, h) for ERA-Interim (panels a, c, e, g) and MERRA (panels b, d, f, h).Panels (a-d) show trends for the period 1980-2011 and panels (e-h) show trends for the period 1980-2011 excluding volcanic periods.Temperature trends are shown in black, dynamically forced temperature trends in brown.Residual trends are shown in blue for cooling, and in red for warming.

Table 5 .
Linear trends (Kelvin/decade) and 1σ statistical uncertainty for tropical average (35 • S-35 • N) temperature, and the dynamically forced temperatures estimated with the eddy heat flux proxy and the momentum balance calculation; for the pressure levels 67 hPa (ERA-Interim) and 70 hPa (MERRA).Results shown for the full period 1980-2011, and the period 1980-2011 excluding vol-

Figure 6 .
Figure 6.(a) Trends for the period 1980-2011 (excluding volcanic periods) in dynamically forced tropical temperatures at 70 hPa for both reanalyses (ERA-Interim and MERRA) and both methods (momentum balance calculation, and eddy heat flux proxy).Data as shown in Fig. 5e-h.(b) The fractional trends in monthly means of the forcing f (t) (Eq.2), i.e. w * for the momentum balance calculation and eddy heat flux for the eddy heat flux proxy.Trends in panel (a) are forced by trends in panel (b), but are phase shifted (due to e-folding with τ ∼ 80 days), and magnitude depends on the empirically determined proportionality coefficients α, K (Eqs. 4 and 5).Colour coding labelled in the figure, conventions as in Fig. 3.

Figure 7 .
Figure 7. (a) Time series of differences between ERA-Interim and MERRA in deseasonalised, monthly means of dynamical proxy for tropical temperatures at 70 hPa and 100 hPa (with ±1 offsets).The proxy time series have been normalised by their standard deviations over the period 1980-2011, with the sign convention as for tropical temperature (i.e. an increase in proxy corresponds to an increase in tropical temperature, and decrease in eddy heat flux).(b) Deseasonalised ratio of e-folded average eddy heat flux at 70 hPa and 100 hPa.The data in (a) and (b) have been smoothed with a running mean of 12 months for visual clarity.

Figure 8 .
Figure 8.(a) Latitudinal structure of ERA-Interim annual mean eddy heat flux trends at 67 hPa as function of month of year, El Chichón and Pinatubo periods excluded.(b) Ditto but for MERRA at 70 hPa.(c) Latitudinal structure of eddy heat flux trends from annual mean data, pressure levels as in panels (a) and (b).Grey shading shows the latitude ranges 75-25 • S and 25 • N-75 • N used to calculate the eddy heat flux temperature proxy.

Figure 9 .
Figure 9. Linear trends for the period 1980-2011 (excluding volcanic periods) in dynamically forced tropical temperature T p based on the ERA-Interim (blue, at 67 hPa) and MERRA (red, at 70 hPa) eddy heat flux calculation as function of equatorward bound of the latitude belt used to average the eddy heat flux, for start latitude 75 • .The error bars show the trend uncertainty in dynamically forced temperature due to the 1 σ uncertainty in the determination of the proportionality coefficient α in Eq. (5).

Table 1 .
Correlation between the dynamical temperature variations predicted by the eddy heat flux proxy (v T ) and the momentum balance calculation (m-bal) based on ERA-Interim and MERRA, for the periods1995-2011 (top) and 1980-2011 (bottom).

Table 2 .
Correlation between temperature and temperature predicted by the eddy heat flux proxy (v T ) and the momentum balance calculation (m-bal), for ERA-Interim and for MERRA: the periods 1995-2011 (no major volcanic eruptions) and the period 1980-2011 with El Chichón and Pinatubo periods excluded (see text).

Table 3 .
Correlation between residuals (as shown in Fig.3d) for the period 1995-2011.Duplicate information not shown.

Table 4 .
Difference in dynamically forced temperature T dyn (in Kelvin) of period after the year 2000 minus period before the year 2000, for T dyn estimated from the momentum balance (mb) calculation and the eddy heat flux proxy (v T ), based on ERA-Interim and MERRA.Also shown (last 2 columns) are the corresponding effective temperature differences ( T ) for ERA-Interim and MERRA.Differences are shown for period lengths of 2-6 years (i.e. the last row shows the difference of2001-2007 minus 1993-1999).