Application of the Singular Spectrum Analysis Technique to Study the Recent Hiatus on the Global Surface Temperature Record

Global surface temperature has been increasing since the beginning of the 20th century but with a highly variable warming rate, and the alternation of rapid warming periods with ‘hiatus’ decades is a constant throughout the series. The superimposition of a secular warming trend with natural multidecadal variability is the most accepted explanation for such a pattern. Since the start of the 21st century, the surface global mean temperature has not risen at the same rate as the top-of-atmosphere radiative energy input or greenhouse gas emissions, provoking scientific and social interest in determining the causes of this apparent discrepancy. Multidecadal natural variability is the most commonly proposed cause for the present hiatus period. Here, we analyze the HadCRUT4 surface temperature database with spectral techniques to separate a multidecadal oscillation (MDV) from a secular trend (ST). Both signals combined account for nearly 88% of the total variability of the temperature series showing the main acceleration/deceleration periods already described elsewhere. Three stalling periods with very little warming could be found within the series, from 1878 to 1907, from 1945 to 1969 and from 2001 to the end of the series, all of them coincided with a cooling phase of the MDV. Henceforth, MDV seems to be the main cause of the different hiatus periods shown by the global surface temperature records. However, and contrary to the two previous events, during the current hiatus period, the ST shows a strong fluctuation on the warming rate, with a large acceleration (0.0085°C year−1 to 0.017°C year−1) during 1992–2001 and a sharp deceleration (0.017°C year−1 to 0.003°C year−1) from 2002 onwards. This is the first time in the observational record that the ST shows such variability, so determining the causes and consequences of this change of behavior needs to be addressed by the scientific community.


Introduction
The record of global mean temperature anomalies (GMTA), available at monthly resolution since 1850, is a noisy time series that is complicated to analyze due to the large number of forcing factors. Among the different contributions, two basic modes of variability stand out from the rest [1,2]: a secular trend (ST) with no obvious oscillation cycles [3][4][5][6], and a multidecadal variability (MDV) resembling natural oscillations such as the Pacific Multidecadal Oscillation (PDO) [7] or the Atlantic Multidecadal Oscillation (AMO) [8,9]. Furthermore, there are numerous internal and external forcings that contribute to climate variability, such as solar cycles [10], stratospheric sulfur aerosols [11] or stratospheric water vapor [12], at a much higher frequency.
The combination of ST + MDV drives the main dynamics of the GMTA series [3], drawing a multidecadal variation over a background trend of increasing temperature. The alternation of periods with accelerated warming rates with others having lower rates, were described almost 30 years ago [1], whereas the actual temperature deceleration period (so called hiatus) had already been predicted 20 years ago [2]. Most of the recent research tends to agree with natural variability being the most probable cause of the last 15-year hiatus in surface temperature [4,5] as the MDV is currently on a negative phase [7]. An increased storage of heat in the deep ocean is one of the most recurrent explanations derived from modeling exercises [8,13] and from measurement records [14,15]; negative tropical sea surface temperature anomalies have also been recently invoked as a major forcing for the recent warming slow-down [16]. However, a global consensus on the causes of the present-day hiatus is still far from being reached.
The identification and attribution of the different modes of variability is a recurrent problem in temperature time-series analysis, and different methods have been employed to discriminate the different signals contributing to the total variability of the series [9,17]. An adequate analysis method should be nonparametric and objective and thus not depend on the choice of the included processes [9]. Singular spectrum analysis (SSA) fulfills these requisites and hence has been successfully employed for analyzing geophysical time series, such as air temperature [2] and sea surface temperature [18].

