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@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. ______________________________________________ R-help@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.
