Applying Spatio-Temporal Model to Crop Yield Data
This example demonstrates the fitting of a spatio-temporal model to predict crop yields in Chatham-Kent, Ontario, and compares the predictions to real data.
``` r
bsts_1 <- fit_bsts(y_train[[medoid_names[1]]], z_train[[medoid_names[1]]], lags = 2, MCMC.iter = 10) # Chatham-Kent
#> =-=-=-=-= Iteration 0 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 1 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 2 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 3 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 4 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 5 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 6 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 7 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 8 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 9 Fri Jul 10 06:44:36 2026 =-=-=-=-=
bsts_2 <- fit_bsts(y_train[[medoid_names[2]]], z_train[[medoid_names[2]]], lags = 2, MCMC.iter = 10) # Wellington
#> =-=-=-=-= Iteration 0 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 1 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 2 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 3 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 4 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 5 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 6 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 7 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 8 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 9 Fri Jul 10 06:44:36 2026 =-=-=-=-=
bsts_3 <- fit_bsts(y_train[[medoid_names[3]]], z_train[[medoid_names[3]]], lags = 1, MCMC.iter = 10) # Lambton
#> =-=-=-=-= Iteration 0 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 1 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 2 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 3 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 4 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 5 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 6 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 7 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 8 Fri Jul 10 06:44:36 2026 =-=-=-=-=
#> =-=-=-=-= Iteration 9 Fri Jul 10 06:44:36 2026 =-=-=-=-=
list_bsts <- list(bsts_1, bsts_2, bsts_3)
GaussianForecast3 <- simul_fun_noGEV_3d(nsim = 10,
n_train = n_train,
n_test = n_test,
copula = copula_list[1],
init_params = init_params_full,
fn = log_likelihood_noGEV_3d,
u1 = u[[1]],
u2 = u[[2]],
u3 = u[[3]],
z1_train = z_train[[1]],
z2_train = z_train[[2]],
z3_train = z_train[[3]],
z1_test = z_test[[1]],
z2_test = z_test[[2]],
z3_test = z_test[[3]],
X_t = (z_train[[1]] + z_train[[2]] + z_train[[3]])/3,
y1_test = y_test[[1]],
y2_test = y_test[[2]],
y3_test = y_test[[3]],
BSTS_1 = bsts_1,
BSTS_2 = bsts_2,
BSTS_3 = bsts_3)
ChathamKent_copula_plot <- plot_forecast(forecast = GaussianForecast3[[2]],
data_train = y_train$`Chatham-Kent`,
data_test = y_test$`Chatham-Kent`,
time = time,
quant_high = 0.95,
quant_low = 0.05,
observed_col = "#96BA66",
forecast_col = "#CF9FFF",
title = "Chatham-Kent - Gaussian copula forecast")
#> Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
#> ℹ Please use `linewidth` instead.
#> ℹ The deprecated feature was likely used in the STCYP package.
#> Please report the issue to the authors.
#> This warning is displayed once per session.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
print(ChathamKent_copula_plot) 