# Punkt b
# Leser inn dataene:
Hawaii=read.table("http://www.uio.no/studier/emner/matnat/math/STK3100/data/Hawaii.txt",header=TRUE)

# Knytter til nlme-bibiloteket:
library(nlme)

# Utfører kommandoene gitt i oppgaven:
Hawaii$Birds=sqrt(Hawaii$Moorhen.Kauai)
M0=gls(Birds~Rainfall+Year,na.action=na.omit,data=Hawaii)
M1=gls(Birds~Rainfall+Year,na.action=na.omit,correlation=corAR1(form=~Year),data=Hawaii)
Hawaii2=Hawaii[!is.na(Hawaii$Birds),]
M2=gls(Birds~Rainfall+Year,data=Hawaii2,correlation=corExp(form=~Year))

# Estimert phi for modell M1 er 0.7734303. Estimert d fra modell M2 er 3.892266. Sjekker at vi har sammenhengen gitt i oppgaveteksten:
-1/log(0.7734303)


# Punkt c

Spruce=read.table("http://www.uio.no/studier/emner/matnat/math/STK3100/data/Spruce.txt", header=T)

# Plukker ut et enkelt tre slik det er gitt i oppgaven:
Spruce1=Spruce[Spruce$Tree=="O1T18",]

# Plotter "logSize" mot "days":
plot(Spruce1$days,Spruce1$logSize)


# Punkt d

# Modell med uavhengige støyledd
M0.1=gls(logSize~days,data=Spruce1)
summary(M0.1)

# Modell med eksponensiell korrelasjonsfunksjon:
M1.1=gls(logSize~days,data=Spruce1,correlation=corExp(form=~days))
summary(M1.1)

# Sammenligner modellene med AIC
AIC(M0.1,M1.1)

# Punkt e

# Ser så på alle trærne og utfører kommandoene gitt i oppgaven (men der vi lar plot være en faktor, ikke en numerisk kovariat)

# Modell med uavhengige støyledd
M0=gls(logSize~days+factor(plot),data=Spruce)
summary(M0)

# Modell med eksponensiell korrelasjonsfunksjon:
cexp=corExp(form=~days|Tree,fixed=FALSE)
M1=gls(logSize~days+factor(plot),Spruce,correlation=cexp)
summary(M1)

# Sammenligner modellene med AIC
AIC(M0,M1)



# Ser først på modellen med uavhengige støyledd. 
# Tilpasser modeller med ulik kompleksitet (bruker ML som estimeringsmetode)

M0.dp=gls(logSize~days+factor(plot),data=Spruce, method="ML")
M0.d=gls(logSize~days,data=Spruce, method="ML")
M0.0=gls(logSize~1,data=Spruce, method="ML")

# Tester modellene
anova(M0.0,M0.d,M0.dp,test=T)


# Gjør så tilsvarende for modellen med eksponensiell korrelasjonsfunksjon:

M1.dp=gls(logSize~days+factor(plot),data=Spruce,correlation=cexp,method="ML")
M1.d=gls(logSize~days,data=Spruce,correlation=cexp,method="ML")
M1.0=gls(logSize~1,data=Spruce, correlation=cexp, method="ML")

# Tester modellene
anova(M1.0,M1.d,M1.dp,test=T)


# Vi ser at plot er en signifikane effekt når vi bruker modellen med uavhengige feilledd, men ikke når vi bruker modellen der feilledden er korrelerte