Simulating Time Series: AR, MA and ARMA

---

Identification stages of Box-Jenkins methodology

	ACF	PACF
AR(p)	Decays exponentially / Damped sine wave pattern / Both	Significant spike through lag p
MA(q)	Significant spike through lag q	Decays exponentially / Damped sine wave pattern / Both
ARMA(p,q)	Exponential decay	Exponential decay

require(astsa)

Simulating AR model

AR(2) with \(\alpha_1 = 0.5\) and \(\alpha_2 = 0.3\) -

set.seed(-99)
ts.sim <- arima.sim(list(order = c(2,0,0), ar = c(0.5, 0.3)), n = 200)

layout(matrix(c(1,1,2,3), 2, 2, byrow = TRUE))
ts.plot(ts.sim)
acf(ts.sim, 12)
pacf(ts.sim, 12)

Generate 100 observations from the AR(1) model -

x <- arima.sim(model = list(order = c(1, 0, 0), ar = .9), n = 100)

# Plot the generated data  
plot(x)

# Plot the sample P/ACF pair
acf2(x)

##      [,1]  [,2]  [,3] [,4]  [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13]
## ACF  0.89  0.78  0.68 0.59  0.52 0.45 0.40 0.36 0.33  0.32  0.27  0.21  0.12
## PACF 0.89 -0.03 -0.05 0.05 -0.02 0.01 0.01 0.05 0.00  0.11 -0.19 -0.10 -0.12
##      [,14] [,15] [,16] [,17] [,18] [,19] [,20]
## ACF   0.06 -0.02 -0.08 -0.16 -0.18 -0.23 -0.25
## PACF -0.01 -0.12  0.02 -0.16  0.13 -0.18  0.06

# Fit an AR(1) to the data and examine the t-table
sarima(xdata = x, p = 1, d = 0, q = 1)$ttable

## initial  value 0.853568 
## iter   2 value 0.549725
## iter   3 value 0.334858
## iter   4 value 0.085441
## iter   5 value 0.062582
## iter   6 value 0.050668
## iter   7 value 0.050503
## iter   8 value 0.050424
## iter   9 value 0.050404
## iter  10 value 0.050286
## iter  11 value 0.050152
## iter  12 value 0.050149
## iter  13 value 0.050149
## iter  13 value 0.050149
## iter  13 value 0.050149
## final  value 0.050149 
## converged
## initial  value 0.058074 
## iter   2 value 0.057969
## iter   3 value 0.057616
## iter   4 value 0.057437
## iter   5 value 0.057383
## iter   6 value 0.057376
## iter   7 value 0.057373
## iter   8 value 0.057371
## iter   9 value 0.057370
## iter  10 value 0.057370
## iter  11 value 0.057370
## iter  12 value 0.057370
## iter  12 value 0.057370
## iter  12 value 0.057370
## final  value 0.057370 
## converged

##       Estimate     SE t.value p.value
## ar1     0.8807 0.0505 17.4271  0.0000
## ma1     0.0555 0.1059  0.5244  0.6012
## xmean  -0.5205 0.8672 -0.6002  0.5497

Generate 100 observations from the AR(2) model -

x <- arima.sim(model = list(order = c(2, 0, 0), ar = c(1.5, -.75)), n = 200)
# Plot x
plot(x)

# Plot the sample P/ACF of x
acf2(x)

##      [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12]
## ACF  0.86  0.55  0.17 -0.19 -0.44 -0.55 -0.50 -0.35 -0.15  0.06  0.24  0.34
## PACF 0.86 -0.78 -0.11 -0.02  0.00  0.02  0.01 -0.08  0.05  0.08  0.07 -0.07
##      [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
## ACF   0.36  0.31  0.19  0.04 -0.11 -0.25 -0.33 -0.36 -0.30 -0.18 -0.03  0.10
## PACF  0.08 -0.01 -0.05 -0.05 -0.04 -0.09 -0.02 -0.06  0.11  0.01 -0.15 -0.04
##      [,25]
## ACF   0.19
## PACF  0.09

# Fit an AR(2) to the data and examine the t-table
sarima(xdata = x, p = 2, d = 0, q = 0)$ttable

## initial  value 1.194170 
## iter   2 value 1.065437
## iter   3 value 0.628967
## iter   4 value 0.377698
## iter   5 value 0.148322
## iter   6 value 0.114817
## iter   7 value 0.035748
## iter   8 value 0.033386
## iter   9 value 0.033154
## iter  10 value 0.032863
## iter  11 value 0.032716
## iter  12 value 0.032716
## iter  13 value 0.032715
## iter  13 value 0.032715
## iter  13 value 0.032715
## final  value 0.032715 
## converged
## initial  value 0.037357 
## iter   2 value 0.037334
## iter   3 value 0.037317
## iter   4 value 0.037316
## iter   5 value 0.037316
## iter   6 value 0.037316
## iter   6 value 0.037316
## iter   6 value 0.037316
## final  value 0.037316 
## converged

