做细菌组装时,一条两万多碱基的序列有很高的深度,相似性也很好,看起来像一条完整质粒。能否把它画成圆环?决定性的证据往往藏在最不显眼的地方,即有没有单条分子真正从序列的末端跨到开头。我们用 Flye 做了两个从头组装,只从第二份数据中拿掉一类读段。剩余的小分子数据仍超过五十倍覆盖,单看局部序列完全正确,组装报告给出的拓扑却变了。两轮的大片段同时充当对照,能够告诉我们改变具体发生在哪里。Flye 项目与重复图组装说明
先定义什么叫一条跨接读段
实验生成两个无重复设计的随机环状模板,长度分别为七万四千和两万六千碱基。随机选择起点,从圆上连续取六千五百至九千五百个碱基,奇偶编号交替使用正向和反向。若起点加读长超过模板长度,说明这条分子在序列文件选定的切点两边都有碱基,无论最后把它写成哪个方向,都记为跨接。这里的切点只是为存储环状序列而任意选择的一处,不能把它理解成细胞中真实存在的物理断口。
大模板抽取四百三十条,小模板抽取二百三十八条,其中七十一条跨过小模板切点。第一轮把全部六百六十八条交给 Flye 2.9.3。第二轮保留大模板的全部读段,只删除小模板的七十一条跨接读段,最后剩五百九十七条。两轮最小重叠长度都设为两千,使用同一高准确纳米孔输入模式与相同参数。这里的读段是根据已知模板生成的无错序列,适合观察局部覆盖、闭环证据及文件表示之间的关系。它们没有模拟真实测序中的同聚物误差和近缘质粒共享重复,不能据此估算实际项目的闭环成功率。
长度、局部覆盖和闭环信号分别去哪里读
第一轮的组装信息表有两个片段,长度恰为两条模板的真值,环状列均为 Y。第二轮大模板仍得到七万四千碱基,仍然标记为环状,小模板则得到二万五千七百四十八个碱基,环状列变为 N。两轮的片段编号还交换了。如果只把第一轮的一号片段和第二轮的一号片段放在同一行,会误以为原来完整的大片段一下缩短了数万碱基。跨轮比较应先按序列比回对应模板,或者按一组可靠的独特标记确认身份,再比较长度与拓扑。
小模板剩余读段的碱基总量除以真实模板长度,第一轮为七十三点三三倍,第二轮仍有五十一点三二倍。Flye 结果表的小片段覆盖列相应为七十七和六十。两个数值由不同方法计算,前者用模拟器记录的原始总产量,后者是软件对组装片段的估计。特别是删除跨接证据后,组装片段本身的有效长度和被覆盖的区间也变了。一个数字大于五十可以和末端完全缺少连续证据同时发生,不能用平均深度替末端投票。Flye 输出文件与覆盖字段
| 进入组装的小模板读段 | 跨切点读段 | 小模板原始平均深度 | 拼出的长度 | 环状列 |
|---|---|---|---|---|
| 全部 238 条 | 71 | 73.33 倍 | 26000 bp | Y |
| 删除跨接后的 167 条 | 0 | 51.32 倍 | 25748 bp | N |

在样本交接时,可以把环状候选单独导出一张表,记录序列名、长度、组装报告的环状标记、首尾连接点和支持连接的独立读名。若读段支持来自同一处很短的重复,表里需要另记接头两边各有多少独特碱基。这样后来即使有人把旋转后的质粒重新编号,也可以回到原分子复核。
实跑经验|让环状序列的起点可以自由旋转
flye --nano-hq reads.fa --out-dir assembly --genome-size 100k --threads 2 --min-overlap 2000
cat assembly/assembly_info.txt
minimap2 -t 1 -x asm5 -c --secondary=no truth_twice.fa assembly/assembly.fasta > check.paf
比对核验时,我们把每条真值首尾连接复制两遍,作为临时参考,查询仍用单份组装序列。环可以从任一点写起,反向互补也表示同一个双链分子。如果只比较两个 FASTA 文件的第一行,完全相同的环可能被切成两个局部比对。重复真值给它留出旋转空间,随后仍然按单个圆的长度计算覆盖。第一轮两个查询在自己的模板上都全长一致,小模板以反向命中,二万六千个对齐列全部相同。两倍参考是为核对方向和切点而构造,不能把重复的参考长度作为生物学分母。
第二轮的小片段同样每个对齐列都匹配,但长度只有二万五千七百四十八。双份参考上的零基区间为 26200 到 51948,折回原模板后是 200 到 25948,意味着文件切点之后二百个碱基、回到切点之前五十二个碱基都没有出现在这一轮组装里。总计缺二百五十二。只报已比对部分百分之百一致,会遮住这个差别。对真实质粒可以先报告已比对长度与环状标记,再从原始读段中找两侧都有足够独特序列的跨接分子,检查其方向、起点分散程度和是否存在另一处同样好的落点。仅有共享重复区中的匹配,仍然难以确定究竟接回了哪条分子。
这次安装还遇到一个可复现的故障。Ubuntu 软件包中的一个解析步骤调用名为 flye-samtools 的程序,系统装有 SAMtools 却没有这个命令名。日志先出现找不到命令,随后报共识文件为空。确认系统版本和调用的子命令后,在实验自己的目录为已安装的程序建一个同名链接并加入该次运行的搜索路径,两轮才都产生了完整 FASTA 和组装信息表。遇到末段报空序列时,向前查第一条具体错误,比反复调整最小重叠更快,也能避免把安装故障当作低覆盖样本的性质。
需要对整批菌株保存组装图、原始读段与环状接头证据时,可了解云生信的生物信息分析服务器。
