The influence of modes of climate variability on the sub-seasonal temporal clustering of extreme precipitation.

Temporally clustered precipitation extremes can have catastrophic impacts. Therefore, understanding their drivers is paramount for risk assessment in current and future climates. Here, we model for each season 3-week extreme precipitation event counts with Poisson Generalized Linear Models and nine major modes of climate variability as covariates. Model goodness-of-fit is highest in the tropics, particularly over the equatorial Pacific, the Maritime Continent, and East Africa, where ENSO, the Indian Ocean Dipole (IOD) and the MJO are the major drivers of sub-seasonal temporal clustering of extreme precipitation. The IOD and MJO also matter over Southwest Asia during boreal fall and winter. In the Northern Hemisphere, the North Atlantic Oscillation impacts clustering west of the Iberian Peninsula and over Scandinavia and Greenland, and the Pacific North American pattern matters over the central/northern Pacific Ocean. Finally, our models show very little skill in the Southern Hemisphere, where temporal clustering is also less frequent.

Introduction

The impacts of extreme events on human and natural systems can be compounded if the extreme events occur in close temporal succession, in which case one speaks of ”temporally compounding” or ”temporally clustered” events (Zscheischler et al., 2020). Precipitation is a case in point: extreme precipitation recurring over weekly to monthly timescales over the same area can significantly increase flood risk (Barton et al., 2016; Tuel and Martius, 2021b). In addition, Tuel and Martius (2021a) recently showed that the temporal clustering of extreme precipitation (TCEP) over sub-seasonal (weekly to monthly) timescales was statistically significant in many regions of the world, particularly the tropics and along the eastern margins of the ocean basins. Taking into account the tendency of extreme precipitation to cluster in time is therefore paramount in climate risk assessment and of interest for sub-seasonal predictions. Forecasting of TCEP remains challenging, notably because we still know little about their physical drivers (Zscheischler et al., 2020). We do know that the likelihood of extreme precipitation frequency is strongly influenced by atmospheric circulation and in particular by large-scale modes of atmospheric variability at sub-seasonal to seasonal timescales (Kenyon and Hegerl, 2010). Among these, the El Niño-Southern Oscillation (ENSO) affects extreme precipitation probability in many parts of the world (e.g., Alexander et al., 2009; Kenyon and Hegerl, 2010). In its positive phase (El Niño), ENSO is associated with more frequent extreme precipitation over southeastern Brazil (Kenyon and Hegerl, 2010) and the western and southern United States (e.g., Schubert et al., 2008; Shang et al., 2011), and less frequent over Southern Africa (Hoell and Cheng, 2018), Southeastern Asia, and Australia (e.g., Min et al., 2013; Jakob and Walland, 2016). Similarly, the Indian Ocean Dipole (IOD), an oscillation of SST anomalies between the western and eastern Indian Ocean, increases extreme precipitation likelihood in Southern Africa (Hoell and Cheng, 2018), East Africa (Wainwright et al., 2021), Southwest Asia (Tuel et al., 2021), and India (Krishnan et al., 2011), and decreases it over eastern Australia (e.g., Min et al., 2013). The occurrence of precipitation extremes in the tropics and subtropics also varies substantially depending on the phase of the Madden-Julian Oscillation (MJO), the dominant mode of sub-seasonal atmospheric variability in the tropics (Shimizu et al., 2017; Schreck, 2021). In Europe, the North Atlantic Oscillation (NAO), the main mode of variability in the location of the North Atlantic storm track, is a major influence on the probability of extreme precipitation during the winter half year (e.g., Haylock and Goodess, 2004; Scaife et al., 2008). In its positive phase, the storm track is displaced polewards, and consequently precipitation extremes are more frequent in Northern Europe, especially Scandinavia, and less frequent in Southern Europe and Northwestern Africa. By contrast, a negative NAO significantly enhances the likelihood of extreme precipitation in Southern Europe, particularly the Iberian Peninsula (Merino et al., 2016). Its hemispheric counterpart, the Arctic Oscillation, or Northern Annular Mode, also impacts extreme precipitation across the Northern Hemisphere (Kenyon and Hegerl, 2010), notably in Asia (Mao et al., 2011; Liu et al., 2019). In North America and Eastern Asia, extreme precipitation probability is also sensitive to Northern Pacific interannual and interdecadal modes of variability like the Pacific Decadal Oscillation (PDO), Western Pacific (WP), or Pacific North American (PNA) pattern (Zhang et al., 2010; Kenyon and Hegerl, 2010; Zarekarizi et al., 2018; Wei et al., 2021). Finally in the Southern Hemisphere, the Southern Annular Mode (SAM) and the Baroclinic Annular Mode (BAM) modulate the incidence of extreme precipitation over southeastern Brazil (Vasconcellos and Cavalcanti, 2010) and Eastern Australia (Min et al., 2013) during the austral summer, as well as in Antarctica (Marshall et al., 2017). These results focus on the occurrence of single extreme events. However, the likelihood of clustered extremes varies with large-scale modes of variability has been a subject of much less attention. It is not straightforward to extend results from individual to clustered extremes, since the temporal persistence of the driver is an important factor in explaining temporal clustering. Villarini et al. (2011) applied Poisson regression to analyze the clustering of extreme precipitation at seasonal scales in the Midwestern US in connection to indices of Atlantic and Pacific variability. Mallakpour et al. (2017) used Cox regression to highlight the influence of the AO and PNA on TCEP in the central United States. Yang and Villarini (2019) took the same approach over Europe, focusing on the AO and NAO. Recently, Barton et al. (2021) fitted Poisson generalized additive models to 3-week extreme precipitation counts over Europe, using principal components of regional 500 hPa geopotential anomalies as covariates. They highlighted the effect on TCEP frequency of major weather regimes, in particular blocking over Greenland. The links between midlatitude cyclone clustering and large-scale climate modes like the NAO and PNA were also investigated in several studies and summarized by Dacre and Pinto (2020). Previous studies were limited to 3–4 modes of variability and took a regional perspective; however, we adopt a global perspective to quantify the role of nine major large-scale modes of climate variability on TCEP likelihood at sub-seasonal timescales. We implement Poisson regressions on three-week extreme event counts with the modes of variability as covariates, with each season analyzed separately (DJF, MAM, JJA, and SON). We also contrast our results with those of Tuel and Martius (2021a) and find that the spatiotemporal patterns of TCEP statistical significance can be mostly explained by our selected large-scale modes of variability.

