Importantly, Y. pestis was not found in other nearby contemporary hunter-gatherers from the Pitted Ware Culture in Sweden (). Altogether, our results support that an epidemic of plague could explain the size and short time frame of the Frälsegården passage grave. The emergence of epidemics in ancient farming villages, which had higher human and animal densities than those of forager groups, had been already proposed based on archaeological records found in different pre-historic settlements, but our findings represent a new and clearer evidence of such events ().

To determine if other individuals from this passage grave or from the same geographic region showed any sign of Y. pestis, we carefully analyzed the reads mapping to the Y. pestis reference genome from individuals that had been also sequenced in previous works (). The analysis revealed the potential existence of Y. pestis in another early farmer individual (Gökhem4, a 20-year-old male, 5,040–4,839 BP) from the TRB culture in the same megalithic burial structure and age ( Figure 1 A; Table S2 ). Here, 2,951 reads mapped to the chromosome and three plasmids. The low number of Y. pestis reads recovered is also consistent with the low coverage of the human genome in this sample (0.04X, versus the 1.33X of Gökhem2) (), indicating a poor overall DNA conservation. To validate their identity, we obtained perfect or near perfect alignment versus non-redundant sequences from NCBI for 1,764 of the reads where 1,116 of these were exclusive to the Y. pestis chromosome and the three plasmids ( Table S3 ). Additionally, the Naive Bayesian classifier classified the reads with a posterior probability of 1 for being Y. pestis ( Figure 1 D). These results confirmed that despite the poor conservation of the sample, there is also an unambiguous presence of Y. pestis in Gökhem4. We can rule out the possibility of bleed-over between lanes and sequencing indexes as there were multiple sequencing libraries from the two individuals, many of which were prepared and sequenced at different times ().

The screening of known human pathogens in public ancient DNA datasets of teeth from individuals found in the Frälsegården passage grave in Sweden () revealed the unambiguous presence of Y. pestis in the “Gökhem2” individual (a.k.a. individual A), a 20-year-old female, dated to 5,040–4,867 BP and belonging to the Funnel Beaker culture (TRB) ( Figure 1 A; Table S1 ). The alignment of all sequencing reads to the Y. pestis genome resulted in 203,733 unique reads recovered, covering 79.3% of the chromosome at 2.7X and 79.0% (14.2X), 83.9% (4.7X), and 58.7% (1.7X) of the pPCP1, pCD1, and pMT1 plasmids, respectively ( Figure S1 ). The lower coverage of the pMT1 is due to the absence of a 20-kb region containing the ymt gene, which was gained after the Bronze Age clade diverged from the main lineage ( Figure S1 ) (). The reads showed the expected degradation signatures for ancient DNA (aDNA), such as, short read length distribution ( Figure 1 B) and deamination at the 5′ and 3′ ends of the sequences ( Figure 1 C). To confirm the identity of the strain (named Gok2 from now), we also re-mapped the reads against the closely related Yersinia pseudotuberculosis IP32953 reference genome and calculated the mapping affinity against both species, as done previously (). Reads showed to be clearly closer to Y. pestis ( Figure 1 D), which was also confirmed using a Naive Bayesian classifier that showed a posterior probability of 0.93 for the reads originating from Y. pestis.

High quality reads mapping across the Y. pestis chromosome and the three plasmids (pMT1, pPCP1 and pCD1). Outer ring is map-ability (gray) and genes are marked as follows: RNA: black, transposon: purple, positive strand: blue, negative strand: red. Read depth of Gok2 (brown), RISE509 and RISE505 (blue), DA101/Tian Shan plague (black), Justinian plague (orange), Black Death (purple), modern Y. pestis D1982001 (green) and Y. pseudotuberculosis IP32881 (red). A red arrow indicates the position of the ymt gene within a missing region (also marked in red) of the pMT1 plasmid in the Bronze Age and Gok2 strains.

(E) Whole genome phylogeny reconstructed using maximum likelihood of modern (163, in gray) and ancient (20, in color) Y. pestis strains. The tree was reconstructed using RAxML from an alignment of 3,558 Y. pestis genes with codon-based partitioning. Bootstrap support is shown at the nodes as a fraction of 100 bootstrap replicates generated using RAxML and hereafter reconstructed using RAxML. The tree is rooted using Y. pseudotuberculosis (not shown) and Gok2 was found to be basal to all known Y. pestis strains.

