Mothur2进阶_Mothur扩增子基因序列_数据预处理

2023-12-21 19:52:28

本人在读研究生，方向环境微生物。之前在学习生物信息分析过程中在网络上四处奔走获取相关学习资料与解决问题，好生麻烦。于是，我就把与同学一起做的一些生物信息分析相关教程与经验总结搬运到这个CSDN这个大平台上来，希望能够与大家一起学习讨论。班门弄斧，大神见文多指教，抱拳抱拳抱拳抱拳！

本节重点讲解使用Mothur软件对扩增子序列进行拼接，切除与唯一化等预处理流程。

01数据的拼接

使用make.contigs完成每一个样品两端测序的拼接以及合并所有样本的数据。以stability.files作为输入文件，首先从fastq文件中提取序列和质量得分数据即生成单独的fasta和qual文件，创建反向读取的反向互补，然后连接序列形成contigs，最后在屏幕上会显示每个样本中的序列数。

命令注释：make.contigs用于读取一个正向fastq文件和一个反向fastq文件并输出一个新的fasta文件和report文件。

屏幕输出：

输出文件：

stability.trim.contigs.fasta

stability.scrap.contigs.fasta

stability.contigs.report

文件注释
Name：序列名称；Length：长度；Overlap_length：重叠长度；Overlap_Start：重叠部分开始；Overlap_End：重叠部分结束；MisMatches：不匹配；Num_NsN：碱基数量；Expected_Error：预期错误率

输出文件：stability.contigs.groups

文件注释：第一列为序列名称，第二列为样本名称。

接下来使用summary.seqs命令查看序列，

命令注释：summary.seqs命令用于总结（概述）一个未排序或已排序的fasta格式文件的序列质量。

屏幕输出：

屏幕输出结果显示：有152360个序列，大部分在248和253个碱基之间变化。数据集中最长的读取时间是502 bp。回想一下，每个读取应该是251 bp。说明这些片段显然没有很好地（或根本没有）组装起来。另外，序列中至少有2.5%不明确的碱基被识别。

输出文件：stability.trim.contigs.summary

文件注释
Seqname：序列标识符； Start：起始位置；End：终止位置；Nbases：总碱基数，可以看做序列长度； ambigs：ambiguous bases 模糊碱基的数目；Polymer: homopolymer 碱基的最大长度；NumSeqs: 序列数，对于每条序列来说，其值总是为1

02数据的切除

运行screen.seqs删除所有碱基不明确和长度超过275 bp的序列，

命令注释：该命令用于保留符合用户自定义标准的序列，也可用于从names、group、contigsreport、alignreport、summary文件中挑选出不符合用户自定义标准的序列。

屏幕输出：

输出文件：

stability.trim.contigs.good.fasta

stability.trim.contigs.bad.accnos

文件注释：第一列序列名称，第二列ambig。

输出文件：stability.contigs.good.groups

文件注释：第一列序列名称，第二列组名。

还有另一种方法可以使用summary.seqs的输出来运行，屏幕输出：

该命令会使结果输出速度更快，因为summary.seqs输出文件已经为序列计算了模糊碱基的数量和序列长度。另外，Mothur记住了make.contigs使用了8个处理器，在当前的程序中也使用8个处理器。

运行get.current查看Mothur还记忆了什么，

命令注释：get.current命令用于显示数据处理过程中目前不同类型的最新版文件。

屏幕输出：

从Mothur保存的当前文档可以看出Mothur记住了最新的fasta文件和group文件以及所拥有的处理器数量。

输出文件：current_files.summary

文件注释：由于文件名很长，对于许多命令使用“current”选项通常是最简单的。为了展示如何具体地使用Mothur，所以下面选择性地使用current选项。

运行命令：mothur > summary.seqs(fasta=stability.trim.contigs.good.fasta)

输出结果如上所示，测序错误率可能下降了一个数量级以上并且有128872个序列。

输出文件：stability.trim.contigs.good.summary

名词解释

Processors：电脑处理器线程

Gigabytes：十亿字节

Fastq文件：是一种包含质量值的序列文件，q为quality，一般用来存储原始测序数据，扩展名为fastq或者fq。目前Illumina测序，BGISEQ，Ion Torrent，pacbio，nanopore都以fastq格式存储测序数据，其中Illumina，BGISEQ一般是双末端测序，一般是一对命名为*_R1.fq.gz与*_R2.fq.gz的文件。

Fasta文件：无质量信息的序列，主要用于存储生物的序列文件，例如基因组，基因的核酸序列以及氨基酸等，是最常见的生物序列格式，一般扩展名为fa，fasta，fna等。该文件中，第一行是由大于号“>”开头的任意文字说明，用于序列标记（为了保证后续分析软件能够区分每条序列，单个序列的标识必须是唯一的，序列ID部分可以包含注释信息）。从第二行开始为序列本身，只允许使用既定的核苷酸或氨基酸编码符号。序列部分可以在一行也可以分成多行，序列中允许空格、换行、空行，直到下一个大于号，表示该序列的结束。

03数据的唯一化

数据中的许多序列是彼此重复的，所以使用unique.seqs命令来唯一化序列。

命令注释：unique.seqs命令只会输出fasta格式测序文件中的唯一化序列和一个表征它们与参考序列一致的文件。

屏幕输出：