ARIMA(2,0,1) with Monte Carlo Simulation in R, always show zero initial values, why?

Question

I have the following R code for using Monte Carlo simulation with ARIMA(2,0,1), but the code is not working, can someone please help out? Many thanks!

I want to simulate 52 weeks ahead for 1000 times.

ar1 = 0.4857
ar2 = 0.0173
ma1 = -0.8054
r0 <-  0.002937432

r.sims <- matrix(rep(0,1000*52),nrow=52, ncol=1000)
e <- rnorm(52*1000, mean =  0.002976873, sd = 0.1056021)

r.sims[1] <- ar1*r0 + e[1] + ma1*e[0]
r.sims[2] <- ar1*r.sims[1] + ar2*r0 + e[2] + ma1*e[1]


for(j in 1:1000){
    for(i in 1:52){

        r.sims[i] <- r.sims[i-1]*ar1 + r.sims[i-2]*ar2 + e[i] + e[i-1]


}
}

Hi everyone, I have modified the R-code, they seem to run, but slight problems with initial values being zeros, can anyone help please?

r.sims <- matrix(rep(0,1000*52),nrow=52, ncol=1000)
e.sims <- matrix( rnorm(52*1000, mean =  0.002976873, sd = 0.1056021), 52, 1000) 

r.sims[1] <- ar1*r0 + e.sims[1] + ma1*e0
r.sims[2] <- ar1*r.sims[1] + ar2*r0 + e.sims[2] + ma1*e.sims[1]

for(j in 1:1000){

              for(i in 3:52)

                r.sims[i,] <- r.sims[i-1,]*ar1 + r.sims[i-2,]*ar2 + e.sims[i,] + e.sims[i-1,]

}


          [,996]       [,997]       [,998]        [,999]      [,1000]
 [1,]  0.000000000  0.000000000  0.000000000  0.0000000000  0.000000000
 [2,]  0.000000000  0.000000000  0.000000000  0.0000000000  0.000000000
 [3,]  0.010016382  0.062227948  0.091050851  0.2800012622  0.066996876
 [4,]  0.176391453  0.209624874 -0.085862618  0.1604297614  0.096992283
 [5,]  0.083494205  0.353051231 -0.292593785 -0.1443630697  0.200975961
 [6,] -0.011723566  0.477519286 -0.415555028 -0.2725228872  0.199274345
 [7,]  0.148164653  0.509141653 -0.202830713 -0.2986390235 -0.013175047
 [8,] -0.022379916  0.431266718 -0.020794588 -0.1296328172 -0.038491656
 [9,] -0.161865835  0.231474893 -0.043960056  0.2651309362  0.053714510
[10,] -0.218850909  0.038706080  0.162595635  0.2791854517 -0.029590991

Answer 1

There were two minor errors in your code.

r.sims[1] <- ar1*r0 + e[1] + ma1*e[1] # I replaced e[0] by e[1]
r.sims[2] <- ar1*r.sims[1] + ar2*r0 + e[2] + ma1*e[1]

for(j in 1:1000){
  for(i in 3:52){ # I replaced 1 by 3
    r.sims[i] <- r.sims[i-1]*ar1 + r.sims[i-2]*ar2 + e[i] + e[i-1]
  }
}

The above code works.