Anthropogenic aerosols mask increases in US rainfall by greenhouse gases

Domain: earth_science
Published: 2024-02-22
Editorial outcome: major_revision

Abstract
A comprehensive understanding of human-induced changes to rainfall is essential for water resource management and infrastructure design. However, at regional scales, existing detection and attribution studies are rarely able to conclusively identify human influence on precipitation. Here we show that anthropogenic aerosol and greenhouse gas (GHG) emissions are the primary drivers of precipitation change over the United States. GHG emissions increase mean and extreme precipitation from rain gauge measurements across all seasons, while the decadal-scale effect of global aerosol emissions decreases precipitation. Local aerosol emissions further offset GHG increases in the winter and spring but enhance rainfall during the summer and fall. Our results show that the conflicting literature on historical precipitation trends can be explained by offsetting aerosol and greenhouse gas signals. At the scale of the United States, individual climate models reproduce observed changes but cannot confidently determine whether a given anthropogenic agent has increased or decreased rainfall.

Introduction
Daily accumulated precipitation, including precipitation associated with extreme events, is an important part of the global water cycle. Precipitation is particularly important considering decreases in natural water storage, including snowpack, glaciers, and groundwater. As a result, a comprehensive understanding of precipitation change is critical to human systems, including agriculture, water resource management, and infrastructure design. Such knowledge can underpin mitigation policies and adaptation in response to changing risks of natural hazards such as flooding and droughts within a nonstationary global climate.

While anthropogenic influence has been identified for many aspects of the Earth system, robust conclusions regarding the human influence on regional (sub-continental) precipitation remain difficult to obtain. Existing studies primarily address changes at the global scale, zonal land-averages, or continental-scale averages. Such large-scale statements about anthropogenic influence on precipitation are highly useful but do not provide the information needed to understand the nature of local climate change, for example, the magnitude and direction (increasing or decreasing) of the change. Attempts to attribute local-scale precipitation trends have proven to be largely inconclusive even over the continental United States (CONUS) where there are well-documented century-length trends in seasonal mean and extreme precipitation. A recent study identifies a statistically significant human influence on regional precipitation over Europe, but only for mean precipitation and only in winter. One of the primary reasons existing studies struggle to robustly assign human influence is due to their reliance on global climate models, and climate model uncertainty is one of the primary factors that limit confidence in regional attribution for precipitation. For example, anthropogenic aerosols have a significant influence on regional precipitation change over the CONUS; however, the signal-to-noise (SNR) ratio for individual ensemble members of single-forcing anthropogenic aerosol climate model runs ranges from –5 to +3 (indicating that aerosols drive both statistically significant increases and decreases in precipitation; see Supplementary Fig. 1). As such, new methods are needed to attribute human influence on regional precipitation, ideally approaches that reduce direct reliance on global climate models, explicitly model natural drivers of precipitation, disentangle the complex causes of regional precipitation change, and account for anthropogenic aerosols.

The purpose of this study is to implement methods developed for regional detection and attribution (D&A) that provide robust conclusions regarding the human influence on seasonal mean and extreme precipitation over the CONUS. Despite dense measurements of long-term rain gauge records, the Sixth Assessment Report of the Intergovernmental Panel on Climate Change indicates low to medium confidence at best and no agreement at worst in the nature of precipitation change over much of North America [Figure SPM.3]. Here, we explicitly decompose the uncertain combined anthropogenic signal into separate contributions from two of the most important forcing agents. This decomposition allows us to conclusively attribute changes to these forcing agents. Importantly, our approach utilizes climate model simulations indirectly to identify an appropriate formula for modeling a timeseries of precipitation [ref. 29, see Eq. (2)] using Pearl-causal inference. Climate model output is then set aside, and we interrogate in situ records from rain gauge measurements which yields Granger-causal attribution statements. Definitive conclusions regarding the spatial patterns and time-to-emergence of human influence on regional precipitation are made possible by simultaneously accounting for both anthropogenic aerosols (globally and locally) and greenhouse gas emissions [much like the two- and three-way analyses in ref. 27]. Complementary analyses of climate simulations assess the degree to which observed relationships can be reproduced by physical models.

Our approach to regional D&A has both similarities and important differences relative to more traditional D&A methods that rely on optimal fingerprinting. A detailed comparison is provided in the Methods (see “Comparison with optimal fingerprinting methods”), but at a high level the important points are as follows. For example, our implementation can be seen as analogous to a two-way regression analysis with single-forcing greenhouse gas-only and anthropogenic aerosol-only experiments. However, unlike optimal fingerprinting, the regressors or so-called “fingerprints” are not model-simulated quantities but instead fixed forcing time series that are reconstructed from observations. The primary motivation of this work is that model uncertainty is a major barrier to regional D&A for precipitation, and hence our methodological choice to not rely on model-simulated responses or inter-model differences is intentional. Also, unlike optimal fingerprinting, we explicitly model part of the internal variability via climate drivers and estimate the magnitude of the internal variability directly from observations, such that uncertainty in this estimate is propagated through to the attribution conclusions. Finally, unlike traditional optimal fingerprinting, we suppose that attribution statements are localized and vary over space (as opposed to a single conclusion for the entire domain), such that D&A statements with uncertainty quantification can be made for either individual grid boxes or spatially-aggregated grid boxes [e.g., all of CONUS or attribution subregions such as in ref. 36] in a single framework. It is also noteworthy that recent research has identified serious problems with the traditional implementation of optimal fingerprinting, namely that it underestimates uncertainty and yields overconfident attribution statements and furthermore produces biased estimates of the scaling factors.

Over the last century, the principal anthropogenic forcing agents for precipitation over the CONUS are well-mixed greenhouse gases (GHGs) and aerosols. Best estimates of the five GHG concentrations (carbon dioxide, methane, nitrous oxide, and chlorofluorocarbon [CFC] 11 and 12) are available from the boundary condition files used in the sixth phase of the Coupled Model Intercomparison Project (CMIP6) and can be converted to their corresponding radiative forcing on the atmosphere; see Supplementary Fig. 2(a). The reconstructed GHG forcing time series used in our analysis involves a time lag to account for the lagged response of sea surface temperatures (SSTs) to GHGs (see the red line in Supplementary Fig. 2d), which dominates the effect of GHGs on precipitation (see “Analysis of GHCN in situ records” in Methods). Based on ref. 45, we assume the GHG forcing is spatially uniform across the CONUS.

Anthropogenic aerosols are more difficult to account for since their effects on precipitation are multi-faceted and relevant century-length observed quantities are significantly limited. For example, while GHG forcing imposes primarily a lagged effect (the “slow” precipitation response) on the climate system, the effects of anthropogenic aerosols on precipitation have non-negligible lagged components due to cooling of SSTs as well as fast components due to aerosol-cloud interactions. In order to account for the slow precipitation response to anthropogenic aerosols, we utilize an observationally constrained time series of historical aerosol effective radiative forcing, denoted “AER-glob”; see Supplementary Fig. 2(c). This forcing time series describes the effects of all non-local anthropogenic aerosols, including remote aerosol emissions from Asia and Europe. As with GHG radiative forcing, we apply a time lag to account for the SST-mediated response on the climate system (see the blue line in Supplementary Fig. 2d). Even though SO2 is the dominant aerosol species for changes in precipitation in the CONUS [Hypothesis 4a of ref. 29, see also Supplementary Table 1], it is nontrivial to explicitly characterize the fast precipitation response to local sulfates in an observational analysis. Without long-term, spatially-resolved observations of, e.g., atmospheric concentrations of SO2, we must rely on climate models. However, the diversity in chemical and physical parameterizations and in atmospheric dynamical formulations across multi-model ensembles yield vastly different concentration and surface deposition rates. Even CONUS-mean SO2-related quantities across simulations from climate models in the Aerosol Chemistry Model Intercomparison Project [AerChemMIP; 50] differ significantly [see Figure G3 of ref. 29]. Aerosol emissions are a prescribed quantity in historical simulations and are hence consistent across climate models, and ref. 29 show that regionally-averaged time series of SO2 emissions can be used to appropriately quantify the fast precipitation response in each season to anthropogenic aerosols over CONUS. Note that our results are insensitive to the specific method used to derive localized estimates of emissions’ influence on precipitation (see Supplementary Fig. 4). Supplementary Fig. 2(b) shows best estimates of CONUS-wide seasonal emissions trajectories from the last century obtained from refs. 51, 52, denoted “AER-local”. These estimates show that SO2 emissions trend upward over the first two-thirds of the 20th century (much like GHG forcing) but then, following the introduction of clean air regulations in the mid-1960s, decline sharply to their low present-day levels.

Ultimately, we use the sum-total lagged GHG and AER-glob forcing time series (see the black line in Supplementary Fig. 2d) to quantify the slow precipitation response to anthropogenic influence, and we employ regionally-averaged local SO2 emissions (AER-local) to quantify the fast precipitation response to anthropogenic aerosol forcing. Each of these forcing time series are actually proportional to the fast and slow precipitation response. Statistical attribution coefficients are estimated from rain gauge data to translate the forcing time series to the corresponding rainfall response (see “Analysis of GHCN in situ records” in Methods).

Results
We analyze in situ measurements of daily precipitation from rain gauges in the Global Historical Climate Network [GHCN; 53, 54], using a set of approximately 2500 high-quality stations with records dating back to 1900. Equipped with the general D&A formula and its simplification for analyzing seasonal mean and extreme daily precipitation over the CONUS presented by [ref. 29, see Eq. (2) in Methods], we generate best estimates of the fast (AER-local) and slow (combined GHG and AER-glob) statistical attribution coefficients, denoted βFast and βSlow, for 0.25°×0.25° longitude/latitude grid boxes (see “Analysis of GHCN in situ records” in Methods). Uncertainty quantification allows us to attribute the observed changes to seasonal mean and extreme precipitation for both individual grid boxes and aggregated subregions while accounting for multiplicity in testing (see Methods, ibid.). We further use the attribution coefficients (and their uncertainties) to examine time-to-emergence for each anthropogenic forcing agent (see Methods, “Summarizing the GHCN analysis”).

