Trends in Winter Warm Spells in the Central England Temperature Record

An important impact of climate change on agriculture and the sustainability of ecosystems is the in- crease of extended warm spells during winter. We apply crossing theory to the central England temperature time series of winter daily maximum temperatures to quantify how increased occurrence of higher temperatures translates into more frequent, longer-lasting, and more intense winter warm spells. We ﬁnd since the late 1800s an overall two- to threefold increase in the frequency and duration of winter warm spells. A winter warm spell of 5 days in duration with daytime maxima above 13 8 C has a return period that was often over 5 years but now is consistently below 4 years. Weeklong warm intervals that return on average every 5 years now consistently exceed ; 13 8 C. The observed changes in the temporal pattern of environmental variability will affect the phenology of ecological processes and the structure and functioning of ecosystems.


Introduction
Alongside the impact (IPCC 2014) of the change in intensity and frequency of the observed summer heat waves and drought, agriculture and the sustainability of ecosystems are affected by increasing occurrence of extended warm spells during winter (Grimm et al. 2013). Experiments in temperate zones show that warming affects soil biochemistry and in particular the nitrogen cycle (Butler et al. 2012). Winter warming can play a distinct role in biochemistry and primary production (Hutchison and Henry 2010) and also impacts ecosystems via changes in the freeze-thaw cycle (Joseph and Henry 2008). In temperate zones, winter warming can lead to ecosystem shifts (Schuerings et al. 2014) and Denotes content that is immediately available upon publication as open access.
Winter warming affects the development of species and the dynamics of populations, resulting in changes in phenology and ecological interactions during spring, with consequences for the structure and functioning of ecosystems. Unusually extended warmer and colder periods can affect development and survival rates of organisms, the composition of communities, and the structure and functioning of ecosystems (Jiguet et al. 2011;Ma et al. 2015). There is, however, little understanding of how climate variability affects species' life cycles and community interactions. Improving understanding of the ecological consequences of extreme climatic events, such as extreme hot and dry, and cold and wet, summers and winters is an increasing research focus (Felton and Smith 2017). A series of recent studies have examined the impacts of these climatic extremes, in the sense of atypical extended warmer and colder periods, on individual plant and animal species, communities, and ecosystems and on agricultural systems. The temporal structure of temperature variability is important in determining the magnitude of any ecological effects. Experimental analyses showed that changes in temporal clustering of warm and cool periods can affect demographic rates and population dynamics of an agricultural pest species (Ma et al. 2018).
Phenological change can be a driver of population trends (Bell et al. 2019). Milder and shorter winters delay the phenology of temperate butterfly species (Stålhandske et al. 2017). Extreme winter warming events in subarctic regions can result in snow and ice melt, exposing terrestrial ecosystems to high air temperatures. This can result in shifts in community composition, affecting ecosystem structure and functioning (Bokhorst et al. 2012). Warm spells during winter can lead to breaks of dormancy of woody species, resulting in advances in bud development (Ladwig et al. 2019). Such development followed by a return to cooler winter conditions can lead to damage of young tissue, affecting their phenology, which can disrupt ecological relationships during spring. Warm spells in winter are generally unfavorable for butterflies. Periods of anomalous heat during overwintering of butterflies have a negative impact on populations (Long et al. 2017). The potential consequences of such changes are not understood, but they are likely to be complex and far reaching, resulting in shifts in phenology affecting the structure and functioning of ecological systems and the services they provide (e.g., affecting water cycles or food security) (Felton and Smith 2017). Persistent warm spells will on average be more frequent and longer lasting when there is an upward drift in the mean and/or the higher moments of a given temperature distribution (IPCC 2012), and to understand them requires quantitative knowledge of how the full distribution is changing. Regional climate patterns drive persistent warm spells and cold snaps; indeed, in the same winter seasons, extreme cold spells have been attributed to the state of regional climate patterns, while increased warm spells were found to be consistent with a long-term overall warming trend (Guirguis et al. 2011). Furthermore, the increased intensity of winter warm spells in one region has also been linked to increased intensity of cold spells in another (Cohen et al. 2018).
Heat waves and warm spells are by definition rare events, so their full distribution will be difficult to determine, requiring ensembles that are only obtainable from model output; however, these also can have significant intermodel variability (Alvarez-Castro et al. 2019). For a single observed time series, provided the overall trend in the high quantiles (the distribution tail) is discernible from shorter-time-scale variability, crossing theory can, however, extract the trends in average warm-spell properties that is model independent. Chapman et al. (2019) recently made what we believe to be the first application of crossing theory to surface temperature time series. Definitions of spells of extreme temperatures can vary (Frich et al. 2002), but if heat waves or warm spells can be usefully defined as runs of consecutive days where the daily maximum temperature exceeds a constant threshold value, we can in principle use crossing theory (Lawrance and Kottegoda 1977;Vanmarcke 2010) to quantify directly from the observations how long-term trends in the higher values in the observed distribution translate into trends in the average properties of heat waves or warm spells. In practice, the utility of this approach depends on the nature of the uncertainties and short-time-scale variability and how these compare with the overall trend in the specific observations of interest. Crucially, this varies with the quantile of the observed distribution; it may provide robust results at some quantiles and not others within the same observed time series. Chapman et al. (2019) focused on heat waves in summer season observations in the central England temperature record (CET). Given the importance of changes in winter warm spells to ecosystem function, we now explore the applicability of crossing theory to the high quantiles of winter season CET observations. We find that a robust trend in winter warm spells does indeed emerge from the analysis. This provides data that could form a direct input into models of ecosystem function and also a consistency check on the outputs of models at high quantiles, which in turn are needed for both attribution and prediction.

