Authors’ reply to the Discussion of ‘Inference for extreme spatial temperature events in a changing climate with application to Ireland’

Firstly, we would like to sincerely thank all the discussants, particularly the proposer and seconder, for their interest in our paper, the variety of questions they have asked, and their suggestions for enhancements to the methods we use. We are especially pleased that the points that the discussants raise are across the entire spectrum of the paper. Below we separate our response into the main themes that were raised, specifically, fundamental modelling of extremes, the use of climate models, the analysis of missing data, marginal modelling choices, dependence modelling choices, accounting for uncertainty, and potential future extensions.

1 Fundamental modelling approaches for extremes

We made the explicit decision for the paper to focus exclusively on the spatial behaviour of temperature data over Ireland and how this changed over time. We achieved this by working ostensibly with parametric models for independent spatial random fields and accounting for the effects of short-term, weather-driven temporal dependence using a novel block bootstrapping scheme. Clearly, a number of potential questions involving dependence within individual extreme events, the effect of clustering of events over time, and any form of accumulation of past behaviour of the process on future events are therefore ruled out in our strategy. We are, therefore, pleased that Opitz’s discussion contribution presents a variety of examples from agriculture, ecology, and epidemiology and that we could add examples of storm sequencing effects on beach management for coastal flood defence or fluvial flooding, for which the temporal behaviour of the process over a much longer time-scale than that of individual extreme events is critical for applications.

Even within the context of independent and identically distributed replicates of a spatial process, Sochaniwsky et al. raise fundamental questions about what an extreme event is. The separation of natural hazards from the subsequent risk to society is key to our definition of an extreme event. Here, the natural hazard is simply the spatial temperature field on a hot day in summer somewhere in Ireland, while the risk is the combination of the natural hazard and the impact caused by exposure to that hazard (e.g. heat stress, loss of crops, etc.). Our focus is on modelling the natural hazard temperature exclusively. So we do not account for how people, agriculture etc. across the world have different capacities to withstand the same temperature. To answer the point of Sochaniwsky et al. more directly, at each site it is natural to describe an extreme as an exceedance of a high quantile (the threshold) by the natural hazard at that location, with the quantile used being selected objectively so that exceedances are well modelled by a generalized Pareto distribution (GPD). We additionally must account for the spatial continuity of the natural hazard. As we do not want to restrict the type of extreme spatial event being modelled, we describe a methodology to statistically model spatial fields which exceed the marginal threshold at least in one site in the region of interest. By design, for any pertinent question about temperature hazards, the probability associated with that event can be derived in a self-consistent way from our single statistical spatial extreme value model. In this way, we avoid the need for different modelling assumptions for any extreme event of interest and need not be concerned with any physical inconsistencies such an approach might produce.

Our focus was to derive a method to produce replicate simulations of potential extreme spatial fields in order to estimate probabilistic properties of any measure of interest in these events, whilst accounting for all marginal and dependence parameters’ uncertainties. As highlighted by Huser, it is possible to use simulations of an extremal process to derive results about hot spots of the process, which are particularly helpful for identifying where the hazard is at its most extreme for non-stationary processes (Hazra & Huser, 2021). Another useful spatial hazard measure which helps quantify how concentrated a spatial extreme event is, at different levels of extremity, is the severity area frequency curve used by Winter et al. (2016). The flexibility in our approach means that all these measures of extreme spatial events can be derived from our framework.

Our view is that one important contribution our paper makes is tackling the issue of how best to present the evidence of, and the uncertainty in, the spatial characteristics of natural hazards to policymakers. Here, the role of the statistician is vital as it requires a critical combination of simplicity of presentation and for the complexity of the various types of extreme events to be captured. Hence, we are particularly pleased that Smith asks about what the paper provides for policymakers.

What is certainly true is that although policymakers have a good handle on extremes of univariate natural hazards, they struggle to account for dependence in their thinking about multivariate or spatial dependence and its implications for risk assessment or even its effect on the frequency of certain types of spatial events. A particular example is illustrated by Figure 10 in Tawn et al. (2018), where the estimate of the probability of a 1 in 100-year event at a site occurring for at least one site (in a river network over England and Wales) in a year is shown to be $0.9$⁠, whereas the policymakers had expected the answer to be much closer to $0.01$—the case of perfect dependence. So to develop a framework which carries over the marginal and dependence estimation uncertainties into such probability estimates and provide a single framework under which probabilities of extreme spatial events can be estimated self-consistently we see it as a major step to help policy makers make reliable and coherent decisions in a critical aspect of risk assessment that has so far been much understudied.

