Trimmomaticの使い方

Trimmomaticの概要
Trimmomaticはマルチスレッド対応のトリミングツールです。

FASTQデータを入力として、アダプターや末端配列の除去に加えて、低品質リード[phredスコア]の除去が行えます。

入力ファイルとして、fastqファイルまたはgzやbz2で圧縮されたファイルがサポートされています。
ペアエンドとシングルエンドの両方で使用できます。

また、アダプター配列としてイルミナ社の汎用アダプター配列がテンプレートして利用できる他、自分で配列を記載したファイルを用意すれば好きな配列を設定することも可能です。テンプレート配列を使用したい場合は、githubのページからダウンロードしてください。

設定が必要な項目が色々あるので、必要な項目を見ていきましょう。

Twitterで記事の更新をお知らせしているので、興味を持たれた方は是非フォローをお願いします。フォローする @harrykun_blog

1. ツールのインストール

$ conda install -c bioconda trimmomatic

インストールが正しく導入されたか確認する。

$ trimmomatic

正しく導入されていれば、下記が表示される。

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. Trimmomaticのコマンド例。

仮にサンプルAについてのペアエンドのWGSデータを持っているとして、コマンドの例を示します。

A_R1.fastq.gz A_R2.fastq.gz

trimmomatic PE \  
 -threads 16  \ 
 A_R1.fastq.gz A_R2.fastq.gz \  
 A_cleaned_1.fastq.gz  A_unpaired_1.fastq.gz    \ 
 A_cleaned_2.fastq.gz  A_unpaired_2.fastq.gz  \ 
 LEADING:30 \  
 TRAILING:30 \  
 SLIDINGWINDOW:4:15 \  
 MINLEN:20 \

オプションの項目について;

PE: PEはペアエンドモードを示す。シングルエンドの場合はSEとする。

-threads 16: パソコンのコア数に合わせて-threadsの後ろの数値を変える。

入力ファイル名と出力ファイル名を続けて入力する。
入力Fastqファイル名 2つ
出力Fastqファイル名４つ (cleanedはフィルターを通過したリード、unpairedは通過しなかったリード）

LEADING:30: しきい値以下の品質の場合、リードの先頭の塩基を削除する。

TRAILING:30: しきい値以下の品質の場合、リードの最後の塩基をカットする。

SLIDINGWINDOW:4:15: 5 "端からスキャンを開始して、ウィンドウ内の平

MINLEN:20: しきい値以下のリード長の場合、リードを削除する。

その他の項目について
ILLUMINACLIP：アダプターやその他指定した配列をリードからカットする。
CROP: 末端の塩基を除去して指定した長さにカッする。
HEADCROP: 読み取りの先頭から指定した数の塩基をカッする。
MINLEN: 指定した長さ以下の場合、読み取りを削除する。
AVGQUAL: 品質の平均値が指定されたレベルを下回った場合にリードをドロップする。
TOPHRED33: 品質スコアをPhred-33に変換する。
TOPHRED64：品質スコアを Phred-64 に変換する。

(なお、本ブログの過去記事で作成したマラリア原虫のFastq.gzファイルを例として用いる。)

3. Trimmomaticを実行する。

以前の記事で、トリミング前のFastqファイルに対して、FastQCによる品質チェックを行った。

FastQCについては別記事を参照してください。

eupatho-bioinfomatics.hatenablog.com

f:id:Harry-kun:20210625111425p:plain — シークエンススコアのヒストグラム

見てわかるように、全区間において平均スコアが30以上と良好なシークエンスを示している。今回は、トリミング前後の変化が分かりやすいように、過剰なトリミングLEADING:35・TRAILING:35を行う。

自身の解析で使う場合には、大体スコアを３０で設定しています。ただし、古いデータを使う際にはクオリティが低い場合があります。そういった場合に厳しい条件を設定しているとほとんどのリードが無くなってしまうので、フィルター条件の緩和を考慮することになります。

    java -jar trimmomatic-0.36.jar PE -threads 16 \
    ERR1081263_1.fastq.gz ERR1081263_2.fastq.gz \
    ERR1081263_cleaned_1.fastq.gz ERR1081263_unpaired_1.fastq.gz \
    ERR1081263_cleaned_2.fastq.gz ERR1081263_unpaired_2.fastq.gz \
    LEADING:35 TRAILING:35 SLIDINGWINDOW:4:15 MINLEN:20

トリミングが正しく行われたことを確認するために、処理後のFastqファイルに対して、再びFastQCを実行する。
cleaned_1.fastq.gzとcleaned_2.fastq.gzがトリミングを条件をPASSしたリードの入ったファイル。
出力された結果に対してmultiQCを使い、１つにまとめる。

MultiQCについてはこちらを参照してください。

eupatho-bioinfomatics.hatenablog.com