Das kleinere Allel ist das weniger häufige Allel in der Population an einem variablen Ort. Wir beschreiben zunächst zwei Hauptansätze zur Schätzung der geringen Allelfrequenz (MAF) an einer bestimmten Stelle im Genom. Der erste Ansatz besteht darin, einzelne Genotypen abzuleiten und diese abgeleiteten Genotypen bei der Schätzung der MAF als vollständig genau zu behandeln. Anschließend untersuchen wir die Leistung eines Likelihood-Frameworks, das die Unsicherheit bei der Zuordnung von Genotypen direkt berücksichtigt. Während unserer gesamten Arbeit, Wir gehen davon aus, dass alle Trennstellen biallelisch sind.

Schätzung der MAF aus Genotypen

Eine Möglichkeit, die MAF aus Next-Generation-Sequenzierungsdaten zu schätzen, besteht darin, zuerst einen Genotyp für jedes Individuum unter Verwendung von Sequenzierungsdaten aufzurufen und diese Genotypen dann so zu verwenden, als wären sie die wahren. Dies war der Ansatz, der traditionell für Genotypdaten und Sanger-Sequenzierungsdaten verwendet wurde. Es ist nicht klar, wie gut es funktionieren wird, wenn es auf Sequenzierungsdaten der nächsten Generation angewendet wird.

Ein Maximum-Likelihood-Ansatz kann verwendet werden, um den Genotyp für jedes Individuum aus den Next-Generation-Sequenzierungsdaten abzuleiten. An jedem Ort j ist für jedes Individuum i die Wahrscheinlichkeit für jeden der drei möglichen Genotypen (vorausgesetzt, wir kennen das kleinere Allel) gegeben als:

wobei D i,j die beobachteten Sequenzierungsdaten in Individuum i an der Stelle j sind, g i , j ∈ {0, 1, 2} die Anzahl der im Genotyp jedes Individuums enthaltenen Nebenallele ist und und Steuerung für Sequenzierungsfehler bzw. Die beobachteten Sequenzierungsdaten für jedes Individuum können als die Ausrichtung der Lesevorgänge an der Stelle j unter Berücksichtigung der Lesequalitätswerte angesehen werden. Dies wird als Genotyp-Likelihood dargestellt und ist in der Genotyp-Likelihood-Datei (GLF) enthalten, die in vielen Programmen zur Analyse von Sequenzierungsdaten der nächsten Generation wie SOAPsnp und MAQ erstellt wird .

Um einem bestimmten Individuum einen Genotyp zuzuordnen, kann die Wahrscheinlichkeit jedes der drei möglichen Genotypen für das Individuum berechnet werden. Der Genotyp mit der höchsten Wahrscheinlichkeit kann dann zugeordnet werden. Forscher bevorzugen jedoch häufig ein strengeres Berufungskriterium und weisen einem Individuum keinen Genotyp zu, es sei denn, der wahrscheinlichste Genotyp ist wesentlich wahrscheinlicher als der zweitwahrscheinlichste. Hier werden die drei möglichen Genotypen nach ihren Wahrscheinlichkeiten sortiert: , wobei g(k)dem Genotyp mit der k-ten größten Wahrscheinlichkeit entspricht. Bei einem gegebenen Schwellenwert f kann man den Genotyp g(1) aufrufen, wenn . Andernfalls wird ein Genotyp nicht aufgerufen und der Genotyp des Individuums wird als fehlt betrachtet. Ein gemeinsamer Schwellenwert von f ist 1, was darauf hinweist, dass der wahrscheinlichste Genotyp mindestens 10-mal wahrscheinlicher ist als der zweitwahrscheinlichste. Beachten Sie, dass diese Art der Filterung zu einem höheren Vertrauen für den „aufgerufenen“ Genotyp führen kann, aber auch zu mehr fehlenden Daten führt.

Maximum likelihood estimator of allele frequency

