Multi-method attribution analysis of extreme precipitation in Boulder, Colorado

Understanding and attributing the characteristics of extreme events that lead to societal impacts is a key challenge in climate science. Detailed analysis of individual case studies is particularly important in assessing how anthropogenic climate change is changing the likelihood of extreme events and their associated risk at relevant spatial scales. Here, we conduct a comprehensive multi-method attribution analysis of the heavy precipitation that led to widespread flooding in Boulder, Colorado in September 2013. We provide clarification on the source regions of moisture associated with this event in order to highlight the difficulty of separating dynamic and thermodynamic contributions. Using extreme value analysis of, first of all, historical observations, we then assess the influence of anthropogenic climate change on the overall likelihood of one- and five-day precipitation events across the Boulder area. The same analysis is extended to the output of two general circulation model ensembles. By combining the results of different methods we deduce an increase in the likelihood of extreme one-day precipitation but of a smaller magnitude than what would be expected in a warming world according to the Clausius–Clapeyron relation. For five-day extremes, we are unable to detect a change in likelihood. Our results demonstrate the benefits of a multi-method approach to making robust statements about the anthropogenic influence on changes in the overall likelihood of such an event irrespective of its cause. We note that, in this example, drawing conclusions solely on the basis of thermodynamics would have overestimated the increase in risk.