(D) Distribution of edit distances of high quality reads from Gok2 and Gok4 mapped to either Y. pestis (dark gray) or Y. pseudotuberculosis (light gray) reference genomes. Reads have a higher affinity to Y. pestis than to the closest related species Y. pseudotuberculosis.

(B and C) The expected ancient DNA degradation patterns of the Gok2 strain showing short read distribution (B), and increased C > T and G > A substitutions in the 5′ and 3′ ends of reads (that are eliminated in downstream analyses) (C).

(A) Archaeological sites and carbon dating (in years before present) of the individuals infected with Gok2/Gok4 and Bronze Age Y. pestis strains. Samples are colored as in the phylogeny (E) or in gray when they had low genome coverage.

To determine the phylogenetic positioning of the Gok2 Y. pestis strain, we reconstructed the phylogeny of this genome together with a collection of Y. pseudotuberculosis (n = 27), and Y. pestis (n = 183) genomes. Using high confidence positions (3.2 million base pairs [Mbp]) and maximum likelihood phylogenetic analysis, we found the Gok2 strain to be basal to all other Y. pestis ( Figures 1 E and S2 ). The branch was supported by 100% bootstrap, independently of using Y. pestis or Y. pseudotuberculosis as reference for genome reconstruction. Importantly, the Gok2 strain was clearly part of the Y. pestis clade and well separated from Y. pseudotuberculosis. Because of the relatively low depth of the genome and the possible confounding effect of DNA damage (i.e., C to T and G to A substitutions), we additionally reconstructed the phylogeny excluding transitions. Here, Gok2 was also clearly basal to the Bronze Age clade ( Figure S2 E). Moreover, the topology of the tree shows that Gok2 is not part of the same monophyletic clade that groups Bronze Age strains ( Figure 1 E), in which the most basal strain is RISE509, that was found in the Altai (Siberia) more than 4,500 km from Gok2 and carbon dated to 4,800–4,700 BP (i.e., slightly younger compared to Gok2).

Full maximum likelihood phylogenetic trees reconstructed from whole-genome information of modern (163) and ancient (20) Y. pestis strains and 27 selected Y. pseudotuberculosis strains. Whole genomes were reconstructed by mapping the reads on either the Y. pestis CO92 genome (A) or the Y. pseudotuberculosis IP32953 (B) genome and phylogenetic reconstructions were made using the same parameters in both cases. Trees excluding Y. pseudotuberculosis strains for better clarity of the topological distribution of Y. pestis is shown in (C) and (D) using Y. pestis and Y. pseudotuberculosis as reference, respectively. The phylogenetic tree was also determined using only transversion (i.e., not using C > T and G > A deamination that can be produced through ancient DNA damage) using Y. pseudotuberculosis IP32953 as reference genome (E). Strains are colored according to major Y. pestis clades ( Table S4 ).

Altogether, our results indicate not only that Gok2 and the other Bronze Age strains belonged to two independent lineages, but also both lineages are extinct and no longer represented among the modern strains of Y. pestis.

Complementary, we used the same strategy to identify SNVs associated with the Bronze Age clade. We identified 154 variants found in at least one of the strains in this clade. Except nine variants that were found in most Y. pestis strains, including Gok2, none of the mutations were found in neither Gok2 nor any other Y. pestis strain ( Figure 2 B). Among these lineage-specific mutations, 44% were transversions and 54% were found in more than one strain, indicating that these are truly specific of the Bronze Age lineage.

