R-sig-mixed-models Digest, Vol 62, Issue 29

9sbnopoi
-----Original Message-----
From:	r-sig-mixed-models-request at r-project.org
Sender:	r-sig-mixed-models-bounces at r-project.org
Date:	Mon, 13 Feb 2012 20:07:11 
To: <r-sig-mixed-models at r-project.org>
Reply-To: r-sig-mixed-models at r-project.org
Subject: R-sig-mixed-models Digest, Vol 62, Issue 29

Send R-sig-mixed-models mailing list submissions to
	r-sig-mixed-models at r-project.org

To subscribe or unsubscribe via the World Wide Web, visit
	https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
or, via email, send a message with subject or body 'help' to
	r-sig-mixed-models-request at r-project.org

You can reach the person managing the list at
	r-sig-mixed-models-owner at r-project.org

When replying, please edit your Subject line so it is more specific
than "Re: Contents of R-sig-mixed-models digest..."

Today's Topics:

   1. Comparing against a negative control in an LMM (Masca, Nick)
   2. MCMCglmm with cross-classified random effects (Agostino Moro)
   3. Re: MCMCglmm with cross-classified random effects
      (Jarrod Hadfield)
   4. Re: Considerable discrepancies between fixed and random
      effect estimates of lme4 (glmer) and glmmADMB (Ben Bolker)
   5. Any package for best subset selection for random effects
      model (Tao Zhang)
   6. Interpretation of nonlinear mixed-effects modeling	results
      (Gang Chen)
   7. Re: Considerable discrepancies between fixed and random
      effect estimates of lme4 (glmer) and glmmADMB (Adam Smith)

----------------------------------------------------------------------

Message: 1
Date: Mon, 13 Feb 2012 12:33:42 +0000
From: "Masca, Nick" <Nick.Masca at effem.com>
To: "r-sig-mixed-models at r-project.org"
	<r-sig-mixed-models at r-project.org>
Subject: [R-sig-ME] Comparing against a negative control in an LMM
Message-ID:
	<8295A4D50D4C644CAC4323DD070D9597106F2D at 034-CH1MPN1-014.034d.mgd.msft.net>

Content-Type: text/plain

Hi all,

I have a problem based on a colleague's experiment that I've been asked to analyse, which is more of a general mixed modelling issue rather than specifically an R issue, and I would be extremely grateful for any help that any readers of this list can provide.

An experiment was conducted in which the aim was to compare 3 concentrations of 2 active treatments (i.e. 6 active treatments in total) to a negative control.  Three batches of each of the actives have been tested, and 3 reps tested for each batch.  In contrast, 20 replicates have been taken of the negative control - but, by definition, there is no "batch" for this treatment.

Here is some code to reproduce the experimental design:

Treat<- factor(c(rep("NC", 20), rep("A", 27), rep("B", 27)))
Conc<-factor(c(rep(1, 20), rep(1:3, each=9), rep(1:3, each=9)))
Batch<-factor(c(rep(1, 20), rep( rep(1:3, each=3), 6)))
Treatment<-factor(Treat:Conc)  #specify new treatment variable (so don't attempt to estimate Conc. 2&3 for NC)

I originally planned to analyses these data in a LMM, with Treat*Conc as a 7 level fixed effect (i.e. 3*2 actives + control), and with Treat:Conc:Batch as random.  The following code simulates my response variable assuming this model:

                Resp<-  rep(9, 74) + #simulate intercept
                                c( rep(rnorm(1, 0, sd=2.5), 20)^2, rep(rnorm(18, 0, sd=2.5), each=3)^2) + #simulate treat.conc.batch variance
                                rep(rnorm(74, 0, sd=.2)^2) + #simulate residual variance
                                c(rep(0,20), rep(c(-4, 0,0,-4, 0,0), each= 9)) #simulate fixed effects
                Data<-data.frame(Treatment, Conc, Batch, Resp)

While this code models the data using lme4:
                Mod<-lmer(Resp ~ Treatment + (1|Treatment:Batch), data=Data)

I can now obtain and plot treatment means/CIs using glht in the multcomp package:
library(multcomp)
                Mean.mat<-diag(rep(1,7))
                                Mean.mat[,1]<-rep(1,7)
                                rownames(Mean.mat)<-levels(Data$Treatment)
                Est.means<-glht(Mod, Mean.mat)
plot(Est.means)

Hopefully from the above plot you can see what my issue is.  The negative control, which I want to compare everything against, has by far the least precision around its estimate, despite the data for the control hardly varying at all.  This happens because the greatest source of variability in the model (by far) is the variability between batches, but different batches of the negative control don't exist.  As such, I'm not sure that this is a fair way to model the data, because the negative control is unfairly penalised by the variability between the batches of the other treatments.

I imagine that this kind of problem isn't particularly uncommon, but it's the first time I've had to deal with something like this myself.   The only potential solution I've come up with so far is to scrap the negative control from the model, and simply subtract the negative control's mean "count" from all other values (either by specifying this mean as an offset or by subtracting it from all data-points).  But this will probably give "anti-conservative" results, as it would assume the mean for the negative control doesn't vary.

I would be extremely grateful if anyone would care to share their thoughts on possible solutions to this problem - and whether anyone has dealt with this kind of issue before.  I feel that I may well be missing something obvious - but can't see at the moment how else to get around it!

Many thanks for any help you can provide.

Cheers,

Nick

------------------------------

Message: 2
Date: Mon, 13 Feb 2012 15:02:20 +0000
From: Agostino Moro <agostino.moro99 at gmail.com>
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] MCMCglmm with cross-classified random effects
Message-ID:
	<CAMS_pxvdSZVhe_qSFgqsnkMofyTGxL6eLAY4g114VzS=k9HFpA at mail.gmail.com>
Content-Type: text/plain; charset=ISO-8859-1

Dear R-users,

I would like to fit ?a glmm with cross-classified random effects with
the function MCMCglmm. Something along the lines:

model1<-MCMCglmm(response~pred1, random=~re1+re2, data=data)

where re1 and re2 should be crossed random effects. I was wondering
whether you could tell me specifying cross-classified random effects
in MCMCglmm requires a particular syntax? Are there any examples
somewhere? I have had a look at the manual and the package vignette,
but I have not been able to find any examples relevant to what I want
to do.

Thanks,

Agostino

------------------------------

Message: 3
Date: Mon, 13 Feb 2012 15:19:07 +0000
From: Jarrod Hadfield <j.hadfield at ed.ac.uk>
To: Agostino Moro <agostino.moro99 at gmail.com>
Cc: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] MCMCglmm with cross-classified random effects
Message-ID: <20120213151907.57957fqy2tmlxcg8 at www.staffmail.ed.ac.uk>
Content-Type: text/plain; charset=ISO-8859-1; DelSp="Yes";
	format="flowed"

Hi,

As long as the levels of re1 and re2 are uniquely labelled any cross  
classification will be dealt with appropriately.

Cheers,

Jarrod

Quoting Agostino Moro <agostino.moro99 at gmail.com> on Mon, 13 Feb 2012  
15:02:20 +0000: