Sensitive Detection of Chromosomal Segments of Distinct Ancestry in Admixed Populations

18 Pages·2009·2.28 MB·English

Sensitive Detection of Chromosomal Segments of
Distinct Ancestry in Admixed Populations
Alkes L. Price1,2,3, Arti Tandon3,4, Nick Patterson3, Kathleen C. Barnes5, Nicholas Rafaels5, Ingo
Ruczinski6, Terri H. Beaty6, Rasika Mathias7, David Reich3,4*, Simon Myers3,8,9*
1DepartmentofEpidemiology,HarvardSchoolofPublicHealth,Boston,Massachusetts,UnitedStatesofAmerica,2DepartmentofBiostatistics,HarvardSchoolofPublic
Health,Boston,Massachusetts,UnitedStatesofAmerica,3BroadInstituteofHarvardandMIT,Cambridge,Massachusetts,UnitedStatesofAmerica,4Departmentof
Genetics, Harvard Medical School, Boston, Massachusetts, United States of America, 5Johns Hopkins Allergy and Asthma Center, Division of Clinical Immunology,
DepartmentofMedicine,SchoolofMedicine,Baltimore,Maryland,United StatesofAmerica,6DepartmentofBiostatistics,JohnsHopkinsSchool ofPublicHealth,
Baltimore,Maryland,UnitedStatesofAmerica,7InheritedDiseaseResearchBranch,NationalHumanGenomeResearchInstitute,NationalInstitutesofHealth,Baltimore,
Maryland,UnitedStatesofAmerica,8DepartmentofStatistics,OxfordUniversity,Oxford,UnitedKingdom,9WellcomeTrustCentreforHumanGenetics,Universityof
Oxford,Oxford,UnitedKingdom
Abstract
Identifying the ancestry of chromosomal segments of distinct ancestry has a wide range of applications from disease
mappingtolearningabouthistory.Mostmethodsrequiretheuseofunlinkedmarkers;but,usingallmarkersfromgenome-
wide scanning arrays, it should in principle be possible to infer the ancestry of even very small segments with exquisite
accuracy. We describe a method, HAPMIX, which employs an explicit population genetic model to perform such local
ancestryinferencebasedonfine-scalevariationdata.WeshowthatHAPMIXoutperformsothermethods,andweexploreits
utilityforinferringancestry,learningaboutancestralpopulations,andinferringdatesofadmixture.Wevalidatethemethod
empiricallybyapplyingittopopulationsthathaveexperiencedrecentandancientadmixture:935AfricanAmericansfrom
the United States and 29 Mozabites from North Africa. HAPMIX will be of particular utility for mapping disease genes in
recently admixed populations, as its accurate estimates of local ancestry permit admixture and case-control association
signalsto becombined, enablingmore powerful testsofassociationthan with eithersignal alone.
Citation:PriceAL,TandonA,PattersonN,BarnesKC,RafaelsN,etal.(2009)SensitiveDetectionofChromosomalSegmentsofDistinctAncestryinAdmixed
Populations.PLoSGenet5(6):e1000519.doi:10.1371/journal.pgen.1000519
Editor:JonathanK.Pritchard,UniversityofChicago,UnitedStatesofAmerica
ReceivedOctober27,2008;AcceptedMay15,2009;PublishedJune19,2009
Thisisanopen-accessarticledistributedunderthetermsoftheCreativeCommonsPublicDomaindeclarationwhichstipulatesthat,onceplacedinthepublic
domain,thisworkmaybefreelyreproduced,distributed,transmitted,modified,builtupon,orotherwiseusedbyanyoneforanylawfulpurpose.
Funding:ThisworkwassupportedbyNIHgrantsU01-HG004168toALP,AT,NP,andDRandbyR01-HL087699(toKCB,THB,IR,RM,andNR).ALPwassupported
byaRuthKirschsteinNationalResearchServiceAwardfromtheNationalInstitutesofHealth;RMwassupportedbytheIntramuralResearchProgramofthe
NationalHumanGenomeResearchInstitute,NationalInstitutesofHealth;KCBwassupportedinpartbytheMaryBerylPatchTurnbullScholarProgram;andDR
wassupportedbyaBurroughsWellcomeCareerDevelopmentAwardintheBiomedicalSciences.Thefundershadnoroleinstudydesign,datacollectionand
analysis,decisiontopublish,orpreparationofthemanuscript.
CompetingInterests:Theauthorshavedeclaredthatnocompetinginterestsexist.
*E-mail:[email protected](DR);[email protected](SM)
Introduction ancestry, leading to false-positive inferences of a deviation in
ancestry at certain loci [13].
The identification of chromosomal segments of distinct Here,wedescribeahaplotype-basedmethod,HAPMIX,which
continental ancestry in admixed populations is an important applies an extension of the population genetic model of Li and
problem,withawiderangeofapplicationsfromdiseasemapping Stephens [14] to the problem of local ancestry inference in
to understanding human history. Early efforts to solve this populationsformedbytwowayadmixture.Weapplythemethod
problem used coarse sets of unlinked markers [1–3] and mostly to simulated mixtures of African and European chromosomes to
focused on populations such as African Americans [4,5] and show that the resulting local ancestry inference is exceedingly
Latinos [6–8] that admixed within the past approximately 10 accurate in comparison to other methods, even in the case of
generations. Applying this approach to more anciently admixed ancient admixture in which the shorter ancestry segments are
populationshasledtoancestrypredictionsthatareambiguousat more difficult to infer. As expected from its use of an explicit
many loci [9]. However, methods based on coarse sets of populationgeneticmodel,HAPMIXmakesmorecompleteuseof
markers do not take advantage of the much richer haplotype dense genome-wide data, producing more accurate results. We
information available in genome-wide data. More recent examinethesensitivityoflocalancestryinferencetoawidearray
methods have been designed to use data from genome-wide of factors. We also explore the utility of HAPMIX for drawing
scanning arrays [10–12], but these methods do not fully model inferences about both the ancestral populations and the date of
linkage disequilibrium (LD) in the ancestral populations. Thus, admixture.
they do not capture all of the available information about We apply HAPMIX to 935 African American individuals
ancestry, and can be far from optimal. Furthermore, unless a genotyped at ,650,000 markers. By studying a large set of
trimming step is applied to remove linked markers [11], individuals from an admixed population of high relevance to
unmodeled LD may cause systematic biases in estimated diseasemapping,wevalidatetheeffectivenessofthismethodina
PLoSGenetics | www.plosgenetics.org 1 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
admixed individual is a better statistical match to one reference
Author Summary
populationortheother.AHiddenMarkovModel(HMM)isused
The genomes of individuals from admixed populations to combine these likelihoods with information from neighboring
consistofchromosomalsegmentsofdistinctancestry.For loci, to provide a probabilistic estimate of ancestry at each locus.
example, the genomes of African American individuals The method allows transition at two scales. The small-scale
contain segmentsof bothAfrican and Europeanancestry, transitions are between haplotypes from within a reference
sothataspecificlocationinthegenomemayinherit0,1, population, typically at a scale of every few tens of thousands of
or 2 copies of European ancestry. Inferring an individual’s bases [14]. The large-scale transitions are between the reference
local ancestry, their number of copies of each ancestry at
populations, at a scale of up to tens of millions of bases for a
eachlocationinthegenome,hasimportantapplicationsin
recentlyadmixedpopulationsuchasAfricanAmericans.Figure1
disease mapping and in understanding human history.
illustrates themethodschematically.
Here we describe HAPMIX, a method that analyzes data
An important strength of HAPMIX is the way it analyzes
from dense genotyping chips to infer local ancestry with
diploid data from admixed individuals. A na¨ıve way to use
very high precision. An important feature of HAPMIX is
population genetic methods to infer ancestry would be to pre-
that it makes use of data from haplotypes (blocks of
nearby markers), which are more informative for ancestry processsuchadatasetusingphasingsoftware,andthentoassume
than individual markers. Our simulations demonstrate the that this guess about the underlying phased haplotype is correct.
utilityofHAPMIXforlocalancestryinference,andempirical However, phase switch errors that arise from this procedure
applications to African American and Mozabite data sets (whicharecommonevenwiththebestphasingalgorithms[17,18])
uncover important aspects of the history of these would inappropriately force the method to infer ancestry
populations. transitions. HAPMIX circumvents this problem by not assuming
thatanyonehaplotypephasesolutioniscorrect.Instead,itusesa
built-inphasingalgorithm,similar tothatof[17],whichallowsit
practical setting and specifically show that the ancestry estimates to average inferences about ancestry over all possible phase
arenotsystematicallybiasedwithinthelimitsofourresolution.To solutions within each admixed individual. We treat the reference
illustrate how the method can provide insights into the history of populations as fully phased, partly because in some cases, e.g.
an anciently admixed population, we also apply HAPMIX to a AfricanandEuropeanchromosomesfromHapMap,thisphasing
data set of 29 individuals from the Mozabite population of uses unambiguous trio information and is therefore highly
northernAfricathatweregenotypedat,650,000markersaspart accurate.Moreimportantly,weexpectourapproachtoberobust
of the Human Genome Diversity Panel (HGDP) [15]. We show toerrorsinphasinginthereferencepopulations,becausetheseare
that the Mozabite have inherited roughly 78% ancestry from a unlikely to force inappropriate ancestry switches, in contrast to
European-relatedpopulationand22%ancestryfromapopulation phasing errorsin theadmixeddata itself.
relatedtosub-SaharanAfricans.Ouranalysisalsoshowsthatthe HAPMIX is also notable in inferring probabilities for whether
Mozabiteadmixturehasoccurredoveraperiodthatbeganatleast anindividualhas0,1,or2allelesofaparticularancestryateach
100 generations ago (,2,800 years ago), and that has continued locus.Asoursimulationsshow,theseestimatesarewell-calibrated.
intothepresentday.Weareabletoinfersmall,ancient,ancestry Thus,whenthemethodgeneratesaprobabilitypforanindividual
segmentsintheMozabite,andwedemonstratethatthesegments being heterozygous for ancestry at a locus, they are in fact
show considerable drift relative to all the other HGDP heterozygous approximately this proportion of the time. A well-
populations,consistentwiththehistoricalisolationoftheMozabite calibrated probability of ancestry at each locus isimportant for a
population. variety of applications, and also allows us to evaluate the
robustness oftheresults.
Materials and Methods HAPMIX is fundamentally different from existing methods
suchasANCESTRYMAPandLAMP[1,11].ANCESTRYMAP
Ethics statement
applies a Hidden Markov Model to unlinked SNPs to model
FortheAfricanAmericandata,informedconsentwasobtained ancestry transitions, while LAMP computes a majority vote of
fromeachstudyparticipant,andthestudyprotocolwasapproved ancestryinformationusingwindowsofunlinkedSNPs,butneither
by the institutional review board at either the Johns Hopkins of those methods makes use of haplotype information. Another
University orHoward University. method for investigating admixture segments, HAPAA, has
recently been published [19]. In common with HAPMIX, the
Overview of haplotype-based inference of local ancestry HAPAAsoftwareusesaHiddenMarkovModeltomodellinkage
HAPMIXassumesthattheadmixedpopulationbeinganalyzed disequilibrium within populations, and infers ancestry segments.
has arisen from the admixture of two ancestral populations, and However, there are also a number of important differences
that phased data are available from unadmixed reference between our model and that used by HAPAA. First, unlike
populations that are closely related to the true ancestral HAPAA, we allow for some rate of miscopying of ancestry
populations (e.g. phased data from HapMap [16]). In theory, segments from the ‘‘wrong’’ population, which we have found
discrepancies between the reference populations and the true greatly improves our ancestry estimation (instead of this, the
ancestral populations may lead to inaccuracies, but in practice HAPAAsoftwareusesapost-hoc‘‘filtering’’ofinferredsegments,
HAPMIX is robust to this concern under a variety of realistic which removes all segments of size below a certain minimum
scenarios (seebelow). threshold).Second,wefullyallowforunphaseddatainourmodel,
The central idea of the method is to view haplotypes of each whiletheHAPAAapproachrequiresapriorphasingofthedata,
admixed individual as being sampled from the reference andthenattemptstoaccountfortheeffectofphase-fliperrorson
populations: for example, haplotypes of an African American ancestryinferenceviaaheuristicprocedure.Webelievethatthese
individual could be sampled from phased African and European features of HAPMIX are likely to be critical in unraveling older
chromosomes from HapMap. At each position in the genome, admixture events, where ancestry segments are much shorter. A
HAPMIX estimates the likelihood that a haplotype from an finaladvantageofHAPMIXoverHAPAAisthatitisdesignedto
PLoSGenetics | www.plosgenetics.org 2 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
Figure1.SchematicoftheMarkovmodelweuseforancestryinference.Theblacklowerlinerepresentsachromosomalsegmentfroman
admixedindividual,carryinganumberoftypedmutations(blackcircles).Theunderlyingancestryisshowninthebottomcolorbar,andrevealsan
ancestrychangefromthefirstpopulation(red)tothesecondpopulation(blue).Theadmixedchromosomeismodeledasamosaicofsegmentsof
DNA from two sets of individuals drawn from different reference populations (red and blue horizontal lines respectively) closely related to the
progenitorpopulationsfortheadmixtureevent.Theyellowlineshowshowtheadmixedchromosomeisconstructedintermsofthismosaic.The
dottedlineabovethebottomcolorbarshowsthereferencepopulationbeingcopiedfromalongthechromosome–notethatatmostpositions,this
isidenticaltothetrueunderlyingancestry,butwithoccasional‘‘miscopying’’fromtheotherpopulation(bluedottedsegmentoccurringwithinred
ancestrysegment).Notealsothatswitchesbetweenchromosomesbeingcopiedfrom,representinghistoricalrecombinations,arerapid(6switches),
whileancestrychanges,representingrecombinationsinceadmixture,aremuchrarer(1switch).Finally,notethatatmostpositionsthetypeofthe
admixedchromosomeisidenticaltothatofthechromosomebeingcopiedfrom,butanexceptiontothisoccursatonesite,shownasagreycircle,
and representing mutationor genotyping error. In our inference framework, we observe only the variation data for the admixed and reference
individuals:theyellowline,andtheunderlyingancestry,mustbeinferredasthehiddenstatesinaHMM.
doi:10.1371/journal.pgen.1000519.g001
produce accurate estimates of uncertainty in inferred segments, fromasingleadmixtureeventbetweentwopopulationswhichare
even foroldadmixture events. geneticallysimilartothetworeferencepopulationsfromwhichP
1
and P are drawn. (The reference populations do not need to
2
Details of haplotype-based inference of local ancestry exactly matchthetrue ancestralpopulations, becauseweallowfor
Modelinggeneticvariationinadmixedpopulations. Our some genetic divergence in our approach.) We will initially
approach to inferring ancestry segments, implemented in consider the case where we have haploid chromosomes from the
HAPMIX, is based on extending a Hidden Markov Model admixed population, and subsequently generalize to the more
(HMM)previouslydevelopedbyLiandStephenstomodellinkage typical case involving unphased genotype data from the admixed
disequilibrium in population genetic data [14]. This model has population.Throughoutthissection,weoperateinunitsofgenetic
been employed in recent years in various population genetic and (not physical) distance.
disease mapping settings [20,21]. Informally, given a previous We begin by modeling the ancestry segments. Assume the
collectionof‘‘parental’’haplotypesfromareferencepopulation,a admixtureeventoccurredatasingletimeTgenerationsago,with
new ‘‘offspring’’ haplotype drawn from the same population is afractionm1ofthehaplotype’sancestrydrawnfrompopulation1,
modeled as a mosaic of these existing haplotypes. This offers a andm2=12m1frompopulation2.Becauserecombinationoccurs
flexible means to account for local linkage disequilibrium (LD), at each generation, it is natural to model ancestry switches as a
because over short distances, the haplotype that an individual Poisson process along the genome [22], at a rate T per unit of
chromosome copiesfromisunlikely tochange. geneticdistance(i.e.TperMorgan).Conditionalonthepositions
We extend the Li and Stephens model to allow inference on of such switches, each segment is independently drawn from
ancestry segments for individuals drawn from an admixed population 1 or 2 with probabilities m1, m2 respectively. In
population. We begin by supposing that we have two previously particular, this implies that not all ancestry switch points will
sampled collections of phased haplotypes, P and P , taken from actually change the underlying ancestry. This model has been
1 2
two reference populations. For example, HapMap provides previouslyusedbyotherauthors[1,22].Sinceancestrycannotbe
phased haplotypes from the CEU, YRI and JPT+CHB popula- directlyobserved,itisnaturaltoviewunderlyingancestrystatusas
tionsgenotypedatover3millionmarkers[16].Wefurtherassume the ‘‘hidden’’ information in an HMM. Our approach probabi-
that P and P have valid data at all sites of interest, with no listically infers this hidden state at each position along a
1 2
missing data. In practice, small amounts of missing data in the chromosome.
reference populations can be filled in by a pre-processing To fully specify our model, we must consider the structure of
imputation step, as has been done for the publicly available variation conditional on these admixture segments. Our model
phased HapMap data. We label P and P as ‘‘parental’’ remains computationally tractable while accommodating impor-
1 2
haplotypes. Next, we sample a new ‘‘offspring’’ haplotype from tantfeaturestypicalofrealdatasuchasmutation,recombination,
anadmixedpopulation.Weassumethatthispopulationiscreated genotyping error, reference populations that are drifted from the
PLoSGenetics | www.plosgenetics.org 3 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
trueancestralpopulations,andincompletesamplingofdiversityin miscopyingparameterisdesignedtoallowforregionsinthegenome
the reference populations reflected in the samples drawn from wheretheoffspringchromosomehasanunusuallydeepcoalescence
thesepopulations.Weassumethatallmutantsitestaketheformof timewiththeothersamplemembers,allowingthe‘‘miscopying’’to
single nucleotide polymorphisms (SNPs) with two alleles that can occur, miscopied regions are likely to have unusually deep
be represented as 0 and 1 (however, our approach could be genealogies. Therefore, we allow a different mutation rate h for
3
extended tomore complexmutation models). suchsegments,whichistypicallyexpectedtobehigherthanh or
1
Wesupposethatsectionsofthegenomewithtrueancestryfrom h .Itmightalsobedesirabletoallowahigherrecombinationratein
2
population 1 are formed as mosaics of the haplotypes in the two suchcases.However,thiswouldresultincomputationalcomplex-
parentalgroups.Specifically,atanygivenpositionwiththisancestry, ities,andwehavechosennottoallowsuchanadditionalparameter.
anindividualfromP1iscopiedwithprobability,andanindividual For a typical application of HAPMIX, we expect to have data
from population P2 is copied with probability p1 (we call this the from a collection of discrete typed sites. Suppose we have S such
‘‘miscopying’’ parameter for population 1). Conditional on the sites, and in addition a map giving the genetic distances
parental group chosen, individuals to copy from are chosen r ,r ,…r between adjacent pairs of sites. In practice, we
1 2 (S21)
uniformly from the n1, n2 respective individuals in that group. interpolate these distances from the genome-wide recombination
SwitchesbetweenindividualsoccurasaPoissonprocesswithrater1, rates estimated using Phase II HapMap [16]. Given the above
the ‘‘recombination’’ parameter, and at each switch point a new parameters,andforahaploidadmixedchromosome,weformalize
copyindividualischosenrandomlyusingtheabovescheme.Finally, thetransitionprobabilitiesasfollows.A(hidden)stateforpositions
atgenotypedSNPs,ifthe‘‘offspring’’copiesa‘‘parent’’fro1{p1m isrepresentedbyatriplet(i,j,k)wherei=1or2representsancestry
population1,theoffspringcarriesanidenticaltypetotheparticular drawn from population 1 or population 2, j=1 or 2 records the
parentitcopiesfromwithprobability(12h1),andcarriestheother population the chromosome copies from at position s (j may be
typewithprobabilityh1,the‘‘mutation’’parameter.Iftheoffspring differentfromiduetomiscopying)andkrepresentstheindividual
instead copies an individual from the other population 2, the from which the chromosomal segment is copied. There are
correspondingmutationparameterish3.Intotalthisapproachleads 2(n1+n2)possiblestates.Letpsði,j,k; l,m,nÞbetheprobabilityof
to4additionalparameters:p1,r1,h1andh3. transitioningfromstate(i,j,k)tostate(l,m,n)betweenadjacentsitess
Forsectionsofthegenomewithancestryfrompopulation2,we and(s+1).Thenwehave thefollowing:
formulate our model in an analogous way, with corresponding
parameters p , r , h and h . We note that h is shared for both
populations, 2a 2cho2ice tha3t is motivated 3by a genealogical psði,j,k;l,m,nÞ~
pbtpitclappimhhornaoaaaeegarrrrSetpyhuaaaalasoleemmmmmaramstmlmchryeeeeeeeeeettnttpneeeewnntletrrrre,at-ssofissabdtf.:stehtaadimmIiTczmnosnitetena,daniisev1donyo,mrennuaeh1aa1rbgaprl,ranyelae1gopspdanua1parbu,dmauternslhsnpalemi2neeame2dm,.2fnt,auipItaoderralsnrl.tknnrei1esm2sad,Arutswfioriugzsarzint2ngaenoeba,eitlcscflfoani,shrtorutotioe1okenotfttm,aen,huhtsehthreatteopot2hhtehfimeeornaterstohhrnsgoheeiueeiidgldesneygptisptetvnohhaehewloaralr3reypinqaihl.pcluutmlaraLtlleadsiploaenitneottat9eairtcepatoiaelrirtunnocrieisnelldonsanyads.dmsu,teeAeiSmwooipbpnrtlfnetiieeboftnhnsopnhue,atodrhdrhgtaumeiethenooghnnoslahdnnyeestfl 8>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><ððeez|11{{{{ððrr11ssTTee{{{{||rreesseð{{TT1{ÞÞrr{ssmmrTrsllrleÞ||Þlmz{|lr1n|pse{m1nrl{{mnlpÞ1mrlp|{nslTmpl1{nmpl iiðiijffff~j=llll==~~mmiiiiaaaoaannnrnndddddkk=mmmm~~=~~nnÞllllaanndd 9>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>= ð0:1Þ
pmegeTatstmbinooi2ffenemhffcuinceecszoeptaieacceeealmsottlusrarmltiittaporvvoisirpomlelweeuyu,ylnielbtwsapasfatpaeoetetseotliocingiirioopkmolaopmdinpneunnunuilolegelsyollaapneayp1eatfrttnugtaitsioosaheroloma.a‘nnaet‘nmlsItfymuecioufaspoeapsirzmspetzsnditlizntte1eizreoinc~2ge,tnyrauoogs.wseenpspHfscrseeihotys2ranp’go1ofo~’nduim,wtfirornet.iohh0hieendoTcin2n,vnesiuceuvothatautossrimtnrhedn,,rreacee.iudmsfbfcneoxeeIaeasncprtrighlrtsltesereri3crasnsnoyocoeottfcgafturcrpmselyloppoereyltomhegmaai1bpnwmerscp=rboetgsoeaeepnipftrappnrnowotupkti2aaatrpesl,soastlrad,uioarbetcocraliamfvohaonoemfnetttrrranehhinddorosoetseeeomemtnishteprnfg[hposrie1esnseo1ssn.ptnw41ron,oe1dF]eo=imdrq.tfsaamiihoirucnnnoenIsriaae2stndnndge,.,,ll >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>:eezz|z{{ððððrr1111ssTT{{{{||eeeeeð{{{{1{rrrr{ssssrTTrTsrleÞÞÞÞlmmmz{|lllr|||senrp{mllÞ1nnppr|{mmnsllTmplnpml ðiijff~j=ll~~mmiiaoaanrnndddkk=mm~==nnÞllaanndd >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>;
theoffspring chromosomeisexpected tohave a deep coalescence Conditional on the underlying hidden state, let e1 ðsÞ denote the
ijk
timewiththeancestorsofthe‘‘correct’’parentalsample,andmay
probabilityoftheoffspringchromosomebeingoftype1atsites,and
insteadcoalescefirstwithanancestoroftheotherparentalsample–
t bethetypeofparentalindividualkinreferencepopulationj.Then
andthereforechooseadescendantofthisancestor,inthe‘‘wrong’’ jk
parentalsample,tocopyfrom.Second,ifourreferencepopulations
aresomewhatinaccuraterelativetothetrueancestralpopulations, e1 (s)~(cid:1)hid(tjk~0)z(1{hi)d(tjk~1) if i~j ð0:2Þ
againitislikelythatincompletelineagesortingwill occur,evenif ijk h d(t ~0)z(1{h )d(t ~1) if i=j
3 jk 3 jk
our‘‘parental’’samplesarebothlarge.Forthesereasons,inpractice
we believe that incorporating non-zero miscopying parameters is
important, and in both real data and simulation we find that it This probability allows us to calculate the likelihood of the
greatly improves our ancestry estimation procedure. Because our observed data in the offspring for each possible underlying state.
PLoSGenetics | www.plosgenetics.org 4 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
At sites with missing data in the offspring chromosome, the in an analogous manner to methods employed in several existing
appropriate likelihoodcontribution issimply 1.0. approaches to imputation for samples drawn from panmictic
Choices of parameter settings. Choices of T and m are populations [20,21]. We will describe this extension, and its
1
specific to each application (see below). However, many of the application todisease mapping,in a separate paper.
remainingparameterswerefixedinallanalysesofbothsimulated Multiple individuals from the admixed
andrealdata.Asdiscussedabove,itisnaturaltoscaler andr , population. Typically, we actually have multiple ‘‘offspring’’
1 2
aswellash andh ,bythenumbersofparentalindividualsn ,n , samples (either haploid chromosomes or diploid genotypes, see
1 2 1 2
respectively.Ourcodeisparameterizedsothisisdoneinternally– below)fromtheadmixedpopulationofinterest.Fortheanalysesin
arbitrarily labeling the European population as ancestral this paper, we used HAPMIX to analyze data from each sample
population 1, we used recombination parameters r =60,000/n independently,usingthesameparentalchromosomesineachcase.
1 1
per Morgan for the European ancestral population and Although in principle improvements to ancestry inference could
r2=90,000/n2 per Morgan for the African ancestral population resultfromconsideringtheprobleminmultiplesamplesjointly,there
(with r2.r1 reflecting the larger effective population size of areformidablecomputationalchallengesinadaptingourapproach
Africans). Further, we set h1=0.2/(0.2+n1) and h2=0.2/(0.2+n2), toallowthis(onepossibilitymightbetoemployMCMC,asusedfor
andh3=0.01(this parameterremains unscaled). Finally,weused unlinked sites [22,23]). To avoid these complications, we simply
miscopying parameters p1=p2=0.05. These values were arrived modeleachadmixedsampleindependently,following[21].Under
atviaaprocessoftrialanderror,basedontheresultsofinferring this scheme, separate HAPMIX runs for each sample enable
parametersviatheEMalgorithm.WehaveimplementedanEM effectiveparallelizationofthesoftware.
algorithm approach to parameter estimation that can infer any Diploid genotype data from the admixed
subset of the HAPMIX input parameters, or all simulataneously population. Typically, real data consists of unphased
(see Text S2). This EM approach to parameter inference is genotypes for individuals drawn from a population, with
currently only implemented for haploid data from the admixed haplotypic phase unknown. Many approaches already exist to
population, but we applied it to haploid data derived from a inferphasefromsuchdata[17,18].However,phaseswitcherrors
phasingofdiploiddata,obtainedbyrunningHAPMIXondiploid that inevitably result from applying such algorithms are likely to
admixed samples and using the software to sample random state result in spurious ancestry switches within regions of the genome
paths. This approach might be applied to diploid samples more whereanindividualisheterozygoteforancestry.Thiswouldlikely
generally, and could be potentially be iterated, by updating lead to considerable overestimation of the time since admixture
phasing based on new parameter sets. However, based on our and a reduction in the accuracy of ancestral inference. To avoid
simulations we believe that for many applications – for example such issues, we have extended our approach to directly analyze
wheneverthesoftwareisappliedtoAfricanAmericandata-itwill diploid genotype data fromtheadmixed population.
besufficienttovaryTandm1andfixtheremainingparametersat The phasing is implemented using a HMM adapted from that
thevalues describedabove.
described above (0.1) and employing a composite hidden state at
Inferring probabilistic ancestry segments and sampling eachlocation,oftheform(i ,j ,k ,i ,j ,k )where(i ,j ,k )represents
1 1 1 2 2 2 1 1 1
from the posterior with HAPMIX. It is easy to see that
the previously defined ‘‘haploid’’ hidden state for the first
equations (0.1)and(0.2) describeaHMM fortheunderlying state
chromosome, and (i ,j ,k ) represents the hidden state for the
2 2 2
(which includes information on ancestry) as we move along the
secondchromosome.Thestatespacethereforenowhasdimension
genome,andthattheunderlyingMarkovprocessisreversible.Given 4(n +n )2.Allowingindependent transitionsbetweenthemarginal
a set of parameters we can exploittheseproperties andHAPMIX 1 2
states for each chromosome, the terms in (0.1) now naturally
implementsstandardHMMtechniquestoefficientlyinferposterior
define an HMM for these composite states (for reasons of space,
probabilities of underlying states, via the forward-backward
we do not explicitly list all of the transition probabilities in the
algorithm, or sample random state paths from the correct joint
modelhere).Thismodelcouldhaveupto18parameters–inour
posterior distribution, using a standard modification of this
implementation, for natural biological reasons we assume all
algorithm. In addition to parameter values, the software takes as
parameters are shared between chromosomes, apart from time
inputarecombinationmapfortheregionstobeanalyzed,phased
since admixture T and admixture proportion m , resulting in 11
‘‘parental’’ chromosomes from the two reference populations, and 1
parameters in total. Further, although our software allows these
‘‘offspring’’datafromtheadmixedpopulationbeinganalyzed.
two parameters to differ, in all applications considered here we
A na¨ıve implementation of the forward/backward algorithm
specifyT andm tobe thesamefor eachchromosome.
wouldrequirecomputationtimeproportionalto4S(n +n )2,inthe 1
1 2 Emission probabilities are also adapted from the haploid case.
above notation. For the original Li and Stephens model, it is
For genotype data, there are 3 possible emissions at typed sites,
possibletoreducecomputationtimesubstantiallybyusingthefact
whichwedenoteasgenotypesg=0,1,or2,withgcountingcopies
that many pairs of transition probabilities between states are
of the ‘‘1’’ allele. Conditional on the underlying hidden state, let
identical, which allows terms to be collapsed in the forward (or
eg ðsÞ denote the probability of observing genotype g given
backward)algorithm,intoexpressionsinvolvingasingletermthat ijklmn
underlyingstate(i,j,k,l,m,n),anddefinet asbeforetobethetypeof
issharedamongalldestinationstates.Calculatingthissharedterm jk
parental individual k inreference population j.Thenusing (0.2)
just once per pair of adjacent sites, and then storing, saves
substantial computational effort [14]. Analogously, in our
somewhat more complicated setting we can exploit a similar e0 ðsÞ~(cid:2)1{e1 ðsÞ(cid:3)|(cid:4)1{e1 ðsÞ(cid:5)
ijklmn ijk lmn
phenomenon,sothatbycalculatingandstoringasomewhatlarger
numberofsharedterms–oneforeachgroupofstatesoftheform e1 ðsÞ~(cid:2)1{e1 ðsÞ(cid:3)|e1 ðsÞze1 ðsÞ|(cid:4)1{e1 ðsÞ(cid:5) ð0:3Þ
(i,j), giving four in total - HAPMIX can complete the forward/ ijklmn ijk lmn ijk lmn
backward algorithm in time proportional to 2S(n1+n2) (with an e2 ðsÞ~e1 ðsÞ|e1 ðsÞ
ijklmn ijk lmn
additional scaling constant).
Itisstraightforwardtoextendourapproachtoallowimputation
ofmissingdata,whilesimultaneouslylabelingunderlyingancestry, where e1 ðsÞand e1 ðsÞ areas defined above.
ijk lmn
PLoSGenetics | www.plosgenetics.org 5 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
Having defined the HMM for this setting, we again use computed using ancestry predictions only) can be compared to
standard techniques to obtain posterior probabilities on (joint) the actual r2 between these quantities (comparing ancestry
ancestry for the two chromosomes, and then sample states from predictions to true ancestries specified in simulations). As we
this posterior distribution. We note that as a by-product of confirmin whatfollows,theestimates ofr2are well calibrated.
sampling complete states jointly for the two chromosomes
together, we are phasing the original data with respect to the Simulations
underlying ancestry.This mayhelpreduce phasingerror ratesin Simulations of local ancestry inference. We simulated
admixed populations compared to methods that ignore local individuals of admixed African and European ancestry by
ancestry, although wedonot pursue thisissuehere. constructing their genomes from a mosaic of real Yoruba and
We can adapt the computational speedups described above to FrenchindividualsgenotypedontheIllumina650Ychipaspartof
the diploid setting, so that while a na¨ıve implementation of the the Human Genome Diversity Panel (HGDP) [15]. We
forwardalgorithmwouldtaketimeproportionalto16S(n +n )4,we downloaded data from 20 Yoruba and 20 French individuals
1 2
can complete the forward/backward algorithm in time propor- from the HGDP data set and jointly phased them using the
tional to 4S(n +n )2. A further speedup for the diploid setting is fastPHASE program [18] to form 40 haploid Yoruba and 40
1 2
described in Text S2. With these speedups implemented, the haploid French genomes.
running time of HAPMIX is roughly 30minutes on a single Weconstructed40haploidadmixedgenomes(n=1to40)from
processorperdiploidgenomeanalyzed(519,248sites).Becausethe the40haploid Yoruba and40haploid Frenchgenomes by using
computationscanbeparallelizedacrossadmixedindividuals(they haploid Yoruba genome n and haploid French genome n to
can also be parallelized across chromosomes), HAPMIX is construct admixed genome n, so that ancestral genomes were
computationallytractableevenforverylargedatasetsifacluster neverreused.Toconstructanadmixedgenome,webeganatthe
ofcomputingnodesisavailable.Forexample,therunningtimefor first marker on each chromosome and sampled French ancestry
adatasetof1,000admixedindividualsonaclusterof100nodesis with probability a and Yoruba ancestry with probability 1-a.
roughly 5hours. Ancestrywasresampledbasedonanexponentialdistributionwith
weightl(thenumberofgenerationssinceadmixture)sothatanew
Measuring the performance of HAPMIX ancestry was sampled with probability 12e2lg when traversing a
Estimate of r2 between predicted and true geneticdistanceofgMorgans.Eachtimeancestrywasresampled,
ancestry. Irrespectiveof whether thetrueancestryisknown (as we sampled French ancestry with probability a and Yoruba
in simulations) or unknown (as in real data), an estimate of the r2 ancestrywithprobability1-a.Foreachindividual,weusedavalue
betweenapredictedancestryvectorYandtrueancestryXcanbe of a to apply to the entire genome by sampling from a beta
computed.Withinanindividual,ateachsites,anaturalmeasureof distribution with mean 0.20 and standard deviation 0.10 (typical
predictedancestryistheexpectednumberY ofhaplotypesfromone for African Americans [4]). We simulated values of l=6 (typical
s
of the two source populations. If HAPMIX provides accurate forAfricanAmericans[4])aswellashighervaluesofl:10,20,40,
ancestry probabilities, the true number of haplotypes from this 60, 100, 200 and 400. Pairs of haploid admixed individuals were
population,X,canbethoughtofasanunknownrandomvariable merged toform20diploid admixed individuals.
s
whichisequalto0,1,or2withprobabilitiesp ,p ,p specifiedbythe It is important to distinguish between the true ancestry
0 1 2
ancestry predictions. We are interested in how correlated the proportion a in a simulated or real admixed individual and the
predictedancestryYandtrueancestryXare,oversamplingsfrom parameterm1usedasinputtoHAPMIX,whichmaydifferfroma(if
thisdistributionofthetrueancestryX.Anaturalwaytoestimatethis aisunknown).Similarly,itisimportanttodistinguishbetweenthe
correlationistocalculatetheexpectedsquaredcorrelationbetween true number l of generations since admixture and the parameter T
XandY,whichwemayapproximateusingaratioofmeans: usedasinputtoHAPMIX.Belowweexploretheconsequencesof
inaccurately specifying theparameters m and T.
1
! ThereferencepopulationsusedasinputtoHAPMIXconsisted
E(cid:4)r2 (cid:5)~E covðX,YÞ2 of 60 YRI individuals (120 haploid chromosomes) and 60 CEU
XY VarðYÞVarðXÞ individuals (120 haploid chromosomes) from the International
HapMap Project [16]. A joint analysis of HGDP and HapMap
EðcovðX,YÞÞ2 data indicated that F (Yoruba,YRI)=0.000 and F (French,-
& , ST ST
VarðYÞEðVarðXÞÞ CEU)=0.001, so that the reference populations used as input to
HAPMIX were extremely accurate. All HAPMIX simulations
where the variances and covariances are taken over loci and were restricted to 519,248 autosomal markers present in HGDP
individuals,andtheexpectationsoversamplingsoftheancestryX. data which were polymorphic in phased YRI and phased CEU
The expected covariance between predicted and true ancestry is data from HapMap. For comparison purposes, we ran the
then the mean value of the covariance between X and Y as we ANCESTRYMAP, and LAMP-ANC programs on the same
sampleancestrypathsatdifferentlociandindifferentindividuals.At simulateddatasets,makinguseofdiploidYRIandCEUgenotype
our single locus, we have E(XsYs)=(p1+2p2)2 and data from HapMap and restricting all input data to subsets of
E(Xs)=E(Ys)=p1+2p2. By separately averaging these three markers that were unlinked in the reference populations, as
expectations across loci and individuals, we can then calculate recommended by those methods [1,11]. For the ANCESTRY-
EðcovðX,YÞÞanalytically.Similarly,wecancalculatethevariance MAPruns,wechoseasubsetof,8,000markerswiththelargest
of Y, and the expected variance of X, across loci and different YRI-CEU differences that were unlinked in both reference
individuals, in a similar way. Combining these variances with the populations. For the LAMP-ANC runs, we set the LD cutoff to
covariancetoestimatecorrelation,andthensquaring,weobtaina 0.10, causing the program to choose a subset of ,260,000
measureofthelevelofcertaintyoftheancestrypredictions. markers. We note that LAMP-ANC differs from the LAMP
Actual r2 between predicted and true ancestry. In program in that LAMP-ANC makes use of input data from
simulated data sets where the true ancestry is known, the reference populations [11], which makes it more comparable to
estimated r2 between predicted and true ancestry (which is HAPMIX. We attempted to run HAPAA on the same data for
PLoSGenetics | www.plosgenetics.org 6 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
comparisonpurposes.However,despiteadvicefromtheauthorsof populationsused(YRIandCEU),andforourinaccurateancestral
the software and extensive effort, we were unable to make the populationsimulationsinwhichtrueancestralpopulations(either
linked applications that form the HAPAA software suite run on Mandenka and Basque, or Yoruba and Druze) are inaccurately
our computers, and hence we were unable to make this modeled by the reference populations (YRI and CEU). We
comparison. restrictedtheseanalysestodatasimulatedusingl=6andl=100
Simulations of local ancestry inference using inaccurate only.
reference populations. We repeated our simulations at l=6 Inference of date of admixture. By comparing the overall
andl=100usingMandenkafromHGDPastheAfricanancestral likelihoodsproducedbyHAPMIXatvariousparametersettings,it
population and Basque from HGDP as the European ancestral ispossibletoevaluatewhichparametersprovidethebestfittothe
population for simulating admixed individuals. We simulated 20 data, irrespective of whether or not the choice of parameter
admixedindividualsusingMandenkaandBasquedata(analogous settings significantly impacts the accuracy of local ancestry
tothesimulationsdescribedaboveusingYorubaandFrenchdata). inference. We investigated how effectively the number of
We continued to use YRIand CEUas thereference populations generations l since admixture can be inferred in this way by
for HAPMIX. A joint analysis of HGDP and HapMap data running HAPMIX at various values of T and computing overall
indicated F (Mandenka,YRI)=F (Basque,CEU)=0.01. We likelihoods, using the data sets simulated at l=6, l=20 and
ST ST
note that these discrepancies between the ancestral populations l=100.Wealsosimulatedadouble-admixturescenarioinwhich
used to construct these simulated data and the reference a50%/50%admixtureofYorubaandFrenchoccurredatl=100
populations used as input to HAPMIX are substantially larger followedbya50%/50%admixtureofthatpopulationandFrench
than the discrepancy between the true African ancestral at l=6 (we call this the l=6+100 run (with a=75%)). We
population of African Americans and YRI, or the true European optimized T at a granularity of 1 for the l=6 and l=20
ancestral population of AfricanAmericansandCEU [4]. simulations and a granularity of 5 for the l=100 and l=6+100
To investigate the scenario of an even more inaccurate simulations.
referencepopulation,aswellastheasymmetricscenarioinwhich
onlyonereferencepopulationisinaccurate,wealsorepeatedour Analysis of 935 African American samples
simulationsatl=6andl=100usingYorubafromHGDPasthe We used HAPMIX to analyze 935 African American samples
African ancestral population and Druze from HGDP as the collected from volunteers living in the Baltimore–Washington,
Europeanancestralpopulationforsimulatingadmixedindividuals. D.C. metropolitan region and genotyped on the Illumina 650Y
We simulated 20 admixed individuals using Yoruba and Druze chip as part of an asthma study. All subjects gave verbal and
data as described above, and continued to use YRI and CEU as written consent. The Johns Hopkins and Howard University
thereferencepopulations.AjointanalysisofHGDPandHapMap Institutional Review Boards (IRBs) determined that the samples
data indicatedthat F (Druze,CEU)=0.02. were consented for genetic research, but not for public release of
ST
Simulations of local ancestry inference as a function of genotype data. Roughly half of these samples were asthma cases
datasizeandparametersettings. We modified our original and half were non-asthmatic controls, but all phenotypic
simulations at l=6 and l=100 to consider different data sizes information was ignored in the current study (disease mapping
and parameter settings. We investigated how the performance of analysesofthesedatawillbedescribedelsewhere;K.Barnesetal.,
HAPMIXvariesasafunctionofdatasizebyvaryingthenumber unpublisheddata).Wenotethatirrespectiveofwhetherasthmatic
ofmarkersfrom5,192randomlyselectedmarkerstothefullsetof cases considered separately exhibit an admixture association
519,248 markers, and by varying the amount of input data from signal, one would not expect to observe such a signal in a
YRI and CEU reference populations from 10 haploid combined analysis of all 935 samples ignoring phenotypic
chromosomes to the full set of 120 haploid chromosomes. We information, due to dilution of the signal. The analyses were
investigated how the performance of HAPMIX varies as a restricted to 510,324 autosomal markers which passed quality
function of parameter settings by incorrectly specifying either the controls in the 935 African Americans and were polymorphic in
Europeanancestryproportionm usedasinputtoHAPMIX(using phased YRI and phased CEU data from HapMap. We ran
1
values different from a<20%) or the number of generations T HAPMIX using YRI and CEU as input reference populations,
sinceadmixtureusedasinputtoHAPMIX(usingvaluesdifferent setting m =20% and running at various values of T to infer the
1
froml=6or l=100,respectively). date of admixture (see above). For comparison purposes, we also
Inferenceofancestralpopulations. ByrunningHAPMIX ran the ANCESTRYMAP and LAMP-ANC programs on this
inthemodethatsamplesrandompaths,whichproducesinteger- data,ineachcaserestrictingallinputdatatoasubsetofmarkers
valued guesses of local ancestry for each individual and each that were unlinked in the reference populations, as described
marker, it is possible to reconstruct chromosomal segments from above.
the ancestral populations. We investigated whether these To draw inferences about the ancestral populations of African
reconstructed segments provide an accurate proxy for the true Americans, we ran HAPMIX in the mode that samples random
ancestral populations by using allele counts to compute values of paths to reconstruct chromosomal segments from the ancestral
F (astandardmeasureofgeneticdistance[24])betweeninferred populations (see above), and used the resulting allele counts to
ST
ancestral segments and true ancestral populations from our compute F values between the inferred ancestral segments and
ST
simulations. Although ancestral individuals are used twice in this the reference populations (YRI and CEU), as well as additional
computation (both to simulate admixed individuals whose populations genotyped as part of the HGDP. To estimate the
ancestral segments are inferred, and in the ancestral populations number of ancestry segment changes in each of the 935 African
themselves), we restricted this analysis to half of the ancestral American individuals, we inferred ancestry using the most likely
individuals for the former and the other half of the ancestral state at each site, and identified ancestry transitions from these
individualsforthelatter,thuspreventinganyduplicationofdatain ancestrystates,assumingzerochangesbetweenpairsofSNPswith
thecomputationofF .Weperformedthiscomputationbothfor identical ancestry states.
ST
our original simulations in which the true ancestral populations To produce an estimator of the number of generations since
(Yoruba and French) are accurately modeled by the reference admixture for each individual with .20 ancestry segments, we
PLoSGenetics | www.plosgenetics.org 7 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
note that the genetic map used as input to the software has total ancestry reconstruction in African Americans: the squared
length35.5Morgans.Foranindividualwithadmixtureproportion correlation between predicted and true number of European
a, we expect to observe a fraction 2a(1-a) of all recombination copies (across all samples) was equal to 0.98, and discernment of
eventsoccurringsinceadmixture(i.e.thosethatresultinachange ancestry transitions was extremely sharp, as seen in a plot of the
in ancestry). Given l generations since admixture, we therefore predicted vs. true number of European copies for an admixed
expecttoseeatotalof142la(1-a)eventsinadiploidindividual. sample on chromosome 1 (Figure 2A). For comparison purposes,
Estimatingausingtheobservedgenome-wideancestryproportion we also computed local ancestry estimates using the
mforthatindividual,ifNancestrytransitionsareobserved,thena ANCESTRYMAP and LAMP-ANC programs [1,11] (see
natural moment estimator of the number of generations since Materials and Methods). (We chose not to explicitly compare
admixture is HAPMIX to additional recently developed methods such as
SABER, LAMP, uSWITCH and uSWITCH-ANC [10–12] ,
becauseinpreviousworktheLAMP-ANCmethod—whichwedo
N
^ll~ : compare HAPMIX to—has been shown to perform
4mð1{mÞ|35:5
approximately as well as each of those methods in a range of
scenarios [11].) The squared correlation between predicted and
true number of European copies was equal to 0.86 for
Weexcluded3clearoutlierindividualswhohadmorethan20
ANCESTRYMAP, 0.83 for LAMP-ANC and discernment of
inferredgenerationsofadmixture,becausewebelievethisislikely
ancestry transitions was less sharp or sometimes missed entirely
toindicatepartialancestryfromathirdsourcepopulationinthese
(Figure 2A).
individuals.
A more challenging setting for ancestry inference is when
admixture occurs further back in time, resulting in smaller
Analysis of 29 Mozabite samples
ancestrysegments.Wethereforerepeatedtheabovecomparisons
Weanalyzed29MozabitesamplesfromtheHGDPdataset.A
with increasing lambda (Figure 3). The results show a uniformly
totalof30Mozabiteindividualswereoriginallygenotypedaspart
better performance by HAPMIX relative to the other two
oftheHGDP,butoneindividual(HGDP01281)wasexcludeddue
methods,withthecomparativeadvantageofHAPMIXincreasing
to cryptic relatedness. We ran HAPMIX on the 29 Mozabite
with timesince admixture.
individuals using YRI and CEU as the input reference
To investigate whether the probabilities of 0, 1, or 2 copies of
populations. We inferred the number of generations since
European ancestry reported by HAPMIX are well-calibrated, we
admixture that provided the best fit to the data, and computed
binned the predicted probabilities into bins of size 0.05 and
F values between the inferred ancestral segments and the
ST compared, for each x=0,1,2 and for each bin, the average
referencepopulations(YRIandCEU),asdescribedaboveforthe
predicted probability vs. the actual frequency in simulations of
African Americandata set.
having x copies of European ancestry. For example, in the l=6
simulation, restricting to instances in which the predicted
Analysis of other HGDP populations
probability of 1 copy of European ancestry was between 0.05
WeranHAPMIXonatotalof13populationsfromtheHGDP and0.10,theaveragepredictedprobabilityof1copyofEuropean
datathatwereofAfrican,European,orMiddleEasternancestry. ancestry was 0.07 and the true frequency of 1 copy of European
Foreachpopulation,weusedYRIandCEUastheinputreference ancestry was 0.08, which is close to 0.07. More generally, we
populations,andestimatedtheEuropean-relatedmixturepropor- observed that HAPMIX predictions from our l=6 and l=100
tion. For populations with European-related ancestry that was simulations were well calibrated for each value of x=0,1,2
estimated to be more than 0% and less than 100%, we also (Figure 4). The calibration of intermediate bins appears visually
estimated thenumber of generationssince mixture. worse for the l=6 simulation; however, the proportion of the
genome that is in the most extreme bins where the method is
Web resources certainis98%,97%,99%,forx=0,1,2inthesesimulations,and
The HAPMIX software is available for downloading at the hencethereliabilityoftheprobabilitiesremainsgoodforrecently
following URL: http://www.stats.ox.ac.uk/,myers/software. admixed populations too.
html. WealsousedtheHAPMIXpredictionstocomputeanestimate
ofthesquaredcorrelationbetweenpredictedandtrue#European
Results copies(seeMaterialsandMethods).Weobtainedestimatesof0.98
forthel=6simulationand0.83forthel=100simulation,which
Simulations are identical to the true r2 values of 0.98 for l=6 and 0.83 for
Simulations of local ancestry inference. We began by l=100,consistentwiththefindingthatHAPMIXpredictionsare
examiningtheperformanceofHAPMIXinasetof20simulated well calibrated.
admixedindividuals,withanaverageof80%Africanancestryand Although most of our simulations focused on individuals of
20%Europeanancestry,andgeneratedwithadmixtureoccurring mixedAfricanandEuropeanancestry,wealsoconsideredamore
6 generations ago (l=6; see Materials and Methods). These general set of two-way mixtures of African, European, Chinese
parameters were chosen to be in the range of typical values for and/or Japanese populations. We again observed that HAPMIX
African Americans. We implemented a simulation framework in outperformedothermethods(seeTextS1).Furthermore,although
which admixed individuals were constructed using genotype data HAPMIX is currently implemented assuming only two reference
from the Human Genome Diversity Project, but modeled using populations, we were able to attain accurate results in a more
referencepopulationsfromHapMap(seeMaterialsandMethods). complex scenario of three-way admixture, by running HAPMIX
WecomparedthelocalancestryestimatesproducedbyHAPMIX inatwo-waymodeusingdifferentchoicesofreferencepopulations
(probabilitiesof0,1,or2copiesofEuropeanancestry)tothetrue (see TextS1).
values of local ancestry that were simulated. These simulation Simulations of local ancestry inference using inaccurate
results suggest that our method is likely to provide near optimal referencepopulations. In many real-world settings, the true
PLoSGenetics | www.plosgenetics.org 8 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
Figure 2. Comparison of ancestry estimates produced by HAPMIX, ANCESTRYMAP, and LAMP-ANC. (A) Results comparison for a
simulatedrecentlyadmixedsampleonchromosome1.Oneachplot,they-axisdenotesthenumberofEuropeanchromosomalcopiespredictedby
eachmethod.Thecentromereofthechromosomeisblankedoutinwhite.ThetopplotshowsthetruenumberofEuropeanchromosomes,whilethe
subsequentlabeledplotsshowtheresultsofapplyingeachrespectivemethod.(B)ResultscomparisonforarealAfricanAmericanindividualacross
chromosome1.Plotsareconstructedasin(A).Wenotethevisiblesimilaritytothesimulationresults.
doi:10.1371/journal.pgen.1000519.g002
reference populations for a particular admixture event may not Simulations of local ancestry inference as a function of
havehadsuitablegeneticdatagathered,ormaynolongerexist. data size and parameter settings. We investigated how the
TotestfortheeffectofthissituationonHAPMIX,werepeated accuracyofHAPMIXvarieswithdatasize,byvaryingeitherthe
our simulations at l=6 and l=100 using 20 admixed samples number of markers or the number of reference chromosomes, in
that were simulated using Mandenka and Basque individuals our l=6 and l=100 simulations (see Materials and Methods).
but modeled using reference populations YRI and CEU, which Accuracy as a function of the number of markers is displayed in
are inaccurate reference populations (see Materials and Figure4A,whichshowsthatasfewas50,000randommarkersare
Methods). For l=6, the squared correlation between close to optimal for l=6 but that hundreds of thousands of
predicted and true #European copies remained high at 0.95, markers are needed to produce optimal results in the more
only marginally worse than the 0.98 obtained using accurate challenging case where l=100. Accuracy as a function of the
reference populations. For l=100, the squared correlation was number of reference chromosomes is displayed in Figure 4B,
0.76, again only slightly worse than the 0.83 obtained using which shows that as few as 40 chromosomes (phased from 20
accurate reference populations. In short, the effects of these diploid samples) from each reference population are close to
levelsofinaccuracyinthereferencepopulations(F =0.01)are optimal.
ST
relatively small. WealsoinvestigatedhowtheaccuracyofHAPMIXisaffected
Wealsorepeatedoursimulationsatl=6andl=100using20 whentheparametersusedasinput areinaccuratelyspecified (see
admixedsamplesthatweresimulatedusingYorubaandDruzebut Materials and Methods). Results of our simulations in which the
modeled using reference populations YRI and CEU, (see genome-wide ancestry proportion m was inaccurately specified
1
Materials and Methods). The squared correlation between (differentfromthevalueausedtosimulatethedata)aredisplayed
predicted and true number of European copies was 0.97 at l=6 inTable1.Weobservedthatevenifm isveryinaccurate(e.g.bya
1
and0.79atl=100,ascomparedto0.98and0.83usingaccurate factor of 4), there is no effect on results for l=6 and only a
reference populations. Thus, HAPMIX is robust to rather minimal effect (which primarily affects the genome-wide average
inaccurate (F =0.02) reference populations, and to the asym- ofHAPMIXancestryestimates,butnottheircorrelationwithtrue
ST
metric casewhere onlyonereference population isinaccurate. ancestry) for l=100. Results of our simulations in which the
PLoSGenetics | www.plosgenetics.org 9 June2009 | Volume 5 | Issue 6 | e1000519
SensitiveAncestrySegmentDetection
Figure3.AccuracyofHAPMIX,ANCESTRYMAP,andLAMP-ANCpredictionsforvariousvaluesofl,thenumberofgenerationssince
admixture.Foreachadmixturetime,resultsarebasedonanalyzing20admixedindividuals,simulatedusinganaveragegenome-wideproportionof
80%Africanand20%Europeanancestry.Foreachmethod,weplotthesquaredcorrelationbetweenpredictedandtruenumberofEuropeancopies
asafunctionofl.
doi:10.1371/journal.pgen.1000519.g003
number of generations T since admixture was inaccurately closely match the YRI and CEU reference populations, the F
ST
specified (different from the value l used to simulate the data) values between segments of inferred ancestry and the
are displayed in Table 2. We observed that even if T is very corresponding ancestral populations were equal to 0.001,
inaccurate(e.g.byafactorof2to5),thereisnoeffectonresultsfor indicating a tight correspondence (Table 3). The l=100
l=6 and only a minimal effect for l=100. Thus, HAPMIX simulations produced a similarly tight correspondence (Table 3),
appears tobe extremely robust toparametermisspecification. even though values of local ancestry could only be inferred with
Inference of ancestral populations. We are interested in moderate accuracy (Figure 3). The correspondence between
applying HAPMIX to improve our understanding of ancestral inferred ancestral segments and true ancestral populations
populations contributing to admixture events. To explore the remained reasonably tight even when the true ancestral
usefulness of the software for this purpose, we analyzed segments populations (either Mandenka and Basque, or Yoruba and
of inferred African or inferred European ancestry from our l=6 Druze) were inaccurately modeled by the reference populations
and l=100 simulations to investigate how closely they (YRI and CEU) used for inference (Table 3). Thus, HAPMIX
corresponded to the true ancestral populations used to simulate shows promise for reconstructing ancestral populations that are
admixedindividuals(seeMaterialsandMethods).Wechosetouse somewhat different fromavailable reference populations.
F ,acommonlyappliedsummarystatistic,toquantifydifferences Although the correspondence between inferred ancestral
ST
betweentheinferredandactualancestralpopulations.Inthel=6 segments and true ancestral populations is reasonably tight, it is
simulationsusingYorubaandFrenchancestralpopulations,which not perfect, with F values as large as 0.007 between inferred
ST
PLoSGenetics | www.plosgenetics.org 10 June2009 | Volume 5 | Issue 6 | e1000519