生信工具Bioinformatics入门实战精通大百科

Trimmomatic:用于Illumina NGS数据的灵活读

2019-07-07  本文已影响11人  Amy_Cui

数据质控软件:
官网:http://www.usadellab.org/cms/index.php?page=trimmomatic
手册:http://www.usadellab.org/cms/uploads/supplementary/Trimmomatic/TrimmomaticManual_V0.32.pdf
个人觉得,直接看上面pdf,然后尝试提交一个代码,理解代码的参数含义,比较质控前后fastqc的结果就可以完成学习了。但是首先你要明白质控在做哪些事情?什么是接头,什么是N。。。

安装

conda安装

$ conda create -n rna3 python=3
$ source activate rna3
$ conda install trimmomatic
$ trimmomatic --help
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
$ trimmomatic -version
0.39

二进制源码安装

$ wget -c http://www.usadellab.org/cms/uploads/supplementary/Trimmomatic/Trimmomatic-0.39.zip
$ unzip Trimmomatic-0.39.zip
$ cd Trimmomatic-0.39/
$ java -jar trimmomatic-0.39.jar
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
$ alias Trimmomatic="java -jar /home/qmcui/Trimmomatic-0.39/trimmomatic-0.39.jar" 
# 可以写入.bashrc里,永久有效
$ Trimmomatic --help

双端任务1

$ Trimmomatic PE  \
/teach/project/1.rna/2.raw_fq/SRR1039510_1.fastq.gz  \
/teach/project/1.rna/2.raw_fq/SRR1039510_2.fastq.gz  \
./output_forward_paired.fq.gz ./output_forward_unpaired.fq.gz  \
./output_reverse_paired.fq.gz ./output_reverse_unpaired.fq.gz  \
ILLUMINACLIP:adapters/TruSeq3-PE.fa:2:30:10:2:keepBothReads \
LEADING:3 \
TRAILING:3 \
MINLEN:36 \
$ ps -ef|grep qmcui

结果

TrimmomaticPE: Started with arguments:
 /teach/project/1.rna/2.raw_fq/SRR1039510_1.fastq.gz /teach/project/1.rna/2.raw_fq/SRR1039510_2.fastq.gz output_forward_paired.fq.gz output_forward_unpaired.fq.gz output_reverse_paired.fq.gz output_reverse_unpaired.fq.gz ILLUMINACLIP:adapters/TruSeq3-PE.fa:2:30:10:2:keepBothReads LEADING:3 TRAILING:3 MINLEN:36
Multiple cores found: Using 4 threads
Using PrefixPair: 'TACACTCTTTCCCTACACGACGCTCTTCCGATCT' and 'GTGACTGGAGTTCAGACGTGTGCTCTTCCGATCT'
ILLUMINACLIP: Using 1 prefix pairs, 0 forward/reverse sequences, 0 forward only sequences, 0 reverse only sequences
Quality encoding detected as phred33
Input Read Pairs: 22852619 Both Surviving: 22393983 (97.99%) Forward Only Surviving: 285634 (1.25%) Reverse Only Surviving: 138619 (0.61%) Dropped: 34383 (0.15%)
TrimmomaticPE: Completed successfully
# 时间8min,输入read1和read2压缩文件为1.3G,均为91410476条reads

双端任务2

# Trimmomatic等于java -jar trimmomatic-0.35.jar 
$ Trimmomatic PE -phred33 \
input_forward.fq.gz input_reverse.fq.gz \
output_forward_paired.fq.gz output_forward_unpaired.fq.gz \
output_reverse_paired.fq.gz output_reverse_unpaired.fq.gz  \
ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 \
LEADING:3 \
TRAILING:3 \
SLIDINGWINDOW:4:15 \
MINLEN:36

单端任务1

java -jar trimmomatic-0.35.jar SE -phred33 input.fq.gz output.fq.gz ILLUMINACLIP:TruSeq3-SE:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36

参数解释

