dataesr/scientific-paragraphs-categorization huggingface dataset:
However, possibly owing to the limited number of patients in both cases results did not reach statistical significance (NL versus LC P = 0. 7600; NL versus HCC P = 0. 6715; LC versus HCC P = 0. 1829). Small interference RNA (siRNA)-mediated NUPR1 gene silencing inhibits cell growth and increases cell response to sorafenib treatment. To investigate the functional role of NUPR1 in the regulation of tumor cell growth and chemoresistance, HCC cells were transfected with specific human NUPR1 siRNA and the effects on cell viability and chemoresistance were evaluated after sorafenib treatment. We found that NUPR1 knockdown decreased cell viability and colony formation and increased tumor cell sensitivity to sorafenib treatment (Figures  4a  and b  and Supplementary Figures  S1A  and B ). These data demonstrate that NUPR1 may have an oncogenic potential in HCC and it may be involved in drug resistance.
Transcriptomic analysis identifies gene expression changes following NUPR1 knockdown. To better understand the molecular mechanism(s) of NUPR1 action in HCC cells, we compared the global gene expression of the Hep3B shNUPR1 clone against the corresponding control clone (pSilencer). For each comparison, a subset of common differentially expressed genes was selected by initial filtering on confidence at P-value ≤ 0. 05, followed by filtering on expression level (≥ 3 fold). Using these stringent selection criteria, we found 273 genes upregulated and 153 genes downregulated in Hep3B cells following NUPR1 knockdown (Supplementary Table  S4 ). Pathway and network analyses generated using the Ingenuity Pathway Analysis (IPA) software showed the major functionally related gene groups, which were differentially expressed in the shNUPR1 clone compared with control cells (Supplementary Table  S5  and Figure  7 ).
,  36  gene expression analysis showed that RUNX2 suppression strongly downregulates NUPR1, as well as NUPR1-regulated genes RELB and IER3 (Figure  8e  and Supplementary Figure  S3C ). This result suggests that RUNX2 might regulate NUPR1 gene expression. Therefore, we performed an in silico analysis by using the open-source Jasper database (http://jaspar. genereg. net/) to predict the binding sites of RUNX2 in the promoter region of NUPR1 gene. This analysis allowed us to identify eight putative   S7  and Supplementary Figure  S4 ). Further studies are needed to clarify whether these binding sites are functional and are involved in the control of NUPR1 expression by RUNX2. However, these data strongly suggest the existence of a NUPR1/RELB/IER3/RUNX2 pathway that may act as an auto-regulatory loop.
Bioinformatics, gene expression profiling and data analysis. Nine HCC public data sets were analyzed for the pattern of NUPR1 gene mRNA expression by the Oncomine software (https://www. oncomine. com). Clinical and pathological features were available in four data sets (Supplementary Table  S1 ). Expression fold change, data normalization and the statistical significance of the differential analysis data were obtained using Oncomine algorithms. Open-source Jasper database (http://jaspar. genereg. net/) was used to predict the binding sites of RUNX2 in the promoter region of the NUPR1 gene. This bioinformatics tool includes a collection of transcription factor DNA-binding preferences, modeled as matrices. To this purpose, the sequence Chr16:28535 905-28541804 (GeneBank: Chromosome 16, NC_000016.
10 -GRCh38. p2 Primary Assembly), including 3003 bp downstream and 2897 bp upstream ATG codon, here considered as NUPR1 promoter region, was selected for the Jaspar prediction analysis. In this analysis, only the predicted sites of RUNX2 factor with a relative profile score threshold 480% was considered. Gene expression was analyzed using Agilent 44K Human Whole Genome Oligonucleotide Microarrays (containing ~44 000 genes) (Agilent, Palo Alto, CA, USA), as previously described. 11, 12, 36 All microarray experiments were performed in duplicate, using dye-swap during labeling. The GeneSpring software (Agilent) was used to generate lists of selected genes for the different statistical and visualization methods. Network and pathway analyses of the microarray data were completed using the IPA software (http://www.
Acknowledgements. This work was supported by grants from the  Italian "Ministero dell'Istruzione, dell'Università e della Ricerca (Ministry for Education, Universities and Research) -MIUR "  FIRB-MERIT  no  RBNE08YYBM  awarded to MC and GM. We thank  Mrs. Antonina Azzolina  for her excellent technical assistance,  Professor Maurizio Soresi  for his contribution in the statistical analysis and  Professor Claudio Tripodo  for his helpful suggestions in the immunohistochemical studies.
Transcriptomics Transcriptome analysis was performed using an Affymetrix Gene 1. 0 ST Array in two biological replicates for each cell line, providing 1 μg of extracted RNA for library production. For comparing BJ, BJEL, and BJELM cells' generated transcriptomes, we normalized all raw CEL files with the Affymetrix software Expression Console version 1. 3. 1 to calculate probe-set signal intensities using RMA algorithms with default settings. High reproducibility between the corresponding biological replicates was evaluated by calculating the Pearson correlation coefficient and skewness parameter between replicates and between BJEL and BJELM relative to BJ (Additional file 2: Figure  S2 ). To identify differentially expressed genes (DEGs), we compared BJEL versus BJ and BJELM versus BJ (in biological replicates).
Inferring transcription factors involved in deregulated gene expression during cell transformation For the selection of transcription factors (TFs) associated with particular co-expression pathways we used the CellNet database of TF-target gene (TG) associations. We first selected TFs that are associated with >10 % of DEGs that constitute a given co-expression pathway. Then we selected TFs with promoter-associated RNA polymerase II (RNA Pol II), which gave rise to 142 TFs. Finally, we assessed the relevance of these TFs in distinct co-expression pathways using a hypergeometric distribution test with subsequent hierarchical clustering (Fig. 2 ).
Quality control of sequence data Sequence-aligned files were then qualified for enrichment performance using the NGS-QC Generator tool  [7]  (http://www. ngs-qc. org/). This methodology provides enrichment quality descriptors in a scale ranging from triple A (best) to triple D (worst). Based on this quantitative method, all ChIP-seq datasets described in this study had at least a triple B quality grade to ensure that only high quality datasets were used for downstream integrative analyses.
Enrichment pattern detection and normalization of ChIP-seq intensity profiles Relevant binding sites in all ChIP-Seq (except the H3K27me3 dataset, which was analyzed with the SICER tool  [8] ) datasets were identified with MeDiChI-Seq  [9]. To enable a comparison of ChIP-seq profiles of the same target between different cell lines, a normalization procedure over profile global amplitudes prior to further analyses was applied using a Quintile-based approach. Briefly, we calculate read count intensity for a non-overlapping window of 100 bp across the genome and then normalize these intensities using quantile normalization from the "limma" package. Quantile normalization is a rankingbased approach where calculated read count intensities are sorted and ranked for each sample. The corresponding ranked values between samples are adjusted into a mean value.
The impact of normalization was assessed using MA plots before and after normalization. First, we normalize all datasets associated with a given target; then normalized target datasets are brought to the same scale via a z-score normalization. A detailed description of this quantile normalization procedure, which is applicable for a variety of ChIP-seq and enrichment-related next-generation sequencing datasets and is available as part of Epimetheus, a userfriendly dedicated tool, is going to be presented in a further publication (in preparation).
Gene regulatory network reconstruction To provide a comparative view of the signal transduction events governing the cell transformation in the stepwise model system, we reconstructed a gene regulatory network (GRN) by combining several layers of information from three different databases: (i) the MiMI, which contains protein-protein, DNA-protein, and other interaction data, querying the interactions only between DEGs  [10] ; (ii) CellNet, a collection of directed TF-TG interactions  [11, 12] , where the TFs listed in the CellNet Gene Regulatory Network (GRN) collection were associated with genes differentially expressed in the BJ/BJEL/ BJELM model system; (iii) several publically available MYC-targeted ChIP-seq datasets (see "Methods"). The integration of DEG-related interactions in the Cytoscape platform (version 2.
Definition of hubs and bottlenecks was performed using topological metrics, such as MNC (Maximum Neighborhood Component), DMNC (Density of Maximum Neighborhood Component) and Bottleneck  [13]. This reduction process generated GRNs composed of 253 nodes and 2657 edges. The organization of the reduced GRN and its visualization were performed with the Cytoscape package Cerebral  [15]. As part of the visualization options in Cytoscape, the differential expression levels in BJEL and BJELM cells per node were presented in a heat map format such that the transcriptome dynamic changes could be visualized. The changes of node color for groups ii, iii, iv, and v in Fig. 4a  and 4b indicate the change in expression of co-regulated genes during the transformation process.
Analysis of publically available ChIP-seq datasets used for targeting MYC enrichment The following ChIP-seq datasets from the Gene Expression Omnibus (GEO) were used to identify MYC enrichment sites: GSM1088663 (HeLa cells), GSM896988 and GSM1000576 (BJ cells), and GSM748557 (NHEK cells). The raw sequencing files were aligned with Bowtie using default parameters and processed with MACS for peak annotation. A threshold of -log10(p value) ≥300 was applied to select peaks with high confidence.
Finally, panel 3 corresponds to the chromatin state cluster that characterizes the group of downregulated genes in BJEL and BJELM cells; the promoters of these genes are in the active state in BJ cells but lose all marks in the BJEL and BJELM cells described in the CellNet database  [11, 20]. This led to the identification of 142 TFs (Additional file 7: Table  S2 ), of which 42 are differentially expressed during stepwise tumorigenesis (Fig. 2a ). Sorting these TFs for their association with particular co-expression pathways led to a clustering of groups of TFs that were preferentially associated with one or more co-expression pathways (Fig. 2b ). According to the hierarchical clustering of TF-specific association with co-expression pathways, we could distinguish at least three subgroups: TFs that are preferentially associated with pathways i and iv, with pathways ii and iii, or with pathways v and vii.
Dashed lines separate the GRN into seven segments corresponding to seven (i to vii) gene co-expression pathways. Functionally related genes are circled under an enriched GO term (DAVID, p < 0. 05) within 1. 5 kb up-and downstream of each TSS for all DEGs, comprising a total of 7616 transcripts (Additional file 10). We detected 26 different combinations of histone marks at DEG promoters and classified them into seven chromatin states (Additional file 11: Figure  S7 ): (a) active (RNA Pol II and at least two histone marks of active transcription are present); (b) weakly active (RNA Pol II and only one histone mark of active transcription); (c) transcription-prone (at least two histone marks of active transcription but no RNA Pol II); (d) bivalent (any of states a to c but accompanied by repressive H3K27me3 marks); (e) ambiguous (only one histone mark or RNA Pol II alone); (f) empty (no signal); and (g) repressed (only H3K27me3).
Reconstruction of GRNs To provide a comprehensive view of the signal transduction events governing the cell transformation in the stepwise human cellular tumorigenesis model, we reconstructed a GRN integrating gene interaction data from publically available databases with gene expression data from our experiments (see "Methods" for details). This integration process resulted in the establishment of a comprehensive GRN of 1265 nodes and 5327 edges. To explore the functionally most relevant aspects of the reconstructed network, we reduced its complexity by applying topological criteria to identify highly connected ("hubs") and key connector nodes ("bottlenecks") that are relevant to the investigated signal transduction processes  [14]. The reduced network of 253 nodes and 2657 edges (Fig. 4a, b ) shows the connectivity between the major nodes, which are possibly functionally involved in the cell transformation process.
In fact, as illustrated in Fig. 5a , normal BJ and immortalized BJEL cells are not able to grow in an anchorage-independent manner, while the overexpression of c-MYC conferred onto BJELM cancer cells the ability to proliferate under these conditions (Fig. 5a, b ). Importantly, CRISPR-Cas9-mediated individual inactivation of all tested "downstream" factors of Myc (CHD7, DHX33, GTF3C4, NOLC1, and PRMT3 genes) resulted in a drastic drop in the ability of BJELM cells to form colonies on soft agar, while BJELM cells, as well as mock-transfected BJELM cells ("siGLO"), showed efficient colony formation in soft agar (Fig. 5a ). Together our data reveal that each of these factors plays a critical role in mediating key oncogenic effects of MYC overexpression in this isogenic model system.
As could be expected, GO analysis of DEGs yielded GO terms indicative of cell transformation. Indeed, a marked increase in rRNA synthesis is a general attribute of many cancers  [40, 41]  and rRNA transcription was shown to be stimulated by c-MYC  [49]. This correlates with our observations showing MYC-mediated upregulation of genes functionally related to RNA biogenesis, such as DHX33, HEATR1, NOLC1, and others (Fig. 4a , b, RNA processing functional group). Notably, disease-oriented functional annotation clustering showed that DEGs in the stepwise BJ transformation system comprise genes that are implicated in different types of cancer, such as breast or bladder cancer. In addition, we used cBioPortal  [50, 51]  to see if genes identified in this study can be correlated with publically available datasets of human cancers (cBioPortal is an exploratory analysis tool that, among others, hosts TCGA (The Cancer Genome Atlas) datasets ready for network analysis).
GRNs of tumorigenesis Cellular phenotypes are determined by the temporal regulation and dynamics of networks of co-regulated genes. Thus, elucidating GRNs is crucial for understanding of normal and cancer cell functioning. Until today only a few studies have addressed this issue-for example, the elucidation of the P53 regulatory network  [53]  or the analysis of locus expression signatures from retroviral insertion-induced tumorigenesis  [54]. To perform a systematic analysis of GRNs underlying tumorigenesis, we used a novel combinatorial approach by (i) integrating transcriptome data during transformation steps with chromatin state dynamics, (ii) complementing this with an analysis of CRMs involved in this process, and (iii) inferring key transformation-related TFs by using a database of established TF-TG associations from multiple human lineages.
Acknowledgements We are very grateful to  Dr William C. Hahn  and  Dr Robert A. Weinberg  for providing plasmids and cells. We thank former members of our laboratory,  Dr Pattabhiraman Shankaranarayanan  and  Dr Jelena Vjetrovic , for providing BJELM transformed cells. We thank  Dr Valeria Pavet-Portal  for help in setting up the stepwise cell system, bioinformaticians  Matthias Blum  and  Pierre-Etienne Cholley  for expert advice,  Benjamin Billore  and  Michele Lieb  for technical support, and all other members of our laboratory for discussions. We thank as well the  IGBMC Microarray and Sequencing platform (France Génomique consortium- ANR-10-INBS-0009 ).
However, our data suggest that the AGM in does not offer the right microenvironment for the development of the second HSC type. Only the AGM explant culture model facilitates study of the two HSC types and their regulation. Others have shown the existence of pre-HSCs in reaggregate AGM cultures or in the early FL  (Taoudi et al. , 2008; Kieusseian et al. , 2012) , and neonatal engrafting HSCs have been described in the E9 yolk sac  (Yoder and Hiatt, 1997). These may represent the precursors to the non-BMP-activated (Hh/VEGF-responsive) HSCs. It will be interesting in future studies to examine these cells and the early hematopoietic tissues in the context of BMP, as well as Hh signaling through fate-mapping approaches. Our RNA sequence datasets show that transcription factors known to be involved in hematopoietic cell development are expressed in both the BMP-activated and non-BMP-activated AGM HPSC fractions, suggesting that these genes are not exclusively regulated by BMP signaling.
In contrast and importantly, Sox17 a known target of the BMP signaling pathway  (Teo et al. , 2012)  is strictly expressed in the BMP-activated HPSC and EC AGM fractions (Figure  5E ). Since Sox17 is required for endothelial to hematopoietic cell transition and since all AGM HSCs in vivo are BMP activated  (Crisan et al. , 2015) , our data suggest that BMP-activated HSCs arise from BMP-activated ECs. This is further reinforced by the dendrogram generated from our RNA sequencing datasets, showing the close relationship of HPSCs to BMP-activated ECs (Figure  5D ). In conclusion, we show here that the BMP and Hh/ VEGF signaling pathways control two HSC types in the Figure  5. Transcriptional Differences between E11 AGM BRE GFP + and GFP À Cell Types (A-C) Expression (FPKMs from RNA sequencing analysis, accession number GEO: GSE76253) detected in sorted BRE GFP + and GFP À hematopoietic progenitor/stem cell (HPSC; CD31 + cKit + ) and endothelial cell (EC; CD31 + cKit À ) fractions of AGM in and AGM ex.
(A) Hematopoietic transcription factor gene. (B) Hematopoietic marker gene. (C) Endothelial marker gene expression in HPSC. (D) Dendrograms showing genetic relationships between HPSCs GFP +/À and ECs GFP +/À in the AGM in and AGM ex as measured by cluster analysis of RNA sequencing datasets. (E) Sox17 expression in GFP +/À ECs in AGM in and AGM ex. AGM explant culture. BMP-activated HPSCs differ in their intrinsic molecular program from the non-BMP-activated Hh responsive HPSCs, and the molecular program of cells within the AGM niches changes upon explant, thus implicating the microenvironment in the control of the distinct HSC types. This knowledge is of high interest to human health, since Hh antagonists are used today as anti-cancer and anti-leukemic drugs  (Von Hoff et al. , 2009; Amakye et al.
Transplantation Assay AGM cells were injected intravenously into female wild-type mice (129SV 3 C57BL/6 or C57Bl/6) irradiated with 900 rad (split dose, 137 Cs source). 2 3 10 5 recipient background spleen cells were coinjected with AGM cells. Unsorted or sorted BM cells from primary recipients were injected into irradiated secondary recipients. For details of AGM embryo equivalents/cell numbers injected and the number of experiments performed, see Table  S2. Chimerism was quantified by DNA PCR for gfp. DNA normalization (myoD) and comparison with gfp contribution controls was performed with Image Quant software. See Table  S1  for PCR primer sequences. Mice showing >5% donor chimerism were considered repopulated.
RNA Preparation, qRT-PCR, and RNA Sequencing For qRT-PCR, total RNA was extracted using TRIzol Reagent (Invitrogen) and quantitated by a NanoDrop 8000 Spectrophotometer (Thermo Scientific, NanoDrop Technologies) or a 2100 Bioanalyzer (Agilent Technologies) with RNA 600 Pico chips (Agilent Technologies). cDNA was generated using SuperScript II Reverse Transcriptase (Invitrogen), from 1 mg of starting material. qRT-PCR primers (Table  S1 ) were used with SYBR Green (Invitrogen) and Platinum Taq DNA Polymerase (Invitrogen), samples run on a CFX96 Real-Time System C1000 Thermal Cycler (Bio-Rad), and analyzed with Bio-Rad CFX Manager v2. 0 software. b-Actin was used as the internal reference. For RNA sequencing, RNA was isolated with the mirVana miRNA Kit and prepared according to SMARTER protocol for the Illumina HiSeq2000 sequencer.
Sequences were mapped to the mouse (NCBI37UCSC/mm910) genome and FPKMs were calculated using Bowtie (v2. 2. 3), TopHat (v2. 0. 12), and Cufflinks (v2. 2. 1). Differential expression was analyzed using Cuffquant with fragment-bias and multi-read corrections, and normalized across all samples using Cuffnorm with geometric library-size normalization  (Trapnell et al. , 2013). Difference in expression was calculated as GFP + FPKM+1/AGMin FPKM+1 with threshold of more than 2-fold change. Genes with more than 2-fold higher or lower expression in each comparison were applied to Enrichr tool. Enrichr output (and calculated p values from chi-squared tests) were imported into R and corrected for multiple testing (FDR), and a threshold of FDR < 0. 05 was used. Hierarchical clusterings were performed using log10-transformed FPKMs of the top 10% most variable genes in the R statistical package using hclust() and dist() commands with default parameters.
Heatmaps were generated using the R heatmap. 2 package with row-scaling. For GSEA, ratios were calculated for AGMex HSPC GFP À versus GFP + FPKMs(+1) and used with GSEA (version 2. 0. 13) pre-ranked method using default options  (Mootha et al. , 2003; Subramanian et al. , 2005). RNAseq data are publicly available at NCBI GEO under accession number GEO: GSE76253.
ACKNOWLEDGMENTS We acknowledge our laboratory colleagues and  Derk ten Berge  for critical discussions of our data. We thank  Polynikis Kaimakis ,  Emma de Pater , and  Catherine Robin  for occasional technical assistance, and our funding organizations:  EMBO Longterm Fellowship  ( ALTF 260-2009 ),  ZonMW Dutch Medical Research Council  ( VENI 916-12-088  and  VICI 911-09-036 ),  FES NIRM  ( Dutch Innovation Grant ),  NIH  ( RO37 DK54077 ),  Erasmus   MC Fellowship  ( 103. 494 ), Landsteiner Society for Blood Research  ( 1109 ), and  KNAW  ( Dutch Royal Akademie Research Master Assistantship ).
Ancillary data Below 20 GHz, several radio maps are available. Reich et al. (2004)  have presented a preliminary map at 1. 4 GHz, but the combination of instrumental artefacts and calibration uncertainties in the survey do not meet our stringent requirements on the error budget, which also caused us to drop a number of Planck datasets as described in Planck Collaboration X (2016). We do use the 408-MHz map of  Haslam et al. (1982) , as re-processed by  Remazeilles et al. (2015)  to reduce instrumental striping and provide better point-source subtraction than previous versions. The noise in the  Haslam et al. (1982)  map is not characterized in detail, but in any case is much less than residual scanning artefacts. We assume an rms error of 1 K, including both effects. We also assume a 10 % uncertainty in the overall intensity scale.
Comparison with WMAP models Here, we discuss the main similarities and differences between the low-frequency foreground component maps from our Commander analysis of the combined Planck and WMAP datasets (Planck Collaboration X 2016), and from the final WMAP analysis  (Bennett et al. 2013a). The WMAP team have published sets of componentseparated maps using a Maximum Entropy Method (MEM) and a Markov Chain Monte Carlo (MCMC) technique on a pixel basis  (Bennett et al. 2013b ). We have compared the synchrotron, free-free and AME maps produced by these methods with the Commander maps through T -T plots with HEALPix N side = 64. Table  3  gives the key model parameters of the WMAP models that we used. Point sources from the PCCS2 catalogue brighter than 1 Jy at 28. 4 GHz have been masked out.