Genetic Structure of the World's First Farmers

new

Admin
Administrator

Posts: 81,732

Genetic Structure of the World's First Farmers Nov 30, 2020 20:42:09 GMT

Quote

Post by Admin on Nov 30, 2020 20:42:09 GMT


298 Advantages of demographic modeling
299 Our sequencing of ancient genomes at >10X, which triples the number of whole genomes
300 available for the Early Holocene period in Europe, has allowed us to perform genetic analyses
301 on an unbiased set of markers minimally impacted by selection, thus ideally suited for
302 reconstructing the Neolithic settlement history of Western Eurasia. Our results fit into the larger
303 picture of an Anatolian origin of the first farmers and a settlement of Europe from the wider
304 Aegean region14. We also confirm the deep structure between early farming communities in
305 the eastern Fertile Crescent and NW Anatolia15,18,51. However, even though early farmers of
306 Central Anatolia are rather similar to those of NW Anatolia17, we show that unlike the
307 populations in the Aegean, they are not part of the Neolithic migration chain to Europe.

308 Two unexpected findings were i) that the group formerly called WHG had already split into
309 two subgroups (west later substructured, and central HGs) approximately 23 kya, and ii) that
310 population with substantial contributions from HGs as they appear later in SE Europe.

311 Even though early farmers were recognized to be genetically intermediate
312 between other Near Eastern groups15, or considered as a mixture of other ancient21 or modern17 
313 populations, this initial admixture signal remained hidden to previous approaches as it was eroded
314 by later genetic drift (Fig. S50-51 and see below).

316 Model validation
317 These results are by no means definitive, and additional high-depth Mesolithic and
318 Epipalaeolithic samples from Greece, SE Anatolia and, importantly, from the Northern Levant
319 are needed to confirm our results. Nevertheless, we can show that the demographic model
320 reported in Fig. 3a-b can reproduce observed population affinities and patterns of admixture.
321 Indeed, genomic data simulated under our best demographic scenario leads to population
322 relationships very similar to those observed on an MDS (Fig. 3c), providing an a posteriori
323 validation of our model-based approach. In addition, we find that the ancestors of NW
324 Anatolian and European farmers simulated just after their defining admixture event (19 kya)
325 were exactly intermediate on a MDS between the two source populations (Fig. S50b), and that
326 their position on the MDS drifts away over time from this initial intermediate position towards
327 that of the early farmers. Similarly, the estimation of admixture on the same simulated data
328 shortly after the admixture event correctly recovers admixture proportions, but this signal
329 declines rapidly and disappears by the time ancient samples are sequenced (Fig. S51). Ten
330 thousand years of genetic drift have erased the initial admixture, thus explaining why the hybrid
331 nature of early farmers had been previously unnoticed. Finally, we used these simulations to
332 construct an admixture graph19 matching the scenario depicted in Fig. 3 (Fig. S53). The
333 estimated qpGraph is compatible with all f-statistics calculated on the real data (Fig. S54-55)
334 and it also recovers an

336 In sum, our population modelling allowed us to extract novel, unforeseen, but complementary
337 and far more detailed information on population affinities and past history than what one can
338 conclude from summary statistics or multivariate analysis alone. In addition, we now have a
339 time frame for the differentiation of the major groups populating SW Asia and Europe from
340 the LGM until the introduction of agriculture, highlighting the crucial role of climatic changes
341 in promoting population fragmentation and secondary contacts52.


412 Population genetics
413 Genetic relationships among individuals were estimated from pairwise average nucleotide
414 divergence !XY67 computed on their neutral genomic portion and represented with a classical
415 multidimensional scaling (MDS) approach implemented in R (cmdscale function).
416 Genomic heterozygosity was computed as the amount of heterozygous sites found in the
417 neutrally evolving portion of each genome divided by the expected number of neutral sites
418 (Suppl. Information - Section 4).

419 Admixture clustering analyses were realized for two subsets of ancient individuals (Suppl.
420 Information - Sections 4&7), focusing on sites with no missing data, with R package LEA68,
421 parameters K, alpha = 100, number of repetitions = 5, and run goodness-of-fit and admixture
422 coefficients being calculated in an unsupervised manner with function snmf.
423 D-statistics were computed on pseudo-haploid data from reference 1240K dataset v42.4 and
424 on majority calls generated with ATLAS at the ~1.2 mil. SNPs from the reference. D-statistics
425 was calculated in a form of D(Individual1, Individual2; Individual3, Outgroup) using
426 ADMIXTOOLS19 and Mbuti as the outgroup. We used qpGraph to check the fit between f
427 statistics.

