[R] linear mixed model question
Daniel Malter
daniel at umd.edu
Mon Sep 7 19:39:53 CEST 2009
The workers=as.factor(workers) codeline in my post dropped below my name. It
should be in the code before the command line for the linear model.
Daniel Malter wrote:
>
> 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-tp25318961p25333981.html
Sent from the R help mailing list archive at Nabble.com.
More information about the R-help
mailing list