dataesr/scientific-paragraphs-categorization huggingface dataset:
For each grouping, we computed the mean LD-corrected pleiotropy score and used two-sample Student's t test to determine whether the mean was significantly different from the baseline. We found (i) that coding regions have higher pleiotropy scores than noncoding regions, (ii) that active promoters and enhancers have the highest pleiotropy scores and quiescent and heterochromatin have the lowest, (iii) that variants that control expression of more genes in more tissues have higher pleiotropy scores, and (iv) that genes associated with more than one model organism phenotype have higher pleiotropy scores Fig. 7  Replication analysis for the genome-wide pleiotropy study. We used 372 UK Biobank heritable medical traits as our discovery dataset, and independent datasets of 73 complex traits and diseases and 430 blood metabolites as replication datasets.
We believe this treatment is appropriate, however counterintuitive. Regardless, these kinds of isolated effects are fairly rare: in our dataset of 372 heritable traits from the UK Biobank, only 6% of variants (42,684 of 767,057) reach genome-wide significance for only a single trait. Indeed, it is unlikely by definition that a variant is associated with only one trait from a set of correlated traits, since we compute our correlations from observed association statistics. Third, we assume all genetic effects are additive and independent, and we do not model epistasis or other more complex genetic architectures. Our findings are in keeping with several recent studies that have found abundant pleiotropy in the genome  [2, 8, 9, 26, 27]. HOPS goes a step further than many of these studies by explicitly removing vertical pleiotropy between traits, which are indicative of fundamental biological relationships between traits  [8, 24, 28].
Collection of genome-wide association (GWA) summary statistics datasets First, we retrieved GWA publicly available summary statistics from 545 continuous traits in 361,194 samples from the UK Biobank  [17] , and 1403 binary traits from the same dataset calculated using SAIGE  [18, 19]. We used LD score regression to calculate heritability for each trait, using the liability scale for binary traits, and restricted the sample to only traits with a significant P value for nonzero heritability after Bonferroni correction. For every pair of traits with correlation coefficient between Z-scores r 2 > 0. 8, we additionally removed the member of the pair with lower heritability. This left a total of 372 traits. Second, we retrieved publicly available genome-wide association (GWA) summary statistics data for 82 complex traits and diseases  (Additional file 1: Table  S3 ).
For each dataset, we retrieved the appropriate variant annotation (build, rsid, chromosome, position, reference, and alternate alleles) and summary statistics (effect size, standard errors, P values, and sample size of the study). All variant coordinates (chr, pos) were lifted over to hg19 using the UCSC Genome Browser LiftOver Tool and aligned to the reference and alternate alleles retrieved from the Ensembl variation database. After performing the same pruning of highly correlated phenotypes, we were left with 73 traits and diseases. Third, we retrieved GWA summary statistics datasets from a GWAS of 453 blood metabolites in 7824 individuals  [25]. After performing the same pruning of highly correlated phenotypes, we were left with 430 metabolites.
Genome-wide pleiotropy study (GWPS) Using the two components of the pleiotropy score, we can run a genome-wide pleiotropy study (GWPS) which consists of computing two P values for each component of the score (P LD m and P LD n ) and for all variants genomewide. We computed the pleiotropy score separately for each of the three datasets described above (372 UK Biobank phenotypes, 73 traits and diseases, and 430 blood metabolites). Additionally, we computed the pleiotropy score on a subset of 372 traits with genomewide significant heritability as calculated by LD Score Regression  [20]  (univariate heritability significant after Bonferroni correction). The 372 UK Biobank heritable traits were used for discovery, and the 73 traits and diseases and 430 blood metabolites datasets were used for replication.
Study of polygenicity on horizontal pleiotropy To study the effect of polygenicity on horizontal pleiotropy, we first estimated the liability scale heritability of all 372 traits in our UK Biobank dataset using LD score regression, and stratified all traits into four equally sized classes of heritability, in order to control for the effect of high heritability separate from the effect of high polygenicity. Next, we estimated the polygenicity of the 372 traits using a corrected version of the genomic inflation factor λ c CG  [20]. The intercept of LD score regression minus one is an estimator of the mean contribution of confounding bias to the inflation in the test statistics. Therefore, we computed a corrected version of the genomic inflation factor by subtracting the quantity (intercept of LD score regression -1) from λ GC.
We divided variants into bins by number of eGenes (below 10, between 10 and 15, between 15 and 20, and over 20) and eTissues (below 30, between 30 and 35, between 35 and 40, and above 40). The reference sets used for these analyses were variants that are not annotated as eQTLs in any gene or tissue. Finally, we used model organism phenotypes measured by the International Mouse Phenotyping Consortium (IMPC)  [67]  and the Saccharomyces cerevisiae Morphological Database (SCMD)  [68]. To map ortholog genes from IMPC and SCMD to human variants, we used orthology annotations of gene orthologs, and eQTLs from GTEx. Thus, variants annotated as associated with a mouse or yeast phenotype are those that are annotated as cis eQTLs in any tissue for a gene whose ortholog in mouse or yeast is associated with that phenotype.
Genome-wide significant pleiotropy loci To detect loci with a genome-wide significant pleiotropy, we used the LD-corrected two-component pleiotropy score (P LD m and P LD n ) computed on the dataset of 372 heritable traits from the UK Biobank described above. We used LD clumping as implemented in PLINK to cluster linked variants, with an r 2 threshold of 0. 1, a distance threshold of 100 kb, and P value thresholds of 5 × 10 -8 for genomewide significance and 0. 05 for nominal significance. The resulting loci were assigned to genes using (1) localization of variants within a gene, as annotated by Ensembl Variant Effect predictor, and (2) annotation as a cis eQTL in any tissue, as annotated by GTEx. We submitted the resulting list of genome-wide significant genes to DAVID for enrichment analysis  [69] [70] [71].
Permutation-based null model for replication analysis In general, we should expect only 5% of loci to replicate by chance in each replication dataset; however, it is possible that this number might increase because of polygenicity in the underlying GWAS statistics and the resulting inflation in our pleiotropy score, which may cause substantially more than 5% of the genome to be assigned P < 0. 05. To correct for this, we performed random permutations of the whitened Z-scores independently for each trait and used these permuted Z-scores to compute our LD-corrected pleiotropy score components ( P LD m and P LD n ). This generates a null expectation that preserves the polygenicity and inflation within each dataset. For both datasets, our null model expected that 5% of loci for P LD m loci and 6% of loci for P LD n should replicate.
Availability of data and materials An R package implementing the HOrizontal Pleiotropy Score (HOPS) method is available on GitHub under the GNU Public License (GPL), at https://github. com/rondolab/HOPS  [72]. The current version of this repository at the time of submission has been deposited in Zenodo, at https://doi. org/10. 5281/ zenodo. 3462163. The dataset of summary statistics for the 372 medical traits from the UK Biobank and the pleiotropy scores computed from these summary statistics are also available in the same GitHub project at https:// github. com/rondolab/HOPS/tree/master/data. The summary statistics for 430 blood metabolites are available from the original publication where this dataset was reported  [61] , and the summary statistics for 73 human traits and diseases are available from the original publications where they were reported, as cited in Additional file 1: Table  S3.
Using a database containing curated resistance/virulence genes is essential for proper detection, identification, and interpretation of data. Databases need to be constantly updated to minimize the possibility of not detecting genes recently discovered. A number of databases are freely available, such as ResFinder, CARD, ARDB, MEGARes, ARG-ANNOT, and SSTAR, or a combination of databases (eg, JMI bioinformatics pipeline)  [44] [45] [46] [47]. In general, genetic markers and their respective annotations have not been standardized, and similar genes can have different annotations from database to database. As alluded to earlier, the organism's biological background is crucial for properly interpreting results. These high-throughput NGS and screening tools are very sensitive in detecting gene determinants, and these systems and databases will not indicate the location (chromosomal, plasmid, or both), relevance, or association of gene(s) with any resistance phenotypes.
Many resistance genes in the databases described here are chromosomally located and ubiquitous to several organisms and, when expressed in low levels, may not necessarily cause resistance phenotypes, or resistance may be limited to certain drugs within a class. Thus, the presence of an intrinsic gene may not be relevant for certain species-drug combinations, but potentially may be important in others. To complicate things further, an isolate can display a resistance phenotype when such intrinsic genes are overexpressed due to alterations in the promoter region, insertion of a promoter upstream, and/or gene mobilization into multicopy plasmids  [48] [49] [50] [51] [52] [53] [54] [55] [56] [57]. Therefore, certain genes detected in some species require careful interpretation. Many organisms express resistance phenotypes due to overexpression of intrinsic systems (efflux pumps), decreased or null expression of other alterations (porin channels needed for antimicrobial cell entry), and/or drug target alterations  [58] [59] [60] [61] [62] [63] [64] [65] [66] [67] [68].
Screen for β-Lactamase Genes by NGS and Bioinformatics Analysis Isolates displaying elevated MIC results (≥2 mg/L) for ceftriaxone, aztreonam, ceftazidime, and/or imipenem/meropenem were selected. A total of 11. 6% (n = 249) and 16. 1% (n = 97) of E. coli-causing UTIs and BSIs, respectively, met the MIC-based criteria, whereas 11. 0% (n = 53) and 13. 7% (n = 40) of K. pneumoniae isolates causing UTIs and BSIs, respectively, met the criteria (data not shown). Total genomic DNA was extracted using the Thermo Scientific KingFisher Flex Magnetic Particle Processor (Cleveland, OH), which was used as input material for library construction. DNA libraries were prepared using the NexteraXT library construction protocol (Illumina, San Diego, CA) following the manufacturer's instructions and were sequenced on a MiSeq Sequencer (JMI Laboratories).
Among KPC-producing K. pneumoniae (2 from BSI and 11 from UTI), 84. 6% (11/13) were ST258 (CC258). Seventeen and 38 unique CCs were noted in E. coli-causing BSIs and UTIs, respectively, and CC131 (or ST131) was the most common CC among BSI (53. 6%) and UTI (58. 2%) E. coli isolates (Figure  2 ). Twenty-four and 26 CCs were noted among K. pneumoniae isolates causing BSIs and UTIs, respectively. CC258 (28. 3%; all but 1 ST1326 were ST258) prevailed in UTI pathogens, whereas CC307 (15. 0%) was the most common CC among BSI isolates (Figure  3 ). Among E. coli isolates, the majority of bla CTX-M-15 (70. 3%) were detected within ST131 isolates, with similar occurrences among those causing BSIs (68. 4%) and UTIs (71. 1%) (data not shown). Most ST131 E. coli isolates causing BSIs carried bla CTX-M genes belonging to CTX-M group 1 (75.
Since phenological observations are prone to larges uncertainties  (Schaber & Badeck, 2002) , and considering the absence of quality flag in the study database, we identified and eliminated anomalous data deviating beyond two standard deviations from the local average. Sites with more than 20 years of phenological time series after outlier removal were used in the study. Finally, we obtained a total of 1107 observations spreading over 38 stations with LU and LC data from 1981 to 2012 (Fig. 1 ). More precisely, there are 226 observations in 9 stations for Hankow willow, 357 observations in 15 stations for Siberian elm, 191 observations in 8 stations for Black locust, and 333 observations in 14 stations for Simon poplar (Table. S1).
Meteorological data The daily mean air temperature and precipitation data used in this study are available from the China Meteorological Data Sharing Service System (http://data. cma. cn/). We tested several formulations for considering the role of water balance on leaf coloration, i. e. seasonal precipitation sums, P-PET (Precipitation minus Potential Evapotranspiration) indexes and a simple soil water balance (Table  S3 ). Finally, the most robust results were obtained with a low precipitation index (LPI), estimated as the number of events with ≥7 consecutive days without precipitation during the leafy season, which had previously been positively correlated to leaf senescence  (Liu et al. , 2018; Xie et al. , 2015). We calculated the daily photoperiod for each station as a function of the latitude of the station and the day of year  (Forsythe et al.
We tested the two hypotheses in all six models using phenological observation sites with more than 20 years of ground data. The original dataset was divided in two subsets of similar distribution using the Wilcoxon rank-sum test (p>0. 5) for each species. One subset including 80% of the data was used for model calibration (Calib Subset) and the remaining 20% were used to validate models (Valid Subset). Since the subsampling of the original data may influence the parameter estimation process and model fitting ability, we calibrated and validated each model 30 times, using different, randomly defined calibration and validation data subsample. The model parameterization was conducted using the simulated annealing method  (Chuine et al. , 1998). The simulated annealing algorithm requires parameter bounds which were specified as follows: Pstart (from minimum photoperiod of all sites to the maximum photoperiod of all sites for each species), Tb (from 10 ºC to maximum temperature of all sites for each species), x and y (from 0 to 2).
Model performance We assessed the quality of the model predictions with three performance indices: the modelling efficiency (ME), = and the corrected Aikaike information criterion (AICc), ME = 1 - ∑ 4 5 - AICc = n × ln A ∑ 4 5 -6 5 7 8 59: = B + 2C + 2C C + 1 = -C -1 Where 4 5 is the observation at data point i, 4 D is the average observation, 6 5 is the prediction at data point i, n is the number of observations, k is the number of model parameters. We considered the average RMSE, ME and AICc calculated from the 30 calibration/validation data splitting as performance metrics for each model. The differences of AICc (∆AICc) between each model and the model with the minimum AICc were calculated. The model with the minimum AICc was selected as the best model for each species.
Site-specific critical cold temperature sum fitting Fitting the critical temperature sums on each site (Table  2 ) decreased the calibration RMSEs as compared to calibration at the species-scale (compare with Table  1 ). The lowest mean sitespecific calibration RMSEs ranged from 5. 65 to 7. 94 days, with minor differences among species (Table  2 ). Yet, the calibration MEs of about 58% models (14/24 models) are negative, which implies these models perform worse than the Null models. Over validation dataset, all MEs are below zero for all species, revealing that all site-specific models applying to the independent dataset perform worse than the null model.
TPDM, Table1) showed the best performances with minimum AICc over calibration data for all species except Black locust. When tested on the validation (i. e. unknown) data, TPDM further displayed the most accurate predictions of LC dates with the lowest RMSE (from 7. 30 days to 9. 93 days) and the highest ME (from 0. 32 days to 0. 62) for all species. Models show consistent performances in calibration and validation for three species (Hankow willow, Siberian elm and Simon poplar), meaning that models giving the best performances in calibration dataset of LC also yielded the most accurate predictions in validation dataset of LC (Table1). For Black locust, we identify TDM as the best model, showing the minimum AICc over calibration data. The TDM predicts LC dates.
Acknowledgements The authors thank the  Meteorological Information Center of the China Meteorological Administration  for providing field phenological data and the phenological group of  Peking University  for digitizing the phenological records. We thank  Isabelle Chuine  for valuable comments and suggestions. We acknowledge support from the  National Natural Science Foundation of China  under grant no. 41471033 ,  41771049  and  31770516. Guohua LIU  also appreciates the  CSC (China Scholarship Council)  for a doctoral fellowship to study in  Université Paris-Sud  for one year and the  Ecologie, Systématique, Evolution lab  for giving her the opportunity to perform this study in  Université Paris-Sud. Y. H. F thanks the support from the  National Key Research and Development Program of China  ( 2017YFA06036001 ), from the  Beijing Normal University  and the  Thousands Talents Program for Young Professionals.
Species Models Calibration subset Validation subset N RMSE ME ∆AICc N RMSE ME TPDM, Table1) showed the best performances with minimum AICc over calibration data for all species except Black locust. When tested on the validation (i. e. unknown) data, TPDM further displayed the most accurate predictions of LC dates with the lowest RMSE (from 7. 30 days to 9. 93 days) and the highest ME (from 0. 32 days to 0. 62) for all species. Models show consistent performances in calibration and validation for three species (Hankow willow, Siberian elm and Simon poplar), meaning that models giving the best performances in calibration dataset of LC also yielded the most accurate predictions in validation dataset of LC (Table1). For Black locust, we identify TDM as the best model, showing the minimum AICc over calibration data.
1 ). This suggested that mRNA response to irradiation depended on its interactions with specific miRNAs. Fig. 1  Levels of Renilla luciferase mRNA in control and irradiated Me45 cells. Cells were transfected with a reporter gene plasmid containing a Renilla luciferase gene with the same promoter and coding sequence but whose primary transcript contained target sequences for different miRNAs versus those that did not contain the targets (empty). Renilla mRNA was assayed by RT-PCR and normalized to the Firefly mRNA to exclude differences in transfection efficiency between experiments. Error bars represent standard deviations of values obtained in three experiments, and p-values show the significance of differences between control and irradiated samples tested by a two sample t-test However, to understand the process in quantitative detail and build a mathematical model, we decided to employ microarrays, which deliver data on multiple mRNAs and miRNA.
Measurement of mRNA and miRNA levels Microarray data for normal cell lines: AG1522, PBMC, HCAEC, and cancer cell lines: DU145, SC3, MOLT4 analyzed in this paper were downloaded from ArrayExpress database  [32]. Additional file 1: Table  S1  includes ID numbers of individual experiments. mRNA and miRNA levels were estimated by Affymetrix (Human Genome U133A) and Agilent (G4870A SurePrint G3 Human v16 miRNA 8x60k) microarrays 1, 12 and 24 h after irradiating the cells, with control being the non-irradiated cells at 1 h, as described in  [30] ; these data are available in the ArrayExpress database (accession numbers E-MEXP-2623, and E-MTAB-5197 for mRNAs and miRNAs, respectively). The mRNA and miRNA datasets for Me45, K562 and HCT116 +/+ cell lines were already explored. In the previous publications authors focused on the relation between reactive oxygen species and miRNA  [1]  and bystander effects  [30, 33].
HCT116 -/-expression levels were not published previously. All datasets were normalized by the standard RMA method  [34]. Only mRNAs and miRNAs with expression above the noise level were selected for analysis using the GaMRed software, which uses Gaussian mixture models  [35]. The default parameter values were applied, i. e. , number of initial conditions 100, number of components 2 to 8, Expectation-Maximization (EM) algorithm precision level 10 -2 , and maximum number of EM interactions 5000. Two methods were used for choosing the number of components to remove, the "top three" rule which finds groups of genes with high, medium or low expression and removes others, or the "k-means" method which uses k-means clustering and statistical information to remove clusters of non-informative genes  [35].
Identification of sequence motifs which potentially influence mRNA levels Additional information, not included directly in the model but useful for interpretation of results, concerned mRNA properties. Adenylate-uridylate-rich (ARE) motifs were identified using NucleoSeq  [40]  based on 3′-untranslated regions (UTRs) from the RefSeq transcript database and the TTATTTAWW consensus sequence  [41]. Messenger RNA turnover times were derived from the experiments of Tani et al. on HeLa cells  [42]. The association of transcription factor (TF) response elements with genes was based on a ChIP-Seq experiment performed by the ENCODE project  [43] ; a TF-gene association table was created using a map of TF-binding site positions in the genome available as a USCS track (wgEncodeRegTfbsClusteredV3). A TF was assumed to regulate a particular gene if its binding site was located between 5000 bp upstream and 1000 bp downstream of the transcription start site, taking into account the strand directionality of the gene.
For example, the motifs targeted by a particular miRNA could vary or have different accessibility in different cell lines. To estimate the influence of the mRNA profile on parameter vector Δk and the quality of model predictions, we divided the dataset of mRNA expressions into two random subsets. These are called the training set and validation set. Parameters estimated using the training set were then used to predict the miRNA-dependent mRNA expressions using the validation set. Resulting miRNA-mRNA scatterplots and histograms are depicted in Fig. 5  (Me45 cell line) and in Additional file 5: Table  Fig. 3  Design of matrix C. Each line joining a mRNA and a miRNA symbolizes the interaction of a miRNA with a target motif in a transcript (one line is one interaction, dashed line for miRNA1, solid line for mRNA2).
Coefficients c ji are obtained based on miRanda3. 3a [34] , TargetScan  [35] , RNAhybrid  [36] , and NucleoSeq  [37]  algorithms (Table  1 ) In the Affymetrix microarrays a significant proportion of transcripts has been present more than once (43% in Me45 cell line, Additional file 6: Figure  S3 ). Since the multiple readings of transcripts of the same gene have highly correlated expression levels, this may result in inflation of the empirical correlation coefficient. Therefore, for the final analysis (Fig. 5  and Additional file 3: Table  S2 ) the superfluous expression levels have been removed. The adjusted (based on genes) and raw (based on transcripts) correlation coefficients are depicted in Additional file 3: Table  S2. As described in the Methods, we also ran our model on different datasets available in public databases including a variety of normal and cancer-related cell lines (see Methods for detail).