Set Options and Load Packages

options(digits = 3)
options(width = 75)
library(forecast)

## Loading required package: zoo
## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## Loading required package: timeDate
## This is forecast 6.1

library(PerformanceAnalytics)

## Loading required package: xts
## 
## Attaching package: 'PerformanceAnalytics'
## 
## The following objects are masked from 'package:timeDate':
## 
##     kurtosis, skewness
## 
## The following object is masked from 'package:graphics':
## 
##     legend

Gaussian White Noise Process

# simulate Gaussian White Noise process
set.seed(123)
y = rnorm(250)
ts.plot(y, main="Gaussian White Noise Process", xlab="time", ylab="y(t)",
        col="blue", lwd=2)
abline(h=0)

GWN Process for Monthly Returns

set.seed(123)
y = rnorm(60, mean=0.01, sd=0.05)
ts.plot(y, main="GWN Process for Monthly Returns",
        xlab="time", ylab="r(t)", col="blue", lwd=2)
abline(h=c(0,0.01,-0.04,0.06), lwd=2, lty=c("solid", "solid", "dotted","dotted"), 
       col=c("black", "green", "red", "red"))

Independent White Noise

## simulate IWN using scaled student's t distribution to have variance 1
set.seed(123)
y = (1/sqrt(3))*rt(250, df=3)
ts.plot(y, main="Independent White Noise Process", xlab="time", ylab="y(t)",
        col="blue", lwd=2)
abline(h=0)

Deterministic Trend + GWN

set.seed(123)
e = rnorm(250)
y.dt = 0.1*seq(1,250) + e
ts.plot(y.dt, lwd=2, col="blue", main="Deterministic Trend + Noise")
abline(a=0, b=0.1)

Random Walk Process

set.seed(321)
e = rnorm(250)
y.rw = cumsum(e)
ts.plot(y.rw, lwd=2, col="blue", main="Random Walk")
abline(h=0)

Simulate MA(1) Process: Loop calculation

First we set the model parameters and run the simulation.

# set simulation parameters
n.obs = 250
mu = 1
theta = 0.9
sigma.e = 1
set.seed(123)
e = rnorm(n.obs, sd = sigma.e)
y = rep(0, n.obs)
# simulate MA(1) using simple loop
y[1] = mu + e[1]
for (i in 2:n.obs) {
  y[i] = mu + e[i] + theta*e[i-1]
}

Next we plot the simulated series.

# plot simulation
ts.plot(y, main="MA(1) Process: mu=1, theta=0.9, sigma.e=1",
        xlab="time", ylab="y(t)", col="blue", lwd=2)
abline(h=c(0,1), col=c("black", "red"), lty = c("solid", "dashed"))

Simulating MA(1) Process

ma1.model = list(ma=0.9)
mu = 1
set.seed(123)
ma1.sim = mu + arima.sim(model=ma1.model, n=250)
ts.plot(ma1.sim, main="MA(1) Process: mu=1, theta=0.9",
        xlab="time", ylab="y(t)", col="blue", lwd=2)
abline(h=c(0,1), col=c("black", "red"))

ACF for MA(1) Process: Theta > 0

ma1.acf = ARMAacf(ar=0, ma=0.9, lag.max=10)
ma1.acf

##     0     1     2     3     4     5     6     7     8     9    10 
## 1.000 0.497 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000

plot(0:10, ma1.acf, type="h", col="blue", lwd=2,
     main="ACF for MA(1): theta=0.9", xlab="lag", ylab="rho(j)")
abline(h=0)

MA(1) Process: Theta < 0

ma1.model = list(ma=-0.75)
mu = 1
set.seed(123)
ma1.sim = mu + arima.sim(model=ma1.model, n=250)
ts.plot(ma1.sim,main="MA(1) Process: mu=1, theta=-0.75",
        xlab="time", ylab="y(t)", col="blue", lwd=2)
abline(h=c(0, 1), col=c("black", "red"))

ACF for MA(1) Process: Theta < 0

ma1.acf = ARMAacf(ar=0, ma=-0.75, lag.max=10)
ma1.acf

##     0     1     2     3     4     5     6     7     8     9    10 
##  1.00 -0.48  0.00  0.00  0.00  0.00  0.00  0.00  0.00  0.00  0.00

plot(0:10, ma1.acf,type="h", col="blue", lwd=2,
main="ACF for MA(1): theta=-0.75", xlab="lag", ylab="rho(j)")
abline(h=0)

AR(1) Process: Phi > 0

ar1.model = list(ar=0.9)
mu = 1
set.seed(123)
ar1.sim = mu + arima.sim(model=ar1.model, n=250)
ts.plot(ar1.sim,main="AR(1) Process: mu=1, phi=0.9",
        xlab="time", ylab="y(t)", col="blue", lwd=2)
abline(h=c(0,1), col=c("black", "red"))

ACF for AR(1) Process: Phi > 0

ar1.acf = ARMAacf(ar=0.9, ma=0, lag.max=10)
ar1.acf

##     0     1     2     3     4     5     6     7     8     9    10 
## 1.000 0.900 0.810 0.729 0.656 0.590 0.531 0.478 0.430 0.387 0.349

plot(0:10, ar1.acf, type="h", col="blue", lwd=2,
     main="ACF for AR(1): phi=0.9", xlab="lag", ylab="rho(j)")
abline(h=0)

AR(1) Process: Phi < 0

ar1.model = list(ar=-0.75)
mu = 1
set.seed(123)
ar1.sim = mu + arima.sim(model=ar1.model, n=250)
ts.plot(ar1.sim,main="AR(1) Process: mu=1, phi=-0.75", col="blue", lwd=2,
        xlab="time",ylab="y(t)")
abline(h=c(0,1), col=c("black", "red"))

ACF for AR(1) Process: Phi < 0

ar1.acf = ARMAacf(ar=-0.75, ma=0, lag.max=10)
ar1.acf

##       0       1       2       3       4       5       6       7       8 
##  1.0000 -0.7500  0.5625 -0.4219  0.3164 -0.2373  0.1780 -0.1335  0.1001 
##       9      10 
## -0.0751  0.0563

plot(0:10, ar1.acf,type="h", col="blue", lwd=2,
     main="ACF for AR(1): phi=-0.75", xlab="lag", ylab="rho(j)")
abline(h=0)

AR(1) Process: Phi = 1 (Random Walk)

set.seed(111)
rw.sim = cumsum(rnorm(250))
ts.plot(rw.sim,main="AR(1) Process: phi=1",
        xlab="time", ylab="y(t)", lwd=2, col="blue")
abline(h=0)

AR(1) Process: Phi > 1 (Explosive)

set.seed(123)
phi = 1.01
e = rnorm(250)
y = rep(0,250)
for (i in 2:250) {
  y[i] = phi*y[i-1] + e[i]
}
ts.plot(y,main="AR(1) Process: phi=1.01",
        xlab="time", ylab="y(t)", lwd=2, col="blue")
abline(h=0)