Results

Figure 1 shows for each season the resulting deviance ratio (DR), a measure of model goodness-of-fit (see STAR Methods). Model goodness-of-fit in all seasons is highest in the tropics (20◦ S-20◦ N), where 7–16% of grid points exhibit DR values larger than 0.2, primarily in the equatorial Pacific east of the date line. Other tropical regions with high model goodness-of-fit include the Maritime Continent in JJA and SON, East Africa in SON and DJF and the subtropical central Pacific Ocean in DJF and MAM. Outside of the tropics, DR values are overall much lower, with only 0.5–4% of grid points where DR 0.2. Model goodness-of-fit is essentially zero southwards of 20◦ S, except over Eastern Australia in SON. In the Northern Hemisphere extratropics, DR rarely exceeds 0.4. DR values between 0.2 and 0.4 can be found around Southwest Asia in SON, in the Sahara in DJF, in Alaska and Northeastern Russia in DJF and SON, west of Morocco and the Iberian Peninsula in DJF, and over Greenland in most seasons. We also indicated on Figure 1 the grid points for which temporal clustering of extreme precipitation was statistically significant at time windows of 10–30 days (results from Tuel and Martius (2021a)). A 0.2 threshold for DR matches quite well with areas of clustering statistical significance; these areas are also those of high dispersion statistics (Figure S1). We note also that almost all grid points pass the Poisson goodness-of-fit test at the 5% confidence level (Figure 2), including in regions of high model goodness-of-fit. Only a handful of grid points along the eastern equatorial Pacific do not pass the test. The reason may be that our 2-day declustering window is too short near the equator and thus we do not remove all the temporal dependence in extreme precipitation counts. Still, these points are a minority, and even when excluding them from the analysis, our conclusions would remain unchanged. At most grid points the optimal model includes two or fewer significant covariates, the exception being East Africa in DJF where we find 3–4 significant covariates (Figure 2); otherwise, most regions where DR 0.2 are characterized by a single significant covariate. To understand the role played by the various covariates, we color on Figure 3 each grid point depending on which of the covariates has the largest absolute coefficient. A clear spatial structure emerges in all seasons. The Niño3.4 index is the dominant (and often only significant, see Figure 2) covariate for most of the tropical Pacific grid points with DR 0.2. The Niño3.4 index is also relevant over large sections of the Maritime Continent and Southeast Asia. In the subtropics, Niño3.4 is also the dominant covariate in Southwest Asia in SON, Western Australia in DJF, Southern Brazil in MAM, and significant portions of the Gulf of Mexico and Caribbean in all seasons except JJA (Figure 2); besides, of these regions exhibit small (0–0.2) DR values. Another important covariate in the tropics is the IOD: it affects clustering in East Africa and the Sahel in DJF, over the tropical Indian Ocean in JJA and over Southwest Asia, East Africa, Indonesia, and Eastern Australia in SON (Figure 2). Finally, the two MJO components also have a major influence on TCEP distribution across much of the tropics, especially over the Arabian Peninsula and the Maritime Continent in DJF and MAM for MJO1 (seeSTAR methods section for details), and from East Africa to the date line and over the tropical Americas in JJA for MJO2. In the Northern Hemisphere extratropics, the PNA appears as the dominant influence on TCEP frequency during DJF and MAM in the North Pacific, particularly over Northeastern Asia, Alaska, and the central Pacific Ocean (Figure 3). The PNA also matters locally in North America and the Caribbean in DJF. The influence of the WP additionally stands out around Alaska and Northeastern Russia in all seasons except JJA. In the North Atlantic, NAO is the prevailing covariate over Greenland and Southwestern Europe, two regions of high model goodness-of-fit (Figure 2), chiefly during DJF. The AO index dominates during DJF and MAM in high-latitude regions of North America and Eurasia, as well as around the Black Sea. We also find the AO influence to be the largest over a large region of the eastern subtropical Atlantic during JJA (Figure 3). Finally, we detect very few coherent signals in the Southern Hemisphere, even for the two Southern Hemisphere modes, SAM and BAM.

