A western United States snow reanalysis dataset over the Landsat era from water years 1985 to 2021

Figure 4 shows a sample of the seasonal cycle and spatial distribution of SWE over HUC2 basins and the entire WUS domain in WY 2019. No SWE or snow depth measurements are assimilated in deriving the WUS–SR dataset. Thus, in situ SWE and snow depth measurements, and ASO SWE and snow depth estimates are used as independent verification datasets. Landsat fSCA measurements are assimilated into the snow reanalysis framework assuming a measurement error (standard deviation) of 10%34. Though Landsat fSCA cannot be used for independent verification, the WUS–SR posterior fSCA estimates, which are fitted to these measurements using a likelihood function, are expected to have comparable bulk error. The snow reanalysis framework has been successfully applied previously to generate datasets over the Sierra Nevada, Andes, and High Mountain Asia33,50,52.

Fig. 4

Illustrative results from the WUS–SR SWE estimates in WY 2019. (a) Seasonal cycle of SWE volume (km3) integrated over HUC2 basins. (b) Spatial distribution of SWE (meters) over part of the Sierra Nevada on March 1st, WY 2019. (c) Spatial distribution of WUS SWE (meters) on March 1st, 2019. The boxed area in (c) represents that shown in (b).

Full size image

Verification with in situ data

In this section, grid-averaged reanalysis SWE and snow depth are compared with point-scale in situ measurements. It should be acknowledged a priori that there are inevitable representativeness issues in the comparison between point-scale in situ data and grid-averaged snow reanalysis data. The WUS–SR estimates are modeled with assumed sub-grid heterogeneity within each ~500 m grid cell (which is modeled via a lognormal distribution) meant to account for the complex sub-grid variations in terrain (elevation, slope, aspect), forest cover, and meteorological forcings. Given that in situ stations are often sited in non-representative regions of a grid cell (i.e., in sheltered flat forest clearings), it is unlikely that the grid-averaged SWE/snow depth (spanning ~ 250,000 m2) should match the point-scale in situ SWE/snow depth (spanning ~10 m2). Nevertheless, in situ measurements, from the SNOTEL and CA Department of Water Resources (CADWR) networks, represent the best available data that covers much of the WUS and extends back several decades. While not expected to match each other, the verification herein is meant to illustrate consistency between the in situ measurements and WUS–SR estimates.

Peak SWE comparison with in situ data

In situ SWE measurements from WY 1985 to 2021 are taken from 1) the SNOTEL network ( managed by the U.S. Natural Resources Conservation Service (NRCS), and 2) CADWR ( from sensor type: “SNO ADJ (82)”), collections of automated snow pillows in the WUS. For in situ verification, we pair each in situ site with the closest snow reanalysis grid based on the geolocation of these two datasets. The precision of in situ coordinate values varies from 0.000001° (<1 m) to 0.01° (>1 km). Considering the potential for geolocation mismatch, the nine nearest pixels32,33,55 are additionally used to compare in situ and WUS–SR peak SWE. In this latter approach, the differences between in situ peak SWE and the neighboring WUS–SR grid cell peak SWE with the smallest difference among the nine nearest snow reanalysis grids are used. To compare the SWE on the same day, peak SWE day determined by in situ SWE is used to extract peak SWE from both datasets throughout the paper.

Figure 5 presents the density scatter plots comparing in situ peak SWE values against collocated grid-cell posterior peak SWE values. Peak SWE values less than 1 cm are screened out from the comparison. In total, 928 in situ sites are used in the comparison with the WUS–SR SWE estimates. To understand the performance of the WUS–SR dataset across different regimes in the WUS, verification is conducted for each HUC2 basin. The comparison is quantified using correlation coefficient (R), mean difference (MD), and root mean square difference (RMSD). Table 5 summarizes the number of total site-years, and statistics for both prior and posterior reanalysis SWE against in situ SWE within each HUC2 basin and over the WUS.

Fig. 5

