The recent dryness in California was unprecedented in the instrumental record. This article employs spatially explicit precipitation reconstructions for California in combination with instrumental data to provide perspective on this event since 1571. The period 2012–15 stands out as particularly extreme in the southern Central Valley and south coast regions. which likely experienced unprecedented precipitation deficits over this time, apart from considerations of increasing temperatures and drought metrics that combine temperature and moisture information. Some areas lost more than two years’ average moisture delivery during these four years, and full recovery to long-term average moisture delivery could typically take up to several decades in the hardest-hit areas. These results highlight the value of the additional centuries of information provided by the paleo record, which indicates the shorter instrumental record may underestimate the statewide recovery time by over 30%. The extreme El Niño that occurred in 2015/16 ameliorated recovery in much of the northern half of the state, and since 1571 very-strong-to-extreme El Niños during the first year after a 2012–15-type event reduce statewide recovery times by approximately half. The southern part of California did not experience the high precipitation anticipated, and the multicentury analysis suggests the north-wet–south-dry pattern for such an El Niño was a low-likelihood anomaly. Recent wetness in California motivated evaluation of recovery times when the first two years are relatively wet, suggesting the state is benefiting from a one-in-five (or lower) likelihood situation: the likelihood of full recovery within two years is ~1% in the instrumental data and even lower in the reconstruction data.
The recent record-breaking dryness in California, considered here for water years 2012–2015, is unprecedented in the instrumental record going back to 1896 (NCEI 2016). [A water year (WY) is defined as running from October of yeart−1 through the following September of yeart, and all years mentioned hereinafter (other than reference citations) are WY periods.] The proximate dynamical causes and impact of this lack of moisture delivery have been much evaluated (e.g., AghaKouchak et al. 2014; Wang et al. 2014; Belmecheri et al. 2015; Diaz and Wahl 2015; Diffenbaugh et al. 2015; Seager et al. 2015; Prein et al. 2016). It has also been analyzed regarding the extent to which associated drought conditions—including temperature impacts that enhance evaporative demand and potential evapotranspiration (Cook et al. 2016)—are unprecedented relative to the past 1–2 millennia (Adams et al. 2015; Robeson 2015; Griffin and Anchukaitis 2014) and can be attributed to anthropogenic influence on regional temperatures and the co-occurrence of dry and warm years (Williams et al. 2015; Diffenbaugh et al. 2015; cf. Cook et al. 2016). In some parts of the state, more than two full years’ moisture delivery was lost during the past four years (Fig. 1).
Here, we evaluate this dry period from a long-term spatial paleoclimate perspective going back to 1571, coupled with instrumental information through 2015. We utilize high-resolution, spatially explicit (0.5° × 0.5°) paleoclimate reconstructions reported and evaluated earlier at the level of the California and western Nevada (CANV) areal average (Diaz and Wahl 2015) to examine the past four years relative to this extended period for the state’s seven climate divisions (Fig. S1 in the supplemental material). From a technical perspective, we also evaluate the extent to which the division-level reconstructions are characterized by nonwhiteness–memory, which was highlighted by Bunde et al. (2013; cf. Zhang et al. 2015) as a potentially significant problem with tree-ring-derived precipitation reconstructions, and which has not been evaluated in related reconstructions (Wise 2016). [Here, we use terms derived from “whiteness” to refer to these characteristics.] Both the reconstructions and instrumental data exhibit whiteness according to this evaluation (supplemental material), and we use this result to motivate evaluating the time for each grid cell to recover back to long-term average moisture delivery from the cumulative 2012–15 precipitation loss. We additionally evaluate this recovery time when the first year in the recovery period is a very-strong-to-extreme El Niño (EN) event, as was the case for 2016.
The truncated EOF–principal components spatial regression (TEOF-PCSR) methodology used to assess preinstrumental climate in western North America was described and evaluated by Wahl and Smerdon (2012; cf. supplemental material therein), and details of its use to reconstruct CANV precipitation are described in Diaz and Wahl (2015). In TEOF-PCSR, least squares regression is applied to the principal component time series (PCs) of a truncated empirical orthogonal function (TOEF) representation of the instrumental climate field and a collection of leading PCs of the predictor data. In this procedure, the n number of reconstructed instrumental PCs n generated by the regression are substituted back into the singular value decomposition of the instrumental field, Precipitationfield–fitted = , to derive the reconstructed field (where the subscript n denotes the rank-reduced matrices). In this formula, n is the matrix of reconstructed instrumental PCs, n is the diagonal matrix of singular values, and is the transposed matrix of EOFs. The reconstruction period is 1571–1977, with a calibration period of 1916–77 and an independent validation period of 1896–1915 (Table 1). The calibration period was set to be as long as possible, given the start year of the instrumental data (1896) and the latest year common to the predictor data (1977), while retaining at least two decades worth of data for independent validation. Along with the high validation performance for reconstruction of the CANV spatial mean reported in Diaz and Wahl (2015), we note the strong spatial validation performance of the reconstruction for California by the standard reduction of error (RE) and coefficient of efficiency (CE) metrics (Figs. S2–S4 and associated captions in the supplemental material).
For the purposes of this paper, the variance of each grid cell’s reconstruction was adjusted to correspond to that of the instrumental data over the calibration period. This adjustment was done to ensure that the amplitude of precipitation variability is properly captured by the reconstructions in light of our foci to: 1) compare the depth of the recent deficits with other dry extremes over the past ~450 yr (section 3a and Figs. 2–4); and 2) to evaluate how long it will take to fully recover back to long-term average moisture delivery (section 3b and Fig. 5). We note that this kind of adjustment, by construction, leads to a loss of optimality relative to the squared error-minimizing character of unadjusted linear regression; however, the loss is relatively small in this case (cf. Figs. S2 and S4) and concentrated in the eastern interior desert and far northwest coastal regions (Fig. S4). The autocorrelation function and persistence length analyses (supplemental material) are not affected by the adjustment.
An important feature of the reconstructions we report and utilize is the large (n = 1000) probabilistic reconstruction ensemble associated with the expected value (EV) reconstructions, which has been newly developed for this paper and is the first of its kind for western North America spatial precipitation reconstructions. The ensemble was generated at the gridcell scale, and the probability ranges and values reported in sections 3a (Figs. 2–4) and 3b (Fig. 5) are developed from it, centered on the corresponding EV reconstructions. This ensemble method is motivated by the U.S. National Research Council (NRC) report of 2006, recommending that the high-resolution paleoclimate reconstruction community provide more systematic estimation of reconstruction uncertainty (NRC 2006). Our implementation utilizes a Monte Carlo generalization of the “bootstrapping the residuals” technique (Dixon 2006), in which the residuals are statistically modeled using the Hosking simulation algorithm to emulate their autoregressive characteristics (Hosking 1984; cf. https://www.rdocumentation.org/packages/waveslim/versions/1.7.5/topics/hosking.sim concerning implementation of Hosking’s method in the R language), and then independent random draws are derived from the modeled distribution to drive generation of the ensemble members (cf. Li et al. 2007; D. Nychka 2010, personal communication). This process produces an ensemble of equally likely precipitation reconstructions conditional on the predictor data, rather than confidence intervals estimated from one EV reconstruction realization (Li et al. 2007). Details of the ensemble generation methodology in the spatially explicit context are provided in Wahl and Smerdon (2012; cf. Wahl et al. 2014).
The instrumental data were regridded to 0.5° scale by area-weighted averages from the 5-km gridded precipitation dataset for the continental United States (CONUS) developed by NOAA’s National Centers for Environmental Information, Center for Weather and Climate (NCEI-CWC) (Vose et al. 2014). There is necessarily a loss of spatial heterogeneity below the 0.5°-scale relative to the 5-km-scale grid; however, this does not affect the interannual variability of the data at the regridded scale.
The whiteness of the precipitation data [supplemental material; cf. Margulis et al. (2016), who find the same whiteness characteristic for Sierra Nevada snowpack data] allows us to examine recovery times from a specified initial precipitation deficit without conditioning the examination on the characteristics of the immediately preceding years. That is, the analysis does not need to depend (in this case) on examining conditions after a run of four particularly dry years but rather can evaluate recovery time from a specified deficit for all possible starting years since the data lack temporal dependence. We define length of recovery to be the number of years it takes for the cumulative precipitation in a grid cell to equal or exceed the cumulative climatological average precipitation for that cell over the same number of years, when the cumulative total for the cell is initialized at the start of the recovery period by the 2012–15 loss (Fig. 1). Cells with no cumulative loss over 2012–15 were assigned recovery times of zero by definition; there are three such cells along the Colorado River in the southeastern portion of the state (Fig. 1). All other cells started with a minimum recovery time of one year. We evaluated recovery periods up to 50 yr in length so that the longest recovery time noted in Fig. 5 is 50 or more years. We thus needed to stop the reconstruction-based analysis at 1928, going forward from the first reconstructed year, 1571. Similarly, we needed to stop the instrumental-based analysis at 1966, going forward from the first WY in the period of record, 1896. Because of the whiteness exhibited by the data, the recovery length analysis can also logically be run in reverse, starting at nearer years and proceeding backward in time. We exploited this capability to fill in the starting years of 1977–29 for the reconstruction-based data and to fill in the starting years of 2015–1967 for the instrumental-based data (Table 1).
For the recovery length evaluation, the means of the gridcell reconstructions over their full period (1571–1977) were set equal to their corresponding instrumental period (1896–2015) means. Analyses we performed showed that the recovery times are related to the baseline period used to determine the climatological precipitation against which recovery is gauged. If the baseline period’s mean is above the full period mean and the data evaluated are not adjusted to the baseline period mean, the recovery times take longer than when the baseline period mean equals the full period mean; similarly, if the baseline period mean is below the full period mean and the evaluated data are not adjusted to the baseline period mean, the recovery times are shorter than when the baseline period mean equals the full period mean. Thus, an analytical decision needs to be made in terms of how recovery length is defined. The entire instrumental period represents the most directly measured information available against which to gauge recovery from the 2012–15 losses, and so we chose this time period as a reasonable baseline for the recovery period evaluation. In turn, we adjusted the means of the EV and ensemble member reconstructions to this norm to ensure full compatibility with the recovery times estimated from the instrumental data per se.
The specification of EN events for recovery length analysis in which the first year is a very-strong-to-extreme EN is taken from the extended multivariate EN–Southern Oscillation (ENSO) index (MEI; Wolter and Timlin 2011) for the period after 1876 (https://www.esrl.noaa.gov/psd/enso/mei.ext/) and a combination of preinstrumental information from Gergis and Fowler (2009, hereafter GF09) and Garcia-Herrera et al. (2008, hereafter G08) prior to 1877. After 1876, a very-strong-to-extreme event is identified by a MEI value ≥ 2; prior to 1877 such events are identified as very strong or extreme by GF09. GF09 is a meta-analysis that employs an objective evaluation methodology to allow combination of quantitative and documentary paleoclimatological information regarding the coupled ocean–atmosphere ENSO system. It includes both event timing along with estimation of relative event strength. G08 is an independent analysis that utilizes historical information from the Spanish archives to identify EN years along the coast of Peru. We use this additional information to help eliminate uncertainty in the GF09 timing chronology relative to the manifestation of ENs in the eastern Pacific that arises from their inclusion of event information across the Pacific basin: information that can include the manifestation of a strong EN in other parts of this region that is not well-represented in the eastern Pacific. The timing of events can be similarly affected by the location in the equatorial Pacific where pre-instrumental events are noted and whether their manifestation is primarily atmospheric or oceanic, all of which is taken into account by GF09. We thus use the G08 chronology of ENs to screen the very strong and extreme events listed by GF09 and employ only those events that agree on their timing in both listings. Together with the MEI events after 1876, this methodology identifies 12 events for analysis (Table 1).
a. Divisional reconstructions
The reconstructions for the seven California climate divisions (cf. Fig. S1) are shown in Figs. 2–4. The 90% reconstruction probability ranges for the annual values (light gray shading in Figs. 2–4) generally encompass the instrumental time series, indicative of the overall spatial skill mentioned in section 2 (cf. Figs. S2–S4). The North Coast drainage (division 1) exhibits the largest number of departures in this regard, as anticipated from Fig. S4, whereas the central coast, south coast, and San Joaquin drainages (divisions 4–6) exhibit nearly zero departures of this kind, also as anticipated (we note the unanticipated few departures in the southeast desert basin, division 7). In general, the divisional results are consistent with the statewide evaluation reported in Diaz and Wahl (2015) and the evaluation for the combined central coast, south coast, and San Joaquin drainages along with the southeast desert basin (divisions 4–7) reported by Griffin and Anchukaitis (2014): that the recent precipitation deficits, while clearly extreme, are not unprecedented in the annually resolved paleo record. It is the combination of precipitation deficits with temperature increases gauged by metrics such as the Palmer drought severity index (PDSI) that suggests uniqueness over the past millennium and longer for the 2012–15 event (Robeson 2015; Williams et al. 2015; Griffin and Anchukaitis 2014). Adams et al. (2015) reach the same conclusion for the past 2000 years from streamflow reconstructions in the southwestern Sierra Nevada driven by PDSI reconstructions (Cook et al. 2007). All the divisions except the San Joaquin and south coast drainages (divisions 5 and 6) show similar characteristics here (i.e., the red lines are above or near the top of the reconstruction 90% probability ranges outlined by the corresponding lighter purple lines in Figs. 2–4).
By contrast, the division-specific results we report allow identification that the precipitation deficits for the San Joaquin and south coast drainages (divisions 5 and 6) are likely unprecedented since the later sixteenth century. These are the regions where the 2012–15 losses were the largest (Fig. 1); in the worst-hit areas representing over two years of average precipitation, and in them no sequential 4-yr period is reconstructed in the EV to have experienced the same loss as 2012–15 (i.e., the heavier purple lines for these divisions are above the corresponding red lines in Fig. 3). The 90% ranges for such extreme losses in these two divisions indicate that there could have been 4-yr sequences prior to the instrumental period that were drier than 2012–15 (i.e., below the corresponding red lines in Fig. 3), but there is a larger probability that 2012–15 was even more extreme compared to the last ~450 yr than indicated by the EV results (i.e., the majority of the driest 4-yr reconstruction sequences are above the corresponding red lines in Fig. 3). From a social geography context, we note that the San Joaquin and south coast drainages encompass the agriculturally rich southern Central Valley and the southern coastal region that includes the greater Los Angeles and San Diego metropolitan areas, which have the third highest regional population (Mexico Consejo Nacional de Población 2016; U.S. Census Bureau 2015) and the second largest gross regional product (U.S. Department of Commerce, Bureau of Economic Analysis, 2016) in North America. In terms of the extreme dry departure of its decadal-smoothed values (heavy black and salmon lines in Fig. 3), the San Joaquin drainage was the hardest hit in terms of precipitation delivery.
Squared correlations between the divisions’ EV reconstructions are presented in Table 2, based on the value of the reconstructions in standardized anomaly units. These highlight the known transition of relatively wetter to drier conditions that occurs from north to south in California (cf. Dettinger et al. 1998), which is particularly notable for the southernmost divisions 6 and 7 versus the northernmost divisions 1–3. The length of the reconstruction period, 407 yr, allows large sample sizes for comparisons of progressively more extreme conditions, and Table 2 indicates that the greater the extreme (either dry or wet), the more the divisions experience similar climatic conditions. However, even in the most extreme situation evaluated (absolute anomaly value ≥ 1.0 for at least one division), the southern desert division 7 shares less than half of its variation with divisions 1 and 3 and just above half with division 2. This difference is particularly notable during the 2012–15 dry period for the Colorado Desert (farthest southeast portion of the state); the only three grid cells in California without a cumulative loss during the four-year period are located in this region (Fig. 1). Table S2 in the supplemental material reports comparable information derived using only the instrumental data. These results exhibit the same general spatial characteristics as the reconstructions (we note that the squared correlation values are generally lower for the shorter time period evaluated).
b. Evaluation of recovery times from recent precipitation deficits
As noted in section 2, the whiteness we find in the reconstructions and instrumental data allows us to evaluate recovery times from a specified initial precipitation deficit without conditioning the evaluation on the characteristics of the immediately preceding years. The spatial results of this evaluation are shown in Figs. 5b–f, along with the characteristics of the 2012–15 deficits at the gridcell scale in Fig. 5a, which replicates Fig. 1 for ease of comparison. The most notable feature is the length of the recovery times, discussed further in section 4. For the instrumental data, the statewide median recovery time across the grid cells is nearly three decades (29.5 yr; Fig. 5c), and nearly 10 yr longer for the reconstruction EV (39.0 yr; Fig. 5e). At the 20th-percentile level (i.e., the chance of a given recovery time or less occurring is 0.2), the corresponding statewide time periods are still at the decadal-plus scale: 8.8 yr for the instrumental data (Fig. 5d) and 14.0 yr for the reconstruction EV (Fig. 5f). We note that the median and 20th-percentile statewide recovery times determined from the instrumental data, while clearly lower than those determined from the reconstruction EV and thus in the lower tail of the corresponding reconstruction ensemble distributions, are still within the 90% probability ranges estimated from the ensemble (Figs. 5e,f). For comparison, the statewide median recovery time for 2014 alone, the single driest year during 2012–15, is 8.0 yr in the instrumental data (Fig. 5b), highlighting the longer-term moisture delivery impact of just one severely dry year. The statewide median recovery time when the first year of the recovery period is screened by very-strong-to-extreme EN conditions is 18.0 yr, indicating a substantial reduction of 54% relative to the full reconstruction statewide median of 39.0 yr (Fig. 5e), or 39% relative to the full instrumental statewide median of 29.5 yr (Fig. 5c). When the analysis is restricted to the reconstruction period only (ending in 1977), the 10 included events show similar results with a statewide median of 21 yr, for a 46% and 29% reduction relative to the full-reconstruction and full-instrumental statewide medians, respectively.
In light of the high precipitation amounts that have occurred to date in California during the current WY (1 October 2016–8 March 2017), we also estimated recovery times and their associated likelihoods under the simple assumption that the first two years of recovery are characterized by generally wet years. We used 1969 and 2005 as examples of such years in recent decades, which exhibited different relative strengths in terms of north–south representation (Figs. 2–4). Each year was used separately to represent the precipitation for the first two recovery years. The results of this analysis indicate much shorter recovery times, with statewide median values of 7.0 and 8.0 yr for the 1969 and 2005 cases, respectively, similar to the 20th-percentile-level case reported in Fig. 5d.
The result that the 2012–15 cumulative precipitation deficits in the San Joaquin and south coast drainages are likely unprecedented since at least the later sixteenth century, separate from considerations of increasing temperatures and related drought metrics, is climatologically important. This result suggests that California precipitation delivery by itself is sufficient under current climate conditions to produce the multiple-year extreme dry conditions recently experienced, at least in these regions, without accounting for the anthropogenic temperature forcing that is partly responsible for amplification of the extreme drought conditions (Williams et al. 2015). This previously unrevealed longer-term perspective underscores the value of paleo-precipitation reconstructions at relatively high spatial resolution. Considering the additional drying impacts of increasing temperatures as anthropogenic warming continues (Cook et al. 2016), the enhanced potential for combinations of severe precipitation deficits with high temperature-driven evaporation and potential evapotranspiration indicates the possibility of “megadroughts” in the North American Southwest, including California, that will exceed the worst such episodes evidenced in the region during medieval times (Cook et al. 2007; Stine 1994) by the second half of the twenty-first century (Cook et al. 2016). The Sierra Nevada snowpack reconstruction of Belmecheri et al. (2015) and related instrumental-period reanalysis of Margulis et al. (2016) suggest that such a signal may already be present in the unprecedented low snowpack that occurred in the Sierra Nevada in early 2015–the intensity of the snowpack loss was greatest at lower montane elevations (<2130 m), where the impact of higher winter temperatures on moisture retention would first be expected to be detectable.
The recovery time results are surprising in two ways. First, and most obviously, the estimated recovery times from the cumulative 2012–15 precipitation deficits are long, both for the instrumental and reconstruction data. To a first order, the longest recovery periods follow the intensity of the cumulative deficits, especially in the midportion of the state (Fig. 5). However, this pattern does not hold as generally in the southeastern deserts and adjacent southern coastal region, particularly in the instrumental data (Figs. 5c,d). In these portions of the state, the recovery times estimated from the instrumental data are relatively shorter in relation to the original 4-yr deficits, although with gridcell medians still on the order of two-plus decades in the heavily populated Los Angeles area (Fig. 5c). We note that surface agricultural amelioration, for example, could occur with even one year of high precipitation and more generally is strongly buffered via storage, transport, and irrigation. The long recovery periods we find are best considered as indicative of full hydrological recovery from the cumulative deficits over 2012–15, including groundwater recharge as it would be represented in the absence of California’s extensive storage, distribution, and usage systems. As mentioned above, the additional impacts of increasing temperatures on moisture balances into the later twenty-first century (Cook et al. 2016) suggest the possibility that the long recovery times our results indicate may actually be conservative relative to what future dry conditions will entail. While the pattern of recovery for 2014 alone (Fig. 5b) continues generally to follow the intensity of the cumulative 4-yr deficits, there is less overall spatial variability associated with this single-year loss.
Second, the differences between the estimated recovery times from the instrumental and reconstruction data (cf. Figs. 5c,d with Figs. 5e,f, respectively) are potentially surprising in that they are not necessarily anticipated from the autocorrelation and persistence length analyses, which show little difference relative to how much of the instrumental data are included (Figs. S5–S10 in the supplemental material). We note that over the short but strictly comparable periods available between the instrumental and reconstruction data (starting years 1896–1928 going forward and 1977–45 in reverse) the recovery times are, in fact, lower for the reconstructions compared to the instrumental data. These results highlight that recovery time estimation is sensitive to the length of data used in the evaluation and argue that the additional centuries of data included in the reconstructions realistically reveal longer recovery times than can be estimated from the 120-plus years of instrumental data. As noted in section 3b, the 90% probability ranges for the reconstruction-based recovery times include those determined from the instrumental data (Figs. 5c,d compared with Figs. 5e,f), but the overall probability distributions for the reconstruction-based recovery times are clearly offset to longer periods compared to the instrumental data.
We note the difference in median recovery times estimated here and the much shorter median recovery time (on the order of 5 yr) estimated for snowpack recovery in the Sierra Nevada region by Margulis et al. (2016). This distinction is primarily due to different ways in which the baseline threshold for deficit is defined. Here, the threshold is defined in relation to the long-term instrumental precipitation mean for each grid cell, as described in section 2. In Margulis et al. (2016) the threshold is represented by a drought level defined to be significantly below the central value of the snow water equivalent time series. As specified by Margulis et al. (2016), nearly 80% of years evaluated (51 out of 65 going back to 1951) are above this threshold and would thus contribute (on average) to recovery. In our analysis only 50% of years would do so (on average), leading to much longer estimated recovery times by construction. The distinction in estimated recovery times is also likely related to the additional information added by the longer period of evaluation we use, mentioned above.
As noted in section 3b, California is experiencing anomalously high precipitation amounts so far during the current WY. Conditional on the assumption described there that the first two years of recovery are characterized by generally wet years, the past half-millennium perspective we utilize suggests the possibility for recovery times of less than a decade at the statewide level, similar to those shown in Fig. 5d. Thus, California appears to be experiencing the benefit of a roughly one-in-five (or lower) likelihood situation for recovery from the 2012–15 precipitation deficits. If 2017 turns out to be an extremely wet WY, such as 1983, that would lead to even shorter, and thus even lower likelihood, estimated recovery times. Our results indicate that the likelihood of a statewide median full recovery within two years is ~1% in the instrumental data and even lower in the reconstruction EV.
Finally, the fact that very-strong-to-extreme EN conditions in the first year of recovery from moisture deficits such as those experienced in 2012–15 lead to a substantial reduction in estimated recovery times is consistent with the predicted winter precipitation probabilities for the extreme EN that occurred during 2016 (NOAA CPC 2015). However, this statewide result disagrees spatially with the distribution of precipitation that accompanied the winter 2015/16 EN conditions, which included above-average moisture delivery for a large portion of California north of ~36°N and generally below-average delivery south of this line (NOAA High Plains Regional Climate Center 2016). We note that only 1 out of the 12 very-strong-to-extreme EN events evaluated (1878) exhibits a clearly similar wetter-north–drier-south pattern, with a more muted similarity in terms of northern wetness in 1973. Thus, conditional on the analysis we report, the 2016 pattern represents a low-likelihood anomaly from an a priori perspective. Reducing uncertainty regarding the timing and event characteristics of pre-instrumental very-strong-to-extreme ENs would greatly help to expand the sample size for this kind of evaluation and thereby enhance its statistical strength, and this represents an important priority for further research.
The authors gratefully acknowledge the evaluations and suggestions of three anonymous reviewers, which substantially helped improve this work. ERW conceived and designed the study and wrote the article with input from HFD and the other authors. ERW performed the analyses with assistance from WSG. All authors contributed to the interpretation of the dataset and discussion. The full ensemble output for each climate division and by grid cell, along with maps of the year-by-year EV reconstructions, are available online (https://www.ncdc.noaa.gov/paleo/study/21793). This is a PAGES 2k (NAM 2k) associated product.
Supplemental information related to this paper is available at the Journals Online website: http://dx.doi.org/10.1175/JCLI-D-16-0423.s1.