schatting van allelfrequentie en associatie mapping met behulp van next-generation sequencing data

Het minor allel is het minder frequente allel in de populatie op een variabele locatie. We beschrijven eerst twee belangrijke benaderingen om de minder belangrijke allelfrequentie (MAF) op een bepaalde plaats in het genoom te schatten. De eerste benadering omvat het afleiden van individuele genotypen en het behandelen van die afgeleide genotypen als zijnde volledig nauwkeurig bij het schatten van de MAF. Vervolgens onderzoeken we de prestaties van een waarschijnlijkheidskader dat direct rekening houdt met de onzekerheid bij het toewijzen van genotypes. Tijdens ons werk gaan we ervan uit dat alle segregerende sites biallelisch zijn.

schatting van MAF uit zogenaamde genotypen

een manier om de MAF te schatten uit sequentiegegevens van de volgende generatie is om eerst een genotype voor elk individu aan te roepen met behulp van sequentiegegevens, en dan deze genotypen te gebruiken alsof ze de ware zijn. Dit was de benadering die traditioneel voor genotype gegevens en sanger wordt gebruikt die gegevens rangschikken. Het is niet duidelijk hoe goed het zal presteren wanneer toegepast op de volgende generatie sequencing gegevens.

een benadering met maximale waarschijnlijkheid kan worden gebruikt om het genotype voor elk individu af te leiden uit de sequentiegegevens van de volgende generatie. Op elke plaats j, voor elk individu i, wordt de waarschijnlijkheid voor elk van de drie mogelijke genotypen (ervan uitgaande dat we het kleine allel kennen) gegeven als:

(1)

waarbij D i,j de geobserveerde sequentiegegevens zijn in individuele i op locatie j, g i , j ∈ {0, 1, 2} is het aantal kleine allelen in het genotype van elk individu, en en controle voor respectievelijk sequencing fouten en lees basiskwaliteiten. De waargenomen het rangschikken gegevens voor elk individu kunnen als de uitlijning van Leest bij plaats j worden gedacht rekening houdend met de gelezen kwaliteitsscores. Dit wordt voorgesteld als de genotype waarschijnlijkheid en wordt gevonden in het genotype waarschijnlijkheid dossier (GLF) dat in vele programma ‘ s wordt geproduceerd die de volgende generatie rangschikkende gegevens, zoals SOAPsnp en MAQ analyseren .

om een genotype aan een bepaald individu toe te wijzen, kan de waarschijnlijkheid van elk van de drie mogelijke genotypen voor het individu worden berekend. Het genotype met de hoogste waarschijnlijkheid kan dan worden toegewezen. Nochtans, verkiezen de onderzoekers vaak een strenger het roepen criterium en zullen geen genotype aan een individu toekennen tenzij het meest waarschijnlijke genotype wezenlijk waarschijnlijker dan de tweede meest waarschijnlijke is. Hier worden de drie mogelijke genotypes gesorteerd naar hun waarschijnlijkheid: , waarbij g(k)overeenkomt met het genotype met de grootste waarschijnlijkheid. Met een bepaalde drempelwaarde f kan men genotype g(1) noemen als . Anders wordt een genotype niet genoemd en wordt het genotype van het individu als vermist beschouwd. Een gemeenschappelijke drempelwaarde van f is 1, wat aangeeft dat het meest waarschijnlijke genotype ten minste tien keer waarschijnlijker is dan het tweede meest waarschijnlijke genotype. Merk op dat dit type van filtering kan resulteren in een hoger vertrouwen voor het “geroepen” genotype, maar het resulteert ook in meer ontbrekende gegevens.

maximale waarschijnlijkheid estimator van allelfrequentie

in plaats van het schatten van de MAF op basis van de genoemde genotypes, een maximale waarschijnlijkheid (ML) methode geïntroduceerd door Kim et al. (zie ook Lynch voor een vergelijkbare benadering) schat MAF ‘ s direct in en houdt rekening met genotype onzekerheid. In het bijzonder, gegeven een klein allel, wordt de waarschijnlijkheid van het observeren van de sequentiegegevens bij elk individu i verkregen door optelling van de waarschijnlijkheden die overeenkomen met alle drie mogelijke genotypes.