Density scatter plot of in situ (snow pillow) peak SWE and collocated posterior (grid-average) peak SWE grouped by HUC2 basins over WYs 1985 to 2021. The solid black line is the 1:1 line. The correlation coefficient (R), mean difference (MD), and root mean square difference (RMSD) are shown for each HUC2 basin. In situ data with peak SWE values greater than 1 cm are included in the comparison.

Full size image
Table 5 Number of in situ sites and comparison metrics between in situ (snow pillow) peak SWE and collocated grid-averaged snow reanalysis prior and posterior (post.) peak SWE grouped by HUC2 basins.
Full size table

Compared with the performance of the prior peak SWE estimates (i.e., not constrained by Landsat fSCA), posterior SWE estimates show a better correlation (higher R) with less bias and random error (lower MD and RMSD) than the prior SWE over most of the HUC2 basins. Posterior SWE in CA has the highest correlation against in situ SWE (R = 0.82). The correlations with in situ SWE over the entire WUS are improved from 0.74 (prior) to 0.77 (posterior). Posterior peak SWE in UCRB has lower bias and uncertainty compared against in situ data with a relatively small MD of 0.06 m in absolute value (reduced by 62% from prior MD) and RMSD of 0.19 m (reduced by 27%). Over the WUS, in situ peak SWE is (on average) larger than the WUS–SR peak SWE (negative MD). Sub-grid topographic variability, snow-forest interactions, and wind-driven snow redistribution may all cause differences seen between grid-averaged peak SWE and point-scale in situ peak SWE. The statistics for PN indicate comparable correlation of in situ and both prior and posterior snow reanalysis, however the MD and RMSD do not get improved from posterior to prior. Fewer cloud-free fSCA measurements are available in PN, which limits the improvement of snow reanalysis SWE via data assimilation.

To acknowledge the potential geolocation mismatch, Fig. 6 provides verification of in situ peak SWE and posterior reanalysis peak SWE using an approach comparing to the best match among the nine nearest pixels. The WUS-wide correlation coefficient (R), MD and RMSD of posterior peak SWE and in situ peak SWE is 0.91, −0.08 m, 0.18 m, respectively. Compared to the approach used in Fig. 5, the posterior reanalysis peak SWE in Fig. 6 (as expected) is more correlated with in situ peak SWE (R values above 0.9), and has lower MD (<0.13 m) and RMSD (<0.24 m) over the WUS and at all HUC2 basins. Posterior reanalysis peak SWE is still lower than the in situ peak SWE at most of the sites, with the largest MD found in the PN. The PN has fewer cloud-free fSCA measurements, which may lead to larger errors than in regions with fewer cloud-contaminated images. The MD in CA is −0.07 m, which is within the range of −0.12 to 0.01 m as reported in Margulis et al.33, where the original 90-m Sierra Nevada SWE reanalysis was compared against in situ peak SWE using the same approach.

Fig. 6

Same as the density scatter plot in Fig. 5 but using posterior (grid-average) peak SWE from the best match among nine closest neighbor pixels.

Full size image

Figure 7a shows that the differences between posterior peak SWE and in situ peak SWE are sensitive to forest fraction exceeding 40%. The median RMSD remains stable at ~ 0.18 m for forest fractions below 40%, and gradually increases to ~ 0.38 m when forest fraction increases to over 60%. The larger RMSD at higher forest fraction pixels might be caused by 1) larger disparities between in situ sites (that tend to be in forest clearings) and collocated pixels with large averaged forest coverage fraction and/or 2) larger estimation errors in WUS–SR peak SWE in areas with large forest coverage. Aside from forest coverage effects, the difference between in situ and posterior peak SWE is impacted by the number of fSCA measurements as illustrated in Fig. 7b. When over 40 fSCA measurements (after cloud screening) are available, the median of absolute difference is as low as ~ 0.11 m. As the number of annual fSCA measurements is reduced, the median and spread of the absolute difference of peak SWE for each year increased. Figure 7c show that the peak SWE days determined by in situ data is highly correlated to peak SWE days determined by posterior WUS–SR SWE (R = 0.73). Overall, in situ SWE peaks later than the WUS–SR SWE with a MD value of −10 days.

Fig. 7

(a) RMSD of peak SWE as a function of averaged forest fraction for each site. RMSD is determined at each site from the 37-year peak SWE from in situ and posterior WUS–SR. (b) Absolute difference of peak SWE over the number of fSCA measurements (after cloud screening) for each year and site. The absolute difference of peak SWE is computed using in situ and posterior peak SWE. (c) Density scatterplot of peak SWE day from in situ and posterior WUS–SR for each year and site.

Full size image

Temporal (daily) SWE comparison with in situ data

Figure 8 shows the spatial distribution of verification statistics at in situ sites by comparing posterior daily SWE against in situ daily SWE greater than 2.54 mm.

Fig. 8

Spatial distribution of evaluation statistics determined via comparison of in situ daily SWE and collocated posterior SWE over WYs 1985 to 2021. Statistics include (a) R, (c) MD in meters, (d) RMSD in meters, (e) MD as percentage of peak SWE, and (f) RMSD as percentage of peak SWE. For reference, the in situ site elevations in meters are shown in (b). Daily SWE values less than 2.54 mm are excluded.

Full size image

Over the entire WUS, posterior daily SWE at in situ sites have high correlations (median of 0.79), small MD (median of −0.08 m) and RMSD (median of 0.17 m) against in situ SWE. The comparison suggests that posterior daily SWE agrees reasonably well with daily in situ SWE, especially in CA and UCRB with higher correlations and relatively lower MD and RMSD. Daily posterior SWE is slightly lower than point-scale in situ SWE (Fig. 6b. negative MD in blue) at most of the sites. At some in situ sites in the western PN, posterior SWE shows higher differences. Figure 8(e,f) show that low MD and RMSD expressed as percent of peak SWE are observed at some sites with high MD and RMSD due to deep snow. For sites with both large absolute and percent of differences, some of these differences may represent larger errors caused by fewer available fSCA measurements after clouds screening. Finer resolutions may be needed to capture large sub-grid SWE values.

Peak snow depth comparison with in situ data

In situ snow depth measurements are taken from the same sources as SWE (i.e., NRCS and CADWR from sensor type: “SNOW DP (18)”). Similar verification steps as with peak SWE (Fig. 5) are conducted for snow depth as shown in Fig. 9. Compared to the SWE measurements, however, in situ snow depth measurements appear to be of lower quality with some station-years showing snow depth with persistently high values throughout the year, non-physical oscillations in the measurements, and other erroneous behavior that are clearly inconsistent with the corresponding SWE measurements. Hence, extra screening is applied to the data before being used for verification. In situ snow depth measurements that changed by more than 1 m in a single day were assumed erroneous and excluded from the analysis. Further, assuming snow density is within the range of 200 to 500 kg/m3 at the peak day, snow depth measurements outside 2–5 times the corresponding SWE measurements were removed. To avoid incorrectly diagnosing peak snow depth day from snow depth measurements with missing data after screening, the in situ peak SWE day was used to determine the in situ snow depth used for comparison with posterior reanalysis estimates. Overall, posterior peak snow depth is correlated with in situ peak snow depth (R = 0.72) and has an MD of −0.36 m and RMSD of 0.66 m over the WUS. Compared to the results from peak SWE verification, the correlation coefficient between in situ and posterior peak snow depth is about the same at all HUC2 basins, with the highest value (R = 0.81) in CA. The MD and RMSD values for peak snow depth are around 2 to 3 times larger than those in peak SWE, partially caused by larger snow depth values than SWE and perhaps the poorer quality of in situ snow depth measurements.

Fig. 9

Same as Fig. 5 but for peak snow depth. Peak day is determined by in situ peak SWE. In situ data with peak snow depth values greater than 5 cm are included in the comparison.

