[R-sig-ME] conditional repeatability from MCMCglmm with random slope

[John Morrongiello]

Dear list
I would like to estimate conditional repeatability of a behavioural trait from a model including random slopes fit with MCMCglmm. Would someone have some tips for how this can be done? The rptR package offers this option using bootstrapping, based on Johnson's 2014 paper<https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4368045/> (equation 11). Here, the average repeatability is estimated across the distribution of a covariate in question. I can readily estimate raw (intercept only) and a conditional (common slope) repeatability from a MCMCglmm model, but I'm not sure how to get the random slopes repeatability. Any help/ advice is much appreciated.

(Johnson, P.C.D. (2014). Extension of Nakagawa & Schielzeth's R2GLMMRGLMM2 to random slopes models. Methods in Ecology and Evolution 5: 944-946).

library(rptR)
library(MCMCglmm)

#######some data#######
behaviour <-
structure(list(fid = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 3L,
3L, 3L, 4L, 4L, 4L, 5L, 5L, 5L, 6L, 6L, 6L, 7L, 7L, 7L, 8L, 8L,
8L, 9L, 9L, 9L, 10L, 10L, 10L, 11L, 11L, 11L, 12L, 12L, 12L,
13L, 13L, 13L, 14L, 14L, 14L, 15L, 15L, 15L, 16L, 16L, 16L, 17L,
17L, 17L, 18L, 18L, 18L), .Label = c("2", "3", "5", "6", "7",
"8", "9", "10", "11", "12", "13", "15", "16", "17", "18", "19",
"20", "21"), class = "factor"), week = c(1L, 2L, 3L, 1L, 2L,
3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L,
1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L,
2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L, 3L, 1L, 2L,
3L), trait1 = c(0.225903614, 0.368017058, 0.805194805, 0.691799861,
0.654561365, 0.815674805, 0.477368175, 0.222362125, 0.09480315,
0.268691229, 0.500809061, 0.262361128, 0.365236524, 0.090701498,
0.344081069, 0.851842461, 0.477558348, 0.324364209, 0.425583266,
0.085414302, 0.311696658, 0.083493708, 0.164316797, 0.104364194,
0.156442474, 0.145193036, 0.237337192, 0.445346658, 0.377818182,
0.256871307, 0.287724075, 0.334800469, 0.324204293, 0.285696594,
0.199901004, 0.227835821, 0.709611452, 0.202676634, 0.084690774,
0.288277218, 0.282670647, 0.408824789, 0.241985203, 0.159210526,
0.0565915, 0.105827529, 0.039312489, 0.049929224, 0.292300806,
0.098316799, 0.15465606, 0.175269499, 0.271719515, 0.230084728
)), .Names = c("fid", "week", "trait1"), row.names = c(NA, -54L
), class = "data.frame")

############# Repeatability from rptR package #####

###intercept only model
rep1 <- rpt(trait1 ~ (1 | fid), grname = "fid", data = behaviour,  datatype = "Gaussian", nboot = 1000, npermut = 0)
###common week slope
rep2 <- rpt(trait1 ~ week+ (1 | fid), grname = "fid", data = behaviour, datatype = "Gaussian", nboot = 1000, npermut = 0)
###individual differences in week slope
rep3 <- rpt(trait1 ~ week+ (week | fid), grname = "fid", data = behaviour,  datatype = "Gaussian", nboot = 1000, npermut = 0)

##raw repeatability: ignore temporal change
print(rep1)
##conditional repeatability 1: shared change over time
print(rep2)
##conditional repeatability 2: change over time that differs among individuals (this is what I want to get from the equivalent MCMCglmm model)
print(rep3)

########### mcmcGLMM code ####
prior1.1 <- list(G = list(G1 = list(V = 1, nu = 0.002)), R = list(V = 1, nu = 0.002))

###intercept only model
MCMC_raw<-MCMCglmm(trait1 ~1,
             random=~fid,family=("gaussian"), prior=prior1.1,
             nitt=110000,thin=25,burnin=10000, data=behaviour, verbose=TRUE)

###common week slope
MCMC_cond1<-MCMCglmm(trait1 ~week,
             random=~fid,family=("gaussian"), prior=prior1.1,
             nitt=110000,thin=25,burnin=10000, data=behaviour, verbose=TRUE)

###common week slope
MCMC_cond2<-MCMCglmm(trait1 ~week,
             random=~us(week):fid,family=("gaussian"), prior=prior1.1,
             nitt=110000,thin=25,burnin=10000, data=behaviour, verbose=TRUE)

####repeatability formula
R= Vg/(Vg + Vr), where R= repeatability, Vg is group level variance and Vr is data level (residual) variance

###MCMCglmm derived raw repeatability
rep4<-(MCMC_raw$VCV[,"fid"]/ (MCMC_raw$VCV[,"fid"]+MCMC_raw$VCV[,"units"]))
posterior.mode(rep4)
HPDinterval(rep4)

###MCMCglmm derived conditional repeatability 1: shared change over time
rep5<-(MCMC_cond$VCV[,"fid"]/ (MCMC_ cond $VCV[,"fid"]+MCMC_ cond $VCV[,"units"]))
posterior.mode(rep5)
HPDinterval(rep5)

####How would I estimate conditional repeatability 2: change over time that differs among individuals using MCMC_cond2?

Cheers
john
--
Dr. John R. Morrongiello
School of BioSciences
University of Melbourne
Victoria 3010, Australia
T: +61 3 8344 8929
M: +61 403 338 554
E: john.morrongiello at unimelb.edu.au<mailto:%20jmorrongiell at unimelb.edu.au>
W: morrongiellolab.com<http://morrongiellolab.com/>


	[[alternative HTML version deleted]]