Latitudinal- and height-dependent long-term climatology of propagating quasi-16-day waves in the troposphere and stratosphere

The global amplitude of the westward propagating quasi-16-day waves (16DW) with wavenumber 1 (Q16W1), the strongest component of 16DW, are derived from the European Centre for Medium-Range Weather Forecasts ERA-Interim reanalysis temperature and zonal wind data sets from February 1979 to January 2018. In terms of temperature and zonal wind, strong climatologically average amplitudes of Q16W1 appear in the upper stratosphere at mid–high latitudes in both hemispheres, and the wave amplitude is stronger in the Northern Hemisphere (NH) than in the Southern Hemisphere (SH). Multivariate linear regression is separately applied to calculate the responses of the Q16W1 temperature and zonal wind amplitudes to the QBO (quasi-biennial oscillation), ENSO (El Niño-Southern Oscillation), solar activity and linear trends of the Q16W1 amplitude. The QBO signatures of the Q16W1 temperature and zonal wind amplitudes are mainly located in the stratosphere. The Q16W1 has significant QBO responses at low latitudes. In addition, only the temperature amplitude presents a larger QBO signature in its strongest climatological amplitude region. No significant responses to ENSO and solar activity are observed in temperature and zonal wind amplitudes. The linear trends of the monthly mean Q16W1 temperature and zonal wind amplitude are generally positive, especially in the mid-upper stratosphere. The trend is asymmetric about the equator and significantly stronger in the NH than in the SH. The seasonal variation in the trend of the temperature amplitude is studied and illustrated to be stronger in winter and weaker in spring and autumn. Further investigation suggests that the background and local instability trends contribute most of the increasing trend of the Q16W1 amplitude. In winter in both hemispheres, a weakening trend of eastward zonal wind provides more favourable background wind for Q16W1 upward propagation, in autumn and winter in the NH and in spring, autumn and winter in the SH, and the increasing trend of local instability may enhance wave excitation.