In answer to Smith’s query about the contribution of the paper more generally, we hope it will be seen as a call for a much stronger combination of skill sets from the climate science and statistical communities. Realistically neither community can sufficiently address issues of natural hazard modelling without taking a collaborative approach which bridges the efforts and knowledge of the two communities. Clearly, there is much to be done in this field, with potential lines of evolution of our proposed methods having been already raised in the discussion. We are excited if our study has raised more questions than it has answered and so it is a trigger point for new collaborative research on this important topic of spatial extreme events.

2 Use of climate and numerical weather prediction models

Our paper relies on climate model output to capture the spatial patterns of the temperature field. However, our choice of a single climate model and our use of it as a covariate has yielded some useful suggestions from the discussants. Similarly, the use of Numerical Weather Prediction Models (NWP), which are constantly ‘nudged’ by observations to match the geophysical model state of the weather state, was also raised. Common themes included how our approach might transfer to other regions or to future time periods; the use of (statistical) down-scaling to obtain higher resolution records; and the construction of the covariates themselves.

Porcu et al. point out that many regions in the world do not have an observational record of similar quantity or quality to Ireland, and so our approach may not translate. This is of course correct, but we do point readers to the considerable efforts being made by NOAA NCEI and the Copernicus Climate Change Service to rescue, collate, and improve global coverage. So in future, our approach may be much more widely applicable. Unfortunately, many regions which have few observations also have climate model outputs with the lowest confidence since we have no means to validate them. A useful future project would be to determine the limits of where our observation-driven approach breaks down. We point the discussants to Noone et al. (2021) to find a discussion on the coverage of observational data.

Cuba et al. note that our approach has much in common with down-scaling, though we note that this term is often applied far more generally to include any approach that uses observations and models to provide information on finer-scaled features. However, in contrast to down-scaling, we do not rely on the spatial dependence structure of the climate model data. Porcu et al. similarly suggest that our approach is somewhat critical of climate models, and we should acknowledge other methods related to down-scaling for extremes, such as those from NWP. We note that NWP is fundamentally different from the climate model projections we use since NWP models are constantly adjusted to align with real-world physics as observed. It is this aspect that makes our approach valid since the two sources of information are entirely independent.

Finally, Huser notes the possibility of our approach being used for extrapolation in time, given that we have found considerable biases in the climate model when looking at extreme temperatures. This was not an approach we took in our paper, as we worked specifically within the time span of the observed data. We believe that it would be feasible to use future projections from a climate model, bias correct them using our observational approach, and thus produce future estimates of return levels.

3 Missing data

Many of the discussants asked about our treatment of missing values that occur in our observational data set. Common themes included the potential interpolation of missing values, the typology of the missing data, and the use of alternative inference regimes (which we defer to Section 6).

Koh suggests using the ERA5 re-analysis to help with event identification, i.e. to help alleviate the problems we encountered with missing data from our study of observational data, extremal spatial dependence, and to better identify large-scale extreme temperature events. However, they note that there is again a fundamental mismatch in spatial resolutions between observations and re-analysis data. This is a very useful suggestion, though we would argue that the main ERA5 data lacks the spatial resolution to capture the fine-grained extremes that may be found at a station level. An alternative would be to use the ERA5-land product, which is a re-running of ERA5 for surface parameters using improved topography and at a higher resolution. However, using this product would entail other changes to our modelling approach, since e.g. HadCRUT5 (which we use as a covariate) would now be redundant. The exploration of changing these data products would be an interesting future exercise.

Richards et al. point out that we assume that the data are Missing (Completely) At Random (MCAR) and propose algorithms such as EM as a useful way to reduce confidence interval widths. We argue that our modelling assumptions correspond to the more general situation of Missing At Random (MAR) and that we have seen very little evidence that more complicated structures, such as Missing Not At Random (MNAR) occur in our data. Indeed the cause of the missingness when we look back through the observational record seems to be related to standard issues concerning the movement of stations and servicing. We see it as exceedingly unlikely that the measuring equipment itself causes data to be missing, though we note that this is only because temperature is relatively easy to monitor. Were we to be looking at extremes of other processes, such as rainfall or flooding, such assumptions may be less valid.

4 Marginal modelling choices