To further understand the genomic variation that differentiated the Gok2 and Bronze Age lineages, we performed a detailed analysis of the single nucleotide variants (SNVs) found in both groups ( Figure 2 ). We first identified all SNVs found in Gok2, relative to the Y. pseudotuberculosis reference, and not fixed in all Y. pestis strains (i.e., identical in all genomes). This resulted in a total of 28 SNVs (11 missense, 14 synonymous, and 3 non-annotated), of which 10 were present in almost all Y. pestis strains, 3 were located at positions that were eventually deleted in most or all Y. pestis genomes, and 15 variants were exclusive to Gok2 compared to all other analyzed strains ( Figure 2 A). Among all 28 SNVs, only 2 were transversions (1 synonymous and 1 non-coding) and 8 were also observed in some Y. pseudotuberculosis strains, which may indicate that they likely did not emerge in Y. pestis. It is worth noting that the 15 unique SNVs from Gok2 (9 missense, 4 synonymous, and 2 non-coding) were supported by 4–8 independent reads, but were all transitions (C > T, G > A), which can potentially result from deamination due to ancient DNA degradation. However, variation predicted to be caused by DNA damage, which are normally observed in the 5′and 3′ends of DNA molecules, were filtered out from reads. All these variants were called from bases located at the center of the reads ( Figure S3 ), strongly suggesting that these SNVs were not due to post-mortem degradation. Interestingly, four of the variants were located in genes with functions related to host-pathogen interactions such as siderophores, iron acquisition, Ton and Tol transport systems, and motility and chemotaxis ( Figure 2 A).

Detailed read mapping supporting each of the 15 SNVs shown in Figure 2 A, that were uniquely found in the Gok2 strain. These results show that each of these mutations (marked with red arrows) were supported by independent reads, and were not placed at the extreme of reads (where deamination due to ancient DNA degradation take place). These results strongly suggest that these are genuine variants that emerged within in the Gok2 strain.

Altogether, (A) and (B) show that the Gok2 and Bronze Age strains belonged to two distinct and independent Y. pestis clades, with genetic variation that is no longer represented in known modern lineages.

(B) Same analysis but showing all variants found in the Bronze Age strains that are not identical in all Y. pestis strains. Most of the variants were only found in the Bronze Age clade and not in the Gok2 strain (highlighted with a darker light blue) or any other strain. Color bar on the right identifies strains as follows: grey, 27 Y. pseudotuberculosis strains; red, Gok2; green, 7 Bronze Age strains; purple, 3 Justinian plague strains; black, 9 strains from or related to the Black Death; cyan, all strains isolated from modern individuals. Samples are placed in the same order as in (A) using also the phylogenetic tree from Figure 1 E.

(A) Single nucleotide variants (SNVs) found in the Gok2 strain compared to the Y. pseudotuberculosis IP32953 reference genome that were not identical in all Y. pestis strains (i.e., already fixed in a common ancestor). Detailed information about each variant is indicated on the right. Positions that are identical to the reference are colored in light blue, unmapped positions in white. The presence of each variant in Y. pseudotuberculosis strains is indicated in gray, for Gok2 in dark-blue, for ancient Y. pestis strains in blue and for all other Y. pestis in dark-red. The ML phylogenetic tree from Figure 1 E was used to sort genomes in the heatmap and is shown at the top of the figure. Ts, transitions; Tv, transversions; N.D, not determined. The results show that 15 of the 28 SNVs identified in Gok2 were unique to this strain.

To explore the history of early Y. pestis strains, we performed a molecular clock analysis to estimate the divergence times between lineages. Our results indicate that Gok2 diverged from all other strains 5,700 years BP (95% HPD interval: 5,250–6,364 BP), the Bronze Age clade 5,300 years BP (95% HPD interval: 4,953–5,731 BP), and the most basal of all known modern clades (most found in China) at 5,100 BP (95% HPD interval: 4,678–5,625 BP) ( Figures 3 and S4 ). Previous genomic and phylogenetic analyses of hundreds of modern and ancient Y. pestis genomes concluded that most modern clades likely originated in East Asia, and more precisely in China (). A recent work has additionally proposed that three basal lineages (from the 0.PE4 clade), two of which persist today, may have split from the main lineage at 4,000 BP, and branched at the Eurasian steppe (). Interestingly, our results indicate that before all these events, a large-scale branching and geographic radiation of Y. pestis occurred between 6,000 and 5,000 years ago, during the period known as the Neolithic decline, and shortly before the massive migrations from the Eurasian Steppe into Europe ().

