科学网 › 标签 › prior

标签: prior

相关帖子	版块	作者	回复/查看	最后发表

没有相关内容

相关日志

BEAST 分析 FAQ （20180302更新）: 热度 5 raindyok 2018-2-6 14:37; Raindy 注：本文首发于本人QQ空间日志（https://user.qzone.qq.com/58001704/blog/1502418934），将不定期更新，希望对初学者有帮助，欢迎关注。如需转载请保留原作者信息。 Q1: BEAST 有1.x 和 2.x 两个系列，如何选择？ A1 : 可以根据分析不同目的进行选择，两个版本的特点简述如下： I、BEAST 1的特点：（1）适用于分析群体或物种级别的数据；（2）可以同时评估树和节点的时间；（3）分子钟模型的范围更广；（4）树先验(Tree Prior) 和 demographic 模型范围更广； II、BEAST 2的特点：（1）重写的程序，更多模块化；（2）可以通过添加包(Packages)来扩展BEAST；（3）BEAST 1中的树先验(Tree Prior) 在该版本中不可用；（4）支持数据模拟；更多的功能比较，可以访问BEAST官网查看（ http://beast2.org/beast-features/ ） Q2: 如何比较BEAST中哪个模型最佳？ A2 : 一般通过贝叶斯因子法（Bayes Factor, BF）来比较两个模型的边际似然值(Marginal Likelihood)，简单办法：通过BF法分析后，哪个模型对应的LnL值越大，相对更优。具体也可以通过以下操作来确定：在Tracer 中同时导入两个模型的xx.log 文件，选中所有要比较的模型后，依次点击菜单Analysis---Model Comparation，Analysis type改为Harmonic mean，Reps=100，最后根据获得的BF值（BF=Pr(D|M1)/Pr(D|M2)至少大于3，则说明M1M2）范围判断模型间的优劣程度。下图所示，第三个数据最佳（P_UCL_Constant_80M.log） Q3: 模型中位点间的速率变异呈gamma分布，如何确定较为准确的Gamma类型数(Number of Gamma Categoies, nCat)？ A3: Gamma分布的形状通过alpha参数来确定，如下所示。当alpha值很小时，大部分位点慢速进化，只有一小部分位点快速进化；但当alpha值非常大时，大部分的位点以接近相同的速度进化。因此，当alpha值比较小时（如：shape=0.32，该值可以通过模型软件计算获得），可以将nCat设置为6（默认值为4）进行分析。 Q4：BEAST 分析过程中，断电了怎么办？ A4: BEAST 1.x 没有断点续行功能，可以更名后重新运行，最后合并前后运行的数据；BEAST 2.x 有断点续行功能，可以追加运行，有点类似于Mrbayes的Checkpoint功能，如下图所示。 Q5：Baysian skyline analysis 时，popsize和groupsize一直不收敛，怎么办？ A5: 这两个参数不收敛，说明BEAUti配置时这两个所使用的默认参数不合适，可以在Operator标签修改（下图）或直接用记事本类工具修改，一般修改为原来值的10倍，如下图所示，将15.0改为150，6.0修改为60即可...此方式也适用于 Exponential growth 模型。 Q6：如何给BEAST及相关的工具增加虚拟内存？ A6: BEAST运算生成的结果文件，特别是tree文件，通常都比较大，使用TreeAnnotator 生成MCMC树，经常因为内存不足报错，其主要原因是虚拟内存设置过低。解决的办法有两种，一种是直接修改Java程序包中的Info.plist文件；另一种是用命令行运行，类似“java -Xms64m -Xmx256m -jar lib/beast.jar”。前者修改后一劳永逸，推荐使用，选择程序对象，如BEAST后，右键选择“显示包内容”，在打开的“Contents”目录内用文本编辑器工具打开Info.plist（mac系统推荐使用免费Pref setter）将JVMOptions下的Item 1的虚拟内存调大，如 -Xmx256m改为-Xmx6g，如下图所示： Q7：如何判断自己的数据适合严格（strict）分子钟还是宽松（relaxed）分子钟模型? A7: 可以先用宽松分子钟结合简单的Tree piror配置个xml文件，MCMC链长不需要很大，运行结束后，将xx.log文件通过Tracer导入，查看xxx.coefficientOfVariation 这个参数，该值介于0和1之间。如果值越接近于0，说明数据更接近严格分子钟模型；反之，更偏向于宽松分子钟模型，如下图，该值为0.511，说明该数据适合宽松分子钟模型。 PS ：如果分析的数据有时间信号，可以通过Phylogenetic dating的treedater包来分析，详见日志： http://user.qzone.qq.com/58001704/blog/1506493780 Q8：如何判断数据中的采样时间是否可用于分子钟校准 ? A8: 一般启用BEAUTi中的Tip dates 功能前，需要做个日期随机化检验(Date-randomizationtests, DRTs)；如果没有通过DRTs，则说明数据中没有时间结构，不可用于分子钟校准。常见的DRTs 分析软件有TipDatingBeast(R包)、least-squares dating (LSD) ，如下图所示，该数据的采样时间具有时间信号，可以用于分子钟校准，相关方法描述详见本人2017年发在EVA上的文章（ http://onlinelibrary.wiley.com/doi/10.1111/eva.12459/abstract ）。 \0 \0 Q9：如何绘制 Extended Bayesian Skyline Plot 图？ A9: 常规的Bayesian Skyline 可以直接通过Tracer 的Analysis菜单中“Bayesian Skyline Reconstruction”来生成，然而Extended Bayesian Skyline Plot (EBSP) 却不能直接应用Tracer操作，早期绘制EBSP图也比较麻烦，通常需要一些脚本来实现。其实一些高版的BEAST已经内置EBSP analyser，只是需要通过命令行来操作，可以通过cd 命令切换到beast.jar 所在的目录后，输入 java -cp beast.jar beast.app.tools.EBSPAnalyser即可调用出EBSP analyser，如下图所示。如果喜欢用R脚本，当然可以下载 RplotEBS.R （ https://sites.google.com/site/santiagosnchezrmirez/home/software/r ）来操作 \0 \0 Q10：如何获得MCC树中的根状态后验概率（ root state posterior probabilities）？ A10：可以打开 TreeAnnotator 对定义Trait的Sample tree进行注释，最后得到MCC树，用Notepad++之类的工具打开，并将光标定位到文件末尾，找到一串类似于“ xxx .set.prob”的字符串，=后紧随的是各定义Trait的后验概率数值，对应的Trait 顺序可以向后找到“ xxx .set=”，比如： Country .set.prob={9.777E-4,0.989,...}, Country .set={China,Turkey,Spain,...} 即中国、土耳其、...对应的root state posterior probabilities分别为9.777E-4,0.989,...。这些数据提取到Excel中通过绘制簇形柱形图即可得到类似下图左上角的效果图。 \0 \0 Q11：如何通过LSD快速获得tMRCA和进化速率？ A11：估算tMRCA和进化速率常见的方法有：（1）Bayesian算法（软件：BEAST）；（2）Least-squares dating（软件LSD）；（3）Root-to-tip (软件TempEst)；（4）Treedater （R包 Treedater ）；(5) TreeTime (Python 包)。其中BEAST最为常用，但其运算量相对较大。相比较之下，LSD以快速著称，可以支持不同算法构建的系统发育分析作为输入源。应用LSD简单三步可以快速估算 tMRCA和进化速率： Step1：获得通过ML法或其他算法构建一个系统发育树； Step2：准备两个文件，一个采样日期（类似格式：taxa1 2001.xx）和序列的Outgroup文件； Step3：命令行运行 ./bin/lsd_mac -i ./sample.tre -d ./sample.date -c -r l -f 100 \0 \0 如上图：估算的某植物病毒CP基因进化速率为 5.972E-5 subs/site/year（95%CI: 2.587E-5-7.321E-5） Q12: 如何应用文献报道的化石记录进行分子钟校准？ Q12: 通常两个步骤可以完成：第一，将序列导入BEAUti后，切换到“Taxon Sets” 标签，点击左下角的“+”添加一个新的sets（示例Taxa sets 名为Test），即有化石记录的Taxa归为一个Taxon Sets，并设置其他参数完毕，切换到“Prior”，双击点开参数“tmrca(Test)”的Prior窗口，选择合适的Piror Distribution（示例为Uniform），文献报道将化石校准点 Min=528.28 Ma 和 Max=636.1 Ma 分别填入 Upper 和 Lower 后的文本框内，“Generate BEAST File” 文件生成xml文件。第二，在生成的xml文件中，找到Starting Tree 位置，把冗余的代码删除，替换上自定义的Newick格式的树，如下图所示： Q13：如何解决设置化石校准后xml 运行过程中Prior出现 zero probability 问题？ A13：主要原因：用户自定义的起始树无法正确约束校准点，例如Q12的校准点为 528.28-636.1 Ma的 uniform prior，也就是说任何 Root age 低于528.28 Myr 或高于636.1 Ma 的树将被指定为概率为零（a zero probability）. 解决办法，将 Prior Distribution 由 Uniform 改为 normalPrior，并设置合适的区间范围，将Root age 限制在528.28-636.1 Ma区间内，如下图所示： Q14：BEAST运行程序过程中，遇到磁盘空间不够或内存不足导致xml文件末尾部分内容未运行完，如何处理？ A14：近期做EBSP分析时就遇到类似的情形，除了用于EBSPplot的xx.csv 未生成外，其他分析均已完成，此时出现内存不足，重新再跑一遍是件非常痛苦的事，尤其对于大数据，有时候时间紧迫也不现实。解决的办法也简单，打开xml文件把已分析完成的部分内容删除，保留尚未分析的事项（下图所示），另存为一个xml，调用BEAST重新运行一下即可，大大节省时间。 Q15：BEAST分析时，链长需要和采样频率设置多少比较合适？ A15：没有普适的标准，一般是10^8（8次方），最主要的一点需要保证最后的样本最少得有10,000，比如：总链长300,000,000，采样频率可以使用30,000，即保证最终样本量达到10,000。 To be continued...; 个人分类: 软件教程|31912 次阅读|8 个评论