stel dat de drie genotype likelihoods gedefinieerd in vergelijking 1 beschikbaar zijn. Met behulp van dezelfde notatie als hierboven, laten D j en p j de waargenomen sequencing gegevens op plaats j en de overeenkomstige MAF, respectievelijk. De genotype waarschijnlijkheid gegeven dat kleine allelfrequentie kan worden berekend door aan te nemen Hardy-Weinberg evenwicht (HWE). Vervolgens, uitgaande van onafhankelijkheid tussen de individuen, de waarschijnlijkheid van de MAF op deze meetkundige plaats is een product van de likelihoods berekend over alle N individuen:

(2)

De ML-schatting van p j kan berekend worden hetzij door rechtstreeks het maximaliseren van de waarschijnlijkheid voor een beperkt parameter ruimte met behulp van de Broyden-Fletcher-Goldfarb-Shanno (BFGS) – methode of met behulp van de verwachting-maximalisatie (EM) algoritme . Bij gebruik van het EM-algoritme wordt voor elk individu de latere verwachting van een genotype berekend en wordt het gemiddelde van die posteriors herhaaldelijk bijgewerkt. Onze implementatie van BFGS was sneller dan het EM algoritme. Bijvoorbeeld, om schattingen te verkrijgen van 100.000 sites, BFG nam ~16 seconden, maar EM nam ~100 seconden. Het verschil in snelheid kan echter specifiek zijn voor de uitvoering. In ons geval, voor beide methoden, stopten we met het bijwerken van parameters toen de toename in de waarschijnlijkheid minder was dan 0,001.

maximale waarschijnlijkheidsschatting met onzeker klein allel

in de praktijk kan vaak het op een na meest voorkomende nucleotide bij individuen als klein allel worden gebruikt. Voor zeldzame SNPs (bijvoorbeeld MAF < 1%) is het echter moeilijk te bepalen welk allel het kleine allel is, aangezien alle vier nucleotiden in sommige lezingen kunnen voorkomen als gevolg van sequentiefouten. Om met deze situatie om te gaan, beschrijven we nu een waarschijnlijkheidskader dat rekening houdt met de onzekerheid bij de bepaling van het kleine allel.

stel dat we voor site j het major allel M. kennen. Merk op dat het beslissen welke van de twee gemeenschappelijke allelen waarschijnlijk de belangrijkste is, niet belangrijk is omdat we ons voornamelijk bezighouden met het schatten van de frequenties bij zeldzame SNPs. Verder is voor allelen met tussenfrequenties (ongeveer 50%) het onderscheid tussen major en minor allel minder belangrijk. Wijs de andere drie niet-belangrijke nucleotiden m1, m2 en m3 toe. De waarschijnlijkheid die in vergelijking 2 wordt geïntroduceerd, gaat uit van een vast major allel M en een vast minor allel m. daarom kan de waarschijnlijkheidsfunctie worden gewijzigd als:

(3)

verder, in de veronderstelling dat een van de drie mogelijke kleine allelen even waarschijnlijk is, verkrijgen we:

(4)

waarbij . Aangezien zeer klein kan zijn met big data sets (bijvoorbeeld met veel individuen), is het nuttig om de waarschijnlijkheid in de log-schaal te berekenen. Orde de drie voorwaardelijke log-likelihoods als naar (l (1), l(2), l(3)), Waar l(1) is de grootste. Vervolgens wordt

G-test met de zogenaamde genotypes voor associatie mapping

In associatiestudies wordt gezegd dat SNP ‘ s die significante verschillen vertonen in allelfrequentie tussen gevallen en controles geassocieerd zijn met het relevante fenotype. Associatie in kaart brengen kan worden uitgevoerd met behulp van gegevens van de volgende generatie sequencing studies. We bespreken eerst benaderingen die individuele genotypes vereisen en voeren dan een test uit voor associatie met behulp van de genoemde genotypes. In deze benadering, wordt een genotype eerst geroepen voor elk individu. De genotypes kunnen worden gefilterd of ongefilterd. Uitgaande van onafhankelijkheid over individuen en HWE, een 2 × 2 contingency tabel kan worden gebouwd door het tellen van het aantal grote en kleine allelen in zowel de gevallen en controles. Dit leidt tot de bekende waarschijnlijkheidsratio-test voor onafhankelijkheid, de G-test:

(5)

