12: Phylogenetic Comparative Methods

General Overview

Dean Adams, Iowa State University

The Accumulation of Biodiversity

Taxonomic Diversity

The Accumulation of Biodiversity

Taxonomic Diversity

Morphological Diversity

The Accumulation of Biodiversity

Taxonomic Diversity

Morphological Diversity

How do we characterize patterns, and hypothesize processes?

Comparative Biology Tradition

Trait correlations often used to study coevolution and adaptation

Species values commonly utilized

The problem?

Species are not independent of one another

Phylogenetic Comparative Methods (PCMs)

Phylogenetic comparative methods condition the data on the phylogeny during the analysis

Allows one to assess trait covariation while accounting for the non-independence due to shared evolutionary history

PCMs: An Incomplete Historical Walk

The conceptual development of PCMs

Outline

PCMs: General Statistical Concepts

PCM: General Model

Most PCMs use GLS (generalized least squares) as a model:

\[\small\mathbf{Y}=\mathbf{X{\hat{\beta}}+\epsilon}\] Here, \(\small\epsilon\) is not iid but is \(\small\sim\mathcal{N}(0,\textbf{V})\): containing expected covariation between taxa as described by the phylogeny (in matrix \(\small\mathbf{V}\))

\(\small\mathbf{V}\) is the phylogenetic covariance matrix

Describes the amount of evolutionary time species share via common ancestry (and thus how similar they are expected to be)

Sometimes V is called C (particularly when derived under Brownian motion)

PCM: General Model

\[\small\mathbf{Y}=\mathbf{X{\hat{\beta}}+\epsilon}\]

Model design (\(\small\mathbf{X}\)) describes the type of analysis

Parameters of model (and model significance) obtained in various ways

see Adams and Collyer. Syst. Biol. (2018); Adams and Collyer. Ann. Rev. Ecol. Evol. Syst. (2019)

1: Phylogenetic Regression

Evaluate \(\small\mathbf{Y}=\mathbf{X}\mathbf{\beta } +\mathbf{E}\) in a phylogenetic context

The workhorse of PCMs

1: Phylogenetic Regression

ANOVA and regression models that account for phylogeny

Requires a model of evolutionary change: typically Brownian motion (BM)

1: Phylogenetic Regression

ANOVA and regression models that account for phylogeny

Requires a model of evolutionary change: typically Brownian motion (BM)

\(\Delta\mu = 0\)

\(\sigma^2\) (variance among taxa) \(\uparrow\) \(\propto\) time

1: Phylogenetic Regression

ANOVA and regression models that account for phylogeny

Requires a model of evolutionary change: typically Brownian motion (BM)

\(\Delta\mu = 0\)

\(\sigma^2\) (variance among taxa) \(\uparrow\) \(\propto\) time

Statistical model: \(\small\mathbf{Y}=\mathbf{X}\mathbf{\beta } +\mathbf{E}\) where \(\small\mathbf{E} \sim \mathcal{N}(0,\textbf{V})\)

Several implementations possible (yield identical results if implemented properly)

Phylogenetically Independent Contrasts

Estimate contrast scores between pairs of taxa (tips or nodes)

Use contrasts for analyses (OLS solution)

see Felsenstein. Am. Nat. (1985)

Phylogenetically Independent Contrasts

Coefficients found as: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T}_{pic} \mathbf{X}_{pic}\right )^{-1}\left ( \mathbf{X}^{T}_{pic} \mathbf{Y}_{pic}\right )\)

Phylogenetically Independent Contrasts: Example

What is the association between Y and X?

Phylogenetically Independent Contrasts: Example

What is the association between Y and X?

anova(lm(pic.y~pic.x + 0))
## Analysis of Variance Table
## 
## Response: pic.y
##           Df  Sum Sq Mean Sq F value  Pr(>F)   
## pic.x      1 14.3519 14.3519  19.285 0.00461 **
## Residuals  6  4.4651  0.7442                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
coef(lm(pic.y~pic.x+0))
##     pic.x 
## 0.8846971