Molecular clock analysis displayed in Figure 3 , but showing all strains used in the analysis and their estimated divergence times. Strains are colored according to major Y. pestis clades and branch lengths are given as years before present (BP).

BEAST2 maximum clade credibility tree from 1,201 million states and 221k trees showing median divergence dates at tree nodes. Dates are given as years before present (BP). Only lineages diverging until the Justinian plague are shown (full version shown in Figure S4 ). The results indicate that at least four deep branching and independent lineages of Y. pestis (the extinct lineages of Gok2 and Bronze Age, the modern 0.PE7 and the root to all other known Y. pestis strains) diverged between 5,700–5,000 BP, a period matching with the decline of Neolithic populations in Europe.

Altogether, our results indicate that the initial and large-scale radiation of plague across Eurasia during the Neolithic ( Figure 3 ) was not contemporary to any major known human migration ( Figure 5 A) and that the genetic background of the infected populations did not significantly change during this period.

(D) The predicted model of early dispersion of Y. pestis during Neolithic and Bronze Age was built by integrating phylogenetic information of Y. pestis strains from this period ( Figure 1 E), their divergence times ( Figure 3 ), the geographic locations, carbon dating and genotypes of the individuals, and the archaeological record. The model suggests that early Y. pestis strains likely emerged and spread from mega-settlements in Eastern Europe (built by the Trypillia Culture) into Europe and the Eurasian steppe, most likely through human interaction networks. This was facilitated by wheeled and animal-powered transports, which are schematized in the map with red lines with arrows pointing in both senses. Our model builds upon a previous model () that proposed the spread of plague to be associated with large-scale human migrations (blue line).

(B) Geographic distribution of the use of animal traction and wheeled transport across Neolithic and Bronze Age populations in Eurasia, which broadly expanded during the period of 5,500 and 5,000 BP. The expansion of these technological innovations overlaps the predicted period for the expansion of the basal Y. pestis strains.

(A) Schematic representation of the trajectories and time periods (thousand years before present, kyr) of major known human migrations in Eurasia during the Neolithic and Bronze Age. The observed geographic distribution and divergence times of Y. pestis strains from the Gok2 and Bronze Age clades cannot be explained by the timings and routes of these human movements.

We then selected only genomes from human remains older than 3,000 BP found near the same geographic regions where Y. pestis was detected (120 genomes in total) and organized them into 6 different groups ( Figure 4 B). These groups were then subdivided when clear distinct populations from different historical periods or archaeological contexts were evidenced e.g., the Gökhem individuals from the TRB culture (group 6B) were clearly different from the preceding hunter-gatherers (group 6A) and the posterior Corded Ware Cultures (group 6C). Importantly, these results showed that the genetic backgrounds of the infected individuals were similar to other individuals within their own regions, even long time before and after the infections happened, with no signs of admixture between the populations of the different infection foci. This is particularly evident for the population of the Gökhem 2 individual (group 6B), which genetically, was significantly different to the individuals from the later Corded Ware Cultures (group 6C) (two-tailed t test, p value = 9 × 10). Similarly, the individuals infected with plague in the Altai region (RISE509 and RISE511, group 1A) was found to be significantly different from every other population where the pathogen was found except compared to group 4 (two-tailed t test, p value = 0.02) ( Figure 4 C).

We used 1,058 ancient human genomes from previously published studies ( Table S2 ) () to investigate whether human migrations may have facilitated the expansion of Y. pestis in the Neolithic and Bronze age across Eurasia. We modeled the ancestry of all these individuals as a mixture of five hypothetical populations (ADMIXTURE analysis with K = 5), which consistently recreated previous descriptions about the genetic landscape and history of the Eurasian populations ( Figures 4 A and S5 ) (). In agreement with these works, the admixture analysis showed clear differences between the populations in the Eurasian Steppe that derived from admixture events between eastern hunter-gatherers (EHG) and populations from Iran/Caucasus. Populations from Eastern and Central Europe, instead, derived from Anatolian populations with some degree of admixture with western hunter-gatherers (WHG). The results also show the major genetic turnover of the European populations after the Neolithic decline in Europe (6,000–5,000 BP), produced by massive human migrations from the Steppe into Eastern and Central Europe around 4,800–4,600 BP (). Another observation was the genetic input from Europe into the Eurasian Steppe after 4,000 BP.