Introduction
Planetary waves (PWs), one of the main components of atmospheric waves, play a key role in transporting energy, momentum and chemical species among different atmospheric regions and are thus important in determining local and global atmospheric climatology and transient structures (Tsuda et al. 1994). PWs are global-scale oscillations predominantly generated by orography and diabetic heating caused by the distribution of land and sea in the troposphere (quasi-stationary PWs) or by irregular thermal or mechanical forcing in the lower atmosphere and/or by instabilities in the middle atmosphere (traveling PWs) with periods near 2, 5, 10 and 16 days (Andrews et al. 1987;Huang et al. 2013). Under certain conditions, PWs can propagate from the troposphere into the mesosphere and lower thermosphere (MLT), and wavenumber 1 and wavenumber 2 components are usually the predominant components (Charney and Drazin 1961).
The quasi-16-day waves (16DW) is one of the PWs and was identified as the second symmetric westward propagating Rossby mode with zonal wavenumber 1 (Salby 1981a, b). In the realistic atmosphere, due to Doppler shifting by the nonzero background flow, the period of 16DWs is from 12 to 20 days (Guharay et al. 2016). 16DWs have been extensively reported in recent decades from ground-based measurements (Mitchell et al. 1999;Das et al. 2010), satellite-borne measurements (McDonald et al. 2011;Alexander and Shepherd 2010) and reanalysis data sets (Vineeth et al. 2010). Combinations of ground-based and satellite-based analyses are also used to reveal the characteristics of 16DWs (Meek and Manson 2009). These findings show that 16DWs are prominent in the MLT from October to April. However, the exploration of 16DWs is relatively insufficient in the lower atmosphere, most likely due to the lack of highquality data sets, i.e., data with a long duration, good continuity, and high resolution. Since local observations cannot provide a global distribution, satellite observations with long-term duration are usually rare. Hence, global data with a long duration are necessary for further study.
Long-term variation is an important topic in atmospheric science. Previous observational and modeling studies have shown some long-term trends in atmospheric parameters, such as temperature (She et al. 2015) and wind fields (Kozubek et al. 2017). These studies revealed some important changes in the atmosphere from the lower atmosphere to the mesosphere. Hu and Tung (2002) determined that there were no obvious stratospheric wavenumber 1 and wavenumber 2 planetary wave (PW) activity trends in early and midwinter from November to January during 1950January during -2000 and 20 hPa along the 60°N latitude circle. However, a significantly negative PW activity trend from January to February 1979-2000 at 100 hPa in the Northern Hemisphere (NH) mid-high latitudes was revealed (Randel et al. 2002). Hu et al. (2019) suggested that the trend for stratospheric wave intensity from 200 to 10 hPa at NH mid-high latitudes strengthened during 1979-2000 and weakened during 2001-2015. However, most studies are related to stationary PWs, and research on traveling PW activity trends, especially Q16W1, is rare. Hence, the global trend in Q16W1 amplitude is far from fully understood.
The background wind can significantly affect the excitation, propagation and dissipation of PWs. On the other hand, PWs impact the background wind by depositing energy and momentum into the background atmosphere through various dissipation processes. Therefore, PWs might be related to the background wind at different time scales. Most previous studies on 16DW-mean flow interactions have focused on dynamic processes at the time scale of wave periods (Huang et al. 2013), seasonal scales (Huang et al. 2017) and intraseasonal variability (Espy et al. 1997;Day et al. 2011) rather than at the climate scale. Thus, to further understand the activity trend of 16DWs, we need more studies on 16DWs and their links with the background wind at the climate time scale.
For the purpose of investigating the global-scale and long-term characteristics of 16DWs in the stratosphere and below, ERA-Interim reanalysis temperature and zonal wind data sets were applied in our study. The responses to the quasi-biannual oscillation ( QBO ), EI Niño-Southern Oscillation ( ENSO ), 11-year solar cycle ( SC ) and linear trend of the strongest wave mode of 16DWs in the troposphere and stratosphere during February 1979-January 2018 were examined. In addition, we attempted to find a possible link between the wave and background wind and instability at the climate time scale. To this end, the rest of the paper is organized as follows. In the following section, we introduce the adopted data, dominant modes of the 16DW and calculation method. Subsequently, the global climatology of wave amplitude is presented. In Sects. "Response to the quasi-biannual oscillation" and "Long-term trend", we present the latitude-and height-dependent responses of the 16DW to the QBO and the linear trend of the strongest wave mode of 16DW, respectively. In the last section, we provide a brief summary of our analyses.

