Re: [R] NaNs in Nested Mixed Model

From: ONKELINX, Thierry <Thierry.ONKELINX_at_inbo.be>
Date: Wed, 09 Mar 2011 09:45:00 +0000

Dear Johan,

A few remarks.

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
>



R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Wed 09 Mar 2011 - 09:48:50 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Thu 17 Mar 2011 - 09:10:23 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive