[R-sig-ME] model selection in MCMCglmm
David Villegas Ríos
chirleu at gmail.com
Fri Oct 23 18:12:24 CEST 2015
Hi all,
I'm running a bivariate model in MCMCglmm to estimate the covariance
between two traits and see if it is significant. For that purpose I run two
models:
In the first model, my random variance-covariance matrix is unconstrained
(us) so the variances and the covariance are estimated. This is the result:
mod1=MCMCglmm(cbind(beh1,beh2)~(trait-1), random=~*us*(trait):id,
rcov=~us(trait):units, data=widew2t, family=c("gaussian","gaussian"),
nitt=1100000, thin=100, burnin=100000, prior=prior1)
G-structure: ~us(trait):id
post.mean l-95% CI u-95% CI
eff.samp
traitbeh1:traitbeh1.id 0.2120 0.1577 0.2715
8872
traitbeh2:traitbeh1.id 0.3437 0.2321 0.4804
10000
traitbeh1:traitbeh2.id 0.3437 0.2321 0.4804
10000
traitbeh2:traitbeh2.id 1.5800 1.2018 2.0011
10000
The DIC of this model is 20363.22. Note that the covariance between both
traits is positive (0.3437) and the confidence interval is far from
including the zero. So I'd have an argument to say that the covariance is
different from zero and therefore is ok to estimate it. If we translate
that covariance into a correlation, that would be 0.59
In my second model, my random variance-covariance matrix is constrained
(idh) so the variances are estimated but the covariance is constrained to
zero. This is the result:
mod2=MCMCglmm(cbind(beh1,beh2)~(trait-1), random=~*idh*(trait):id,
rcov=~us(trait):units, data=widew2t, family=c("gaussian","gaussian"),
nitt=1100000, thin=100, burnin=100000, prior=prior1)
G-structure: ~idh(trait):id
post.mean l-95% CI u-95% CI eff.samp
traitbeh1.id 0.2062 0.1549 0.265 10000
traitbeh2.id 1.5310 1.1574 1.940 10000
The DIC of this model is 20364.96 . This is deltaDIC (unconstrained minus
constrained) is -1.74 (which is less than -2, the generally agreed thresold
to support the unconstrained model). So I should prefer the second model.
So if I look at the confidence interval in model 1, I would assume that
both traits covary/are correlated. However if I look at the DICs, the
second model seems better and supports a covariance of zero.
Is there any preferred or better method to test the significance of the
covariance between beh1 and beh2 in this case?
Many thanks.
Regards
David
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list