Trimmomatic:用于Illumina NGS数据的灵活读取修剪工具

数据质控软件:
官网: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文件会创建所有读取修剪的日志,指示以下详细信息:

  • the read name
  • the surviving sequence length
  • the location of the first surviving base, aka. the amount trimmed from the start the location of the last surviving base in the original read
  • the amount trimmed from the end

ILLUMINACLIP:

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

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

  • fastaWithAdaptersEtc:指定包含所有adapters,PCR序列等的fasta文件的路径。(the path to a fasta file containing all the adapters, PCR sequences etc.)
  • seedMismatches:指定仍允许执行完全匹配的最大不匹配数目,如2
  • palindromeClipThreshold:指定两个“适配器连接”读取之间的匹配对于PE回文读取对齐的准确程度。
  • simpleClipThreshold:指定任何适配器等序列之间的匹配必须与读取的准确程度。

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.

  • mySampleFiltered_1P.fq.gz - for paired forward reads
  • mySampleFiltered_1U.fq.gz - for unpaired forward reads
  • mySampleFiltered_2P.fq.gz - for paired reverse reads
  • mySampleFiltered_2U.fq.gz - for unpaired reverse reads
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.


结束

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

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 218,284评论 6 506
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 93,115评论 3 395
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 164,614评论 0 354
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 58,671评论 1 293
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 67,699评论 6 392
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 51,562评论 1 305
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 40,309评论 3 418
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 39,223评论 0 276
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 45,668评论 1 314
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 37,859评论 3 336
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,981评论 1 348
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 35,705评论 5 347
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 41,310评论 3 330
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,904评论 0 22
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 33,023评论 1 270
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 48,146评论 3 370
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 44,933评论 2 355

推荐阅读更多精彩内容

  • Introduction What is Bowtie 2? Bowtie 2 is an ultrafast a...
    wzz阅读 5,660评论 0 5
  • 一、测序数据的格式转换 sra文件下载好后,使用fastq-dump转换数据格式: fastq-dump --sp...
    BioLearner阅读 4,206评论 0 1
  • 昨天到脚本虽然顺利运行,但是路径太麻烦,没有检测输出文件是否qc成功。 下面是师兄脚本里的语句。 PE模式,HiS...
    KK_f2d5阅读 8,055评论 0 6
  • ​很多人都知道人脉重要,都知道人脉即是钱脉。 但并不明白什么是真正的人脉!更不懂得如何获得人脉。 所谓的人脉,简单...
    蒋轩霖阅读 380评论 0 0
  • 莫说人寰好无情, 恐别离怨难无欲! 暗殇落痕摧无眠, ...
    顾旸阅读 175评论 0 1