Skip to content

Instantly share code, notes, and snippets.

@davharris

davharris/geiger.R

Created January 30, 2015 00:52
Show Gist options
  • Save davharris/d268c42eb9502542f7c7 to your computer and use it in GitHub Desktop.
Save davharris/d268c42eb9502542f7c7 to your computer and use it in GitHub Desktop.
Download ZIP
geiger versus multivariate normal
Raw
geiger.R
# Confirming that I understand how geiger's variance-covariance matrices work
library(geiger)
library(mvtnorm)
geo = get(data("geospiza"))
# Simulate one trait using geiger::sim.char
sim = t(
sim.char(
drop.tip(geo$phy, "olivacea"),
matrix(1),
nsim = 1E5,
root = 0
)[ , 1, ]
)
# Simulate one trait using multivariate normal
sim2 = rmvnorm(1E5, sigma = vcv(drop.tip(geo$phy, "olivacea")))
# Simulate again to assess sampling error
sim3 = rmvnorm(1E5, sigma = vcv(drop.tip(geo$phy, "olivacea")))
plot(
round(cov(sim), 3),
round(cov(sim2), 3)
)
abline(0,1)
plot(
round(cov(sim2), 3),
round(cov(sim3), 3)
)
abline(0,1)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment