# Some Code Examples from Day 3
# First, and example of housekeeping in S-Plus
# Then a calculation of rho hat, with some discussion of
# the "right" formula.

########################################################
#Housekeeping

ls()

remove(ls())
######################################################

#Generating a stationary AR(1) model sigma=1, rho=.3, T=100

y<-rep(0,100)

y[1]<- rnorm(n=1, sd=1/sqrt(1-.09))

for (i in 2:100) { y[i]= .3*y[i-1] +rnorm(1)}

####################################################
#calculate rho hat
#We'll use the method of moments estimator here
###################################################

C<-sum(y[1:99]*y[2:100])/99

V<-sum(y[1:100]^2)/100

#So rho hat is given by

C/V

#This is the estimated covariance (assuming mean zero) divided by the
#estimated variance (assuming mean zero).

####################################################
#Note: you could just as well use the estimator

sum(y[1:99]*y[2:100])/sum(y[1:100]^2)

# which is one percent different
# or used 

sum(y[1:99]*y[2:100])/sum(y[2:100]^2)

# which is also about one percent different
###################################################

# We'll discuss the various T versus T-1 issues as the class progresses.
# The difference is not meaningful in any real-world problem, and the
# AR(1) model is about the only case where we really get to write down
# a formula.