科学网

 找回密码
  注册
科学网 标签 LRT

tag 标签: LRT

相关帖子

版块 作者 回复/查看 最后发表

没有相关内容

相关日志

混合线性模型如何检测固定因子和随机因子的显著性以及计算R2
热度 2 yijiaobai 2019-5-11 19:46
很多朋友写信问我, 像要知道固定因子的显著性和随机因子的显著性如何计算,他们使用的是 lme4 这个R包, 但是这个包使用anova时没有P值,还要手动计算, 随机因子也需要自己计算loglikehood值, 然后使用LRT的卡方检验进行显著性检验, 其实lme4包有扩展的包可以非常友好的做这件事情. 1. 载入数据和软件包 ###载入软件包和数据 library(lme4) library(lmerTest) library(sjstats) library(learnasreml) data(fm) 2. 软件包介绍 lme4 R语言中最流行的混合线性包 结果不太友好, 所以才有下面两个包作为辅助 安装方法 install.packages(lme4) lmerTest 主要是用于检测lme4对象的固定因子和随机因子,它有两个函数: lmerTest::anova.lmerModLmerTest 用于检测固定因子的显著性, 方差分析表采用III平方和的形式. lmerTest::ranova 用于检测随机因子的显著性, 使用的是LRT检验, 给出的是卡方结果. 安装方法 install.packages(lmerTest) sjstats 可以计算R2 可以提取方差组分 安装方法 install.packages(lmerTest) 3. 使用lme4进行混合线性分析 模型介绍 固定因子: Spacing + Rep 随机因子: Fam 建模 ###固定因子:Spacing+Rep,随机因子:Fam fm1-lmer(h1~Spacing+Rep+(1|Fam),fm) 固定因子检验 anova(fm1)#固定因子显著性检验 可以看到Spacing 和Rep都达到极显著 随机因子显著性检验 ranova(fm1)#随机因子显著性检验,LRT 可以看到Fam达到极显著 计算R2 r2(fm1)#计算R2 R-SquaredforGeneralizedLinearMixedModel  #方差组分 还有一个包叫 MuMIn ,也可以计算R2 5. 关于混合线性模型计算R2 library(MuMIn) r.squaredLR(fm1)#计算R2 0.217233511687581 6. 完整代码分享 #混合线性模型,如何检测固定因子和随机因子 ###载入数据 library(lme4) library(lmerTest) library(sjstats) library(learnasreml) data(fm) str(fm) ###固定因子:Spacing+Rep,随机因子:Fam fm1-lmer(h1~Spacing+Rep+(1|Fam),fm) summary(fm1) anova(fm1)#固定因子显著性检验 ranova(fm1)#随机因子显著性检验,LRT r2(fm1)#计算R2 p_value(fm1)#计算每个水平的显著性 re_var(fm1)#计算方差组分 ###对比asreml fm2=asreml(h1~Spacing+Rep,random=~Fam,data=fm) anova(fm2)#固定因子显著性检验,这里anova是anova.asreml fm_Null=asreml(h1~Spacing+Rep,data=fm) lrt.asreml(fm2,fm_Null)#随机因子显著性检验LRT summary(fm2)$varcomp #方差组分 library(MuMIn) r.squaredLR(fm1)#计算R2 ​ 如果您对于数据分析,对于软件操作,对于数据整理,对于结果理解,有任何问题,欢迎联系我。 邮箱: dengfei_2013@163.com 微信公众号:R-breeding
个人分类: 农学统计|9772 次阅读|2 个评论
选择压力分析之EasyCodeML完整篇(By Raindy)
raindyok 2019-2-19 13:58
絮语: 自然选择是五种遗传力(突变、重组、选择、基因流、漂变)之一,选择压力分析更是进化分析中不可或缺的一项重要内容。虽然本人也整理过包括《 PAML FAQ 》、《正选择分析之 Branch site model 篇( By Raindy )》和《 EasyCodeML 使用指南》等多篇相关的日志教程和一份课件《 EasyCodeML 选择压力分析 _ 北林》。由于相关内容比较零碎,值EasyCodeML软件文章接收在版之际,特整理此日志以飨有需要的科研工作者,并祝大家元宵节快乐。 如果您的文章使用 EasyCodeML 进行选择压力分析,欢迎参考以下格式进行引用: Gao F, Chen C, Arab DA, Du Z, He Y, Ho SYW. 2019. EasyCodeML: A visual tool for analysis of selection using CodeML. Ecology and Evolution .DOI: 10.1002/ece3.5015. 软件下载 :https://github.com/BioEasy/EasyCodeML 一、自然选择的检测 非同义替代与同义替换的比值,即: ω 值,也就是通常所说的 dN/dS (或 Ka/Ks )。 ( 1 )当 ω =1 时,中性进化 ( Neutral selection ),即不受选择: ( 2 )当 ω 1 时,正选择( Positive selection ); ( 3 )当 0 ω 1 时,负选择( Negative selection ,也叫净化选择或纯化选择 Purifying selection ) 二、 CodeML 中四个常见的模型假设及其工作流 在 CodeML 中,考虑不同序列间(考虑位点)或系统发育上的支系间(考虑支系 Lineage )的 ω 值不同,主要有以下四类常见的模型: 1. 枝模型( Branch model ) 主要用于对系统发育树中不同支系 ω 值差异性进行界定,主要有三个模型: ( 1 ) One-ratiomodel :假设系统发育树中所有支系的 ω 值相等; ( 2 ) Free-ratiomodel :假设系统发育树中所有支系的 ω 值不相等; ( 3 ) Two-ratiomodel :假设前景枝和背景枝的 ω 值不同; 2. 位点模型 ( Sitemodel ) 主要假设数据集中不同氨基酸位点受的选择压力不同(而不考虑不同支系间受的选择压力差异)。 该模型主要用于检测正选择( ω 1 )作用,共有 8 个不同假设的模型: ( 1 ) M0 (单一比率),即: One-ratio model ,假设所有位点具有相同的 ω 值; ( 2 ) M1a (近中性),假设仅有保守位点( 0 ω 1 )和中性位点( ω = 1 )而没有正选择位点( ω 1 )存在,这两类位点的比率分别为 p0 和 p1 ,其对应的 ω 值分别为 ω 0 、 ω 1 ; ( 3 ) M2a (正选择),该模型在 M1 基础上增加了第三类 ω 值 ,即假设除了保守位点和中性位点外,还存在处于正选择压力下的位点 ( ω 1 ) ,这三类位点的比率分别为 p0 、 p1 和 p2 ,其对应的 ω 值分别为 ω 0 、 ω 1 和 ω 2 ; ( 4 ) M3 (离散),假设所有的位点 ω 值呈简单的离散分布趋势; ( 5 ) M7 ( beta ),假设所有位点的 ω 属于矩阵( 0, 1 )并呈 beta 分布; ( 6 ) M8 ( beta ω ) ,该模型在 M7 基础上增加另一类 ω 值( ω 1 ); ( 7 ) M8a ( beta ω = 1 ),与 M8 模型类似,但将 ω 值固定为 1 ( ω =0 ); 3. 枝位点模型 ( Branch site model ) :主要假设不同氨基酸位点的和不同支系间受的选择压力均存在差异(既考虑位点间也考虑支系间的 ω 值存在差异),共有四个模型 Model A 、 Model B 、 Model C 和 ModelD ,主要参数如下: ( 1 ) Model A (Model 2 , NSites= 2 , ncatG=ignored) ( 2 ) ModelB(Model 2 , NSites= 3 , ncatG=ignored) ( 3 ) ModelC(Model 3 , NSites= 2 , ncatG=ignored) ( 4 ) ModelD(Model 3 , NSites= 3 , ncatG= 2 or 3 ) 4. 进化枝模型 ( Clade Model ):与枝位点模型类型,能同时检测多个进化枝( Clade ),共有 CmC 和 CmD 两种模型,主要参数如下: ( 1 ) CmC(Model 3 , NSites= 2 , ncatG= 2 or 3 ) ( 2 ) CmD(Model 3 , NSites= 3 , ncatG= ignored ) CodeML 分析大致工作流简要如下: ( 1 )配置零假设模型和备选模型的参数( codeml.ctl ); ( 2 )运行 CodeML 程序进行分析获得对应的 LnL 和 np 值; ( 3 )通过似然率检验( LRT )(卡方检验)获得 p 值确定零假设模型和备选模型之间是否存在差异; ( 4 )根据结果进行解读。 三、 EasyCodeML 的主要功能及其改进 EasyCodeML 是一款以 CodeML 为内核的、通过可视化交互式的操作进行选择压力分析的工具,并整合了 CodeML 中主要的基于密码子的模型( Codon-based models )。该软件由 Java 程序语言编写,并已编译好适用于 Windows 、 Mac OS X 和 Linux 三个操作系统预编译版本(其他系统操作未测试,不能保证其正常运行)。在运行 EasyCodeML 前,请检查对应的操作系统中是否已经安装最近版本的 Java 运行环境( Java Runtime Environment, JRE 1.6 或更高版本 ) 。与其他的选择压力分析工具比较, EasyCodeML 主要功能如下: 1. 两种运行模式 EasyCodeML 中提供了两种运行模式, 第一种是预置模式( PresetMode ),让用户彻底告别原 CodeML 晦涩高深的操作,只是轻松点击即可完成,尤其适合新手使用;第二种是定制模式( Custom Mode ),让用户根据需要随时调整或修改相关参数,适合于 CodeML 老手使用; ( 1 ) 预置模式 ( PresetMode ),即在该运行模式下,已经整合了用于选择压力分析的常见成对模型的主要参数,并通过流水线式自动完成分析。这些成对模型包括枝模型( M0 vs. Two-ratio model )、位点模型( M0 vs. M3 、 M1a vs. M2a 、 M7 vs. M8, M8a vs. M8 )、枝位点模型( Model Anull vs Model A )和进化枝模型( M2a_rel vs. CmC )四大类共 7 对。用户只需要根据分析目的在该运行模式下,添加自己的数据(格式化后的序列文件和树文件),选择其中一类模型分析即可; ( 2 ) 定制模式 ( Custom Mode ) , 有点类似于 PamlX ,用户可以任何修改模型中的所有参数以满足分析目的的需要。为了用户调用方便, EasyCodeML 在定制模式下 “Load” 模块中特别整合了一个 “ControlFileViewer” ,可以预加载上述的四大类 14 个常见模型的参数。用户可以快速载入任何一个模型的预先优化的参数,并在此基础调整部分参数。 2. 可视化标记前景枝(或进化枝) 在 CodeML 中,除了位点模型( Site model )外的其他三类模型在分析前都需要先定义前景枝( Foreground branch ),但早期的工具一般是通过文本编辑器在树文件的进行手动标记,这对于新手是个大挑战,一是操作繁琐不直观,二是容易标记效率低。 EasyCodeML 中提供可视化的前景枝标记功能,所见即所得,具有高效、准确等特点,还可以避免由手动修改造成误标问题。 用户只需要在载入的树文件,在需要标记的分支上通过鼠标点击即可完成前景枝的定义。 注意 :定义前景枝一般是基于一定的生物学证据。如果暂无生物学信息可用,可以将所有可能的分支分别定义为前景枝,并通过不同假设进行检验; 3. 其他附带功能 ( 1 )序列格式转换器,除了可以将多种序列格式转为 CodeML 分析所需要的 PML 格式,还可以实现这些序列格式之间的任意转换; ( 2 )自动完成似然率检验( LRT )分析,让用户省去繁琐处理过程(预置模式下); ( 3 ) “Export” 模块可生成准发表级的表格,用户只需要在导出主要结果数据的表格基础上简单修改即可用于发表; ( 4 ) “Check” 模块,用于检查序列名称和树文件中的 Taxa 名称相一致,可以为顺利进行 CodeML 分析作为铺垫; ( 5 )支持文件或文件夹拖曳功能,工作目录、序列文件和树文件,可以直接拖入对应的输入框内; ( 6 )支持多线程操作,主要应用于位点模型( Site model ),多线程运行可以充分利用计算机资源; 四、 EasyCodeML 选择压力分析简明图解 EasyCodeML 运行流程如上图所示 1. 数据准备(序列文件、树文件及其标记) ( 1 )用于选择压力的序列文件必须是比对后的 PML ,如果格式尚未转换,可以通过 EasyCodeML 工具栏下的 “Sequence Format Convertor” 进行转换,支持拖曳操作,自动识别序列格式类型;如果核苷酸序列,则序列长度必须是 3 的倍数( Codon 方式比对)。 ( 2 )用于选择压力分析的树文件必须是 Newick format ( 如: Examples/Example1.tree) ,可以通过 Figtree 等软件导出,如图 2 所示。注意:树中的类别名称( Taxon name )不能带有空格、逗号等非法字符。 ( 3 )除 位点模型( Site models ) 外,与枝相关的模型均需要提前进行前景枝的标记。枝模型( Branch model )和进化枝模型( Clade model )可以标记多个前景枝(或进化枝),但枝位点模型( Branch-site model )只能一次标记一个前景枝,图 3 。 简明操作 Tip :载入 Newick 格式的树文件后,点击 “Label” 按钮,自动弹出操作窗口,选择一个分支后,该分支将以橙色显示,并自动标记上 #1 符号( Clade model 使用 $ ),如果要多个分支需要标记,点击左侧菜单点击 2nd 、 3rd 、 4th 和 5th 以同样方式操作。注意: Branch Model 和 CladeModel 最多支持标记五个前景枝(或进化枝)。 2. 运行模式 2.1 预置模式 ( 1 )选择一个本地文件夹作为工作目录,可以直接通过拖曳操作; ( 2 )选择一个选择分析的模型类型(示例为 CladeModel ); ( 3 )载入 PAML 格式的比对序列(非 PAML 的序列也可以直接拖入,程序会自动转换为 PAML 格式); ( 4 )载入 Newick 格式的树文件(需要事先转换,注意树中的类别名称与比对中的序列名称一致); ( 5 )通过 “Check” 检查树文件中类别名称与比对序列的名称一致性; ( 6 )通过 “Lable” 模块定义前景枝( Site Model 不需要 ) ; ( 7 )保存当前参数配置信息( 强制操作,否则不能运行 ); ( 8 )启动 CodeML 分析(如果需要静默运行,请勾选 In Slient Mode 选项,即可不弹出相关提示信息); ( 9 )自动完成成对模型的似然法检验( LRTs ); ( 10 )通过 “Export” 将主要结果及参数导出并生成为准发表级的表格; ( 11 )启动 Microsoft Excel 查看生成的表格文件; 2.2 定制模式 ( 1 )选择一个本地文件夹作为工作目录; ( 2 )载入 PAML 格式的比对序列; ( 3 )载入 Newick 格式的树文件; ( 4 )选择当前的数据类型(示例数据为 Codon ); ( 5 )检查树文件中类别名称与比对序列的名称一致性; ( 6 )模块定义前景枝( Site Model 不需要 ) ; ( 7 )通过 “Load” 模块调出 Control File Viewer (Fig. 4) ; ( 8 )保存当前参数配置信息(强制项); ( 9 )启动 CodeML 分析; ( 10 )查看运行结果,获得 LnL1 和 np1 值。同样方式获得备选模型的 LnL2 和 np2 值,用于 LRT 分析; ( 11 )通过工具栏的 “LRTCalculator” 进行成对模型的 LRTs 分析 (Fig. 5) ; 3. 结果解读 相关的结果解读,详见 EasyCodeML 的软件文章以及参考文献后的两篇示例数据文献 。 五、 常见问题 1.CodeML 中的四类模型应如何选择? ( 1 )枝模型( Branch model )一般用于检测支系间的选择约束强度( selective constraints , 0 ω 1 )。当比较不同支系间的 ω 值是否显著差异时,一般使用成对模型 One-ratio model vs. Free-ratio model ;当比较前景枝与背景枝的 ω 值是否显著差异时,一般使用成对模型 One-ratio model vs. Two-ratio model ; ( 2 )位点模型 ( Sitemodel )一般用于检测正选择位点(不考虑支系间的 ω 差异),常用的成对模型一般有 M2a vs. M1a 、 M8 vs.M7 和 M8 vs. M8a ,后面两对更为常用; ( 3 )枝位点模型 ( Branch site model )一般用于基因复制事件发生后,检测前景枝中正选择作用对部分部分的影响。常用的成对模型为 Model A vs. Model A null ; ( 4 )进化枝模型 ( Clade Model )与枝位点模型相似,但不限于检测正选择作用,还可以整个进化枝或部分支系上特异位点的选择约束性。常用的成对模型是 CmC model vs. M2a_rel model ; 2. 为什么我的 EasyCodeML 被识别为压缩包? EasyCodeML.jar ,是主程序文件,正常情况下,直接双击即可运行。但如果 *.jar 文件被识别为压缩包,则程序会被系统关联的解压缩软件进行解压缩处理。解决办法是去除 *.jar 与解压缩软件的关联,比如 WinRAR 之类的,在参数设置去除 jar 的关联即可 3. 为什么找到的正选择位点在原始比对序列中找不到? 最主要原因可能是比对序列中带有 gap , codeml 默认忽略 gap 所带的一列数据 ( 即启用 Clean data =1) ,从而导致位置发生偏移。因此,分析前最好将带 gap 的同一列序列全部手工删除,找到正选择位点的氨基酸位置后,再还原对应到原始比对序列上。当然也可以在 EasyCodeML 主界面中把 Clean data 前面的选项勾选取消掉,重新运行分析即可。 4. 前景枝(进化枝)标记问题 与 CodeML 相一致, EasyCodeMLv1.2 支持枝模型( Branch model )和进化枝模型( Clade model ) 一次同时标记多个前景枝(或进化枝) , 但枝位点模型( Branch-site model )只能一次标记一个前景枝 ,而位点模型( Site model )不需要标记。 六、 参考文献及推荐阅读材料 1. 综述文献 ( 1 ) Vitti, J.J., Grossman, S.R., Sabeti, P.C., 2013. Detecting natural selection in genomic data. Annual Reviewsof Genetics47, 97-120. ( 2 ) Sironi, M., Cagliani, R., Forni, D., Clerici, M., 2015. Evolutionary insights into host-pathogen interactions from mammalian sequence data. Nature Reviews Genetics 16, 224-236. 2. 相关软件: ( 1 ) Yang, Z., 2007. PAML 4 : Phylogenetic analysis by maximum likelihood. Molecular Biology and Evolution 24, 1586-1591 ( 2 ) Egan, A., Mahurkar, A., Crabtree, J., Badger, J., Carlton, J., Silva, J., 2008. IDEA : Interactive Display for Evolutionary Analyses. BMC Bioinformatics 9, 524. ( 3 ) Busset, J., Cabau, C., Meslin, C., Pascal, G., 2011. PhyleasProg : a user-oriented web server for wide evolutionary analyses. Nucleic Acids Research 39, W479-W485. ( 4 ) Xu, B., Yang, Z., 2013. pamlX : A graphical user interface for PAML. Molecular Biology and Evolution 30, 2723–2724. ( 5 ) Maldonado, E., Sunagar, K., Almeida, D., Vasconcelos, V., Antunes, A., 2014. IMPACT_S : integrated multiprogram platform to analyze and combine tests of selection. PLOS ONE 9, e96243. ( 6 ) Maldonado, E., Almeida, D., Escalona, T., Khan, I., Vasconcelos, V., Antunes, A., 2016. LMAP : Lightweight multigene analyses in PAML. BMC Bioinformatics 17, 354. ( 7 ) Schott, R.K., Gow, D., Chang, B.S.W., 2016. BlastPhyMe : A toolkit for rapid generation and analysis of protein-coding sequence datasets. BioRxiv. 3. 示例数据 : ( 1 ) Bielawski, J.P., Yang, Z., 2003. Maximum likelihood methods for detecting adaptive evolution after gene duplication. Journal of Structural and Functional Genomics 3, 201-212. ( 2 ) Padhi, A., Verghese, B., Otta, S.K., 2009. Detecting the form of selection in the outer membrane protein C of Enterobacter aerogenes strains and Salmonella species. Microbiological Research 164, 282-289. 七、特别致谢 EasyCodeML 测试版自 2015 年暑假首次推出,至 2019 年 2 月 11 日软件文章被正式接收。四年间 EasyCodeML 历经无数次的调试和修改,感谢陈程杰博士的辛勤付出,也特别感谢陕西博瑞德生物科技有限公司陈振玺、西南大学家蚕基因组生物学国家重点实验室李寒博士、南京师范大学生命科学学院张麟博士和四川农业大学园艺学院陈清老师等的意见和建议,使得程序日臻完善。 Raindy: 本文首发于本人的QQ空间(http://user.qzone.qq.com/58001704/blog/1549031707),欢迎转载,但请保留作者原信息。
个人分类: 软件教程|31184 次阅读|0 个评论

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-5-17 08:13

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部