[R-sig-ME] Doug Bates' new book? (Douglas Bates

Carl Von Ende cvonende at niu.edu
Tue Oct 18 20:01:28 CEST 2011


Hi Leo,
 
I don't know if has been published, yet, but preprints of chapters are available at R-forge.
 
 http://lme4.r-forge.r-project.org/book/
 
Enjoy.
 
Carl

 
******************************************
Carl N. von Ende
Department of Biological Sciences
Northern Illinois University
DeKalb, IL 60115-2861

Office:            (815) 753-7826
Fax:                (815) 753-0461
email:              cvonende at niu.edu
******************************************
>>> <r-sig-mixed-models-request at r-project.org> 10/18/2011 12:35 PM >>>
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. To interpret the interactions when the predictors are
      treatment factors in lme4 (Likan Zhan)
   2. Re: parameter restrictions (dave fournier)
   3. Doug Bates' new book? (lopez toledo)
   4. Re: Doug Bates' new book? (Douglas Bates)


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

Message: 1
Date: Tue, 18 Oct 2011 10:05:52 +1100
From: Likan Zhan <likan.zhan at mq.edu.au>
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] To interpret the interactions when the predictors
aretreatment factors in lme4
Message-ID: <011385BA-E66C-4EAD-B9A6-F6A185321E06 at mq.edu.au>
Content-Type: text/plain


> Dear Dr. Maechler,

> I am writing to ask a clarification when I am trying to use "lme4" package in R.

> If in our experiment, both A and B are three levels treatment factors:
> A:1,2,3;
> B:1,2,3

> And we use the following model to fit our data:
> model=lmer(response~A*B+(1|subject)+(1|trial),data=xx)

> and the following is the result of the fixed-effects:

> Fixed-effects

> (Intercept)    t=100, p<.000
>  A2           t=1.0, p<.9
>  A3           t=1.0, p<.9
>  B2           t=10,  p<.008
>  B3           t=1.0, p<.9
>  A2:B3        t=100, p<.000   <===
>  A3:B3        t=1.8, p<.6
>  A3:B2        t=1.0, p<.7
>  A3:B3        t=1.0, p<.6

> Could you tell me what is the meaning of the significant effect of "A2:B3",
> what is the null hypothesis of it?

> Thank you very much.

> Likan

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

Likan Zhan

PhD Candidate

Macquarie Centre for Cognitive Science (MACCS)

Macquarie University

Sydney NSW 2109, Australia

Ph:  +61 2 9850 4440

Cellphone: +61 (0)425 961 818

Fax: +61 2 9850 6059

Email: likan.zhan at mq.edu.au

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



[[alternative HTML version deleted]]



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

Message: 2
Date: Mon, 17 Oct 2011 17:20:14 -0700
From: dave fournier <davef at otter-rsch.com>
To: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] parameter restrictions
Message-ID: <4E9CC63E.8040908 at otter-rsch.com>
Content-Type: text/plain; charset=ISO-8859-1; format=flowed

I don't think it is very difficult. Something like this should do the job.
It depends on how you want to model the correlation of the random effects.
Perhaps the ADMB list would be the place to discuss it.



DATA_SECTION
   int nfixed
  !! nfixed=7;  // for this example
   init_int nobs   // number of observations
   init_int nrand   // number of observations
   init_vector Y(1,nobs)   // observations
   init_matrix X(1,nobs,1,nfixed)  // fixed effects design matrix
   init_matrix Z(1,nobs,1,nrand)   // random effects design matrix
PARAMETER_SECTION
   init_number a
   init_bounded_vector alpha(2,nfixed,0.0,1.0)
   objective_function_value f
  !! int nr1=nrand*(nrand-1)/2; // number of parameters needed for 
correlation
   init_vector c(1,nr1)   // parameters for choleski decomp
                                        // of the correlation matrix
   init_vector log_sigma(1,nrand)           // log of variances
   random_effects_vector u(1,nrand)
PROCEDURE_SECTION
   // parameterize the b's to satisfy constraints
   dvar_vector b(1,nfixed);
   b(1)=-exp(a);
   dvariable ssum=0.0;
   for (int i=2;i<=nfixed;i++) {
     ssum-=b(i-1);
     b(i)=alpha(i)*ssum;
   }
   // choleski decomp of correlation matrix
   dvar_matrix C(1,nrand,1,nrand);
   C.initialize();
   int ioffset=1;
   for (int i=1;i<=nrand;i++) {
     C(i,i)=1;
     if (i>1)
     {
       C(i)(1,i-1)=c(ioffset,ioffset+i-1).shift(1);
       C(i)/=norm(C(i));
       ioffset+=i-1;
     }
     C(i)*=exp(log_sigma(i));   // turn it into Choleski decomp of 
covariance
                                // matrix
   }
   dvar_vector v=C*u;           // now v has right covariance
   dvar_vector pred_Y=X*b+Z*v;  // predicted Y
   dvar_vector r=Y-pred_Y;
   dvariable w=norm2(r)/nobs;   // mean square error
   f=0.5*nobs*log(w);    // concemntrated likelihood
   f+=0.5*norm2(u);    // log prior  for u



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

