dataesr/scientific-paragraphs-categorization huggingface dataset:
Gene expression analysis. We used t-distributed stochastic neighbor embedding (t-SNE) to classify HCC based on their gene expression profiles. We selected the 1000 most variably expressed genes, and we used 1 minus the weighted Pearson correlation coefficient as the distance measure. Pairwise Pearson correlation was calculated using the wtd. cors function of the weights R package. We used standard deviation substracted by 0. 2 as the weight, giving more variable genes greater influence. The resulting distance matrix was used to perform the t-SNE analysis using the R package Rtsne 61 with default parameters except the following: theta = 0, is_distance = T, pca = F, max_iter = 2000. We used the Bioconductor limma package  62  to test for differential expression between CCN-HCC and other HCC of all genes expressed in at least five samples (FPKM > 0).
We applied a q-value threshold of ≤ 0. 05 to define differentially expressed genes. We used an in-house adaptation of the GSEA method  63  to identify gene sets from the MSigDB v6 database overrepresented among upregulated and downregulated genes. Viral insertion screening. AAV2 insertions had previously been screened by viral capture and whole exome sequencing in 83 tumors from the LICA-FR cohort  16. We extended this screen to AAV2 and HBV insertions in all HCC from the LICA-FR cohort using RNA-seq and WES data. In the ICGC-JP cohort, AAV2 and HBV insertions had already been screened using WGS data and were provided by Fujimoto et al. 11  In the TCGA cohort, we screened AAV2 and HBV insertions using RNA-seq data from all tumors and WES data from 37 tumors showing viral reads or overexpression of CCNA2 or CCNE1 in RNA-seq data.
For these 10 tumors, we were able to extract reads covering the full length of the inserted viral genome and to reconstruct the complete human-virus-human chimeric sequence. Consequences of cyclin A2 alterations on protein structure. All tumors from the LICA-FR series harboring AAV2 or HBV insertions in CCNA2 were analyzed by WGS or viral capture  16  to determine the precise boundaries of viral insertion breakpoints. RNA-seq reads were then aligned on the reconstructed chimeric sequence with TopHat2 58 , and we used Cufflinks v2. 2. 1  66  to identify and quantify the different transcripts. We used ElemeNT 67 to predict transcription initiation sites and Alamut Visual software (Interactive Biosoftware) to identify splicing signals on the chimeric DNA sequence. We used ATGpr  68  to identify translation initiation sites on abnormal transcripts resulting from viral insertion or gene fusions.
Signal detection was performed using the ChemiDoc XRS system and the Image Lab software (Bio-Rad). All antibodies were used at 1:1000 dilution except secondary antibody, which was used at 1:2000. Mutational and rearrangement signature analysis. We used the Palimpsest R package  69  to extract mutational and rearrangement signatures from WGS data. For point mutations, we quantified the contribution of the 10 mutational signatures referenced on the COSMIC website (https://cancer. sanger. ac. uk/cosmic/signatures) and described as operative in liver cancers (signatures 1, 4, 5, 6, 12, 16, 17, 22, 23,  24)   22  to each tumor genome. For structural rearrangements, we performed a de novo signature analysis across the 350 HCC genomes from the LICA-FR, TCGA and ICGC-JP datasets. We 6 rearrangement signatures that were very similar to the 6 signatures we previously obtained on a smaller dataset  22  , except that the two initially described deletion signatures were now merged into signature RS5, and that a new signature emerged (RS6, dominated by inversions < 10 kb).
We used Palimpsest to quantify the contribution of each signature to each tumor genome and to estimate the probability of each structural rearrangement being due to each process. Identification of rearrangement hotspots. We identified 8466 breakpoints attributed to signature RS1 (probability > 0. 5) across the 350 HCC genomes from the LICA-FR, TCGA and ICGC-JP datasets. To account for the uneven distribution of rearrangements in the genome, we then modeled the background distribution of breakpoints considering various genomic features as described by Glodzik et al. 27  , with some modifications. In short, we divided the genome into 500 kb bins, and we characterized for each bin 17 genomic features likely to influence the density of rearrangements: replication timing in HepG2 cell line (ENCODE  70  ), highly expressed (top 25%) and low-expressed (remaining 75%) genes in normal liver, average copy-number in the cohort, repetitive sequences (segmental duplications, ALU elements and other repeats), number of N bases in the reference genome, known fragile sites  71  , chromatin staining, DNAse hyper-sensitive sites and 6 histone marks (H3K4me1, H3K4me3, H3K9me3, H3K27me3, H3K36me3, H3K27ac) in adult liver (ROADMAP  72  ).
Clinical associations. We tested the association of CCN-HCC in the LICA-FR cohort with gender, age, etiology, liver fibrosis, Edmonson grade, and vascular invasion using Wilcoxon rank sum test for continuous variables, Fisher's exact test for binary variables and Chi square test for trend for categorical variables. We used log-rank test and Kaplan-Meier method to compare overall survival between CCN-HCC and others, considering only HCC with curative resection (R0) and excluding patients who died within 3 months after surgery. Computing codes. The functions used to perform the signatures analysis and associated figures are available as an open-source R package, Palimpsest, available on Github: https://github. com/FunGeST/Palimpsest. URLs. ICGC data portal, https://dcc. icgc. org/; COSMIC database, https://cancer.
42regulatory N-terminal domain. Apart from liver cancers, none of the 2606 tumor genomes from the ICGC PCAWG dataset displayed a rearrangement breakpoint in CCNA2 intron 2. Consistently, a recent RNA-seq analysis of 9,624 TCGA samples from 33 cancer types42did not reveal any CCNA2 fusion in other cancer types. CCNA2 fusions thus appear to be specific of liver cancers. Rearrangements affecting CCNE1 promoter region result in the overexpression of cyclin E1 by bringing active enhancer regions upstream the transcription start site, mirroring the effect of viral enhancers. This mechanism was more frequent than CCNE1 amplification in the liver cancer series we analyzed.
Acknowledgements We thank  Hidewaki Nakagawa  for fruitful discussions and providing TERT promoter mutation data for the ICGC-JP series. We thank  Rameen Beroukhim  and  Joachim Weischenfeldt  for helping us access ICGC structural variant tables, and  Aristotelis Tsirigos  for providing topologically associated domain boundaries. We thank  Tatiana Popova  and  Céline Vallot  for critical discussion of the results. We thank the principal investigators of the liver cancer TCGA ( Lewis Roberts ,  David Wheeler ) and  ICGC-JP (Tatsuhiro Shibata ,  Hidewaki Nakagawa ) projects, and the  ICGC consortium  as a whole for providing the high quality data sets used in this study. We thank all the clinician surgeons and pathologists who have participated to this work:  Jean Saric ,  Christophe Laurent ,  Laurence Chiche ,  Brigitte Le Bail ,  Claire Castain  ( CHU Bordeaux ),  Alexis Laurent ,  Daniel Cherqui ,  Daniel Azoulay  ( CHU Henri Mondor, Créteil ,.
Data availability The sequencing data reported in this paper have been deposited to the EGA (European Genome-phenome Archive) database (RNA-seq accession [EGAS00001002879]; WES accessions [EGAS00001000217], [EGAS00001001002] and [EGAS00001003063]; WGS accessions [EGAS00001002408], [EGAS00001000706] and [EGAS00001002888]) and the Inter-national Cancer Genome Consortium (ICGC) data portal (http://dcc. icgc. org/; release 27, April 2018).
42regulatory N-terminal domain. Apart from liver cancers, none of the 2606 tumor genomes from the ICGC PCAWG dataset displayed a rearrangement breakpoint in CCNA2 intron 2. Consistently, a recent RNA-seq analysis of 9,624 TCGA samples from 33 cancer types42did not reveal any CCNA2 fusion in other cancer types. CCNA2 fusions thus appear to be specific of liver cancers. Rearrangements affecting CCNE1 promoter region result in the overexpression of cyclin E1 by bringing active enhancer regions upstream the transcription start site, mirroring the effect of viral enhancers. This mechanism was more frequent than CCNE1 amplification in the liver cancer series we analyzed. NATURE COMMUNICATIONS | DOI: 10. 1038/s41467-018-07552-9 CLL THYCA CMD LGG AML PEDBCA KIRC DLBCL GBM BILCA SARC OCA ML BOCA PACA SKIN HNSCC PROST LUAD LIHC COLORECT BLCA LUSC CSCC UTCA ESOCA GASTR BRCA OVCA ARTICLE NATURE COMMUNICATIONS | (2018) 9:5235 | DOI: 10.
unige. ch/montoya-currat/scripts/. It is easily modifiable to fit particular case studies with specific parameter values or to investigate other theoretical questions about hybridization. Further studies may, for instance, consider hybridization success rate as a variable instead of a fixed parameter, which may depend on demographic frequencies or given environmental conditions. More complex demographic models, incorporating environmental and demographic stochasticity, can also be envisaged in order to maximize the likelihood to explain historical datasets  35, 36. New hybridization opportunities induced by humans are expected for the near future either due to the invasion of translocated species  37  , domesticated animals  16, 38  and genetically modified organisms  9  , or due to the range expansion of species induced by climate change  4, 6.
Acknowledgements This study was financed by a fellowship from the  Center for Advanced Modeling Science (CADMOS)  and by grants from the  Swiss National Science Foundation  ( n°31003A_141233  to J. I. M. B. and n°  31003A_156853  to M. C. ). C. S. Q. acknowledges support from  CONICYT-Becas Chile , from the iGE3 student salary award and from the  Swiss National Science Foundation  ( P2GEP3_168973 ). We thank  Alicia Sanchez-Mazas ,  Bastien Chopard ,  William Amos  and  Tim Coulson  for their constructive comments on earlier versions of this manuscript. All computations were performed in the High Performance Computing (HPC) cluster baobab. unige. ch.
According to the archaeological data the oldest are 5 samples from Uybat V, kurgan 1, and Uybat III, kurgan 1  (86). Radiocarbon data on two of them (RISE 675 и RISE 677: 2600-2400 BCE) support their early dating. Other samples belong to the Chernovaya period: Okunev Ulus, Verkhniy Askiz, kurgans 1 and 2, Uybat V, kurgan 4, and Syda V, kurgan 3  (86, 87). 8 radiocarbon dates of these samples are within 2300-1900 BCE. The only deviant dating of 2600-2400 BCE was obtained for samples of Syda V, kurgan 3.
We merged this novel data with high coverage genomes from previous studies  (1, 115). For exhaustive description of genotyping parameters see Section S2. 6. All saliva samples used for generating high coverage genomes were collected by a close collaborator of the Eske Willerslev research group complying with legal requirements. All donors provided informed written consent stored in Copenhagen. Permission for undertaking the study in the country of the corresponding author in Denmark was obtained according to the Danish National Committee who deemed the study non-notifiable according to the Committee Law paragraph 14. The samples were all anonymized and remain identifiable only by the first author. In addition, we genotyped 140 individuals from 5 populations in Pakistan (Gujar, Kohistani, Tarkalani, Uthmankhel, and Yusufzai), using the Infinium OmniExpressExome-8 v.
1. 3 BeadChip array platform. All samples were collected by a member of the Eske Willerslev research group for demographic analyses in the districts of Swat and Dir. All donors provided informed written consent, and permission for undertaking the study in the country of the corresponding author was obtained according to the Danish National Committee who deemed the study non-notifiable according to the Committee Law paragraph 14. We merged this novel data with genotype data from present-day Indian populations (  43 ) and with the merged dataset from (3), which is enriched in individuals with Eurasian ancestry from various time periods ranging between the Mesolithic to the present. The merged dataset consisted of 236811 SNP sites for 1805 individuals from 165 populations.
github. io/picard/), and we added read groups to reads with AddOrReplaceRG. We merged bam files belonging to the same sample, which we then processed with the Genome Analysis Toolkit (GATK) Target Creator  (119) , providing known indels from the 1000 Genomes followed by Indel Realignment. Finally, we used samtools calmd to generate the MD tag with extended BAQ calculation. Genomic coverage was calculated using qualimap with default parameters  (120). We present basic sequencing statistics and post-mortem DNA damage in Table  S1.
S2. 3 Contamination estimates We estimated contamination using two approaches: first, using contamMix  (121) , an approach that compares the mapping affinities of each mitochondrial read to the consensus sequence of the individual with the mapping affinities to worldwide dataset of putative contaminants assembled in  (122). This approach can be used successfully on all individuals with a mitogenomic coverage > 10X. Secondly, we estimated contamination using a method developed for males in  (123)  implemented in ANGSD, taking advantage of variation at the X-chromosome to assess contamination. We show estimated contamination values in Table  S1.
S2. 5 Relatedness Including relatives in population frequency-based statistics could lead to incorrect assessments. Secondly, related individuals may be informative for interpretation on social organization. For these reasons, we estimated relatedness between all pairs of individuals using a two-step approach. We first calculate all the outgroup-f3 statistics of the form f3(Individual X, Individual Y; Mbuti) in order to identify and flag pairs of individuals with inflated levels of shared ancestry (Fig. S6 ). To follow up on this method, we estimated biological relatedness between pairs of individuals using LCMLKIN (  125 ) (https://github. com/COMBINE-lab/maximum-likelihood-relatedness-estimation). An advantage of LCMLKIN is to use genotype likelihoods instead of genotypes and therefore not assuming that genotypes are ascertained without error.
This is of particular importance in ancient DNA studies, where low coverage data is abundant. First, we selected 300,000 SNPs at random from the Human Origins dataset  (42). Next, we called genotype likelihoods at these SNP positions using 'SNPbam2vcf. py' provided with LCMLKIN. Finally, we estimated biological relatedness between pairs of individuals using LCMLKIN. Individuals with high relatedness are shown on Table  S5. Having verified that a large number of Okunevo_EMBA pairs present high levels of relatedness and given that we sampled individuals from 4 distinct burial sites (Syda 5, Uybat, Okunev Olus, and Verkhni Askiz), we wanted to investigate whether these values represented mobility across different communities or instead were the result of temporal and geographic proximity within communities.
S2. 9 D-statistics We computed allele frequency-based D-statistics (with AdmixTools) to formally test hypotheses about the ancestry composition of different groups in the merged dataset. In brief, D-statistics of the form D(H1, H2; H3, H4) are expected to be consistent with 0 if H1 and H2 form a clade in the unrooted tree (((H1, H2), H3), H4). Significant deviations from this expectation may arise due to the proposed tree being wrong, gene flow between the lineages in the tree, or differential error rates between H1 and H2. In order to assess the statistical significance of the deviation, we estimated the standard error for each statistic using a weighted block jackknife approach over 5Mb blocks and computed Z-scores for each value of D. We consider D-statistics for which |Z|>3. 3 (p-value < 0.
S2. 11 qpGraph shows no evidence of Botai-Yamnaya gene flow To validate our finding of no Botai-Yamnaya admixture, we used qpGraph (Admixtools https://github. com/DReichLab/AdmixTools) to fit a simple admixture graph on Yamnaya, Botai, EHG, CHG, Xiongnu (representing East Asian ancestry), and Mbuti (outgroup), using transversion SNPs and a jackknife block size of. 05 Morgans. This graph (Fig. S28 ) had no direct Botai-Yamnaya gene flow and fit all f4 statistics (|Z| <= 1. 77), agreeing with other results that show no evidence of direct gene flow between Yamnaya and Botai.
S2. 12 Chromopainter We extracted from our call set 621,799 positions genotyped in the Human Origins dataset  (42). We merged variants in our call set with the Human Origins genotype dataset using PLINK 1. 9, and filtering for missingness per individual (--mind 0. 51) and missingness per marker (--geno 0. 05), resulting in a total of 1,250 individuals genotyped for 581,755 SNPs, including the newly sequenced ancient samples BOT2016 (Botai), Sholpan (Central Steppe EMBA), and Yamnaya Karagash, and the previously published Ust-Ishim  (128). We then used SHAPEIT v2. r790  (129)  in mode "check" to detect variant alignment errors in our data, which we excluded from the dataset, resulting in 540,070 SNPs. We subsequently phased these genotypes using SHAPEIT with default parameters, providing the 1000 Genomes Phase 3 haplotypes and recombination map as a reference (http://mathgen.
In the remainder of this supplemental section we will usually refer to these populations by shortened names, so that they fit more easily in the figures. These shortened names are "Yamnaya", "Sidelkino", "Botai", "Sholpan", "Okunevo", and "ShamEN". Sholpan is the site of the 9x-coverage CentralSteppe_EMBA sample; the other shortened names are self-explanatory. We used the method momi (28) to compute expected SFS values under the multipopulation coalescent, which were then combined into a composite likelihood, where the observed SFS was modeled to be drawn from a multinomial distribution, while the total number of heterozygotes per individual were modeled as independent Poisson variables (we used heterozygotes per individual, rather than the total number of SNPs in the dataset, because it is easier to account for the effect of missing data).
However, instead of using the usual f 4 statistic, which is based on the difference of BABA and ABBA counts, we used a modified version of it, which we denote by f 4 * , and define as f 4 * = log(BABA) -log(ABBA) = log( ). That is, instead of using the difference of BABA and ABBA counts, we use the difference of their logarithms. f 4 * is robust to certain biases that may affect f 4 = through the normalization constant N (the total number of observed SNPs). In particular, missing data or reference bias may cause a decrease in observed singletons, especially in lower-coverage individuals, leading to a decrease in the total number of SNPs. By contrast, f 4 * only depends on BABA and ABBA counts, which require 2 copies of each allele and thus are not affected by singleton counts. To compute the empirical f 4 * , we counted the number of BABA and ABBA SNPs in every subsample of 4 alleles and took the log-ratio.
We called Y-chromosomal variants in 45 ancient and 103 modern samples (Section S2. 1) using bcftools (http://www. htslib. org/doc/bcftools. html) (134)  mpileup and bcftools call emitting all sites within mappable Y-chromosomal regions  (135). Haplogroup determination was done with the script callHaplogroups. py distributed with Yhaplo  (136) , with the parameter --ancDer, which outputs the allele counts for ancestral and derived SNPs along a path of branches of the Y-chromosome tree. In total, approximately 20,000 phylogenetically informative SNPs from the ISOGG 2016 database (http://isogg. org/tree/ISOGG_YDNA_SNP_Index. html) were used for haplogroup determination. Given the low coverage of the ancient DNA samples and the effect of deamination on lineage determination, we manually inspected ancestral and derived alleles to evaluate their phylogenetic consistency, ensuring that the lineages identified were the most likely considering the data observed.
S2. 14. 1. 2 Y-chromosome phylogeny We investigated Y-chromosomes in our dataset by first selecting 103 present-day individuals from Africa, Eurasia and the Americas, including the ones newly sequenced in the present work and 6 additional high-coverage ancient samples: Yamnaya (present study), Clovis  (137) , Ust-Ishim  (128) , Saqqaq  (138) , KK1  (7) , and BR2  (139). We filtered heterozygous SNPs from this dataset to remove potential deamination and errors and selected variants with exactly 2 alleles, minimum allele count of 1, depth of coverage >=5 and genotyping rate 0. 97, and restricted variants within callable regions of the Y-chromosome. This resulted in a VCF file with 19534 SNPs, which we converted to tab format with vcf-to-tab  (9) , and then to multi-fasta with vcf_tab_to_fasta_alignment.
S2. 14. 1. 3 Adding low-coverage ancient branches to a tree estimated with high-coverage Ychromosomal data The ancient DNA (aDNA) field is abundant in low/medium-coverage data, but considering the difficulties inherent to estimating accurate phylogenies from it, datasets with large number of ancient samples are rarely represented in the form of a tree. Therefore, we aimed to incorporate low-coverage ancient samples into a pre-computed Y-chromosome tree with high-coverage modern and ancient samples. The main idea behind this approach is that haplogroup names identified in aDNA samples can be informative about their relative position on the tree. For example, a given sample carrying the M269-R1b1a2 lineage should be placed within the same clade as other R1b1a2-derived individuals. In the case where further downstream markers are not available for that particular sample, which would allow placing it at a more specific branch, the upper bound of confidence for placing the ancient sample would be at the node of all clade(s) containing R1b1a2-derived individuals.
S2. 14. 1. 4 Visualizing ancestral and derived SNPs Given the incompleteness typical of low coverage ancient DNA data, full Y-chromosome haplotype resolution was not possible for the majority of our samples. With this in mind, we generated a visual representation of allele status and missing data at important branches of the tree for ancient and modern samples in our dataset. Yhaplo's default behavior uses a decision table, which specifies the number of ancestral and derived SNPs required to continue traversing the tree and which nodes to visit. In this mode, the output only includes derived and ancestral alleles observed in the tree path travelled for lineage assignment. We altered the code of Yhaplo so that positions with missing data (no alleles observed) were also outputted in addition to derived and ancestral alleles.
In the phylogenetic tree, this sample was placed more closely to one R1b1a2a2 Avar and 1 Okunevo individual. The Upper Paleolithic MA1, whose ancestry is present in both Yamnaya and Botai, carries derived alleles at markers defining the basal R haplogroup, and, therefore, it is placed at the root of all R clades. The geographical distribution of R clades found in our dataset can be seen in Fig. S25.
One Okunevo sample (RISE683) belongs to Q1a1b1 (xQ1a1b1a), also identified in 1 Karasuk individual  (1)  and is extremely rare in present-day populations. In our modern dataset, 1 sample from Uzbekistan carrying Q1a1b1a is the closest match to Q1a1b1. We note that these lineages are distinct than the one presented by Saqqaq Q1a1a-F746, which is prevalent in Inuviats from the Canadian Western Territories  (143). The Okunevo individual RISE670 belongs to Q1a2b-L940 (xQ1a2b1,Q1a2b2), which has a mostly Central Asian distribution. In our modern dataset, 1 Dungan is the closest match. 2 Okunevo and 1 UstidaLN and UstidaBA individuals belong to Q1a2-M346. In  (147)  this lineage appeared only in 2 individuals, one from the South Asian Brahmin population and the other from European Croats. In our modern dataset, Q1a2 has been identified in a Tajik individual.