Yield change scenarios: based on theoretical models

In this script, we compute global yields change based on theoretical model provided by Agnolucci. We use the following variables:

acumulated growing season precipitation: the growing season depends of each product;
mean growing season temperature.

Yield changes are computed for product:

wheat,
barley,
potatoes,
oats,
rapeseed.

Load packages, read data and source custom scripts

rm(list = ls())
library(stars)

#> Loading required package: abind

#> Loading required package: sf

#> Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1

path_proj <- day2day::git_path()
path_data <- file.path(path_proj, "data")
path_cleaned <- file.path(path_data, "cleaned")
path_processed <- file.path(path_data, "processed")
path_modelled <- file.path(path_data, "modelled")

source(file.path(path_proj, "src", "72-yield-change.R"))

speed_meta <- readr::read_csv(file.path(path_cleaned, "ukcp18-speed_rcp85_bias_corrected_01_metadata.csv"))

#> Parsed with column specification:
#> cols(
#>   file = col_character(),
#>   variable = col_character(),
#>   summary = col_character(),
#>   temp_resolution = col_character(),
#>   from = col_double(),
#>   to = col_double(),
#>   from_year = col_double(),
#>   from_month = col_double(),
#>   to_year = col_double(),
#>   to_month = col_double(),
#>   description = col_character()
#> )

yields <- readRDS(file.path(path_cleaned, "yield-params-and-season.rds"))
products <- names(yields$season)

Compute yield change

periods <- unique(subset(speed_meta, select = c("from_year", "to_year")))

for (product in products) {
    custom_fun <- function(x, y, season) {
        yield_change_metadata(speed_meta, x, y, yields$params[[product]],
                              yields$season[[product]], product,
                              vars = c("pr", "tas"),
                              path_cleaned, path_modelled, path_processed)
    }
    purrr::walk2(periods$from_year, periods$to_year, ~ custom_fun(.x, .y, season))
}

Time to execute the task

Only useful when executed with Rscript.

proc.time()

#>    user  system elapsed 
#>  36.284   2.033  38.344