3.5: Multivariate Brownian motion

The Brownian motion model we described above was for a single character. However, we often want to consider more than one character at once. This requires the use of multivariate models. The situation is more complex than the univariate case – but not much! In this section I will derive the expectation for a set of (potentially correlated) traits evolving together under a multivariate Brownian motion model.

Character values across species can covary because of phylogenetic relationships, because different characters tend to evolve together, or both. Fortunately, we can generalize the model described above to deal with both of these types of covariation. To do this, we must combine two variance-covariance matrices. The first one, C, we have already seen; it describes the variances and covariances across species for single traits due to shared evolutionary history along the branches of a phylogentic tree. The second variance-covariance matrix, which we can call R, describes the variances and covariances across traits due to their tendencies to evolve together. For example, if a species of lizard gets larger due to the action of natural selection, then many of its other traits, like head and limb size, will get larger too due to allometry. The diagonal entries of the matrix R will provide our estimates of σi2, the net rate of evolution, for each trait, while off-diagonal elements, σij, represent evolutionary covariances between pairs of traits. We will denote number of species as n and the number of traits as m, so that C is n × n and R is m × m.

Our multivariate model of evolution has parameters that can be described by an m × 1 vector, a, containing the starting values for each trait – $\bar{z}_1(0)$, $\bar{z}_2(0)$, and so on, up to $\bar{z}_m(0)$, and an m × m matrix, R, described above. This model has m parameters for a and m ⋅ (m + 1)/2 parameters for R, for a total of m ⋅ (m + 3)/2 parameters.

Under our multivariate Brownian motion model, the joint distribution of all traits across all species still follows a multivariate normal distribution. We find the variance-covariance matrix that describes all characters across all species by combining the two matrices R and C into a single large matrix using the Kroeneker product:

$\textbf{V} = \textbf{R} ⊗ \textbf{C} \label{3.23}$

This matrix V is n ⋅ m × n ⋅ m, and describes the variances and covariances of all traits across all species.

We can return to our example of evolution along a single branch (Figure 3.4a). Imagine that we have two characters that are evolving under a multivariate Brownian motion model. We state the parameters of the model as:

$$\begin{array}{lcr} \mathbf{a} = \begin{bmatrix} \bar{z}_1(0) \\ \bar{z}_2(0) \\ \end{bmatrix} \\ \mathbf{R} = \begin{bmatrix} \sigma_1^2 & \sigma_{12} \\ \sigma_{12} & \sigma_2^2 \\ \end{bmatrix} \\ \end{array} \label{3.24}$$

For a single branch, C = [t1], so:

$$\mathbf{V} = \mathbf{R} \otimes \mathbf{C} = \begin{bmatrix} \sigma_1^2 & \sigma_{12} \\ \sigma_{12} & \sigma_2^2 \\ \end{bmatrix} \otimes [t_1] = \begin{bmatrix} \sigma_1^2 t_1 & \sigma_{12} t_1 \\ \sigma_{12} t_1 & \sigma_2^2 t_1 \\ \end{bmatrix} \label{3.25}$$

The two traits follow a multivariate normal distribution with mean a and variance-covariance matrix V.

For the simple tree in figure 3.4b,

 \begin{align} \mathbf{V} &= \mathbf{R} \otimes \mathbf{C} \5pt] &= \begin{bmatrix} \sigma_1^2 & \sigma_{12} \\ \sigma_{12} & \sigma_2^2 \\ \end{bmatrix} \otimes \begin{bmatrix} t_1+t_2 & t_1 \\ t_1 & t_1+t_3 \\ \end{bmatrix} \\[5pt] &= \begin{bmatrix} \sigma_1^2 (t_1+t_2) & \sigma_{12} (t_1+t_2) & \sigma_1^2 t_1 & \sigma_{12} t_1 \\ \sigma_{12} (t_1+t_2) & \sigma_2^2 (t_1+t_2) & \sigma_{12} t_1 & \sigma_2^2 t_1 \\ \sigma_1^2 t_1 & \sigma_{12} t_1 & \sigma_1^2 (t_1+t_3) & \sigma_{12} (t_1+t_3) \\ \sigma_{12} t_1 & \sigma_2^2 t_1 & \sigma_{12} (t_1+t_3) & \sigma_2^2 (t_1+t_3) \\ \end{bmatrix} \\ \end{align} \label{3.26} Thus, the four trait values (two traits for two species) are drawn from a multivariate normal distribution with mean \[a=[\bar{z}_1(0), \bar{z}_1(0), \bar{z}_2(0), \bar{z}_2(0)]

and the variance-covariance matrix shown above.

Both univariate and multivariate Brownian motion models result in traits that follow multivariate normal distributions. This is statistically convenient, and in part explains the popularity of Brownian models in comparative biology.