The cd package builds on ERA5-Land hourly reanalysis
(1950–present, ~9 km native grid) fetched from the DestinE Earth Data
Hub, aggregated to monthly, seasonal, and annual periods over British
Columbia, and published as Cloud-Optimized GeoTIFFs alongside a static
STAC catalog in a public S3 bucket. R consumer functions read those
GeoTIFFs directly via GDAL, crop to any area of interest, and compute
baselines, anomalies, and Mann-Kendall / Theil-Sen trend statistics.
This vignette runs the consumer pipeline on a regional administrative boundary — the Fish & Wildlife Compensation Program (FWCP) Peace Region, ~73,000 km² in northeastern British Columbia covering Williston Reservoir and surrounding watersheds.
Area of Interest
The bundled area of interest is the FWCP Peace Region polygon (single
multi-polygon, EPSG:4326). Any sf polygon works. The
watersheds within this region all drain to the Mackenzie River — the
Peace flows east into the Slave, which flows into the Mackenzie, which
empties into the Beaufort Sea on the Arctic Ocean.
We also carry British Columbia ecoregions through the analysis. Ecoregions partition the province into broad climate-physiography zones based on latitude, elevation, and dominant vegetation. Climate departure can vary across them in ways the regional average hides, and we use them later as the sub-region for the per-ecoregion breakdown.
FWCP Peace Region area of interest (~73,000 km²) in northeastern British Columbia, coloured by ecoregion. Williston Reservoir dominates the basin.
Connect to the Data Catalog
The producer pipeline in this repository fetches ERA5-Land hourly reanalysis data, derives variables (vapour pressure deficit, relative humidity, soil moisture), aggregates to monthly / seasonal / annual, and writes Cloud-Optimized GeoTIFFs to S3. Alongside the GeoTIFFs we publish a static SpatioTemporal Asset Catalog (STAC) — a small set of JSON files that index the assets. Both the GeoTIFFs and the catalog JSON live at https://stac-era5-land.s3.us-west-2.amazonaws.com/catalog.json.
cd::cd_catalog() reads that URL by default. The
Cloud-Optimized GeoTIFFs are also usable directly outside R — for
example, in QGIS via the STAC plugin, in gdalcubes, or with
any STAC-aware client.
catalog <- cd::cd_catalog()
kableExtra::kable_styling(
knitr::kable(catalog, label = NA,
caption = "Available climate variables and periods in the STAC catalog."),
bootstrap_options = c("striped", "hover", "condensed")
) |>
kableExtra::scroll_box(height = "320px")Extract Climate Time Series
cd::cd_extract() crops each cloud-hosted GeoTIFF to the
area of interest and computes the spatial mean per year. For an area of
interest of this size a live extraction takes a few seconds per
variable. To keep the vignette fast and reproducible we load a
pre-computed result of exactly that call.
ts <- cd::cd_extract(catalog, aoi)
head(ts, 10)| variable | period | year | value |
|---|---|---|---|
| prcp | annual | 1950 | 645.5450 |
| prcp | annual | 1951 | 665.1381 |
| prcp | annual | 1952 | 673.8270 |
| prcp | annual | 1953 | 828.5204 |
| prcp | annual | 1954 | 942.6607 |
| prcp | annual | 1955 | 729.9778 |
| prcp | annual | 1956 | 752.0231 |
| prcp | annual | 1957 | 783.3298 |
| prcp | annual | 1958 | 698.6779 |
| prcp | annual | 1959 | 882.8473 |
Trends
Anomalies are computed against a pre-warming reference period — 1951–1980, the three decades before climate change accelerated. This is the same base period Hansen et al. (2012) use to detect the emergence of 3-sigma summertime-temperature outliers globally. Saying a year is “+1.5 °C” means it was 1.5 °C warmer than the average year between 1951 and 1980.
The trend table that follows has two rows per variable. We compute trends from two different start years:
- 1951–present (75 years) — the long view. Captures the full magnitude of warming since the pre-warming reference.
- 1981–present (45 years) — starts at the beginning of the World Meteorological Organization’s most recent 30-year “climate normal” (1981–2010). This is the reference period used in most published climate products, so it makes results easy to compare against Intergovernmental Panel on Climate Change reports and government climate summaries.
Comparing the two slopes is informative. If the 45-year slope is steeper than the 75-year slope, warming has accelerated — recent decades are heating faster than the long-term average. If the 45-year slope is shallower, warming has slowed (though “slower” almost never means “stopped”). When the two slopes are similar, the rate of change has been roughly steady across the full record.
Total Change is the slope multiplied by the number of years — the cumulative shift over the trend window.
bl_early <- cd::cd_baseline(ts, baseline_years = 1951:1980)
ano <- cd::cd_anomaly(ts, bl_early)
trn <- cd::cd_trend(ano, trend_start = c(1951, 1981))
cd::cd_summary(trn)| Parameter | Period | Slope | Years | Total Change | Unit | p-value |
|---|---|---|---|---|---|---|
| Precipitation | Annual | 0.085 | 75 | 6.4 | % | 0.1971 |
| Relative humidity | Annual | 0.001 | 75 | 0.1 | % | 0.8620 |
| Snow cover | Annual | -0.053 | 75 | -4.0 | % | 0.0026 |
| Snowfall | Annual | -0.088 | 75 | -6.6 | % | 0.2528 |
| Snowfall fraction | Annual | -0.079 | 75 | -5.9 | % | 0.0078 |
| Snowmelt | Annual | -0.039 | 75 | -2.9 | % | 0.5279 |
| Day of 50% melt | Annual | -0.150 | 75 | -11.3 | day | 0.0018 |
| Peak weekly melt rate | Annual | 0.101 | 75 | 7.6 | mm/wk | 0.3848 |
| Soil moisture | Annual | -0.002 | 75 | -0.1 | % | 0.8981 |
| Snow water equivalent | Annual | -0.098 | 75 | -7.3 | % | 0.2644 |
| Annual peak snow water equivalent | Annual | -0.075 | 75 | -5.7 | mm | 0.7697 |
| Maximum temperature | Annual | 0.027 | 75 | 2.0 | °C | 0.0000 |
| Mean temperature | Annual | 0.030 | 75 | 2.2 | °C | 0.0000 |
| Minimum temperature | Annual | 0.032 | 75 | 2.4 | °C | 0.0000 |
| Vapour pressure deficit | Annual | 0.005 | 75 | 0.3 | Pa | 0.0008 |
| Precipitation | Fall | 0.065 | 75 | 4.9 | % | 0.4926 |
| Relative humidity | Fall | -0.016 | 75 | -1.2 | % | 0.0906 |
| Snow cover | Fall | -0.032 | 75 | -2.4 | % | 0.4051 |
| Snowfall | Fall | 0.019 | 75 | 1.4 | % | 0.9344 |
| Snowmelt | Fall | -0.149 | 75 | -11.2 | % | 0.4208 |
| Soil moisture | Fall | 0.000 | 75 | 0.0 | % | 0.9964 |
| Snow water equivalent | Fall | 0.006 | 75 | 0.4 | % | 0.9854 |
| Maximum temperature | Fall | 0.013 | 75 | 1.0 | °C | 0.0659 |
| Mean temperature | Fall | 0.020 | 75 | 1.5 | °C | 0.0099 |
| Minimum temperature | Fall | 0.025 | 75 | 1.9 | °C | 0.0027 |
| Vapour pressure deficit | Fall | 0.002 | 75 | 0.2 | Pa | 0.0323 |
| Precipitation | Spring | 0.184 | 75 | 13.8 | % | 0.1644 |
| Relative humidity | Spring | -0.009 | 75 | -0.7 | % | 0.3848 |
| Snow cover | Spring | -0.047 | 75 | -3.5 | % | 0.0014 |
| Snowfall | Spring | -0.102 | 75 | -7.6 | % | 0.4587 |
| Snowmelt | Spring | 0.559 | 75 | 41.9 | % | 0.0006 |
| Soil moisture | Spring | 0.032 | 75 | 2.4 | % | 0.0151 |
| Snow water equivalent | Spring | -0.119 | 75 | -8.9 | % | 0.2723 |
| Maximum temperature | Spring | 0.025 | 75 | 1.9 | °C | 0.0007 |
| Mean temperature | Spring | 0.028 | 75 | 2.1 | °C | 0.0001 |
| Minimum temperature | Spring | 0.027 | 75 | 2.0 | °C | 0.0001 |
| Vapour pressure deficit | Spring | 0.005 | 75 | 0.4 | Pa | 0.0001 |
| Precipitation | Summer | 0.175 | 75 | 13.1 | % | 0.1847 |
| Relative humidity | Summer | -0.007 | 75 | -0.6 | % | 0.7558 |
| Snow cover | Summer | -0.121 | 75 | -9.1 | % | 0.0008 |
| Snowfall | Summer | -0.735 | 75 | -55.1 | % | 0.0015 |
| Snowmelt | Summer | -0.807 | 75 | -60.5 | % | 0.0066 |
| Soil moisture | Summer | -0.036 | 75 | -2.7 | % | 0.1728 |
| Snow water equivalent | Summer | -0.792 | 75 | -59.4 | % | 0.0029 |
| Maximum temperature | Summer | 0.026 | 75 | 2.0 | °C | 0.0004 |
| Mean temperature | Summer | 0.031 | 75 | 2.4 | °C | 0.0000 |
| Minimum temperature | Summer | 0.035 | 75 | 2.6 | °C | 0.0000 |
| Vapour pressure deficit | Summer | 0.009 | 75 | 0.7 | Pa | 0.0338 |
| Precipitation | Winter | -0.050 | 75 | -3.8 | % | 0.6606 |
| Relative humidity | Winter | 0.035 | 75 | 2.7 | % | 0.0026 |
| Snow cover | Winter | 0.000 | 75 | 0.0 | % | 0.2816 |
| Snowfall | Winter | -0.044 | 75 | -3.3 | % | 0.7076 |
| Snowmelt | Winter | 0.056 | 75 | 4.2 | % | 0.5987 |
| Soil moisture | Winter | -0.005 | 75 | -0.4 | % | 0.4983 |
| Snow water equivalent | Winter | 0.009 | 75 | 0.7 | % | 0.9126 |
| Maximum temperature | Winter | 0.045 | 75 | 3.4 | °C | 0.0000 |
| Mean temperature | Winter | 0.044 | 75 | 3.3 | °C | 0.0002 |
| Minimum temperature | Winter | 0.043 | 75 | 3.3 | °C | 0.0003 |
| Vapour pressure deficit | Winter | 0.001 | 75 | 0.1 | Pa | 0.0007 |
| Precipitation | Annual | -0.010 | 45 | -0.4 | % | 0.9298 |
| Relative humidity | Annual | -0.016 | 45 | -0.7 | % | 0.3947 |
| Snow cover | Annual | -0.063 | 45 | -2.8 | % | 0.0516 |
| Snowfall | Annual | -0.073 | 45 | -3.3 | % | 0.6317 |
| Snowfall fraction | Annual | -0.003 | 45 | -0.1 | % | 0.9143 |
| Snowmelt | Annual | -0.107 | 45 | -4.8 | % | 0.3947 |
| Day of 50% melt | Annual | -0.119 | 45 | -5.4 | day | 0.2141 |
| Peak weekly melt rate | Annual | 0.090 | 45 | 4.1 | mm/wk | 0.6598 |
| Soil moisture | Annual | -0.017 | 45 | -0.7 | % | 0.6041 |
| Snow water equivalent | Annual | -0.116 | 45 | -5.2 | % | 0.6598 |
| Annual peak snow water equivalent | Annual | -0.251 | 45 | -11.3 | mm | 0.6884 |
| Maximum temperature | Annual | 0.021 | 45 | 0.9 | °C | 0.0449 |
| Mean temperature | Annual | 0.023 | 45 | 1.0 | °C | 0.0264 |
| Minimum temperature | Annual | 0.024 | 45 | 1.1 | °C | 0.0116 |
| Vapour pressure deficit | Annual | 0.005 | 45 | 0.2 | Pa | 0.0963 |
| Precipitation | Fall | 0.083 | 45 | 3.7 | % | 0.6884 |
| Relative humidity | Fall | -0.022 | 45 | -1.0 | % | 0.1678 |
| Snow cover | Fall | -0.078 | 45 | -3.5 | % | 0.3630 |
| Snowfall | Fall | -0.046 | 45 | -2.1 | % | 0.7767 |
| Snowmelt | Fall | 0.177 | 45 | 8.0 | % | 0.6740 |
| Soil moisture | Fall | -0.003 | 45 | -0.2 | % | 0.9376 |
| Snow water equivalent | Fall | -0.443 | 45 | -19.9 | % | 0.3328 |
| Maximum temperature | Fall | 0.029 | 45 | 1.3 | °C | 0.0338 |
| Mean temperature | Fall | 0.038 | 45 | 1.7 | °C | 0.0053 |
| Minimum temperature | Fall | 0.046 | 45 | 2.1 | °C | 0.0014 |
| Vapour pressure deficit | Fall | 0.005 | 45 | 0.2 | Pa | 0.0834 |
| Precipitation | Spring | -0.231 | 45 | -10.4 | % | 0.3044 |
| Relative humidity | Spring | -0.054 | 45 | -2.4 | % | 0.0227 |
| Snow cover | Spring | -0.047 | 45 | -2.1 | % | 0.2444 |
| Snowfall | Spring | -0.267 | 45 | -12.0 | % | 0.4631 |
| Snowmelt | Spring | 0.550 | 45 | 24.7 | % | 0.0944 |
| Soil moisture | Spring | -0.010 | 45 | -0.5 | % | 0.7100 |
| Snow water equivalent | Spring | -0.127 | 45 | -5.7 | % | 0.6041 |
| Maximum temperature | Spring | 0.011 | 45 | 0.5 | °C | 0.6179 |
| Mean temperature | Spring | 0.006 | 45 | 0.3 | °C | 0.7468 |
| Minimum temperature | Spring | 0.003 | 45 | 0.1 | °C | 0.7917 |
| Vapour pressure deficit | Spring | 0.008 | 45 | 0.4 | Pa | 0.0141 |
| Precipitation | Summer | -0.121 | 45 | -5.4 | % | 0.7468 |
| Relative humidity | Summer | 0.001 | 45 | 0.0 | % | 0.9922 |
| Snow cover | Summer | -0.118 | 45 | -5.3 | % | 0.0766 |
| Snowfall | Summer | -0.544 | 45 | -24.5 | % | 0.1295 |
| Snowmelt | Summer | -0.822 | 45 | -37.0 | % | 0.0944 |
| Soil moisture | Summer | -0.069 | 45 | -3.1 | % | 0.2952 |
| Snow water equivalent | Summer | -0.659 | 45 | -29.7 | % | 0.0799 |
| Maximum temperature | Summer | 0.030 | 45 | 1.3 | °C | 0.0617 |
| Mean temperature | Summer | 0.034 | 45 | 1.5 | °C | 0.0113 |
| Minimum temperature | Summer | 0.036 | 45 | 1.6 | °C | 0.0015 |
| Vapour pressure deficit | Summer | 0.009 | 45 | 0.4 | Pa | 0.4057 |
| Precipitation | Winter | 0.187 | 45 | 8.4 | % | 0.4281 |
| Relative humidity | Winter | 0.021 | 45 | 1.0 | % | 0.3527 |
| Snow cover | Winter | 0.000 | 45 | 0.0 | % | 0.3945 |
| Snowfall | Winter | 0.145 | 45 | 6.5 | % | 0.4631 |
| Snowmelt | Winter | 0.419 | 45 | 18.9 | % | 0.3376 |
| Soil moisture | Winter | -0.003 | 45 | -0.1 | % | 0.8911 |
| Snow water equivalent | Winter | 0.068 | 45 | 3.1 | % | 0.8220 |
| Maximum temperature | Winter | 0.015 | 45 | 0.7 | °C | 0.4752 |
| Mean temperature | Winter | 0.007 | 45 | 0.3 | °C | 0.7321 |
| Minimum temperature | Winter | 0.005 | 45 | 0.2 | °C | 0.7917 |
| Vapour pressure deficit | Winter | 0.000 | 45 | 0.0 | Pa | 0.8220 |
cd::cd_plot_timeseries(
ano, variable = "tmean", period = "annual", trend = trn,
title = "Annual Mean Temperature Anomaly — FWCP Peace Region"
)
Annual mean temperature anomaly for the FWCP Peace Region relative to 1951-1980 baseline.
cd::cd_plot_timeseries(
ano, variable = "prcp", period = "annual", trend = trn,
title = "Annual Precipitation Anomaly — FWCP Peace Region"
)
Annual precipitation anomaly (% of 1951-1980 baseline) for the FWCP Peace Region.
Recent Decade vs Pre-Warming Reference
The table below compares two windows directly — the recent decade (2015–2025) against the pre-warming reference (1951–1980). Δ absolute is the difference of the two means in the variable’s native units. Δ % is shown only for variables where it is meaningful (precipitation, soil moisture, vapour pressure deficit, relative humidity). Two p-values are reported, answering different questions. Δ p (windows) is the Welch t-test on the annual values within each window — does the recent decade differ from the pre-warming reference? Trend p (75-yr) is the Mann-Kendall test on the full 1951–present series — is there a steady year-on-year ramp? Step changes show up as significant Δ p with non-significant trend p; gradual ramps show up as both significant.
The recent decade was 1.7 to 2.4 °C warmer than the pre-warming reference for annual mean, daytime maximum, and overnight minimum, with both window and trend p-values below 0.001. Vapour pressure deficit is up significantly on both tests. Annual precipitation was about 3 to 4 % higher in the recent decade; neither test confirms the shift (Δ p ≈ 0.4, trend p ≈ 0.20). Soil moisture and relative humidity show no meaningful change.
cmp <- cd::cd_compare(ts) # defaults: 2015–2025 vs 1951–1980, Welch t-test| Variable | Period | Recent (2015–2025) | Pre-warming (1951–1980) | Δ absolute | Δ % | Δ p (windows) | Trend p (75-yr) |
|---|---|---|---|---|---|---|---|
| prcp | annual | 835.20 | 803.11 | 32.08 | 4.0 | 0.308 | 0.197 |
| prcp | fall | 239.18 | 230.76 | 8.42 | 3.7 | 0.622 | 0.493 |
| prcp | spring | 155.72 | 150.62 | 5.10 | 3.4 | 0.643 | 0.164 |
| prcp | summer | 266.04 | 242.15 | 23.89 | 9.9 | 0.283 | 0.185 |
| prcp | winter | 174.26 | 179.59 | -5.33 | -3.0 | 0.674 | 0.661 |
| rh | annual | 74.08 | 74.08 | 0.00 | 0.0 | 0.992 | 0.862 |
| rh | fall | 79.38 | 80.16 | -0.78 | -1.0 | 0.279 | 0.091 |
| rh | spring | 68.43 | 69.49 | -1.06 | -1.5 | 0.159 | 0.385 |
| rh | summer | 67.25 | 67.63 | -0.38 | -0.6 | 0.760 | 0.756 |
| rh | winter | 81.25 | 79.05 | 2.20 | 2.8 | 0.004 | 0.003 |
| snow_cover | annual | 59.89 | 63.84 | -3.95 | NA | 0.001 | 0.003 |
| snow_cover | fall | 49.76 | 52.03 | -2.27 | NA | 0.424 | 0.405 |
| snow_cover | spring | 88.85 | 93.80 | -4.95 | NA | 0.004 | 0.001 |
| snow_cover | summer | 4.22 | 12.80 | -8.58 | NA | 0.000 | 0.001 |
| snow_cover | winter | 96.73 | 96.71 | 0.02 | NA | 0.238 | 0.282 |
| snowfall | annual | 405.53 | 432.37 | -26.84 | -6.2 | 0.089 | 0.253 |
| snowfall | fall | 136.62 | 135.00 | 1.61 | 1.2 | 0.895 | 0.934 |
| snowfall | spring | 94.08 | 109.02 | -14.94 | -13.7 | 0.074 | 0.459 |
| snowfall | summer | 4.87 | 12.02 | -7.15 | -59.5 | 0.002 | 0.002 |
| snowfall | winter | 169.96 | 176.32 | -6.36 | -3.6 | 0.619 | 0.708 |
| snowfall_fraction | annual | 48.45 | 53.56 | -5.10 | NA | 0.005 | 0.008 |
| snowmelt | annual | 381.07 | 409.27 | -28.20 | -6.9 | 0.139 | 0.528 |
| snowmelt | fall | 25.92 | 30.93 | -5.02 | -16.2 | 0.214 | 0.421 |
| snowmelt | spring | 289.89 | 212.24 | 77.65 | 36.6 | 0.002 | 0.001 |
| snowmelt | summer | 63.82 | 165.05 | -101.23 | -61.3 | 0.001 | 0.007 |
| snowmelt | winter | 1.44 | 1.05 | 0.39 | 37.2 | 0.599 | 0.599 |
| snowmelt_doy_50 | annual | 137.74 | 148.50 | -10.75 | NA | 0.000 | 0.002 |
| snowmelt_rate_peak | annual | 118.44 | 116.52 | 1.93 | 1.7 | 0.805 | 0.385 |
| soil_moisture | annual | 0.32 | 0.32 | 0.00 | -0.3 | 0.749 | 0.898 |
| soil_moisture | fall | 0.32 | 0.32 | 0.00 | -0.4 | 0.821 | 0.996 |
| soil_moisture | spring | 0.32 | 0.32 | 0.01 | 1.9 | 0.029 | 0.015 |
| soil_moisture | summer | 0.33 | 0.33 | -0.01 | -2.5 | 0.132 | 0.173 |
| soil_moisture | winter | 0.31 | 0.31 | 0.00 | -0.2 | 0.819 | 0.498 |
| swe | annual | 121.65 | 135.11 | -13.46 | -10.0 | 0.104 | 0.264 |
| swe | fall | 29.15 | 30.24 | -1.09 | -3.6 | 0.751 | 0.985 |
| swe | spring | 259.56 | 292.80 | -33.24 | -11.4 | 0.100 | 0.272 |
| swe | summer | 5.27 | 21.50 | -16.23 | -75.5 | 0.000 | 0.003 |
| swe | winter | 192.64 | 195.91 | -3.27 | -1.7 | 0.754 | 0.913 |
| swe_max | annual | 333.64 | 348.42 | -14.79 | -4.2 | 0.433 | 0.770 |
| tmax | annual | 3.66 | 1.99 | 1.67 | NA | 0.000 | 0.000 |
| tmax | fall | 3.49 | 2.55 | 0.93 | NA | 0.039 | 0.066 |
| tmax | spring | 3.93 | 2.16 | 1.76 | NA | 0.000 | 0.001 |
| tmax | summer | 16.62 | 15.00 | 1.62 | NA | 0.001 | 0.000 |
| tmax | winter | -9.40 | -11.77 | 2.37 | NA | 0.001 | 0.000 |
| tmean | annual | -0.59 | -2.41 | 1.82 | NA | 0.000 | 0.000 |
| tmean | fall | -0.32 | -1.58 | 1.26 | NA | 0.005 | 0.010 |
| tmean | spring | -0.75 | -2.57 | 1.82 | NA | 0.000 | 0.000 |
| tmean | summer | 11.64 | 9.76 | 1.88 | NA | 0.000 | 0.000 |
| tmean | winter | -12.92 | -15.27 | 2.34 | NA | 0.003 | 0.000 |
| tmin | annual | -4.12 | -6.05 | 1.93 | NA | 0.000 | 0.000 |
| tmin | fall | -3.19 | -4.71 | 1.52 | NA | 0.001 | 0.003 |
| tmin | spring | -4.76 | -6.49 | 1.74 | NA | 0.000 | 0.000 |
| tmin | summer | 6.93 | 4.81 | 2.12 | NA | 0.000 | 0.000 |
| tmin | winter | -15.47 | -17.81 | 2.33 | NA | 0.004 | 0.000 |
| vpd | annual | 2.14 | 1.86 | 0.28 | 15.1 | 0.002 | 0.001 |
| vpd | fall | 1.48 | 1.32 | 0.16 | 12.3 | 0.065 | 0.032 |
| vpd | spring | 2.05 | 1.67 | 0.38 | 22.8 | 0.000 | 0.000 |
| vpd | summer | 4.62 | 4.06 | 0.55 | 13.6 | 0.045 | 0.034 |
| vpd | winter | 0.43 | 0.40 | 0.03 | 6.5 | 0.023 | 0.001 |
Daytime Highs and Overnight Lows
The cd package ships daytime maximum (tmax) and overnight minimum (tmin) temperatures alongside the daily mean. They carry distinct information. Overnight minimums warming faster than daytime maximums — the “day-night asymmetry” — is one of the textbook fingerprints of greenhouse warming. Karl et al. (1993) documented this first at the global scale, finding overnight minimums rose at roughly three times the rate of daytime maximums between 1951 and 1990 (0.84 vs 0.28 °C). Whether a watershed or region shows that signal depends on local geography (valley inversions, snow cover, slope-aspect mix).
For the FWCP Peace Region, overnight minimums are warming faster than daytime maximums — the textbook day-night asymmetry. Daytime maximums warmed about +0.027 °C per year since 1951 (+2.0 °C cumulative), while overnight minimums warmed about +0.032 °C per year (+2.4 °C cumulative). The overnight side warmed roughly 0.4 °C more than the daytime side over the full record, narrowing the diurnal temperature range by the same amount. The three figures below show the tmax, tmin and diurnal-range time series that yield those numbers.
trn_tmax <- cd::cd_trend(
ano[ano$variable == "tmax" & ano$period == "annual", ],
trend_start = c(1951, 1981)
)
cd::cd_plot_timeseries(ano, variable = "tmax", period = "annual", trend = trn_tmax,
title = "Daytime Maximum (tmax) — Annual Anomaly")
Annual daytime maximum temperature (tmax) anomaly for the FWCP Peace Region relative to the 1951-1980 baseline.
trn_tmin <- cd::cd_trend(
ano[ano$variable == "tmin" & ano$period == "annual", ],
trend_start = c(1951, 1981)
)
cd::cd_plot_timeseries(ano, variable = "tmin", period = "annual", trend = trn_tmin,
title = "Overnight Minimum (tmin) — Annual Anomaly")
Annual overnight minimum temperature (tmin) anomaly for the FWCP Peace Region relative to the 1951-1980 baseline.
Diurnal temperature range (daytime maximum minus overnight minimum) annual mean for the FWCP Peace Region. The downward trend indicates overnight lows are warming faster than daytime highs — the textbook day-night asymmetry shows up here.
Snowpack
Snowpack is the hinge of BC hydrology: winter precipitation falls as snow, accumulates on the ground, and releases as meltwater across spring and summer. That seasonal storage is the difference between a late-summer creek that still flows and one that doesn’t. It also sets the timing of the spring freshet — the annual flood pulse that salmonids time their up-river migrations to. Cordillera-wide, snowpack has been declining for decades (Mote et al. 2018; Pederson et al. 2011). For four BC river basins (Fraser, Peace, Columbia, Campbell), Najafi et al. (2017) attribute observed spring SWE decline to anthropogenic forcing using formal detection-attribution. For the Fraser specifically, Kang et al. (2016) document a ~10-day advance of the spring freshet over 1949-2006, with declining summer flows during the salmon migration window.
For the FWCP Peace Region, the snowpack signal in our data is unambiguous and concentrated in the warm shoulders of the year. Annual snow water equivalent (SWE) declined ~10% (135 → 122 mm), but the seasonal breakdown is sharper: summer SWE collapsed by 75% (21.5 → 5.3 mm) and spring snowmelt rose 37% (212 → 290 mm) as the freshet shifted earlier into the calendar. Total annual snowfall barely changed (-6%) — the SWE decline is mostly about warming removing snow earlier, not less snow falling.
A few notes on how to read the snow numbers below.
ERA5-Land represents snow on a roughly 9 km grid — each grid cell is a single number summarizing snow averaged over about 80 km² of mixed terrain (saddles, slopes, valleys, forest, exposed alpine, all combined). When we compare the model against a snow station inside that cell, two distinct errors can stack: a scale mismatch (a single point measurement isn’t the same thing as an 80 km² average, especially in mountain terrain where snow accumulation varies sharply over short distances), and a cell-mean bias that the model has at the Northern Hemisphere scale (ERA5-Land overestimates mountain SWE by 150-200% even when compared against area-averaged satellite estimates (Kouki et al. 2023), traced to a simplified snow layer in the underlying atmospheric model).
Our QA cross-check at four British Columbia automated snow stations inside the FWCP Peace Region (1985-2025, 95 paired station-years) sees both kinds of error stacked and cannot fully separate them. At Pine Pass — a Coast Mountain saddle that catches orographic precip from systems coming inland from the Pacific — the model is 61% too low at the station, likely mostly scale mismatch: the station sits at a high-snow microsite within an averaged cell. At Aiken Lake — a drier interior valley basin — the model is 54% too high, likely a mix of scale mismatch (averaging snowier surrounding terrain into the same cell) and the Kouki-style cell-mean bias for interior-continental cells.
What does survive both kinds of error: the gap between the model and the stations is the same size in 2020 as it was in 1990. The bias is stable over time at every site (regression of model-minus-station on year is non-significant, p > 0.2). That means the changes over time shown below — peak snowpack dropping, freshet shifting earlier — are real even though the absolute mm numbers shouldn’t be quoted as ground truth. For regional aggregates (the numbers in this section are spatial means over hundreds of cells), random point-vs-cell mismatches partially cancel, and stable cell-mean biases preserve the trend signal.
The trend tests below use raw Mann-Kendall plus Theil-Sen — no pre-whitening — which is the right call for our 76-year series with strong climate trends per Yue and Wang (2002) (pre-whitening underestimates slope when a real trend exists).
Seasonal snowpack curve
The four monthly-native snow variables — SWE, snowfall, snowmelt, and snow cover — show when snow accumulates and melts. Aggregated to the four standard meteorological seasons — winter (December–February), spring (March–May), summer (June–August), and fall (September–November) — the table below compares the recent decade (2015-2025) against the pre-warming reference (1951-1980) directly. The headline numbers above (summer SWE collapse, spring snowmelt rise) are in the summer and spring rows.
| Variable | Season | Pre-warming (1951–1980) | Recent (2015–2025) | Δ % |
|---|---|---|---|---|
| SWE (mm) | winter | 195.9 | 192.6 | -1.7 |
| SWE (mm) | spring | 292.8 | 259.6 | -11.4 |
| SWE (mm) | summer | 21.5 | 5.3 | -75.5 |
| SWE (mm) | fall | 30.2 | 29.1 | -3.6 |
| Snowfall (mm) | winter | 176.3 | 170.0 | -3.6 |
| Snowfall (mm) | spring | 109.0 | 94.1 | -13.7 |
| Snowfall (mm) | summer | 12.0 | 4.9 | -59.5 |
| Snowfall (mm) | fall | 135.0 | 136.6 | 1.2 |
| Snowmelt (mm) | winter | 1.1 | 1.4 | 37.2 |
| Snowmelt (mm) | spring | 212.2 | 289.9 | 36.6 |
| Snowmelt (mm) | summer | 165.0 | 63.8 | -61.3 |
| Snowmelt (mm) | fall | 30.9 | 25.9 | -16.2 |
Annual climate-departure signals
Four derived annual scalars capture the climate-departure signals the literature treats as headline metrics for snow hydrology: peak snowpack (Mote et al. 2018, 2005), the date of melt midpoint (Stewart et al. 2005; Cayan et al. 2001), freshet flashiness, and the fraction of precipitation falling as snow (Knowles et al. 2006). Two notes on our specific implementation. We use the actual annual maximum of daily SWE rather than the April-1 SWE canon (Pederson et al. 2011) — equivalent in effect for BC pixels (peak is at or near April 1) but date-insensitive. And our freshet-flashiness metric (annual maximum of 7-day rolling daily snowmelt) is upstream of the streamflow-based flashiness measures in the literature (Stewart et al. 2005; Kang et al. 2016) — diagnostic of snowpack-side intensity before routing through soil and channel storage.
trn_swe_max <- cd::cd_trend(
ano[ano$variable == "swe_max" & ano$period == "annual", ],
trend_start = c(1951, 1981)
)
cd::cd_plot_timeseries(ano, variable = "swe_max", period = "annual",
trend = trn_swe_max,
title = "Annual peak SWE — Anomaly")
Annual peak snow water equivalent (swe_max) for the FWCP Peace Region. ERA5-Land mm SWE (regional spatial mean).
trn_doy <- cd::cd_trend(
ano[ano$variable == "snowmelt_doy_50" & ano$period == "annual", ],
trend_start = c(1951, 1981)
)
cd::cd_plot_timeseries(ano, variable = "snowmelt_doy_50", period = "annual",
trend = trn_doy,
title = "Snowmelt 50% DOY — Anomaly")
Day of year when half the annual snowmelt has accumulated (snowmelt_doy_50). Earlier dates indicate an earlier freshet centroid.
trn_rate <- cd::cd_trend(
ano[ano$variable == "snowmelt_rate_peak" & ano$period == "annual", ],
trend_start = c(1951, 1981)
)
cd::cd_plot_timeseries(ano, variable = "snowmelt_rate_peak", period = "annual",
trend = trn_rate,
title = "Peak weekly melt rate — Anomaly")
Annual maximum of 7-day rolling daily snowmelt (snowmelt_rate_peak). Higher values indicate more concentrated freshet pulses.
trn_frac <- cd::cd_trend(
ano[ano$variable == "snowfall_fraction" & ano$period == "annual", ],
trend_start = c(1951, 1981)
)
cd::cd_plot_timeseries(ano, variable = "snowfall_fraction", period = "annual",
trend = trn_frac,
title = "Snowfall fraction — Anomaly")
Annual snowfall fraction (snowfall_fraction): the percent of annual precipitation that fell as snow. Anomaly is in percentage points.
What this means for the FWCP Peace Region
Three findings carry the snowpack story for the Peace.
Snow is leaving the region earlier each year, not falling less. Total annual snowfall is roughly stable (-6%); the snowpack decline is dominated by earlier melt, not by less snow. The seasonal data make this concrete: spring snowmelt is up 37% while summer snowmelt is down 61% — same total melt, redistributed earlier into the year. This matches the broader Pacific-Northwest signal documented since Mote et al. (2005).
The freshet is shifting into spring. Spring snowmelt up 37% is direct evidence of an earlier freshet centroid — consistent in direction and rough magnitude with Stewart et al.’s (2005) 1-4 week earlier streamflow timing across western North America, and with the 10-day Fraser freshet advance documented by Kang et al. (2016) for a basin whose southern boundary is just south of the FWCP Peace Region.
Summer is becoming snow-free. Summer SWE has collapsed by 75% and summer snowmelt is down 61%. The high-elevation snowpack that historically lingered into summer no longer does in the recent decade. For aquatic ecosystems downstream, this is a loss of late- season cold-water input to streams during the warmest part of the year.
Spatial Pattern
The zonal mean reduces a region this size to a single number; the spatial pattern carries the rest of the story.
Warming is not spatially uniform. Total departures range from about +1.5 °C in the southeast to over +2.1 °C in the west and centre. The dominant gradient runs east-to-west — the western half of the region averages roughly 0.2 °C more warming than the eastern half — with a smaller south-to-north component. The west-of-Rockies pattern is consistent with windward-slope amplification along the Continental Divide.
Spatial pattern of annual mean temperature departure across the FWCP Peace Region (2015-2025 mean minus 1951-1980 mean, degrees C).
Per-Ecoregion Variation
The regional zonal mean averages over five ecoregions with different elevations and exposures. To check whether the regional story holds within each ecoregion — and where it does not — we run the same pipeline on each ecoregion polygon individually.
All five ecoregions warmed at roughly the same rate — about +0.3 °C per decade since 1951, with cumulative changes within 0.2 °C of each other. Precipitation is the variable that diverges: the two northernmost ecoregions (Boreal Mountains and Plateaus, Northern Canadian Rocky Mountains) show statistically significant precipitation increases, while the southern and central ecoregions do not. Vapour pressure deficit is up significantly in all five.
# For each ecoregion polygon, run the same cd pipeline as the regional one:
results <- lapply(seq_len(nrow(ecoregions)), function(i) {
poly <- ecoregions[i, ]
ts_i <- cd::cd_extract(catalog, poly)
bl_i <- cd::cd_baseline(ts_i, baseline_years = 1951:1980)
ano_i <- cd::cd_anomaly(ts_i, bl_i)
trn_i <- cd::cd_trend(ano_i, trend_start = c(1951, 1981))
cmp_i <- cd::cd_compare(ts_i, window_a = 2015:2025, window_b = 1951:1980)
list(ano = ano_i, trn = trn_i, cmp = cmp_i)
})
Annual mean temperature anomaly relative to the 1951-1980 baseline, by ecoregion. Dashed line is the 75-year Theil-Sen trend (1951-present); solid line is the 45-year trend (1981-present). A solid line steeper than the dashed line indicates accelerating warming.
Annual precipitation anomaly (% of the 1951-1980 baseline) by ecoregion. Dashed line is the 75-year trend (1951-present); solid line is the 45-year trend (1981-present). The two northernmost ecoregions (BMP, NRM) show statistically significant precipitation increases over the full record; the others do not.
Snow per ecoregion
Two snow metrics are worth viewing per-ecoregion: peak snow water equivalent (annual snowpack magnitude) and DOY-50 (the day-of-year by which half the year’s snowmelt has already happened — earlier values mean an earlier freshet). The two figures below show how each varies ecoregion to ecoregion. Pair them with the Watershed Groups Across Ecoregions section below to read each watershed group’s dominant ecoregion’s signal.
Annual peak snow water equivalent (SWE) anomaly by ecoregion, relative to the 1951-1980 baseline. Bars are mm of water-equivalent snowpack departure from baseline. Dashed line is the 75-year Theil-Sen trend (1951-present); solid line is the 45-year trend (1981-present).
Snowmelt 50% day-of-year (DOY-50) anomaly by ecoregion, relative to the 1951-1980 baseline. Negative values (red) mean the freshet midpoint shifted earlier in the year. Dashed line is the 75-year Theil-Sen trend; solid line is the 45-year trend.
| Ecoregion | tmean degC/dec | tmax degC/dec | tmin degC/dec | prcp mm/yr | prcp p | vpd hPa/dec | vpd p | prcp pct change | soil moisture pct change |
|---|---|---|---|---|---|---|---|---|---|
| FAB | 0.29 | 0.26 | 0.30 | 0.032 | 0.634 | 0.061 | 0.002 | 0.8 | -1.3 |
| CRM | 0.29 | 0.25 | 0.30 | 0.077 | 0.253 | 0.045 | 0.001 | 4.1 | -0.5 |
| OMM | 0.32 | 0.29 | 0.34 | 0.052 | 0.421 | 0.054 | 0.000 | 2.5 | -0.7 |
| BMP | 0.30 | 0.27 | 0.33 | 0.144 | 0.023 | 0.033 | 0.003 | 6.3 | 0.3 |
| NRM | 0.29 | 0.26 | 0.31 | 0.156 | 0.015 | 0.030 | 0.004 | 6.6 | 0.6 |
Ecoregion codes used in the table above: FAB — Fraser Basin; CRM — Central Canadian Rocky Mountains; OMM — Omineca Mountains; BMP — Boreal Mountains and Plateaus; NRM — Northern Canadian Rocky Mountains.
Watershed Groups Across Ecoregions
The FWCP Peace Region is reported and managed at the watershed-group scale — the British Columbia Freshwater Atlas hydrological reporting unit, and the unit the Fish & Wildlife Compensation Program funds work in. Sixteen watershed groups make up the canonical FWCP Peace list. They span the five ecoregions in different ways: some sit entirely within one, others split across two or three. The map below shows the watershed groups labelled with their codes, on top of the ecoregion fills. The table that follows gives, for each watershed group, the share of its area falling in each ecoregion.
The Upper Peace River watershed group is also intersected by the FWCP boundary, but only at its upstream end — its main drainage runs east beyond the region, so it is excluded here to keep the reporting unit clean.
Watershed groups intersecting the FWCP Peace Region (full extent, including the parts spilling outside the FWCP boundary), labelled with their codes, on top of ecoregion fills. The Fraser Basin, Central Canadian Rockies, and Omineca Mountains ecoregions occupy the southern and central portions; Boreal Mountains and Plateaus and Northern Canadian Rocky Mountains sit in the north.
| Code | Name | km² in FWCP | FAB % | CRM % | OMM % | BMP % | NRM % |
|---|---|---|---|---|---|---|---|
| CARP | Carp Lake | 1794 | 100.0 | 0.0 | 0.0 | 0.0 | 0.0 |
| CRKD | Crooked River | 2172 | 100.0 | 0.0 | 0.0 | 0.0 | 0.0 |
| FINA | Finlay Arm | 7383 | 0.3 | 20.2 | 61.8 | 10.3 | 7.4 |
| FINL | Finlay River | 5479 | 0.0 | 0.0 | 0.0 | 30.4 | 69.6 |
| FIRE | Firesteel River | 4376 | 0.0 | 0.0 | 1.0 | 99.0 | 0.0 |
| FOXR | Fox River | 4279 | 0.0 | 0.0 | 0.0 | 16.8 | 83.2 |
| INGR | Ingenika River | 5328 | 0.0 | 0.0 | 0.7 | 99.3 | 0.0 |
| LOMI | Lower Omineca River | 4013 | 0.0 | 0.0 | 100.0 | 0.0 | 0.0 |
| MESI | Mesilinka River | 3294 | 0.0 | 0.0 | 94.8 | 5.2 | 0.0 |
| NATR | Nation River | 6889 | 49.3 | 0.0 | 50.7 | 0.0 | 0.0 |
| OSPK | Ospika River | 2973 | 0.0 | 60.5 | 0.0 | 0.0 | 39.5 |
| PARA | Parsnip Arm | 3729 | 15.1 | 31.3 | 53.5 | 0.0 | 0.0 |
| PARS | Parsnip River | 5583 | 28.6 | 71.4 | 0.0 | 0.0 | 0.0 |
| PCEA | Peace Arm | 5861 | 0.0 | 98.8 | 1.2 | 0.0 | 0.0 |
| TOOD | Toodoggone River | 4849 | 0.0 | 0.0 | 0.0 | 100.0 | 0.0 |
| UOMI | Upper Omineca River | 3905 | 0.0 | 0.0 | 100.0 | 0.0 | 0.0 |
Interpretation
Three findings emerge from the climate departures shown above.
Warming is broad, fast, and uniform. All five ecoregions warmed between +1.7 °C and +1.9 °C cumulative since 1951, at a rate near +0.3 °C per decade, with Mann-Kendall p-values below 0.001 in every ecoregion. The trend is statistically significant beyond reasonable doubt of being random noise. Daily maximum, daily minimum, and daily mean temperatures all tell the same story. There is no thermal hot spot. The whole region has moved together. A subtle gradient does emerge — the northern, higher-elevation ecoregions (Omineca Mountains, Boreal Mountains and Plateaus) warmed about 0.2 °C more than the southern Fraser Basin, consistent with the elevation-dependent warming signal documented at mid-latitude mountain sites (Pepin et al. 2015), though the regional evidence base remains heterogeneous and not every mountain region shows the same pattern (Rangwala and Miller 2012). The gradient is small relative to the regional signal.
Precipitation is increasing significantly only in the two northernmost ecoregions. Across the FWCP Peace as a whole, the long-term trend in annual precipitation is not statistically significant. Broken out by ecoregion, Boreal Mountains and Plateaus and Northern Canadian Rocky Mountains — both in the northern half of the region — show real increases in annual precipitation (p = 0.02 and 0.01). The southern and central ecoregions show no significant change. This is a north-versus-south contrast that the regional average hides, and the place where breaking the analysis down by ecoregion was worth doing.
Mapping that contrast onto the watershed-group reporting unit: five watershed groups sit in or are dominated by the wetter ecoregions — Toodoggone, Firesteel and Ingenika sit almost entirely within Boreal Mountains and Plateaus, and Finlay and Fox are mostly in Northern Canadian Rocky Mountains. For these, the small upward precipitation trend is statistically defensible. The remaining 13 watershed groups sit elsewhere — Fraser Basin, Central Canadian Rockies, Omineca Mountains — where annual precipitation has not changed in a way distinguishable from natural variability.
The atmosphere is drying despite, in places, more precipitation. Vapour pressure deficit — the gap between how much water the air could hold and how much it actually does — rose significantly in every ecoregion (p < 0.005 across all five), mirroring the continental-scale drying that Ficklin and Novick (2017) documented for the United States as a whole, driven by combined air-temperature increases and relative-humidity changes. Warmer air holds more water before saturating, and that pulls moisture out of soil and vegetation through evaporation and transpiration. Soil moisture is essentially flat in every ecoregion — even where precipitation increased — because the warmer atmosphere is drinking the extra water back. This is the headline ecological finding for the region: water inputs may be rising in places, but soil and vegetation are not seeing more available moisture.
Snow is leaving the region earlier, not falling less. Annual snowfall across the FWCP Peace is essentially unchanged (-6%); the snowpack story is about timing, not quantity. Spring snowmelt rose 37% (212 → 290 mm) and summer snowmelt fell 61% (165 → 64 mm) — same total annual melt, redistributed earlier in the calendar. Annual peak snow water equivalent (SWE) is down about 10% (135 → 122 mm), and most strikingly, summer SWE has collapsed by 75%: the high-elevation snowpack that historically lingered into June and July no longer does. This is the direct mechanism behind the freshet shift described in the snow-trend literature: warmer winters and springs do not stop snow from falling — they shorten how long it stays on the ground.
The freshet-timing signal is broadly uniform across the region — DOY-50 (day-of-year by which half the year’s snowmelt has happened) is shifting earlier at roughly 1 day per decade in every one of the five ecoregions, with Mann-Kendall p-values below 0.01 in each (median p ≈ 0.005). All sixteen FWCP Peace watershed groups, regardless of which ecoregion they sit in, are seeing the same order-of-magnitude earlier melt. The peak-SWE signal is noisier — small recent-decade declines show up in every ecoregion, but the year-to-year variability is large enough that the long-term Theil-Sen slopes do not reach statistical significance at the ecoregion scale. The 75% summer-SWE collapse is detectable in the regional aggregate because spatial averaging suppresses that variability.
For the cold-water resident salmonids the FWCP Peace supports — bull trout, Arctic grayling, mountain whitefish, rainbow trout, kokanee — these signals compound. Stream temperatures are likely rising in step with warmer ambient air temperatures. Mantua et al. (2010) modelled this for Washington State watersheds and found that warming summer stream temperatures combined with altered low flows would reduce reproductive success for many salmon populations, and Eaton and Scheller (1996) showed across 57 US fish species that cold- and cool-water species lose substantially more thermal habitat than warm-water species under the same forcing. Both findings carry directly into the BC interior cold-water cohort the FWCP Peace supports. The evapotranspiration imbalance means low-flow conditions in late summer are not being relieved by the precipitation increase that did occur; and the cold-water input that high-elevation snowpack provides to streams during the warmest, most thermally stressful weeks of summer is dropping in parallel with summer SWE. The spring freshet — the dominant high-flow event that shapes channel morphology, mobilizes spawning gravels, and refills off-channel rearing habitat — is shifting weeks earlier; the neighbouring Fraser Basin shows the same shift — Kang et al. (2016) documented a ~10-day advance in the spring freshet at Hope between 1949 and 2006, with persistent declines through autumn recession just when salmon are migrating up the Fraser.