Dear R-users,
let's assume a simple linear model in which to each original observation one
attaches random response and covariate, but with weights that are practically
zero. If one runs a linear model to both original and augmented data, one
obtains the same estimates, but obviously different standard errors.
Could we reproduce the lm() outcomes in a mixed effect model setting using, for
instance, lme(nlme)? Is there any way to get the lm() estimates and standard
deviations by accounting for correlation within each new pseudo-observation?
Below a reproducible and commented instance of the problem .
Thanks in advance for your help,
GC
library(nlme)
## simulating data
m <- 1000
x <- runif(m)
y <- 3 + 4*x + rnorm(m, sd=0.5)
id <- 1:m
dati <- cbind(id=id, x=x, y=y, w=1)
## new data: repeat n times original data
n <- 10
datiA <- kronecker(dati, matrix(1,n,1))
## trasform in dataframes with suitable names
dati <- as.data.frame(dati)
datiA <- as.data.frame(datiA)
names(datiA) <- names(dati)
## assign ~0 weight to everyone, but the original observations
fw <- 10^-30
datiA$w <- fw
datiA$w[seq(1,m*n,n)] <- 1 - n*fw
## assign random numbers to all response/covariate, but for the original
observations
datiA$x[datiA$w==fw] <- runif(m*n-m)
datiA$y[datiA$w==fw] <- runif(m*n-m)
## fits with lm()
fit <- lm(y~x, weights=w, data=dati)
summary(fit)$coef[,1:2]
fitA <- lm(y~x, weights=w, data=datiA)
summary(fitA)$coef[,1:2]
## mixed model approach, lme()
## attempting to reproduce results from "fit"
fitMM <- lme(y ~ x, weights=varFixed(~1/w),
random=~1|id,
data=datiA)
summary(fitMM)$tTable[,1:2]
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