Anstatt die MAF aus den genannten Genotypen zu schätzen, wurde eine von Kim et al. (siehe auch Lynch für einen ähnlichen Ansatz) schätzt MAFs direkt und berücksichtigt die Genotypunsicherheit. Insbesondere wird bei einem kleinen Allel die Wahrscheinlichkeit, die Sequenzdaten an jedem einzelnen i zu beobachten, durch Summieren der Wahrscheinlichkeiten erhalten, die allen drei möglichen Genotypen entsprechen.

Angenommen, die drei in Gleichung 1 definierten Genotyp-Wahrscheinlichkeiten sind verfügbar. Unter Verwendung derselben Notation wie oben seien D j und p j die beobachteten Sequenzierungsdaten an der Stelle j bzw. die entsprechende MAF. Die Genotypwahrscheinlichkeit angesichts dieser geringen Allelhäufigkeit kann unter der Annahme des Hardy-Weinberg-Gleichgewichts (HWE) berechnet werden. Dann, unter der Annahme der Unabhängigkeit unter Individuen, ist die Wahrscheinlichkeit des MAF an diesem Ort ein Produkt aller Wahrscheinlichkeiten, die über alle N Individuen berechnet werden:

Die ML-Schätzung von pj kann entweder durch direkte Maximierung der Wahrscheinlichkeit für einen begrenzten Parameterraum unter Verwendung des Broyden-Fletcher-Goldfarb-Shanno BFGS) Methode oder unter Verwendung des Erwartungsmaximierungsalgorithmus (EM) . Bei Verwendung des EM-Algorithmus wird die posteriore Erwartung eines Genotyps für jedes Individuum berechnet, und der Mittelwert dieser Posterioren wird wiederholt aktualisiert. Unsere Implementierung von BFGS war schneller als der EM-Algorithmus. Um beispielsweise Schätzungen von 100.000 Standorten zu erhalten, dauerte BFGS ~ 16 Sekunden, EM jedoch ~ 100 Sekunden. Der Geschwindigkeitsunterschied kann jedoch implementierungsspezifisch sein. In unserem Fall haben wir für beide Methoden die Aktualisierung der Parameter gestoppt, als die Erhöhung der Wahrscheinlichkeit weniger als 0,001 betrug.

Maximum-Likelihood-Schätzer mit unsicherem Minor-Allel

In der Praxis kann oft das zweithäufigste Nukleotid über Individuen hinweg als Minor-Allel verwendet werden. Bei seltenen SNPs (z. B. MAF < 1%) ist es jedoch schwierig zu bestimmen, welches Allel das kleinere Allel ist, da alle vier Nukleotide aufgrund von Sequenzierungsfehlern in einigen Lesevorgängen auftreten können. Um dieser Situation zu begegnen, beschreiben wir nun einen Likelihood-Rahmen, der die Unsicherheit bei der Bestimmung des Minor-Allels berücksichtigt.

Angenommen, wir kennen für den Ort j das Hauptallel M. Beachten Sie, dass die Entscheidung, welches von zwei gemeinsamen Allelen wahrscheinlich das wichtigste ist, nicht wichtig ist, da wir uns hauptsächlich mit der Schätzung der Frequenzen bei seltenen SNPs befassen. Für Allele mit Zwischenfrequenzen (etwa 50%) ist die Unterscheidung zwischen Haupt- und Nebenallel weniger wichtig. Weisen Sie die anderen drei Nicht-Hauptnukleotide m1, m2 und m3 zu. Die in Gleichung 2 eingeführte Wahrscheinlichkeit geht von einem festen Hauptallel M und einem festen Nebenallel m aus. Um Unsicherheit bei der Bezeichnung des Nebenallels zu berücksichtigen, kann die Wahrscheinlichkeitsfunktion daher wie folgt modifiziert werden:

Unter der Annahme, dass eines der drei möglichen Nebenallele gleich wahrscheinlich ist, erhalten wir:

wobei. Da bei großen Datensätzen (z. B. bei vielen Personen) sehr klein sein kann, ist es nützlich, die Wahrscheinlichkeit in der Log-Skala zu berechnen. Ordnen Sie die drei bedingten Log-Wahrscheinlichkeiten nach (l (1), l (2), l(3)), wobei l (1) die größte ist. Dann

