次世代シークエンシングデータを用いた対立遺伝子の頻度と関連マッピングの推定

マイナー対立遺伝子は、可変部位の集団における頻度の低い対立遺伝子である。 まず,ゲノム中の特定の部位におけるマイナー対立遺伝子頻度(MAF)を推定するための二つの主要なアプローチについて述べた。 最初のアプローチは、個々の遺伝子型を推定し、MAFを推定する際にそれらの推定された遺伝子型を完全に正確であるとして扱うことを含む。 次に、遺伝子型を割り当てる際の不確実性を直接考慮する尤度フレームワークのパフォーマンスを調べます。 私たちの仕事を通して、私たちはすべての分離サイトが双方向であると仮定します。

呼び出された遺伝子型からMAFの推定

次世代の配列データからMAFを推定する一つの方法は、最初に配列データを使用して個々の遺伝子型を呼 これは、伝統的に遺伝子型データとサンガー配列決定データのために使用されるアプローチでした。 次世代シークエンシングデータに適用した場合、それがどれだけうまく機能するかは明らかではありません。最尤アプローチは、次世代配列決定データから各個体の遺伝子型を推測するために使用することができる。

最尤アプローチは、次世代配列決定データから各個体の遺伝子型を推測するために使用することができる。 各部位jにおいて、各個体iについて、3つの可能性のある遺伝子型のそれぞれの尤度(我々が小対立遺伝子を知っていると仮定して)は、次のように与え:

(1)

ここで、d i,jはサイトjにおける個々のiにおける観察された配列決定データであり、g i,j≤{0,1,2}は各個人の遺伝子型に含まれるマイナー対立遺伝子の数であり、シーケンスエラーの制御と読み取りベース品質をそれぞれ制御します。 各個体について観察された配列決定データは、読み取り品質スコアを考慮したサイトjでの読み取りの整列と考えることができる。 これは遺伝子型尤度として表され、SOAPsnpやMAQなどの次世代配列決定データを分析する多くのプログラムで生成される遺伝子型尤度ファイル(GLF)に見特定の個体に遺伝子型を割り当てるために、3つの可能性のある遺伝子型のそれぞれの尤度をその個体について計算することができる。

次に、最も高い可能性を有する遺伝子型を割り当てることができる。 しかし、研究者はしばしばより厳格な呼び出し基準を好み、最も可能性の高い遺伝子型が第二の最も可能性の高いものよりも実質的に可能性が高い ここで、3つの可能な遺伝子型は、その尤度によってソートされます:、ここでg(k)はk番目の最大尤度を持つ遺伝子型に対応します。 与えられたしきい値fで、の場合、遺伝子型g(1)を呼び出すことができます。 それ以外の場合、遺伝子型は呼び出されず、個人の遺伝子型は欠落しているとみなされます。 Fの共通した閾値は1であり、最も可能性の高い遺伝子型が第二の最も可能性の高い遺伝子型よりも少なくとも10倍可能性が高いことを示す。 このタイプのフィルタリングは、「呼び出された」遺伝子型の信頼性が高くなる可能性があるが、より多くの欠損データをもたらすことにも注意してく

対立遺伝子頻度の最尤推定器

呼び出された遺伝子型からMAFを推定する代わりに、Kimらによって導入された最尤(ML)法。 (同様のアプローチについてはLynchも参照)Mafを直接推定し、遺伝子型の不確実性を考慮に入れる。 具体的には、小対立遺伝子が与えられた場合、各個体iで配列データを観察する確率は、3つの可能な遺伝子型すべてに対応する確率を合計することに

式1で定義された三つの遺伝子型尤度が利用可能であると仮定します。 上記と同じ表記法を使用して、djおよびpjを、それぞれ、サイトjおよび対応するMAFでの観察された配列決定データとする。 マイナー対立遺伝子頻度がHardy-Weinberg平衡(HWE)を仮定することによって計算できることを与えられた遺伝子型確率。 次に、個人間の独立性を仮定すると、この軌跡におけるMAFの尤度は、すべてのN個の個人にわたって計算されるすべての尤度の積です。

(2)

