Whole-Exome Sequencing in the Isolated Populations of Cilento from South Italy

The present study describes the genetic architecture of the isolated populations of Cilento, through the analysis of exome sequence data of 245 representative individuals of these populations. By annotating the exome variants and cataloguing them according to their frequency and functional effects, we identified 347,684 variants, 67.4% of which are rare and low frequency variants, and 1% of them (corresponding to 319 variants per person) are classified as high functional impact variants; also, 39,946 (11.5% of the total) are novel variants, for which we determined a significant enrichment for deleterious effects. By comparing the allele frequencies in Cilento with those from the Tuscan population from the 1000 Genomes Project Phase 3, we highlighted an increase in allele frequency in Cilento especially for variants which map to genes involved in extracellular matrix formation and organization. Furthermore, among the variants showing increased frequency we identified several known rare disease-causing variants. By different population genetics analyses, we corroborated the status of the Cilento populations as genetic isolates. Finally, we showed that exome data of Cilento represents a useful local reference panel capable of improving the accuracy of genetic imputation, thus adding power to genetic studies of human traits in these populations.

Genomic studies are playing a crucial role in enlarging knowledge of genetic variation between and within populations 1 and in identifying the genetic bases of complex traits/diseases. In the last few years, a huge amount of genetic data have been collected by whole-exome sequencing (WES) and whole-genome sequencing (WGS) studies 2-10 , performed on different general and isolated populations worldwide, including isolates from Northern Italy 11 . These studies revealed an enrichment for rare and low frequency variants, and identified variants private to each examined population. However, none of these studies has included populations from Southern Italy, leading to an underrepresentation of their genetic diversity in the genomic data collected to date. To try to fill this gap, we have performed a whole-exome sequencing study on individuals from the populations of three small villages (Campora, Gioi, and Cardile) located in an inland area of the "National Park of Cilento and Vallo di Diano", in Southern Italy, included in the "Genetic Park of Cilento and Vallo di Diano Project". Campora and Gioi were settled at the beginning of the 11 th century by Greek and Byzantine monks, while Cardile was founded much later, in the mid-18th century, by a group of individuals coming from Gioi. Afterwards these populations suffered severe bottlenecks due to famine and the epidemic of bubonic plague in the 17 th century. Then, they were affected by different strong waves of emigration that drastically reduced the population size and experienced a high level of reproductive isolation until mid-20th century. As a consequence of the high percentage of endogamous marriages among their inhabitants, in each isolate the majority of the current population are connected in a unique huge genealogy, and descend from a few number of founding lineages. Also, the presence of inbreeding and a higher level of linkage disequilibrium characterize the Cilento isolates when compared to outbred populations 12,13 . Data of the Cilento isolates have contributed to many Genome-wide association studies (GWAS) for several complex traits of clinical interest [14][15][16][17][18][19][20] .

