Calculate measures of latent VAR community stability

Compute reactivity, return rates and contributions of interactions to stationary forecast variance from mvgam models with Vector Autoregressive dynamics

Usage

stability(object, ...)

# S3 method for mvgam
stability(object, ...)

Arguments

object: list object of class mvgam resulting from a call to mvgam() that used a Vector Autoregressive latent process model (either as VAR(cor = FALSE) or VAR(cor = TRUE))
...: ignored

Value

A data.frame containing posterior draws for each stability metric.

Details

These measures of stability can be used to assess how important inter-series dependencies are to the variability of a multivariate system and to ask how systems are expected to respond to environmental perturbations. Using the formula for a latent VAR(1) as: $$ \mu_t \sim \text{MVNormal}(A(\mu_{t - 1}), \Sigma) \quad $$ this function will calculate the long-term stationary forecast distribution of the system, which has mean $\mu_{\infty}$ and variance $\Sigma_{\infty}$, to then calculate the following quantities:

prop_int: Proportion of the volume of the stationary forecast distribution that is attributable to lagged interactions (i.e. how important are the autoregressive interaction coefficients in $A$ for explaining the shape of the stationary forecast distribution?): $$ det(A)^2 \quad $$
prop_int_adj: Same as prop_int but scaled by the number of series $p$ to facilitate direct comparisons among systems with different numbers of interacting variables: $$ det(A)^{2/p} \quad $$
prop_int_offdiag: Sensitivity of prop_int to inter-series interactions (i.e. how important are the off-diagonals of the autoregressive coefficient matrix $A$ for shaping prop_int?), calculated as the relative magnitude of the off-diagonals in the partial derivative matrix: $$ [2~det(A) (A^{-1})^T] \quad $$
prop_int_diag: Sensitivity of prop_int to intra-series interactions (i.e. how important are the diagonals of the autoregressive coefficient matrix $A$ for shaping prop_int?), calculated as the relative magnitude of the diagonals in the partial derivative matrix: $$ [2~det(A) (A^{-1})^T] \quad $$
prop_cov_offdiag: Sensitivity of $\Sigma_{\infty}$ to inter-series error correlations (i.e. how important are off-diagonal covariances in $\Sigma$ for shaping $\Sigma_{\infty}$?), calculated as the relative magnitude of the off-diagonals in the partial derivative matrix: $$ [2~det(\Sigma_{\infty}) (\Sigma_{\infty}^{-1})^T] \quad $$
prop_cov_diag: Sensitivity of $\Sigma_{\infty}$ to error variances (i.e. how important are diagonal variances in $\Sigma$ for shaping $\Sigma_{\infty}$?), calculated as the relative magnitude of the diagonals in the partial derivative matrix: $$ [2~det(\Sigma_{\infty}) (\Sigma_{\infty}^{-1})^T] \quad $$
reactivity: A measure of the degree to which the system moves away from a stable equilibrium following a perturbation. Values > 0 suggest the system is reactive, whereby a perturbation of the system in one period can be amplified in the next period. If $\sigma_{max}(A)$ is the largest singular value of $A$, then reactivity is defined as: $$ log\sigma_{max}(A) \quad $$
mean_return_rate: Asymptotic (long-term) return rate of the mean of the transition distribution to the stationary mean, calculated using the largest eigenvalue of the matrix $A$: $$ max(\lambda_{A}) \quad $$ Lower values suggest greater stability
var_return_rate: Asymptotic (long-term) return rate of the variance of the transition distribution to the stationary variance: $$ max(\lambda_{A \otimes{A}}) \quad $$ Again, lower values suggest greater stability

Major advantages of using mvgam to compute these metrics are that well-calibrated uncertainties are available and that VAR processes are forced to be stationary. These properties make it simple and insightful to calculate and inspect aspects of both long-term and short-term stability. But it is also possible to more directly inspect possible interactions among the time series in a latent VAR process. To do so, you can calculate Generalized or Orthogonalized Impulse Response Functions using the irf function, or you can calculate Forecast Error Variance Decompositions using the fevd function.

References

AR Ives, B Dennis, KL Cottingham & SR Carpenter (2003). Estimating community stability and ecological interactions from time-series data. Ecological Monographs. 73, 301-330.

Author

Nicholas J Clark

Examples

# \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)

# }