Skip to content

Help with coding an intervention analysis with lmer()

7 messages · Dennis Murphy, Dave Schoeman, Ista Zahn +2 more

Original
Ista Zahn
# 
Hi all,
I hesitate to weigh in on this, because I still consider myself a
newbie with lme4. But I really think you're on the wrong track here,
and I trust that someone will correct me if I'm wrong.
On Sat, Aug 28, 2010 at 2:58 AM, Dave Schoeman <d.schoeman at ulster.ac.uk> wrote:
I think this part is reasonable.
I think this is where you go off the rails. Nesting refers to
_grouping variables_, not to fixed effects. If I'm understanding
correctly, group is not a grouping variable, it is a fixed effect with
2 repeatable levels: control and experimental.

 In this case, and generally with my real data, the estimated variance
for Group is very small (zero here), so could probably be eliminated,
reducing the random term to (1|Subject). This would make sense,
anyway, if lmer() automatically accounts for nesting (uniquely-coded
Subjects).
I don't think this is right. I don't think you need to specify nesting
this way, especially in the fixed-effects part of the model. The way I
see it, you have only one grouping factor, and that is Subject.

If we go back to your original model specification, you had

mod <- lmer(Attack ~ Group*BA + BA/Year + BA/Year:Group +
Group/Subject:BA + (1|Group/Subject), data = dat, family = poisson)

in my view that is gobbledygook. You have Group as a fixed effect
predictor, but also as a grouping factor, and I don't think that makes
any sense. I also don't think that those nesting statements in the
fixed-effects specification make any sense, but it's possible I'm the
one who doesn't understand lmer syntax... Anyway, I would do something
like

Subject <- c(rep(c(1:20), 8))
Group <- c(rep(c(rep("C", 10), rep("E", 10)), 8))
BA <- c(rep("Bef", 80), rep("After", 80))
Year <- c(rep(1990, 20), rep(1991, 20), rep(1992, 20), rep(1993, 20),
rep(1994, 20), rep(1995, 20), rep(1996, 20), rep(1997, 20))

Attack <- c(rpois(10, lambda = 10), rpois(10, lambda = 12), rpois(10,
lambda = 10), rpois(10, lambda = 12), rpois(10, lambda = 10),
rpois(10, lambda = 12), rpois(10, lambda = 10), rpois(10, lambda =
12), rpois(10, lambda = 15), rpois(10, lambda = 3), rpois(10, lambda =
15), rpois(10, lambda = 3), rpois(10, lambda = 15), rpois(10, lambda =
3), rpois(10, lambda = 15), rpois(10, lambda = 3))

dat <- data.frame(Attack, Subject, Group, BA, Year)
dat$Subject <- factor(dat$Subject)
dat$Year <- dat$Year - 1994
dat$BA <- relevel(dat$BA, ref="Bef")


library(lme4)

mod <- lmer(Attack ~ Group*Year + Group*BA + (1|Subject), data=dat,
family=poisson)

This tests whether the difference in Attack before and after differs
by Group, and whether the trajectory of Attack over time (Year)
differs by Group. Note that year is numeric and re-centered at 1994
(the year the intervention started).

Hope this was helpful. And please, those of you who know more than I
do, please jump in and correct me if I've given bad advice here...

Ista
--
Ista Zahn
Graduate student
University of Rochester
Department of Clinical and Social Psychology
http://yourpsyche.org
Dave Schoeman
# 
Many thanks for the input! This is very interesting.
On 28 Aug 2010, at 17:00, Ista Zahn wrote:

            

      
      
      
    

        
<snip>

            

      
      
      
    

        
OK, I think this is a fundamental question. In my understanding, nesting arises when levels of one factor (here Subject) are present in only one level of another factor (here Group). In other words, Subjects 1-10 occur only in the Control Group, and Subjects 11-20 occur only in the Experiment Group. If each Subject occurred in each Group, Subject and Group would be crossed. The same goes for BA and Year: Years 1990-1993 are Before and Years 1994-1997 are After, so again Year is nested within BA. If this nesting is correct, then should I not try to estimate the effects of each? And if I should, the question is HOW? I agree with you that my model specification looks like gobbledygook, but it is just an attempt to model all components of variance.

<snip>

            

      
      
      
    

        
Yes, I agree that your model tests the hypotheses I'm interested in, and the model-building process ends up in the same place, irrespective if whether you start with my model or yours (both with the simulated data provided, or with the real data). BUT, should the model not start by including all terms?

            

      
      
      
    

        
Yes, very helpful, indeed. Thanks. And thanks, also, to Dennis who wrote earlier.

I do think that it would be useful, though, to get further input on the nesting concepts here and how to specify them using the fixed and random effects...

- Dave

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

      
      
    

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Andrew Robinson
    

    

      
      #
    

  


  

          

      
      
      
    

        
Just as a general point, there is a case when fixed effect predictor can also appear as a grouping effect: in split-plot designs. Indeed, this occurs in one of the examples in Pinheiro & Bates (2000).  It is a bit unclear in coding, and it would be cleaner if the coder were to create separate factors for the fixed and the random effects, even if they are the same. The algorithms allow the shortcut, though.    

Cheers

Andrew

            

      
      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Douglas Bates
    

    

      
      #
    

  


  

      
On Sun, Aug 29, 2010 at 5:26 AM, Andrew Robinson

        
<mensurationist at gmail.com> wrote:

            

      
      
      
    

        
There is such an example but it doesn't apply in this case.  Even
though the whole plot factor was being used as a grouping factor for
the random effects it was used in the form of an interaction with the
blocking factor.  In other words, if F1, F2, etc are fixed-effects
factors and B1, B2, etc. are environmental factors (such as location)
then a formula of

Y ~ F1 + (1|B1/F1)

which expands to

Y ~ F1 + (1|B1) + (1|B1:F1)

is sensible but a term of the form (1|F1) is not.