G-Test unter Verwendung genannter Genotypen zur Assoziationszuordnung

In Assoziationsstudien sollen SNPs, die signifikante Unterschiede in der Allelhäufigkeit zwischen Fällen und Kontrollen zeigen, mit dem interessierenden Phänotyp assoziiert sein. Assoziationsmapping kann unter Verwendung von Daten aus Next-Generation-Sequenzierungsstudien durchgeführt werden. Wir diskutieren zunächst Ansätze, die den Aufruf einzelner Genotypen erfordern, und führen dann einen Test zur Assoziation mit den genannten Genotypen durch. Bei diesem Ansatz wird zunächst für jedes Individuum ein Genotyp aufgerufen. Die Genotypen können gefiltert oder ungefiltert sein. Unter der Annahme der Unabhängigkeit zwischen Individuen und HWE kann eine 2 × 2-Kontingenztabelle erstellt werden, indem die Anzahl der Haupt- und Nebenallele sowohl in den Fällen als auch in den Kontrollen gezählt wird. Dies führt zu dem bekannten Likelihood-Ratio-Test für die Unabhängigkeit, dem G-Test:

wobei O k, h die in einer Zelle beobachtete Frequenz und E k,h die unter der Nullhypothese erwartete Frequenz ist, bei der die Allelfrequenz zwischen Fällen und Kontrollen gleich ist. Der bekannte Pearson-Chi-Quadrat-Test ist asymptotisch äquivalent zum G-Test. Wenn die Tabelle aus echten Genotypen generiert wird, folgt die G-Statistik asymptotisch einer Chi-Quadrat-Verteilung mit 1 Freiheitsgrad (χ2(1)). In unseren Studien konstruieren wir jedoch die G-Statistik unter Verwendung von „genannten“ Genotypen, so dass HWE aufgrund von Über- und Unterbezeichnung von Heterozygoten nicht halten kann. Darüber hinaus führt die Konstruktion der Teststatistik durch Zählen von „aufgerufenen“ Genotypen anstelle von „beobachteten“ Genotypen wahrscheinlich zu einer zusätzlichen Variabilität. Daher ist die statistische Theorie möglicherweise nicht mehr gültig. Beachten Sie, dass, wenn ein Genotyp für eine bestimmte Person nicht aufgerufen wird, die Daten als fehlen gelten und nicht in der 2 × 2-Tabelle enthalten sind.

Likelihood Ratio-Test unter Berücksichtigung der Unsicherheit in den beobachteten Genotypen für die Assoziationszuordnung

Anstatt Genotypen aufzurufen, ermöglicht das Likelihood Framework die Unsicherheit in den Genotypen und testet an jeder Stelle, ob die Allelhäufigkeit zwischen Fällen und Kontrollen gleich ist. Formal berechnen wir die Wahrscheinlichkeit der Hypothesen H O : p j ,1 = p j ,2(= p j ,0) und H A : p j ,1 ≠ p j ,2, wobei p j,1 und p j ,2 die MAFs in Fällen bzw.

Unter der Annahme, dass Minor (m) – und Major (M) -Allele bekannt sind, kann die Wahrscheinlichkeit der Minor-Allelhäufigkeit wie in Gleichung 2 beschrieben berechnet werden, und die Likelihood Ratio-Teststatistik wird berechnet als:

wobei und die beobachteten Daten für Fälle bzw. Kontrollen sind und und sind die mlEs der MAFs in Fällen bzw.

Wenn das Minor-Allel unbekannt ist, wird die Wahrscheinlichkeit unter der Nullhypothese wie in Gleichung 3 berechnet, und die LRT-Statistik wird wie folgt modifiziert:

wobei DJ die beobachteten Daten für beide Fälle und Kontrollen und die Allelhäufigkeit unter der Nullhypothese ist. Andere Notationen sind die gleichen wie in Gleichung 6.

