Värdgenetisk variation och dess mikrobiominteraktioner inom Human Microbiome Project

högkvalitativ DNA-sekvensering av HMP-kohorten

HMP-kohortdesignen och provsamlingen har beskrivits ingående . Kortfattat omfattar HMP-kohorten 300 givare rekryterade på två platser i USA. Majoriteten av givarna (71%) var av kaukasiskt ursprung; de återstående givarna var av afrikansk (6%), Asiatisk (9%), Latino (11%) eller blandad (3%) anor. Förhållandet mellan man och kvinna var ungefär lika med 151 kvinnor och 149 män. Målet med att välja givare var att hitta friska individer utan ny medicinsk användning eller sjukdomshistoria, som tillhörde en liknande åldersgrupp (19-40 år) och som hade ett relativt hälsosamt kroppsmassindex (BMI på 19-34 kg/m2).

för att erhålla värdgenominformation sekvenserades genomiskt DNA från blodet hos 298 av de 300 individerna med användning av PCR-fri sekvensering. Den genomsnittliga sekvenseringstäckningen var 32,77 x, med ett intervall på 23,9 till 56,7 till 56,7 (Fig. 1). Kontaminering och andelen chimära läsningar var båda väl under standardavsnittet på 5% i alla prover (ytterligare fil 1: figur S1A). Fördelningen av andra kvalitetsmått som insatsstorlek och procentandel av läsningar som justerades i par markerade inte några tydliga outlierprover; därför inkluderades alla prover i ytterligare analys. Variantnumret var också anmärkningsvärt stabilt vid ~ 2 M singelnukleotidpolymorfismer (SNP) och 200 K indels per person (Fig. 1), med undantag för de afroamerikanska givarna, som hade högre genetisk mångfald (ytterligare fil 1: Figur S1B). Det fanns ingen detekterbar korrelation mellan sekvenseringsdjup och antal varianter som återhämtades, vilket indikerar att djupet i alla prover var tillräckligt.

översikt över Human Microbiome Project värd genom och metagenom täckning. Sekvenseringsdjup för varje värdgenom (vänster) och antal läsningar för alla tillgängliga prover med hel metagenomsekvensering

dessa data ger en nästan fullständig parning av humant genomsekvensering till mikrobiella amplikoner och metagenomer över hela HMP-kohorten. Genetisk variation i denna kohort härleddes tidigare med hjälp av ”förorenande” mänskliga läsningar från 93 försökspersoners WMS-data . Även om detta gav en genomsnittlig mänsklig genomtäckning på ~ 10 kg, varierade det mycket mellan proverna och för många nådde endast 5 kg. Detta var tillräckligt för att detektera 13 M genetiska varianter totalt, 5,5 M med MAF > 0,05. I jämförelse har vår studie mer än tredubblat antalet givare, och genom att direkt rikta in sig på värd-DNA identifierade vi dubbelt så många varianter totalt sett med jämn täckning (minst 25 kcal) mellan proverna. Detta ökade kvaliteten över alla prover och genererade en komplett dataset som kan brytas i detta arbete såväl som i framtida studier.

Sekvenseringsresultat överensstämmer med andra jämförbara populationer

efter filtrering enligt kvalitet och plats i regioner med låg komplexitet kvarstod 29 m varianter, bestående av 26,7 m SNV och 2,3 m insertions-deletions (indels) (ytterligare fil 1: Tabell S1). Jämfört med GoNL-och 1000-Genomprojekten var 5.1 M SNV och 856 K indels nya, men majoriteten av dessa var sällsynta (Fig. 2a). Däremot delades de vanliga varianterna vi identifierade (MAF > 5%) nästan universellt mellan de tre kohorterna. Sammantaget identifierade vi 7.8 M fler varianter jämfört med GoNL-konsortiet med samma storlek. Eftersom många av varianterna också var närvarande i 1000-genomprojektet, tillskriver vi skillnaden till den större etniska mångfalden i HMP-kohorten. Även om ett stort antal SNV var unika för varje kohort, var andelen varianter som faller i introniska, exoniska och intergena regioner i genomet nästan identiska mellan kohorter (ytterligare fil 1: figur S2).

fördelning av genetiska varianter och jämförelse med andra kohorter. en upptäckt varianter kategoriseras efter frekvens och överlappar med andra kohorter. AC-allelräkning, MAF-mindre allelfrekvens. B Fördelning av antalet kodande mutationer efter frekvens och beräknad påverkan

vi kommenterade nästa kodningsvarianter med LOFTEE-plugin för VEP-verktyget , som kategoriserar varianter i klasser baserat på deras inverkan på kodningssekvensen. Antalet varianter med hög effekt, definierade som de som skulle leda till förlust av funktion hos en viss gen, var 2670 (ytterligare fil 1: Tabell S2); detta resultat överensstämmer med aktivt negativt urval mot dessa varianter. Negativt urval var också tydligt från allelfrekvensfördelningen, eftersom svårighetsgraden av en allels påverkan var starkt relaterad till dess frekvens i populationen. Till exempel berikades varianter med hög effekt kraftigt i varianter som endast observerades en gång i vår dataset (Fig. 2B, AC1). Fördelningen av kodande mutationer bland gener var inte heller enhetlig, med ett litet antal gener som fångar ett stort antal varianter. Trettio gener visade mer än fem potentiella varianter med hög effekt, och sex gener hade mer än tio varianter. Det lilla antalet gener med kodande mutationer med hög effekt föreslog att denna kohort var för liten för belastningstestning för att dra korrelationer mellan mutationsfrekvens inom en gen och mikrobiella egenskaper. Istället fokuserade vi vår analys på att identifiera samband mellan vanliga varianter och mikrobiell taxa eller funktionell potential.

mikrobiell taxa och funktionell potential vid sex kroppsställen

för HMP samlades mikrobiomprover från 18 kroppsställen och föll i fem huvudområden: gastrointestinal (GI) kanal, munhålan, hud, nares och vagina. I vissa fall samlades replikatprover över tid för att bedöma tidsmässig stabilitet hos mikrobiomen. Totalt karakteriserades mer än 5000 prover med användning av 16S rRNA-gensekvensering och mer än 2000 med användning av Hagelgevär WMS. Den tidigare metoden ger en översikt över taxonomisk sammansättning på hög nivå, medan den senare tillåter artnivåidentifiering och profilering av mikrobiomets funktionella potential. Vi använde därför WMS-data i efterföljande analyser. Fördelningen av prover med WMS var inte ekvivalent mellan kroppsställen, med de flesta prover från sex platser som representerar fyra av de viktigaste områdena som beskrivs ovan: tarm (avföring), oral (buckal slemhinna, supragingival plack, tunga dorsum), nares (anterior nares) och vaginal (posterior fornix); inga WMS-data fanns tillgängliga från hudproverna. Inom de sex kroppsställena varierade antalet givare från 80 för vaginal posterior fornix till 209 för tarmprover; det genomsnittliga antalet avläsningar per prov varierade från 34 m i posterior fornix till 86 M i tunga dorsum. Med hjälp av WMS-data från dessa prover identifierade vi taxonomisk komposition med MetaPhlan2 och funktionell potential med HUMAnN2 . Dessa utgångar analyserades sedan för föreningar med värdgenetisk variation.

humana genomiska huvudkomponenter korrelerar med mikrobiomkomposition

för att jämföra värdgenetisk variation med mikrobiell variation bedömde vi först graden till vilken genetiska mönster på hög nivå kunde korreleras med mikrobiomkomposition. PCA på de gemensamma SNV: erna (MAF > 0.05) visade att de första fem huvudkomponenterna huvudsakligen representerade givarnas etniska och rasliga anor. Till exempel visade Värdgenetiken hos afroamerikan, en av två grupper av asiatisk-amerikanska och de kaukasiska ämnena den starkaste effekten (Fig. 3a). För att ytterligare jämföra övergripande genetisk variation till andra kohorter ordinerade vi också gemensamt en kombinerad dataset av hmp300 och 1000 genomdeltagare, med SNV med MAF > 0.05 i båda kohorterna. Individer från båda kohorterna fördelade i det resulterande huvudkomponentutrymmet nästan identiskt enligt anor (ytterligare fil 1: figur S3).

