# \donttest{
# Simulate some time series that follow a latent VAR(1) process
simdat <- sim_mvgam(family = gaussian(),
n_series = 4,
trend_model = VAR(cor = TRUE),
prop_trend = 1)
plot_mvgam_series(data = simdat$data_train, series = 'all')
# Fit a model that uses a latent VAR(1)
mod <- mvgam(y ~ -1,
trend_formula = ~ 1,
trend_model = VAR(cor = TRUE),
family = gaussian(),
data = simdat$data_train,
chains = 2,
silent = 2)
# Calulate stability metrics for this system
metrics <- stability(mod)
# Proportion of stationary forecast distribution
# attributable to lagged interactions
hist(metrics$prop_int,
xlim = c(0, 1),
xlab = 'Prop_int',
main = '',
col = '#B97C7C',
border = 'white')
# Within this contribution of interactions, how important
# are inter-series interactions (offdiagonals of the A matrix) vs
# intra-series density dependence (diagonals of the A matrix)?
layout(matrix(1:2, nrow = 2))
hist(metrics$prop_int_offdiag,
xlim = c(0, 1),
xlab = '',
main = 'Inter-series interactions',
col = '#B97C7C',
border = 'white')
hist(metrics$prop_int_diag,
xlim = c(0, 1),
xlab = 'Contribution to interaction effect',
main = 'Intra-series interactions (density dependence)',
col = 'darkblue',
border = 'white')
layout(1)
# How important are inter-series error covariances
# (offdiagonals of the Sigma matrix) vs
# intra-series variances (diagonals of the Sigma matrix) for explaining
# the variance of the stationary forecast distribution?
layout(matrix(1:2, nrow = 2))
hist(metrics$prop_cov_offdiag,
xlim = c(0, 1),
xlab = '',
main = 'Inter-series covariances',
col = '#B97C7C',
border = 'white')
hist(metrics$prop_cov_diag,
xlim = c(0, 1),
xlab = 'Contribution to forecast variance',
main = 'Intra-series variances',
col = 'darkblue',
border = 'white')
layout(1)
# Reactivity, i.e. degree to which the system moves
# away from a stable equilibrium following a perturbation
# (values > 1 suggest a more reactive, less stable system)
hist(metrics$reactivity,
main = '',
xlab = 'Reactivity',
col = '#B97C7C',
border = 'white',
xlim = c(-1*max(abs(metrics$reactivity)),
max(abs(metrics$reactivity))))
abline(v = 0, lwd = 2.5)
# }
Run the code above in your browser using DataLab