寫在前面

Trimmomatic工具是用于illumina二代測(cè)序數(shù)據(jù)的reads處理，主要對(duì)接頭(adapter)序列和低質(zhì)量序列進(jìn)行過濾。下面是使用該工具處理雙端測(cè)序(PE)數(shù)據(jù)時(shí)，常用參數(shù)的一些說明。

參考文檔

• Trimmomatic工具的參考文獻(xiàn)
• Trimmomatic工具官網(wǎng)
• Trimmomatic工具使用手冊(cè)

軟件使用

• 執(zhí)行命令

  ## 雙端測(cè)序數(shù)據(jù)使用方法 # 使用v0.32版本: 1. java -jar trimmomatic-0.32.jar PE \ 2. [-threads <threads>] \ 3. [-phred33|-phred64] \ 4. [-trimlog <trimLogFile>] \ 5. [<inputFile1> <inputFile2>] \ 6. [<outputFile1P> <outputFile1U> <outputFile2P> <outputFile2U>] \ 7. [ILLUMINACLIP:<fastaWithAdaptersEtc>:<seed mismatches>:<palindrome clip threshold>:<simple clip threshold>:<minAdapterLength>:<keepBothReads>] \ 8. [SLIDINGWINDOW:<windowSize>:<requiredQuality>] \ 9. [LEADING:<quality>] \ 10. [TRAILING:<quality>] \ 11. [MINLEN:<length>]# 舉例: 1. java -jar trimmomatic-0.32.jar PE \ 2. -threads 16 \ 3. -phred33 \ 4. -trimlog trim.log \ 5. input_forward_R1.fq.gz input_reverse_R2.fq.gz \ 6. output_forward_paired_R1.fq.gz output_forward_unpaired_R1.fq.gz output_reverse_paired_R2.fq.gz output_reverse_unpaired_R2.fq.gz \ 7. ILLUMINACLIP:TruSeq3-PE.fa:2:30:10:8:true \ 8. SLIDINGWINDOW:5:20 \ 9. LEADING:3 \ 10. TRAILING:3 \ 11. MINLEN:36

• 參數(shù)說明

  PE 設(shè)置使用trimmomatic處理雙端數(shù)據(jù)，單端數(shù)據(jù)用(‘SE’)
  -thread 16 設(shè)置線程數(shù)為16
  -phred33 設(shè)置堿基的質(zhì)量格式(默認(rèn)-phred64，自v0.32版本之后可自動(dòng)識(shí)別是phred33還是phred64)
  -trimlog trim.log 設(shè)置trimmommatic工具處理的日志文件為’trim.log’，每兩行為一對(duì)reads信息

  • 生成的日志文件’trim.log’包含5列信息：

    1) read名稱2) 切除后剩余的read長度3) 切除后第一個(gè)堿基所在位置(0-base)，也就是起始位置切除了幾個(gè)堿基4) 切除后最后一個(gè)堿基所在位置5) read末尾切除了幾個(gè)堿基# 由于該生成文件較大，如后續(xù)步驟不需該文件信息，可考慮不設(shè)置

  input_forward_R1.fq.gz 輸入的forward正向鏈R1對(duì)應(yīng)的fastq文件
  input_reverse_R2.fq.gz 輸入的reverse反向鏈R2對(duì)應(yīng)的fastq文件
  output_forward_paired_R1.fq.gz 處理后輸出的R1對(duì)應(yīng)reads成對(duì)fastq文件
  output_forward_unpaired_R1.fq.gz 處理后輸出的R1對(duì)應(yīng)reads不成對(duì)的fastq文件
  output_reverse_paired_R2.fq.gz 處理后輸出的R2對(duì)應(yīng)reads成對(duì)fastq文件
  output_reverse_unpaired_R2.fq.gz 處理后輸出的R2對(duì)應(yīng)reads不成對(duì)的fastq文件
  ILLUMINACLIP:TruSeq3-PE.fa:2:30:10:8:true 切除illumina接頭參數(shù)設(shè)置。說明ILLUMINACLIP的各參數(shù)之前，先解釋一下要使用的兩種切除接頭的模式和比對(duì)分值計(jì)算方法

• Simple mode, simple模式，它可用于切除任何序列(any technical sequence, 暫且稱之為污染序列)。該模式的方法是將污染序列直接與一條read進(jìn)行比對(duì)（局部比對(duì))，將進(jìn)行下面4個(gè)步驟(情形)：
  

  Simple mode是根據(jù)與污染序列的比對(duì)情況判斷切除序列A) 污染序列3'端的一部分與read5'端有overlap，去除整條read；B) 污染序列從read的5'端開始能與read完全overlap, 去除整條read；C) 污染序列在read中間有overlap，5'端一部分沒有overlap的保留，overlap的部分及之后的進(jìn)行切除；D) 污染序列3'端的一部分與read有overlap，切除overlap部分.

