macでインフォマティクス

macでインフォマティクス

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

review article要約 SNPs callingビギナーズガイド

 8/24 誤字修正

A beginners guide to SNP calling from high-throughput DNA-sequencing data

(Andre ́ Altman et al.,  2012)より

 

 ハイスループットDNAシークエンシング(HTS)は、ライフサイエンスにおいてますます重要になっている。その最も顕著な用途の1つは、全ゲノムまたは全エキソン領域(exome)のようなターゲット領域のシーケンシング解析である。ここでの目的は、一塩基多型(SNP)などの遺伝的バリアントの同定である。rawシーケンシングリードからのSNP抽出には、多くの処理工程および多様のツールの適用が必要である。この論文では未処理のHTSデータからSNPを呼び出すパイプラインの不可欠なビルディングブロックを確認する。パイプラインには、クオリティ管理、リファレンスゲノムへのショートリードのマッピング、ベースクオリティ再較正を含むアラインメントの視覚化および後処理が含まれる。パイプラインの最後のステップには、SNPコールプロシージャと、SNP候補のフィルタリングが含まれる。

 

http://folk.uio.no/jonkl/StuffForMBV-INFx410/Articles/AAltmann.pdf

f:id:kazumaxneo:20180822230838p:plain

SNPsコールパイプライン。論文より転載。

 

本文要約

Introduction

2001年に公開された最初のヒトゲノム全シーケンス配列決定は、大規模な研究コンソーシアムによってのみ達成できるeffortであったが、10年の時間と大きな財源を必要とした(Consortium 2004; Lander et al、2001; Venter et al、2001)。得られたヒトゲノムの青写真は、ゲノムワイド関連研究およびマイクロアレイを用いたゲノムワイドな遺伝子発現プロファイリングのようないくつかのフォローアップ技術を促進した。このラインに沿った最新の技術進歩、すなわち次世代シークエンシング(NGS、HTSとも呼ぶ)は、単一の研究室で2〜3週間以内に、比較的低コストで、単一個体の全ゲノム配列をシーケンスすることを可能にする。

  • ヒトゲノムプロジェクトのシーケンスデータは、伝統的なキャピラリーベースのサンガーシーケンシング技術を用いて生成され、500-1,000ヌクレオチドのリードを生成した。今日最も広く適用されているHTSプラットフォームのワークフローは、より小さいセグメントのDNAを使う。これらの断片のヌクレオチド配列を、合成またはライゲーションのいずれかによって決定する。使用されるヌクレオチドは、塩基の同一性を示す光シグナルがそれらの組込みの際に放出されるように修飾され、これは、例えば、化学発光または蛍光を利用している。画像キャプチャ装置は、成長中の第2のストランドに沿って生成された光信号を記録する。しかし、DNAの一本鎖での合成またはライゲーションは、記録することができる十分に強い光シグナルを放出するには不十分である。結果として、単一の断片は、合成またはライゲーション段階前に増幅され、培地上に固定されてコロニーを形成する。これらのコロニーから放出された光シグナルは、画像捕捉装置によって記録され、ヌクレオチド配列決定のために分析される。現在利用可能なプラットフォームは、どのようにコロニーが形成され増幅されるかと、最終的にどのようにヌクレオチド配列が決定されるかが異なる(HTS技術のレビューについてはMetzker 2010を参照 pubmed)。
  • 高いスループットは、これらのコロニー数百万を並行して配列決定することによって達成される。伝統的なサンガーシークエンシングと比較して、プラットフォームによって生成されるreadouts(しばしばリードと呼ばれる)はかなり短く[例えば、SOLiDプラットホームについては75ヌクレオチド(Shendure et al、2005)、illuminaナプラットフォームについては150ヌクレオチド(Federco et al、2006)、454パイロシーケンシングプラットフォーム(Margulies et al、2005)500ヌクレオチド]、さらに、各プラットフォームはシーケンシングワークフローの特徴を反映したシーケンシングエラーを持っている。したがって、サンガーシーケンシングと比較して、HTSははるかに短く、品質が劣ることも多い。このことは、得られた配列を下流の分析でどのように処理するかに大きな影響を与える。HTSの技術はon goingで開発が進行中であり、現代世代の技術は、技術的問題の解決を目指し、より近代的なアプローチ(すわわちthe ‘‘next–next generation’’)に取って代えられる 。最近の方法は、光学システムから、ナノテクノロジー(Clarke et al、2009)、半導体(Rothberg et al、2011)、および顕微鏡(Tanaka and Kawai、2009)に依存したシステムに向かっている(論文執筆時点)。
  • この原稿では、従来のアレイベースのゲノムワイド関連研究と密接に関連するアプリケーションであるSNPの決定に焦点を当てている。 Nielsenらのレビューですでに詳細に説明されているSNPコール方法の確かな背景には集中せず、HTSデータからSNPをコールするための実地ガイドを提供することを目指している。この目的のために、1,000ゲノムプロジェクト(コンソーシアム2010)(サンプルID:NA12287; runID:ERR034546; ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/data/ NA12287/sequence_read/)から得られた全エキソームシーケンシングデータセットを用いてパイプラインに沿ったデータ処理を説明する。

 
