Simulating the error term with an AR(2) model I'm working on some simulations in R where I have a time series in the format
$$y_t=\mu_t+e_t$$
where $\mu_t$ is a trend function that I have defined as $\mu_t=5*t+4*sin(2t)$ and where $t$ goes from $ \{ 0,4\pi \}$ with $n=100$ datapoints. Generally $e_t \sim N(0, \sigma^2)$ however in this case of these simulations I wish to replace $e_t$ with an AR(2) process. The reasoning behind this is that normally if $e_t \sim N(0, \sigma^2)$ then we assume that each $e_t$ is i.i.d however by limiting the error term to a known process the i.i.d is dropped from our assumptions, which is what I am trying to simulate here and I would like to know how I would go about changing the error term to be AR(2).
I've posted a short R-code snippet that shows the function defined above with a normal $e_t$ term which would be replaced with an AR(2) process.
n_sim <- 100 # number of data points
t_sim <- seq(0,4*pi,,100)
c_norm_sim <- rnorm(n_sim)
y2_sim <-5*t_sim+4*sin(2*t_sim)+2*c_norm_sim # Gaussian/normal error
plot(t_sim, y2_sim, t="l", ylim=range(y1_sim,y2_sim)*c(1,1.2))


Question: How do I simulate a vector of error values from the AR(2) process?


Use the rGARMA function in the ts.extend package
You can simulate an $\text{AR}(2)$ model using the rGARMA function in the ts.extend package (see O'Neill 2021 for details).  This package allows you to simulate from the marginal or conditional distribution from any stationary Gaussian ARMA process.$^\dagger$  To use the function you need to stipulate the parameters of your model (and any conditioning values you want to use).
Below I give an example of some code where you first stipulate the parameters of the model and then simulate $n=12$ vectors of $m=100$ values from the process.  (By default the function uses a process with zero mean and unit error variance, but you can change these values if you want.)  The generated object ERRORS is a $12 \times 100$ matrix showing the simulated time-series vectors.  If you just want to simulate a single time-series vector you can set $n=1$ in the simulation.
#Load the package
library(ts.extend)

#Set the model parameters
AR    <- c(0.8, -0.1)
SIGMA <- 5

#Simulate and plot a vector of 100 values from an AR(2) process
set.seed(1)
ERRORS <- rGARMA(n = 12, m = 100, ar = AR, errorvar = SIGMA)
plot(ERRORS)


The plot generated by the package shows each of the simulated time-series vectors; the points in the background are the points from all $n=12$ vectors.  (If you want to get rid of these then just set background = FALSE when you call the plot.)

$^\dagger$ The "GARMA" in rGARMA refers to the Gaussian Auto-Regressive Moving Average model.