Spatial scales of attribution: fast versus slow response
The detection and attribution of anthropogenic climate change is inherently a signal-to-noise problem, and a common approach for increasing the signal-to-noise ratio (SNR) is spatial aggregation. Given the general challenges associated with attributing changes to precipitation for individual grid boxes and limited success with attribution for sub-continental scales, our first result examines the spatial scales for which we can confidently attribute changes to regional precipitation using our D&A framework. Starting with the entire CONUS, we subsequently divide the CONUS into two, four, 13, and 75 nested subregions [using the attribution regions defined in ref. 36, see Supplementary Fig. 3 and the right-side panel of Fig. 1] while also assessing individual 0.25°×0.25° grid boxes. The attribution regions correspond to spatial scales of ≈8 Mm2 (all of CONUS), ≈4 Mm2 (two subregions), ≈2 Mm2 (four subregions), ≈0.5 Mm2 (13 subregions), and ≈0.1 Mm2 (75 subregions), where 1 Mm2 = 1 million km2; the grid boxes are ≈600 km2. For each set of subregions, we area-average the statistical attribution coefficients, test the two subregion-specific null hypotheses H0,f : βf = 0 for f ∈ {Slow, Fast} (one for each response), and after applying a multiple testing adjustment, identify the subregions for which we can confidently attribute changes in both the slow and fast precipitation response to external forcing.

Figure 1 tallies the fraction of CONUS for which we can attribute changes with both moderate and strong significance as well as the sign of the attributed change (i.e., if the forcing agent or agents drives increases or decreases in precipitation) across the spatial scales considered. When considering the entire CONUS, we can attribute changes to the fast and slow response with confidence across many seasons, although in some cases attribution can only be made for moderate significance. Note that there is a near monotonic decrease in the fraction of CONUS with significant attribution as we move from large to small scales across forcings, seasons, and precipitation type. For individual grid boxes, the slow response is significant for at least some grid boxes for spring, summer, and autumn for both mean and extreme precipitation. The fast response has a highly significant drying effect for CONUS-average mean precipitation in the winter; in summer and autumn, the fast response instead results in (primarily) enhancements to mean and extreme precipitation, as found in prior multi-model experiments concerning the response of precipitation to sulfates (see below for further detail). In JJA and SON, the enhancements in extreme precipitation are significant for ≈0.1 Mm2 (and larger) regions for extreme precipitation and ≈0.5 Mm2 (and larger) regions for mean precipitation. These attribution conclusions, which are based on the relative comparison of the signal (β) and the noise (comprised of uncertainty from precipitation vs. f relationships and from internal variability of the climate system), imply that the SNR for the slow precipitation response (AER-glob plus GHG forcing) remains large for very small spatial scales, while the SNR for the fast response (AER-local) is detectable down to spatial scales of ≈ 0.1-0.5 Mm2. Furthermore, since uncertainty due to internal variability of the climate system is the same for all forcing agents and the magnitude of the signal is comparable for the fast and slow precipitation response (see Fig. 2), this implies that we have higher certainty for quantifying the slow precipitation response relative to the fast precipitation response at the finest spatial scales considered here.

Grid-box attribution for precipitation change
Because both GHG and AER forcing have an attributable human influence on mean and extreme precipitation response for individual 0.25°×0.25° grid boxes in at least one season (Fig. 1), we first explore the statistical attribution coefficients and their significance for GHG, AER-glob, and AER-local forcing at these very small spatial scales. Figure 2 shows the product of the attribution coefficients and the range of each forcing agent such that the plotted units describe the effect of each forcing agent on precipitation (see “Summarizing the GHCN analysis” in Methods); hatching indicates a statistically significant attribution (i.e., that the null hypothesis of no anthropogenic influence is rejected; determined with both moderate and strong significance). The anthropogenic signal is strongest for the GHG influence on 20-year return values of extreme daily precipitation (Fig. 2b), where there is evidence that human-induced GHG forcing causes changes in extremes in all seasons except winter, ranging from 21% of CONUS in summer to 36% of CONUS in spring. GHG forcing primarily causes present-day extreme values to exceed their early-1900 totals by as much as 10 mm day−1, with the largest increases in the central US (in winter, spring, and summer), the northern Great Plains (in spring), the southeast (in autumn), and the northeast (in spring and autumn). GHG-driven changes in mean precipitation (Fig. 2a) have similar patterns to those present in extreme precipitation, as evidenced by the high pattern correlation between the effect on mean and extreme precipitation ranging from 0.65 in winter to almost 0.9 in autumn. These relatively strong correlations are noteworthy since changes in means and changes in extremes are driven by different dynamical and thermodynamic mechanisms, such that there is no guarantee of any correspondence in these spatial patterns. Finally, with at least moderate significance we can also conclude that GHGs cause average daily precipitation totals to increase by as much as 1 mm day−1 in many places (ranging from 6% of CONUS in spring to 18% of CONUS in autumn). We cannot ascribe statistical significance to grid box mean changes in the winter.

The slow AER response has the opposite sign to that of the slow GHG response by construction (see Methods) since the AER-glob and GHG forcings have opposite signs and the magnitudes of both monotonically increase over most or all of the GHCN record. Therefore significance statements for the slow precipitation response to aerosols are identical to the GHG hatching. As expected from atmospheric theory, the slow precipitation response to aerosols is negative almost everywhere, with only limited areas (and generally non-significant) showing increases to precipitation.

While GHG- and AER-glob-driven increases in precipitation are in line with model-based analyses, the fast precipitation response to local aerosols (as quantified by SO2 emissions) is much more nuanced. In the winter and spring, increases in SO2 emissions result in drying for both mean and extreme precipitation over 60–75% of the CONUS (see Table 1). This is consistent with the global-mean drying seen in single-aerosol-forcing historical model runs. On the other hand, increased SO2 emissions enhance both mean and extreme precipitation during the summer and autumn. Multi-model mean estimates of the seasonal fast and slow precipitation response to sulfate aerosols (see “Fast versus slow precipitation response to aerosols” in Methods) derived from experiments in the Precipitation Driver and Response Model Intercomparison Project [PDRMIP; 43, 57] reveal that the best estimates of the fast response to aerosols averaged over the CONUS are increases to precipitation, especially in the summer (see Supplementary Fig. 5a). In summer, the spatial patterns of the fast precipitation response to aerosols estimated from PDRMIP (Supplementary Fig. 5b) show strong correspondence with our GHCN-based estimates in Fig. 2: decreases in the eastern U.S. with large and statistically significant increases in the central and northwest U.S. (Supplementary Fig. 6 shows a side-by-side comparison of the spatial patterns).

It is important to note that while our best estimates sometimes indicate SO2 enhancements to both mean and extreme precipitation, in nearly all cases these local changes are often not statistically significant. This lack of statistical significance is because, as mentioned previously, the fast precipitation response to local aerosols is simply much more uncertain than the slow precipitation response.

Lastly, it is clear from Fig. 2 that the the fast and slow precipitation responses are of comparable magnitude, up to ±1 mm day−1 for mean precipitation and ±10 mm day−1 for 20-year return values. Table 1 summarizes the joint distribution of the statistical attribution coefficients for the fast versus slow precipitation response: in winter and spring, it is most common for the slow response (Slow+) to be offset by the fast response (Fast−; for both mean and extreme precipitation), although for roughly 20–30% of the CONUS increases in the slow response are enhanced by additional local AER-driven increases in the fast response (Fast+). For summer and fall, the dominant category is slow response increases further enhanced by fast response increases (again for both mean and extreme precipitation), although the Slow+/Fast− category maintains 20–26% of the domain. The fact that these two anthropogenic agents have an equal effect on precipitation while both their joint behavior and the fast vs. slow precipitation response to aerosols vary by season reiterates that D&A studies must account for both forcings when attempting to interpret the historical record [as discussed in 28, 30].

GHG signal emergence masked by aerosols
We next show CONUS-average anomalies for each forcing agent relative to the early 20th century climate in Fig. 3 to assess anthropogenic influence as a function of time. While spatial averaging disguises important heterogeneity in the statistical attribution coefficients (shown in Fig. 2), assessing CONUS-wide averages allows us to evaluate the overall trajectories of precipitation change over the last century. Furthermore, while we are best able to identify anthropogenically-induced changes at the ≈ 8 Mm2 scale (i.e., all of CONUS), we also explore smaller spatial scales later in the section. An important conclusion of ref. 29 was that the total anthropogenic response (denoted ANT) can be represented by summing the individual effects of GHG, AER-glob, and AER-local forcing; hence, we can also compare the CONUS-wide trajectory of the combined anthropogenic influence on precipitation. For both GHG and ANT trajectories we can then identify the year in which each signal emerges (see “Summarizing the GHCN analysis” in Methods; “emergence” is shown by dashed vertical lines in Fig. 3). We do not identify an emergence time for AER-glob or AER-local spatial scales smaller than ≈8 Mm2 (all of CONUS)? For comparison, we show corresponding time-to-emergence plots for the ≈4 Mm2 (dividing CONUS into two subregions) and ≈2 Mm2 (dividing CONUS into four subregions) scales in Supplementary Figs. 7 and 8. To summarize these and even smaller spatial scales, where CONUS is divided into 13 subregions, for each season and precipitation type, we identify the emergence time for the isolated GHG signal and sum-total ANT signal as in Fig. 3 for each CONUS subregion and calculate the difference between the two emergence times (ANT minus GHG). These differences are shown in Fig. 4; we also note cases where only one or neither of the GHG/ANT signals emerge. For all spatial scales where either signal emerges, local aerosols mask the GHG signal. This outcome is actually more common for extreme precipitation relative to mean precipitation, indicating that masking from the fast precipitation response to local aerosols is more prominent for extremes. There are of course some cases where only the ANT signal emerges, but in all cases this occurs at the very end of the record (2010 or later).