a. High quantile time dynamics in the CET
We will analyze an extensively curated dataset, the CET (Parker et al. 1992). The CET daily extrema are reported to a precision estimated to be better than a degree (Parker and Horton 2005) for the record since 1878. The CET has warming trends both in its mean and seasonal extremes (Benner 1999;Brabson and Palutikov 2002). To form a distribution at a specific time, we aggregate the daily observed maximum temperature observations for the winter season, that is, December-February (DJF) over several consecutive years. We will then use distributions constructed at different times to determine how warm-spell properties are changing. Clearly there will be a trade-off between the statistical precision of the estimate of the distribution (which improves with sample size) and the time resolution of any trends estimated by comparing an earlier distribution with a later one (which deteriorates as we aggregate over more successive seasons). The choice made in this analysis is such that the statistical precision of the estimate of the distribution approximately matches the precision of the observations at the distribution quantiles of interest. For the winter season CET, aggregating nine successive years gives a cumulative density function (cdf) with uncertainty [estimated as 95% confidence bounds using the Greenwood (1926) formula] at the 0.95-0.99 quantiles of about a degree. Figure 1 plots the estimates of the cdf and probability density function (pdf) aggregated over 9 successive winter seasons, one for each year that is at the middle of each 9-yr sample from 1882 to 2015. [The online supplemental material contains the same figure constructed for distributions aggregated over different numbers of seasons (Figs. S1-S3) and a plot showing the 0.95-0.99 quantiles and Greenwood uncertainties for the 9-yr aggregates (see Fig. S4).] In Fig. 1 we can see that the shorter-time-scale variability is high at low quantiles, below about 0.25. At high quantiles, the figure shows that above 0.75, this variability is smaller, and at the highest quantiles is within the CET precision (Parker and Horton 2005). This confirms that we can resolve winter season trends in the CET on a multiyear time scale at the 0.9-0.95 quantiles, that is, daily maximum temperatures of ;128-138C, which are typical warm winter day temperatures. This will not necessarily be the case in general, high variability of winter temperatures tend to mask the amplitude of regional temperature anomalies more than is the case in summer (Hansen et al. 2012). The aggregated DJF daily maximum temperatures in the CET from 1878 to 2019 then give samples centered on years 1882-2015.
We will use kernel density estimates to sample the cdfs at every 0.18, the same resolution as the CET time series. Estimates of the cdf uncertainty can then be transformed to give a corresponding uncertainty in the temperature threshold for a given average return period and run length.

b. Crossing theory
Crossing theory (Lawrance and Kottegoda 1977;Cramér and Leadbetter 1967;Vanmarcke 2010) quantifies the properties of runs above or below a threshold directly from the distribution of observations. It has previously found application in hydrology (Nordin and Rosbjerg 1970;Bras and Rodríguez-Iturbe 1993). Most applications of the theory following Rice (1944) require a Gaussian distribution; however, there is one result  (Lawrance and Kottegoda 1977;Vanmarcke 2010) for random variables that does not require any particular form for the observed distribution, as long as we can meaningfully take averages. Within a given time interval the maximum daily temperature record T k on day k is just one of a sequence of samples {T k , k 5 1, 2, . . .} that, following Lawrance and Kottegoda (1977), we treat as a stationary discrete time series. We do not require the T k to follow any particular distribution, and we do not consider any particular structure or correlation in time. A warm spell is a sequence of consecutive days with T . u, and the beginning of a warm spell occurs when the time series crosses a threshold u from below to above, which occurs on day k if T k21 , u and T k . u with probability P(T k21 , u, T k . u). The number of times that the time series goes from below to above the threshold is just the number of warm spells in time interval of M daily observed records [0, M], it is N M (u), and it has expectation value E[N M (u)] 5 P(T k21 , u, T k . u)M. A given data record has probability P(T k . u) of being above the threshold u, and so the expected number of records that are above u in time interval [0, M] is P(T k . u)M. The mean duration of sequences of daily records above the threshold u is just the mean duration of a warm spell, and in days this is t(u) 5 P(T k . u) P(T k21 , u, T k . u) . (1) The expected number of runs or warm spells that occur in time [0, M] also determines their average return period, which is Here M 0 determines the units in which the return period is expressed; we will use years, so that for a winter seasonal average M/M 0 5 90. Now P(T k . u) 5 1 2 C(u), where C(T) is the cdf of daily temperature observations T k . Combining this with Eqs. (1) and (2) gives Warm-spell return periods and durations are not determined independently by this expression because it does not capture the distributions of either of these quantities. It simply translates observed changes in the distribution of daily temperatures to changes in average warm-spell (sequences of consecutive days above a threshold) properties. We previously (Chapman et al. 2019) verified Eq.
(3) using model test time series with different time-correlation properties.

