In daijiang/phyr: Model Based Phylogenetic Analysis

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "##",
  fig.path = "man/figures/README-",
  out.width = "100%"
)

phyr

Installation

To install this package:

options(repos = c(
  phyr = 'https://daijiang.r-universe.dev',
  CRAN = 'https://cloud.r-project.org'))
install.packages('phyr')

# or
devtools::install_github("daijiang/phyr")

Main functions

The phyr package has three groups of functions:

community phylogenetic diversity metrics (alpha: psv, psr, pse, etc. and beta: pcd), which were included in the picante package originally. They were updated with c++ to improve speed.
models to estimate correlation between functional traits while accounting for phylogenetic relationships (cor_phylo), which was included in the ape package originally. It has new syntax, much improved performance (c++), and bootstrapping option.
phylogenetic generalized linear mixed models (pglmm), which was originally included in the pez package. It has new model formula syntax that allows straightforward model set up, a faster version of maximum likelihood implementation via c++, and a Bayesian model fitting framework based on INLA.
- We hope the model formula proposed here can be used to standardize PGLMMs set up across different tools (e.g. brms for Stan).
- PGLMM for comparative data (pglmm.compare), which was originally from ape::binaryPGLMM() but has more features.

Usage examples of `pglmm()`

pglmm use similar syntax as lme4::lmer to specify random terms: add __ (two underscores) at the end of grouping variable (e.g. sp) to specify both phylogenetic and non-phylogenetic random terms; use (1|sp__@site) to specify nested term (i.e. species phylogenetic matrix V_sp nested within the diagonal of site matrix I_site) to test phylogenetic overdispersion or underdispersion. This should be the most commonly used one and is equal to kronecker(I_site, V_sp).

We can also use a second phylogeny for bipartite questions. For example, (1|parasite@host__) will be converted to kronecker(V_host, I_parasite); (1|parasite__@host__) will be converted to kronecker(V_host, V_parasite).

For details about model formula, see documentation ?phyr::pglmm. More application examples can be found in Ives 2018 Chapter 4.

library(phyr)

library(dplyr)
comm = comm_a
comm$site = row.names(comm)
dat = tidyr::gather(comm, key = "sp", value = "freq", -site) %>% 
  left_join(envi, by = "site") %>% 
  left_join(traits, by = "sp")
dat$pa = as.numeric(dat$freq > 0)
head(dat)
# phy-LMM
test1 = phyr::pglmm(freq ~ 1 + shade + (1|sp__) + (1|site) + (1|sp__@site), 
                    data = dat, family = "gaussian", REML = FALSE,
                    cov_ranef = list(sp = phylotree))
test1
# phy-GLMM
test2 = phyr::pglmm(pa ~ 1 + shade + (1|sp__) + (1|site) + (1|sp__@site), 
                    data = dat, family = "binomial", REML = FALSE,
                    cov_ranef = list(sp = phylotree))
test2
# bipartite
tree_site = ape::rtree(n = n_distinct(dat$site), tip.label = sort(unique(dat$site)))
z_bipartite = phyr::pglmm(freq ~ 1 + shade + (1|sp__) + (1|site__) + 
                            (1|sp__@site) + (1|sp@site__) + (1|sp__@site__), 
                          data = dat, family = "gaussian",REML = TRUE,
                          cov_ranef = list(sp = phylotree, site = tree_site))
z_bipartite