waarbij O k,h de frequentie is die in een cel wordt waargenomen, en E k,h de frequentie die wordt verwacht onder de nulhypothese waarin de allelfrequentie tussen gevallen en controles gelijk is. De bekende Chi-kwadraattest van Pearson is asymptotisch equivalent aan de G-test. Als de tabel wordt gegenereerd uit echte genotypen, dan volgt de G-statistiek asymptotisch een chi-kwadraat verdeling met 1 vrijheidsgraad (χ2 (1)). Nochtans, in onze studies, construeren wij de G-statistiek gebruikend “geroepen” genotypes, dus kan HWE wegens over-en onder-roepen van heterozygotes niet houden. Bovendien, het construeren van de teststatistiek door het tellen van “genaamd” genotypen in plaats van “waargenomen” genotypen waarschijnlijk leidt tot extra variabiliteit. Daarom is de statistische theorie misschien niet meer geldig. Merk op dat wanneer een genotype niet wordt gevraagd voor een bepaalde persoon, de gegevens worden geacht ontbreken en is niet opgenomen in de 2 × 2 tabel.

Waarschijnlijkheidsratio test rekening houdend met onzekerheid in de waargenomen genotypes voor associatie mapping

in plaats van genotypes aan te roepen, maakt het waarschijnlijkheidskader onzekerheid mogelijk in de genotypes en tests op elke locatie j of de allelfrequentie gelijk is tussen gevallen en controles. Formeel berekenen we de waarschijnlijkheid van de hypothesen H o: p j, 1 = p j, 2 (=p j ,0) en H A : p j ,1 ≠ p j, 2, waarbij p j, 1 en p j ,2 de MAF ‘ s zijn in respectievelijk gevallen en controles.

aangenomen dat kleine (m) en grote (M) allelen bekend zijn, kan de waarschijnlijkheid van de kleine allelfrequentie worden berekend zoals beschreven in vergelijking 2, en wordt de waarschijnlijkheid ratio teststatistiek berekend als:

(6)

waar en zijn de waargenomen gegevens voor gevallen en controles, respectievelijk, en en zijn de MLEs van de MAFs in de gevallen en controles, respectievelijk.

indien het kleine allel onbekend is, wordt de waarschijnlijkheid onder de nulhypothese berekend zoals in vergelijking 3, en wordt de LRT-statistiek gewijzigd als:

(7)

waarbij D j de waargenomen gegevens is voor zowel gevallen als controles, en is de allelfrequentie onder de nulhypothese. Andere notaties zijn hetzelfde als in vergelijking 6.

schatten van MAF in gesimuleerde gegevens

vergelijken we de schattingen van de allelfrequentie op gesimuleerde gegevens met behulp van echte genotypen (True), genotypen genoemd zonder filtering (Call NF), genotypen genoemd met filtering (f = 1; Call F), en de maximum likelihood method (ML). Voor zeldzame SNPs is het kleine alleltype vaak niet duidelijk. Wanneer het roepen van genotypes, wordt de tweede gemeenschappelijkste nucleotide verondersteld om het minder belangrijke allel te zijn. De ML-methode bevat direct onzekerheid bij het bepalen van het minor allel en tenzij anders vermeld, worden resultaten met behulp van de onbekende minor allel-methode (vergelijking 3) weergegeven. Merk op dat de onbekende minor allel ML methode presteert vergelijkbaar met de bekende minor allel ML methode, maar de eerste beter voor zeer zeldzame SNPs (extra bestand 1).

we hebben eerst geëvalueerd hoe goed de verschillende benaderingen in staat waren om de MAF te schatten in 200 individuen over een reeks van sequentiediepten voor 1.000 SNPs met een echte MAF van 5%. Figuur 1 toont boxplots van de distributies van geschatte MAF ‘ s met behulp van de vier verschillende benaderingen. Zoals verwacht, voor gegevens met een hogere dekking, zoals een individuele diepte van 12×, presteren alle methoden evenals wanneer de genotypes met zekerheid bekend zijn (True). Echter, wanneer de diepte afneemt, worden de schattingen van de MAF verkregen door eerst genotypes te roepen bevooroordeeld. Bijvoorbeeld, de mediane MAF geschat met behulp van de Call F methode is 5,3% bij 6× dekking en is 12,5% bij 2×. De reden voor de opwaartse vooringenomenheid is dat het moeilijker wordt om heterozygoten te noemen omdat echte heterozygoten er vaak uitzien als sequentiefouten. Daarom hebben meer heterozygoten dan kleine homozygoten de neiging om ontbrekende genotypes te hebben. Echter, de Algemene bias in MAF schattingen van de zogenaamde genotypes is niet altijd in één richting (gegevens niet weergegeven). Interessant is dat de bias slechter lijkt te zijn voor de Call F methode dan de Call NF methode. Dit patroon kan contra-intuïtief lijken omdat het filteren van de genotype-oproepen de kans lijkt te verminderen om een sequentiefout een heterozygoot te noemen. Nochtans, resulteert de aanroep F methode ook in een grotere hoeveelheid ontbrekende gegevens aangezien vele homozygotes voor het belangrijke allel niet wegens het rangschikken van fouten zullen worden geroepen. Dus, in dit geval, het aanroepen van genotypes zonder filtering lijkt de betere strategie dan het filteren van genotypes bij het proberen om de MAF schatten.

