Trinity - JI

アッセンブル：Trinity

2020 年 10 月 11日　改訂

Trinity はトランスクリプトーム解析専用の k-mer アセンブラです．ここでは， Trinity とその他のプログラムを用いて NGS で得られたトランスクリプトームデータを de novo assembly (参照配列を用いずに配列張り合わせ) する手順を紹介します．得られたリードが 150bp より長い場合は，Newbler の方が良いみたいです．

SGE (Sun Grid Engine) スケジューラーで解析を行う場合は，こちらのジョブスクリプトを参照してください．

種の系統解析を行う場合：詳しい同僚は、Trimmomatic, Trinity, CD-HIT, TransDecoder による処理が必要、と言っていました。
　CD-HIT による解析で得られる配列を絞り込みたい場合は、塩基配列だけでなく、さらに、これを翻訳したアミノ酸配列でも解析した方が良いそうです。塩基およびアミノ酸配列の解析、それぞれ 95% 以上の相同性を示す配列は捨てて良いそうです。5 万配列で 10 分程度で終了します (2020 年 10 月)。

1. Trinity

ここでは Trinity を用いて，参照配列を用いない配列張り合わせ (de novo assembly) を行います．Trinity は fastq ファイルを読み込み，一つの fasta file (すべてのトランスクリプトーム配列を含む) を算出します．

種の系統解析が目的ならば、largest variant を抽出する処置を行うと良いそうです。K mer coverage でエラーが出ることがあります。その場合は、デフォルトでは 2 なので、得られる配列数は少なくなりますが、3, 5, 10 など変化させて解析して良いそうです (2020 年 10 月)。

解析を行う前に，Solexa QA 2 や condetri などによって (上記の手順)，quality が低いリードを除外してからアセンブルを実行した方が良いらしいです [<= これは古い技術のようで、2018 年の段階では、おすすめされませんでした]．

こちらは英語の資料ですが，Trinity の全体像，つまり，Inchworm[シャクトリ虫], Chrysalis[サナギ], Butterfly[チョウ] の役割，をつかむには良いです (とくに HOW IT WORKS というビデオ)．日本語の資料としてはこちら (1, 2) を参考にしました．

インストール
こちらをご覧ください．

テストデータの解析

ダウンロード & コンパイルして得られたファイル，

trinityrnaseq_r2013-02-25/sample_data/test_Trinity_Assembly

に入ってください．そこで，

sh runMe.sh

と入力すれば，解析が始まるはずです．trinity_out_dir にアウトファイルが保存されます．Trinity.fasta が最終的に必要なファスタファイルです．
実際にこのスクリプト内部で解析を行っている行は，

../../Trinity.pl --seqType fq --JM 2G --left reads.left.fq --right reads.right.fq --SS_lib_type RF --CPU 4 --no_cleanup --monitoring

です．もう一行，

../../util/RSEM_util/run_RSEM_align_n_estimate.pl --transcripts trinity_out_dir/Trinity.fasta --seqType fq --left reads.left.fq --right reads.right.fq --SS_lib_type RF

がありますが，これは発現量解析などで使うコマンドのようです．こちらをご覧下さい．私は使わないです．

．

Slurm の job file

テストファイルの解析には、以下の jobfile を使いました。

#!/bin/bash
#SBATCH --job-name=Trinity
#SBATCH --mail-user="jun.inoue@oist.jp"
#SBATCH --partition=compute
#SBATCH --mem=2G
#SBATCH --cpus-per-task=4
#SBATCH --ntasks=1 # 1 task

software_dir='/bucket/SatohU/inoue/softwares'
module load cmake/3.18.4
export PATH=$software_dir/salmon-latest_linux_x86_64/bin:$PATH
export PATH=$software_dir/bowtie2-2.4.2-sra-linux-x86_64:$PATH
export PATH=$software_dir/jellyfish-2.3.0/bin:$PATH
export PATH=$software_dir/trinityrnaseq-Trinity-v2.8.4:$PATH
export PATH=$software_dir/samtools-1.11/bin:$PATH

