[R] linear mixed model question
Daniel Malter
daniel at umd.edu
Mon Sep 7 19:36:49 CEST 2009
Wen, to follow up on Thierry, your workers are nested in machines (since each
worker only works one machine). Consider fitting a nested model. Though,
with few observations, you might run into the same problem. Further, if you
have observation triplets, and you expect systematic differences between
each trial, then you would have to include a trial effect in some way.
Have you plotted the data? This can be very informative.
Finally, you may want to consider a fixed-effects approach. Throw in worker
dummies only (without intercept) and test the linear hypothesis that the sum
of the worker dummy coefficients is equal between two machines. The
following example does that for two machines (if you have n machines you
would have binomial coefficient (n,2) linear hypotheses):
#simulate data
machines=rep(c(0,1),each=9)
workers=rep(c(1:6),each=3)
workers.re=rep(rnorm(6),each=3)
error=rnorm(18,0.3)
output=2*machines+workers.re+error
#code machines and workers as factors/dummies
machines=as.factor(machines)
#run OLS (fixed effects) of output on worker dummies without intercept
fm4=lm(output~-1+workers)
summary(fm4)
#test the linear hypothesis that coefficients on the worker dummies are
equal for both machines
#the equivalent formulation used is: is the sum of the coefficients across
machines equal to zero?
library(car)
linear.hypothesis(fm4,hypothesis.matrix=c(1,1,1,-1,-1,-1),rhs=c(0))
hope that helps,
daniel
workers=as.factor(workers)
Wen Huang-3 wrote:
>
> Hello,
>
> I wanted to fit a linear mixed model to a data that is similar in
> terms of design to the 'Machines' data in 'nlme' package except that
> each worker (with triplicates) only operates one machine. I created a
> subset of observations from 'Machines' data such that it looks the
> same as the data I wanted to fit the model with (see code below).
>
> I fitted a model in which 'Machine' was a fixed effect and 'Worker'
> was random (intercept), which ran perfectly. Then I decided to
> complicate the model a little bit by fitting 'Worker' within
> 'Machine', which was saying variation among workers was nested within
> each machine. The model could be fitted by 'lme', but when I tried to
> get
> confidence intervals by 'intervals(fm2)' it gave me an error:
>
> Error in intervals.lme(fm2) :
> Cannot get confidence intervals on var-cov components: Non-positive
> definite approximate variance-covariance
>
> I am wondering if this is because it is impossible to fit a model like
> 'fm2' or there is some other reasons?
>
> Thanks a lot!
>
> Wen
>
> #################
>
> library(nlme)
> data(Machines)
> new.data = Machines[c(1:6, 25:30, 49:54), ]
> fm1 = lme(score ~ Machine, random = ~1|Worker, data = new.data)
> fm1
> fm2 = lme(score ~ Machine, random = ~Machine-1|Worker, data = new.data)
> fm2
> intervals(fm2)
>
> ______________________________________________
> R-help at r-project.org mailing list
> 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.
>
>
--
View this message in context: http://www.nabble.com/linear-mixed-model-question-tp25318961p25333956.html
Sent from the R help mailing list archive at Nabble.com.
More information about the R-help
mailing list