macでインフォマティクス

macでインフォマティクス

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

de novo transcriptome解析のクラスタリングとclosely rellatedな種の情報を用いたアノテーションを行う Grouper

2018 Bioinformatics annotation human genome large genome RNA seq RNAseqの定量 de novo transcriptome

 

　シーケンシング技術の進歩により、モデル生物の範囲を超えてトランスクリプトームを効率的かつ正確に探索することが可能になった（Ekblom and Galindo、2011; Marioni et al、2008）。トランスクリプトームシークエンシングは、高品質のリファレンスゲノムを持たない種でも遺伝子発現を理解するための扉を開く。これは既知の種の大多数を構成している（Martin and Wang、2011）。トランスクリプトームシークエンシングのいくつかの一般的な用途には、変異検出、フュージョンおよび選択的スプライシングイベントの発見、ならびにdifferential expression analysisが含まれる（Soumana et al、2015； Stubben et al、2014）。転写物の配列が先験的に知られていない生物のショートリードシーケンスでは、決定的な最初のステップは、ショートシークエンシングリードを使ってトランスクリプトームをアセンブリすることである。
　トランスクリプトームは典型的にはディープにシーケンシングされ、その結果、数千から数億のリードが生じ、それを元の転写産物シーケンスを再構築するために de novo assemblyというプロセスを通してアセンブリする必要がある。これらの転写産物は後の分析のための「リファレンス」として後で機能するであろう。例えば、differential expression分析パイプラインでは、リードが最初にアセンブリされたトランスクリプトにマッピングされ、トランスクリプトの量が推測される（Li and Dewey、2011）。
Trinity（Haas et al、2013）、Oasis（Schulz et al、2012）、Trans-ABySS（Robertson et al、2010）など、冗長性​​とリード間の重複を利用してショートリードから完全長のトランスクリプトをアセンブリするための多くの一般的なツールがある。同様に、これらの方法で生成されたアセンブリの品質を向上させるために使用できるツールがある（Cabau et al、2017; Durai and Schulz、2016）。ただし、強力なアルゴリズムとヒューリスティックが採用されているにもかかわらず、最終的な出力シーケンス（つまりコンティグ）は、完全長のトランスクリプトを表していないことがよくある。言い換えれば、多くの転写産物は重複しないコンティグのセットに断片化されている。これは、シーケンシングされたリードにおけるエラー、オルタナティブスプライシングおよびパラロガス遺伝子から生じるトランスクリプトーム複雑性、シーケンシングされた分子の長さにわたる不均一なまたは不十分なカバレージ、アセンブリ（例えば、解決される計算上の問題が根底にある生物学とうまく一致しない）ならびに以下に使用される基礎的アプローチの欠陥を含む。これらのすべてが、アッセイされた実際の転写産物のセットよりもはるかに大きい出力コンティグ集団をもたらし得る。
　潜在的に低い品質のde novoアセンブリアウトプットのために、転写産物レベル（すなわち、コンティグレベル）分析の代わりに遺伝子レベルのdifferential expression分析を追求することがより有望でありそして頑強である。遺伝子レベルの情報を得るためには、同じ転写産物の一部を表すコンティグと同じ遺伝子のアイソフォームを表すコンティグを推測し、それらをグループ化する必要がある。これは、Davidson and Oshlack（2014）、Ptitsyn et al（2015）、Srivastava et al（2016）が提案しているものである。 Corsetは、コンティグの共有カウントに基づいて階層的クラスタリングを計算し、RapClustは、コンティグの存在量に基づいて重み付けされた、リードを共有するエッジとコンティグを表すmapping ambiguity graphと呼ばれるスパースグラフをクラスタリングした。さらに、これらのクラスターの意味のある解釈と様々な分析のために、de novoアセンブリのコンティグが何を表しているかについてのいくつかの概念を持つことが重要である（Garber et al、2011）。多くのケースで、我々は研究中の非モデル生物とclosely relatedな種からのゲノムまたはトランスクリプトームを使いアノテーションを付ける。この情報を活用して、de novoアセンブリのコンティグに正確にアノテーションを付け、differential expressionな遺伝子（DEGs）の機能に関する理解を深めることができる。伝統的に、BLASTのバリエーション（Altschul et al、1990）を用いてこのアノテーション付けを行い、次いで遺伝子オントロジー（GO）analysisを完成させる（Ji et al、2012； Parchman et al、2010）。
　現在、下流分析を改善し、そのような分析からの結果の解釈を容易にするために、新たなアセンブリを処理するいくつかの方法が存在する。しかし、それぞれの方法には限界があり、Corset以外には、クラスタリングや新規のコンティグへのアノテーション付けのための完全なパイプラインを提供するツールは存在しない。 RapClustの枠組みを包み込んで改良し、de novoトランスクリプトームアセンブリの正確で効率的な処理に使用できるGrouperを紹介する。Grouperの基礎となるアルゴリズムはde novoアセンブリからグラフに情報を転送するため、アセンブリされたコンティグの直接解析を超え、将来のde novoトランスクリプトーム解析の有用な方向性を示唆している。完全なGrouperパイプラインは、differential expression analysesに使用できる定量情報も提供する。 Grouperソフトウェアはhttps://github.com/COMBINE-lab/grouperにて自由に利用できる。

 

 

