博文
什么是系统发育独立差（Phylogenetic Independent Contrast）？
||
张金龙、徐姗
系统发育独立差（Phylogenetic Independent Contrast, PIC）是去除性状分析中物种系统发育关系的一种方法，是美国进化生物学家 J. Felsenstein于1985年提出的。
举例来说，如果要研究哺乳动物脑容量和体重的关系，就要先获得这些种之间的进化关系，也就是经过分子钟校订的进化树。然后再计算脑容量的PIC和体重的PIC，再用脑容量的PIC和体重的PIC建立线性模型，进行相关性检验。
如果假设每个性状的进化服从布朗运动模型，也就是每个子代的性状的值是从祖先节点的性状随机变化而来的，那么进化树上相邻节点的两个种性状偏离于最近共同祖先性状的期望应该是0，且二者之间的差异正比于分化时间，也就是进化树的枝长。
假设相邻的分类单元为i和j，其性状分别为X_i和X_j，i和j的共同祖先是节点为k。
先计算X_i - X_j，作为X和Y的contrast。该contrast的期望为0，方差与枝长v_i+v_j成正比
从进化树中去掉节点i和j，保留节点k，将k当做末端节点。具体来说，就是以X_i和X_j的加权平均值作为X_k，公式如下：
4. 将k节点下的枝长从v_k延长为v_k + v_i*v_j/(v_i + v_j)。这样做的原因是，第3步中，k节点的X_k并非真实值，而是用X_i 、 X_j计算到的。
图 2. 拥有5个分类单元的进化树，分类单元名称(1-5为末端节点，6-8为内部节点)、枝长（V）
图3. PIC的计算公式
library(ape)
### The example in Phylip 3.5c (originally from Lynch 1991)
cat("((((Homo:0.21,Pongo:0.21):0.28,",
"Macaca:0.49):0.13,Ateles:0.62):0.38,Galago:1.00);",
file = "ex.tre", sep = "\n")
tree.primates <- read.tree("ex.tre")
X <- c(4.09434, 3.61092, 2.37024, 2.02815, -1.46968)
Y <- c(4.74493, 3.33220, 3.36730, 2.89037, 2.30259)
names(X) <- names(Y) <- c("Homo", "Pongo", "Macaca", "Ateles", "Galago")
pic.X <- pic(X, tree.primates)
pic.Y <- pic(Y, tree.primates)
cor.test(X, Y)#### Pearson's product-moment correlation #### data: X and Y ## t = 2.5736, df = 3, p-value = 0.08224 ## alternative hypothesis: true correlation is not equal to 0 ## 95 percent confidence interval: ## -0.1964310 0.9884173 ## sample estimates: ## cor ## 0.8296107
cor.test(pic.X, pic.Y)
#### Pearson's product-moment correlation #### data: pic.X and pic.Y ## t = -0.85623, df = 2, p-value = 0.4821 ## alternative hypothesis: true correlation is not equal to 0 ## 95 percent confidence interval: ## -0.9874751 0.8823934 ## sample estimates: ## cor ## -0.5179156
lm(pic.Y ~ pic.X - 1) # both regressions
#### Call: ## lm(formula = pic.Y ~ pic.X - 1) #### Coefficients: ## pic.X ## 0.4319
lm(pic.X ~ pic.Y - 1) # through the origin
#### Call: ## lm(formula = pic.X ~ pic.Y - 1) #### Coefficients: ## pic.Y ## 0.998
setwd("...")library(plantlist)## This is plantlist 0.6.1.
# devtools::install_github("helixcn/plantlist")
library(V.PhyloMaker)
# devtools::install_github("jinyizju/V.PhyloMaker")
library(picante)library(ape)
dat <- read.csv("trait_herb_114.csv", header = TRUE)
head(dat)## species height ## 1 Stipa_albasiensis 48.5 ## 2 Heteropappus_altaicus 15.6 ## 3 Roegneria_alashanica 43.7 ## 4 Saussurea_alaschanica 21.5 ## 5 Astragalus_alaschanus 18.0 ## 6 Panzerina_lanata 12.5
height <- dat$height
names(height) <- dat$species
# Build a phylogeny using phylo.maker
species <- gsub("_", " ", as.character(dat$species))
tab2 <- subset(TPL(gsub("_", " ", species)),
select = c("YOUR_SEARCH", "POSSIBLE_GENUS", "FAMILY"))
colnames(tab2) <- c("species", "genus", "family")
result1 <- phylo.maker(tab2, scenarios = c("S1"))## [1] "Taxonomic classification not consistent between sp.list and tree."## genus family_in_sp.list family_in_tree## 77 Hemerocallis Xanthorrhoeaceae Asphodelaceae
# write.tree(result1[[1]], "tree1.newick") # If you want to examine the tree tree <- multi2di(result1[[1]]) # plot(tree) # "If x has names, its values are matched to the tip labels of phy, otherwise its values are taken to be in the same order than the tip labels of phy." data_matched <- match.phylo.data(phy = tree, data = height ) # Phylogenetic Independent Contrastheight_pic <- pic(data_matched$data, data_matched$phy)
参考资料
Felsenstein, J. (1985) Phylogenies and the comparative method. American Naturalist, 125, 1–15. Jin, Y., & Qian, H. (2019). V. PhyloMaker: an R package that can generate very large phylogenies for vascular plants. Ecography, 42(8), 1353-1359. Kembel, S. W., Cowan, P. D., Helmus, M. R., Cornwell, W. K., Morlon, H., Ackerly, D. D., ... & Webb, C. O. (2010). Picante: R tools for integrating phylogenies and ecology. Bioinformatics, 26(11), 1463-1464. Kembel S. 编， 张金龙 译， 在R中分析物种和系统发育多样性 http://blog.sciencenet.cn/home.php?mod=space&uid=255662&do=blog&id=1116212 Paradis, E., & Schliep, K. (2019). ape 5.0: an environment for modern phylogenetics and evolutionary analyses in R. Bioinformatics, 35(3), 526-528. Zhang, J.L. (2018). plantlist: Looking Up the Status of Plant Scientific Names based on The Plant List Database. R package version 0.6.1. https://github.com/helixcn/plantlist/
https://blog.sciencenet.cn/blog-255662-1219388.html
上一篇：再忆铁峰
下一篇：用plantlist程序包查询和处理植物学名 (plantlist 0.6.5)