Message: 3
Date: Tue, 18 Oct 2011 09:39:38 -0700 (PDT)
From: lopez toledo <llopezt2000 at yahoo.com.mx>
To: "r-sig-mixed-models at r-project.org"
<r-sig-mixed-models at r-project.org>
Subject: [R-sig-ME] Doug Bates' new book?
Message-ID:
<1318955978.8077.YahooMailNeo at web161209.mail.bf1.yahoo.com>
Content-Type: text/plain

Hi guys:
Some months ago, I did see a message from Dr Bates, announcing a book he was writing on mixed models using lme4. Have any of you got an idea whether this book is already published?

Cheers
Leo




________________________________
De: "r-sig-mixed-models-request at r-project.org" <r-sig-mixed-models-request at r-project.org>
Para: r-sig-mixed-models at r-project.org
Enviado: Martes, 18 de octubre, 2011 5:00:02
Asunto: R-sig-mixed-models Digest, Vol 58, Issue 20

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. Re: longitudinal analysis of nested samples (Christopher Eckert)
   2. Re: parameter restrictions (Ben Bolker)
   3. Re: Observer-level random effects & the varianc (Ross Culloch)
   4. to interpret the interactions of treatment factors in    lme4
      (Likan Zhan)
   5. mcmcsamps for binomial models in lmer (Gillian Raab)
   6. Re: longitudinal analysis of nested samples (ONKELINX, Thierry)


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

Message: 1
Date: Mon, 17 Oct 2011 09:37:31 -0400
From: Christopher Eckert <chris.eckert at queensu.ca>
To: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] longitudinal analysis of nested samples
Message-ID: <8DC43B42-85D9-4F10-9F0E-B44C6194F76D at queensu.ca>
Content-Type: text/plain; charset=us-ascii


Chris Eckert
Department of Biology
Queen's University
Kingston, Ontario K7L 3N6 Canada
chris.eckert at queensu.ca


> 
> Date: Mon, 17 Oct 2011 08:20:19 +0000
> From: "ONKELINX, Thierry" <Thierry.ONKELINX at inbo.be>
> To: Giancarlo Sadoti <gcsadoti at yahoo.com>,
>     "r-sig-mixed-models at r-project.org"    <r-sig-mixed-models at r-project.org>
> Subject: Re: [R-sig-ME] longitudinal analysis of nested samples
> Message-ID:
>     <AA818EAD2576BC488B4F623941DA74273B8A26EA at inbomail.inbo.be>
> Content-Type: text/plain; charset="iso-8859-1"
> 
> Dear Giancarlo,
> 
> You will need to have a look at your data. Since you claim that TR has 4 levels but the model output indicates only one (Number of obs: 560, groups: TR, 1).
> 
> The trees are nested within point. Therefore a more likely model is lmer(V ~ YR + (1|PT/TR), data=TREES) or lmer(V ~ YR + (1|PT) + (1|PT:TR), data=TREES) Both models are identical, the second is a bit more verbose but more clear. If the points are nested within transects then you can simply add it to the model: lmer(V ~ YR + (1|Transect) + (1|Transect:PT) + (1|Transect:PT:TR), data=TREES)
> 
> These are models with only a random intercept. You can add random slopes as well. E.g. if you want a random slope along year at the tree level:
> lmer(V ~ YR + (1|Transect) + (1|Transect:PT) + (1 + YR|Transect:PT:TR), data=TREES)
> 
> Best regards,
> 
> Thierry
> 
I am interested in this analysis and have a question for Dr. Onlelinx:
He suggested that for a random slopes and intercepts model, one would use:
> lmer(V ~ YR + (1|Transect) + (1|Transect:PT) + (1 + YR|Transect:PT:TR), data=TREES)
Why is the 4th term (1+ YR|Transect:PT:TR) why not just (YR|Transect:PT:TR)?
Thanks very much for any light you can shed on this.
I'm still coming to terms with lmer syntax.

Chris Eckert
Department of Biology
Queen's University
Kingston, ONtario, Canada



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