In summary, our results show that uncertainties related to the emergence of a detectable and attributable human influence on regional precipitation over the CONUS over the historical record [Figure SPM.3; 32] can in most cases be explained by decomposing the combined ANT signal into the fast precipitation response to local aerosol forcing and the slow precipitation response to both GHG and AER-glob forcing. Furthermore, we find that the emergence of the isolated GHG signal is masked by local aerosols for spatial scales as small as ≈0.5 Mm2.

Challenges in model-based attribution
Turning to an assessment of global climate models (GCMs), we apply our general D&A formula to individual ensemble members from the CMIP6-historical experiment and apply a weighting scheme that emphasizes internal consistency of historical results relative to pre-industrial and single-forcing runs (see “Analysis of the CMIP6 historical multimodel ensemble” in Methods). Figure 5 summarizes weighted averages of CONUS-wide estimates (including uncertainties) of fast precipitation response to aerosols and the GHG and AER-glob contributions to the slow precipitation response from 316 individual ensemble members representing 25 distinct GCMs. Figure 5 also shows corresponding CONUS-wide estimates from the GHCN rain gauge analysis. As in the previous section, we again focus on CONUS-wide averages for both to simplify presentation and to focus on the outermost scale at which GCMs should have the maximum skill in attributing human-induced changes to precipitation. A large fraction of the multimodel ensemble (MME) captures the expected GHG-driven increases in mean and extreme precipitation, particularly in winter, spring, and autumn. However, it is noteworthy that a non-negligible fraction of the MME cannot rule out GHG-driven drying, particularly for mean and extreme precipitation in the summer. Furthermore, for extremes in the summer, even the best estimate from the average across the MME suggests that GHG forcing results in drying. The slow precipitation response to aerosols (summarized by AER-glob forcing) is generally drying, although the effect on precipitation is generally smaller for Slow-AER versus Slow-GHG since the maximum forcing range is reduced (−0.89 Wm−2 for AER-glob versus +2.61 Wm−2 for GHG, see Supplementary Fig. 2).

A similar result holds for the MME estimates of the fast precipitation response to aerosols, where large numbers of ensemble members suggest that increased SO2 emissions lead to either increases or decreases in mean and extreme precipitation. The station-data-based 90% confidence intervals shown in Fig. 5 have at least some overlap with the central 90% of the weighted CMIP6 MME across all seasons, precipitation type (mean and extreme), and anthropogenic forcing agent. This implies that climate models are consistent with GHCN-based estimates, although the model uncertainty is so high that the sign of the trend cannot be discerned (whereas it can for GHCN due to lower uncertainty). The degree of consistency between climate models and observations is of course differentiated, with higher consistency for GHG-driven changes in extreme precipitation in winter and spring than local AER-driven changes in extreme summer precipitation.

Our analysis of the large CMIP6 MME reiterates that if one is to pursue D&A of regional precipitation with climate models (a pursuit we advise against with current models), it is important to use multiple climate models. The weighted summaries of individual ensemble members in Fig. 5 illustrate the broad range in GCM-simulated responses to anthropogenic forcing, meaning that individual climate models cannot confidently attribute even the sign of the effect of AER and GHG forcing on CONUS precipitation up through present day. This is particularly true for aerosols, although (as previously mentioned) in some cases individual ensemble members indicate that increases to GHG forcing result in drying. The implication is that D&A for regional precipitation based on a single-model can yield opposite results relative to identical assessments using other climate models in the MME, including GHG-driven drying and AER-driven moistening (sometimes “conclusively”; see Supplementary Fig. 9). Reassuringly, the weighted CMIP6 multimodel ensemble means are much more consistent with the observations (also shown in Fig. 5), although the ensemble means appear to be overconfident in their assessment of uncertainty [as has been observed elsewhere, e.g., 59]. A notable exception is in forced changes to summer precipitation, particularly for GHG-driven changes to extremes, where the multimodel mean indicates, seemingly erroneously, high confidence in GHG-induced CONUS-wide drying rather than moistening.

Methods
Analysis of GHCN in situ records
The central results of this paper are based on an analysis of in situ measurements of daily precipitation from the GHCN-daily database [GHCN; 53, 54]. Specifically, we use measurements from a set of n = 2480 high-quality weather station records in the contiguous United States (CONUS) that have at least 66.7% non-missing daily values over December 1, 1899 to November 30, 2020. This set of stations is a subset of the more than 21,000 total GHCN stations in the CONUS; the geographic distribution of the high-quality stations is shown in Supplementary Fig. 11. For each station, we calculate seasonal mean precipitation rates (mm day−1) and seasonal maximum daily precipitation rates (mm day−1, often referred to as Rx1Day) in each year so long as there are no more than 33.3% missing measurements in the season-year 29. Developed a general framework for modeling a timeseries of seasonal mean or maximum daily precipitation {P(t)}, where the temporal units refer to a year. For simplicity, we suppress notation for geographic location; we later apply this formula to seasonal precipitation from each of the GHCN gauged locations. Our approach statistically models precipitation as

P(t) = P0 + PF(t) + PD(t) + PW(t) (1)

where the time series is decomposed into a pre-industrial component P0 that represents an overall average uninfluenced by external forcing, an externally forced component PF(t), and an internal variability component comprised of a “driven” term PD(t) that describes year-to-year variability due to known modes of large-scale oceanic and atmospheric drivers together with a term PW(t) associated with short-term weather variability. One strength of this framework is that it explicitly characterizes part of the internal variability via PD(t) [as advocated by ref. 28]. The various analyses in ref. 29 (summarized in Supplementary Table 1) verified that the additive framework of Eq. (1) is appropriate and can safely be simplified when considering seasonal precipitation over the CONUS in the historical record (1900-present) as follows:

PF(t) ≈ βSlow FSlow(t, τSlow) + βFast FFast(t, τFast),
PD(t) ≈ βD d(t),
Var(PW(t)) = σW2(t) ≈ σW2 × exp{V1t}, (2)
Var(PW(t)) / Var(P(t)) ≈ Constant:

Historically, the principal anthropogenic forcing agents for CONUS precipitation are well-mixed greenhouse gases (GHGs) and aerosols (AER), which collectively define a fast and slow precipitation response:

FSlow(t, τSlow) = FGHG(t, τGHG) + FAER-glob(t, τAER-glob),
FFast(t, τFast) = FAER-local(t, τAER-local), (3)

Here, FGHG(t, τGHG) is the lagged GHG forcing time series shown in Supplementary Fig. 2(a), FAER-glob(t, τAER-glob) is the lagged aerosol effective radiative forcing [which applies globally, hence “AER-glob”; 47] shown in Supplementary Fig. 2(c), and FAER-local(t, τAER-local) is the regionalized (local, hence “AER-local”) SO2 emissions. All forcing time series are reconstructed from observations and considered fixed. Equation (2) specifies that the forced component can then be described by a linear sum of a coefficient βf(⋅) multiplied by the forcing time series Ff(⋅)(t). We henceforth denote βSlow and βFast as the “statistical attribution coefficients,” since in our approach these are the analog of scaling factors in an optimal fingerprinting analysis, e.g., ref. 35 (see Methods below for additional details). Furthermore, note that the attribution coefficients (which are estimated from data) translate the fixed forcing time series into the corresponding fast and slow precipitation response. Following ref. 29, we set τFast = 0 years, and for the GHCN analysis we use τSlow = 14 years [the CMIP6 multimodel ensemble average of the lagged response to GHG forcing; see Appendix A of 29]. Note that we use stochastically regionalized (i.e., spatially-varying) SO2 emissions data to characterize FFast(t, τFast); again see ref. 29.

Following ref. 67, the driven component PD(t) can be well-approximated by a linear function of climate drivers D, comprised of the El Niño Southern Oscillation (ENSO) Longitude Index (ELI), the Arctic Oscillation (AO), the North Atlantic Oscillation (NAO), the Pacific-North American pattern (PNA), and the Atlantic Multidecadal Oscillation (AMO). Note that, as described in ref. 67, the AO is excluded from the December/January/February (DJF) analysis due to its strong coupling (and high correlation) with the NAO in this season. The residual weather variability term PW(t) is modeled statistically as following either a Gaussian distribution for the seasonal mean (with mean zero and variance σ2(t)) or a Generalized Extreme Value distribution for the seasonal maxima (centered on zero with variability described by σ2(t) and time-invariant shape parameter). The seasonal time series of the drivers for the GHCN analysis are calculated from observations; see “Analysis of GHCN in situ records” in ref. 67. The three-way variability variance V0 and a time-varying quantity exp{V1t}, but the P(⋅) are otherwise statistically independent (i.e., we assume that all autocorrelation in the time series is fully captured by PF(⋅) and PD(⋅)). Finally ref. 29, shows that the signal-to-noise ratio of the time series is approximately constant.

