[R-sig-ME] comparing random effect variation between data sets
John Morrongiello
jrmorrongiello at gmail.com
Thu Jun 25 15:07:47 CEST 2015
Hi all
I have movement data for two species of fish (say species X and species Y).
The response variable is 'distance from home' (distance) and multiple
measurements (5-17) are made for 20-30 individuals per species over a
period of 200 days. Distance is continuous and strictly positive so I am
using a gamma distribution with a log link. My collegues would like me to
fit seperate models for each species to explore temporal patterns in
movement (time). We believe there will be non-linear patterns in movement
through time so have I have fitted GAMMs
I'm happy with the two models and their interpretation. I would, however,
like to compare the level of among-individual variation in movement between
species. The species-average trends are pretty similar but it looks like
there is more variance among individuals for species X than species Y.
I was thinking of extracting the individual level random effects from the
two models and performing a variance test on these. My logic is that given
the response variables are the same and the model structures are the same,
I could make a meaningful comparison of random effects. The overall model
intercepts are different but I thought this is not an issue here as I'm
only interested in the variance of the random effects, not their average.
Is this approach valid, or am I violating assumptions/ these data are not
directly comparable ? If not, what would be a potential way to perform this
comparision? Below is the code I'm thinking of using.
Thanks very much for your time
John
#####GAMMs
X1<-gamm4(distance ~ s(time,k=4), random=~(1|CODE),data=speciesX,
family=Gamma(link=log))
Y1<-gamm4(distance ~ s(time,k=4), random=~(1|CODE),data=speciesY,
family=Gamma(link=log))
######extract random effects from the two models
X1ranef<-ranef(X1$mer)###extract random effects from model
X1ranef<-data.frame(matrix(unlist(X1ranef$'CODE'),ncol=1))##convert list of
random effects to data frame
names(X1ranef)<-c('CODE')
Y1ranef<-ranef(X1$mer)###extract random effects from model
Y1ranef<-data.frame(matrix(unlist(Y1ranef$'CODE'),ncol=1))##convert list of
random effects to data frame
names(Y1ranef)<-c('CODE')
######perform variance test
var.test(X1ranef$CODE,Y1ranef$CODE)
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list