Heat extremes in Western Europe increasing faster than simulated due to atmospheric circulation trends

Domain: earth_science
Published: 2023-10-26
Editorial outcome: major_revision

Abstract
Over the last 70 years, extreme heat has been increasing at a disproportionate rate in Western Europe, compared to climate model simulations. This mismatch is not well understood. Here, we show that a substantial fraction (0.8°C [0.2°−1.4°C] of 3.4°C per global warming degree) of the heat extremes trend is induced by atmospheric circulation changes, through more frequent southerly flows over Western Europe. In the 170 available simulations from 32 different models that we analyzed, including 3 large model ensembles, none have a circulation-induced heat trend as large as observed. This can be due to underestimated circulation response to external forcing, or to a systematic underestimation of low-frequency variability, or both. The former implies that future projections are too conservative, the latter that we are left with deep uncertainty regarding the pace of future summer heat in Europe. This calls for caution when interpreting climate projections of heat extremes over Western Europe, in view of adaptation to heat waves.

Introduction
Extreme heat has been increasing at global scale1,2, with a rapid rate in several regions. In Western Europe3, summer temperatures and heat extremes have warmed much faster than elsewhere in the mid-latitudes over the last two decades3,4. As a consequence, several unprecedented heat waves took place in the last 20 years. In 2003, the full summer season mean temperature was unprecedented in Europe5. Northwestern Europe was hit by record temperatures in 20186,7. In 2019, two short (3-day) but intense heat waves saw all-time temperature records broken in many places, associated with a rapid northward advection of Saharan air6. All-time records were broken again in 2022, with temperatures above 40°C reaching far north (eg. Brittany, U.K.) (https://www.worldweatherattribution.org/without-human-caused-climate-change-temperatures-of-40c-in-the-uk-would-have-been-extremely-unlikely/, (2022)). Unprecedented, and even record-shattering extremes are plausible in climate projections8, but the pace of their increasing magnitude in Western Europe is generally not predicted by these climate models, as well as trends in mean summer temperatures4,9–12.

Here we focus on summer (JJA) maximum and mean of daily maximal temperatures (resp. denoted hereafter TXx and TXm for simplicity), and the regional amplification of their trends relative to the global temperature trend. Trends in TXx and TXm are calculated over the 73-year 1950–2022 period using a linear regression with the Global mean Surface Air Temperature (GSAT, see methods section) from ERA5, and are expressed in °C per global warming degree (GWD). As shown in Fig. 1 and Supplementary Fig. 1, both ERA5 reanalyses13 and E-OBS interpolated observations14 exhibit trends reaching more than 5°C/GWD for TXx in northern France and Benelux. Over the limited areas spanning 5W-15E; 45N-55N (blue box, called hereafter “Western Europe”), the land area-average TXx trend is 3.4°C/GWD for ERA5 and E-OBS [2.4–4.3°C/GWD]. It exceeds the more moderate TXm trends by about 40% for ERA5 (2.4°C/GWD [1.7–3.0°C/GWD] and 30% for E-OBS (2.6°C/GWD [1.9–3.3°C/GWD]). These rapid warming trends are exceptional on a global scale: The 20°×10° Western Europe region has the highest TXx (all year round) trend of all regions of the same size around the globe between 75°S and 75°N shifted by steps of 5° (including sea points).

A variety of processes have been proposed for explaining these over proportional warming trends with respect to global temperature change. For mean summer temperatures, changes in mean atmospheric circulation15,16, changes in aerosol17 and changes in early summer soil moisture18 and related feedbacks were considered for explaining (part of) the trends. For extreme heat, the increase in the frequency and persistence of split midlatitude jet states over the last 40 years, possibly associated with the reported weakening of the mean summer zonal circulation19, can explain about a third of the amplified trend in heatwave intensity3. Changes in atmospheric circulations around Europe that favor heat were also emphasized20,21, in particular a positive trend in a dipole structure with a low pressure over the Eastern Atlantic22,23 and a high pressure over the Mediterranean extended towards central Europe24. Yet, no increasing trend was found in blocking over Scandinavia that has led to the 2018 heat wave6,25. Moreover, reported changes in Rossby waves are not robust and are sensitive to their exact definition26. In addition, variability of summer temperatures has been shown to be large in Central Europe27. Thus, while several studies have hinted at a potential role of dynamical changes in amplifying European heatwaves, a systematic analysis is lacking, including also how models simulate these changes.

Results
Role of dynamical changes in the temperature trends
We used a method based on circulation analogues to assess the role of dynamical changes in the TXx and TXm trends (see the methods section for a full description). Regional atmospheric circulation patterns are characterized by their 500 hPa stream function over the domain shown in Fig. 1a (black box). We identify circulation analogues for a given day by searching for other summer dates (JJA months) with similar anomaly structures, measured by the spatial anomaly correlation coefficient (ACC). A set of dates with circulation analogues allows us to calculate statistics conditionally to a given circulation28–31, or to assess the role of dynamical changes in circulation-conditioned variables32,33.

In order to estimate the contribution of dynamical changes to TXx and TXm trends (called hereafter the “dynamical TXx and TXm trends”), we replace each daily temperature field by the temperature field from a different day that had the best analogue circulation. In the absence of long-term trends in circulation, this is equivalent to shuffling the temperature time series while keeping the dynamics, thereby creating a trend-free “analogue temperature time series”. In the presence of long-term circulation trends, the trend in the analogue temperature time series comes from the changes in circulations (e.g. an increase in circulations favorable to heat, or vice versa). Replacement by analogues should in principle remove thermodynamical effects from global warming. As global warming is not homogeneous across the time period, and to ensure analogue regional temperatures represent a given global warming level, we further apply a correction by scaling all analogue temperatures to a reference year for global warming (2022) (see Methods). We verified that results were similar in both cases (with and without scaling).

The dynamical TXx trend (Fig. 1b) is generally positive over Western Europe and reaches about 1.5°C/GWD in several areas. The dynamical TXm trend is found to exceed 1°C/GWD over Southwestern Europe (Fig. 1d). Over Western Europe, the average TXm and TXx dynamical trends are respectively 0.74°C/GWD [0.26–1.21°C/GWD] and 0.79°C/GWD [0.24–1.35°C/GWD]. For E-OBS the dynamical trends are 0.78°C/GWD [0.27–1.29°C/GWD] and 0.86°C/GWD [0.29–1.43°C/GWD] for TXm and TXx respectively.

We verify these findings on the dynamical contributions to extreme temperatures trends with a second method, called “dynamical adjustment”34: The method uses a spatial circulation field (here: z500 for consistency with previous studies) as a proxy in order to estimate the contribution of circulation to temperature variability. Here, we use ridge regression, a linear regression technique that regularizes the coefficients of the high-dimensional circulation predictors35, and we subsequently evaluate the dynamical contribution of z500 to the Western Europe TXx trends and averaged results over Western Europe (see method details in the Methods section). Results are consistent with the analogue approach (Supplementary Fig. 2), although with a slightly weaker dynamical TXx trend of 0.56°C/GWD.

To test the sensitivity of our results to the analogue domain, we performed sensitivity experiments by extending and reducing the domain by 10° longitude and 5° latitude (leaving about 2/3 or more of the domain common with the reference one). The dynamical trend is significant and within 0.5°C/GWD and 0.9°C/GWD, except when reducing the domain towards the North-Eastern part (20W-20E; 35N-60N), (dynamical tendency reduced to 0.38°C/GWD) a probable consequence of the key role of the upstream part of the pattern.

Further, we investigate the specific streamfunction patterns associated with summer maximum extreme temperatures over central France [1.5E; 46.5N]—i.e., a region where the TXx dynamical trend is large (see Fig. 1). We select the referenced date (29/06/2019) for which the streamfunction pattern (Fig. 2a) has a maximal average ACC (0.59) with other streamfunction patterns occurring each year when maximal temperature (TXx) is reached at this grid point, so it is most representative of those “TXx days”. We find that about 15% of the summer days in total have an ACC larger than 0.5 with the 29/06/2019 pattern, and that 53 out of 72 other TXx patterns also correlate by more than 0.5. For the sake of simplification, we will refer this class of patterns as the “Southerly Flow” patterns (SF), since almost all of the patterns bear a positive west-east streamfunction gradient (eg. 99% of patterns when considering the gradient between 15°W and 5°E at 50°N), inducing southerly flows over the Western margin of Europe. This pattern also includes a strong anticyclonic component over Central Europe, which induces increased radiation and potential land-atmosphere feedbacks if persistent. As another example, the outstanding temperatures in London on 19/07/2022 were also obtained with a similar circulation pattern (ACC=0.81 with 29/06/2019). To assess sensitivity to the reference pattern we also repeat all calculations with the 10 most representative TXx patterns (Supplementary Fig. 3) in the above sense. In these other cases, the frequency of associated correlated flows is within the 10–20% range.

To check how the SF days contribute to the dynamical trend, we recalculated the dynamical trend excluding the SF days: we removed SF days from the time series, calculated the analogue temperatures of remaining days, the resulting yearly TXx, and recalculated the dynamical trend. We also did the opposite operation by keeping only SF days in the time series. On average over Western Europe (Fig. 2b), the dynamical TXx trend without SF patterns becomes insignificant over Western Europe (0.08°C/GWD on average over Western Europe), while the SF-only TXx dynamical trend is both high and statistically significant (1.3°C/GWD). Similar results are found when using a different reference date among the 10 most representative patterns.

Dynamical TXx trends over Western Europe can therefore largely be explained by changes in the characteristics of SF patterns. First, their frequency has increased by 43% [10%; 76%] per GWD (52% with time between 1950 and 2022) (see Supplementary Table 1). Second, the number of “events” (one event is defined here as a set of consecutive days) per year and their mean persistence have increased (see Supplementary Fig. 4). The persistence of SF patterns has increased by about 24% along the period [-1%, +50%] as a function of GWD. Such changes all give more chance, within a season, to reach the high end of the conditional temperature distribution. Other characteristics may also have changed (eg. amplitude) but were not investigated here.

Significant frequency increases are also found for at least the 10 most representative patterns of Supplementary Fig. 3, with rates in the range of 35% to 55%.

Note that SF is not the only flow pattern changing, and not all patterns associated with TXx days have an increasing frequency or persistence. For instance, the 23/07/2021 pattern, corresponding with summer TXx in central France for 2021, shows no particular evolution (Supplementary Fig. 4). Our results are also consistent with the increase in occurrence and persistence of the specific class of double jet circulation systems explaining a large fraction of European heat extremes3, and about half (i.e., much more than the mean probability, 15%) of double-jet days are found within the SF days.

Simulated temperature trends and their dynamical contributions
The representation of summer TXx and TXm trends has also been analyzed for a large number of CMIP6 model simulations (273 simulations in total for 36 models) (see Methods section for data processing). Over Western Europe, almost all CMIP6 simulations fail to simulate the observed strong TXx trends, as seen in Fig. 3a, plotting the percentage of simulations with larger trends than observed, for each grid point. These differences are less pronounced for TXm (Fig. 3b) but the number of runs reaching the ERA5 trend remains small here too (10-20% in large parts of South-Western Europe). There are also other land areas outside Western Europe where the CMIP6 simulations are mostly above the observed warming TXx trend (i.e. Sahara, Northern Scandinavia, Southern Balkans). This suggests that there is no general underestimation of extreme heat trends over all regions (or land regions). However, understanding these regional discrepancies across the globe is beyond the scope of this article.

When averaging TXx trends over the Western Europe region above defined, only 4 of the 273 individual runs analyzed (members of 3 models out of 36, ACCESS-ESM1, NorESM2-LM and KIOST-ESM) have a larger trend than the observations. The strong TXx trends observed correspond to the ~98-99th percentile of the overall CMIP6 distribution and could, from a statistical standpoint, be interpreted as consistent with Western Europe witnessing a very unlikely phase of low-frequency internal variability. However, in the five large model ensembles that were at our disposal (eg. ACCESS-ESM1-5, CanESM5, IPSL-CM6-LR, MIROC6, MPI-ESM1-2-LR), only ACCESS-ESM1-5 has a few members for which TXx warms as rapidly as observed (Fig. 3c), but this ensemble strongly overestimates the TXm trend (Fig. 3d). Hence, this ensemble does not correctly estimate the daily maximum temperature distribution as observed in ERA5.

Our results are qualitatively robust to the way trends are calculated. We estimated trends relative to time instead of GWD, and to each model initial-condition ensemble mean GWD instead of individual member GWD. In the first (resp. second) case, 9 (resp. 5) simulations (from 4 different models) slightly exceed the ERA5 TXx trend. Trends relative to time allowed in particular two members of CanESM5 to reach observations thanks to the strong global warming (about 1.7°C since 1950), while the regional response to global warming (the regional trend as a function of GWD is about twice weaker than in ERA5.

We also implemented a multiple testing procedure, the False Discovery Rate36–38, to test the significance of the result in Western Europe. Under the hypothesis that “models are indistinguishable from reality”, the rank of the observed TXx and TXm trends in the distribution of members is uniform and there can be regions over which the observation falls outside the model range only by chance. Supplementary Fig. 5 shows that even taking into account the multiple nature of the test, Western Europe is among the regions where the mismatch between observed and simulated TXx trends is significant at the 95% confidence level in the sense of the FDR procedure, while no significant mismatch is found in this region for TXm trends.

Climate simulations do not capture the dynamical changes underlying these temperature extreme changes. We applied the analogue analysis to all available realizations for each model for which 500 hPa wind fields were available (170 simulations in total). This set was found to be rather representative of the overall simulation distributions, albeit with more weight on faster-warming simulations (see Fig. 3a, b histograms) regarding TXx trends. None of their dynamical TXx trends reach the amplitude of the observed one over Western Europe (Fig. 4a). This shows that there is less than 1% chance that the observed trend estimate is drawn from the same population as simulation estimates, accounting for all uncertainties. Remarkably, all members of the three available large ensembles (ACCESS-ESM1-5 [40 members], IPSL-CM6A-LR [31 members] and MPI-ESM1-LR [30 members]) exhibit values lower than observed, despite a few members exceeding the overall TXx trend. Also, on average over Western Europe, for TXm, a handful of models do have dynamical trends comparable to or larger than observations, but all others exhibit lower trends (Supplementary Fig. 6).

We also calculated the thermodynamical trend obtained as a residual by subtracting the dynamical trend from the total trend and reported the result in Fig. 4b. This shows that climate models exhibit thermodynamical contributions that are broadly consistent with ERA5, but there is a tendency for an underestimation of TXx thermodynamical trends, and a general agreement for TXm trends (see Supplementary Fig. 6). This analysis clearly shows that dynamical changes are largely responsible for the mismatch between modeled and observed TXx trends.

All 170 climate simulations realistically simulate the climatological mean frequency of the SF patterns (range from 12.5% to 18%). However, the rapid observed increase in frequency of this flow field (+43%/GWD [10–76%]) is only roughly captured by one among the 170 simulations (NorESM2-LM, and weaker in the others (Supplementary Table 1).

Methods
Calculation of dynamical contributions to mean and extreme summer temperature trends
The method used to estimate dynamical contribution to the change in one variable follows the conceptual framework developed in Vautard et al. (2016), with a different implementation here. It is based on the estimation of the change in the variable solely due to the changes in regional upper-air circulations. For instance, even without extra heating from radiative and diabatic processes, an increase in the frequency of southerly flows in Western Europe would induce a mean regional warming. An increase in anticyclonic conditions would similarly lead to increased radiation and thus temperature. This can also lead to a cooling if increasingly frequent circulations are linked to cooler temperatures (eg. in Northerly winds). To estimate this dynamical effect of changing circulations on temperatures, we need to carefully remove any thermodynamical effect of climate change.

We assume that daily temperature T (which can be mean, minimum or maximum daily temperature, and in the current article will be maximum temperature) has a distribution at a given location or grid point which depends on the atmospheric circulation and on other processes, including global warming. We then assume a decomposition into:

T = hT|X i GWD + T0 (1)

where X is the 500 hPa stream function anomaly, characterizing the atmospheric circulation (simultaneous to the temperature), GWD stands for the global warming degree, <T|X>GWD is the average daily maximum temperature conditioned to the circulation, assumed to be dependent on GWD, and T’ is a fluctuation. This circulation-conditioned temperature includes not only advection effects (i.e. from cooler/warmer regions), but also all processes linked to the circulation (subsidence in anticyclone, increased radiation, surface-atmosphere feedbacks, …), so the overall dynamical trend includes all underlying processes tied to the dynamical conditions. In order to remove thermodynamical effects due to climate change, we scale all temperatures to a reference warming level. For this, we assume that the circulation-conditioned mean temperature depends linearly on the global warming level, so the decomposition can be written:

T = hT|X i ref + b(X) · (GWD - GWD ref) + T0 (2)

where ref refers to a reference global warming level, taken here as that of 2022, so all changes are expressed relatively to 2022. The coefficient b(X) represents the mean warming rate conditioned to the circulation X, which includes thermodynamical effects of the climate change response—it is therefore assumed that the amount of warming depends on the circulation type. Assuming one can calculate b(X) and GWD, all daily temperatures are then scaled to the reference level with the following thermodynamical correction:

T s = T - b(X) · (GWD - GWD ref) (3)

The dynamical contribution to any temperature trend constructed from daily temperatures (e.g. here TXm, TXx) can then be calculated from the Ts time series, because changes with GWD are only through the changes in the frequency of occurrences of X for given GWDs. Trends should also not depend on the particular time Ts values are drawn as long as they occur simultaneously to a stream function anomaly which is similar to that encountered in the same sequence order as that of the series. Hence to increase statistical robustness and remove any residual link to the specific order of temperatures, we replace Ts temperatures by those occurring in circulations X along the time series. This has the advantage of “randomizing” the timing of analogues and providing multiple realizations to calculate dynamical trends. A new temperature analogue series is created by replacing each daily with that of the best circulation analogue, then another new series is made with the second best analogue, etc… (see below for practical analogue calculation). From each of these analogue time series, TXm and TXx are recalculated for each year, then averaged across analogues, and a regression with GWD is calculated at each grid point, together with its confidence interval, (plus or minus twice the standard error of the regression coefficient). To keep analogue quality high, we limit the number of time series to 3. To calculate time series of averages over Western Europe land, we apply the 0.5°x0.5° land mask of E-OBS and average over the grid points included in [-5W - 15E; 45N - 55N].

Estimation of yearly GWD
In practice, GWD is calculated as a moving centered 5-year average of the global temperature with available data, for reanalyses and models, accounting for series ends in ERA5 (i.e. for 1950, taking into account an average only over 1950 to 1952, and for 2022 an average over 2020 and 2021). The 2022 value is then subtracted to all values, so GWD is 0 in 2022, and generally negative before.

Selection of circulation analogues
In practice, circulations are characterized by the 500 hPa stream function over the [-30 + 20°E; 30 60°N] domain. Analogs of a given circulation are characterized by anomaly correlation coefficient (ACC) between stream function fields. For each summer day, we collect the best analogues (highest ACCs), and impose that they remain spaced by 6 days or more within a season, and self-analogues are not considered. This is done by successively testing fields in descending order of the ACC, and skipping days not respecting the separation with previously selected fields.

Calculation of the circulation-conditioned thermodynamical trend b(X)
To calculate b(X), we also use analogue circulations, in a different way than above: For each summer day d of the 1950-2022 period, we estimate b(X(d)) using a regression of each raw temperature T(d) (before thermodynamical correction) associated with a large set of best analogue circulations of X(d) found between 1950 and 2022 with the GWD values of their respective year. We use the best 1% summer analogues (67 days) with the same spacing of at least 6 days. 99% of the worst of these 67 analogues across all summer days have ACC > 0.5, 65% have ACC > 0.7. Imposing a quality criterion on analogues such as ACC > 0.7 or more would leave days with an insufficient number of analogues for regression.

Dynamical adjustment
Dynamical adjustment is used as a second, alternative technique to estimate the influence of circulation-induced temperature trends. This method relies on the idea that temperature variability can be decomposed into a component that is driven by circulation-induced variability, and a residual, thermodynamical component. The “thermodynamical” component is expect to contain a forced signal as well as any other unexplained variability or feedbacks43. Most applications of this technique characterize circulation-induced temperature variability using a proxy variable such as geopotential height34,35,44,45. Dynamical adjustment techniques typically rely on linear methods such as variants of linear regression or circulation analogue techniques.

Here, we use the spatial pattern of z500 in a relatively large circulation domain over Europe and the North Atlantic (-30 to 20°E, 30 to 60°N, similar to Fig. 1), following the method outlined in 46. However, we introduce some modifications and additional details. We use a regularized regression technique, called “ridge regression”, which is well-suited to deal with the large number of circulation predictor grid cells and a relatively short observed record. For TXx, we train our ridge regression model on the 15 warmest days in each summer during 1950-2021 at each grid cell in the ERA5 reanalysis, resulting in a total of 1080 observations (72 summers and 15 days per summer). Since the z500 field contains information about the lower troposphere, and is affected by temperature change via thermal expansion, we detrend the spatial z500 field by subtracting the global average z500 at each timestep and over each grid cell in the circulation domain. Hence, the analysis is based only on relative changes within the z500 field. To obtain regional estimates of the circulation-induced component of TXx, we performed an area-weighted average across the grid cells within the study domain.

Data availability
All analyzes have been conducted using 3 main datasets. The ERA5 reanalysis and the E-OBS datasets (processed from the https://climate.copernicus.eu) has been downloaded, and are available from the Climate Explorer https://climexp.knmi.nl. CMIP6 model simulations are available from the IPSL ESGF node https://esgf-node.ipsl.upmc.fr/.

Code availability
Codes used in this article develop classical statistical algorithms, and are available upon request. Application codes are provided in the archive: https://zenodo.org/record/8310140.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "This study examines trends in extreme heat over Western Europe in models and observations. The\nauthors show that a substantial fraction of the observed trends in extreme heat are linked with\nchanges in atmospheric circulation. The observed trends are on the extreme end of the model\nensemble spread, but those few simulations that do capture the magnitude of the observed trends,\ndo so for the wrong reasons. No simulation shows an increase in circulation-linked heat extremes\nas high is as seen in the models. The authors conclude that models are either missing a circulation\nresponse to external forcing or the models systematically underestimate internal variability –\neither would have important societal consequences.\nOverall, I think that this is well done study on an important topic. The results are novel and, for\nthe most part, convincing. This work will help motivate further research into the causes of the\nmodel-observations discrepancy. The paper is, in my opinion, potentially suitable for publication in\nthe Nature Communications after the following comments addressed:\nThere are two issues that could point to authors potentially underestimating the chance that the\nmodel-observation difference is from internal variability that is accurately captured by the models.\n-First, the authors are examining very small region that was not chosen a priori. It was chosen\nbecause the trends are very extreme. Even if models were perfect, if you looked at every small\nregion there is a good chance that somewhere the trends will be very extreme due to internal\nvariability, and potentially outside model spread. While the observed trends do look very extreme\nrelative to the models at this location, the chances of trends this extreme occurring somewhere by\nchance are substantially higher. This is essentially the ‘multiple hypothesis testing problem’, where\nthe null hypothesis is that the models and observations agree.\n-The second is issue is that the methods the authors use will likely artificially reduce the ensemble\nspread from internal in some models. As the authors are likely aware, a subset of the CMIP6\nmodels have very high climate sensitivity and show too much global warming over the historical\nperiod relative to observations. By calculating the trends with respect to global mean surface\ntemperature, the authors are essentially accounting for these differences in the ensemble means.\nHowever, the spread is artificially altered if the global mean trends are different to that of the\nobservations. It is ratio of the global mean trends to the spread from internal variability that is\nmaintained when calculating trends this way. If two models have the same spread due to internal\nvariability in the raw trends with respect to time, but one model has twice as much global\nwarming, when calculating the trends with respect global mean temperatures, this resulting model\nspread will be reduced by a factor of 2 relative to the other model. This is particularly an issue for\nthe CanESM5 model which warms substantially more than observations and is one of the large\nensembles used to conclude that the observed trends are outside the model spread. Not\nsurprisingly, the ensemble spread in this model also appears to be quite small (Fig 3). Of course,\nmodels that warm less than observations will show too much spread, but I suspect that these are\nless common.\nI don’t think either of these are sufficient to completely invalidate the conclusions, but I think the\nauthors should either account for these issues and/or make it clear that these issues exist.\nOther comments/questions:\nFigure 1: It would help if the different panels had the same colour scales.\nL138-140: Do you have any insights into why the dynamical adjustment methods result in smaller\ntrends? Although the methods are different, they are in principle trying to accomplish the same\nthing.\nL216: It looks like UK-ESM1-0-LL only has 4 members (Fig 3). I wouldn’t consider this a large\nensemble.\nL251-258/Figure 4: I am really surprised that most of the model spread in the trends in extremes\ncomes from the thermodynamic component and that there is very little spread in the dynamical\ncomponent. We always hear about how the thermodynamic responses to climate change are better\nknown and that the dynamics and associated impacts are much more uncertain ( e.g. Shepherd\n2014, 10.1038/ngeo2253). Any comment/speculation on why that is different here?\nL295-298: It should be mentioned that there is evidence that climate models underestimate the\nmultidecadal variability of the atmospheric circulation in the North Atlantic (see e.g. O’Reilly et\nal 2021", "reviewer_doi": "10.1038/s43247-021-00268-7, and references therein)", "reviewer_2": "Review of NCOMMS-23-02718-T: Heat extremes in Western Europe are increasing faster than\nsimulated due to missed atmospheric circulation trends by R. Vautard et al.\nSummary:\nSummer temperatures and the amplitude of heat extremes in Western Europe have increased\nmuch faster than elsewhere in the mid-latitudes during the past decades, leading to frequent\nunprecedented heatwaves and threatening human lives. It is crucial to project how these events\nwill change in the future so that mitigation and adaptation measures can be taken. The climate\nmodel is an essential tool to make such projections and possibly the only tool we can rely on.\nHowever, the models may not sufficiently capture the observed trends, and the decisions based on\nthis biased projection may harm the adaptation and resilience plans. Based on the ERA5 and E-\nOBS datasets and 273 runs from the CMIP6, this study finds that a substantial fraction of the\nobserved trends of the summer temperatures and heat extremes in Western Europe arises from\nthe dynamical component due to changes in the atmospheric circulation, with a southerly flow\npattern. It also suggests that the paces of the increasing magnitude of the summer temperatures\nand heat extremes in Western Europe are severely underestimated by the current state-of-the-art\nclimate models. The reason is that the models underestimated or erroneously represented the\nsoutherly flow-like circulation response to external forcing. Moreover, it suggests that the few\nsimulations reaching the observed warming trends in extreme heat are accompanied by a decrease\nin the occurrence of southerly flows, meaning that the models capture the warming trend for the\nwrong reason. These results are important because they imply a conservative projection and high\nuncertainties for future changes in the summer heat over Europe. They guarantee caution for\nstakeholders and decision-makers. The topic is relevant, the methods are reasonable, and the\nresults are overall supported by the evidence provided. I recommend the authors revise their\nmanuscript to further clarify several places. The following comments are for the authors'\nconsideration.\nMajor comments:\n1. The authors did not describe clearly what benchmark was used to calculate the analogue. It\nseems to be the 500hPa streamfunction on 29 June 2019 (line 147), but it might also be the\n500hPa streamfunction on other dates or its trends (lines 106-108). Without an explicit\ndescription, I cannot judge the reliability of the subsequent analyses and results. Using the\nstreamfunction on 29 June 2019 is OK. Still, the authors need to justify that the circulation in the\n2019 heat extreme is typical enough to represent most of the heat extremes. If it is not typical,\nthe subsequent analysis may be biased. I would also like the authors to discuss the sensitivity of\nthe results to the choice of domain for the 500hPa streamfunction.\n2. Following on from the previous comment, what is the most typical circulation configuration for\nthe occurrence of heat extremes in Western Europe, the southerly flow pattern, or others? The\nanalysis in this study indicates that the southerly flow pattern is important, but it does deny the\nimportance of other patterns. There is a possibility that a second pattern is more dominant than\nthe southerly flow pattern to cause heat extremes and that this second pattern is captured by\nmodels that reproduced the observed TXx trends. I strongly suggest the authors examine this\npossibility because it is directly relevant to the conclusion of this study.\nMinor comments:\n1. The authors interpreted that the modeled TXx trends are likely associated with internal\ninterdecadal climate variability. I believe this is very likely the case. I can easily understand\nthe link between the trend of the TXm and the internal interdecadal climate variability, but I\nhave no idea about the TXx case. Can the authors provide some explanations or interpretations?\n2. Line 58. TXx and TXm-> TXm and TXx.\n3. Line 122. What are the 'both cases'?", "reviewer_3": "This paper explores the mechanisms behind the strong positive trend of heat extremes in Western\nEurope over the past decades. The work separates the trend into dynamical and thermodynamical\ncomponents, using two methods separate methods, showing that a dynamical trend has\ncontributed around 25% of the total trend in annual summertime maximum temperatures (TXx).\nThe research then goes on to explore this trend in the historical simulations of the CMIP6, finding\nthat such a trend for this Western Europe region is produced in very few (3 out of 273 simulations,\ni.e. around 1%) of the CMIP6 simulations. The authors explore the mechanisms behind the trends\nin the CMIP6 simulations using the thermodynamic-dynamic separation, arguing that, even for\nmodels that exhibit a strong trend, the contribution from the dynamical component is too weak.\nThey conclude that the strong trend in extreme heat seen over Western Europe is therefore either\nI. unforced internal variability, which is stronger than seen in the model simulations; or II. a\nforced trend due to anthropogenic warming that is missing in the models.\nThis work is original, provides interesting results, using recently developed methodologies to\nexplore a topic of interest to the many people living in western Europe who have experienced\nthese heat extremes in recent decades, as well as to others who are looking to this region as a\npossible sign of the potential impacts of climate change. The work is certainly of interest within the\nfield of climate change, extreme events, and dynamical contributions; although I don’t see any\nimmediate significance to related fields.\nMy main concern is whether the findings are highly sensitive to the relatively small study region\nchosen. If natural (internal) variability does provide an unusually large contribution to the\nincreasing trend, then, by design, the study has picked the one region that is more likely to be an\noutlier. This region has been studied precisely because there is a strong trend in extremes here,\nand it is important to understand whether this is a forced trend, or internal variability; however, I\nam not yet convinced by the presented results that the finding is as unexpected as the manuscript\nsuggests.\nThe area of the study region (‘Western Europe’) is 20x10 degrees. There are 18x18 = 324 such\nareas on the globe. Given the strong degree of natural variability in the variables of interest,\nparticularly TXx, some regions of the globe are likely to experience particularly anomalous trends,\nand I would expect at least one of these 324 regions is likely to experience a trend that is far into\nthe tail of the distribution from the 273 CMIP6 simulations. Indeed, in Figs 3a,b, there is a lot of\nvariability in this field (% of CMIP6 models with a trend greater than observed) with some regions\nof the Atlantic/northern Africa/Russia with over 90% or even over 95% of simulations exceeding\nthe ERA5 TXm or TXx trend, and other regions with very low values, as mentioned in the text. I\nagree that understanding the global distribution of these values seems outside the scope of this\narticle, but does this large variability not imply that perhaps 273 model simulations are not\nsufficient to sample internal variability for regions this small when looking at a highly variable\nmetric such as TXx? The authors do mention this in line 214: ‘could, from a statistical standpoint,\nbe interpreted as consistent with Western Europe witnessing a very unlikely phase of internal\ninterdecadal variability’; however, my feeling is that this is understated throughout the\nmanuscript. I agree that it is very unlikely for any one specific region to experience this, but it is\nmuch less unlikely that one region somewhere on Earth is experiencing this internal variability,\nand that, by selecting the region based on the trend in TXx, you may have preferentially selected\nthis exact region. I think this needs to be addressed throughout the manuscript. If you could show\nthat no region of Earth saw as high a TXx trend as was seen over Western Europe, this would be a\nmore convincing argument of a discrepancy between observations and the models.\nThat said, it is interesting that the models that do meet or exceed the observed trend do so for\ndifferent reasons, and do not reproduce the observed dynamical trend. However, I again caution\nabout over-interpretation of a rare result in the context of looking only at one small region. Similar\nto my question above: does any model simulate the observed Western Europe dynamical trend for\nany region on Earth? I realise that we don’t expect internal variability to be the same everywhere,\nbut if you could show that the models don’t simulate as strong a dynamical trend anywhere in the\nmid-latitude regions, that would be much more convincing that something unusual, and out of the\nrange of the models, has been happening over Western Europe. I understand the logistical reasons\nfor not running the circulation analogue analysis for all the ensemble members, but when you’ve\npicked one region out of 324 to study, using only 32 models simulations does not feel conclusive\nenough that the observations are so far outside what the models might suggest. Also, why not at\nleast pick the ensemble member closest to the observed TXx trend for the ensemble simulations,\nsince you’re then comparing the dynamical trend to the observed dynamical trend? From Fig 4, I\ndon’t imagine this will make much difference, but it does seem a more reasonable approach.\nI think the wording of the text also understates the possible role of internal variability\n(although it does acknowledge it in places). For example, line 246-247 states: ‘Among those\nsimulations with relatively similar total TXx trends as ERA5, none exhibit a statistically\nsignificant dynamical trend.” What about KIOST-ESM? This model exceeds the ERA5 trend in TXx,\nand the dynamical trend appears to be about half of the observed value – is this not large enough to be\nsignificant? How do you test the significance of the trends in the models and where is this shown?\nOverall, I think this is an interesting study, and (with the exceptions of the points noted above)\nuses valid methodologies. To be published in Nature Communications, I think more concrete\nevidence (as suggested above) that the observations really are so far outside the modeled internal\nvariability plus trend is needed.\nMinor comments\nYou refer frequently to summer mean trends, when I think you mean summer mean daily\nmaximum trends? Also ‘warming trend’ when you mean increase in TXm or TXx. This could be\nclearer – I think most readers will automatically think summer mean trends or ‘warming trend’\nrefers to mean temperature, not mean of maximum temperatures, even if you have defined it\nearly in the paper.\nFor the analogue analysis carried out on the CMIP models, are the patterns you are looking for\nfound in the models, or you find them in the observations.\nFigures 1, 3 etc. You show TXm before TXx, but then discuss TXx first in the text – I suggest\nrewording to make your point more clear, and perhaps showing TXx first if that is your main\npoint? It could also be more clear in the text that ERA5 is shown in Fig 1 and E-Obs is in S.\nFig 1.\nLine 63. Correction needed: Suppl. Fig. 1.\nLine 75 ‘split mid-latitude jet states’\nLine 82 ‘not robust and are sensitive to ….’ (as written it’s ambiguous whether they are, or are no\nsensitive to definition)\nLine 108 ‘allows us’\nLine 114. Could be re-phrased to be clearer, e.g. “we replace each daily temperature field by the\ntemperature field from a different day that had an analogue circulation.”\nLine 119-121. This sounds like you don’t trust that the method does actually remove all\nthermodynamic effects? Or is the lack of difference between the methods proof that it does?\nIf so, could you rephrase as ‘this was tested by re-calculating results with a scaling applied to all\nanalogue temperatures….’.\nLine 140. Is it correct to say that the 0.58C/GWD is equivalent to a linear change of 0.89C over\n1950-2021? If so, this might be a better way of saying it. If not, please clarify what this\nlinear change is, as you didn’t report this for your analogues method.\nFig. 2b. A bit confusing to work out what the red and blue colours are – is it just that blue is\nfor decreasing analogue trends and red is increasing analogue trends? This is a little\nconfusing given that the decreasing one isn’t significant (as inferred from the text). Is blue\ntherefore ‘not significant’ trends. I think a different color scheme could make this plot easier to\ninterpret.\nLine 165. What do you mean by ‘the contribution of SF to yearly TXx becomes almost systematic’?\nLine 193. Please specify what the 5% in parentheses represents? Is this the new (tripled)\noccurrence of double-jet days?\nExtended Fig. 2. ‘Circulation’ is mis-spelled in the caption\nLine 212. Extended data figure 3 doesn’t seem like the correct figure to reference here.\nLine 217, strongly overestimates\nLine 217. The phrase ‘warm as rapidly as observed’ is not clear that you are referring to TXx,\nthis sounds more like you’re referring to mean temperatures. Suggest to rephrase for clarity.\nLine 218 ‘overestimates the mean summer trend’ do you mean average summer temperatures,\nwhich the phrase implies, but which you don’t show, or are you referring to the average daily\nmaximum temperatures as shown in Figure 3c? Please clarify.\nLine 243. What is ‘resp.’ standing for in resp. maximum? Is it necessary?"}, "author_response": "REVIEWER COMMENTS\nReviewer #1 (Remarks to the Author):\nThis study examines trends in extreme heat over Western Europe in models and\nobservations. The authors show that a substantial fraction of the observed trends in extreme\nheat are linked with changes in atmospheric circulation. The observed trends are on the\nextreme end of the model ensemble spread, but those few simulations that do capture the\nmagnitude of the observed trends, do so for the wrong reasons. No simulation shows an\nincrease in circulation-linked heat extremes as high is as seen in the models. The authors\nconclude that models are either missing a circulation response to external forcing or the\nmodels systematically underestimate internal variability – either would have important\nsocietal consequences.\nOverall, I think that this is well done study on an important topic. The results are novel and,\nfor the most part, convincing. This work will help motivate further research into the causes of\nthe model-observations discrepancy. The paper is, in my opinion, potentially suitable for\npublication in the Nature Communications after the following comments addressed:\nWe thank the reviewer for the time spent and careful reading and comments. We hope to\nhave addressed all comments. Please note that a few errors were found in the calculation of\nstreamfunction: (the cdo calculation of vorticity requires a decreasing latitude grid, which was\nerroneously not what we used, and corrupted file problem for ERA5). The results are now\nreplaced and the errors impacted essentially Figures 1 and Figure 4. In Figure 1 the\ndynamical trend is now about 20% smaller. In Figure 4, the mean summer temperature\nresults (TXm) are significantly modified, with higher amplitude in TXm trends for some\nmodels. For TXx the results also show less difference between simulations and\nobservations, but yet none of the analyzed simulation reaches the ERA5 dynamical trend\namplitude. These changes do not alter the overall conclusions, and the abstract is reworded\naccordingly.\nWe then addressed all the questions raised by the reviewer, in particular the two major ones\nbelow.\nThere are two issues that could point to authors potentially underestimating the chance that\nthe model-observation difference is from internal variability that is accurately captured by the\nmodels.\n-First, the authors are examining very small region that was not chosen a priori. It was\nchosen because the trends are very extreme. Even if models were perfect, if you looked at\nevery small region there is a good chance that somewhere the trends will be very extreme\ndue to internal variability, and potentially outside model spread. While the observed trends\ndo look very extreme relative to the models at this location, the chances of trends this\nextreme occurring somewhere by chance are substantially higher. This is essentially the\n‘multiple hypothesis testing problem’,wherethe nullhypothesis is thatthe models and\nobservations agree.\nThe initial motivation to investigate this region was to focus on a region where regional\nattribution studies of extreme temperatures systematically exhibit a strong observed trend\nand a systematic underestimation in model trends. The region is highly populated, and we\nneed to understand these trends for a proper adaptation. It was previously suggested that\nthere are dynamical drivers causing this increase in Tx (e.g. double jet paper by Rousi et al.\n2022), which motivated us to look further into this region as well.\nThat being said, based on comments from 2 of the reviewers on the choice of region and its\nrelevance, we added two analyses.\nThe first one is an investigation of trends in TXx (averaging over land and sea) in each of the\n2088 regions of the world of 20° longitude x 10° latitude between 75°S and 75°N shifted by\nsteps of 5° longitude and 5° latitude (then allowing overlaps). Despite the initial region\nselection was not based on this fact, it turns out that the region selected has the highest\ntrend of all the 2088 regions investigated. This is now mentioned in the discussion of Figure\n1 :“Interestingly,the20°x10°Western Europeregionhas the highest TXx (all year round)\ntrend of all the 2088 20°x10° regions around the globe between 75°S and 75°N shifted by\nsteps of5°,eventhough the region was notadjusted forthis purpose.”\nThe second one is the implementation of a multiple test: the False Discovery Rate\n(Benjamini & Hochberg, 1995; Wilks 2006, 2016), as suggested by the reviewer. The fraction\nof members with a trend higher than that of ERA5 behaves like a traditional p-value, since\nunder H0 \"members are indistinguishable from reality\", its distribution is uniform over [0,1].\nFor TXx trends, with alpha=0.1 (equivalent to a 95% confidence level according to Wilks),\nwe reject H0 on several grid points, including some in our study region (see the Figure below\nwhich has been included as Extended Data Figure 5). This suggests that the model / obs\ndiscrepancy is significant even in the sense of a multiple testing procedure, and gives\nconfidence that the result in Europe is not obtained by chance.\nWilks,D.S.(2006).On“field significance”and thefalsediscoveryrate.Journalofapplied\nmeteorology and climatology, 45(9), 1181-1189.\nWilks,D.(2016).“The stippling shows statisticallysignificantgrid points”:Howresearch\nresults are routinely overstated and overinterpreted, and what to do about it. Bulletin of the\nAmerican Meteorological Society, 97(12), 2263-2273\nBenjamini, Y., & Hochberg, Y. (1995). Controlling the false discovery rate: a practical and\npowerful approach to multiple testing. Journal of the Royal statistical society: series B\n(Methodological), 57(1), 289-300.\n-The second is issue is that the methods the authors use will likely artificially reduce the\nensemble spread from internal in some models. As the authors are likely aware, a subset of\nthe CMIP6 models have very high climate sensitivity and show too much global warming\nover the historical period relative to observations. By calculating the trends with respect to\nglobal mean surface temperature, the authors are essentially accounting for these\ndifferences in the ensemble means. However, the spread is artificially altered if the global\nmean trends are different to that of the observations. It is ratio of the global mean trends to\nthe spread from internal variability that is maintained when calculating trends this way. If two\nmodels have the same spread due to internal variability in the raw trends with respect to\ntime, but one model has twice as much global warming, when calculating the trends with\nrespect global mean temperatures, this resulting model spread will be reduced by a factor of\n2 relative to the other model. This is particularly an issue for the CanESM5 model which\nwarms substantially more than observations and is one of the large ensembles used to\nconclude that the observed trends are outside the model spread. Not surprisingly, the\nensemble spread in this model also appears to be quite small (Fig 3). Of course, models that\nwarm less than observations will show too much spread, but I suspect that these are less\ncommon.\nOn this point, the three ways of calculating trends and the argument in comment\nI don’tthink eitherof thesearesufficientto completelyinvalidate the conclusions,butIthink\nthe authors should either account for these issues and/or make it clear that these issues\nexist.\nWe understand the argument proposed by the reviewer but we think the different ways of\ncalculating trends mean different things. The regional trend per degree of global warming\n(called GWD) is used here to represent the regional response to global warming whatever\nthe global model warming rate is relative to time. As in many studies, and the IPCC report,\nregional consequences of climate change are studied per degree of warming (see eg. all\nfigures of the WGI report). We agree then that, if the model spread does not vary with the\nclimate sensitivity, its spread is reduced.Butthe word“artificially” is debatable.Ifa model\nhas a too strong global warming with time and a correct regional TXx time trend, one may\nalso expect strong regional responses and say that this good performance is artificial\nbecause of some error compensation.\nIn any case we have recalculated TXx trends relative to time and comparisons between\nmodel and observations provide qualitatively similar conclusions. In this case, 9(8)\nsimulations have trends above ERA5(E-OBS) instead of 4(4). Now indeed 2 simulations\nfrom CanESM5 do have trends exceeding observations. However they have global\ntemperatures that have warmed by about 1.7°C since 1950 while ERA5 witnessed a\nwarming of 1°C. We also calculated trends in a third way, relative to the ensemble mean\ntrend for each model, and found again only a small number of simulations 5(5) with trends\nexceeding observation. Our conclusions are therefore not qualitatively modified by the way\ntrends are calculated. We added a paragraph:\n“Ourresults are qualitatively robust to the way trends are calculated. We estimated trends\nrelative to time instead of GWD, and to each model initial-condition ensemble mean GWD\ninstead of individual member GWD. In the first (resp. second) case, 9 (resp. 5) simulations\n(from 4 different models) slightly exceed the ERA5 TXx trend. Trends relative to time allowed\nin particular two members of CanESM5 to reach observations thanks to the strong global\nwarming (about 1.7°C since 1950), probably compensating another bias as the trends per\nGWD, assumed to represent the regional response to global warming are about twice\nweakerthan in ERA5.”\nOther comments/questions:\nFigure 1: It would help if the different panels had the same colour scales.\nWe did not change as the two trends (overall and dynamical) have very different upper\nbounds (+6°C/GWD and +1.5°C/GWD)\nL138-140: Do you have any insights into why the dynamical adjustment methods result in\nsmaller trends? Although the methods are different, they are in principle trying to accomplish\nthe same thing.\nThis was at least partly due to the ERA5 streamfunction dataset, which was partly corrupted\nin the submitted version. Now the dynamical trend is smaller, 0.79°C/GWD, while the\ndynamical adjustment method, which now calculates gridpoint adjustment before averaging\nover Western Europe, instead of the other way around, gives a dynamical trend of 0.56\n0C/GWD. In the revised version, the dynamical trends from the dynamical adjustment and\nanalogue approach are relatively closer than in the previous version. Dynamical adjustment\nshows relatively smaller trends, which may be because this approach is based on the\nprinciples of statistical learning, and trained regression models may not fully capture the\ncirculation-induced variability. This may be caused by several factors, such as a limited\ntraining data or non-linear relationship between the predictors and target variables.\nL216: It looks like UK-ESM1-0-LL onlyhas 4members(Fig 3).I wouldn’t considerthis a\nlarge ensemble.\nThe reviewer is right. The full sentence had a problem and was reformed.\nL251-258/Figure 4: I am really surprised that most of the model spread in the trends in\nextremes comes from the thermodynamic component and that there is very little spread in\nthe dynamical component. We always hear about how the thermodynamic responses to\nclimate change are better known and that the dynamics and associated impacts are much\nmore uncertain ( e.g. Shepherd 2014, 10.1038/ngeo2253). Any comment/speculation on why\nthat is different here?\nThis was due to the error in the streamfunction calculation which is now corrected (see\nabove). Now the results show a larger variability in the dynamical trends, despite\nunderestimated.\nL295-298: It should be mentioned that there is evidence that climate models underestimate\nthe multidecadal variability of the atmospheric circulation in the North Atlantic (see e.g.\nO’Reilly et al 2021"}]

Decision letter