Modelling univariate extremes above a high threshold using asymptotically justified parametric models while accounting for covariates is a well-established approach that goes back to a previous RSS discussion paper by Davison and Smith (1990). The novelty in the marginal modelling methodology in our paper involves dealing with spatiotemporal variation in the threshold itself, and the data below the threshold. In multivariate extremes problems (Coles & Tawn, 1994), an extremal process need not be extreme at all sites at once (relative to each site’s respective marginal threshold), so it is important to account for the below-threshold component of the marginal distribution. Unlike in the multivariate case, the spatial context requires the model to apply across all of the sites in a continuous spatial region, and hence, it requires a structure to be placed on the marginal parameters so that they are well-defined at all sites in the region, not just for those with data. Our contribution here was based on the choice of a form of semi-parametric model, built around a parametric quantile regression, using covariates based on the equivalent quantile estimates from the climate model at the associated sites and the annual mean temperature covariate $MI(t)$⁠. The discussants raise several queries about our strategy, concerning the specific features of our modelling choices and parameterization.

Most fundamentally, Opitz questions our use of a threshold above which we assume that the GPD holds exactly and below which we adopt a different model formulation. In particular, Opitz suggests that we could adopt a fully parametric model for the entire distribution, using a class of sub-asymptotic models which have as much flexibility in the tail as the entire GPD family, but smoothly transition away from that family into the body of the distribution, a class of models termed the extended GPD (Naveau et al., 2016; Papastathopoulos & Tawn, 2013). As one of the authors of this discussion paper was involved in the original construction of the extended GPD models, we are of course familiar with this approach. However, our principal concern with using extended GPD models to describe the whole distribution was that any failure of that model to capture the behaviour of the body of the distribution could lead to bias in the upper tail, which was the region of most concern for Papastathopoulos and Tawn (2013). Secondly, although there is an increasingly large number of published parametric extended GPD models, we were not convinced that a single model form would be appropriate over the whole of Ireland, and such an assumption would lead to alternative modelling challenges. Thirdly, for parsimony, such approaches require that the inclusion of covariates in the model must be common across the whole distribution, which may not be appropriate. However, we are intrigued by the possibilities and look forward to seeing evidence of the successful application of the extended GPD models in spatiotemporal contexts. Of course, there are potential issues with our quantile regression approach, where we avoid the well-established crossing problem, by applying the regression at well-separated percentiles, with Castillo-Mateo et al. suggesting an alternative fitting by using joint quantile modelling for the entire distribution. It would be interesting to see if this approach made tangible differences in our inferences.

Optiz raises questions about the possibility of incorporating additional information into the GPD fitting. Opitz particularly focuses on sources of information about the upper endpoint at each site of the spatially varying marginal distribution. Domain experts in different fields often raise similar queries, at the very least due to the bounds on the energy in any environmental system imposing that the GPD has a finite upper endpoint, i.e. implying that the shape parameter $ξ<0$⁠. Of course, this can be imposed in likelihood inference by restricting the parameter space for ξ. Alternatively, when the endpoint is known exactly, as in Dryden and Zempléni (2006), one can fix the shape parameter through a reparametrization of the scale and shape parameters of the GPD. However, it is more suitable to impose such knowledge (when the endpoint is uncertain) through Bayesian inference. In our experience, experts realistically have a clearer prior knowledge about the values of high quantiles of the distribution, and this is easily incorporated, as in Coles and Tawn (1996). We are grateful for the interesting recent references that Opitz provides on endpoint inferences, and how they are likely to change due to climate change. Such collaborations of climate scientists and extreme value statisticians is an approach we would strongly support in furthering the analysis of environmental extremes.

In a similar vein, Smith wonders if we could exploit information from the climate model shape parameter estimates to better model the observation data than our proposed constant over space and time shape parameter estimate, i.e. as we did for the GPD scale parameters. In Section 3.2 of the supplementary material, we found no evidence from the observation data against having a constant shape parameter, but rather than showing that this implies it’s very likely that the shape parameter $ξo(s)=ξo$ we do appreciate that it is more likely to show that there is very little information in the observed data for such a subtle parameter.

