Trimmomatic安装和使用
Trimmomatic安装和使用
本文作者：Sunny-King
发布时间：2022-06-14
本文链接：https://www.cnblogs.com/Sunny-King/p/Bioinformatics-Trimmomatic.html
2022-06-14
Trimmomatic: A flexible read trimming tool for Illumina NGS data
Bolger, A. M., Lohse, M., & Usadel, B. (2014). Trimmomatic: A flexible trimmer for Illumina Sequence Data. Bioinformatics, btu170.
一、Trimmomatic简介
官网：http://www.usadellab.org/cms/?page=trimmomatic
官方文档：TrimmomaticManual_V0.32.pdf
Github网址：https://github.com/usadellab/Trimmomatic
Trimmomatic是一个快速的多线程命令行工具，可以用来处理Illumina (FASTQ)数据以及删除adapter。程序有两种主要模式:Paired-End模式和Single-End模式。Paired-End模式将保持read-pairs的对应关系，并利用read-pairs所包含的附加信息，更好地找到文库制备过程中引入的adapter或PCR引物片段。支持使用gzip或bzip2压缩的文件，并通过使用.gz或.bz2文件扩展名来识别。
Trimmomatic是基于Java开发的，因此需要提前安装Java，才能使用Trimmomatic。目前的Trimmomatic主要包含下面步骤
- ILLUMINACLIP: 从read中切除adapter和Illumina的特殊序列。
- SLIDINGWINDOW: 从read的5'端开始滑动窗口，当窗口内的average quality小于阈值时，切除后边的碱基.。
- MAXINFO: 自适应的trim策略，可以平衡read的长度和错误率保证read的长度最长。
- LEADING: 切除read的开头质量低于阈值的碱基。
- TRAILING: 切除read的末尾质量低于阈值的碱基。
- CROP: 从read的末尾切除特定长度。
- HEADCROP: 从read的开头切除特定长度。
- MINLEN: 舍弃小于最小长度的reads 。
- AVGQUAL: 舍弃平均质量 小于阈值的reads。
- TOPHRED33: 将质量分数转换为Phred-33
- TOPHRED64: 将质量分数转换为Phred-64
二、Trimmomatic安装(trimmomatic-0.39)
1、下载安装
下载编译好的二进制文件
#下载
wget http://www.usadellab.org/cms/uploads/supplementary/Trimmomatic/Trimmomatic-0.39.zip
#解压
unzip Trimmomatic-0.39.zip
解压后的文件包括以下几个：
Trimmomatic-0.39
├── adapters
│ ├── NexteraPE-PE.fa
│ ├── TruSeq2-PE.fa
│ ├── TruSeq2-SE.fa
│ ├── TruSeq3-PE-2.fa
│ ├── TruSeq3-PE.fa
│ └── TruSeq3-SE.fa
├── LICENSE
└── trimmomatic-0.39.jar
1 directory, 8 files
- adapters路径下的文件是ILLUMINA常规的接头文件，有单端和双端。
- trimmomatic-0.39.jar是编译之后的文件，需用java执行。
$ java -jar Trimmomatic-0.39/trimmomatic-0.39.jar -h
Usage:
PE [-version] [-threads <threads>] [-phred33|-phred64] [-trimlog <trimLogFile>] [-summary <statsSummaryFile>] [-quiet] [-validatePairs] [-basein <inputBase> | <inputFile1> <inputFile2>] [-baseout <outputBase> | <outputFile1P> <outputFile1U> <outputFile2P> <outputFile2U>] <trimmer1>...
or:
SE [-version] [-threads <threads>] [-phred33|-phred64] [-trimlog <trimLogFile>] [-summary <statsSummaryFile>] [-quiet] <inputFile> <outputFile> <trimmer1>...
or:
-version
2、conda安装
conda install trimmomatic
三、Trimmomatic使用
1、SE模式
对于单端数据，指定一个输入文件和一个输出文件。所需的处理步骤(裁剪、适配器裁剪等)作为附加参数指定在输入/输出文件之后。
java -jar <path to trimmomatic jar> SE [-threads <threads>] [-phred33 | -phred64] [-trimlog
<logFile>] <input> <output> <step 1> ...
指定trimlog文件将创建一个记录所有reads 的日志，指出以下细节:
- read的名字，即QNAME
- 保留下的序列长度
- 从起始位置，第一个保留碱基的位置。
- 原始read中最后一个碱基的位置
- 从末尾切除掉的碱基个数
可以根据需要指定多个步骤，方法是在部分处理步骤中描述的在最后使用附加参数。
2、PE模式
对于Pair-End数据，需要两个输入文件，会输出4个文件。其中两个文件对应于"paired"数据，即read1和read2都保留。另外两个对应于"unpaired"数据，在处理的过程中会过滤掉其中一端的reads。
java -jar <path to trimmomatic.jar> PE [-threads <threads] [-phred33 | -phred64] [-trimlog <logFile>] >] [-basein <inputBase> | <input 1> <input 2>] [-baseout <outputBase> | <unpaired output 1> <paired output 2> <unpaired output 2> <step 1> ...
- -phred33|-phred64：指定碱基质量编码规则。在0.32版本以后如果缺省，则会根据数据自动选择。默认值是-phred64。
- -threads：表示要使用的线程数，用于提高多核计算机的性能。
## 示例
java -jar trimmomatic-0.30.jar PE \
# 输入文件
Sample_R1.fq.gz Sample_R2.fq.gz \
# 输出文件
Sample_Trim_paired_R1.fq.gz Sample_Trim_unpaired_R1.fq.gz \
Sample_Trim_paired_R2.fq.gz Sample_Trim_unpaired_R2.fq.gz \
# 去接头参数
ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 \
# 切除TruSeq3-PE中提供的Illumina适配器。最初Trimmomatic将寻找种子匹配(16个碱基)，最多允许2个不匹配。如果在配对端读取的情况下达到30分(约50个碱基，50*0.6)，或在单端读取的情况下达到10分(约17个碱基, 17*0.6)，这些"seed"序列将被修剪。
LEADING:3 \
# 删除前低质量碱基(低于质量3)
TRAILING:3 \
# 删除后低质量碱基(低于质量3)
SLIDINGWINDOW:4:15 \
# 扫描4个碱基宽滑动窗口，当每个碱基的平均质量下降到15以下时进行剪切
MINLEN:36
# 上述步骤完成后，删除小于36个碱基的reads
四、详细的处理细节
1、ILUMINACLIP
这一步用于查找和切除Illumina adapters。理论上虽然adapters和其他实验产生的序列可能位于在reads的任何位置，到目前为止，adapters污染最常见的原因是测序的DNA片段短于读取长度，即segment < read_length。在这种情况下，reads的开头包含有效的数据，但是reads的末尾会读取到到adapters。这将导致在reads的3'端出现部分或全部adapters序列。虽然识别完整的adapters序列比较容易，但识别短的部分adapters序列本质上是困难的。
有趣的是，在Pair-End的数据中，由于reads完全从两端开始，前向和反向reads的非adapters部分将是反向互补。对于Illumina的Pair-End的数据，Trimmomatic包含了第二种adapters识别策略，专门用于adapters测通的情况，并利用Pair-End数据中的证据。这种策略被称为回文(palindrome)模式。
下面的图说明了两种不同的策略。
- A：reads中包含读数据中的整个技术产生的序列，因此标准比对方法足以准确识别并切除。
- B：只有部分的技术序列包含在reads的3'端，因此只能使用短的对齐。这种情况会平衡假阳性和假阴性来进行取舍，对于一些长度低于阈值的的adapters序列片段可能会残留下来。
- C：回文模式检测长的reads。在本例中，根本没有任何有用的片段，并且两条reads都从adapters的序列开始。
- D：使用回文方法来，在reads的3'端有一个非目标序列。然而，由于利用了在配对端模式下可用的额外数据，作为校准的一部分进行比较的区域要长得多。不仅一次比较了两个adapters序列，而且还检查每条reads序列。因此，这种对齐比B中的短对齐要可靠得多，可以去除仅剩有1个碱基的adapters序列。
Trimmomatic使用两步来识别adapters和reads。首先，在reads的每个可能位置测试每个adapters的短段(最大16bp)。这种匹配成为"seed",由seedMismatch参数调控。如果"seed"接近于完全匹配，则adapters和reads之间将计算一个full alignment score。这种两步策略带来了相当大的效率增益，因为"seed"计算较快，同时减少full alignment score的计算。
full alignment score计算方法如下。每匹配一个碱基，评分增加0.6，而每错配一次，比对评分减少Q/10。通过考虑碱基的的质量，读取错误导致的miss-match的影响较小。一个完美匹配的12个碱基序列的score在7个以上，而需要25个碱基才能得到15。因此，我们建议在7 - 15之间的值作为简单对齐模式的阈值。对于回文匹配，可以采用更长的对齐方式，如上所述。因此，在30的范围内这个阈值可以更高。即使这个阈值非常高(需要匹配近50个碱基)，Trimmomatic仍然能够识别非常、非常短的adapters片段。
使用参数
ILLUMINACLIP: <fastaWithAdaptersEtc> :<seed mismatches>:<palindrome clip threshold>:<simple clip threshold>
ILLUMINACLIP:<fastaWithAdaptersEtc>:<seed mismatches>:<palindrome clip threshold>:<simple clip threshold>:<minAdapterLength>
ILLUMINACLIP:<fastaWithAdaptersEtc>:<seed mismatches>:<palindrome clip threshold>:<simple clip threshold>:<minAdapterLength>:<keepBothReads>
- fastaWithAdaptersEtc：指定包含所有adapters、PCR序列等fasta文件的路径。
- seedMismatches：指定"seed"匹配允许的最大miss-match。
- palindromeClipThreshold：指定PE回文读取对齐时，两个“adapter ligated”读取之间的full alignment score阈值。
- simpleClipThreshold:指定任何adapter序列之间的full alignment score阈值。
- minAdapterLength:除了对齐评分外，回文模式还可以验证已检测到adapter的最小长度。如果未指定，则由于历史原因，默认为8个碱基。然而，由于回文模式有一个非常低的假阳性率，这可以安全地减少，甚至减少到1，以允许更短的adapter片段被删除。
- keepBothReads:在通过回文模式检测到读透后，删除适配器序列，反向读取包含与正向读取相同的序列信息，尽管是反向补码。因此，默认行为是完全删除反向读取。通过为该参数指定true，反向读取也将被保留，这可能会很有用，例如，如果下游工具不能处理成对和非配对读取的组合。一般设置为true。
2、SLIDINGWINDOW
执行滑动窗口修剪，一旦窗口内的平均质量低于阈值就进行切割。通过考虑多个碱基，避免因为一个低质量的碱基导致删除后面的reads中高质量的数据。
使用参数
SLIDINGWINDOW:<windowSize>:<requiredQuality>
- windowSize：滑动窗口的大小
- requiredQuality：窗口内的平均质量。
3、MAXINFO
执行自适应的高质量调整，平衡保持较长reads的好处与保留错误的基础的成本。对于许多应用程序，一次读取的值是三个因素之间的平衡
使用参数
MAXINFO:<targetLength>:<strictness>
- targetLength: This specifies the read length which is likely to allow the location of the read within the target sequence to be determined.
- strictness: 这个值应该设置在0到1之间，它指定尽可能多地保留读长度与删除不正确的碱基之间的平衡。该参数的低值(<0.2)有利于读取较长的数据，而高值(>0.8)有利于读取的正确性。
4、LEADING
从起始位置移除低质量的碱基。当碱基的值低于这个阈值，这个碱基就会被移除，然后研究下一个碱基。
使用参数
LEADING:<quality>
- quality:最低质量。
5、TRAILING
从末端移除低质量的碱基。只要一个碱基的值低于这个阈值，这个碱基就会被移除，并且从下一个碱基开始处理(由于从3 '端开始，下一个碱基其实是在刚刚移除的碱基之前)。这种方法可以用来去除特殊的illumina低质量段区域(质量分数为2)，但我们推荐滑动窗口或MaxInfo代替。
使用参数
TRAILING:<quality>
- quality:最低质量。
6、CROP
无论读取的质量如何，都将从读取的末尾移除碱基，以便在执行此步骤后读取最大限度地达到指定的长度。在CROP之后执行的步骤可能会进一步缩短reads。
使用参数
CROP:<length>
- length:从reads开始开始所要保留的碱基数。
7、HEADCROP
无论质量如何，从读取开始移除指定数量的碱基。
使用参数
HEADCROP:<length>
- length:从reads起始所要移除的碱基数。
8、MINLEN
此模块删除低于指定最小长度的reads。它通常在所有其他处理步骤之后。这一步删除的读数将被计数，并包括在trimmomatic摘要中。
使用参数
MINLEN:<length>
- length:指定要保留的读取的最小长度。
本文作者：Sunny-King
本文链接：文章来源于博客园 https://www.cnblogs.com/Sunny-King/p/Bioinformatics-Trimmomatic.html
转载要求：欢迎转载，转载之后请务必在文章明显位置标出原文链接和作者
错误修复：如有错误或疑问请联系博主
版权声明：本作品采用署名－非商业使用-禁止演绎 (by-nc-nd)许可协议进行许可
如果本文对您有帮助，请点个赞吧！志同道合的朋友可以点个关注。
浙公网安备 33010602011771号