Genome-wide ADMIXTURE analysis, with K = 5, with the predicted ancestry proportions of the 1,058 individuals used to build Figure 4 . The ID, lifestyle (H.G: hunter-gatherer, Farm.: Farmer), period, country where remains were found and carbon dating of each individual is indicated in the y axis labels separated by a pipe character.

(C) Boxplots showing the symmetric Kullback-Leibler (KL) divergence between the ADMIXTURE components of all individuals from groups G6B (Gökhem individuals, left) or G1A (RISE509/RISE511, right) with each other (bottom boxes) or versus the individuals from other groups (all other boxes). Lower and upper hinges correspond to the first and third quartiles and median is indicated with a red perpendicular line. Upper and lower whisker extends up to 1.5 ∗ IQR, with a red dot indicating the mean of all values. A two-sided t test was used as significance test (p values are indicated; NS, not significant; asterisk, significant). These results indicate that the ancestry of the populations infected with the most basal, and thus, initially spread Y. pestis strains, were significantly different from all other infected populations.

(B) Individuals infected with Y. pestis were divided into six different geographic regions (G1–G6, written in red) and the ADMIXTURE proportions of all available genomes (n = 120) within each region is indicated with their corresponding time periods below (in kyr). Subgroups (e.g., G1A-G1B, G6A-G6B-G6C) were defined when discontinuous ancestries/time-periods of individuals were identified. Skulls indicate the excavation sites of the infected individuals, black arrows indicate their ADMIXTURE proportions, and their carbon dates are written in green when the human genomes were available or in blue when they were not. Groups G7 (Globular Amphora Culture) and G8 (Trypillia Culture) are defined to support the model presented in Figure 5

(A) Predicted ancestry of the Eurasian and Near East populations showing an overview of their genetic makeup during the Neolithic and Bronze Age also reflecting the major genetic turnovers due to massive migrations. For example, there was a drastic genetic turnover in the populations from Ukraine around 6 kya. Similarly, massive migrations of the Eurasian Steppe populations significantly changed the genetic background of Central Europe populations after the Neolithic decline, around 4.8 kya, which was later followed by a new migration from Central Europe back to the Eurasian Steppe at 4 kya.

A total of 1,058 ancient human genomes were re-analyzed using an ADMIXTURE analysis (K = 5) to evaluate if the Neolithic dispersion of Y. pestis could be explained by human migrations or admixture events between the infected populations.

A Model of Plague Dispersion Associated with the Neolithic Decline

