de novo アセンブリー

Trinity

Trinity は de novo アセンブリーを行うプログラムの一つである。シーケンサーから出力されるリードからトランスクリプトームを作成するのがメイン機能である。

Trinity によるアセンブリ

single-end read

single-end リードの場合は Trinity コマンドを実行するときにオプションとして --single を与える。

Trinity --seqType fq \
        --single SRR504111.fastq \
        --CPU 8 \
        --max_memory 8G \
        --output SRR504111_ref

リードが starnd-specific RNA-seq により出力されているならば、その方向を --SS_lib_type オプションで指定する。single-end リードならば F または R を指定する。また、アセンブルされるコンティグの長さの最小値は --min_contig_length で指定する。

Trinity --seqType fq \
        --single SRR504111.fastq \
        --SS_lib_type F \
        --min_contig_length 300 \
        --CPU 8 \
        --max_memory 8G \
        --output SRR504111_ref

paired-end read

paired-end リードの場合は Trinity コマンドを実行する際にオプション --left--right にそれぞれ *_1.fastq と *_2.fastq ファイルを指定する。

Trinity --seqType fq \
        --left SRR504904_1.fastq \
        --right SRR504904_2.fastq \
        --CPU 8 \
        --max_memory 8G \
        --output SRR504904_ref

FASTQ ファイルが複数セットある場合、カンマ区切りで与える。

Trinity --seqType fq \
        --left SRR504904_1.fastq,SRR504905_1.fastq,SRR504906_1.fastq \
        --right SRR504904_2.fastq,SRR504905_2.fastq,SRR504906_2.fastq \
        --CPU 8 \
        --max_memory 8G \
        --output SRR504904_ref

この他に特に注意すべきオプションとしては、--SS_lib_type--min_contig_length などがある。リードが starnd-specific RNA-seq により出力されているならば、その方向を --SS_lib_type オプションで指定する。paired-end リードの場合は RF または FR を指定する。また、アセンブルされるコンティグの長さの最小値は --min_contig_length で指定する。

Trinity --seqType fq \
        --left SRR504904_1.fastq \
        --right SRR504904_2.fastq \
        --SS_lib_type FR \
        --min_contig_length 300 \
        --CPU 8 \
        --max_memory 8G \
        --output SRR504904_ref

転写物アノテーション

Trinity を実行し終えると、ディレクトリ下に、Trinity.fasta ファイルが得られる。この FASTA ファイルに、アセンブリ後の転写物の配列が保存されている。アセンブリされた転写物に、アノテーションをつけるには、Trinotate を利用する。

まず、次のコマンドを利用して、必要なデータをインターネットからダウンロードし、ファイルを準備する。このコマンドが正しく実行されると、ディレクトリには、Trinotate.sqlite、uniprot_sprot.pep、Pfam-A.hmm.gz の 3 つのファイルが作成される。

Trinotate-3.0.2/admin/Build_Trinotate_Boilerplate_SQLite_db.pl Trinotate

次に、アセンブリされた塩基配列を、アミノ酸配列に翻訳する。これを実行すると、Trinity.fasta.transdecoder.pep ファイルが作成される。

TransDecoder-3.0.1/TransDecoder.LongOrfs -t Trinity.fasta
TransDecoder-3.0.1/TransDecoder.Predict -t Trinity.fasta

準備が整いたら、さきほどダウンロードした UniProt のデータに対して BLAST を実行し、ホモログを探す。BLAST は数日かかる場合がある。

makeblastdb -in uniprot_sprot.pep -dbtype prot

blastx -query Trinity.fasta -db uniprot_sprot.pep -num_threads 8 -max_target_seqs 1 -outfmt 6 > blastx.outfmt6
blastp -query Trinity.fasta.transdecoder.pep -db uniprot_sprot.pep -num_threads 8 -max_target_seqs 1 -outfmt 6 > blastp.outfmt6

BLAST の実行が終了すると、BLAST の実行結果を収集して、結果をエクセルファイルに保存する。

trinityrnaseq-Trinity-v2.4.0/util/support_scripts/get_Trinity_gene_to_trans_map.pl Trinity.fasta > Trinity.fasta.gene_trans_map

Trinotate-3.0.2/Trinotate Trinotate.sqlite init --gene_trans_map Trinity.fasta.gene_trans_map \
                                                --transcript_fasta Trinity.fasta \
                                                --transdecoder_pep Trinity.fasta.transdecoder.pep

Trinotate-3.0.2/Trinotate Trinotate.sqlite LOAD_swissprot_blastp blastp.outfmt6
Trinotate-3.0.2/Trinotate Trinotate.sqlite LOAD_swissprot_blastx blastx.outfmt6

Trinotate-3.0.2/Trinotate Trinotate.sqlite report > trinotate_annotation_report.xls

エラーが発生した場合

同じデータに対して、Trinity を連続して実行した場合、2 回目で以下のエラーになる。これは、1 回目の実行ですでに出力ディレクトリが作成されているためである。対処法としては、出力ディレクトリ --output を変更するか、既存のディレクトリを削除する。

terminate called after throwing an instance of 'jellyfish::file_parser::FileParserError'
  what():  'single.fa': Invalid input file. Expected fasta or fastq

References

  • Grabherr MG, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol. 2011, 29(7):644-52. PubMed Abstract
  • Haas BJ, et al. De novo transcript sequence reconstruction from RNA-seq using the Trinity platform for reference generation and analysis. Nat Protoc. 2013, 8(8):1494-512. PubMed Abstract