##       Estimate     SE  t.value p.value
## ar1     1.5343 0.0436  35.1623  0.0000
## ar2    -0.7745 0.0435 -17.8185  0.0000
## xmean   0.2707 0.3032   0.8931  0.3729

ADF test

Simulating a I(1) data -

set.seed(-99)
ts.sim <- arima.sim(list(order = c(2,1,0), ar = c(0.5, 0.3)), n = 200)
plot(ts.sim)

Testing stationarity using ADF test -

tseries::adf.test(ts.sim, alternative="stationary")

## 
##  Augmented Dickey-Fuller Test
## 
## data:  ts.sim
## Dickey-Fuller = -2.6241, Lag order = 5, p-value = 0.3147
## alternative hypothesis: stationary

Hence the time series is not stationary.

If we take first difference -

diff.ts.sim <- diff(ts.sim, 1)
tseries::adf.test(diff.ts.sim, alternative="stationary")

## 
##  Augmented Dickey-Fuller Test
## 
## data:  diff.ts.sim
## Dickey-Fuller = -3.4616, Lag order = 5, p-value = 0.04768
## alternative hypothesis: stationary

So after first difference, the time series is stationary.

ts.plot(diff.ts.sim)

Simulating MA model

MA(2) with \(\beta_1 = 0.5\) and \(\beta_2 = 0.3\) -

set.seed(-99)
ts.sim <- arima.sim(list(order = c(0,0,2), ma = c(0.5, 0.3)), n = 200)

layout(matrix(c(1,1,2,3), 2, 2, byrow = TRUE))
ts.plot(ts.sim)
acf(ts.sim, 12)
pacf(ts.sim, 12)

x <- arima.sim(model = list(order = c(0, 0, 1), ma = -.8), n = 100)

# Plot x
plot(x)

# Plot the sample P/ACF of x
acf2(x)

##       [,1]  [,2]  [,3]  [,4]  [,5] [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12]
## ACF  -0.46  0.05 -0.15  0.11 -0.02 0.07 -0.07 -0.01 -0.05  0.15 -0.18  0.14
## PACF -0.46 -0.19 -0.28 -0.13 -0.07 0.03  0.00 -0.04 -0.11  0.07 -0.13  0.00
##      [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20]
## ACF  -0.06 -0.02 -0.08  0.20 -0.22  0.19  -0.2  0.17
## PACF  0.03 -0.08 -0.16  0.09 -0.15  0.05  -0.1  0.00

# Fit an MA(1) to the data and examine the t-table
sarima(x, p = 0, d = 0, q = 1)

## initial  value 0.272440 
## iter   2 value 0.122906
## iter   3 value 0.108026
## iter   4 value 0.101302
## iter   5 value 0.096222
## iter   6 value 0.095735
## iter   7 value 0.095695
## iter   8 value 0.095626
## iter   9 value 0.095614
## iter  10 value 0.095614
## iter  10 value 0.095614
## iter  10 value 0.095614
## final  value 0.095614 
## converged
## initial  value 0.070607 
## iter   2 value 0.066854
## iter   3 value 0.065328
## iter   4 value 0.062169
## iter   5 value 0.054403
## iter   6 value 0.053703
## iter   7 value 0.053655
## iter   8 value 0.052978
## iter   9 value 0.048972
## iter  10 value 0.048675
## iter  11 value 0.045764
## iter  12 value 0.044893
## iter  13 value 0.044125
## iter  14 value 0.044045
## iter  15 value 0.044044
## iter  16 value 0.044044
## iter  17 value 0.044044
## iter  18 value 0.044044
## iter  18 value 0.044044
## iter  18 value 0.044044
## final  value 0.044044 
## converged

## $fit
## 
## Call:
## arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
##     xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
##     optim.control = list(trace = trc, REPORT = 1, reltol = tol))
## 
## Coefficients:
##           ma1   xmean
##       -1.0000  0.0025
## s.e.   0.0321  0.0035
## 
## sigma^2 estimated as 1.043:  log likelihood = -146.3,  aic = 298.6
## 
## $degrees_of_freedom
## [1] 98
## 
## $ttable
##       Estimate     SE  t.value p.value
## ma1    -1.0000 0.0321 -31.1077  0.0000
## xmean   0.0025 0.0035   0.7198  0.4733
## 
## $AIC
## [1] 2.985965
## 
## $AICc
## [1] 2.987202
## 
## $BIC
## [1] 3.06412

Simulating ARMA model

ARMA(1,1) with \(\alpha_1 = 0.5\) and \(\beta_1 = 0.5\) -

set.seed(-99)
ts.sim <- arima.sim(list(order = c(1,0,1), ar = 0.5, ma = 0.5), n = 200)

layout(matrix(c(1,1,2,3), 2, 2, byrow = TRUE))
ts.plot(ts.sim)
acf(ts.sim, 12)
pacf(ts.sim, 12)