We should have investigated the evidence for spatial variation in the shape parameter for the climate model data more thoroughly. To this end, Figure 1 shows the results of some subsequent analyses that we have conducted. Figure 1 (left) shows estimates of the GPD shape parameter $ξ^c(s)$ over Ireland that are obtained when fitting a GPD separately to each site, i.e. with no spatial parametric constraints on either the scale or shape parameters. The point estimates appear to show some spatial structure with an increase in value from west to east, localized differences particularly on the south coast, and with the pattern inland that is not obviously explainable from standard geographical or climatology features. We, therefore, test site-by-site whether there is evidence that $ξ^c(s)$ differs significantly from the constant shape parameter estimate $ξ^c=−0.189$ that we used in the paper. Based on an assumption of the $ξ^c(s)$ being much more variable than $ξ^c$ and the normality of $ξ^c(s)$⁠, we used the p-value $p(s)=2[1−Φ(|ξ^c(s)−ξ^c|/SE(s))]$⁠, where SE$(s)$ is the standard error of the estimate $ξ^c(s)$ and Φ is the standard normal distribution function. In Figure 1 (right), we present p-values in colour for the grid points where $p(s)<0.05$ and give a black cross otherwise. The plot shows that the majority of point estimates do not differ statistically significantly from $ξ^c$⁠, and indicates that the primary difference from a constant shape parameter arises on the southwest and south coasts, with the former the more spatially contiguous evidence and with smaller p-values, i.e. for the counties Clare, Limerick, and Kerry (north to south). Of course, our testing here does not account for multiple testing or associated spatial dependence, but the $p(s)$-values and their spatial pattern indicate that there is probably a need for some additional spatial variation to be accounted for in both the climate and observation data modelling of the shape parameters in these coastal regions. The differences are probably reasonably small and so we do not believe that incorporating information from the climate model to inform $ξo(s)$ would have made any substantial practical difference to our results. However, we see that there is potential for such an approach, as suggested by Smith, to be important when modelling more topographically diverse regions than Ireland.

In univariate extreme value analysis, there are two ways to define and subsequently model extreme values, namely, the block maxima approach, using the generalised extreme value distribution and the threshold exceedance approach, using the GPD (Coles, 2001). It is helpful to see that both these methods can be derived from a common Poisson process result, as shown by Smith (1989), revealing the interlinked nature of the two parameterizations. Smith questions whether we would be better off using the parameterization of the block maxima for the GPD as the parameters are then independent of the threshold choice (provided a high enough threshold is selected for the GPD to be the appropriate model for excesses of the threshold), unlike in our choice where the scale parameter varies with the threshold choice. We were swayed to our formulation based on a combination of the approach of Chavez-Demoulin and Davison (2005) which shows that with our approach the threshold exceedance rate can be independently estimated from the threshold excess values, and the evidence in Sharkey and Tawn (2017), where it is shown that the choice of the block size to specify the parametrization of the block maxima has as a substantial effect on the inference, with the optimal choice being found to correspond to our GPD parametrization.

Porcu et al. raise a query about the interpretation of our finding of more rapid temporal changes in extreme temperatures than mean temperatures, and in particular notes that we do not provide a ‘comparative analysis of variance and distribution shifts’. We are confident in our findings as the fact that extremes change faster than the mean has been well rehearsed in the literature and is amply documented in IPCC (2021, Chapter 11). Furthermore, we are not modelling the data for a trend in the mean only, but explicitly modelling changing features of the body and tails of the distribution as a whole, and hence capture variance as well as mean changes. For the extremes, this was explicitly found to be relevant to the threshold exceedance rate and the GPD scale parameter. The results we report about changes in the 100-year level integrate all these temporal changes (through the value of $MI(t)$⁠) and their associated uncertainties.

5 Dependence modelling choices

Our modelling of the spatial dependence structure of Ireland’s extreme temperature events through r-Pareto processes enables the estimation of probabilities of a range of features of spatial temperature hazards. Though the paper contained a range of exploratory analyses to support the assumptions we made, a number of discussants have presented interesting comments about our approach, which correspond to the appropriate univariate scale on which to assume an r-Pareto process, the suitability of the r-Pareto process, and on aspects of our spatial covariance of the latent Gaussian process of the r-Pareto process.

Koh points out that de Fondeville and Davison (2022) express two variants of the r-Pareto process, the form that we use in Section 4.3 of the paper, and one on the original margins, which allows the scale parameter to vary over space but not the shape parameter. The key difference between the two approaches is the selection of what constitutes the set of identified extreme events. Although the latter approach allows a clear focus on a set of events of interest and maintains a physical interpretation throughout the modelling procedure, our approach is consistent with the vast majority of the literature on modelling dependence, via copulas and in multivariate extremes, of standardizing marginals before investigating dependence structure and it avoids picking a dependence threshold for a single purpose.