figuur 1
figure1

schattingen van de allelfrequentie op locaties met een echte MAF van 5% voor verschillende dekkingsdiepte. Op elke diepte werden 1.000 locaties gesimuleerd met 200 individuen, en op elke locatie wordt een schatting van de allelfrequentie berekend met behulp van: (1) echte genotypen (True); (2) genotypen zonder filtering genoemd (Call NF); (3) genotypen met filtering genoemd (Call F); en (4) de maximale waarschijnlijkheidsmethode (ML). Zie methoden voor meer informatie over de schattingsmethoden.

de resultaten zijn dramatisch verschillend voor de nieuwe ML-methode. Deze methode verstrekt onbevooroordeelde schattingen van de MAF (mediaan van ~4,9%) over een waaier van diepten. Zelfs bij 2× tonen de schattingen slechts een iets grotere variantie dan die op basis van de echte genotypes.

We vergeleken ook de geschatte gemiddelde kwadraatfout (MSE; verwachting () van de verschillende schattingen van de MAF over een reeks sequentiediepten (Figuur 2). De ML methode heeft een lagere MSE dan de aanroepende methoden met 50 of 200 individuen. In het bijzonder, de MSE berekend op basis van de aanroep F methode is veel hoger dan die van de andere methoden vooral wanneer de diepte afneemt. De MSE van de schattingen van de MAF op basis van de werkelijke genotypes weerspiegelt de ondergrens van de MSE en is niet constant over dieptes toe te schrijven aan steekproefvariantie en een eindige steekproefgrootte. Gebruikend 50 individuen, benadert MSE 0.0005 met stijgende diepte en wanneer het gebruiken van een steekproefgrootte van 200 individuen, benadert het 0.0013 met stijgende diepte.

Figuur 2
figure2

gemiddelde squred error (MSE; verwacht ) van vier verschillende typen allelfrequentie-schatters voor verschillende steekproefgrootten (linker-en rechterpaneel) en dekkingsdiepte (X-as). Op elke diepte werd MSE berekend op basis van de allelfrequentie-schattingen die werden gemaakt met behulp van vier verschillende methoden: True, aanroep NF, aanroep F en ML (voor details van de methoden, zie het bijschrift van Figuur 1).

over het geheel genomen presteert de nieuwe ML-methode beter dan genotype aanroepende methoden.

schatten van een verdeling van MAFs op basis van gesimuleerde gegevens

vervolgens onderzoeken we hoe de verschillende schattingsbenaderingen werden uitgevoerd bij het schatten van het aandeel SNP ‘ s op verschillende frequenties in de populatie (vergelijkbaar met het frequentiespectrum van de locatie, maar gebaseerd op populatieallelfrequentie in plaats van steekproeffrequentie). Hier simuleerden we 20.000 SNP ’s waarbij de verdeling van de echte MAF’ s de standaard stationaire verdeling volgde voor een effectieve bevolkingsgrootte van 10.000 (zie methoden). Merk op dat het in de praktijk echter zeer moeilijk is om een zeer zeldzame SNP te onderscheiden van een sequentiefout. Daarom hebben we SNP ‘ s met een geschatte MAF van minder dan 2% weggegooid voor vergelijkingsdoeleinden met echte gegevens. Figuur 3 toont het aandeel SNP ’s dat in elke verschillende frequentiebak valt na uitsluiting van SNP’ s met een geschatte MAF<2%.

Figuur 3
figure3

verdeling van allelfrequenties van SNP ‘ s gesimuleerd uitgaande van de standaard stationaire verdeling van allelfrequenties. Op elke diepte (elk paneel) werden 20.000 SNPs gesimuleerd, en voor elke SNP werden schattingen van de MAF verkregen met behulp van vier verschillende methoden (zie de titel van Figuur 1). Vervolgens worden voor elke methode (elke kleur) alleen sites met geschatte allelfrequenties > 2% gebruikt om elk histogram (x-as) te genereren.

zoals verwacht, met een hoge dekkingsdiepte, zoals 10× per individu, bieden alle methoden geschatte MAF-distributies die vergelijkbaar zijn met de verwachte distributie op basis van de werkelijke genotypes (Figuur 3). Bij een ondiepere dekkingsdiepte, zoals minder dan 4× per individu, wijken de distributies van Maf ‘ s die zijn verkregen met genotype calling-methoden significant af van de verwachte MAF-distributie op basis van echte genotypes (Figuur 3). Met name deze methoden overschatten het aandeel van laagfrequente SNP ‘ s. Bijvoorbeeld, het verwachte aandeel van SNP ‘ s in de tweede bin (geschatte MAF tussen 2-4%) is 18%. Het overeenkomstige aandeel op basis van de Call NF-methode op een diepte van 4× is 26%, wat 1,4 maal hoger is dan verwacht. De overschatting van het aandeel van laagfrequente SNPs gebeurt als gevolg van verwarring van het rangschikken van fouten met echte heterozygoten, wat resulteert in overcalling heterozygote genotypes. De omvang van deze inflatie verschilt tussen de verschillende filtering cut-offs, maar een grotere cut-off hoeft niet per se de inflatie te verhogen of te verlagen.

de afbeelding is geheel anders voor de ML-methode. De met de nieuwe ML-methode verkregen geschatte verdeling van de MAF volgt de werkelijke verdeling nauwkeurig, zelfs bij ondiepe dieptes van de dekking. Hier is er bijna geen overmaat aan laagfrequente SNPs. Bij een diepte van 4× is het aandeel SNP ‘ s in de tweede bak van het histogram 18,4%, wat zeer dicht bij het verwachte aandeel (18%) ligt. Zo kunnen betrouwbaardere schattingen van het frequentiespectrum worden gemaakt op basis van gegevens met een lage dekking door gebruik te maken van onze waarschijnlijkheidsbenadering dan door gebruik te maken van de genotype-aanroepbenaderingen.

associatie mapping in gesimuleerde gegevens

We vergelijken de prestaties van methoden die afgeleide genotypen behandelen als echte genotypen in associatietests (met behulp van een G-test) met onze waarschijnlijkheidsratio-test (LRT) die de onzekerheid in de genotypen verklaart. We onderzoeken de verdeling van de teststatistiek onder de nulhypothese van geen allelfrequentie verschil tussen cases en controles. We vergelijken ook de kracht van de verschillende benaderingen.

bij redelijk grote steekproefgroottes suggereert de standaard asymptotische theorie dat onder de nulhypothese zowel de G-statistiek als de LRT-statistiek een chi-kwadraatverdeling volgen met één graad van vrijheid (χ2(1)). Daarom hebben we de nulverdeling van de G-statistiek berekend op basis van aanroepende methoden en de LRT-statistiek vergeleken met de χ2(1) – verdeling met behulp van QQ-plots (Figuur 4). We simuleerden 5.000 SNP ‘ s over verschillende sequentiediepten in 500 cases en controles waarbij de MAF gebruikt om genotypes te simuleren 5% was in beide gevallen en controles. De verdeling van de G-statistiek, berekend aan de hand van de werkelijke genotypen, vertoont een zeer goede overeenkomst met een χ2(1) – verdeling. De verdeling van de G-statistiek berekend op basis van de genoemde genotypes wijkt echter aanzienlijk af van een χ2(1) verdeling. Het aanroepen van genotypen en vervolgens het behandelen van deze genotypen Als accuraat produceert een enorme overmaat aan vals-positieve signalen als de p-waarden worden berekend met behulp van een χ2(1) distributie. Op een diepte van 2× had 11% van de SNP ‘ s bijvoorbeeld een p-waarde van minder dan 5%, vergeleken met de verwachte 5%. Het effect wordt veroorzaakt door een toename van de variantie, als gevolg van overcalling homozygoten als heterozygoten, in de allelische test die hier wordt gebruikt voor het detecteren van associatie. Genotypische tests zoals Armitage trendtest, die robuust zijn voor afwijkingen van het Hardy-Weinberg-evenwicht, tonen geen vergelijkbare toename van de vals-positieve snelheid (aanvullend dossier 2). In overeenstemming met deze observatie, het filteren van de genoemde genotypes resulteert in een afname van de fractie van significante tests bij het gebruik van de G-test, hoewel het filteren het probleem niet volledig oplost. Aan de andere kant laat de LRT-statistiek slechts een zeer geringe afwijking zien van een χ2(1) – verdeling voor 2× of 5× dieptes van de dekking.

Figuur 4
figure4

QQ-plots die de nulverdeling van de teststatistiek vergelijken met een χ2(1) – verdeling. Elke kolom komt overeen met een andere teststatistiek: (1) g-statistiek berekend met behulp van de echte genotypen (True); (2) G-statistiek berekend met behulp van de zogenaamde genotypen zonder filtering (aanroep NF); (3) G-statistiek berekend met behulp van zogenaamde genotypes met filtering (Call F); en (4) de waarschijnlijkheid ratio test statistiek met onbekend minor allel (LRT). Uitgaande van 500 gevallen en 500 controles, onder de nulhypothese, werd een reeks van 5.000 plaatsen gesimuleerd met een MAF van 5% met een sequentiediepte van 2× (bovenste panelen) en 5× (onderste panelen). De” inflatie ” factor wordt weergegeven in de linkerbovenhoek van elk cijfer.

we hebben ook curves (receiver operating characteristic) (ROC) gegenereerd voor elk van de verschillende associatietesten. Deze curven tonen het vermogen van de test bij verschillende fout-positieve snelheden. Aangezien de verdelingen van sommige van de teststatistieken niet volgens de nulhypothese de χ2(1) – verdeling volgen, verkregen wij, om een eerlijke vergelijking te maken, de kritische waarde voor elke fout-positieve rate op basis van de empirische nulverdeling. Het vermogen wordt berekend als de fractie van gesimuleerde ziekteloci die een statistiek hebben die de kritische waarde overschrijdt. In het algemeen vinden we dat de LRT beter presteert dan de G-test op basis van beide genotype calling methode (Figuur 5). Bij een foutpositief percentage van 5% en met een sequentiediepte van 5× is het vermogen om een ziektelocus met een MAF van 1% en een relatief risico (RR) van 2 op te sporen bijvoorbeeld 51% met de LRT, maar het vermogen daalt tot 33% met behulp van de oproepende methode zonder filtering en tot 34% met behulp van de oproepende methode met filtering. In het bijzonder, op lage diepte, presteert de G-test toegepast op de zogenaamde genotypes met filtering zeer slecht (links meeste kolom in Figuur 5). Als we de kracht van de LRT vergelijken met de Armitage trendtest met de zogenaamde genotypes, zien we dat de LRT ook een hoger vermogen heeft dan de Armitage trendtest (aanvullend dossier 3). Dit suggereert dat als men wil gebruiken genaamd genotypes, filteren ze op basis van call vertrouwen kan resulteren in een verlies van macht.

Figuur 5
figure5

Curves van vier associatietests. Voor de definitie van de vier statistieken, zie het bijschrift van Figuur 4. Uitgaande van 500 gevallen en 500 controles, werden een reeks van 20.000 plaatsen gesimuleerd onder de nul en onder het alternatief bij individuele sequentiediepte van 2×, 5× en 10× (drie kolommen). Bij elke fout-positieve snelheid (x-as) werd de overeenkomstige kritische waarde berekend met behulp van de empirische nulverdeling. De werkelijke positieve snelheid (vermogen; y-as) werd verkregen door het berekenen van de fractie van oorzakelijke sites met teststatistieken die de kritische waarde overschrijden.

Application to real data

we analyseerden 200 exomen uit de controlegroep voor een ziekte associatie studie die zijn gesequenced met behulp van Illumina technologie op een per-individuele diepte van 8× . We gebruikten het genotype likelihoods gegenereerd door het” SOAPsnp ” programma voor onze gevolgtrekking. Zie methoden voor meer informatie.

eerst hebben we de nauwkeurigheid van de schattingen van de MAF van de volgende generatie sequencing data voor 50 SNP ’s onderzocht door ze te vergelijken met de geschatte MAF’ s van sequenom genotype data. Zowel de schattingen met behulp van de ML-methode als de genotype calling-methode zonder filtering zijn sterk gecorreleerd met de schattingen gemaakt van de Sequenom genotype data (d.w.z., een klein gestandaardiseerd verschil tussen de twee schattingen in Figuur 6). Echter, schattingen gebaseerd op genotype bellen met filtering tonen slechte overeenstemming met de frequenties geschat op basis van de Sequenom genotype data, vooral wanneer sequencing diepte laag is. Interessant is dat er één SNP is waar de geschatte MAF uit de resequencing gegevens zeer verschillend is van de schatting verkregen uit de Sequenom genotype gegevens, hoewel de sequencing diepte zeer hoog is (14×). Specifiek is de geschatte MAF van de sequenom genotype gegevens 22,5%, maar is 17,2% wanneer geschat met behulp van de ML benadering. Het individuele onderzoek toont aan dat in vele individuen, het hoogst gesteunde die genotype op de rangschikkende gegevens wordt gebaseerd van de Sequenomgenotypes verschilt. Gezien het feit dat deze SNP wordt gedekt door veel reads in deze individuen en dat de waargenomen leesbases hebben hoge kwaliteit scores (>Q20), is het waarschijnlijk dat het verschil is te wijten aan Sequenom genotyping fouten. Merk op dat er een paar SNP ‘ s zijn waarin de geschatte MAFs van de genotype aanroepende benadering zonder filtering beter lijken te corresponderen met de Mafs geschat van de Sequenom genotyping dan de schattingen van de ML benadering doen. Bijvoorbeeld, bij één SNP is de geschatte MAF 25,7% van de sequenom genotype gegevens, 25.9% van de genotype calling-methode zonder filtering, en 27,2% van de ML-methode. Nochtans, onthult de individuele inspectie er een paar individuen zijn waarvoor het geroepen genotype van de het rangschikken gegevens van het Sequenom genotype verschilt. In deze gevallen, annuleerden de fouten in de geroepen genotypes, die de verschijning van betere correspondentie met de sequenom genotype gegevens geven. Daarom is het voor deze SNP ‘ s moeilijk om te zeggen welke methode het beste presteert.

Figuur 6
figure6

schattingen van allelfrequentie berekend uit 200 personen met behulp van sequenom-gegevens van de volgende generatie. Op elke locatie werden alleen individuen met zowel sequenom genotype data als sequencing data gebruikt voor de schatting van allelfrequentie. Voor het rangschikken van gegevens, werden schattingen van MAF verkregen gebruikend drie verschillende methodes (aanroep NF; aanroep F; en ML). Het gestandaardiseerde verschil voor elke schatting werd berekend als , waarbij en de geschatte MAF ‘ s zijn op basis van respectievelijk de sequencinggegevens en de sequenom-genotype-gegevens, en n het aantal personen is dat Voor de schatting is gebruikt. Elke locatie wordt ingedeeld in een van de vier bakken op basis van de gemiddelde individuele dekkingsdiepte (kleur): minder dan 4×, hoger dan 4× maar minder dan 8×, hoger dan 8× maar minder dan 16×, en hoger dan 16×.

vervolgens onderzochten we de verdeling van MAFs die werd berekend met behulp van verschillende benaderingen over een reeks sequentiedieptes van onze volgende generatie exome sequentiegegevens (Figuur 7). We hebben SNP ‘ s met geschatte MAF <2% weggegooid omdat het moeilijk is om deze SNP ‘ s met zeer lage frequentie te onderscheiden van sequentiefouten in deze dataset. We verwijderden verder sites waarin er een significant verschil was (p-waarde minder dan 10-5 met behulp van een rank-sum-test ) in de kwaliteitsscore van leesbases tussen de kleine en grote allelen. Deze plaatsen zullen waarschijnlijk kunstmatige SNPs zijn die wegens onjuiste in kaart brengen of Onbekende vooroordelen kunnen voorkomen die tijdens de experimentele procedure worden geà ntroduceerd. Dan hebben we elke site ingedeeld in bakken op basis van de diepte van de dekking. Het aantal SNP ‘ s in elke bak is weergegeven in Tabel 1. Wanneer de gemiddelde diepte minder dan 9× is, zijn de distributies van geschatte MAFs op basis van genotype aanroepende methoden zeer verschillend van die op basis van de ML-methode. In het bijzonder, de genotype roepen benaderingen leiden tot een grote overmaat van laag-frequente SNPs (MAF tussen 2% en 4%). Dit patroon weerspiegelt wat werd gezien in onze simulatie studies (Figuur 3). Ook voor het genotype roepende methoden, verandert de verdeling van de allelfrequentie dramatisch als het rangschikken van diepte verandert. Daarom, zoals eerder besproken, wanneer de diepte niet zeer hoog is, zullen de genotyperende roepende methodes waarschijnlijk heel wat valse SNPs omvatten die fouten rangschikken. Deze fouten verschijnen als een overmaat van laagfrequente SNPs in de frequentieverdeling. De verdeling op basis van de ML-methode is stabieler over dieptes, maar er is nog steeds een overmaat aan SNPs met een lage allelfrequentie met een diepte van minder dan 9×in vergelijking met het aandeel van SNPs met een lage frequentie op grotere dieptes.

Figuur 7
figure7

verdeling van de minimale allelfrequentie geschat op basis van de exomes van 200 gesequenced individuen. Voor elke plaats werd de minimale allelfrequentie geschat met behulp van vier verschillende methoden.: (1) De ML-methode met onbekend klein allel, (2) De ML-methode met een bekend of vast klein allel, (3) aanroepgenotypen zonder filtering (aanroep NF) en (4) aanroepgenotypen met filtering (aanroep F). Elke site is ingedeeld in bakken op basis van de diepte van de dekking. Bovendien worden in elk histogram locaties met een geschatte MAF van minder dan 2% niet in aanmerking genomen. Voor het aantal SNP ‘ s dat Voor deze analyse werd gebruikt, zie Tabel 1.

Tabel 1 Aantal SNP ‘ s met een geschatte MAF groter dan 2% met behulp van een bepaalde methode (rij) binnen elke bin (kolom) gedefinieerd door gemiddelde sequentiediepte tussen individuen.

ten slotte hebben we deze exome-resequencing gegevens gebruikt om een case-control associatiestudie te simuleren. Om de verdeling van de associatieteststatistieken onder de nulhypothese te onderzoeken, hebben we willekeurig 100 individuen toegewezen aan een case-groep en de andere 100 aan de controlegroep. Voor alle SNP ‘ s op chromosoom 2 met MAF schattingen > 2% (gebaseerd op de onbekende Minor allel ml methode), hebben we getest op allelfrequentie verschillen tussen gevallen en controles door de G-statistiek te berekenen met behulp van zogenaamde genotypes zowel met als zonder filtering, evenals de LRT statistiek. Figuur 8 toont de QQ-percelen die de verdelingen van de teststatistieken vergelijken met de standaard χ2(1) – verdeling. Zoals blijkt uit simulatiestudies, wijkt de nulverdeling van de G-statistiek, berekend bij het aanroepen van genotypen zonder filtering, aanzienlijk af van de χ2(1) – verdeling. De null-verdeling van de LRT-statistiek volgt echter nauw de verdeling van χ2(1). De inflatiefactor is 1,01, wat betekent dat de GRF-statistiek goed presteert wanneer deze wordt toegepast op reële gegevens.

Figuur 8
figure8

QQ-plots die de associatieteststatistieken voor allelfrequentieverschillen tussen 100 gevallen en 100 controles vergelijken met een χ2(1) – verdeling. Fenotypen werden willekeurig toegewezen aan individuele individuen in de exome resequencing dataset zodat er 100 gevallen en 100 controles. Voor elke site werden drie statistieken berekend: de G-statistiek die genotypes zonder filtering gebruikt (aanroep NF), de G-statische die genotypes met filtering gebruikt (aanroep F), en de LRT-statistiek. Om de opname van valse SNP ‘ s te minimaliseren, worden sites met ML MAF-schattingen minder dan 2% weggegooid. Voor weergavedoeleinden worden resultaten van sites op chromosoom 2 getoond. Merk op dat de inflatiefactor wordt weergegeven in de linkerbovenhoek van elke QQ-grafiek.

Related Posts

Geef een antwoord

Het e-mailadres wordt niet gepubliceerd. Vereiste velden zijn gemarkeerd met *