Figure 1

GLM deviance ratio (DR) for each season: DJF (topleft), MAM (topright), JJA (bottomleft), and SON (bottomright)

Grid points where temporal clustering of extreme precipitation at the 10–30 days timescale is statistically significant (results of Tuel and Martius (2021a)) are indicated by black dots.

Figure 2

Number of significant covariates in the optimal GLM model for each season: DJF (topleft), MAM (topright), JJA (bottomleft), and SON (bottomright)

Grid points where the Poisson goodness-of-fit test is rejected at the 5% level of confidence are shown in green.

Figure 3

Covariate with the largest coefficient, in absolute terms, among all optimal model covariates, for each season

DJF (topleft), MAM (topright), JJA (bottomleft), and SON (bottomright). Grid points where DR 0.15 are masked.

Discussion

Physical perspectives on TCEP drivers

How the selected modes of variability correlate to TCEP frequency is generally consistent with their known effect on the likelihood of extreme precipitation. In the Northern Hemisphere, the influence of the NAO is mainly felt during winter over the North Atlantic sector, in particular over the Atlantic Ocean directly west of the Iberian Peninsula and above Greenland and Scandinavia (Figures 2 and 4). Persistent positive NAO states divert extratropical cyclones toward Scandinavia and away from western Greenland and Southwestern Europe and toward Northern Europe, whereas persistent negative NAO states lead to enhanced Rossby wave breaking and cyclone development around the Iberian Peninsula (Dacre and Pinto, 2020). In the North Pacific, the effects of PNA variability are mainly negative (Figure 4). The positive PNA phase, characterized by a pronounced trough over the central North Pacific, is associated with an equatorward displaced stormtrack and enhanced cyclone clustering in the Gulf of Alaska (Bunker, 2018). The negative PNA phase, by contrast, exhibits a strong ridge over the central North Pacific resembling a Rex-type block and a poleward shifted stormtrack, leading to more frequent TCEP at high latitudes. The significant increase in TCEP probability over the central subtropical Pacific during negative PNA events may result from enhanced Rossby wave breaking downstream of the anomalous North Pacific ridge (Martius et al., 2007), as well as from positive subtropical SST anomalies associated with a negative PNA (Sheng, 1999). The weak model goodness-of-fit in the Southern Hemisphere extratropics and the lack of a significant influence from either the SAM or the BAM are in agreement with the quasi-absence of significant clustering in precipitation extremes at sub-seasonal timescales highlighted by Tuel and Martius (2021a) (Figure 1). It does not mean that the SAM and BAM indices have no impact on extreme precipitation frequency, but that they are not relevant for TCEP likelihood on sub-seasonal timescales. Moving on to the tropics, we find that ENSO exerts a strong influence on tropical TCEP frequency, increasing it in the equatorial Pacific and decreasing it over the Maritime Continent during an El Niño phase, mainly during SON and DJF (Curtis et al., 2007) (Figure 4). In DJF, a positive ENSO similarly tends to suppress clustered precipitation extremes over southwestern Australia and enhance them southern North America (Figure 4, topleft), in keeping with previous findings focusing on extreme precipitation frequency only (e.g., Kenyon and Hegerl, 2010; Shang et al., 2011; Sun et al., 2015). ENSO’s weak positive influence in southern Brazil during MAM and Southeastern China in DJF were also noted by Sun et al. (2015). The IOD primarily affects TCEP frequency in the Indo-Pacific region in SON and DJF (Figures 3 and 5). In SON, a positive IOD increases TCEP frequency over Southwest Asia and decreases it in Eastern Australia, consistent with IOD impact on single extremes described by Tuel et al. (2021) and Min et al. (2013). It also strongly enhanced extreme precipitation likelihood from September to February over East Africa, at a time when prevailing winds over the Indian Ocean are usually easterly and advect moisture away from the African continent (Wainwright et al., 2021). A positive IOD, characterized by warmer than average SSTs in the Western Indian Ocean, thus leads to anomalous westerlies that concentrate moisture over East Africa and to extreme precipitation (Wainwright et al., 2021). Finally, we turn to the effects of MJO variability. Taken together, the two MJO proxies (MJO1 and MJO2) significantly modulate TCEP frequency across the tropics for much of the year, but especially over the Maritime Continent (Figure 7). Their influence is coherent with their large-scale climate footprint discussed in Wheeler and Hendon (2004). The positive MJO1 phase (MJO phases 4–5) is linked to enhanced convection, low-level moisture convergence and TCEP frequency over the Maritime Continent, and the opposite over the equatorial Pacific and Indian Oceans (MJO phases 8-1). The MJO2 index, by contrast, correlates positively to convection, low-level moisture convergence and TCEP frequency in the tropical Western Pacific, and negatively over the Maritime Continent. A positive MJO2 index corresponds to MJO phases 6–7, and a negative index to MJO phases 2–3 (Wheeler and Hendon, 2004). Overall, MJO influence is very limited in the extratropics. The main exception is the weak, though significant, impact of MJO1 over Southwest Asia and North Africa in SON and DJF (Figures 2 and 7). Therefore, negative MJO1 values (MJO phases 8-1) are associated with more frequent TCEP and increased extreme precipitation frequency also Schreck (2021), possibly related to enhanced upper-level divergence over Northern Africa in the corresponding MJO phases. Interestingly, most MJO-dominated regions do not exhibit clustering significance at the chosen sub-seasonal timescale. It could be explained by weak MJO persistence over sub-seasonal timescales in space (at the grid point level), even though the MJO phase may persist in time (Lafleur et al., 2015).