Our analysis of the GHCN records proceeds in three steps, following the methodology developed in refs. 68 and 26. First, independently at each station, we obtain maximum likelihood estimates of all unknown statistical parameters in Eq. (2), including the statistical attribution coefficients βSlow and βFast. Second, utilizing all stations, we apply a spatial statistical approach using second-order nonstationary Gaussian processes to interpolate each statistical parameter to obtain best estimates of each field for a high-resolution 0.25°×0.25° longitude-latitude grid over CONUS. Third, we quantify uncertainty via resampling methods, specifically a block bootstrap for estimating standard errors and confidence intervals [as in ref. 68] and a permutation/reshuffling approach for assessing statistical significance and spatial patterns [as in ref. 26]. In describing the statistical significance of the statistical attribution coefficients, the attribution exercise proceeds by testing the null hypothesis H0,f : βf = 0, f ∈ {Slow, Fast} for either individual grid boxes or spatially-aggregated grid boxes (e.g., the nested attribution regions defined in 36). Rejecting H0,f implies that there is a significant relationship between forcing agent(s) f and seasonal precipitation, i.e., we can conclusively attribute changes in precipitation to the human activity described by agent f. Note that for a given spatial location we can only attribute the sum-total GHG and AER-glob slow precipitation response (via statistical inference on βSlow); however, we can separate the individual effect of GHG and AER-glob forcing on precipitation over time when assessing emergence times (see “Summarizing the GHCN analysis” below). As in refs. 26, 67, the statistical significance of each H0,f : βf = 0 is determined with both “moderate” and “strong” significance; furthermore, when a set of simultaneous tests is conducted (e.g., for a set of CONUS subregions) we apply a multiple testing adjustment that bounds the proportion of type I errors at 0.33 (for moderate significance) and 0.1 (for strong significance).

Comparison with optimal fingerprinting methods
In contrast with our approach, traditional D&A methods rely on optimal fingerprinting [OF;] 35], which regresses an observed quantity onto a linear combination of model-simulated responses to a set of external forcings. OF-based analyses comprise the large majority of attribution statements cited in the IPCC but, as mentioned in the introduction, have yielded largely inconclusive results for attributing changes to regional precipitation. The OF approach proposes the following statistical model:

Y = Σf∈F βf(Ff + ef) + e0 (4)

where Y represents an observed quantity, F is a set of relevant experiments involving one or more external forcings, the Ff are corresponding model-simulated quantity from experiment f (often called the “fingerprints” and typically representing some sort of model ensemble average), ef represents sampling uncertainty in the model-simulated F, e0 quantifies the effect of internal variability on the observations, and the βf are scaling factors used to make causal statements. Specifically, a fingerprint Ff is “detected” in the observations if the uncertainty limits on βf do not include zero, and Ff is furthermore attributable if the uncertainty limits on βf are consistent with one. In Eq. (4), Y and the Ff are typically vectors representing values indexed over space and/or time as anomalies from some baseline period. Typical applications explore either so-called “one-way” regression analysis wherein kFk = 1 (considering only an all-forcings experiment or only a single-forcing experiment) or “two-way” regression analysis with kFk = 2 where one attempts to separate the total anthropogenic and natural responses [using two single- or multiple-forcing experiments; see, e.g., ref. 27].

Our approach (Eq. (2)) can be seen as a special case of Eq. (4). First, we implement a two-way regression where we consider two individual forcing agents, namely GHGs and AER, separated into fast and slow response components; i.e., F={Slow, Fast}, where Slow = GHG+AER-glob and Fast = AER-local. Recall that 29 explicitly showed that a two-way regression is appropriate by verifying that (1) including only GHG and AER is sufficient for D&A over CONUS and all other external forcings are negligible; and (2) the total anthropogenic response can be represented by the sum of GHG and AER forcing (global and local). The observed vector in Eq. (4) represents yearly measurements of seasonal mean and maximum precipitation over time from a single geospatial location, i.e., Y={P(t): t=1,…,T}. Unlike Eq. (4), in our implementation the Ff are not model-simulated precipitation but instead fixed forcing time series that are reconstructed from observations. The mechanisms for how GHG forcing impacts precipitation are relatively well understood and the forcing time series FGHG(⋅) has relatively low uncertainty; however, unlike GHG forcing, the AER-glob forcing time series has non-negligible uncertainty, such that we must account for non-zero ef. Similarly, while the SO2 emissions are well-observed, their influence on precipitation is less certain, and so we must account for non-zero ef. Both the probability-weighted trajectories of AER-glob forcing for [the 5th, 16th, 50th, 84th, and 95th percentile trajectories; see Supplementary Fig. 2 and 47, 48] and (2) an ensemble of regionalized emissions for e29. Since the units of the Ff are no longer mm day−1, the βf are no longer unitless (as they are in Eq. (4)) and instead have units mm day−1 per unit increase in the forcing time series. Therefore, in our case the D&A exercise is reduced to assessing a single null hypothesis test (i.e., H0,f : βf = 0) for each forcing agent, while explicitly modeling internal variability via the driven component, as advocated in ref. 28: while the elements of e0 = [e0(1),…,e0(T)] describe all non-externally-forced variability, in our approach the non-forced variability is explained by both the driven term PD(t) and the weather variability term PW(t). In other words, one can relate e0(t) = PD(t) + PW(t), such that by definition we have Var(PW(t)) ≤ Var(e0(t)). Since uncertainty in βf is largely a function of the magnitude of internal variability (here Var(PW(t)) vs. Var(e0(t))), our approach will increase the signal-to-noise ratio relative to a more traditional approach. Unlike Eq. (4), where the variance-covariance matrix of e0 is estimated offline from a set of pre-industrial control runs, the variance-covariance matrix of the weather variability PW(t) is estimated using only the observations and simultaneously with the βf. Finally, our methodology proposes that the βf is not a scalar quantity but instead a spatial process, i.e., varying across the spatial domain of interest. Critically, this allows us to estimate different β values for a high-resolution 0.25°×0.25° grid, such that D&A statements can be made for either individual grid boxes or aggregated grid boxes in a single framework. In other words, it is trivial to generate spatially aggregated D&A statements across various spatial scales.

Summarizing the GHCN analysis
While the results presented in this paper focus on the statistical attribution coefficients from the forced component, i.e., βSlow and βFast, for the sake of interpretability and cross-comparison we convert each of these values to an effect on precipitation. These summaries directly account for uncertainty in the aerosol forcing time series FAER-glob and FAER-local. For each season and precipitation type, the relevant output of our statistical analysis are estimates of the attribution coefficients

{βSlowa,r(g), βFasta,r(g): a=1,…,5; r=1,…,100} (5)

at each grid cell g, for AER-glob trajectory a=1,…,5 and regionalized SO2 emissions trajectory r=1,…,100. (Note that there are corresponding estimates for bootstrap and permutation resampling; these are used to quantify uncertainty.) For the “Slow-AER” (AER-glob) maps shown in Fig. 2, each grid cell shows the best estimates

ΔP^GHCN(g) = 1/100 Σa=1 5 wa Σr=1 100 βSlowa,r(g), (6)

where ΔAER-glob is the change in lagged AER-glob forcing for 1900 vs. 2010 with τSlow = 14 years for each of the five trajectories shown in Supplementary Fig. 2, the wa = ϕ(pa) are weights derived from the standard Normal probability density function (where {pa} = (0.05, 0.16, 0.5, 0.84, 0.95)), and βSlowa,r(g) are from Eq. (5). The reference years 1900 and 2010 are chosen since these are when AER-glob forcing is at its minimum and maximum, respectively, in the time period we are analyzing. For the “Slow-GHG” maps shown in Fig. 2, each grid cell shows the best estimates

ΔP^GHCN(g) = ΔGHG 1/100 Σa=1 5 wa Σr=1 100 βSlowa,r(g), (7)

where ΔGHG = 2.61 Wm−2 is the change in lagged GHG forcing for 1900 vs. 2020 with τSlow = 14 years and βSlowa,r(g) are from Eq. (5). The reference years 1900 and 2020 are chosen since these are when GHG forcing is at its minimum and maximum, respectively, in the time period we are analyzing. For the “Fast-AER” (AER-local) maps, each grid cell shows

ΔP^GHCN(g) = 1/100 Σa=1 5 wa Σr=1 100 ΔrAER-local(g) βFasta,r(g), (8)

where ΔrAER-local(g) is the change in the rth stochastically-regionalized SO2 emissions trajectory in grid cell g for 1900 vs. 1966 and βFasta,r(g) are from Eq. (5). The reference years 1900 and 1966 are chosen since these are the years in which SO2 emissions were at their minimum and maximum, respectively, in the time period we are analyzing. Uncertainty assessments and stippling are based on applying Eqs. (6)–(8) to resampling-based estimates of these quantities; see ref. 26.

The signal emergence plots in Fig. 3 also show the effect of each forcing agent on precipitation but now over time and aggregated spatially. When assessing temporal changes, we can now explicitly separate trajectories over time due to GHGs, AER-glob, and AER-local. For a given collection of 0.25°×0.25° grid cells A, we compute area-averaged anomalies in year t using

P^Obs(t) = 1/|A| Σg∈A a(g)βbSlow(g) * [FSlow(t, τSlow) - FSlow(pi-clim, τSlow)]
P^Obs(t) = 1/|A| Σg∈A a(g)βbSlow(g) * [FSlow(t, τSlow) - FSlow(pi-clim, τSlow)]
P^Obs(t) = 1/|A| Σg∈A a(g)βbFast(g) * [FFast(t, τFast) - FFast(pi-clim, τFast)], (9)

where a(g) is the area of grid cell g, |A| = Σg∈A a(g), and “pi-clim” refers to a pre-industrial climate (represented by the 1900-1929 average). For brevity, we omit the dependence of the attribution coefficient estimates βb(⋅) on forcing trajectory, but note that the above calculations involve (weighted) averages of trajectory-specific estimates and forcing anomalies similar to, e.g., Eq. (6). The sum-total anthropogenic forcing shown in Fig. 3, denoted ANT, is simply