　 Grouperには、クラスタリングモジュールとラベリングモジュールの2つの主要モジュールがある。 前者はツールRapClustに基づいており、迅速な転写産物レベル定量のためSailfishまたはSalmonツールの下流で実行するように設計されている。 Grouperのラベリングモジュールは、closely relatedな種のゲノムのアノテーション情報を使い、de novoアセンブリ内のコンティグにアノテーションを付けるために使用できる。 

 

インストール

ubuntu16.04にて、condaの仮想環境を作ってテストした（conda create -n Grouper python=2.7.15）。

依存

  • python2
  • Click
  • PyYAML
  • Pandas
  • NumPy
  • Networkx v1.11
  • mcl

本体 Github

#依存も含めpipで導入可能
pip install biogrouper

#MCLがなかったので別途導入した
conda install -c bioconda -y mcl

> Grouper -h

$ Grouper -h    

Usage: Grouper [OPTIONS]

 

Options:

  --config TEXT  Config file describing the experimental setup  [required]

  -h, --help     Show this message and exit.

 dockerイメージ（未テスト）

docker pull combinelab/grouper

 

実行方法

１、Githubの解説の通り、ヒトRNA seqのバイオプロジェクトPRJNA162905 (link) の６つのSRAファイルをダウンロードする。

ここではfasterq-dumpを使っている（紹介）。forで回す。

sample=(SSRR493366 SRR493367 SRR493368 SRR493369 SRR493370 SRR493371) 

for name in ${sample[@]}; do
 fasterq-dump $name -e 8 -p
done

 

 

２、リファレンスtranscriptsのindexing（NCBIよりGRCh38のRefSeq Transcriptsダウンロードして使用 => link）

salmon index -p 20 -t GRCh38_latest_rna.fna.gz -i ref_index
  • -p  [ --threads ] arg (=2) Number of threads to use (only used for computing bias features)
  • -t  [ --transcripts ] arg Transcript fasta file.
  • -k  [ --kmerLen ] arg (=31) The size of k-mers that should be used for the quasi index.
  • -i  [ --index ] arg Salmon index.

ref_index/が出力される。これをstep3で指定する。 

 

 

３、salmonによる定量.-dumpEqと--writeOrphanLinksオプション付きで実行する。ここではGNU parallelで６ジョブを並列処理している。

parallel -j 6 "samp={}; salmon quant -i ref_index -l a \
 -1 reads/{$samp}_1.fastq -2 reads/{$samp}_2.fastq \
 -o {$samp}_quant --dumpEq --writeOrphanLinks -p 4" \
 ::: SRR493366 SRR493367 SRR493368 SRR493369 SRR493370 SRR493371

 > ls -alth

$ ls -alth

合計 111M

drwxrwxr-x  5 kazu kazu 4.0K  7月 10 18:29 SRR493371_quant

drwxrwxr-x  5 kazu kazu 4.0K  7月 10 18:29 SRR493368_quant

drwxrwxr-x  5 kazu kazu 4.0K  7月 10 18:29 SRR493370_quant

drwxrwxr-x  5 kazu kazu 4.0K  7月 10 18:28 SRR493367_quant

drwxrwxr-x  5 kazu kazu 4.0K  7月 10 18:28 SRR493369_quant

drwxrwxr-x  5 kazu kazu 4.0K  7月 10 18:28 SRR493366_quant

 

それぞれのディレクトリ（{$samp}_quant）に各サンプルの定量結果が出力される。

 

 

４、Grouperによるクラスタリング。

 yamlファイルを準備して、そこにグループ、step３で得た処理群とコントロール群それぞれの定量ファイルパス、出力ディレクトリなどの情報を記載する。Githubの通りだと、以下のようになる。

conditions:

    - Control

    - HOXA1 Knockdown

samples:

   Control:

       - SRR493366_quant

       - SRR493367_quant

       - SRR493368_quant

   HOXA1 Knockdown:

       - SRR493369_quant

       - SRR493370_quant

       - SRR493371_quant

outdir: human_grouper

orphan: True

mincut: True

 

ファイル名"config.yaml"として保存し、Grouperを実行。

Grouper --config config.yaml

このステップにはかなりの時間がかかる（*1）。ランが終わると、指定したhuman_grouper/に複数ファイルが出力される。mag.flat.clustがクラスター情報のファイル。 tximportのようなツールと互換性のある "transcript-to-gene"マッピング形式で計算されたクラスター情報が含まれている。

 

近縁種のtranscriptsのアノテーションは、transcripts.fastaを用意し、それをyamlファイルに記載することで機能します。詳細はGithubで確認してください。

 

引用

Grouper: graph-based clustering and annotation for improved de novo transcriptome analysis
Laraib Malik, Fatemeh Almodaresi, Rob Patro
Bioinformatics, Volume 34, Issue 19, 01 October 2018, Pages 3265–3272

 

関連


*1

テスト時は数日かかった。