428 Runs of Homozygosity (ROHs) were identified in genomes of modern and ancient Western
429 Eurasian individuals imputed using IBDSeq v. r120669 with default parameters but errormax =
430 0.005 and ibdlod = 2, and after artificially long tracts spanning assembly gaps or centromeres
431 were split into shorter tracks excluding the gap stretch.


437 Phenotype predictions
438 Pigmentation phenotypes of hair, skin and eyes were predicted for each of the newly sequenced
439 samples with the HirisPlexS webtool72,73 on genotypes or BAM files directly when genotypes
440 were missing (Suppl. Information - Section 5).

441 Genotypes for SNPs associated with additional phenotypes of interest were inspected manually
442 for each sample: rs4988235 variant in MCM6 gene and associated with lactase-persistence in
443 Eurasia; rs3827760 in EDAR gene; rs17822931 in ABCC11; seven SNPs located in the
444 FADS1/2 gene complex.

445 For predicted standing height, polygenic scores (PS) were computed based on a set of 670
446 SNPs74.

447 Demographic Analyses
448 MSMC2 analyses75 were performed on the phased-imputed dataset following the author’s
449 recommendation, including masking for chromosome mappability on the hs37d5 human
450 reference genome and for sample-specific sequencing quality (Supp. Information - Section 6).

451 fastsimcoal2 analyses25 were carried out on seven different data sets of newly sequenced
452 individuals, on the neutral SFS and with neutral mutation rate defined for each dataset (see Supp.
453 Information - Section 6). We used 50 independent runs with 100 expectation conditional
454 maximization (ECM; 150 in one model) cycles per run and 500,000 coalescent simulations per
455 estimation of the expected SFS. Confidence intervals for the maximum-likelihood parameters
456 point estimates were computed with a parametric bootstrap approach.

Last Edit: Dec 1, 2020 6:50:20 GMT by Admin

Admin
Administrator

Posts: 81,732

Genetic Structure of the World's First Farmers Dec 1, 2020 6:43:20 GMT

Quote

Post by Admin on Dec 1, 2020 6:43:20 GMT

Section 4: Genetic structure analysis
Multi-Dimensional Scaling (MDS)
For each pair of individuals X and Y, we identified the sites showing no missing data for the two
individuals, and we computed their average nucleotide divergence %XY162 over these sites. We then
represented the relationships between the ancient genomes (Fig. S19a) or modern and ancient
genomes from Europe and SW Asia (Fig. S19b), using a classical multidimensional scaling (MDS)
approach implemented in R (cmdscale function). We also performed the same analyses on a
restricted set of sites defined as neutral according to Pouyet et al.’s163 definition (Fig. 2a, Fig. S50a).

The analysis shown in Fig. S19a reveals a slightly different picture than Fig. 3 performed only on
neutrally evolving sites, as the Upper Palaeolithic sample of Bichon (shown as an inverted triangle)
is here well within the cluster of the other four hunter-gatherers (HGs; triangles with the tip up).
Furthermore, Fig. S19b suggests strongest affinities of European Neolithic farmers with modern
individuals from Southern Europe (Crete, Greece, Italy, Albania, Spain) rather than with Sardinians,
at odds with Fig. 2a and previous analyses based on the projections of ancient individuals on
Principal Components computed from modern individuals only164. Note also that NW Anatolian
Early Neolithic individuals (AKT16 and Bar8) slightly diverge from the Neolithic cluster and are
closer to modern individuals from Turkey and Bulgaria, thus showing some surprising genetic
continuity across millennia within this region. CarsPas1 from England is also an outlier found in
between Southern Europeans and a Bulgarian individual.

Genetic diversity
We estimated the level of genetic diversity of each individual as the proportion of heterozygous
sites found in the neutrally evolving portion of its genome. Thus, for every genome, we divided the
amount of observed neutral heterozygous sites by the expected number of neutral sites. This number
is obtained by considering all those sites that i) have the same reference allele for both the
chimpanzee and gorilla reference genomes (for polarization), ii) are not CpG sites and out of CpG
islands, iii) are in regions with recombination rate > 1.5 cM/Mb based on 163. Finally, in order to
consider only sites unaffected by BGC, we divided the number of these previously defined sites by
3, as only one third of all mutations should be BGC free (i.e. A↔T and G↔C polymorphisms). We
represented this genomic neutral heterozygosity as a function of the geographical air distance to
Boncuklu, Central Anatolia (Fig. 2b).

