序言： 准确推断快速进化病原（如 RNA病毒等）的进化速率或最近共祖时间，可以通过对 含有时间信号（Temporal signal）或时间结构（Time-structured）的序列进行系统发育分析获得。目前判断数据集中是否有足够的时间信号的方法主要有Root-to-tip （RTT）线性回归分析和日期随机化检验（DRT，详见日志 https://user.qzone.qq.com/58001704/blog/1537536180 ）。然而，这两种方法都有各自的优缺点，比如RTT主要基于严格分子钟假设，而DRT需要比较原始数据和随机化数据，运算量非常大。 本文介绍一个完全基于贝叶斯模型比较的一种检测方法BETS，即 Beysian evaluation of Temporal Signal。该方法主要引入一组成对的模型：异时模型（Heterochronous model, M het )和等时模型（Isochronous model, M iso ），采用广义垫脚石采样法（Generalized stepping-stone sampling, GSS）获得两个模型的边际似然估值（marginal likelihood estimation, MLE），最后应用贝叶斯因子法（Bayes Factor, BF）确定哪个模型更适合分析的数据集。其中，异时模型（ M het )是应用实际采样时间校准分子钟（Tip date）计算替换速率，而等时模型（ M iso ）不启用采样时间校准功能（即：without tip date）而假定以一个随意的替换速率估值。 是否有 时间信号的判断依据 ： M het 优于 M iso 时，说明数据集含有时间信号。 当logBF=Log(P(Y|Mhet))-Log(P(Y|Miso)) 1时，说明 M het 略优于 M iso ； 当logBF=Log(P(Y|Mhet))-Log(P(Y|Miso)) 1时3时， M het 明显优于 M iso ； 当logBF=Log(P(Y|Mhet))-Log(P(Y|Miso)) 1时5时， M het 绝对优于 M iso 。 准备工作 1. nexus格式的比对序列 为方便提取序列的时间信息，建议每条序列名称命名为类似 Accession_date格式，如：AY99893_2001.67211 2.BEAST 1.x系列 需要提前计算好序列数据的最适替代模型和树先验，本实例为简化直接采用GTR+F+I和Constant size的树先验进行配置。 分析流程： 1. M het 模型的设置 Step 1 ：打开BEAUti配置XML文件，在Paritions标签下载入Nexus格式的序列，如下图所示 Step 2 ：Tips标签下，勾选“Use tip dates”，点击“Parse Dates”并设置规则提取每个序列中年份信息（示例数据年份前有 _符号间隔） Step 3: 在Site标签下，设置数据集对应的替换模型，示例数据为GTR+F+I，这里+F为Empirical 碱基频率，+I为InvariantSites Step 5: 为演示方便，分子钟头模型用Strict模型，Tree Prior用Coalescent的ConstantSize模型，这里不作截图，直接切换MCMC标签，在MLE下拉菜单选择GSS算法后，点击“Settings”进行设置： GSS的step数、链长还有采样频率可以参考PS/SS一文（ https://user.qzone.qq.com/58001704/blog/1529390689 ）设置，注意选择Tree working prior，与Tree Prior的对应，示例为Constant size，故这里选择Matching coalescent model，如下图所示： 设置完成，可以通过界面右下角的“GenerateBEAST File”生成XML，保存为xxxx_tipdate.xml后用 BEAST打开并运行，直至完成。 2. M iso 模型的设置 同上，将标签切换至Tips，取消勾选“Use tip dates”，如下所示： 由于没启用tip功能，需要通过Prior标签设置一下固定的速率值，选择Clock.rate参数，Prior Distribution选择CMTC Rate Reference输入一个固定的速率值，如2.291E-3 （理论上可以任意值，不过建议参考 Root-to-tip等的结果进行设置）。 设置完毕，将生成并保存的XML通过BEAST运行，直至结果，在BEAST软件目录会生成对应的文件 xxxxx.mle.result.log。 3. BF法比较 M het 模型与 M iso 模型 打开运行结果文件xxxxx.mle.result.log，找到尾部复制log marginal likelihood后的数值进行比较。 M het 的MLE值为 -4133.385， M iso 的MLE值为 -4411.406，logBF= 278.021，远大于5 。 4.结果解读 根据计算获得的BF值可知， M het 绝对 优于 M iso 表明数据集具有足够的时间信号，可用于后续的Bayesiandating分析。 参考： 1. BETS简明教程： http://beast.community/bets_tutorial 2.BETS软件文章 (Preprint版)： https://www.biorxiv.org/content/10.1101/810697v1 推荐阅读 BET 引用文章 （1）Livia et al., 2020, Nature Microbiology, https://doi.org/10.1038/s41564-020-0706-0 （2）Düx et al., 2020, Science, https://doi:10.1126/science.aba9411 示例数据 请访问我的网盘 http://raindy.ys168.com/ 示例数据目录内
标签: Tip date
