modeltest用法.docx
- 文档编号:18138032
- 上传时间:2023-08-13
- 格式:DOCX
- 页数:18
- 大小:33.13KB
modeltest用法.docx
《modeltest用法.docx》由会员分享,可在线阅读,更多相关《modeltest用法.docx(18页珍藏版)》请在冰点文库上搜索。
modeltest用法
Medeltest3.7:
准备工作:
下载modeltest软件windows版本,解压到系统盘,这里为C盘,我们这里将其置于根目录。
从http:
//www.rhizobia.co.nz/phylogenetics/modeltest.html该网站下载文件modelblockPAUPb10.txt和文件ML-search.txt。
下载Paup4.0
1.打开paup,打开你自己的filename.nex文件,然后再打开并执行所下载的modelblockPAUPb10.txt文件。
Paup将运行并产生一个名为model.scores的文件。
2.将产生的model.scores文件与modeltest执行程序放置在相同的文件夹内,这里为C:
\modeltest3.7\Modeltest3.7folder\bin\modeltest3.7.win。
3.在开始-运行中输入cmd命令,确认,弹出一个dos界面的窗口,将其默认目录改为model.scores文件所在文件夹,修改方式为在dos窗口输入:
cd盘符:
\filename\filename,确认,这里我们输入
cdC:
\modeltest3.7\”Modeltest3.7folder”\bin\。
带有空格的文件夹名要用双引号括起来。
4.修改好后,在窗口内输入:
Modeltest3.7.win.exe
在model.scores所在文件夹内产生test.outfile文件,即我们所需要的文件。
文件test.outfile的文件名可以修改。
5.用记事本打开test.outfile文件,里面列出了两种检验标准LRT和AIC分别选出的最优DNA进化模型。
6.如果你要在paup中利用最大似然法建树,那么只需将其中一种标准的对应命令拷贝到下载的ML-serch.txt文件中标记的位置,然后用paup打开你的nex文件之后执行该文件即可。
或者将命令拷贝到你自己的ML批处理文件内执行。
命令既是beginpaup;和end;之间的部分。
7.如果是用于mrbayes分析,不能直接拷贝上述命令。
其中base即Basefrequencies,代表各碱基出现的频率,顺序为A-C-G-T,其中T的频率没有列出来,不过这些值都可以在test.outfile中找到。
在mrbayes中,对应于Statefreqpr选项,该选项默认状态为Dirichlet(1.0,1.0,1.0,1.0),依次为A-C-G-T的频率,即各碱基的出现频率相同。
在此需将该选项状态修改为fixed(),括号中依次填入base中对应的值,中间逗号隔开,如Statefreqpr=fixed(0.1649,0.3340,0.3209,0.1801)。
其中各值相加为1,但在mrbayes中,各值范围不定,各选项的值可同时乘或除以某值,如100,1000等,只要比例不变。
如在此可修改为Statefreqpr=fixed(1649,3340,3209,1801)。
其中rmat即ratematrix,是碱基位点的变异频率或者说异质率(rateheterogeneity),共6个值,R(a)[A-C],R(b)[A-G],R(c)[A-T],R(d)[C-G],R(e)[C-T],R(f)[G-T]。
rmat命令中只列出前五个,最后一个R(f)[G-T]未列出,不过这些值都可以在test.outfile中找到。
在mrbayes中,该命令对应于revmatpr选项,该选项默认状态为Revmatpr=Dirichlet(1.0,1.0,1.0,1.0,1.0,1.0),括号中对应于rmat的各值,顺序不变,依次为[A-C],[A-G],[A-T],[C-G],[C-T],[G-T],各数值间用逗号隔开;在此,我们需将该选项状态修改为fixed,即recmatpr=fixed(),括号中依次填入rmat中对应的各值,如revmatpr=fixed(1.0000,3.7410,1.0000,1.0000,2.0672,1.0000)。
其中pinvar即Proportionofinvariablesites,是不变位点的比率或称among-siteratevariation,在mrbayes中对应于Pinvarpr选项,该命令默认为uniform状态,即pinvarpr=uniform(0,1),取0-1内的一个值,0-1也是该选项的有效范围,在此,我们将该选项修改为fixed固定值选项,即pinvarpr=fixed(),括号内填入pinvar中对应的数值,如pinvarpr=fixed(0.3076)。
其中Shape即Gammadistributionshapeparameter,确定伽马分布的形状参数α,在mrbayes中对应于Shapepr选项,该选项默认状态为α范围内的均一分布,shapepr=uniform(0,α),在此修改为shapepr=fixed(1.5422)。
现在假设原命令是:
lsetBase=(0.16490.33400.3209)Nst=6Rmat=(1.00003.74101.00001.00002.0672)Rates=gammaShape=1.5422Pinvar=0.3076;
修改后对应的mrbayes命令为:
Lsetnst=6rates=gamma;
PrsetStatefreqpr=fixed(0.1649,0.3340,0.3209,0.1801)
revmatpr=fixed(1.0000,3.7410,1.0000,1.0000,2.0672,1.0000)shapepr=fixed(1.5422)pinvarpr=fixed(0.3076);
在statefreqpr和revmatpr中,fixed一般设置为dirichlet,不知道两者有没有明显的差别?
另,statefreqpr以及shapepr、pinvarpr参数在mrbayes中似乎没见别人用过,不晓得为何?
ModelTest检测(转)
2009-08-3110:
37
建树最首先的一个条件是,检测你的序列是否存在重排,是否替换已经饱和,如果是的话,说明你获得的这些基因的可靠度很有问题,用这些基因得到的进化树肯定是不可靠的。
那么确定基因可以用之后呢,最重要的就是使用ModelTest检测这些序列最适合的是那个Model。
现行的ModelTest是2005的版本,ModelTest3.7。
但是这个版本,是DOS环境下运行的,对于电脑水平一般的,使用起来不是很方便。
有人做了一个MrMTgui软件,可以把ModelTest转化成界面操作。
另外,ModelTest的运行,还需要PAUP来计算这组序列在所用Model下的似然值,然后ModelTest才能根据这些似然值来判断,那个Model是最合适的。
1、把序列使用ClustalX1.83进行比对,然后保存成PAUP识别的格式Nexus,然后导入至PAUP。
2、导入序列之后,在Modeltest的文件夹中,有一个paupblock的文件夹,里面有一个文件叫
“modelblockPAUPb10”,然后在paup软件中,file,open,就可以执行这个命令了。
就是paup用ML
计算不同model下的似然值。
如果序列比较多,且长的话,需要计算一阵呢!
3、运行完这步之后,在paupblock文件夹中,会生成一个名为“model.scores”的文件。
4、安装MrMTgui,完成之后,打开界面,在下面会有MTpath和PAUPpath两个按钮,然后打开,找到
Modeltest和paup所在的文件夹,双击他们的安装程序,就可以连接上两个软件了。
5、然后点击左边的按钮“selectfile”,找到“model.scores”,之后,点击按钮
“Modeltest!
!
”,然后就会出现运行结果了。
6、会出现两个运行结果,一个是hLRT得出的结果,另一个是AIC给出的结果。
[!
Likelihoodsettingsfrombest-fitmodel(TrN+I+G)selectedbyhLRTinModeltest3.7onFriJul2415:
23:
492009
]
BEGINPAUP;
LsetBase=(0.29240.27190.2306)Nst=6Rmat=(1.00002.50531.00001.00003.5293)Rates=gammaShape=1.1541Pinvar=0.0131;
END;
[!
Likelihoodsettingsfrombest-fitmodel(GTR+I+G)selectedbyAICinModeltest3.7onFriJul2415:
23:
492009
]
BEGINPAUP;
LsetBase=(0.27890.26370.2411)Nst=6Rmat=(2.15624.49091.97511.87026.0884)Rates=gammaShape=1.2557Pinvar=0.0132;
END;
Beginpaup后面的一段命令是在paup中运行,来执行这个模型的命令。
PAML软件的一些简单的具体的使用操作(转)
2009-08-3110:
45
PAML软件的一些简单的具体的使用操作(2008-11-14153642)标签:
杂谈
1.首先用ClustalX进行序列比对:
要保证:
保证核苷酸序列是三的倍数,没有终止密码子,核苷酸序列的第一位是密码子的第一位。
假设序列名为cox1.fas
2.使用DAMBE软件进行转换成PML格式。
打开要转换换的文件,然后“file”“saveandconvertsequenceformat”,在保存类型中选择“Yang’sPAML”。
那么此时的序列名为“cox1.PML”
3这样就可以得到文件“.PML”,然后就直接把后缀改成“.nuc”。
那么此时的序列名为“cox1.nuc”这样就完成了文件格式的转换。
(好像可以不用改)
4打开PAML软件的文件夹,找到文件名是“bin”的文件夹,打开之后,找到程序“codeml.exe”,把该程序复制到D盘的根目录下。
(这一步并不是必要的,只是要把用到的几个程序放在同一个目录下)
5在你使用ClustalX进行序列比对的时候,会生成一棵进化树,适用treeview软件可以打开,你需要的是把文件的后缀名改称“.trees”。
即树的文件名是“cox1.trees”,这就完成了树的格式的转换。
(这个有待证实,目前PAML识别PAUP和PHYLIP建的进化树)
6然后再PAML4的文件夹中找到一个是“codeml.ctl”的文件,复制到和“codeml.exe”相同的地方。
7要对codeml.ctl文件中的各个选项的值进行修改,具体内容如下:
8seqfile=cox1.nuc按你自己的文件名进行修改,就可以了,
9treefile=cox1.trees
outfile=mlcmainresultfilename,
noisy=90,1,2,3,9howmuchrubbishonthescreen,
verbose=00concise;1detailed,2toomuch
runmode=0
seqtype=11codons;2AAs;3codons--AAs
CodonFreq=20161each,1F1X4,2F3X4,3codontable
clock=0
aaDist=00equal,+geometric;-linear,1-6G1974,Miyata,c,p,v,a
aaRatefile=wag.datonlyusedforaaseqswithmodel=empirical(_F) dayhoff.dat,jones.dat,wag.dat,mtmam.dat,oryourown
model=0,这是使用的最简单的模型,modelsforcodons 0one,1b,22ormoredNdSratiosforbranches
modelsforAAsorcodon-translatedAAs
0poisson,1proportional,2Empirical,3Empirical+F286FromCodon,7AAClasses,8REVaa_0,9REVaa(nr=189)
NSsites=031278,依次选取了6个模型。
也可以选其中的两个,但必须是0和3,1和2,7和8。
相互配对
icode=40universalcode;1mammalianmt;2-10seebelow如果是核基因的话就选0。
fix_kappa=0
kappa=5
fix_omega=0
omega=0.2
getSE=0
RateAncestor=0
Small_Diff=.5e-6
cleandata=1removesiteswithambiguitydata(1yes,0no)
method=00simultaneous;1onebranchatatime
10最后保证在同一个文件夹内同时具有:
三个文件“codeml.exe”,“cox1.nuc”“cox1.trees”,这时候你双击codeml.exe,就可以运行程序。
如果不能正确运行的话,你可以通过运行cmd,在dos情况下,运行codeml.exe,这样会有错误提示,知道你错在哪里了。
常见命令解释:
1.Baseml.ctl的命令说明:
2.noisy用来控制输出结果的多少,如果模型适用的运算比较多的话,noisy的值可以选择的比较大,verbose可以控制结果文件中结果的多少。
3.runmode=0表明在树的结构文件中估算树的拓扑结构。
这个选项是我们通常情况下选择的,基本上可以满足我们的需要。
4.Runmode=1or2表明通过星状-分解算法来进行启发式搜索树。
Runmode=2这种算法是从星状树开始搜索,而runmode=1则表明软件读取多歧树是从树的结构文件中,并且同过比较去估计最佳二歧树。
5.runmode=3表明是逐步增加的。
6.runmode=4通过简约法来搜索具有NNIperturbation的起始树。
7.runmode=5表明从树的结构文件中来读取NNIperturbationwith起始树。
8.Model0,1,…,8分别代表以下模型:
JC69,K80,F81,F84,HKY85,T92,TN93,REV(alsoknowasGTR),andUNREST。
9.Mgene用于和序列数据文件中的optionG进行联合,用于多个基因和多个位点的联合分析。
如果不使用optionG的话,则选择0。
10.ndata用于指定文件中的分隔的数据集的数目。
它的变化被用于模拟,你可以使用evolver来产生200个复制数据集,这是设置ndata=200,然后用baseml进行分析。
11.clock用于指定谱系之间速率恒定或变化的模型。
Clock=0,意味着整棵树中,不同分支之间不存在clock现象;Clock=1,意味着globalclock,所有的分支具有相同的进化速率;clock=2意味着localclock,所有分支之间的进化速率被分成几个部分;clock=3意味着多个基因或多重分隔数据,允许分支的进化速率以不同的方式变化。
;;
12.Codeml.ctl的使用说明:
13.CodonFreq用于平衡密码子替换模型中的密码子使用频率。
Codonfreq=0说明每种密码子的使用频率是相同的;codonfreq=1说明是从平均核苷酸频率中计算出来的;codonfreq=2说明是从三个密码子位置的平均核苷酸频率得来的;codonfreq=3则使用了三个参数。
Codonfreq=0,1,2和3所代表的模型中使用的参数的数目分别为:
0,3,9,和60。
14.aadist用于指定氨基酸距离是否是相同的(=0),还是属于Grantham’smatrix(=1)。
15.runmode=-2执行ML方法来推测蛋白序列两两之间的dn和ds。
16.model用于估计各个分支之间的w值。
Model=0,表明所有的谱系具有一个w比率(onewratio);model=1,每一支具有一个速率(free-ratio);model=2表明速率的任意数字。
17.NSsites主要是用于指定模型允许dnds(w)在不同的位点之间变化。
NSsites=m表明对应于model=m。
变化的ncatG被用来指定在一些特定的模型下的w分布的类型的数目。
NcatG的值被用于执行一下分析:
paperare3forM3(discrete),5forM4(freq),10forthecontinuousdistributions(M5gamma,M62gamma,M7beta,M8beta&w,M9beta&gamma,M10beta&gamma+1,M11beta&normal1,andM120&2normal1,M133normal0).ThismeansM8willhave11siteclasses(10fromthebetadistributionplus1additionalclass)。
通过NSsites可以同时执行多个模型,例如:
NSsites=012378,的意思就是同时执行M0,M1,M2a,M3,M7,和M8。
作者建议:
使用M1a和M2a来重建LRT,使用M7和M8来重建LRT,使用M2a和M8来鉴别受到正选择的位点。
18.icode用来更改所选序列的遗传密码子,以期得到更加准确的结果。
19.RateAncestor=1表明你想重建原始序列,如果RateAncestor=0说明你将避免这个计算。
不过使用效果并不明显,还需要进一步研究如何使用。
PAML使用中最重要的就是模型的选择:
PAML中所有的模型都在baseml和codeml这两个程序中使用。
这两个程序是最大似然程序,它们使用数值优化算法来最大化对数似然值。
这些模型最大的用途就是适用likelihoodratiotest(似然比率检验)来检测有趣的生物学假设。
这些模型是在Baseml中使用的,软件中常用的数学模型有:
JC69(JukesandCantor1969),K80(Kimura1980),F81(Felsenstein1981),F84(Felsenstein1984),HKY85(Hasegawa1984,1985),Tamura(1992),TamuraandNei(1993),andREV,alsoknowasGTRforgeneral-time-reversible(Yang1994)。
模型的一般遵循以下假设:
1.在不同的谱系中替换是独立发生的。
2.在不同的位点中替换也是独立发生的。
3.替换的过程我们通过时间均匀马尔科夫过程(time-homogeneousMarkovprocess)。
常用的两种检测方法:
1.Maximumlikelihoodestimates(MLEs):
观测到的数据X的概率(probability),当做为一个未知参数θ的函数的时候,就叫做似然函数(likelihoodfunction):
L(θ:
∣X)=f(θ∣X)。
根据似然规则(likelihoodprinciple),似然函数包括数据中关于参数θ所有的信息。
参数θ的最佳点估计(optimalpointestimate)可以通过最大化似然L的θ值或l(θ;X)的似然对数进行估计。
并且,似然曲线可以为未确定的点估计提供信息。
2.Likelihoodratiotests(LRTs):
假设一个简单模型或无效模型(simplerornullmodel)有一个参数p0,更通用的模型或可选择的模型(generaloralternativemodel)有一个参数p1,两个模型的最佳似然值分别为l0和l1。
那么对数似然值差异(loglikelihooddifference)的两倍是:
2△l=2(l1-l0),如果无效模型(nullmodel)成立的话,那么对数似然值差异的二倍将与自由度是d.f.=p1-p0的卡方分布具有渐进关系(asymptotically)。
因此,对数似然值差异的二倍的检验统计可以通过比较卡方分布来检验无效模型(nullmodel)是否拒绝备择模型(alternativemodel)。
所谓Likelihoodratiotest(似然比率检测)是用来检验两个模型的。
离散伽玛模型(discrete-gammamodel)允许不同位点具有不同的变化速率。
Baseml中有核苷酸替换模型,Codeml中有不同位点替换速率变化的模型。
1.作者在Codeml中进行比较的两个模型比较有:
M1a(NearlyNeutral)和M2a(PositiveSelection);M7(beta)和M8(beta&ω)。
2.作者认为M3对于正选择的LRT检测并不是十分适合,并不推荐适用M3模型。
使用似然比率检测可以验证正选择(Testingpositiveselectionusingthelikelihoodratiotest)。
作者推荐使用二到三种LRT来验证正选择。
第一个检测是比较M1a和M2a,
第二个检测是比较M7和M8。
Gamma分布中形状参数所表示的含义:
1.α>1,大多数位点的替换速率在1附近,但有少数位点具有比较高或比较低的替换速率。
曲线形状为bell-shaped2.α→∞,表明所有的位点具有一个相同速率。
3.α≤1,表明大部分位点的替换速率比较低,或接近于不变,可是有一些位点具有比较高的替换速率。
曲线形状为L-shape。
PAML的一个重要功能就是检测基因是否受到正选择,即适应性选择。
但是现在用于估计适应性选择的方法,忽略了氨基酸的化学性质,这样得出的结果是不准确的,作者表示,直接通过dn和ds的比较来确定受到什么样的选择压力,是不准确的。
PAML中的无效模型(nullmodel)是指不允许任何位点的ω值大于1。
PAML中nullmodel是不允许w值大于1,如果nullmodel成立,则w小于1,基因受到负选择;如果nullm
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- modeltest 用法
![提示](https://static.bingdoc.com/images/bang_tan.gif)