基因序列中的密碼子隨著時(shí)間的推移以幾乎恒定的比例進(jìn)行替換,具有恒定的演化速率,因此兩個(gè)物種的遺傳距離與物種分歧時(shí)間成正比榔至。序列比對(duì)時(shí),我們通常會(huì)利用單拷貝直系同源基因多序列比對(duì)文件用于構(gòu)建進(jìn)化樹欺劳。為了獲取更多的遺傳信息唧取,我們會(huì)將蛋白質(zhì)比對(duì)轉(zhuǎn)換成密碼子比對(duì),這時(shí)我們可以將這些密碼子比對(duì)文件串聯(lián)起來形成超級(jí)比對(duì)文件划提。而我們繪制時(shí)間樹枫弟,需要穩(wěn)定的核酸替換率,而三位密碼子所受到的選擇壓力各不相同鹏往,因此有必要使用 PartitionFinder2 將串聯(lián)而成的多序列比對(duì)文件進(jìn)行分區(qū)(可以分成三個(gè)部分)淡诗。將分區(qū)的比對(duì)文件可以分別計(jì)算進(jìn)化模型,RaxML伊履,MrBayes韩容,Beast,iqtree 都支持分區(qū)建樹唐瀑。而推導(dǎo)時(shí)間樹的程序 mcmctree 和 BEAST 都支持分區(qū)信息群凶。
PROGRAM | SUPPORT MULTITHREADING | INPUT FILE(Fossil calibration information) | TIME | EASY INSTALLATION |
---|---|---|---|---|
Mega | YES | A tree without branching information,Multi-sequence comparison file | About 1 hour | YES |
R8s | NO | A tree with branching information(root) | Within 10 minutes | NO |
MCMCtree | NO | A tree without branching information(root),Multi-sequence comparison file | More than four days | NO |
BEAST | YES | Multi-sequence comparison file | More than four days | YES |
一.mcmctree
三個(gè)輸入文件(最好是核酸比對(duì)文件):
- 1.序列比對(duì)文件
- 2.帶有校準(zhǔn)點(diǎn)的有根樹
- 3.配置文件mcmctree.ctl
1.序列比對(duì)文件
2.帶有校準(zhǔn)點(diǎn)的有根樹
- 1.獲取樹文件(樹文件通過核酸比對(duì)建立更加可靠)
- 2.去除樹文件中分支長度信息
- 3.訪問TimeTree網(wǎng)站或 古生物化石數(shù)據(jù)庫獲取物種分歧時(shí)間信息
TimeTree網(wǎng)站
Cl:置信區(qū)間(Confidence Interval)
MYA:百萬年以前(Million Years Ago)的縮寫
CI: (6.0 - 6.5 MYA)表示置信區(qū)間為從600萬年前到650萬年前的時(shí)間段。換句話說介褥,這個(gè)置信區(qū)間的長度是0.5百萬年座掘,也就是50萬年。
7 1
((((human, (chimpanzee, bonobo)) '>.06<.08', gorilla), (orangutan, sumatran)) '>.12<.16', gibbon);
這棵帶有校準(zhǔn)點(diǎn)的有根樹中時(shí)間單位是:一億年 --100個(gè)百萬年
上圖這棵樹中:人類(human)和黑猩猩(chimpanzee)之間的分歧時(shí)間在600萬年-800萬年之間
要注意的點(diǎn)是:樹在根上沒有化石校正點(diǎn)
3.配置文件mcmctree.ctl
使用PAML進(jìn)行分歧時(shí)間計(jì)算 | 陳連福的生信博客 (chenlianfu.com)
估算系統(tǒng)樹分歧時(shí)間 —— paml.mcmctree,r8s | 生信技工 (yanzhongsino.github.io)
【比較基因組】從進(jìn)化樹到分化時(shí)間 - 簡書 (jianshu.com)
***輸入輸出參數(shù):**
seed = -1 *設(shè)置一個(gè)隨機(jī)數(shù)來運(yùn)行程序柔滔。若設(shè)置為-1溢陪,表示使用系統(tǒng)當(dāng)前時(shí)間為隨機(jī)數(shù)。
seqfile = input.txt *設(shè)置輸入的多序列比對(duì)文件路徑
treefile = input.trees *設(shè)置輸入的帶校準(zhǔn)點(diǎn)信息有根樹的文件路徑
mcmcfile = mcmc.txt *設(shè)置輸出的mcmc信息文本文件睛廊,可以用Tracer軟件查看
outfile = out.txt *設(shè)置輸出文件路徑形真,該文件中記錄了超度量樹和進(jìn)化速率等參數(shù)信息。
***數(shù)據(jù)使用說明參數(shù):**
ndata = 3 *設(shè)置輸入的多序列比對(duì)的數(shù)據(jù)個(gè)數(shù)。這里指密碼子3個(gè)位置的數(shù)據(jù)咆霜。
seqtype = 0 *設(shè)置多序列比對(duì)的數(shù)據(jù)類型:0邓馒,表示核酸數(shù)據(jù);1蛾坯,表示密碼子比對(duì)數(shù)據(jù)光酣;2,表示氨基酸數(shù)據(jù)脉课。根據(jù)不同的數(shù)據(jù)類型救军,程序調(diào)用相應(yīng)的baseml或codeml進(jìn)行相應(yīng)的參數(shù)計(jì)算。
usedata = 1 *設(shè)置是否利用多序列比對(duì)的數(shù)據(jù):0倘零,表示不使用多序列比對(duì)數(shù)據(jù)唱遭,則不會(huì)進(jìn)行l(wèi)ikelihood計(jì)算,雖然能得到mcmc樹且計(jì)算速度飛快呈驶,但是其分歧時(shí)間結(jié)果是有問題的拷泽;1,表示使用多序列比對(duì)數(shù)據(jù)進(jìn)行l(wèi)ikelihood計(jì)算袖瞻,正常進(jìn)行MCMC司致,是一般使用的參數(shù); 2,進(jìn)行正常的approximation likelihood分析虏辫,此時(shí)不需要讀取多序列比對(duì)數(shù)據(jù)蚌吸,直接讀取當(dāng)前目錄中的in.BV文件锈拨。該文件是使用usedata = 3參數(shù)生成的out.BV文件重命名而來的砌庄。此外,由于程序BUG奕枢,當(dāng)設(shè)置usedata = 2時(shí)娄昆,一定要在改行參數(shù)后加 *,否則程序報(bào)錯(cuò) Error: file name empty.. 3缝彬,程序利用多序列比對(duì)數(shù)據(jù)調(diào)用baseml/codeml命令對(duì)數(shù)據(jù)進(jìn)行分析萌焰,生成out.BV文件。由于mcmctree調(diào)用baseml/codeml進(jìn)行計(jì)算的參數(shù)設(shè)置可能不太好(特別時(shí)對(duì)蛋白序列進(jìn)行計(jì)算時(shí))谷浅,推薦自己修該軟件自動(dòng)生成的baseml/codeml配置文件扒俯,然后再手動(dòng)運(yùn)行baseml/codeml命令,再整合其結(jié)果文件為out.BV文件一疯。
cleandata = 0 *設(shè)置是否移除不明確字符(N撼玄、?墩邀、W掌猛、R和Y等)或含以后gap的列后再進(jìn)行數(shù)據(jù)分析:0,不需要眉睹,但在序列兩兩比較的時(shí)候荔茬,還是會(huì)去除后進(jìn)行比較废膘;1,需要慕蔚。
clock = 2 *設(shè)置分子種方法:1丐黄,global clock方法,表示所有分支進(jìn)化速率一致孔飒;2孵稽,independent rates方法,各分支的進(jìn)化速率獨(dú)立且進(jìn)化速率的對(duì)數(shù)log(r)符合正態(tài)分布; 3十偶,correlated rates方法菩鲜,和方法2類似,但是log(r)的方差和時(shí)間t相關(guān)惦积。
* TipDate = 1 100 *當(dāng)外部節(jié)點(diǎn)由取樣時(shí)間時(shí)使用該參數(shù)進(jìn)行設(shè)置接校,同時(shí)該參數(shù)也設(shè)置了時(shí)間單位。具體數(shù)據(jù)示例請(qǐng)見examples/TipData文件夾狮崩。
RootAge = '<1.0' *設(shè)置root節(jié)點(diǎn)的分歧時(shí)間蛛勉,一般設(shè)置一個(gè)最大值。
***位點(diǎn)替換模型參數(shù):**
model = 4 *設(shè)置堿基替換模型:0睦柴,JC69诽凌;1,K80坦敌;2侣诵,F(xiàn)81;3狱窘,F(xiàn)84杜顺;4,HKY85蘸炸。當(dāng)設(shè)置usedata = 1時(shí)躬络,不能使用其它堿基替換模型。若需要選擇其它堿基替換模型搭儒,則考慮使用usedata = 3參數(shù)穷当,就可以設(shè)置model參數(shù)值為其它的堿基替換模型。此時(shí)淹禾,程序會(huì)將數(shù)據(jù)進(jìn)行分割馁菜,例如,ndata = 3稀拐,則程序分別對(duì)3個(gè)數(shù)據(jù)進(jìn)行分析火邓,生成baseml的輸入文件,并調(diào)用baseml進(jìn)行分析(若數(shù)據(jù)量較大,分析較慢铲咨,推薦按ctrl + c 強(qiáng)行終止躲胳,程序則終止baseml的運(yùn)行,進(jìn)行下一個(gè)數(shù)據(jù)的輸入文件生成并運(yùn)行下個(gè)數(shù)據(jù)的baseml命令纤勒,再按 ctrl + c 強(qiáng)行終止坯苹,這樣可以讓程序生成3各數(shù)據(jù)的輸入文件和baseml配置文件,然后手動(dòng)并行化運(yùn)行摇天,加快時(shí)間)粹湃;每個(gè)分析結(jié)果中會(huì)得到rst2結(jié)果文件;可以在此處手動(dòng)修改tmp0001.ctl配置文件泉坐,例如为鳄,修改其model參數(shù)值,推薦修改outfile參數(shù)值腕让,從而可以手動(dòng)并行化運(yùn)行baseml命令孤钦;然后將所有數(shù)據(jù)的rst2結(jié)果使用cat命令合并成文件in.BV,用于下一步設(shè)置參數(shù)usedata = 2纯丸,進(jìn)行MCMC分析偏形。
alpha = 0.5 *核酸序列中不同位點(diǎn),其進(jìn)化速率不一致觉鼻,其變異速率服從GAMMA分布俊扭。一般設(shè)置GAMMA分布的alpha值為0.5。若該參數(shù)值設(shè)置為0坠陈,則表示所有位點(diǎn)的進(jìn)化速率一致萨惑。此外,當(dāng)userdata = 2時(shí)畅姊,alpha咒钟、ncatG吹由、alpha_gamma若未、kappa_gamma這4個(gè)GAMMA參數(shù)無效。因?yàn)閡serdata = 2時(shí)倾鲫,不會(huì)利用到多序列比對(duì)的數(shù)據(jù)粗合。
ncatG = 5 *設(shè)置離散型GAMMA分布的categories值。
BDparas = 1 1 0.1 *設(shè)置出生率乌昔、死亡率和取樣比例隙疚。若輸入有根樹文件中的時(shí)間單位發(fā)生改變,則需要相應(yīng)修改出生率和死亡率的值磕道。例如供屉,時(shí)間單位由100Myr變換為1Myr,則要設(shè)置成".01 .01 0.1"。
kappa_gamma = 6 2 *設(shè)置kappa(轉(zhuǎn)換/顛換比率)的GAMMA分布參數(shù)伶丐。
alpha_gamma = 1 1 *設(shè)置GAMMA形狀參數(shù)alpha的GAMMA分布參數(shù)悼做。
***進(jìn)化速率參數(shù):**
rgene_gamma = 2 20 1 *設(shè)置序列中所所有位點(diǎn)平均[堿基/密碼子/氨基酸]替換率的Dirichlet-GAMMA分布參數(shù):alpha=2、beta=20哗魂、初始平均替換率為每100 million 年(取決于輸入有根樹文件中的時(shí)間單位)1個(gè)替換肛走。若時(shí)間單位由100Myr變換為1Myr,則要設(shè)置成"2 2000 1"录别⌒嗌總體上的平均進(jìn)化速率為:2 / 20 = 0.1 個(gè)替換 / 每100Myr,即每個(gè)位點(diǎn)每年的替換數(shù)為 1e-9组题。
sigma2_gamma = 1 10 1 *設(shè)置所有位點(diǎn)進(jìn)化速率取對(duì)數(shù)后方差(sigma的平方)的Dirichlet-GAMMA分布參數(shù):alpha=1葫男、beta=10、初始方差值為1崔列。當(dāng)clock參數(shù)值為1時(shí)腾誉,表示使用全局的進(jìn)化速率,各分枝的進(jìn)化速率沒有差異峻呕,即方差為0利职,該參數(shù)無效;當(dāng)clock參數(shù)值為2時(shí)瘦癌,若修改了時(shí)間單位猪贪,該參數(shù)值不需要改變;當(dāng)clock參數(shù)值為3時(shí)讯私,若修改了時(shí)間單位热押,該參數(shù)值需要改變。
finetune = 1: .1 .1 .1 .1 .1 .1 *冒號(hào)前的值設(shè)置是否自動(dòng)進(jìn)行finetune斤寇,一般設(shè)置成1桶癣,然程序自動(dòng)進(jìn)行優(yōu)化分析;冒號(hào)后面設(shè)置各個(gè)參數(shù)的步進(jìn)值:times, musigma2, rates, mixing, paras, FossilErr娘锁。由于有了自動(dòng)設(shè)置牙寞,該參數(shù)不像以前版本那么重要了,可能以后會(huì)取消該參數(shù)莫秆。
***MCMC參數(shù):**
print = 1 *設(shè)置打印mcmc的取樣信息:0间雀,不打印mcmc結(jié)果;1镊屎,打印除了分支進(jìn)化速率的其它信息(各內(nèi)部節(jié)點(diǎn)分歧時(shí)間惹挟、平均進(jìn)化速率、sigma2值)缝驳;2连锯,打印所有信息归苍。
burnin = 2000 *將前2000次迭代burnin后,再進(jìn)行取樣(即打印出該次迭代計(jì)算的結(jié)果信息运怖,各內(nèi)部節(jié)點(diǎn)分歧時(shí)間霜医、平均進(jìn)化速率、sigma2值和各分支進(jìn)化速率等)驳规。
sampfreq = 10 *每10次迭代則取樣一次肴敛。
nsample = 20000 *當(dāng)取樣次數(shù)達(dá)到該次數(shù)時(shí),則取樣結(jié)束(程序也將運(yùn)行結(jié)束)吗购。至少要一萬個(gè)
二.MAGA
- CLOCKS -> Compute Timetree -> RelTime-ML
- LOAD SEQUENCE DATA : alignment data (fasta format)
- LOAD TREE FILE : tree data (nwk format,no branch)
- Select Taxa : add outgroup data(click OK)
- Select Branch : (click a branch then click OK)
-
CALIBRATE NOTE -> (click clock icon ) -> SELECT CALIBRATE TYPE(Normal Distribution) -> (click OK)
可以多線程運(yùn)行趴俘,但是由于電腦的內(nèi)存不足可能無法運(yùn)行(8條序列蚤吹,1101462 * 8個(gè)堿基,54GB內(nèi)存)
也可以下載運(yùn)行l(wèi)inux版本的maga。只不過下載的都是預(yù)編譯的版本令野,需要系統(tǒng)內(nèi)部特定版本的glibc 纳击,為此還得專門去編譯 glibc俏讹。
三.BEAST
- 將多序列比對(duì)文件(fasta文件)轉(zhuǎn)換成nexus文件(可以使用mega的另存為)
- 使用BEAST下面的子程序 beauti (可視化頁面)將nex文件轉(zhuǎn)換成xml文件蝴罪。
# 導(dǎo)入beast 所需要的庫文件,否則不能正常運(yùn)行
git clone https://github.com/beagle-dev/beagle-lib.git
cd beagle-lib/
mkdir build
cd build/
cmake -DCMAKE_INSTALL_PREFIX=/home/bqxiao/program/beagle ..
make
make install
export LD_LIBRARY_PATH=/home/bqxiao/program/beagle/lib:$LD_LIBRARY_PATH
3.1 BEAST 1
BEAST 運(yùn)行的本身命令比較簡單:
nohup beast file.xml -threads N &
# beast(linux版) 本身是一個(gè)shell腳本埠偿,通過調(diào)用java來運(yùn)行beast.jar程序透罢,在beast腳本內(nèi)部,可以設(shè)置給java虛擬機(jī)的內(nèi)存大小冠蒋。當(dāng)選擇核算替換模型和線程數(shù)時(shí)羽圃,應(yīng)當(dāng)考慮堆內(nèi)存大小。
beauti 主要參數(shù):
patitions : 在這個(gè)頁面顯示導(dǎo)入的數(shù)據(jù)(xml抖剿,nex)的分區(qū)信息朽寞。在這里可以設(shè)置分區(qū)之間是否使用相同的模型,先驗(yàn)樹等信息
Taxa : 添加校準(zhǔn)節(jié)點(diǎn)斩郎,可以在序列數(shù)據(jù)中設(shè)置分類單元子集(相同子集拖入 Included Taxa 內(nèi))脑融,在后續(xù)可以添加校準(zhǔn)時(shí)間
Tip : 為每一條序列添加注釋(如:采集時(shí)間,采集地點(diǎn)等)
Sites : 進(jìn)化模型選擇
Clocks : 分子鐘模型設(shè)置
Tree : 種群或者物種規(guī)模隨時(shí)間的變化模型缩宜,主要考慮序列來自物種內(nèi)還是物種間肘迎。
Priors : 為前面的所有模型添加先驗(yàn)信息,在這里可以設(shè)置前面的 Taxa 校準(zhǔn)時(shí)間脓恕。
MCMC :
- Length of chain :MCMC分析運(yùn)行的總代數(shù)膜宋。可以通過以增加代數(shù)來提高ESS值(Effective Sample Size)使參數(shù)收斂炼幔。增加采樣頻率也可以提高ESS值,但是過于采樣頻繁也不會(huì)提高ESS值史简。
- Log parameters every : 采樣頻率乃秀,即多少代生成一個(gè)樣本肛著。這個(gè)值的設(shè)置受Length of chain影響,最主要的一點(diǎn)需要保證最后的樣本最少得有10,000個(gè)跺讯。
-
Echo state to screen : 將狀態(tài)摘要輸出到BEAST控制臺(tái)窗口的頻率枢贿。
最后點(diǎn)擊生成xml文件
3.2 BEAST 2
BEAST 2中beast的使用方法和BEAST 1 中的beast使用方法一致。而子程序 beauti 有些許偏差刀脏。
四.r8s
r8s 分析速度非尘旨裕快,對(duì)于一棵樹只需要幾秒鐘可以得出結(jié)果愈污。其中 pyr8s 是 r8s 的重寫版本耀态。其中 r8s 中的很多參數(shù)的(模型)算法只有部分在 pyr8s 得到的。
4.1 r8s v1.8.1
原本的r8s目前已經(jīng)不支持了暂雹,僅有 r8s download | SourceForge.net 這個(gè)網(wǎng)站還有源碼下載首装,但是依然有一些問題。
比如在編譯 fortran 一些文件時(shí)杭跪,有些代碼過于古老仙逻,特性已被刪除,不能通過編譯涧尿。
因此可以將make報(bào)錯(cuò)部分手動(dòng)編譯
# 原始make部分
gfortran -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /home/bqxiao/miniconda3/include -c -o tn.o tn.f
# 手動(dòng)編譯系奉,保留已刪除的特性: -std=legacy
# 其他標(biāo)準(zhǔn)還有 -std=f95 -std=f2003 -std=f2008
gfortran -std=legacy -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /home/bqxiao/miniconda3/include -c -o tn.o tn.f
另一部分編譯報(bào)錯(cuò)于全局變量重定義 -> 只需要在重定義部分前置 extern
修飾即可。
使用 :
r8s -b -f file.nex > file.out
# -b 使用命令行而不是交互式
配置文件如下:
#NEXUS
begin trees;
tree PAUP_9 = (Marchantia:0.033817,(Lycopodium:0.040281,((Equisetum:0.048533,(Osmunda:0.033640,Asplenium:0.036526):0.000425):0.011806,((((Cycas:0.009460,Zamia:0.018847):0.005021,Ginkgo:0.014702):1.687e-86,((Pinus:0.021500,(Podocarpac:0.015649,Taxus:0.021081):0.006473):0.002448,(Ephedra:0.029965,(Welwitsch:0.011298,Gnetum:0.014165):0.006883):0.016663):0.006309):0.010855,((Nymphaea:0.016835,(((((Saururus:0.019902,Chloranth:0.020151):1.687e-86,((Araceae:0.020003,(Palmae:0.006005,Oryza:0.031555):0.002933):0.007654,Acorus:0.038488):0.007844):1.777e-83,(Calycanth:0.013524,Lauraceae:0.035902):0.004656):1.687e-86,((Magnolia:0.015119,Drimys:0.010172):0.005117,(Ranunculus:0.029027,((Nelumbo:0.006180,Platanus:0.002347):0.003958,(Buxaceae:0.013294,((Pisum:0.035675,(Fagus:0.009848,Carya:0.008236):0.001459):0.001994,(Ericaceae:0.019136,Solanaceae:0.041396):0.002619):1.687e-86):0.004803):1.687e-86):0.006457):0.002918):0.007348,Austrobail:0.019265):1.687e-86):1.687e-86,Amborella:0.019263):0.003527):0.021625):0.012469):0.019372);
End;
begin rates;
blformat nsites=952 lengths=persite;
collapse;
mrca LP marchantia pisum;
mrca ANGIO amborella pisum;
mrca VP lycopodium pisum;
mrca V1 Equisetum pisum;
mrca FERN Osmunda Asplenium;
mrca SP Ginkgo pisum;
mrca GYMNO Ginkgo Gnetum;
mrca CYC Cycas Zamia;
mrca G2 Taxus Podocarpac;
mrca GNET Ephedra Gnetum;
mrca CL Calycanth Lauraceae;
mrca ChS Chloranth Saururus;
mrca MON Acorus Oryza;
mrca PO Palmae Oryza;
mrca MAG Magnolia Drimys;
mrca EUD Ranunculus Carya;
mrca NEL Nelumbo Platanus;
mrca BX Buxaceae Carya;
mrca ERIC Ericaceae Solanaceae;
mrca FC Fagus Carya;
fixage taxon=LP age=450;
constrain taxon=CL min_age=110;
constrain taxon=ChS min_age=125;
constrain taxon=MON min_age=120;
constrain taxon=PO min_age=85;
constrain taxon=MAG min_age=125;
constrain taxon=NEL min_age=110;
constrain taxon=BX min_age=105;
constrain taxon=ERIC min_age=91;
constrain taxon=FC min_age=96;
constrain taxon=VP min_age=420;
constrain taxon=V1 min_age=390;
constrain taxon=FERN min_age=255;
constrain taxon=GYMNO min_age=310;
constrain taxon=CYC min_age=220;
constrain taxon=GNET min_age=115;
constrain taxon=G2 min_age=210;
constrain taxon=EUD max_age=125;
constrain taxon=SP max_age=330;
unfixage taxon=lp;
set num_time_guesses=5;
unfixage taxon=platanus;
constrain taxon=platanus maxage=90;
divtime method=LF algorithm=tn;
showage;
describe plot=chronogram;
describe plot=tree_description;
end;
4.2 pyr8s
安裝非常簡單
pip install git+https://github.com/iTaxoTools/pyr8s.git
4.2.1 python 內(nèi)部使用
pyr8s 原文檔 使用方法有問題
import pyr8s.parse as p
a = p.parse('tests/legacy_1')
# 其中 parse.py 中根本沒有 parse 函數(shù)
通過查看 run.py 可以發(fā)現(xiàn) 姑廉,應(yīng)當(dāng)是以下代碼:
import pyr8s.parse as p
a=p.from_file("tests/legacy_1")
result=a.run()
4.2.2 shell 環(huán)境中使用
pyr8s file.nex > result.txt
配置文件如下:
#NEXUS
begin trees;
tree Apiaceae = [&U] (((Angelica_sinensis:0.0263287636,((Apium_graveolens:0.0460894249,Heracleum_sosnowskyi:0.0322190162):0.0024624336,(Coriandrum_sativum:0.0283896098,Peucedanum_praeruptorum_Dunn:0.0307326278):0.0032486520):0.0020207251):0.0171794647,Oenanthe_sinensis:0.0537308346):0.0110179388,Daucus_carota:0.0638858678,cacao:0.6170394091);
End;
begin rates;
[
nsites : 核苷酸的個(gè)數(shù)
persite : 最大似然樹
]
blformat nsites=1101462 lengths=persite;
collapse;
mrca DO Daucus_carota Oenanthe_sinensis;
mrca CH Coriandrum_sativum Heracleum_sosnowskyi;
constrain taxon=DO min_age=24 max_age=44;
constrain taxon=CH min_age=5 max_age=40;
set num_time_guesses=5;
divtime method=NP algorithm=PL;
[
showage : 展示各個(gè)分支的節(jié)點(diǎn)的年齡
]
showage;
[
describe plot=chronogram : 在終端畫出樹的形狀
]
describe plot=chronogram;
[
describe plot=tree_description : 在終端導(dǎo)出這棵樹的nwk格式
]
describe plot=tree_description;
end;