[R-sig-ME] multivariate mixed nested model

Claudio oppela at gmail.com
Fri Feb 3 16:22:28 CET 2017


Dear Thierry,
yes, I recorded only once each trait for each individual (i.e., it was
not a longitudinal study).
And, yes, the idea is that each individual has a different effect.
thanks again
Claudio


Il giorno gio, 02/02/2017 alle 10.12 +0100, Thierry Onkelinx ha
scritto:
> Dear Claudio,
> 
> It looks like you have only one observation (of each trait) per
> individual. That is not sufficient to fit (0 + trait | individual).
> Hence the error message.
> 
> Using (1 | individuals) implies that each individual has a effect
> which is common for all traits. Whether that is a relevant assumption
> depends on your experiment. Tip: write the model as an equation for
> each trait. Then see if the random intercept for individual makes
> sense.
> 
> I'm not proficient with MCMCglmm. So I can help you with that.
> 
> Best regards,
> 
> 
> 
> ir. Thierry Onkelinx
> Instituut voor natuur- en bosonderzoek / Research Institute for
> Nature and Forest 
> team Biometrie & Kwaliteitszorg / team Biometrics & Quality
> Assurance 
> Kliniekstraat 25
> 1070 Anderlecht
> Belgium
> 
> 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
> 
> 2017-02-01 20:01 GMT+01:00 Claudio <oppela at gmail.com>:
> > Dear Thierry,
> > thanks a lot, it is exactly the kind of suggestion I was looking
> > for!
> > However, when using (0 + traits|individuals) for the random effect
> > of
> > individuals I got the message:
> > "Errore: number of observations (=1431) <= number of random effects
> > (=1431) for term (0 + traits | individuals); the random-effects
> > parameters and the residual variance (or scale parameter) are
> > probably
> > unidentifiable"
> > The models work when I use for the individual part (1|individuals).
> > 
> > A second further question: in order to include in the model a
> > continuous fixed covariate, am I doing the right thing when
> > scripting:
> > 
> > mumo2 <- lmer(value~0 + traits + traits:species + traits:covariate
> > (0 +
> > traits|species:populations) + (1|individuals), mdata, REML=FALSE)
> > 
> > About cbind and MCMCglmm, below there is an example which causes
> > the
> > message:
> > "Error in `[<-.data.frame`(`*tmp*`, , response.names, value = c(1,
> > 2,
> > 3,  : missing values are not allowed in subscripted assignments of
> > data
> > frames"
> > 
> > a = c(1, 2, 3, 5, 5, 7, 8, 9, 4, 5, 7, 8, 9, 4, 1, 2, 3, 5, 5, 4,
> > 5, 7,
> > 8, 9, 4, 1) 
> > b = c(8, 9, 4, 5, 7, 8, 9, 4, 1, 2, 3, 1, 2, 3, 5, 5, 7, 4, 1, 2,
> > 3, 1,
> > 2, 3, 6, 6)
> > d = c(8, 9, 4, 5, 5, 7, 8, 9, 4, 5, 7, 8, 9, 4, 1, 2, 3, 2, 3, 1,
> > 2, 3,
> > 5, 5, 7, 6)  
> > s = c("m", "m", "f", "f", "m", "m", "m", "f", "m", "f", "f", "f",
> > "m",
> > "f", "f", "f", "m", "m", "m", "m", "m", "f", "m", "f", "f", "f") 
> > df = data.frame(a, b, d, s)
> > 
> > y<-cbind(a, b, d)
> > prior <- list(R = list(V = 1, nu = 0.002))
> > m <- MCMCglmm(y ~ s, family = "gaussian" , data = df, prior =
> > prior,
> > verbose = FALSE, pl = TRUE)
> > summary(m)
> > 
> > all the best (and thanks again)
> > Claudio
> > 
> > Il giorno lun, 30/01/2017 alle 16.02 +0100, Thierry Onkelinx ha
> > scritto:
> > > Dear Claudio,
> > >
> > > I this you need to add the interaction with traits to all the
> > fixed
> > > and random effects. Otherwise you assume that these have the same
> > > effect for each trait. Note that 0 + traits is identical to
> > traits -
> > > 1.
> > >
> > > mumo1 <- lmer(value~0 + traits + traits:species + (0 +
> > > traits|species:populations) + (0 + traits|individuals), mdata,
> > > REML=FALSE)
> > >
> > > Your second question needs a reproducible example.
> > >
> > > Best regards,
> > >
> > > ir. Thierry Onkelinx
> > > Instituut voor natuur- en bosonderzoek / Research Institute for
> > > Nature and Forest 
> > > team Biometrie & Kwaliteitszorg / team Biometrics & Quality
> > > Assurance 
> > > Kliniekstraat 25
> > > 1070 Anderlecht
> > > Belgium
> > >
> > > 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
> > >
> > > 2017-01-28 16:43 GMT+01:00 Claudio <oppela at gmail.com>:
> > > > Hi all.
> > > > I collected six body features (bf1-bf6)from three populations
> > of a
> > > > salamander and from two populations of another sister species
> > of
> > > > salamander.
> > > > I would evaluate how the species (fixed) and population
> > belonging
> > > > (random) affect the body features, by comparing models built
> > with
> > > > lme4.
> > > > For some models, I also want to include bf6 as covariate. Thus,
> > in
> > > > case
> > > > of univariate analyses, some models, for example, could be:
> > > > mo1<-lmer(bf1~species+(1|species:population), data, REML=FALSE)
> > > > mo2<-lmer(bf1~species+bf6+(1|species:population), data,
> > REML=FALSE)
> > > >
> > > > However, I want to fit multivariate models, and my post is
> > about
> > > > this.
> > > > First, I melted the data:
> > > > mdata<-melt(data, id.vars = c("species", "population", "bf6"),
> > > > measure.vars = c("bf1", "bf2","bf3","bf4","bf5"), variable.name
> > =
> > > > "traits)
> > > >
> > > > Now the question.
> > > > 1) Are the multivariate versions of the models mo1 and mo2
> > above
> > > > mumo1<-lmer(value~traits -1 + species + (1|species:populations)
> > +
> > > > (1|individuals), mdata, REML=FALSE)
> > > > mumo1<-lmer(value~traits -1 + species + bf6 +
> > > > (1|species:populations) +
> > > > (1|individuals), mdata, REML=FALSE)
> > > >
> > > > A secondary question, which in case I will move to a new post:
> > > > it seemed to me that building multivariate models with MCMCglmm
> > is
> > > > easier. However, cbind did not work, even without missing
> > values:
> > > > to
> > > > your knowledge, is there any issue?
> > > >
> > > > thanks in advance
> > > > Claudio  
> > > >
> > > > _______________________________________________
> > > > R-sig-mixed-models at r-project.org mailing list
> > > > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >



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