Hi folks!
I’m a 2nd year PhD student studying root ecological strategies and mycorrhizal outsourcing. I’m interested in understanding how the evolution of three response variables (three root traits) is influenced by a categorical trait (mycorrhizal states). I’m relatively new to this space (both phylogenetics and MR-PMM) and would appreciate your guidance on the following aspects:
After reading through some vignettes and papers, I’ve put together the following formula:
brms::brmsformula(mvbind(rt1, rt2, rt3) ~ state + (1|gr(species, cov = corrmat)) + (1|p|taxa)) + set_rescor(TRUE)
The rt variables are the root traits and state is the categorical trait (with 6 states). The first random intercept is the phylogenetically structured trait correlation between species and the second random intercept is the phylogenetically unstructured trait correlations within and between species (taxa is identical to species). My data has records for 1,300 unique species (5,200 total records - many species have repeated measurements).
First question: (1|p|taxa) fits a random intercept over species and accounts for correlations between the three root traits and the random intercepts based on species, as mentioned here. Is that correct?
Second question: The variance covariance matrix used here is generated by ape::vcv.phylo which uses a BM process while the rest of my analysis is based on OU models. This discussion is basically the problem I have at the moment. Is there a way to convert my phylogeny into a vcv matrix based on an OU process without knowing how to accurately parametrize alpha and theta? I have three different continuous traits and how do I decide these two parameters which could vary depending on the traits?
Third question: If I were to also include higher taxonomic ranks in addition to the species names (say genus names and family names) so I could see the (potentially) hierarchical nature of random effects, would the following formula be appropriate? I intend to compare these using the LOO method eventually to see how their fits compare to the first one above?
-
brms::brmsformula(mvbind(rt1, rt2, rt3) ~ state + (1|gr(species, cov = corrmat)) + (1|p|taxa)) * (1|q|genus) + set_rescor(TRUE) -
brms::brmsformula(mvbind(rt1, rt2, rt3) ~ state + (1|gr(species, cov = corrmat)) + (1|p|taxa)) * (1|q|genus) * (1|r|family) + set_rescor(TRUE)
Allowing interaction between these new random intercepts with taxa since these are taxonomically hierarchically strucured. Any help will be greatly appreciated!. Thanks!