Introduction
Episodes of extreme precipitation have recently received significant attention, both in terms of scientific studies and media exposure. Following such events, many stakeholders are interested in understanding whether and to what extent anthropogenic climate change played a role in the likelihood of such an event occurring. Given their huge societal cost, the perceived increase in large-scale floods over the course of the last few decades is of particular interest. On a global scale, the atmosphere has become warmer and moister over the course of the last century and observations suggest that the extent of this increase in moisture is largely in line with the Clausius-Clapeyron relation of approximately 7% per K [1]. However, the relationship between atmospheric moisture content and heavy precipitation is known to exhibit substantial seasonal and regional dependencies, with atmospheric circulation, vertical stability, and actual moisture availability regularly playing a more important role than the moisture-holding capacity of the atmosphere [2,3]. Addressing the question of attribution must thus be placed in a region-specific context that allows for results to differ from case-to-case [4]. Previous work has demonstrated such dependencies, with the likelihood of heavy precipitation associated with flooding in Thailand [5] and Central Europe [6] shown to exhibit no detectable change, while heavy rain events in Southern France [7] and Northern England in December 2014 [8] have become more likely. By contrast, there are a number of counterexamples: a Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence.
Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. decrease in the likelihood of floods in Northern England in spring [9] and an increase in the likelihood of pressure patterns associated with low rainfall in Southern Australia [10].
The episode of extreme precipitation in and around Boulder, Colorado that led to widespread floods in September 2013 is one such example of a high-impact event that established a new record for flooding damages in the state of Colorado. A daily precipitation total of 230.6 mm was recorded on 12 September 2013 in Boulder, almost double the previous highest annual daily maximum from over a century of observations, while less than 150 km away, the town of Fort Carson set a new state record for daily rainfall with 301.0 mm (http://weather.gov/pub/new 24HourRecordColoradoRainfall). Combined with sustained heavy precipitation over a week-long period and subsequent flooding this resulted in ten fatalities and property damages currently estimated at almost $4 billion (http://denverpost.com/2015/09/12/twoyears-later-2013-colorado-floods-remain-anightmare-for-some/ [3].
In a previous study focusing on a class of events closely resembling that of September 2013, Hoerling et al [3] found an increase in the intensity of high fiveday average precipitable water across the Boulder area to an extent that would be expected from a greater atmospheric moisture capacity in a warming climate. As a result of these thermodynamic consequences, an increase in the intensity of heavy precipitation over the region would also be anticipated. However, Hoerling et al [3] assert that the likelihood of an extreme fiveday precipitation event of this nature has likely decreased as a result of climate change, suggesting a counteracting role of changes in other features of the climate system including atmospheric circulation.
Recently, Trenberth et al [11] proposed an alternative framing of the attribution question that seeks to separate the thermodynamic and dynamic contributions to the likelihood of a given extreme event. The two important differences to previous attribution studies for events similar to the observed case are (1) the focus on only the thermodynamic component and (2) the consideration of the individual event as observed, not a class of events similar to the observed one. To illustrate their point, the authors also discuss the Boulder event, strongly suggesting that the attribution question should consider the anthropogenic contribution to the role of anomalous Pacific sea surface temperature that led to this episode of heavy precipitation but excluding any potential influence on atmospheric circulation. While Hoerling et al [3] identified the importance of large moisture quantities over the region, Trenberth et al [11] suggested that a more comprehensive analysis of the physical mechanisms of moisture transport is required in order to fully understand the role of climate change in the individual event.
In order to provide further clarification on some of the issues raised in previous work, we conduct a comprehensive probabilistic event attribution analysis of the Boulder heavy precipitation event of September 2013. We first of all provide a full diagnosis of the source and transport of moisture associated with this event using back trajectory analysis applied to reanalysis data. We then apply statistical methods using extreme value analysis to both historical observations and the output of two climate model ensembles in order to assess the influence of anthropogenic climate change on the likelihood of one-and five-day precipitation events across the Boulder area. In our conclusions we synthesise the results into attribution statements.

Moisture sources and transport
An extensive report of the heavy rainfall and widespread flooding across the Colorado Front Range during September 2013 is given in Gochis et al [12], to which the reader is directed for a detailed overview of the underpinning meteorological mechanisms. As identified by the authors, while episodes of precipitation leading to flash flooding in this region are relatively common, it was the persistence and spatial extent of heavy precipitation over the course of a week that contributed to the extreme nature of the September 2013 event. Gochis et al [12] also describe the prevailing synoptic conditions, characterized by a slow-moving cyclone over the southwestern United States and a blocking ridge to the north, centred over the Canadian Rockies. This north-south contrast fostered the development of large-scale atmospheric flow and the drawing of moisture northwards toward Colorado. While Hoerling et al [3] acknowledge the abundance of precipitable water across the region, they infer its origin from the southern Great Plains and the Gulf of Mexico via 700 hPa flow anomalies (their figure 5.1(b)). Trenberth et al [11], in contrast, claim that the region off the west coast of Mexico was the most important source of moisture associated with this event based on satellite imagery. Trenberth et al [11] also speculate that anomalous SST in this region and subsequent record amounts of water vapour would most likely not have occurred without climate change. This assertion is key to the argument presented by Trenberth et al [11] that suggests a separation of dynamic and thermodynamic contributions to the impact of a particular event. Gochis et al [12] also highlighted the importance of large-scale tropical SST in generating anomalous moisture. However, in this case, the authors state that such source regions include not only the eastern tropical Pacific but also the Gulf of Mexico.
Here, we seek to provide evidence of the source and transport of moisture specifically related to this heavy precipitation episode using back trajectory analysis following a two-dimensional Lagrangian method [13,14]. For each 6 hourly interval for the period 10-16 September 2013 and each grid cell with precipitation within the target region (107°W-103°W; 42°N-38°N), the product of the horizontal wind and atmospheric moisture content was vertically integrated to construct a series of upstream atmospheric columns up to a period of ten days. All meteorological fields were taken from the ERA-interim reanalysis [15]. Figure 1(a) shows the cumulative evaporation at each grid point that is passed by a trajectory as a contribution to the total diagnosed precipitation. Where a grid point is crossed by two or more trajectories, evaporation is summed to give a total. For comparison, the source regions for all September precipitation between 1979 and 2012 are shown in figure 1(b). There is strong consistency in the spatial spread of the 2013 trajectories and the most prominent evaporation sources appear within the Gulf of Mexico. Net moisture gain (evaporation minus precipitation) is greatest further east across the Straits of Florida and parts of the western North Atlantic (not shown). These findings contradict the suggestion of Trenberth et al [11] based on satellite imagery that anthropogenically-driven SST anomalies in the eastern tropical Pacific were the major source of moisture, a conclusion that would also have been expected on the basis of long term diagnosed source regions for September precipitation ( figure 1(b)).
The results shown here, in addition to the findings of previous work, demonstrate one aspect of the difficulty in separating the contributions of dynamic and thermodynamic processes to the likelihood of an event of this nature. Indeed, limiting the analysis to anticipated thermodynamic consequences might present a misleading impression of the role of anthropogenic climate change. A robust attribution statement requires clear communication of the question that is answered and an experimental set up that does not prescribe the answer. While conditioning on aspects of the circulation will improve the understanding of the event itself, a thorough assessment of a change in risk due to external drivers requires more than one source of evidence drawn from both observations and model simulations [16]. Such an approach is undertaken with regard to the Boulder September 2013 event in the following sections.

Trends in return times of extremes
Extreme precipitation in the Colorado area is known to exhibit large seasonal and regional variability, which was reported on by Mahoney et al [17] in a recent paper. The authors showed that the largest daily precipitation episodes east of the Continental Divide involve the unimpeded movement of moisture from the Gulf of Mexico and occur most likely during the summer months. However, the study concluded that extreme precipitation has historically occurred throughout the year in a state with highly variable topography and that risk assessment should consider potential impacts during all seasons. It follows that attribution analysis should also extend to all times of year and cover all areas with extreme precipitation risks.
Trends in precipitation maxima within the observational record and climate model simulations were investigated. In line with the most damaging aspects of the Boulder floods [12] and in order to avoid limiting the analysis to extreme precipitation events over short periods that are most likely caused by convective activity [17], we consider both one-(RX1day) and five-day (RX5day) maxima. The generalized extreme value (GEV) distribution fitted with annual maxima is a common approach to statistically modelling the distribution of precipitation extremes [19]. When it is anticipated that the distribution of extremes exhibits a degree of dependence on an external environmental process, it is necessary for the distribution parameters to reflect such a dependency. Here, the GEV fit is assumed to scale with a measure of anthropogenic climate change, in this case taken as the global mean temperature smoothed with a 4 year running mean to dampen the influence of interannual variability of the climate system, particularly the El Niño Southern Oscillation. The adjusted location μ′ and scale σ′ parameters are determined such that the ratio between the two remains constant, with an exponential dependence on temperature inspired by the Clausius-Clapeyron relation: where μ and σ are the location and scale parameters of the original distribution and α is the linear trend in precipitation maxima as a function of global mean temperature T. The remaining GEV parameter, ξ, was assumed to be stationary. Such an approach has been applied in previous work in the context of event attribution (e.g. [6,8]). For the observational analysis, all stations within the spatial domain 42°-37°N and 105.5°-103°W with at least 30 years of data were chosen from the GHCN-D v2 dataset [18]. The north-south dimension was chosen in order to include notable Front Range heavy precipitation events earlier in the historical record, such as in June 1965 south of Denver. The western border more or less follows the Rocky Mountain foothills as the mountain stations further to the west have different characteristics [17]. We thus do not study extremes in the Rocky Mountain area, as these are too diverse in the observations and depend on orography that is not resolved in the models. The eastern border is not too far away as annual maximum rainfall increases to the east. Our definition allows the maximum precipitation to be considered similarly distributed. A GEV was fitted separately to RX1day and RX5day aggregated over the set of 116 stations, providing 7452 station years for analysis. In each case, we compare the GEV fit when scaled with the global mean temperature of 1915 and 2013. The exceptional events of 2013 are not included in the GEV fit. The uncertainty margins were estimated using bootstrapping with a sample size of 1000. Spatial dependencies between the stations were taken into consideration using a moving block technique in which, for a given sample station, all remaining stations with RX1day and RX5day precipitation exhibiting a correlation greater than 1/e with the station in question were included in the sample alongside the original bootstrap. This gives about 21 degrees of freedom in the 116 stations for RX1day. For RX5day there are approximately 12 degrees of freedom. This does not affect the central value, but increases the uncertainty margins to more realistic values compared to the wrong assumption that all stations are independent.
As a crosscheck, figure 2(a) shows spatial maxima of RX1day and RX5day between 1879 and 2015, highlighting the exceptional magnitude of the September 2013 event with RX1day and RX5day totals of 266.7 mm and 391.0 mm respectively in our set of stations. The number of stations in this set is roughly constant from 1945 to now, with a slight decline over the last twenty years. The changes in number of stations imply that we cannot use a fit to establish how rare the event was and whether the probability has changed. Figures 2(b), (c) use all data to show differences in the return levels associated with 1915 and 2013 fits indicated by the solid blue and red lines respectively. The upper and lower lines indicate the range of the 95% confidence interval associated with a given return period. The likelihood of a one-day precipitation annual maximum of the magnitude of the 2013 event, indicated by the horizontal line in figure 2(b), has increased by a factor of 1.4 between 1915 and 2013, within a 95% confidence interval range of 1.0-2.0 (p<0.05). For the 2013 five-day precipitation maximum, a larger increase in likelihood is found: 1.9 within a range of 1.1-4.1 (figure 2(c)) (p<0.01). These results correspond to an estimated increase in extreme one-and five-day precipitation of 3.8% (within a confidence interval range of −0.6% to 10.5%) and 4.8% (within a range of 1.0%-13.7%) respectively. This suggests only a small increase in the likelihood of extreme precipitation events in this region. However, the uncertainty bands make the findings compatible with both the null result found by Hoerling et al [3] and an increase in accordance with the Clausius-Clapeyron relation as suggested by Trenberth et al [11].
The same analysis was applied to ensemble output from two general circulation models (GCMs): EC-Earth2.3 and HadGEM3-A. As noted by van Oldenborgh et al [8], extending the analysis to model data allows us to partially remove the statistical uncertainty associated with observations at the expense of greater uncertainty associated with systematic model bias. The spatial resolution of each model is too low to sufficiently resolve complex topography and its influence on daily (and sub-daily) precipitation characteristics. However, it may give an indication of the probabilities on the plains, even though the grid-box averages of a climate model have different characteristics than the point observations of the station data. Indeed, the GEV fit parameters show that the coarser-resolution EC-Earth model has a smaller scale parameter σ relative to the location parameter μ compared to the observations, with the higher-resolution HadGEM3-A model in between these (table 1). The shape parameters do not differ significantly. In the following, we assume that trends in grid box averages and point values are the same. These results correspond to a change of between −1.0% and 6.4% (−1.2% to 5.9%) in one-day (five-day) extreme precipitation.
The analysis was repeated for two 15-member ensembles of HadGEM3-A at a N216 (approximately 0.833°×0.555°) resolution for the period 1960-2013 [21]. The first was driven with observed forcings (including anthropogenic forcings) and sea-surface temperatures (hereafter 'historical' or ANT). The second was a 'counterfactual' ensemble driven with preindustrial forcings and sea-surface temperatures ('historicalNat' or NAT), and thus representative of the evolution of a climate system in a world without anthropogenic climate drivers. Again, comparison is made between the GEV fit scaled to the climates of 1915 and 2013. For the historical HadGEM3-A ensemble we find a trend in the likelihood of a 2013type one-day precipitation event (p<0.1), which is estimated to have changed by a factor of between 0.88 and 2.34 ( figure 2(f)). This is equivalent to a change in the magnitude of such events of between −2.6% and 16.2%. The five-day precipitation maximum is found to be between 0.62 and 2.23 times more likely in 2013 than 1915 (figure 2(g), although the runs start only in 1960 we can set the global mean temperature to the value of 1915). This trend is not significantly different from zero, even at p<0.1. This corresponds to a change in the expected magnitude of −6.9% to 12.8%. No trends are found in the historicalNat ensemble either, allowing us to conclude that the natural forcings have not significant changed the likelihood of either one-or five-day precipitation extremes (figures 2(h), (i)).
The differing setups used permit additional analysis in which we calculated inter-scenario differences and deduce any change in probability as a result of anthropogenic climate change. By comparing the 2013 return times in the historical and historicalNat Had-GEM3-A ensembles, we find no change in the likelihood of a one-day precipitation event similar in magnitude to that observed in 2013 due to anthropogenic emissions. For the five-day precipitation event, the risk ratio is 0.72 (within a range of 0.25-1.22). The large uncertainties do not permit us to conclude that this trend is statistically significant. Unlike in the comparison of the GEV fit scaled to different climate in single scenario analysis, the percentage change in precipitation magnitude is not independent of the return time for which it is evaluated. In order to produce uncertainty quantities that are comparable with the previous observation-and model-based analyses, we evaluate the ANT-NAT differences in the magnitude of precipitation events associated with a 100 year return time. For one-day precipitation, the degree of change is almost zero (−0.2%) with a large uncertainty range (−10.7% to 11.1%) but the magnitude of five-day precipitation totals appears to have decreased (−9.0% within a range of −18.8% to 1.0%). The suggestion in this model that the magnitude of five-day precipitation extremes has decreased in response to climate change is broadly consistent with previous work [3].

Summary and conclusion
A multi-method approach in attribution analysis has the potential to increase confidence in results but the interpretation of findings that are either contradictory in nature or associated with large uncertainties is problematic. Here, both observation-and modelbased methods were used to conduct a probabilistic attribution analysis of the extreme precipitation event that led to widespread flooding in the Boulder region of Colorado in September 2013. Trends in one-and Table 1. Summary of the GEV fits for observation-and model-based data: the location μ, the scale σ and the shape ξ parameters of the GEV fit for 2013. Confidence intervals (95%) shown in brackets. A summary of the risk ratios calculated for each method is given in figure 3. Figure 4 summarizes the equivalent change in the magnitude of extremes in percentage terms, allowing a direct comparison with the change that would be expected in a warming world in line with the Clausius-Clapeyron increase of the moisture-holding capacity of the atmosphere (6%-7%/K) for a warming since the pre-industrial era of 1 K, which is roughly what is observed world-wide. Regional Jun-Sep trends differ from 0.8 times the global mean in the source region of this event, the Gulf of Mexico, to 1.5 times the global mean in the coastal Pacific that corresponds to the climatological source and in Colorado itself. In both figures 3 and 4, we include the multi-method average weighted against the degree of uncertainty associated with each method. Note that this only includes the uncertainty due to natural variability and not the model uncertainty, so the true uncertainty range is larger. However, it was not possible to estimate the model uncertainty with the data available. A chi-squared goodness of fit test was used to assess consistency among methods. As the model spread is comparable to the natural variability in RX1day (χ 2 =1.07) and not much larger in RX5day (χ 2 =1.38) the model uncertainty is estimated not to be larger than the uncertainty due to natural variability. This also validates the assumption that the  trends in grid box averages are similar to the trends in station data, at least to the level discernible in the natural variability.
For one-day precipitation, all results are compatible with a modest increase in the likelihood of extremes ( figure 3(a)). The magnitude is less than what we would expect due to thermodynamics alone. On the low side, the 95% uncertainty range among methods includes no change in extremes ( figure 4(a)). For five-day precipitation, while a significant trend is found in the observed record, no trend is apparent in any of the model-based methods. We conclude that we cannot detect a change in the likelihood of five-day precipitation extremes as a result of anthropogenic climate change ( figure 3(b)) and that it is less than expected by the Clausius-Clapeyron relation alone ( figure 4(b)).
Clearly, while there is considerable uncertainty in the results shown, their collation suggests that dynamic effects and/or local forcings may exert a counteracting influence on the changing likelihood of extremes in a warming world beyond simple thermodynamics. It should not be unexpected that anthropogenic climate change has altered relationships between large-scale drivers and local events. A socalled conditional approach to attribution may seek to address this issue specifically and remains an important means by which to understand the mechanisms involved. However, the diagnosis of moisture sources and transport detailed in section 2 illustrates the Boulder event as an example where there are inconsistencies in the conclusions of different approaches seeking to understand the large-and local-scale linkages. In the meantime, most stakeholders are ultimately interested in defining an event in terms of the damage caused and attribution analysis aiming to answer questions of most relevance to stakeholders should instead follow a holistic approach in order to prioritise an understanding of the changes of overall risk rather than the contribution of different factors [16].
While increases in the frequency and magnitude of precipitation extremes as a result of a warming world are anticipated on a global scale, the understanding of region-and case-specific changes remains a key challenge in climate science. The extent to which the evaluation of the risk of such changes should be conditional on different contributory factors has been the topic of recent debate (e.g. [11,16,23]). In the case of the extreme precipitation event in Boulder, it was possible to draw reasonably robust conclusions about the influence of anthropogenic climate change on the overall likelihood of an event of such a magnitude irrespective of its cause, showing that thermodynamics alone would have overestimated the increase in risk. Ultimately, to be of maximum benefit to the stakeholder, attribution studies should clearly state the assumptions made at their outset, particularly if the attribution is conditional upon rather than independent of a prescribed set of meteorological circumstances.