-
Population Genetic Structure of Traditional Populations in the Peruvian Central Andes and Implications for South American Population History
Molecular-based characterizations of Andean peoples are traditionally conducted in the service of elucidating continent-level evolutionary processes in South America. Consequently, genetic variation among “western” Andean populations is often represented in relation to variation among “eastern” Amazon and Orinoco River Basin populations. This west-east contrast in patterns of population genetic variation is typically attributed to large-scale phenomena, such as dual founder colonization events or differing long-term microevolutionary histories. However, alternative explanations that consider the nature and causes of population genetic diversity within the Andean region remain underexplored. Here we examine population genetic diversity in the Peruvian Central Andes using data from the mtDNA first hypervariable region and Y-chromosome short tandem repeats among 17 newly sampled populations and 15 published samples. Using this geographically comprehensive data set, we first reassessed the currently accepted pattern of western versus eastern population genetic structure, which our results ultimately reject: mtDNA population diversities were lower, rather than higher, within Andean versus eastern populations, and only highland Y-chromosomes exhibited significantly higher within-population diversities [End Page 147] compared with eastern groups. Multiple populations, including several highland samples, exhibited low genetic diversities for both genetic systems. Second, we explored whether the implementation of Inca state and Spanish colonial policies starting at about ad 1400 could have substantially restructured population genetic variation and consequently constitute a primary explanation for the extant pattern of population diversity in the Peruvian Central Andes. Our results suggest that Peruvian Central Andean population structure cannot be parsimoniously explained as the sole outcome of combined Inca and Spanish policies on the region’s population demography: highland populations differed from coastal and lowland populations in mtDNA genetic structure only; highland groups also showed strong evidence of female-biased gene flow and/or effective sizes relative to other Peruvian ecozones. Taken together, these findings indicate that population genetic structure in the Peruvian Central Andes is considerably more complex than previously reported and that characterizations of and explanations for genetic variation may be best pursued within more localized regions and defined time periods.
Mtdna, Y-Chromosome, Native Americans, Andean Populations, Amazonian Populations, Peru, Peopling of South America
[Errata]
Since at least the 1960s, characterizations of Andean peoples from a molecular perspective have been conducted peripherally, as part of larger investigative efforts in South America as a whole. These efforts have had two main goals: first, to inform debates on the peopling of South America and of the Americas writ large, and second, to contribute to broad understandings of human micro- to macroevolutionary change. Both goals were vitalized by the work of James V. Neel, Francisco M. Salzano, and their colleagues among tribal groups in Venezuela and Brazil. Thus, from the 1960s into the 1990s, researchers working in South America conducted numerous regional studies based on “classical” markers (blood groups and proteins), primarily among peoples from the Amazon and Orinoco River Basins (e.g., Matson et al. 1966; Neel 1970, 1978a, 1978b; Salzano et al. 1973, 1974, 1977; Ward et al. 1975; O’Rourke and Suarez 1985), and only occasionally among peoples in the Andean region (Best et al. 1966; Matson et al. 1966; Modiano et al. 1972).
It was not until the late 1980s that researchers began to focus explicitly on the genetic history and structure of Andean peoples. Rothhammer and Silva (1989) were the first to argue for a distinct Andean pattern of biological variation as ascertained via craniometric isoline maps; this interpretation was supported by follow-up studies of craniometric and classical marker data (Rothhammer and Silva 1990, 1992; Cavalli-Sforza et al. 1994). Beginning in the early 2000s, spatial analyses of molecular data revealed that Andean populations contain a pattern of high within- and low between-population genetic variation compared with populations in the Amazon and Orinoco River Basins, which appear to harbor much higher levels of between-group variation (Luiselli et al. 2000; Simoni et al. 2000; Tarazona-Santos et al. 2001; Fuselli et al. 2003).
Subsequently, a distinction between “west” (i.e., the Andean mountain range and surroundings regions) and “east” (i.e., Amazonian lowlands and surrounding regions) has prevailed in discussions of continental South American population genetic variation (e.g., Pucciarelli et al. 2006; Lewis et al. 2007b; Wang et al. 2007; Lewis and Long 2008; Lewis 2009; Yang et al. 2010). This dichotomized representation of genetic variation in South American populations, in turn, has had important implications for the reconstruction of the continent’s early population history. Whereas some researchers view the contrasting patterns as the outcome of two distinct initial colonizing events (Luiselli et al. 2000; Rodríguez-Delfín et al. 2001), others posit that the two regions experienced differing population histories only after the initial entry into South America through the Panamanian isthmus (Tarazona-Santos et al. 2001; Fuselli et al. 2003; Rothhammer and Dillehay 2009; Yang et al. 2010; Bodner et al. 2012).
However, the above characterization of Andean diversity may be incomplete, for three reasons. First, it derives from a small number of population samples, most of which contain few individuals. Until recently, only 14 western Andean populations have been represented in the scholarly literature; these populations are spread over more than 4,000 km, ranging from Ecuador to Peru, northwest Argentina, and northern Chile. Many of the post-1990s molecular studies are based on repeated genotyping efforts of an even smaller set of population samples (e.g., Luiselli et al. 2000; [End Page 148] Simoni et al. 2000; Tarazona-Santos et al. 2001; Fuselli et al. 2003; Battilana et al. 2006; Scliar et al. 2012; Roewer et al. 2013). Second, most relevant studies are based on either mtDNA or Y-chromosome polymorphisms, which, in isolation, contain limited information about amounts and patterns of diversity (see Lewis and Long 2008).
A third weakness of the current perspective is that views on the nature and causes of population genetic diversity within the western Andean region of South America remain underdeveloped. This is despite the fact that archaeological and ethno-historical studies attest to complex and varying human-environment interactions through time and space in the Andes (Sandweiss et al. 2001; Wernke 2007; Gosling and Williams 2013). Specifically, early patterns may have been most recently overlaid by Inca state practices and, later, by Spanish colonial-driven policies. For example, Inca rulers routinely moved individuals (males and females) and entire families, often throughout locales in the mountainous highlands [e.g., Matienzo 1967 (1567); Cieza de León 1985 (c. 1550); Wightman 1990], but also throughout the coast and areas to the south of the Lake Titicaca basin (Covey and Elson 2007). Subsequent European invasions commencing in the 1530s precipitated the rapid decline of many indigenous coastal communities, whereas highland communities had higher initial rates of survival and were used for long-distance service (tributary) activities (Cook 1980). The Spanish conquest reinforced a frontier between the highlands and the Amazon Basin, while colonial policies indirectly encouraged population mixture and migration within the Andean region.
The expected cumulative outcome of these phenomena on present-day indigenous population diversity is an overall maintenance of genetic variation combined with a pronounced lack of population structure throughout the region, particularly in the highland regions. It is therefore possible that the region’s extant pattern of population genetic variation is best explained by events from the last half millennium, instead of or in addition to early migratory phenomena.
Given the current paucity of genetic data representing Andean indigenous peoples, we may have an unclear picture of the structure of Andean population genetic diversity, including the role of postpeopling events in shaping this structure (Pucciarelli et al. 2006; Dillehay 2009). The present study reexamines Andean genetic diversity using both mtDNA and Y-chromosome data from a large sample of traditional populations from within the Peruvian Central Andes, a region that includes the coast and highlands of Peru, as well as the tropical forest that begins on the eastern slopes of the Andes and descends to the Amazon Basin (Von Hagen and Morris 1998: 14).
Our study contributes 17 population samples from multiple ecological zones (highlands, coast, and lowlands) associated with broad cultural patterns revealed by archaeological, ethnohistorical, and ethnographic studies (see D’Altroy 2015). We investigate the region’s population genetic structure to test two hypotheses: (1) In response to previous characterizations of population genetic diversity, we tested the hypothesis that populations throughout the Peruvian Central Andes in general, and in highland areas in particular, show little or no population genetic structure in the form of high within- and low between-population genetic variation relative to eastern regions of South America. (2) We tested whether extant patterns of variation within the Peruvian Central Andes can be parsimoniously explained by events from the last half millennium. Given that both Inca and Spanish policies were strongly implemented among all individuals (i.e., both males and females), primarily within highland communities, we explored (a) whether the pattern of genetic structure differs in highland versus coastal and lowland communities in the Peruvian Central Andes and (b) whether any sex biases in patterns of variation can be detected.
Materials and Methods
Population Samples
This study is based on multilocus genetic data from samples collected through an international collaborative effort within the Republic of Peru between 2001 and 2005. The senior author (A.C.S.) moved from the University of New Mexico (UNM) to Arizona State University (ASU) during the course of the project. Thus, this project was approved by institutional review boards from both UNM and ASU, as well as from the Universidad Nacional Mayor de San Marcos and the Universidad Ricardo Palma, Lima, Peru. Biological samples in the form of buccal [End Page 149] cells or blood spots were collected with informed consent from 611 individuals from 17 populations distributed throughout 9 of 25 regiones (first-level political and administrative subdivisions) and three ecological contexts in Peru (Figure 1).
Map of Peru and eastern South America (inset) with sampling locations. For the Peruvian Central Andes, coastal samples are designated with stars, highland samples with circles, and lowland samples with squares. Eastern samples from the Amazon and Orinoco River Basins are designated with Xs.
The Andes range defines a succession of ecozones running from the Pacific coastal desert, across areas of montane valleys and high tundra, and onto the humid eastern slope that descends to the Amazon Basin. Although mindful of the considerable local environmental diversity of the Andean region, we have found it useful to group our samples into three regions: coastal (areas lying to the west of the Andes), highland (areas lying in the mountain valleys), and lowland (areas lying in the humid tropics of the eastern escarpment).
Populations were assigned sample names according to the city, village, or locale from which individual biological samples were taken (Table 1, Figure 1). We adopt a geographically based nomenclature given that the samples reflect relatively vulnerable populations, therefore warranting a level of conservatism in community identification (see Obregón-Tito 2013). Each individual in the sampled populations was from one of three ecological contexts: coast (n = 78), highlands (n = 430), and lowlands (n = 103). Individuals in the coastal sample were from Caleta Santa Rosa (n = 31), Catacaos (n = 16), Islilla (n = 18), and Trujillo (n = 13). Individuals in the highlands sample were from Ancash (n = 73), Andahuaylas (n = 56), Cara Cara (n = 26), Cusco North (n = 34), Cusco South (n = 95), Huancapi (n = 13), Otuzco (n = 26), Puno (n = 69), Santiago de Chuco (n = 21), and Tupe (n = 17). The lowlands sample includes individuals from Picota (n = 21), Picota Centro (n = 22), and San Juan del Oro (n = 60). All population samples were rural indigenous communities with the exception of An-cash, which comprises recent local rural migrants to the cities of Huaraz, Chancay, and Lima or from the small town of Yungay in the región of Ancash. Inclusion criteria included unrelated individuals [End Page 150] over 18 years of age with all four grandparents from the designated sample area.
Peruvian Central Andean Populations: Sample Locations and Data Sources
To conduct our tests, we sequenced the mtDNA first hypervariable region (HVI) and determined the mtDNA haplogroup of 611 individuals, and we typed 10 Y-chromosome short tandem repeats (Y-STRs) and determined the Y-chromosome haplogroup of 272 males. These data were combined with published data for an additional 374 individuals distributed throughout 15 populations in the Peruvian Central Andes (Table 1, Figure 1). A comparative sample of 366 individuals from 10 populations from the Amazon and Orinoco River Basin regions (henceforth also designated as “eastern”) was compiled from published sources (Table 2, Figure 1).
Laboratory Methods
Each laboratory conducted DNA extractions on a subset of the total sample using a standard phenol/chloroform method (Green and Sambrook 2012). mtDNA haplogroup assignments were conducted at Universidad Nacional Mayor de San Marcos and the Universidad Ricardo Palma on nonoverlapping population samples. mtDNA sequencing of all population samples took place at UNM/ASU, and Y haplogroup assignments and STR fragment analysis took place at ASU. [End Page 151]
Amazon and Orinoco River Basin Populations: Sample Information and Data Sources
mtDNA
HVI was amplified using one of four PCR primer sets: L15996 and H16401 (Vigilant et al. 1989), L15926 (Kocher et al. 1989) and H16498 (DiRienzo and Wilson 1991), L16055 and H16410 (Handt et al. 1996), or L15936 (5’-CTTGTAAACCGGAGATGAAA-3’) and H16498 (DiRienzo and Wilson 1991). PCR products were purified using shrimp alkaline phosphatase (Affymetrix/USB) and Exonuclease I (Affymetrix/USB) and sequenced in forward and reverse directions using the BigDye protocol (Applied Biosystems) on an Applied Biosystems 3730 capillary sequencer. We typed four major Native American founder haplogroups (A–D) through the identification of diagnostic HVI polymorphisms. Haplogroup assignments were verified through PCR amplification of characteristic mtDNA markers followed by restriction enzyme typing (HaeIII, HincII, or AluI), or by the presence/absence of the COII/tRNALys nine-base-pair deletion, as described previously (Stone and Stoneking 1998). When mtDNA sequences were not diagnostic of Native American haplogroups A–D, we sequenced additional sections of the mtDNA coding region using primers as described in Ramos et al. (2009, table 1, primer sets 2, 8, and 9). Haplogroups for these samples were defined using PhylotreeMT (mtDNA tree build 16; van Oven and Kayser 2009).
Y-Chromosome
Samples were typed for 10 Y-STRs: DYS19, DYS389I, DYS389II, DYS390, DYS391, DYS392, DYS393, and DYS439 using the Powerplex Y kit (Promega), and DYS388 and DYS426 typed using primers and PCR conditions from Kayser et al. (1997) and Jobling et al. (1997), respectively. STR fragment analysis was conducted either on an ABI 377 gel sequencer or an Applied Biosystems 3730 capillary sequencer and processed using Applied Biosystem’s GeneScan and Genotyper software.
All samples were initially screened for M242, M3, and M19 SNPs that are diagnostic of Y haplogroups Q*, Q1a2a1a1, and Q1a2a1a1a, respectively. Any sample within the Q* paragroup was subsequently typed for the M346 SNP that is diagnostic of haplogroup Q1a2. If a sample did not fall into Y-chromosomal haplogroup Q, haplogroups were predicted via a Bayesian approach (Athey 2006) implemented in the Haplogroup Predictor web interface (www.hprg.com/hapest5). Haplogroup predictions were verified through direct sequencing of diagnostic SNPs. All SNP primer sets were derived from Karafet et al. (2008), except M242, M3, M415, and M420. M242 was PCR-amplified using the primer set 5’-CCTTGCTGTCTAGTTCACTAG-3’ and 5’-AATACCTTACCTAGAACAACTC-3’ and typed using the BsiHKAI restriction enzyme (New England BioLabs). M3 was amplified using mismatch primers from Santos et al. (1999) and typed using the MfeI restriction enzyme (New England BioLabs). Primer sets for M415 and M420 are from Myres et al. (2011) and Hinds et al. (2005), respectively.
Analytical Methods
Analyses of mtDNA and Y-chromosome variation of traditional populations in the Peruvian Central Andes were limited to Native American founder lineages. Y-chromosome lineages use a notation of lineage-based name followed by the diagnostic SNP (e.g., “Q1a2a1a1-M3”). Y haplogroup designations are based on the International Society of Genetic Genealogy 2014 Y-DNA Haplogroup Tree (version 9.7).
mtDNA sequences were aligned to the revised Cambridge Reference Sequence (Andrews et al. 1999) and edited using SeqMan Pro (DNASTAR Lasergene 11). We excluded nucleotide positions (np) 16182 and 16183 because these positions are dependent on the presence of a C mutation at np 16189 (Pfeiffer et al. 1999). Analyses also excluded any cytosine inserts between np 16193–and np 16194, as these may be present due to individual sequence length heteroplasmy (see Irwin et al. 2009) [End Page 152] or to differences in sequencing methods among labs. To facilitate comparisons with other published data sets, we cropped all mtDNA HVI sequences from np 16084 to 16364 and limited our Y-chromosome analyses to six STRs (DYS19, DYS389I, DYS389B, DYS390, DYS391, DYS393) within Q1a2a1a1-M3. DYS389B was calculated by subtracting DYS389I from DYS389II. Additionally, Y-chromosomes from three eastern population samples (Gavião, Suruí, Zoró) were combined into a single sample to replicate the analysis in Tarazona-Santos et al. (2001); however, we found that our results were robust compared with alternative groupings, such Gavião only, Gavião with either Suruí or Zoró, or all three separately.
To test the hypothesis that populations throughout the Peruvian Central Andes contain patterns of high within- and low between-population genetic diversity, we assessed distributions of population diversity by estimating average gene diversity within mtDNA HVI and Y-STR haplotypes (π and D estimates, respectively; Nei 1987). We further assessed levels of population genetic differentiation within ecozones through the estimation of ΦST between pairs of populations (Excoffier et al. 1992). Because previous claims of “high” within- and “low” between-population variation in the Peruvian Central Andes are based on comparisons with the same or similar measures in eastern South American populations, we conducted analyses on available population samples from the Amazon and Orinoco River Basins (see Table 2). Furthermore, because previous analyses focused exclusively on highland populations, we also conducted these tests using only Peruvian highland samples versus eastern samples. These estimates were generated using Arlequin, version 3.5.1.3 (Excoffier and Lischer 2010). Nonparametric randomization tests, in which the observed difference in average π or D estimates of populations within regions (“west” or “east”) was compared with a distribution of estimates obtained from 5,000 random permutations of the original data with replacement, were executed in Microsoft Excel, version 14.4.3.
We used two approaches to ascertain the extent to which events from the last millennium could account for extant genetic structure in the Peruvian Central Andes: (a) we explored relative patterns of the apportionment of variation among populations in different ecozones using analyses of molecular variance (AMOVA; Excoffier et al. 1992), and (b) we investigated potential effects of sex-biased gene flow on patterns of genetic variation, in two ways: genetic and geographic distances, and differences in genetic structure of males and females. We assessed correlations between genetic and geographic distances via Mantel tests (Mantel 1967; Smouse et al. 1986) using PASSaGE 2 (Rosenberg and Anderson 2011). Mantel tests were statistically evaluated using Bonferroni-corrected probabilities. For mtDNA, ΦST genetic distances were estimated using the Tamura-Nei evolutionary model (Tamura and Nei 1993) with the gamma distribution parameter alpha (α) set at 0.26 (Meyer et al. 1999). For Y-STRs, distances were estimated using a stepwise mutation evolutionary model (RST; Slatkin 1995). Genetic distance analyses were performed using Arlequin, version 3.5.1.3 (Excoffier and Lischer 2010). Great circle geographic distances were calculated from the longitude and latitude coordinates for each sample site using PASSaGE 2. We assessed differences in the genetic structure of males and females using the parameter Nν, calculated as (1/ΦST) – 1 (Cavalli-Sforza and Bodmer 1971). Different Nν values between the two genetic systems reflect differences in migration rates and/or effective population sizes, as mutation rate may be considered negligible (see Destro-Bisol et al. 2004; Bolnick et al. 2006).
Last, to visualize patterns captured by the above analyses, we used principal coordinate (PCO) plots of ΦST and RST genetic distances, as above. PCO analyses were carried out using QIIME Pipeline python scripts (Caporaso et al. 2010) and visualized using KiNG, version 2.21 (Chen et al. 2009).
Results
mtDNA and Y-Chromosome Haplogroup/Haplotype Distributions
mtDNA
The vast majority (99.3%) of individuals in the study sample belonged to Native American founder haplogroups A–D (Table 3). Three population samples (Huancapi, Picota, and Picota Centro) exhibited a 100% frequency of haplogroup B, and several others contained high frequencies of haplogroup B, such as Cara Cara (96%) and Cusco [End Page 153] North (78%). The rest of the population samples contained varying proportions of lineages A–D. Four individuals exhibited nonfounder mtDNA lineages, specifically, African mtDNA haplogroup L haplotypes (Salas et al. 2004; Behar et al. 2008; van Oven and Kayser 2009; Soares et al. 2012), and likely represent recent (post-ad 1492) admixture. These four individuals were distributed in three population samples (Andahuaylas, Catacaos, and Caleta Santa Rosa; see Supplemental Appendix S1).
mtDNA and Y-Chromosome Native American Founder Haplogroup Counts (Frequencies) in Peruvian Central Andean Population Samples
Distribution of mtDNA B Haplotype (16189C, 16217C, 16242T, 16324C) in the Peruvian Central Andes
We additionally identified a geographically uneven distribution of a single mtDNA haplotype (16189C, 16217C, 16242T, and 16324C) among our sampled populations (and in one published sample, the Tayacaja; see Table 4). Several of the lowland samples demonstrated either 100% (Picota and Picota Centro) or a very high percentage (e.g., San Juan del Oro, 43%) of this haplotype; one of our highland samples (Huancapi) also demonstrated 100% of this haplotype, and another (Cara Cara) was at 96%.
Y-Chromosome
Sixteen biallelic polymorphisms on the Ychromosome were characterized in nine population samples, resulting in 12 Y-chromosome haplogroup lineages (see Supplemental Appendix S2). Approximately 92% of males characterized exhibited one of two Native American founder haplogroups: Q1a2-M346, comprising 1.4% of the total sample, and Q1a2a1a1-M3, comprising 90% of the total sample. Other, presumably nonfounder, Y haplogroups were R1b1-M415 (9.2%), E-M96 (1.1%), [End Page 154]
mtDNA (Haplogroups A–D) HVI Summary Statistics
[End Page 155]
J-M304 (2.6%), G-M201 (1.8%), R1b-M343 (0.7%), I-M170 (0.7%), and R1a1-SRY10831.2 (0.4%). Many of these nonfounder haplogroups were distributed within the Ancash sample, followed by the Andahuaylas sample, and then were sparingly present in the rest of the samples, except for Picota and Picota Centro, which contained none (see Supplemental Appendix S2). We additionally screened R1b-M415 samples for R1b1a1-M73 (Malyarchuk et al. 2011) to investigate the possibility that these R1b lineages, typically attributed to European male admixture, might instead be of Asian origin (Bortolini et al. 2003; Kemp and Schurr 2010; Lell et al. 2002), but no sample had the diagnostic M73 SNP.
The Q1a2a1a1-M3 “modal haplotype,” 13.10.17.24.10.13, identified in previous Y-chromosome studies among Native South Americans (Bianchi et al. 1998; Tarazona-Santos et al. 2001) was absent in our sample. The most frequent (22%) Q1a2a1a1-M3 haplotype, 13.14.18.23.10.14, occurred mostly in lowland samples (Picota, Picota Centro, and San Juan del Oro) and, furthermore, was the only haplotype found in the Picota and Picota Centro samples (excluding individuals with missing data).
Patterns of Within- and Between-Population Genetic Variation
Prior characterizations of western Andean population genetic variation assume a strong distinction from eastern South American populations. In our reanalysis, the contrast between “west” and “east” is not as stark as reported previously.
In the case of mtDNA, the relative comparison between west and east is inconsistent with previous reports (i.e., Fuselli et al. 2003; Cabana et al. 2006). The western Andes showed lower, not higher, estimates of within-population diversity. In our study, average within-population π estimates for the entire Peruvian Central Andes and the Peruvian highlands are 0.0186 and 0.0192, respectively, relative to 0.0320 in the eastern region (Table 5). Comparisons of the Peruvian highland and coast with eastern populations are statistically significant at P ≤ 0.05, but the Peruvian lowland–eastern South America comparison is not (P = 0.37). In terms of population differentiation, the Peruvian Central Andes as a whole (ΦST = 0.1833) and its highland and coastal areas (ΦST = 0.1522 and 0.1170, respectively) have lower estimates than do the Amazon/Orinoco River Basins (ΦST = 0.2155). In other words, throughout the Peruvian Central Andes, and especially within the highlands, approximately 18% of the variation is apportioned among populations, versus 21% in eastern South America; these values are not substantially different from one another. The Peruvian lowland ecozone is exceptional, exhibiting relatively high population differentiation estimates (ΦST = 0.3189).
A handful of highland populations exhibited low mtDNA diversity levels. For example, Huancapi and Cara Cara demonstrated extraordinarily low diversity estimates (π = 0 and 0.0021, respectively) due to their high frequencies of the mtDNA haplogroup B haplotype (16189C, 16217C, 16242T, and 16324C) described above. To date, no other reported highland Andean population sample has exhibited such low mtDNA diversity.
In contrast to the mtDNA results, Peruvian Central Andean Q1a2a1a1-M3 Y-chromosomes exhibited higher, but statistically nonsignificant, average gene diversity estimates relative to the eastern region (D = 0.3601 vs. 0.3007, respectively; P = 0.16; Table 6); the Peruvian highland area showed a much higher and significant relative value (D = 0.3936; P = 0.03). Population differentiation was lower in the Peruvian Central Andes overall (ΦST = 0.4147), as well as in the highlands (ΦST = 0.4042) and lowlands (ΦST = 0.2500), compared with the eastern estimate (ΦST = 0.6265). This indicates that the eastern region contains about 20% more variation apportioned among populations relative to the Peruvian Central Andes. Interestingly, however, highland population Y-chromosomes were much more differentiated than found previously (Tarazona-Santos et al. 2001; Yang et al. 2010). For example, Tarazona-Santos et al. (2001) report a pooled Central Andean ΦST estimate of 0.024, but our comparative analyses demonstrated a much higher estimate for highland populations, as well as those within the entire Peruvian Central Andes.
Population Genetic Structure by Ecozone
To address the second hypothesis, that variation within the Peruvian Central Andes can be explained by events from the last half millennium, we first sought to assess the degree to which highland populations were differentially affected by Inca and Spanish policies relative to lowland and coastal populations, if at all. The hypothesis is [End Page 156] supported for mtDNA but not for Y-chromosome data. AMOVA analyses for both genetic systems show minimal, often nonsignificant, differences in the apportionment of variation among ecozones. When within- and between-population apportionments of variation were analyzed independently for each ecozone, highland groups had higher amounts of within-group variation for mtDNA but not for Y-chromosome data (Table 7).
Y-Chromosome Q1a2a1a1-M3 Summary Statistics
More specifically, mtDNA-based AMOVA analyses of multiple comparisons among ecological zones showed most (~67–82%) genetic variation distributed within populations; the little genetic variation that was apportioned among zones rarely reached statistical significance after 10,000 permutations (Table 7). This general lack of structure by ecological zone was apparent in the mtDNA PCO plot depicting the first two PCOs (Figure 2A) in which the first and second PCOs captured 46% and 24% of the total variance in the data, respectively. Highland, coastal, and lowland samples were intermixed. Four samples that form a somewhat [End Page 157] distant cluster along the first PCO axis—Cara Cara, Huancapi, Picota, and Picota Centro—all share very high frequencies of the mtDNA haplogroup B haplotype (16189C, 16217C, 16242T, and 16324C) mentioned above. A comparison of the apportionment of variation within and between populations by ecozone shows that a large proportion of genetic variation in coastal and highland zones is found within (~85–88%) rather than among (12–15%) populations, whereas in the lowlands, only 68% was distributed within, and up to 31% among populations.
Analyses of Molecular Variance (AMOVA)
Y-chromosome-based AMOVA analyses showed genetic variation also as weakly and insignificantly apportioned between highland and lowland ecozones, with a little more than half (~55%) of genetic variation distributed within populations and the rest among populations within ecozones. The Y-chromosome PCO plot visually depicts the general lack of structure between highland and lowland samples (Figure 2B). The first and second PCO captured 47% and 24% of the total variance in the data, with highland and lowland samples for the most part intermixed along both axes, with the exception of two clusters of samples. One cluster, consisting of Cara Cara, Picota, and Picota Centro, shares high frequencies of a single Y haplotype (13.14.18.23.10.14) within Q1a2a1a1-M3. A second cluster, consisting of Tayacaja and Arequipa, shares many haplotypes in common to the exclusion of most or all other Peruvian Central Andean populations. Finally, as pointed out above, much more variation is apportioned among, versus within, groups in the highlands relative to those in the lowlands.
Sex-Specific Patterns of Population Genetic Structure
A secondary hypothesis related to the potential effect of Inca and Spanish policies was an expected lack of sex bias in the extant pattern of population genetic variation throughout the Peruvian Central Andes, and in the highlands especially. This secondary hypothesis is rejected: though Mantel tests lacked evidence of distance-based patterning of either maternal or paternal genetic variation, Nν values between the two genetic systems reveal higher female effective sizes and/or gene flow in the Peruvian highlands; the converse is true for the Peruvian lowlands.
Mantel tests yielded low and statistically insignificant correlations between genetic and geographic distances for both genetic systems. For mtDNA, the Mantel statistic was low (r = 0.1555; critical P-value = 0.001). For the Y-chromosome, the Mantel statistic was slightly negative (r = –0.5142; critical P-value = 0.003). Both PCO plots depict this lack of correspondence: although pairs of spatially proximate population samples cluster together somewhat along either the first or second PCO axis, relationships among population samples generally do not readily relate to geographic location (Figure 2).
Maternally and paternally inherited loci showed opposite patterns of genetic differentiation (ΦST) among populations (Tables 5, 6). mtDNA indicates high differentiation in the lowlands (ΦST = 0.32) and approximately equal levels in the highlands (ΦST = 0.15) and the coast (ΦST = 0.12), whereas Y-chromosomes indicate low differentiation in the lowlands (ΦST = 0.25) and high differentiation in the highlands (ΦST = 0.40). The ΦST values for the Peruvian Central Andes produce a ratio of mtDNA to Y-chromosome Nν of 3.05, indicating a female migration rate and/or effective size that is [End Page 158] more than three times that of males. In the Peruvian highlands, the ratio is even larger (Nν = 3.78), whereas in the lowlands it is lower (Nν = 0.71), suggesting the male migration rate and/or effective size is almost 1.5 times higher than that of females.
Discussion
This study presents a comprehensive investigation of population genetic variation in the Peruvian Central Andes with the purpose of contributing to the reconstruction of the region’s population history. To this end, the study analyzes maternally inherited mtDNA HVI and paternally inherited Y-chromosome STRs of 17 novel population samples combined with 15 published population samples. With this large number of samples covering multiple ecological areas throughout the Peruvian Central Andes, we revisit the long-standing characterization of Andean population genetic variation as relatively unstructured compared with eastern South American populations, and we assess the extent to which state expansion and colonialism during the last half millennium might have differentially affected highland communities and/or generated a lack of sex bias in population genetic structure.
As noted in the introductory remarks, an overall pattern of high within- and low between-population variation in the Andes relative to eastern South America has been a consistent finding across multiple genetic-based studies. Our results reject a consistent pattern of relatively high within-and low between-population variation for both genetic systems. For the Peruvian Central Andes, mtDNA HVI within-population diversities were low, rather than high, relative to values for eastern populations, and only highland Y-chromosomes exhibited significantly higher within-population diversities compared with eastern Y-chromosomes. Also in contrast to previous studies, we found several population samples within the Peruvian highlands, as well as throughout the Peruvian Central Andes, with relatively low mtDNA genetic diversity estimates. These results are partly due to the fact that the characterization of eastern population structure differs from previous reports because of our inclusion of novel samples. Eastern population diversities for both genetic systems are higher than, and eastern Y-chromosomes are twice as differentiated as, previously reported (Tarazona-Santos et al. 2001; Fuselli et al. 2003; Cabana et al. 2006).
First two principal coordinates based on mtDNA HVI pairwise ΦST (A) and Y haplogroup Q1a2a1a1-M3 pairwise RST (B) distance estimates. Axes are labeled with the percentage of variance explained. Coastal population samples are designated with stars, highland samples with circles, and lowland samples with squares. Abbreviations: AMT, Amantaní; ANC, Ancash; AND, Andahuaylas; ARQ, Arequipa; CAR, Cara Cara; CAT, Catacaos; CHM, Chimú; CJM, Cajamarca; CPC, Capachica; CSN, Cusco North; CSR, Caleta Santa Rosa; CSS, Cusco South; CUS, Cusco; HCP, Huancapi; HNC, Huancavelica; ISL, Islilla; MAT, Matsiguenga; OTZ, Otuzco; OXP, Oxapampa; PCC, Picota Centro; PCT, Picota; PUN, Puno; SCH, Santiago de Chuco; SJO, San Juan del Oro; SMP, San Martín de Pangoa; TAY, Tayacaja; TQL, Taquile; TRJ, Trujillo; TUP, Tupe; URO, Los Uros.
Up until this point, the characterization of a relative lack of population genetic structure in the Andean region versus eastern South America has gone unchallenged. Instead, recent work has built [End Page 159] on the presumed weak structure to argue against a dual founder model for the peopling of South America (Lewis and Long 2008; Lewis 2009) and to argue for the possibility of range expansion by a single founding population into the Central Andes (Yang et al. 2010; Batai and Williams 2014). However, by simply increasing the number and geographical coverage of populations sampled, we show that populations on either side of the Andean range do not have contrasting—or even opposing—patterns of genetic variation and thus cannot be easily characterized as two distinct “meta-” (Lanata and García 2005) or “mega-” (Pucciarelli et al. 2006) populations. This is the case even if we confine our definition of “western Andes” to highland regions only.
Explanations for the presumed dichotomous pattern have focused on migratory events in early South American prehistory, presumably because large-scale migrations of distinct founder populations could reasonably explain why such clear and distinctly opposing patterns of population structure were observed over such large, well-defined geographic areas several millennia later. The fact that this explanation has dominated discussions, however, means that the role of other or additional (pre)historical phenomena in shaping extant patterns of population genetic variation in the Andes has been underexplored. In order to forward that exploration, we chose to focus on the Inca imperial and Spanish colonial policies from the last half millennium because (1) they are recent phenomena relative to the region’s ~15,000 years of prehistory (Dillehay 2009), and (2) they generate clear expectations of the structure of population genetic diversity.
Thus, our second hypothesis investigates population structure solely within the Peruvian Central Andean region, by ecozone and/or by genetic system. Archaeological and ethnohistorical accounts led us to expect weak signals of population genetic structure and a lack of sex bias in gene flow and/or effective size estimates, particularly in the Peruvian highlands relative to the coast and lowlands. Our results show that Peruvian highland communities do exhibit weak population genetic structure relative to other ecozones only for mtDNA data, and not for Y-chromosome data. Second, highland populations exhibit higher levels of female- versus male-mediated gene flow and/or high effective sizes relative to other regions, whereas the lowlands show the opposite pattern. Though our results provide strong evidence for high levels of gene flow and/or effective sizes among females in Peruvian highland groups, they do not fully support our expectations for both genetic systems. These results therefore provide, at best, equivocal evidence for a strong influence of combined Inca and Spanish policies on population genetic variation throughout the Andes or in highland regions only.
Population genetic structure in the Andes does not appear to be dictated by large-scale phenomena and may instead be understood via an exploration of more spatiotemporally defined processes. Archaeological, ethnohistorical, and recent paleogenetic studies of the Central Andes encourage a finer-grained analysis. These studies show early geographic and diachronic variation in the area: from the time of initial colonization, hunter-gatherers engaged in diverse subsistence practices with varying patterns of mobility, including sedentism, transhumance, and long-distance foraging—practices that probably isolated some human populations regionally (e.g., Aldenderfer 2008; Sandweiss 2008). As food production became central to the subsistence economies and social organization of most Andean societies (particularly after 3000 bc), groups increasingly engaged in contacts that spread crops, craft goods, and cultural practices within and among regions (Perry et al. 2006; Haas et al. 2013). Andean urban centers emerged around ad 400–600, signaling the presence of centralized and hierarchical societies with subsistence strategies favoring rapid population growth rates (Read and LeBlanc 2003). State expansion, in which central political control was exercised over large territories, occurred in multiple locations. In the case of Wari and Tiwanaku, state populations also established colonies. Wari and Tiwanaku disintegrated around ad 1000, replaced by agropastoral societies that used broad kin-based networks to access diverse resources. Populations grew markedly across the highlands until around ad 1400, when imperial expansion of the Inca state likely altered the region’s demography. Spanish documents indicate that by the mid-sixteenth century, the region displayed diverse ethnic groups, languages, and dialects (Mannheim 1991).
Moreover, paleogenetic studies suggest that [End Page 160] local dynamics may have deviated from broader-scale processes over the last 500 years in the Andean region. Studies by Baca et al. (2014), Kemp et al. (2009), and Lewis et al. (2007a) present evidence for spatiotemporal continuity of populations, whereas Fehren-Schmitz et al. (2011) found evidence for both continuity and discontinuity of populations that differ between maternally and paternally inherited loci. Notably, our study’s identification of an uneven geographic distribution of a single mtDNA B haplotype hints at localized instances of gene flow and drift.
Conclusion
Over the last two decades, studies of continent-level population genetic variation and structure have been executed under the presumption that either a single- or dual-founder colonization model would best explain observed patterns. Our study upends this research trend by rejecting the presumed dichotomous pattern of South American population genetic structure.
Importantly, our study contributes greater nuance to the “western” pattern. We have shown that with denser population and genomic-level sampling, the “west” reveals finer-scale patterns of population genetic variation than formerly recognized. These patterns may result from relatively localized and historically contingent processes, rather than from single major trends or causal variables. We suggest that future endeavors in the Peruvian Central Andes strive for denser geographic and genomic coverage using both ancient and modern DNA techniques to better access local (pre)histories at multiple spatial and temporal scales.
These associate editors of Human Biology were excluded from the editorial process for and any decisions by the journal about this article, to minimize conflict of interest.
These associate editors of Human Biology were excluded from the editorial process for and any decisions by the journal about this article, to minimize conflict of interest.
Supplemental Appendices
Complete population-level data for mtDNA HVI and Y-chromosome STRs are provided in Supplemental Appendixes S1 and S2, respectively. Supplemental Appendix S1 (available from http://goo.gl/j2bgf8) provides mtDNA haplogroup and full (uncropped) sequence data, along with four spreadsheets summarizing sequence variants by founder haplogroups A–D. Supplemental Appendix S2 (available from: http://goo.gl/Sa345l) provides Y-chromosome haplogroup and haplotype profiles based on 10 STRs. mtDNA HVI sequences are also deposited in GenBank (KM888878–KM889488).
acknowledgments
We are grateful to the peoples in the Peruvian Central Andes who agreed to participate in this study, and to Mariana Contreras, Leslie Castillo Pampas, and Julian Ccahuana Quispe for help with sample collection. We thank Daniel Mendoza and Gabriel Ríos of the Instituto de Ciencia y Tecnología at Universidad Ricardo Palma, Lima, Peru, for their professional support, and Keith L. Hunley, Elizabeth A. Sobel, and three anonymous referees for insightful comments. The research reported here was supported by the National Science Foundation through grants BCS-0401434 (A.C.S.), BCS-0242958 (A.C.S./C.M.L.), and SES-0309904 (G.S.C.).
