library(weathergenr)
library(ggplot2)
ncfile <- system.file("extdata", "ntoum_era5_data.nc", package = "weathergenr")
ncdata <- read_netcdf(ncfile)1 Purpose
This vignette demonstrates how to apply climate perturbations to daily weather series using weathergenr. The workflow modifies precipitation and temperature in a controlled way to reflect user-defined stress scenarios while preserving realistic daily variability. This is useful for hydrological stress testing, climate sensitivity analysis, and impact modeling.
2 Setup
3 Define climate scenario
We define a stress-test scenario representing a drier and warmer climate: 30% reduction in mean precipitation and +2 degrees Celsius temperature increase.
obs_date <- ncdata$date
obs_month <- as.integer(format(obs_date, "%m"))
obs_year <- as.integer(format(obs_date, "%Y"))
obs_year_idx <- obs_year - min(obs_year) + 1L
n_years <- max(obs_year_idx)
# Monthly perturbation vectors (length 12)
precip_mean_change <- rep(0.7, 12)
precip_var_change <- rep(1.0, 12)
temp_mean_change <- rep(2, 12)
# Expand to matrices (n_years x 12) for quantile mapping
precip_mean_change_m <- matrix(precip_mean_change, nrow = n_years, ncol = 12, byrow = TRUE)
precip_var_change_m <- matrix(precip_var_change, nrow = n_years, ncol = 12, byrow = TRUE)4 Preview quantile mapping
Before applying perturbations across all grid cells, we preview the quantile-mapping diagnostics on a single cell to verify target changes are achieved.
precip_qm <- adjust_precipitation_qm(
precip = ncdata$data[[1]]$precip,
mean_factor = precip_mean_change_m,
var_factor = precip_var_change_m,
scale_var_with_mean = TRUE,
exaggerate_extremes = FALSE,
enforce_target_mean = TRUE,
month = obs_month,
year = obs_year_idx,
min_events = 10,
validate_output = TRUE,
diagnostics = TRUE,
verbose = FALSE
)
print(precip_qm$diagnostics)
=== Precipitation QM Diagnostics ===
Mean and variance changes:
Statistic Original Adjusted Target Achieved Error
Mean 7.64 5.35 -30% -30% +0%
Variance 87.85 32.48 -51% -63% -12%
Wet/dry day counts:
Category Original Adjusted Difference
Dry days 161 183 22
Wet days 7510 7488 -22
Percentiles (wet days):
Percentile Original Adjusted Change
P50 (Median) 5.62 4.03 -28%
P90 16.13 11.41 -29%
P95 20.52 14.55 -29%
P99 39.35 25.23 -36%
Extreme tails:
Tail Threshold change Mean change
Top 10% (P90) -29% -34%
Top 5% (P95) -29% -36%
Top 1% (P99) -36% -43%
Top 0.1% (P99.9) -47% -51%
Spell lengths:
Type Original (days) Adjusted (days) Change
Dry spells 1.3 1.3 +3%
Wet spells 60.6 54.7 -10%
Monthly mean changes:
Month Target Achieved Error
Jan -30% -30% +0%
Feb -30% -30% +0%
Mar -30% -30% +0%
Apr -30% -30% +0%
May -30% -30% +0%
Jun -30% -30% +0%
Jul -30% -30% -0%
Aug -30% -30% +0%
Sep -30% -30% +0%
Oct -30% -30% +0%
Nov -30% -30% +0%
Dec -30% -30% +0%5 Apply climate perturbations
The apply_climate_perturbations() function processes all grid cells. It accepts a list of data frames (data), grid metadata (grid), and a date vector (date). Monthly perturbation factors control precipitation mean/variance scaling and temperature shifts. When diagnostic = TRUE, the function returns both the perturbed data and validation diagnostics.
rlz_future <- apply_climate_perturbations(
data = ncdata$data,
grid = ncdata$grid,
date = ncdata$date,
precip_mean_factor = precip_mean_change,
precip_var_factor = precip_var_change,
temp_delta = temp_mean_change,
precip_occurrence_factor = NULL,
precip_intensity_threshold = 0,
temp_transient = TRUE,
precip_transient = TRUE,
precip_occurrence_transient = TRUE,
compute_pet = FALSE,
pet_method = "hargreaves",
qm_fit_method = "mme",
scale_var_with_mean = TRUE,
exaggerate_extremes = FALSE,
extreme_prob_threshold = 0.95,
extreme_k = 1.2,
enforce_target_mean = TRUE,
precip_cap_mm_day = NULL,
precip_floor_mm_day = NULL,
precip_cap_quantile = NULL,
seed = 42,
diagnostic = TRUE,
verbose = FALSE
)The output contains data (list of perturbed data frames) and diagnostic (list of per-grid validation objects, precipitation only).
6 Visualizations
6.1 Monthly precipitation
obs_cell <- ncdata$data[[1]]
sim_cell <- rlz_future$data[[1]]
obs_monthly <- tapply(obs_cell$precip, obs_month, mean, na.rm = TRUE)
sim_monthly <- tapply(sim_cell$precip, obs_month, mean, na.rm = TRUE)
plot_data <- data.frame(
Month = factor(month.abb[1:12], levels = month.abb),
Observed = as.numeric(obs_monthly),
Perturbed = as.numeric(sim_monthly)
)
plot_data_long <- tidyr::pivot_longer(
plot_data, cols = c("Observed", "Perturbed"),
names_to = "Series", values_to = "Precipitation"
)
ggplot(plot_data_long, aes(x = Month, y = Precipitation, fill = Series)) +
geom_col(position = "dodge", width = 0.7) +
scale_fill_manual(values = c(Observed = "steelblue", Perturbed = "coral")) +
labs(
title = "Monthly Mean Precipitation",
subtitle = "Observed vs. perturbed (30% reduction scenario)",
y = "Precipitation (mm/day)", x = NULL
) +
theme_bw() +
theme(legend.position = "bottom")
6.2 Distribution comparison
wet_obs <- obs_cell$precip[obs_cell$precip > 0]
wet_sim <- sim_cell$precip[sim_cell$precip > 0]
dist_data <- data.frame(
Precipitation = c(wet_obs, wet_sim),
Series = factor(
c(rep("Observed", length(wet_obs)), rep("Perturbed", length(wet_sim))),
levels = c("Observed", "Perturbed")
)
)
ggplot(dist_data, aes(x = Precipitation, color = Series, fill = Series)) +
geom_density(alpha = 0.3, linewidth = 1) +
scale_color_manual(values = c(Observed = "steelblue", Perturbed = "coral")) +
scale_fill_manual(values = c(Observed = "steelblue", Perturbed = "coral")) +
scale_x_continuous(limits = c(0, quantile(c(wet_obs, wet_sim), 0.98))) +
labs(
title = "Wet-Day Precipitation Distribution",
x = "Precipitation (mm/day)", y = "Density"
) +
theme_bw() +
theme(legend.position = "bottom")
6.3 Temperature time series
plot_idx <- seq_len(nrow(obs_cell))
temp_ts <- data.frame(
Day = plot_idx,
Observed = obs_cell$temp[plot_idx],
Perturbed = sim_cell$temp[plot_idx]
)
temp_ts_long <- tidyr::pivot_longer(
temp_ts, cols = c("Observed", "Perturbed"),
names_to = "Series", values_to = "Temperature"
)
ggplot(temp_ts_long, aes(x = Day, y = Temperature, color = Series)) +
geom_line(alpha = 0.7) +
scale_color_manual(values = c(Observed = "steelblue", Perturbed = "coral")) +
labs(
title = "Daily Temperature Time Series",
subtitle = "Observed vs. +2 deg C shift scenario",
x = "Day", y = "Temperature (deg C)"
) +
theme_minimal() +
theme(legend.position = "bottom")
7 Use cases and extensions
7.1 Step changes
For immediate perturbations instead of transient ramps:
rlz_step <- apply_climate_perturbations(
data = ncdata$data, grid = ncdata$grid, date = ncdata$date,
precip_mean_factor = precip_mean_change,
precip_var_factor = precip_var_change,
temp_delta = temp_mean_change,
temp_transient = FALSE,
precip_transient = FALSE,
diagnostic = TRUE
)7.2 Occurrence perturbations
To modify wet-day frequency:
occurrence_factor <- rep(1.0, 12)
occurrence_factor[6:8] <- 0.8 # 20% fewer wet days in summer
rlz_dry_summer <- apply_climate_perturbations(
data = ncdata$data, grid = ncdata$grid, date = ncdata$date,
precip_mean_factor = precip_mean_change,
precip_var_factor = precip_var_change,
precip_occurrence_factor = occurrence_factor,
temp_delta = temp_mean_change,
diagnostic = TRUE
)7.3 PET recalculation
To propagate temperature changes to PET:
rlz_with_pet <- apply_climate_perturbations(
data = ncdata$data, grid = ncdata$grid, date = ncdata$date,
precip_mean_factor = precip_mean_change,
precip_var_factor = precip_var_change,
temp_delta = temp_mean_change,
compute_pet = TRUE,
pet_method = "hargreaves",
diagnostic = TRUE
)7.4 Extreme amplification
To amplify upper-tail changes:
rlz_extreme <- apply_climate_perturbations(
data = ncdata$data, grid = ncdata$grid, date = ncdata$date,
precip_mean_factor = precip_mean_change,
precip_var_factor = precip_var_change,
temp_delta = temp_mean_change,
exaggerate_extremes = TRUE,
extreme_prob_threshold = 0.95,
extreme_k = 1.2,
diagnostic = TRUE
)8 Session info
R version 4.5.2 (2025-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64
Running under: Windows 11 x64 (build 26200)
Matrix products: default
LAPACK version 3.12.1
locale:
[1] LC_COLLATE=English_United States.utf8
[2] LC_CTYPE=English_United States.utf8
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C
[5] LC_TIME=English_United States.utf8
time zone: Europe/Amsterdam
tzcode source: internal
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] ggplot2_4.0.1 weathergenr_1.1.0
loaded via a namespace (and not attached):
[1] generics_0.1.4 tidyr_1.3.2 class_7.3-23 lattice_0.22-7
[5] digest_0.6.39 magrittr_2.0.4 evaluate_1.0.5 grid_4.5.2
[9] RColorBrewer_1.1-3 waveslim_1.8.5 multitaper_1.0-17 iterators_1.0.14
[13] fastmap_1.2.0 foreach_1.5.2 jsonlite_2.0.0 Matrix_1.7-4
[17] e1071_1.7-17 survival_3.8-6 purrr_1.2.1 scales_1.4.0
[21] codetools_0.2-20 cli_3.6.5 rlang_1.1.7 splines_4.5.2
[25] withr_3.0.2 yaml_2.3.12 otel_0.2.0 tools_4.5.2
[29] dplyr_1.1.4 ncdf4_1.24 fitdistrplus_1.2-4 vctrs_0.6.5
[33] R6_2.6.1 proxy_0.4-29 lifecycle_1.0.5 MASS_7.3-65
[37] pkgconfig_2.0.3 pillar_1.11.1 gtable_0.3.6 glue_1.8.0
[41] xfun_0.55 tibble_3.3.1 tidyselect_1.2.1 rstudioapi_0.18.0
[45] knitr_1.51 farver_2.1.2 htmltools_0.5.9 patchwork_1.3.2
[49] labeling_0.4.3 rmarkdown_2.30 compiler_4.5.2 S7_0.2.1