Schätzung der MAF in simulierten Daten

Wir vergleichen die Schätzungen der Allelhäufigkeit in simulierten Daten mit echten Genotypen (True), Genotypen ohne Filterung (Call NF), Genotypen mit Filterung (f = 1; Call F) und der Maximum-Likelihood-Methode (ML). Bei seltenen SNPs ist der Minor-Alleltyp oft nicht ersichtlich. Beim Aufruf von Genotypen wird angenommen, dass das zweithäufigste Nukleotid das kleinere Allel ist. Die ML-Methode beinhaltet direkt Unsicherheit bei der Bestimmung des Minor-Allels, und sofern nicht anders angegeben, werden Ergebnisse unter Verwendung der unbekannten Minor-Allel-Methode (Gleichung 3) gezeigt. Beachten Sie, dass die unbekannte Minor-Allel-ML-Methode ähnlich wie die bekannte Minor-Allel-ML-Methode funktioniert, die erstere jedoch für sehr seltene SNPs besser ist (Zusätzliche Datei 1).

Wir haben zuerst ausgewertet, wie gut die verschiedenen Ansätze in der Lage waren, den MAF bei 200 Individuen über einen Bereich von Sequenzierungstiefen für 1.000 SNPs mit einem echten MAF von 5% abzuschätzen. Abbildung 1 zeigt Boxplots der Verteilungen der geschätzten MAFs unter Verwendung der vier verschiedenen Ansätze. Wie erwartet, für höhere Abdeckungsdaten, wie eine individuelle Tiefe von 12 ×, Alle Methoden funktionieren so gut wie wenn die Genotypen mit Sicherheit bekannt sind (Wahr). Wenn jedoch die Tiefe abnimmt, werden die Schätzungen der MAF, die durch das erste Aufrufen von Genotypen erhalten werden, voreingenommen. Beispielsweise beträgt die mediane MAF, die unter Verwendung der Call F-Methode geschätzt wird, 5,3% bei 6 × Abdeckung und 12,5% bei 2 ×. Der Grund für die Tendenz nach oben ist, dass es schwieriger wird, Heterozygoten zu nennen, da echte Heterozygoten oft wie Sequenzierungsfehler aussehen. Daher neigen mehr Heterozygoten als kleinere Homozygoten dazu, fehlende Genotypen zu haben. Die Gesamtverzerrung in MAF-Schätzungen von Genotypen ist jedoch nicht immer in eine Richtung (Daten nicht gezeigt). Interessanterweise scheint die Verzerrung für die Call F-Methode schlechter zu sein als für die Call NF-Methode. Dieses Muster mag kontraintuitiv erscheinen, da das Filtern der Genotypaufrufe die Wahrscheinlichkeit zu verringern scheint, einen Sequenzierungsfehler als Heterozygot zu bezeichnen. Die Call F-Methode führt jedoch auch zu einer größeren Menge fehlender Daten, da viele Homozygoten für das Hauptallel aufgrund von Sequenzierungsfehlern nicht aufgerufen werden. Daher scheint in diesem Fall das Aufrufen von Genotypen ohne Filterung die bessere Strategie zu sein als das Filtern von Genotypen, wenn versucht wird, die MAF abzuschätzen.

Die Ergebnisse unterscheiden sich dramatisch für die neue ML-Methode. Diese Methode liefert unvoreingenommene Schätzungen der MAF (Median von ~ 4,9%) über einen Bereich von Tiefen. Selbst bei 2 × zeigen die Schätzungen nur eine geringfügig größere Varianz als diejenigen, die auf den wahren Genotypen basieren.