Phylogenetic Generalized Least Squares (PGLS)

Use GLS model with non-independent error structure

Statistical model: \(\small\mathbf{Y}=\mathbf{X}\mathbf{\beta } +\mathbf{E}\) where \(\small\mathbf{E} \sim \mathcal{N}(0,\textbf{V})\)

Phylogenetic Generalized Least Squares (PGLS)

Use GLS model with non-independent error structure

Statistical model: \(\small\mathbf{Y}=\mathbf{X}\mathbf{\beta } +\mathbf{E}\) where \(\small\mathbf{E} \sim \mathcal{N}(0,\textbf{V})\)

\(\textbf{V}\) is the phylogenetic covariance matrix

Describes expected covariance among taxa due to shared evolutionary history (typically under BM)

Coefficients found as: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T} \mathbf{V}^{-1} \mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{V}^{-1}\mathbf{Y}\right )\)

see Grafen. Phil. Trans. Roy. Soc. (1989)

PGLS: Example

bm.gls<-gls(Y~X, correlation=corBrownian(phy=tree,form = ~species), data=df)
anova(bm.gls)
## Denom. DF: 6 
##             numDF  F-value p-value
## (Intercept)     1 12.87792  0.0115
## X               1 19.28544  0.0046
coef(bm.gls)  
## (Intercept)           X 
##  -2.3296557   0.8846971

Identical results to PICs!

Statistical Digression: OLS vs GLS

To understand what PGLS is doing, consider the following methods for obtaining model parameters

OLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T}\mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{Y}\right )\)

Statistical Digression: OLS vs GLS

To understand what PGLS is doing, consider the following methods for obtaining model parameters

OLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T}\mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{Y}\right )\)

Unweighted GLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T} \mathbf{V}_{id}^{-1} \mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{V}_{id}^{-1}\mathbf{Y}\right )\) where

\[\tiny\mathbf{V}_{id}= \left( \begin{array}{ccc} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{array} \right) \]

Statistical Digression: OLS vs GLS

To understand what PGLS is doing, consider the following methods for obtaining model parameters

OLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T}\mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{Y}\right )\)

Unweighted GLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T} \mathbf{V}_{id}^{-1} \mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{V}_{id}^{-1}\mathbf{Y}\right )\) where

\[\tiny\mathbf{V}_{id}= \left( \begin{array}{ccc} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{array} \right) \]

Weighted GLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T} \mathbf{V}_{phy}^{-1} \mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{V}_{phy}^{-1}\mathbf{Y}\right )\) where

\[\tiny\mathbf{V}_{id}= \left( \begin{array}{ccc} (v_1+v_2) & v_{12} & 0 \\ v_{12} & (v_1+v_2) & 0 \\ 0 & 0 & 1 \end{array} \right) \]

In PGLS, the weights are the phylogenetic distances, which describe the phylogenetic non-independence

Statistical Digression: OLS vs GLS

To understand what PGLS is doing, consider the following methods for obtaining model parameters

OLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T}\mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{Y}\right )\)

Unweighted GLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T} \mathbf{V}_{id}^{-1} \mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{V}_{id}^{-1}\mathbf{Y}\right )\) where

\[\tiny\mathbf{V}_{id}= \left( \begin{array}{ccc} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{array} \right) \]

Weighted GLS model: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T} \mathbf{V}_{phy}^{-1} \mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{V}_{phy}^{-1}\mathbf{Y}\right )\) where

\[\tiny\mathbf{V}_{id}= \left( \begin{array}{ccc} (v_1+v_2) & v_{12} & 0 \\ v_{12} & (v_1+v_2) & 0 \\ 0 & 0 & 1 \end{array} \right) \]

In PGLS, the weights are the phylogenetic distances, which describe the phylogenetic non-independence

Phylogenetic Transformation