• Palindrome mode，palindrome模式，它僅用于read測(cè)穿(read-through)的情形， 該模式的方法是：在識(shí)別接頭序列之前，先將接頭序列加在成對(duì)reads起始（‘a(chǎn)dapter ligated’ read，read-with-adapter, 暫且稱加adapter后的這對(duì)reads為a-R1和a-R2），再將a-R1和a-R2進(jìn)行比對(duì)（全局比對(duì))。 該模式將進(jìn)行如下4個(gè)步驟(情形)：
  

  Palindrome mode: 是根據(jù)a-R1, a-R2的overlap判斷切除的序列A) overlap是在接頭序列與reads的起始之間，即overlap只有接頭序列，沒有有用的信息，這種情形將這兩條read都整條去除；B) overlap有一部分是DNA插入片段，一部分是整個(gè)接頭序列，該情形會(huì)將這對(duì)reads從overlap序列及剩余序列(3'端)切除；C) overlap中有接頭序列的一小部分，與B)相同方式切除；D) overlap中不存在接頭序列.

  再說明一下比對(duì)分值的計(jì)算方法：

  match：匹配一個(gè)堿基加分，分值為$log_10{4}$(約等于0.602)；mismatch: 錯(cuò)配一個(gè)堿基減分，分值為Q/10，其中，Q為對(duì)應(yīng)堿基的堿基質(zhì)量值，由于一般堿基質(zhì)量值區(qū)間為: [0-40]，即錯(cuò)配的罰分值區(qū)間[0-4]。該分值公式說明堿基質(zhì)量越高的堿基出現(xiàn)錯(cuò)配時(shí)罰分越多.## 例如：# 12bp的接頭序列能夠完全比對(duì)到read上的分值為：12 * 0.602 ≈ 7；# 25bp的接頭序列完全比對(duì)到read上的分值為：25 * 0.602 ≈ 15；# a-R1和a-R2有50bp完全匹配的分值為50 * 0.602 ≈ 30.

  下面是ILLUMINACLIP的各參數(shù)說明：

  ":<fastaWithAdaptersEtc>:<seed mismatches>:<palindrome clip threshold>:<simple clip threshold>:<minAdapterLength>:<keepBothReads>"- ":TruSeq3-PE.fa" # <fastaWithAdaptersEtc>是fasta格式的接頭序列文件路徑- ":2" # <seed mismatches>是將接頭序列的一段(不超過16bp)作為seed，與reads進(jìn)行比對(duì)，能夠容忍的最大錯(cuò)配(mismatch)數(shù)，這里是最多2個(gè)錯(cuò)配- ":30" # <palindrome clip threshold>是 a-R1, a-R2的比對(duì)分值閾值，達(dá)到閾值，才進(jìn)行切除，這里設(shè)置閾值為30（大約比對(duì)50bp堿基）- ":10" # <simple clip threshold>是任意(接頭)序列與read比對(duì)最低分閾值，大于這個(gè)閾值，才進(jìn)行切除，這里設(shè)置閾值為10（大約比對(duì)17bp堿基）- ":8" # <minAdapterLength>只作用于Palindrome模式，是設(shè)置檢測(cè)到接頭序列的最小長度（默認(rèn)為8，甚至可設(shè)置為1）- ":true" # <keepBothReads>只作用于palindrome模式，是設(shè)置是否保留反向鏈，這里是說去除接頭序列后，由于正反鏈包含相同的序列信息(盡管序列是反向互補(bǔ)的)，默認(rèn)情況(":false")下會(huì)去除反向鏈，設(shè)置為":True"則保留反向鏈

SLIDINGWINDOW:5:20 設(shè)置滑動(dòng)窗口閾值，以為5bp為窗口，這5bp堿基的平均質(zhì)量值低于20，要進(jìn)行切除
LEADING:3 設(shè)置從reads起始開始，去除質(zhì)量低于閾值或?yàn)?#39;N'的堿基，直到達(dá)到閾值不再去除，這里設(shè)置閾值為3
TRAILING:3 設(shè)置從reads末尾開始，去除質(zhì)量低于閾值的堿基或?yàn)?#39;N'的堿基直到達(dá)到閾值不再去除，這里設(shè)置閾值為3，這種方法是去除特定的illumina平臺(tái)低質(zhì)量區(qū)域(由于illumina會(huì)將低質(zhì)量堿基標(biāo)記為2)，官方操作文檔中建議使用 SLIDINGWINDOW 或 MAXINFO代替[這里未給出MAXINFO參數(shù)說明] ）
MINLEN:36 設(shè)置read切除后的最短長度，這里設(shè)置長度至少為36bp，長度小于36bp的reads被去除

總結(jié)

以上是生活随笔為你收集整理的学习使用一款数据质控软件（Trimmomatic）的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。