Simulate Unit Root Plus Drift Process in R

[Amarjit_Chandhial]

Hi,

How do I simulate 200 realizations of the following Unit Root Plus Drift process:

Δy(t) = 0.5 + 0.5Δy(t-1) + 0.2Δy(t−3) + eps(t)

in either base R, library(forecast), library(fable) or library(gratis), and estimate?

It's from W. ENDERS (2014) Applied Econometric Time Series, 4thEdn, CH4, p.217.

Amarjit

[startz]

This isn't very elegant, but more or less works.

T <- 200
deltaY <- rnorm(T+3)
for (t in 4:(T+4)) deltaY[t]=0.5+0.5*deltaY[t-1] + .2*deltaY[t-3] + rnorm(1)
y <- cumsum(deltaY)[4:203]

[Amarjit_Chandhial]

Thanks. I have retrieved/estimated the parameters using base R, library(forecast) and library(fable). The labelling is a bit confusing in R: 0.5 is the constant in the model retrieved from fit_tsbl |> report(), and mean = constant/(1-phi1-phi3), I think that's correct, please let me know if I have any errors.

#===============================
set.seed(1)

#===============================
c1 <- +0.5
phi1 <- +0.5
phi3 <- +0.2
mean1 <- (c1/(1.0 - phi1 - phi3))
mean1

T <- 300
deltaY <- rnorm(T+3)

for (t in 4:(T+4)) deltaY[t] = c1 + phi1*deltaY[t-1] + phi3*deltaY[t-3] + rnorm(1)

y <- cumsum(deltaY)[104:303]

plot.ts(y)

y.ts <- as.ts(y)

#===============================
# base R

regs <- cbind(const=rep(1,200))
fit.b <- arima(y.ts, order=c(3,1,0), fixed=c(NA, 0, NA, NA), transform.pars = FALSE, xreg = apply(regs, 2, cumsum))
fit.b

#===============================
library(forecast)

fit <- Arima(y.ts, order=c(3,1,0), fixed=c(NA, 0, NA, NA), transform.pars = FALSE, include.constant = TRUE)
fit

#===============================
library(fable)

fit_tsbl <- 
  y.ts |>
  as_tsibble() |>
  model(arima = ARIMA( value ~ 1 + pdq(3, 1, 0), fixed=c(NA, 0, NA, NA), transform.pars = FALSE) )

fit_tsbl |> report()  
fit_tsbl[[1]][[1]][["fit"]][["par"]][["estimate"]]

drift <- 0.5519750/(1-0.4565424-0.1854919)
drift

fit_tsbl[[1]][[1]]$fit$par # constant
fit_tsbl[[1]][[1]]$fit$model$coef # drift or mean