IQtree 默认构建的是 无根 系统发育树,但 v2.1.3 版本后 IQtree 提供了使用 Lie Markov Model 算法对系统发育树实现 无外群生根 功能(Rooting phylogenetic trees)。但在使用过程中,由于功能较新,IQtree 在 处理 SNP 数据 的过程中包含较多 Documentation 未说明的 Bug 。
本文主要包含两部分:系统发育树生根操作及踩坑内容,不包含生根原理。其中踩坑内容了解即可,可能随着 IQtree 版本的更新而得到修复。
PS:SNP 数据的处理步骤参见 IQtree:使用 SNP 数据(vcf file)构建玉米群体的 无根 系统发育树
本文使用的数据集矩阵大小为 577 × 10000。样本取自 577 个玉米品系,包含 70 个大刍草(par)、208 个热带品系(tst)、172 个温带品系(temp)、127 个混合品系(mixed)。SNP 取自 3 号染色体上 743w SNP,经过过滤(maf > 0,miss = 0)与随机抽样后,最终得到 1w SNP。数据已上传(链接如下,0 积分下载)。
IQtree 使用 不可逆(non-reversible)模型中的 Lie Markov 模型进行无外群系统发育树根的位置的推断。
1)使用 IQtree 对所有 Lie Markov 模型进行 遍历,寻找 最适 模型
iqtree -s SNP.varsites.phy -mset liemarkov -nt 100 -fast --prefix SNP.varsites.liemarkov.fast
-mset 遍历所有指定集合(本命令为 liemarkov)的替换模型
计算每种模型下建树的似然值,并筛选出最优模型
本命令会逐个检索 Lie Markov 模型,当遍历到 6.7a 之后时,IQtree 会提示 2 个 Warning:
WARNING: Numerical underflow for non-rev lh-branch Noname
WARNING: evec not invertible
当遍历到 9.2b 后,IQtree 会 报错。上述未报错的各种模型里 6.7a+ASC+R5 的 log-likelihood 值最大(-78103)。经过不同数据集的多次实验,作者发现使用 SNP.varsites.phy 这种 不含 任何 constant sites 的序列作为输入时,IQtree 大概率会在遍历到某个模型时 报错 。
2)使用最适模型 6.7a+ASC+R5 构建系统发育树
iqtree -s SNP.varsites.phy -m 6.7a+R5+ASC -nt 100 -fast --prefix SNP.varsites.6.7_R5_ASC.fast
系统发育树构建结果如上图所示,其中红色为大刍草品系,蓝色为热带品系,黑色为温带品系。
本文的生根步骤与 IQtree 教程差异较大,IQtree 建议包含 2 步,先局部再合并生根:
1. 使用常规模型对序列进行局部最优建树(-p)。
iqtree2 -s SNP.varsites.phy -p SNP.varsites.nex -T AUTO --prefix rev_snp
2. 使用 Lie Markov 模型将多个局部系统发育树合并且生根(--model-joint 12.12)。
iqtree2 -s SNP.varsites.phy -p rev_snp.best_scheme.nex --model-joint 12.12 -T AUTO --prefix nonrev_snp
先局部再合并生根策略的目的可能是提高建树与生根的准确性,但 SNP.varsites.phy 文件运行第二步合并时会报错。除此之外,Lie Markov 模型与 许多功能连用 都会 报错,现记录如下:
1. --model-joint
iqtree -s SNP.varsites.phy -p rev_snp.best_scheme.nex --model-joint 6.7a+ASC+R5 -nt 100 -fast --prefix nonrev_snp
2. --lmap
iqtree -s SNP.varsites.phy -m 6.7a+ASC+R5 -lmap 10000 -fast -nt 100 --prefix nonrev_snp.lmap1w
3. --root-test
iqtree -s SNP.varsites.phy -m 6.7a+ASC+R5 --root-test -zb 1000 -te SNP.varsites.6.7_R5_ASC.fast.treefile -fast -nt 100 --prefix nonrev_snp.root_test
PS:除 IQtree 教程推荐的先局部再合并生根策略外,作者还考虑过 限制发育树结构 ,以提高建树与生根的准确性。
1. 先构建无根进化树
iqtree -s SNP.varsites.phy -m GTR+F+R6+ASC -fast -nt 100 --prefix rev_snp
2. 再以无根进化树为模板(-te)生根
iqtree -s SNP.varsites.phy -m 6.7a+ASC+R5 -te rev_snp.treefile -fast -nt 100 --prefix nonrev_snp
但实验结果显示,IQtree 不修改 -te 所提供的树的拓扑结构,只修改 整体 的 分枝长度。对比 log 日志可以发现,-te 参数下 iqtree 只进行 FINALIZING TREE SEARCH 步骤,因为 Estimate model parameters (epsilon = 0.050)。下表罗列了 3 次重复实验的结果,其中 treefile 输入 iTOL 后可以发现拓扑结构完全一致,分枝之间的长度关系基本完全一致,区别 仅为 比例尺(tree scale)。
Log-likelihood | Total tree length | Sum of internal branch lengths ratio | |
---|---|---|---|
repeat 1 | -77441 | 1.136 | 19.94 % |
repeat 2 | -77442 | 1.077 | 19.94 % |
repeat 3 | -77397 | 1.248 | 19.91 % |