Full size image

Verification with airborne snow observatory (ASO) data

The WUS–SR estimates are further verified against gridded SWE and snow depth estimates from ASO11,56,57,58. The lidar-based ASO measures snow depth via an airborne laser scanner (ALS) based on the differences in elevations between a snow-off day and snow-on days. ASO SWE is estimated from the high-resolution snow depth measurements and modeled snow density11. For comparison, the 50-m ASO SWE and snow depth snapshots are aggregated to the WUS–SR SWE model resolution. ASO data is available over select sites in California, Colorado, and Washington starting from 2013. While abundant snapshots are available in the Tuolumne River Basin in California, limited snapshots (commonly once per year) were taken at most of the ASO sites. ASO snow depth is a relatively accurate measurement with measurement error less than 0.02 m at a 50 m × 50 m grid. Model error (5%–8%11) could exist in modeled snow density, which is expected to propagate to ASO SWE estimates.

ASO SWE and snow depth estimates are compared with prior and posterior ensemble median SWE and snow depth maps on coincident days (Figs. 10 and 11). Tables 6 and 7 reports the statistical metrics for comparisons closest to April 1st at sampled ASO basins: USCATB (Tuolumne River Basin, California), USWAOL (Olympic Mountains, Washington), and USCOCM (Aspen/Castle-Maroon, Colorado).

SWE map comparison

For the California domain (USCATB, Fig. 10 left column, Table 6), posterior SWE is highly correlated with ASO SWE (ranging from 0.81 to 0.91) compared against prior SWE (ranging from 0.50 to 0.71). A negative MD indicates that the WUS–SR SWE is less than ASO SWE (on average) in Tuolumne. The difference significantly decreases from prior to posterior estimates in most years, along with decreased RMSD. WY 2015 was a historically dry year, in which posterior SWE shows no bias compared with ASO SWE, with a small RMSD of 0.07 m. Posterior SWE in WY 2017 has the highest correlation (0.91) with ASO SWE compared with a lower correlation (0.56) in prior SWE. MD drops from −0.13 m to −0.04 m, and RMSD decreases by half from prior to posterior in WY 2017. Figure 10 (left column) illustrates that Tuolumne-averaged posterior SWE (1.23 m) is comparable with ASO SWE (1.27 m), suggesting that the posterior WUS–SR SWE and ASO are in good agreement with respect to the basin-wide mean SWE. The prior underestimates SWE at high elevations in the northern parts and southern edges of Tuolumne basins whereas it overestimates shallow SWE near the basin outlet. The performance of the spatially distributed posterior SWE is considerably improved over the prior compared with ASO SWE. Though MD in WY 2019 increases from −0.06 m to −0.14 m (from prior to posterior), RMSD in that year decreases from 0.34 m to 0.27 m. The differences between prior SWE and ASO SWE are large in absolute values, while large positive differences are offset by negative differences causing a low MD for prior SWE in WY 2019.

Fig. 10

Comparison of ASO SWE with prior and posterior SWE at three ASO sites (top four rows): Tuolumne River Basin, California, (USCATB) in WY 2017 (left column); Olympic Mountains, Washington, (USWAOL) in WY 2016 (middle column); Aspen/Castle-Maroon, Colorado (USCOCM) in WY 2019 (right column). The prior maps are not shown, but instead included implicitly via the difference maps. The bottom row shows the relative RMSD between ASO and WUS–SR SWE as a function of forest fraction. RMSD (from pixels with both ASO and WUS–SR SWE greater than 1 cm) is computed for each forest fraction bin and then normalized by bin-averaged ASO SWE to get relative RMSD.

Full size image
Table 6 SWE comparison statistics between ASO SWE estimates and prior and posterior (post.) snow reanalysis SWE on ASO measurement days (Day of Water Year; DOWY) closest to April 1st.
Full size table