Data description and analytical approach
The global 6-hourly temperature and zonal wind data from the ERA-Interim data sets were downloaded from the European Centre for Medium-Range Weather Forecasts (ECMWF) for the period from February 1979 to January 2018 (39 years in total) at 37 pressure levels from 1000 to 1 hPa and were utilized in this study (Dee et al. 2011). ERA-Interim uses cycle 31r2 of the ECMWF'S Integrated Forecast System, which was introduced operationally in September 2006, with a reduced Gaussian grid with approximately 79 km spacing from the surface and other grid-point fields and a 60-level vertical resolution extending to 0.1 hPa. The vertical resolution of the product decreases with altitude, with a range of 0.1-1 km in the troposphere, which is reduced to 1-4.5 km in the stratosphere. For each day, the temperatures and zonal wind are interpolated at four UTCs with a 6-h interval: 0000 UTC, 0600 UTC, 1200 UTC and 1800 UTC. We chose the products on a grid of 72 × 143 points with 2.5° longitude and 2.5° latitude resolution.
At each pressure level, we first calculated the background temperature using a linear fitting within a 60-day sliding window with a 1-day shift interval. Then, the temperature fluctuation could be obtained by subtracting the background temperature and zonal mean temperature from the raw data. To demonstrate the dominant modes of 16DWs, we calculated the frequency-wavenumber spectra by performing a twodimensional fast Fourier transform on the temperature disturbance in each sliding 60-day window at each height (Huang et al. 2013;Gong et al. 2019). Furthermore, at each height, all spectra at different days were averaged to determine the temporally averaged spectra. Finally, a mean spectrum was obtained by averaging the temporally averaged spectra at all pressure levels. The calculation of the zonal wind spectrum is exactly the same as that of the temperature spectrum. The mean frequency-wavenumber spectra of temperature and zonal wind averaged from 38,475,745 spectra are shown in Fig. 1a, b, respectively, which illustrate that the most prominent spectral peak has a period of 15 days and a wavenumber of W1, and the largest amplitude is 0.46 K in temperature and 1.08 m/s in zonal wind. The secondary wave mode with a period of 20 days and a wavenumber of E1 can also be recognized in terms of temperature and zonal wind. Here, we focus only on the strongest PW, i.e., 16DW with wavenumber W1, which is named Q16W1 for simplicity.
To extract the monthly average amplitudes of Q16W1 and to determine the amplitude at the centre day of the window in each sliding 60-day window (Wu et al. 1995), we employed a harmonic fitting to temperature or zonal wind perturbations at each latitude bin and each pressure level, according to the following equation: where T ′ is the time series of the temperature or zonal wind perturbation; p j , t i , and i are the jth wave period, ith time, and longitude, respectively; and s = 1 is the wavenumber. In each window, p j varies from 12 to 20 days with an interval of 0.25 days, which accords with the 6-hourly data interval in ECMWF. B j and C j are the two coefficients to be fitted for the jth period, and then the wave amplitude T j can be specified by T j = B 2 j + C 2 j . In each 60-day window, we can specify 33 wave amplitudes corresponding to 33 periods p j in total, among which, the maximum amplitude and the corresponding period are regarded as the Q16W1 amplitude and period, at the centre day of the 60-day window, respectively. Finally, all Q16W1 amplitudes at different centre days within each month were averaged to derive the monthly mean amplitude of the Q16W1 amplitude, i.e., T A . The 1σ uncertainties of the sinusoidal fitted amplitude were employed to evaluate the derived Q16W1 amplitudes in the present study. The results demonstrated that, regardless of in temperature or zonal wind, 1σ uncertainties in amplitude for each fit were approximately two orders of magnitude smaller than fitted amplitudes, which indicated that the extracted global Q16W1 amplitude was reliable.
Multivariate linear regression (MLR) analysis was extensively utilized to isolate specific signals for zonalmean anomalies in temperature, zonal wind, ozone, and gravity wave energy from simulation and observation data (Gan et al. 2017;Weber et al. 2018). PWs could be modulated by QBO , ENSO and SC , which display interannual variation in PWs. There are also seasonal changes in PWs, which are mainly reflected in annual, semiannual, tri-annual and quarter-annual oscillations. Therefore, when performing MLR analysis on the amplitude of Q16W1, the influences of all these factors should be considered. Therefore, we chose a particular set of indices for the regression.
Before implementing MLR, inflection point analysis based on piecewise fittings of the Q16W1 amplitude was performed during the 39-year period to confirm that the linear trend fit over the entire time range was appropriate. In this study, MLR analysis was separately performed on the time series of the monthly mean temperature and zonal wind amplitudes of Q16W1, i.e., T A (t). The fitting equation is written as where t is time in months (478 months over 39 years), ω = 2π/12months , and B is the coefficient of the linear trend. The third to fourteenth terms on the right side of Eq. (2) are the linear correlations between the T A (t) and the SC , ENSO , two QBO components, annual oscillation ( AO ), semiannual oscillation ( SAO ), ter-annual oscillation ( TAO ) and quarter-annual oscillation ( FAO ), respectively, which are thought to be the dominant contributors to the variations in Q16W1. Among them, the third to sixth terms denote the interannual variations in Q16W1; and the seventh to fourteenth terms are the intraannual variations in Q16W1. G = G 2 1 + G 2 2 , H = H 2 1 + H 2 2 , I = I 2 1 + I 2 2 and J = J 2 1 + J 2 2 represent the annual, semiannual, ter-annual and quarterannual oscillation components, respectively. Here, the indices of F 10.7 cm ( unit : sfu,1sfu = 10 −22 Wm −2 Hz −1 ) and multivariate ENSO index (MEI) (downloaded from http:// www. esrl. noaa. gov/ psd/), shown in Fig. 2a, b, respectively, are used as proxies of the SC and ENSO activities. The time series QBO1 and QBO2 , shown in Fig. 2c are the two QBO components with a quarter-cycle phase difference (Gan et al. 2017). These two QBO time series correspond to the first two components of the empirical orthogonal functions (EOFs) extracted from
The delay autocorrelation coefficients of all independent variables, including SC , ENSO , QBO1 , QBO2 , AO , SAO , TAO and FAO were calculated, and all were close to zero. The delay cross-correlation coefficients between all variables were also obtained, except that the correlation between the SOLAR and ENSO index was 0.27, and all others were lower than 0.20. These small correlations indicated that they were all independent. In the MLR calculation, the confidence levels of the regression coefficients were estimated according to the variance-covariance matrix and Student's t test (Liu et al. 2017). Figure 3 shows the latitude-pressure distribution of the monthly mean Q16W1 amplitude in temperature and zonal wind fields averaged over 39 years. It is notable that the latitude-pressure distribution in both temperature and zonal wind exhibit a well-defined hemispheric symmetry. Moreover, both the temperature and zonal wind amplitude have maximum amplitudes in the upper stratosphere due to the decreasing atmospheric density decreasing with height. In addition, remarkable differences between the latitudinal variations in the temperature and zonal wind amplitudes can also be observed. The temperature wave amplitudes (Fig. 4a) in the NH are larger than those in the Southern Hemisphere (SH), with a maximum of 1.98 K (at 60°N, 2 hPa) in the NH and 1.32 K (at 45°S, 2 hPa) in the SH, respectively, which are quantitatively consistent with the results in Gong et al. (2019). Referring to Fig. 4b for zonal wind, there are 3 peak amplitudes at 1 hPa altitude that reach ~ 4.5 m/s, ~ 4 m/s, and ~ 3.5 m/s at ~ 85°N, ~ 45°N and ~ 15°S, respectively.

Global Q16W1 climatology
The month-pressure sections of the monthly zonalmean Q16W1 amplitudes in temperature and zonal wind, at six latitudes, representing high, middle, and low latitudes in the two hemispheres, are shown in Figs. 4 and 5, respectively. At 65°N/S, in both temperature and zonal wind, the AO of the Q16W1 amplitude are dominant in the entire stratosphere (100-1 hPa), with peak amplitudes appearing in winter months, and the amplitudes are larger in the NH than in the SH. A good illustration of this in regard to temperature is the maximum amplitude in the NH, which is 7.8 K in December 2000 at 1 hPa, which is much larger than the maximum of 4.3 K in June 2002 at 3 hPa in the SH, while in regard to zonal wind, the peak amplitude is ~ 12.5 m/s in the NH, which appears in December 2000 at ~ 1 hPa, and relatively weak amplitude of ~ 9.2 m/s in the SH is found in July 2002 at 2 hPa. At 35°N/S, significant AO of the Q16W1 amplitude in both temperature and zonal wind are concentrated only in a height range of 20-1 hPa, with peak amplitude also appearing in winter months. It should be noted that, different from the temperature amplitude, the zonal amplitude exhibits maximum amplitudes at middle latitudes in both hemispheres, one of which is ~ 16 m/s in January 2008 at 1 hPa in the NH, and the other is approximately 11 m/s in June 2016 at 1 hPa in the SH. At 10° N/S, except for the zonal wind in the SH, the amplitudes are smaller than those at middle and high latitudes, with a maximum peak smaller than 1.2 K in temperature and 12 m/s in zonal wind, and the annual variation of wave amplitudes in temperature still exists but is not as clear as those in zonal wind. In addition to the annual variation, the temperature amplitude also displays SAO above 5 hPa and QBO above 20 hPa.

Response to the quasi-biannual oscillation
The latitude-pressure patterns of the response of the monthly mean Q16W1 temperature amplitude to the QBO1 and QBO2 components, and the regression coefficients E and F calculated by Eq. (2) are shown in Fig. 6a, b, respectively. Regarding the zonal wind amplitude, its QBO1 and QBO2 signals are illustrated in Fig. 6c, d, respectively. It should be noted that QBO1 and QBO2 are not time series at the same height (Yamashita et al. 2018;Crooks et al. 2005). We confirmed that the QBO1 component can best describe the evolutions of the equatorial zonal wind at 10-30 hPa, and the major portion of the QBO1 is the zonal wind at 20 hPa. While the QBO2 component mainly represents a composite zonal wind at the heights of 50 hPa (major portion) and 70 hPa (minor portion). Therefore, the positive response to QBO1 or QBO2 signifies roughly in-phase variation between the Q16W1 and the zonal wind at the corresponding heights. As expected, at low latitudes, Q16W1 in temperature has a significantly positive response to QBO1 at ~ 100-30 hPa, and a positive and negative response located at ~ 70-30 hPa and ~ 30-7 hPa, respectively, and these responses are almost symmetrical about the equator. In addition, a positive response of zonal wind amplitude to QBO1 emerges at 5°N-25°S, 10-3 hPa, with a peak of ~ 0.3-0.4 m/s at ~ 10°S, 5 hPa. In addition to the low latitude responses, a strong response in temperature amplitude can also be clearly seen in the higher stratosphere at the middle and high latitudes. In the NH, strong positive and negative responses to QBO1 occur at 40°N-75°N, 10-1 hPa and 50°N-85°N, 100-30 hPa, with a positive maximum of 0.12 K at 60°N, 3 hPa and a negative maximum magnitude of ~ 0.06 K at 72.5°N, 50 hPa, respectively. In the SH, the response to QBO1 illustrates a different pattern. The positive value region is mainly concentrated in mid-latitudes, that is, in 25°S-55°S, 20-1 hPa, with a positive peak of ~ 0.11 K at 35°S, 7 hPa.
The response to the QBO2 illustrates a different pattern. At low latitudes, the responses of temperature amplitude to QBO2 have weak negative values at 100-50 hPa at low latitudes, which are nearly symmetrical about the equator. A similar symmetrical structure also exists in a height range of 70-30 hPa in the zonal wind. Furthermore, in zonal wind amplitude, a significantly negative response and positive response can be observed at 10°S-25°N above 7 hPa with peak amplitude The Same as Fig. 4 but for the Q16W1 zonal wind amplitude of ~ − 0.40 m/s at ~ 2.5°N, 3 hPa, and at the mid-latitude in the SH above 10 hPa with an amplitude between 0.1 and 0.2 m/s. Similar to the response to QBO1, the response of the temperature amplitude to QBO2 has the strongest value at higher latitudes and heights but with a negative peak of -0.16 K at 62.5°N at 2 hPa. Moreover, a stronger negative response of the temperature amplitude at high latitudes appears in the NH, exhibiting significant hemispheric asymmetry, which is evidently different from the response to QBO1. Our multiple regression analysis results do not show statistically significant responses of Q16W1 amplitude to the ENSO and solar cycle. Figure 7 diaplays the long-term trends in the monthly mean Q16W1 temperature and zonal wind amplitude over 39 years as a function of pressure and latitude. Regardless of the temperature or zonal wind, the spatial regions of the significantly positive trend regions in Fig. 7 are almost coincident with those of the strong Q16W1 climatologically average amplitude, as highlighted in Fig. 3, which indicates that the wave activities in the mid-upper stratosphere are generally strengthened. In addition, in temperature and zonal wind, the preponderant positive trends are asymmetric about the equator, and there is a superior positive trend in the NH. There are also distinctions between their distributions. In terms of temperature, a prominent positive trend appears at 20-1 hPa in the NH and 7-1 hPa in the SH, and the trend value increases with height from 30 to 3 hPa and then decreases with height in the NH. Moreover, with increasing latitude in both hemispheres, the positive trend increases poleward, reaches a maximum in the mid-high latitudes and then decreases with latitude. In the NH, the trend has a peak of 0.21 K/decade at 3 hPa around a latitude of 67.5°N. In the SH, a positive peak of 0.10 K/decade occurs at 3-2 hPa around a latitude of 60°S. In zonal wind, the overriding positive trend is located at 5-1 hPa in both hemispheres. In the NH, the trend has three peaks, 0.25-0.30 m/s/decade, 0.45-0.50 m/s/decade, and 0.30-0.35 m/s/decade, which are Fig. 6 Latitude-pressure sections of the responses of the monthly-mean Q16W1 temperature amplitude to (a) QBO1 and (b) QBO2. Similar to the temperature amplitude, (c) and (d) correspond to the responses of the zonal wind amplitude to QBO1 and QBO2, respectively. Only the results with a confidence level at/above 95% are plotted as contours. The solid and dashed contours denote the positive and negative responses, respectively located at ~ 15°N, ~ 37.5°N, and 85°N in the uppermost height, respectively. In the SH, there is one trend peak, 0.25-0.30 m/s/decade, which emerges at ~ 10°S in the uppermost height. Furthermore, in temperature, a negative trend occurs below 300 hPa, notably at high latitudes in the NH, with a peak negative trend of -0.08 K/decade at 950 hPa around a latitude of 75°N, while in zonal wind, a negative trend with a peak value of ~ 0.08 m/s/decade can be clearly discerned closely below 30 hPa, particularly in the low and mid-latitudes in the NH. Some previous studies revealed significant seasonal variation in Q16W1 (Williams and Avery 1992;Luo et al. 2002b;Hibbins et al. 2009;Day and Mitchell 2010). Here, the Q16W1 temperature amplitude trends in four seasons are investigated. We take December, January, and February to represent the NH winter and SH summer, and March to May, June to August, and September to November to represent the NH (SH) spring (autumn), summer (winter), and autumn (spring), respectively. In each season of each year, there are 3 monthly mean temperature amplitudes of Q16W1. Thus, in the calculation (2) to calculate the trend in each season. Then, in each season, all monthdependent trends are averaged to specify the trend in this season. From Fig. 8, it can be seen that the positive trends are in 40°N-85°N, 20-1 hPa and 67.5°S-82.5°S, 5-2 hPa in spring; 50°N-82.5°N, 3-1 hPa and 20°S-65°S, 5-1 hPa in autumn; and 30°N-87.5°N, 30-1 hPa and 20°S-85°S, 7-1 hPa in winter. To better display the relation between the Q16W1 amplitude trend and the climatological mean Q16W1 temperature amplitude, we present the climatological mean of seasonal Q16W1 temperature amplitude in Fig. 8. The strong enhancement trend is generally consistent with the latitude, height, and seasonal changes of the climatologically average peak amplitude and generally appears in the mid-high latitudes and high stratosphere in winter. The strongest trend appears in the NH winter (Fig. 8d) with a peak of 0.54 K/decade at 3 hPa around a latitude of 67.5°N, which is significantly larger than the trend in the SH winter (Fig. 8b), in which a positive peak of 0.19 K/decade occurs at 3 hPa around a latitude of 50°S, indicating obvious hemispheric asymmetry in the Q16W1 amplitude trend.

Long-term trend
From the above analyses, it can be seen that the Q16W1 climatologically average temperature amplitudes and their trends appear in middle and high latitudes. Dickinson (1968) proposed a polar waveguide theory, which predicted PW propagation from the troposphere to the stratosphere at high latitudes. Mastuno (1970) obtained similar results using a quasi-geostrophic model. Luo et al. (2002a, b) confirmed that westward-traveling 16-day waves with a small phase velocity can propagate upward through the winter polar night jet to reach the MLT region only in an eastward background flow of moderate speed which is present in the winter hemisphere. The zonal wind was confirmed by Smith (1997) to be a key factor for the vertical propagation and dissipation of PWs. For simplification, a PW can propagate through a region only when the PW zonal phase velocity c and background zonal wind u o satisfy the following conditions (McDonald et al. 2011): ( Fig. 8 Contours represent the seasonal variation of the long-term trend (in K per decade) of the monthly-mean Q16W1 temperature amplitude from the 39 year ERA-Interim temperature data set. The solid and dashed contours denote the positive and negative trends, respectively. Only the trends with a confidence level at/above 95% are plotted. The colors represent the climatological mean of seasonal Q16W1 temperature amplitude where k = 2π/ x and 2π/ y are the zonal and meridional wavenumbers of Q16W1, respectively, in which x and y are the zonal and meridional wavelengths of Q16W1, respectively. f o and β represent the Coriolis parameter and Rossby parameter, respectively. The zonal phase velocity of Q16W1 c = ω/k , where ω = 2π/T is the frequency of Q16W1. The positive (negative) values of u o and c indicate the eastward (westward) direction, and U c is the critical Rossby velocity determined by the PW. Equation (3) manifests that weak eastward zonal wind is favorable for Q16W1 propagation; otherwise, energy is trapped/reflected in regions, where the zonal winds are westward or strongly eastward. Therefore, the zonal mean wind plays an important role in PW propagation, and its long-term variation may be related to the longterm trend of Q16W1 amplitude. Figure 9 displays the trends and climatological mean of seasonal zonal wind. Combined with Fig. 8, it is revealed that the amplitude of Q16W1 is usually strong in regions, where the background zonal wind is weakly eastward (Day et al. 2011). The zonal wind trend has significantly positive values in a height range from 1000 to 7 hPa in mid-high latitudes during SH summer. More interestingly, the zonal wind trend has significantly negative values in a height range from 1000 to 50 hPa at high latitudes in the NH winter, with a negative peak of − 1.4 m/s/decade at 50 hPa around a latitude of 67.5°N. Another negative trend of the zonal wind appears at 30°S-50°S, at 30-20 hPa. It should be noted that these negative trend regions correspond to the climatological eastward wind region, implying the gradual weakening of the eastward wind u o . This weakened eastward wind leads to a more favorable background for the upward propagation of Q16W1, which in turn may cause an increasing trend of wave amplitude. Obviously, the negative zonal wind trend could contribute to the positive trend of the Q16W1 amplitude.
It is known that the mean flow barotropic and baroclinic instabilities are important in local excitation for PWs (Hartmann 1983;Huang et al. 2021;William and Leslie 1991). Therefore, we will investigate the possible influences of the local excitation mechanism on the trend of the Q16W1 amplitude by analyzing the mean flow instability. A necessary condition for instability of The solid and dashed contours denote the positive and negative trends, respectively. The thick black contour represents the 0 value. The stippled regions represent the trends at/above the 95% confidence level. The colors represent the climatological mean of seasonal zonal wind the mean flow is that the quasi-geostrophic potential vorticity gradient ( q ) must change sign somewhere in the flow domain (Andrews et al. 1987). Thus, the overturning of q can be regarded as local instability, which indicates that relatively small q implies a high probability of the overturning of q , as well as local instability (Lu et al. 2020). q can be expressed as where overbar, prime and subscript denote zonal average and derivative; u , ρ o , a , , ∅ , and f represent the zonal wind, background density, Earth radius, angular velocity of the Earth, latitude and Coriolis parameter, respectively; z is the altitude above the Earth's surface; and N is the buoyancy frequency, which is specified by calculating N 2 = ḡ T ( ∂T ∂z + g c p ) , where g , T , and c p are the heightdependent gravitational acceleration, background temperature, and specific heat at constant pressure, respectively. (4) Considering the validity limitations of the quasi-geostrophic approximation near the equator and surface, and numerical singularities near the poles, the climatological distributions and negative trends of monthly mean q values in two latitude zones, 15°S-85°S and 15°N-85°N, and altitudes above 500 hPa, are shown in Fig. 10. In the NH, it is intriguing to note that in all seasons, the monthly mean q values are mostly positive in almost all regions except at approximately 300 hPa in the low and middle latitudes and in the high latitudes of summer. Large values of q are found in the troposphere. We concentrate on the regions with relatively small q and negative q trends, where an increasing instability trend is most likely to occur. These regions are located at 70°N-80°N from 5 to 2 hPa in summer, 75°N-85°N from 200 to 100 hPa in autumn and 55°N-75°N from 250 to 70 hPa. Combined with the wave trend in the NH in Fig. 8 and the probability of upward wave propagation, the increased possibility of instability could partly explain the Q16W1 amplitude trends in autumn and winter.
Similar to that in the NH, the climatological mean of seasonal Q16W1 in the SH in Fig. 10 is also mainly positive; interestingly, it has significantly large positive value Fig. 10 Trend (contours, units: 10 -5 m -1 s -1 /decade) of the monthly mean q in four seasons. The dashed contours denote the negative trend. The stippled regions represent the trends at/above 95% confidence level. The colors represent the climatological mean of seasonal q in the middle and high latitudes of the stratosphere in both spring and winter and a negative value at high latitudes in all seasons except autumn. The regions are located at 57.5°S-62.5°S at 2 hPa in spring, high latitudes from 200 to 100 hPa in autumn and 30°S-40°S at 50-20 hPa in winter, where q represents negative trends, indicating an enhancement of instability. In general, the increasing trend of instability located in 50°N-82.5°N from 3 to 1 hPa in autumn and 30°N-87.5°N from 30 to 1 hPa in winter in the NH and in 70°S-80°S, from 3 to 2 hPa in spring, 20°S-62.5°S from 5 to 1 hPa in autumn and 20°S-85°S from 7 to 1 hPa in winter in the SH, may contribute to the increasing trend of Q16W1 amplitude.

Summary
Q16W1 was demonstrated to be the dominant mode of 16DW globally from the surface to the stratopause by analyzing the ERA-Interim reanalysis data sets from February 1979 to January 2018. The climatological characteristics of the Q16W1 temperature and zonal wind amplitudes were presented. The responses of the monthly mean Q16W1 temperature and zonal wind amplitudes to QBO, ENSO and solar activity and their linear trends were investigated by multiple regression analysis. Possible mechanisms of seasonal variation in the Q16W1 temperature amplitude trend were discussed. The primary results are summarized as follows.
For the 39 year climatologically average temperature and zonal wind amplitudes, a strong Q16W1 occurs in the upper stratosphere at mid-high latitudes in both hemispheres, and the wave amplitude is stronger in the NH than in the SH. There are differences in the latitudinal variations of the amplitude between zonal wind and temperature, mainly reflected in the distribution of the peak amplitude. Q16W1 presents obvious seasonal variation, that is, strong in winter and weak in summer. The obvious AO in the entire stratosphere (100-1 hPa) in the high latitudes and at 20-1 hPa in the midlatitudes can also be observed in both temperature and zonal wind. However, at low latitudes, the zonal wind amplitude has a more significant AO signature than the temperature amplitude, followed by all latitudes in the two hemispheres. Furthermore, stronger peak amplitudes in zonal wind exist in the upper stratosphere at middle latitudes, while stronger peak amplitudes in temperature are found at high latitudes.
In temperature and zonal wind, the responses of the Q16W1 amplitude to QBO1 and QBO2 are mainly concentrated in the stratosphere and exhibit distinct signatures at low latitudes. In the temperature amplitude, two strong positive responses to QBO1 appear at mid-high latitudes in the NH and middle latitudes in the SH, while an evidently negative response to QBO2 emerges in the high latitudes in the NH. The zonal wind amplitude demonstrates a prominent positive response for QBO2 in the middle latitudes in the SH. No obvious responses to solar activity and ENSO are found in temperature and zonal wind amplitudes.
In the temperature and zonal wind, the long-term trends of the monthly mean Q16W1 amplitude are generally positive and mainly concentrated in the stratosphere. Whether in temperature or in zonal wind, the spatial regions of the significantly positive trend regions are almost coincident with those of the strong climatological amplitude, the preponderant positive trend is asymmetric about the equator, and there is a higher positive trend in the NH. There are also distinctions in the trends between temperature and zonal wind amplitude. These are mainly reflected in the areas and locations, where obvious positive trends and their peak amplitudes appear and the distribution of negative trends, which are described in detail in the previous sections.
The trend of the monthly mean Q16W1 temperature amplitude has evident seasonal variation; that is, stronger in winter and weaker in spring and autumn in the two hemispheres. To investigate the possible causes of the Q16W1 trend, we calculated the long-term trend of the monthly mean zonal wind and q in each of the four seasons. In winter in both hemispheres, the weakening trend of eastward zonal wind contributes to the increasing trend of the Q16W1 amplitude, because the weak eastward zonal wind is favorable for Q16W1 propagation. Moreover, in autumn and winter in the NH and in spring, autumn and winter in the SH, the increasing trend of instability could enhance the wave excitation, and thus lead to the increasing trend of the Q16W1 amplitude.