Message: 2
Date: Mon, 17 Oct 2011 14:48:17 +0000 (UTC)
From: Ben Bolker <bbolker at gmail.com>
To: r-sig-mixed-models at r-project.org
Subject: Re: [R-sig-ME] parameter restrictions
Message-ID: <loom.20111017T155049-287 at post.gmane.org>
Content-Type: text/plain; charset=utf-8

Shu Ding <dingx099 at ...> writes:

> *is there a simple way in R to add parameter restrictions to linear mixed
> effects models. Specifically I would like to restrict slope coefficients to
> be*
> 
> *b1<0 and 0<b2< -b1 and 0<b3< -b1-b2 and ? and 0<b7< -b1-b2-b3-b4-b5-b6.*
> 
> * *
> 
> *Can I use lmeControl, optim, optimx? And How?*
> 

  I don't think there's an easy way to do this.
  You could use AD Model Builder, which allows box constraints 
(lower < parameter < upper ) on individual parameters; reparameterization
(e.g. make b2' = -b2/b1, constrain it between 0 and 1, then b2 will
be constrained between 0 and -b1); and generalized penalties.  But
there's a pretty steep learning curve.


  Ben Bolker



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

Message: 3
Date: Mon, 17 Oct 2011 08:27:47 +0100
From: Ross Culloch <ross.culloch at durham.ac.uk>
To: r-sig-mixed-models at r-project.org
Cc: highstat at highstat.com
Subject: Re: [R-sig-ME] Observer-level random effects & the varianc
Message-ID: <6136B577-57A8-4C53-AB48-3AA1EA86CBEA at dur.ac.uk>
Content-Type: text/plain


> Subject: Re: [R-sig-ME] Observer-level random effects & the varianc
> Message-ID:
>     <1b2697e9233051b4681cc1534e6e183f10ecf803 at sitemail.hostway.com>
> Content-Type: text/plain
> 
> 

Hi Alain, 

thanks for the reply - much appreciated..

> Dear Ross,
> 
> In our upcoming book "Zero Inflated Models and Generalized Linear
> Mixed Models with R" (sorry for the advertising here)

No doubt, like your other books that one will also be on the shelf in our lab once it is published (although sadly that will be after the thesis is due) :)

> we have a
> section in which we simulate data and apply GLMMs with and without an
> observation level random effect. Just a plain simple linear mixed
> model with a random intercept, and some extra noise added for the
> OLRE.
> 
> Y_ij = alpha + beta * X_ij + a_i + eps_ij
> 
> a_i      ~ N(0, sigma1^2)
> eps_ij ~ N(0, sigma2^2)
> 
> In most scenarios (different values for sigma1 and sigma2) both
> variances are estimated close to what they should be. Occasionally I
> see (real) data sets in which the OLRE consumes all the variation,

Indeed, for a data set that I did not present that is the case - which I must admit, I really don't understand. 

> and
> the variance of the random effect a_i is nearly 0. This is confusing
> as the steps are often along the lines of:
> 
> Y_ij = alpha + beta * X_ij + a_i
> a_i      ~ N(0, sigma1^2)
> 
> Sigma1 bigger than 0....but..he.....there is overdispersion....add
> the OLRE.....and then the sigma1 is close to 0. I guess that is a
> scenario where the a_i are all close to 0 anyway.

That seems to be true for my data sets, the values are very low in the first instance, add the OLRE and then it just = 0

> 
> To fully understand the mechanism of the OLRE I suggest that you
> simulate data according to the model above..use different valeus for
> sigma1 and sigma2...and then estimate them. It helps understanding
> what such models do!

I'll give that a go, thanks for the advice!

> 
> In my opinion it is perhaps better to ask yourself the question why
> there is overdispersion (zero inflation, correlation, outliers,
> missing covariates, large variation in the data, non-linear patterns
> modeled as linear). Instead of fixing the problem with a big canon the
> like observation level REs I prefer to dig deeper and fix the problem
> in the appropriate way.