Non-seasonal SWE in portions of the PN (USWAOL) site is a potential error source in both ASO SWE and WUS–SR SWE. Snow depth retrieved from ASO may be erroneous at glacier pixels due to the lack of snow-off flights. The snow reanalysis framework does not include explicit modeling of glaciers. Therefore, non-seasonal snow pixels are removed when comparing the ASO SWE with WUS–SR SWE. This paper generates the WUS–SR non-seasonal snow mask following the method described in Liu et al.50. To summarize the method herein, a pixel is considered as a non-seasonal snow pixel if the annual minimum SWE exceeds 10% of the annual maximum SWE at least once over the dataset period. After applying the non-seasonal snow mask, the mean posterior SWE is 0.51 m which is slightly lower than 0.55 m in ASO SWE. Though the correlation coefficient is high (over 0.8) between prior snow reanalysis SWE and ASO SWE, the MD and RMSD in absolute value is over 0.50 m and 0.60 m respectively, which are both reduced significantly (by 94% and 44% respectively) in the posterior.

In Colorado (USCOCM), the mean of posterior SWE (0.55 m) is comparable with ASO SWE (0.54 m). The MD is reduced by 98% (to 0.01 m) and RMSD is reduced by 64% (to 0.17 m) from prior to posterior estimates. Although the posterior correlation coefficient is significantly improved over the prior, it is lower than the values seen at the USCATB and USWAOL sites. In Colorado, snow albedo has been shown to be affected by dust, black carbon, and other light-absorbing particles in recent decades59. In the current snow reanalysis framework, the impact of dust on snow albedo is modeled through an unconstrained uncertainty parameter. Future work could be done to apply a more explicit treatment of dust impacts on snow albedo to yield potentially improved results.

The effect of forest fraction on the performance of reanalysis SWE estimates is further illustrated using ASO SWE in Fig. 10. The Olympics basin has denser forest fraction with a mean of 58%, while the Tuolumne and Aspen/Castle-Maroon basins have mean forest fractions of 17% and 20%, respectively. At all three ASO basins, the relative RMSD of posterior SWE increases with the forest fraction. This is expected since Landsat-derived fSCA is only available over bare areas and/or forest gaps within a pixel. As forest cover increases, less useful information is available, while information is maximized at 0% of forest cover. However, the improvement in prior to posterior SWE estimates increases with forest coverage. This is likely related to the increased complexity of modeling SWE in dense forest areas where the larger uncertainty in forest areas is still reduced with the assimilation of fSCA.

Snow depth map comparison

Similar to the SWE comparison, posterior snow depth is verified against the ASO snow depth measurements (Fig. 11, Table 7). The spatial distribution of snow depth differences is comparable to the SWE differences with a correlation coefficient (R) of 0.85 and 0.76 in Washington and Colorado, respectively, and a value above 0.82 in California. In California, the MD of posterior snow depth is reduced by over 30% and RMSD is reduced by over 20% compared to the statistics of prior snow depth over WY 2015 to 2018, and WY 2021. In WY 2019 and 2020, while the posterior MD values are larger than the prior MD (positive and negative differences cancel each other out), the R values are as high as 0.9, and the RMSD values are reduced by 28% and 30%, respectively. In Washington, the posterior MD is close to 0 with RMSD significantly reduced by over 50% from the prior to the posterior estimates. In Colorado, despite the absolute values of MD and RMSD for posterior snow depth being more than twice the values of posterior SWE statistics (due to the larger dynamic range), the estimation of posterior snow depth is significantly improved from the prior snow depth with MD and RMSD reduced by 60% and 40%, respectively.

Fig. 11

Same as Fig. 10 (top four rows) but for snow depth.

Full size image
Table 7 Same as Table 6 but for snow depth.
Full size table

Source: Resources -

A short exposure to a semi-natural habitat alleviates the honey bee hive microbial imbalance caused by agricultural stress

Phosphate limitation intensifies negative effects of ocean acidification on globally important nitrogen fixing cyanobacterium