korrelation mellan genetisk variation på hög nivå och mikrobiomkomposition. a de två första komponenterna i den genetiska huvudkomponentanalysen visas, baserat på vanliga singelnukleotidvarianter, överlagrade av självrapporterad givar etnicitet. AA Afroamerikansk. B visas är hur mycket varians i mikrobiomdata i genomsnitt kan förklaras av de genetiska huvudkomponenterna, jämfört med permutation på samma data. Värden som visas är Z-poäng baserat på permutationer, som också användes för att beräkna empiriska p-värden. C fördelning av genetiska huvudkomponent R2-värden för olika arter och vägar i avföring. Y-axeln visar variansen förklaras, och X-axeln visar permutationsbaserade empiriska p-värden för vart och ett av dessa siffror. Endast namnen på arter med false discovery rate (FDR) < 0,05 och pathways’ FDR < 0,01 visas. Histogrammet nedan visar fördelningen av empiriska p-värden, och Y-axeln visar antalet arter i ett fack. Gröna staplar under väghistogrammet visar hur vägarna som är associerade med jäsning rankas av R2

därefter beräknade vi för HMP300 vilken procentandel av mikrobiell variation i de sex kroppsställena som kunde förklaras (R2) av de första fem värdgenomets huvudkomponenter. I avföringsprover var procentandelen av artnivåvariationen förklarad av värdens huvudkomponenter 3.8%, högre än förväntat av en slump ensam (empirisk p = 0.0001; Fig. 3b). Fördelningen av empiriska p-värden för R2-värdena för den enskilda arten skiftades starkt mot noll (Fig. 3c), vilket indikerar att starka korrelationer inte var begränsade till några få arter, men den genetiska befolkningsstrukturen påverkade övergripande mikrobiella konfigurationer. Vi observerade en liknande effekt på artnivån på orala platser. I buckal slemhinna beskrivs de genomiska huvudkomponenterna i genomsnitt 5,2% (empirisk p = 0,0008) av artnivåvariationen; i tungdorsum var denna siffra 4,1% (empirisk p = 0,0034). I en identisk analys av metacyc metabolisk väg överflöd, fann vi bara vägarna i tarmmikrobiomerna att vara signifikant korrelerade med vanliga varianthuvudkomponenter. Sammanfattningsvis var sambandet mellan genetiska egenskaper på hög nivå och mikrobiomegenskaper signifikant vid flera kroppsställen.

när man undersökte korrelationen mellan enskilda mikrobiella egenskaper med värdgenetik som bidrog till dessa medelvärden visade vissa egenskaper mycket starkare enskilda föreningar. I avföring, där den genetiska korrelationen var den starkaste, var fem arter av 118 signifikant associerade (false discovery rate (FDR) < 0, 05 genom permutationstest), med R2-värden som nådde nästan 10% (Fig. 3c). Av dessa fem arter, lachnospiraceae bakterie, Roseburia intestinalis och Subdoligranulum (oklassificerad) var alla positivt korrelerade med den första genomiska huvudkomponenten, vilket visar att dessa arter har högre överflöd hos givare av kaukasiskt ursprung. En annan betydande Art, Sutterella wadsworthensis, var associerad med PC4, som separerar givare av asiatiskt ursprung i två grupper. Genom att undersöka andra kroppsställen fann vi att Porphyromonas catoniae, Propionibacterium propionicum och oklassificerad Gemella var signifikant associerade med värdgenomisk variation i buckal slemhinna (ytterligare fil 1: figur S4 och tabell S1).

en liknande vägnivåanalys avslöjade ett stort antal (82 av 541) vägar signifikant (FDR < 0,05 som ovan) korrelerade med genetiska huvudkomponenter i avföring (Fig. 3c). Flera vägar var relaterade till aminosyra och kortkedjig fettsyrabiosyntes och nedbrytning. I en mer systematisk vy fann vi att medlemmarna i jäsningssuperklassen i MetaCyc-databasen berikades signifikant i toppbanans rankningar (Fig. 3c). De flesta av dessa vägar var associerade med den första genetiska huvudkomponenten som skiljer vita givare från andra ras-eller etniska förfäder. Sådana funktionella berikningar kan peka på etniska skillnader i kost, men också på genetisk variation i förmågan att metabolisera vissa näringsämnen.