Figure 4

Regression coefficients in the Poisson GLM model during DJF for the (top) NAO and (bottom) PNA indices

Only grid points where the coefficients are statistically significant at the 5% confidence level (modified to take the false discovery rate into account) are shown.

Figure 5

Regression coefficients in the Poisson GLM model for the Niño 3.4 index in all seasons

Only grid points where the coefficients are statistically significant at the 5% confidence level (modified to take the false discovery rate into account) are shown.

Figure 7

Regression coefficients in the Poisson GLM model for the MJO1 and MJO2 indices in all seasons

Only grid points where the coefficients are statistically significant at the 5% confidence level (modified to take the false discovery rate into account) are shown.

Regression coefficients in the Poisson GLM model for the IOD index in all seasons

Link to statistical significance of clustering

It is interesting to contrast the results of the present study with the seasonal TCEP statistical significance maps shown by Tuel and Martius (2021a) (added here in Figure 1). They used Ripley’s K function to assess the statistical significance of TCEP, i.e., the tendency of extreme precipitation events to be more clustered in time than what a random Poisson process would imply. The Ripley’s K approach of Tuel and Martius (2021a) and the Poisson GLM we adopt here are complementary. Ripley’s K is a purely statistical approach; it quantifies how frequent TCEP is and to what extent TCEP frequency deviates from that of a time-independent process, regardless of any potential drivers. By contrast, the Poisson GLM relates TCEP frequency to selected physical drivers, regardless of how often TCEP occurs in the precipitation time series. Tuel and Martius (2021a) found TCEP to be statistically significant over many regions of the globe that often coincide with those of high model goodness-of-fit (Figure 1). They identified clustering significance in particular over the equatorial Pacific, East Africa, the Maritime Continent, Southwest Asia, Southwestern Europe, Australia, Alaska, and Northeastern Russia in seasons where the DR values also tend to be high (Figure 2). The spatiotemporal pattern of clustering significance in the equatorial Pacific (Figure 1) notably coincides well with that of the Niño 3.4-dominated region of Figure 3. A handful of regions exhibiting clustering significance do not stand out in the present analysis–chiefly the tropical Atlantic in DJF and MAM, and Western Australia in SON–where local processes may be more important, or where we may have missed the relevant modes of variability. Nevertheless, it appears that our reduced selection of large-scale modes of variability is sufficient to explain much of the clustering significance at the 3-week timescale. Clustering significance seems to be primarily linked to four modes of variability: the IOD and Niño 3.4 indices in the tropics, and the PNA and NAO in the Northern Hemisphere. The IOD and Niño 3.4 both represent long-lived SST anomalies; therefore, it is not surprising that they would influence weather on sub-seasonal timescales. The NAO and PNA also exhibit persistence at sub-seasonal timescales (Barnes and Hartmann, 2010; Li et al., 2019). By contrast, few regions where TCEP frequency is dominated by the MJO exhibit clustering significance.