The SNP calling pipeline

SNPコールパイプラインは、論文図1のフローチャートで視覚化された7つのステップから構成されている(PDF link)。パイプラインの最初のステップはbase callingと呼ばれていて、シーケンシングプロセス中に取得されたイメージを評価し、短いリードを生成するステップである。これに続いて、リードの初期クオリティ管理、基準配列へのアライメント(*1)、およびアライメントの後処理が行われる。これらの手順は、ほぼすべてのHTSアプリケーションで共有されている(*2)。残りの3ステップは、SNPコールパイプラインに特有で、すなわちクオリティスコア再較正、SNPコールとSNP候補のフィルタリングである。ここでは基本的に単一個人のExomeシーケンシングのデータ解析の流れを概説する。

 

ステップ0: ベースコール(base calling)

  • 上述したように、撮像装置は、合成またはライゲーションプロセスによって生成された光信号を、新たに生成されたストランドに記録する。画像データの取得後、これらの記録されたシグナルはヌクレオチド塩基に変換されなければならない。 SOLiDプラットフォームの場合、光シグナルは単にcolorsと呼ばれる隣接ジヌクレオチドをコードする。ヌクレオチド配列と色とを区別するために、base spaceとcolor spaceという用語が利用されている。
  • さらに、統計的モデルにより、ヌクレオチド自体に加えて各ベースコールの確実性の尺度が提供される。これらの統計モデルは、記録された画像からのシグナル強度、シーケンシングサイクル数および他の配列のコロニーまでの距離などの情報に関する誤差推定値に基づいている。これらの確実性は、通常、Phredライクなクオリティスコア、すなわち、ベースコールの予想されるエラー確率の常用対数として表される。

f:id:kazumaxneo:20180822010015p:plain

  • この式を使用すると、例えばエラー率5%はPhredスコア13に変換される。このステップはbase callingと呼ばれ、通常はシーケンシングプラットフォーム自体によって自動的に実行される。ここでも、各シーケンシングプラットフォームは、基礎となるシーケンシング手法に固有の課題を解決する必要がある。したがって、base callingステップは各プラットフォームに特化している。しかし、近年、いくつかのアルゴリズムが提案されている [Illumina (Kao and Song 2011; Kircher et al. 2009), Roche 454 (Quinlan et al. 2008), SOLiD (Wu et al. 2010)]。しかし残念なことに、メーカーのベースコールアルゴリズムを変更または改良するには、簡単には理解できないほど深い知識が必要で、操作方法も知っている必要がある。したがって、ほとんどのユーザは、依然としてシーケンシングプラットフォームによって提供されるbase callingアルゴリズムに完全に頼っている。

 

