dataesr/scientific-paragraphs-categorization huggingface dataset:
, 2021). Physically based models can represent pesticide dynamics at an hourly or daily time step  (Holvoet et al. , 2007; Gevaert et al. , 2008; Boithias et al. , 2014; Ammann et al. , 2020). Current models can represent various processes influencing pesticide transport and fate in soils and surface waters  (Holvoet et al. , 2007; Ippolito and Fait, 2019; Wang et al. , 2019). However, these models need continuous input and detailed landscape data and are time-consuming to set up. Moreover, they are usually applied at the scale of headwater catchments, and they are challenging to up-scale to larger areas  (Gassmann et al. , 2014). Another method to predict pesticide risks relies on using indicators, combining multiple large-scale datasets related to pesticide application and transport, such as slope, land use, soil types, or distance to the stream.
Database of monitoring sites The NAWA TREND program consists of multiple monitoring sites draining watersheds of different sizes and various anthropogenic pressures. We excluded large watersheds with stream orders (the number representing the branching level in a river system) greater than five for most analyses, as they make up only a small part of the Swiss river network (Appendix 1). In this study, we used the Strahler stream order attributing the first order to each headwater reach  (Strahler, 1957 (Strahler, , 1964)). We also excluded sites not operational before May 2019. The final set contains 21 independent, non-nested stations across the Swiss Plateau and the Rhône valley (Fig. 1 ). These sites represent various hydrological and agricultural characteristics with an average elevation ranging from 415 to 683 m.
Thus, 23 % of the total mass is applied to corn, giving f load,corn,S À Metolachlor,4 = 0. 24. For monitoring site 10, there are no beets or other uses, resulting in f load,corn,S À Metolachlor,10 = 1, as the total S-Metolachlor found in-stream should come exclusively from corn fields. Relation with corn fields. We used land use data and information about average application rates to estimate the average concentration from application to corn: C i,corn,k,av ¼ f load,corn,i,k :C i,k,av Assuming spatially uniform application rates for any given herbicide i on a given crop across catchments, we may also expect the cumulative concentration to be proportional to the areal fraction of corn: C i,corn,k,av $ f corn,k As for single compounds, we tested the hypothesis that the sum of the normalised, corn-specific cumulative herbicide concentrations varied linearly with the fraction of corn.
As an alternative hypothesis, we combined f corn,k with the connectivity map (see Table  1 ) to obtain the fraction of corn connected to the stream in each catchment (f corn connected,k ), which has been compared to ∑ i C i,corn,k,av. Nevertheless, this comparison was performed in 18 sites instead of 21 as the high-resolution land use data was not available for sites 7, 8, and 9 (Appendix 2): ∑ i C i,corn,k,norm $ f corn,k or f corn connected,k 2. 3. 1. 2. Herbicides selection. To calculate average concentrations, we had to consider which corn herbicides were included in the NAWA TREND program. Corn herbicides were characterised by three factors: i) the fraction of corn fields treated with these herbicides ( area of corn fields receiving the herbicide total area of corn fields ; x-axis in Fig.
Approach 2: Concentrations of multiple pesticides from multiple crops The NAWA TREND program monitors 48 mandatory pesticides and ten optional ones. As of 2021, the Waters Protection Ordinance (Gewässerschutzverordnung; GSchV) defines Environmental Quality Standards (EQS; Appendix 6) for 17 of the mandatory analytes, selected to cover a significant part of the ecotoxicological risk, which we selected in our study (Appendix 6). We excluded two substances (Diazinon and Terbutryn) from the analysis because they are no longer registered as Plant Protection Products (PPP) in Switzerland. The final list comprises three fungicides, eight herbicides, and six insecticides. We computed average expected concentrations at the NAWA TREND sites for the pesticides listed in Appendix 6. We calculated an average concentration C k for each NAWA TREND site during the growing season, i.
These findings demonstrate that the resulting herbicide patterns in the streams vary substantially, even for a crop with a simple herbicide application pattern (one application per year). Based on the available data and information, predicting the respective herbicide concentrations seems challenging. There are many possible reasons for this high uncertainty. Firstly, only a fraction of all corn herbicides is measured within the NAWA TREND program (Fig. 2 ). If the application of corn herbicides was spatially uniform, this fact would not be problematic. However, comparing corn herbicides across the monitoring sites (Appendix 7) suggests a very heterogeneous application pattern in space, especially for Foramsulfuron and Dimethenamid-P. Hence, one major limitation for the successful prediction of corn herbicides is the poor knowledge of the actual input in each catchment.
Methodology This paper explores relationships to predict in-stream pesticide concentrations with explanatory variables. We tried various approaches with a large dataset concerning pesticide concentrations, i. e. , high frequency and numerous stations, and detailed explanatory variables (e. g. , detailed land use data at the parcel scale). Our analyses indicate that moderate relationships between aquatic pesticide pollution levels and land use data only appear with a combination of high-resolution land use and connectivity data sets for specific crops, e. g. , corn fields, in this study. The high uncertainty shown by our spatial extrapolation models can be attributed to many sources, where the limited knowledge about the pesticides inputs plays a major role. Although different models showed high uncertainty when extrapolating in space, interannual variability showed a more consistent pattern.
Missing samples, substances, and catchments The number of samples in each catchment is not always the same due to differences in sampling protocols. As was done in this paper, normalising the number of samples may attenuate the problem. However, lacking samples may induce missing parts of the signal important to reflect the total amount of pesticides exported to the outlet, as concentration patterns are highly dynamic. Therefore, missing samples will impact the average concentrations we calculated for 2019 and may influence all the correlations with the explanatory variables we tried to highlight in this paper. This statement is even more true if missing samples occur during the growing season, which typically have higher concentrations than in the non-growing season. To avoid this issue, future research may focus on multiyear datasets in order to attenuate the impact of missing samples on the final signal.
In Approach 2, we explore multiple compounds from the Water Protection Ordinance of Switzerland, selected to cover a significant part of the ecotoxicological risk. In this paper, we selected substances based on their impact on freshwater ecosystems but also depending on the sampling frequencies in the NAWA TREND program. It was important to choose compounds measured at all the sites selected in our study to cover Switzerland and keep as many watersheds of the NAWA TREND program as possible. However, these selections imply that we miss signals from unselected compounds, as in Approach 1 with the unselected corn herbicides in Fig. 2. These missing substances may be applied in substantial amounts in some watersheds and not others, which may influence the total signal received at the outlet. Therefore, the correlations could be improved by analysing more compounds and representing a broader picture of the pesticides used on the national scale.
This pesticide is mainly used on potatoes (Appendix 6), which could explain the outlier as the fraction of potatoes in site 8 is much higher than in other sites (12 % against 1 % on average). The models provided in this paper depend on the pesticide concentration and the quality of the explanatory datasets. As mentioned in Section 2, land use data are available at a high-resolution in many cantons but only at the municipality scale in others. This difference may impact the spatial extrapolation depending on the land use datasets. For example, the unavailability of high-resolution land use data in watersheds 7, 8, and 9 did not allow calculating a fraction of connected corn fields in these areas. In the same way, soil characteristics are qualitative indices that limit the statistical analysis.
Acknowledgments We acknowledge the  Federal Office for the Environment  and the cantons involved in the NAWA micropollutant monitoring program for sharing their pesticide data sets. We thank the organisations of the cantons for the sharing of their high-resolution crop data: Daten des Kantons Aargau;  Geodaten Kanton/Bezirke Appenzell I. Rh ;  Daten des Kanton Berns ;  Geodaten des Kantons Basel-Landschaft; Système d'information du territoire à Genève (SITG) ;  Kanton Glarus ;  Géodonnées de la République et Canton du Jura ;  Geoinformation Kanton Luzern ;  Kanton Nidwalden ;  Kanton Obwalden ;  St. Gallen ;  Amt für Landwirtschaft , Kanton Schwyz ; Thurgau;  Uri: Lisag AG ;  GIS Kanton Zug.
The majority of gene expression studies focus on the search for genes whose mean expression is different between two or more populations of samples in the so-called "differential expression analysis" approach. However, a difference in variance in gene expression may also be biologically and physiologically relevant. In the classical statistical model used to analyze RNA-sequencing (RNA-seq) data, the dispersion, which defines the variance, is only considered as a parameter to be estimated prior to identifying a difference in mean expression between conditions of interest. Here, we propose to evaluate four recently published methods, which detect differences in both the mean and dispersion in RNA-seq data. We thoroughly investigated the performance of these methods on simulated datasets and characterized parameter settings to reliably detect genes with a differential expression dispersion.
We applied these methods to The Cancer Genome Atlas datasets. Interestingly, among the genes with an increased expression dispersion in tumors and without a change in mean expression, we identified some key cellular functions, most of which were related to catabolism and were overrepresented in most of the analyzed cancers. In particular, our results highlight autophagy, whose role in cancerogenesis is context-dependent, illustrating the potential of the differential dispersion approach to gain new insights into biological processes and to discover new biomarkers.
Therefore, dedicated methods based on discrete probability distributions were developed to analyze these data  [32]. The negative binomial (NB) distribution has become the ubiquitous distribution to model RNAseq read count data by providing the best fit for the extra-variance commonly observed in datasets composed of biological replicates  [33]. In this model, the random variable describing the count of reads mapped to gene i in sample j is denoted as Y ij � N Bðm ij ; � i Þ, where μ ij is the expected value and ϕ i is the dispersion parameter. The variance is given by VarðY ij Þ ¼ m ij þ � i m 2 ij. Analyzing the variance independently with respect to the mean expression can therefore be achieved by analyzing the dispersion parameter ϕ i. In the classical RNA-seq data analysis workflow, differential expression detection methods based on the NB distribution consider the dispersion as a noise parameter to be estimated prior to identifying a difference of mean expression  [34].
The generally low sample sizes of RNA-seq datasets at the time when the first versions of these methods were published made per-gene dispersion estimation unreliable. In addition, the very high number of genes made estimation difficult. Thus, Robinson et al. proposed an accurate shared estimator based on the expression of sets of genes across all samples, independent of biological condition  [33]. Pergene estimators were then shrunk towards this shared estimator using different levels of shrinkage  [35] [36] [37] [38]. Aggregating all the samples that compose the dataset implies that no difference of dispersion in the expression of genes between the conditions of interest can be modeled, which is not biologically realistic.
Objectives The performances of methods identifying differences in mean expression in the so-called "differential expression analysis" using RNA-seq data have been extensively studied  [43] [44] [45] [46]. The large amount of publicly available RNA-seq data opens new perspectives for researchers in the search for genes whose expression exhibits a difference of dispersion between samples from different conditions and a new space for the discovery of biomarkers. To that end we first evaluate the performances of four NB-based methods, MDSeq, DiPhiSeq, GAMLSS and DiffDist, and a statistical test of difference of variances, Levene's test  [47] , to identify differentially dispersed (DD) genes using simulated RNA-seq datasets. Based on our simulation study results, we reliably applied these methods to The Cancer Genome Atlas datasets and identified DD genes that could not be identified by the classical differential expression analysis.
Evaluation of differential dispersion detection performances Differential dispersion detection for genes with unconstrained differences in mean expression. We simulated RNA-seq datasets to evaluate the performances of MDSeq, DiPhi-Seq, GAMLSS and DiffDist to identify differential dispersion between two sets of samples of equal size that represent two conditions of interest. We also included Levene's test to identify differences of variances after log 2 -transformation of the data. Differences in the mean and dispersion between the two sets of samples were introduced and defined for DE and DD genes, respectively (see the Methods section for more details). Overall, GAMLSS and DiffDist outperform all the other evaluated methods for all samples sizes (green boxplots in  Fig 1). It is noteworthy that Levene's test performs only a bit worse than DiPhiSeq and better than MDSeq.
Given the results in The trends observed with datasets with high mean fold change values are confirmed by the results obtained with simulated datasets with moderated mean fold change values. GAMLSS and DiffDist outperform again all the other evaluated methods. They are indeed much sensitive to detect DD genes than all the other evalutated methods while controlling the FDR, in particular for low sample sizes. DiPhiSeq is the more robust method with very low FDR values for all the evaluated sample sizes while being also the less sensitive method, in particular for low sample sizes. DiPhiSeq is indeed unable to detect any DD genes with fewer than 40 samples per population. However, increasing the sample size to 100 samples per population enabled DiPhi-Seq to reach a sensitivity similar to the one observed with Levene's test and close to those of GAMLSS and DiffDist.
Regarding MDSeq, increasing the maximum tolerated mean fold change value increases the FDR for the detection of differential dispersion. However, the FDR remained below 0. 05 for datasets composed of 30 to 50 samples per population with maximum tolerated mean fold changes up to 1. 5. When only 20 samples per population are available, the maximum tolerated mean fold change value must be at most 1. 3 to keep the FDR below 0. 05. Similarly to simulated datasets with high values of mean fold change, we evaluated the Benjamini-Hochberg procedure and observed that it was again unable to maintain the FDR below 0. 05 for almost all datasets, reaching values above 0. 1 for large mean fold change values (S2 Fig).
Nevertheless, we observed some differences with these simulated datasets. Contrary to datasets with high mean fold change values, the performances of Levene's test, GAMLSS and DiffDist are not affected by the presence of the low mean fold change values. Besides, Levene's test appears to be even more sensitive than MDSeq and DiPhiSeq while maintaining the FDR below 0. 05 for most replicates. The limitation in the application of MDSeq to lowly DE genes is not prohibitive since the purpose of our approach is to identify genes that would not be detected by the classical differential expression analysis or, at least, that would not appear in the top results of these analyses. Thus, in our approach, these genes represent the set of genes of primary interest among which to search for differential dispersion in expression.
Differential dispersion in gene expression in cancer We applied the five evaluated methods to The Cancer Genome Atlas (TCGA) datasets  [49]  to identify DD genes when comparing normal and tumor samples. We used RNA-seq data from patients for whom tumor tissue and adjacent normal tissue samples were available, to limit individual variability. In agreement with the results of our simulation study, only the datasets with more than 30 samples for both conditions were analyzed, in order to maintain the FDR below 0. 05 with MDSeq and to ensure sufficient power with DiPhiSeq. We list these datasets in Table  1. Identification of DD genes among non-DE genes. Taking benefit of the GLM implemented in MDSeq, we identified DE genes while taking into account batch effects and identified DD genes among non-DE genes using Levene's test, MDSeq, DiPhiSeq, GAMLSS and DiffDist.
We used a fold change threshold of 1 to identify both DE genes and DD genes among non-DE genes (Fig  4 ). The identification of DE genes by MDSeq leaves several thousand genes among which DD genes can be searched for in any dataset. Among non-DE genes, the majority of DD genes are overdispersed in tumors (DD+). Overall, all the methods generate consistent results: some cancers are characterized by a high number of DD+ genes (breast, colon, kidney, liver and lung), while others contain much less DD genes (head and neck, prostate and thyroid). Consistent with our simulation study, most of the DD+ genes are identified by all the five evaluated methods and GAMLSS and DiffDist are the methods that specifically identifies the most DD+ genes for most datasets, like in the kidney renal clear cell carcinoma dataset (TCGA-KIRC, Fig  5A).
Hundreds of DD+ genes are indeed only detected by these two methods (S2 File). These sets of specific genes may represent a very significant part in the overall set of DD+ genes detected by all the methods like for the prostate cancer (S2 File), highlighting the breakthrough of GAMLSS and DiffDist among the methods identifying differential dispersion in RNA-seq data. Nevertheless, the two other differential dispersion methods specifically identify significant sets of genes for some datasets, such as the kidney renal papillary cell carcinoma (TCGA-KIRP, Fig  5B ) and prostate adenocarcinoma datasets (TCGA-PRAD, S2 File) for MDSeq, the thyroid carcinoma dataset (TCGA-THCA) for DiPhiSeq (S2 File) and, to a lesser extent, the head and neck squamous cell carcinoma dataset (TCGA-HNSC) for both (S2 File).
These sets of genes are all the more important that they deal with the three datasets for which the less DD + genes are detected over all the analyzed datasets. Finally, Levene's test is the method which specifically identifies the less DD+ for most datasets. Gene Ontology term enrichment analysis. To gain biological insight, we conducted an analysis of enrichment in Gene Ontology (GO) terms among the previously identified DD + genes for each TCGA dataset. Since DD+ genes are identified with an FDR below 0. 05 for all methods according to our simulation study, the entire set of DD+ genes identified by at least one method is taken into account to gain the most biological knowledge in the GO term enrichment analysis for each dataset. As we showed in our simulation study, GAMLSS and DiffDist may incorrectly estimate the sign of dispersion log 2 -fold change.
Since this type of error is detrimental for the subsequent GO term enrichment analysis, we only retained the DD genes identified by these two methods whose dispersion log 2 -fold change sign is in agreement with the one predicted by DiPhiSeq (S4 Fig). We used redundancy reduction methods to ease the comparison of enriched GO terms across all the analyzed datasets (see Methods for more details). The top  40    Interestingly, among DD+ genes, the most significantly enriched GO terms were the most widespread across all the analyzed tissues and focused on some key cellular functions, such as catabolism. In contrast, GO terms that were found to be significantly enriched for only a few datasets tended to have higher p-values than the most widely enriched GO terms (S3 File). This striking result suggests some common features in tumoral development and progression, regardless of the tissue of origin, whose involved gene expression is characterized more by an increase in dispersion than by a change in the mean in tumors.
It is worth noting that many of these common biological processes are related to catabolism, e. g. "mRNA catabolic process", "protein targeting" or "proteasomal protein catabolic process", as previously shown by Han et al. [27]. In particular, several processes related to the ubiquitinproteasome system, which is a major controller of the protein degradation process and is highly involved in cancer  [50] , are found among the most significant results ("protein polyubiquitination", "proteasomal protein catabolic process"). In contrast, no process related to anabolism was found among the most frequently enriched processes among DD+ genes, suggesting that catabolic processes are much more affected by the dysregulation of gene expression than are anabolic processes. Autophagy was also found as the most widespread and significantly enriched biological process among DD+ genes for all the analyzed datasets.
Similar to the proteasome, it is a main recycling system for biological molecules that enables cells to survive critical situations such as nutrient starvation and the degradation of damaged organelles or pathogens. In pre-malignant cells, autophagy actively acts to preserve the physiological homeostasis of multiple functions, e. g. elimination of mutagenic entities, decrease local inflammation, and thus aid the struggle against tumor development. In malignant cells, autophagy affects the tumor progession and the response to treatment in multiple ways, some of which act in opposition. Autophagy Enriched GO terms among overdispersed genes in tumors identified by the evaluated methods. Top 40 representative enriched Gene Ontology (GO) terms among overdispersed genes in tumors (DD+) among non-differentially expressed (non-DE) genes, ordered first by the number of datasets for which they are enriched (decreasing order) and second by the mean p-values of enrichment across all datasets (increasing order).
Gene expression dispersion in cancer Overall increase of dispersion and robustness. By applying the five evaluated methods to TCGA datasets, we identified an overall increase in the dispersion in the expression of many non-DE genes in tumors in comparison with normal tissues. Also analyzing TCGA datasets, Han et al. have already revealed an increase in the coefficient of variation of gene expression in tumors of breast, colon, lung and liver cancers  [27]. In addition, using microarray data, Ho et al. [8]  also noticed that an increase in gene expression variance in a disease condition such as cancer is more common than a decrease. Our work confirms these results, extends them to other cancers and increases their reliability by using RNA-seq data and methods based on a more appropriate statistical framework, and rigorously validates them in a simulation study.