[R] Which Durbin-Watson is correct? (weights involved) - using durbinWatsonTest and dwtest (packages car and lmtest)
Dimitri Liakhovitski
dimitri.liakhovitski at gmail.com
Fri Aug 12 17:43:00 CEST 2011
Hello!
I have a data frame mysample (sorry for a long way of creating it
below - but I need it in this form, and it works). I regress Y onto X1
through X11 - first without weights, then with weights:
regtest1<-lm(Y~., data=mysample[-13]))
regtest2<-lm(Y~., data=mysample[-13]),weights=mysample$weight)
summary(regtest1)
summary(regtest2)
Then I calculate Durbin-Watson for both regressions using 2 different packages:
library(car)
library(lmtest)
durbinWatsonTest(regtest1)[2]
dwtest(regtest1)$stat
durbinWatsonTest(regtest2)[2]
dwtest(regtest2)$stat
When there are no weights, the Durbin-Watson statistic is the same.
But when there are weights, 2 packages give Durbin-Watson different
statistics. Anyone knows why?
Also, it's interesting that both of them are also different from what
SPSS spits out...
Thank you!
Dimitri
############################################
### Run the whole code below to create mysample:
intercor<-0.3 # intercorrelation among all predictors
k<-10 # number of predictors
sigma<-matrix(intercor,nrow=k,ncol=k) # matrix of intercorrelations
among predictors
diag(sigma)<-1
require(mvtnorm)
set.seed(123)
mypop<-as.data.frame(rmvnorm(n=100000, mean=rep(0,k), sigma=sigma,
method="chol"))
names(mypop)<-paste("x",1:k,sep="")
set.seed(123)
mypop$x11<-sample(c(0,1),100000,replace=T)
set.seed(123)
betas<-round(abs(rnorm(k+1)),2) # desired betas
Y<-as.matrix(mypop) %*% betas
mypop<-cbind(mypop, Y)
rSQR<-.5
VARofY<- mean(apply(as.data.frame(mypop$Y),2,function(x){x^2})) -
mean(mypop$Y)^2
mypop$Y<-mypop$Y + rnorm(100000, mean=0, sd=sqrt(VARofY/rSQR-VARofY))
n<-200
set.seed(123)
cases.for.sample<-sample(100000,n,replace=F)
mysample<-mypop[cases.for.sample,]
mysample<-cbind(mysample[k+2],mysample[1:(k+1)]) #dim(sample)
weight<-rep(1:10,20);weight<-weight[order(weight)]
mysample$weight<-weight
--
Dimitri Liakhovitski
marketfusionanalytics.com
More information about the R-help
mailing list