Utilize OLS transformation of GLS models

Phylogenetic transformation matrix: \(\small\mathbf{T}=\left ( \mathbf{U}\mathbf{W}^{-1/2} \mathbf{U}^{T}\right )^{-1}\)

U and W are eigenvectors and eigenvalues of V

Transformed data: \(\small\tilde{\mathbf{Y}}=\mathbf{TY}\)

Transformed design matrix: \(\small\tilde{\mathbf{X}}=\mathbf{TX}\)

Coefficients solved as: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{\tilde{X}}^{T} \mathbf{\tilde{X}}\right )^{-1}\left ( \mathbf{\tilde{X}}^{T} \mathbf{\tilde{Y}}\right )\)

see Garland and Ives. Am. Nat. (2000); Adams. Evol. (2014); Adams and Collyer. Evol. (2018)

Phylogenetic Transformation: Example

pgls.res<-lm.rrpp(Y~X, Cov = df$Cov,data=df, print.progress = FALSE) 
anova(pgls.res)
## 
## Analysis of Variance, using Residual Randomization
## Permutation procedure: Randomization of null model residuals 
## Number of permutations: 1000 
## Estimation method: Generalized Least-Squares (via OLS projection) 
## Sums of Squares and Cross-products: Type I 
## Effect sizes (Z) based on F distributions
## 
##           Df      SS      MS     Rsq      F      Z Pr(>F)   
## X          1 14.3519 14.3519 0.76271 19.285 2.9475  0.002 **
## Residuals  6  4.4651  0.7442 0.23729                        
## Total      7 18.8170                                        
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Call: lm.rrpp(f1 = Y ~ X, data = df, Cov = df$Cov, print.progress = FALSE)
pgls.res$LM$gls.coefficients
## (Intercept)           X 
##  -2.3296557   0.8846971

Identical results to PICs & PGLS!

Phylogenetic Regression Via LMM

One can also use linear mixed models (PGLMM: phylogenetic generalized linear mixed models) to evaluate trait associations

Here the phylogeny is considered a random effect (sense Lynch 1991)

Results from MCMC runs are analogous (though not always identical) to the linear model solutions

pgls.res$LM$gls.coefficients
## (Intercept)           X 
##  -2.3296557   0.8846971
summary(model_simple)$solutions
##              post.mean   l-95% CI u-95% CI eff.samp      pMCMC
## (Intercept) -1.4721864 -5.1694902 1.813651       98 0.38775510
## X            0.8097609  0.4361165 1.224440       98 0.01020408

Assessing Significance

PIC, PGLS, and Phylogenetic transform yield identical parameter estimates

\[\tiny \hat{\mathbf{\beta }}=\left ( \mathbf{X}^{T} \mathbf{V}^{-1} \mathbf{X}\right )^{-1}\left ( \mathbf{X}^{T} \mathbf{V}^{-1}\mathbf{Y}\right )\]

\[\tiny =\left ( \mathbf{X}^{T}_{pic} \mathbf{X}_{pic}\right )^{-1}\left ( \mathbf{X}^{T}_{pic} \mathbf{Y}_{pic}\right )\]

\[\tiny =\left ( \mathbf{\tilde{X}}^{T} \mathbf{\tilde{X}}\right )^{-1}\left ( \mathbf{\tilde{X}}^{T} \mathbf{\tilde{Y}}\right )\]

\[\tiny \log{L}=\log{ \left(\frac{exp{\left({-\frac{1}{2}{\left({\mathbf{Y}-E(\mathbf{Y})}\right)^{T}}} \mathbf{V}^{-1}\left({\mathbf{Y}-E(\mathbf{Y})}\right)\right)}} {\sqrt{\left({2\pi}\right)^{Np}\times\begin{vmatrix} \mathbf{V} \end{vmatrix}}}\right)}\]

Assessing Significance

PROBLEM: Parametric significance testing suffers from Rao’s paradox

