macでインフォマティクス

macでインフォマティクス

HTS (NGS) 関連のインフォマティクス情報についてまとめています。

de novo transcriptomeのアセンブル評価ツール RSEM-EVAL(DETONATEパッケージ)

2018、8/9 誤字脱字修正

2021 12/24 タイトル変更

 

 RNAシーケンシング(RNA-Seq)技術は、トランスクリプトームの大規模分析を可能にすることによって、ゲノムの配列がまだ決定されていない種の研究に革命をもたらしている。このようなトランスクリプトームを研究するためには、ゲノム配列情報の助けを借りずに、RNA-Seqのリードからde novoトランスクリプトームアセンブリを介してtranscripts配列を再構成し、転写物配列のセットを決定しなければならない。 Roche 454 Life Scienceのプラットフォームデータ[論文より ref.9]〜[12]を対象としたイルミナのプラットフォームデータ[1]〜[8]やその他のものについて、多くのde novo transcriptome assemblersが現在利用可能である。これらのアセンブラは、しばしばかなりのユーザ調整可能パラメータのセットと組み合わされて、単一データセットの候補アセンブリの大きなスペースの生成を可能にしている。しかし、特にgrand truthがわからないときに、アセンブリの精度を適切に評価することは困難である。

 de novo transcriptome assemblyの評価には多くの研究が行われている[ref.13]〜[20]。このような研究で使用されるアセンブリ評価尺度は、リファレンスベースとリファレンスなしの2つのクラスに分類することができる。基準に基づく尺度は、既知の系列を用いて計算される尺度である。例えば、アセンブリとリファレンス転写産物セットとの間の対応を確立した後、リファレンス転写産物とアセンブリのマッチした割合(recall)、または転写産物とこれらの2つの組み合わせ(例えば、F1 measure)[5]、[16]、[17]。転写産物セットに加えて、ゲノムおよびタンパク質配列もまた、アセンブリ評価のためのリファレンスとして使用することができる[ref.2]、[ref.4]、[ref.8]、[ref.13]、[ref.15]、[ref.20]。

 しかし、デノボアセンブリに関心のあるほとんどのケースでは、リファレンスシーケンスは利用できないか、不完全であるか、または関心のあるサンプルのgrand truthからかなり離れているため、アセンブリ評価作業が著しく困難になる。そのような場合には、リファレンスフリーの手段に頼らなければならない。一般的に使用されるリファレンスフリー測定値には、コンティグ長の中央値、コンティグ数およびN50が含まれる[13]、[16]、[17]。残念なことに、これらの措置は原始的であり、しばしば誤解を招く[20]。たとえば、最も一般的なリファレンスフリー測定の1つであるN50は、簡単なアセンブリで最大化できる。 N50は、少なくともその長さのすべてのコンティグがアセンブリの塩基の少なくとも50%を構成するような最長コンティグの長さと定義される[21]。この対策の動機は、より良いアセンブリは、入力リード間の識別されたオーバーラップの数が多いことに起因するため、より長いコンティグに組み合わされたリードが多くなることである。しかし、すべての入力を1つのコンティグに連結することによって構築された簡単なアセンブリは、この尺度を最大化してしまうことは容易に理解できる。要するに、N50はコンティグの連続性を測定するが、精度は測定しない[22]。アセンブリがシングルトン(すなわち、シングルリードから得られるコンティグ)を含む場合に潜在的に有益であることが示されているが、他の単純化されたリファレンスフリー測定は同様に誤解を招く可能性がある[20]。

 著者らは、DETONATE方法論とソフトウェアパッケージを提示することにより、トランスクリプトームアセンブリ評価における最先端技術を向上させる。 DETONATEは、RSEM-EVALとREF-EVALの2つのコンポーネントで構成されている。 DETONATEの主な寄稿者であるRSEM-EVALは、アセンブリにのみ依存する新規確率モデルに基づくリファレンスフリーの評価方法である。 RSEM-EVALは、統計モデルを用いてゲノムおよびメタゲノム[24]、[25]アセンブリを評価または構築する最近のアプローチに類似しているが、 著者らが論じるように、転写産物のabundanceおよびオルタナティブスプライシング、N50のような単純化されたリファレンスフリー測定とは異なり、RSEM-EVALは、アセンブリのコンパクトさやRNA-Seqデータからのアセンブリのサポートなど、複数の要因を単一の統計的に原理的な評価スコアに組み合わせる。このスコアは、最適なアセンブラを選択し、アセンブラのパラメータを最適化し、新しいアセンブラの設計を目的関数として導くために使用できる。さらに、アセンブリ内の各コンティグについて、RSEM-EVALは、コンティグがRNA-Seqデータによってどの程度うまくサポートされているかを評価するスコアを提供し、不要なコンティグをフィルタリングするために使用することができる。 DETONATEの第2のコンポーネントであるREF-EVALは、リファレンスベースの対策のツールキットである。これは、既存のリファレンスベースの測定値よりも、アセンブリの精度をより洗練された視点で提供する。

 著者らは、RSEM-EVALスコアの価値を実証するために、実際のデータとシミュレートされたデータの両方について多数の実験を行った(以下略)。

 