building modeling: 热度 2 justinzhao 2012-1-10 15:27; 好久没写博文了,新学期里, 计划有空就多写写, 好记性不如烂笔头, 看懂一篇paper, 和能把这篇paper 写出来是两码事. 看的时候总是从宏观层面来理解, 而文中中的一些细节会忽略掉. 而写出来, 为了保证准确无误,及时对于一个公式, 往往都要查阅很多资料. PCA, SVM, SVD, LDA 等等都是机器学习中随处可见的概念, 相信只要对机器学习一点点了解的同学, 对它们都不陌生, 如果要你不查阅任何资料, 写一份tutorial, 估计好多人都傻了眼.所以决定以后勤写读书笔记, 勤写survey. 昨天看了两篇关于building modeling from LIDAR 的paper, 一篇是法国INRIA的Florent Lafarge的2011ICCV文章, 另一篇是sarnoff公司08年CVPR文章, after initial classification, 接下来核心的一步是labeling, 两篇文章都用到了MRF, 他们的不同在于pairwise energy term, 前者定义了breakline, 位于breakline两侧的points如果有不同的labels, penalty会小;反之会大. 后者的employ了potts model. 在去年自己投的一篇CVPR文章中, 同样是用了MRF优化 initial labels.不像Florent侧重pairwise energy, 我把重点放在了data energy上面.还有, 前2篇文章的高明之处是, 一个object 就是一个点, 而我则是在grid depth image上面操作的. 之所以用MRF, 是为了embed prior knowledge into unorganized points. MRF是一种很好的加入先验知识的模型, 有严格的数学基础. 印象很深刻的是sarnoff公司的那篇paper,他们用了EM来拟合平面, 其间计算出每个点属于平面的概率. NOTE: 做3D building modeling, 法国的inria 在全球可谓遥遥领先, 做相关研究的学者多多关注他们. 对于机器学习, 特别是data mining, 关注怎样从 Unstructured data中提取information.对于自动的信息提取, 人们往往限定一些constraints,提取出符合constraints的信息, 其实这些constraints就是prior knowledge. 如果能用数学模型对prior knowledge描述, 至少增添了方法的理论基础. 其实MRF便是有数学模型的prior knowledge. 90年代末期出现的boosting,给模式分类带来了新的生机, 不过boosting是基于在training data 完全可信的前提下的,并且 typical boosting是没有考虑sample之间的关联的, 如是schapire在99年提出了incorporate prior knowledge into boosting, 提高了 newsgroups, AP-titles等bench-mark datasets的分类精度. 现在有两点值得研究: 1) 如果training data中有少量noise怎么办? 虽然training data 一般情况下是manual sampling的结果; 2) 是否可以考虑基于prior knowledge的boosting对 LIDAR data 分类? 下面说说看书与动手的关系. 看懂了, 有很大的可能你实现不了这个算法; 而从实现算法的过程中,你既理解了一些概念,也许还能发现一些新的问题. 以后做一个课题, 如果有很多中方法, 你又没有自己的idea, 还不如实现几种经典的方法. 总之, 勤于动手! 接下来的几天, 要把Kullback-leibler divergence between two matrix, Levenberg-marquardr 搞懂. 总之, 要多记录, 彷徨时多动手. 参考文献: florent. Building large urban environments from unstructured point data, 2011, ICCV BC Matei, HS Sawhney. Building segmentation for densely built urban regions using aerial lidar data , 2008, CVPR John Secord and Avideh Zakhor. Tree Detection in Urban Regions Using Aerial LiDAR and Image Data, berkeley Robert E. Schapire. The Boosting Approach to Machine Learning: an overview. 2003, ATT lab Robert E. Schapire.incorporate prior knowledge into boosting, 2002.; 个人分类: 读书日记|3807 次阅读|2 个评论

更多...

帐号		自动登录	找回密码
密码			注册

关闭 安全验证

标签: prior

相关帖子

相关日志

关闭安全验证