Limitations of the study

Because we regress (3-week extreme precipitation counts against concurrent values of large-scale modes of variability, some of the influence of these modes on TCEP frequency may end up diluted in our approach, in case their effect takes some time to be felt remotely. The extratropical influence of ENSO or the MJO, for instance, results from Rossby wave trains that typically take 1–2 weeks to travel from their source (Cassou, 2008). In addition, our GLM model does not take into account potential nonlinearities in covariate effects, which have been documented e.g., in the case of ENSO (Sun et al., 2015). Poisson GLMs are nevertheless useful to understand, at first order, the global patterns of clustering significance and its major climate drivers. Still, detailed regional analyses could rely on more complex approaches like Generalized Additive Models (Barton et al., 2021).

Conclusions

Based on a simple regression approach of 3-week extreme precipitation counts, we investigate the seasonal influence of eight major modes of climate variability on the frequency of temporally clustered precipitation extremes (TCEP) across the globe. Removing the seasonal cycle in extreme precipitation magnitude allows us to concentrate on sub-seasonal TCEP drivers. Model goodness-of-fit is highest in the tropics, particularly over the equatorial Pacific, the Maritime Continent, and East Africa, but also in a handful of extratropical regions like Southwest Asia, Alaska, Northeastern Russia, and Southwest Europe. These regions largely correspond to those where temporal clustering in extreme precipitation time series is statistically significant. In the tropics, ENSO and the IOD have the most influence on TCEP frequency, and their persistence in time likely explains most of the clustering significance in the tropics. In the tropics still, the MJO also modulates TCEP frequency, especially over the Maritime Continent. There is little significant TCEP in the Southern Hemisphere extratropics. In the Northern Hemisphere, the NAO, AO, WP, and PNA all significantly affect TCEP frequency, particularly over the North Atlantic and Pacific basins. Although several of the selected modes of variability remain challenging to forecast even at sub-seasonal timescales (e.g., Feng et al., 2021), others such as ENSO or the IOD, are routinely well predicted by weather models (e.g., Jin et al., 2008). Our results thus open the door for improved forecasts of sub-seasonal TCEP frequency and associated impacts.

STAR★methods

Key resources table

Resource availability

Lead contact

Further information and requests for resources should be directed to the lead contact, A. Tuel (alexandre.tuel@giub.unibe.ch).

Materials availability

This study did not generate new datasets.

Data and code availability

This paper analyzes existing, publicly available data. These accession numbers for the datasets are listed in the key resources table

All original code has been deposited on Github and is publicly available as of the date of publication. Links are listed in the key resources table.

Any additional information required to reanalyze the data reported in this paper is available from the lead contact upon request.

Method details

Precipitation

We use daily precipitation data at one ◦ resolution from the ERA5 reanalysis (Hersbach et al., 2020) over the 1979–2019 period. While precipitation in ERA5 is a prognostic variable and precipitation observations are not directly assimilated into the model, Rivoire et al. (2021) showed that ERA5 was able to correctly capture the timing of extreme precipitation events in the extratropics, though less so in the tropics. Tuel and Martius (2021a) additionally showed that ERA5 compared well to several satellite-and observation-based products in terms of TCEP statistics.

Modes of climate variability

We consider here nine canonical modes of climate variability: the El Niño-Southern Oscillation (ENSO), the Indian Ocean Dipole (IOD), the Madden-Julian Oscillation (MJO), the North Atlantic Oscillation (NAO) and Arctic Oscillation (AO), the Pacific North American pattern (PNA), the Western Pacific pattern (WP), the Southern Annular Mode (SAM) and the Baroclinic Annular Mode (BAM). Daily time series of the NAO, AO, SAM and PNA indices are obtained from NOAA’s Climate Prediction Center at https://www.cpc.ncep.noaa.gov/products/precip/CWlink/daily_ao_index/history/history.shtml. Details of the calculations are provided at https://www.cpc.ncep.noaa.gov/products/precip/CWlink/daily_ao_index/history/method.shtml. We define the WP index following Wallace and Gutzler (1981) as the normalised difference in 500 hPa geopotential height between 60◦ N/155◦ E and 30◦ N/155◦ E. Normalisation is performed with the mean and SD of a 30-day moving window centered on each day of the year. The geopotential data is taken from ERA5. We calculate daily Niño 3.4 and IOD indices using SST data from ERA5. The Niño 3.4 index is defined as the mean SST anomaly in the 5◦ S-5◦ N/170◦ W120◦ W region, and the IOD index as the anomaly in SST gradient between the western equatorial (10◦ S-10◦ N/50◦ E-70◦ E) and southeastern equatorial (10◦ S-0◦/90◦ E-110◦ E) Indian Ocean. Daily Niño 3.4 and IOD indices are not especially relevant, but we subsequently average them over sub-seasonal timescales (see section Precipitation extremes and sub-seasonal averaging). For the MJO, we use the two real-time multivariate indices of Wheeler and Hendon (2004), denoted MJO1 and MJO2, available from the Australian Bureau of Meteorology at http://www.bom.gov.au/climate/mjo/. Finally, we compute a daily BAM index following Thompson and Woodworth (2014) as the time series of the leading principal component of eddy kinetic energy between 70 ◦ and 20◦ S and 1000–200 hPa (after weighting the data according to latitude and pressure). Eddy kinetic energy is calculated from ERA5 6-hourly wind on all available pressure levels; the long-term monthly-means are subtracted from the wind components at all grid points and pressure levels.

Precipitation extremes and sub-seasonal averaging

Following Tuel and Martius (2021a), we define precipitation extremes at each grid point as days when daily precipitation exceeds its 99th all-days percentile of the corresponding month over the full 1979–2019 period. This approach allows to remove the seasonality in precipitation magnitude to concentrate on sub-seasonal timescales, and to have similar extreme precipitation frequencies in all seasons. Note however that the absolute extreme precipitation values vary considerably across location and season (see Tuel and Martius (2021a) for more details). The choice of the 99th percentile as threshold is a compromise to select reasonably extreme events (3–4 per year) while retaining a high enough frequency of TCEP for the statistical analysis. Because individual weather systems responsible for extreme precipitation can linger several days over a given location, we implement at each grid point a runs declustering method (Coles, 2001), in which extreme precipitation events separated by less than 2 days are grouped together and assigned to the same event. The resulting extreme event frequency is shown on Figure S2. Extreme event counts are then summed over successive, non-overlapping 3-week windows. This choice of summing window is appropriate to quantify TCEP at sub-seasonal timescales (Barton et al., 2016; Tuel and Martius, 2021a). Our results remain the same for two or 4-week windows (compare Figures 3 and S3). We also average daily climate mode indices over the same 3-week windows, and assign each window to the season with which it overlaps most: winter (DJF), spring (MAM), summer (JJA) and fall (SON). We use non-overlapping windows to ensure that we fit the model on independent time periods.

Poisson regression

We model 3-week extreme event counts n, at each grid point and for each season separately, with Poisson Generalized Linear Models (GLMs). n is assumed to follow a Poisson distribution with rate of occurrence λ function of covariates as follows:

The covariates are the rescaled 3-week averaged time series of climate modes introduced in section Modes of climate variability. The rescaling consists in transforming each covariate series into a series with mean 0 and SD 1. We carry out a covariate selection with a forward stepwise approach. Starting from a model with no covariates, we add covariates, one at a time, to a list of ”best” covariates leading to an ”optimal” model. Model performance is determined by the Akaike Information Criterion (AIC), equal to the difference between the number of model covariates and the model log likelihood. A smaller AIC indicates a better model. At each step, each covariate, among all those not yet included in the model, is added to the model. The covariate leading to the largest drop in AIC is added to the list of ”best” covariates, and the sequence continues as long as AIC can decrease by adding further covariates. We also checked that a random shuffling of the covariates at the beginning of the stepwise covariate selection did not impact the results. We assess the goodness-of-fit of the optimal model with its deviance ratio (DR) (McCullagh and Nelder, 1989). Model deviance D is equal to , where is the log likelihood of the ”perfect” (or saturated) model (one free parameter per observation) and the log likelihood of the fitted model. Letting D be the deviance of the ”null model” (in which λt = λ0), DR is then defined as:

DR is therefore a measure of model goodness-of-fit, since it quantifies the relative improvement (in terms of deviance) of the model with respect to the null model. DR close to 0 means that the model performs hardly better than the null one, while DR close to one implies an almost perfect model. For example, DR = 0.2 means that the model offers a 20% improvement on the null model. The validity of the Poisson distribution assumption is assessed with a chi-square test on the fitted model deviance D. Under the null hypothesis that the assumed Poisson model is correct, D follows a distribution with n degrees of freedom, where n equals the number of observations minus the number of model parameters. The Poisson model (1) implies that for a unit increase in covariate , λt (the average expected number of extremes, ) is multiplied by . We determine the significance of the regression coefficients (Equation 1) by filtering the corresponding p values with a 5% significance threshold that takes the false discovery rate into account (Wilks, 2016). We retain significant coefficients only before determining, at each grid point, the number of selected model covariates and the covariate with the largest coefficient (in absolute terms).

REAGENT or RESOURCE	SOURCE	IDENTIFIER
Deposited data
ERA5 reanalysis data	ECMWF	https://dx.doi.org/10.24381/cds. adbb2d47
Daily NAO, AO, SAM and PNA indices	National Centers for Environmental Prediction	https://www.cpc.ncep.noaa.gov/products/precip/CWlink/daily_ao_index/history/history.shtml
Software and algorithms
Code to reproduce the results of this study	Authors	https://github.com/Quriosity129/TCEP_modes_variability