Simulation of AR(2)-process with \( \phi_{1} \) = 0.6 and \( \phi_{2} \) = -0.28 (In ggplot2)

Setup

I am working my way through Bernhard Pfaff's Analysis of Integrated and Cointegrated Time Series with R. The book contains some excellent examples. In this RPub I'll show you how to change the code of Example 1.2 to use ggplot2.

You can also see my work for Example 1.1

The Existing Example

As it is the example demonstrates how to simulate an AR(2)-process. Additionally the code calculates the roots of the characteristic function.

set.seed(123456) # I added this part
series <- rnorm(1000)
y.st <- filter(series, filter = c(0.6, -0.28), method='recursive')
ar2.st <- arima(y.st, c(2, 0, 0), include.mean=FALSE, transform.pars=FALSE,
                method="ML")
ar2.st$coef

##    ar1    ar2 
##  0.656 -0.309 

polyroot(c(1, -ar2.st$coef))

## [1] 1.062+1.452i 1.062-1.452i

Mod(polyroot(c(1, -ar2.st$coef)))

## [1] 1.799 1.799

root.comp <- Im(polyroot(c(1, -ar2.st$coef)))
root.real <- Re(polyroot(c(1, -ar2.st$coef)))

# Plotting the roots in a unit circle
x <- seq(-1,1,length=1000)
y1 <- sqrt(1 - x^2)
y2 <- -sqrt(1 - x^2)
plot(c(x, x), c(y1, y2), xlab = 'Real Part', ylab = 'Complex Part',
     type = 'l', main='Unit Circle', ylim=c(-2,2), xlim=c(-2, 2))
abline(h=0)
abline(v=0)
points(Re(polyroot(c(1, -ar2.st$coef))), 
       Im(polyroot(c(1, -ar2.st$coef))), pch=19)
legend( -1.5, -1.5, legend="Root of AR(2)", pch=19)

Converting To ggplot2

The Time Series and Autocorrelation Plots

Pfaff does not actually include the code for plotting the figures. It isn't too hard to figure out how to plot things by looking at Example 1.1 and adapting it to the code in Example 1.2.

You'll end up with a plot like this:

And if you look at my previous RPub, you can similarly adopt it to produce a plot like this:

## Loading required package: ggplot2

## Loading required package: gridExtra

## Loading required package: grid

The Unit Circle Plot

We can use the same y1 and y2 from the original example. So we have,

require(ggplot2)

unitCircle <- ggplot(data.frame(x=x,y1=y1,y2=y2)) + #Using the results from earlier
  geom_hline(aes(yintercept=0)) +
  geom_vline(aes(xintercept=0)) +
  geom_line(aes(x=x,y=y1)) + # top of circle
  geom_line(aes(x=x,y=y2)) + # bottom of circle
  geom_point(aes(x=Real, y=Imaginary), 
             data.frame(Real=Re(polyroot(c(1, -ar2.st$coef))), # Using the same code as above
                        Imaginary=Im(polyroot(c(1, -ar2.st$coef))))) +
                          scale_y_continuous("Imaginary Part", limits=c(-2,2)) +
                          scale_x_continuous("Real Part", limits=c(-2,2))

unitCircle

There you have it! I hope you enjoyed the RPub! Check out my blog at http://www.quantifyingQuality.com