Following your GLMM book and your paper (of which I don't have to hand to properly reference) among other sources, I have checked for all of these potential concerns and accounted for them all (well, missing covariates are hard to account for, so that may be the only exception). 

I really appreciate that point of view, ('fixing with cannon') and I will have a look at simulating the data - but as a skeptic (or perhaps better put as a field research pretending to play statistician) my main concern is (especially considering your long list of potential 'problems') is that I don't know what this cannon ball is trying to fix. But if the OLRE does affect the result would it be safe to say that there is 'a problem'? 

Cheers,

Ross

> 
> Alain
> www.highstat.com [1]
> 
>     Hello all, I am comparing two GLMMs (please see below), one that
> includes the observer-level random effect (OLRE) and one that does
> not. As I understand it, without the OLRE the variance is absorbed in
> the FINAL_ID_FAC (variance = ca. 0.14). By comparing the two models
> there is clearly a considerable difference in the variance of the
> FINAL_ID_FAC, and more importantly the OLRE, which now explains
> considerably more of the variance (ca. 0.4) as compared to the
> FINAL_ID_FAC (ca. 0.12). Ultimately, I am interested in using the
> random effect (FINAL_ID_FAC) to make inference on the importance of
> considering individual in the model. My concern is that the OLRE is
> accounting for the variation within and between individuals, which is
> what I am trying to (somehow) quantify. Previous researchers have
> published without using OLRE (as this field of ecology is constantly
> evolving this very understandable!). However, my results lead me to
> different conclusions (by using the BLUPs for each ID in a subsequent
> analysis) depending on whether or not I include or exclude the OLRE. I
> have posted about this before (many thanks to David Duffy for his
> reply) and I have read numerous posts from other members that
> constantly post on this forum, perhaps I have missed the point, but I
> am still not entirely comfortable with exactly what variation the OLRE
> is accounting for. Many thanks for your time, and any advice you are
> willing to share! Ross
> ########################################################################################
> Generalized linear mixed model fit by the Laplace approximation
> Formula: cbind(PC_C_T, PC_C_F) ~ S_ACT_PUP_PROP + S_POOL_DIST + (1 |
> FINAL_ID_FAC) Data: MODEL_DATA_ID_2009 AIC BIC logLik deviance 1149
> 1170 -570.4 1141 Random effects: Groups Name Variance Std.Dev.
> FINAL_ID_FAC (Intercept) 0.14343 0.37872 Number of obs: 1525, groups:
> FINAL_ID_FAC, 15 Fixed effects: Estimate Std. Error z value Pr(>|z|)
> (Intercept) -3.95641 0.11727 -3374 < 2e-16 *** S_ACT_PUP_PROP 0.42450
> 0.05467 7.76 8.21e-15 *** S_POOL_DIST 0.15400 0.07245 2.13 0.0336 *
> --- Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1
> Correlation of Fixed Effects: (Intr) S_ACT_ S_ACT_PUP_P -0.170
> S_POOL_DIST -0.079 -0.045
> ########################################################################################
> Generalized linear mixed model fit by the Laplace approximation
> Formula: cbind(PC_C_T, PC_C_F) ~ S_ACT_PUP_PROP + S_POOL_DIST + (1 |
> FINAL_ID_FAC) + (1 | OBS_INDEX) Data: MODEL_DATA_ID_2009 AIC BIC
> logLik deviance 1144 1171 -567.2 1134 Random effects: Groups Name
> Variance StdDev. OBS_INDEX (Intercept) 0.40142 0.63358 FINAL_ID_FAC
> (Intercept) 011623 0.34092 Number of obs: 1525, groups: OBS_INDEX,
> 1525; FINAL_ID_FAC, 15 Fixed effects: Estimate Std. Error z value
> Pr(>|z|) (Intercept) -4.13969 0.11132 -37.19 < 2e-16 ***
> S_ACT_PUP_PROP 0.43440 0.05960 729 3.12e-13 *** S_POOL_DIST 0.15137
> 0.07498 2.02 0.0435 * --- Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*?
> 0.05 ?.? 0.1 ? ? 1 Correlation of Fixed Effects: (Intr) S_ACT_
> S_ACT_PUP_P -0.179 S_POOL_DIST -0.081 -0038
> ########################################################################################
> 
> 
> Links:
> ------
> [1] http://www.highstat.com/
> 
> 
>     [[alternative HTML version deleted]]
> 
> 
> 
> ------------------------------


    [[alternative HTML version deleted]]



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

Message: 4
Date: Tue, 18 Oct 2011 10:09:15 +1100
From: Likan Zhan <zhanlikan at hotmail.com>
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] to interpret the interactions of treatment factors
    in    lme4
Message-ID: <BLU0-SMTP2152436453390AB4DC4379FB3E40 at phx.gbl>
Content-Type: text/plain; charset="us-ascii"


> Dear Dr. Maechler,

> I am writing to ask a clarification when I am trying to use "lme4" package in R.

> If in our experiment, both A and B are three levels treatment factors:
> A:1,2,3;
> B:1,2,3

> And we use the following model to fit our data:
> model=lmer(response~A*B+(1|subject)+(1|trial),data=xx)

> and the following is the result of the fixed-effects:

> Fixed-effects

