Genome-wide data across the Bronze/Iron Age transition in Ashkelon

We extracted and sequenced DNA from 108 skeletal elements excavated in Ashkelon. In line with the low DNA preservation previously reported for the southern Levant (11–14), only 10 yielded sufficient amounts of human DNA (data file S1). Sequencing libraries for these 10 individuals were enriched for human DNA using an established in-solution DNA capture targeting 1.24 million genome-wide single-nucleotide polymorphisms (SNPs) (“1240 k capture”) (15–17). The earliest three individuals were excavated from a Bronze Age necropolis (18–20) and dated by the archeological context to the Middle Bronze IIC–Late Bronze Age II (labeled ASH_LBA; two of the individuals were directly carbon-dated to 1746 to 1542 cal BCE). Four early Iron Age infants were recovered from burials beneath late 12th century Iron Age I houses (labeled ASH_IA1; directly carbon-dated to 1379 to 1126 cal BCE), and three later Iron Age individuals were recovered from a cemetery adjacent to the city wall of ancient Ashkelon, which was estimated to have been used between the 10th and the 9th century BCE (labeled ASH_IA2; one individual directly carbon-dated to 1257 to 1042 cal BCE) (Fig. 1, A and B, Table 1, table S1, and text S1). While the chronological ranges estimated by the archeological context are approximately within the range of the radiocarbon assays conducted directly on the petrous bones sampled for DNA (table S1), we caution that they should be taken as a rough estimate, since the presence of marine carbon in the diet could make the calibrated radiocarbon dates older than the true age.

To control for the quality of the dataset, we estimated exogenous DNA levels and relatedness. Mitochondrial contamination for all 10 individuals ranged between 0 and 9%, with only two of them above 5% (table S2). For males (n = 4), nuclear contamination was estimated to a maximum of 4.3% (table S2). None of the individuals were either first- or second-degree related to each other (fig. S1). We performed population genetic analysis on a merged dataset, including the genotype data of the newly reported individuals from this study and previously published datasets from 638 ancient individuals and 4943 individuals belonging to 298 present-day populations. (11, 12, 21). (An overview of the key ancient individuals analyzed is given in data file S2).

Persistence of the local gene pool in the Bronze Age Levant

To understand the genetic profile of the ancient Ashkelon individuals, we began by projecting them onto the first two axes of variation (PC1 and PC2) of present-day west Eurasians, inferred from principal component analysis (PCA) (Fig. 2A and fig. S2) (12). The ASH_LBA individuals overlap with the cline of present-day Near Easterners and are close to earlier Bronze Age Levantines and Anatolians [Early Bronze Age individuals from ‘Ain-Ghazal, Jordan labeled “Jordan_EBA” (12); Middle Bronze Age individuals from Sidon, Lebanon labeled “Lebanon_MBA” (13); and Early to Late Bronze Age individuals from Central Anatolia labeled “Anatolia_EBA” (21) and “Anatolia_MLBA” (21), respectively]. Compared to earlier Levantines [Levantine early farmers from present-day Jordan and Israel labeled “Levant_N” (12) and Chalcolithic individuals from Peqi’in, Israel labeled “Levant_ChL” (11)], the Bronze Age individuals including ASH_LBA are all shifted along PC2 toward ancient Iranian and Caucasus individuals [e.g., the Caucasus Mesolithic hunter-gatherers labeled “CHG” (22); early farmers excavated from present-day Iran labeled “Iran_N” (12); and the Chalcolithic individuals excavated from present-day western Iran labeled “Iran_ChL” (12)], in agreement with previous observations (11–13). Unsupervised genetic clustering using ADMIXTURE (23) shows a similar pattern [with K = 9 clusters; Fig. 2B and fig. S3]; ASH_LBA is assigned the same major ancestral component as all earlier Levantine populations (shown in orange). Consistent with their PCA positions, a second component (shown in green) that is maximized in Iran_N is, on average, higher in ASH_LBA and in each of the earlier Bronze Age Levantines compared to all earlier Levantines (20 to 30% and 3 to 8%, respectively).