Wir verglichen auch den geschätzten mittleren quadratischen Fehler (MSE; estimated mean squared error) () der verschiedenen Schätzungen des MAF über einen Bereich von Sequenzierungstiefen (Abbildung 2). Die ML-Methode hat eine niedrigere MSE als die aufrufenden Methoden mit 50 oder 200 Personen. Insbesondere ist die MSE, die basierend auf der Call F-Methode berechnet wird, viel höher als die der anderen Methoden, insbesondere wenn die Tiefe abnimmt. Die MSE der Schätzungen der MAF basierend auf den wahren Genotypen spiegelt die untere Grenze der MSE wider und ist aufgrund der Stichprobenvarianz und einer endlichen Stichprobengröße nicht über die Tiefe konstant. Bei Verwendung von 50 Personen nähert sich die MSE mit zunehmender Tiefe 0,0005 und bei Verwendung einer Stichprobengröße von 200 Personen nähert sie sich mit zunehmender Tiefe 0,0013.

Insgesamt übertrifft die neue ML-Methode Genotypaufrufmethoden.

Schätzen einer Verteilung von MAFs aus simulierten Daten

Als nächstes untersuchen wir, wie die verschiedenen Schätzansätze bei der Schätzung des Anteils von SNPs bei verschiedenen Frequenzen in der Population (ähnlich dem Standortfrequenzspektrum, jedoch basierend auf der Populationsallelfrequenz anstelle der Stichprobenfrequenz). Hier simulierten wir 20.000 SNPs, wobei die Verteilung der wahren MAFs der stationären Standardverteilung für eine effektive Populationsgröße von 10.000 folgte (siehe Methoden). Beachten Sie, dass es in der Praxis jedoch sehr schwierig ist, einen sehr seltenen SNP von einem Sequenzierungsfehler zu unterscheiden. Daher haben wir zum Vergleich mit realen Daten SNPs mit einem geschätzten MAF von weniger als 2% verworfen. Abbildung 3 zeigt den Anteil der SNPs, die in die verschiedenen Frequenzbereiche fallen, nachdem diese SNPs mit geschätzten MAF<2% ausgeschlossen wurden.

Wie erwartet liefern alle Methoden mit einer hohen Erfassungstiefe von 10 × pro Individuum geschätzte MAF-Verteilungen, die der erwarteten Verteilung auf der Grundlage der wahren Genotypen ähnlich sind (Abbildung 3). Bei einer geringeren Erfassungstiefe, wie z. B. weniger als 4 × pro Individuum, weichen die durch Genotypisierungsmethoden erhaltenen Verteilungen der MAFs signifikant von der erwarteten MAF-Verteilung ab, die auf echten Genotypen basiert (Abbildung 3). Insbesondere überschätzen diese Methoden den Anteil niederfrequenter SNPs. Beispielsweise beträgt der erwartete Anteil an SNPs im zweiten Bin (geschätzter MAF zwischen 2-4%) 18%. Der entsprechende Anteil basierend auf der Call-NF-Methode in einer Tiefe von 4× beträgt 26%, was 1,4-fach höher ist als erwartet. Die Überschätzung des Anteils an niederfrequenten SNPs erfolgt aufgrund von Verwechslungen von Sequenzierungsfehlern mit echten Heterozygoten, was zu einer Überbewertung heterozygoter Genotypen führt. Die Höhe dieser Inflation unterscheidet sich über verschiedene Filtergrenzen hinweg, aber ein größerer Grenzwert erhöht oder verringert nicht unbedingt die Inflation.

Das Bild ist für die ML-Methode völlig anders. Die geschätzte MAF-Verteilung, die mit der neuen ML-Methode erhalten wird, folgt der wahren Verteilung auch bei geringen Bedeckungstiefen genau. Hier gibt es fast keinen Überschuss an niederfrequenten SNPs. In einer Tiefe von 4× beträgt der Anteil der SNPs im zweiten Teil des Histogramms 18,4%, was dem erwarteten Anteil (18%) sehr nahe kommt. Somit, Zuverlässigere Schätzungen des Frequenzspektrums können aus Daten mit geringer Abdeckung unter Verwendung unseres Likelihood-Ansatzes vorgenommen werden als unter Verwendung der Genotypaufrufansätze.

Assoziationsmapping in simulierten Daten

