In Deltares/weathergenr: Gridded semi-parametric weather generator

knitr::opts_chunk$set(collapse = TRUE, comment = "#>", message=FALSE)

Overview

This vignette shows how the climate forcing data for a climate stress test can be prepared with weathergenr.

library(weathergenr)
ncfile <- system.file("extdata", "ntoum_era5_data.nc", package = "weathergenr")
ncdata <- readNetcdf(ncfile)

Define stochastic weather realizations

# Set path to store weather generator results
output_path <- "C:/testrun/"
variables <- c("precip", "temp", "temp_min", "temp_max")
realization_num <- 3

stochastic_weather <- generateWeatherSeries(
     weather.data = ncdata$data,
     weather.grid = ncdata$grid,
     weather.date = ncdata$date,
     variable.names = variables,
     variable.labels = variables,
     variable.units = NULL,
     sim.year.num = 20,
     sim.year.start = 2020,
     month.start = 1,
     realization.num = realization_num,
     warm.variable = "precip",
     warm.signif.level = 0.90,
     warm.sample.num = 10000,
     warm.subset.criteria = NULL,
     knn.sample.num = 120,
     mc.wet.quantile= 0.2,
     mc.extreme.quantile = 0.8,
     evaluate.model = FALSE,
     evaluate.grid.num = 20,
     output.path = output_path,
     seed = 123)

Define climate change domain

Metereological input to a climate stress test can be generated by obtaining a wide range of plausible weather conditions. \

In this tutorial, we will generate 6 climate change scenarios by combining 3 monthly precipitation changes and 2 monthly temperature changes. We will then apply these scenarios to 3 natural variability realizations, resulting in a total of 6 x 3 = 18 scenarios.\

Defining the stress testing matrix

The first step is define a data table to store all information regarding the scenarios, e.g., how weather statistics are being changed. To do this, we first define a bandwith range (e.g., minimum and maximum values) to define the delta factors for each climate statistic being perturbed. In the example below, we vary three statistics: mean of precipitation, variance of precipitation, and mean of temperature.

# Temp mean changes Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
delta_temp_mean_min <- c(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)
delta_temp_mean_max <- c(3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0)

# Precip mean changes   Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
delta_precip_mean_min <- c(0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7, 0.7)
delta_precip_mean_max <- c(1.3, 1.3, 1.3, 1.3, 1.3, 1.3, 1.3, 1.3, 1.3, 1.3, 1.3, 1.3)

# Precip variance changes   Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec
delta_precip_variance_min <- c(1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0)
delta_precip_variance_max <- c(1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0)

# Number of incremental step changes for precip and temp variables
precip_step_num <- 3
temp_step_num <- 2

precip_mean_steps <- sapply(1:12, function(m)
         seq(delta_precip_mean_min[m], delta_precip_mean_max[m],
             length.out = precip_step_num))

precip_variance_steps <- sapply(1:12, function(m)
         seq(delta_precip_variance_min[m], delta_precip_variance_max[m],
             length.out = precip_step_num))

temp_mean_steps <- sapply(1:12, function(m)
         seq(delta_temp_mean_min[m], delta_temp_mean_max[m],
             length.out = temp_step_num))

 df1 <- as.data.frame(precip_mean_steps) %>% mutate(level = 1:n(),
   variable = "precip_mean", .before = 1)
 df2 <- as.data.frame(precip_variance_steps) %>% mutate(level = 1:n(),
   variable = "precip_variance", .before = 1)
 df3 <- as.data.frame(temp_mean_steps) %>% mutate(level = 1:n(),
   variable = "temp_mean", .before = 1)
 df <- bind_rows(df1, df2, df3) %>% gather(month, value, V1:V12) %>%
   mutate(month = factor(month, levels = paste0("V",1:12), labels = 1:12))

 p <- ggplot2::ggplot(df, aes(x = month, y = value, group = level, color = level)) +
   facet_wrap(. ~ variable, scales = "free_y", ncol = 2) +
   geom_line() +
   labs(x="month", y = "delta factor") +
   scale_color_distiller(palette = "Set1") +
   guides(color = "none")

 p

Now lets create the scenario matrix using the bandwiths and incremental step sizes specified for each variable

# Stress test matrix
 strtest_matrix <- tidyr::expand_grid(stoc_ind = 1:realization_num,
   precip_ind = 1:precip_step_num, temp_ind = 1:temp_step_num)

 # Total number of scenarios
 smax <- nrow(strtest_matrix)

 # Stress test delta factors for each variable/climate statistic
 strtest_matrix_precip_mean <- precip_mean_steps[strtest_matrix$precip_ind, ]
 strtest_matrix_precip_variance <- precip_variance_steps[strtest_matrix$precip_ind, ]
 strtest_matrix_temp_mean <- temp_mean_steps[strtest_matrix$temp_ind, ]

# Write stress test matrices to file (optional)

 # write.csv(strtest_matrix,
 #   paste0(output_path, "strtest_matrix.csv"), row.names = FALSE)
 # write.csv(strtest_matrix_precip_mean,
 #   paste0(output_path, "strtest_matrix_precip_mean.csv"), row.names = FALSE)
 # write.csv(strtest_matrix_precip_variance,
 #   paste0(output_path, "strtest_matrix_precip_variance.csv"), row.names = FALSE)
 # write.csv(strtest_matrix_temp_mean,
 #   paste0(output_path, "strtest_matrix_temp_mean.csv"), row.names = FALSE)

Finally, lets generate the stress test input data

 # Read-in resampled dates & date series (from csv files included with the package)
 resampled_dates <- read.csv(system.file("extdata", "resampled_dates.csv", package = "weathergenr"),
   colClasses = "Date")
 sim_dates <- read.csv(system.file("extdata", "sim_dates.csv", package = "weathergenr"),
   colClasses = "Date")[[1]]

 # Use results from generateWeatherSeries function output
 # resampled_dates <- stochastic_weather$resampled
 # sim_dates <- stochastic_weather$dates

# progress bar (optional)
pb = txtProgressBar(min = 1, max = smax, initial = 0, style = 3)
 for (s in 1:smax) {

   setTxtProgressBar(pb,s)

   # Find the current scenario indices for the stochastic realization and delta factors
   stoc_ind <- strtest_matrix$stoc_ind[s]

   # Obtain stochastic series by re-ordering historical data
   day_order <- match(resampled_dates[[stoc_ind]], ncdata$date)
   rlz_historical <- lapply(ncdata$data, function(x) x[day_order,])

   # Apply climate changes to climate data
   rlz_future <- imposeClimateChanges(
     climate.data = rlz_historical,
     climate.grid = ncdata$grid,
     sim.dates = sim_dates,
     change.factor.precip.mean = strtest_matrix_precip_mean[s,],
     change.factor.precip.variance = strtest_matrix_precip_variance[s,],
     change.factor.temp.mean = strtest_matrix_temp_mean[s,],
     change.type.temp = "transient",
     change.type.precip = "transient")

     # Save to netcdf file
     # writeNetcdf(
     #   data = rlz_future,
     #   coord.grid = ncdata$grid,
     #   output.path = output_path,
     #   origin.date =  stochastic_weather$dates[1],
     #   calendar.type = "noleap",
     #   nc.template.file = ncfile,
     #   nc.compression = 4,
     #   nc.spatial.ref = "spatial_ref",
     #   nc.file.prefix = "climx",
     #   nc.file.suffix = s)
 }
 close(pb)