Power reduces as p-dimensions increase

Another solution required

Adams. Evol. (2014); Adams and Collyer. Syst. Biol. (2018)

Assessing Significance: RRPP

Transform data: \(\small\tilde{\mathbf{Y}}=\mathbf{TY}\), \(\small\tilde{\mathbf{X}}_{F}=\mathbf{TX}_{F}\), \(\small\tilde{\mathbf{X}}_{R}=\mathbf{TX}_{R}\)

Obtain Parameter estimates: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{\tilde{X_{F}}}^{T} \mathbf{\tilde{X_{F}}}\right )^{-1}\left ( \mathbf{\tilde{X_{F}}}^{T} \mathbf{\tilde{Y}}\right )\)

Permute residuals from reduced model (\(\small\tilde{\mathbf{E}}_{R}\)), add to predicted values (\(\small\hat{\mathbf{Y}}_{R}\)), assess model, repeat

Assessing Significance: RRPP

Transform data: \(\small\tilde{\mathbf{Y}}=\mathbf{TY}\), \(\small\tilde{\mathbf{X}}_{F}=\mathbf{TX}_{F}\), \(\small\tilde{\mathbf{X}}_{R}=\mathbf{TX}_{R}\)

Obtain Parameter estimates: \(\small\hat{\mathbf{\beta }}=\left ( \mathbf{\tilde{X_{F}}}^{T} \mathbf{\tilde{X_{F}}}\right )^{-1}\left ( \mathbf{\tilde{X_{F}}}^{T} \mathbf{\tilde{Y}}\right )\)

Permute residuals from reduced model (\(\small\tilde{\mathbf{E}}_{R}\)), add to predicted values (\(\small\hat{\mathbf{Y}}_{R}\)), assess model, repeat

Appropriate for ANOVA, regression and more complex PGLS models with high-dimensional data (Adams & Collyer. Evol. 2018)

RULE: Transform first, permute second!

For discussion on permutation in PGLS see: Adams and Collyer. Evol. (2015); Adams and Collyer. Evol. (2018)

Phylogenetic Regression: Example

Are body dimensions associated with overall size across Plethodon salamander species?

Phylogenetic Regression: Example

Are body dimensions associated with overall size across Plethodon salamander species?

PGLS.reg<-procD.pgls(Y ~ size, phy = phy, data = gdf, iter = 999, print.progress = FALSE)
PGLS.reg$aov.table
##           Df     SS      MS     Rsq      F      Z Pr(>F)   
## size       1 214.21 214.214 0.35852 19.561 3.0528  0.001 **
## Residuals 35 383.29  10.951 0.64148                        
## Total     36 597.50                                        
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Phylogenetic ANOVA: Example

Do body dimensions differ by salamander ecotype?

PGLS.reg<-procD.pgls(Y ~ gp, phy = phy, data = gdf, iter = 999, print.progress = FALSE)
PGLS.reg$aov.table
##           Df     SS     MS     Rsq      F       Z Pr(>F)
## gp         1  23.44 23.435 0.03922 1.4288 0.77209  0.236
## Residuals 35 574.07 16.402 0.96078                      
## Total     36 597.50

Phylogenetic ANOVA: Group Aggregation

How groups are distributed on the phylogeny can make a difference: Adams and Collyer. Evol. (2018)

Too few group ‘state’ changes across phylogeny results in lower power

2: Phylogenetic PLS

Assess correlation between two blocks of variables in a phylogenetic context

PLS of evolutionary covariance (rate) matrix

\[\small\mathbf{R}=\frac{\left(\mathbf{Y}-E(\mathbf{Y})\right)^{T}\mathbf{V}^{-1} \left(\mathbf{Y}-E(\mathbf{Y})\right)} {N-1}\]

SVD of \(\small\mathbf{R}_{12}=\mathbf{U}_{12}\mathbf{D}\mathbf{V}^{T}_{12}\)

