SPAdes 只拼出一条序列，就算组装完整了吗？用重复区、双端跨度和环形坐标做一次验算

一份细菌组装结果只有一条序列，里面没有未知碱基，长度也接近预期，很多人会直接在项目表里勾选完整。但这些信息还无法回答两个重要问题。两段相同的重复有没有放在各自正确的邻居之间，序列头尾是否只是输出时人为选择的断点。我们安装环境中实际使用 SPAdes 3.15.5，用同一个含重复区的已知模板跑了两轮组装。两轮都只有一条序列，完整性却不同。下面把读对怎样提供跨重复区的证据，以及怎样检查最终序列，写成可以独立重跑的步骤。

单条读段读到什么，双端片段又连接了什么

双端测序的一对读段来自同一个 DNA 片段的两端。每端一百五十碱基描述的是实际测到的序列长度，片段长三百或一千零五十碱基描述的是两端之间的几何关系。即使总的已测碱基数相同，较长的片段也可能让两端落在一段重复的左右两个独特区域。此时读对把两个独特邻居关联起来。中间未直接读出的区域仍要由其他读段支持，片段跨度能够帮助选择路径，却不会凭空给出每个未测碱基的身份。

考虑两份完全一样的重复，分别夹在甲乙与乙丙的独特序列之间。只读到重复内部时，读段可以来自任一拷贝。在一端落于重复内部、另一端落于独特区域的情况下，信息会增加，但仍可能不足以同时确定重复两边的关系。如果同一读对的两端完整落在重复左右的独特区，就更容易判断这一次重复连接的是哪两个邻居。真正有价值的是这种能够区分候选路径的跨度。仅仅提高重复内部的覆盖，会让已有重复片段被更多次观察，未必补足其两端的连接关系。SPAdes 原始论文 · 作者使用手册

实跑经验 · 读段长度和对数相同，只改变片段跨度

实验先生成一个一万二千零四十碱基的圆形模板，其中两个六百二十碱基的重复拷贝完全一致，其他三个独特区依次长三千六百、三千五百和三千七百碱基。用固定种子产生一千八百个片段起点，每一对读段各读取一百五十碱基。第一组片段长三百，第二组长一千零五十。两组的第一端序列和起点完全相同，第二端随着片段跨度改变。没有人为加入碱基错误，也没有混入其他基因组，因此可以把输出直接对照设计真值。重跑脚本和逐项实测结果

两轮均从双端 FASTQ 进入普通 SPAdes 的仅组装模式，使用词长二十一、三十三和五十五，固定两个线程和相同内存参数。本次运行研究的是一个基因组内的重复，不涉及群落混合输入或宏基因组模式的表现。去掉纠错步骤，是因为模拟读段本身没有加入错误。真实测序若包含错误，应另行评估适合该数据的纠错与模式选择。把这些前提写进记录，其他人才能知道这组结果能回答哪一个具体问题。

python3 106-spades/run.py "$(command -v spades.py)" work106
# 复查生成的一份输入时，请另用输出目录
spades.py --only-assembler -1 work106/short/r1.fq \
  -2 work106/short/r2.fq -k 21,33,55 -t 2 -m 2 -o short_again

脚本除了启动组装，还逐条检查输入读对的已知位置。要求两端都完全位于重复两侧的独特区，短片段组在两个重复位置各有零对满足条件，长片段组分别有十五对和十八对。这个计数描述的是模拟输入确实具备的跨区证据。现实数据的相同统计需要对齐后的坐标、方向和比对质量，还需要判断多处比对与光学或扩增重复。不能把软件估计的平均片段长度直接当作已经验证的跨区读对数量。

两轮的 scaffolds 与 contigs 在这些统计上相同，没有通过插入未知碱基将片段勉强连接。短片段结果少了五百一十个参考位置，长片段结果恢复了全部参考位置。后者的输出长度比模板多五十五碱基，正好对应当前最后一轮词长的末端重叠。因为有已知真值，可以检查这些额外碱基是否确实和首端重复。若只有天然样本结果，看到一个相同长度的末端重叠仍应继续确认头尾连接证据，不能凭长度巧合宣布一个新环状分子已经闭合。

一条长序列和一次完整的参考覆盖，为什么需要分开统计

汇总序列长度时，同一参考位置可能因重复输出或头尾重叠而被计算多次。我们的脚本先找每条序列是否能完整精确落到设计参考的某个方向，再将覆盖的坐标按模板长度取模，统计不同参考位置的并集。这样长片段结果虽然输出一万二千零九十五个碱基，覆盖量仍是一万二千零四十。对圆形模板，这种计数能剥离输出起点及末端重叠的影响。对真实参考的近缘样本，则要采用容许变异的比对，并单独记录错配、插入缺失和多重放置。