Results
Description of Cilento exome sequencing. The analysis of the whole-exome data on the Cilento populations revealed high-performance sequencing, with an aligned read mean depth on target regions of 75X.
We found 347,684 variants, 43.8% of which are rare (Minor Allele Frequency, MAF ≤ 1%) and 23.6% have MAF between 1% and 5%, giving a total of 67.4% having low frequency (MAF ≤ 5%). 20.7% of all variants identified in Cilento are singletons (Minor Allele Count, MAC = 1). We grouped sequence variants into four functional impact categories, in order of decreasing severity: 1) HIGH (Loss of Function), including stop-gain or -loss variants, frameshift variants, splice donor or acceptor variants and initiator codon variants; 2) MODERATE, including missense variants and in-frame insertions and deletions; 3) LOW, including synonymous variants, stop retained variants, incomplete terminal codon variants and splice-region variants; and 4) MODIFIER, including intronic and intergenic variants, 5′ and 3′ UTR variants, regulatory region and transcription factor binding site (TFBS) variants, miRNA variants and non-coding exon variants. We identified 3,562 HIGH impact variants, representing approximately 1% of all the variants. On average, each individual carries 319 (range 279-365) HIGH impact variants, of which 50 (range 31-68) are in a homozygous state. Rare variants are the most represented in the HIGH impact category (58.6%), while frequent variants (MAF > 1%) are more prevalent in both the LOW and MODIFIER categories (58.1% and 58.2% respectively) ( Fig. 1

Cilento novel variants.
Comparing the Cilento exome data with five largely used reference datasets (see Methods), we found that 39,946 variants, corresponding to 11.5% of all the identified variants, are novel, while 88.5% were found in at least one reference database and 10.9% were shared with all the five databases (Supplementary Table 2). The percentage of novel variants in each frequency class decreases as MAF increases (as presented in Fig. 2 Table 2). In accordance with the capacity of WES in identifying new rare variants 2,10 , 87.3% of Cilento novel variants have a MAF below 1% while the remaining 12.7% is distributed in the MAF classes above 1%. However, we also identified 19 novel variants (0.047%) with a MAF > 5%.

and in Supplementary
According to their impact, 70.2% of the novel variants were included in the MODIFIER class, 9.9% in the LOW, 18.1% in the MODERATE, and 1.8% (716 variants) in the HIGH class (Supplementary Table 2). We observed a statistical significant enrichment of Cilento novel variants in the HIGH and MODIFIER categories and a depletion in the MODERATE and LOW classes compared to variants shared with the other databases (Supplementary Table 3 and Fig. 3). Overall, these differences were more significant for those variants belonging to the lower MAF classes. In detail, for all variants with MAF < 2% we observed a significant difference in   . Functional enrichment of Cilento novel variants. The analysis was performed comparing novel variants with those shared with at least one reference database. The size of the circles represents the significance level of the two-sided test based on asymptotic normal distribution (Fisher exact test for the following categories: HIGH/MAF > 5%, MODERATE/MAF > 5%, and LOW/MAF > 5%). The x-axis indicates the fold enrichment, the vertical line indicates no enrichment (the proportions in the NOVEL and the SHARED set are equal). NS = not significant. Gene-enrichment analysis. To discern whether the genes containing the variants showing AF increases in each isolate are located in specific pathways, an over-representation analysis was performed using the ConsensusPathDB (CPDB) program 22 . Because of the small sample size of the Cilento isolates (Cardile in particular), we only considered variants with a fold increase in frequency >5 respect to TSI or variants monomorphic in TSI that have a MAF ≥ 0.0223 in Cilento villages (corresponding to a fold increase of 5 for a singleton in TSI). This MAF threshold ensures that variants analysed have at least MAC = 3 in Cardile and MAC = 5 in Gioi and Campora.
According to these criteria, 9,126 variants in Campora, 6,882 in Gioi and 11,512 in Cardile (only 354 in common between the three sets), were used for the analysis. These variants mapped to 8,442, 6,710, and 9,979 genes respectively (of which 2,087 are in common). Out of these genes, 4,087 in Campora, 3,292 in Gioi, and 4,766 in Cardile were found in at least one pathway.
The over-representation analysis identified 23 significantly enriched pathways in Campora, 69 in Gioi, and 73 in Cardile (Supplementary Table 5). We consider pathways to be significant if they have a q-value (multiple-testing adjusted p-value) < 0.05. Of these, 13 pathways were found enriched in all the three analyses, suggesting that, despite the very low number of increased allele frequency variants in common between the three villages, many are located in genes belonging to the same pathways. Furthermore, additional pathways were found in common between two villages: 13 were found both in Gioi and Cardile, 3 in Campora and Cardile and 1 in Campora and Gioi. (Table 1). Interestingly, the majority the pathways significantly enriched in common among the three villages (8 out of 13) are related to the extracellular matrix (ECM) formation and organization -in particular to the collagen component -and to its interaction with cell receptors. The remaining enriched pathways are related to Vesicle-mediated transport and Transmembrane transport, Protein digestion and absorption, Axon guidance, and Stimuli-sensing channels.
As sensitivity analysis we restricted the CPDB over-representation analysis to those variants that showed a statistically significant increase (p < 0.05) in AF compared to TSI. Although, in this case, about two-thirds of increased allele frequency variants were retained, the identified pathways remain significant at least in one of the three villages: among the 13 over-represented pathways, all of them showed a p-value < 0.05 in the three villages, and 10 pathways in Gioi and 11 in Cardile showed also a q-value < 0.05. Results of this analysis are reported in the Supplementary Table 6.
Rare genetic disease-causing variants. Taking into account the same variants used for the analysis described in the previous section, and following a search for them in the ClinVar archive of NCBI (https://www.  Table 1. List of the 13 over-represented pathways in common between the ConsensusPathDB analyses performed on the three villages. p-values are calculated according to the hypergeometric test based on the number of genes present in both the pathway-based set and input list of genes. q-values represent the p-values corrected for multiple testing using the false discovery rate method. The last two columns in the table represent the genes used as input for the analyses and the increased allele frequency variants located in those genes, that are in common between the three isolates. www.nature.com/scientificreports www.nature.com/scientificreports/ ncbi.nlm.nih.gov/clinvar) 23 , we found 20 variants increased in AF described as pathogenic for rare genetic diseases in at least one village: 9 variants were identified in Campora, 6 in Gioi, and 12 in Cardile. Also, 10 out of the 20 variants were monomorphic in TSI but they showed a MAF ≥ 0.0223 in at least one isolate. The list of the identified variants, together with the rare genetic pathologies for which they are responsible is reported in Table 2. For each variant, a description of its characteristics as well as those of the related disease is also reported in Supplementary Data. Although the overall reported diseases affect a large number of organs/systems, different conditions affecting the metabolism, the neurological system, and the eye emerge.
To confirm the increase in allele frequency of the identified variants, we directly genotyped two of them (rs72653706 and rs201680145). In particular, the variant rs72653706, selected as a representative of increased frequency in all three villages, was directly genotyped in all available samples in Cilento and showed the following allele frequencies: Campora AF = 0.027 (40 carriers on 750 genotyped individuals), Gioi AF = 0.008 ( www.nature.com/scientificreports www.nature.com/scientificreports/ on 720), and Cardile AF = 0.046 (50 carriers on 544). No individual was homozygous for the deleterious allele in the Cilento villages. One individual, not included in the WES study, was found to be heterozygous for the rs72653706 variant and was affected by Pseudoxanthoma Elasticum, as reported in the clinical data collected. The hypothesis that this individual is a compound heterozygous should be further verified. The variantrs201680145, responsible for CADASIL disease, and chosen due to its presence in the Cardile population, was directly genotyped on the overall Cardile sample, confirming the presence of the disease allele in this population with an AF = 0.073 (76 heterozygous and 3 homozygous individuals on 564 genotyped individuals).

Population Genetics analyses.
To characterize the genetic differentiation of the Cilento isolates, a Principal Component Analysis (PCA) was performed including data of each village and TSI population. The analysis, both using common (MAF > 5%) and rare and low-frequency (MAF ≤ 5%) variants revealed three defined clusters corresponding to the three isolates, with the first principal component discriminating Campora and Cardile from Gioi and the second principal component separating Campora and Cardile from TSI. We observed that Gioi partially overlaps with TSI, in particular when considering rare and low frequency variants (Fig. 4a,b). Effective population size (Ne) history based on LD estimations reveals that the three isolates have substantially smaller Ne compared to the Italian general population (TSI). The analysis showed that 30,000 years ago there was a comparable size among the three isolates (about 3,500 for isolates vs 4,500 for TSI). After, while the Ne of Campora and Cardile approximately remains unchanged over time, the Ne of Gioi tends to increase but always remaining smaller than TSI. At 5,000 years ago; the isolate population sizes diverge considerably with Ne~3,000 in Campora and Cardile, Ne~5,000 in Gioi vs Ne~10,000 in TSI (Supplementary Fig. 1).
The analysis also revealed that, while in TSI the majority of ROH show a length between 1 and 2 Mb (~60%), in Cilento populations these short segments represent less than 30% of the total, while a higher prevalence of long segments (>2 Mb) can be observed in these populations compared to TSI. This is especially evident for segments with a length >10 Mb, which occur with a frequency of 1.27% in TSI, while they reach a frequency of 7.51% in Campora, 8.98% in Cardile and 6.38% in Gioi (Supplementary Fig. 2). These results are in agreement with isolation and the presence of consanguineous unions in the history of these populations. Improvement of Imputation. Two imputation strategies were compared on the Cilento dataset on chromosome 10: firstly using the latest release of the 1KG_Ph3 and secondly using a combination of the 1KG_Ph3 and a local reference panel of haplotypes created by phasing the Cilento WES data (we denote these two imputation strategies as '1KG_Ph3' and '1KG_Ph3 + WES'). Over the entirety of chromosome 10, 3,921,562 and 3,923,571 variants were imputed by the strategies 1KG_Ph3 and 1KG_Ph3 + WES respectively. Indeed, including this WES local reference panel will allow the imputation of the previously discussed novel variants of Cilento. We compared imputation quality scores for variants that were neither imputed as monomorphic or present on either of the two genotyping arrays used in Cilento; our comparison involved 2,321,569 variants including 7,545 exonic variants present in the WES panel. In agreement with a recent study on simulated data in Campora 24 and other studies examining population specific imputation panels [25][26][27][28] , our analysis revealed the addition of the WES panel improved the quality of imputation, particularly for low frequency exonic variants. Indeed, the inclusion of the WES panel led to an average increase in the 'info' score by 2.1% (0.96 to 0.98) across exonic variants with a MAF > 0.05 and by 15.9% (0.77 to 0.89) over exonic variants with MAF ≤ 0.05 (Fig. 6). From the set of variants used for the comparison of imputation strategies, 1KG_Ph3 + WES imputed variants had higher 'info' scores: the number of poorly imputed variants ('info' <0.4) was 1,259,194 under 1KG_Ph3 + WES compared to 1,300,113 under 1KG_Ph3 and the number of confidently imputed variants ('info' >0.7) was 872,425 under 1KG_Ph3 + WES compared to 756,520 under 1KG_Ph3 (+15.3% of confidently imputed variants). www.nature.com/scientificreports www.nature.com/scientificreports/

Discussion
Next-generation sequencing approaches have produced a more comprehensive survey of the genetic diversity in human population improving knowledge of population history and disease mapping studies. These studies are particularly interesting in founder populations which constitute a special resource due to the effects of bottlenecks and drift on their genetic variation. In fact, several examples of complex trait associations in isolates involve variants that are rare or at low frequency in the general population yet are shown to occur with greater frequency in the studied population isolates [29][30][31][32] . This work provides a deep description of the exonic architecture of the isolated populations of Cilento, and represents the largest exome variant analysis in a population-based study in Southern Italy. The study identifies several variants with an increased allele frequency that could be prioritized in genetic association studies and for which fewer individuals would be required to achieve sufficient power to detect association. In line with a founder effect, some mutations, responsible for rare genetic conditions, are at high frequency in Cilento isolates.
In this study we have identified 347,684 variants, the vast majority of which (over two thirds) are rare or at low frequency. This is in accordance with other previous studies on whole-exome sequencing analysis 2,4 . According to the effect of variants, each person in our study carries, on average, 319 HIGH impact variations; a lower burden of deleterious variations per person was reported in a similar study on the isolated population of VIS in Croatia 2 . However, in that case a slightly different classification of loss of function variants was used including only the stop-gained (nonsense) or splice site-disrupting (splice donor or acceptor) variants.
By comparing the Cilento sequences with data from five public databases we have identified many novel, Cilento specific variants, representing more than 11% of the total set of variations found in Cilento. This was similar to the percentage observed in an exome-sequencing study on an isolated population from the Dalmatian island of Vis (9%; where five databases were explored) 2 and lower than the percentage reported in a whole-exome sequencing study of a French-Canadian founder population (20%; but here only the 1000 genome database was explored) 4 . As already found in other studies on isolated populations 2,10 , novel variants were enriched in rarity and functional relevance, which is expected as these variants have not fully undergone purifying selection due to the relatively recent foundation of these isolates.
The comparison of allele frequencies of the Cilento isolate variants with those from the Tuscan population, a general Italian population, highlighted the existence of variations that have drifted positively in frequency in Cilento. Interestingly, a pattern of common enriched pathways, mainly related to extracellular matrix formation and organization, was found in all the three villages. Genes related to collagen proteins and the extracellular matrix, and genes associated with the anatomy and physiology of cilium were reported to be potentially targets of positive selection in a study on mammal adaptation to environment 33 . Also, it has been reported in other two human populations that extracellular matrix genes show strong signs of balancing selection 34 . A work by Fumagalli et al. showed that genes coding for extracellular matrix component are implicated in an anti-viral response and can be thought of as targets of a virus-driven selective pressure 35 .
Genetic variations with increased AF in Cilento also included several rare disease causing variants responsible especially for pathologies compatible with life and reproduction for which high frequencies of carriers were observed. As these disease-variants were responsible for a number of diseases, it is unlikely that all these different conditions have been the subject of a carrier selective advantage. Therefore, although that hypothesis cannot be excluded, the results suggest that a stochastic effect of drift is most likely to represent the main force generating high carrier frequencies in Cilento. Recent studies on exome data have been focused on a survey of disease-allele www.nature.com/scientificreports www.nature.com/scientificreports/ carrier frequency at population level in both isolate communities and general populations revealing strong population specificity and variability of disease alleles 36,37 . In this context, our study provides a contribution to the frequency of rare disease causing variants in the Cilento region. Studies as such applied to larger exome/whole sequencing data from general population in Southern Italy, might contribute to create comprehensive panels of genetic screening of clinical interest.
We have added confirmation of the genetic isolation of the Cilento populations. Indeed, autozygosity, detected through ROH analysis on WES variants, revealed longer homozygous segments, especially in Campora and Cardile, compared to the general population, in accordance with similar findings obtained in two studies on Greek and other European isolated populations, in which ROH were called from SNP array data and whole-genome sequencing data, respectively 11,21 . The exome analysis confirmed previous findings from SNP chip data 38 and places Cilento among the other European and Italian isolates for which WGS data are available. 11 . Furthermore, the results obtained from the PCA and the Effective Population Size analyses, confirmed a genetic differentiation of the Cilento isolates from the general Italian population. This differentiation was greater for Campora and Cardile compared to Gioi, as previously reported 12,13 . Effective Population Size analysis shows that isolates from Southern Italy differentiated from Italian general population similarly to isolates from Northern Italy 11 .
Finally, in accordance with previous studies demonstrating that the addition of a local reference panel in the genome-wide imputation procedure substantially improves the imputation quality, especially for rare and low frequency variants 24-28 , we reported a greater accuracy in the imputation especially for rare variants in Cilento when using a local reference panel. This will open up the opportunity for future GWAS for complex traits in Cilento including well imputed low frequency variants and indeed novel Cilento variants. This may potentially provide, in turn, a more complete description of the genetic architecture for such traits.
Methods participant recruitment and sample collection. We conducted exome-wide sequencing of 247 individuals from the isolated villages (Campora (N = 93), Gioi (N = 94) and Cardile (N = 60)) in the Cilento area, South Italy. The sample for this study was based on the initial cohort of 2,137 participants, which were recruited in the Genetic Park of Cilento and Vallo di Diano Project between 2003 and 2008. Sequenced samples were selected on the basis of three criteria: i) they (and both their parents) were born in one of the three Cilento villages; ii) they had been already genotyped with Illumina 370 K or OmniExpress arrays; iii) they maximize the whole sample haplotype representation, as they were selected using the procedure proposed by Uricchio 39 . We adapted Uricchio's procedure from pedigree to exome data as follow: for the 1,617 individuals genotyped with one of the two Illumina arrays, the Identical by Descent (IBD) segments were reconstructed with GERMLINE using the 192,092 SNPs in common between the two chips. Then the segments located on exome regions were selected to calculate the exome kinship using the estimated Jacquard identity coefficients 40 . Finally, the person with the highest average kinship to the other individuals in the study was chosen and subsequently, in an iterative manner, the individuals who had the highest average kinship to the remaining unselected individuals, but who had kinship not exceeding 0.1 with any other individual already selected for sequencing, were selected.
Mean genomic kinship of the 247 individual sample was 0.007, as estimated by GenABEL library 41 .
The study design was approved by the ethics committee of Azienda Sanitaria Locale Napoli 1 and the ethics committee "Carlo Romano" University of Naples "Federico II". The study was conducted according to the criteria set by the declaration of Helsinki and each subject signed an informed consent form before participating in the study. exome capture and sequencing. Sequencing was performed at the Wellcome Trust Sanger Institute, Hinxton, Cambridge, UK. Genomic DNA (approximately 1 μg) was fragmented to an average size of 150 bp and subjected to DNA library creation using established Illumina paired-end protocols. Adapter-ligated libraries were amplified and indexed via Polymerase Chain Reaction (PCR). A portion of each library was used to create an equimolar pool comprising eight indexed libraries. Each pool was hybridised to SureSelect RNA baits (Human All Exon V5; Agilent Technologies) and sequence targets were captured and amplified in accordance with the manufacturer's recommendations. Enriched libraries were subjected to 75 base paired-end sequencing (HiSeq 2000; Illumina) following the manufacturer's instructions.
After data pre-processing, a multi-sample variant-calling was performed using Samtools mpileup (v0.1.19). GATK Unified Genotyper (v1.5-32-g2761da9) was used to recall the sites discovered by samtools and to perform the Variant Quality Score Recalibration (VQSR). Single Nucleotide Variants (SNVs) and indels were filtered by minimum VQSLOD score corresponding to the truth sensitivity threshold of 99.5% and 95% respectively. Genotypes with a depth (DP) < 4 and/or a genotype quality (GQ) < 60 for that particular call, were set as null genotypes.
Called variants were annotated with dbSNP137 rsIDs and 1000 Genomes super population allele frequencies that were extracted from the final 1000 Genomes Phase 1 integrated (v3) callset.
Performance of whole-exome sequencing. Performance of whole-exome sequencing was evaluated using GATK 42 . We collected high-coverage exome sequencing data for 247 samples with a mean of 83.3 million reads per subject. When PCR duplicates were removed, 69.8 million reads (83.9%) per individual were retained.
The mean value of aligned-read depth on the target regions (~50 Mb) was ~75x (range 57-103x) with more than 85% of target regions covered by at least 30x, and 66% covered by at least 50x. Considering also the calling regions outside the targets (for a total of ~90 Mb), the mean value of depth per individual was 50 (range 39-72), with 75% of positions covered by at least at 30x.
Per-sample and per-variant Quality control. The sample-level Quality Control (QC) criteria were: (a) genotype concordance for the overlapping set of individuals/SNPs between the Cilento genome-wide data (genotyped using Illumina 370 K or OmniExpress arrays) and whole-exome sequences. Two individuals (from Cardile population) with more than 50 discordant sites were excluded. A very high genotype concordance rate (mean 99.8%, SD 0.05%) was observed for the remaining individuals; (b) a sample call rate exclusion threshold of 90% was applied with no sample exclusion; (c) the heterozygosity exclusion threshold of ±3 standard deviations (SD) from the mean was used again with no sample exclusion; (d) outliers were searched for using multidimensional scaling analysis (MDS), plots were visually inspected and no subjects were excluded. Per-sample QC steps were performed using PLINK 1.9 44 . Samples failing sample-level QC were removed prior to performing variant-level QC.
The variant-level QC step was carried out under the following criteria: (a) only variants that were autosomal, biallelic, and polymorphic in the Cilento sample (MAC ≥ 1) were included; (b) variant call rate threshold of 95% was applied; (c) deviation from Hardy-Weinberg equilibrium (HWE) was calculated and variants falling below the threshold for HWE, exact p-value < 1 × 10 −4 , were excluded. All variants failing QC steps were removed. All subsequent analyses were performed using clean post-QCed datasets. Per-variant QC steps were performed using VCFtools (v0.1.12b).
Ti/Tv ratio. To assess the probability of false positives we evaluated aggregate transition-transversion (Ti/Tv) ratio for all variants that passed QC. Average Ti/Tv ratio for targeted regions over all subjects and all locations was 2.57, in accordance with the value expected in such studies.

