The vast majority of variability in the instrumental surface temperature record is at annual frequencies. Systematic changes in the yearly Fourier component of surface temperature have been observed since the midtwentieth century, including a shift toward earlier seasonal transitions over land. Here it is shown that the variability in the amplitude and phase of the annual cycle of surface temperature in the northern extratropics is related to Northern Hemisphere atmospheric circulation as represented by the northern annular mode (NAM) and the Pacific–North America mode (PNA). The phase of the seasonal cycle is most strongly influenced by changes in spring atmospheric circulation, whereas amplitude is most strongly influenced by winter circulation. A statistical model is developed based on the NAM and PNA values in these seasons and it successfully predicts the interdecadal trends in the seasonal cycle using parameters diagnosed only at interannual time scales. In particular, 70% of the observed amplitude trends and 68% of the observed phase trends are predicted over land, and the residual trends are consistent with internal variability. The strong relationship between atmospheric circulation and the structure of the seasonal cycle indicates that physical explanations for changes in atmospheric circulation also extend to explaining changes in the structure of the seasonal cycle.
The annual cycle of surface temperature has changed over the last half century with systematic shifts toward earlier seasonal phasing on land, later seasonal phasing over the ocean, and smaller annual amplitude over land (Thomson 1995; Thompson 1995; Mann and Park 1996; Wallace and Osborn 2002; Stine et al. 2009). Many mechanisms have been suggested to explain variability in the phase of the annual cycle of surface temperature, and we briefly review four of these and focus on a fifth.
First, Thomson (1995) proposed that changes in temperature seasonality exist due to changes in the relative sensitivity of surface temperature to annual insolation forcing associated with two different processes—the annual variability in the Earth–Sun distance (the anomalistic year) and the annual variability in the tilt of the Earth’s rotation axis relative to the sun (the tropical year). Second, loss of sea ice was found to result in a shift toward later seasons, making this possibility unlikely to be responsible for terrestrial shifts toward earlier seasons (Mann and Park 1996), though this mechanism may be a factor contributing to phase delays observed in the high-latitude North Atlantic (Dwyer et al. 2012). Third, large-scale decreases in soil moisture were indicated by model calculations (Seager et al. 2007; Sheffield and Wood 2008; Dai 2011) and these could decrease the lag time between insolation and temperature by lowering the thermal inertia of the surface (Stine et al. 2009, hereafter SHF2009). However, such moisture trends are contradicted by the sparse available observations (Robock et al. 2000; Vinnikov and Yeserkepova 1991), Fourth, changes in shortwave optical properties alter the lag of temperature behind insolation in an energy budget model (SHF2009), but modeling studies have failed to explain the observed twentieth-century phase shift as a response to either direct CO2 forcing (Mann and Park 1996), to combined greenhouse gas and aerosol forcing (Wallace and Osborn 2002), or to the complete twentieth-century forcing history (SHF2009). For instance, the shift toward earlier seasons is not predicted by any of 72 simulations of twentieth-century climate made using 24 different general circulation models included in the CMIP3 archive (SHF2009).
A fifth, and much older idea, is that variability in the seasons arises from the changes in the relative influences of land and ocean at different locations (Brooks 1917, 1918, 1919; Gorczyński 1920; Brunt 1924; Simpson et al. 1924; Conrad 1946). More recently, Thompson (1995) suggested that trends in the phase and amplitude of the annual cycle of surface temperature found in six long European temperature records resulted from variability in ocean and atmosphere circulation. The presence of phase trends in a 100-yr control simulation of the Geophysical Fluid Dynamics Laboratory climate model (Manabe et al. 1991) also lead Mann and Park (1996) to suggest that observed trends in seasonality may be a manifestation of century-scale natural climate variability. Furthermore, SHF2009 found a significant relationship between Northern Hemisphere averaged phase and the northern annular mode (NAM) index, though that initial analysis was not conclusive.
Here we expand upon the analysis of SHF2009 to include the Pacific–North America (PNA) pattern, spatially resolve the relationship between seasonality and atmospheric circulation, and account for the seasonal variability in atmospheric circulation. Variability in solar insolation and its relationship with the calendar year are also better accounted for. In the context of this hypothesis, it is worth noting that general circulation models forced with the observed twentieth-century forcing not only fail to capture observed trends in temperature seasonality, as mentioned, but also generally fail to reproduce the observed trends in atmospheric circulation (Gillett et al. 2003; Gillett 2005; Cohen et al. 2005; though also see Deser and Phillips 2009). Thus, the hypothesis that changes in atmospheric circulation are responsible for changes in the structure of the seasonal cycle is consistent with the failure of general circulation models to reproduce the trends in either, but this affords rather circumstantial evidence, and here we pursue testing for a direct relationship using an observational analysis approach.
2. Estimation of temperature seasonality
Variations in the structure of the seasonal cycle can be quantified in many ways including threshold crossing statistics (Schwartz et al. 2006; Sparks and Menzel 2002; Cayan et al. 2001), changes in variance (Thompson and Wallace 2000), and other seasonal indicators (Westerling et al. 2006; Abatzoglou and Redmond 2007; Luterbacher et al. 2004). Here, following SHF2009, we rely upon fitting an annual sinusoid to the observations. In particular, given 12 monthly values of surface temperature we estimate the yearly component as
where X(t) represents 12 monthly values of surface temperature or solar insolation. From Eq. (1) we recover the amplitude,
of the annual sinusoid. This represents a parsimonious way to describe the observed surface temperature record, explaining on average 97% of the variance of monthly departures from the annual mean temperature over Northern Hemisphere extratropical land and 95% over ocean. The annual sinusoid is well resolved in monthly average data (Thomson 1995). For example, our tests using the long Northern Hemisphere extratropical records in the Global Historical Climatological Network (Peterson and Vose 1997) indicate an average correlation coefficient (i.e., Pearson’s R) between estimates of phase calculated using daily versus monthly average temperature of >0.99. Although this approach will not capture seasonally asymmetric changes that project, for example, onto the two cycles per year Fourier component, empirical examination has shown that the structure of large-scale temperature changes over the past century are well described as a combination of changes in annual mean temperature and transformations of the 1/yr component (SHF2009).
The amplitude and phase of surface temperature at each point for which we have observations are referenced to the amplitude and phase of local incoming solar radiation at the top-of-atmosphere (insolation). We define the gain (G) as the ratio of the amplitude of local surface temperature variability to the amplitude of insolation,
and define the lag (λ) as the difference between the phase of the annual cycle in surface temperature and the phase of the annual cycle in insolation, expressed in days,
In SHF2009 we defined φsun and Asun as simple functions of latitude, but here we also apply time-varying corrections associated with the calendar and changes in Earth’s orbital configuration. In particular, we account for the mismatch between the length of the Gregorian year and the tropical year, the change in the amplitude of total insolation forcing due to changes in the phase relationship between tropical and anomalistic year insolation, and changes in the phase of the total year in insolation forcing relative to the tropical year due to the secular trend in the phase of anomalistic year forcing relative to the tropical year. In general, these lead to a decrease in the magnitude of phase trends over Northern Hemisphere land, though they are still statistically significant. The details and magnitudes of each correction are discussed in detail in the appendix.
We use the University of East Anglia Climate Research Unit’s 5° × 5° gridded surface temperature anomalies (Brohan et al. 2006), plus gridded climatology (New et al. 1999) to identify surface phase and amplitude variability. Over land we use the Climatic Research Unit gridded land surface air temperature version 3 (CRUTEM3) dataset, and over the ocean we use the the Second Hadley Centre Sea Surface Temperature (HadSST2) dataset. Grid boxes are classified as land or ocean depending on which accounts for the larger proportion of the grid box, as determined by the Clark U.S. Navy Fleet Numerical Oceanographic Center Land/Ocean Mask (Cuming and Hawkins 1981). Here G/λ estimates are made only when temperature estimates exist for all 12 months in a given year. We consider grid boxes where G/λ estimates exist for at least 50 of 60 years, and where the 1/yr Fourier component explains, on average, >85% of within-year variance. In section 3 we also examine the vertical distribution of G and λ in the atmosphere at 45°N, which we estimate using the National Centers for Environmental Prediction (NCEP) reanalysis product (Kalnay et al. 1996).
The analysis is restricted to Northern Hemisphere extratropical locations because this region has the best data coverage, is most likely to be influenced by Northern Hemisphere atmospheric dynamics, and has variability that is better described by annual period variability than that of the tropics. Twice-per-year variability in insolation forcing becomes important in the tropics. To the extent that dynamical modulation of the annual cycle exists in the Southern Hemisphere, it is more likely to be related to Southern Hemisphere atmospheric variability, which is not considered here.
Winter atmospheric circulation is stronger and more variable than that of any other season, and we expect the dynamics in this seasons to most influence the seasonal cycle. Thus, we adopt a winter-centered year, extending from August to July for the purposes of applying Eq. (1), denoted by the calendar year in which January occurs. This choice of calendar year has essentially no effect on our estimate of the trends in G and λ but does reduce the interannual variability in both of these quantities. Furthermore, the variations in both quantities also becomes spatially smoother, and this has the effect of slightly decreasing the effective spatial degrees of freedom (ESDOF), which we estimate using the moment matching method of Bretherton et al. (1999). When using winter-centered years, there are 11 ESDOF for Gland, 24 for Gocean, 17 for λland, and 15 for λocean, which is 8 ESDOF smaller for λ and 1 ESDOF smaller for G than the values recovered using summer-centered years for both land and ocean. The reduction in degrees of freedom supports using the winter-centered calendar and has the effect of making it more difficult to identify changes in the structure of the seasonal cycle as being statistically significant. As a check of the robustness of our approach, we also performed the analysis shown in this paper using three-year increments of data tapered using a hamming window, rather than the 12 month boxcar window implied by Eq. (1), and obtained qualitatively identical results.
3. The seasonal cycle and trends in its gain and lag
We consider G and λ trends for 1951–2010, a period for which there is good spatial coverage for surface temperature variability and for which reliable estimates of atmospheric circulation are available. It is useful to break up the analysis between land and ocean because variability in the structure of the seasonal cycle is proportional to G, which is much lower over the ocean than the land. In particular, the variability in G is directly proportional to the local climatological average G (Fig. 1), and the variability in λ is inversely proportional to the climatological G, where the latter relationship arises because randomly oriented noise more easily alters the phase of a sinusoid whose amplitude is small (SHF2009, supplementary material). The pointwise significance of trends (Fig. 2) are determined by testing local trends against a Monte Carlo distribution of trends, generated at each grid box, by randomizing the local G and λ time series and calculating the trend of the resulting time series 10 000 times (Schreiber and Schmitz 2000). To estimate the ocean-average and land-average trends, we use the mean of the distribution of 1951–2010 local trends over each domain, and evaluate the significance of the distribution means using a t test with the ESDOF described above.
Land gain: The mean trend in gain over land, Gland, is −2.8°C kW−1 m2 (60 yr)−1 (p < 0.05)—a trend toward smaller amplitude seasons. Specifically, gain trends are negative in the interior of Eurasia, western Canada, and the central United States, although they are positive near the Mediterranean, Labrador, and in parts of the western United States (Fig. 2a). Unlike in SHF2009, the time variability of the amplitude of solar insolation, which varies because of the precession of the equinoxes and changes in Earth’s obliquity, are accounted for, although this effect is ~100 times smaller than the observed trends (see appendix).
Ocean gain: The mean trend in gain over the ocean, Gocean, is 0.43°C kW−1 m2 (60 yr)−1 (p = 0.2)—an insignificant trend toward larger amplitude seasons. The strongest Gocean trends are toward larger amplitude and are found in the North Atlantic near the Labrador coast and in the central North Pacific. Weaker trends toward smaller Gocean are seen in the western Pacific Ocean, off the coast of British Columbia, and in a tongue extending from the eastern coast of the United States. Note that the heterogeneous pattern of changes in ocean gain and other components of the seasonal cycle suggests that the seasonal cycle over land and ocean can be further differentiated, which we will return to in the next section.
Land lag: The mean trend in lag over land, λland, is −1.3 days (60 yr)−1 (p < 0.02)—a significant trend toward earlier seasons. Trends toward earlier seasonal transitions are seen over most of North America, over northern Asia, and over Europe where these trends are particularly strong. Trends toward later seasonal transitions are seen in Labrador, in the southeastern United States, and in the vicinity of the Indian subcontinent. The effects of calendric corrections, in the form of leap years and changes in the phase of the seasonal cycle because of the precession of the equinoxes, is to shift the phase of insolation forcing toward earlier seasonality by, on average, 0.6 days, making the magnitude of the average land phase trend reported here smaller than those identified in SHF2009.
Ocean lag: The mean trend in lag over the ocean, λocean, is 0.83 days (60 yr)−1 (p = 0.3). This is larger than the trend reported in SHF2009 because of a 0.6-day correction toward later seasons. Strong trends toward earlier seasons are present in the eastern Pacific Ocean, and the midlatitude western Atlantic has weaker trends in the same direction. Trends toward later seasons are seen in the central and western Pacific and throughout the eastern Atlantic.
The contrast between the trends observed on land and the trends observed over the ocean is worth emphasizing. Land trends are toward smaller gain and earlier seasons, whereas ocean trends are toward larger gain and later seasons. This pattern implies that mechanisms that would tend to produce the same sign of trends on land and ocean, such as solar modulation or direct greenhouse forcing, are less likely to be responsible. However, a viable way to produce opposite sign trends on the land compared to the ocean is through systematic changes in the seasonal atmospheric advection of heat associated with variability in circulation.
The influence of the atmospheric circulation upon the mean structure of the seasonal cycle is apparent in the increase of amplitude and decrease in lag associated with increasingly continental conditions. Indeed, a strong relationship exists at any given surface location in the Northern Hemisphere extratropics between the long-term-mean λ the long-term-mean G and distance of that location from the coast, where distance is measured in the westward direction—climatologically “upwind” in the latitude range where circulation is dominated by midlatitude westerlies (SHF2009). This effect can also be discerned in the vertical dimension. For example, moving from the Atlantic to Eurasia at 45°N, the Eurasian G and λ values are more oceanlike in the middle and upper troposphere than at the surface (Fig. 3). Likewise, moving from Eurasia to the Pacific, landlike G values extend over the western Pacific in the middle and upper troposphere.
Overall, the climatological G and λ share a pattern analogous to that of passive tracers advected by prevailing westerly winds, though there are notable exceptions attributable to the dynamic influence of temperature upon circulation, meridional structures in circulation, and the nonlinear relationship between amplitude and phase. To give three examples—cyclonic flow about the Icelandic low produces onshore wintertime flow in northeastern North America, which may account for the appearance of oceanlike G/λ values in eastern North America (Fig. 3). Second, Eurasian land extends eastward as far as the date line at 65°N, which likely contributes to the high G values seen aloft at 45°N in the western Pacific. Third, Rossby waves forced by relatively warm ocean temperatures in winter propagate with westward group velocities (Kaspi and Schneider 2011), amplifying the high G values seen in eastern Eurasia. Because the temperature contrast between land and ocean is larger in early winter than in late winter, these should also force earlier λ in the same region, consistent with an overall Rossby wave influence toward more continental conditions along eastern continental boundaries. The obvious role of atmospheric circulation in producing the climatological G and λ, both at the surface (Fig. 1) and throughout the extratropical troposphere (Fig. 3), provides additional support for the hypothesis that variability in atmospheric circulation plays an important role in modulating seasonality.
4. The relationship between circulation and interannual changes in seasonality
Variability in wintertime Northern Hemisphere extratropical circulation is dominated by two modes of variability: the northern annular mode (NAM; Thompson and Wallace 1998, 2000; Walker and Bliss 1932) and the Pacific–North American mode (PNA; Wallace and Gutzler 1981; Lorenz 1951). Together, these two modes of variability explain 36% of month-to-month variability and virtually all of the long-term trends in wintertime sea level pressure (Quadrelli and Wallace 2004). The NAM and PNA then represent the natural modes of variability to test for control of trends in G and λ by variations in Northern Hemisphere circulation. Given the association of variations in the NAM and PNA with anomalies in seasonal temperatures (Thompson and Wallace 2001; Hurrell 1996; Leathers et al. 1991), they must also influence the amplitude and phase of the seasonal cycle, and our focus here is upon the degree of control and whether these modes are also associated with the trends in the seasonal structure described in the foregoing section.
Estimates of the NAM and PNA are obtained from the National Oceanic and Atmospheric Administration Climate Prediction Center (Zhou et al. 2001; Barnston and Livezey 1987). Note that the PNA records used here are not orthogonal to the NAM, but that we have computed all results using an orthogonalized version, found no difference, and therefore opted to use the standard PNA index due to its greater familiarity. For comparison with atmospheric circulation, we use surface temperature data from the 1951–2010 winter-centered years.
To illustrate the temporal relationship between atmospheric circulation and G and λ, we initially focus on northern Europe because the influence of the NAM on wintertime temperature is quite clear in this region (Thompson and Wallace 2001). In particular, we consider the 1951–2010 annual cycle averaged over land from 50° to 60°N, 5°W to 30°E. Variability in the NAM index is significantly correlated with local gain (R = −0.65, p < 0.001) and with local lag (R = −0.44, p < 0.001). Both correlations are substantially stronger after 1979 [R(G, NAM) = −0.76, R(λ,NAM) = −0.62] when satellite observations are available (Figs. 4a,b). Comparing the annual cycle averaged over the 10% of years with the highest December–March (DJFM) NAM index values against the average over the 10% of lowest DJFM NAM values (Fig. 4c) shows that high NAM years are characterized by a seasonal gain of temperature that is 20% smaller and a lag 5.1 days earlier than the low NAM years.
Suggestively, the maps of trends in both G and λ (Fig. 2) show spatial structures similar to the wintertime surface temperature anomalies associated with the NAM and PNA (Fig. 5, left column). To further examine the degree to which these spatial patterns resemble one another, the NAM and PNA indices are directly correlated against G and λ at each grid point after detrending all time series in order to focus on higher frequency variations. The variability and trends in the NAM and PNA indices are a function of season (Abatzoglou and Redmond 2007), as is their influence upon changes in G and λ. Correlations were performed for every season, and here we focus on those using winter [December–February (DJF)] for G and spring [March–May (MAM)] for λ and later discuss why these seasons are appropriate.
The correlation maps between G and the DJF NAM index is very similar to the correlation maps of temperature and the NAM index, and the same holds when using the PNA index, with the maps having a cross-correlation of −0.93 and −0.89 (p < 0.001) for the NAM and PNA, respectively (Fig. 5). High NAM years in regions such as the eastern United States, northern Eurasia, and the Atlantic off the east coast of the United States show both winter warmth and low G values, and the converse holds in regions such as Quebec and the Middle East. Similar results are obtained for the PNA, with western Canada and the western margin of the United States experiencing winter warming and low gain in high PNA years, and the eastern United States and the central North Pacific experiencing winter cooling and high gain. Of course the relationship between atmospheric circulation and gain is no surprise: changes in wintertime temperature change the amplitude almost by definition and winter temperatures vary more than summer temperatures. More to the point, however, is that the strong correspondence between the correlation maps suggests that atmospheric circulation is a primary cause of G variability.
Spring (MAM) λ correlation maps between atmospheric circulation and λ also resemble the pattern of circulation influence on temperature (Fig. 6), with the relationship stronger for the NAM (R = −0.84, p < 0.001) than for the PNA (R = −0.77, p < 0.01). High NAM index values are associated with warm surface temperatures and earlier λ in northern Eurasia, the eastern United States, and the region of the Atlantic off the east coast of the United States. Conversely, they are associated with cooler surface temperatures and later λ over the Middle East and North Africa, Labrador, and the high-latitude North Atlantic. High PNA values are associated with warmer surface temperatures and earlier λ over northwestern Canada and Alaska, and with cooler surface temperatures and later λ over the southeastern United States and the central North Pacific. Correlations for fall (September–November) are very weak, and we will discuss the reason for this later in this paper.
The foregoing results indicate that circulation influences seasonality in the locations and direction expected from its influence upon temperature. It is also useful to explore what fraction of changes in seasonality can be directly explained by the NAM and PNA. A multiple linear regression against wintertime NAM and PNA indices shows that a significant (p < 0.05) amount of G variability can be explained at 65% of land grid boxes and 33% of ocean grid boxes and that these modes of circulation account for, on average, 22% and 9% of G variability over land and ocean respectively (Fig. 7a). Spring indices of the NAM and PNA explain a significant amount of λ variability at 54% of of land grid boxes and 30% of ocean grid boxes, accounting for an average of 11% and 8% of total λ variability over land and ocean, respectively (Fig. 7b).
5. Predicting trends in seasonality
The existence of a relationship between atmospheric circulation and variability in seasonality at interannual time scales is related to, but does not resolve, the question of whether changes in these modes of atmospheric circulation are responsible for the trends in G and λ. Both the NAM and PNA shifted toward more positive index values from 1951 to 2010. We posit that the trends in seasonality are driven by shifts in atmospheric circulation and that a self-similar relationship exists between these quantities at interannual and longer time scales. This hypothesis need not be true. For instance, the strength of the relationship could differ between high and low frequencies, perhaps because year-to-year heat storage in the ocean mixed layer would lead to a larger sensitivity to persistent changes in atmospheric circulation. Furthermore, if seasonality changes respond to slow changes, for example, in the thermal mass of the surface due to drying (SHF2009) or due to slow astronomical changes (Thomson 1996), then the trends in seasonality could be decoupled from atmospheric dynamics despite the relationship at interannual time scales.
To explore whether the trends in G and λ can be explained by changes in atmospheric circulation, we call upon the multiple linear regression coefficients that we obtained in the foregoing section using the detrended time series. Specifically, we multiply these coefficients by the full (not detrended) NAM and PNA time series, add the two scaled time series, and calculate the resulting 1950–2010 linear trend. This gives the trend in G expected at each point from the NAM and PNA trends, assuming an equal response between the high and low frequency seasonal variability. An analogous procedure is used to predict the λ trends.
The mean G trend over land is explained by the trends in the NAM and PNA. Specifically, the mean of the distribution of observed 1951–2010 Gland trends is −2.8°C kW−1 m2 (60 yr)−1 (p < 0.05), a trend toward a lower amplitude seasonal cycle (Fig. 8). Removal of the dynamic correction decreases the magnitude of the mean Gland trend by 70% to −0.8°C kW−1 m2 (60 yr)−1, which is statistically indistinguishable from zero (p = 0.5). The dynamic correction also narrows the distribution of trends across extratropical Northern Hemisphere grid points, decreasing the variance by 29% and indicating that the detrended relationships between Gland and the first two modes of atmospheric circulation predict the observed trend in both regions with positive and negative phase lag trends.
Removal of the dynamically predicted λ trends likewise decreases the mean trends over land by 68% from −1.3 (p < 0.02) to −0.42 (p = 0.3) days (60 yr)−1, again rendering the residual insignificant. As with Gland, the variance of the λland trend distribution decreases by 22%. Note that the above approach likely underestimates the magnitude of G and λ trends that can be explained by atmospheric circulation because the regression coefficients are subject to regression dilution, which tends to decrease the magnitude of the predictor variable in the presence of noise in the independent variable (Frost and Thompson 2000); we do not pursue a correction for this effect.
Oceanic G and λ do not show a mean trend of any significance (Table 1). However, subtraction of the ocean trends predicted by atmospheric variability decreases the variance of the Gocean trend distribution by 14% and of the λocean trend distribution by 41%, again indicating shifts in atmospheric variability account for regional trends in the seasonal cycle. Because we derive the relationship between atmospheric circulation and the seasonal cycle using detrended data, there is no expectation that the predicted trends will decrease the magnitude of the mean seasonality trends or the width of the trend distributions absent a true relationship, and we take these results as a strong demonstration that trends in the seasonal cycle are attributable to the trends in atmospheric circulation.
The foregoing demonstrates that trends toward earlier λland are associated with circulation-related spring warming, and the question arises as to whether other seasons are also important for determining changes in seasonality? To address this question, we decompose the change in λland attributable to typical NAM atmospheric circulation variability on a given month, m, into three contributions,
That is, changes in λland(m) are determined by multiplying a one standard deviation change in the NAM index on a given month by two different sensitivities: 1) the sensitivity of extratropical land surface temperature to changes in atmospheric circulation, which is determined by regression, and 2) the sensitivity of λland to changes in land temperature, which is determined numerically by perturbing the average seasonal cycle (Fig. 9). The same procedure is also applied to determine how variations in the PNA contribute to changes in λland.
Application of Eq. (6) shows that March NAM variability and February–March PNA variability are most responsible for the variability in λland (Fig. 9). Although the sensitivity of λland to temperature is of equal magnitude in autumn and spring (Fig. 9b), the NAM and PNA are substantially less active in autumn (Fig. 9a). Also notable is that the sensitivity of λ over the ocean is shifted about one month later relative to that of land, as follows from the difference in the mean values of λocean and λland. The annual cycle in λocean variability associated with the NAM and PNA is, therefore, more symmetric between spring and fall over the ocean, and it is more likely in any given year that spring and fall circulation variability will produce offsetting perturbations to λocean, consistent with the more limited explanatory power we find over the oceans.
6. Discussion and conclusions
Observed average trends in Gland and λland can be understood as arising from changes in atmospheric circulation as can, to a lesser extent, regional trends in both G and λ. For instance, λ trends over Europe, Asia, Labrador, and the North Atlantic are consistent with positive spring NAM trends, whereas λ trends over western Canada and the North Pacific Ocean are consistent with positive trends in the spring PNA, and trends over the continental United States can be understood as arising from the competing influence of NAM and PNA trends. Winter trends in the PNA are associated with negative G trends in western Canada, and more positive values of the winter NAM coincide with negative G trends over northern Asia and negative G trends over the Middle East, North Africa and Labrador. One notable exception, however, is that trends in G are positive in much of Europe, despite a winter NAM trend suggesting the opposite response, because of larger warming trends during summer.
The observed 1951–2010 trends in gain and lag can be explained in relationship with changes in the NAM and PNA. The simplest interpretation of how the NAM and PNA influence seasonality is through advection of heat that changes local temperature. Atmospheric circulation could also affect surface temperature through changing rain and snowfall patterns that can, in turn, alter the surface wind field, atmospheric wave field, and setup blocking conditions (Thompson and Wallace 2001, 2000; Hurrell et al. 2003; Thompson et al. 2003).
One could also envision causality propagating in the opposite direction. For example, an earlier spring warming could, in principle, drive a springtime circulation that favors the positive NAM index. Under this hypothesis seasonal temperature anomalies would themselves give rise to NAM- and PNA-like circulation patterns. However, the broad spatial extent of the relationship between local λ and circulation-modulated wintertime temperature (Fig. 5) and the fact that this relationship is observed over land and ocean, and over both arid and wet continental regions, suggests that changes in atmospheric circulation are the predominant source of the change in the seasonality of temperature. Indeed, circulation-induced temperature changes are well documented in observations and in model results (Hurrell 1996; Thompson et al. 2000; Abatzoglou and Redmond 2007), and the primary contribution made here is to demonstrate that these changes in temperature are adequate to explain a significant fraction of the interannual gain and lag variability and virtually all of their multidecadal hemispheric-scale trends. Our conceptual model for the coupling between circulation and seasonality is thus one in which pressure gradients give rise to anomalous circulation with this anomalous circulation then causing convergence or divergence of heat by acting across either the land–ocean temperature gradient (a mechanism which is particularly important in Northern Europe) or across the pole-to-equator temperature gradient (an important mechanism in eastern North America and in East Asia).
The discrepancy in explanatory power between interannual variability and multidecadal trends bears further comment. The NAM and PNA explain 52% of average G variability (using wintertime NAM and PNA index values) and 29% of average λ variability (using spring NAM and PNA index values), where averages are taken over Northern Hemisphere land, and explain a substantially smaller fraction of the pointwise variability. However, the NAM and PNA explain the preponderance of the spatially averaged trends in these quantities. There is a parallel here to the explanatory power of the leading dynamical modes for sea level pressure. The NAM and PNA together explain almost all hemispheric-scale trends in wintertime sea level pressure, but only explain 50% of the interannual variations in sea level pressure variability and only 36% of the monthly variability (Quadrelli and Wallace 2004). Apparently, the NAM and PNA explain relatively more of the low- than high-frequency variability for both sea level pressure and the structure of the seasonal cycle.
SHF2009 found that the shift toward earlier seasons over land during the last half century was anomalous relative to variability between 1850 and 1953, but the present results weaken this conclusion in two ways. First, the demonstrations that the leading two modes of atmospheric circulation explain the majority of the land and ocean trends in seasonal structure suggests that a smaller number of degrees of freedom exist in these quantities than was estimated using the full interannual variability, making it more difficult to reject any null hypothesis that recent mean trends differ from earlier ones (Table 1). Second, SHF2009 found that those regions with continuous monthly thermometer data since 1850 showed that phase trends during the last half century were unprecedented. However, repeating this analysis for winter-centered years and making the necessary temporal calendrical corrections, we find that the 1863–1922 trend, and the 60-yr trends ending one year earlier and one year later, are now the strongest, though recent trends would remain the strongest if summer-centered years were used. More generally, these results indicate that the separation of trend analysis into land and ocean is less fundamental than grouping regions by the sign of their response to NAM or PNA perturbations. The prior finding that seasons are shifting earlier over land appears to be a consequence of the fact that positive springtime excursions of the NAM and PNA cause more of the northern extratropical land surface to warm than to cool and of the fact that the United States and northern Europe are more densely represented in the instrumental surface temperature record than the Arctic and the Middle East.
If observed changes in seasonality are caused by changes in atmospheric circulation, then what produces these changes in atmospheric circulation? Despite an extensive literature on the origin of variability in atmospheric circulation, the source of trends in Northern Hemisphere atmospheric circulation remains controversial.
Natural variability in the midlatitude jet arises from wave–mean flow interactions (Lorenz and Hartmann 2003; Limpasuvan and Hartmann 2000; Feldstein and Lee 1998) and involves interactions between the troposphere and the stratosphere (Baldwin et al. 1994; Baldwin and Dunkerton 1999; Baldwin et al. 2003; Thompson and Wallace 1998; Thompson et al. 2003; Charney et al. 1961; Cohen et al. 2001, 2002, 2010a). Forced variability may be associated with ozone depletion (Volodin and Galin 1998, 1999), greenhouse gas forcing (Shindell et al. 1999; Kuzmina et al. 2005; Rind et al. 2005; Cohen 2010; Cohen et al. 2010b), or some combination of both (Hartmann et al. 2000; Shindell et al. 2001). Models driven with greenhouse gas forcing alone, however, do not capture the magnitude of the NAM trend of the past half-century (Gillett et al. 2005; Deser and Phillips 2009). Simulations forced with observed sea surface temperatures, however, have captured the 1950–2000 deepening of the Aleutian low associated with observed PNA trends (Deser and Phillips 2009), and a more in-depth discussion of the origin of the trends in Northern Hemisphere atmospheric circulation can be found in Deser and Phillips (2009).
Numerous hypotheses have been proposed to explain the existence of trends in the amplitude and phase of the annual cycle of surface temperature, including changes in the effective thermal mass of land surface associated with changes in soil moisture (SHF2009), changes in optical properties of the atmosphere (Wallace and Osborn 2002; SHF2009), and changes in the relative sensitivity of the surface to anomalistic- and tropical-year insolation forcing (Thomson 1995). However, observational evidence for these hypotheses has not been forthcoming. Here we demonstrate a strong relationship between interannual variations in the seasonal cycle and atmospheric circulation, and show that the observed trends in seasonality can be linearly predicted as a response to variations in atmospheric circulation. Thus, it appears that a complete explanation for trends in atmospheric circulation over the last 60 years will also represent a sufficient explanation for the observed changes in the seasonality of surface temperature.
ARS was funded by the Harvard Center for the Environment and by NSF ATM-0628678, “Carbon-Climate Interactions with Increasing Water Demand,” and NOAANA05OAR4311167, “Carbon Dynamics in Changing Hydrologic Regimes.” PH was funded by NSF Award 0902374. Comments from Thomas Laepple, Jeffrey Park, and two anonymous reviewers improved this manuscript.
Variability in the Annual Cycle of Insolation Forcing
We determine G and λ by referencing the phase and amplitude of the annual cycle in surface temperature to the annual cycle in solar insolation [Eqs. (4) and (5)], taking into account temporal variability in insolation due to orbital and calendrical drift. These represent geometrical calculations of top-of-atmosphere insolation as a function of latitude and time. In this appendix, we describe how the temporal variability in Asun and φsun is calculated. Accounting for time variability in insolation forcing also decreases the average variance in λ by 1% and has negligible influence on the amplitude variability.
a. Orbital effects
Referencing the annual cycle in surface temperature to the annual cycle in solar insolation is complicated by the fact that there are two distinct years in solar insolation. Annual forcing arises due to the tilt of the Earth’s rotation axis relative to the plane of the Earth’s orbit around the Sun. This is the tropical year, which is measured as the time interval between which the Earth’s rotation axis returns to a particular tilt relative to the Sun and which has a mean frequency of ftropical = 1/365.242 day−1 (Seidelmann 2005). The anomalistic year in solar insolation, in contrast, arises due to the fact that the Earth rotates about the Sun in an ellipse, rather than a circle. At perihelion, the point of the Earth’s orbit when it is closest to the Sun, the intensity of solar insolation at the top of Earth’s atmosphere in the plane perpendicular to the Sun’s rays is 1412 W m−2, as compared to 1321 W m−2 at aphelion. Due to the precession of the Earth’s rotation axis (at ~1 revolution per 21 kyr), the anomalistic year has a slightly higher frequency of fanomalistic = (365.259 day)−1 (Seidelmann 2005), than the tropical year, which results in a gradual drift in the phase of the anomalistic year relative to the tropical year. As the anomalistic year shifts toward later phases relative to the tropical year, the phase of incoming solar radiation at the annual period also shifts relative to the tropical year, and the amplitude of solar insolation changes due to constructive and destructive interference between the two annual components of insolation. In addition, changes in the tilt of the spin axis change the amplitude of seasonal insolation and the relative strength of the tropical and anomalistic year contributions.
We estimate the time variability of Asun and φsun using the orbital calculations of Berger (1978, 1991) as implemented by Huybers and Eisenman (2006). The calculations of Berger (1991) give the solar insolation as a function of year, latitude, and day of the (tropical) year. We calculate time series of insolation for each desired calendar year, at each latitude, with a time step of 6 h and average these to monthly resolution insolation time series. From these time series we calculate the phase and amplitude of solar insolation using Eq. (1) for each calendar year and at each latitude. These amplitude time series then represent the Asun values used in Eq.( 4). The φsun values must be further corrected for calendric effects. Trends in daily insolation from 1951 to 2010 range from −0.5 to 0.3 W m−2 per 60 years (Fig. A1).
Orbitally induced phase trends over the length of the instrumental record are toward earlier seasonal forcing in the Northern Hemisphere and toward later seasonal forcing in the Southern Hemisphere (Fig. A2a). Trends are strongest closer to the equator and decrease in magnitude as we approach the poles. Orbitally induced Northern Hemisphere extratropical phase trends range from 0.15 days per 60 years at the edge of the tropics to 0.04 days per 60 years at the North Pole. Orbitally induced amplitude trends in the Southern Hemisphere currently decrease the amplitude of solar insolation by about 0.1 W m−2 per 60 years (Fig. A2b), a value which is relatively uniform with latitude, though it represents a larger fractional change in amplitude as one moves toward the equator (Fig. A2c). In the Northern Hemisphere, changes in the orientation of the spin axis are decreasing the amplitude of solar insolation north of 42°N by up to 0.06 W m−2 per 60 years and increasing the amplitude of solar insolation south of 42°N.
b. Calendrical effects
The preference for calendar years that contain an integer number of days necessitates some mismatch between the length of the calendar year and the length of the mean tropical year. The Gregorian calendar damps long-term phase drift between the calendar and solar insolation by introducing leap days that leads to a drift relative to the tropical year of only 1 day in 3268 years over the entire 400-yr leap day cycle, though on shorter time scales the Gregorian calendar follows the tropical year less closely.
To quantify the effects of the Gregorian calendar on phase estimates, we first generate a sinusoid with the frequency and phase of the tropical year,
The Fsun represents solar insolation for a location with sinusoidal tropical year insolation. Here Fsun is sampled daily using the rules of the Gregorian calendar. Specifically, February has 29 days in leap years and 28 days in nonleap years. After taking monthly averages, we compute the phase of insolation, φsun, for each year using Eq. (1), and use the anomalies from the 1950–2010 mean value as an additive correction to the orbital phase correction. The phase drifts later in the year relative to the Gregorian calendar with leap years introducing corrective jumps (Fig. A2d). The sign and magnitude of the induced trends depends on which interval is selected. The net effect of calendrical misfit to the tropical year over the 1951–2010 period considered here is that the phase of solar insolation has shifted 0.5 days earlier on the calendar.