病原微生物基因检测的两大核心任务是物种组成和功能组成的鉴定，而扩增子测序的首要目的是找到致病的细菌或者病毒，即鉴定物种组成。
经过质量控制，测序数据中已经不包含非生物的异常序列了，下面我们用vsearch软件完全后续分析。
1.去重（dereplication）
同一对引物的扩增产物，理论上应该是完全一样的，这些冗余的信息会造成比较大的运算负担，因此需要先去冗余，相同的序列只保留一条就好了。
vsearch去冗余有以下三种模式：
-
全长模式：两条序列长度完全相同，碱基完全一样，总之，序列要完全一模一样才认为是重复
-
一致模式：在全长模式的基础上，进一步要求序列的名称要完全相同
-
前缀模式：只要一条序列是另一条序列的前缀，也就是一条序列与另一条序列的前面部分相同，就可以认为它们是重复。
可以看出，这三种模式的严格程度依次是：一致模式 > 全长模式 > 前缀模式。
我们用全长模式去除冗余：
|
|
--derep_fulllength，后面跟通过质控的fastq文件；
--output，输出文件，fasta格式；
--minuniquesize，最低丰度值，低于该丰度的序列会被过滤掉；
--minseqlength，最低长度值，低于该长度的序列会被过滤掉；
--strand，当判断两条序列是否一致时，默认只考虑正链plus，both表示考虑正反两个方向都考虑；
--sizeout，在结果文件中序列名称后面添加丰度信息；
--fasta_width，限定fasta结果文件中每条序列在一行中最多显示的字符数，默认是80，0表示不做限制；
2.降噪（denoise）
按97%相似度对序列进行聚类曾经是扩增子序列分析的金标准，但这有一个问题，就是物种只能鉴定到属或种，因为有些物种之间的相似度接近99%，甚至有的菌株之间只有SNP的差别。因此现在流行的做法是进行降噪处理，通俗地说就是按100%相似度进行聚类，从而使鉴定达到种或株的水平。
|
|
--cluster_unoise，上一步去重后的fasta文件；
--centroids，fasta结果文件，包含每一个聚类中的种子序列；
--consout，fasta结果文件，包含每一个聚类的一致性序列；
--minsize，降噪最低的序列丰度要求，默认是8条；
--strand，当判断两条序列是否一致时，默认只考虑正链plus，both表示考虑正反两个方向都考虑；
--sizeout，在结果文件中序列名称后面添加丰度信息；
--fasta_width，限定fasta结果文件中每条序列在一行中最多显示的字符数，默认是80，0表示不做限制；
--clusterout_sort，结果文件中序列的顺序默认是按其在输入文件中的顺序，设定该参数则是按照降噪后序列的丰度排序；
3.去嵌合体
测序数据中可能存在嵌合体序列，每个嵌合体至少来源于两个或两个以上的扩增子模板，因而需要事先去除。
去嵌合体有两种方法：一是基于参考数据库，一是基于无参考数据库的denovo方法。
vsearch基于denovo方法去除嵌合体的命令：
|
|
--uchime3_denovo，上一步去重后获得的一致性序列文件；
--nonchimeras，不含有嵌合体的结果文件；
--chimeras，含有嵌合体的结果文件；
--uchimealns，以人类易于阅读的形式呈现嵌合体与其两个亲本进行比对的结果文件；
--sizeout，在结果文件中序列名称后面添加丰度信息；
--fasta_width，限定fasta结果文件中每条序列在一行中最多显示的字符数，默认是80，0表示不做限制；
4.创建OTU表
OTU（operational taxonomic units），即可操作分类单元，通过vsearch创建OTU表，可以看到每个OTU的丰度情况。
a.拷贝OTU文件：经过去重、降噪、去嵌合体之后，获得的无嵌合体的fasta文件可作为最终OTU结果
|
|
b.格式转换：将fastq格式的质控文件转换成fasta格式，以下命令只做格式转换，没有进行任何过滤操作
|
|
c.创建OTU表：本质上是将经过质控的序列比对到OTU序列上，从而可以得出每一个OTU序列的丰度
|
|
--usearch_global，后面跟经过质控的fastq文件，fasta格式；
--db，OTU文件，fasta格式；
--id，相似度阈值：当查询序列与目标序列之间的相似度达到多少时，才算比对上；
--query_cov，覆盖度：满足相似度的情况下，同时要求查询序列的覆盖度达到多少；
--strand，当判断两条序列是否一致时，默认只考虑正链plus，both表示考虑正反两个方向都考虑；
--outtabout，输出OTU表。
--samheader，当要保存SAM文件时，默认不包含SAM文件头，设定此参数即会包含；
--samout，输出比对的SAM格式结果文件，可用于后续IGV软件可视化查看比对情况。
5.物种鉴定
有了样本的代表性序列后，想知道样本的物种组成，可将这些序列比对到参考数据库，通过相似性及覆盖度进行判断。
|
|
--usearch_global，OTU文件，fasta格式；
--db，参考序列库，fasta格式；
--id，相似度阈值：当查询序列与目标序列之间的相似度达到多少时，才算比对上；
--query_cov，覆盖度：满足相似度的情况下，同时要求查询序列的覆盖度达到多少；
--strand，当判断两条序列是否一致时，默认只考虑正链plus，both表示考虑正反两个方向都考虑；
--biomout ，输出biom格式的结果文件；
--alnout，输出比对结果文件；
--blast6out，输出blast 6格式的比对结果文件；
--fastapairs，fasta格式文件，保存了查询序列以及对应的目标序列；
--notmatched，没有搜索到目标序列的查询序列集合；
--userfields，自定义输出结果的列名，从左至右分别为：查询序列id，目标序列id，相似度，查询序列覆盖度，目标序列覆盖度；
--userout，按–userfields定义的表头输出自定义的结果文件。
上述命令得到一个自定义的结果文件${outdir}/${sample}_stat.xls，根据其中的相似度、覆盖度以及前面的OTU表${outdir}/${sample}_otutab.txt，就可以知到物种的组成情况，以及它们的丰度信息。
现在，通过3篇文章，一套病原微生物扩增子物种鉴定的流程就建立起来了。
往期精彩