Last Edit: Dec 1, 2020 6:43:50 GMT by Admin

Admin
Administrator

Posts: 81,732

Genetic Structure of the World's First Farmers Dec 2, 2020 3:44:30 GMT

Quote

Post by Admin on Dec 2, 2020 3:44:30 GMT

Runs of Homozygosity (ROHs)
We inferred the segments of homozygosity-by-descent or Runs of Homozygosity (ROHs) from the
imputed genomes (step described in Section 3) of the 89 European and SW Asian modern and
ancient individuals by using IBDSeq v. r1206 165 with default parameters but errormax = 0.005 and
ibdlod = 2. We further processed the artificially long tracts spanning assembly gaps or centromeres
(their genomic locations were obtained from UCSC genome browser:
hgdownload.cse.ucsc.edu/goldenPath/hg19/database/gap.txt.gz) and split them into shorter
tracks excluding the gap stretch, inspired by 166. As advised when genomes do not come from an
homogeneous population165, we focused on short (2-10 Mb) and long (>10 Mb) ROHs similarly to
167. We show their number and total length in ancient genomes (Fig. 2c), as well as in moderns for
comparison (Fig. S20).

In Fig. S20a, we expect the number and total genome size made up of the short ROHs to be
negatively correlated with population size, so that a larger number of short ROHs indicate smaller
population sizes. Consequently, the HG individuals seem to be drawn from small populations, as
well as Bon002, WC1 and LEPE52 among Early Neolithic individuals. In Fig. S20b, long ROHs
are indicative of recent inbreeding between close relatives, and therefore WC1 and LEPE52 are
inferred to be the most inbred ancient individuals, followed by Stuttgart, Bichon and Loschbour.
Note however that these ancient individuals appear less inbred than some modern individuals
(Jordanian-2 or BedouinB-1 for instance).

Last Edit: Dec 2, 2020 3:45:39 GMT by Admin

Admin
Administrator

Posts: 81,732

Genetic Structure of the World's First Farmers Dec 2, 2020 20:51:57 GMT

Quote

Post by Admin on Dec 2, 2020 20:51:57 GMT

Admixture analysis
We calculated levels of admixture for 20 ancient individuals: the Late Palaeolithic individual
(Bichon), 3 Mesolithic individuals (VLASA7, VLASA32 and Loschbour but not SF12 from
Sweden, as it is a geographical outlier), 16 Early Neolithic (the lowest quality Neolithic genome
Bon002 was not considered for this analysis in order to maximize the number of genomic sites to
be used, as well as Lepenski Vir genomes as it is unclear if their cultural attribution corresponds to
their genetic relatedness, Table 1, Fig. 2a). We restricted our analyses to sites without missing data
among all 20 genomes, i.e. in total 742,403. Admixture coefficients of each genome were estimated
using the R package LEA168 with parameters K=2 or 3, alpha = 100, and number of repetitions = 5.
The function snmf calculates the fit (entropy) of each run and the admixture coefficients in an
unsupervised manner.

This admixture analysis reveals for K=3, a cluster corresponding to hunter-gatherers (HGs); a
second with the Neolithic individuals; and a last cluster including the Iranian Neolithic individual.
A few Neolithic individuals (CarsPas1, NE1, AKT16 and Bar8) show some HG background
admixture. NW Anatolia in particular shows some heterogeneity: AKT16 and Bar8 are found
admixed between the Neolithic and the HG component while Bar25 is found with only the Neolithic
component, and to lower extents with Iranian Neolithic. Those results are in agreement with the
MDS analysis (Fig. 2a) and the demographic inferences (Fig. 3a, Fig. S37).

Sex-specific markers and haplogroups
We used phy-mer169 to determine mitochondrial haplogroups from BAM files for the 15 newly
sequenced genomes, with the minimal number of occurrences of a K-mer set to 10 (Table 1; Supp.
Table 3). Y-chromosomal haplogroups were determined from BAM files using Yleaf170, using the
recommended minimal base-quality of 20 (-q 20) and base-majority to determine an allele of 90%
(-b 90).