To formally test the qualitative observations based on the PCA and ADMIXTURE analyses, we used f4-statistics (24). Consistent with the ASH_LBA positions on PC2, the f4-statistic of the form f4 (ASH_LBA, Levant_N/Levant_ChL; test, Mbuti) (figs. S4 and S5) estimated excess allele sharing between ASH_LBA and ancient Iranian/Caucasus-related populations (such as CHG, Iran_N, and Iran_ChL) compared to the Neolithic/Chalcolithic Levantines (Levant_N and Levant_ChL), confirming previous reports of post-Neolithic gene flows from prehistoric populations related to Iran or the Caucasus into the Levant (11–13). Accordingly, modeling ASH_LBA as a two-way admixture using qpAdm (12, 17) produced a fitting model (χ²P = 0.445), in which ASH_LBA derives around 60% of their ancestry from Levant_ChL (60.0 ± 6.5%; estimate ±1 SE) and the rest from Iran_ChL. The spatiotemporal origins of this post-Neolithic gene flow remain broadly defined since alternative two-way models also fit, albeit with smaller P values. In these models, Anatolia_EBA/Anatolia_BA/Armenia_ChL are used as the second source, replacing Iran_ChL (χ²P = 0.053 to 0.060; table S3).

We then examined whether there are noticeable genetic differences between ASH_LBA and the earlier Bronze Age Levant populations by measuring f4 (ASH_LBA, Jordan_EBA/Lebanon_MBA; test, Mbuti (figs. S6 and S7). We find that ASH_LBA and both earlier Bronze Age populations are symmetrically related to all test populations within our data’s resolution (|f₄| < 3 SE). However, we observe a marginal excess of affinity between ASH_LBA and populations genetically related to ancient Iran and the Caucasus (such as CHG and Steppe_Eneolithic) when compared to Jordan_EBA (fig. S7). Within our current resolution, we cannot distinguish whether this affinity is due to a small-scale gene flow from an Iranian/Caucasus-related population entering the Levant between the Early and Late Bronze Age or whether it reflects a certain population structure during this period. Nonetheless, the apparent symmetry suggests a high degree of genetic continuity throughout at least a millennium between culturally and geographically distinct Bronze Age Levantine groups occupying a region that stretches from the inland southern Levant where present-day Jordan is located and along the coastal regions of present-day Israel and Lebanon.

Genetic discontinuity between the Bronze Age and the early Iron Age people of Ashkelon

In comparison to ASH_LBA, the four ASH_IA1 individuals from the following Iron Age I period are, on average, shifted along PC1 toward the European cline and are more spread out along PC1, overlapping with ASH_LBA on one extreme and with the Greek Late Bronze Age “S_Greece_LBA” (25) on the other (Fig. 2A). Similarly, genetic clustering assigns ASH_IA1 with an average of 14% contribution from a cluster maximized in the Mesolithic European hunter-gatherers labeled “WHG” (shown in blue in Fig. 2B) (15, 22, 26). This component is inferred only in small proportions in earlier Bronze Age Levantine populations (2 to 9%). In light of these observations, we formally tested the variation within each of the three Ashkelon populations by computing the variance of the square rooted statistic f4 (Ashkelon individual 1, Ashkelon individual 2; test, Mbuti) (fig. S8). Consistent with the wide spread of ASH_IA1 in the west Eurasian PCA, the variance is significantly higher in ASH_IA1 than in either ASH_LBA or ASH_IA2 (P < 2.2 × 10⁻¹⁶ for both by a Fligner-Killeen test) (table S4 and fig. S8), suggesting that the ASH_IA1 individuals were more heterogeneous in their genetic affinities to global populations in respect to both ASH_LBA and ASH_IA2.

We next formally quantified the genetic differences between ASH_IA1 and the earlier ASH_LBA by calculating f4 (ASH_IA1, ASH_LBA; test, Mbuti) (Fig. 3A). In agreement with the PCA and ADMIXTURE results, only European hunter-gatherers (including WHG) and populations sharing a history of genetic admixture with European hunter-gatherers (e.g., as European Neolithic and post-Neolithic populations) produced significantly positive f4-statistics (Z ≥ 3), suggesting that, compared to ASH_LBA, ASH_IA1 has additional European-related ancestry. To estimate the levels of the European-related ancestry in all Ashkelon populations, we compared qpAdm models of two-way mixtures (Levant_ChL and Iran_ChL; i.e., 0% contribution from WHG) to three-way ones, in which we add WHG as the third source (Fig. 3C and table S5). ASH_LBA and ASH_IA2 fit well with the two-way model (χ²P = 0.445 and χ²P = 0.313, respectively; table S5), whereas the three-way one infers small nonsignificant proportions (−2.3 ± 2.2 and 2.0 ± 2.2% for ASH_LBA and ASH_IA2, respectively; table S5). In contrast, for ASH_IA1, the two-way model is inadequate (χ²P = 3.80 × 10⁶; table S5) and the three-way one fits (χ²P = 0.765). Thus, of the three, additional WHG-like ancestry is only necessary to model ASH_IA1.

We find that the PC1 coordinates positively correlate with the proportion of WHG ancestry modeled in the Ashkelon individuals (fig. S9 and table S6), suggesting that WHG reasonably tag a European-related ancestral component within the ASH_IA1 individuals. However, these Mesolithic individuals are unlikely to be a good proxy for the true source in the much later early Iron Age. To examine more proximate sources, we compiled a set of chronologically and geographically relevant candidate populations, including populations that shared higher affinity with ASH_IA1 compared to ASH_LBA in the above f4-statistic. Subsequently, we modeled ASH_IA1 as two-way and three-way mixtures of ASH_LBA and combinations of the candidate populations (table S7). Of the 51 tested models, we find four plausible ones (χ²P > 0.05), all are two-way mixtures. The best supported one (χ²P = 0.675) infers that ASH_IA1 derives around 43% of ancestry from the Greek Bronze Age “Crete_Odigitria_BA” (43.1 ± 19.2%) and the rest from the ASH_LBA population. ASH_IA1 could also be modeled with either the modern “Sardinian” (35.2 ± 17.4%; χ²P = 0.070), the Bronze Age “Iberia_BA” (21.8 ± 21.1%; χ²P = 0.205), or the Bronze Age “Steppe_MLBA” (15.7 ± 9.1%; χ²P = 0.050) as the second source population to ASH_LBA. To check whether these results are due to the low coverage of ASH_LBA, we repeated this analysis, but this time, we modeled ASH_IA1 as a three-way mixture of each of the candidate populations, Levant_ChL and Iran_ChL. The two latter populations have higher genome coverage and can model ASH_LBA well in combination (table S3). In this analysis, only the models including “Sardinian,” “Crete_Odigitria_BA,” or “Iberia_BA” as the candidate population provided a good fit (χ²P = 0.715, 49.3 ± 8.5%; χ²P = 0.972, 38.0 ± 22.0%; and χ²P = 0.964, 25.8 ± 9.3%, respectively). We note that, because of geographical and temporal sampling gaps, populations that potentially contributed the “European-related” admixture in ASH_IA1 could be missing from the dataset. Therefore, better proxies might be found in the future when more data is available. Nonetheless, the tested candidate populations from Anatolia, Egypt, and the Levant that did not produce well-fitting models can be excluded as potential sources of the admixture observed in ASH_IA1.

The transient impact of the “European-related” gene flow on the Ashkelon gene pool

The ASH_IA2 individuals are intermediate along PC1 between the ASH_LBA ones and the earlier Bronze Age Levantines (Jordan_EBA/Lebanon_MBA) in the west Eurasian PCA (Fig. 2A). Notably, despite being chronologically closer to ASH_IA1, the ASH_IA2 individuals position closer, on average, to the earlier Bronze Age individuals. Such a reduced affinity to the European populations is also apparent in the genetic clustering results, showing that the European hunter-gatherer–related ancestry contributes considerably less to ASH_IA2 than to ASH_IA1 (8 and 14%, respectively; Fig. 2B). To test for differences in affinities between ASH_IA2 and ASH_IA1, we measured the f4-statistic of the form f4 (ASH_IA2, ASH_IA1; test, Mbuti). This statistic produced no significantly positive results, and all significantly negative statistics (Z ≤ −3) are with test populations from either Europe or Anatolia (fig. S10). This, again, highlights the increased European-related affinity of ASH_IA1.

The transient excess of European-related genetic affinity in ASH_IA1 can be explained by two scenarios. The early Iron Age European-related genetic component could have been diluted by either the local Ashkelon population to the undetectable level at the time of the later Iron Age individuals or by a gene flow from a population outside of Ashkelon introduced during the final stages of the early Iron Age or the beginning of the later Iron Age. Considering the symmetry measured between ASH_IA2 and ASH_LBA in the f4-statistic of the form f4 (ASH_IA2, ASH_LBA; test, Mbuti) (|Z| < 2.8; fig. S11), the replacing population is likely to have stemmed from the Late Bronze Age Levantine gene pool, whether or not from Ashkelon. By modeling ASH_IA2 as a mixture of ASH_IA1 and earlier Bronze Age Levantines/Late Period Egyptian, we infer a range of 7 to 38% of contribution from ASH_IA1, although no contribution cannot be rejected because of the limited resolution to differentiate between Bronze Age and early Iron Age ancestries in this model (table S8).