Results and Discussion
By applying SSA to the latest surface temperature database (HadCRUT4) at annual resolution, it is found that nearly 88% of the total variability of GMTA is due to the combination of ST (78.8%) and MDV (8.6%) [4,6] (Fig. 1). ST is obtained from the combination of eigenvectors (EVs) #1 + #2, while MDV is the result of EVs #3 + #4 (Fig. 2a).
The remaining 12% of the series contain very little information, with no significant trend and almost no autocorrelation (Fig. 2), so it could be considered as white noise. Higher-frequency oscillations (such as those associated with solar activity) [10] are contained in the residuals; e.g., EVs #7 and #8 define a 13year period signal, while EVs #9 and #10 define an 18-year period signal (also identified through the slightly significant autocorrelation values at time lags 13 and 18 years in Fig. 2c). However, as their contribution to the total signal is rather weak (below 1%, Fig. 2a) and their oscillation period is rather short, it is acceptable not to consider them in subsequent long-term analyses [17]. The presence of these higher-frequency cycles in the residuals is, however, evidence of the capability of SSA to separate the signals existing in the original time-series.
The reconstructed surface temperature series obtained by adding ST+MDV (black line in Fig. 1) presents three hiatus periods with no warming or even with small cooling trends, 1878-1907; 1945-1969 and 2001-2013, roughly one every 62 years [19]. This is more noticeable if the annual gradient of each signal is plotted against time (Fig. 3). Obviously in this figure, the stalling periods of the reconstructed series are coincident with cooling trends in the MDV. Results of our analysis are, henceforth, clearly suggesting that natural climate variability is a key player in global temperature hiatus periods, as previously suggested [4,5]. On the other hand, when MDV is on a positive phase (e.g., from 1970 to 1995) the total warming rate greatly accelerates (Fig. 3) as also suggested in previous analysis [19,20].
On the background, the ST has shown positive warming rates since the beginning of the 20 th century (Fig. 3a, red line) with an increasing trend up to the final years of the 1990's. It is also clear that the ST continued to increase during the two first hiatus events ( Fig. 1), showing small variations of its warming trend (Fig. 3a). However, the current hiatus is somehow different (details in Fig. 3b). After the maximum warming rate associated with MDV was reached by approximately 1990, ST showed a distinct peak from 1992-2001, with an unprecedented increase of its warming rate from 0.0085uC year 21 to 0.017uC year 21 , almost doubling in one decade. After this warming rate peak, the ST shows a pronounced decline, 0.017uC year 21 in 2001 to 0.003uC year 21 , in 2013. This type of quick fluctuations in the ST warming rate has no precedent in the observational record (Fig. 3a). Thus, our analysis appears to indicate that the recent hiatus in global surface air temperature is mostly attributable to a negative phase of the MDV (as the two previous ones) plus a reduction of the warming rate associated with the ST (but still showing positive values).
The fluctuations observed in the warming rate associated with the ST during the last 20 years (1992-2013) could not, in principle, be due to a statistical artifact of the analysis method, the so-called 'end-effect', because of the underlying mathematics in such method. However, and in order to assess the existence of such a problem, we applied the same SSA analysis described below (see Material and Methods) to a fraction of the total GMTA from 1850 to 1960 when the previous hiatus was present (see Fig. 1). In this analysis (not shown) the ST did not display any variations at the end of the series so the recent fluctuation could not be attributed to a methodological problem with the used statistics but rather seems to correspond to a distinct dynamics of the observed temperature series. However, the relatively short time-period in which such fluctuation happened (,20 years) prevents more solid conclusions to be derived.
Analyzing surface air temperature data from the different hemispheres with the same technique gives results that are quite  Therefore, the very recent strong changes observed in the warming rate associated with the ST appear to be a global phenomenon that had not occurred before (at least not during the last 160 years). It could not be attributable to MDV or any other form of climatic variability (such as solar cycles), as the different contributions are effectively separated by the SSA analysis (Fig. 2). This unprecedented modification of the ST behavior should be more deeply studied by the scientific community in order to address whether a change in the global climate sensitivity [21] has recently occurred.
The unique fluctuation in ST warming rate during the recent decades could have different origins, such as the proposed recent shift of the tropical Pacific Ocean to a quasi-permanent La Nina condition [16], could be related to the enhanced melting of Arctic ice in recent decades [22] or to the increasing stratospheric aerosol content [11]. A spatially-explicit analysis of the characteristics of the ST and MDV could, potentially, serve to better identify the likely mechanisms involved. Unfortunately, a 2D analysis could not be performed on the gridded HadCRUT4 dataset due to the sensibility of SSA to the presence of gaps in the analyzed data.   Such gaps become more common when the studied area reduces preventing, thus, the analysis of the spatial distribution of the different signals.
Such an analysis could have also helped in clarifying the relation of the MDV to some of the previously described natural modes of surface temperature variability (e.g., PDO and AMO). However, a simple phase analysis (not shown) seems to indicate that the MDV is closer to AMO than to the PDO.
While it is beyond the scope of the present paper to elucidate the causes of the observed variability in the ST of surface temperatures, its unique presence in the entire observational timeseries deserves some more detailed analysis in the coming future. Further studies on the roles and interplays of the atmosphere, the cryosphere and the hydrosphere are necessary to understand the causes of the observed changes. Such studies are crucial to being able to develop models that allow a realistic representation of the present day conditions [13,20] and hence to create plausible scenarios of future climate evolution [5].

Surface temperature data
Mean surface air temperature anomalies were obtained from the HadCRUT4 dataset [23] (available at http://www.cru.uea.ac. uk/cru/data/temperature/, data downloaded on 01/2014). Different datasets spanning the period 1850-2013 were obtained: global mean temperature, temperatures of the Southern and Northern Hemispheres and the gridded dataset at 5u65u resolution.

Singular Spectrum Analysis (SSA)
Singular spectrum analysis (SSA) is a methodology specifically designed to extract information from noisy time series [24] and is analogous to applying an extended empirical orthogonal function (EEOF) analysis to successive lags of a univariate time series [25]. The particular SSA applied here was performed with the Rssa package of the statistical software R, freely available from the Comprehensive R Archive Network (CRAN, http://cran.rproject.org/). Using such a freely available code supports the reproducibility of the analysis performed in the present manuscript [26]. SSA allows decomposition of the time series into a sequence of elementary patterns of behavior that are classified as either trends (ST in this work) or oscillatory patterns. SSA was applied to the annual temperature anomalies in the different studied regions and to the different time series obtained for individual months. The latter analysis allows for the determination of whether the observed ST presents some type of seasonality and also allows for the computation of confidence intervals around the mean, as shown in the corresponding figures. In total, 81 eigenvectors were determined from each time series (of 161 years) using the Broomhead/King method to calculate covariance matrices [27].