Normal-inverse-Wishart Matrix

We provide functions to generate matrix-variate Normal and inverse-Wishart.

sim_mnormal(num_sim, mu, sig): num_sim of \(\mathbf{X}_i \stackrel{iid}{\sim} N(\boldsymbol{\mu}, \Sigma)\).
sim_matgaussian(mat_mean, mat_scale_u, mat_scale_v): One \(X_{m \times n} \sim MN(M_{m \times n}, U_{m \times m}, V_{n \times n})\) which means that \(vec(X) \sim N(vec(M), V \otimes U)\).
sim_iw(mat_scale, shape): One \(\Sigma \sim IW(\Psi, \nu)\).
sim_mniw(num_sim, mat_mean, mat_scale_u, mat_scale, shape): num_sim of \((X_i, \Sigma_i) \stackrel{iid}{\sim} MNIW(M, U, V, \nu)\).

Multivariate Normal generation gives num_sim x dim matrix. For example, generating 3 vector from Normal(\(\boldsymbol{\mu} = \mathbf{0}_2\), \(\Sigma = diag(\mathbf{1}_2)\)):

sim_mnormal(3, rep(0, 2), diag(2))
#>        [,1]   [,2]
#> [1,] -0.626  0.184
#> [2,] -0.836  1.595
#> [3,]  0.330 -0.820

The output of sim_matgaussian() is a matrix.

sim_matgaussian(matrix(1:20, nrow = 4), diag(4), diag(5))
#>      [,1] [,2]  [,3] [,4] [,5]
#> [1,] 1.49 5.74  9.58 12.7 18.5
#> [2,] 2.39 5.38  7.79 15.1 18.0
#> [3,] 2.98 7.94 11.82 15.6 19.9
#> [4,] 4.78 8.07 10.01 16.6 19.9

When generating IW, violating \(\nu > dim - 1\) gives error. But we ignore \(\nu > dim + 1\) (condition for mean existence) in this function. Nonetheless, we recommend you to keep \(\nu > dim + 1\) condition. As mentioned, it guarantees the existence of the mean.

sim_iw(diag(5), 7)
#>         [,1]    [,2]    [,3]    [,4]    [,5]
#> [1,]  0.1827  0.0894 -0.0411 -0.0924 -0.1305
#> [2,]  0.0894  0.6110  0.0860 -0.3754 -0.1015
#> [3,] -0.0411  0.0860  0.2189 -0.1649 -0.0988
#> [4,] -0.0924 -0.3754 -0.1649  0.4577  0.2136
#> [5,] -0.1305 -0.1015 -0.0988  0.2136  0.7444

In case of sim_mniw(), it returns list with mn (stacked MN matrices) and iw (stacked IW matrices):

sim_mniw(2, matrix(1:20, nrow = 4), diag(4), diag(5), 7)
#> $mn
#>      [,1] [,2]  [,3] [,4] [,5]  [,6] [,7]  [,8] [,9] [,10]
#> [1,] 1.19 4.96  8.03 12.7 14.6 0.323 5.07  8.24 13.1  17.4
#> [2,] 2.13 5.17 10.41 14.1 19.3 2.022 5.56  9.96 13.9  19.4
#> [3,] 3.01 7.15 10.29 14.9 18.1 3.380 7.24 11.34 15.5  17.9
#> [4,] 4.59 7.99 12.12 15.0 18.5 3.755 8.79 11.84 15.6  19.6
#> 
#> $iw
#>         [,1]     [,2]     [,3]    [,4]   [,5]     [,6]    [,7]     [,8]    [,9]
#> [1,]  0.2887 -0.04393  0.11892 -0.2959  0.110  0.28118 -0.1397  0.00529  0.0283
#> [2,] -0.0439  0.33357  0.00197  0.0106 -0.167 -0.13965  0.2876  0.06156 -0.1575
#> [3,]  0.1189  0.00197  1.30074 -0.0291  2.210  0.00529  0.0616  0.31282 -0.1372
#> [4,] -0.2959  0.01061 -0.02913  0.5072  0.250  0.02831 -0.1575 -0.13724  0.3441
#> [5,]  0.1104 -0.16736  2.20957  0.2504  4.623  0.03803 -0.1856 -0.29860  0.1217
#>       [,10]
#> [1,]  0.038
#> [2,] -0.186
#> [3,] -0.299
#> [4,]  0.122
#> [5,]  0.738

