dataesr/scientific-paragraphs-categorization huggingface dataset:
2009). Only data sets containing a minimum of four homologous sequences were retained for phylogenetic reconstruction. Trees were inspected manually to correct or refine the core/ shell classification and to determine the origin of the transfers. The phylogenetic origin of a cluster or a single gene was determined from the closest neighbors in the tree at two levels of precision depending on the quality of the phylogenetic signal retained: Domains (Eucarya, Archaea, and Bacteria) or phyla/ classes (Euryarchaeota, Crenarchaeota, and Thaumarchaeota for the Archaea; Alpha-, Beta-, Gamma-, Delta-Proteobacteria, Firmicutes, Chlamydia, CFB, Actinobacteria, Acidobacteria, and Cyanobacteria for the Bacteria). We inferred bacteria-to-archaea HGT events for a given gene 1) when phylogenetic trees reproduced with reasonable support the overall monophyly of recognized bacterial phyla (though local zones of low resolution might sometimes occur, as well as limited HGT among bacteria) and Thaumarchaeota or GII/III-Euryarchaeota genes formed a monophyletic group branching within a bacterial phylum or group of monophyletic phyla and/ or 2) when Thaumarchaeota or GII/III-Euryarchaeota appeared alone to the exclusion of other archaea forming a monophyletic group among many bacterial phyla including many sequences.
Poorly resolved trees (when the phylogenetic signal of a given gene was too low and bacterial and archaeal sequences were intermixed), trees showing a very high level of HGT among bacteria or trees with few archaeal and bacterial members were (conservatively) excluded from our analysis. The manual analysis of our trees led to a more reliable list of clusters, which were classified according to their phylogenetic origin and distribution (table 1 and fig. 3 ). A proportion of genes had no similarity in the database (orphans) (table 1). To verify the prediction that HT-genes present in marine Thaumarchaeota were also shared by soil Thaumarchaeota, we included the genome of Nitrososphaera gargensis  (Spang et al. 2012 ) in our database, looked for homologs to the HT-genes that we had identified in our thaumarchaeal FIG.
Synteny Analysis and Horizontally Transferred Genes Synteny blocks in archaeal fosmids were defined as arrays of one of more contiguous genes of same origin class (archaeal core, lineage-specific core, early HT-genes, late HT-genes, and others, including the remaining predicted genes without homologs in the database). Each gene or synteny block is flanked by blocks of one or two different origins. Because we consider five possible origin gene classes for Thaumarchaeota or GII/III-Euryarchaeota, there are 15 (5 + 4 + 3 + 2 + 1) different possibilities for any synteny block (or gene) to be bounded. Bounding-couple occurrence was compiled for each synteny block in Thaumarchaeota and GII/III-Euryarchaeota fosmids separately and the corresponding data matrix subjected to hierarchical clustering analysis  (Eisen 1998 ) using the MeV package  (Saeed et al.
Codon Usage and Codon Adaptation Index Analysis of HT-Genes The codon adaptation index (CAI) for a total of 26,678 genes acquired through HGT by Thaumarchaeota and GII/III-Euryarchaeota was calculated as follows: 1,244 ribosomal protein genes (459 from GII-Euryarchaeota and 785 from Thaumarchaeota fosmids) were first selected as a reference pool of highly expressed genes, either together or in groups of similar origin, and their codon usage table calculated using the cusp program from the EMBOSS suite, version 6. 5. 7  (Rice et al. 2000). The CAI for all genes was then calculated with the three sets of ribosomal genes codon usage tables (Thaumarchaeota, Euryarchaeota, or Thaum + Euryarchaeota) serving as reference with the CAI program (EMBOSS). Codon usage values were then submitted to PCA analysis  (Raychaudhuri et al.
2012). Urea degradation seems to be a metabolic feature of deep-sea Thaumarchaeota thriving in highly oligotrophic conditions, irrespective of the geographic region or local temperature, because deep-Mediterranean waters are relatively warm (14 C on average)  (Martin-Cuadrado et al. 2007 ). Also, all the genes that have been proposed to take part in the 3-hydroxypropionate/4-hydroxybutyrate cycle for autotrophic carbon fixation in N. maritimus  (Walker et al. 2010)  were present in the thaumarchaeal pangenome, reinforcing the idea that deep-sea planktonic Thaumarchaeota have the potential for chemolithoautotrophic growth. In contrast, despite a minimum of nine complete genomes were represented in the GII/III-Euryarchaeota data set, genes encoding proteorhodopsin homologs were not detected.
To define Thaumarchaeota-specific core genes, we considered OGs shared by our deep-Mediterranean fosmids and the genomes of their closest phylogenetic relatives from aquatic environments, namely, N. maritimus SCM1 (NC_010085), C. symbiosum A (NC_014820), and N. limnia SFB1  (Blainey et al. 2011 ) (fig. 1 ) to the exclusion of other archaea, which resulted in a total of 706 Thaumarchaeotaspecific core genes. Likewise, for GII/III-Euryarchaeota, we used the composite genome built from surface seawater metagenome (CM001443. 1) (Iverson et al. 2012) , resulting in a remarkably similar number, 704, of GII/III-Euryarchaeota-specific core genes (table 1 and fig. 3 ). The remaining OGs present in only a subset of Thaumarchaeotal or GII/III-Euryarchaeota fosmids, having only one to three hits in the database or lacking homologs in archaea but not in other life domains were classified as shell genes.
Functional Classes of Transferred Genes We looked for potential functional differences between late and early HT-genes to Thaumarchaeota and GII/III-Euryarchaeota pangenomes, and between these and the corresponding archaeal core and lineage-specific core genes. Shell non-HGT genes were not included in this analysis, because they correspond to predicted genes with one to a few homologs only in fosmids and lacking clear homologs in the database (hence, nonannotated). Overall differences were already seen at a very general level of functional classification in COG classes and KEGG superclasses between gene origin classes. However, there were remarkable similarities in the functional patterns observed for the different gene origin classes between Thaumarchaeota and GII/III-Euryarchaeota (fig. 7 ), suggesting similar underlying processes and/or mechanisms of adaptation by gene acquisition.
Acknowledgments We are extremely grateful to the patients who contributed to this study and we thank all the participating investigators:  Abitbol Vered ,  Service de Gastroentérologie, Hôpital Cochin, Paris ;  Allez Matthieu ,  Service d'Hépato-Gastroentérologie, Hôpital Saint-Louis, Paris ;  Bonaz Bruno  and  Mathieu Nicolas ,  Service d'Hépato-Gastroentérologie ,  CHU Hôpital Michallon, Grenoble; Bonnaud Guillaume, Service de Gastroentérologie, Clinique des Cèdres, Cornebarrieu; Bouhnik Yoram, Service de Gastroentérologie et Assistance Nutritive, Hôpital Beaujon, Clichy ;  Bourreille Arnaud  and  Flamant Mathurin ,  Service d'Hépato-Gastroentérologie -IMAD ,  CHU Hôtel Dieu, Nantes; Cadiot Guillaume and Brixi Hedia, Service d'Hépato-Gastroentérologie, CHU Hôpital Robert Debré, Reims; Carbonnel Franck, Service d'Hépato-Gastroentérologie-Nutrition, Centre Hospitalier de Bicêtre, Le Kremlin Bicêtre ;  Colombel Jean-Frédéric  and  Nachury Maria ,  Hôpital Claude Huriez, Lille ;  Cuillerier Emmanuel ,  Service d'Hépato-Gastroentérologie ,  Centre Hospitalier Victor Jousselin ,  Dreux ;  Dupas Jean-Louis, Service d'Hépato-Gastroentérologie, CHU Amiens ;  Faure Patrick ,  Clinique Saint-Jean du Languedoc , Toulouse;  Flourié Bernard  and  Nancey Stéphane ,  Service d'Hépato-Gastroentérologie ,  CHU Hôpital Lyon Sud, Lyon; Grimaud Jean-Charles, Service d'Hépato-Gastroentérologie, CHU Hôpital Nord, Marseille; Hagège Hervé, Centre Hospitalier Intercommunal, Créteil; Hebuterne Xavier, Service de Gastroentérologie et Nutrition, CHU Hôpital de l'Archet 2, Nice; Laharie David, Service d'Hépato-Gastroentérologie et Nutrition, CHU Hôpital Haut-Lévêque, Pessac ;  Lerebours Eric  and  Savoye Guillaume ,  Service d'Hépato-Gastroentérologie et Nutrition, CHU Hôpital Charles-Nicolle, Rouen; Levy Patrick, Strasbourg; Moreau Jacques, Service de Gastroentérologie et Nutrition, CHU de Rangueil, Toulouse; Nahon Stéphane, Centre Hospitalier Intercommunal, Le Raincy Montfermeil; Peyrin-Biroulet Laurent, Service d'Hépato-Gastroentérologie, CHU Hôpital du Brabois, Nancy ;  Pingannaud Marie-Pierre ,  Marseille ;  Reimund Jean-Marie ,  Service d'Hépato-Gastroentérologie et Nutrition, CHU Caen ;  Roblin Xavier ,  Service de Gastroentérologie ,  CHU Saint Etienne; Tucat Gilbert, Paris ;  Vandromme Luc , Reims; and  Veyrac Michel ,  Service d'Hépato-Gastroentérologie ,  CHU Hôpitaux Saint-Eloi ,  Montpellier.
Transcriptome studies Microarray analysis was carried out at the Unité de Recherche en Génomique Végétale (Evry, France), using the CATMAv6. 2 array based on Roche-NimbleGen technology. CATMAv6. 2 microarray slides contain 12 chambers, each containing 219,684 primers representing all the Arabidopsis thaliana genes: 37,309 probes corresponding to TAIRv8 annotation (including 476 probes of mitochondrial and chloroplast genes) and 1,796 probes corresponding to EUGENE software predictions. The slides also include 5,328 probes corresponding to repeat elements, 1,322 probes for miRNA/MIR, 329 probes for other RNAs (rRNA,tRNA, snRNA, soRNA) and several controls. In each chamber, probes are present in triplicates and in both strands. Three independent biological replicates of the microarray analysis were produced.
For each biological repetition and each point, 14-day-old seedlings grown in long day conditions were collected and RNA samples were obtained by pooling more than 50 plants. Total RNA was extracted using Qiagen RNAeasy according to the supplier's instructions. For each comparison, one technical replicate with fluorochrome reversal was performed for each biological replicate (that is, six hybridisations per comparison). The labelling of cRNAs with Cy3-dUTP or Cy5-dUTP (Perkin-Elmer-NEN Life Science Products) and the hybridisation to the slides were performed as previously described  [83]. Two micron scanning was performed with InnoScan900 scanner (Innopsys R , Carbonne, FRANCE) and raw data were extracted using Mapix R software (Innopsys R , Carbonne, FRANCE).
We considered as being differentially expressed the probes with an adjusted P value ≤0. 05. Analysis was done with the R software. The function SqueezeVar of the library limma has been used to smooth the specific variances by computing empirical Bayes posterior means. The library kerfdr has been used to calculate the adjusted P values. The overlap between different sets of genes was generated by the Venn diagram generator Venny  [87]. The analysis to find overrepresented categories in the gene sets was obtained with AmiGO  [88] , which is based on a hypergeometric test. Co-expression analysis was performed with ATTEDII version 6. 1 using the Network Drawer tool and 'add a few genes' settings for co-expression and Protein Protein Interaction options  [89, 90]. The thickness is representative of the rank of correlation between two genes of interest via the calculation of a geometric averaged rank (MR).
Clustering of microarray data The dataset for the co-expression analysis was built from the results of the differential analyses. Probes with at least one Bonferroni pvalue lower than 0. 05 were considered. It leads to a dataset of 4,378 probes described by seven expression differences, each one being the average of the three biological replicates. The clustering was performed with a multidimensional Gaussian mixture with unequal proportions and a component number varying from 2 to 40. Covariance matrices are constrained so that their volumes differ and their orientation and shape are equal. Estimations were done with the MIXMOD software  [92]  and a mixture of 29 components was selected according to the BIC criterion. Probes were assigned in the cluster for which the conditional probability is the highest and interpretation was done only for probes for which this probability is greater than 0.
Detection of the cis-elements We analysed the presence of conserved motifs in the 5′ region of genes, also known as cis-elements. The Arabidopsis promoter dataset was downloaded from FLAGdb++ based on TAIRv8  [94]. The dataset includes 27,025 promoters containing 1,000 base pairs upstream known transcription starting sites (TSSs) or upstream the ATG start otherwise. A list of 140 motifs known to be involved in stress responses was extracted from the databases PLACE  [95]  and AGRIS  [96]. For each cluster, the presence of these motifs was identified by the Preferentially Located Motifs (PLMs) method  [97]. This method determines the preferential location of each motif relative to the TSS and a functional window derived from the peak boundaries of the region in which the transcription factor binding site is over-represented.
Gene interaction network construction A total of 12,741 protein-protein interactions (PPI) data were extracted from: (1) Arabidopsis Interactome Consortium  [98] , where a matrix of 9 k × 9 k full length protein encoding ORFs were tested by yeast 2-hybrid assay and a total of 6,475 positive interactions were detected; (2) public databases: BioGRId, IntAct, TAIR and BIND (6365 experimental PPI data). Concerning the TF-target data, 769 confirmed interactions were downloaded from AtRegNet database  [96]. These interaction information on protein-protein interactions (PPI) and transcription factor-target interactions were combined to the coexpression clusters using a home-made Perl program leading to a gene interaction network of 839 genes linked with 983 edges. Network visualisation and analysis were done with Cytoscape  [99].
Evaluation of transcriptional diversity metrics in simulated datasets We used two different strategies to assess transcriptional diversity among a group of profiles. In the first strategy, we considered all possible pairwise distances among the profiles the group using Pearson or cosine correlation for pairwise similarity of transcriptional profiles. In the second strategy, we used a dispersion distance metric to summarize the diversity as distance from the population centroid  [19]  (details in Methods). We considered either the mean or the median of the distribution of pairwise or dispersion distances as a summary of the central tendency in the group. To evaluate the performance of the different metrics, we generated artificial gene expression datasets representing scenarios that differed in the number of latent subgroups present and in the within (σ g ) and between (σ p ) tumor sample variance using the R package Umpire  [20].
Transcriptional and genomic diversity within TCGA breast cancer subtypes We assessed the transcriptional heterogeneity within the four breast cancer subtypes as assigned by the PAM50 classifier using gene expression data from The Cancer Genome Atlas (TCGA) breast cancer project (N = 547; Table  1 )  [4]. Gene expression data were used as provided and probe sets with low expression or low variation across the dataset were filtered as described in the Methods. Mean dispersion distance was computed in groups of 100 cases from each subtype, repeatedly sampled with replacement over 500 bootstrap iterations. Overall, tumors of the luminal A subtype were significantly less heterogeneous than other subtypes (permutation test P < 10 -5 ), while basallike tumors had the highest transcriptionally diversity Figure  1  Assessment of different transcriptional diversity metrics in simulated datasets.
A) Simulated gene expression profiles generated using a hierarchical model to independently control within sample (σ g ) and between samples (σ p ) transcriptional variation and the number of latent subgroups. Each profile consists of 50 genes (rows) and 40 samples (columns). Profiles for 1, 2, 4 and 40 latent subgroups are shown for low (σ p /σ g = 0. 5/1. 5) and high (σ p /σ g = 0. 5/0. 5) relative between-to-within sample variation. B) Transcriptional diversity within the simulated profiles assessed using the mean pairwise Pearson distance, the mean pairwise cosine distance, or the mean dispersion distance. The boxplots represent the distributions of these metrics obtained from 500 independent simulations of each dataset (blue: low σ p /σ g , green high σ p /σ g ). C) Same metrics as in B assessed in a two latent subgroup dataset of 40 samples with increasing proportion of the smaller subgroup.
Validation of transcriptional diversity patterns within breast cancer subtypes using microarray datasets We evaluated the generalizability of the results obtained on the TCGA data using publicly available gene expression datasets (Table  2 ) to assess the transcriptional diversity of basal-like (n = 403), HER2-enriched (n = 171), luminal A (n = 515) and luminal B (n = 210) breast cancers. We calculated the mean dispersion distance of transcriptional profiles within each subtype over 500 bootstrap iterations in groups of 100 cases per subtype. The pattern of transcriptional diversity was very similar to that observed in the TCGA dataset (Figure  2E ), with basal-like cancers being the most diverse and luminal A cancers being the most uniform (permutation test P < 10 -5 ), except that HER2enriched cancers in this cohort appeared less diverse.
Basal-like breast cancers that respond to chemotherapy are transcriptionally homogeneous Among the 403 basal-like cases in the Affymetrix dataset, 255 cases received preoperative combination chemotherapy (paclitaxel, 5-fluoracil, cyclophosphamide and doxorubicin), which allowed us to compare the transcriptional diversity of exceptionally chemotherapy sensitive cancers that achieved pathological complete response (pCR, n = 96) to those who had partial or no response as evidenced by viable amounts or residual cancer after preoperative chemotherapy (RD, n = 159)(Table  2 ). We compared the transcriptional diversity of cancers from the two response groups by calculating the mean dispersion distance in sets of 100 cases per group over 500 bootstrap iterations. Chemotherapy resistant basal-like cancers showed significantly greater transcriptional diversity compared to chemotherapy sensitive cancers (permutation test P < 10 -5 ; Figure  2 F ).
Breast cancer datasets and processing The Cancer Genome Atlas (TCGA) Dataset Datasets were downloaded from the TCGA breast carcinoma web site (https://tcga-data. nci. nih. gov/docs/publications/brca_2012/) [4]  (Table  1 ). Gene expression data by Agilent 244 K arrays were available for 547 cases. DNA copy number alteration data by Affymetrix 6. 0 SNP arrays were available for 466 cases. Somatic mutation data by whole exome sequencing were available on 463 cases. The TCGA PAM50 subtype classification was downloaded also from the same site. A sample list is provided in Additional file 2: (TCGA sample ID. xlsx).
Affymetrix U133A dataset A total of 1460 breast cancer profiles from eight breast cancer datasets that were uniformly profiled on the Affymetrix U133A microarray were obtained form the Gene Expression Omnibus (http://www. ncbi. nlm. nih. gov/ geo/). Because of partial overlap in several of these datasets, the selected unique samples from each GEO dataset are provided in Additional file 3: (Affy_sample_info1460. txt). A summary of the combined data set is given in Table  2.
Data preprocessing The TCGA gene expression data were used as provided without any additional processing. The segmented log2 DNA copy number values were trichotomized as deletion (score < -0. 3), amplification (score > 0. 3) or no effect. We aggregated the segment-level alterations at the gene level using RefSeq gene models. A gene was labeled as deleted or amplified if the entire transcribed region was in a deleted or amplified segment. Of the 45,918 RefSeq genes, 22,271 were amplified or deleted in at least one of the patients in the cohort and were used in the calculation of the Hamming distance. The somatic mutation data were provided at the gene level in the TCGA dataset (any mutation vs no mutations) and did not require any further processing. The Affymetrix U133A data files were compiled in a single dataset, normalized by MAS5, scaled to a target intensity of 600 and log2 transformed.
Principal component analysis (PCA) was applied to detect potential batch effects across the datasets. We used the empirical Bayes framework procedure ComBat  [28]  as implemented in the R package sva  [29]  to remove batch effects due to different datasets. We used PCA after batch effect removal to verify lack of batch effects after correction (Additional file 1: Figure  S4 ). A previously validated PAM50 algorithm for Affymetrix-based expression data was used to generate breast cancer subtype classification for each cancer in the dataset  [30]. Table  2  provides a summary of all the datasets used in this study and the subtype distribution within each dataset. The 161 normal-like cases were excluded from the analysis. Prior to analysis, the combined, batchadjusted dataset was filtered to remove probe sets in the lowest expression quartile and those in the lowest median absolute deviation (MAD) decile.
Collection and processing of KEGG pathways The 186 pathways used in this study were collected from the Molecular Signatures Database  [31]  (MSigDB; http:// www. broadinstitute. org/gsea/msigdb/index. jsp) and correspond to gene sets derived from the KEGG pathway database  [23]  (http://www. genome. jp/kegg/pathway. html). We processed the original list of 186 pathways to produce a reduced set of pathways that had minimal overlap in gene membership. We broadly followed the steps described previously  [32]. First, we calculated a matrix of pairwise pathway similarities derived from the hypergeometric distribution as follows. If N, N A and N B represent the genes in all pathways, in pathway A and in pathway B respectively, and these two pathways have n genes in common, then the degree of enrichment of pathway B in genes also included in pathway A is given by the hypergeometric distribution h(N, N A , N B )  [33].
Generation of simulated datasets We used the Umpire R language package to simulate realistic gene expression datasets with controlled within and between patient variation in gene expression  [20]. In each simulated scenario we considered 50 genes per case and 40 cases. We assumed all the genes are expressed in all patients and that the expression of each gene across patients follows a log-normal distribution with parameters μ g and σ P , i. e. log(X i ) ∼ Normal(μ g , σ P ). The mean log expression μ g is assumed to have a normal distribution across genes with mean μ 0 and standard deviation σ g , i. e. μ g ∼ Normal(μ 0 , σ g ). It was further assumed that the between-sample gene expression variance σ P follows an inverse Gamma distribution with hyperparameters α = 15 and β = 7, which correspond a mean σ P equal to 0.
Figure3Transcriptional diversity of 50 KEGG biological pathways within breast cancer subtypes from the Affymetrix dataset. A) Heatmap of mean dispersion distance within breast cancer subtypes considering genes in each of the 50 KEGG pathways. Pathways (rows) are ranked from the least diverse at the top to the most diverse at the bottom. Blue represents low and red high mean dispersion. B) Detailed expression heat maps for basal-like cancers showing the heterogeneity of gene expression for genes in the least heterogeneous (ribosome metabolism; top) and the most heterogeneous pathways (linoleic acid metabolism; bottom). Blue represents low and red high expression level. C) Distribution of pathway-based transcriptional diversity within each subtype. Pathway-level mean dispersion distances were calculated by bootstrap as described in Supplementary Methods.