P^Obs(t) = P^Obs(t) + P^Obs(t) + P^Obs(t): (10)

Basic bootstrap confidence intervals are calculated by applying Eq. (9) to bootstrap estimates of the forcing coefficients.

Analysis of the CMIP6 historical multimodel ensemble
To compare our GHCN analysis with corresponding estimates from Global Climate Models (GCM), we apply the D&A formula described by Eq. (2) to individual ensemble members of each GCM in the CMIP6 historical experiment, matching the period covered by the in situ records (i.e., 1900-present). Note that Eq. (2) is applicable to the GCMs precisely because all hypotheses enumerated in Supplementary Table 1 were tested using the CMIP6 multimodel ensemble. The experimental protocol for these simulations prescribes external forcing agents that correspond to the historical period and hence the forcing time series FSlow and FFast are as in Supplementary Fig. 2 (but with GCM-specific values of τ). Since the historical runs are fully coupled, each ensemble member has its own set of driver time series: these are calculated via the Climate Variability Diagnostics Package [for everything except ELI; 69] and the Toolkit for Extreme Climate Analysis [for ELI; 70]. Our analysis of the historical ensemble members mirrors the GHCN analysis but without the use of spatial statistical methods: we simply (1) obtain maximum likelihood estimates of all unknown statistical parameters at each model grid cell and (2) utilize resampling methods to quantify uncertainty. After conducting these analyses, we are left with coefficient estimates βSlow and βFast as well as measures of uncertainty at each model grid cell and each ensemble member. To explore CONUS-wide changes, we then calculate area-weighted averages of all coefficient estimates and their uncertainties, denoted

{βf(j,m), σf(j,m): f ∈ {Slow, Fast}, j=1,…,nm}, (11)

where m=1,…,M indexes climate models and j=1,…,nm indexes the ensemble members from model m (for brevity we again omit the dependence of the attribution coefficient estimates βb(⋅) on forcing trajectory). Recall that these quantities are calculated separately for each season and precipitation type (mean and extreme). As with the GHCN analysis, for plotting we convert the coefficient estimates and uncertainties to a precipitation response following Eqs. (6)–(8).

Model weighting. Our philosophy for combining estimates of βSlow and βFast and their effect on seasonal precipitation across the CMIP6 historical multimodel ensemble emphasizes internal consistency of climate models for deriving weights:

1. Down-weight ensemble members j from a given model m with values of βSlow and βFast that could arise purely by chance, due to internal variability [as quantified by comparing with estimates from CMIP6 Diagnosis, Evaluation, and Characterization of Klima (DECK) pre-industrial control runs; 58], and
2. Up-weight ensemble members j from a given model m with values of βSlow that are consistent with estimates from corresponding runs of transient CO2-only forcing [i.e., the CMIP6 DECK 1pctCO2 runs; 58].

All comparisons are made within-model, and note that we do not attempt to specify which climate models are “better” or “worse” with respect to observations. Furthermore, this approach to model weighting emphasizes trends (as opposed to mean climatologies) which is most relevant for this exercise.

Pre-industrial control (piControl) fitting. To deprecate “false positive” detection of trends due to anthropogenic forcing agents, we apply Eq. (2) to i=1,…,N overlapping 121-year segments of pre-industrial control runs (i.e., esm-piControl and/or piControl experiments) from each model (to mirror the length of the GHCN observational record). In each of these fits, we use the observed Slow and Fast forcing time series (from the historical 1900-2020 period; we again use stochastically-regionalized SO2 emissions) but drivers corresponding to the realized conditions in each ensemble member. As with the historical estimates, we apply the formula in each model grid cell and then obtain area-weighted CONUS averages. This procedure detects false positives since the anthropogenic forcing agents are not actually present in the simulations; the piControl estimates reveal the magnitude of secular 121-year trends that can arise from internal climate variability. Let {βPI,f(i,m), σPI,f(i,m)} represent pairs of area-weighted, CONUS-averaged coefficient estimates and bootstrap standard errors for forcing f and segment i=1,…,N of model m. Then, we calculate a single effect (with uncertainty) for each model as:

βPI,f(m) = 1/N Σi=1N βPI,f(i,m), σPI,f(m) = sqrt(1/N2 Σi=1N [σPI,f(i,m)]2): (12)

Note that our inclusion criterion is that the esm-piControl and/or piControl runs have at least 500 years of data; we then select the last 500 years of each run, which allows us to fit N=19 overlapping segments of 121 years.

The first component of the model ensemble weights is then one minus the inverse squared exponential of the standardized difference between the piControl estimates (which have arisen purely by chance, from Eq. (12)) and the historical estimates (from Eq. (11)):

wPI,f(j,m) = 1 - exp{-ZPI,f(j,m)2}, (13)

where

ZPI,f(j,m) = (βf(j,m) - βPI,f(m)) / sqrt(σf(j,m)2 + σPI,f(m)2): (14)

Transient 1pctCO2 fitting. Next, we prioritize climate models with historical ensemble members that are internally consistent with Slow coefficients estimated from transient CO2-only runs. Here, we again apply Eq. (2) to 1pctCO2 runs (again applied to each grid cell followed by calculating area-weighted averages over CONUS), now considering the initialized ensemble members up through the year in which the lagged forcing from increased CO2 equals the lagged forcing from the collection of GHGs in the historical run, using the lag time-constant τSlow specific to each GCM. In the fitting formula, the Slow forcing time series FSlow(t, τGHG) is set to be the lagged CO2-only forcing (i.e., the AER-glob forcing is omitted since these runs have no anthropogenic aerosols), the Fast forcing FFast(t, τFast) and βFast are set to 0, and the drivers correspond to the real conditions within each ensemble member. This fitting procedure yields {β1%f(m), σ1%f(m)} for each model. For climate models that provide data for more than one ensemble member, we average over the ensemble members as in Eq. (12). The second component of the model ensemble weights is then

w1%f(j,m) = exp{-Z1%f(j,m)2}, (15)

where

Z1%f(j,m) = (βf(j,m) - β1%f(m)) / sqrt(σf(j,m)2 + σ1%f(m)2): (16)

Model ensemble weighting. Combining all of the above, we obtain a weight for each ensemble member j=1,…,nm of the model m=1,…,M:

w(j,m) = ν w1%f(j,m) wPI,f(j,m) s.t. 1/nm Σj w(j,m) = 1 (17)

where ν is a normalization factor that enforces the right-hand condition (recall f ∈ {Fast, Slow}). Note that including the 1/nm in the summations gives each model equal weight; again, recall that these weights are calculated separately for each season and precipitation type (mean and extreme). The number of ensemble members we have for each of the historical, 1pctCO2, and piControl experiments are shown in Supplementary Table 2. A total of M=25 climate models have the requisite data for these three experiments; across these models, we have 316 historical ensemble members (i.e., Σm=1M nm = 316). Individual weights (wPI,f(j,m) and w1%f(j,m)) are plotted against corresponding estimates of each forcing’s influence on precipitation (see Eqs. (18) and (19)) in Supplementary Fig. 9. The overall (non-normalized) weights w(j,m) are shown in Supplementary Fig. 10.

Comparing weighted climate models with observations. To summarize fits from the individual historical ensemble members, we convert the historical coefficient estimates into their effect on precipitation, averaged over CONUS, as follows. First, in each season and for each precipitation type, the slow precipitation response for ensemble member j of model m is

PSlow(j,m) = ΔSlow(m) Σg=1Gm a(g)βSlow(j,m,g) / Σg=1Gm a(g), (18)

where ΔSlow(m) is the change in lagged GHG and AER-glob forcing for 1900 vs. 2014 for model m (this quantity is model-specific since each model has its own lag τSlow), a(g) is the area of grid cell g=1,…,Gm of model m, and βSlow(j,m,g) is the Slow coefficient estimate in each grid cell g. The reference years 1900 and 2014 are chosen since these are when the lagged GHG and AER-glob forcing is at its minimum and maximum, respectively, in the time period bounded by the start of the GHCN record and end of the CMIP6 historical simulations. Next, again in each season and for each precipitation type,

PFast(j,m) = Σg=1Gm ΔFast(m,g) a(g)βFast(j,m,g) / Σg=1Gm a(g), (19)

where ΔFast(m,g) is the change in stochastically-regionalized SO2 emissions in the grid cell g of the model m for 1900 vs. 1966. Note that ΔFast(m,g) also depends on the season. The reference years 1900 and 1966 are chosen since these are the years in which SO2 emissions were at their minimum and maximum, respectively, in the period of interest. For each season and precipitation type across all ensemble members, we obtain the best estimates, a lower 90% confidence bound, and an upper 90% confidence bound, denoted

{P^f(j,m), blf(j,m), ubf(j,m)}, (20)

of the effect on precipitation for f ∈ {Slow, Fast} after applying Eqs. (18) and (19) to best estimates and bootstrap estimates of each coefficient. In light of the large number of ensemble members we need to summarize (across seasons, precipitation type, and forcing), we define

hf(x) = Σm=1M Σj=1nm w(j,m) × I(x - blf(j,m)) × I(ubf(j,m) - x) (21)

(notation for season and precipitation type are suppressed), where

I(z) = 1 if z ≥ 0, 0 otherwise. (22)

The function hf(x) summarizes the weighted proportion of ensemble members for which x is included in their 90% confidence interval and ranges between 0 and 1 by construction. Here x is the anthropogenically forced change in precipitation rate caused by the forcing agent(s) f. Given this concise summary of the model results across seasons for each precipitation type, we can plot hf(x) for a range of x values and compare the resulting curve with observational estimates corresponding to Eq. (20) obtained from the GHCN analysis, denoted