The two newly sequenced HGs from the Danube Gorges have mitochondrial haplogroup U5
(VLASA32: U5a2a, VLASA7: U5a2a), the most common haplogroup in individuals from the
Danube Gorges published to date127,171. Y-chromosomal haplogroups of these two individuals
(VLASA32: R1b1, VLASA7: I2) were also consistent, as all Mesolithic individuals from the
Danube Gorges published to date127,171 exclusively belong to haplogroups I or R (R1b*).
In contrast, Y-chromosomal haplogroups of the non-Mesolithic males in our dataset either belong
to haplogroup G (Asp6: G2a2b2a3, Bar25: G2a2b2a1, Ess7: G2a2b2a1a1, LEPE52: G2a2b2a1a1c,
VC3-2: G2a2a1a3~) or C (LEPE48: C1a2b, Dil16: C1a2b). While haplotypes from haplogroup G
(G2a*) were common among Early Neolithic farmers127, haplogroup C was less frequent in
comparison (G: 56.52%, C: 17.39%, estimated based on version v42.4 of
reichdata.hms.harvard.edu/pub/datasets/amh_repo/curated_releases/V42/V42.4/SHARE/pu
blic.dir/v42.4.1240K.anno).

Mitochondrial haplogroup frequencies observed in the Neolithic individuals were consistent with
those previously reported from a Neolithic context. We found five individuals with haplogroup K
(Nea2: K1a, Nea3: K1a2c, LEPE48: K1a1, Herx: K1a4a1i, AKT16: K1a3), as well as one individual
each with haplotypes N1a1a1 (Bar25), J1c6 (Dil16), W1-119 (Klein7), T2e2 (STAR1) and HV16311 (VC3-2).

U5 haplotypes, common among HGs127 and increasing in frequencies during the
middle and later Neolithic periods, were found in Early Neolithic individuals from Southern
Germany and Lower Austria (Ess7: U5b2c1, Asp6: U5a1c1). Additionally, haplogroup H3 was
inferred for an individual from the Early-Middle Neolithic in the Danube Gorges (LEPE52). H3*
haplotypes are rare or even absent in Early Neolithic individuals, but are found more frequently in
the Middle Neolithic in Germany172,173 as well as the Iberian Peninsula78,174,175. While it was
previously suggested that H3* haplogroups were associated with a glacial Iberian refugium and
have spread throughout Europe from there176, our data indicates that H3* haplotypes were also
already present in Neolithic Serbia.

Last Edit: Dec 2, 2020 20:53:02 GMT by Admin

Admin
Administrator

Posts: 81,732

Genetic Structure of the World's First Farmers Dec 3, 2020 5:18:58 GMT

Quote

Post by Admin on Dec 3, 2020 5:18:58 GMT

Neanderthal introgression
As a result of interbreeding between Neanderthals and modern humans some time after the out-ofAfrica migration, all non-African individuals derive a small portion of their genome from
Neanderthals.

We quantified proportions of Vindija Neanderthal177 introgression by computing f4 ratio statistics
of the form /((0123!,56!78; :,;!<=3)
/((0123!,56!78; >!<?!@3,;!<=3) as previously suggested178 with qpF4ratio from the
ADMIXTOOLS package179. To do so, we extracted the diploid HOIll sites180 covered by the 1240k
capture array and merged them with the reference dataset (version 37.2) available from David
Reich’s lab webpage (https://reich.hms.harvard.edu/downloadable-genotypes-present-day-andancient-dna-data-compiled-published-papers).

For the Neanderthal admixture proportions inferred for the ancient individuals used in fastsimcoal2
demographic modelling (Fig. S22), besides a notable outlier (Bon002), proportions roughly lie in
an interval between between 2% and 3%, without any clear pattern correlated to the population
history inferred in this study. However, comparisons of Neanderthal ancestry levels between ancient
hunter-gatherers (HGs), Neolithic (NEO) and modern European and SW Asian (MOD) individuals
(Fig. S23) reveal a significantly lower proportion of Neanderthal ancestry in modern individuals
(Kolmogorov-Smirnov tests: HGs vs. MOD, p-value=0.0184, NEO vs. MOD, p-value=0.0002), yet
no significant difference between HG and Neolithic groups (p-value=0.5114).

Our observations support previous results by Petr et al.178 (Fig. 2 and Fig. S7). The high ancestry
proportions inferred in Bichon and especially in Bon002 are likely artifacts from lower data quality