Wir vergleichen die Leistung von Methoden, die abgeleitete Genotypen als echte Genotypen in Assoziationstests (unter Verwendung eines G-Tests) behandeln, mit unserem Likelihood Ratio Test (LRT), der die Unsicherheit in den Genotypen berücksichtigt. Wir untersuchen die Verteilung der Teststatistik unter der Nullhypothese, dass keine Allelhäufigkeitsdifferenz zwischen Fällen und Kontrollen besteht. Wir vergleichen auch die Kraft der verschiedenen Ansätze.

Bei relativ großen Stichprobengrößen legt die asymptotische Standardtheorie nahe, dass unter der Nullhypothese sowohl die G-Statistik als auch die LRT-Statistik einer Chi-Quadrat-Verteilung mit einem Freiheitsgrad folgen (χ2 (1)). Daher haben wir die Nullverteilung der G-Statistik, die basierend auf aufrufenden Methoden berechnet wurde, sowie die LRT-Statistik mit der χ2 (1) -Verteilung unter Verwendung von QQ-Plots verglichen (Abbildung 4). Wir simulierten 5.000 SNPs über eine Vielzahl von Sequenzierungstiefen in 500 Fällen und Kontrollen, wobei die zur Simulation von Genotypen verwendete MAF in beiden Fällen und Kontrollen 5% betrug. Die unter Verwendung der wahren Genotypen berechnete Verteilung der G-Statistik zeigt eine sehr gute Übereinstimmung mit einer χ2(1)-Verteilung. Die Verteilung der G-Statistik, die auf der Grundlage der genannten Genotypen berechnet wird, weicht jedoch wesentlich von einer χ2 (1) -Verteilung ab. Das Aufrufen von Genotypen und das anschließende Behandeln dieser Genotypen als genau erzeugt einen großen Überschuss an falsch positiven Signalen, wenn die p-Werte unter Verwendung einer χ2 (1) -Verteilung berechnet werden. Zum Beispiel hatten 11% der SNPs in einer Tiefe von 2× einen p-Wert von weniger als 5%, verglichen mit den erwarteten 5%. Der Effekt wird durch eine Erhöhung der Varianz verursacht, die durch eine Überbezeichnung von Homozygoten als Heterozygoten im hier zum Nachweis der Assoziation verwendeten Alleltest verursacht wird. Genotypische Tests wie der Armitage-Trendtest, die robust gegenüber Abweichungen vom Hardy-Weinberg-Gleichgewicht sind, zeigen keinen ähnlichen Anstieg der Falsch-Positiv-Rate (Zusätzliche Datei 2). In Übereinstimmung mit dieser Beobachtung führt das Filtern der genannten Genotypen zu einer Verringerung des Anteils signifikanter Tests bei Verwendung des G-Tests, obwohl das Filtern das Problem nicht vollständig löst. Auf der anderen Seite zeigt die LRT-Statistik nur eine sehr geringe Abweichung von einer χ2 (1) -Verteilung für 2 ×- oder 5 × -Deckungstiefen.

Wir haben auch ROC-Kurven (Receiver Operating Characteristic) für jeden der verschiedenen Assoziationstests generiert. Diese Kurven zeigen die Leistung des Tests bei unterschiedlichen falsch-positiven Raten. Da die Verteilungen einiger Teststatistiken nicht der χ2 (1) -Verteilung unter der Nullhypothese folgen, haben wir für einen fairen Vergleich den kritischen Wert für jede falsch positive Rate basierend auf der empirischen Nullverteilung erhalten. Die Potenz wird als Anteil simulierter Krankheitsorte berechnet, deren Statistik den kritischen Wert überschreitet. Insgesamt stellen wir fest, dass der LRT eine bessere Leistung erbringt als der G-Test, basierend auf beiden Genotypen oder Methoden (Abbildung 5). Beispielsweise beträgt bei einer falsch-positiven Rate von 5% und einer Sequenzierungstiefe von 5× die Leistung zum Erkennen eines Krankheitsorts mit einer MAF von 1% und einem relativen Risiko (RR) von 2 51% mit der LRT, die Leistung sinkt jedoch auf 33% unter Verwendung der aufrufenden Methode ohne Filterung und auf 34% unter Verwendung der aufrufenden Methode mit Filterung. Insbesondere in geringer Tiefe schneidet der G-Test, der auf Genotypen mit Filterung angewendet wird, sehr schlecht ab (Spalte ganz links in Abbildung 5). Wenn wir die Leistung des LRT mit dem Armitage-Trendtest unter Verwendung sogenannter Genotypen vergleichen, stellen wir fest, dass der LRT auch eine höhere Leistung als der Armitage-Trendtest aufweist (Zusätzliche Datei 3). Dies deutet darauf hin, dass, wenn man Genotypen verwenden möchte, das Filtern auf der Grundlage des Anrufvertrauens zu einem Leistungsverlust führen kann.