Huser presents some very interesting findings, based on our observational and climate model data, raising issues of the suitability of the r-Pareto process. In particular, that class of processes requires the underlying process to exhibit a property known as asymptotic dependence, which in the notation of the paper corresponds to the limit

being positive, i.e. the process can be larger than its $(1−p)$th marginal quantile at every site in the region $S$⁠, given that the process at one site is larger than its $(1−p)$th marginal quantile, as p tends to 0. In contrast, the process is termed to be spatially asymptotically independent when $χ(S)=0$⁠.

In Figure 4 of the paper, we show empirical evidence that a non-zero extremal occurrence probability (estimated pairwise over different inter-site distances) appears to be a reasonable approximation for the data for both observation and climate model data. Using our transformed to Pareto margins data, assessing higher-order properties of these data, Huser’s extended analysis casts some doubt on the evidence for asymptotic dependence. Specifically, his Figure 1 shows that for the observational data in the very southwest and northwest regions of Ireland (his regions 1 and 8, respectively), asymptotic dependence is probably not a justified assumption, given the gradual decay of his χ estimates to 0 as $p→0$⁠. For the other regions in Figure 1, the estimates are not tending to zero smoothly, or, the estimates seem to be highly variable as the limit is approached. Here confidence intervals will be distorted by missing data as well. In these situations, where there is no additional information source, using graphical tools to separate between asymptotic dependence and asymptotic independence is typically unreliable. Measures of asymptotic dependence are found to be more helpful in drawing conclusions in such unclear settings, see Coles et al. (1999) for a discussion of this in the bivariate context and Eastoe and Tawn (2012) for measures of asymptotic independence that are higher order.

Craigmile et al. suggest that it may be possible to extract more key information from the climate model data in our spatial analysis. Throughout the paper, we emphasized the importance of the climate model data to capture core physical properties of the data and any analysis of observation data needs to embody that knowledge. So, we are intrigued to see that Huser’s Figure 2, reporting similar analyses for the climate model data, has already conducted the analysis we would have done to investigate further the evidence from Huser’s Figure 1. In particular, the results indicate that it is only the two regions 1 and 8 that appear to show a discrepancy from asymptotic dependence, as it appears the limit $χ(Si)>0$ for all regions $Si$⁠, for $i=2−7$⁠, where the numbering matches that of Huser and ${S1,…S8}$ form a partition of $S$⁠. Given these findings, there is some evidence that our modelling in the paper could be improved to allow for parts of Ireland to be asymptotically independent, which suggests that models of the form and flexibility of Wadsworth and Tawn (2022) are required. We are not convinced that allowing for some small regions of Ireland to be asymptotically independent, while the majority is asymptotically dependent, will make much difference to our resulting spatial hazard assessments, but we look forward to seeing future research which explores the sensitivity to this aspect.

The above discussion indicated that the extremal spatial dependence structure over Ireland is spatially varying, a feature not captured via our modelling choice of an r-Pareto process latent stationary Gaussian process with Matérn correlation function. Huser and Porcu et al. also raise the issue of global stationarity imposed by our choice of the Matérn correlation function. They suggest different ways forward, either using a deformation of space, using the approach of Richards and Wadsworth (2021) or a form of parametric local non-stationarity respectively. Whilst we agree such approaches are worth undertaking we are not sure they will have a large effect on our results. Although we have reservations about the ability of climate models to fully capture the precise level of the spatial variability of the observational data, we trust them to contain the core physical properties that determine the drivers of spatial extremal dependence. So, we believe it is best to first explore the climate model data for evidence of factors other than distance, i.e. topography or coastal proximity, which may affect these latent correlations, given the richness in space of these data and their absence of missing data or spatial bias in data locations. Only if evidence is present from such an analysis is it worth investigating further in the observational data. In Section 7.2 of the supplementary material, we suggest a possible way forward to model dependence which changes in a more complex way than simply with distance. Specifically, we suggested an alternative distance metric which also incorporates distance from the coast, which appears to be the secondary contributing factor to dependence over Ireland, though the differences found for regions 1 and 8 by Huser indicate other factors need to be accounted for.

6 Uncertainty

