#Examples Based on Chapter 2 material
mort = scan("cmort.dat")
temp = scan("temp.dat")
part = scan("part.dat")
par(mfrow=c(3,1))
ts.plot(mort)
ts.plot(temp)
ts.plot(part)
pairs(cbind(mort, temp, part)) 
temp = temp - mean(temp)
temp2 = temp^2
temp3 = temp^3
t = 1:length(mort)
t = 1:length(mort)
summary(fit1)
(AIC(fit1)/508) -log(2*pi)
fit2 = lm(mort~t + temp + temp2 + part)
# Results
AIC(fit2)/508 -log(2*pi)
# R gives n*AIC
fit3=lm(mort ~ t + temp + temp2 + temp3 + part)
(AIC(fit3)/508) -log(2*pi)
temp4=temp^4
fit4=lm(mort ~ t + temp + temp2 + temp3 + temp4 + part)
(AIC(fit4)/508) -log(2*pi) 
# Differencing global temp series
gtemp=scan("globtemp.dat")
x=gtemp[45:142]
t=1900:1997
par(mfrow=c(2,1))
plot(t[-1],diff(x), type="l", xlab="year", ylab="first diff")
acf(diff(x),lag=20)
plot(t[-(1:2)],diff(x,lag=2),type="l",xlab="year",ylab="2nd diff")
acf(diff(x,lag=2),lag=20)
# Moving Averages mortality data
mort=scan("cmort.dat")
ma1=filter(mort,sides=2,rep(1,53)/53)
ma2=filter(mort,sides=2,rep(1,11)/11)
plot(mort,pch="o",xlab="time",ylab=" ")
lines(ma1,col=5)
lines(ma2,col=6)
# lowess plot
l1=lowess(mort,f=0.1)
l2=lowess(mort,f=0.3)
l3=lowess(mort,f=0.6)
plot(mort,pch="o",xlab="time",ylab=" ")
lines(l1$y,col="blue")
lines(l2$y,col="green")
lines(l3$y,col="black")
# supsmu plot
t=seq(1:length(mort))
s1=supsmu(t,mort,span=0.01)
s2=supsmu(t,mort,span=0.05)
plot(mort,pch="o",xlab="time",ylab=" ")
lines(s1$y,col=3)
lines(s2$y,col=4)
#
t=1:500
xt=2*cos(2*pi*(1/50)*t + .6*pi ) + rnorm(500,0,2)
plot(xt,type='l',col=2)
z1t=cos(2*pi*(1/50)*t); z2t=sin(2*pi*(1/50)*t)
fit=lm(xt ~ z1t + z2t -1 )
lines(fit$fitted.values,col=1)
#  Periodogram
P=(4/500)*abs(fft(xt)/sqrt(500))^2
# Values of j/n
f=0:250/500
plot(f,P[1:251],type="l",xlab="frequency",ylab=" ")
#
J=1:249
beta1=rep(NA,249)
beta2=rep(NA,249)
for(j in 1:249)
{
beta1[j]=(2/500)*(sum(xt*cos(2*pi*t*(j/500))))
beta2[j]=(2/500)*(sum(xt*sin(2*pi*t*(j/500))))
}
P=beta1^2 + beta2^2
plot(J/500,P,type="l",xlab="frequency",ylab=" ")