仅对比总长度还会掩盖另一类问题。缺失一小段独特区，同时把某处重复多放了一份，二者长度可能相互抵消。若只报告总长接近预期、最大片段足够长或连续性指标升高，错误路径仍然可能混在结果中。具体检查时应同时保留覆盖了哪些参考位置、同一位置出现几份，以及邻居连接是否正确。参考是近缘对象时，还应考虑真实结构差异，让新增的连接接受原始读段复查。把差异全部强行修成参考样子，也会掩盖待研究样本真正的结构。

本例中所有计入恢复量的序列，都要求沿一个连续方向精确匹配设计模板。这个要求比“每个小片段都能在参考上找到相似区”更严格。后者即便把两段错误交换顺序，也可能得到很高的总匹配比例。检查一个读者真正关心的连接，可以选择重复两边各四十个碱基作为独特标记，确认同一条输出序列中二者之间的长度和实际序列。我们同时寻找预期的两种连接和两种错误交叉连接，因此有办法区分恢复了正确邻居与仅仅找到了重复本身。

环形结果的一个容易误判之处

长片段组得到完整输出后，直接在原始线性记录上寻找两个预期连接，只找到其中一个。另一处正好跨过输出文件选择的头尾断点，简单字符串搜索就会漏掉。若在这里停住，会把一个完整结果判成少了一条边。正确顺序是先独立确认整条输出能够沿设计圆环连续走满一周，再按这个已经验证的圆周坐标检查跨断点的连接。在这一步，两处预期连接都出现，两个错误交叉连接都没有出现。

这个顺序非常重要。若对任何一条未闭合序列都先把头尾拼起来，程序当然能制造出额外的候选连接，但那个连接没有新增实验依据。短片段结果在脚本中没有通过完整圆周匹配，所以没有获得这种跨断点解释的资格。分析天然样本时，可先检查组装图中的环路、原始读段跨头尾的支持，以及首尾是否有可解释的重叠。若这些证据尚未齐全，就保存线性坐标下确实观察到的结果，并明确哪些邻接关系还在等待确认。

对于已经验证的圆形分子，文件从哪里开始通常还关系到注释展示。相同序列可以旋转到不同起点，导致基因坐标全面变化，甚至让某个开放阅读框在文件头尾被拆开。比较两个版本时应先统一方向与旋转起点，或者使用支持环形坐标的比对，再判断基因是否真的缺失。否则一次纯粹的坐标变换，也会在逐行比较注释表时制造出大量看似新增、删除或截短的记录。记录转换关系能让后续实验定位与可视化使用同一套坐标。

还有一项容易被图形布局影响的细节。组装图查看器为了排版可能把相邻节点画得很远，把无连接的节点画得很近。判断邻接时应读取图中实际记录的边及方向，同时查看节点覆盖和长度，不能只根据屏幕距离推断。如果同一个重复节点被不同路径经过多次，路径记录的次数也应与最终序列中的拷贝数相互核对。保存用于出图的图文件、选定路径和对应序列，可以让其他人从相同图结构重建图中的局部结论。

怎样把这套检查用在自己的数据上

第一步确认文库描述是否准确。保留每对读段的配对标识，质控时成对同步处理，检查第二端是否按正常测序方向提供。不要用一段任意改写的插入长度替代实际读对的分布。片段跨度较长也会带来其他文库性质的变化，因此这次固定跨度的模拟只能帮助理解几何约束。真实文库的跨度分布、错误、嵌合片段和覆盖均匀性，都需要从具体数据重新检查，尤其应关注决定连接的那一小群读对。

第二步先挑研究结论依赖的重复区。若文章要讨论某个耐药基因处于哪个邻近区域，就优先检查连接这个基因与左右独特序列的证据。把每条支持读对的独立起点、方向和比对质量导出，再与候选错误路径比较。读对可以同样合理地放到两条路径时，应记录这种歧义。有足够跨区读对且两端位置明确时，可以收紧解释范围。把全部重复统一标为“已解决”，会把证据丰富的区域与证据不足的区域混到一起。

第三步同时检查序列文件和组装图。前者便于注释与回贴，后者能保留未被输出为唯一长路径的分支。某个目标序列出现在图的两条候选路径上时，记录两种邻接关系，再查看原始数据能否选择其中一种。对必须确认的连接，增加能跨越整个重复和两侧独特区的测序证据，是一个可以明确验收的后续实验目标。补测后仍然使用相同坐标表复核，比较的是同一条连接获得了多少新支持。

最后交付一个简短的审查表，列出输入模式、两端长度、实测跨度分布、最终序列数、未知碱基、不同参考位置覆盖量和重点连接的证据。对于已经闭合的分子，再附方向、旋转起点和末端重叠处理方式。同事收到这样的表后，即使只看到一条序列，也知道它在完整性和结构正确性上各通过了哪些检查，后续注释或实验设计便有明确的坐标依据。

若要对更多分离株逐批执行组装、回贴和连接复核，可据试跑规模与云生信生物分析服务器团队安排批量分析环境。