> (Intercept)    t=100, p<.000
>  A2           t=1.0, p<.9
>  A3           t=1.0, p<.9
>  B2           t=10,  p<.008
>  B3           t=1.0, p<.9
>  A2:B3        t=100, p<.000   <===
>  A3:B3        t=1.8, p<.6
>  A3:B2        t=1.0, p<.7
>  A3:B3        t=1.0, p<.6

> Could you tell me what is the meaning of the significant effect of "A2:B3",
> what is the null hypothesis of it?

> Thank you very much.

> Likan



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

Message: 5
Date: Tue, 18 Oct 2011 09:32:10 +0100
From: Gillian Raab <gillian.raab at gmail.com>
To: r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] mcmcsamps for binomial models in lmer
Message-ID:
    <CAJCiz0xAYZUQw-JJiW-a1BkUwfd4=bgxVzt6i+Vsy2ph=NU8PA at mail.gmail.com>
Content-Type: text/plain

I have some functions that uses an old version of the lme4 package where
mcmcsamp was implemented for binomial models. WHen I try to run it on the
latest lmer implementation I get messages saying "Update not yet written".

Is there any chance it may be written some time soon? The old version used
to work fine for getting reasonable confidence intervals for variance
components. If not are there alternatives that would not take too much
rewriting. Would MCMCsamp have similar syntax?

Here is some sample models from my functions.

#A univariate model:
    fitres1<-lmer(cbind(pos,neg)~1+(1|STUDY),
    family=binomial(link = "logit"),data=TPR )
    spl1<-mcmcsamp(fitres1,Nmcmc+5000)[-1:-5000,]

# a bivariate model
fitres<-lmer(cbind(pos,neg)~(GROUP)*(logit_-1)+(logit_-1|STUDY),
        family=binomial(link = "logit"),data=dd)
         spl<-mcmcsamp(fitres,Nmcmc+5000)[-1:-5000,]



With many thanks to all you R folk for making things available.

-- 
Gillian M Raab
10 Ainslie Place EH3 6AS
tel 0131 226 6234
mobile 07748 678 551

    [[alternative HTML version deleted]]



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

Message: 6
Date: Tue, 18 Oct 2011 08:43:07 +0000
From: "ONKELINX, Thierry" <Thierry.ONKELINX at inbo.be>
To: Christopher Eckert <chris.eckert at queensu.ca>,
    "r-sig-mixed-models at r-project.org" <r-sig-mixed-models at r-project.org>
Subject: Re: [R-sig-ME] longitudinal analysis of nested samples
Message-ID:
    <AA818EAD2576BC488B4F623941DA74273B8A312D at inbomail.inbo.be>
Content-Type: text/plain; charset="us-ascii"

Dear Chris,

(1 + A|B)  and (A|B) are identical. (1 + A|B) just is more explicit about the presence of the random intercept.

Best regards,

Thierry

----------------------------------------------------------------------------
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek
team Biometrie & Kwaliteitszorg
Gaverstraat 4
9500 Geraardsbergen
Belgium

Research Institute for Nature and Forest
team Biometrics & Quality Assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium

tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be

To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey



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

_______________________________________________
R-sig-mixed-models mailing list
R-sig-mixed-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


End of R-sig-mixed-models Digest, Vol 58, Issue 20
**************************************************
[[alternative HTML version deleted]]



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

Message: 4
Date: Tue, 18 Oct 2011 12:35:34 -0500
From: Douglas Bates <bates at stat.wisc.edu>
To: lopez toledo <llopezt2000 at yahoo.com.mx>
Cc: "r-sig-mixed-models at r-project.org"
<r-sig-mixed-models at r-project.org>
Subject: Re: [R-sig-ME] Doug Bates' new book?
Message-ID:
<CAO7JsnSf-c84ufSXeDe7t8fwPbFoQcZjQE+2PzXpBDE-KxAzmQ at mail.gmail.com>
Content-Type: text/plain; charset=ISO-8859-1

On Tue, Oct 18, 2011 at 11:39 AM, lopez toledo <llopezt2000 at yahoo.com.mx> wrote:
> Hi guys:
> Some months ago, I did see a message from Dr Bates, announcing a book he was writing on mixed models using lme4. Have any of you got an idea whether this book is already published?

After he finishes writing it.  The problem with writing a book about
code that you are also writing is that you are your own worst enemy.
Every time you change the code, you need to revise the book and every
time you write about how to use the software you decide that it should
be enhanced.



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

_______________________________________________
R-sig-mixed-models mailing list
R-sig-mixed-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


End of R-sig-mixed-models Digest, Vol 58, Issue 21
**************************************************


More information about the R-sig-mixed-models mailing list