Capturing uncertainty in extremes was a theme that appeared in a number of comments. These useful suggestions tended to focus on the role of deterministic models in our data analysis, either concerning the information they contain with regard to the data, or the uncertainty (or lack of it) that they themselves might exhibit. A second theme was the use of alternative inference approaches, particularly that of Bayesian inference.

Porcu et al. focused on the differing resolution between the gridded regional climate model (RCM) data (approx 10 km resolution) and our point-based station observations. This is of course part of the well-known spatial misalignment problem (Alizadeh, 2022; Klaver et al., 2020). We believe we have partially ameliorated this problem by using the gridded data as covariates, taking the closest grid point to each observation station in each case. More generally, we argue that climate model data are not represented even at a single grid box scale, being usually more representative when using several grid boxes together. However, this problem is far more acute in other variables than temperature, with precipitation being the obvious example.

A number of discussants (Smith, Huser, and Craigmile et al.) raised points about the use of ensembles of climate models, many of which are valid and point to interesting directions of future work. We used just one RCM/earth system model (ESM) combination and so have not attempted to quantify what might occur had we used alternative models that might differ on the temperature field over Ireland, as Smith suggests. Finding which (if any) climate models accurately capture extreme climate behaviour would be of great value, not only for our model but much more widely. There are interesting choices to be made in terms of RCM/ESM combinations that may allow for designed experiments based on a range of features, such as climate sensitivity, jet stream location, etc. A possible starting point for such a data set might include the TRANSLATE ensemble (O’Brien & Nolan, 2023).

Huser raised the related point that we use the HadCRUT5 data as fixed covariates, and so are neglecting the known uncertainty in these values. This could potentially be included in our model using an errors-in-variables approach. Since HadCRUT5 comes with a 100-member parametric uncertainty estimate, these could be used to calculate a joint uncertainty range for use as the errors in variables structure. However, it is unclear how much difference it would make to our output given the volume of data we have, the relatively small uncertainty in HadCRUT5, and the general smoothness of the temperature field. We note that many alternatives to HadCRUT5 have become available since the start of our analysis (e.g. the high-resolution Berkeley Earth dataset¹) though many of these also lack clear uncertainty quantifications.

Rizelli and Castillo-Mateo et al. point to alternative methods for capturing the uncertainty beyond our block-based bootstrapping approach. The most obvious is that of Bayesian inference, where a superior model would capture marginals and spatial patterns simultaneously in a joint model. We did start on such an approach early on in our exploration but found the computational cost of fitting such models using e.g. MCMC, far too demanding. It remains an ambition for future work. We found the paper cited by Rizelli particularly interesting as it uses t-distributions for the extremes. We wonder how this would match with our findings that Ireland’s temperature data had tails that were lighter than the exponential. However, their computation is tractable and so points to an interesting furrow of future research.

7 Possible extensions of our modelling framework

Probably the most over-simplifying element of our approach is its restriction to the analysis of independent spatial random fields, thus suppressing the temporal dependence within events and the carryover effects of previous realizations of the spatial process on the risk and impacts and even on the behaviours of the current events. This was obviously something we were explicit about in the formulation of our paper but is clearly an area of much opportunity for future work. Hence we were not surprised to see a range of valuable suggestions, from across the discussants, for extensions of the work to include in this direction. We certainly see our work as a stepping stone towards a fully spatial–temporal analysis of extremes, which is a much-understudied area for research, despite interesting papers (Davis et al., 2013; Simpson & Wadsworth, 2021) and an increasing focus on compound events and event sequencing from the statistical extreme value community. We agree with Opitz’s statement that ‘for assessing climate-change impacts, it appears equally important to improve how temporal dependence of weather extremes is taken into account, especially regarding the persistence of moderate and extreme anomalies, and compounding effects’. Also, we agree with Porcu et al. that, when considering this more general scenario, special care is needed in addressing how spatial–temporal covariance is modelled in a non-stationary context.

de Fondeville and Davison raise some important points about using mixture models to potentially account for large-scale climate features, such as blocking or major/sustained variations in the jet stream or using suitable covariates for these, and similarly, Koh asks about the effects of local land surface feedbacks and regional dynamics related to prior weather conditions. Porcu et al. questions whether regional climate indices may be useful covariates. There are, of course, significant covariates in the physical mechanisms which lead to extreme weather events. The relevant covariates differ depending on the physical context. For Ireland, the primary covariates are the build-up of a soil moisture deficit and the associated rebalancing of heat exchanges and advection of warm dry air from the continent followed by stilling and clear skies. Our view is that for our analysis at least, including such covariates would lead to a substantial additional complication in modelling due to long-term temporal dependence and could be computationally implausible, and critically it risks massive overfitting issues with great complexity as to interpretation of results. Despite the potential for improvements in inference by incorporating such additional physical covariates, our finding from exploratory analysis was that even taking the well-known NAO index as a covariate was not particularly helpful for our analyses for Ireland as it was very weakly related to the extreme events of interest.

