dataesr/scientific-paragraphs-categorization huggingface dataset:
Distribution of real structures in function of the number of solutions returned This section is devoted to the distribution of real structures found by BiokoP in function of the number of solutions returned on our dataset of pseudoknotted RNAs, and to the comparison with Mod1 and Mod2, in order to show the pertinence of combining these two models on one hand, and to return several solutions on another hand. We extended Mod1 and Mod2 so that they return the kbest solutions, using the constraint  [31]  presented in the "Algorithm for finding the k-best Pareto sets" section. We refer to these extensions as Mod1 so and Mod2 so (so stands for sub-optimal). The results are reported in Fig. 5. BiokoP is made to return sets of solutions and all the solutions belonging to one Pareto set are not comparable.
Then, this experiment requires to rank the solutions of the Pareto sets returned by BiokoP in order to compare the solutions one against the others. The solutions of each Pareto set are ranked in the following manner: the solutions optimizing equally the two objectives, i. e. , the solutions closer to the diagonal, are better ranked. The results on the dataset of pseudoknotted RNA show that, as expected, BiokoP predicts more real structures than Mod1 and Mod2 (corresponding respectively to Mod1 so and Mod2 so for one solution returned). Indeed, The more there are solutions, the more BiokoP is likely to find the real structure, and with a fast increase in probability. We observe that after about 20 solutions returned (for about 2 or 3 Pareto sets), the number of real structures found seems to be stable, which supports the results of the previous section.
Considered software To evaluate the performances of BiokoP, we compare it with other methods predicting pseudoknotted RNA secondary structures that are able to return several solutions. To our knowledge, only two methods are available in the literature, namely pKiss  [4]  and McGenus  [5]. The principle of pKiss is to decompose the RNA sequence into every possible sub-words and to compute the MFE secondary structure of the decompositions. To reduce the search space, pKiss is based on the canonical rules which reduce the number of possible predicted pseudoknots (only certain canonical and kissing pseudoknots) and the redundancy thanks to a non-ambiguous dynamic programming algorithm. McGenus is based on a Monte Carlo algorithm which search for a minimum score which includes the energy and the genus of the secondary structure.
The genus expresses the complexity of a pseudoknot. McGenus performs a stochastic search that allows to find various types of pseudoknots. We also compare BiokoP with IPknot  [10]  and RNAsubopt from the ViennaRNA package  [22]. RNAsubopt predicts pseudoknot-free RNA secondary structures using an MFE algorithm to compute all the sub-optimal structures in an energy range. For the evaluation, we consider the first solution returned by IPknot and the 30 first solutions returned by BiokoP, pKiss, McGenus and RNAsubopt. IPknot (version 0. 0. 4) was executed with the Dirks and Pierce set of thermodynamic parameters and with the options -g 2 and -g 4. pKiss (version 2. 2. 12) was executed with the default parameters. We used the option -relativeDeviation to obtain up to 30 solutions for each RNA.
Overall results In this section are presented the results obtained on the dataset of 198 pseudoknotted RNAs. Table  2  reports the weighted means of sensitivities and PPVs in function of the number of solutions returned for BiokoP, pKiss and McGenus. We observe that BiokoP has better sensitivities than pKiss and McGenus and that, when the number of returned solutions increases, the gap between the sensitivity of BiokoP and the one of the other tools increases. Regarding the weighted means of PPVs, we observe that BiokoP outperforms McGenus. In Fig. 6  we present the weighted means of F 1 -scores obtained by each tool, in function of the number of solutions returned. BiokoP has higher F 1 -scores than pKiss and McGenus. The F 1 -scores of BiokoP are quite stable. There is only a decrease of 10 points going from 1 to 30 returned solutions, whereas there is a decrease of 15 and 18 points for pKiss and McGenus.
Results over optimal solutions The purpose of this section is to complete and to precise the results given by the previous statistics. It is not obvious to compare the optimal solution returned either by pKiss, McGenus or IPknot with only one solution obtained by an arbitrary ranking of the solutions of the optimal Pareto set given by BiokoP. Indeed, the solutions of a Pareto set are not comparable. We thus focus here on the comparison of the first solutions returned, i. e. the optimal solutions of BiokoP (the first Pareto set) and the optimal (one) solution returned by the other tools. Finally, we observe that the gap between the minimum and the maximum F 1 -scores of BiokoP can be important. This shows that BiokoP returns a diversified set of optimal solutions. The results for the H-type pseudoknots are very similar to the results of the entire dataset, which is not surprising since the H-type is largely represented in it (154 among 198 RNAs).
Results by pseudoknot types The BOIP of BiokoP has been modeled to be able to predict any kind of pseudoknots. This is confirmed by the results obtained that are very homogeneous. Indeed, the F 1 -scores of BiokoP are never lower than 70% for any number of solutions returned. This is not the case for pKiss and McGenus, for which we can observe that the results depend greatly on the pseudoknot type. In particular, they obtain F 1 -scores around 50% for the HLin type. Since the datasets of some pseudoknot types are small  (3 HHH, 26 HLout, 4 HLin, 11 LL and 154 H) , further experiments need to be done to confirm the results. Finally, when one wants to predict a secondary structure of an RNA, there is generally no information about the pseudoknot type. Therefore, it is better to use BiokoP, which allows to predict structures as close as possible to the real ones, regardless of the pseudoknot type.
Results on pseudoknot-free RNAs In order to determine if BiokoP is able to predict pseudoknot-free structures of RNAs, we evaluated it on our dataset of pseudoknot-free RNAs and compared the obtained results with those obtained by pKiss, McGenus, IPknot and also RNAsubopt which is a tool dedicated to pseudoknot-free RNA prediction. In Fig. 9 , we present the weighted means of F 1 -scores obtained by BiokoP, pKiss, McGenus, RNAsubopt and IPknot on the pseudoknot-free RNAs, in function of the number of returned solutions. BiokoP, pKiss and RNAsubopt give comparable results, showing that BiokoP and pKiss are both able to predict pseudoknot-free structures, unlike McGenus. We expected a bigger difference with RNAsubopt since this tool is made to generate only pseudoknot-free structures, but finally, it seems that BiokoP and pKiss do not suffer from a bias due to their purpose of returning pseudoknotted structures.
RNAsubopt finds the two first helices and some base pairs of the third one. The structure predicted by IPknot is identical to the one predicted by BiokoP, but its pseudoknot possesses an additional base pair. This example supports the previous results on the pseudoknot-free RNA dataset. As RNAsubopt predicts structures without pseudoknots, its results are better than the ones obtained by the other tools which can predict extra pseudoknots. However, the best solution of RNAsubopt is predicted at rank 25 which is far compared to the one of BiokoP.
Conclusion In this paper, we provide an original approach for predicting RNA secondary structures with pseudoknots that is based on the computation of several optimal and several sub-optimal solutions with respect to two different models. Our method is based on integer programming approach, that presents the advantage to be more flexible compared to the dynamic programming approach usually used in the RNA secondary structure prediction. In order to combine the two models, we developed a biobjective integer program (BOIP), and proposed a generic novel algorithm to compute the sets of optimal and suboptimal solutions for any BOIP. We applied it to the prediction of RNA secondary structures with pseudoknots, obtaining the software BiokoP. The experimental tests performed have confirmed the importance of considering sub-optimal solutions in addition to optimal ones.
Even if BiokoP is made to predict pseudoknotted RNAs, the F 1 -score results for the dataset of pseudoknot-free RNAs are close to the results of RNAsubopt. The results of BiokoP are widely homogeneous regardless of the pseudoknot type or the presence or not of pseudoknots. Indeed, the F 1 -scores are always higher than 70% and regardless of the number of solutions returned, unlike pKiss and McGenus. A drawback of our approach is the time complexity due to the need to optimally solve several integer programs. As stated before, this complexity is exponential, with respect to the length of the sequence. A perspective to decrease the execution time is to use parallelism. Indeed, due to its recursive shape, the generic algorithm we propose to solve BOIP is suitable to parallelization. Another idea to decrease the complexity of the BOIP would be to add constraints when some information on the structure are known: base pairs, presence or not of pseudoknots, type of pseudoknots.
reports the F 1 -score results for the optimal solutions of BiokoP versus pKiss, McGenus and IPknot, for each RNA of the dataset of pseudoknotted RNAs. The RNAs are sorted according to the ascending order of the maximum F 1 -score of BiokoP. For BiokoP, we report the maximum and minimum F 1 -scores of the set of solutions for each RNA. BiokoP finds a better solution than pKiss for 84 RNAs (among 198) and than McGenus for 103 RNAs while the optimal solutions found by pKiss and McGenus are better than the optimal solutions of the set generated by BiokoP for respectively 54 and 39 RNAs. The results show that BiokoP returns 61 better solutions compared to IPknot, while IPknot does not return better solutions compared to BiokoP. Returning several optimal solutions allows BiokoP to obtain the best solution more times than the other tools.
The BOIP of BiokoP has been modeled to be able to predict any kind of pseudoknots. This is confirmed by the results obtained that are very homogeneous. Indeed, the F 1 -scores of BiokoP are never lower than 70% for any number of solutions returned. This is not the case for pKiss and McGenus, for which we can observe that the results depend greatly on the pseudoknot type. In particular, they obtain F 1 -scores around 50% for the HLin type. Since the datasets of some pseudoknot types are small (3 HHH, 26 HLout, 4 HLin, 11 LL and 154 H), further experiments need to be done to confirm the results. Finally, when one wants to predict a secondary structure of an RNA, there is generally no information about.
A commonly used approach to establish a link between metabolic fingerprints and endogenous metabolism consists in applying enrichment analysis  (Chagoyen and Pazos, 2011; Persicke et al. , 2012)  on biochemical pathway collections provided by publicly available databases such as KEGG  (Kanehisa et al. , 2014). This technique aims at finding which metabolic pathways contain a significant number of metabolites that belong to the fingerprint. The notion of metabolic pathway is informative since it assigns functions to a set of reactions. Nevertheless, the definition of their functional boundaries (input and output compounds) often varies from one database to another  (Altman et al. , 2013)  potentially leading to different interpretations when analyzing the same fingerprint in various databases.
Adapting human genome scale metabolic network and defining transition probabilities The last two decades have seen an exponential growth of metabolic network reconstructions, which are made available through public databases  (Wimalaratne et al. , 2014)  or alongside articles  (Hucka et al. , 2003). However, for practical reasons, most methods designed to analyse metabolomic results in the context of those networks are database-dependent, restricting for example their use to KEGG  (Kanehisa et al. , 2014). In order to apply our analysis to networks coming from various sources as well as home-brewed networks, we propose a generic method that can be applied to any network described in the standard SBML format  (Hucka et al. , 2003)  with sufficient information on metabolites. We applied our method to Recon2 human genome scale metabolic network  (Thiele et al.
Since this information was not available for all metabolites in the Recon2 network, we automatically retrieved this knowledge form chemical databases using web services [PubChem, ChEBI  (Davies et al. , 2015) , HMDB  (Wishart et al. , 2013) ]. When only the InChI identifiers  (Heller et al. , 2015)  were available, we converted them into SMILES using the Chemistry Development Toolkit  (Steinbeck et al. , 2003). Finally, all substrate-to-product transitions that do not involve the preservation of at least one carbon atom between the source and the target were removed. The previously described steps allow building the graph on which PPR and PCR will be computed. The algorithm used to compute PPR considers every outgoing edges of a node to have the same probability to be traversed in the random walk.
Metabolic fingerprint of HE Metabolomics experiments have been conducted on CSF samples from 14 patients suffering from HE against samples from 27 control patients without any proven neurological disease  (Weiss et al. , 2016). The CSF metabolome was analysed by LC/MS [Orbitrap-Exactive and Q-Exactive Plus: Positive and negative ESI Scanning from m/z 75 to m/z 1000. Mass resolution: 100 000 FWHM. The relevant LC/HRMS features were kept by applying analytical filters from XCMS data matrix as followed: correlation between dilution factor of diluted QC (Quality Control, mix of all sample) and their area corresponded and area. r 2 > 0:7, Mean QC/mean Blank > 3, CV (QC) < 30%; Feature annotation using public databases, CAMERA R package and ESI-MS and HCD MSMS spectral inhouse database] and the discriminating fingerprint was built using multivariate [SIMCA-P software (version 11.
Mining literature to corroborate and enrich suggestions Social-network recommendation system efficiency is commonly evaluated by measuring the number of suggestions followed by a user during future web browsing. Unfortunately, this methodology cannot be applied to metabolic recommendation system assessment since recommendations in our case is not part of a decision process but is involved in data interpretation. More generally, assessing the quality of methods providing support to biological interpretation is still a key challenge in the field. In fact, since the disease mechanisms are still partially known, we do not have gold standard datasets (other validated biomarkers) to compare with our recommendation system suggestion. In order to establish a link between metabolites of interest and HE, we used the Metab2MeSH tool  (Sartor et al.
, 2012). MeSH (Medical Subject Headings) is a controlled vocabulary thesaurus Fig. 1. Effect of the number-of-product bias on compound graph transition probabilities. By overshadowing the reaction levels, seen on the bipartite graph (A), the use of compound graph will favour reactions that involve more products than other consuming reactions (B). The hybrid weighting policy (C) allows suppressing that bias hierarchically structured used to index scientific publications from the MedLine and PubMed databases. Metab2MeSH performs enrichment analysis to identify MeSH terms significantly associated with metabolite names, based on their occurrences in scientific publications. Compound names from the Recon2 model were converted to PubChem entry names using the Chemical Translation Service web service  (Wohlgemuth et al.
Discussion PR based methods have recently gained much interest for proteinprotein interaction networks. For instance, the Protein Rank  (Freschi, 2007)  is designed for protein functional annotation, using PPR favouring protein with a selected function as random walk seeds. Another example is the SubNet approach, providing a subnetwork extraction from interaction network based on PR scoring  (Zhang and Zhang, 2013). It uses a 'global' PR implementation where random walks can start from any compounds, but favour starts from nodes of interest by adding a constant parameter to bias Fig. 6. Kynurenic acid concentrations (arbitrary units) in CSF samples from HE and control patients. Data were obtained by LC/MS using a HILIC column and ESI mass spectrometry detection in the negative mode. Kynurenic acid identification level 1 according to  (Sumner et al.
, 2007)  (i. e. the same chromatographic retention time, accurate measured mass and MS/MS spectrum as those of the reference compound) Fig. 5. Association between brain and liver diseases MeSH terms and suggested compounds. Black and grey cells represent an unexpected number of co-occurrences between the compound name and the MeSH annotation in PubMed, defined accordingly to Smalheiser and Bonifield's metric, with an odds ratio threshold of 3. Only suggested compounds that are found by Metab2Mesh tool are represented the damping. Previously, Ivan and Grolmusz also proposed to use PPR for protein-protein interaction networks (Ivan and Grolmusz, 2011) and successfully retrieved cancer related proteins from proteomics data of melanoma patients. A similar approach, known as GeneRank  (Morrison et al.
Acknowledgements 'PageRank' is a registered trademark of Google. Inc. C. F would like to thank  Dr Dima Shepelyansky  ( Paul Sabatier University, Toulouse, France ) and  Dr Leonardo Ermann  ( CNEA, Argentina ) for their insightful courses about  Google Matrix  and fruitful discussions at the  Luchon Summer School. The authors also thank  Dr Sara Fontanella  for proof-reading the manuscript.
Several recent studies have portrayed DNA methylation as a new player in the recruitment of transcription factors (TF) within chromatin, highlighting a need to connect TF binding sites (TFBS) with their respective DNA methylation profiles. However, current TFBS databases are restricted to DNA binding motif sequences. Here, we present MethMotif, a two-dimensional TFBS database that records TFBS position weight matrices along with cell type specific CpG methylation information computed from a combination of ChIP-seq and whole genome bisulfite sequencing datasets. Integrating TFBS motifs with TFBS DNA methylation better portrays the features of DNA loci recognised by TFs. In particular, we found that DNA methylation patterns within TFBS can be cell specific (e. g. MAFF). Furthermore, for a given TF, different DNA methylation profiles are associated with different DNA binding motifs (e.
These systematic in vitro studies have unveiled the disparate mechanisms involved in TFs/DNA interactions with DNA methylation being a key leverage in this process. However, it has been shown that only a small fraction of TF binding sites (TFBS) found in vitro can be validated in vivo  (10) , due to a more complex and crowded nuclear environment involving the chromatin structure and non-DNA molecules such as RNAs, histones D146 Nucleic Acids Research, 2019, Vol. 47, Database issue and non-histone proteins  (11). Therefore, further efforts are required to characterise the methylation status of TFBS in the cell context. Most of the current TFBS databases such as JASPAR  (12) , HOCOMOCO  (13)  and GTRD  (14)  provide TFBS motifs derived from in vivo ChIP-seq datasets. Nevertheless, none of these databases have combined TFBS motifs with a precise depiction of their respective DNA methylation profiles yet.
Data collection and pre-processing Raw ChIP-seq datasets (fastq files) in different cell types were downloaded from ENCODE consortium and GEO databases (Figure  1A , Table  1 ). After read quality checking with FastQC (https://www. bioinformatics. babraham. ac. uk/ projects/fastqc/) and read trimming when required, using Trimmomatic (15), short reads were aligned to the human genome assembly hg38, using STAR  (16)  with the spliced function blocked, for its speedy and accurate aligning performance (BioRxiv: https://doi. org/10. 1101/053686). Relying upon the Irreproducibility Discovery Rate (IDR) pipeline  (17) , TF genomic occupancy regions with high consistency in replicates were called by MACS2 (18) using default parameters. We took advantage of Whole Genome Bisulfite Sequencing (WGBS) datasets to measure the genome-wide DNA methylation levels in each cell type at the CpG resolution.
WGBS datasets were obtained from ENCODE consortium and GEO databases (Table  2 ). Raw datasets underwent read quality check using FastQC, as well as quality and adapter trimming using Trim Galore. High quality short reads were aligned to hg38 using Bowtie2 with Bismark (v0. 16. 3) and default parameters  (19, 20). Methylation states in CpG context were extracted with the Bismark methylation extractor module. Whenever their correlation coefficients were confirmed over 0. 8, as recommended by ENCODE standards, WGBS datasets with biological replicates were merged to increase the sequencing coverage using methylKit  (21).
In house WGBS data generation We complemented the existing published datasets from EN-CODE and GEO by generating a WGBS of the HCT116 cell line (human colorectal carcinoma). HCT116 cells were grown in RPMI 1640 medium with 10% FBS. The DNA was extracted using Qiagen DNeasy Blood & Tissue Kits. After sample DNA testing, negative control DNAs were added into the initial DNA pool, and then fragmented into a 200-300 bp range using Covaris S220. Next, terminal repairing, A-ligation, methylation sequencing adapters ligation were performed on DNA fragments. Then, DNA fragments underwent a bisulfite treatment using the EZ DNA Methylation Gold Kit from Zymo Research. Library concentration was firstly quantified by Qubit 2. 0, before a dilution to 1 ng/ul followed by insert size verification using the Agilent Bioanalyzer 2100 and qPCR.
Profiling the DNA methylation landscapes surrounding EN-CODE ChIP-ed protein binding regions Integrative analyses of ChIP-seq and WGBS datasets were performed to profile DNA methylation landscapes surrounding the genomic occupancy regions across all EN-CODE ChIP-ed proteins. For each DNA binding protein, methylation scores (beta scores)  (22)  of CpGs within ±100 bp surrounding peak summits were collected using the in-tersectBed module of bedtools  (23). For each cell type, the distributions of the collected CpG methylation scores across all ChIP-ed proteins were shown in a heatmap, and proteins were clustered into different groups based on their DNA methylation landscapes surrounding their binding regions, using hierarchical clustering with the Euclidean distance. This representation allows us to classify ENCODE ChIP-ed proteins according to the DNA methylation levels of their respective binding loci (Figure  1B ).
Characterisation of DNA methylation levels within TFBS We employed the MEME-ChIP package with the default parameters to identify TF motifs in ±100 bp surrounding TF peak summits  (24). Since the actual binding sites for the ChIP-ed TFs are prone to central enrichment and maximum probability in occurrence at peak centres  (25) , each TF motif was determined based on its centrally enriched propensity profiled by CentriMo from the MEME-ChIP package. Then, the exact TF binding sites were localised using the FIMO module, while the methylation states of CpGs within the binding motif were assessed based on the WGBS dataset using the intersectBed module from bedtools  (23). In order to intuitively exhibit the DNA methylation levels within all TFBSs, we adopted a novel MethMotif logo combining the classical motif logo commonly used, with a cumulative bar chart describing the methylation level of each CpG present in the motif (Figure  1C, D ).