> mdata <-read.table("ch2001.csv", header=T, sep=",")
> year<-mdata[,1]
> ch  <-mdata[,2]
> ydh <-mdata[,3]
> wh  <-mdata[,4]
> m1 <-lm(ch~ydh+wh)
> logLik(m1)
'log Lik.' -137.7718 (df=4)
> AIC(m1)
[1] 283.5436
> BIC(m1)
[1] 286.3758
> yhat <-fitted.values(m1)
> sres <-rstandard(m1)
 plot(yhat,sres, ylim=c(-3,3))
> abline(2,0)
> abline(-2,0)
> abline(0,0)

> library(car)
> vif(m1)
     ydh       wh 
2.698702 2.698702

> #Breusch-Pagan test for heteroskedasticity 
> ncvTest(m1)
Non-constant Variance Score Test 
Variance formula: ~ fitted.values 
Chisquare = 1.986688, Df = 1, p = 0.15869

># Durbin-Watson test (not Durbin's alternative h)
> dwt(m1)
 lag Autocorrelation D-W Statistic p-value
   1     -0.02119744      2.013588   0.526
 Alternative hypothesis: rho != 0