How to generate a VAR(1) model? [closed]

Question

I already written the code but something went wrong and I don't know why...

here is the code.

    > set.seed(777)
    > 
    > n = 200
    > k = 2
    > p = 1
    > 
    > phi = matrix(c(.2, -.6, .3, 1.1), k)
    > 
    > rho = 0.5*sqrt(2)
    > mu1 = 0
    > s1 = sqrt(2)
    > mu2 = 0
    > s2 = 1
    > 
    > mu = c(mu1,mu2)
    > sigma = matrix(c(s1^2, s1*s2*rho, s1*s2*rho, s2^2), 2)
    > 
    > library(MASS)
    > at = mvrnorm(n, mu = mu, Sigma = sigma )
    > 
    > var1 = matrix(0, k, n + 2 * p) 
    > for (i in (p + 1):(n + 2*p)){ 
    +   var1[, i] = phi %*% var1[, i-1] + as.vector(at)
    + }
    Error in phi %*% var1[, i - 1] + as.vector(at) : 
      dims [product 2] do not match the length of object [400]
    > 
    > var1 = ts(t(var1[, -(1:p)]))

seems like the length didn't match. But considered I know nothing about VAR model, I don't know how to fix it. I got the code about VAR online and because my homework **require a bivariate normal distribution as a(t)**, I made a `at` like that and replace the `rnorm` in the example code, then the error pop up.

can anyone help me with this? thank you!

update---------------------------------------------------------

now I change the code to below:

> set.seed(777)
> 
> n = 200
> k = 2
> p = 1
> 
> phi = matrix(c(.2, -.6, .3, 1.1), k)
> 
> n_at = 202
> rho = 0.5*sqrt(2)
> mu1 = 0
> s1 = sqrt(2)
> mu2 = 0
> s2 = 1
> 
> mu = c(mu1,mu2)
> sigma = matrix(c(s1^2, s1*s2*rho, s1*s2*rho, s2^2), 2)
> 
> library(MASS)
> at = mvrnorm(n_at, mu = mu, Sigma = sigma )
> 
> var1 = matrix(0, k, n + 2 * p) 
> for (i in (p + 1):(n + 2*p)){ 
+   var1[, i] = phi %*% var1[, i-1] + at[, i]
+ }
Error in phi %*% var1[, i - 1] + at[, i] : 
  dims [product 2] do not match the length of object [202]
> 
> var1 = ts(t(var1[, -(1:p)]))
> 
> plot.ts(var1)

still, there is error, but if I change at[, i] to as.vector(at)[i], there will not be any error but by the plot.ts, 2 series are exactly the same...

how to fix this?

Answer 1

I found out by myself

here is the simplified version:

> set.seed(777)
> 
> phi = matrix(c(.2, -.6, .3, 1.1), 2)
> 
> library(MASS)
> at = mvrnorm(201, mu = c(0,0), Sigma = matrix(c(2, 1, 1, 1), 2))
> 
> var1 = matrix(0, 2, 201) 
> for (i in 2:201){ 
+   var1[, i] = phi %*% var1[, i-1] + at[i,]
+ }
> 
> var1 = ts(t(var1[, -1]))

It should be at[i,] but not at[,i] if I want to iterate it by row.