Dear Johan,
A few remarks.
- R-sig-mixed models is a better list for asking questions about mixed model.
- I presume that Nymphs is the number of insects? In that case you need a
generalised linear (mixed) model with poisson family
- What are you interessed in? The variability among genotypes or the effect of
each genotype.
You can achieve the first with a glmm like glmer(Nymphs ~ Species +
(1|Genotype), family = "poisson"). Genotype will be implicitly nested
in Species. Note that since you have only 4 genotypes, you will not get very
reliable estimates of the genotype variance.
For the latter you cannot use a mixed model so you need a simple glm(Nymphs ~
Species/Genotype, family = "poisson"). Note that several coefficients
will be NaN, because you cannot estimate them.
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
> -----Oorspronkelijk bericht-----
> Van: r-help-bounces at r-project.org
> [mailto:r-help-bounces at r-project.org] Namens Johan Stenberg
> Verzonden: dinsdag 8 maart 2011 16:52
> Aan: r-help at r-project.org
> Onderwerp: [R] NaNs in Nested Mixed Model
>
> Dear R users,
>
> I have a problem with something called "NaNs" in a nested mixed
model.
>
> The background is that I have studied the number of insect
> nymphs emerging from replicated Willow genotypes in the
> field. I have 15 replicates each of 4 Willow genotypes
> belonging two 2 Willow species.
> Now I want to elucidate the effect of Willow genotype on the
> number of emerging nymphs. Previously I performed a simple
> one-way anova with "genotype" as explanatory factor and
> "number of nymphs emerging" as dependent variable, but the
> editor of the journal I've submitted this piece to wants me
> to nest Willow genotype within Willow species before he
> accepts the paper for publication [Species*Genotype(Species)].
>
> The fact that I didn't include "Willow species" as a factor
> in my initial analysis reflects that I am not very interested
> in the species factor per se - I am just interested in if
> genetic variation in the host plant is important, but
> "species" is of course a factor that structures genetic
diversity.
>
> I thought the below model would be appropriate:
>
> > model<-lme(Nymphs~Species*Genotype,random=~1|Species/Genotype)
>
> ...but I then get the error message "Error in MEEM(object, conLin,
> control$niterEM) : Singularity in backsolve at level 0, block 1"
>
> I then tried to remove "Genotype" from the fixed factors, but
> then I get the error message "NaNs produced".
>
> > model<-lme(Nymphs~Species,random=~1|Species/Genotype)
> > summary(model)
> Linear mixed-effects model fit by REML
> Data: NULL
> AIC BIC logLik
> 259.5054 269.8077 -124.7527
>
> Random effects:
> Formula: ~1 | Species
> (Intercept)
> StdDev: 0.9481812
>
> Formula: ~1 | Genotype %in% Species
> (Intercept) Residual
> StdDev: 0.3486937 1.947526
>
> Fixed effects: Nymphs ~ Species
> Value Std.Error DF t-value p-value
> (Intercept) 2.666667 1.042243 56 2.558585 0.0132
> Speciesviminalis -2.033333 1.473954 0 -1.379510 NaN
> Correlation:
> (Intr)
> Speciesviminalis -0.707
>
> Standardized Within-Group Residuals:
> Min Q1 Med Q3 Max
> -1.4581821 -0.3892233 -0.2751795 0.3439871 3.1630658
>
> Number of Observations: 60
> Number of Groups:
> Species Genotype %in% Species
> 2 4
> Warning message:
> In pt(q, df, lower.tail, log.p) : NaNs produced
> ***********
>
> Do you have any idea what these error messages mean in my
> case and how I can get around them?
>
> Thank you on beforehand! (data set attached).
>
> Johan
>