# Commands for the examples from the lectures in week 34
#=======================================================

# Birthweight example
#---------------------

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

# Plot data
plot(firstborn$weeks,firstborn$weight,col=c(rep('blue',28),rep('red',28)),
     pch=c(rep(15,28),rep(16,28)),xlab="Weeks",ylab="Weight")

# Fit linear regression model and plot regression lines:
fit=lm(weight~weeks+factor(sex),data=firstborn)
summary(fit)
abline(fit$coef[1],fit$coef[2],col='blue',lwd=2)
abline(fit$coef[1]+fit$coef[3],fit$coef[2],col='red',lwd=2)

#=====================================

# Beetle example
#---------------

# Read data:
data="http://www.stat.ufl.edu/~aa/glm/data/Beetles2.dat"
beetle=read.table(data,header=T)

# Plot data:
plot(beetle$logdose,beetle$dead/beetle$n,pch=16,col='blue',
     ylim=c(0,1),xlab="Dose",ylab="Proportion dead")

# Fit logistic regression model and plot fitted probabilites:
fit=glm(cbind(dead,n-dead)~logdose,family=binomial,data=beetle)
summary(fit)
x=seq(1.65,1.95,0.01)
fitx=exp(fit$coef[1]+fit$coef[2]*x)/(1+exp(fit$coef[1]+fit$coef[2]*x))
lines(x,fitx,lwd=2,col='red')

#==========================================

# Example: number of children
#----------------------------

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

# Plot data (and add a small amount of noise)
plot(jitter(birth$age,0.25),jitter(birth$children,0.25),pch=1,col='blue',
     ylim=c(0,7),xlab="Age",ylab="Number of children")

# Fit Poisson regression model and plot fitted values:
fit=glm(children~age,family=poisson,data=birth)
summary(fit)
x=seq(17,43,1)
fitx=exp(fit$coef[1]+fit$coef[2]*x)
lines(x,fitx,lwd=2,col='red')