Create AR-Process with a for-loop?-CodePudding

I tried to create a AR-process using a for-loop. Starting values should be 0 (be varied later), AR-process should have e defined constant (=3), a random error and a time trend (=t):

t=1:1000
x=c(rep(0,1000))
b=c(rnorm(1000))
for (i in 1:length(t)){
x[i] = 0.01*t[i] 3  0.4*x[i-1] b
}

When I run this code, somewhow x is NULL. I can not comprehend why. Any ideas on how to fix that?

CodePudding user response：

I don't get NULL; I get "replacement" error. There are 2 issues at code level:

loop counter should start from i = 2, because you are to access x[i - 1] in the loop;
the error is b[i], not b.

t = 1:1000
x = rep(0, 1000)
b = rnorm(1000)
for (i in 2:length(t)){
  x[i] = 3   0.01 * t[i]   0.4 * x[i - 1]   b[i]
}

In addition, at model level, I wonder if you should really do:

t = 1:1000
x = rep(0, 1000)
b = rnorm(1000)
## AR(1)
for (i in 2:length(t)){
  x[i] = 0.4 * x[i - 1]   b[i]
}
## linear trend
x <- 3   0.01 * t   x

In this way, the auto-correlation is not messed up with the linear trend.

I highly encourage you to understand how a loop works. On one hand, this basic syntax is useful. On the other hand, it is the canonical way to express math algorithms. But once you are comfortable with loops, you are encouraged to

Here is the problem. diff is much stronger than you think. It is able to remove piecewise random linear trend (the characteristic of a Brownian motion). When we have a simple, global, deterministic linear trend, diff tends to overwork and introduce negative auto-correlation. The correct way to estimate this lag-1 auto-correlation is:

fit <- lm(x ~ t)
acf(fit$residuals)