Results
We now apply the above to winter warm spells seen in the CET time series. An example of how the change in the properties of winter warm spells can impact on ecosystem function is desynchronization, where foraging insects such as bees emerge in advance of flowering of their food resources. Negative effects on fitness are seen at a mismatch of 3 days, and these become catastrophic at 6 days (Schenk et al. 2018). In Fig. 2 we use Eq. (3) and the cdfs formed from nine consecutive winter (DJF) seasons plotted in Fig. 1 to show how the average return period of a run of 5 warm days above a given temperature threshold has changed over the last 140 years. Each curve is an estimate based on a different 9-yr cdf. The curves are color coded to indicate the time interval from which each 9-yr sample is drawn. Blue indicates the epoch up to about 1910, green indicates from about 1930 to 1970, and pink indicates about 1980 onward. To generate each of these curves in Fig. 2, we fix the average warm-spell duration t(u) 5 5 days in Eq. (3) and plot the value of return period R(u) versus temperature threshold u. Each curve then just depends on how the cdf of daily maximum temperatures C(u) varies with temperature u. The different curves correspond to the observed cdfs from different 9-yr winter aggregated samples.
In Fig. 2, we can see that the curves move progressively to the right with increasing time. In reading vertically down these plots, that is, at a fixed threshold temperature, it is seen that 5-day warm spells of daily maxima over 128-138C had a return period of ;2-7 years before 1910, and now this is ;0.5-3 years. Reading horizontally across these plots, that is, at a fixed quantile of the (nine winters' aggregate) cdf, one sees that about 140 years ago a 5-day run would have a return period of about a year (the 0.95 quantile) at temperatures around 118C; this has now increased to over 12.58C for some recent 9-yr intervals. The sequence of curves captures how the drift in the temperature distribution with time translates into a change in average warm-spell properties. Warm-spell durations become longer as the distribution drifts toward higher temperatures.
One can instead plot the change in run length at a given return period, this is shown in Fig. 3. To generate each of these curves in Fig. 3 we fix the average return period R(u) 5 5 years in Eq. (3) and plot the value of warm-spell duration t(u) versus temperature threshold u. Again, each curve then just depends on how the cdf of daily maximum temperatures C(u) varies with temperature u. The different curves correspond to the observed cdfs from different 9-yr winter aggregated samples. Reading horizontally across these plots (at the 0.99 quantile), one sees that 5-day-long warm intervals that return on average every 5 years now consistently exceed ;138C; in reading vertically upward, it is seen that a winter warm spell with daytime maxima above 138C with 5-yr return period has a duration that was consistently less than 8 days and is now consistently more than 8 days. Warm-spell return periods become shorter as the distribution drifts toward higher temperatures.
The region of the cdf that we need to resolve, and the accuracy to which we need to resolve it, are determined by the specific threshold(s) that are relevant to how an ecosystem is operating, and the return period and duration of warm spells that are likely to impact that ecosystem. This sets a constraint on how small the sample interval can be that we use to form the cdfs as both the number of daily observations in the sample, and the functional form of the cdf, directly translate to an uncertainty that we can estimate empirically. The plotted uncertainties on our results then derive from those of the observations (Parker and Horton 2005) and those of the estimated cdf. They linearly combine to ;618-28C [supplemental Fig. S4 plots sample cdfs and uncertainty estimated using the Greenwood (1926) formula]. This level of uncertainty is such that we can indeed discern long-term trends in the properties of winter warm spells. To see the detailed time dynamics, we have taken cuts through Figs. 2 and 3 at 128 (blue lines) and 138C (green lines), and we plot these in Figs. 4a and 4b. This approach plots how the return period (Fig. 4a) and duration (Fig. 4b) of winter warm spells have changed, and we can see that both their frequency and duration have increased. To relate these results to a specific ecosystem resilience or agricultural impact, we can instead plot the change in overall winter warm-spell intensity obtained for a specific return period and duration. Figure 4c fixes the return period at 5 years and the warm-spell duration at 5 days and plots the change in the threshold that the maximum daily temperature will exceed in each day of the warm spell. A 5-day warm spell with an average return period of 5 years has a threshold temperature that has increased from being typically below 138C to typically above it. The essential features of Fig. 4 are not strongly sensitive to the interval chosen to aggregate the cdf (see Figs. S5 and S6 in the online supplemental material). Previous studies (King et al. 2015) that use Kolmogorov-Smirnov statistical significance to test against similarity with the quasi-natural temperature distributions, while being successful in models on regional to continental scales, have been inconclusive for the CET. Here, crossing theory reveals a persistent multidecadal enhancement in the warm-spell threshold, which after about 1960 remains above 138C in Fig. 4c, corroborating the overall findings of King et al. (2015).

Conclusions
We have shown that crossing theory can be used to quantify directly from the CET record how long-term trends of the higher values in the observed distribution translate into trends in the average properties of winter warm spells. Our approach depends on an overall longtime-scale trend in the specific observations of interest being discernable from short-time-scale variability as well as statistical uncertainties. We therefore performed these first applications of crossing theory to the CET (from 1878), which is one of the longest temperature records available, initially to summer heat waves in Chapman et al. (2019), and here in this paper to winter warm spells. Multiple-station and spatially gridded daily temperature observations are available since the 1950s, and these can be used to obtain maps of changes in the distribution at high quantiles , which then translate directly into maps of changes in average warm-spell properties. Whether robust signals of change in warm spells can be resolved will, however, have both geographical and quantile dependence (Chapman et al. , 2015. FIG. 4. Winter warm-spell changes over the last 140 years: (a) average return periods for runs of 5 consecutive days with maximum winter daily temperatures above 128 (blue diamonds) and 13 8C (green squares); (b) the average duration of runs of consecutive days with maximum winter daily temperatures above 128 (blue diamonds) and 13 8C (green squares) with average return period of 5 years; (c) the threshold of maximum daily temperature that is exceeded for 5 consecutive days on average every 5 years. Color indicates the sample central year in the time sequence as in Figs. 2 and 3. In (a)-(c), data sampled over nine consecutive winter (DJF) seasons centered on each of the years 1882-2015 are used to form cdfs. Equation (3) relates these cdfs to average return periods and run lengths where the daily maximum temperature is above a given threshold. Gray shading in (c) indicates uncertainties estimated as the larger of that from 618C in the temperature time series (Parker and Horton 2005) and the 95% confidence bounds in the underlying cdf estimated using the Greenwood (1926) formula.
Crossing theory, as applied here to an individual time series, could also in principle be used to compare the performance of model outputs with the observations. With any such comparison it is important to optimize how aggregation across spatial grid points and over a time window is performed to obtain samples of the timevarying distribution that are both statistically significant at the quantiles of interest and are not overly coarse grained such that detailed trends are ''washed out.'' This optimization will vary geographically since individual spatially localized distributions are known from the observations to vary substantially both in their time dynamics and its uncertainty Stainforth et al. 2013).
Ecological studies have highlighted the need for improved resolution of long-term monitoring and understanding of ecosystem-level effects of previously atypical climatic events (Ladwig et al. 2019;Long et al. 2017;Friedl et al. 2014). The importance of developing appropriate definitions and measures of untypical runs of high temperatures for ecological analyses has also been highlighted (Friedl et al. 2014;Bailey and van de Pol 2016). The current study provides a useful basis for quantifying temporal environmental variability and defining these events. This study has also shown that it is not just the frequency but also the duration of winter warm spells that has increased. Improving understanding of the ecological responses to such changes in the temporal structure of environmental variability will be crucial for generating understanding of the potential impacts of future climate-driven change (Ma et al. 2018).
Warm spells are defined in a climatic sense as sustained excursions above a fixed temperature threshold. Ecological and agricultural systems do not in general respond in a binary manner to the sustained crossing of a temperature threshold; their response is far more varied and complex (Williams et al. 2015). Winter warm spells can contribute to ecosystem responses that are nonlinear, such as outbreaks. These can be triggered when a winter warm spell favors reproduction that is just sufficient to move a population above its outbreak threshold. Increase in outbreak frequency has occurred, for example, in bark beetles across European forests with both ecosystemwide and economic impact (Hlásny et al. 2019).