PjのML推定値は、Broyden-Fletcher-Goldfarb-Shanno(Bfgs)を使用して制限されたパラメータ空間の尤度を直接最大化することによって計算することができます。方法または期待最大化(em)アルゴリズムを使用して。 EMアルゴリズムを使用する場合、遺伝子型の事後期待値が個々ごとに計算され、それらの事後の平均が繰り返し更新されます。 BFGSの実装はEMアルゴリズムよりも高速でした。 たとえば、100,000サイトから推定値を取得するには、BFGには約16秒かかりましたが、EMには約100秒かかりました。 ただし、速度の違いは実装固有のものである可能性があります。 私たちの場合、両方の方法で、尤度の増加が0.001未満のときにパラメータの更新を停止しました。

不確かなマイナー対立遺伝子を持つ最尤推定量

実際には、多くの場合、個人間で二番目に一般的なヌクレオチドは、マイナー対立遺伝子とし しかし、まれなSnp(例えば、MAF<1%)では、配列決定エラーのために4つのヌクレオチドすべてが一部の読み取りに現れる可能性があるため、どの対立遺伝子がマイナー対立遺伝子であるかを決定することは困難である。 この状況に対処するために、我々は今、アカウントにマイナー対立遺伝子の決定の不確実性を取る尤度の枠組みについて説明します。サイトjについて主要な対立遺伝子Mを知っていると仮定します。

サイトjについて主要な対立遺伝子Mを知っていると仮定します。 まれなSnpでの頻度の推定に主に関心があるため、2つの共通の対立遺伝子のどれが主要な対立遺伝子である可能性が高いかを決定することは重要ではないことに注意してください。 さらに、中間頻度(約50%)を有する対立遺伝子については、主要対立遺伝子とマイナー対立遺伝子との区別はあまり重要ではない。 他の三つの非主要ヌクレオチドm1、m2、およびm3を割り当てます。 したがって、マイナー対立遺伝子の指定の不確実性を可能にするために、尤度関数は次のように変更することができます。:

(3)

さらに、三つの可能なマイナー対立遺伝子のいずれかが等しく可能性が高いと仮定すると、次のようになります。

(4)

ここで、

ここで、

は大きなデータセット(例えば、多くの個人)では非常に小さくなる可能性があるため、対数スケールで尤度を計算すると便利です。 3つの条件付き対数尤度を(l(1),l(2),l(3))に順序付けします(l(1)は最大の対数尤度です)。 その後、

関連マッピングのための呼び出された遺伝子型を使用したGテスト

関連研究では、症例と対照の間の対立遺伝子頻度に有意な差を示すSnpは、関心の表現型と関連していると言われている。 関連マッピングは、次世代配列決定研究からのデータを使用して実行することができます。 我々は最初に個々の遺伝子型を呼び出す必要があるアプローチを議論し、呼び出された遺伝子型を使用して関連付けのテストを実行します。 このアプローチでは、遺伝子型は、最初に各個人のために呼び出されます。 遺伝子型はろ過するか、またはろ過しないことができます。 個体間の独立性とHWEを仮定すると、2×2分割表は、症例と対照の両方における主要対立遺伝子と副対立遺伝子の数を数えることによって構築することができます。 これは、独立性のためのよく知られた尤度比検定、G検定につながります:

(5)

ここで、O k,hはセルで観察される頻度であり、E k,hは対立遺伝子の頻度がケースとコントロールの間で同じである帰無仮説の下で期待される頻度である。 よく知られているピアソンのカイ二乗検定は、g検定と漸近的に等価です。 表が真の遺伝子型から生成された場合、G統計量は1自由度(σ2(1))のカイ二乗分布に漸近的に従います。 しかし、我々の研究では、我々は”と呼ばれる”遺伝子型を使用してG統計を構築するため、HWEは、ヘテロ接合体のオーバーとアンダーコールのために保持しない可能性が さらに、「観察された」遺伝子型の代わりに「呼ばれた」遺伝子型を数えることによって試験統計を構築することは、おそらく余分な変動性を導入する。 したがって、統計理論はもはや有効ではないかもしれません。 遺伝子型が特定の個体のために呼び出されない場合、データは欠落しているとみなされ、2×2表には含まれないことに留意されたい。