Since each matrix is stacked rowwise, it is actually not that readable and practical. This function has been defined for the next simulation functions.

Minnesota Prior

BVAR

Consider BVAR Minnesota prior setting,

\[A \sim MN(A_0, \Omega_0, \Sigma_e)\]

\[\Sigma_e \sim IW(S_0, \alpha_0)\]

From Litterman (1986) and Bańbura et al. (2010)
Each \(A_0, \Omega_0, S_0, \alpha_0\) is defined by adding dummy observations
- build_xdummy()
- build_ydummy()
sigma: Vector \(\sigma_1, \ldots, \sigma_m\)
- \(\Sigma_e = diag(\sigma_1^2, \ldots, \sigma_m^2)\)
- \(\sigma_i^2 / \sigma_j^2\): different scale and variability of the data
lambda
- Controls the overall tightness of the prior distribution around the RW or WN
- Governs the relative importance of the prior beliefs w.r.t. the information contained in the data
  - If \(\lambda = 0\), then posterior = prior and the data do not influence the estimates.
  - If \(\lambda = \infty\), then posterior expectations = OLS.
- Choose in relation to the size of the system (Bańbura et al. (2010))
  - As m increases, \(\lambda\) should be smaller to avoid overfitting (De Mol et al. (2008))
delta: Persistence
- Litterman (1986) originally sets high persistence \(\delta_i = 1\)
- For Non-stationary variables: random walk prior \(\delta_i = 1\)
- For stationary variables: white noise prior \(\delta_i = 0\)
eps: Very small number to make matrix invertible

bvar_lag <- 5
(spec_to_sim <- set_bvar(
  sigma = c(3.25, 11.1, 2.2, 6.8), # sigma vector
  lambda = .2, # lambda
  delta = rep(1, 4), # 4-dim delta vector
  eps = 1e-04 # very small number
))
#> Model Specification for BVAR
#> 
#> Parameters: Coefficent matrice and Covariance matrix
#> Prior: Minnesota
#> # Type '?bvar_minnesota' in the console for some help.
#> ========================================================
#> 
#> Setting for 'sigma':
#> [1]   3.25  11.10   2.20   6.80
#> 
#> Setting for 'lambda':
#> [1]  0.2
#> 
#> Setting for 'delta':
#> [1]  1  1  1  1
#> 
#> Setting for 'eps':
#> [1]  1e-04

sim_mncoef(p, bayes_spec, full = TRUE) can generate both \(A\) and \(\Sigma\) matrices.
In bayes_spec, only set_bvar() works.
If full = FALSE, \(\Sigma\) is not random. It is same as diag(sigma) from the bayes_spec.
full = TRUE is the default.

(sim_mncoef(bvar_lag, spec_to_sim))
#> $coefficients
#>            [,1]     [,2]      [,3]     [,4]
#>  [1,]  0.196255 -0.19084 -0.026363  0.07600
#>  [2,] -0.017158  0.27873 -0.033404 -0.00816
#>  [3,] -0.225869  0.31165  0.127705 -0.48148
#>  [4,] -0.002706 -0.26068  0.038566  0.04105
#>  [5,] -0.023064 -0.11717 -0.030881  0.35496
#>  [6,]  0.000253 -0.05523 -0.018351  0.03580
#>  [7,] -0.001364 -0.06563 -0.051615 -0.33330
#>  [8,] -0.025569  0.11938 -0.011720 -0.15248
#>  [9,]  0.062006 -0.07985 -0.001830  0.15440
#> [10,]  0.008609 -0.03437  0.016242  0.00866
#> [11,] -0.069804  0.11559  0.003275  0.30043
#> [12,]  0.001677  0.01789 -0.000582 -0.02224
#> [13,] -0.000864  0.05893  0.038960  0.06788
#> [14,]  0.008799 -0.04337  0.005557  0.02273
#> [15,] -0.053924  0.16120 -0.006870  0.10085
#> [16,] -0.009109  0.00284 -0.008880 -0.01168
#> [17,]  0.008001 -0.05816  0.012293 -0.05256
#> [18,] -0.006107  0.03645 -0.004491 -0.00687
#> [19,] -0.011204  0.05703  0.036651  0.13185
#> [20,] -0.003148 -0.05735  0.016794  0.04086
#> 
#> $covmat
#>        [,1]     [,2]    [,3]   [,4]
#> [1,]  2.615  -5.1044  0.1522   2.44
#> [2,] -5.104  32.2768 -0.0549 -12.04
#> [3,]  0.152  -0.0549  1.4573   0.31
#> [4,]  2.440 -12.0446  0.3101  30.83