ILLUMINACLIP 去除adapter,后面有6个参数;如adapters/TruSeq3-PE.fa:2:30:10:2:keepBothReads (index.fa、错配个数2、比对分值30、最低比对分值10、切除接头最短序列2(default8)
、keepBothReads设置为true还是false(双端序列最好加上,再比如bowtie下游处理软件必需))
LEADING:3 去除领先的低质量或N碱基(低于质量3)
TRAILING:3 去除拖尾的低质量或N个碱基(低于质量3)(TRAILING:3)
SLIDINGWINDOW:4:15 使用4个底宽的滑动窗口扫描读数,当每个基础的平均质量降至15以下时切割
MINLEN:36 丢掉长度小于36nt长的reads

全部参数

不理解滑窗的可以适当google了解一点

# --help出现的参数
[-threads <threads>] 控制线程数,如果不提供,自动选择
[-phred33 | -phred64] illumina的锅,目前数据测得都是-phred33
# 你可以去网上搜索一段代码,怎么判断自己的数据是phred33还是phred64
# 运行任务的时候,提示信息里会出现程序判断的数据类型,显示信息如下
# Quality encoding detected as phred33
[-trimlog <trimLogFile>] 对每条reads出一条log,一般不用这个参数,日志条数太多了
[-basein <inputBase> | <inputFile1> <inputFile2>] 输入文件
[-baseout <outputBase> | <outputFile1P> <outputFile1U> <outputFile2P> <outputFile2U>] 输出文件,注意顺序

## 其他调节参数
ILLUMINACLIP:从读取中剪切适配器和其他Illumina特定序列。
SLIDINGWINDOW:执行滑动窗口修剪,一旦窗口内的平均质量低于阈值,就切割。
LEADING:如果低于阈值质量,则在读取开始时切断碱基
TRAILING:如果低于阈值质量,则在读取结束时剪切碱基
CROP:将读取切割为指定的长度
HEADCROP:从读取开始剪切指定的碱基数
MINLEN:如果读取低于指定长度,则删除读取
TOPHRED33:将质量得分转换为Phred-33
TOPHRED64:将质量得分转换为Phred-64

-trimlog

指定trimlog文件会创建所有读取修剪的日志,指示以下详细信息:

ILLUMINACLIP:

为TruSeq2(用于GAII机器)和TruSeq3(用于HiSeq和MiSeq机器)提供,用于单端和双端模式。这个有仪器决定,一般你用的PCR的接头的试剂盒。

ILLUMINACLIP:<fastaWithAdaptersEtc>:<种子不匹配>:<palindrome clip threshold>:<简单剪辑阈值>

SLIDINGWINDOW:<windowSize>:<requiredQuality>

滑动窗口:<windowSize>:<requiredQuality>
windowSize:指定平均值的基数
requiredQuality:指定所需的平均质量。

LEADING:<quality>

LEADING:<质量>
quality:指定保持基础所需的最低质量。

TRAILING:<quality>

TRAILING:<质量>
quality:指定保持基础所需的最低质量。

CROP:<length>

CROP:<长度>
length:从读取开始时要保留的碱基数。

HEADCROP:<length>

HEADCROP:<长度>
length:从读取开始时删除的碱基数。

MINLEN:<length>

MINLEN:<长度>
length:指定要保留的最小读取长度。

输入文件要求

For input and output files adding .gz/.bz2 to an extension tells Trimmomatic that the file is
provided in gzip/bzip2 format or that Trimmomatic should gzip/bzip2 the file, respectively.

双端输出文件要求

对照着--help的4个文件理解输出文件
For paired-end data, two input files, and 4 output files are specified, 2 for the 'paired' output
where both reads survived the processing, and 2 for corresponding 'unpaired' output where a
read survived, but the partner read did not.

image.png

质量值类型:自动选择

以前版本默认64,现在版本自动根据数据来做判断
-phred33 or -phred64 specifies the base quality encoding. If no quality encoding is specified,it will be determined automatically (since version 0.32). The prior default was -phred64.

线程数:自动选择

-threads indicates the number of threads to use, which improves performance on multi-core
computers. If not specified, it will be chosen automatically.


结束

修剪按照在命令行上指定步骤的顺序进行。在大多数情况下,建议尽可能早地完成适配器剪辑。

上一篇下一篇

猜你喜欢

热点阅读