ステップ1: クオリティコントロール

 ほとんどのプラットフォームは、FASTQ(Cock et al。2009)のようなフラットファイル形式で直接読み取ったデータを提供するか、またはネイティブ出力形式を準標準FASTQに変換するためのツールを提供する。生成されたシーケンスデータの品質をチェックすることは、ベースラインまたはカラースペース内の実際のシーケンスデータを処理するパイプラインの最初のステップである。

  • 各シーケンス位置でのクオリティスコアの分布は、ランの全体的な品質にとって最も興味深い品質パラメータの1つである。典型的には、manufacturerのベースコールソフトウェアは、既に クオリティに関する最初のオーバービューを提供している。より徹底したクオリティオーバービューについては、SolexaQA(Cox et al。2010)やFastQC(http://www.bioinformatics。babraham.ac.uk/projects/fastqc/)などの自由に利用できるツールを適用することができる。
  • FastQCは現在のすべてのHTSプラットフォームのデータを扱うことができ(論文執筆時点)、ユーザーフレンドリなグラフィカルインターフェイスを提供し、過剰発現された配列、期待されるGC含量からのずれ、読み取り位置あたりのヌクレオチドの分布を確認することができ、サンプルの調製およびシーケンシング中に発生する可能性がある問題の特定するのに役立つ。
  • 論文図2のSolexaQA出力例は、リード長が長くなるとエラー確率が増加することを示していおり、HTSプラットフォームでは一般的なことである。リードのトリミングは、シーケンシングエラーを含む可能性のあるリードの最後の塩基を除去することによって、結果として、マッピング可能なリード数を増加させる。したがって、リードトリミングはすべてのアライメントが重要な分析に価値がある。
  • イルミナは、例えば、製造業者のソフトウェアによるクオリティ管理を受ける。 SOLiDプラットフォームの場合、クオリティ管理は提供されない。このプラットフォームは、不十分なクオリティりリードが参照の配列とアライメントしないという事実に依拠している。ただし、マッピング時間を短縮するために、SOLiDプラットフォームの平均クオリティスコアが10未満のリードを破棄することは理にかなっている。この作業はPRINSEQ(Schmieder and Edwards、 2011)のようなツールを介して、またはSHRiMP2(David et al、2011)などのアライナーを介して明示的に行うこともできる。

 

ステップ2:アラインメント/マッピング

  • ほとんどすべてのパイプラインにおける次のステップは、リファレンス配列、すなわち本稿ではヒトゲノムに対するリードのアライメントである。少しの偏差(例えば、SNP、indels)およびシーケンスエラーを含む数百万の短いリードを、リファレンス配列またはシーケンスデータベースにアライメントさせるため、多くの効率的なアルゴリズムがもたらされた。
  • 簡単に言えば、タスクを解決するために2つのアプローチが一般的に使用される。最初のものは、効率的なデータ圧縮のために、無損失のBurrows-Wheeler変換(BWT)(Burrows and Wheeler 1994)(PDF)(*3)を適用する。2つ目のアルゴリズムは、アライメントステップを加速するためハッシングに依存している。ハッシングの使用は、リファレンス配列内の部分配列位置情報への迅速なアクセスを可能にする。ハッシュベースのアライナは、シーケンスリードをハッシングするEland(IlluminaのCASAVAスイートの一部)、またはSOAP(Li et al。2009b)などのリファレンス配列をハッシングするものがある。
  • 一般に、アライメントツールの選択とそれに対応する設定は、結果に大きな影響を与える。これは、誤って位置合わせされたリードがアライメントからのアーティフィシャルな偏差を招く可能性があるため、SNPコールでは特に当てはまる。これらの偏差は、下流の処理において誤ってSNPとして分類されることがある。
  • ミスアライメントリードは、他のアプリケーションではバリアントコールより重要性が低くなる。これには、主に遺伝子発現プロファイリングのような定量的解析が含まれる。ここではリードのシーケンス内容は2次的関心事にすぎない。それにもかかわらず、依然としてミスアライメントは推測される発現レベルを歪める可能性がある。
  • 自由に利用できるソフトウェアに焦点を当て、Bowtie(Langmead et al、2009)とその後継Bowtie2、BWA(Li and Durbin、2009)、MAQ(Li et al、2008)、およびstampy(Lunter and Goodson、2011)を使用して、サンプルリードをマッピングした。UCSCによって提供されたヒトゲノムバージョン18に対してアライメントされた(Kent et al、2002)。 Bowtieはギャップのあるアライメントを実行できないため、短いindelsを見つけることはできない。Bowtie2の後継Bowtie2には、ギャップアライメントなどの多くの便利な機能が含まれている。しかし、Bowtie2はSOLiDのcolor spaceのリードに適用することはできない。 
  • 表2はペアエンドリードをヒトゲノムにアライメントさせるために単一のCPUに必要な時間を示す。アライメントステップをさらに加速するために、多くのアルゴリズムは、単にプログラムの実行時にパラメータを追加するだけで、より多数のCPU(例えばBowtieやBWA)上に容易にタスクを分配することができる。
  • BWTベースのアルゴリズム(ここではBWA、Bowtie、Bowtie2)のハッシュベースのアルゴリズム(ここではMAQ、stampy、stampyBWA)に対する処理上の利点は、処理速度である(表2を参照)。しかし、BWTアルゴリズムはハッシュベースのアライナほど感度が高くないため、変動性の高い領域でマッピングバイアスを導入する可能性がある(詳細な感度/特異性分析については、Lunter and Goodson 2011を参照)。著者らの意見では、BWTアプローチで導入されたスピードアップゲインと、ハッシュベースのアラインメントの感度を組み合わせたstampy のようなアプローチは良い妥協点である。
  • リードとリファレンス間の完全一致のみが許される場合、参照とシーケンスされたリードとの間に何らの相違も見出さないので、SNPは検出されない。逆に、リファレンスとリード間に多くのミスマッチが許容されると、下流の解析で偽陽性SNPが高くなる。したがってアライメントされたリードの数を最大にすることは、最良の戦略ではない。受け入れられるミスマッチの最良の数を選択することはまた、種に大きく依存する。例えば、Mus musculus(ハツカネズミ)の系統は利用可能なリファレンスからかなり大きく外れる可能性がある(Keane et al、2011)。一方、ヒトのサンプルは変化が少ない傾向がある(コンソーシアム2010)。残念ながら、このステートメントはゲノム全体では成立しない。例えば、主要組織適合複合体(MHC)は、ヒト個体間で高い変動性を示す。したがって、この領域で良好なアラインメントを実行することは、一般的に非常に困難になる。
  • リファレンスへのリードのアライメント後は、多くのアルゴリズムは結果をsequence align- ment/map (SAM) format (Li et al、2009a)に格納する。SAMフォーマットは、アライメントされた各リードに関する情報、特に、リードのリファレンス上の位置、リードの向き、アライメントクオリティ、およびリードのさらなるアライメントの可能性を格納している(下のPDF参照)。アライナの出力がSAMとは異なるフォーマットである場合、SAMフォーマットへの変換のために第三者のツールが利用可能であるかもしれない。
  • SAMフォーマットとそのバイナリバージョンであるBAMフォーマットは、現在ではアライメント結果を格納するための準標準となっている。したがって、多くの下流処理ツールはSAM / BAM形式に依存している。さらに、SAMおよびBAM形式で保存されたマッピングを効率的に操作するツールが公開されている。最も広く適用されるツールキットはSAMtools(Li et al、2009a)、GATK(McKenna et al、2010)、Picard(表1参照)である。
  • マッピング後、再度アラインメントをチェックすることを勧める。例えば、SAMtoolsのflagstatコマンドまたはPicardのCollectAlignmentSummaryMetricsモジュールを使いマッピング統計値、すなわちリファレンスにうまくマッピングされたリードの割合を計算することができる。さらに、ペアエンドリードで作業している場合、正常にペアでアライメントされたリードの割合とインサート・サイズの分布は関心のあるパラメータである(論文図3 インサートサイズ)。著者らのケースでは、以下の処理のために1,000 bpより大きなインサートサイズを選択する。
  • ターゲットリシーケンシングからのデータ分析の場合、ターゲット領域由来リードのoff-target領域からのエンリッチメントは重要である。 PicardのCalculateHsMetricsモジュールは、この比率を計算するのに有用である。サンプルデータでは、アラインメント方法にかかわらず、34倍をエンリッチメントを達成していた。これはエンリッチメントプロセスが成功したことを示す。これに加えて、CalculateHsMetricsモジュールは、いかなるリードもカバーしていなかったターゲット領域の割合に関する情報(データの約15%)、ターゲット領域中で最小カバレッジ10以上の割合(73%)、また、平均カバレッジ情報(約130)を提供する。
  • 通常、全ゲノムシーケンシング実験の目視検査は現実的ではない。しかし、ターゲット領域内のアラインメントを例えば、SAMtoolsのviewコマンドを使用して抽出し、ゲノムブラウザ(例えばIntegrative Genomics Viewer : IGV(Robinson et al、2011 ))でその特定の領域のみを視覚化することができる。論文図4は、IGVを用いたSLC6A15遺伝子座における読み取りのアライメント視覚化を示す。
  • 視覚化は、アライメントされたリードがエキソン領域に集中していることを明らかにする。したがって、少なくともこの目視検査された領域では完全なエクソンシークエンシングが成功しているが(論文図4を参照)、隣接するイントロン領域も高度にエンリッチされているため、エンリッチメントプロセス全体があまり正確ではないことに注目することは興味深い。しかし、隣接しているイントロン領域のカバレッジフォワードリードまたはリバースリードのいずれかによって支配されており、それゆえ、その領域からコールさえるSNPにバイアスを導入する可能性がある。
  • IGV以外にも、アラインメントの視覚的検査を可能にするツール、例えばGenomeView(Abeel et al、2012)、SAVANT(Fiume et al。2010)などのツールがある。後者はユーザーが提供するソフトウェアモジュールによっても拡張することができる。
  • 一般に、アラインメントパラメータを微調整するには、ある程度の労力が必要となる。ほとんどのプロジェクトには大量のデータが含まれているため、無作為に選択されたデータのサブセット、たとえば10%のデータでのみ作業することで時間とハードドライブのスペースを節約できる(*4)。これは、ステップ1のクオリティ管理にも当てはまる。

 

ステップ3:アライメント後のポストプロセシング

  • 実際のバリアントコールの前に、アルゴリズムはアライメントをその染色体位置に対してソートする必要がある。これは、SAMtoolsやPicardのようなツールを使って簡単に行うことができる。次に、ライブラリーのPCR増幅ステップで(同じ位置で開始し、同じインサート長を有するリードまたはペアエンドリードを導入する)PCR duplicationが発生することがあるので、そのようなPCRアーチファクトをマークする。SAMtoolsとPicardはこのタスクを解決する手段も提供している。
  • 次の後処理ステップは、すべての非ユニークなアライメントのリード、すなわち最適アライメントを複数有するリードの除去である。このような場合、どのサイトから実際にシーケンスされたのかを判断することができないからである。
  • 最後に、小さなindelsの周りのリードのリアライメントをするのが一般的である。たとえば、GATKソフトウェアはこのリアライメント手段を提供している(*5)。


ステップ4:クオリティスコアの再較正

  • これまでの研究では、シーケンシングプラットフォームによって発行されたPhredのようなクオリティスコアが真のエラーレートから逸脱することがあることが示されている(Li et al、2009b link)。現代のSNPコールアルゴリズムでは検査対象のサイトをカバーするベースのPhredスコアをスコアリング関数に統合するため、正確なクオリティスコアを持つことが必須である(ステップ5を参照)。
  • クオリティスコアの再計算を行う最初のソフトウェアはSOAPsnpであった(Li et al、2009b)。このアプローチは、既知SNP情報なしでリファレンスゲノム中の部位を探索する。SOAPsnpは、探索したサイトの真のクオリティ推定値としてempirical(経験的)なミスマッチ率を計算する。より正確には、所与のシーケンシングマシンによって提供されたクオリティスコア、シーケンシングサイクル(すなわちリードにおけるベースの位置)および置換タイプについて、所与されたリファレンスに対するミスマッチ率を計算する。このミスマッチ率は、再調整したクオリティスコアとして使用される。
  • 類似の概念に基づいて、GATKソフトウェアはまた再較正機能も提供する:第1に、塩基は、いくつかの特徴(例えば、rawクオリティ、ジヌクレオチドの内容)に関してグループ化される。第2に、そのようなカテゴリーごとに、empirical(経験的)なミスマッチ率が計算され、rawクオリティスコアを補正するために使用される。 GATKの再較正機能は、さまざまなプラットフォームのシーケンシングデータに適用できる。図5は、BWAでアライメントされ、GATKで再較正されたクオリティスコアとオリジナルクオリティスコアを示す。
  • わずかな差異しか予想されないので、ステップ3とステップ4を入れ替えてもよい。

 

ステップ5:バリアントコールおよびジェノタイプコール

  • 初期のバリアントまたはSNPコールアプローチは、単一部位(例えば、Wang et al。2008)における高品質ヌクレオチド存在量のカウントに全面的に頼っていたが、最近のアプローチは、確率論的枠組みの中でいくつかの情報源を統合する。この手順は中程度 ~ 低カバレッジの領域でSNPコールを容易にする。たとえば、わずか5リード分のカバレッジしかない部位のSNP可能性など。さらに、これらの確率論的アプローチは、バリアントコールに関する不確実性を定量化するための自然な方法を提供する。使用される統計モデルの詳細は、オンラインでアクセス可能なSupplementary Materialに記載されている。
  • 統計的フレームワークの1つの大きな利点は、所与の位置におけるSNPの事前確率の使用である。これらの事前確率は、dbSNP(Sherry et al、2001)などの報告済みSNPデータベースを利用するか、同時に複数個人からSNPコールを行うことによって導き出せる。 SAMtoolsとGATKで実装されたSNPコールルーチンは、複数サンプルのSNPコールをサポートしている。
  • 連鎖不平衡(LD)情報を組み込むことにより、さらなる改善が達成され得る。この機能は、例えばBeagleソフトウェアに実装されている(Browning and Yu、 2009)。しかし、exomeデータ全体を扱う場合、LD構造の大部分が欠落しているため、このステップを適用することでの改善は期待できない。
  • 利用可能なSNPコールアルゴリズムに関する徹底的なレビューは、Nielsen et al、(2011)(link)およびその中の参考文献を参照されたい。サンプルデータでは、SAMtoolsとGATKという2つのSNPコールプログラムを適用した。この例では単一個人のデータしか扱ってないため、マルチサンプルSNP callingによる精度向上は行えなかった。

 

ステップ6:SNP候補のフィルタリング

  • フィルタリングは、偽陽性SNPコールの数を減らすための不可欠なステップである。一般的には、Hardy–Weinberg equilibrium(HWE)(wiki)、最小および最大シーケンスデプス、indelsへの隣接度合い、strandバイアスなどからの偏差をチェックしてフィルタリングを適用する。フィルタリングは候補リストから実際のSNPも削除する可能性はゼロではないが、 SNPコールのアーティファクトを最小にするために重要なステップである。
  • フィルタリングは、GATK、SAMtools(スクリプト 'vcfutils.pl')とVCFツール(Danecek et al、2011)によって提供される。特にGATKは、SNP候補のフィルタリングを含むバリアントコールパイプラインの「ベストプラクティス」セッティングを提供している。
  • サンプルデータでは、SAMtoolsのデフォルトフィルタ(表1を参照)とGATK VariantFiltrationを使用した。ここでは、GATKベスト・プラクティス・ガイドラインバージョン3の推奨設定に従い、サンプル数が少ない、つまり単一の個人のため、ハード・フィルターを使用し、30.0未満のクオリティ、5.0未満のカバレッジのクオリティを有するSNP、または6-bp以上のホモポリマー内のSNPを破棄した。
  • ほとんどのSNPコールツールには、VCF形式のデータを生成するオプションがある(Danecek et al、2011)。 VCF形式は、染色体位置、リファレンス塩基、同定されたオルタナティブの塩基、または3対立遺伝子の塩基のような各SNP候補の基本情報を記録する。さらに、SNPコールのクオリティに関する情報ならびにコールに利用可能なシーケンスデータの量が格納される。 VCFtools(紹介)は、複数のVCFファイルをマージしたり、定義された領域内のSNPを抽出するなど、VCFファイルを簡単に操作できるようにする。
  • 論文表4は、初期フィルタリングステップの後に、アライメントアルゴリズムとSNPコーラーで呼び出されたSNP数の比較となる。 1つのアライメントアルゴリズムのみ与えられると、2つの異なるツールで呼び出されたSNPは主に重なり合い、すなわち、SNPの約85%が共有された。しかし、コールされたSNP数は予想された数を1桁上回った。これは、ターゲット領域の外側のリード由来の可能性が高い。
  • したがって、VCFtoolsを使用してさらにフィルタリングした。SNPがエンリッチされたターゲット領域±50bpに存在し、SNPを有する位置は、最低10のシーケンスデプスを示すフィルタリングを行った。表5は、第2のフィルタリング工程の後のSNP数を示す。同じアライメント由来であれば、再び、約80%のSNPは両方のSNPコーラーによって重なり合った。より正確には、GATKはSAMtoolsと比較して約5,000の追加SNP候補を生成した。
  • 表6は、SNPコーラーが固定で、アラインメントアルゴリズムが異なる場合の影響を示している。ここで、大部分のSNP(約85%)は、4つの利用されたアライナー由来アライメント全てにおいて見出された。したがって、アライナーおよびSNPコーラーの両方が、結果に有意な影響を示した。しかし、SNPの大部分は、ここで使用されたツールの任意の組み合わせで発見された。

 

SNPデータの意味を理解する

  • ExomeデータのSNPコールプロセスでは、約24,000のバリアントを生成した。したがって、各バリアントを手動でフォローアップすることは明らかに範囲外になる。全ゲノムリシーケンシング研究では尚更そうである。この要件のために、自動化されたバリアントアノテーションツールが開発されている。
  • 例えば、ANNOVAR(Wang et al。2010)は異なる種のバリアントに注釈を付けるコマンドラインインタフェースを提供している。このツールは、UCSCゲノムブラウザを介して入手可能な情報に依存し、非同義アミノ酸置換の可能性のある機能的結果を予測するSIFTスコア(NgおよびHenikoff 2003)、などの事前計算スコアをさらに提供する。ANNOVARを使用すると、さらに、dbSNPまたは1,000ゲノムプロジェクトの情報を使用して、既知のSNPを簡単にフィルタリングすることができる。
  • オンラインSupplementary materialの表S2は、アライナーとSNPコーラーの組み合わせで呼び出されたSNPの、ANNOVARによるfunctional annotationの結果を示している。発見されたSNPの約50%は実際にエキソン上のバリアントである。残り50%は、ターゲット領域の周囲の±50bpの領域に由来する。
  • SNPデータを分析する際のさらなる関心のある値は、transitionsとtransversionsの比率である。この統計は、VCFtoolsを使用して簡単に計算できる。著者らのケースでは、GATKでコールされたSNPのTs / Tv比は2.1〜2.4で、SAMtoolでは2.6 ~ 2.7であった。したがって、GATKでコールされたSNPは期待される2に近く、SAMtoolsでコールされたSNPと比べてtransitionsが強化されているように見える。
  • 生成されたSNPデータの後処理および解釈は大きな課題であり、この重要な課題に関連する努力は過小評価されるべきではない。 ANNOVAR以外にも、SNPデータの解釈を支援するツールがある。たとえば、http://www.gen2phen.org/wiki/tools-predicting- overal-functional-consequences-snpsに掲載されているツールを参照。グラフィカルユーザインタフェースを提供するアノテーションソフトウェアの例は、sequence variant analyzerである(Ge et al、2011)。さらに、これらのフリーツールに加えて、多くの商用スイートが存在する。
  • バリアントアノテーションだけでなく、重要なassociationsを見つけるための統計も適合させる必要がある。全ゲノムおよび全exomeシークエンシング研究がますますレアなSNP、すなわち1%未満の小さな対立遺伝子頻度を有するSNP(wiki)を生成するので、標準的な統計的アプローチは、現在利用可能で現実的なサンプルサイズとの有意な関連性を見出すのに十分な能力を有さない。したがって、同じ機能的注釈を有するSNP、同じ生物学的パスウェイ内のSNP、または互いに近いSNP(Cohen et al。2004)、のような類似の特性を有するバリアントは、しばしば、力を増加させるために代理の変数に分類される。これらの代理変数は、統計的有意性を分析する対象となる。
  • Bansal et al (2010)は、レアバリアントを分析するための現在の統計的アプローチに関するレビューを提供している(link)。

 

結論

  • 多くの他のHTSパイプラインと同様に、SNPコールパイプラインは、多くの異なるステップからなる。通常、これらの手順は、Sun Grid Engine(Gentzsch 2001)などのキューイングシステムと組み合わせてシェルスクリプトを使用して統合できる。
  • コンピュータ体験の少ないユーザーが多くのHTSタスクを実行できるもう1つのオプションは、Galaxy Webサービスである(Goecks et al、2010)。 Galaxyは、データ集約型生物医学研究のためのWebベースのプラットフォームであり、開発者は自由にアクセス可能なWebサービスとして提供する。
  • 新しいパイプラインのセットアップと統合には、徹底的なコンピュータ知識が必要になる。このレビューで取り上げなかったオプションは、市販のソフトウェア製品の使用である。多くの企業が、HTSデータの基本的かつ高度な分析を可能にするソフトウェアパッケージを開発してきている。これらのソフトウェアスイートは、一般に、グラフィカルユーザインタフェースに基づいている。しかし、使いやすさは、柔軟性の欠如によって損なわれる可能性がある。どのような場合でも、ソフトウェアライセンスに関連するコストは無視できない。

この論文では、シーケンシングリードからバリアントアノテーションに至る、SNPコールパイプラインを説明した。この例では、ツールおよびパラメータの選択が最終結果に重要な影響を及ぼすことを示した。したがって、SNPコールパイプラインの結果はまだ確定されていない。例えば、推奨される戦略は、SNP候補を生成するために異なるアライナーとSNPコーラーを使用することである。信頼性の高い候補は、複数の組み合わせで共通して表示される候補である。HTSシステムの作業は学際的な取り組みであり、データ生成は主に研究所中心で行われるが、シーケンシングデータの初期処理はバイオインフォマティクスの領域に入り、ある程度自動化されている。しかし、結果の解釈は、データからの最大限の洞察を導き出すために、生物学とバイオインフォマティクス間での緊密な相互作用を必要とする。これらの3つの領域(データの生成、処理、解釈)に関連する作業は無視できないものであり、研究プロジェクトの成功を保証するための専用リソース(人的資源を含む)を必要とするものである。


Further reading

 

引用

A beginners guide to SNP calling from high-throughput DNA-sequencing data

Altmann A, Weber P, Bader D, Preuss M, Binder EB, Müller-Myhsok B.

Hum Genet. 2012 Oct;131(10):1541-54.

http://folk.uio.no/jonkl/StuffForMBV-INFx410/Articles/AAltmann.pdf

 

参考

*1

このブログではよくアライメントと書いてしまっているが、方式の性質上、多数のリードをリファレンスに当たる場合はマッピングと呼ぶのが妥当で、1つのリードをリファレンスに当たる場合はアライメントと呼ぶのが妥当と思われる。

 

 *2

例外もある。k-merスペクトラム比較による、アライメントフリーのバリアントコールなど。

 

*3

explaining BWA and Bowtie algorithm


*4 

ランダムサンプリングにはseqkit(紹介)が利用できる。

 

*5

最近のマッパーやバリアントコールツール(例えばGATK4)は改良が進んでおり、リアライメントが必ずしも必要というわけでもなくなってきている(データ、手法による)。

 

Hands-on Tutorial on SNP Calling 

http://www.transplantdb.eu/sites/transplantdb.eu/files/SNPcallingTutorial-Haberer_010715.pdf