Consider a system process given by $$x_t = -0.9x_{t-2} + w_t, \hspace{2mm} t=1,\dots n$$ where $x_0 \sim N(0,\sigma^2_0), x_{-1} \sim N(0,\sigma_1^2)$, and $w_t$ is Gaussian white noise with variance $\sigma^2_w$.
My question is how can I simulate this series and include $x_0$ and $x_{-1}$? I know how to use the arima.sim function to simulate general AR processes, but I don't know how to include the initial conditions $x_0$ and $x_{-1}$.
Programming the simulation from scratch requires just a few lines of $\textsf{R}$ code:
n <- 500
sigma_sq_minus1 <- 3
sigma_sq_0 <- 2
sigma_sq_w <- 1
x_minus1 <- rnorm(1, sd = sqrt(sigma_sq_minus1))
x_0 <- rnorm(1, sd = sqrt(sigma_sq_0))
w <- rnorm(n, sd = sqrt(sigma_sq_w))
x <- c(x_minus1, x_0, rep(NA, n))
for (t in seq_len(n)+2) {
x[t] <- -0.9 * x[t-2] + w[t-2]
}
plot(x, type = "l")
arima.sim to do it but it's not really its purpose: it assumes that you always want the long-run distribution, so it forces an initial burn-in period that is discarded. You can put in extra zeros for $x_{-2}$ and $x_{-3}$ so that $x_{-1} = w_{-1}$ and $x_0 = w_0$, and then you can supply the initial values as initial innovations, like this: arima.sim(model = list(ar = c(0, -0.9)), innov = c(x_minus1, x_0, w), start.innov = c(0, 0), n = n+2, n.start = 2). I wouldn't recommend it; what you have here is a lot easier to read, and has similar performance.
$\endgroup$
Commented
Sep 23, 2022 at 0:08
arima.sim. Tbh, I'm not really familiar with $\textsf{R}$'s ARIMA functionality. Since the OP explicitly mentions the arima.sim function your comment (which I encourage you to repost as an answer) might actually be a better answer than mine.
$\endgroup$
Commented
Sep 23, 2022 at 10:10