Skip to content

Advice for analysis of biological data - Mixed model or NESTED-Anova?

4 messages · Evan Palmer-Young, Savani Anbalagan

#
Dear all,

I was suggested in the stack exchange.com to consult in this maling list.

I have data from image analysis of zebrafish brain structures. I will
discuss our data below with some analogy to make my explanation clear.

   1. Data model: Group>Animal 1..2...3....10>Volume 1..2..3.....1000
   2. Data model: Group>Drug treatment..1..2>Animal 1..2...3....10>Volume
   1..2..3.....1000
   3. I am studying axonal synapses in Brain.
   4. I have 3 or more groups (Genotypes: Wild type, Hetero, Homozygous
   mutant)
   5. Animals are sacrificied to image them.
   6. I have 10+ animals from each group.
   7. The number and volume of the synapses are variable.
   8. Within the group, some animals have 300 synapses, some have 450
   synapses.
   9. The volume of the synapses range from 0.2 to 50. The histrogram is
   highly skewed towards lower values. A log transformation makes it look more
   normal.
   10. Some times, we also treat the different groups to a drug. So, it
   makes another level.
   11.

Analogy:

   1. > (Imagine a tree with fruits of different sizes. And I am interested
   in the size of the fruits)
   2. >(Lets say, I have trees of different species. example Indian Mango
   vs Brazilian Mango vs another Mango)
   3. >(To collect fruits, The trees are cut. )
   4. >(10+ trees in each groups)
   5. >(The number of fruits vary depending on tree to tree even within
   same group. The size of the fruit varies. There are relatively too many
   small fruits).
   6. >(Some times, fertilizers are added to tree, and then effect of fruit
   count/size is also checked)

My questions:
Could you please let me know,


   1. Should I perform Nested ANOVA or Mixed model analysis?
   2. If mixed model design, should I run the analysis on log transformed
   data or raw data? Is the distribution important for mixed model analysis?
   3. If drug treatment is added, Is it Nested or Mixed model design?
   4. For mixed model analysis how can I calculate p-value? Could you
   please let me know for both the cases. For experiments, without any drug.
   And for experiments with the drug treated vs control.
   5. These are the codes that I use to analyze my data: Could you check if
   it is correct?


My nested anova code I use:
logGFPVol.anova = aov(logVolume ~ Group + Error(Animal_ID/Group), data=data)
summary(logGFPHBVol.anova)


Mixed model code:
model2=lmer(logVolume ~ Group + (1|Animal_ID/Group ), data=data, REML =
FALSE)
summary(model2)


Please feel free to ask if I am unclear.

Many thanks,
Savani

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

*Savani Anbalagan, Ph.D*

*Dept. of Mol. Cell Biology*


*Weizmann Institute of Science234 Herzl St., Rehovot 76100,*


*ISRAELPhone: +972-8934-6158*
#
Dear Savani,
I think you are on the right track. If you use function nlme, you can get
your p-values straightaway.
With lme4, you have to employ another function (Likelihood ratio test on
full and reduced models, or Wald tests with Anova in car) to extract them:
see:
http://www.inside-r.org/packages/cran/lme4/docs/pvalues

For your model coding, make sure that the biggest group is listed FIRST.
So for you:
model2=lmer(logVolume ~ Group + (1|Animal_ID/Group ), data=data, REML =
FALSE)
Instead use
brainmodel<-nlme(logVolume ~ Group , random= ~Group/Animal_ID)
See some examples under "model specification" on this very helpful page:
http://glmm.wikidot.com/faq

Here are some nlme examples:
http://www.stat.ubc.ca/~lang/Stat527a/ex4.r

Good luck!

On Wed, May 25, 2016 at 3:32 PM, Savani Anbalagan <savani1987 at gmail.com>
wrote:

  
    
#
Hi Evan,

Thanks a lot.
When you say the biggest group. I assume it is just about the number of
animals in the group. Not about about the total number of observations for
a group.

And, when I run the code, I get the error

?
brainmodel<-nlme(logVolume ~ Group , random= ~Group/Animal_ID,
?
data=dat)
Error in model[[3]][[1]] : object of type 'symbol' is not subsettable

And

In the full or reduced models with lme4: What might be an reduced model in
my case?
Because, in expt design 1 : I only have 3 groups.


thanks,
Savani
On 25 May 2016 at 21:58, Evan Palmer-Young <epalmery at cns.umass.edu> wrote:

            

  
  
5 days later
#
Savani,
I think your problem is that your model is using group as both fixed and
random effect.
Try this instead:
So for your experiments with NO treatment, you need a model to assess the
effect of group, if you want to say, "Did the groups differ":
simple_model<-nlme(logVolume ~ Group , random= ~Animal_ID, data=savanidata)

For your experiments WITH a treatment,
treatment_model<-nlme(logVolume ~ Treatment , random= ~Group/Animal_ID,
data=savanidata)

If you actually want to know the interaction between treatment and group,
then you must have "Group" as fixed effect:
interaction_model<-nlme(logVolume ~ Treatment*Group , random= ~Animal_ID,
data=savanidata)
Then you can test the terms like this:

model1<- update(interaction_model, ~. - Treatment:Group) #excludes
interaction

anova(interaction_model, model1) #likelihood ratio test; if p<0.05, the
interaction term is doing a good job in your model, so you should keep the
interaction term
Otherwise, you can keep simplifying, like this:
model2<-update(model1, ~. - Group)
anova(model1, model2)
And so you can test the importance of each term.

Happy modeling,
Evan






On Thu, May 26, 2016 at 5:17 AM, Savani Anbalagan <savani1987 at gmail.com>
wrote: