Trimmomatic的一些使用说明整理

一、前言

关于这篇博客只是网上调研的结果，博主本人并没有对这个软件进行过实际操作，纯属纸上谈兵，更多的只是对官网使用工具书的翻译，且不全，没有涉及具体算法，仅涉及用法，如有出错或侵权，请联系abel.k@qq.com，谢谢。

后期博主有使用的话，会持续更新本文，感谢关注。

二、软件介绍

NGS 原始数据过滤对后续分析至关重要，去除一些无用的序列也可以提高后续分析的准确率和效率。

Trimmomatic是一个快速的多线程命令行工具，可用于修剪和裁剪Illumina（FASTQ）数据以及移除接头。根据文库制备和下游应用程序的不同，这些接头可能会带来真正的问题。

程序有两种主要模式：双端模式和单端模式。双端模式将保持read成对的对应关系，并且还使用成对reads中包含的附加信息来更好地找到由文库制备过程引入的接头或PCR引物片段。

Trimmomatic使用FASTQ文件（根据使用的Illumina通道，使用phred-33或phred-64质量分数）。支持使用gzip或bzip2压缩的文件，并通过使用.gz或.bz2文件扩展名进行标识。

下面是使用该工具时的一些说明。

三、相关链接

1、Trimmomatic软件官网
2、Trimmomatic使用手册
3、Trimmomatic参考文献

四、软件下载

点击Trimmomatic软件官网，在Downloading Trimmomatic目录下，选择一个版本，点击binary，即可下载一个名trimmomatic-x.xx(x.xx为版本号，本文编写时，已有0.39版本)的ZIP压缩包，解压后得到trimmomatic-x.xx.jar，随后各种操作都是用java调用这个jar包。

五、软件使用

1、修剪步骤说明

Trimmomatic为illumina对端和单端数据执行各种有用的修剪任务。修剪步骤及其相关参数的选择在命令行中提供。

正确的修剪步骤:

1)ILLUMINACLIP: 切除read中的接头以及Illumina特异序列；
2)SLIDINGWINDOW: 划窗修剪方法。它从5'端开始扫描，当窗口内的平均质量低于阈值时，它会剔除该窗口内的所有碱基；
3)MAXINFO: 自适应质量微调器，它平衡读取长度和错误率，最大化每条read的价值；
4)LEADING: 切除read起始端低于阈值的碱基；
5)TRAILING: 切除read末端低于阈值的碱基；
6)CROP: 切除read末端指定数量的碱基；
7)HEADCROP: 切除read起始端指定数量的碱基；
8)MINLEN: 丢弃低于指定长度的read；
9)AVGQUAL: 丢弃平均质量低于指定质量的read；
10)TOPHRED33: 转换质量分数为Phred-33；
11)TOPHRED64: 转换质量分数为Phred-64。

2、运行程序

不同的处理步骤按照在命令行上指定步骤的顺序进行。在大多数情况下，如果需要的话，建议尽早进行接头剪切，因为使用修剪后部分来正确识别接头会比较困难。

（1）单端模式

对于单端数据，指定一个输入和一个输出文件。所需的处理步骤（修剪、剔除、接头剪切等）是被指定为输入/输出文件之后的附加参数。命令模式如下：

java -jar <path to trimmomatic jar> SE 
	[-threads <threads>] 
	[-phred33 | -phred64] 
	[-trimlog <logFile>] 
	<input> <output> 
	<step 1> 
	<step 2>
	...

或者

java -classpath <path to trimmomatic jar> org.usadellab.trimmomatic.TrimmomaticSE 
	[- threads <threads>] 
	[-phred33 | -phred64] 
	[-trimlog <logFile>] 
	<input> <output> 
	<step 1> 
	<step 2>
	...

（2）双端模式

对于双端测序数据，指定了2个输入文件和4个输出文件，其中2个用于“成对”输出（数据处理后，read1，read2同时存在），2个用于“未成对”输出（数据处理后，read1或read2缺失）。命令模式如下：

java -jar <path to trimmomatic.jar> PE 
	[-threads <threads] 
	[-phred33 | -phred64] 
	[-trimlog <logFile>] >] 
	[-basein <inputBase> | <input 1> <input 2>] 
	[-baseout <outputBase> | <unpaired output 1> <unpaired output 1> <paired output 2> <unpaired output 2> 
	<step 1> 
	<step 2>
	...

或者

java -classpath <path to trimmomatic jar> org.usadellab.trimmomatic.TrimmomaticPE 
	[-threads <threads>] 
	[-phred33 | -phred64] 
	[-trimlog <logFile>] 
	[-basein <inputBase> | <input 1> <input 2>] 
	[-baseout <outputBase> | <paired output 1> <unpaired output 1> <paired output 2> <unpaired output 2> 
	<step 1> 
	<step 2>
	...

（3）输入/输出文件

双端模式需要2个输入文件（用于正向和反向reads）和4个输出文件（用于正向配对、正向未配对、反向配对和反向未配对reads）。由于这些文件通常具有相似的名称，因此用户可以选择提供单独的文件名，或者只提供一个可以从中派生文件名的名称。

对于输入文件，以下两种情况任选一种：
1)明确2个输入文件的命名；
2)使用-basein标志命名正向文件，这样可以自动确定反向文件名称。第二个文件是通过查找文件命名的常见模式，并更改适当的字符来命名反向文件。正确处理的示例包括：

a)Sample_Name_R1_001.fq.gz —— Sample_Name_R2_001.fq.gz
b)Sample_Name.f.fastq —— Sample_Name.r.fastq 
c)Sample_Name.1.sequence.txt —— Sample_Name.2.sequence.txt

对于输出文件，以下两种情况任选一种：
1)明确4个输出文件的命名；
2)使用-baseout标志提供基本文件名，从中派生四个输出文件。如果提供的基本文件名是 mySampleFiltered.fq.gz，则将使用以下文件名:

a)mySampleFiltered_1P.fq.gz - 正向配对reads；
b)mySampleFiltered_1U.fq.gz - 正向未配对reads；
c)mySampleFiltered_2P.fq.gz - 反向配对reads；
d)mySampleFiltered_2U.fq.gz - 反向未配对reads。

（4）注解

对于两个模式，都有：
-phred33或-phred64指定碱基质量编码。如果未指定质量编码，则将自动识别（从0.32版起）。之前版本默认的是-phred64。

命令中[-trimlog <logFile>]会生成一个指定的trimlog，其中含有以下信息：
1)read名称；
2)剪切后剩余的read长度；
3)剪切后第一个碱基的所在位置，也就是从起始位置开始剪切了几个碱基；
4)剪切后最后一个碱基在原始read中位置；
5)末端剪切碱基的数量。
注：由于该文件较大，如后续步骤不需该文件信息，可考虑不设置。

可以根据需要指定多个步骤，方法是在命令行末尾使用附加参数，参数如“修剪步骤”一节中所述，附加参数直接添加在 <step 1> <step 2>处，一个参数一个step。
对于输入和输出文件，将.gz/.bz2添加到扩展名，让程序知道文件是以gzip/bzip2格式提供的，或者让程序分别处理不同格式的文件。

3、参数详情

（1）详情说明

大多数处理步骤采用一个或多个设置，由“:”分隔。

1)ILUMINACLIP:<fastaWithAdaptersEtc>:<seed mismatches>:<palindrome clip threshold>:<simple clip threshold>:<minAdapterLength>:<keepBothReads>

该步骤用于寻找并去除Illumina接头。
a)fastaWithAdaptersEtc1: 指定包含所有接头、PCR序列等的fasta文件的路径。此文件中各种序列的命名决定了它们的使用方式；
b)seed mismatches：指定仍允许执行完全匹配的最大不匹配数；
c)palindrome clip threshold：指定两个成对接头reads之间的匹配对于双端回文read对齐的精度；
d)Simple clip threshold：指定任何接头等序列与read之间的匹配精度；
e)minAdapterLength：除了对齐分数之外，回文模式还可以验证检测到接头的最小长度。如果未指定，出于历史原因，默认为8个bp。但是，由于回文模式的假阳性率非常低，因此可以安全地减少，甚至减少到1，以允许删除较短的接头片段；
f)keepBothReads：在回文模式检测到read测穿并删除接头序列后，反向读取包含与正向读取相同的序列信息。因此，默认行为是完全删除反向读取。通过为该参数指定true，反向读取也将被保留，这可能是有用的，例如，如果下游工具无法处理成对和非成对reads的组合。

2)SLIDINGWINDOW:<windowSize>:<requiredQuality>

当窗口内的平均质量低于阈值时，执行滑动窗口剔除。通过考虑多个碱基，单个质量差的碱基不会导致删除了在read后期高质量的数据。
a)windowSize: 设定窗口涵盖碱基数量；
b)requiredQuality: 设定需要的平均质量。

3)MAXINFO:<targetLength>:<strictness>

该步骤进行自适应质量微调，它平衡读取长度和错误率，根据需求获取长读取或者高质量，最大化节约成本。
a)targetLength: 指定目标位置的读取长度；
b)strictness: 这个值应该设置在0和1之间，它制定了保持尽可能多的读取长度与删除不正确的碱基之间的平衡。此参数值设置偏低<0.2有利于较长的读取，而偏高>0.8有利于读取正确性。

4)LEADING:<quality>

从起始端开始去除低质量的碱基。只要一个碱基的质量值低于阈值，就会切除该碱基，并调查下一个碱基。
quality: 指定保留碱基所需的最低质量。

5)TRAILING:<quality>

从末端移除低质量的碱基。只要碱基的质量值低于阈值，则切除该碱基，并调查下一个碱基（因为Trimmomatic从3'prime end开始，将是位于刚切除碱基之前的碱基）。此方法可用于去除Illumina低质量段区域（质量分数标记为2），但官方建议使用SLIDINGWINDOW或MAXINFO代替。
quality: 指定保留碱基所需的最低质量。

6)CROP:<length>

不管质量如何，从read末端删除碱基，指定read最大长度。
length: 从read起始端开始要保留的长度。

7)HEADCROP:<length>

不管质量如何，从read起始端端删除碱基。
length: 从read起始端开始要切除的长度。

8)MINLEN:<length>

此模块删除低于指定最小长度的读取。通常应该在所有其他处理步骤之后。此步骤删除的read将被计数，并包含在Trimmomatic摘要中显示的dropped reads计数中。
length: 设置保留reads的最小长度。

9)TOPHRED33/TOPHRED64
质量分数的编码转换。该指令没有其他参数。

（2）示例

1)完整运行示例

java -jar trimmomatic-0.32.jar PE \
    [-threads <threads>] \
    [-phred33|-phred64] \
    [-trimlog <trimLogFile>] \
    [<inputFile1> <inputFile2>] \
    [<outputFile1P> <outputFile1U> <outputFile2P> <outputFile2U>] \
    [ILLUMINACLIP:<fastaWithAdaptersEtc>:<seed mismatches>:<palindrome clip threshold>:<simple clip threshold>:<minAdapterLength>:<keepBothReads>] \
    [SLIDINGWINDOW:<windowSize>:<requiredQuality>] \
    [LEADING:<quality>] \
    [TRAILING:<quality>] \
    [MINLEN:<length>]

2)具体参数设置示例

java -jar trimmomatic-0.32.jar PE \
    -threads 16 \
    -phred33 \
    -trimlog trim.log \
    input_forward_R1.fq.gz input_reverse_R2.fq.gz \
    output_forward_paired_R1.fq.gz output_forward_unpaired_R1.fq.gz output_reverse_paired_R2.fq.gz output_reverse_unpaired_R2.fq.gz \
    ILLUMINACLIP:TruSeq3-PE.fa:2:30:10:8:true \
    SLIDINGWINDOW:5:20 \
    LEADING:3 \
    TRAILING:3 \
    MINLEN:36