2: Phylogenetic PLS

Assess correlation between two blocks of variables in a phylogenetic context

PLS of evolutionary covariance (rate) matrix

\[\small\mathbf{R}=\frac{\left(\mathbf{Y}-E(\mathbf{Y})\right)^{T}\mathbf{V}^{-1} \left(\mathbf{Y}-E(\mathbf{Y})\right)} {N-1}\]

SVD of \(\small\mathbf{R}_{12}=\mathbf{U}_{12}\mathbf{D}\mathbf{V}^{T}_{12}\)

Equivalently found from PLS of \(\small\Sigma\) from \(\small\tilde{\mathbf{Y}}=\mathbf{TY}\)

Significance from RRPP of \(\small\tilde{\mathbf{Y}}\)

Adams and Felice. PLoS One (2014); Adams and Collyer. Syst. Biol. (2018)
Use of PICs equivalent (Klingenberg and Marugan-Loban. Syst. Biol. 2013)

Phylogenetic PLS: Example

Do head traits correlate with body traits (same data)?

PLS.Y <- phylo.integration(A = Y[,2:3], A2 = Y[,c(1,4:6)], phy= plethtree, print.progress = FALSE)
PLS.Y$r.pls; PLS.Y$P.value
## [1] 0.6133857
##   obs 
## 0.001

Phylogenetic PLS: Example 2

Is there evolutionary integration across species between cranial and mandible morphology?

Yes

Phylogenetic Ordination

  1. Phylomorphspace: project phylogeny into PCA space (using ancestral states)
    • Rotates data to directions of maximal variation
  2. Phylogenetic PCA (pPCA): account for phylogeny in PCA computations:
    • SVD of evolutionary rate matrix (the covariance matrix ‘standardized’ by phylogeny): \(\tiny\mathbf{R}=\frac{(\mathbf{Y}-E(\mathbf{Y}))^{T}\mathbf{C}^{-1}(\mathbf{Y}-E(\mathbf{Y}))}{N-1}\)
    • Rotates data such that pPC1 is independent of phylogenetic signal
  3. Phylogenetically-aligned component analysis (PACA)
    • Rotates data to directions that maximize phylogenetic signal (Collyer and Adams 2021)

Phylo-PCA: Example

Phylomorphospace: Examples

Visual inspection of phylomorphospaces useful for hypothesis generation

PACA

4: Phylogenetic Signal

The degree to which phenotypic similarity associates with phylogenetic relatedness

Phylogenetic Signal

Kappa statistic: Blomberg et al. Evol. (2003)

\[\small K=\frac{\left(\mathbf{Y}-E(\mathbf{Y})\right)^{T} \left(\mathbf{Y}-E(\mathbf{Y})\right)} {\left(\mathbf{Y}-E(\mathbf{Y})\right)^{T}\mathbf{V}^{-1} \left(\mathbf{Y}-E(\mathbf{Y})\right)} / \frac{tr(\mathbf{V})-N(\mathbf{1}^{T}\mathbf{V1})^{-1}}{N-1}\]

Univariate response data only!

Multivariate Phylogenetic Signal

\(\small{K}_{mult}\): A multivariate generalization of Kappa (Adams. Syst. Biol. 2014)

Based on distances between objects in original and phylo-transformed dataspaces

\[\small\ K_{mult}= \frac{\mathbf{D}_{\mathbf{Y},E(\mathbf{Y})}^{T}\mathbf{D}_{\mathbf{Y},E(\mathbf{Y})}} {\mathbf{PD}_{\tilde{\mathbf{Y}},0}^{T}\mathbf{PD}_{\tilde{\mathbf{Y}},0}} / \frac{tr(\mathbf{V})-N(\mathbf{1}^{T}\mathbf{V1})^{-1}}{N-1}\]

Multivariate Phylogenetic Signal

\(\small{K}_{mult}\): A multivariate generalization of Kappa (Adams. Syst. Biol. 2014)

