Genetic Ancestry and Natural Selection Drive Population Differences in Immune Responses to ...

35 Pages·2016·4.55 MB·English

by  Yohann Nédélec

Article
Genetic Ancestry and Natural Selection Drive
Population Differences in Immune Responses to
Pathogens
Graphical Abstract Authors
YohannNe´de´lec,Joaquı´nSanz,
GolshidBaharian,...,JennyTung,
VaniaYotova,LuisB.Barreiro
Correspondence
[email protected]
In Brief
Differencesinthetranscriptional
responsetoinfectioninhuman
populationsareunderstronggenetic
influence,dictatedbytheirancestryand
byrecentnaturalselectionevents.
Highlights Data Resources
d Thousandsofgenesshowpopulationdifferencesin GSE81046
transcriptionalresponsetoinfection
d Africanancestryisassociatedwithastrongerinflammatory
response
d Populationdifferencesinimmuneresponseareoften
geneticallycontrolled
d Naturalselectioncontributedtoancestry-associated
differencesingeneregulation
Ne´de´lecetal.,2016,Cell167,657–669
October20,2016ª2016ElsevierInc.
http://dx.doi.org/10.1016/j.cell.2016.09.025
Article
Genetic Ancestry and Natural Selection
Drive Population Differences
in Immune Responses to Pathogens
YohannNe´de´lec,1,2,11Joaquı´nSanz,1,2,11GolshidBaharian,1,2,11ZacharyA.Szpiech,3AlainPacis,1,2AnneDumaine,2
Jean-ChristopheGrenier,2AndrewFreiman,4AaronJ.Sams,5StevenHebert,2ArianePage´ Sabourin,2FrancescaLuca,4,6
RanBlekhman,7RyanD.Hernandez,3,8RogerPique-Regi,4,6JennyTung,9VaniaYotova,2andLuisB.Barreiro2,10,12,*
1DepartmentofBiochemistry,FacultyofMedicine,Universite´ deMontre´al,Montreal,QCH3T1J4,Canada
2DepartmentofGenetics,CHUSainte-JustineResearchCenter,Montreal,QCH3T1C5,Canada
3DepartmentofBioengineeringandTherapeuticSciences,UniversityofCalifornia,SanFrancisco,SanFrancisco,CA94143,USA
4DepartmentofObstetricsandGynecology,WayneStateUniversity,Detroit,MI48202,USA
5DepartmentofBiologicalStatistics&ComputationalBiology,CornellUniversity,Ithaca,NY14850,USA
6CenterforMolecularMedicineandGenetics,WayneStateUniversity,Detroit,MI48202,USA
7DepartmentofGenetics,CellBiology,andDevelopmentandDepartmentofEcology,UniversityofMinnesota,TwinCities,MN55108,USA
8InstituteforHumanGeneticsandQuantitativeBiosciencesInstitute,UniversityofCalifornia,SanFrancisco,SanFrancisco,CA94143,USA
9DepartmentsofEvolutionaryAnthropologyandBiologyandDukePopulationResearchInstitute,DukeUniversity,Durham,NC27708,USA
10DepartmentofPediatrics,FacultyofMedicine,Universite´ deMontre´al,Montreal,QCH3T1J4,Canada
11Co-firstauthor
12LeadContact
*Correspondence:[email protected]
http://dx.doi.org/10.1016/j.cell.2016.09.025
SUMMARY frompathogens(BarreiroandQuintana-Murci,2010;Fumagalli
et al., 2011; Karlsson et al., 2014). When human populations
Individualsfromdifferentpopulationsvaryconsider- migrated out of Africa, they encountered markedly different
ably in their susceptibility to immune-related dis- pathogenicenvironments,likelyresultinginpopulation-specific
eases.Tounderstandhowgeneticvariationandnat- selectionontheimmuneresponse(BarreiroandQuintana-Murci,
ural selection contribute to these differences, we 2010;Fumagallietal.,2011;Karlssonetal.,2014).Substantial
tested for the effects of African versus European evidencesupportsthishypothesisatthegeneticlevel.However,
westillknowlittleabouttheextenttowhichneutraloradaptive
ancestry on the transcriptional response of primary
inter-population genetic differences affect the actual immune
macrophages to live bacterial pathogens. A total
responsetopathogens.
of 9.3% of macrophage-expressed genes show
Addressingthisgapisnotonlyimportantforunderstanding
ancestry-associated differences in the gene regula-
recenthumanevolution,butmayalsohelprevealthemolecu-
toryresponsetoinfection,andAfricanancestryspe-
larbasisofancestry-relateddifferencesindiseasesusceptibil-
cifically predicts a stronger inflammatory response ity. Individuals from different populations vary considerably in
and reduced intracellular bacterial growth. A large their susceptibility to many infectious diseases, chronic in-
proportion of these differences are under genetic flammatory disorders, and autoimmune disorders. For tuber-
control:for804genes,morethan75%ofancestryef- culosis, systemic lupus erythematosus, systemic sclerosis,
fectsontheimmuneresponsecanbeexplainedbya psoriasis, and septicemia, African American (AA) and Euro-
single cis- or trans-acting expression quantitative pean American (EA) individuals exhibit an up to 3-fold differ-
ence in prevalence (reviewed in Brinkworth and Barreiro,
trait locus (eQTL). Finally, we show that genetic ef-
2014; Pennington et al., 2009; Richardus and Kunst, 2001).
fectsontheimmuneresponsearestronglyenriched
These observations argue in favor of significant ancestry-
for recent, population-specific signatures of adap-
relateddifferencesinimmuneresponse,especiallyinsuscep-
tation. Together, our results demonstrate how his-
tibilitytoinflammation(Penningtonetal.,2009;Richardusand
torical selective events continue to shape human
Kunst, 2001).
phenotypic diversity today, including for traits that Such differences almost certainly involve major contribu-
arekeytocontrollinginfection. tionsfromtheenvironment.However,genome-wideassociation
studies (GWAS) also support a keyrole for genetic factors, as
manyoftheGWAS-variantsassociatedwithinfectious,autoim-
INTRODUCTION mune,andinflammatorydiseasespresentextremedifferencesin
allele frequency (F > 0.4) between human populations, again
st
Asourprimaryinterfacewiththeenvironment,theimmunesys- supporting a possible history of population-specific selection
temisthoughttohaveevolvedunderstrongselectivepressure (BrinkworthandBarreiro,2014).
Cell167,657–669,October20,2016ª2016ElsevierInc. 657
GWASresultsalsoindicatethatsusceptibilitytomanycom- witheitherListeriaorSalmonella(Figure1A).Wefoundextensive
mon immune-related diseases is primarily controlled by non- differencesingeneexpressionlevelsbetweeninfectedandnon-
codingvariants(Gusevetal.,2014;Hindorffetal.,2009;Schaub infected cells, with 5,201 (44%) and 6,701 (56%) differentially
etal.,2012).Thus,manyancestry-relateddifferencesindisease expressed genes after infection with Listeria and Salmonella,
susceptibility may result from genetically controlled transcrip- respectively (see the STAR Methods, false discovery rate
tionaldifferencesinimmuneresponsestoinflammatorysignals. [FDR] < 0.01 and jlog (fold change)j > 0.5; Table S2A). As ex-
2
Thisideaisconsistentwithrecentexpressionquantitativetrait pected, the sets of genes that responded to either infection
locus (eQTL) mapping studiesin innateimmunecellsexposed were strongly enriched (FDR < 0.01) for gene sets involved in
to immune antigens or live infectious agents (Barreiro et al., immune function, including the regulation of inflammatory re-
2012; C¸alısxkan et al., 2015; Fairfax et al., 2014; Lee et al., sponses, cytokine production, T cell activation, and apoptosis
2014). Such immune ‘‘response eQTL’’ studies have identified (TableS3).
hundreds of genetic variants that both explain variation in the
host immune response and are significantly enriched among Ancestry-RelatedDifferencesintheInnateImmune
GWAS-associated loci. However, because studies to date ResponsetoInfection
have mostly focused on individuals of European ancestry, the We first aimed to characterize European versus African
degree to which such variants contribute to population differ- ancestry-relatedtranscriptionaldifferencesinnon-infectedand
encesintheimmuneresponseremainsunclear. infected macrophages. Because self-identified ethnicity is an
Here, we report an RNA-sequencing (RNA-seq)-based im- impreciseproxyfortheactualgeneticancestryofanindividual,
mune response eQTL study to test for the effects of African we used the genotype data to estimate genome-wide levels
versus European ancestry on the transcriptional response to ofEuropeanandAfricanancestryineachsampleusingthepro-
several live bacterial pathogens. We integrate statistical and gramADMIXTURE(Alexanderetal.,2009).Consistentwithpre-
evolutionary genetic analyses with primary macrophage gene viousreports(Brycetal.,2010;Tishkoffetal.,2009),wefound
expression levels, before and after infection, to characterize thatmanyself-identified AA individualshave ahigh proportion
ancestry-relateddifferencesintheimmuneresponse.Ourana- of European ancestry (mean = 30%, range 0.9%–100%; Fig-
lysesaddressthreefundamentalquestionsaboutrecentevolu- ure S1B). In contrast, self-identified EA showed more limited
tioninthehumanimmunesystem:(1)thedegreetowhichinnate levels of African admixture (mean = 0.4%, range 0%–18%;
immuneresponsesaredifferentiatedbyEuropeanversusAfrican Figure S1B). Thus, we used these continuous estimates (as
ancestry, (2) the genetic variants that account for such differ- opposedtoabinaryclassificationofindividualsintoAfricanor
ences,and(3)theevolutionarymechanisms(neutralgeneticdrift Europeanancestry)toidentifyancestry-associateddifferentially
versus positive selection) that led to their establishment in expressed genes (i.e., pop-DE genes: genes for which gene
modern human populations.Finally, to facilitate theuse of our expression levels are linearly correlated with ancestry levels;
databytheresearchcommunity,wehavedevelopedanacces- see the STAR Methods for details on the nested linear model
sible, publicly available browser for exploring our results: the usedforthisanalysis).
ImmunPopQTLbrowser(http://www.immunpop.com). Of the 11,914 genes we tested, we identified 3,563 pop-DE
genes(30%)inatleastoneoftheexperimentalconditions,ex-
RESULTS plaining a mean 8.2% of expression variance (range 1.8%–
44%)(FDR<0.05:1,745innon-infected[NI],1,336inListeria-
TranscriptionalResponseofMacrophagestoListeria infected[L],and2,417inSalmonella-infected[S]macrophages)
andSalmonella (Figures 1B and 1C; Table S2B). These differences pri-
Weinfectedmonocyte-derivedmacrophages—phagocyticcells marilyinfluencemeangeneexpressionlevelsacrosstranscript
that are essential for fighting foreign invaders, tissue develop- isoforms,asopposedtotheproportionofisoformusagewithin
ment,andhomeostasis(OkabeandMedzhitov,2016)—derived genes. Specifically, among genes with at least two annotated
from80AAand95EAindividuals(TableS1)witheitherListeria isoforms (n = 10,223), only 62, 39, and 48 genes exhibited
monocytogenes (a Gram-positive bacterium) or Salmonella evidence for ancestry-associated differential isoform usage,
typhimurium (a Gram-negative bacterium). Following 2 hr of in the non-infected, Listeria-infected, and Salmonella-infected
infection, we collected RNA-seq data from matched non-in- conditions, respectively (multivariate generalization of the
fectedandinfectedsamples,foratotalof525RNA-seqprofiles Welch’s t test; FDR < 0.05) (Figures 1D and S2A; Table S2D).
across individual-treatment combinations (mean = 36 million Theseresultswereunalteredbyusinganalternativeidentifica-
readspersample;seetheSTARMethods;FigureS1A).Eachin- tion approach (Wilcoxon rank sum test, as in Lappalainen
dividual was genotyped for over 4.6 million single nucleotide etal.,2013;seetheSTARMethodsfordetails)orwhenrelaxing
polymorphisms(SNPs),withadditionalimputationto(cid:1)13million the FDR threshold used to define significance (Figure S2B).
SNP genotypes (see the STAR Methods). After quality control Despitethelownumberofgenesshowingancestry-associated
(Figure S1A), we were able to study 171 individuals with high- differences in isoform usage, many of these genes are key
qualityRNA-seqdata,amongwhich168werealsosuccessfully regulators of innate immunity, including OAS1 that encodes
genotyped. isoformswithvaryingenzymaticactivityagainstviralinfections
Thefirstprincipalcomponentoftheresultinggeneexpression (Bonnevie-Nielsenetal.,2005).
dataaccountedfor85%ofthevarianceinourdatasetandsepa- Next we sought to identify genes for which the response to
rated non-infected macrophages from macrophages infected infection(i.e.,foldchangeingeneexpressionininfectedversus
658 Cell167,657–669,October20,2016
A B
C D
E G
F
Figure1. EuropeanandAfricanAncestry-AssociatedDifferencesinImmuneResponse
(A)Principalcomponentanalysisofgeneexpressiondatafromallsamples.PC1(xaxis)andPC2(yaxis)clearlyseparatenon-infectedmacrophagesfrom
macrophagesinfectedwitheitherListeriaorSalmonella.
(legendcontinuedonnextpage)
Cell167,657–669,October20,2016 659
non-infectedmacrophages,culturedinparallel)significantlycor- immune-relatedandtightlyconnectedtoadysregulatedinflam-
relates with ancestry (see the STAR Methods). We term these matory response. Further, among the diseases most signifi-
genes ‘‘population differentially responsive’’ (pop-DR) genes. cantly enriched for pop-DE genes were rheumatoid arthritis,
We detected 1,005 and 206 pop-DR genes (FDR < 0.05) in systemic sclerosis, and ulcerative colitis, all of which have
response to Salmonella and Listeria, respectively (Figure 1E; beenreportedtodifferinincidenceordiseaseseveritybetween
Table S2C) (the increased power for Salmonella likely results AAandEAindividuals(BrinkworthandBarreiro,2014;Penning-
from the stronger transcriptional response induced by Salmo- ton et al., 2009). Thus, ancestry-associated gene regulatory
nellarelativetoListeria,seeFigure1A).Thesegenesexplaina differences likely contribute to known ethnic disparities in
mean 7.4% (range 2.6%–24%) of variance in transcriptional inflammatory and autoimmune disease susceptibility, in part
response to infection. Overall, we found that macrophages throughaffectingtheabilityofmacrophagestocontrolbacterial
from individuals of African ancestry produced a markedly infections.
stronger transcriptional response to both bacterial infections
(Mann-Whitneytest,p<1310(cid:3)15,Figure1F).GOtermenrich- GeneExpressionQTLinNon-infectedandInfected
mentanalysesfurtherrevealedthatgenesrelatedtoinflamma- Macrophages
tory processes were the most enriched among pop-DR genes Toidentifywhetherpop-DEandpop-DRgenesareexplainedby
showingastrongerresponsetoinfectioninAfrican-descentindi- geneticdifferencesbetweenEuropeanandAfricanpopulations,
viduals (Figures 1G and S2C). Together, these results indicate wefirstmappedgeneticvariantsthatareassociatedwithgene
thatincreasedAfricanancestrypredictsastrongerinflammatory expressionlevels(i.e.,eQTL)ortranscriptisoformusage(alter-
responsetoinfection. nativesplicingQTL[asQTL])inthecompletesample.Todoso,
We hypothesized that ancestry-associated differences in weusedalinearregressionmodelthataccountsforpopulation
the transcriptional response to infection could translate into structure and principal components of the expression data,
ancestry-associated differences in the ability of macrophages thus limiting the effect of unknown confounding factors (see
tocleartheinfection.Wetestedthishypothesisinasubsetof the STAR Methods for details). Given that our sample size is
89individualsbyquantifyingthenumberofbacteriaremaining toosmalltorobustlydetecttrans-actingeQTL,wefocusedour
inside the macrophages rightafter the infectionstep (T0),2hr analyses on local associations that, for simplicity, we refer to
(T2),and24hr(T24)post-infection.Forbothbacteria,increased ascis-eQTL.Wedefinecis-eQTLandcis-asQTLhereasSNPs
Africanancestrypredictedimprovedcontrolofintracellularbac- locatedinthegenebodyorinthe100kbflankingthegeneof
terialgrowth.Thiseffectwasparticularlynoticeableinourinfec- interest.
tionexperimentswithListeria.Despitenosignificantdifference We identified cis-eQTL for 1,647 genes (14% of all genes
in the initial number of bacteria infecting macrophages (Fig- tested;FDR<0.01)inatleastoneoftheexperimentalconditions
ure 2A, p = 0.95), the number of bacteria inside the macro- (875innon-infectedmacrophages,1,087intheListeria-infected
phages of individuals with high levels of African ancestry at condition, and 983 in the Salmonella-infected condition; Fig-
T24was3.2-foldlowerthanthatofEuropeans(Figure2A,p= ure 3A; Table S4A; Figure S3A for number of eQTL found at
2.0310(cid:3)4). morerelaxedcutoffs).Similarly,wedetectedalargenumberof
Finally, we tested if pop-DE genes were enriched among cis-asQTL affecting the ratio of alternative isoforms used for
GWAS-associated genes. We found seven diseases for which the same gene (1,120 genes, 10% of all genes tested; FDR <
susceptibility genes reported by GWAS were significantly en- 0.01[Figure3A;TableS4C]:886innon-infectedmacrophages,
riched among genes classified as pop-DE, in at least one 746inListeria-infectedsamples,and615inSalmonella-infected
experimental condition (Figure 2B). Contributing to these en- samples).
richments are several HLA genes (HLA-DQA1, HLA-DPA1, Outofallgeneswithcis-eQTL,alargefraction(21.8%)were
HLA-DRB1, HLA-DPB1, HLA-DRA), known to be the main associated with an eQTL only in infected macrophages. In
genetic risk factors for several immune disorders. Strikingly, contrast, only 7.3% of genes showed an infection-specific
six of these seven diseases (all but Parkinson’s disease) are cis-asQTL (Figures 3A and 3B). Infection-specific cis-eQTL
(B)Venndiagramillustratingthenumberofpop-DEgenes(FDR<0.05)innon-infected(black),Listeria-infected(yellow),andSalmonella-infected(orange)
macrophages.
(C)Exampleofagene,thechemokineCCL15,forwhichexpressionlevelsinallconditionsaresignificantlyassociatedwithlevelsofAfricanversusEuropean
ancestry.Theaveragesequencingdepthforeachbase(normalizedpermillionmappedreads)isshownontheyaxis.
(D)Exampleofthreegenes(POLR1A,NDUFS5,andOAS1)forwhichancestrypredictsdifferencesinisoformusage.
(E)Exampleofthreeimmune-relatedpop-DRgenes.Theyaxisshowsthelog2foldchangesingeneexpressionlevelsinresponsetoListeriaandSalmonella,asa
functionofcontinuousdifferencesinAfricanancestry(xaxis).
(F)Absolutedifferenceinthelog2foldchangeresponsetoSalmonella(toppanel)andListeriainfection(bottompanel)betweenEuropeanandAfricanindividuals
(xaxis),amongallpop-DRgenes(red)andpop-DRgenesassociatedwiththeinflammatoryresponse(blue).Thenullexpectationfrompermutingadmixturelevels
acrossindividualsisshowninlightgrayforcomparison.AshiftinthedistributiontotherightreflectsastrongerresponsetoinfectioninAfrican-ancestry
individuals.
(G)GOenrichmentanalysisforgenesshowingasignificantinteractionbetweenancestryandtheresponsetoSalmonella.OnlyGOtermswithanenrichmentat
FDR<0.1aredisplayed,andGOtermsarecolor-codedintofunctionallyrelatedtermsbasedontheoverlapamonggenesets(Bindeaetal.,2009).ForeachGO
term,theaverageinteractioneffectisplottedonthexaxisandthemeanlog2foldchangeingeneexpressionlevelsinresponsetoinfectionisplottedontheyaxis.
SeealsoTablesS2andS3.
660 Cell167,657–669,October20,2016
Listeria Salmonella Figure 2. Increased African Ancestry Pre-
A dictsImprovedControlofBacterialGrowth
3
2 insideMacrophages
d ratio 1 d ratio 12 (n(Ayu)amxBbisoe)xrimpolfomtbesadcisatheteroilawyiiannfgsteidrtehineifnefqcetucioatenndt(ilTme0a),nc2ororhmprhapaliogzseetds-
e 0 e
z z 0 infection (T2), and 24 hr post-infection (T24)
mali−1 mali−1 (xaxis).Wequantilenormalizeddataacrossin-
r r dividuals and time points. Analyses were con-
o o
N−2 N−2 ductedusingacontinuousmeasureofancestry;
P = 0.95 P = 0.05 P = 2e-4 P = 0.69 P = 0.16 P = 1.1e-2 however,forvisualizationpurposes,Africanand
−3 −3 Europeansamplesweredefinedasthosewithan
T0 T2 T24 T0 T2 T24 estimated African ancestry >75% (green) and
Timepoints Timepoints <25%(pink),respectively.
(B) Fold enrichment of pop-DE genes (y axis)
European-American African-American amonggenesidentifiedindiseasesusceptibility
GWAS,atprogressivelystringentpvaluethresh-
B olds(xaxis).Grey,yellowandorangelinesshow
25 FDR < 0.1 significantenrichments(filledcirclesatanFDR<
0.1 and open circles at anominal p <0.05) for
Pvalue < 0.05
pop-DE genes identified in non-infected and
Listeria- and Salmonella-infected macrophages,
20 rheumatoid arthrittiiss
respectively.Lightgraylinesshownon-significant
nt multiple
e sclerosis diseases/traits.
m SeealsoFigureS2.
h 15
nric ssycsleteromsiics ulcecraotliivties
E
old 10
F Parkinson's disease tle interaction effects: eQTL can be
Parkinson's sharedacrossconditionsaslongastheir
disease
5 effect size differs between infected and
non-infectedsamples.Wedetected244
Inflammatory bowel disease and 503 genes with a cis-reQTL (FDR <
0 0.01,TableS4B)fortheresponsetoLis-
teriaandSalmonella,respectively.Inter-
10 15 20 40
estingly,amonggenesassociatedwitha
GWAS -log (P)
10 cis-reQTL, we found several key regu-
latorsoftheimmuneresponse,including
the transcription factors STAT4 and
werefurthersupportedbyanalysisofallele-specificexpression IRF2 (Figure 3D). We also found cis-reQTL for known suscep-
(ASE) levels, which provides independent but complementary tibility loci for ulcerative colitis (e.g., HLA-A, HLA-DQA2,
evidence for functional cis-regulatory variation. As expected, PMPCA), systemic lupus erythematosus (ITGAX, HLA-DQA1),
geneswithcis-eQTLweresignificantlyenrichedforgeneswith and the infectious diseases hepatitis B and leprosy (e.g.,
ASE,comparedtothebackgroundofall9,588genestested(Fig- HLA-C,NOD2).
ure S3B, Fisher’s exact test, p < 13 10(cid:3)15 for all conditions). To investigate the regulatory mechanisms that account for
Further, genes harboring infection-specific eQTL also tended immune reQTL, we next profiled the genome-wide chromatin
toexhibitinfection-specificASEinthesamecondition(Listeria accessibility landscape of non-infected and Listeria and Sal-
or Salmonella) in which the eQTL was identified (Figure 3C, monella-infected cells using assay for transposase-accessible
(cid:1)27 fold-enrichment of infection-specific ASE among infec- chromatin using sequencing (ATAC-seq) (Buenrostro et al.,
tion-specific eQTL, relative to shared eQTL; p < 1 3 10(cid:3)15). 2013).Thisapproachallowedustoidentifytranscriptionfactor
Thus, in agreement with previous studies (Fairfax et al., 2014; (TF) binding motifs likely to be occupied by their respective
Lee et al., 2014), genotype-environment (G 3 E) interactions TFs,inbothconditions(seetheSTARMethods).Wefoundthat
are common in the context of immune responses to infection, SNPswithinaccessibleTFbindingsitesweregreaterthanfour
albeitmoresoformeanexpressionlevelsthanalternativeiso- timesmorelikelytobeidentifiedasreQTL(Figure3E).Further,
formusage. reQTL in our analyses were strongly enriched (>20-fold) for
AcomplementaryapproachtoidentifyingG3Einteractions PU.1bindingsites(apioneerTFinvolvedinregulatingenhancer
forexpressionlevelsistodirectlymapresponseeQTL(reQTL): activityinmacrophages)(Garberetal.,2012)andforvirtuallyall
QTL associated with the magnitude of change in expression TFsthatorchestrateinnateimmuneresponsestoinfection(Fig-
levelsafterinfection(Barreiroetal.,2012;C¸alısxkanetal.,2015; ure3E)(e.g.,nuclearfactorkB[NF-kB]:>50-fold;AP1:>55-fold;
Lee et al., 2014). In contrast to condition-specific eQTL (an and IRFs: 14-fold for Salmonella only). In striking contrast, we
extremecaseofG3Einteraction),reQTLcancapturemoresub- found no such enrichment for eQTL identified in non-infected
Cell167,657–669,October20,2016 661
A C
B
D
E
Figure3. eQTLandASEAnalysesRevealExtensivecis-RegulationofGeneExpressionResponsestoPathogensinMacrophages
(A)Schematicrepresentationofthenumberofcis-eQTLandcis-asQTLsharedacrossallconditions,oronlyfoundinnon-infectedmacrophagesorListeriaand/or
Salmonellainfectedmacrophages.Infection-specificeQTLweredefinedasthoseshowingverystrongevidenceofeQTLintheinfectedstate(FDRalwayslower
than0.01),andverylimitedinthenon-infectedstate(FDRalwayshigherthan0.3).
(B)Examplesofacis-eQTLobservedinallconditions(HLA-DQB1),acis-eQTLobservedonlyininfectedmacrophages(ADSS),andacis-eQTLobservedonlyin
Salmonella-infectedmacrophages(HLA-C).PinkandgreendotsinsidetheboxplotsdistinguishAfrican(>75%Africanancestry)andEuropean(<25%African
ancestry)samples,respectively.
(C)PlotcontrastingtheevidenceforASE(-log10pvalues)innon-infectedmacrophages(yaxis)andinmacrophagesinfectedwithSalmonella(xaxis),forgenes
whereweidentifiedcis-eQTLinbothconditions(purple),genesforwhichcis-eQTLwereonlyfoundinnon-infectedmacrophages(gray),andgenesforwhichcis-
eQTLwereonlyfoundinSalmonella-infectedmacrophages(orange).Qualitativelysimilarresultswereobtainedwhencontrastingnon-infectedandListeria-
infectedcells(FigureS3C).
(D)Examplesoftwocis-reQTLwheregenotype(xaxis)hasasignificanteffectontheresponseofSTAT4(left)andIRF2(right)toSalmonellainfection.
(E)reQTLenrichments(xaxis)inactivelyregulatedTFbindingsitesannotatedbyATAC-seqfootprinting.Errorbarsshow95%confidenceintervals.Bindingsites
weregroupedintofunctionallyoverlapping‘‘TFclusters’’usingsequencesimilarityandco-localizationinthegenome(TableS6;STARMethods).
SeealsoTableS4.
macrophages(p>0.05forNF-kB,AP1,andIRFs)(FigureS3D). GeneticBasisofAncestry-AssociatedDifferencesinthe
TheseresultsshowthatreQTLvariantsareoftenconditionallysi- ImmuneResponsetoPathogens
lent in resting macrophages but become functionally relevant Wehypothesizedthatdifferencesinallelefrequenciesforsome
post-infection, and this transition is explained by disruption of of the eQTL identified above could explain the observed
bindingsitesforimmuneresponse-activatedTFs. ancestry-associateddifferencesinthetranscriptionalresponse
662 Cell167,657–669,October20,2016
A B
C D
Figure4. ContributionofcisandtransGeneticVariationtopop-DEandpop-DRGenes
(A)Theproportionofpop-DE,pop-DR,andgenesthatexhibitancestry-associatedisoformusagethatareassociatedwithacis-eQTL,cis-reQTL,orcis-asQTL,
respectively(FDR<0.01).Nullexpectations(basedonthegenome-wideproportionofgenesassociatedwitheachQTLclass)areshowningray.Similarresults
areobtainedwhenfocusingontranscriptionalQTLidentifiedatanFDRof0.05(FigureS4B).
(B)AverageDPSTobtained(±SE)whenregressingoutthegenotypeeffectoftheleadcis-ortheleadtrans-SNPforpop-DEandpop-DRgenes(yaxis),defined
usingprogressivelystringentFDRcutoffs(xaxis).ColoredlinesshowaverageDPSTvaluesbasedontherealdata;graylinesshowthesamevalueswhenre-
gressingoutthegenotypeeffectoftheleadSNPidentifiedbasedonpermutedgenotypes.
(C)Numberofgenesidentifiedaspop-DEandpop-DRatanFDR<0.05(yaxis)before(dashedbars)andafter(filledbars)regressingouttheeffectoftheleadcis-
SNPassociatedwiththesegenes.
(D)Examplesofgenesforwhichtheleadcis-(blue)ortheleadtrans-SNP(green),explainsatleast75%ofthedifferencesingeneexpressionassociatedwith
AfricanversusEuropeanancestry.
to infection. In support of this hypothesis, we found that pop- proportionofoverallgeneexpressionvarianceexplainedbybe-
DEgeneswereenrichedupto3.3-foldforgeneswithcis-eQTL tween-populationphenotypicdivergence(asopposedtowithin-
(p < 1 3 10(cid:3)10), and pop-DR genes were enriched up to 5.8- populationdiversity).P valuesrangefrom0to1,withvalues
ST
fold for genes with cis-reQTL (p < 1 3 10(cid:3)10) (Figures 4A and closeto1implyingthatthemajorityofagene’sexpressionvari-
S4A).Additionally,(cid:1)60%ofgenesthatexhibitedancestry-asso- ance is due to differences between populations. Our score,
ciated isoform usage were associated with an asQTL (up to deltaP (DP ),isdefinedasthedifferencebetweenP values
ST ST ST
24-fold enrichment, p < 1 3 10(cid:3)10). Thus, although rare, beforeandafterregressingouttheeffectofthecis-SNPthatwas
ancestry-associatedchangesinisoformusagearelargelygenet- moststronglyassociatedwiththetargetgene’sexpressionlevel
icallydriven. (regardless of significance level), divided by the P value
ST
To explicitly quantify the contribution of our eQTL set to observedbeforeremovingthegenotypeeffect.DP therefore
ST
transcriptional differences detected between populations, we quantifiestheproportionofancestry-associatedexpressionlevel
devised a new score based on P estimates (Leinonen etal., differencesthatstemfromthestrongestcis-associatedvariant.
ST
2013; Pujol et al., 2008). P is the phenotypic analog of the Among all pop-DE genes, we found that cis-regulatory
ST
population genetic parameter F , providing a measure of the variants explained an average of 31%, 31%, and 26% of
ST
Cell167,657–669,October20,2016 663
ancestry-related differences in expression observed in non- (MAF) >5% in the CEU and YRI samples. In contrast to F ,
ST
infected, Listeria-infected and Salmonella-infected samples, iHSisawithin-populationmeasureofrecentpositiveselection
respectively(Figure4B).Further,thelargertheeffectofancestry that is not affected by the levels of population differentiation
in the original pop-DE analysis, the larger the contribution of (Voightetal.,2006).
cis-regulatoryvariationtothesedifferences:forpop-DEgenes Our analyses identified significantly higher mean F values
ST
identifiedatastringentFDRof1310(cid:3)4,cis-regulatoryvariation among genes that were pop-DE, pop-DR, or showing differ-
explained close to 50% (on average) of ancestry effects (Fig- ences in isoform usage between populations (p % 1 3 10(cid:3)3;
ure4B).Weobservedasimilarpatternforpop-DRgenesafter Figures5AandS5Aforsimilarresultswhenusingalternativewin-
regressing out the genotype effect of the lead cis-reQTL SNP dowsizes).Further,variantsidentifiedascis-eQTLweresignifi-
(Figure 4B). In support of the substantial role of cis-regulatory cantly enriched ((cid:1)2-fold) for high iHS values (i.e., iHS > 99th
variationinexplainingpop-DEandpop-DRgenes,geneexpres- percentile of genome-wide distribution, Figure 5B, p < 1 3
sion values for 30% and 45% of pop-DE and pop-DR genes, 10(cid:3)8),consistentwiththeimportanceofregulatorygeneticvari-
respectively, were no longer significantly associated with ation in recent human evolution (Fraser, 2013). cis-reQTL and
ancestryonceweregressedoutcis-geneticeffects(Figure4C). cis-asQTL were even more strongly enriched among high iHS
Importantly, DP values never exceeded 5% when we re- values(upto3.6-fold;Figure5B,p<1310(cid:3)5).
ST
gressedouteither(1)thegenotypeeffectofrandomlyselected Overall,withinthesetofcis-eQTL-,cis-reQTL-,orcis-asQTL-
SNPsmatchedfortheallelefrequencyoftheleadcis-SNP,or associated genes, 258 carried a signature of recent positive
(2)leadcis-SNPsidentifiedafterpermutingthegenotypedata. selectionineitherCEUorYRIsamples(jiHSjR99thpercentile
Thus, our results cannot be simply explained by population ofthegenome-widedistribution)(Figure5C;TableS5A).These
structure(Figure4B). variantswerealsosignificantlyenrichedforhighXP-EHHvalues
Based on their known importance in the genetic control of (Sabetietal.,2007)((cid:1)6-fold,p<1310(cid:3)10,FigureS5C),further
generegulationandbecauseofpowerlimitations,ourmainanal- supporting that these variants have been important in recent,
ysis of ancestry-associated gene expression patterns focused population-specific human adaptation. However, because
ontheroleofcis-eQTL.However,inaseparateanalysis,were- outlier methods for detecting selection can be susceptible
calculatedDP usingtheleadtrans-SNPforeachgeneinplace to false positives (Kelley et al., 2006), we complemented our
ST
of the lead cis-SNP (although only 51, 21, and 22 trans- iHS analysis with a model-based approach. Specifically, we
eQTL genes survived genome-wide multiple testing correction comparedtheobservediHSvalueforeachputativelyselected
(FDR<0.1)innon-infected,Listeria-infectedandSalmonella-in- allele to those observed under neutral coalescent simulations
fected samples, respectively). Intriguingly, we found that lead matchedtotheknowndemographichistoryofAfricanandEuro-
trans-SNPs accounted for an average of (cid:1)23% and (cid:1)20% of peanpopulations(Gutenkunstetal.,2009),thecandidateallele’s
ancestry effects on gene expression levels for pop-DE and observedfrequency,andthelocalrecombinationrate.Thevast
pop-DR genes, respectively (Figure 4B; at least 2-fold more majority (92%) of all sites tested exhibited significantly larger
thanestimatesbasedonpermuteddata,p<1310(cid:3)10).These observed iHS statistics than expected under a neutral model
resultssuggestthatleadtrans-SNPs,althoughdifficulttodetect (p<0.01,TableS5B),providingstrongconvergentsupportfor
atagenome-widesignificancelevel,areenrichedfortruetrans- recentpositiveselectionattheseloci.Farmoreofthesegenes
associations that could be resolved with larger sample sizes. are pop-DE or pop-DR than expected by chance (47% and
Together, a single cis- or trans-acting variant was sufficient to 23%,respectively:Figure 5D,p <0.001),showing thatnatural
explainalmostallancestryeffects(DP >75%)ongeneexpres- selectionhascontributedtopresent-dayinter-populationdiffer-
ST
sionlevelsfor804pop-DEgenesandpop-DRgenes(Figure4D), encesininnateimmuneresponsestoinfection.
includingformasterregulatorsoftheimmuneresponsesuchas Neanderthalancestrymakesup(cid:1)2%oftheancestryofliving
CASP1,STAT4,andMICA.Ourresultsthusprovideacompre- humans found outside of Africa (Kelso and Pru¨fer, 2014). It
hensive genome-wide map of cis- and trans-genetic variants is therefore plausible that interbreeding between Neanderthal
associated with African and European ancestry-related differ- and modern human populations could also contribute to
encesintheimmuneresponsetoinfection. some of the ancestry-related differences in gene expression
we observed, especially if it enabled the ancestors of modern
NaturalSelectionandGeneticAncestryEffectsonGene Europeans to more rapidly adapt to a new pathogen environ-
ExpressionDivergence ment(Se´gurelandQuintana-Murci,2014).Totestthishypothe-
Finally, we sought to determine the impact of recent local sis,weidentifiedsiteswherethederivedalleleissharedbetween
positive selection ineither African orEuropean populations on Neanderthalsandnon-Africanpopulations,butisabsentinsub-
ancestry-related divergence in gene expression levels. To do SaharanAfricanssamplesconsidered.Thisclassofsites,which
so, we first calculated F values between the Yoruba African wecall‘‘Neanderthal-likesites’’(NLS),isaconservativeindicator
ST
(YRI) and the western European population (CEU) in Phase 3 of Neanderthal introgression (Sankararaman et al., 2014).
data from the 1000 Genomes Project (Auton et al., 2015). To Among the 18,862 NLS tested in our cis-QTL analyses, 297
generate gene-specific estimates, we averaged F values for were significantly associated with transcriptional variation of
ST
variantswithinawindowof10kbaroundthetranscriptionstart 145 genes (NLS-QTL). Among these 145 genes, 46% (FDR <
site (TSS) of each gene we analyzed (11,914 genes). As a 0.05) were differentially regulated in at least one experimental
complementaryapproach,wealsocalculatedintegratedhaplo- condition (non-infected, Listeria-infected, Salmonella-infected,
type scores (iHS) for all SNPs with a minor allele frequency orintheresponsetoeithertypeofinfection)betweenEuropeans
664 Cell167,657–669,October20,2016
A B C
D
Figure5. NaturalSelectiononeQTLandItsContributiontoAncestry-AssociatedRegulatoryDifferences
(A)MeanF valuesinawindowof±10kbaroundtheTSSofallgenes,pop-DEgenes,pop-DRgenes,andgenesshowingdifferencesinisoformusagebetween
ST
populations(toppanel).
(B)ProportionofallSNPs,cis-eQTL,cis-reQTL,andcis-asQTLwithaniHSvalueabovethe99thpercentileofthegenome-widedistributionintheCEU(jiHSj>
2.70)andtheYRI(jiHSj>2.68)populations.SeeFigureS5BforsimilarresultsconsideringQTLidentifiedatanFDR<0.05(insteadof0.01).
(C)Manhattanplotofagenome-widescanforselectioninCEU(top)andYRI(bottom)forSNPsidentifiedasregulatoryQTLinmacrophages.Thedashedline
representsthe99thpercentileofthegenome-widedistribution.DarkershadesofbluerepresentlargerFSTvaluesforSNPswithelevatedjiHSjvalues;bluecircled
dotshighlightgenesthatshowoneormoretranscriptionalassociationswithAfricanversusEuropeanancestry.GenesinredareregulatedbyNLSwithelevated
jiHSjvaluesinCEU(jiHSj>2.7),supportingadaptiveintrogressionfromNeanderthalsintotheancestorsofmodernEuropeans.
(D)ProportionofgenesregulatedbyeQTLandtargetedbyrecentpositiveselection(amongthe258representedbythebluecirclesinC)thatarepop-DE,
pop-DR,orshowpopulationdifferencesinisoformusage(bluetriangles),comparedtorandomexpectationswhensamplingthesametotalnumberofgenes
10,000timesfromallgenestested(violinplots).
SeealsoTableS5.
and Africans (63% at a more relaxed FDR < 0.1). Thus, a ofLCL lines, which were generated more than 20 years apart
non-negligible proportion of ancestry-related gene expression (Daussetet al., 1990).
divergence probably results from introgression of functional Oneofthemoststrikingobservationsfromourstudywasthe
Neanderthalvariantsintotheancestors ofmodernEuropeans. markedly stronger response to infection induced in macro-
Interestingly,someofthesevariants(n=16)alsohaveelevated phagesfromindividualsofAfricandescent,particularlyamong
iHS values (jiHSj R 2) (Figure 5C; Table S5A) and therefore inflammatoryresponsegenes.Thisresultagreeswithprevious
representnewcandidatesforadaptiveintrogressioninhumans. reports showing that AAs have higher frequencies of alleles
associatedwithanincreasedpro-inflammatoryresponse(Ness
DISCUSSION et al., 2004), increased levels of circulating C-reactive protein
(Kelley-Hedgepethetal.,2008),andamuchhigherrateofinflam-
Together, our results provide a comprehensive characteriza- matory diseases than EA individuals (Pennington et al., 2009).
tion of genes for which the transcriptional responses of pri- Althoughtheexactcausallinkbetweenancestryandthepro-in-
mary cells to live pathogenic bacteria differs depending on flammatory response has yet to be established, we speculate
European versus African ancestry. We show that 34% of thatthestrongerinflammatoryresponseassociatedwithAfrican
genes expressed in macrophages show at least one type of ancestry accounts for the increased ability of macrophages in
ancestry-related transcriptional divergence, whether in the African ancestry individuals to control bacterial growth post-
form of differences in gene expression (30%), the transcrip- infection.
tional response to infection (9.3%), or, less commonly, differ- Nevertheless, the evolutionary pressures that explain these
encesinisoformusage(1%).Notably,themodestcontribution differencesremainanopenquestion.Onepossibilityisthat,after
ofdifferencesinisoformusagetoancestry-relatedexpression humanpopulationsmigratedoutofAfrica,theywereexposedto
levelsdiffersfrompreviousresultsinlymphoblastoidcelllines lowerpathogenlevels(Guernieretal.,2004),whichreducedthe
(LCLs),wheretheywerefoundtobequitecommon(Lappalai- needforstrong,costlypro-inflammatorysignals.Changeinthis
nen et al., 2013). The discrepancy between our results and directionmayhavebeenfavoredduetothedetrimentalconse-
those reported for LCLs may be related to differences in quencesofacuteorchronicinflammation,whicharekeycontrib-
the experimental procedures used to produce the two sets utorstothedevelopmentofautoinflammatoryandautoimmune
Cell167,657–669,October20,2016 665