{P^f(GHCN), blf(GHCN), ubf(GHCN)}: (23)

Figure 5 shows our results, comparing fits obtained from the CMIP6-historical multimodel ensemble with those from the in situ records in the GHCN. Note that we have applied a spline-smoothing to the hf(⋅) curves for visual appeal.

Fast versus slow precipitation response to aerosols
The analysis in ref. 43 presents a multimodel assessment of the fast and slow precipitation response to individual climate forcings, including sulfate aerosols (SO4). In the main text, their Fig. 3 shows geographical patterns of multimodel mean precipitation change (fast, slow, and total), while Fig. 4 shows the same quantities for various spatial aggregations, including continental averages. However, their analysis only looked at annual mean precipitation, and furthermore from their Fig. 3 it appears that the North America land average summarized in Fig. 4 is heavily influenced by the precipitation responses in Canada, which is of course not included in our CONUS-specific analysis. To explicitly evaluate the geographical patterns and CONUS-wide means of seasonal mean and extreme precipitation responses, we repeat the analysis described in ref. 43. Output from four modeling experiments is utilized, all of which are part of the Precipitation Driver and Response Model Intercomparison Project [PDRMIP; 57]:

1. base-fsst: the external forcings correspond to present-day conditions and the sea-surface temperatures are fixed.
2. sulx5-fsst: the external forcings correspond to present-day conditions except for sulfate aerosol concentrations, which are multiplied by a factor of 5, and the sea-surface temperatures are fixed.
3. base-coupled: the external forcings correspond to present-day conditions and the sea-surface temperatures are prognostic.
4. sulx5-coupled: the external forcings correspond to present-day conditions except for sulfate aerosol concentrations, which are multiplied by a factor of 5, and the sea-surface temperatures are prognostic.

