前回に引き続き、SeqCap Epi連載の第2回です。
今回は、SeqCap Epiによる実験データをどのように解析していくのか、ご紹介いたします。用いるソフトはすべてオープンソースのフリーウェアです。
このパイプラインでは、SeqCap Epi CpGiant Enrichment Kit（ロシュ・ダイアグノスティック）を使用して、ライブラリー調整を行い、その後、NGSによるシーケンスを実行して得られたFASTQファイルから解析を行います。
解析ワークフロー
解析のワークフローは次の通りです。
| 順序 | 項目 | 使用ソフト |
|---|---|---|
| 1 | QC（クオリティコントロール） | FastQC, Trimmomatic |
| 2 | マッピング | BSMAP |
| 3 | フィルタリング | bamtoolsなど |
| 4 | 基本データの算出 | picardなど |
| 5 | メチル化解析 | BSMAP(methratio.py) |
| 6 | メチル化程度の比較 | methylKit |
| 7 | アノテーションの付加 | methylKit |
| 8 | SNPを考慮したメチル化解析 | BisSNP |
- FastQCを用いて、クオリティをチェックし、Trimmomatic*1による低クオリティリードのトリミングを行います。
- リファレンスゲノムに対して、リードをマッピングします。
- BSMAP*2でマッピングするとSAMファイルを得られますので、Picardを使用してBAMに変換します。その後、bamtoolsやBamUtilなどのソフトを使用して、適切な整形を行います。
- 次に、整形後のBAMファイルから、リード数、カバレージ、ハイクオリティリードの割合、ターゲット領域にマッピングされたリードの数などの基本的なデータを算出します。この際には、PicardやGATK toolkitを用います。
また、この際に実験が上手く行えたかどうかを確認するため、バイサルファイト変換効率を算出します。 - BSMAPのmethratio.py*3というスクリプトを用いて、BAMファイル内のメチル化を検出します。
- methylKit*4を用いてサンプル間のCpGサイトにおけるメチル化程度を比較します。
- metylKitを使用して近傍の遺伝子の転写開始点との距離や、遺伝子名などのアノテーション情報を付加します。
- BisSNP*5によりSNPを検出します。
以上が全体的な流れになっています。
このパイプラインの特徴
最大の特徴は、メチル化解析にBSMAPとBisSNPの2つのソフトを使用している点にあります。本パイプラインにおけるこれらの違いは以下です。
BSMAP：すべてのシトシンのメチル化率を算出
BisSNP：すべてのSNPおよびCpGサイトのメチル化率を算出
すなわち、これらの使い分けをすることにより、BSMAPでは、より網羅的にメチル化率を検出し、BisSNPではSNPの可能性のあるシトシンを検出することが可能になっています。このように、多くのCpGサイトをターゲットとし、SNPも考慮することを実現したメチル化解析パイプラインです。
更に、アメリエフでは、この他にも機能を追加して、メチル化程度をIGVで可視化ができるよう、解析パイプラインをブラッシュアップ中です。どうぞご期待ください。
*1:usadellab - Trimmomatic: A flexible read trimming tool for Illumina NGS dataにて公開中。トリミングの際にパラメータを設定することで、指定したアダプタ配列のトリミングが可能。
*2:bsmap - Bisulfite Sequence Mapping Programにて公開中。バイサルファイトシーケンスデータ用のマッピングソフト。
*3:BSMAPに含まれているPythonスクリプト。BAMファイルからメチル化シトシンを検出する。
*4:methylkit - R package for DNA methylation analysisにて公開中。メチル化解析データをサンプル間で比較するためのRパッケージ。
*5:Bis-SNPにて公開中。GATK toolkitに基づくSNPの検出とメチル化シトシンの検出が可能。
