dataesr/scientific-paragraphs-categorization huggingface dataset:
The first approach is Mutual Information (MI) which measures a statistical dependence between two variables  8. It is based on density function estimates and has been shown to perform well with non linear relationships  9. The second approach which relies on integrating multiple transcriptional associations is Partial Correlation (PC). First-order PCs may be simply calculated using the following equation with R xy R yz and R xz the simple correlation coefficients between genes x and y, y and z, x and z respectively  10  , but such computations may be very time consuming for large datasets. In this case, PCs should rather be calculated by multiple linear regression including a feature selection step  11. In this way, PCs aim at explaining a gene's expression profile by a small number of strongly correlated genes after eliminating those less correlated that do not significantly explain this gene's expression profile.
The two last methods are Correlation Coefficients (CCs), either Pearson CC (PCC) or Spearman CC (SCC), which are the classical estimators of linear transcriptional relationship among genes  11, 12. CCs are 2-dimensional distance measurements because a CC between two genes does not take into account the expression of the remaining transcripts in the whole transcriptome. To improve CCs, it has been proposed to recursively correct expression values according to a neural network algorithm with a weighted CC matrix until an equilibrium is reached between input and output expression  13. Used on a small dataset, this method was useful to uncover interaction networks even for genes included in background cellular activities. R xy z R R R R R ( )( ) ( 1 Neural networks have however been reported to be computationally expensive for larger datasets  14.
MR and HRR are thought to be more integrative than unranked CCs because they depend on other CC values around that of a gene pair. It might also be expected that ranking CCs partially correct for the "range restriction effect" observed for CCs leading to robust correlations for high variance genes only  18. Although not as robust as supervised methods, unsupervised methods can efficiently capture relevant gene associations as previously shown  9. These authors have shown that non parametric CC and MI calculations were more efficient than PCC on a small dataset. Among other unsupervised methods, SCC calculations have been similarly shown to outperform other distance measurements in Human expression data  19. In this case, SCC were calculated from RNA-seq or microarray data in order to construct several smaller networks subsequently aggregated to yield the final network.
We firmly believe that genome-scale networks inferred with CCs, especially when combined with a ranking procedure, are helpful to find new associations between genes. Although CCs are not efficient in detecting non linear associations  9  , gene-to-gene relationships have been predicted to be essentially linear  20  suggesting that CCs are valuable distance measurements. To date, there is no clear evaluation of how ranked CCs affect genome-scale network reconstruction with RNA-seq data in comparison with other unsupervised methods. We evaluated ranked CC, raw CC, MI, and PC performance in global and targeted network construction using Arabidopsis microarray and differentially processed RNA-seq expression data (Fig. 1 ). By using unsupervised methods, our aim was to highlight transcriptional associations between genes from an expression dataset.
Querying genome-scale datasets is a useful preliminary step in identifying genes involved in biological pathways. Because correlation does not imply causality, the study of regulatory interactions or sequential orders through supervised methods is beyond the scope of our work. Although specific dedicated conditions are useful to find groups of intercorrelated genes  21  , large scale expression matrices were used here in order to increase the number of different experimental conditions and obtain robust averaged correlations. This approach should be useful for non-model species which global transcriptional landscapes have to be studied. Performance was measured as network ability to capture biologically relevant gene associations found in a Gene Ontology (GO) annotation reference set but also to correctly cluster guide genes in PLC.
Global network quality was first evaluated according to the different dataset × distance measurement combinations. The resulting global networks were next interrogated in PLC analyses with five different guide gene sets corresponding to four different metabolic pathways and one signaling pathway. Whereas metabolic pathways have relatively clearly defined and partially linear partitions, signaling pathways usually involve post transcriptional regulations and a more intricate organization, which might render gene transcriptional associations less evident. We looked at the dataset × distance measurement combinations optimizing pathway reconstruction and maximizing co-occurrence quality between microarray and RNA-seq networks. Our results show that, of the six methods evaluated, PCC ranked with HRR generated the best biologically relevant networks according to initial guide gene representation and clustering in distinct modules.
Results Inferring global co-expression networks and comparing correlation measurements. Large co-expression networks were obtained by varying the confidence threshold (correlation value above or rank below) within lists containing the 10 million best gene pairs from eight different datasets and six data measurement combinations (Fig. 1 ). Each of the 10 million best pair lists was filtered at different confidence thresholds (1, 5, 10, 20, 40, 60 or 80% best pairs from these lists) to evaluate the effect of network size on performance. Expression datasets included a microarray-based expression matrix and seven RNA-seq based expression matrices normalized with different methods to evaluate their effect on network inference: transcript per Million (TPM), log2 TPM, sample scaled (ss) TPM, ss log2 TPM, raw counts, variance stabilized transformed (VST) raw counts and VST-TPM.
The six distance measurements were: raw PCC, raw SCC, PCC-HRR, SCC-HRR, PC and MI. Each network performance was considered as a network ability to capture edges corresponding to functional associations found in the GO reference dataset and was evaluated in 4 different ways (Fig. 2 ): GO term enrichment (GO terms that are significantly enriched with gene pairs from the co-expression network), a ROC curve constructed with TPR and FPR calculated for each confidence threshold and two ROC analyses based on the GBA concept, an average 3-fold cross validated neighbor voting (NV) AUROC and a global AUROC. AUROCs correspond to Area Under Receiver Operating Characteristic curves calculated for every network either from each GO (with three test sets obtained after hiding part of the gene labels, NV AUROC corresponding to the average of AUROCs for all GO terms) or the whole annotation dataset (global AUROC).
AUROCs are used as global indicators of a dataset performance, a value of 0. 5 indicating a random attribution of labels in the network and a value of 1 indicating a perfect match with the reference dataset. AUROC > 0. 6 may be considered as moderate  19. In global Each of these matrices was thresholded to obtain global networks at different confidence thresholds. Global networks were evaluated and also queried with specific guide gene sets reflecting 5 different pathways in a process named Pathway Level Correlation (PLC). The resulting subnetworks were evaluated and used to construct co-ocurrence networks between microarray and RNA-seq datasets. In white are indicated the figures corresponding to the different steps analyzed. Dataset × distance combinations are indicated in blue and characteristics that are improved by these combinations.
In addition to these 3 evaluations (FPR vs TPR, global AUROC and NV AUROC), GO term significant enrichment was statistically tested with a hypergeometric distribution (not shown in this example). Figure  3  displays TPM network evaluation at different confidence thresholds and Fig. 4  shows networks having 1 million of edges across all dataset × distance combinations. Metrics for all other dataset × distance measurement combinations are presented in Supplementary Fig. 1. All networks combined, pairwise correlations between enriched GO counts, global and NV AUROC performance metrics were moderate (Spearman's rho > 0. 4) but significant (p < 0. 001) indicating these three performance metrics evaluated networks in different ways. The highest correlation was observed between NV AUROC and enriched GO counts (rho = 0.
The global AUROC displayed a very low variation (min = 0. 55, average = 0. 61, max = 0. 68) and was significantly correlated to vertex number (rho = 0. 43, p < 0. 001) only. This observation suggests that the global AUROC is not an appropriated measure in our case. At equivalent edge numbers, different distance measurements generated networks varying considerably in vertex number (Fig. 4A , Supplementary Fig. 1 ). Considering all datasets and distance measurements, raw PCC, raw SCC and MI resulted on average in fewer vertices and higher node degree (vertex number/node degree: 13,164/511, 9,986/465 and 14,074/468 respectively) than PCC-HRR, SCC-HRR or PC (26,645/116, 24,731/124 and 23,966/166 respectively). This trend was clearly observed when setting an edge number to 1 million (Fig. 4A ).
Expression networks constructed from microarrays, TPM, TPM log2, and counts displayed very similar ROC curves: PC based networks followed random predictions (NV AUROC = 0. 5) and the other distance measurements were above the random prediction with similar AUC (Supplementary Fig. 1 ). This was confirmed for PC by NV AUROC and enriched GO term counts. Performance of the other distance measurements in the global TPR/FPR curves did not exactly match that measured with AUROCs. Taking the TPM dataset as an illustration (Fig. 3 ), the MI ROC curve was above the others while NV AUROC for similar edge numbers was slightly below that measured for SCC. This was probably due to differences in network topologies (see above) and the procedures underpinning the two evaluations. The global TPR/FPR curve does not measure a network predictability per se as NV AUROC does and considering any gene pair sharing a same GO term as valid could have overestimated TP (Fig.
Area 1 corresponds to RNA-seq networks having few genes in common with PC networks and microarrays networks and area 2 to combinations maximizing edge cooccurrence between microarray and RNA-seq. between NV AUROC and degree AUROC (rho = 0. 47, p < 2e-16) indicating that highly predictive networks (NV AUROC > 0. 7) also had a higher node degree AUROC. Node degree AUROC was generally under 0. 55. We therefore considered that in our conditions, this bias was only limited. Concerning edge co-occurrence between the different dataset × distance combinations, the lowest conservation was observed with raw (MI, PCC and SCC) RNA-seq datasets and PC networks and microarrays networks (Fig. 4B , area 1). Interestingly, many edges were conserved between VST or sample scaled normalized datasets when correlations were calculated with raw PCC (4B; below block 1).
5A ). For PLC, we systematically removed all genes showing a degree value of 1 (i. e. , those connected to only one guide gene). However we included edges between associated genes if they were found among edges retained at the selected threshold. Using five pathways (Table  1 , Fig. 5B , Supplementary Table  2  and Supplementary Fig. 3 ), we extracted five PLC from the global networks generated above to determine the best suitable dataset × distance measurement combinations. All pathways have modular structures with gene sets forming specific sub-pathways (also called partitions or modules). We expected that PLC would be able to reconstruct such a partitioning, by connecting guide genes with associated genes. Pathways investigated here are typically reaction chains involving several enzymes working in coordination and having no direct effect on the transcription of their gene mates.
As a consequence, specific sub-pathways are expected to form coordinated and specific correlation units highlighted within specific experimental conditions  21. Concerning the signaling pathway, some of the query gene products could directly bind the promoter of other guide genes to control their expression  22. Its study may clearly require a regulatory network approach but we investigated how co-expression network performed with such a signaling pathway. The phenylpropanoid pathway contains a core module composed of 3 genes leading to a precursor used by 3 other distinct subpathways  [23] [24] [25] [26]  (Fig. 5B , Table  1 , Supplementary Table  2 ). The three other metabolic pathways, carbohydrates, fatty acids and terpenoids, were structured in modules as described on the KEGG database  27  (Table  1 , Supplementary Table  2 ).
The fatty acid pathway contains 97 genes divided into 6 modules. The central carbohydrate metabolism contains 202 genes partitioned in 8 modules. Finally, the terpene pathway has 64 genes partitioned into 6 modules. Pathway organizations were used as indicated in the KEGG database (apart from phenylpropanoid pathway which was manually curated from our previous work) and compared to PLC subnetworks. The plant cytokinin (CK) pathway is known to regulate many processes in plant physiology and is hierarchically organized in three levels: a histidine kinase receptor, a transducer (histidine phosphotransfer proteins) and a response regulator (type A/B/C) which may act as a transcription factor  28  (Table  1 , Supplementary Table  2 ). Although CK pathway members are relatively well known, each level is represented by several members which may have specific roles and it is still unclear how they biologically interact with each other to drive a specific physiological response.
Each point in the lower panels (scatterplots) represent one network for which 2 characteristics (eg NV AUROC and modularity) are compared. Data are presented for each pathway separately with a specific color. (C) The expected partitioning of phenylpropanoid related guide genes was compared to two PLC: (D) higher predictability and lower modularity (microarrays raw PCC) and (E) lower predictability and higher modularity (microarrays PCC-HRR). In (D and E), colored vertices correspond to genes encoding enzymes catalyzing steps of similar color in (C). Community (surrounded by grey polygons) numbers in (E) are indicated in deep blue and can be used to access Supplementary Table  3. normalized Chi 2 (rho > 0. 59, p < 0. 001). We found that PLC performance (NV AUROC) was almost negatively correlated with normalized Chi 2 (rho < -0.
Taking the phenylpropanoid pathway as an example, the PCC-HRR based TPM network (Fig. 5E , with a higher modularity) correctly clustered genes from the core phenylpropanoid (PP) and the flavonoid modules while the raw PCC network did not (Fig. 5D , with a higher NV AUROC). Similar results were observed with the four other pathways with either microarray or RNA-seq datasets (Supplementary Fig. 3 ). Modularity and normalized Chi 2 could therefore be considered as consistent quality metrics for PLC. NV AUROC should also be considered to ensure that subgraphs had a minimum predictability (>0. 55).
HRR-CCs optimize recovery and clustering of guide genes in PLC. The best performing dataset × distance measurement combinations were searched by analyzing NV AUROC, modularity and normalized Chi 2 among networks with a Chi 2 p < 0. 05. Statistical effects of dataset, distance, ranking and their interactions on subgraph characteristics were analyzed by ANOVA for each pathway. Ranking and distance measurements had generally the strongest effects on modularity and normalized Chi 2 (p < 2e-5) (Fig. 6A, B ). Ranking had a significant effect on NV AUROC (p < 0. 01) but was weaker than distance measurement (p < 1e-4). Datasets only had a significant effect on modularity (p < 0. 002). Significant interactions were rarely observed between these three factors (i. e. in few pathways and with a weak effect).
This revealed that the different RNA-seq normalizations had only minor effects on these PLCs. Taken as a whole, networks obtained with raw datasets had a significant higher NV AUROC (t-test, mean in raw = 0. 58, mean in HRR = 0. 57, p < 0. 01) but significant lower modularity (mean in raw = 0. 35, mean in HRR = 0. 68, p < 2. 2e-16) and lower normalized Chi 2 value (mean in raw = 0. 20, mean in HRR = 0. 35, p < 2. 2e-16) (Fig. 6A ). It therefore appeared that clustering guide genes correctly was improved with CC ranked with HRR at the expense of performance. NV AUROCs in HRR-based networks were generally higher than 0. 55, indicating an average low performance in GO capture (Fig. 6A ). In non-ranked distances, PC resulted in the weakest NV AUROC, while MI and raw SCC based networks displayed the highest NV AUROC (Fig.
6B ). This weakness in PC based networks was compensated neither by a higher modularity nor by a higher normalized Chi 2 statistic. A more detailed examination of best PLC subgraphs maximizing either modularity or NV AUROC, revealed that each of the five pathways involved specific dataset × distance measurement combinations. PCC-HRR based networks were always found to maximize modularity (Fig. 6C ) and normalized Chi 2 (Fig. 6D ) with almost all datasets. Raw distance based PLCs had a higher NV AUROC and some of them also had a good modularity but they also had a lower normalized Chi 2 statistic indicating they contained fewer guide genes (e. g. raw RNA-seq counts with raw SCC in the terpene PLC). The results suggest that PCC-HRR could be used as a reliable distance measurement whatever the dataset.
Vertex and edge co-occurence in microarray and RNA-seq based PLC subgraphs. Edge co-occurrence in networks constructed from expression datasets obtained by different technologies may be considered as a further validation. Quantifying gene expression with microarrays relies on probe hybridization by sequence complementary while with RNA-seq, short reads are mapped back in silico to the reference transcriptome. The two main differences Table  1. Pathway description. ? Indicates that partition in sub-pathway is not known. between these technologies are (i) the number of quantified transcripts (due to the completion of genome annotation) and (ii) the dynamic range (fluorescent probe intensities for microarrays, in silico read counts for RNA-seq). Because microarrays and RNA-seq technologies differ, edges co-occurring in networks obtained from these two technologies are probably more relevant.
In Fig. 4C , we analyzed co-occurrence in global networks and found that HRR ranked CCs apparently increased the number of co-occurring edges between microarrays and RNA-seq. To get more insights into co-occurrence in PLCs, common edges and vertices were counted in pairwise intersections of networks (RNA-seq vs microarrays) obtained with the six distance measurements and set at a 1,000 vertices. The resulting intersection networks were further characterized by the number of represented guide genes, their normalized Chi 2 statistic, modularity and NV AUROC. This evaluation was performed with the RNA-seq dataset expressed as TPM only because we showed in the previous section that normalization methods had a minor impact on PLC. In addition, TPM networks with raw distance methods had enough vertices to correctly extract PLC (it was not the case with raw distances, e.
6 can be considered as moderate). Different distance measurements displayed different efficiencies according to the dataset but as a general trend, performance of the different combinations were similar (e. g. between microarrays and RNA-seq TPM in Fig. 3B ). The same performance was obtained for different topologies: high node degree (more edges and fewer vertices) for MI and raw CC networks vs lower node degree (fewer edges and more vertices) for CC-HRR networks. PC networks displayed a high performance with microarray data only, suggesting that PCs calculated with 'corpcor' R package may not be recommended for RNA-seq data. A recent study has focused on metabolic pathways in plants using mutual ranks, another CC ranking method  16. Complementary to this previous work, we found that ranking CCs increases vertex number without penalizing absolute network performance.
Contrastingly, an opposite trend was observed in another study  31  , where larger networks displayed a lower Matthew Coefficient when compared to protein-protein interactions or regulatory networks. This indicates that different absolute performance measurements lead to different results and interpretations but this might also be due to our datasets which were larger than theirs. Another advantage of CC-HRR was that it clearly homogenized network characteristics from differently normalized RNA-seq datasets in addition to increase the number of co-occurring edges between microarrays and RNA-seq. The improved performance of HRR may be explained by its ability to integrate the whole transcriptional landscape of each gene thanks to the reciprocal rank calculation. In addition to normalize biases observed in raw CC values when using different data normalization procedures  15  , it probably added robustness to correlations because a low HRR value for a gene pair is observed only if the two genes had no higher correlations with other genes and this probably correct for noise in the expression dataset.
In addition, we saw clear associations between ETR1 and AHK3 in individual PLC subgraphs. Such association highlights crosstalk already known between CK and ethylene signaling pathways  29. The co-occurrence pathway was relatively sparse in contrast to the metabolic pathways (Supplementary Fig. 5 ). It is possible that vertex number for this analysis (1,000) might have been too small to capture complex associations within this signaling pathway. Using VST normalized TPM increased edge and vertex number in the co-occurrence network (Supplementary Fig. 6  and Table  5 ). The above-described associations were also found in this co-occurrence network. While effective in revealing strong gene associations, merging PLC from microarray and RNA-seq data could miss other relevant associations. First, experimental conditions represented by each starting dataset are not completely overlapping.
Together with inherent differences due to dynamic range, this leads to networks with very different edge compositions and node degrees  19  , explaining the relative weak overlap between networks. Second, RNA-seq expression data include genes that are not included in the GPL198 microarray. As an example, some important genes in aliphatic glucosinolate biosynthesis were not represented in a previous Arabidopsis microarray dataset but found in RNA-seq expression matrices from other related species  16. To capture transcriptional environment of a query gene list, distance calculations have to be performed on the whole transcriptome. Calculating partial correlations was particularly challenging but using a covariance shrinkage estimator worked well in terms of computing performance. It took less than 2 h for RNA-seq expression matrices but more than 12 h for the microarray dataset.