#==========================================================
# Simulation of an AR(1) Model using arima.sim()
# Code is from Z&W (slightly modified)
# Also illustrates par(mfrow(2,2)) and graph labeling
#=======================================================
#The simulation using start.innov burn-in. We simply add "something" (here it is 1)
#to get a series that has a mean differing from zero.
e<-rnorm(100,sd=1)
e.start=rnorm(25, sd=1)
y.ar1<- 1+arima.sim(model=list(ar=0.75), n=100,innov=e, start.innov=e.start)
#Here the values in start.innov are used for the "burn in." For a nice AR(1) model
#a burn-in period of 25 probably suffices, but for more complicated models you may want a
#longer burn-in period.
#By the way, you "need" just as much burn-in to generate a series of length 100 as to
#generate one of length 5000.
# Let's check on the mean to see if it makes sense.
# We may run our simulation several times.
mean(y.ar1)
# If we want to have the simulation begin with the same random seed each time
# we can explicitly set the pseudo-random number generator seed
set.seed(101)
# When we use this seed we get a larger mean than we might have expected, but this
# but due to the reasonably large variance of this series (calculate tau!), the mean
# will bounce a lot from seed to seed.
#Now lets do some pictures --- and put three on just one page:
#
#(1) the simulated series,
#(2) the theoretical ACF
#(3) the empirical ACF for our simulation
par(mfrow=c(2,2))
tsplot(y.ar1, main="Simulated AR(1)")
abline(h=1)
gamma.j<-rep(.75,10)^(0:9)
tsplot(gamma.j, type="h", main="Theoretical ACF for AR(1)", ylab="Autocorreleation", xlab="Lag+1")
temp<-acf(y.ar1, lag.max=10, plot=F)
tsplot(temp$acf,type="h", main="Empirical ACF for Simulated AR(1)",
ylab="Autocorreleation", xlab="Lag+1")
#One line to make it prettier:
abline(h=0)
#How do the series compare? To me the typical run looks pretty reasonable,
#though you may notice that at higher lags the sample ACF tends to "turn up."
#Thus, there is a sampling effect that is not present in the theoretical ACF
#You might want to look at the help screens for the acf function.