[R-sig-ME] Multi-response MCMCglmm and multivariate location effects

Rafael Maia rmaialists at gmail.com
Tue Aug 29 20:22:52 CEST 2017


Ah I see, thanks for pointing that out and that thread (11 years old!). I really appreciate everyone’s help.

Thanks, 
-R


> On Aug 29, 2017, at 11:10 AM, Ben Bolker <bbolker at gmail.com> wrote:
> 
> 
> 
> On 17-08-29 02:04 PM, Rafael Maia wrote:
>> Thank you so much Jarrod, that seems to be exactly what I need, cheers. Would that be equivalent to looking at the joint posterior distribution of coefficient estimates is greater/smaller than zero given their means? Crudely, something like 
>> 
>> sum(m1$Sol[,4] >= 0 & m1$Sol[,5] <= 0 & m1$Sol[,6] >= 0)
>> 
>> Given the coefficients resulting from using set.seed(1982) ? 
>> 
>> Thanks,
>> -Rafael
> 
>   Jarrod's response accounts for the correlation of the estimates as
> well, effectively giving an elliptical boundary.  Your solution gives a
> rectangular/hypercube boundary, which will be very far off in some
> cases.  (About a million years ago Doug Bates offered a similar solution
> for output from the now-defunct mcmcsamp(), which did post-hoc MCMC
> sampling from an lmer fit:
> https://stat.ethz.ch/pipermail/r-help/2006-September/113184.html <https://stat.ethz.ch/pipermail/r-help/2006-September/113184.html> )
>> 
>> 
>>> On Aug 29, 2017, at 4:00 AM, Jarrod Hadfield <j.hadfield at ed.ac.uk> wrote:
>>> 
>>> Hi,
>>> 
>>> I think the appropriate model is 2) but you want to know whether the 3 predictor effects could *all* be zero or not. I usually mix up my paradigms and perform a Wald test using the posterior means and posterior (co)variances. Simulation suggests it works well: the p-values have a uniform distribution under the null hypothesis. For example,
>>> 
>>> 
>>> n<-100
>>> 
>>> V<-matrix(0.25,3,3)
>>> diag(V)<-1
>>> 
>>> mu<-rep(0,3)
>>> beta<-c(-0.2,0.2,0)
>>> 
>>> predictor<-rnorm(n)
>>> 
>>> y<-predictor%*%t(beta)+MASS::mvrnorm(n, mu, V)
>>> 
>>> dat<-data.frame(x=y[,1], y=y[,2], z=y[,3], predictor=predictor)
>>> 
>>> m1<-MCMCglmm(cbind(x,y,z)~trait+trait:predictor-1, rcov=~us(trait):units, data=dat, family=rep("gaussian",3))
>>> 
>>> aod::wald.test(cov(m1$Sol[,4:6]), colMeans(m1$Sol[,4:6]), Terms=1:3)$result$chi2["P"]
>>> 
>>> Cheers,
>>> 
>>> Jarrod
>>> 
>>> 
>>> 
>>> On 22/08/2017 01:29, Rafael Maia wrote:
>>>> Dear list,
>>>> 
>>>> Suppose I have a dataset consisting of 3 continuous, Gaussian variables and 1 categorical variable. I would like to run a multi-response model using the 3 continuous variables as response variables and the categorical variable as a predictor. These 3 variables are somewhat arbitrary, simply representing XYZ Cartesian coordinates, and it is not particularly easy to interpret them individually. They are also in different scales and the effect of the predictor can have different signs and magnitudes on the variables. Therefore, I’m mostly interested in understanding if the categorical variable influences the overall location of the sample in this multivariate cartesian space - as a simple MANOVA would.
>>>> 
>>>> I’ve looked into MCMCglmm as what seems to be the best approach to this question, and one I’m a bit familiar with. From my understanding, there are two different ways I could parametrize the fixed part of my model:
>>>> 
>>>> 1) fixed = cbind(x,y,z) ~ trait:predictor-1, which would estimate the effects of the predictor across each of the 3 responses individually, and is like an interaction between the predictor and the responses;
>>>> 
>>>> 2) fixed = cbind(x,y,z) ~ trait+predictor-1, which would estimate a single effect across all responses, and from what I understand, fundamentally assumes that the effect of the predictor is of the same magnitude and in the same direction across all responses.
>>>> 
>>>> So from my understanding, I _think_ I would want something similar to (2); however, because of the structure of the data as described above, the effect of the categorical predictor will often not be of the same magnitude/direction across xyz, so that doesn’t seem to be an appropriate model. I would rather not make my inference based on the trait:predictor fixed effects structure because, from my naive understanding, I guess this would be akin to running 3 ANOVAs instead of a MANOVA, and beat the purpose of using a multi-response model? And I think it would also be possible to not have above-zero effect on any of the 3 variables independently, yet still have a measurable difference in the multivariate location of the groups?
>>>> 
>>>> I guess this question is somewhat related to model selection, for which there doesn’t seem to have a satisfactory approach for when using  MCMCglmm? I remember reading somewhere that DIC model selection for models with different fixed effects wasn’t really recommended.
>>>> 
>>>> (As you’d expect, this is a simplified version of the dataset I have, where there are multiple predictors and random effects, to demo these issues, which is why I’m looking for a solution using MCMCglmm!)
>>>> 
>>>> I’ve uploaded some code to reproduce what I’m trying to describe here: http://rpubs.com/rmaia/multimcmc <http://rpubs.com/rmaia/multimcmc><http://rpubs.com/rmaia/multimcmc <http://rpubs.com/rmaia/multimcmc>><http://rpubs.com/rmaia/multimcmc <http://rpubs.com/rmaia/multimcmc> <http://rpubs.com/rmaia/multimcmc <http://rpubs.com/rmaia/multimcmc>>>
>>>> 
>>>> Any help would be greatly appreciated. Thanks!
>>>> 
>>>> Best,
>>>> Rafael Maia
>>>>>>>> Simons Junior Fellow, Columbia University
>>>> Rubenstein Lab
>>>> Department of Ecology, Evolution and Environmental Biology
>>>> New York, NY
>>>> http://www.rafaelmaia.net <http://www.rafaelmaia.net/> <http://www.rafaelmaia.net/ <http://www.rafaelmaia.net/>>
>>>> 	[[alternative HTML version deleted]]
>>>> 
>>>> _______________________________________________
>>>> R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org> <mailto:R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>> mailing list
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models> <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>>
>>> 
>>> 
>>> -- 
>>> The University of Edinburgh is a charitable body, registered in
>>> Scotland, with registration number SC005336.
>>> 
>>> _______________________________________________
>>> R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org> <mailto:R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>> mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models> <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>>
>> 
>> 	[[alternative HTML version deleted]]
>> 
>> _______________________________________________
>> R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org> mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>
>> 
> 
> _______________________________________________
> R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org> mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>

	[[alternative HTML version deleted]]



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