The eight GCMs and number of simulated years used are given in Supplementary Table 3. Following ref. 43, we utilize the last 10 years of the fsst experiments and the last 50 years of the coupled experiments; we calculate the seasonal mean daily precipitation and seasonal maximum daily precipitation (Rx1Day) in each grid box of each model-experiment. The fast precipitation response is calculated as the mean of sulx5-fsst minus the mean of base-fsst; the total precipitation response is calculated as the mean of sulx5-coupled minus the mean of base-coupled; and the slow precipitation response is the total precipitation response minus the fast precipitation response. To assess geographical patterns of change, all model results are conservatively remapped to a 1°×1° longitude-latitude grid. For all quantities we calculate a 90% bootstrap confidence interval to summarize uncertainty. Our results are shown in Supplementary Fig. 5.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Review of the revised version of “Anthropogenic aerosols mask increases in US rainfall by greenhouse gases” for Nature Communications (initial submission to Nature)\n#reviewer 1\nI appreciate the authors’ responding to my comments, including both major ones on interpretation and minor ones on language etc. The revisions made together with the explanations in the response to reviewers alleviates the points I raised in the previous review. Besides the new text in 318ff and 499ff I also appreciate some of the changes made in response to the other reviewers, for instance, the “Comparison with optimal fingerprinting methods” section in the appendix and what is now Fig. 1, and think they improve the manuscript. The following few minor points remain.\nComment previous review\n- Fig. 2 what does the linewidth of the dashed lines show? > All dashed lines are the same width | Yes, sorry, that seems to be the way okular shows pdf figures, depending on how much I zoom in, some lines are thinner and thicker.\nMinor points new text\n- 32 I find “equal and opposite [aerosol and greenhouse gas signals]” super confusing, at least as non-native English speaker, does it mean “equal (in magnitude) and opposite (in sign)” (aerosols vs ghg) or “same sign (in some season) and different sign (in other seasons)” aerosols vs ghg or a third meaning I’m not getting? If there is a way to clarify this in the text (I’m aware of the abstract’s word limit), that would be good.\n- 35 these emission sources >each< increase or decrease rainfall?\n- 60-61 “meaningful”, unclear what that means (having impacts? statistically significant? detectable?)\n- 66-72 This new argument needs to be better linked. Firstly, I think this “even” is opaque, especially since in line 67 “aerosol” is the first time mentioned in the text (discounting the abstract), so how ant. aerosol and human influence are linked has not been introduced. Secondly I disagree with the logic when comparing the SNR from individual simulations with that for the multi-model mean, since both internal variability and model uncertainty play a role in that comparison, so from this alone you can’t conclude on the model uncertainty contribution I think.\n- 111 what does this “otherwise” mean?\n- 153 I’d remove “here” since the reader might not be in the conus…\n- 154 characterise so2 – is very broad a term, I think you mean characterising some specific aspects of so2 that are relevant for its influence on climate, perhaps consider specifying.\n- 162 sounds a bit weird since something that is prescribed IMO doesn’t count as “showing agreement across” models, perhaps rephrase\n- 232, 237, 238 “vs” here I find quite confusing\n- Tab. 1. Good idea to synthesise this. Though if I understand this correctly, this now includes any grid boxes, regardless of whether the statistical attribution is significant. Must be, since the significance hatching differs for aerosols and GHG grid boxes. I suggest it would be more meaningful if the authors could add a third option for each AER and GHG so it’s “- | + | nonsignificant” or “- | + | <0.3mm/day” (or 0.1 at least, what’s white in Fig. 2. Then I also suggest to make the diagonal entries (-/-, +/+, n/n) in bold or perhaps rather the ones where the signals clearly counteract (-/+ and +/-) to make the table easier to read again…\n- 354 typo precipitaiton\n- 424-428 similar point as above to 66-72: does the broad range of individual ensemble members really illustrate the broad range in responses to anthropogenic forcing meaning that individual climate models [not: climate model runs?] cannot confidently attribute [...]?\n431 “wrong” seems a strong word here, if “right” is the result derived w/ your forcing-only (rather than model) based method? perhaps specify instead what you are comparing to, or use a more neutral term like “different”\nKind regards, Sabine Undorf", "reviewer_2": "Thanks for the responses from the authors, which considerably improved the manuscript. But I don’t think the authors completely addressed my previous major comment #1 on the remote aerosol effect and comment #2 on the fast vs slow responses. These two hypotheses basically constitute the foundation of their framework. I don’t think this manuscript can be accepted in Nature Communication until they are well justified.\n1. To demonstrate the marginal effect of remote aerosol on the CONUS precipitation, the author generated Fig. B2 showing that the CONUS-average responses to the Asian aerosols are statistically indistinguishable from zero. Thanks for the work here but I still have a couple of remaining questions here.\n(a) The authors mainly use the bootstrap method and permutation approach (for spatial pattern) to estimate the uncertainty based on their published paper. I think the authors should briefly summarize how they conducted these two methods and why they fit to make the reader understand the contents easily rather than refer to their paper to learn the complicated statistical contents in much detail. For instance, how the bootstrap is conducted for Fig. B.2.\n(b) The authors just show the continental-scale precipitation mean changes. What about the spatial distribution of precipitation responses to aerosols? There could be strong increases and decreases on a regional scale. Averaging these anomalies could smooth out the aerosol signals.\n(c) As the authors pointed out in this study, the aerosol influences on precipitation are in comparable magnitude to that of GHGs. The changes indistinguishable for zero don’t necessarily mean that the effects of remote aerosols are not important. It may be only because of the large uncertainty in the presentation of aerosol impact in climate models. For example, the changes in HadGEM3 are not marginal to me, as well as in some seasons in other models. Discussion on this should also be made.\n(d) I understand that the Asian aerosol level is low during 1900–1966, which may not affect the detection and attribution strongly in this work. Instead, the European aerosols could be important in affecting the CONUS precipitation. Looking at the influence of European aerosols (e.g., black carbon and sulfate) on the CONUS precipitation is essential.\n2. For the fast vs slow responses, the authors mentioned that the argument in Richardson et al. (2018) supports the hypothesis that the sulfate aerosol-induced fast precipitation changes dominate on the global land scale, which is however not cited in the revised manuscript. I agree that fast responses are more important than slow responses in land-mean precipitation changes supported by Fig. 4 in Samset et al. (2016) and Fig. 8a Richardson et al. (2018). But this argument may vary regionally as seen in Fig. 4 in Samset et al. (2016) (as demonstrated by the authors also), particularly in the CONUS, where is the focus of this study. The slow responses are far more important than fast responses in the North America (NAM) region. Again, this means the authors probably have to consider the lagged effects of sulfate aerosols rather than omit them.\nSpecific comments:\nL242: please provide details on how you define 20-year return values of extreme daily precipitation.\nReference\nRichardson, T. B., Forster, P. M., Andrews, T., Boucher, O., Faluvegi, G., Fläschner, D., et al. (2018). Drivers of precipitation change: An energetic understanding. Journal of Climate. https://doi.org/10.1175/JCLI-D-17-0240.1\nSamset, B. H., Myhre, G., Forster, P. M., Hodnebrog, Andrews, T., Faluvegi, G., et al. (2016). Fast and slow precipitation responses to individual climate forcers: A PDRMIP multimodel study. Geophysical Research Letters. https://doi.org/10.1002/2016GL068064", "reviewer_3": "I appreciate the authors’ effort to consider my previous comments. Although the revised manuscript provides more explanations and arguments, unfortunately, I find my major concerns not addressed clearly.\nFirst, authors’ intentional choice of not using GCMs to estimate fingerprints would help reduce uncertainties in AER fingerprints, but directly linking SO2 emissions with local precipitation is not supported well by physical processes, implying the large uncertainty remaining in their ‘attribution’ results to AER. In addition, authors describe similarities and differences of their d&a method compared to traditional d&a methods (OF). This helps to understand the authors’ method, but the main question about how ‘noise’ of scaling factors is estimated and boosts S/N remains unclear. As I understand from their responses, authors seem to use weather variability (P_W) estimated from the observations as a major component of noise but OF usually considers P_D (its long-term variations). This indicates that the proposed d&a method will have larger noise than OF and therefore increases S/N by ignoring signal uncertainties due to model differences (i.e. using forcing data to estimate signal) rather than by reducing noise.\nSecond, main conclusions “AER masks rainfall increase by GHG in US” remain not robust. Responding to reviewer #1’s concerns, authors removed “resolve outstanding uncertainties regarding the human influence on regional precipitation” or similar interpretations throughout their manuscript (except Abstract?). This clarifies that this study does not advance our understanding of ‘regional scale’ d&a of precipitation changes, which is suggested as the aim. Also, authors still apply different periods when attributing the observed changes to GHG and AER (Figures 1, 2, and 4, and Table 1 as well) such that 1900-1966 trends are used for AER and 1900-2020 trends for GHG. Before and after 1966, we can expect difference contribution by AER as reviewer #1 pointed. Under this different timing, how can they conclude that AER mask GHG influences? Their important motivation is the large uncertainty in Figure SPM.3 of IPCC AR6 but that figure considered the observed trends since 1950. Thus, results from this study cannot be compared properly with the IPCC results.\nLastly, another important finding was suggested as attribution at local scales, particularly AER influences as the title says. However, major conclusions are obtained based on US averages and local scale D&A is found very limitedly in space and time. They even say that “The AER signal is essentially non-existent outside JJA when considering spatial scales smaller than CONUS.” Given the limited impact of AER (aside from the different period issue), how can we conclude that “AER signal masks GHG signal in US”? In terms of physical mechanisms for AER attribution in JJA, authors only speculate that “convective invigoration by SO2” may have worked, without showing clear evidences.\nSome specific points:\nL28-30: “Aerosol emissions offset these increases in the winter and spring but appear to enhance rainfall during the summer and fall.” This will be dependent on periods considered. I think authors need to consider different periods and quantify GHD and AER contributions to the observed trends to clarify this issue.\nL74-76: Proposed method does not “disentangle the complex causes of regional precipitation change”.\nL84-87: “this decomposition allows us to conclusively attribute changes to these forcing agents” This decomposition simply estimates contribution from each forcing, not conclusively attributing \"observed\" change to \"each\" forcing.\nL381-385: “In summary, our results show that uncertainties related to the emergence of a detectable and attributable human influence on regional precipitation over the CONUS over the historical record [Figure SPM.3; 31] can in most cases be explained by decomposing the combined ANT signal into its two components, GHG and AER forcing.” Results do not support this solid conclusion for regional attribution. Also, Figure SPM3 considers heavy precipitation trends since 1950 but the authors’ analysis covers different periods, particularly for AER.\nLine 617-620: “we quantify uncertainty via resampling methods, specifically a block bootstrap for estimating standard errors and confidence intervals [as in 69] and a permutation/reshuffling approach for ascribing statistical significance to spatial patterns [as in 26].” This is the important procedure for estimating uncertainty in scaling factors so authors need to explain details in comparison with traditional d&a methods. How does this consider internal variability, particularly P_D and P_W?\nL630-631: “the proportion of type I errors at 0.33 (for moderate significance) and 0.1 (for strong significance)”. Figure 2 shows that grid-box attribution results for AER have ‘moderate’ significance only in all cases (almost no grid shows ‘strong’ significance), indicating AER signal is not detectable when applying 5% level as used in OF. Is there any specific reason to use the lower significance level?\nLine 676-679: “Therefore, in our case the D&A exercise is reduced to assessing a single null hypothesis test (i.e., rejecting H0,f : βf = 0 both detects and attributes forcing f), such that larger βf (in absolute value) implies a stronger influence of forcing f.” I think this is equivalent to the ‘detection’ of signals in OF. Can authors think about another hypothesis test for ‘attribution’ by which they check agreement of the forced trend with the observed trend. I think this sort of quantification of forcing contribution is needed for robust ‘attribution’."}, "author_response": "Comments from Reviewer #1\n97\nones on language etc. The revisions made together with the explanations in the response to reviewers allevi-\n98\nate the points I raised in the previous review. Besides the new text in 318ff and 499ff I also appreciate some\n99\nof the changes made in response to the other reviewers, for instance, the “Comparison with optimal finger-\n100\nprinting methods” section in the appendix and what is now Fig. 1, and think they improve the manuscript.\n104\n• Fig. 2 what does the line width of the dashed lines show? “All dashed lines are the same width” –\n106\nsome lines are thinner and thicker.\n109\n• 32 I find “equal and opposite [aerosol and greenhouse gas signals]” super confusing, at least as a\n110\nnon-native English speaker, does it mean “equal (in magnitude) and opposite (in sign)” (aerosols vs\n112\nthird meaning I’m not getting? If there is a way to clarify this in the text (I’m aware of the abstract’s\n113\nword limit), that would be good.\n115\nequal in magnitude and opposite in sign. We have replaced “equal and opposite” with “offsetting,”\n116\nwhich we hope is more clear.\n117\n• 35 these emission sources each increase or decrease rainfall?\n119\nhave changed this sentence to “whether each emissions source increases or decreases rainfall.”\n120\n• 60-61 “meaningful”, unclear what that means (having impacts? statistically significant? detectable?)\n122\ncant.”\n123\n• 66-72 This new argument needs to be better linked. Firstly, I think this “even” is opaque, especially\n124\nsince in line 67 “aerosol” is the first time mentioned in the text (discounting the abstract), so how ant.\n125\naerosol and human influence are linked has not been introduced. Secondly, I disagree with the logic\n126\nwhen comparing the SNR from individual simulations with that for the multimodel mean, since both\n127\ninternal variability and model uncertainty play a role in that comparison, so from this alone you can’t\n128\nconclude on the model uncertainty contribution, I think.\n129\n130\nfocused the argument on the behavior of individual ensemble members; see lines 69-74.\n131\n• 111 what does this “otherwise” mean?\n133\n• 153 I’d remove “here” since the reader might not be in the CONUS...\n135\n• 154 characterize so2 – is very broad a term, I think you mean characterizing some specific aspects of\n136\nso2 that are relevant to its influence on climate, perhaps consider specifying.\n137\n• 162 sounds a bit weird since something that is prescribed IMO doesn’t count as “showing agreement\n139\nacross” models, perhaps rephrase\n141\nare hence consistent across climate models, ...”\n142\n• 232, 237, 238 “vs” here I find quite confusing\n144\nresponses.\n145\n• Tab. 1. Good idea to synthesize this. Though if I understand this correctly, this now includes any grid\n146\nboxes, regardless of whether the statistical attribution is significant. Must be, since the significance\n147\nof hatching differs for aerosols and GHG grid boxes. I suggest it would be more meaningful if\n148\nthe authors could add a third option for each AER and GHG, so it’s “ +  non-significant” or\n149\n“+ 0.3mm/day” (or 0.1 at least, what’s white in Fig. 2. Then I also suggest making the\n150\ndiagonal\nentries (-/-, +/+, n/n) in bold or perhaps rather the ones where the signals clearly counteract (-/+ and\n151\n+/-) to make the table easier to read again...\n153\nbusy (e.g., we’d need to summarize where one of Fast vs. Slow is significant, as well as both), we’ve\n154\nelected to keep the table as is. The new Figure 4 summarizes the overall effect of GHGs and aerosols\n155\nas well as their combined significance.\n156\n• 354 typo precipitaiton\n158\n• 424-428 similar point as above to 66-72: does the broad range of individual ensemble members\n159\nreally illustrate the broad range in responses to anthropogenic forcing meaning that individual climate\n160\nmodels [not: climate model runs?] cannot confidently attribute [...]?\n162\n• 431 “wrong” seems a strong word here, if “right” is the result derived w/ your forcing-only (rather\n163\nthan model) based method? perhaps specify instead what you are comparing to, or use a more neutral\n164\nterm like “different”\n167\n170\nbeen made much better as a result of your comments.\n173\nthe authors completely addressed my previous major comment #1 on the remote aerosol effect and com-\n174\nment #2 on the fast vs slow responses. These two hypotheses basically constitute the foundation of their\n175\nframework. I don’t think this manuscript can be accepted in Nature Communication until they are well\n176\njustified.\n178\nto address your feedback from the previous review. We refer you to summary of the large changes made\n179\ndescribed in our response to the AE above. We feel that the manuscript and results are greatly improved\n180\nthanks to your suggestion to model the slow precipitation response to aerosols and hope you do as well.\n181\n1. To demonstrate the marginal effect of remote aerosol on the CONUS precipitation, the author gen-\n182\nerated Fig. B2 showing that the CONUS-average responses to the Asian aerosols are statistically\n183\nindistinguishable from zero. Thanks for the work here, but I still have a couple of remaining ques-\n184\ntions here.\n185\n(a) The authors mainly use the bootstrap method and permutation approach (for spatial pattern) to\n186\nestimate the uncertainty based on their published paper. I think the authors should briefly sum-\n187\nmarize how they conducted these two methods and why they fit to make the reader understand\n188\nthe contents easily rather than refer to their paper to learn the complicated statistical contents in\n189\nmuch detail. For instance, how the bootstrap is conducted for Fig. B.2.\n190\n(b) The authors just show the continental-scale precipitation mean changes. What about the spatial\n191\ndistribution of precipitation responses to aerosols? There could be strong increases and decreases\n192\non a regional scale. Averaging these anomalies could smooth out the aerosol signals.\n193\n(c) As the authors pointed out in this study, the aerosol influences on precipitation are in comparable\n194\nmagnitude to that of GHGs. The changes indistinguishable for zero don’t necessarily mean that\n195\nthe effects of remote aerosols are not important. It may be only because of the large uncertainty\n196\nin the presentation of aerosol impact in climate models. For example, the changes in HadGEM3\n197\nare not marginal to me, as well as in some seasons in other models. Discussion on this should\n198\nalso be made.\n199\n(d) I understand that the Asian aerosol level was low during 1900-1966, which may not affect the\n200\ndetection and attribution strongly in this work. Instead, the European aerosols could be important\n201\nin affecting the CONUS precipitation. Looking at the influence of European aerosols (e.g., black\n202\ncarbon and sulfate) on the CONUS precipitation is essential.\n204\ntion of the Asian aerosol analysis, including Supplemental Figure B.2.\n205\n2. For the fast vs slow responses, the authors mentioned that the argument in Richardson et al. (2018)\n206\nsupports the hypothesis that the sulfate aerosol-induced fast precipitation changes dominate on the\n207\nglobal land scale, which is however not cited in the revised manuscript. I agree that fast responses are\n208\n6\nmore important than slow responses in land-mean precipitation changes supported by Fig. 4 in Samset\n209\net al. (2016) and Fig. 8a Richardson et al. (2018). But this argument may vary regionally as seen\n210\nin Fig. 4 in Samset et al. (2016) (as demonstrated by the authors also), particularly in the CONUS,\n211\nwhere is the focus of this study. The slow responses are far more important than fast responses in the\n213\neffects of sulfate aerosols rather than omit them.\n215\naerosols. A version of the Samset et al. (2016) annual mean precipitation result (reproduced in our\n216\npaper, see “Fast versus slow precipitation response to aerosols” in the Methods and also Supplemental\n218\nespecially in the summer, which was something we previously had trouble interpreting.\n220\n• L242: please provide details on how you define 20-year return values of extreme daily precipitation.\n222\nof Figure 2.\n225\nmore explanations and arguments, unfortunately, I find my major concerns not addressed clearly.\n228\nin AER fingerprints, but directly linking SO2 emissions with local precipitation is not supported well by\n229\nphysical processes, implying the large uncertainty remaining in their ‘attribution’ results to AER. In addition,\n230\nauthors describe similarities and differences of their d&a method compared to traditional d&a methods (OF).\n232\nis estimated and boosts S/N remains unclear. As I understand from their responses, authors seem to use\n233\n234\nweather variability (PW) estimated from the observations as a major component of noise but OF usually\n235\nconsiders PD (its long-term variations). This indicates that the proposed d&a method will have larger noise\nthan OF and therefore increases S/N by ignoring signal uncertainties due to model differences (i.e. using\n236\nforcing data to estimate signal) rather than by reducing noise.\n238\nresponse to aerosols (and favorable comparison with PDRMIP simulations) resolves some of your concerns\n239\nwith the uncertainty assessments in our attribution statements.\n241\nto reviewer #1’s concerns, the authors removed “resolve outstanding uncertainties regarding the human influ-\n242\nence on regional precipitation” or similar interpretations throughout their manuscript (except the Abstract?).\n244\nchanges, which is suggested as the aim. Also, authors still apply different periods when attributing the ob-\n245\nserved changes to GHG and AER (Figures 1, 2, and 4, and Table 1 as well) such that 1900–1966 trends are\n246\nused for AER and 1900–2020 trends for GHG. Before and after 1966, we can expect difference contribution\n247\nby AER, as reviewer #1 pointed. Under this different timing, how can they conclude that AER mask GHG\n248\ninfluences? Their important motivation is the large uncertainty in Figure SPM.3 of IPCC AR6, but that fig-\n249\nure considered the observed trends since 1950. Thus, results from this study cannot be compared properly\n250\nwith the IPCC results.\n251\n252\npretation of the attribution coefficients, and All of our results (including Figures 1, 2, and 4, and\n253 Fast Slow.\nTable 1) are based on estimates of the attribution coefficients obtained from analyzing the entire 121-year\n254\nrecord. As such, there is no difference in time period for the AER trends vs. the GHG trends. The 1900-\n255\n1966 for AER and 1900-2020 for GHG time periods are simply a way to rescale the attribution coefficients\n256\nso that (1) they are in sensible/interpretable units (mm/day) and (2) we can compare the maximum effect\n257\nof each forcing agent in units that make sense. We have changed various things to drive this point home,\n258\ne.g., the label on the color bar in Figure 2 is “Change in precipitation (mm/day).” Also, the new Figure 4\n259\nshows how the aerosol masking impacts GHG trends at a variety of spatial scales, from CONUS down to\n260\n261 05Mm2.\nLastly, another important finding was suggested as attribution at local scales, particularly AER influences as\n262\nthe title says. However, major conclusions are obtained based on US averages and local scale D&A is found\n263\nvery limited in space and time. They even say that “The AER signal is essentially non-existent outside JJA\n264\nwhen considering spatial scales smaller than CONUS.” Given the limited impact of AER (aside from the\n265\ndifferent period issue), how can we conclude that “AER signal masks GHG signal in the US”? In terms of\n266\nphysical mechanisms for AER attribution in JJA, authors only speculate that “convective invigoration by\n269\nat local scales. Additionally, the PDRMIP comparison of fast precipitation responses gives credence to our\n270\nobservational results.\n272\n• L28-30: “Aerosol emissions offset these increases in the winter and spring but appear to enhance\n273\nrainfall during the summer and fall.” This will be dependent on periods considered. I think authors\n274\nneed to consider different periods and quantify GHG and AER contributions to the observed trends to\n275\nclarify this issue.\n277\none in the abstract address the isolated effect of individual forcing agents (answering questions like\n278\n“when aerosol emissions increase, does precipitation increase or decrease?”).\n279\n• L74-76: Proposed method does not “disentangle the complex causes of regional precipitation change”.\n281\nresponses and the localized conclusions that can be made.\n282\n• L84-87: “this decomposition allows us to conclusively attribute changes to these forcing agents”\n284\nserved” change to ”each” forcing.\n286\n• L381-385: “In summary, our results show that uncertainties related to the emergence of a detectable\n287\nand attributable human influence on regional precipitation over the CONUS over the historical record\n288\n[Figure SPM.3; 31] can in most cases be explained by decomposing the combined ANT signal into its\n289\ntwo components, GHG and AER forcing.” Results do not support this solid conclusion for regional\n290\nattribution. Also, Figure SPM3 considers heavy precipitation trends since 1950, but the authors’\n291\nanalysis covers different periods, particularly for AER.\n293\nthe IPCC result was more used as evidence that large uncertainties remain in attributing changes to\n294\nregional precipitation.\n295\n8\n• Line 617-620: “we quantify uncertainty via resampling methods, specifically a block bootstrap for\n296\nestimating standard errors and confidence intervals [as in 69] and a permutation/reshuffling approach\n297\nfor ascribing statistical significance to spatial patterns [as in 26].” This is the important procedure for\n298\nestimating uncertainty in scaling factors, so the authors need to explain details in comparison with\n299\n300\ntraditional d&a methods. How does this consider internal variability, particularly PD and PW?\nThank you for raising this point, which is well-taken. In this case, we obtain best estimates of the\n301\nreshuffling) allow us to build a sampling distribution (via bootstrap) and null distribution for testing\n303\n(via reshuffling) that explicitly account for internal variability by quantifying the relative magnitude\n304\n305\nof year-to-year variability beyond what is described by PF.\n• L630-631: “the proportion of type I errors at 0.33 (for moderate significance) and 0.1 (for strong\n306\nsignificance)”. Figure 2 shows that grid-box attribution results for AER have ‘moderate’ significance\n307\nonly in all cases (almost no grid shows ‘strong’ significance), indicating AER signal is not detectable\n308\nwhen applying 5% level as used in OF. Is there any specific reason to use the lower significance level?\n310\nthe IPCC and elsewhere. Ultimately, any specific significance level (including, e.g., the ubiquitous\n311\nα = 0.05) is a subjective choice, for which the evidence against falsifiable hypotheses should always\n312\nbe judged in light of. Here, we have found it helpful to present multiple significance statements to\n313\nhelp the reader understand the strength of evidence for the claims made throughout the paper.\n314\n• Line 676-679: “Therefore, in our case the D&A exercise is reduced to assessing a single null hypoth-\n315\n316\nesis test (i.e., rejecting H0,f : βf = 0 both detects and attributes forcing f), such that larger βf (in\nabsolute value) implies a stronger influence of forcing f.” I think this is equivalent to the ‘detection’\n317\nof signals in OF. Can authors think about another hypothesis test for ‘attribution’ by which they check\n318\nagreement of the forced trend with the observed trend. I think this sort of quantification of forcing\n319\ncontribution is needed for robust ‘attribution’.\n321\nour case, our estimates of the forced trend are the observed trends. We would like to point out that\n322\nwe have a brief mention in the discussion of how one could compare overall trends from observations\n323\nwith corresponding quantities from the GCMs – see our Supplemental Figure B.14."}, {"round": 2, "reviewer_comments": {"reviewer_2": "I suggest accepfing this manuscript published in Nature Communicafions.", "reviewer_3": "I am safisfied with the authors' responses to my previous comments and would like to recommend acceptance in the present form. In parficular, the separate considerafion of fast and slow precipitafion response to aerosols as well as the improved explanafions of their methods have clarified the major issues."}, "author_response": ""}]

Decision letter