dataesr/scientific-paragraphs-categorization huggingface dataset:
OptiCLIP employs benchmarked tools in order to extract reproducible data to gain insights into the miRNA mode of binding. Within optiCLIP, we developed a novel method that addresses the reproducibility of the peaks and integrated a recently developed method to identify miRNAbinding sites that do not exclusively use perfect seed-match. To investigate the validity and efficiency of our new framework, we systematically applied optiCLIP to eleven Ago2 CLIP-seq datasets from mouse and human samples and identified 98,999 miRNA-binding sites. Briefly, the new method we use in optiCLIP consists in performing the peak calling step in each replicate separately and then quantify the overlap of peaks among multiple replicates through the Jaccard index (J-index) calculation. The consensus peaks are provided by the union of all the peaks that fulfill the reproducibility constrains.
In summary, here we proposed an optimized bioinformatic framework, called optiCLIP, to analyze Ago2 CLIPseq datasets, and showed by comparative studies its reliability and robustness in comprehensively identifying miRNAbinding sites. We used a systematic approach based on an adaptation of the game theory, to propose a quantitative model for stratification of Ago2 CLIP-seq peaks, which revealed that the synergy among miRNA-binding sites in each Ago2 peak finely controls target repression. This model can be also used for in silico stratification of predicted miRNAbinding sites. In conclusion, with increasing number of high-throughput data for miRNA targeting, our method is providing a powerful tool to interpret the dynamics and the complexity of miRNA-dependent post-transcriptional gene expression.
FUNDING  ANR  through the ' Investments for the Future ' [ ANR-11-LABX-0028-01  to LABEX SIGNALIFE]; FRM [ DEQ20140329551  to M. T. ]; F. S. was supported by ' UCA JEDI: Excellence Initiative ' - Université Côte d'Azur ; S. B. was supported by ' Fondation Lefoulon -Delalande -Institut de France '. Funding for open access charge: Inserm.
Results In this study, we propose a novel tool (named DeepGRN) for transcription factors binding site prediction based on the combination of two components: single attention module and pairwise attention module. The performance of our methods is evaluated on the ENCODE-DREAM in vivo Transcription Factor Binding Site Prediction Challenge datasets. The results show that DeepGRN achieves higher unified scores in 6 of 13 targets than any of the top four methods in the DREAM challenge. We also demonstrate that the attention weights learned by the model are correlated with potential informative inputs, such as DNase-Seq coverage and motifs, which provide possible explanations for the predictive improvements in DeepGRN.
Several popular methods for binding prediction, such as DanQ  [15] , DeeperBind  [16] , and FactorNet  [17] , are built on such model architecture. Recently, the concept of attention mechanism has achieved great success in neural machine translation  [18]  and sentiment analysis  [19]. It enhances the ability of DNNs by focusing on the information that is highly valuable to successful prediction. Combining with RNNs, it allows models to learn the high-level representations of input sequences with long-range dependencies. For example, long short-term memory (LSTM) models with attention mechanism have been proposed in relation classification  [20]  and sentence compression  [21]. Because of the input context similarities between language processing (sentences) and the DNA binding site prediction (sequences and results from massively parallel sequencing), similar approaches can be applied improve the performance of existing methods  [22] [23] [24].
Implementation Datasets from ENCODE-DREAM Challenge The datasets used for model training and benchmarking are from the 2016 ENCODE-DREAM in vivo Transcription Factor Binding Site Prediction Challenge. The detailed description of the pre-processing of the data can be found at https://www. synapse. org/#!Synapse:syn6131484/. For all TF and cell-types provided in the challenge datasets, the label of the binding status of the TFs is generated from ChIP-Seq experiments and used as ground truth. Chromatin accessibility information (DNase-Seq data), and RNA-Seq data are provided as input features for model training. For model training, we follow the rules and restrictions of the DREAM challenge: the models are trained on all chromosomes except 1, 8, and 21, and chromosome 11 is used as validation.
The model with the best performance in validation data is used for final prediction if no "leaderboard" dataset is provided by the challenge. The leaderboard data are available for some TFs for benchmarking, and each participant can test the performance on these TFs with up to ten submissions. Thus, if such data are provided, we pick the top 10 best models from the first step as an optional model selection step. The final performance of our models is reported based on the final test data that are used to determine the rank of the submissions in the challenge (Figure  S1  and Table  S1 , see Additional file 1). We use the similar organization of input features introduced by FactorNet  [17] : DNA Primary sequence, Chromatin accessibility information (DNase-Seq data) are transformed into sequential features and become the input of the convolution layers at the first part of the models.
Gene expression and annotations are transformed into non-sequential features and feed into the intermediate dense layers of the model (Details are described in the "Deep neural network models with attention modules" section). We also collected DNase and ChIP profiles for additional cell lines from the Encode Project (https://www. encodeproject. org) and Roadmap Epigenomics databases (http://www. roadmapepigenomics. org/data/) to improve the capability of generalization of our model. The performance of models trained with and without external datasets are evaluated separately.
DNase-Seq data Chromatin accessibility refers to the accessibility of regions on a chromosome and is highly correlated with TF binding events  [4]. DNase-Seq experiment can be used to obtain genome-wide maps of chromatin accessibility information as chromatin accessible regions are usually more sensitive to the endonuclease DNase-I than non-accessible regions  [31]. DNase-Seq results for all cell-types are provided in the Challenge datasets in the BigWig format. Normalized 1x coverage score is generated from the BAM files using deepTools  [32]  with bin size = 1 and is represented by a 1000×1 vector for each input bin.
Gene expression and annotation The annotation feature for each bin is encoded as a binary vector of length 6, with each value represent if there is an overlap between the input bin and each of the six genomic features (coding regions, intron, promoter, 5'/3'-UTR, and CpG island). We also include RNA-Seq data since they can be used to characterize the differences in gene expression levels among different cell-types. Principal Component Analysis (PCA) is performed on the Transcripts per Million (TPM) normalized counts from RNA-Seq data of all cell-types provided by the Challenge. The first eight principal components of a cell-type are used as expression scores for all inputs from that celltype, generating a vector of length 8. The processed data files for these features are provided in the FactorNet Repository (https://github.
Challenge. For each TF and cell-type combination, our attention model has better performance on 69% (9/13) of the prediction targets than Anchor  [40] , 85% (11/13) than FactorNet  [17] , 85% (11/13) than Cheburashka  [7] , and 77% (10/13) than Catchitt  [41]. Among all methods benchmarked, our method has the highest ranking in 7 out of 13 targets (CTCF/iPSC, FOXA1/liver, FOXA2/liver, GABPA/liver, HNF4A/liver, NANOG/iPSC, and REST/liver), with the best average score (0. 31) across all TF/ cell-types pairs (Table  2 ). To precisely evaluate the capability of deepGRN under the restrictions of the ENCODE DREAM Challenge, we also compared the performance of deepGRN trained using datasets provided by the challenge with four available features: Genomic sequence features, DNase-Seq and RNA-Seq data.
We found that the attention scores on the two DNA strands have a higher correlation (ρ=0. 90, σ=0. 79) than the saliency scores (ρ=0. 78, σ=0. 51) (Figure  5a, 5b ). Across all TF and cell-type pairs, we found that there is a positive correlation between the attention weights and normalized ChIP-Seq Fold (Figure  5c ), and such relationship is not detected globally in saliency scores (Figure  5d ). For all TF and cell-types in the benchmark datasets, we select at least four different genomic regions that have a clear ChIP-Seq peak signal in each target for demonstration. We show that the averaged attention weights put more focus on the actual binding region for each cell-type and these focusing points shift along with the shift of TF binding signals (see Additional file 3). Since the accessibility of the genome plays an important role in TF binding, it is expected to find high DNase coverage for those openly accessible areas that can explain the binding event detected by the ChIP-Seq experiment.
Motif detection over high attention scores regions For those positive samples without distinct DNase-Seq peaks, the patterns of genomic sequences are critical information for successful prediction. To test the ability of attention weights to recognize motifs that contribute to binding events from the genomic sequences, we use an approach similar to DeepBind  [13]. For the model trained for each TF, we first acquire the coordinates on the relative positions of maximum column sum of the attention weights from all positive bins in test datasets and extract a subsequence with a length of 20bp around each coordinate. To exclude samples that can be easily classified from patterns of DNase-Seq signal, we only select positive bins that have no significant coverage peaks (ratio between the highest score and average scores < 15).
Then we run FIMO  [42]  to detect known motifs relevant to the TF of the model in the JASPAR database  [43]. From the extracted subsequences, we discover motif MA0139. 1 (CTCF) in the prediction for CTCF/induced pluripotent cell and MA0148. 4 (FOXA1) in the prediction for FOXA1/liver cell. Figure  7a  and 7b show the comparison between the sequence logo of the motif rebuilt from the subsequences and the actual known motifs. We also plot the attention scores of the samples that contain these subsequences (Figure  7c , 7f) and the relative location of the regions with detected motifs in FIMO (Figure  7d , 7g). Furthermore, we show that these maximum attention weights do not come from the DNase-Seq peaks near the motif regions by coincidence since no similar pattern is detected from the normalized DNase scores in the same regions (Figure  7e , 7h).
Discussion The attention mechanism is attractive in various machine learning studies and has achieved superior performance in image caption generation and natural language processing tasks  [38, 44]. Recurrent neural network models with attention mechanism are particularly good at tasks with long-range dependency in input data. Inspired by these works, we introduce the attention mechanism to DNN models for TF binding site prediction. The benchmark result using ENCODE-DREAM Challenge datasets shows that the performances of our model are competitive with the current state-of-the-art methods. It is worth mentioning that the DNase-Seq scores are the most critical feature in the attention mechanism from our experiments according to the feature importance analysis. Many prediction tools for binding site prediction before the challenge, such as DeepBind or TFImpute, are not able to utilize the DNase-Seq data and are not as suitable as the four methods that we used for benchmarking in this study.
Selection of HA sequences and multiple alignment All available HA coding sequences (complementary DNA (cDNA)) have been uploaded from the NIAID Influenza Research Database (IRD) through http://www. fludb. org on 15 April 2020. Duplicates and laboratory-originating sequences were removed from the data set. In the case of H1-H16 subtypes, only sequences from avian hosts were selected. Moreover, HPAIV sequences descending from the unique H5 emergence in Guangdong in 1996 were also removed from the analysis to keep only independent H5 HPAIV emergence events. Only the HPAIV sequences with the accession numbers AF144305 and AF148678 were retained to represent this emergence event. Multiple alignments of HA cDNA sequences per subtype were performed by  using MUSCLE program (3. 8. 31 version). More precisely, HA nucleotide sequences were downloaded, translated, and then subjected to alignments based on the amino acid sequence.
Model selection, HA5/HA7 DNA tree inferences The phylogenetic analysis was conducted with the IQ-Tree-2. 0. 6 software using a maximum-likelihood (ML) approach  (Nguyen et al. 2015; Minh et al. 2020). The best DNA models were estimated and selected via Akaike Information Criterion  (Kalyaanamoorthy et al. 2017 ). HA5 and HA7 phylogenetic trees were constructed using a general timereversible nucleotide substitution model with a relaxed-Gamma Figure  1. Analysis pipeline to the study of the possible relation between cSL stability evolution and the HPAIV emergence event. The analytic procedure is constituted of five steps. First, HA5 and HA7 encoding cDNA sequences were downloaded from the NIAID IRD (Step 1). Sequences were subjected to multiple alignments, classified as LPAIV or HPAIV-HA sequence, thanks to data from literature/cleavage site motif.
Finally, clustering of different RNA structures families was performed on every H5/H7 LPAIV-HA sequences and their estimated ancestral sequences by using RNAclust software (Step 5). dGO: free energy (kcal/mol); MFE: minimal free energy structure; D: genetic distance; cSL: previously described conserved stem-loop at the HA1/HA2 boundary region. distribution: GTR + F + R6 and GTR + F + R5 for HA5 and HA7 trees, respectively (GTR-Free Rate Model). Topological support was estimated by the standard non-parametric bootstrap with 100 replicates.
Molecular clock Molecular clock analysis using the least square dating (LSD2) method was performed  (To et al. 2016). The collection date of every HA5 and HA7 DNA sequence was retrieved from the IRD website (uploaded on 15 April 2020). A home-designed Python script edited the date format to fulfil the IQ-Tree-command requirements. The best-used DNA model and previous ML-DNA trees were used to speed up calculation and fix the new time-tree global topology. The HA5 and HA7 time-trees with estimated date of every internal node (complemented with the confidence interval at 0. 95) were generated and compared to previous trees and data from literature. The pairwise genetic distances from generated time-trees were computed and used as data for the analysis of RNA structure stability evolution based on the genetic distance to HPAIV emergence events.
Analysis of the link between the RNA structure stability and the genetic distance to emergence event We used an analytic pipeline to address the question of how structural features of the RNA secondary stem-loop evolve along the evolutionary pathway of HA sequences. Most of the codes were realized with the R-software 4. 0. 2 (R Core Team 2020) and Python 3. 7 (Von Rossum and Drake 2009). First, phylogenetic trees of HA5 and HA7 sequences were analysed with the Ape 5. 4 package  (Paradis and Schliep 2019). A matrix containing the pairwise genetic distances between all terminal and estimated internal nodes was computed. By using a slightly modified version of the Ape function 'seq_root2tip', all lists of tree-nodes from the root to every tip were calculated, referred as 'evolutionary pathways'.
Data from other variables were also extracted: the list of every node from the tree-root to a specified terminal node (HPAIV or LPAIV-HA sequence or taxon), the HA subtype (H5 or H7), the evolutionary group membership (mechanism of genetic evolution involved for the evolutionary pathways leading to HPAIV emergence event or LPAIV group when the last sequence from the list of internal nodes was an LPAIV), the MFE-dG0 values (kcal/mol) from RNA structure predictions, and finally the genetic distance between every node to the specified terminal node (HPAIV or LPAIV-HA sequence). Line plots of RNA structure stabilities based on the genetic distance to tree-tip and regressions using the 'lowess' method were generated with ggplot2 3. 3. 2 and tidyverse 1. 3. 0 packages  (Wickham 2009).
Clustering analysis of predicted RNA structures from ancestors A complementary approach was performed to evaluate structural differences or similarities in the predicted RNA structures encompassing the conserved HA-CS encoding region. We performed a clustering of RNA structures from all LPAIV-HA sequences uploaded from the NIAID IRD and from the predicted ancestors (both on the vRNA and cRNA) by RNAclust perl-script  (Will et al. 2007). In brief, groups of sequences that share a common secondary RNA structure motif were identified. The base pair probability matrix of each sequence-secondary structure distribution was computed (by using RNAfold from Vienna RNA package). A sequence-structure alignment was subsequently calculated using LocARNA. Eventually, a hierarchical clustertree derived by WPGMA clustering of the pairwise alignment distances was generated, as previously described  (Will et al.
01 in terms of genetic distance to the tree-tips). When data from all evolutionary pathways were compiled and fitted using the loess regression method, we observed that the pattern of cSL structure stability evolved differently depending on the sequence sense (vRNA or cRNA) and as a function of the major lineage Table  2. HA sequences from known HPAIV emergence events retrieved from our phylogenetic analysis. For each HPAIV emergence, we supplied the country and the date (year) of isolation, the HA subtype, the HA-HPAIV sequence accession number, the sequence id, the phylogenetic-tree lineage membership, the closest LP HA sequence identified, the estimated LPAIV-HA sequence ancestor, and the genetic mechanism involved into the LPAIV-to-HPAIV transition.
Location (country) Year of HP emergence For some emergence events, the closest LPAI strain and/or the estimated ancestor (internal node harboring monobasic HA-CS motif) were impossible to be identified. The phylogenetic clade that includes the corresponding HPAI strain was only constituted of HPAI viruses. b At the time of the data collection, the HPAI sequence was not present in the used public database for this work. The corresponding HPAI sequence for this emergence was added to our initial data set, and only the closest LPAI strain was identified for the RNA structure clustering analysis. membership (Fig. 4 ). The pattern of vRNA cSL stability was similar between the Eurasian 'LPAIV-only' group and 'LPAIV-to-HPAIV' group by insertion. Similarly, the pattern of vRNA cSL stability was also comparable between the American 'LPAIV-only' group and 'LPAIV-to-HPAIV' group by substitutions.
Clustering of predicted RNA structures from LPAIV-HA sequences and predicted ancestors RNA molecules with similar thermodynamic stabilities can have different secondary conformations. Thus, to determine if specific RNA structure topologies could be identified among HA sequences evolving towards HPAIV, we analysed RNA structures using the RNA structure clustering software RNAclust. The analysis was performed on all LPAIV-HA sequences and the HA sequences estimated at each internal node of the timetrees, including the HPAIV progenitors. RNA structure clusters that contain ancestors specifically leading towards HPAIV Downloaded from https://academic. oup. com/ve/article/7/2/veab093/6415083 by guest on 11 September 2022 emergence events were highlighted by a blue rectangle on the computed cluster trees  (Figs 7, 9 and 10).
Two HPAIV-HA sequences that were not present in the NIAID IRD at the time of initial data retrieval were subsequently downloaded from the IRD or GISAID databases and included in the RNA clustering analysis. For A/chicken/Spain/6279-2/2009 (H7N7) (Spain 2009) and A/chicken/Netherlands/03008927/2003 (H7N3)  (NL 2003) , the phylogenetically and genetically closest LPAIV-HA nucleotide sequences were considered as their direct ancestors (Table  2 ). In the case of HA5 sequences, we identified eight different clusters containing HPAIV ancestors on the vRNA and eight clusters containing HPAIV ancestors on the cRNA (Table  3 ). All the ancestors directed towards HPAIV events via insertions were grouped in three and two different clusters on the vRNA The computed groups from H5 sequence clustering on vRNA (A) and cRNA (B) with additional information are shown: -The unique evolutionary pathway that the HP-ancestor present in the cluster is following.
The American lineage 'LPAIV-only' sequences were evenly distributed among the different clusters (Fig. 10 ). The ancestor of the H7N3-HP emergence in Chile (accession number AY303631) that did not cluster with the other ancestors of recombination events was located in cluster # 1181. We noted that the MFE consensus structure of this cluster appeared similar in shape to cluster #802 (Fig. 11 ). This analysis thus revealed that most of the American H7 ancestors (four over five events of recombination belonging to the American lineage) leading towards HPAIV emergence events via recombination shared similar cSL RNA structure topologies. The computed groups from H7 sequences clustering on vRNA (A) and cRNA (B) with additional information are shown: -The unique evolutionary pathway that the HP ancestor present in the cluster is following.
Funding This work was funded by a grant from the  Agence Nationale de la Recherche  ( ANR-16-CE35-0005-01 ) to  Romain Volmer. Pierre Bessière  was supported by a Ph. D. scholarship funded by the  Region Occitanie (France)  and by the  Chaire de Biosécurité at the École Nationale Vétérinaire de Toulouse (French Ministry of Agriculture). Gabriel Dupré  is supported by  Ph. D. scholarship funded by the  French Ministry of Research and Education.
Numerous satellite SSS datasets and validation approaches have been produced by the different SMOS, Aquarius, and SMAP ground segments, associated missions' expertise groups and users from the scientific community. This diversity is explained by the exploratory nature of the satellite SSS missions, their new instrumental concepts and the outcomes of the first spaceborne SSS measurements. As an example of this variety, there are about 10 satellite SSS data production centers for the three missions: the ESA/Data Production Ground Segments (DPGS) and Expert Support Laboratories (ESLs), the Centre National d'Etudes Spatial/Centre Aval de Traitement des Données SMOS (CNES/CATDS), the Institute of Marine Sciences/Barcelona Expertise Center (ICM-CSIC/BEC), the NASA/Jet Propulsion Laboratory (JPL), the NASA/Goddard Space Flight Center (GSFC), Earth and Space Research (ESR), the International Pacific Research Center (IPRC), Remote Sensing Systems (RSS), the Integrated Climate Data Center (ICDC), and finally, the ESA Climate Change Initiative Sea Surface Salinity project.