[R-sig-ME] R-sig-mixed-models Digest, Vol 161, Issue 1

Szymek Drobniak ger@|ttee @end|ng |rom gm@||@com
Fri May 1 12:20:42 CEST 2020

Hi Carolina, your second model assumes homogenous tips variance across variables (random=~tips) - but in your prior you structure both tips and residual variances as 2x2 matrices. In the second model you have to use us/idh variance function for tips effects (e.g. us(trait):tips), and add rcov term to model similar structure for units (e.g. rcov=~us(trait):units).


Dr hab. Szymon Drobniak

Institute of Environmental Sciences
Jagiellonian University, Kraków, Poland

School of Biological, Environmental and Earth Sciences
University of New South Wales, Sydney, Australia

Google Scholar profile

> Message: 1
> Date: Thu, 30 Apr 2020 08:30:55 +0000
> From: Carolina Karoullas
> <carolina.karoullas using postgrad.manchester.ac.uk>
> To: "r-sig-mixed-models using r-project.org"
> <r-sig-mixed-models using r-project.org>
> Subject: [R-sig-ME] Creating a phylogenetically corrected multivariate
> linear model using MCMCglmm
> Message-ID:
> <56C3CB6634CED94898F1E5FC8E71A09CB602BF47 using MBXP03.ds.man.ac.uk>
> Content-Type: text/plain; charset="utf-8"
> Hi all,
> I'm trying to use the package MCMCglmm to run a multivariate linear model that is phylogenetically corrected. Here is a subset of my data (there are 210 entries in total for 67 species and 6 clusters):
> Names PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8 PC9 PC10 PC11 tips Clusters
> Accipiter_gentilis1 -1.34E-01 1.98E-02 -1.72E-02 4.00E-02 -1.93E-03 2.45E-02 -1.28E-04 1.03E-02 9.09E-03 -5.33E-03 -1.30E-02 Accipiter_gentilis "Soaring"
> Accipiter_gentilis2 -1.26E-01 1.22E-02 1.66E-02 9.22E-03 1.15E-02 1.68E-02 -1.50E-02 1.02E-02 7.93E-03 -8.47E-03 -1.02E-02 Accipiter_gentilis "Soaring"
> Accipiter_nisus1 -1.81E-01 5.76E-02 -1.82E-02 6.15E-02 -9.25E-03 3.40E-02 -1.77E-02 5.45E-03 7.01E-03 -2.07E-02 -8.78E-03 Accipiter_nisus "Soaring"
> Accipiter_nisus2 -2.00E-01 7.05E-02 -1.12E-02 5.94E-02 3.49E-03 3.10E-02 -1.58E-03 -1.55E-03 6.92E-03 -3.54E-02 -1.80E-02 Accipiter_nisus "Soaring"
> Accipiter_nisus3 -8.14E-02 -3.39E-04 -8.88E-03 4.25E-02 -5.48E-04 -8.51E-03 5.07E-03 4.56E-03 1.97E-02 -1.46E-02 -1.43E-02 Accipiter_nisus "Soaring"
> Accipiter_nisus4 -2.06E-01 7.05E-02 -2.17E-02 6.38E-02 -1.61E-02 2.80E-02 8.70E-03 -5.96E-03 6.15E-03 -5.29E-02 -2.05E-02 Accipiter_nisus "Soaring"
> Actitis_hypoleucos1 2.27E-02 -2.74E-03 4.79E-02 -2.30E-02 -2.76E-02 -2.36E-02 1.70E-02 2.43E-03 3.82E-03 1.15E-02 -9.87E-03 Actitis_hypoleucos "Continuous flapping"
> Actitis_hypoleucos2 6.67E-02 -1.05E-02 5.12E-02 -2.65E-02 -3.21E-02 -2.61E-02 3.21E-03 7.46E-03 7.29E-03 4.70E-03 -1.37E-02 Actitis_hypoleucos "Continuous flapping"
> Aix_sponsa1 -3.70E-02 -1.41E-02 1.13E-02 3.16E-02 2.32E-02 -1.70E-02 2.32E-02 1.91E-03 2.91E-02 -7.71E-03 7.40E-03 Aix_sponsa "Continuous flapping"
> Aix_sponsa2 1.03E-02 -4.08E-02 -6.62E-03 1.19E-02 2.83E-02 -1.49E-02 3.78E-02 6.98E-03 2.91E-02 -4.32E-03 2.54E-03 Aix_sponsa "Continuous flapping"
> Aix_sponsa3 1.19E-02 -3.48E-02 -1.53E-02 3.76E-03 2.17E-02 1.47E-02 8.84E-03 2.39E-02 -9.20E-03 -1.78E-02 8.76E-04 Aix_sponsa "Continuous flapping"
> Aix_sponsa4 -3.37E-02 -1.75E-02 -8.06E-03 3.64E-02 -5.50E-03 1.03E-02 2.37E-02 3.33E-03 -1.04E-03 -2.00E-02 5.89E-03 Aix_sponsa "Continuous flapping"
> Aix_sponsa5 -2.30E-02 -9.59E-03 1.06E-02 3.01E-02 7.10E-03 -1.23E-02 2.08E-02 1.17E-02 1.59E-03 2.83E-03 8.75E-03 Aix_sponsa "Continuous flapping"
> Aix_sponsa6 -1.70E-02 -2.98E-02 -1.96E-02 1.76E-02 1.23E-02 4.92E-03 5.45E-03 1.99E-02 -6.43E-03 -9.63E-04 1.99E-03 Aix_sponsa "Continuous flapping"
> Aix_sponsa7 2.57E-02 -4.22E-02 -1.60E-02 1.75E-02 3.41E-03 5.80E-03 2.89E-02 6.10E-03 7.12E-03 2.75E-03 4.99E-03 Aix_sponsa "Continuous flapping"
> Aix_sponsa8 4.09E-02 -4.46E-02 -4.49E-03 2.24E-02 2.37E-03 -5.90E-03 2.78E-02 -8.26E-04 1.17E-02 -5.71E-03 -1.77E-03 Aix_sponsa "Continuous flapping"
> I created a univariate model taking inspiration from this link:
> https://github.com/JonBrommer/Multivariate-Mixed-Models-in-R/wiki/MCMCglmm-examples#organisational-level-4https://github.com/JonBrommer/Multivariate-Mixed-Models-in-R/wiki/MCMCglmm-examples#organisational-level-4
> And when I try to run it, it works (code below, phylo refers to my tree):
> Ainv<-inverseA(phylo,nodes="TIPS",scale=F)$Ainv
> p.var=var(data[,c("PC1")])
> prior1<-list(R=list(V=(p.var),nu=0.002),G=list(G1=list(V=(p.var),nu=0.002)))
> m7.phylo<-MCMCglmm(PC1~Clusters,
> random=~tips,
> family=rep("gaussian",1),
> ginverse=list(tips=Ainv),
> data=data,
> prior=prior1)
> However, as soon as I try to make a multivariate model, I get an error:
> Ainv<-inverseA(phylo,nodes="TIPS",scale=F)$Ainv
> p.var=var(data[,c("PC1","PC2")])
> prior1<-list(R=list(V=(diag(2)*p.var),nu=0.002),G=list(G1=list(V=(diag(2)*p.var),nu=0.002)))
> m7.phylo<-MCMCglmm(cbind(PC1,PC2)~Clusters,
> random=~tips,
> family=rep("gaussian",2),
> ginverse=list(tips=Ainv),
> data=data,
> prior=prior1)
> Error in priorformat(if (NOpriorG) { :
> V is the wrong dimension for some prior$G/prior$R elements
> I do want to use all 11 PCs in the model so it's not encouraging that I can't seem to get it to work with just 2 of them... Does anyone have any ideas of what could have gone wrong? I would like to create other models with the same structure using different data and trees so it would be good to understand what's going on and how to create a prior properly for next time.
> Thanks,
> Carolina

	[[alternative HTML version deleted]]

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