MitoZ实战：动物线粒体基因组组装注释全流程与避坑指南

1. 做线粒体基因组之前先想清楚这几件事我第一次真正被线粒体基因组组装折磨是在一个涉及上百个昆虫样本的项目里。当时我还在用老思路先拿SOAPdenovo乱跑一版再把contigs拿去BLAST比对近缘物种然后手动拼接、手动补gap、手动把序列卷成环。遇到控制区重复段多一点的物种一整天就耗在一段不到200 bp的序列上。直到我换了MitoZ整个流程从纯手工熬人变成半自动流水线。这篇文章就把我跑MitoZ做动物线粒体基因组组装注释的完整经验写下来从软件逻辑、命令参数到输出文件解读再到各种让人头皮发麻的坑一次性说清楚。MitoZ项目地址在GitHub名字由mitochondrial和Z组成Z代表Z制终点是一款专门面向动物线粒体基因组的组装和注释软件。它最大的特点是一体化把二代测序原始数据喂进去经过组装、优化、成环、注释几个环节最后直接给你一个带完整基因注释的GenBank文件、GFF注释文件和圈图数据。对于做系统发育、种群遗传、DNA条形码研究的人来说这套流程能省掉大量拼接和手工整理的时间。不过我得先泼一盆冷水MitoZ不是魔法。它能把90%的流程自动化但剩下10%的判断力仍然需要你亲自上场。比如线粒体基因组是否混入了核基因组同源片段NUMTs注释出的tRNA边界是否准确控制区重复是否被压缩——这些都得结合覆盖度、BAM比对结果和近缘物种参考序列做人工判断。所以这篇文章也会花大量篇幅讲怎么排查异常因为这才是真正拉开新手和老手差距的地方。1.1 线粒体基因组为什么好组装又不好组装动物线粒体基因组通常是一个12~20 kb的环状DNA分子包含13个蛋白编码基因PCG、22个tRNA、2个rRNA以及一个控制区D-loop / AT-rich region。和核基因组动辄几Gb的规模相比它的体量极小加上单个细胞里有几百上千个线粒体拷贝所以从理论上看它是全基因组里最容易组装出来的部分。但实际跑起来远没那么轻松。第一道坎是测序数据里的杂质我们做线粒体基因组测序时通常提取的是总DNA里面99%以上是核基因组片段线粒体reads占比可能只有0.5%~5%左右。MitoZ这类工具需要先从海量reads里把线粒体reads钓出来这依赖参考种子序列和k-mer频次分析。如果某个物种的线粒体序列和核基因组中插入的线粒体假基因NUMTs相似度太高组装器就可能被干扰最后拼出个核线粒体嵌合体。第二道坎是环状逻辑。绝大多数二代测序平台产出的读长只有150 bp左右软件需要利用k-mer之间的小重叠逐步延伸拼到首尾重叠时判断成环。一旦控制区含有长串联重复单元或极其AT富集比如很多昆虫的AT富集区AT含量超过90%延伸过程容易断掉、错配或者把重复单元数拼错导致最终得到的环长度和真实长度有明显差异。第三道坎是注释。线粒体遗传密码与核基因组不同又因物种分类群而异脊椎动物普遍用遗传密码表2绝大多数无脊椎动物用表5无脊椎动物线粒体密码表或表4后生动物线粒体密码表。选错密码表蛋白编码基因的翻译就全乱套了。MitoZ把密码表验证做进了注释流程但前提是你在跑之前就得选对clade。1.2 传统流程和MitoZ的差异到底在哪在MitoZ出现之前一个标准的线粒体基因组组装注释流程大概是这样用NOVOPlasty或MITObim基于近缘种子进行组装然后用MITOS做初步注释再用tRNAscan-SE单独验证tRNA最后把14个PCG逐一BLAST到近缘物种上人工确认边界和翻译起点。这一套下来顺利的话两三天能出一个样本不顺利的话一两周也正常。MitoZ把上面这套流程写成了流水线。它内部会用SOAPdenovo-Trans和ABySS等组装器做多组k-mer的组装尝试再基于覆盖度和连接关系筛选最优候选接着用片段覆盖图read graph对候选进行延伸、修正、环化。注释阶段则整合了同源比对BLAST against reference database、蛋白域搜索HMMER以及tRNAscan-SE等工具最终自动输出标准GFF3、GenBank、tbl文件。对一个课题组来说这套流程带来的改变不是省了几个小时而是一个人可以从只能处理少量样本变成批量处理几十上百个样本。我后来跑的几百个样本基本全是MitoZ all命令一条龙跑出来的每个样本从原始fastq到最终注释文件常规计算节点上大约需要1~3小时。这个吞吐量在以前完全不可想象。2. MitoZ的安装与三个核心模块的调用逻辑2.1 环境依赖与三种安装方式MitoZ对计算环境的要求不算苛刻但对依赖版本有讲究。它是在Python 2时代开发的虽然也适配了Python 32.4.0以后版本做了适配但如果你直接pip install不指定版本很容易撞上Python 3.10以上版本里某些依赖库的兼容问题。我个人最推荐的是用conda安装简单粗暴环境隔离不容易把系统搞乱。用conda的话一条命令搞定conda create -n mitoz python2.7 conda activate mitoz conda install -c bioconda mitoz如果你需要最新版本或者要改源码调试也可以直接从GitHub拉仓库git clone https://github.com/linzhi2013/MitoZ.git cd MitoZ python setup.py install这里我明确一个容易踩的坑conda安装时一定要把channel优先级配好最好加上bioconda和conda-forge否则某些核心依赖比如perl模块、abyss会安装失败。另外MitoZ部分依赖工具是Perl写的比如tRNAscan-SE的wrapper所以系统里至少要有Perl 5并且确保perl -v能直接调用。安装完成后用mitoz命令检查是否可用。如果你是用conda建的独立环境每次新开终端都要先conda activate mitoz。我第一次没注意这个在base环境里直接敲mitoz报了module not found排查了半天才发现是环境没切换。2.2 模块拆解assembly、annotate、extract和all的关系MitoZ的本质是一套命令集合它按照经典的分步流程把功能拆成了几个子命令。理解每一个子命令在干什么比单纯会敲mitoz all重要得多因为后续排查问题、自定义流程时你得知道卡点在哪一个环节。子命令作用核心产物assembly从原始reads组装出线粒体候选序列组装结果文件夹内的fasta序列annotate对组装好的circle或linear序列做基因注释GFF3、GenBank、tbl、cds、protein序列extract从已有GenBank等注释文件中提取特定序列cds、tRNA、rRNA、基因区间的fastaall一站式完成组装注释上述所有产物汇总extract这个模块容易被忽略但它在实际使用中的价值很高。比如你已经有一个旧的GenBank文件只想重新提取CDS序列做密码子偏好性分析就不用重新注释一遍直接mitoz extract搞定。再比如你手动修正了注释结果也可以通过extract模块重新生成所有导出序列保持文件间的一致性。值得多说一句的是MitoZ在all模式下的组装并不是简单调用单一组装器跑一遍而是内部做了多k-mer。它会对不同k-mer长度的组装结果进行综合评估选出覆盖度最均匀、结构最完整的候选再进入延伸环节。这个设计的隐蔽好处是当你面对不同测序深度、不同GC含量的样本时不需要手动调参就能获得稳定结果。但代价是运行时间会增加尤其当输入数据里有大量核基因组reads时多k-mer的比对和评估阶段会比较吃内存。3. 从二代测序原始数据到完整环状基因组的实操流程3.1 输入文件准备和必设参数MitoZ的assembly模块最少只需要一对双端fastq文件但你得给它提供几个关键参数才能真正跑出好东西。下面是我常用的命令模板以一条脊椎动物样本为例mitoz all \ --outprefix vertebrate_sample \ --thread_number 16 \ --clade Chordata \ --genetic_code 2 \ --data_type Illumina \ --insert_size 350 \ --reads_p1 ./clean_1.fq.gz \ --reads_p2 ./clean_2.fq.gz \ --fastq_read_length 150这些参数看着简单但每一个都值得展开讲。先是--clade。MitoZ的注释阶段会根据clade选择近缘参考数据库注释的灵敏度和准确性直接受这个参数影响。可选的clade包括Arthropoda、Mollusca、Chordata、Echinodermata、Nematoda、Platyhelminthes等基本覆盖了常见动物门类。选大方向的时候就选到门级不要试图填一个非常小众的目级名称反而容易因为没有参考而报错或者注释不完整。然后是--genetic_code遗传密码表编号。这个参数在annotate阶段是决定性的。我要强调一下虽然MitoZ的all流程能自动识别一些情况但不要把命运完全交给默认值。不同分类群的线粒体遗传密码差异很大像棘皮动物用表9、扁形动物用表9或表14选错一个数字蛋白注释的准确率会断崖式下降。跑之前一定先去NCBI Genetic Codes页面查清自己的物种该用哪个表。--insert_size是插入片段平均长度。这个参数会用在reads比对和延伸阶段。很多人图省事随便填一个300但其实最稳妥的做法是用fastqc或者qualimap看清楚插入片段分布或者直接看文库构建记录的估算值。如果填的insert_size和真实文库偏差太远MitoZ在计算成环时的连接关系容易产生误判。最后提一句--data_type。如果你给MitoZ喂的是已经做过trimming的clean reads就填Illumina如果你喂的是PacBio长readsMitoZ也能处理但要填PacBio且参数逻辑完全不同。绝大多数人用的二代短reads所以本文默认围绕Illumina展开。3.2 组装过程发生了什么当命令敲下去之后MitoZ会依次做几件事。首先是对reads进行质量统计和k-mer频率分析这一步会生成一些统计图和一个read_support的中间目录。接着进入多k-mer组装阶段每个k-mer值都会尝试独立组装出contigs。在组装阶段产生的日志里你会看到类似下面的输出kmer41, graph_sizexxx, contigsxxx kmer51, graph_sizexxx, contigsxxx kmer61, graph_sizexxx, contigsxxx不要急着关掉日志——这里藏着一个非常实用的判断依据如果所有k-mer梯度下的contigs都非常破碎碎片数量超过几百个那么大概率是线粒体reads占比太低或者测序质量太差。反之如果某个k-mer下能明显看到大部分contigs被连接成一个大结构就说明这个k-mer比较匹配当前数据的复杂度。组装完成后MitoZ会在输出目录下生成*_assembly文件夹里面包含candidate_1.fasta、candidate_2.fasta等候选序列。命名里的candidate意味着软件自己也没法100%确定哪一条就是你样本真实的线粒体基因组。判断标准有几个序列长度是否符合目标类群线粒体基因组的正常范围、覆盖度是否均匀、是否能环化。MitoZ会在circle.fasta文件里给出被视为成功环化的序列但这个结果仍然值得人工复核。3.3 为什么遗传密码表不能选错遗传密码表的选择是MitoZ流程中试错成本最高的参数。因为组装阶段基本不会用到它很多人在初步组装出完整环之后才跑注释结果注释出大量提前终止的蛋白序列或者出现一堆根本不存在的ORF此时再回去重新注释时间和算力都浪费了。线粒体遗传密码表与核基因组标准密码表差异巨大。举个经典的例子在标准密码表表1里AGA和AGG编码精氨酸但在脊椎动物线粒体遗传密码表表2里这两个密码子是终止密码子而表5无脊椎动物线粒体遗传密码表里AGA/AGG又会被翻译为丝氨酸。如果一只昆虫样本错选了表2那么本应编码丝氨酸的位置会被识别成终止信号导致cox1等基因被拦腰截断。这里我给你一个可以抄作业的选择规则脊椎动物、头索动物、被囊动物表2绝大多数原口动物昆虫、软体动物、环节动物等表5棘皮动物、扁形动物、部分刺胞动物表9后生动物通用兜底表4不确定的时候最简单的办法是先跑一次注释打开蛋白fasta文件检查几个高保守基因比如COX1、CYTB。如果所有基因都有完整的、无明显提前终止的开放阅读框那么密码表大概率选对了如果大量蛋白只有一半长度基本可以确定是密码表的问题。4. 注释模块的产物解读从BED文件到圈图4.1 注释流程内部做了什么MitoZ的annotate模块不是简单地把序列丢给MITOS跑一下。它内部有一个多步骤注释流程主要包括四个层次第一rRNA和tRNA注释。tRNAscan-SE负责识别大部分tRNA同时ARAGORN会做交叉验证。这一步是注释结果里比较稳定的部分绝大多数tRNA都能被准确识别出来。但昆虫等类群中个别tRNA比如tRNA-Ser(AGN)缺少典型的三叶草结构tRNAscan-SE可能报低分甚至漏检这时就需要结合近缘物种参考序列检查该区域是否有遗漏。第二蛋白编码基因注释。MitoZ首先会把组装序列翻译成六个阅读框然后分别与近缘参考蛋白序列库做BLAST比对找到同源区域后再根据遗传密码表和序列保守性逐步修剪出基因的精确边界。对于线粒体基因组这种基因排列极其紧密、甚至存在基因重叠的现象比如ATP8和ATP6经常共用或紧邻核苷酸这个边界确定过程非常关键MitoZ会用HMMER的结构域搜索来辅助确定起始和终止位置。第三控制区和其他非编码区检测。MitoZ会找出一段序列中最长的非编码区域将其默认标注为控制区D-loop或AT-rich region。这个判断在大多数情况下是对的但如果线粒体基因组中存在插入序列或者重复扩张控制区识别可能出现偏差。这个环节我个人建议对照近缘物种的圈图检查一下控制区位置是否合理。第四格式统一输出。做完上述遗留注释MitoZ会生成标准GFF3和GenBank格式文件同时导出cds.fasta、protein.fasta、tRNA.fasta、rRNA.fasta等常用文件。这些标准格式可以无缝进入后续的系统发育分析、选择压力分析和序列提交流程。4.2 输出文件里哪些值得重点关注默认输出目录里文件很多很多刚接触的人会觉得眼花缭乱。下面是我认为最值得关注的几个文件和目录文件内容用途*_annotate.gff3全部基因的坐标和结构注释进IGV/JBrowse可视化、自定义统计*_annotate.gbGenBank格式完整注释提交NCBI、导入Geneious查看*_annotate.tbl表格化特征表提交NCBI的必需格式*_cds.fasta所有蛋白编码基因的CDS序列系统发育和密码子分析*_protein.fasta所有蛋白编码基因的氨基酸序列蛋白树构建*_tRNA.fasta所有tRNA基因序列tRNA结构验证*_circle.fasta环化后的完整线粒体基因组序列作为参考序列*_bam/*_bam.bai原始reads比对到组装结果的bam文件检查覆盖度和嵌合体其中bam文件是我每次必看的。不要只看组装结果漂亮就觉得万事大吉真正能反映组装质量的是reads回帖的覆盖度曲线。一个健康的线粒体基因组组装覆盖度在全环上应该相对一致上下波动通常在一个数量级以内。如果你发现某一段几千bp的覆盖度突然掉到平均值的十分之一那这段很可能在组装中被错误连接或者来自核基因组污染。4.3 怎么把这些产物变成圈图MitoZ本身不直接画圈图但它的注释结果可以直接交给常见绘图工具。我个人用得最顺的是用mitoz extract先导出每个基因的坐标和序列然后用R的gggenes包或circlize包画圈图。如果你更习惯图形界面也可以把gb文件直接拖进Geneious在Generate Circular Map功能里一键生成漂亮的全景圈图。这里有个细节要注意MitoZ输出的GB文件中序列方向是按照组装时的随机起点写的。如果你是本着发表的目的画圈图通常需要把序列的起点旋转到某个保守基因的位置比如COX1的起始密码子或者脊椎动物的tRNA-Phe这会让所有样本的圈图在视觉上具有可比性。做法是在Geneious或序列处理脚本中找到目标基因的起始坐标然后将整个环状序列重新以该坐标作为起点展开最后再重新生成GB文件。注意做完旋转后一定要同步更新GB文件里所有基因feature的坐标否则提交NCBI时会报错。5. 我在数百个样本中踩过的坑与排查思路5.1 组装不出完整circle时的排查链路实际上并不是每个样本都能一次跑出完整环。在我经手的样本中成功率大概在八成左右剩下两成多多少少都有点问题。如果你发现MitoZ最终没有得到circle.fasta或者组装结果是一堆碎片我建议按下面的链路排查不要一上来就盲调参数。第一步检查数据质量。用fastp或fastqc看一下reads的碱基质量分布、接头残留和duplication比例。低质量样本的线粒体reads占比会进一步稀释导致k-mer覆盖度不足组装自然成不了环。第二步检查线粒体reads占比。可以用BWA把reads比对到近缘物种的完整线粒体基因组上看比对率。如果比对率小于0.1%那说明测序深度或富集方式有问题优先考虑补充数据或者富集线粒体。如果比对率偏高但组装仍失败往往说明近缘参考和你目标物种差异太大需要换一个更近的参考。第三步调整k-mer策略。MitoZ默认的k-mer范围覆盖大多情况但你可以通过参数手动指定。当样本GC含量极低或者极高比如AT-rich控制区特别长的昆虫样本适当降低k-mer值有助于提高含重复区域的延伸能力。我一般在默认跑不圆的时候会加一次--kmer 31/41/51/61的手动尝试。第四步检查是否有其他物种污染。双端测序数据里偶尔会混入共生菌、寄生蜂或者饲料生物的reads。如果这些污染reads高到一定程度组装出的候选序列可能不是目标物种的线粒体。此时用NCBI BLAST比对组装结果看最佳匹配是否真的是同属或同科物种是快速验证的手段。5.2 NUMTs和异质性导致的嵌合体这个坑比较隐蔽也是我认为MitoZ这类自动流程最需要人工介入的地方。NUMTsNuclear Mitochondrial DNA segments核线粒体假基因是核基因组中插入的线粒体片段。在部分物种中这些假基因可能很长甚至达到数kb。由于它们和真实线粒体序列相似性很高如果测序数据里NUMTs reads占比较高组装器可能把线粒体序列和NUMT拼接在一起产出一个看似完整、实则嵌合的环。怎么识别这种嵌合体我的经验是看覆盖度曲线。真实线粒体基因组在所有位置都会维持一个稳定的高覆盖度但NUMT区域的reads在核基因组里是单拷贝的覆盖度会明显低于旁边来自线粒体的区域。打开bam文件在IGV里看如果某段覆盖度出现断崖式下跌就要高度怀疑嵌合。另一种办法是拿组装完成的序列重新BLAST到近缘物种的核基因组上检查是否有大片段的同源区域。此外线粒体异质性heteroplasmy也会造成类似现象。异质性指的是同一个体内存在多个版本的线粒体基因组本该在某个位点上出现的一致性碱基变成了混合信号。轻微异质性影响不大但如果两个版本差异太大组装器可能产生气泡bubble导致错误连接。MitoZ在评估时会倾向于选择覆盖度更高的分支但这不代表低频分支完全没有生物学意义——如果你的研究课题恰好关注异质性就不能完全依赖自动流程而是要专门对变异位点做深度验证。5.3 近缘物种污染与数据筛选做线粒体基因组批量项目时最让人崩溃的不是单个样本跑不出来而是跑出来的样本序列和数据库里某个近缘物种一模一样。这种情况通常是样本间交叉污染或者物种鉴定本身出了问题。交叉污染的一个常见来源是建库时使用了同一个index的接头或者测序仪flowcell上的lane发生泄漏。判断方法很简单把疑似样本的原始reads比对到近缘物种A和B的线粒体上看reads是否随机分布。如果大部分reads只比对上其中一个物种的特定区间而不是完整覆盖整个环那大概率是污染reads混入而不是真实杂交个体。MitoZ里没有专门的去污染模块但你可以通过上游过滤来减少影响。最直接的办法是用fastp按overrepresented sequence过滤高丰度非目标序列或者用Kraken2等分类工具对reads做物种分类保留目标物种reads后重新组装。有一点要清楚如果污染reads占比超过30%MitoZ的k-mer频率分析可能会把污染序列当作高支持度的候选这时要先清洗数据再跑。5.4 注释缺失tRNA或假阳性的应对MitoZ自动注释的tRNA准确率总体不错但总会碰到个别案例某个tRNA没被识别出来或者识别出的tRNA在反密码子位置明显异常。这主要有两个原因一是该类群某个tRNA确实存在结构变异自动工具训练集覆盖不到二是该区域序列太短、太AT富集导致score达不到阈值。遇到这种情况我的处理方法是先看近缘物种的注释文件尤其是同属或同科的GenBank确认该物种是否也有相同的tRNA缺失或结构变异。如果近缘物种正常那就需要针对该区域手动检查。你可以用ARAGORN单独跑一次或者把该区域序列丢到MITOS的Web Server做交叉验证。如果三种工具都不认而你又确信该区域是tRNA比如反密码子清晰、能折叠出三叶草结构可以在GB文件中手动添加feature然后重新导出GFF3。这里要特别提醒手动更改注释后一定要把起始和终止坐标的补链complement方向核对清楚。线粒体基因组两条链都编码基因tRNA或蛋白基因可能在轻链或重链上坐标方向一旦搞反后续所有基于该注释的分析都会出错。6. 结果可靠性检验与后续必备操作6.1 组装注释完成后必须做的几项验证即使MitoZ跑得再顺我也不会直接拿结果去提交数据库或做分析。每个样本我都会做一遍系统性的验证这几项加起来大概半小时到一小时但能省下后面几个月返工的麻烦。首先是reads回帖覆盖度检查这一步参照前面说的bam文件分析。其次是PCR验证在环上随机选2~3个跨区域位置其中一个尽量跨厨房和编码区边界设计引物做常规PCR并Sanger测序和二代组装结果比对。对于线粒体这种多拷贝结构跨区域PCR的一致性能有效排除NUMTs带来的嵌合风险。然后是与近缘物种的共线性分析。线粒体基因排列在昆虫等类群中极为保守同一科内的物种基本共享同一套基因排列顺序。如果你的注释结果里出现了基因相对位置颠倒、基因缺失或重复的现象首先应该怀疑组装或者注释有问题其次才考虑真实的基因重排事件。用MAUVE或Geneious的Annotate Predict和近缘物种做共线性比较是判断注释可靠性的好办法。最后是蛋白编码基因完整性检查。把所有CDS翻译成蛋白逐一确认起始密码子、终止密码子和氨基酸长度。比如COX1在多数动物中长度在510~520个氨基酸左右CYTB在380个氨基酸左右。如果某个样本的COX1只有200个氨基酸除非有真实生物学变异否则大概率是注释出了问题。MitoZ的蛋白文件里如果出现大量*提前终止的序列基本就是密码表选错了别犹豫回去改密码表重跑注释。6.2 把线粒体基因组数据提交到公共数据库前的准备工作很多人的最终目的是发表文章并把线粒体基因组提交到NCBI GenBank或ENA。MitoZ输出的GB文件和tbl文件已经比较接近提交格式了但直接用自动产物提交经常会收到数据库审核员的修改意见。我总结了几条高频退回原因一是起始位点选择不符合要求。NCBI允许你以任意位点为起点但要保证注释中所有feature都完整且坐标正确。很多审核员建议把CDS的起始密码子作为序列起点这样显得规范。具体做法前面提过就是旋转序列再重新生成GB文件。二是gene和product字段填写不规范。比如线粒体基因的gene名应按标准命名cox1、nad2、atp6等product字段要写清楚是cytochrome c oxidase subunit I而不是缩写。三是来源信息source modifier不完整。提交时需要填写isolate、collection_date、country、host等来源信息。MitoZ不会自动帮你填这些需要你通过BankIt或表格录入。四是重复序列或控制区的注释。控制区有时会标注为D-loop或AT-rich region具体用哪个词汇要看你所在的分类群习惯。对昆虫来说AT-rich region更常见对哺乳动物D-loop更通用。提交前打开GB文件检查一下避免低级错误。6.3 从MitoZ结果到下游分析的无缝衔接MitoZ真正让我觉得香的地方是它可以无缝衔接下游分析。有了标准格式的CDS序列和蛋白序列我可以直接跑IQ-TREE建树、用PAML算选择压力、用DnaSP统计多态性不需要再花时间做格式转换。对于批量项目我建议用脚本把样本的cds.fasta集中起来拼接成concatenated matrix。需要注意的是不同物种之间基因顺序可能不同拼接前先用同源基因序列做ortholog配对再按固定顺序拼接避免出现样品间基因错位。MitoZ导出的基因名基本遵循标准命名这为脚本处理提供了极大便利。另外如果你做的是比较线粒体基因组学建议在MitoZ输出后用GBDP等工具先做一次全基因组相似度计算确认各样本之间的遗传距离合理再决定后续的分析策略。一次大批量的组装注释之后用一个简单的QC脚本把所有样本的序列长度、GC含量、基因数量汇总成一张表是我强烈推荐的好习惯——这张表能在你后期汇总数据、写文章方法部分的时候省下大量时间。说到底MitoZ是一个把重复劳动自动化程度拉满的工具但它取代不了研究者对每个样本数据的审视。我的习惯是批量跑完流程后抽10%的样本做人工全流程复核包括打开bam文件看覆盖度、逐个基因检查蛋白长度、用近缘物种比对验证基因排列。至今这个习惯帮我拦截了至少三五个会毁掉整篇文章的嵌合体样本。刚开始接触MitoZ的朋友我建议你手里留一份已发表近缘物种的线粒体基因组做参考答案不管是组装还是注释阶段遇到拿不准的情况就比一比。这套工具人工复核的打法目前在我处理过的大量动物样本中一直很稳定。