全文链接:http://tecdat.cn/?p=31434
我们围绕进化树技术进行一些咨询,分析生物类群在时间上的多样性是如何变化的(点击文末“阅读原文”获取完整代码数据)。
我们将用到分类单元数-时间图(Lineages-through-time plot),该图可以用来描述物种多样化的总体趋势。
数据
3500trees.nexus是nexus格式的文件,里面有3500棵树。
besttree.nexus也是nexus格式的文件,里面有1颗树,是从3500颗树中筛选出来的一致树。
各支系图示
这棵树总共有4大支系(Lineage),现在我需要做的分析都是需要分别做总的,以及4个支系的,也就是说同样的分析要做5次,针对5组不同的对象。
分析方法
方法主要是物种多样化速率(diversification rate) 相关的内容。
trees=read.nexus("3500trees.nexus") besttree=read.nexus("besttree.nexus")
点击标题查阅往期内容
生态学建模:增强回归树(BRT)预测短鳍鳗生存分布和影响因素
01
02
03
04
1.mltt plot (multiple lineage through time)
分类单元数-时间图
lingeage的数值取log的,95%置信区间的ltt plot,中间黑色线的是besttree的,要显示出来。分别做总的,以及4个支系的,共5个图。
plot(trees)
,log='y')
besttree
# 95% ltt置信区间 ltt.ci<-function(tree.all){ ntip=length(tree.all[[1]]$tip.label) ntree=length(tree.all)
2.gamma statistic
检验分化速率的变化趋势,看γ的值是正的还是负的。结果需要得到每组的γ值及P值。
mmaStat(besttree) ## [1] -3.693285
3. Monte Carlo constant rates test
检验样品不全是否对分化速率的结果有显著的影响,应该也是每组都要做的。
mc.out <- mcmc.pop line(tree.hiv) plot(sk, l
4.对每个组做几个模型的检验,主要包括Pure-birth, birth-death, Yule 2-rate,density-dependent logistic,density-dependent exponential模型。
tdAICr ## --------------Model Summary---------------- ## ## MODEL pureBirth ## ## Parameters: r1 ## ## LH 535.1086 ## ## AIC -1068.217 ## ## r1 0.1817879 ## ## a -1068.217 ## ## ## -------------------------- ## MODEL bd ## ## Parameters: r1, a ## ## LH 535.1086 ## ## AIC -1066.217 ## ## r1 0.1817879 ## ## a 0 ## ## ## -------------------------- ## MODEL DDL ## ## Parameters: r1, k ## ## LH 542.2213 ## ## AIC -1080.443 ## ## r1 0.2537928 ## ## a -1080.443 ## ## k 554 ## ## ## -------------------------- ## MODEL DDX ## ## Parameters: r1, X ## ## LH 536.991 ## ## AIC -1069.982 ## ## r1 0.3098342 ## ## a -1069.982 ## ## x 0.1131752 ## ## ## -------------------------- ## ## Best Constant Rate Model = pureBirth AIC -1068.217 ## ## Best Rate Variable Model = DDL AIC -1080.443 ## ## delta AICrc = 12.2254 ## model params np mtype LH r1 r2 a xp k ## 1 pureBirth r1 1 RC 535.1086 0.1817879 NA -1068.217 NA NA ## 2 bd r1, a 2 RC 535.1086 0.1817879 NA 0.000 NA NA ## 3 DDL r1, k 2 RV 542.2213 0.2537928 NA -1080.443 NA 554 ## 4 DDX r1, X 2 RV 536.9910 0.3098342 NA -1069.982 0.113175