DETONATE:

http://deweylab.biostat.wisc.edu/detonate/

 

インストール

公式サイトからダウンロードする。

DETONATE: DE novo TranscriptOme rNa-seq Assembly with or without the Truth Evaluation

解凍してビルド。

cd detonate-1.11/rsem-eval
make

#またはcondaを使う
#bioconda(link)
conda install -c bioconda -y detonate

 

ラン

1、初めに近縁なモデル生物のtranscripts.fastaの長さを出力する。

cd detonate-1.11/rsem-eval/
./rsem-eval-estimate-transcript-length-distribution transcripts.fa parameter_file

 parameter_fileが出力される。

 

2、アセンブルしたFASTAとそれに使ったfastqを指定して RSEM-EVALスコアを計算する。

Paired-end

./rsem-eval-calculate-score -p 8 --transcript-length-parameters parameter_file --paired-end pair1.fq pair2.fq assembly.fasta L

 L For single-end data, L represents the average read length. For paired-end data, L represents the average fragment length. It should be a positive integer (real value will be rounded to the nearest integer).

 

Single-end

./rsem-eval-calculate-score -p 8 --transcript-length-parameters parameter_file single.fq assembly.fasta L

 

 解析が終わるといくつかファイルができるが、~socreが評価ファイルである。下はテストデータを使った時の~socreファイル出力。

$ cat detonate-1.11/rsem-eval/toy_assembly_1.score 

Score -87426.14

BIC_penalty -8.25

Prior_score_on_contig_lengths_(f_function_canceled) -7.91

Prior_score_on_contig_sequences -867.82

Data_likelihood_in_log_space_without_correction -86542.17

Correction_term_(f_function_canceled) -0.00

Number_of_contigs 1

Expected_number_of_aligned_reads_given_the_data 3812.00

Number_of_contigs_smaller_than_expected_read/fragment_length 0

Number_of_contigs_with_no_read_aligned_to 0

Maximum_data_likelihood_in_log_space -86541.97

Number_of_alignable_reads 3812

Number_of_alignments_in_total 3812

1行目はスコア。負の値であるが、大きい方が優れている(-20万より-8万がより良い)

スコアが続き、7行目の Number_of_contigはアセンブリに含まれるコンティグの数。

8行目のNumber_of_contigs_smaller_than_expected_read/fragment_lengthは、長さが期待されるフラグメント長よりも短いコンティグの数。Number_of_alignable_readsは少なくとも1つのアラインメントが検出されたリードの数(テストデータなので3800少ししかない)で、Number_of_alignments_in_totalがその総数。

 

計算の詳細は論文と公式ページで確認してください。

 

引用
Evaluation of de novo transcriptome assemblies from RNA-Seq data

Bo Li, Nathanael Fillmore, Yongsheng Bai, Mike Collins, James A Thomson, Ron Stewart and Colin N Dewey

Genome Biol. 2014 Dec 21;15(12):553

 

RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome
Li B, Dewey CN.

BMC Bioinformatics. 2011 Aug 4;12:323.