The detection of climate change signals in rather short satellite datasets is a challenging task in climate research and requires high-quality data with good error characterization. Global Navigation Satellite System (GNSS) radio occultation (RO) provides a novel record of high-quality measurements of atmospheric parameters of the upper-troposphere–lower-stratosphere (UTLS) region. Because of characteristics such as long-term stability, self calibration, and a very good height resolution, RO data are well suited to investigate atmospheric climate change. This study describes the signals of ENSO and the quasi-biennial oscillation (QBO) in the data and investigates whether the data already show evidence of a forced climate change signal, using an optimal-fingerprint technique. RO refractivity, geopotential height, and temperature within two trend periods (1995–2010 intermittently and 2001–10 continuously) are investigated. The data show that an emerging climate change signal consistent with the projections of three global climate models from the Coupled Model Intercomparison Project cycle 3 (CMIP3) archive is detected for geopotential height of pressure levels at a 90% confidence level both for the intermittent and continuous period, for the latter so far in a broad 50°S–50°N band only. Such UTLS geopotential height changes reflect an overall tropospheric warming. 90% confidence is not achieved for the temperature record when only large-scale aspects of the pattern are resolved. When resolving smaller-scale aspects, RO temperature trends appear stronger than GCM-projected trends, the difference stemming mainly from the tropical lower stratosphere, allowing for climate change detection at a 95% confidence level. Overall, an emerging trend signal is thus detected in the RO climate record, which is expected to increase further in significance as the record grows over the coming years. Small natural changes during the period suggest that the detected change is mainly caused by anthropogenic influence on climate.
1. Introduction and background
With the Global Positioning System/Meteorology (GPS/Met) experiment on board the MicroLab-1 satellite, launched in 1995 (Ware et al. 1996), a new active limb-sounding technique to scan the earth’s atmosphere was put into practice. The Global Navigation Satellite System (GNSS) radio occultation (RO) technique allows the retrieval of information in the upper-troposphere–lower-stratosphere (UTLS) region. This part of the atmosphere is governed by a complex interaction between dynamics, transport, radiation, chemistry, and microphysics (Mohanakumar 2008) and reacts particularly sensitive to climate change. This study assesses the usability of RO data in climate change science. It investigates the response to El Niño–Southern Oscillation (ENSO) and the quasi-biennial oscillation (QBO) in data and whether a climate change signal is emerging in the RO record.
Most upper-air datasets, such as radiosondes or microwave satellite data, were not primarily intended for climate monitoring. They are affected by uncertainties, for example, in cross calibration between satellite data, inhomogeneities due to introduction of improved instrumentation, and other systematic changes, all of which yield long-term trends uncertainties. Requirements for data suitable for climate monitoring are outlined by, for example, Ohring et al. (2007), for a discussion of drawbacks see, for example, Karl et al. (2006) or Randel et al. (2009). The RO method shows promise to overcome these deficiencies. It provides vertically well-resolved profiles of bending angle, refractivity (N), geopotential height (Z), pressure (p), and temperature (T) for the UTLS, which we define as the region between 5 and 35 km. The RO method uses GNSS radio signals (practically so far from the U.S. GPS), which are received by a low earth orbit satellite (LEO). Scanning through an atmospheric limb results from the relative motion between the transmitter and receiver satellite placed on opposite sides of the limb. An occultation event occurs when a GNSS satellite sets down or rises from behind the Earth as seen from the LEO satellite. Detailed descriptions of the method are given by, for example, Kursinski et al. (1997), Steiner et al. (2001), and Hajj et al. (2002).
Very briefly summarized, the measured core quantity is the signal phase at two GPS frequencies, which is proportional to the distance between the transmitter and the receiver. These phase measurements are based on precise atomic clocks, which enable long-term stability and International System of Units (S.I.) traceability of RO measurements (Leroy et al. 2006b). The so-called excess phase is the part of the signal phase caused by the atmospheric refractivity field, adding to the vacuum phase path. The latter can be derived from precise orbit position and velocity information of the GPS and LEO satellites and is subtracted. Based on the excess phase and orbit arc data, acquired as a function of time during the occultation event, atmospheric variables as a function of altitude are then retrieved that are found closely consistent among different satellites and over multiple years.
Remarkable consistency between data from the two satellites Challenging Minisatellite Payload (CHAMP) and Satélite de Aplicaciones Científicas-C (SAC-C) was shown by Hajj et al. (2004). Similar excellent consistency results were gained more recently by Foelsche et al. (2008b, 2009) for different Taiwan/United States Formosa Satellite Mission 3/Constellation Observing System for Meteorology, Ionosphere and Climate (hereafter F3C) based RO products, where seasonal RO temperature climatologies of different satellites were found to be in agreement to within <0.1 K in the UTLS. Also, for multisatellite monthly climatologies at zonal resolutions as used in this study the deviation of temperatures from the satellite mean was found <0.1 K between 8 and 25 km, in line with Pirscher (2010). For global means, consistency is found even within 0.05 K (Foelsche et al. 2011). Furthermore, for the two months with GPS/Met data used in this study (October 1995 and February 1997) it has been shown that the error characteristics are not worse than those of the more recent F3C data (Pirscher 2010). Thus, RO data from different satellites or sensors can be combined without intercalibration or temporal overlapping as long as the same processing scheme is employed and data quality is tracked (Foelsche et al. 2008b; Pirscher 2010). Nearly all-weather capability is achieved by the use of radio signals, which penetrate the atmosphere in nearly all conditions and are almost insensitive to clouds and aerosols. Near-polar orbits ensure global coverage of measurements with already a single LEO satellite. Regarding typical numbers of occultation events, the CHAMP satellite with one receiving antenna provided about 250 events per day while the six-satellite F3C constellation provided up to 2500 events per day.
The utility of RO data for climate monitoring has been analyzed in various studies, using 1) simulated data or proxies such as climate model data (e.g., Steiner et al. 2001; Foelsche et al. 2008c; Ringer and Healy 2008; Lackner et al. 2009; Lackner 2010); 2) real RO data of various satellite missions (e.g., Schmidt et al. 2004; Leroy et al. 2006b; Foelsche et al. 2009; Steiner et al. 2009; Schmidt et al. 2010); or 3) comparing RO data to other satellite datasets (e.g., Gobiet et al. 2007; Steiner et al. 2007), reanalyses or analyses data from Numerical Weather Prediction (Borsche et al. 2007; Foelsche et al. 2008a). Structural uncertainty, defined as unintentional bias arising from the chosen methodological approaches (Thorne et al. 2005a), was analyzed for RO data from four different processing centers by Ho et al. (2009). Low structural uncertainty estimates together with full error characterizations (e.g., Steiner and Kirchengast 2005; Pirscher et al. 2007; Foelsche et al. 2008a) underline the value of RO data for climate applications.
2. Setup of detection study
The optimal-fingerprinting technique for signal detection relies on data from observations and climate model simulations. Beside those datasets, reanalyses and radiosonde records are employed for trend comparisons and for better understanding of the RO data response to the ENSO and the QBO. Linear regression is applied to estimate the variability of these large-scale atmospheric patterns governing variability in the ULTS. The first part of this section presents the spatiotemporal study setup, the data, and their characteristics. The second part gives an overview of the methods used.
The study is based on single refractivity, geopotential height, and temperature (all as functions of pressure) trend patterns of UTLS RO multisatellite climatologies. The anticipated climate change signal is estimated from an ensemble of anthropogenically forced scenario runs of three global climate models (GCM) of the Fourth Assessment Report (AR4) of the Intergovernmental Panel on Climate Change (IPCC; Meehl et al. 2007). Natural variability is based on preindustrial control simulations of the same models.
The focus region of our study is the UTLS within 50°S–50°N, where the best RO data quality is provided (Steiner et al. 2009). A horizontal resolution based on five zonal regions was employed, consisting of a tropical band between 10°S and 10°N and 4 bands at 30°–10°S/10°–30°N and 50°–30°S/30°–50°N. The data are used on 8 pressure levels spanning from 300 hPa (approximately 9 km) to 30 hPa (approximately 25 km). For all investigations, we consider monthly mean time series for the 1995 to 2010 period from which the average seasonal cycle over the analysis period was removed. Table 1 gives an overview on the datasets and their resolutions.
1) RO data
The setup of the analysis is closely tied to the availability of RO measurements. First RO measurements of the earth’s atmosphere are available from the GPS/Met experiment for several periods from 1995 to 1997. We only use the best quality data given for October 1995 and February 1997 (Schreiner et al. 1998), when sufficient data were gained to calculate monthly means. A continuous multiyear RO record became available with CHAMP (e.g., Wickert et al. 2004), which was launched in 2001 and delivered measurements until 2008. Because of the sparse data caused by satellite instrument problems in July and August 2006, these two months are excluded from the analysis. Besides GPS/Met and CHAMP, data from the following satellites are used for the RO record until mid-2010: SAC-C (e.g., Hajj et al. 2004), Gravity Recovery and Climate Experiment (GRACE-A) (e.g., Beyerle et al. 2005; Wickert et al. 2005), and the six satellites of the F3C mission (e.g., Anthes et al. 2008). We employ monthly-mean zonal-mean multisatellite climatologies within 1995 to 2010 (see Table 1 for details). Together with the two GPS/Met months, a total of 107 monthly-mean climatologies spanning across about 15 years are at hand. Two time periods are considered for the detection study: 1) the intermittent period 10/1995–07/2010 including GPS/Met data, and 2) the continuous period 9/2001–7/2010 (only gap 7–8/2006).
The generation of profiles of atmospheric RO parameters and climatologies is based on the Wegener Center occultation processing system OPSv5.4 (Gobiet et al. 2007; Kirchengast et al. 2007; Fritzer et al. 2009). The data are available online at www.globclim.org (December 2010). A short overview on OPSv5.4 is given by Steiner et al. (2009), and a detailed description can be found in Pirscher (2010). Briefly, RO phase and orbit data from the University Corporation for Atmospheric Research (UCAR; Boulder, Colorado) are used to derive bending angle profiles, which in turn lead to refractivity. Vertical integration of refractivity, which is closely equivalent to air density, yields pressure. Geopotential height is computed by integrating the latitude- and height-dependent acceleration of gravity over the RO-derived altitude (divided by the standard acceleration of gravity) and interpolating it to pressure levels. The equation of state is used to compute dry temperature. The relation between refractivity (N), pressure (p), temperature (T), and humidity is given for each measurement point by the Smith–Weintraub formula (Smith and Weintraub 1953):
Here, n is the refractive index, pw the partial pressure of water vapor, and k1 = 77.6 K hPa−1 and k2 = 3.73 × 105 K2 hPa−1 are empirically determined constants. Above the midtroposphere, atmospheric humidity is small. Thus, the second right-hand-side term in Eq. (1) can be neglected, meaning that a “dry” atmosphere is assumed (valid for the UTLS from about 300 hPa upward). If water vapor is available, dry temperature is always lower than physical temperature (Foelsche et al. 2008c). As atmospheric water vapor is very likely increasing when the climate warms (Held and Soden 2006), dry temperature trends associated with warming cannot exceed physical temperature trends, unless there would indeed be a decrease in absolute humidity.
RO climatologies derived from binning and averaging the RO profiles (Foelsche et al. 2008a, 2009; Pirscher 2010) feature well-defined error characteristics. The overall total error for single-satellite temperature climatologies of 10° zonal means below 30 km was estimated as <0.3–0.5 K, where the greatest contribution is due to uneven sampling (Pirscher et al. 2007; Foelsche et al. 2008a) and getting smaller when more satellites are used (Foelsche et al. 2009). The sampling error (SE) can be estimated via reference atmosphere fields, such as the European Centre for Medium-Range Weather Forecasts (ECMWF) analysis fields, which can be assumed to statistically represent the “true” atmospheric variability. ECMWF profiles, collocated in space and time to RO profiles, are computed and averaged to monthly mean zonal bins. The SE is then estimated as the difference between the monthly mean field from these collocated profiles and the monthly mean from the full ECMWF analysis fields (Pirscher et al. 2007; Foelsche et al. 2008a; Pirscher 2010). It is subtracted from the RO climatologies to mitigate the influence from uneven sampling, further reducing the error of the climatologies (as discussed in Foelsche et al. 2008a; Steiner et al. 2009; Foelsche et al. 2011). For processing scheme induced structural uncertainty an upper bound of 0.03% (5 yr)−1 was estimated for global refractivity trend monitoring, which corresponds to about 0.06 K (5 yr)−1 for temperature trends (Ho et al. 2009).
2) Model data
To estimate the expected climate signal and the natural climate variability for optimal detection a multimodel, multirealization dataset is used. It is based on three GCMs of the IPCC AR4 archive: 1) the Community Climate System Model, version 3 (CCSM3) (Collins et al. 2006) from the National Centers for Environmental Prediction (NCEP) and the National Center for Atmospheric Research (NCAR), United States; 2) ECHAM5 (Roeckner et al. 2003) from the Max Planck Institute for Meteorology (MPI-M), Hamburg, Germany; and 3) the third climate configuration of the Met Office Unified Model (HadCM3) (Gordon et al. 2000; Pope et al. 2000) from the Met Office Hadley Centre (UKMO), United Kingdom. The model output is available via the World Climate Research Programme’s (WCRP) Working Group on Coupled Modeling [Coupled Model Intercomparison Project cycle 3 (CMIP3)] multimodel database, along with model information (available online at http://www-pcmdi.llnl.gov/ipcc/model_documentation/ipcc_model_documentation.php, December 2010). The three GCMs perform without flux corrections and include stratospheric ozone depletion and recovery forcings (Roeckner et al. 2005; J. Meehl et al. 2007, personal communication), which is particularly important for an adequate simulation of stratospheric temperatures (Forster et al. 2007). According to an analysis of Reichler and Kim (2008), the models belong to the five best-performing models without flux correction, meaning that they show good agreement with observations in their time-mean state of the climate.
The expected anthropogenically forced climate change signal (hereinafter fGCM) is estimated from climate simulations of the twentieth century (experiments run with greenhouse gas increases as observed in the twentieth century), which are concatenated with Special Report on Emission Scenarios (SRES; Nakićenović et al. 2000) A2 and B1 simulations beyond 1999. For these simulations, direct (all GCMs) as well as indirect (HadCM3) or semidirect (CCSM3, ECHAM5) sulfate aerosol forcings are used in addition to greenhouse gas and ozone forcings. The period to calculate the trend signals is chosen consistently to the RO data. Natural variability is assessed with preindustrial control simulations (hereinafter PICTRL) of the models. Details about the available number of forced simulations and used PICTRL years are given in Table 1.
Climate model and reanalyses data (see section 3) were first regridded to a common 2.5° × 2.5° grid in latitude and longitude, then area-weighted zonal means were calculated. The model and reanalysis data are used after interpolating them to the same spatial resolution as the RO data and are limited to the same spatial and temporal coverage. Geopotential height, temperature, and specific humidity (all provided by the data centers as functions of pressure) are used to calculate refractivity based on Eq. (1) for GCMs and the reanalyses.
3) Reanalyses and radiosonde data
For the analysis of large-scale patterns, such as ENSO or the QBO, 40-yr ECMWF Reanalysis (ERA-40) output (Simmons and Gibson 2000; Uppala et al. 2005) is used. The reanalysis was chosen because it provides a record with continuous spatial and temporal coverage. Furthermore, ERA-40 was found to match best the in situ observations in the climate mean state (Reichler and Kim 2008). To complete the tropospheric picture of the patterns, we employ additionally all available pressure height levels from 300 hPa to the surface. The ECMWF ERA-Interim reanalysis (Simmons et al. 2007a,b; Uppala et al. 2008; Poli et al. 2010; Simmons et al. 2010) and the radiosonde temperature records Radiosonde Observation using Reanalysis (RAOBCORE, version RAOBCOREv1.4; Haimberger 2007), Radiosonde Innovation Composite Homogenization (RICH; Haimberger et al. 2008) from the University of Vienna, and the Hadley Centre radiosonde temperature dataset (HadAT; Thorne et al. 2005b) were used to compare with RO results, especially for assessing the consistency of trends based on the full temporal coverage compared to those based on the incomplete RO record between 1995 and 2001.
The detection study is based on an ordinary least squares fingerprinting approach following Hegerl et al. (1996). Additional analysis of large-scale patterns and their share in the total variability are implemented with multiple linear regressions.
1) Optimal fingerprinting
Optimal fingerprinting (see, e.g., Hegerl et al. 1997; Barnett et al. 2005; Hegerl et al. 2007) can be considered as generalized multivariate regression. An observed record y is represented by a scaled (β) estimated externally forced climate change pattern (guess pattern, X), which contains the response to anthropogenic forcings and an estimate for internal climate variability (ε):
As observed climate change signals y we use the spatial trend patterns of RO parameters, representing latitude–height slabs.
The guess pattern is presented by the average of the 20-member fGCM ensemble. Averaging over 20 ensemble members will reduce the influence of internal climate variability on the guess pattern by a factor of more than 4. While this noise may still lead to a slight underestimate of the forced signal, it will not lead to spurious detection (Allen and Stott 2003; see discussion in Hegerl et al. 2007). Therefore, the use of an ordinary least squares fit is appropriate here.
The PICTRL data are on the one hand used to represent the data (RO, fGCM) in a dimension-reduced EOF space (where the covariance matrices of the noise are diagonal, truncation levels are discussed below) and to estimate the scaling factors , which modulate the amplitudes of the guess patterns. On the other hand, a second statistically independent sample of PICTRL data (referred to as control) is employed to access the uncertainty in the detection variables. To obtain the two samples, the PICTRL data of each GCM were divided into halves. The shortest PICTRL simulation (HadCM3) provides 341 yr. To avoid disproportional weighting, 170 yr of each of the three models were used to compose one sample (of 510 yr) and likewise for the other sample. Based on the two samples with 510 yr of data each, the PICTRL trend pattern matrices equivalent to the used RO periods were then generated. The individual trend patterns were calculated without temporal data overlap.
The eigenvectors (pooled in F) and eigenvalues (pooled in the diagonal matrix Λ) of the decomposed covariance matrix of the first PICTRL trend sample were used to represent the trend patterns in the dimension-reduced EOF space (bdata) and to estimate in the following the respective scaling factors :
The subscript “data” stands for RO, fGCM, and control trend patterns for their presentation in the EOF space in Eq. (3) and for RO and control concerning the estimation of the scaling factors in Eq. (4). The null hypothesis that the scaling factors of the observations are zero is tested by means of the distribution of the control scaling factors.
To assess whether the climate variability of the PICTRL simulations adequately represents the variability of the observations in the truncated space, a residual consistency test following Allen and Tett (1999) was performed. The aim is to have no explicit reason to distrust the uncertainty estimates in the analysis. This is true if the regression residuals (r) are consistent with the noise estimate from the control simulations. If the noise in the truncated EOF space is normally distributed, its variance will be Chi-square (χ2) distributed,
where k is the number of retained EOFs. If the residuals are outside 5%–95% confidence bands, there is evidence that the optimization has focused on poorly sampled or poorly simulated features of the response, yielding unreliable results (Allen and Tett 1999).
2) Influence from modes of variability on RO data
The UTLS variability is mainly influenced by two atmospheric patterns, the QBO and ENSO. Originating in the tropics, they also impact higher latitudes’ conditions. A multiple linear regression model was employed 1) to estimate the share of the patterns in UTLS variability for RO data and 2) to remove the ENSO signal (from GCMs and RO) and the QBO signal (from RO). It is given for each gridpoint by
The regression coefficients are the constant term a0, the trend coefficient a1, the QBO coefficient a2, and the ENSO coefficient a3. Here, ε denotes the residual error term. For the QBO the monthly mean 50-hPa zonal wind index (available online at www.cpc.ncep.noaa.gov/data/indices/qbo.u50.index, November 2010) of the Climate Prediction Center (CPC)/National Oceanic and Atmospheric Administration (NOAA) was employed. For ENSO we used for the observational records the seasonally smoothed monthly N3.4 index (available online at www.esrl.noaa.gov/psd/forecasts/sstlim/global/indices_global, November 2010), made available from the Physical Sciences Division of the NOAA. For the GCMs we estimated the N3.4 index as standardized 5-month moving average of the 1000-hPa temperature mean in the N3.4 region (5°S–5°N, 170°–120°W) of the respective trend period (tests confirmed a correlation of at least 98% between 1000 hPa and sea surface temperature fields, the latter were not available for all models). We identified and considered a four-month atmospheric lag for ENSO, which is consistent with Seidel et al. (2004) and Steiner et al. (2009). The coefficients of determination (R2) were used for an estimate of the model-explained variability due to the individual patterns.
The IPCC AR4 models do not produce meaningful QBO variability. Figure 1 shows tropical time series of the RO parameters N, Z, and T from observations (left) and the fGCMs (right). While four QBO periods are clearly visible in the RO records (best in temperature above 100 hPa), the signal is not only lacking in the multimodel mean fGCM (right) but also in the variability in individual GCM simulations (not shown). To get a cleaner climate signal and comparable internal variability for the datasets we removed the ENSO signal from the GCMs and the ENSO and QBO signal from the RO data by applying Eq. (6).
In this section, we first briefly discuss the solar influence and then in particular the influences of the ENSO and the QBO on the RO trend patterns. Detection results are presented in the second part of this section.
The contribution of solar forcing to global warming was broadly discussed by Hegerl et al. (2007) and more recently reviewed by Gray et al. (2010). Stott et al. (2003) concluded that in the second part of the twentieth century anthropogenic forcings explained most of the warming, with solar signals being likely much smaller. Over the about last 20 yr solar trends are likely to have been in opposite direction in regard to expected global tropospheric temperature changes (Lockwood and Fröhlich 2007; Kopp and Lean 2011). The solar influence on the RO record comes in via residual ionospheric errors and was estimated as a temperature bias of 0.1 to 0.2 K around 20-km height for solar maximum versus solar minimum conditions (e.g., Gobiet and Kirchengast 2004). This effect should become negligible if more than a single solar cycle is covered. GPS/Met measurements took place around the beginning of the last solar cycle, when solar activity was at a minimum. The beginning of the CHAMP measurements coincided with the maximum of the cycle and lasted until the ensuing minimum, so that the solar signal during the CHAMP period would rather counteract an anthropogenic climate change signal.
a. ENSO and QBO during the RO period
Since the RO record is rather short, ERA-40 was used as proxy to learn about ENSO and QBO in the UTLS for adequate interpretation of these patterns in the RO data. ERA-40 offers a sufficient record length to calculate mean patterns of the different modes.
The N3.4 index for ENSO and the evolution of the 50-hPa zonal wind index for the QBO are shown in Fig. 2 for the RO period. A positive QBO index (QBO+) refers to a phase governed by westerly winds and positive temperature anomalies around 50 hPa. The easterly phase is hereafter referred to as QBO−. The QBO index shows a rather steady progress of the oscillation, covering almost four periods between fall 2001 and summer 2010. The two months of the GPS/Met measurement period in 1995 and 1997 were each influenced by a different QBO phase of medium amplitude. ENSO exhibited only weak to moderate events during the RO period and, most importantly, the first two measurement months fall in a period with almost no ENSO signal. During the fall 2001 to winter 2007 period and after mid-2009, weak to moderate El Niño conditions prevailed, while the months in mid-2006 as well as mid-2007 to mid-2009 were governed by La Niña conditions.
Figure 3 illustrates zonal mean average El Niño, La Niña, and neutral conditions based on ERA-40 data between 1980 and 2001. The RO patterns (not shown) are based on few and not very distinct events but are similar to the ERA-40 patterns; they feature smaller amplitudes and are more confined to the tropics, so that the change in sign of the anomalies occurs around 20°S/N. The mean El Niño patterns were calculated (including strong events which are not present in the RO period) from the 1986/87, 1987/88, 1994/95, and 1997/98 events; the mean La Niña patterns from the 1984/85, 1988/89, 1998/99, and 1999/2000 events; and the neutral patterns are based on the years 1980/81, 1981/82, 1989/90, and 1992/93. All patterns were gained by averaging August to ensuing July values. The 1982/83 and 1991/92 El Niño events, which coincided with volcanic eruptions (El Chichón in 1982 and Pinatubo in 1991), were disregarded. For refractivity and temperature, the strongest and spatially most confined UTLS ENSO characteristics (positive tropospheric and negative stratospheric temperature anomalies for El Niño phases) are found in the tropics and subtropics around 300 to 200 hPa along with significant differences between the two phases, as shown in Fig. 3 (right panels). The significance of the difference of the four El Niño and four La Niña patterns was determined with a Mann–Whitney U test (Milton 1964). For geopotential height, the strongest ENSO pattern signal emerges around 100 hPa, thus significant differences appear also at higher levels than for refractivity and temperature and also for the tropical troposphere. Yet the significant differences are more restricted to the tropical bin (±10°S/N).
The zonal mean average ERA-40 QBO signal is depicted in Fig. 4 and shows a well-defined stratospheric signature of the pattern. The RO patterns (not shown) feature stronger amplitudes of the phases and a more pronounced contrast between the low and midlatitudes. The significance test (Fig. 4, right panel) is based on 8 westerly and 8 easterly QBO phases. They were gained by averaging over all months during a phase with index values >5 for QBO+ and index values <−5 for QBO−. The 1988/1989 QBO+ phase was ignored as only 3 successive months exceeded this threshold value. Significant differences between the phases are restricted to the tropics above 100 hPa and polewards of 30° latitude around 200 and 100 hPa.
The proportion of QBO and ENSO variability on the total variability is finally based on the RO data and illustrated in Fig. 5. It was assessed with R2 in the multilinear regression model (see section 2). The result is consistent with the ERA-40 mean patterns, showing for the monthly RO record of 5 zonal mean bands the greatest QBO-explained variability in the tropical stratosphere (≈30% QBO-explained variability). ENSO governs the tropical UT and the LS above 50 hPa for refractivity and temperature and the tropical troposphere and tropopause region for geopotential height (≈40%–50% ENSO-explained variability). QBO and ENSO both impact the stratospheric refractivity and temperature variability at northern midlatitudes (≈15%–25% explained variability).
b. Characteristics of trend patterns
The RO and fGCM refractivity, geopotential height, and temperature trend patterns are presented in Fig. 6 for the intermittent 1995–2010 (Fig. 6a) and the 2001–10 period (Fig. 6b). The fGCM ensemble shows for temperature across the latitudes a rather smooth tropospheric warming (up to around 50 hPa) and a stratospheric cooling above. Refractivity features a reversed trend pattern, and the geopotential height field exhibits a general increase, following the thermal expansion of the troposphere. The RO record, in contrast, presents more distinct patterns, which are more strongly affected by atmospheric variability.
The trend patterns do not depend much on taking the GPS/Met measurements into account, only the pattern amplitudes increase for the 2001 to 2010 period (cf. left panels of Figs. 6a and 6b). The strong RO temperature trends near 100 hPa in the tropics and above 70 hPa in the northern subtropics and midlatitudes are the most striking differences to the fGCM pattern.
Since the GCMs do not reproduce the QBO signal, it was removed from the observations for all further calculations, using linear regression [Eq. (6)]. The ENSO signal was removed from both observations and models (see section 2). The RO trend patterns without QBO and ENSO signals are shown in Fig. 6, middle panels. Elimination of the QBO and ENSO signal primarily influences the trend amplitudes, but it hardly affects the trend patterns and also does not considerably change the pattern correlations with the fGCMs. Below 100 hPa the tropical and subtropical trends are slightly stronger after removal of QBO and ENSO. The reason for this is that mid-2002–end-2005 was a period of rather pronounced El Niño conditions and regression of ENSO leads to cooler temperature anomalies. On the other hand, mid-2007–mid-2009 was a period of La Niña conditions and regression of ENSO results in warmer temperature anomalies, leading to overall somewhat stronger trends. Above 100 hPa the combination of a positive ENSO and negative QBO index, as occurring in mid-2003 and in the second half of 2005, leads to warmer anomaly values after removal of QBO and ENSO; and vice-versa, anomalies turn cooler when negative ENSO and positive QBO indices coincide as is the case past 2005 (e.g., mid-2006, mid-2008, and mid-2009). For the fGCM patterns, the ENSO signal was removed for each simulation separately before calculating the ensemble mean. Because of the large ensemble, this has only a very small influence on the guess pattern. The ENSO signal has also been removed from the control simulations used to estimate the uncertainty in the signal estimates to reflect the observations and their processing.
To cross check the RO patterns and to make a rough estimate of the influence of the missing months between 1995 and 2001 on the trend signal, we compared them to ERA-Interim reanalysis and RICH, RAOBCORE, and HadAT radiosonde trend patterns. Figure 7 depicts the RO temperature trend patterns for the two investigated periods (left panels) and the equivalent trend patterns of the other datasets. Above about 150 hPa the patterns of the datasets are fairly similar and exhibit comparable amplitudes (ERA-Interim shows somewhat stronger trends and has a characteristic warming/cooling feature above/below 150 hPa, which is a peculiarity of this reanalysis; Poli et al. 2010). Below about 150 hPa all comparison datasets differ more from RO and also amongst themselves, indicating more long-term stability weaknesses of the non-RO data in the troposphere (significant evidence exists for a high long-term stability of the RO data as discussed in section 1).
The ERA-Interim and radiosonde trend patterns for the continuous period 1995–2010 are also inspected (lower panels). Apart from the strength of the trends, they are similar to the intermittent and shorter patterns of the considered RO periods above about 200 hPa; below especially the radiosondes show more tropospheric warming, similar to RO for its more limited periods (RICH and RAOBCORE had shown tropospheric cooling in the tropics for the limited periods).
The potential dependence of the RO 1995–2010 trend patterns on the GPS/Met months can be considered as the main weak point regarding the RO pattern. But the close match of the intermittent 1995–2010 and the 2001–10 RO pattern as well as the basic similarity of the intermittent and continuous patterns of ERA-Interim and radiosondes (with their own limitations in long-term stability) underpin the sufficient quality and robustness of the RO trend patterns and their adequacy for the detection study. Also when comparing continuous 1995–2010 anomaly time series of the other datasets to the RO data good agreement was found throughout the UTLS, including the two GPS/Met months; the latter also show good representativity when compared to annual averages of the years 1995 and 1997 of the other datasets. Furthermore, we recall that the sampling of the GCM data to calculate trends was done fully matching the “observation mask” of the available RO data. Hence, we imposed the intermittent coverage in time from the RO data on both the model fingerprints and the model-based samples of internal climate variability. Sampling model data in the same way as observations is standard practice in detection and attribution and avoids introducing bias to results due to missing data in observations (see Hegerl et al. 1996, 1997, 2007; Stott and Tett 1998). This subsampling can be applied in space or time, and unless the sampling characteristics of data and models at the sampled points are substantially different, this will yield reliable results. Estimating the submonthly sampling error in the RO, which is found very small anyway (see section 2), and then correcting for it (Foelsche et al. 2008a; Steiner et al. 2009; Foelsche et al. 2011), avoids any remaining possibility of the spatial and temporal sampling influencing our detection results.
The similarity of trends for the long and short periods (see below) and the consistency of detection results between both periods (with less significant, but similar signals for the shorter period) further increases confidence that the sampling characteristics in the RO data does not introduce a problem. It is worth remembering that since RO data need no cross-satellite calibration, the early data are fully comparable with the later data (this is also demonstrated by similarity of trend patterns with other dataproducts, discussed below). In addition, we coanalyze the shorter continuous 2001–10 period, which contains no time gap, as well as coanalyze the broad 50°S–50°N single-band region, as a complement to the latitude-resolved analysis.
c. Detection requirements
For the detection analysis, the RO, fGCM, and control trend patterns are transformed into a truncated PICTRL EOF space (see section 2). It is required that these truncated patterns are still able to represent most of the anticipated signal of the data. The fGCM ensembles with their mostly dipole trend patterns (see Fig. 6) only need few EOFs to be reasonably well rebuilt. Figure 8 shows the original RO trend patterns (QBO and ENSO removed), the rebuilt patterns using the first k = 5 eigenvectors, and the pattern correlations between original and rebuilt patterns for 1–20 EOFs. Generally, the pattern correlations increase quickly with the number of retained EOFs, yielding around 60% for 5 retained EOFs. Slightly more EOFs are needed to rebuild the temperature and geopotential height patterns, though the displayed rebuilt patterns in Fig. 8 already capture most of the expected signals and are very similar to the fGCM patterns. For the 2001–10 period (not shown), the pattern correlations between RO and rebuild patterns are slightly lower but similar to the 1995–2010 period.
Optimal fingerprinting is based on matching observed and model-simulated patterns and thus relies on reasonably realistic simulated variability from the climate models at the analyzed space and time scales. The three climate models used feature different variability characteristics (shown with respect to RO variability in Fig. 9). The CCSM3 temperature and refractivity variability patterns show best agreement to the RO patterns, even though the amplitudes are slightly smaller, particularly in the stratosphere. HadCM3 exhibits a less distinct variability between UT and LS, but slightly higher amplitudes than CCSM3 in the troposphere. ECHAM5 shows for refractivity and temperature very strong UT variability (3 times stronger values than the RO record) and, compared to the other models, average LS variability. For geopotential height, the tropical and subtropical ECHAM5-simulated variability is stronger than the RO variability and about two times as strong as in the other two models. All three models show for this parameter slightly less midlatitudinal variability than the RO record does. Even though variance estimates are noisy in time, particularly for the climate system where large-scale variations tend to show autocorrelation, these model variability features are based on several hundreds of years of data, suggesting real differences between the models. Checking RO against ERA-40 and ERA-Interim variability (not shown) exhibits very similar variability amplitudes with only minor differences (below a factor of ≈1.5 for anomaly standard deviations) in the tropical tropopause region for temperature and refractivity. Figure 9 depicts the ratios between the RO and model-simulated variability. The ratios are based on standard deviations of detrended data from 1995 to 2010. The fGCM variability was calculated by stringing the 20 different simulations of 15 yr each together. Results based on PICTRL data (not shown) are consistent to the fGCM ones. The observed and model variance vary at most within a factor of 2 when all model simulations are considered (Fig. 9, left column). Deviations in UT refractivity and temperature are mainly caused by ECHAM5, which contributes with 6 out of 20 simulations (besides CCSM3 with 12 out of 20 simulations) considerably to the fGCM ensemble variability. For geopotential height, the observed and model variability agree very well over the study region. At large, the ensembles of the three models do cover a representative range of observed variability, as required for optimal fingerprinting, particularly for geopotential height. Variability in the upper part of the stratosphere seems underestimated, particularly by ECHAM5 and HadCM3. However, our key results indicate no inconsistency between model and residual variability in the truncated space used in our analysis (see below).
4. Discussion of trend detection results
A residual consistency test from Allen and Tett (1999), see section 2, was used to test the null hypothesis that the climate variability of the control simulations is a realistic representation of the observed variability in the truncated EOF space. The regression residuals in the EOF space are related to a χ2 distribution and can be used to determine the number of retained EOFs, which give a realistic estimate of climate noise. Figure 10 depicts the model to observations ratio of cumulative residual variances for 2 to 16 retained EOFs, and the respective χ2 values corresponding to a 5%–95% confidence limit. The test values are required to fluctuate around one and to remain within the confidence limits for a low number of EOFs. Compliance with these requirements shows that observations and models exhibit comparable variance for the number of retained EOFs and corresponds to being unable to reject the null hypothesis.
For the 1995–2010 period (Fig. 10, left), the refractivity values for less than 17 retained EOFs are arranged close to the upper confidence limit, indicating a relatively high model variability compared to the observations but within consistency. For geopotential height and temperature, reliable test results are obtained for about up to 10 retained EOFs. Similar results are gained for the 2001–10 period, where residual consistency for refractivity is given for up to 14 retained EOFs, for geopotential height and for temperature for up to 6 EOFs. These numbers of EOFs retain in case of refractivity 99%, in case of geopotential height 98%, and in case of temperature between 91% (6 EOFs) and 97% (10 EOFs) of the variance of the signal. Beyond that, the regression residual increases sharply. This may either indicate sampling problems of noise, or, more likely, focus on small-scale features of the signal that are not reliably modeled in the GCMs. We exclude these features by limiting truncation as described above. Different truncations for different data products are reasonable, as the EOF truncations resolve different data characteristics.
Detection of a climate change signal can be claimed if the null hypothesis that the observed climate change is part of the natural climate variability can be rejected. In this case, the uncertainty range (10 percentiles or 5 percentiles) excludes zero, that is, there is a less than 10% or 5% chance of a type-1 error of the null hypothesis being true. Based on the distribution of the control scaling factors , which can be assumed to be Gaussian (Hegerl et al. 1996), the uncertainty in the RO-based scaling factors can be determined in a one-tailed t test, with the results indicating that the observed trend pattern is significantly distinct from internal climate variability and thus likely represents a changing climate, for example, in response to greenhouse gas increases. Consistency between the observed and forced climate signal is given, when the RO scaling factors and uncertainty ranges include unity.
Figure 11 shows the RO scaling factors and their 10%–90% and 5%–95% uncertainty ranges (error bars), which are estimated as ±1.3 control standard deviations and ±1.7 control standard deviations, respectively. For the shorter 2001–10 trend period (Fig. 11, right), the uncertainty ranges are broader than for the 1995 to 2010 period (Fig. 11, left). The RO scaling factors for refractivity and geopotential height, which correspond to numbers of retained EOFs that passed the residual consistency test, are close to unity. For temperature, the scaling factors vary around 0.4 if the large-scale aspects of the pattern are considered. When 7 or more EOFs are retained, the scaling factors “jump” by about 1.5. This jump reflects the strong RO temperature trends in the tropics, which start to be resolved when 7 or more EOFs are considered (Fig. 8, bottom left and right).
At the chosen α = 10% significance level, climate change based on the intermittent 1995–2010 period cannot be detected for refractivity, where the α values range between about 13%–30% for up to 10 retained EOFs. For geopotential height, detection can be affirmed for 5–9 retained EOFs at the 10% significance level. For temperature detection cannot be claimed at a 10% significance level (α values around 30%) for up to 6 EOFs. For 7–10 EOFs, which also pass the residual consistency test, a scaling of about 2 is needed to adjust the fGCM trend pattern to the RO pattern. This indicates that RO temperature trends, especially in the tropical lower stratosphere (see Fig. 8, bottom left), are significantly stronger than fGCM trends. Steiner et al. (2009) found, in a classical time series analysis for trends in February 1997 and February 2002–08, tropical lower-stratosphere trends in RO data also somewhat stronger than in GCMs (while at the same time consistent with radiosonde trends). For the shorter 2001–10 period, detection at the 10% significance level is not yet possible except for the full 50°S to 50°N band (see below), but the scaling factors for low-order EOFs vary only little around one. Depending on the scale of resolved patterns, confidence levels of around 70%–75% are achieved for refractivity and of around 80% for temperature. Geopotential height results are close to detection at a 10% significance level, with confidence levels of about 85% for more for 3–6 retained EOFs.
Averaging over all latitudes between 50°S and 50°N (not separately shown) yields detection at a 10% significance level or better for geopotential height and temperature, for both the intermittent period and the continuous 2001–10 period. Only for refractivity, the detection results are less clearly pronounced. For the continuous period, detection depends on the number of retained EOFs, while for the intermittent period residual consistency is not achieved for refractivity. Our results are in addition robust to not subtracting the sampling error (at least 90% confidence for geopotential height and refractivity). If ENSO and QBO are not subtracted, confidence levels decrease but results are broadly consistent. Furthermore, the RO detection results are consistent with results of a stability test (also called perfect model study) employing a GCM simulation for the RO record (not shown either). Using RICH and RAOBCORE temperature fields at the same observation mask as RO data (not shown) leads to detection results consistent with RO in both periods, similarly showing an emerging climate change signal, confirming that the RO measurement of changing atmospheric temperatures is robust. Using ERA-Interim instead, detection results are comparable to RO results for geopotential height in the intermittent period. For temperature and refractivity, signals of climate change cannot yet be identified from the ERA-Interim product. As the ERA-Interim is a blended data and model product, this difference is not a concern. For the shorter continuous period, all data products, including ERA-Interim, show consistent detection results compared to those we obtained using the RO data.
Summarizing, the multisatellite RO records show an emerging climate change signal, which is—most clearly in geopotential height fields that reflect overall tropospheric warming—consistent with the GCM projections for the given intermittent 15-yr period and the continuous 9-yr period. The results are consistent with theoretical RO detection time estimates from two simulation studies. Leroy et al. (2006a) used 12 IPCC AR4 GCMs and an optimal detection approach to estimate the detection time in log-dry pressure trends (which can be related via the local scale height to geopotential height trends as function of pressure as we used here) with a 95% confidence as 6–13 yr. These authors also stated that geopotential height may be the quantity detected earliest, which is found confirmed. Simulations of UTLS RO bending angle profiles (from which refractivity derives by weighted integration of bending angles over height) were used by Ringer and Healy (2008) for detection time estimates of 10–16 yr based on an autoregressive model approach. The results are furthermore consistent with the estimated trend patterns in refractivity, relative pressure/geopotential height, and temperature from GCM simulations discussed by Foelsche et al. (2008c); these authors also discussed why refractivity is expected to have less detection capability, which is confirmed here as well.
Naturally forced climate change is probably small and of opposite direction over our study period. No major volcanic eruption took place during the investigated period. Changes in solar radiation probably provides only a minor forcing for the period 1995–2010 (see Kopp and Lean 2011, who discuss natural and anthropogenic changes over 1980–2010 in detail), and the solar cycle would have most likely caused tropospheric cooling over the shorter period 2001–10. This suggests that the detected change, most clearly the overall tropospheric warming, has been mainly caused by anthropogenic influence on climate.
5. Summary and conclusions
This study aimed to investigate the usability of the RO record in climate change science. An optimal-fingerprinting approach, which enhances the signal (trend) to noise (natural variability) ratio, was applied to test whether RO observations within the period 1995 to 2010 exhibit an UTLS climate change pattern in refractivity, geopotential height, and temperature fields, which is consistent with the expected climate change signal as projected in GCMs. Former studies (Leroy et al. 2006a; Ringer and Healy 2008; Foelsche et al. 2008c), based on GCM or simulated RO data, have shown that a climate change signal should become detectable in RO parameters within a 6–16-yr record.
The influence of QBO and ENSO, the main patterns of UTLS variability, on trend estimates was assessed via multilinear regression. The RO data show a clear response to both modes of variability, with significant increases in geopotential height throughout the atmospheric column in the tropics, and significant tropospheric warming in response to El Niño. The response to the QBO is strongest in the stratosphere, with increased geopotential height and stratospheric temperature in the tropics. The QBO was removed from the observations for further analysis only since the considered GCMs do not include QBO variability. The ENSO signal was removed for both observations and GCMs to get a cleaner climate change signal. Two trend periods were considered, 10/1995 to 07/2010 intermittently and 09/2001 to 07/2010 continuously. Though only an intermittent RO record is available, the trend pattern is consistent with the intermittent 15-yr as well as the continuous 9-yr radiosonde and reanalysis pattern above around 150 hPa.
An emerging climate change signal, consistent with model estimates, can be detected at a 90% confidence level both for the intermittent 15-yr and the continuous 9-yr RO geopotential height record; the 9-yr record so far in the broad 50°S to 50°N band only. The latitude-resolved trend signals of the 9-yr period are still masked by natural variability, indicated by a broad distribution of control trends, which makes climate change detection difficult, even though the scaling factors indicate a consistency between RO and GCM trend patterns. Low confidence is achieved for temperature when only large-scale patterns are considered. Retaining more than 6 EOFs, RO temperature trends are stronger than the GCM projected trends and allow for climate change detection at a 95% confidence level. For refractivity, detection is not yet possible at a 90% confidence level. These results are in line with the expected emergence of climate change signals as discussed in the previous simulation-based studies cited above.
Since natural changes (solar, volcanic, and internal long-term variability modes) have been found of minor relevance during the period, our results suggest that the changes detected in geopotential height—reflecting overall tropospheric warming—and temperature have been mainly caused by anthropogenic influence on climate. Concluding, our study showed that the GNSS RO method offers a high-quality climate record for UTLS climate variability and change studies, and that these data capture an emerging trend signal, which is expected to gain further in significance over the coming years.
We thank B. Scherllin-Pirscher (UCAR, Boulder, CO, USA and Wegener Center) and J. Fritzer, F. Ladstädter, and U. Foelsche (Wegener Center) for their efforts in RO processing and OPS system development; we are especially thankful to B. S-P. for climatology preparations. The authors acknowledge UCAR/CDAAC (Boulder, CO, USA) for the provision of RO phase delay and orbit data. ECMWF (Reading, UK) is acknowledged for access to their global operational analysis, reanalyses and forecast data, Hadley Centre/MetOffice (Exeter, UK) for the provision of HadAT2 radiosonde data, and L. Haimberger, University of Vienna, for the provision of RICH and RAOBCORE radiosonde data. We also thank two anonymous reviewers for valuable comments that helped to improve the paper. We acknowledge the modeling groups, the Program for Climate Model Diagnosis and Intercomparison (PCMDI) and the WRCP’s Working Group on Coupled Modelling (WGCM) for their roles in making available the WRCP CMIP3 multimodel dataset. Support of this dataset is provided by the Office of Science, U.S. Department of Energy. The Climate Prediction Center (United States) is acknowledged for making ENSO and QBO indices available. The work was financed by the Austrian Science Fund (FWF) under research Grants P18733-N10 (INDICATE), P21642-N21 (TRENDEVAL), P22293-N21 (BENCHCLIM), and by NSF Grant ATM-0296007, and regarding OPS system development by ESA/ESTEC Noordwijik and FFG/ALR Austria.