Trinity --seqType fq --max_memory 2G \
--left /.../reads.left.fq.gz \
--right /.../reads.right.fq.gz \
--SS_lib_type RF \
--CPU 4

実際のデータは以下の jobfile です。

#!/bin/bash
#SBATCH --job-name=TriTrimed
#SBATCH --mail-user="jun.inoue@oist.jp"
#SBATCH --partition=compute
#SBATCH --mem=100G
#SBATCH --cpus-per-task=20
#SBATCH --ntasks=1 # 1 task

software_dir='/work/SatohU/0inoue/software'
module load cmake/3.11.4
export PATH=$software_dir/salmon-0.11.3-linux_x86_64/bin:$PATH
export PATH=$software_dir/bowtie2-2.3.4.3-linux-x86_64:$PATH
export PATH=$software_dir/jellyfish-2.2.10/bin:$PATH
export PATH=$software_dir/trinityrnaseq-Trinity-v2.8.4:$PATH
export PATH=$software_dir/samtools-1.9/bin:$PATH

TRINITY_HOME='/work/SatohU/0inoue/software/trinityrnaseq-Trinity-v2.8.4'

${TRINITY_HOME}/Trinity --seqType fq --max_memory 100G \
--left paired_output_sub1.fq \
--right paired_output_sub2.fq\
--CPU 20

--max_memory
スケジューラーで割り当てたメモリー数を書きます．

．

結果

Trinity.pl がうまく走ると，以下のようなスクリーンアウトが得られるはずです

[jun-inoue:test_Trinity_Assembly]$ perl /home/j/jun-inoue/bin/trinityrnaseq_r2013-02-25/Trinity.pl --seqType fq --JM 2G --left reads.left.fq --right reads.right.fq --SS_lib_type RF --CPU 4 --no_cleanup --monitoring
Use of uninitialized value in pattern match (m//) at /home/j/jun-inoue/bin/trinityrnaseq_r2013-02-25/Trinity.pl line 474.
WARNING, --monitoring can only be used on linux. Turning it off.

Current settings:
core file size (blocks, -c) 0
data seg size (kbytes, -d) unlimited
scheduling priority (-e) 0
file size (blocks, -f) unlimited
pending signals (-i) 399360
max locked memory (kbytes, -l) unlimited
max memory size (kbytes, -m) unlimited
open files (-n) 1024
pipe size (512 bytes, -p) 8
POSIX message queues (bytes, -q) 819200
real-time priority (-r) 0
stack size (kbytes, -s) unlimited
cpu time (seconds, -t) unlimited
max user processes (-u) 399360
virtual memory (kbytes, -v) unlimited
file locks (-x) unlimited

Paired mode requires bowtie. Found bowtie at: /usr/bin/bowtie

-since butterfly will eventually be run, lets test for proper execution of java
#######################################
Running Java Tests
CMD: java -Xmx64m -jar /imports/home/j/jun-inoue/bin/trinityrnaseq_r2013-02-25/util/ExitTester.jar 0
CMD finished (0 seconds)
CMD: java -Xmx64m -jar /imports/home/j/jun-inoue/bin/trinityrnaseq_r2013-02-25/util/ExitTester.jar 1

-we properly captured the java failure status, as needed. Looking good.
Java tests succeeded.
###################################

..........

---------------------------------------------------------------
-------------------- Butterfly --------------------------------
-- (Reconstruct transcripts from reads and de Bruijn graphs) --
---------------------------------------------------------------

CMD: /imports/home/j/jun-inoue/bin/trinity........
Number of Commands: 52
succeeded(52) 100% completed.

All commands completed successfully. :-)

CMD finished (16 seconds)
CMD: /imports/home/j/jun-inoue/bin/trinityrnaseq_r2013-02-25/util/print_butterfly_assemblies.pl /work/........
CMD finished (0 seconds)

###################################################################
Butterfly assemblies are written to /work/SinclairU/inoue/tr
###################################################################