på andra kroppsställen var variationen på vägnivå i genomsnitt inte korrelerad med genetiska huvudkomponenter, även om vissa enskilda korrelerade vägar hittades (ytterligare fil 1: figur S5 och tabell S1). Till exempel korrelerade ett antal vägar i tunga dorsum-mikrobiomer starkt med genetiska huvudkomponenter. Intressant nog var nästan alla associerade vägar relaterade till andning och trikarboxylsyra (TCA) – cykeln, vilket indikerar en syregradient och skillnader i aerob andning av orala organismer i tungans dorsum mellan givare. Anrikningen av TCA-cykeln i oral mikrobiom och jäsning i tarmmikrobiom återspeglar de dominerande metaboliska egenskaperna hos motsvarande mikrobiomer och visar hur dessa kan påverkas av värdgenetik och miljöfaktorer korrelerade med genetisk anor.

relaterade givare har liknande mikrobiomer

även om HMP-kohorten inkluderade givare relaterade till varandra, var denna information inte tillgänglig i de insamlade metadata. Genomisk sekvensering av givarna gjorde det möjligt för oss att dra slutsatsen om förhållandet mellan alla givarpar och identifiera upp till tredje graders släktingar bland dem. Med hjälp av vanliga SNV (MAF > 0.05) för analysen identifierade vi 11 par första graders släktingar och ett par tredje graders släktingar.

vi försökte därefter bestämma om graden av relation återspeglades i likheten hos deras mikrobiomer. För denna analys beräknade vi Bray-Curtis-avståndet mellan alla givarpar och delade paren i tre grupper: samma etnicitet, olika etnicitet och släktingar (Fig. 4a). Som kan förväntas från PCA, proverna inom etniska grupper var i genomsnitt något mer lika än prover från olika etniska grupper, men mikrobiom likhet mellan besläktade givare var mer uttalad. Med undantag för tarmen, på alla testade kroppsställen, var mikrobiomsamhällets sammansättning mellan släktingar mer lika än mellan slumpmässiga givarpar; i främre nares och buckal slemhinna var effekten också statistiskt signifikant genom t-test mellan orelaterade och relaterade likhetspoäng. För vaginala prover uttalades också effekten, men vi hade inte tillräckligt med kvinnliga kvinnliga par för att uppnå statistisk signifikans.

släktskap och mikrobiom likhet och replikering av kända föreningar. en Bray-Curtis likhet mellan 12 par nära släktingar (tredje graden eller närmare) identifieras från genetiska data jämfört med likheter mellan andra par. P-värdena motsvarar resultaten av T-test mellan likhetspoäng för släktingar, mot alla andra par. b förening mellan FUT2 secretor variant och B. longum. C-samband mellan genetisk variant rs4988235 nära LCT-genen och B. longum. I både b och c visar vi log10 transformerad relativ överflöd

FUT2-och LCT-genotyper är associerade med Bifidobacterium longum

för att studera påverkan av enskilda genetiska varianter på mikrobiomkomposition började vi med kända föreningar av FUT2 och LCT. FUT2 ansvarar för överföringen av terminala fukosrester till mukosala glykaner . Bifidobakterier använder också slem-härledd fukos som en kolkälla, och överflöd och mångfald av B. longum är signifikant lägre i icke-sekretorer (individer med ett för tidigt stoppkodon i FUT2, rs601338) jämfört med sekretorer .

för att avgöra om denna förening kunde verifieras i HMP-kohorten sökte vi efter mikrobiella arter korrelerade med värdfucosyltransferas sekretorgenotyp (MAF = 0,41). Vi använde linjär regression för att förutsäga det relativa överflödet av varje enskild art baserat på sekretorns genotypdosering. B. longum hade den starkaste korrelationen mellan de 118 testade arterna (FDR = 0.018; Fig. 4b), med ökad relativ överflöd i sekretorgenotypen i förhållande till icke-sekretorgenotypen. Detta resultat överensstämmer med tidigare experimentella observationer och visar att kohorten är tillräckligt driven för att validera riktade mikrobiella värdföreningshypoteser.

på liknande sätt har LCT associerats med ökat överflöd av laktosmetaboliserande bifidobakterier i tarmen . LCT kodar laktas, enzymet som är ansvarigt för att bryta ner laktos i övre GI-kanalen; i kombination med ökade bifidobakterier tyder detta på att mer dietlaktos samlas i tjocktarmen. Förmågan att producera laktas vid vuxen ålder eller laktosintolerans (hypolaktasi) styrs av närvaron av en homozygot g-allel i rs4988235 SNV nära LCT . En ny upptäckt att mjölkkonsumtion och bifidobakterier överflöd är positivt korrelerad endast hos personer med hypolactasia genvarianten stöder denna mekanism.

