# Example R code for STAT 520 HW 2



#2

library(TSA); data(wages); 

#(a)
plot(wages, type='o')

#(b)

wage.lin.fit <- lm(wages ~ time(wages))
summary(wage.lin.fit)

plot(wages, type='o'); abline(wage.lin.fit)

#(c)
plot(rstudent(wage.lin.fit),type='o')

#(d)
time.sq = (time(wages))^2
wage.quad.fit <- lm(wages ~ time(wages) + time.sq)
summary(wage.quad.fit)

plot(wages, type='o'); lines(as.vector(time(wages)), fitted(wage.quad.fit))

#(e)
plot(rstudent(wage.quad.fit),type='o')

#(f)
runs(rstudent(wage.quad.fit))

#(g)
acf(rstudent(wage.quad.fit))

#(h)
hist(rstudent(wage.quad.fit))
qqnorm(rstudent(wage.quad.fit))
shapiro.test(rstudent(wage.quad.fit))

#3
library(TSA); data(prescrip)

#(a)
plot(y=prescrip,x=as.vector(time(prescrip)),type='l',xlab='year',ylab='costs')
points(y=prescrip,x=as.vector(time(prescrip)),pch=as.vector(season(prescrip)))

#(b)
har.=harmonic(prescrip,m=1) 

model.har.pre=lm(prescrip~har.+time(prescrip))
summary(model.har.pre)


plot(as.vector(time(prescrip)),fitted(model.har.pre),ylab='Costs',type='l',ylim=range(c(fitted(model.har.pre),prescrip))) 
# the ylim option ensures that the 
# y axis has a range that fits the raw data and the fitted values
points(prescrip,pch=as.vector(season(prescrip)) )

#(c)

plot(rstudent(model.har.pre), type='o')

#(d)

runs(rstudent(model.har.pre))

#(e)

acf(rstudent(model.har.pre))

#(f)
hist(rstudent(model.har.pre))
qqnorm(rstudent(model.har.pre))
shapiro.test(rstudent(model.har.pre))