Based on distances between objects in original and phylo-transformed dataspaces

\[\small\ K_{mult}= \frac{\mathbf{D}_{\mathbf{Y},E(\mathbf{Y})}^{T}\mathbf{D}_{\mathbf{Y},E(\mathbf{Y})}} {\mathbf{PD}_{\tilde{\mathbf{Y}},0}^{T}\mathbf{PD}_{\tilde{\mathbf{Y}},0}} / \frac{tr(\mathbf{V})-N(\mathbf{1}^{T}\mathbf{V1})^{-1}}{N-1}\]

Where \(\small\mathbf{D}_{\mathbf{Y},E(\mathbf{Y})}\) is the distance from each object to the phylogenetic mean, and \(\small\mathbf{PD}_{\tilde{\mathbf{Y}},0}\) is the distance between each object in the phylogenetically-transformed space (found as: \(\small\tilde{\mathbf{Y}}=\mathbf{P}(\mathbf{Y-E(\mathbf{Y})})\)) and the origin

Multivariate Phylogenetic Signal: Example

Do body dimensions in Plethodon display phylogenetic signal?

PS.shape <- physignal(A=Y,phy=plethtree,iter=999, print.progress = FALSE)
PS.shape
## 
## Call:
## physignal(A = Y, phy = plethtree, iter = 999, print.progress = FALSE) 
## 
## 
## 
## Observed Phylogenetic Signal (K): 0.4342
## 
## P-value: 0.002
## 
## Effect Size: 2.6861
## 
## Based on 1000 random permutations

Yes they do.

5: Comparing Evolutionary Models

What evolutionary model best describes trait variation?

Fit data to phylogeny under differing evolutionary models

5: Comparing Evolutionary Models

What evolutionary model best describes trait variation?

Fit data to phylogeny under differing evolutionary models

Model comparisons of:

1: Evolutionary rate models: BM1, BMM, etc.

2: Evolutionary ‘mode’ models: BM1, OU1, OUM, etc.

Evolutionary Rates

“How fast, as a matter of fact, do animals evolve in nature?” (Simpson, 1944)

Comparing Rates Among Traits

Is there evidence that one trait evolves faster than another?

Obtain \(\small\mathbf{R}_0\) and \(\small{logL}\):

\[\tiny\mathbf{R}=\frac{(\mathbf{Y}-E(\mathbf{Y}))^{T}\mathbf{C}^{-1}(\mathbf{Y}-E(\mathbf{Y}))}{N-1}=\left( \begin{array}{ccc} \sigma^2_{1} & \sigma^2_{1,2} & \sigma^2_{1,3} \\ \sigma^2_{2,1} & \sigma^2_{2} & \sigma^2_{2,3} \\ \sigma^2_{3,1} & \sigma^2_{3,2} & \sigma^2_{3} \end{array} \right)\]

\[\tiny\log{L}=-\frac{np}{2}ln{2\pi}-\frac{n}{2}ln{\begin{vmatrix}{(\mathbf{R\otimes{C}})}\end{vmatrix}}-\frac{1}{2}tr{(\mathbf{Y}-E(\mathbf{Y}))^{T}(\mathbf{R\otimes{C}})^{-1}(\mathbf{Y}-E(\mathbf{Y}))}\]

Comparing Rates Among Traits

Is there evidence that one trait evolves faster than another?

Obtain \(\small\mathbf{R}_0\) and \(\small{logL}\):

\[\tiny\mathbf{R}=\frac{(\mathbf{Y}-E(\mathbf{Y}))^{T}\mathbf{C}^{-1}(\mathbf{Y}-E(\mathbf{Y}))}{N-1}=\left( \begin{array}{ccc} \sigma^2_{1} & \sigma^2_{1,2} & \sigma^2_{1,3} \\ \sigma^2_{2,1} & \sigma^2_{2} & \sigma^2_{2,3} \\ \sigma^2_{3,1} & \sigma^2_{3,2} & \sigma^2_{3} \end{array} \right)\]