i HMP-kohorten jämförde vi bakteriearter överflöd i avföring mellan givare med hypolaktasi och alternativa varianter. Eftersom hypolaktasi är ett recessivt drag använde vi ett T-test för att jämföra de 64 givarna med den homozygota g-genotypen med resten av de 145 givarna. Efter att ha testat var och en av de 118 enskilda arternas överflöd mot närvaron av denna variant fann vi att B. longum hade den starkaste effekten (FDR = 0,095), vilket bekräftar den tidigare hittade föreningen (Fig. 4c).

mikrobiella föreningar med värdgenomvarianter är kroppsplatsspecifika

slutligen bedömde vi föreningar mellan värdgenom och mikrobiomvariation på ett icke-målinriktat sätt direkt genom en genomövergripande associeringsstudie. Vi utförde analysen separat för varje kroppsställe och koncentrerade oss på SNV med MAF > 0.1 och jämförde dem med både mikrobiella arter – och metaboliska vägnivåer. Vi använde vanliga linjära regressionsmodeller, med hänsyn till effekterna av kön, etnicitet och provsamlingsplats. Efter filtrering av de mikrobiella egenskaperna (se metoder) identifierade vi 120-160 arter i GI-tarmkanalen och orala prover och cirka 25 arter i hud-och vaginalprover. Antalet metaboliska vägar som passerade filtrering var betydligt högre, mellan 350 och 530 huvudvägar per plats. Tillsammans ställer det stora antalet SNV: er, kroppsställen och mikrobiella egenskaper i analysen ett strikt signifikanskriterium (p < 3 kcal 10-12 enligt Bonferroni-korrigering för flera tester), vilket i kombination med vår blygsamma provstorlek begränsar vår upptäcktspotential till föreningar med mycket stora effektstorlekar.

av denna anledning begränsade vi först vår analys till SNV: er som finns i katalogen National Human Genome Research Institute (NHGRI) Genome-Wide Association Studies (GWAS) och antog att dessa SNV: er berikades med genomiska varianter som har potentiell inverkan på mikrobiomegenskaper. Denna uppsättning inkluderade SNV associerade med en mångsidig uppsättning kvantitativa egenskaper som sträcker sig från komplexa sjukdomar till antropometriska mätningar. Totalt 16,869 av dessa SNV hittades i våra data, men vi upptäckte inte några signifikanta föreningar med hjälp av denna delmängd av SNV. Vidare, enligt kvantil-kvantildiagrammet för jämförelserna, fanns det ingen systematisk anrikning av mindre p-värden bland jämförelserna (ytterligare fil 1: figurerna S6 och S7). Vi fick inte signifikanta resultat med ännu mer begränsade variantuppsättningar associerade med inflammatorisk tarmsjukdom eller med någon av GWAS-katalogens underkategorier på hög nivå (t.ex. ”immunförsvar”, ”matsmältningssystem”).

vi körde nästa föreningsanalys på alla vanliga SNV. Vi såg inga föreningar med p-värden mindre än flera testkorrigerade signifikansgränser. Det fanns emellertid ett antal föreningar (ytterligare fil 1: figurerna S8 och S9 och tabell S2) med relativt små p-värden. Dessa föreningar utgör en rik informationskälla för framtida mikrobiomomfattande föreningsstudier och utredare intresserade av specifika gener eller mikrobiella arter.

på grund av den unika utformningen av HMP-kohorten kunde vi också undersöka arten av föreningarna mellan kroppsställen. Studier i uttryck kvantitativ egenskap loci har rapporterat anmärkningsvärd stabilitet av genuttryck-genotypföreningar över vävnader . Här observerade vi inte några SNV-mikrobpar i toppen av föreningsrankingen för flera kroppsplatser, vilket visar den mer indirekta karaktären hos dessa föreningar samt de unika samhällsbyggnadsprinciperna för varje kroppsplats. När vi överväger fördelningen av de bästa p-värdena från varje kroppsställe oberoende av den associerade organismen, kunde vi fortfarande inte hitta några SNV som var starkare associerade med mikrober än förväntat av en slump. Detta resultat tyder vidare på att de genetiska mekanismerna bakom mikrobiomkompositionen är kroppsspecifika.