Anwendung auf reale Daten

Wir analysierten 200 Exome von Kontrollen für eine Krankheitsassoziationsstudie, die mit Illumina-Technologie in einer Tiefe von 8 × pro Individuum sequenziert wurden . Wir haben die vom „SOAPsnp“ -Programm generierten Genotyp-Wahrscheinlichkeiten für unsere Inferenz verwendet. Weitere Informationen finden Sie unter Methoden.

Zunächst untersuchten wir die Genauigkeit der Schätzungen des MAF aus Sequenzierungsdaten der nächsten Generation für 50 SNPs, indem wir sie mit den geschätzten MAF aus Sequenzgenotypdaten verglichen. Sowohl die Schätzungen mit der ML-Methode als auch die Genotyp-Aufrufmethode ohne Filterung korrelieren stark mit den Schätzungen aus den Sequenom-Genotypdaten (d. h. ein kleiner standardisierter Unterschied zwischen den beiden Schätzungen in Abbildung 6). Schätzungen, die auf Genotypdaten mit Filterung basieren, zeigen jedoch eine schlechte Übereinstimmung mit den aus den Sequenzgenotypdaten geschätzten Häufigkeiten, insbesondere wenn die Sequenzierungstiefe gering ist. Interessanterweise gibt es einen SNP, bei dem sich die geschätzte MAF aus den Resequenzierungsdaten stark von der Schätzung aus den Sequenzgenotypdaten unterscheidet, obwohl die Sequenzierungstiefe sehr hoch ist (14 ×). Insbesondere beträgt die geschätzte MAF aus den Sequenom-Genotypdaten 22,5%, beträgt jedoch 17,2%, wenn sie mit dem ML-Ansatz geschätzt wird. Die individuelle Untersuchung zeigt, dass sich bei vielen Individuen der auf den Sequenzierungsdaten basierende stark unterstützte Genotyp von den Sequenzgenotypen unterscheidet. Angesichts der Tatsache, dass dieser SNP von vielen Lesevorgängen bei diesen Personen abgedeckt wird und dass die beobachteten Lesebasen hohe Qualitätswerte aufweisen (>Q20), ist es wahrscheinlich, dass der Unterschied auf Sequenzgenotypisierungsfehler zurückzuführen ist. Beachten Sie, dass es einige SNPs gibt, bei denen die geschätzten MAFs aus dem Genotyp-Aufruf-Ansatz ohne Filterung besser den aus der Sequenzgenotypisierung geschätzten MAFs zu entsprechen scheinen als die Schätzungen aus dem ML-Ansatz. Zum Beispiel beträgt bei einem SNP der geschätzte MAF 25,7% aus den Sequenom-Genotyp-Daten, 25.9% aus der Genotyp-Ii-Methode ohne Filterung und 27,2% aus der ML-Methode. Die individuelle Untersuchung zeigt jedoch, dass es einige Individuen gibt, bei denen sich der Genotyp aus den Sequenzierungsdaten vom Sequenom-Genotyp unterscheidet. In diesen Fällen werden die Fehler in den genannten Genotypen aufgehoben, was den Anschein einer besseren Übereinstimmung mit den Sequenom-Genotypdaten erweckt. Daher ist es für diese SNPs schwer zu sagen, welche Methode am besten funktioniert.