[jun-inoue:test_Trinity_Assembly]$

リードの方向

こちらのページに，まとめられています．
(以下，作成中です)
--SS_lib_type オプションで設定します．わからない時はこのオプションは使わない方が良いです．ライブラリの作り方がわかっていない場合におかしな設定を行うと，結果が無茶苦茶になるそうです．以下は mouse SRA のページを例としています．

青のラインがある部分で，NEBnext Ultra RNA Library Prep Kit を使った，と書いてあります．この製品試薬のページから判定します．

アウトファイル
trinity_out_dir というディレクトリが自動的に作成され，ここに結果が保存されます．このなかにある，

Trinity.fasta

というファイルに，アッセンブルの結果得られた contig が fasta 形式で保存されています．Trinity.fasta が出力されない場合は，何か解析に問題があったことになります．Butterfly はメモリをたくさん使うようで，私の場合は Butterfly が終了せず Trinity.fasta が得られないことがよくあります．
　他にたくさんのアウトファイルが得られますが，普通は使いません．

#####

TrinityStats.pl を使うと，得られた結果の内容を知ることができます．

[jun-inoue:TEST]$ perl $TRINITY_HOME/util/TrinityStats.pl trinity_out_dir/Trinity.fasta
Total trinity transcripts: 90
Total trinity components: 47
Contig N50: 4954

データ量

リード数が 30M～60M が常識的なデータ量だそうです．これぐらいなら，Trinity でそのまま解析して良いです．これより多い場合は，リード数を間引くなどの処理をしないと，結果の質が悪くなることがあるそうです．良いアッセンブルを行うのに 100M 必要，という人もいますが，30M で十分だそうです。

人によっては、250MB 必要、という場合もあります　(2020 年 10 月)。

エラー
ヘッダー：ヘッダーを変更するように，エラーメッセージが出る．

[cluster:Symsagittifera-roscoffensis]$ cat slurm-4622056.out
Trinity-v2.8.4
....
CMD: seqtk-trinity seq -A /work/SatohU/0inoue/Symsagittifera-roscoffensis/paired_output_sub2.fq >> right.fa
Error, not recognizing read name formatting: [SRR5760179.40838091]

If your data come from SRA, be sure to dump the fastq file like so:

SRA_TOOLKIT/fastq-dump --defline-seq '@$sn[_$rn]/$ri' --split-files file.sra

Thread 1 terminated abnormally: Error, cmd: seqtk-trinity seq -A /work/SatohU/0inoue/Symsagittifera-roscoffensis/paired_output_sub1.fq >> left.fa died with ret 512 at /work/SatohU/0inoue/software/trinityrnaseq-Trinity-v2.8.4/util/insilico_read_normalization.pl line 762.

この例題は SRA データを解析中に出たものです．bold の行は，ヘッダーに問題があるから，変更する必要がある，というメッセージが出ています． NCBI が配っている SRA tool kit に入っている fastq-dump を用いて SRA ファイルを二つの fastq ファイルに分割する際に，上記 bold のオプションをつけるように支持されています．fastq-dump の使い方はこちらをご覧ください．

bowtie：bowtie が見つからないというメッセージが出ました．

[i2:test_Trinity_Assembly]$ sh runMe.sh
....
Error, cannot find path to bowtie, which is now needed as part of Chrysalis' read scaffolding step at ../../Trinity.pl line 622.

こちらのサイトにあるダウンロードサイトから bowtie-1.0.0-macos-i386.zip をダウンロードしました (uname -p で CPU の種類を確認)．解凍して得られた bowtie を ~/bin にコピーしたら問題なく動きました．bowtie-build もない場合は，~/bin にコピーしましょう．

データの前処理は必要？
Solexa QA 2 と filterPCRdupl.pl によって前処理を行った場合，Trinity の結果がどの程度少なくなるか Trinity のテストデータ (reads.left.fq と reads.right.fq) を用いて比較してみました．解析には TrinityStats.pl を使いました．

前処理をやらない場合