Bemonstering
Domesticatie schapen (O. aries) en geiten (C. hircus) werden bemonsterd in Iran (IROA en IRCH groepen, respectievelijk) en Marokko (MOOA en MOCH groepen, respectievelijk) voor een totaal van 20 dieren per groep (Supplementary Fig. 6). Deze monsters werden verzameld tussen januari 2008 en maart 2012 in het noordelijke deel van Marokko en tussen augustus 2011 en juli 2012 in het noordwesten van Iran, in het kader van het Europese Nextgen-project (subsidieovereenkomst nr. 244356), overeenkomstig de ethische voorschriften van Richtlijn 86/609/EEG van de Europese Unie. Oorknips werden verzameld uit het distale deel van het oor van willekeurig gekozen dieren, en onmiddellijk opgeslagen in 96% ethanol gedurende een dag alvorens te worden overgebracht in silica-gel korrels tot DNA-extractie.
De wilde soorten Aziatische moeflon (O. orientalis) en Bezoar steenbok (C. aegagrus) werden bemonsterd in het noordwesten van Iran binnen de domesticatie wieg21,22. Dertien Aziatische moeflons en 18 Bezoar steenbok weefsels (respectievelijk IROO en IRCA groepen, supplementaire fig. 6) werden verzameld van in gevangenschap gehouden of recent gejaagde dieren, en van ingevroren monsters die beschikbaar waren op het Iraanse Ministerie van Milieu. Deze op het individu gebaseerde bemonsteringsaanpak is ontworpen om potentiële bias te minimaliseren door de oververtegenwoordiging van lokale effecten te vermijden (b.v,
Aanvullende gegevens
Extra werd een wereldwijd rassenpanel samengesteld voor schapen en geiten (respectievelijk wpOA en wpCH). wpOA omvatte 20 “whole-genome re-sequencing” (WGS) monsters met een 12x-dekking die 20 verschillende wereldwijde rassen vertegenwoordigden, verstrekt door het International Sheep Genome Consortium. wpCH bestond uit 14 WGS-monsters met een 12x-dekking die 9 Europese individuen vertegenwoordigden, d.w.z, 2 Franse Alpiene, en 2 Franse Saanen monsters gesequenteerd door INRA, 5 Italiaanse Saanen monsters geleverd door Parco Tecnologico Padano, en 5 Australische individuen, d.w.z, 2 Boer, 2 Rangeland, en 1 Kasjmier monsters geleverd door de CSIRO (Supplementary Data 5).
Productie van WGS data
Genomisch DNA werd met succes geëxtraheerd uit alle weefselmonsters met behulp van de Macherey Nagel NucleoSpin 96 Tissue kit, het aanpassen van de fabrikant protocol. Weefselmonsters werden genomen in MN square-well blokken om 25 mg fragmenten per monster te verkrijgen. Drie en een half MN vierkant-96 blokken werden bereid, en extractie werd uitgevoerd met behulp van een Tecan Freedom Evo Liquid handler volgens het protocol van de fabrikant. Een pre-lyse stap werd uitgevoerd om monsters te homogeniseren met 180 µl van T1 Buffer en 25 µl van proteinase K overnacht bij 56 ° C. Om de bindingscondities aan te passen, werd 200 µl BQ1-buffer toegevoegd en werd de monsterplaat 1 uur bij 70 °C geïncubeerd; vervolgens werd 200 µl 100% ethanol toegevoegd. De lysaten werden overgebracht naar de Nucleospin Tissue binding plaat en een vacuüm (-0.2 bar, 5 min) werd toegepast om de doorstroom te verwijderen. Er werden drie wasstappen uitgevoerd met BW- en B5-buffers, respectievelijk, en er werd opnieuw vacuüm gezogen om het doorstroomde materiaal te verwijderen. Vóór de elutie van het genomisch DNA werd een silicamembraan van de Nucleospin-weefselbindingsplaat gedurende 10 minuten gedroogd onder vacuüm bij ten minste -0,6 bar. De elutiestap werd uitgevoerd met 100 µl voorverwarmde BE-buffer (70 °C) en een centrifugatiestap bij 3700 rcf gedurende 5 min in 96-PCR-wells. Genomisch DNA werd bewaard bij 4 ° C om vriesdooi te voorkomen en getest op concentratie (als ng / ul) met behulp van de Picogreen-methode en met behulp van een Nanodrop.
Whole genomen werden geresequenced van 500 ng genomisch DNA die werden geschoren tot een 150-700 bp bereik met behulp van de Covaris ® E210 instrument voor elk monster en gebruikt voor Illumina ® bibliotheek voorbereiding door een semi-geautomatiseerde protocol. Eindreparatie, A-tailing, en Illumina ® compatibele adapters (BioScientific) ligatie werden uitgevoerd met behulp van de SPRIWorks Bibliotheek Voorbereiding Systeem en SPRI TE instrument (Beckmann Coulter) volgens het protocol van de fabrikant. Een 300-600 bp grootteselectie werd toegepast om de meeste fragmenten te recupereren. DNA-fragmenten werden geamplificeerd door 12 PCR-cycli met Platinum Pfx Taq Polymerase Kit (Life® Technologies) en Illumina® adapter-specifieke primers. De bibliotheken werden gezuiverd met 0,8x AMPure XP-korrels (Beckmann Coulter), en geanalyseerd met de Agilent 2100 Bioanalyzer (Agilent® Technologies) en qPCR-kwantificering. Bibliotheken werden gesequenced met behulp van 100 base-lengte lees chemie in paired-end flow cell op de Illumina ® HiSeq2000.
Illumina paired-end leest voor Ovis werden in kaart gebracht om de schapen referentie-genoom (OAR v3.1, GenBank assemblage GCA_000298735.146), en voor Capra aan de geit referentie-genoom (CHIR v1.0, GenBank assemblage GCA_000317765.147), met behulp van BWA-MEM48. Het voor elk individu geproduceerde BAM-bestand werd gesorteerd met Picard SortSam en verbeterd met sequentieel Picard MarkDuplicates (http://picard.sourceforge.net), GATK RealignerTargetCreator en GATK IndelRealigner49, en SAMtools calmd50.
Variant ontdekking werd uitgevoerd met behulp van drie verschillende algoritmen: Samtools mpileup50, GATK UnifiedGenotyper51, en Freebayes52. Variant sites werden onafhankelijk geïdentificeerd voor elk van zes groepen, met behulp van de multi-sample modes van de aanroepende algoritmen: (i) 162 monsters van MOOA; (ii) 20 monsters van IROA; (iii) 14 monsters van IROO; (iv) 162 monsters van MOCH; (v) 20 monsters van IRCH; (vi) 19 monsters van IRCA. Voor sommige groepen waren de WGS van meer individuen beschikbaar in het kader van het NextGen-project (zie hierboven). De in deze studie gebruikte monsters werden geselecteerd om waar mogelijk evenwichtige groepen van 20 individuen te verkrijgen. Voor IRCA- en IROO-groepen kwamen in een later stadium extra monsters beschikbaar, die werden toegevoegd voor downstream-analyses. Dieren met een lage alignment en calling kwaliteit werden verwijderd om de uiteindelijke data set te verkrijgen (Supplementary Data 5).
In elke groep waren er twee opeenvolgende rondes van variant site kwaliteit filtering. Filteringstap 1 voegde de oproepen van de drie algoritmen samen, terwijl de oproepen met de laagste betrouwbaarheid werden uitgefilterd. Een variantplaats werd goedgekeurd als hij door ten minste twee verschillende aanroepende algoritmen werd opgeroepen met een fred-variantkwaliteit >30. Een alternerend allel op een site passeerde als het door een van de aanroepende algoritmen werd aangeroepen, en het genotypegetal >0 was. Filterstap 2 gebruikte Variant Quality Score Recalibration door GATK. Eerst genereerden we een trainingsset van de meest betrouwbare variantensites binnen de groep waarbij (i) de site door alle drie de variantcallers werd aangeroepen met een phred-scaled variantkwaliteit >100, (ii) de site biallelisch is, (iii) het minor allele count ten minste 3 is terwijl alleen monsters met genotype phred-scaled kwaliteit >30 werden geteld. De trainingsset werd gebruikt om een Gaussiaans model te bouwen met behulp van de tool GATK VariantRecalibrator, waarbij gebruik werd gemaakt van de volgende variantannotaties van UnifiedGenotyper: QD, HaplotypeScore, MQRankSum, ReadPosRankSum, FS, DP, InbreedingCoefficient. Een Gaussiaans model werd toegepast op de volledige dataset, wat een VQSLOD (log odds ratio van het zijn van een echte variant) opleverde. Sites werden gefilterd als VQSLOD <cutoff-waarde. De afkapwaarde werd voor elke groep als volgt bepaald: Minimum VQSLOD = {de mediane waarde van VQSLOD voor trainingsset varianten}-3 × {de mediane absolute afwijking VQSLOD van trainingsset varianten}. De overgang/transversie SNP-verhouding suggereerde dat het gekozen afkapcriterium het beste evenwicht gaf tussen selectiviteit en gevoeligheid.
SNPs call sets voor zes groepen Ovis en Capra dieren werden gegenereerd (d.w.z., Iraanse en Marokkaanse domestics, en wilds voor elk genus). Omdat de in deze studie uitgevoerde analyses vergelijkingen tussen groepen vereisten, creëerden we genotype-call sets op een consistente set SNP-locaties voor alle dieren uit elke groep. Voor elk genus voegden we de variante-call-sites van zijn drie groepen samen, en behielden alleen biallelische posities zonder ontbrekende gegevens. Genotypes werden opnieuw opgeroepen op elke biallelische SNP-locatie voor alle individuen van belang door GATK UnifiedGenotyper met behulp van de optie GENOTYPE_GIVEN_ALLELES. In dit stadium werd de lijst van individuen uitgebreid met de dieren die behoren tot de panels van wereldrassen van schapen en geiten (wpOA en wpCH) en bijkomende wilde monsters die in dit stadium beschikbaar kwamen (4 O. orientalis en 4 C. aegagrus). Genotypes werden binnen groepen verbeterd en gefaseerd door Beagle 453, en vervolgens uitgefilterd wanneer de genotypekans kleiner was dan 0,95. Tenslotte hebben we sites uitgefilterd die monomorf waren in de verschillende subsets van individuen die in deze studie zijn gebruikt (zie hieronder).
Om de signalen van selectie die tussen Ovis en Capra zijn gedetecteerd te vergelijken, hebben we een kruisvergelijking tussen de twee referentie-genomen uitgevoerd. Eerst hebben we de paarsgewijze uitlijningspijplijn uit de Ensembl release 69 codebase54 gebruikt om de referentie-genomen van schapen (OARv3.1) en geiten (CHIR1.0) uit te lijnen. Deze pijplijn maakt gebruik van LastZ55 om op DNA-niveau uit te lijnen, gevolgd door post-processing waarbij uitgelijnde blokken aan elkaar worden geketend volgens hun plaats in beide genomen. De LastZ pijplijn voor paarsgewijze uitlijning wordt routinematig door Ensembl uitgevoerd voor alle ondersteunde soorten, maar de geit is nog niet in Ensembl opgenomen. Om een vertekening naar een van beide soorten te vermijden, produceerden we twee verschillende interspecifieke uitlijningen. De ene gebruikte schapen als referentiegenoom en geiten als niet-referentie, terwijl de andere geit als referentiegenoom gebruikte en schapen als niet-referentie. Het verschil is dat genomische regio’s van de referentiesoort gedwongen worden zich op unieke wijze te laten verbinden met afzonderlijke loci van de niet-referentiesoort, terwijl genomische regio’s van niet-referentiesoorten zich op meerdere plaatsen van de referentiesoort mogen laten verbinden. Wij verkregen voor segmenten van chromosomen van een referentiegenoom de coördinaten op het niet-referentiegenoom. Tenslotte hebben we voor de SNP’s die in het ene genus werden ontdekt, de volledige genoomuitlijning met het referentiegenoom van het andere genus gebruikt om de overeenkomstige posities te identificeren (supplementaire tabel 6).
Genetische structuur
Om de genetische diversiteit binnen groepen te beschrijven, hebben we VCFtools56 gebruikt om samenvattende genetische variatiestatistieken te berekenen voor de 73 individuen voor Ovis (d.w.z, 13 IROO, 20 IROA, 20 MOOA, en 20 wpOA) en 72 individuen voor Capra (d.w.z. 18 IRCA, 20 IRCH, 20 MOCH, en 14 wpCH). De gemeten statistieken waren het totale aantal polymorfe varianten (S) voor de hele reeks individuen in elk genus en binnen elke groep, de gemiddelde nucleotidediversiteit (π) binnen elke groep en de inteeltcoëfficiënt (F) voor elk individu. Binnen elk genus werden de verschillen tussen de wilde groep en elke gedomesticeerde groep getest met behulp van een eenzijdige t-toets voor individuele inteelt- en genetische belastingwaarden, en een tweezijdige Mann-Whitney toets voor nucleotidediversiteit per site.
De totale divergentie tussen de vier groepen binnen elk genus (d.w.z, wilde, Iraanse en Marokkaanse domestics, en wereldpanel) werd geschat met behulp van alle biallelische SNPs en de gemiddelde gewogen paarsgewijze Fst volgens Weir en Cockerham57 zoals geïmplementeerd in VCFtools56. De genetische structuur tussen groepen werd beoordeeld met de clustermethode sNMF26, na het snoeien van de gegevensreeks om SNP’s met een linkage disequilibrium (r²) van meer dan 0,2 met VCFtools te verwijderen. Linkage disequilibrium (r²) werd berekend tussen paren SNPs binnen glijdende vensters van 50 SNPs, waarbij één SNP per paar willekeurig werd verwijderd wanneer r² groter was dan 0,2. Voor elke sNMF analyse werden vijf runs van hetzelfde aantal clusters (K) uitgevoerd met waarden van K van 1 tot 10. We gebruikten het cross-entropie criterium om de meest waarschijnlijke clustering te identificeren, maar alternatieve partities voor verschillende aantallen K werden ook onderzocht om te beoordelen hoe individuen werden verdeeld over clusters.
Om een onderscheid te maken tussen gedeelde voorouders en vermenging, voerden we TreeMix27 uit om gezamenlijk populatiesplitsingen en daaropvolgende vermengingsgebeurtenissen te schatten met behulp van de gesnoeide gegevensverzameling gebruikt voor sNMF. We voerden TreeMix uit met de -globale optie om onze maximale likelihood conclusies te verfijnen. We hebben de TreeMix boom geworteld met de splitsing tussen wilde en gedomesticeerde individuen. De blokgrootte voor jackknifing was -k 500 SNPs, wat ongeveer overeenkomt met 150 kb, meer dan de gemiddelde blokken van LD gevonden in zowel schapen als geiten. We genereerden een boom met maximale waarschijnlijkheid zonder migratie en voegden vervolgens migratiegebeurtenissen toe en onderzochten de incrementele verandering in de variantie die door het model werd verklaard en de residuele waarden tussen individuen. Het doel was om een potentieel hoge restwaarde of migratierand tussen wilde en gedomesticeerde individuen op te sporen. Om de statistische relevantie van mogelijke vermengingsvectoren, geïdentificeerd door TreeMix (supplementaire tabel 3), verder te onderzoeken, berekenden wij de drie-populatie test f328 als een formele test van genetische introgressie, met behulp van het qp3Pop programma van de ADMIXTOOLS suite58 voor elke combinatie van groepen. Voor Capra, werd de wpCH groep verdeeld tussen Australische rassen, Franse rassen, en Italiaanse rassen. De resultaten worden gerapporteerd in de aanvullende gegevens 2.
Demografische inferentie
Voor elk geslacht hebben we analyses uitgevoerd voor de voorouderlijke demografische inferentie met behulp van het MSMC model, geïmplementeerd in de MSMC2 software25. MSMC is gebaseerd op de paarsgewijze sequentiële Markoviaanse coalescentie59 ; het gebruikt echter haplotypes van gefaseerde genoomsequentiegegevens als invoer. Voor elke analyse gebruikten we twee individuen uit één groep, dus 4 haplotypes. Elke analyse werd herhaald voor een andere willekeurige set van twee individuen, d.w.z. een replicaat van de analyse per groep. Invoer- en uitvoerbestanden werden gegenereerd en geanalyseerd met de python-scripts die bij de MSMC-software werden geleverd en die te vinden zijn op https://github.com/stschiff/msmc-tools. De analyseparameters werden op de standaardwaarden gehouden, behalve de mutatiesnelheid die op 2,5×10-8 werd gesteld en de generatielengte die op 2 jaar werd gesteld. Om de onzekerheid op de tijdschattingen te schatten, varieerden we deze parameters (mutatiesnelheid van 2,5×10-8 en 1,0×10-8 in combinatie met generatielengte van 2 en 4 jaar) en gaven we een ruwe schatting van de domesticatieperiode (zie supplementair fig. 2).
Genetische belasting
Genetische belasting werd op twee manieren geschat. Ten eerste door de genetische belasting voor elk individu te berekenen als de som van de schadelijke fitheidseffecten over alle eiwitcoderende genomische posities volgens de methode van Librado et al.60. Kort samengevat gebruikten wij als proxy voor evolutionaire beperking de PhyloP-scores van de 46-weg zoogdieruitlijning (http://hgdownload.cse.ucsc.edu/goldenPath/hg19/phyloP46way/placentalMammals/). Van deze uitlijning, identificeerden we eiwit-coderende sites evolueren onder functionele beperkingen (phyloP score ≥1.5). Voor elk Ovis of Capra genoom, onderzochten we vervolgens of deze sites gemuteerd waren. Zo ja, dan hebben we de fyloP-scores over alle gemuteerde sites bij elkaar opgeteld, zodat mutaties in sterk beperkte sites proportioneel meer bijdragen aan de totale belastingschatting. Dit leverde een schatting van de belasting voor elk schapen/geitengenoom op. Tenslotte, om een gemiddelde belasting per site te verkrijgen, hebben we gedeeld door het totaal aantal geanalyseerde posities. Op te merken valt dat we de voorwaarden op homozygote plaatsen hebben toegepast om de dominantiecoëfficiënt van mutaties op heterozygote plaatsen niet te modelleren (b.v. recessief, intermediair, dominant). Ten tweede vergeleken we gen-voor-gen de genetische schadelijke belasting in wilde en gedomesticeerde Ovis groepen door het uitvoeren van een Wilcoxon test, met als alternatieve hypothese dat de gedomesticeerde dieren meer belasting hebben dan wilde verwanten. p-waarden werden gecorrigeerd voor meervoudige testen61 en we pasten een drempel van aangepaste p-waarden < 0.05 toe. We voerden een gene ontology enrichment analyse uit op de set van genen die een significante toename in genetische belasting lieten zien met WebGestalt62,63. Aangezien de referentie-genomen slecht geannoteerd zijn voor genen, hebben we ons gebaseerd op single-copy orthologs tussen onze soort en de mens en muis. Genen van het X-chromosoom werden uitgesloten van de achtergrond set. We voerden deze analyse niet uit op Capra vanwege de hogere inteelt die werd waargenomen in de wilde monsters.
Opsporing van selectiehandtekeningen
Voor het opsporen van handtekeningen van selectie in verband met domesticatie gebruikten we alle biallelische SNPs met een minor allelfrequentie van meer dan 0,10 in ten minste één van de drie geteste groepen (d.w.z. de Iraanse en de Marokkaanse domesticatiegroepen, en de wilde groep voor elk genus). Omdat we verwachtten dat tekenen van selectie gerelateerd aan het domesticatieproces aanwezig zouden zijn in alle gedomesticeerde dieren, hebben we de volgende algemene strategie gevolgd: we hebben met hapFLK29 (zie aanvullende noot 5 en aanvullende figuren 9, 10 en 11) voor elk genus de wilde groep getest tegen elk van de traditioneel beheerde gedomesticeerde groepen (d.w.z. Iraans en Marokkaans) en ons geconcentreerd op die gemeenschappelijke regio’s die vermoedelijk onder selectie stonden en die in beide gevallen werden gedetecteerd. De omvang van de steekproeven (n = 13-20) was in overeenstemming met de eisen van de methode29. We hebben visueel gecontroleerd of de consistente tekenen van selectie gevonden met hapFLK ook aanwezig waren in de corresponderende wereld panel set van elk genus, maar hebben deze groepen niet meegenomen in de statistische test vanwege hun multi-ras samenstelling. Tenslotte hebben we gezocht naar gedeelde signalen van selectie tussen Ovis en Capra met behulp van een gestratificeerde FDR benadering. De strategie is weergegeven in supplementaire Fig. 4.
We voerden hapFLK tests uit voor het contrasteren van de wilde groep met elk van de Iraanse en Marokkaanse groepen in elk genus. De verwantschapsmatrix werd berekend uit de genetische afstanden64 van Reynold tussen groepsparen, met gebruikmaking van een willekeurige subset van één procent van de varianten. De afgeleide populatieboom werd gebouwd met behulp van het neighbor-joining algoritme. Voor elke SNP voerden we de hapFLK-test uit die haplotypische informatie integreert om het vermogen om selectieve vegen te detecteren te verhogen. Voor elk getest SNP berekende de hapFLK statistiek de afwijking van de haplotypische frequenties ten opzichte van het neutrale model geschat door de verwantschapsmatrix65. Om gebruik te maken van linkage disequilibrium informatie, maakt hapFLK gebruik van het Scheet en Stephens’66 multipoint model voor multilocus genotypen dat kan worden toegepast op ongefaseerde gegevens. Een van de belangrijkste toepassingen van dit model is het uitvoeren van faseschatting (fastPHASE software66). In onze analyse werd het model getraind op nietgefaseerde gegevens, en daarom wordt in onze analyse rekening gehouden met fase-onzekerheid. De methode werd gebruikt om lokale haplotypen langs chromosomen te hergroeperen in een gespecificeerd aantal clusters K ingesteld op 25, met behulp van een verborgen Markov Model.
Om de gemeenschappelijke regio’s te identificeren die vermoedelijk onder selectie staan in de twee traditioneel beheerde huisgroepen voor elk genus, combineerden we de twee voorgaande hapFLK analyses. Voor elke analyse werden de hapFLK scores aangepast aan een χ2 verdeling om p-waarden te verkrijgen (script beschikbaar bij https://forge-dga.jouy.inra.fr/projects/hapflk/documents). De resultaten van de twee contrasten tussen de wilde groep en elk van de gedomesticeerde groepen werden gecombineerd met behulp van Stouffer’s methode67 om enkele p-waarden te verkrijgen voor de vergelijking van wilde vs. gedomesticeerde dieren. Tenslotte werd het FDR-kader68 toegepast op de gehele SNP-set om de gecombineerde p-waarden om te zetten in q-waarden. SNPs met q-waarden < 10-2 werden behouden en gegroepeerd in genomische regio’s wanneer zij minder dan 50 kb van elkaar verwijderd waren.
Om te onderzoeken of het signaal van selectie gedeeld werd tussen Ovis en Capra, gebruikten wij eerst de kruisuitlijning van de twee referentie-genomen om homologe segmenten te identificeren. Vervolgens hebben wij een gestratificeerd FDR raamwerk69 toegepast. Deze aanpak is gebaseerd op het feit dat er een inherente stratificatie in de tests is gegeven de voorafgaande informatie in de genetische gegevens69, omdat de onderliggende verdeling van ware alternatieve hypothesen verschillend kan zijn naar gelang van de verschillende dynamiek van verschillende genoomgebieden, wat leidt tot verschillende verdelingen van p-waarden. Dit vereist dat FDR-gecorrigeerde p-waarden (d.w.z. q-waarden) afzonderlijk voor de verschillende strata worden verkregen. We hebben gezocht naar convergenties in elk genus door de regio’s die homoloog zijn met de in het andere genus ontdekte regio’s (het gedeelde stratum genoemd) te scheiden van de rest van het genoom (het algemene stratum genoemd). Wij extraheerden de p-waarden afzonderlijk voor elk van de twee gedefinieerde strata en berekenden vervolgens q-waarden via het FDR-kader. Deze gestratificeerde q-waarden waren de uiteindelijke grootheden die voor statistische significantie (<10-2) in aanmerking werden genomen om SNP’s onder selectie op te sporen en in de overeenkomstige genomische regio’s samen te voegen.
Om te testen op convergente handtekeningen van selectie die wilde van gedomesticeerde dieren in beide geslachten onderscheidt, onderzochten we de relatie tussen de significantiedrempel toegepast op q-waarden (die we lieten variëren van 0,2 tot 0,002) in het ene geslacht en de geschatte waarschijnlijkheid dat een SNP is geselecteerd in het gedeelde stratum van het andere geslacht met behulp van Storey et al.70 benadering. Een toename van de afgeleide waarschijnlijkheid met een verlaging van de drempelwaarde die op de q-waarde wordt toegepast (toename van de striktheid) geeft aan dat hoe significanter de regio in het ene geslacht is, hoe waarschijnlijker het is dat we significante SNP’s in het andere geslacht zouden vinden.
We filterden de selectiesignalen die niet consistent waren tussen de drie inheemse groepen. Voor elke gedetecteerde regio gebruikten we de gefaseerde haplotypes van individuen die werden geclusterd met behulp van Neighbor-Joining bomen op basis van het percentage van identiteit tussen de sequenties. Alleen regio’s met consistente signalen werden behouden (aanvullende Fig. 5).
Om af te leiden of de met hapFLK gedetecteerde selectiesignalen duidden op ontspanning van selectie of positieve selectie bij de inheemse bevolking, schatten wij het verschil in nucleotide diversiteit (π) op elke vermeende regio onder selectie tussen de wilde en inheemse groepen. We drukten dit verschil uit als de Δπ index, die voor elke genomische regio werd berekend als het verschil tussen π berekend voor de wilde groep en het gemiddelde van π voor de Iraanse en Marokkaanse domesticatiegroepen, min het verschil in π tussen deze twee groepen berekend over het hele genoom:
Een negatieve waarde zou erop wijzen dat de nucleotidediversiteit lager is in de wilde groep vergeleken met het gemiddelde van de twee huiselijke groepen, en zou worden beschouwd als een teken van verslapping van selectie in deze laatste groepen, diversifiërende selectie in de huiselijke groepen of positieve selectie in de wilde groepen. Een positieve waarde daarentegen zou wijzen op directionele positieve of stabiliserende selectie die in de gedomesticeerde groepen heeft plaatsgevonden. We gebruikten ook de haplotype clustering om in elke regio manueel na te gaan of de gedetecteerde selectieve sweep de door de Δπ index gegeven indicaties bevestigde.
We voerden als volgt functionele interpretaties uit. Voor elke regio onder selectie, beschouwden we de regio plus 50 kb aan elke kant om functionele rollen te identificeren en 5 kb stroomopwaarts en stroomafwaarts van genen en we beoordeelden de overlap tussen deze coördinaten om de genen van belang te behouden. Tenslotte gingen we ervan uit dat een gen gerelateerd was aan een bepaalde gedetecteerde regio wanneer de posities van de regio en het gen elkaar overlapten. Vervolgens beoordeelden we welk gen het meest waarschijnlijk het doelwit van selectie was door het gen te beschouwen dat zich het dichtst bij het topsignaal bevond, d.w.z. de positie van de laagste q-waarde binnen de regio. Genen werden functioneel geannoteerd met behulp van Uniprot (http://www.uniprot.org/), door hun betrokkenheid bij 30 kind-termen (d.w.z. de directe afstammelingen van de termen) van de categorie “Biologisch proces” (d.w.z., GO:0008150) in overweging te nemen. We hebben alle GO-termen die corresponderen met elk gen opgehaald (aanvullende gegevens 4) voor 30 van de 33 categorieën, omdat we drie termen die niet betrokken waren bij zoogdierfuncties (d.w.z. GO:0006791 zwavelgebruik, GO:0006794 fosforgebruik, GO:0015976 koolstofgebruik) buiten beschouwing hebben gelaten. We hebben twee χ2-toetsen uitgevoerd om de verdeling van genen in de GO-categorieën te vergelijken, d.w.z. (i) genen onder selectie uit genus-specifieke regio’s versus die uit homologe regio’s, en (ii) alle genen onder selectie versus de 18.689 menselijke genen geassocieerd met GO-termen in Swiss-Prot. Om de functies van genen in een veeteeltcontext te interpreteren, hebben we ook de beschikbare informatie uit de literatuur over hun fenotypische effecten opgehaald.
Finitief, om de SNPs binnen de eerder gedetecteerde regio’s te vinden die het meest gedifferentieerd waren tussen wilde en gedomesticeerde groepen, hebben we de FLK statistiek gebruikt. Zoals voor hapFLK, vertegenwoordigt het de afwijking van single-marker allelic frequenties ten opzichte van het neutrale model geschat door de verwantschapsmatrix65. Dezelfde procedure werd gebruikt om de scores van de twee analyses aan een χ2 verdeling aan te passen en de verkregen p-waarden te combineren als werd gebruikt voor de hapFLK test. De niet-uniforme verdeling van de p-waarden maakte het echter onmogelijk het FDR-kader toe te passen en wij selecteerden SNP’s binnen de met hapFLK gedetecteerde regio’s met p-waarden <10-4. Voor deze SNP’s gebruikten wij de Variant Effect Predictor (VEP) annotaties71 die werden gegenereerd uit de Ensembl v74 schapen OARv3.1 genoomannotatie voor Ovis (http://www.ensembl.org/Ovis_aries/Tools/VEP) en uit de geiten CHIR1.0 genoomannotatie geproduceerd door de NCBI eukaryotic genome annotation pipeline voor Capra (https://www.ncbi.nlm.nih.gov/genome/annotation_euk/process/). SNPs werden geclassificeerd als intergenic, upstream en downstream (inclusief UTRs), en intronic en exonic posities. De verschillen tussen de distributies van SNP’s met FLK p-waarden <10-4 en alle SNP’s gebruikt voor het detecteren van selectiehandtekeningen werden onderzocht met een χ2-test.
Beschikbaarheid van gegevens
Sequenties en metadata gegevens gegenereerd voor de 73 Ovis en 72 Capra monsters gebruikt in deze analyses zijn publiekelijk beschikbaar. Algemene informatie en alle vcf bestanden zijn te vinden op de Ensembl website (http://projects.ensembl.org/nextgen/). Alle Fastq bestanden, Bam bestanden, en de novo assemblies van O. orientalis en C. aegagrus zijn te vinden op het European Nucleotide Archive (https://www.ebi.ac.uk/ena) onder de toetredingscode van het Nextgen project (PRJEB7436).