th millennium BP, the Globular Amphora Culture (GAC) started forming at the north of the Carpathians ( Woidich, 2014 Woidich M. The Western Globular Amphora Culture. A new model for its emergence and expansion. Ebbesen, 1977 Ebbesen K. Die Jungere. Trichterbecherkultur auf den Dänischen Inseln. Johannsen and Laursen, 2010 Johannsen N.

Laursen S. Routes and wheeled transport in late 4th–early 3rd millennium funerary customs of the Jutland Peninsula: regional evidence and European context. Sjögren, 2015 Sjögren K.G. News from Frälsegården. Aspects on Neolithic burial practices. Meller and Friederich, 2017 Meller H.

Friederich S. Salzmünde – rule or exception?. By the end of the 6millennium BP, the Globular Amphora Culture (GAC) started forming at the north of the Carpathians (), which would soon expand to encompass a region from the forest steppe of the former Trypillia settlements to Poland and Germany. There are also signs of GAC influence in the TRB culture (6,000–4,800 BP) of southwestern and southern Denmark (), which was in turn related to the culture of the Gökhem 2 and Gökhem 4 individuals in Sweden (). The results from human genome analyses also showed that the Gökhem2 population (G6B) was not significantly different from the Globular Amphora (G7) and Trypillia (G8) Cultures (two-tailed t test, p value = 0.21 and 0.44) ( Figure 4 C), which is likely indicative of a genetic continuity and a connection that could potentially explain the dispersion route of the pathogen. Moreover, although the TRB settlements in Scandinavia were rather small, they also had large ceremonial gathering places (causewayed camps), which were also characteristic in other cultures in Central and Western Europe and would have been instrumental in creating social bonds and favoring frequent movements of people, sometimes over large distances (). These places could also have favored the spread of diseases between distant regions. Altogether, the archaeological record also indicates that there was contact, although indirect, between the Trypillia Culture and the TRB populations infected with Y. pestis in Sweden. Therefore, we propose that the Y. pestis found in the TRB populations of Sweden derived from a lineage that originally emerged in the Trypillia Culture, which spread and diverged with the Globular Amphora Cultures.

Zimbler et al., 2015 Zimbler D.L.

Schroeder J.A.

Eddy J.L.

Lathem W.W. Early emergence of Yersinia pestis as a severe respiratory pathogen. Rasmussen et al., 2015 Rasmussen S.

Allentoft M.E.

Nielsen K.

Orlando L.

Sikora M.

Sjögren K.G.

Pedersen A.G.

Schubert M.

Van Dam A.

Kapel C.M.

et al. Early divergent strains of Yersinia pestis in Eurasia 5,000 years ago. Ahlström, 2009 Ahlström T. Underjordiska dödsriken - humanosteologiska undersökningar av neolitiska kollektivgravar. Coast to Coast Books no 18. Sjögren, 2015 Sjögren K.G. News from Frälsegården. Aspects on Neolithic burial practices. Rasmussen et al., 2015 Rasmussen S.

Allentoft M.E.

Nielsen K.

Orlando L.

Sikora M.

Sjögren K.G.

Pedersen A.G.

Schubert M.

Van Dam A.

Kapel C.M.

et al. Early divergent strains of Yersinia pestis in Eurasia 5,000 years ago. Andrades Valtuena et al., 2017 Andrades Valtuena A.

Mittnik A.

Key F.M.

Haak W.

Allmae R.

Belinskij A.

Daubaras M.

Feldman M.

Jankauskas R.

Jankovic I.

et al. The Stone Age plague and its persistence in Eurasia. Spyrou et al., 2018 Spyrou M.A.

Tukhbatova R.I.

Wang C.-C.

Valtueña A.A.

Lankapalli A.K.

Kondrashin V.V.

Tsybin V.A.

Khokhlov A.

Kühnert D.

Herbig A.

et al. Analysis of 3800-year-old Yersinia pestis genomes suggests Bronze Age origin for bubonic plague. Although we do not have experimental evidence of the pathogenicity of early Y. pestis, its presence in ancient DNA obtained from human teeth is a strong sign of a high titer in the bloodstream at death. This most likely indicates that it was deadly. At the genomic level, these strains contain the plasminogen activator gene that is sufficient to cause pneumonic plague—the deadliest form of historic and modern plague (). Furthermore, the RISE509 and RISE511 individuals from Afanasievo Gora were found in a mass grave containing remains of multiple individuals and likely produced by an epidemic (), while the Gökhem 2 and Gökhem 4 individuals from Frälsegården were found in a Neolithic collective grave used for an unusually short period of time (). Finally, we and others have frequently identified plague in individuals across Eurasia during the Neolithic and Bronze Age (), suggesting it was fairly prevalent and virulent. We therefore hypothesize that these early forms of plague may have contributed to the decline of Neolithic Cultures in Europe through the emergence and expansion of a prehistoric plague pandemic.

In this work, we report the discovery of plague infecting Neolithic farmers in Scandinavia, which not only pre-dates all known cases of plague, but is also basal to all known modern and ancient strains of Y. pestis. We identified a remarkable overlap between the estimated radiation times of early lineages of Y. pestis, toward Europe and the Eurasian Steppe, and the collapse of Trypillia mega-settlements in the Balkans/Eastern Europe. Our results are consistent with the existence of a prehistoric plague pandemic, spreading mainly through early trade networks, rather than massive human migrations, which allowed a rapid and large-scale expansion of the pathogen, that persisted through the Bronze Age with lineages that got eventually extinct. We propose that plague may have contributed to the Neolithic decline, which paved the way for the later steppe migrations into Europe.