\[\tiny\log{L}=-\frac{np}{2}ln{2\pi}-\frac{n}{2}ln{\begin{vmatrix}{(\mathbf{R\otimes{C}})}\end{vmatrix}}-\frac{1}{2}tr{(\mathbf{Y}-E(\mathbf{Y}))^{T}(\mathbf{R\otimes{C}})^{-1}(\mathbf{Y}-E(\mathbf{Y}))}\]

Estimate \(\small\mathbf{R}_C\) and \(\small{logL}_C\) where diagonal of \(\mathbf{R}\) is contrained to be equal:

\[\tiny\mathbf{R}_C=\left( \begin{array}{ccc} \sigma^2_{p} & & \\ \sigma^2_{2,1} & \sigma^2_{p} & \\ \sigma^2_{3,1} & \sigma^2_{3,2} & \sigma^2_{p} \end{array} \right)\]

Compare \(\small{logL}_C\) against \(\small{logL}\)

Adams. Syst. Biol. (2013); Denton and Adams. Evol. (2015)

Compare Rates Among Traits: Example

Compare Rates Among Clades

Is there evidence for multiple evolutionary rates on the phylogeny?

Procedure: Fit data to phylogeny under single-rate and multi-rate models. Estimate \(\small\mathbf{R}\) and \(\small\log{L}\) & compare \(\small\log{L}\) between models (LRT, AIC, simulation, etc.)

Compare Rates Among Clades

Is there evidence for multiple evolutionary rates on the phylogeny?

Procedure: Fit data to phylogeny under single-rate and multi-rate models. Estimate \(\small\mathbf{R}\) and \(\small\log{L}\) & compare \(\small\log{L}\) between models (LRT, AIC, simulation, etc.)

Works well for univariate \(\mathbf{Y}\), but Type I error \(\small\uparrow\) with p-dimensions

\(\small\log{L_{mult}}\) not a general solution for high-dimensional data

see: Adams. Syst. Biol. (2014); Adams and Collyer. Syst. Biol. (2018); Adams and Collyer. Ann. Rev. Ecol. Evol. Syst. (2019)

Compare Multivariate Rates: \(\small\sigma^{2}_{mult}\)

Generalizes \(\small\sigma^{2}\) to estimate the net evolutionary rate in multivariate data

Procedure: Phylogenetic-transform data \(\small\tilde{\mathbf{Y}}=\mathbf{TY}\)

Estimate \(\small\sigma^{2}_{mult}=\frac{\mathbf{PD}_{\tilde{\mathbf{Y}},0}^{T}\mathbf{PD}_{\tilde{\mathbf{Y}},0}}{N}\) for regimes in BMM

Obtain rate ratio: \(\small{RR}=\frac{\sigma^{2}_{mult.max}}{\sigma^{2}_{mult.min}}\); evaluate via simulation or permutation (RRPP).

Displays appropriate type I error and power (breaks Rao’s paradox)

Adams. Syst. Biol. (2014); Denton and Adams. Evol. (2015)

Comparing Evolutionary ‘Modes’

Compare models that describe the manner in which trait variation accumulates over evolutionary time

BM1, BMM, OU (Ornstein-Uhlenbeck), Early Burst (EB), ACDC, etc.

Fit data to phylogeny under differing evolutionary models and compare fit (LRT, AIC, etc.)

Comparing Evolutionary ‘Modes’: Example

Comparing Multivariate Evolutionary ‘Modes’

Unfortunately, ALL current implementations for multivariate OU models display high misspecification rates

More work is needed in this area

see: Adams. Syst. Biol. (2014); Adams and Collyer. Syst. Biol. (2018); Adams and Collyer. Ann. Rev. Ecol. Evol. Syst. (2019)

Univariate PCMs: Conclusions

Multivariate PCMs: Conclusions