The covariate or latent mixture features which de Fondeville and Davison, Koh, and Porcu et al. raise can potentially affect both the marginal distributions (at a site given the covariates or mixture component) and the spatial dependence structure. The recent paper by Richards et al. (2023) illustrates the values of this type of approach when dealing with rainfall fields with different frontal and convective behaviours. The authors show that the use of high-resolution climate models enables a suitable dichotomous covariate to be derived. Extracting such covariate information from observational data at a limited number of sites is not practical. Even if clear mixture structures are found then this presents two further problems. Firstly, to derive overall results there is a need to integrate over the mixture components (Eastoe & Tawn, 2009) which requires the estimation of the mixture probabilities, and how the distributional properties of each mixture term and the mixture probabilities each evolve under climate change. Secondly, and possibly more importantly, is the question of how best to incorporate this mixture structure knowledge into advice about spatial hazards for public policymakers. We do agree that having a coherent modelling strategy for dealing with mixtures in extreme value is an important and missing element from the current extreme value modeller’s toolbox, which we hope will be addressed soon. As in the spirit of the core message of our paper, in the context of climate extremes, we do believe that information from climate models will play a valuable input in this development.

We thank Richards et al. for drawing our attention to the interesting paper of Turkman et al. (2021) which shows encouraging results from the calibration of numerical model output to observational data within a GPD setting. We agree that this paper is tackling a very connected issue to ours. However, we see no direct applicability of that work to our application/setting, as unlike in that case our climate model data has not been run using the same weather events as for the observed data, so other than distributional links between $Xc(t,s)$ and $Xo(t,s)$ the processes are independent of each other. Furthermore, their focus is purely on univariate calibration, with no consideration of the spatial dependence of the processes $Xc(t,s)$ and $Xo(t,s)$⁠. Their use of a spatially dependent latent process to help in the calibration is a concept we believe could improve our approach. Specifically, we have used simply the nearest site on the climate model grid to give $σc(s)$ as a covariate for $σo(s)$⁠, but it may be more sensible to use locally weighted averages of these values. A further difference in our approach to the work in Turkman et al. (2021) is that they have a fully parametric marginal distribution (the extended GPD), which for the reasons we expanded on in Section 4, of this response, in relation to Opitz’ point about sub-asymptotic models, we have reservations about using for our context.

Craigmile et al. ask about whether using an approach in Aberg and Guttorp (2008) to derive the distribution of the maximum in a region, could be helpful for our investigation of the probability of $At,S(T)$⁠. In that paper, they use a variant of Rice’s formula for deriving the distribution of the maximum over a region of a continuous random field. There are two aspects of this approach which we do not feel is best suited to our needs. Firstly, evaluation using Rice’s formula is most tractable for Gaussian processes, as considered in Aberg and Guttorp (2008) and secondly, we have no value in deriving the distribution of the maximum over continuous space, as our climate model is already on a fine enough grid relative to the smooth spatial variation in the temperature process over Ireland. Trying to derive features of the process at a finer scale than this using statistical methods risks introducing noise rather than capturing any additional physical properties available in the data.

We thank Opitz for providing some helpful references as an introduction to the rapidly evolving area of generative adversarial networks and their adaption to extreme events which we hope to be influential to others, who, like us, are new to this intersectional literature. Clearly, it is early days to evidence the value of these approaches relative to established approaches of extreme value theory and climate science encapsulated in climate models. We struggle to see that any entirely data-based generation can provide better inference at levels that are more extreme marginally than any event previously observed nor in new unseen climate regimes than inferences from some combination of extreme value methods in conjunction with syntheses of climate model-produced data. Simply, as the climate models capture physical constraints and can explore future emissions scenarios outside the range of the observed data and extreme value methods have a long history of providing a parsimonious method for tail extrapolation. Time will tell, and we look forward to seeing future developments.

Data availability

Not applicable.

References

Footnotes

Author notes