【kallistoの使い方】アラインメントフリーで遺伝子発現量を計算する

RNA-seq

より高速な発現定量データ解析手法として、発現定量を主眼とした、ゲノムにアラインメントしない方法（アラインメントフリー）が広まりつつあります。

アラインメントフリーの手法では、これまでに知られたすべてのトランスクリプトーム配列を使って、そのどれに解読した配列が相当するかをカウントします。

この記事では、アラインメントフリーでRNA-seqのデータを定量化するkallistoの使い方を解説します。

目次

kallistoとは

kallistoとは、高速で正確なRNA-seqデータの定量化を行うためのツールで、k-merインデックスと呼ばれるデータ構造を用いて、リードの配列をリファレンスのトランスクリプトームにマッピングします。

アラインメントを必要としないため、リファレンスゲノムにマッピングする手法と比較して高速です。

手法 特徴 ツール
ゲノムマッピング ゲノムにマッピングするため、新規の転写産物を発見できる可能性があるが、その分データ解析に時間がかかる。 RSEM + STAR
アラインメントフリー 既知のトランスクリプトーム配列と照合して各遺伝子の発現量を計算する。非常に高速であるが、新規の転写産物は発見できない。 kallisto, salmon

 

kallistoのインストール

Mac

macの場合、homebrewで簡単にインストールすることができます。

brew install kallisto

Windows

Windowsの場合、以下の手順でkallistoをインストールします。

  1. 公式のページからWindows用のzipファイルをダウンロード
  2. ダウンロードしたzipを展開
  3. 展開したkallistoフォルダをC¥ProgramFilesに移動
  4. C¥ProgramFiles¥kallistoフォルダをパスに追加（環境変数Pathに追加）

 

データの準備

低酸素状態のヒトの細胞株RCC4-EV (DRR100656)と、低酸素状態をレスキューしたRCC4-VHL (DRR100657)RNA-seqデータを使用します。

日本のシーケンスデータアーカイブであるDDBJからデータを取得すると、FASTQファイルの圧縮ファイル(.bz2)をブラウザから直接ダウンロードすることができます。

DRR100656

DRR100657

ダウンロードした圧縮ファイルをpbzip2を使って、以下のコマンドで解凍します。bpzip2はMacの場合、brew install pbzip2でインストールできます。

pbzip2 -d DRR100656_1.fastq.bz2
pbzip2 -d DRR100656_2.fastq.bz2
pbzip2 -d DRR100657_1.fastq.bz2
pbzip2 -d DRR100657_2.fastq.bz2

kallistoの使い方

kallistoの使い方は、大きく分けて以下の2つのステップになります。

  • kallisto index: トランスクリプトームの配列からk-merインデックスを作成する
  • kallisto quant: リードの配列からトランスクリプトの発現量を推定する

トランスクリプトーム配列の取得

インデックス作成に必要なヒトトランスクリプトーム配列をEnsemblのFTPサイトから取得しておきます。

解析の目的に応じて、必要なトランスクリプトーム配列を用意しましょう。

https://ftp.ensembl.org/pub/release-110/fasta/homo_sapiens/cdna/Homo_sapiens.GRCh38.cdna.all.fa.gz

インデックスの作成

kallisto indexは、トランスクリプトームの配列からk-merインデックスを作成するコマンドです。

k-merインデックスは、リードのマッピングに必要なデータ構造です。k-merインデックスを作成するには、以下のコマンドを実行します。

ダウンロードしたトランスクリプトーム配列はgzip圧縮したまま使えます。

kallisto index -i human_transcripts_index Homo_sapience.GRCh38.cdna.all.fa.gz

ここで、human_transcripts_indexは作成されるk-merインデックスのファイル名、Homo_sapience.GRCh38.cdna.all.fa.gzはトランスクリプトームの配列が含まれるFASTAファイルの名前です。

k-merインデックスを作成する際には、以下のようなオプションを指定できます。

  • -k: k-merの長さを指定する（デフォルトは31）

発現量の推定

kallisto quantは、リードの配列からトランスクリプトの発現量を推定するコマンドです。発現量の単位は、TPM（transcripts per million）やest_counts（expected counts）などがあります。

TPMは、トランスクリプトの長さを考慮した相対的な発現量で、est_countsは、リードの数に基づく絶対的な発現量です。

発現量を推定するには、以下のようにコマンドを実行します。

kallisto quant -i human_transcripts_index -o output_dir DRR100656_1.fastq DRR100657_2.fastq

ここで、human_transcripts_indexは作成したk-merインデックスのファイル名、output_dirは出力フォルダの名前です。

発現量を推定する際には、以下のようなオプションを指定できます。

  • -b: ブートストラップサンプリングを行って、発現量の分散を推定する（デフォルトはオフ）
  • -t: 並列処理に用いるスレッド数を指定する（デフォルトは1）

私の環境では、2200万リードのサンプルをシングルスレッドで約6分で処理できました。圧倒的な速さですね！

kallisto quantの実行が終了すると、出力フォルダに以下のファイルが生成されます。

  • abundance.h5: HDF5形式のバイナリファイルで、発現量や分散などの情報が含まれる
  • abundance.tsv: TSV形式のテキストファイルで、トランスクリプトごとの発現量が含まれる
  • run_info.json: JSON形式のテキストファイルで、実行時のパラメータや統計情報が含まれる

サンプルごとの発現量の違いを解析するには、これらのファイルを解析していくことになります。

おわりに

kallistoは、高速で正確なRNA-seqデータの定量化を行うことができる優れたツールです。

ぜひ、kallistoを使って、RNA-seqデータの解析に挑戦してみてください。

RNA-seqの最新記事8件