尤度比検定は、関連マッピングのために観察された遺伝子型の不確実性を考慮する

遺伝子型を呼び出す代わりに、尤度フレームワークは、遺伝子型の不確実性を可能にし、各部位jで対立遺伝子の頻度がケースとコントロールの間で同じであるかどうかを検定する。 形式的には、仮説H Oの尤度を計算します:pj,1=pj,2(=pj,0)およびH A : ここで、pj,1およびpj,2は、それぞれケースおよびコントロールにおけるMafである。マイナー(m)およびメジャー(M)対立遺伝子が既知であると仮定すると、マイナー対立遺伝子頻度の尤度は式2に記載されているように計算することができ、尤度比検定統計量は次のように計算されます。

:p>

(6)

ここで、は、それぞれケースとコントロールのmafのmleです。小対立遺伝子が不明な場合、帰無仮説の尤度は式3のように計算され、LRT統計量は次のように修正されます:P>

(7)

ここで、djはケースとコントロールの両方の観測データであり、は帰無仮説の下での対立遺伝子頻度です。 その他の表記は、式6と同じである。

シミュレートされたデータでMAFを推定する

我々は、真の遺伝子型(True)、任意のフィルタリングなしの遺伝子型(Nfを呼び出す)と呼ばれる、フィルタリング(f=1; 稀なSnpの場合、小対立遺伝子型はしばしば明らかではない。 遺伝子型を呼び出すとき、第二の最も一般的なヌクレオチドは、マイナー対立遺伝子であると仮定されます。 ML法は、マイナー対立遺伝子を決定する際の不確実性を直接組み込み、特に明記しない限り、未知のマイナー対立遺伝子法(式3)を用いた結果が示される。 未知の小対立遺伝子ML法は、既知の小対立遺伝子ML法と同様に機能するが、非常にまれなSnp(追加ファイル1)に対しては前者の方が優れていることに注

我々は最初に、異なるアプローチが200人の個体において、真のMAFが5%の1,000Snpのシーケンス深度の範囲にわたってMAFを推定することができたかを評価した。 図1は、4つの異なるアプローチを使用して推定されたMafの分布の箱ひげ図を示しています。 予想されるように、1 2×の個々の深さなどのより高いカバレッジデータについては、遺伝子型が確実に知られている場合と同様に、すべての方法が実行さ しかし、深さが減少すると、最初に遺伝子型を呼び出すことによって得られたMAFの推定値は偏ったものになる。 たとえば、Call F法を使用して推定されたMAFの中央値は、6倍のカバレッジで5.3%、2倍で12.5%です。 上方バイアスの理由は、真のヘテロ接合体はしばしば配列決定エラーのように見えるので、ヘテロ接合体を呼び出すことが困難になるからである。 したがって、マイナーなホモ接合体よりも多くのヘテロ接合体は、遺伝子型が欠落している傾向があります。 しかしながら、呼び出された遺伝子型からのMAF推定の全体的な偏りは、常に一方向にあるとは限らない(データは示されていない)。 興味深いことに、バイアスはCall Fメソッドの方がCALL NFメソッドよりも悪いように見えます。 このパターンは、遺伝子型呼び出しをフィルタリングすると、配列決定エラーをヘテロ接合体と呼ぶ確率が低下するように見えるため、直感的ではないように見えるかもしれません。 しかし、Call Fメソッドは、主要な対立遺伝子の多くのホモ接合体が配列決定エラーのために呼び出されないため、大量の欠損データをもたらす。 したがって、この例では、フィルタリングせずに遺伝子型を呼び出すことは、MAFを推定しようとするときに遺伝子型をフィルタリングするよりも

図1
図1

カバレッジの異なる深さのための5%の真のMAFを持つサイトでの対立遺伝子頻度の推定。 各深さで、1,000サイトは、200人を使用してシミュレートされ、各サイトで、対立遺伝子頻度の推定値は、(1)真の遺伝子型(True)を使用して計算され、(2)フィルタリング 推定方法の詳細については、”方法”を参照してください。新しいMLメソッドの結果は劇的に異なります。

新しいMLメソッドの結果は劇的に異なります。 この方法は、深さの範囲にわたってMAF(-4.9%の中央値)の公平な推定値を提供します。 2倍であっても、推定値は真の遺伝子型に基づくものよりもわずかに大きな分散しか示さない。また、シーケンス深度の範囲にわたってMAFの異なる推定値の推定平均二乗誤差(MSE;Expectation()を比較しました(図2)。 MLメソッドのMSEは、50人または200人の呼び出し元メソッドよりも低くなります。 特に、Call Fメソッドに基づいて計算されたMSEは、特に深度が減少する場合に、他のメソッドからのMSEよりもはるかに高くなります。 真の遺伝子型に基づくMAFの推定値のMSEはMSEの下限を反映し、サンプリング分散と有限サンプルサイズのために深さにわたって一定ではありません。 50人を使用すると、MSEは深さの増加に伴って0.0005に近づき、200人のサンプルサイズを使用すると、深さの増加に伴って0.0013に近づきます。

図2
図2
図2

異なるサンプルサイズの対立遺伝子頻度推定量の四つの異なるタイプの平均スクレードエラー(MSE;期待)左右のパネル)および適用範囲(x軸)の深さ。 各深さにおいて、MSEは、4つの異なる方法:真、呼気NF、呼気F、およびML(方法の詳細については、図1のキャプションを参照)を使用して行われた対立遺伝子頻度推定値から計算された。P>

全体的に、新しいMLメソッドは、遺伝子型呼び出しメソッドをアウト実行します。

シミュレートされたデータからMafの分布を推定する

次に、母集団の異なる周波数でのSnpの割合を推定する際に、異なる推定アプローチがどのよ ここでは、真のMafの分布が10,000の有効な母集団サイズの標準定常分布に従った20,000Snpをシミュレートしました(方法を参照)。 しかしながら、実際には、非常にまれなSNPを配列決定誤差と区別することは非常に困難であることに留意されたい。 したがって、実際のデータとの比較目的のために、推定MAFが2%未満のSnpを破棄しました。 図3は、推定MAF<2%を持つSnpを除外した後、各異なる周波数ビンに落ちるSnpの割合を示しています。

図3
図3

対立遺伝子頻度の標準的な定常分布を仮定してシミュレートされたSnpの対立遺伝子頻度の分布。 各深さ(各パネル)で、20,000Snpをシミュレートし、各SNPについて、MAFの推定値を四つの異なる方法を用いて得た(図1のキャプションを参照)。 次に、各方法(各色)について、推定された対立遺伝子頻度>2%を有する部位のみを使用して、各ヒストグラム(x軸)を生成する。

予想されるように、個々の10倍などのカバレッジの深さが高いため、すべてのメソッドは、真の遺伝子型に基づいて期待される分布に類似した推定MAF分布を提供します(図3)。 個体あたり4倍未満のようなカバレッジの浅い深さでは、遺伝子型呼び出し方法によって得られたMafの分布は、真の遺伝子型に基づいて予想されるMAF分布から有意に逸脱する(図3)。 特に、これらの方法は、低周波Snpの割合を過剰推定する。 例えば、第2のビン中のSnpの予想される割合(2〜4%の間の推定MAF)は1 8%である。 4×の深さでのCall NFメソッドに基づく対応する割合は26%であり、これは予想よりも1.4倍高い。 低周波Snpの割合の過剰推定は、ヘテロ接合遺伝子型を過負荷になる真のヘテロ接合体との配列決定エラーの混乱のために発生します。 このインフレの大きさは、フィルタリングカットオフによって異なりますが、カットオフを大きくしても、必ずしもインフレを増減するとは限りません。ML法では画像がまったく異なります。

新しいML法から得られた推定MAF分布は,カバレッジの浅い深さでも真の分布に密接に従う。 ここでは、過剰な低周波Snpはほとんどありません。 深さが4×の場合、ヒストグラムの2番目のビン内のSnpの割合は18.4%であり、これは予想される割合(18%)に非常に近い。 したがって、周波数スペクトルのより信頼性の高い推定値は、遺伝子型呼び出しアプローチを使用するよりも、私たちの尤度アプローチを使用して低カバレージデータから行うことができます。

シミュレートされたデータにおける関連マッピング

私たちは、(G検定を使用して)関連のテストで真の遺伝子型として推定された遺伝子型を扱う方法のパフォーマンスを、遺伝子型の不確実性を考慮した尤度比検定(LRT)と比較します。 症例と対照の間に対立遺伝子頻度差がないという帰無仮説の下での検定統計量の分布を調べた。 また、さまざまなアプローチのパワーを比較します。

サンプルサイズがかなり大きい場合、標準漸近理論は、帰無仮説の下では、G統計量とLRT統計量の両方が1自由度(σ2(1))のカイ二乗分布に従うことを示 したがって、呼び出しメソッドとLRT統計に基づいて計算されたG統計量のnull分布を、QQプロットを使用したσ2(1)分布と比較しました(図4)。 我々は、遺伝子型をシミュレートするために使用されるMAFは、ケースとコントロールの両方で5%であった500例とコントロールのシーケンスの深さの様々な5,000Snpをシミュ 真の遺伝子型を使用して計算されたG統計量の分布は、σ2(1)分布と非常に良好な対応を示しています。 しかし、呼び出された遺伝子型に基づいて計算されたG統計量の分布は、実質的にσ2(1)分布から逸脱する。 遺伝子型を呼び出し、それらの遺伝子型を正確であると処理すると、p値がσ2(1)分布を使用して計算される場合、偽陽性信号が膨大に過剰に生成さ 例えば、2×の深さでは、Snpの1 1%は、予想される5%と比較して、5%未満のp値を有していた。 この効果は、会合を検出するためにここで使用される対立遺伝子試験において、ヘテロ接合体としてホモ接合体を過負荷にすることによる分散の増加 Hardy-Weinberg平衡からの逸脱に対して頑健であるArmitage trend検定のような遺伝子型検定は、偽陽性率の同様の増加を示さない(追加ファイル2)。 この観察と一致して、呼ばれた遺伝子型をフィルタリングすると、G検定を使用するときに有意な検定の割合が減少しますが、フィルタリングは問題を完全に解決するわけではありません。 一方、LRT統計は、カバレッジの2×または5×深さのいずれかのσ2(1)分布からの非常にわずかな逸脱のみを示しています。

図4
図4

QQ-関心のある検定統計量のヌル分布とσ2(1)分布を比較するプロット。 各列は、(1)真の遺伝子型を使用して計算されたG統計量(True)、(2)フィルタリングなしで呼び出された遺伝子型を使用して計算されたG統計量(NFの呼び出し); (3)フィルタリングを用いた呼出された遺伝子型を用いて計算されたG統計量(Call F);および(4)未知の小対立遺伝子(LRT)を用いた尤度比検定統計量。 500のケースと500のコントロールを仮定すると、帰無仮説の下で、5,000のサイトのセットは、5%のMAFで2×(上のパネル)と5×(下のパネル)の配列決定の深さでシミュ “インフレ”係数は、各図の左上隅に示されています。また、異なる関連テストごとに受信機動作特性(ROC)曲線も生成しました。 これらの曲線は、異なる偽陽性率での検定の検出力を示しています。 いくつかの検定統計量の分布は帰無仮説の下でσ2(1)分布に従わないので、公正な比較を行うために、経験的帰無分布に基づいて各偽陽性率の臨界値を得た。 検出力は、臨界値を超える統計量を有するシミュレートされた疾患遺伝子座の割合として計算される。 全体として、我々は、LRTは、いずれかの遺伝子型呼び出し方法に基づいてG-テストよりも優れて実行することがわかります(図5)。 例えば、5%の偽陽性率で、配列決定の深さが5×である場合、MAFが1%、相対リスク(rr)が2の疾患遺伝子座を検出するパワーは、LRTでは5 1%であるが、フィルタリング 特に、低深度では、フィルタリングを用いて呼び出された遺伝子型に適用されるG検定は非常に不十分に実行される(図5の最も左の列)。 呼び出された遺伝子型を使用してLrtの検出力をArmitage trend検定と比較すると、LRTもArmitage trend検定よりも高い検出力を持っていることがわかります(追加ファイル3)。 これは、呼び出された遺伝子型を使用したい場合、呼び出しの信頼性に基づいてそれらをフィルタリングすると、電力が失われる可能性があること

図5
図5

レシーバ動作特性(ROC)関連の四つのテストの曲線。 4つの統計の定義については、図4のキャプションを参照してください。 500ケースと500コントロールを仮定すると、20,000サイトのセットは、ヌルの下で、2×、5×、および10×(三列)の個々の配列決定深さで代替の下でシミュレートされました。 各偽陽性率(x軸)において、対応する臨界値を経験的ヌル分布を用いて計算した。 真の陽性率(パワー;y軸)は、臨界値を超えるテスト統計量を有する原因部位の割合を計算することによって得られた。

実際のデータへの適用

Illumina技術を使用して配列決定された疾患関連研究のための対照から200のエクソームを8×の個々の深さで分析しました。 “SOAPsnp”プログラムによって生成された遺伝子型尤度を推論に使用しました。 詳細については、”方法”を参照してください。

まず、Sequenom遺伝子型データから推定されたMafと比較することにより、50Snpの次世代シーケンシングデータからMAFの推定の精度を検討しました。 ML法を使用した推定値とフィルタリングなしの遺伝子型呼び出し法の両方は、シークノム遺伝子型データから行われた推定値(すなわち、図6の2つの推定値の間の小さな標準化された差)と高い相関を示しています。 しかし、フィルタリングを用いた遺伝子型呼び出しに基づく推定値は、特に配列の深さが低い場合、配列の遺伝子型データから推定された頻度への対応が悪いことを示している。 興味深いことに、配列決定の深さが非常に高い(14×)にもかかわらず、再配列データから推定されたMAFが配列遺伝子型データから得られた推定値とは非常に異 具体的には、Sequenom遺伝子型データから推定されたMAFは22.5%であるが、MLアプローチを用いて推定された場合は17.2%である。 個々の検査は多くの個人で、配列決定データに基づく非常に支えられた遺伝子型がSequenomの遺伝子型と異なることを示します。 このSNPがこれらの個体の多くの読み取りによってカバーされ、観察された読み取りベースが高品質スコアを有することを考えると(>Q20)、その差はシークエノムジェノタイピングエラーによるものである可能性が高い。 フィルタリングなしの遺伝子型呼び出しアプローチから推定されたMafは、MLアプローチからの推定よりもSequenom genotypingから推定されたMafによく対応するように見えるSnpのカップルがあることに注意してください。 例えば、1つのSNPでは、推定されるMAFは、シークエノム遺伝子型データから2 5.7%であり、2 5である。フィルタリングなしの遺伝子型呼び出し法から9%、およびML法から27.2%。 しかし、個々の検査は、配列決定データから呼び出された遺伝子型がSequenom遺伝子型とは異なるいくつかの個体があることを明らかにする。 これらのケースでは、呼び出された遺伝子型のエラーがキャンセルされ、Sequenom遺伝子型データとのより良い対応の外観が与えられる。 したがって、これらのSnpでは、どの方法が最適かを判断することは困難です。

図6
図6

次世代配列決定データとシークエノム遺伝子型データを用いて200人から計算された対立遺伝子頻度の推定値。 各部位において、配列遺伝子型データおよび配列決定データの両方を有する個体のみが、対立遺伝子頻度の推定のために使用された。 配列決定データについて、MAFの推定値は、3つの異なる方法(Call NF;Call F;およびML)を使用して得られた。 各推定値の標準化された差は、は、それぞれ配列決定データおよびSequenom遺伝子型データからの推定Mafであり、nは推定に使用される個体数である。 各サイトは、平均的な個々のカバレッジの深さ(色)に基づいて、4×未満、4×より高いが8×未満、8×より高いが16×未満、および16×より高い四つのビンのいず

次に、次世代のエクソームシーケンスデータからシーケンス深さの範囲にわたっていくつかのアプローチを使用して計算されたMafの分布を調べました(図7)。 このデータセットのシーケンスエラーからこれらの非常に低周波Snpを区別することは困難であるため、推定MAF<2%を持つSnpを破棄しました。 我々はさらに、マイナーとメジャー対立遺伝子の間の読み取り塩基の品質スコアに有意な差(p値未満10-5ランク合計テストを使用して)があったサイトを取 これらの部位は、誤ったマッピングまたは実験手順中に導入された未知のバイアスのために発生する可能性のある人工Snpである可能性が高い。 次に、各サイトをカバレッジの深さに基づいてビンに分類しました。 各ビン内のSnpの数を表1に示す。 平均深度が9×未満の場合、遺伝子型呼び出し法に基づく推定Mafの分布は、ML法に基づくものとは非常に異なる。 具体的には、遺伝子型呼出アプローチは、低周波Snp(2%〜4%のMAF)の大過剰を生じさせる。 このパターンは、シミュレーション研究で見られたものを反映しています(図3)。 また、遺伝子型呼び出し方法では、配列決定の深さが変化するにつれて対立遺伝子の頻度分布が劇的に変化する。 したがって、前述したように、深さがあまり高くない場合、遺伝子型決定呼び出しメソッドは、配列決定エラーである多くの偽のSnpを含む可能性が高い。 これらの誤差は、周波数分布に過剰な低周波Snpとして表示されます。 ML法に基づく分布は、深さにわたってより安定であるが、深さが9×未満の低い対立遺伝子頻度を有するSnpの過剰は、より深い深さでの低周波Snpの比

図7
図7

配列された200人の個体のエクソームから推定されたマイナー対立遺伝子頻度の分布。 各サイトについて、マイナー対立遺伝子の頻度は、四つの異なる方法を使用して推定されました: (1)未知のマイナー対立遺伝子を持つML法、(2)既知または固定のマイナー対立遺伝子を持つML法、(3)フィルタリングせずに遺伝子型を呼び出す(Nfを呼び出す)、および(4)フィルタリングリングを持つ遺伝子型を呼び出す(Fを呼び出す)。 各サイトは、カバレッジの深さに基づいてビンに分類されます。 さらに、各ヒストグラムでは、推定MAFが2%未満のサイトは考慮されません。 この分析に使用されたSnpの数については、表1を参照してください。

表1個人間の平均シーケンス深度によって定義された各ビン(列)内の特定の方法(行)を使用して、推定MAFが2%を超えるSnpの数。

最後に、このエクソーム再シーケンスデータを使用して、ケースコントロール関連研究をシミュレートしました。 帰無仮説の下での関連検定統計量の分布を調べるために、100人をケースグループに無作為に割り当て、他の100人を対照グループに割り当てました。 MAF推定>2%(未知のマイナー対立遺伝子ML法に基づいて)と染色体2上のすべてのSnpのために、我々は、フィルタリングの有無だけでなく、LRT統計の両方と呼ばれる遺伝子型を使用してg統計量を計算することにより、ケースとコントロールの間の対立遺伝子頻度差をテストした。 図8は、検定統計量の分布を標準のσ2(1)分布と比較したQQプロットを示しています。 シミュレーション研究に見られるように、実質的にフィルタリングせずに遺伝子型を呼び出すときに計算されたG統計のヌル分布は、σ2(1)分布から ただし、LRT統計量のヌル分布は、σ2(1)分布に密接に従います。 インフレーション係数は1.01であり、LRT統計が実際のデータに適用された場合に良好に機能することを意味します。

図8
図8

QQ-100ケースと100コントロール間の対立遺伝子頻度差の関連テスト統計量をσ2(1)分布に比較するプロット。 表現型は、100のケースと100のコントロールがあるように、エクソームリシーケンスデータセット内のindivdidualsにランダムに割り当てられました。 各サイトについて、三つの統計が計算されました: フィルタリングなしで呼び出された遺伝子型を使用するG統計量(Nfの呼び出し)、フィルタリング付きで呼び出された遺伝子型を使用するG静的統計量(Fの呼び出し)、およびLRT統計量。 偽Snpの包含を最小化するために、2%未満のML MAF推定値を有する部位を廃棄する。 表示目的のために、染色体2上のサイトからの結果が示されています。 インフレーション係数は、各QQプロットの左上隅に表示されていることに注意してください。

Related Posts

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です