Als nächstes untersuchten wir die Verteilung von MAFs, die mit verschiedenen Ansätzen über eine Reihe von Sequenzierungstiefen aus unseren Exomsequenzierungsdaten der nächsten Generation berechnet wurden (Abbildung 7). Wir haben SNPs mit geschätzten MAF <2%verworfen, da es schwierig ist, diese sehr niederfrequenten SNPs von Sequenzierungsfehlern in diesem Datensatz zu unterscheiden. Wir entfernten weiter Stellen, an denen ein signifikanter Unterschied (p-Wert weniger als 10-5 unter Verwendung eines Rang-Summen-Tests) im Qualitätsfaktor der Lesebasen zwischen den Neben- und Hauptallelen bestand. Diese Stellen sind wahrscheinlich künstliche SNPs, die aufgrund falscher Zuordnung oder unbekannter Verzerrungen auftreten können, die während des experimentellen Verfahrens eingeführt wurden. Dann haben wir jeden Standort basierend auf der Deckungstiefe in Behälter eingeteilt. Die Anzahl der SNPs in jedem Behälter ist in Tabelle 1 angegeben. Wenn die durchschnittliche Tiefe weniger als 9 × beträgt, unterscheiden sich die Verteilungen der geschätzten MAFs basierend auf Genotyp-2-Methoden stark von denen, die auf der ML-Methode basieren. Insbesondere führen die Genotyp-Aufrufansätze zu einem großen Überschuss an niederfrequenten SNPs (MAF zwischen 2% und 4%). Dieses Muster spiegelt das wider, was wir in unseren Simulationsstudien gesehen haben (Abbildung 3). Auch für die Genotyp aufrufenden Methoden ändert sich die Allelhäufigkeitsverteilung dramatisch, wenn sich die Sequenzierungstiefe ändert. Wenn die Tiefe nicht sehr hoch ist, enthalten die Genotypisierungsaufrufmethoden daher, wie zuvor erörtert, wahrscheinlich viele falsche SNPs, die Sequenzierungsfehler darstellen. Diese Fehler treten als Überschuss an niederfrequenten SNPs in der Häufigkeitsverteilung auf. Die auf der ML-Methode basierende Verteilung ist über Tiefen stabiler, aber es gibt immer noch einen Überschuss an SNPs mit niedriger Allelfrequenz mit einer Tiefe von weniger als 9 × im Vergleich zum Anteil an niederfrequenten SNPs in größeren Tiefen.

Schließlich haben wir diese Exom-Resequenzierungsdaten verwendet, um eine Fall-Kontroll-Assoziationsstudie zu simulieren. Um die Verteilung der Assoziationsteststatistiken unter der Nullhypothese zu untersuchen, haben wir zufällig 100 Personen einer Fallgruppe und die anderen 100 Personen der Kontrollgruppe zugeordnet. Für alle SNPs auf Chromosom 2 mit MAF-Schätzungen > 2% (basierend auf der unknown Minor Allele ML-Methode) testeten wir auf Allelhäufigkeitsunterschiede zwischen Fällen und Kontrollen, indem wir die G-Statistik unter Verwendung von Genotypen mit und ohne Filterung sowie die LRT-Statistik berechneten. Abbildung 8 zeigt die QQ-Diagramme, in denen die Verteilungen der Teststatistiken mit der Standard-χ2 (1) -Verteilung verglichen werden. Wie aus Simulationsstudien hervorgeht, weicht die Nullverteilung der G-Statistik, die beim Aufruf von Genotypen ohne Filterung berechnet wird, wesentlich von der χ2 (1) -Verteilung ab. Die Nullverteilung der LRT-Statistik folgt jedoch eng der χ2(1) -Verteilung. Der Inflationsfaktor beträgt 1,01, was bedeutet, dass die LRT-Statistik bei Anwendung auf reale Daten eine gute Leistung erbringt.