Cilento variant analysis.
As the sequenced individuals were selected to be not closely related to each other, their relationship should not influence the allele frequency estimation. For this reason, the MAF of the variants identified in the Cilento sample was calculated by naive counting. According to their MAF, we classified the 347,684 variants in 6 categories: singletons (MAC = 1), doubletons (MAC = 2), ≤1%, 1-2%, 2-5%, >5%. The accuracy of singletons was verified by looking at allele transmission in a subset of 6 out of the 245 individuals who were part of parent-child couples for which additional whole-genome sequences were also available. On average, we observed a transmission rate of 50.5% for these variants, as expected. This result indicates that the identified singletons are unlikely to be false positive variants.

Functional enrichment of novel variants.
To test for functional enrichment of Cilento novel variants, we divided all variants in two groups: those found in at least one of the five reference database (SHARED) and those private to Cilento (NOVEL). Variants were then grouped by MAF and impact categories, as described above. Two sided tests based on asymptotic normal distributions were used to compare the proportion between the SHARED and NOVEL datasets for each MAF and Impact category. A Fisher exact test was used instead for the categories containing variant counts <5. Fold enrichment was calculated with the following equation: A Principal Component Analysis (PCA) was performed on Campora, Gioi, Cardile, and TSI with PLINK 1.9 software 44 , using the 106,828 common (MAF > 5%) and the 59,668 rare and low-frequency (MAF ≤ 5%) variants separately.
LD-based demographic inference for the Cilento and TSI populations was performed in the NeON R-package 45 , reconstructing the long term effective population size (Ne) from 5000 to 30000 years ago.
To detect the ROH, 60,837 among the 106,828 common variants were selected to be independent using the LD pruning option of PLINK 1.9: In windows of 50 SNPs, LD is calculated between each SNP pair, if for a pair of SNPs, the LD r 2 exceeds 0.8 then one SNP is removed. This process is then repeated by sliding the window down the chromosome in increments of five SNPs. The ROH were reconstructed in Campora, Gioi, Cardile, and TSI applying the parameters optimized for WES analysis. That is, ROH was required to comprise a minimum of 50 consecutive homozygous SNPs without any heterozygous sites (options in Plink were as follows: -homozyg-snp 50, -homozyg-window-het 0) 46 . Only ROH with a length >1 Mb were selected. For each individual included in the analysis, the total number of segments, the total length (corresponding to the sum of ROH), and the mean length (the total length divided by the total number) of ROH were calculated. T-tests were used to test the statistical significance of the differences in means between populations. Imputation quality. We sought to verify if the use of the Cilento WES data as a local reference panel (LRP) in addition to the 1KG_Ph3 1 could improve the accuracy of genome-wide imputation. To do this, we performed two imputations: the first using only the 1KG_Ph3 as reference panel, the second using a combination of both the 1KG_Ph3 and a Cilento LRP of WES data (1KG_Ph3 + WES). Variants from both Illumina 370 K and OmniExpress arrays and WES on chromosome 10 were filtered for genotyping call rate ≥95% and MAF ≥ 0.01. Mendelian errors in the genotyping data were located with PLINK 1.9 and set to missing. Two individuals were excluded from this analysis for having array genotyping rates below 70%. Haplotypes were reconstructed in two groups, one for each genotyping platform. To create the LRP for Cilento, WES data and SNPs present on both genotyping arrays were merged for the 245 exome sequenced individuals with high quality genotypes 47 . Discordant genotypes were set to missing (a mean concordance of 99.6% between WES and array genotypes was observed on chromosome 10). SHAPEIT2 was used for phasing both the 370 K and OmniExpress groups, as well as the LRP. We ran SHAPEIT2 with default settings apart from the following specifications: 15 burn-in iterations, 15 pruning iterations, 35 main iterations, and the 'duohmm' option 48 . We also supplied the 1KG_Ph3 as a reference panel for phasing. Imputation was performed on the 370 K and OmniExpress groups using IMPUTE2 49 with default settings, windows of size 5 Mb, and buffer regions of size 250 Kb. The 'merge-ref-panels' option was used to perform imputation with a combination of the 1KG_Ph3 and the phased WES LRP. For each variant imputed in both groups, the IMPUTE2 quality score ('info') was calculated using QCTOOL (http://www.well.ox.ac.uk/~gav/qctool_v2/) from the imputed dosages of all individuals except those in the LRP. The 'info' score is an "information measure" that assesses the imputation quality of a SNP. The 'info' takes values between 0 and 1, where higher values indicating that a SNP has been imputed with higher certainty. This serves as a proxy for a measurement of imputation accuracy 50 . A threshold of 0.4 for the info score is generally used to select well-imputed variants to be used in GWAS. We then compared 'info' scores from different imputation pipelines against MAFs calculated on the European populations in the 1KG_Ph3.