BVHAR

sim_mnvhar_coef(bayes_spec, full = TRUE) generates BVHAR model setting:

\[\Phi \mid \Sigma_e \sim MN(M_0, \Omega_0, \Sigma_e)\]

\[\Sigma_e \sim IW(\Psi_0, \nu_0)\]

similar to BVAR, bayes_spec option wants bvharspec. But
- set_bvhar()
- set_weight_bvhar()
There is full = TRUE, too.

BVHAR-S

(bvhar_var_spec <- set_bvhar(
  sigma = c(1.2, 2.3), # sigma vector
  lambda = .2, # lambda
  delta = c(.3, 1), # 2-dim delta vector
  eps = 1e-04 # very small number
))
#> Model Specification for BVHAR
#> 
#> Parameters: Coefficent matrice and Covariance matrix
#> Prior: MN_VAR
#> # Type '?bvhar_minnesota' in the console for some help.
#> ========================================================
#> 
#> Setting for 'sigma':
#> [1]  1.2  2.3
#> 
#> Setting for 'lambda':
#> [1]  0.2
#> 
#> Setting for 'delta':
#> [1]  0.3  1.0
#> 
#> Setting for 'eps':
#> [1]  1e-04

(sim_mnvhar_coef(bvhar_var_spec))
#> $coefficients
#>          [,1]    [,2]
#> [1,]  0.06349 -0.2954
#> [2,] -0.01177  0.4277
#> [3,]  0.05543 -0.2689
#> [4,]  0.00235  0.0574
#> [5,] -0.05121  0.1928
#> [6,] -0.00891  0.0245
#> 
#> $covmat
#>        [,1]   [,2]
#> [1,]  0.322 -0.407
#> [2,] -0.407  3.195

BVHAR-L

(bvhar_vhar_spec <- set_weight_bvhar(
  sigma = c(1.2, 2.3), # sigma vector
  lambda = .2, # lambda
  eps = 1e-04, # very small number
  daily = c(.5, 1), # 2-dim daily weight vector
  weekly = c(.2, .3), # 2-dim weekly weight vector
  monthly = c(.1, .1) # 2-dim monthly weight vector
))
#> Model Specification for BVHAR
#> 
#> Parameters: Coefficent matrice and Covariance matrix
#> Prior: MN_VHAR
#> # Type '?bvhar_minnesota' in the console for some help.
#> ========================================================
#> 
#> Setting for 'sigma':
#> [1]  1.2  2.3
#> 
#> Setting for 'lambda':
#> [1]  0.2
#> 
#> Setting for 'eps':
#> [1]  1e-04
#> 
#> Setting for 'daily':
#> [1]  0.5  1.0
#> 
#> Setting for 'weekly':
#> [1]  0.2  0.3
#> 
#> Setting for 'monthly':
#> [1]  0.1  0.1

(sim_mnvhar_coef(bvhar_vhar_spec))
#> $coefficients
#>          [,1]    [,2]
#> [1,] -0.01953  0.2058
#> [2,]  0.14585 -0.2401
#> [3,]  0.01613 -0.0491
#> [4,]  0.05043 -0.1112
#> [5,]  0.02742  0.3732
#> [6,] -0.00569 -0.0749
#> 
#> $covmat
#>       [,1]   [,2]
#> [1,] 0.566  0.263
#> [2,] 0.263 14.995