dataesr/scientific-paragraphs-categorization huggingface dataset:
Putative chimeras were detected by UCHIME  [277]. The combining of these cleaning procedures with a read length (2 × 250 bp) covering the entire amplicon length (~220 bp without the primers) allowed the sequencing of each base twice and the drastic minimization of the sequencing errors. After this cleaning step and demultiplexing process, the remaining DNA sequences were clustered at a 95% similarity threshold with UPARSE 7. 0  [276]  to obtain the seed MOTUs (option-cluster_fast). The taxonomic affiliation was performed by BLAST against the SSURef SILVA database  [278]  after application of the following selection criteria: length >1200 bp, quality score >75% and a pintail value >50; the SILVA database was enriched by lacustrine DNA sequences originating from various studies on lacustrine systems  [279] [280] [281] [282].
The taxonomy of a MOTU corresponded to the best hit given by the similarity search. If an OTU was associated with several best hits (hits with the same identity), then the taxonomy was the common taxonomy of these hits named also the lowest common ancestor. MOTUs affiliated to embryophytes, cercozoa nuclear and unclassified Eukaryota were removed from these analyzes in order to consider only unicellular eukaryotes, fungi and metazoa in this analysis. The two molecular inventories obtained from one SDH DNA extract obtained using the protocol PB-0. 75g showed low amounts of read abundances and were thus removed for further analysis of the dataset. For each sample, the molecular inventories were merged and MOTUs detected in only one of the two duplicates were discarded. In addition, all MOTUs with a total abundance lower than 3 DNA read sequences were removed.
We applied negative controls during both DNA extraction and PCR setup. After PCR, we pooled and purified PCR products into two pools following Clarke et al. [26]. One pool (JIE4) was shipped to FASTERIS, SA (Switzerland) and converted into a single-indexed PCRfree MetaFast library (following  [26] ). The second pool (AOHL-3-8) was converted into a unique dual-indexed library in Tromsø following the Illumina TruSeq DNA PCR-free protocol (Illumina, Inc. , CA, USA), with the bead-cleanup steps modified to retain our short amplicon inserts. Both libraries were sequenced on ~10% of separate 2x 150 cycle mid-output flow cells on the Illumina NextSeq platform at FASTERIS. Bioinformatics and data analysis-We followed the bioinformatics pipeline presented by  Rijal et al. (2020) , which uses the ObiTools software package  [261]  and custom R scripts (available at https://github.
com/Y-Lammers/MergeAndFilter (accessed on 1 February 2021)). Briefly, we merged overlapping paired-end reads and retained only those that merged. We then demultiplexed our tagged amplicons using the tag-PCR replicate lookup presented in. We then collapsed identical sequences, removed putative artifactual sequences from our data that may have derived from Illumina library index-swaps or PCR/sequencing errors, and removed sequences that had fewer than 3 reads in all PCR replicates. We next identified amplicon sequence variants (ASVs) that had 100% identity agreement with either the ArctBorBryo and/or EMBL nucleotide (rl133 release) databases, following Volstad et al. [49]. We further removed identified ASVs that 100% matched two "blacklists" consisting of known contaminants/exotics or synthetic sequences  [16].
The final taxonomic assignment of the 71 retained ASVs was determined by one of us (Alsos) using regional botanical taxonomic expertise. We calculated the proportion of weighted replicates (wtRep) for each retained ASV in each sample using the method presented by Rijal et al. [16] , which weights PCR replicate detections based on relative sequencing depth. For each of the six taxa that were represented by more than one ASV, we selected the ASV that had the greatest dataset-combined wtRep. This resulted in a final data set of 60 taxa that were used for all taxonomic richness and PCR replicability statistics. For visualization, we selected the top 40 taxa that had the greatest dataset-combined wtRep. We visualized the results in R v3. 4. 4 (R Core Team 2018) using the rioja v0. 9-21 package  [299].
The relationship between different levels of integration is a key feature for understanding the genotype-phenotype map. Here, we describe a novel method of integrated data analysis that incorporates protein abundance data into constraint-based modeling to elucidate the biological mechanisms underlying phenotypic variation. Specifically, we studied yeast genetic diversity at three levels of phenotypic complexity in a population of yeast obtained by pairwise crosses of eleven strains belonging to two species, Saccharomyces cerevisiae and S. uvarum. The data included protein abundances, integrated traits (life-history/fermentation) and computational estimates of metabolic fluxes. Results highlighted that the negative correlation between production traits such as population carrying capacity (K) and traits associated with growth and fermentation rates (J max ) is explained by a differential usage of energy production pathways: a high K was associated with high TCA fluxes, while a high J max was associated with high glycolytic fluxes.
However, previous studies on the HeterosYeast data set failed to find a clear link between variation at the trait level and variation at the proteomic level. Given the important genomic resources in yeast  [23] , a number of curated genomescale metabolic models are now available  [24]. Among these, the DynamoYeast model  [25]  describes central carbon metabolism in yeast. It is small enough (70 reactions and 60 metabolites) to remain tractable, and has been tested against experimental data  [26]. The availability of the HeterosYeast dataset, combined with a curated metabolic model of central carbon metabolism in yeast and a probabilistic approach to explore the solution space, encouraged us to integrate experimental proteomic data into the metabolic model in order to predict unobserved metabolic fluxes.
Results The HeterosYeast dataset provided valuable observations on the genetic diversity of yeast strains involved in the wine-making process at different levels of cellular organization, i. e. phenotypic traits related to life-history or fermentation  [21] , and quantitative proteomic data  [20]. All traits were estimated or measured at 18˚C and 26˚C on a half-diallel design comprising 7 strains of S. cerevisiae and 4 strains of S. uvarum, with a total of 127 strain × temperature combinations. In order to access an intermediate level of integration between protein abundances and traits, we used the DynamoYeast model, a curated Constraint-Based Model (CBM) of central carbon metabolism in yeast consisting of 70 reactions and 60 metabolites  [25]  (Fig  1 ). In the DynamoYeast model, the only entry is glucose, and the model does not take into account the complexity of metabolism like the recycling of amino-acids or extra-cellular TCA supplements.
However, it was shown to accurately predict growth on complex medium like grape must  [26]. Using gene-protein-reaction associations, enzymatic proteins and protein complexes were linked to the reactions of the CBM. Among the 70 enzymatic proteins and protein complexes, the abundances of 33 of them were retrieved from the dataset of 615 protein abundances quantified in the HeterosYeast project. Thus, the metabolic fluxes that best matched the observed patterns of variation of enzymatic protein abundance were retained (see Material and methods). In brief, the strategy we proposed was: (i) to characterize the feasible solution space of the DynamoYeast model, L, through the posterior density distribution of the fluxes given by the Expectation Propagation algorithm (hereafter denoted EP,  [17] ); (ii) to select a unique solution through minimization of the objective function that measures the Euclidean distance between observed enzyme abundances and reaction rates (Z, Eq 11 in Material and methods).
e. the number of degrees of freedom of the DynamoYeast metabolic model. When considering the actual set of enzyme abundances that were matched between the HeterosYeast proteomic data and the DynamoYeast CBM (N obs = 33), we observed a high correlation between v initial and v predicted. This correlation increased when increasing the number of sampled points (see for example Fig 2C  and 2D ). We also checked that the 33 reactions were distributed over the modules of the metabolic reaction tree  [28]  (see Testing the CBM coverage with the observed proteomic dataset), as shown in S4 Fig. Altogether , we considered that our algorithm was efficient for predicting unobserved fluxes from enzyme abundances, given the correlation values obtained through simulations. More generally, simulations provide a method to check for the good coverage of the metabolic model with observations, without needing to predict elementary modes, which becomes computationally heavy as the size of the metabolic model increases.
Predicting unobserved fluxes from the observed variation in protein abundances The HeterosYeast proteomic data were used in the context of the DynamoYeast model of central carbon metabolism in yeast. In addition, for each strain × temperature combination, the observed CO 2 release rate was used as an additional constraint in the form of a priori knowledge to delimit the feasible solution space L. We sampled N s = 10 6 points in L to select a unique solution that minimizes the Euclidean distance between fluxes and enzyme abundances. As a result, we predicted 69 unobserved fluxes in the CBM for each of the 127 strain × temperature combinations. Statistical approaches were then used to investigate the variation components and the structure of the new dataset D consisting of 615 protein abundances (E), 70 metabolic fluxes (V) and 28 fermentation and life-history traits (T): D ¼ ðE; V; TÞ.
Patterns of variation depend on the integration level The 127 observations in the new dataset D had a specific structure. There were 7 parental strains (S. c. ) and 21 intraspecific hybrids from S. cerevisiae (S. c. × S. c. ), 4 parental strains (S. u. ) and 6 intraspecific hybrids from S. uvarum (S. u. × S. u. ), and 28 interspecific hybrids (S. c. × S. u. ). All strains were observed during alcoholic fermentation of wine grape juice at two temperatures, 18˚C and 26˚C  [21]. To better understand the patterns of variation at each integration level, Principal Component Analyses (PCA) were performed for each type of trait separately. Results are presented in Fig  3 , where strains are identified by species, type of cross (intraspecific hybrid, interspecific hybrid or parental strain) and temperature.
Discussion We applied cutting-edge methods of data integration to an original yeast dataset. The Hetero-sYeast dataset comprised quantitative proteomic data as well as fermentation and life-history traits measured during wine fermentation on a range of strains from two yeast species. The objective was to integrate information at different levels of cellular organization (proteomic and metabolic fluxes) to better understand the metabolic bases of phenotypic variation in yeast, in particular life-history traits related to fitness. The key point of this study was to incorporate proteomic data in a constraint-based metabolic model to estimate the values of unobserved metabolic fluxes. Using a combination of multivariate analyses of heterogeneous highdimensional datasets, we were able to show that the metabolic flux level retains information that is not directly interpretable at the proteomic or trait level.
Despite the fact that the DynamoYeast metabolic model is an oversimplified model of central carbon metabolism with glucose as the only external carbon source, we show that protein abundance variations were sufficient to capture quantitative changes in the orientation of central carbon metabolism that occurred between strains and between growing temperatures in our dataset. Even though our flux predictions may not be very accurate, we are confident that we captured the main patterns of flux variation. Predicting unobserved fluxes from observed protein abundances overall adds information about the functioning of the actual metabolic network. The method described here relies on a probabilistic approach. Following  [17] , we chose to characterize the feasible solution space L by means of its posterior density distribution calculated with the Expectation Propagation (EP) algorithm.
The final dataset consisted of three matrices of 127 × 615 protein abundances, 127 × 70 central carbon fluxes, and 127 × 28 fermentation/life-history traits. As the total number of characters (713) greatly exceeded the number of samples (127), we used regularization techniques for the multivariate analyses  [33]. In order to relate the variation patterns we observed at different integration levels, we used a top-down strategy from the most integrated to the least integrated level. First, we explored the correlation between traits and metabolic fluxes. Second, we identified the proteins which were not included in the metabolic model that best explained the correlation between traits and fluxes. We showed that the integration of the flux level allowed us to better understand the patterns of variation at the trait level.
In our dataset, we found a negative correlation between traits associated with growth and CO 2 fluxes, and traits associated with population size and fermentation length. These negative correlations reflected different life-history strategies, as has been observed previously in different yeast collections from either industrial  [34]  or natural sources  [35, 36]. This broadly corresponds to the well-known r-K trade-off in ecology  [37]. More recently,  [38]  suggested that such a trade-off could arise from eco-evolutionary feedback loops because competing strains also modify their environment through the production of different sets of metabolites. The HeterosYeast dataset shows that the choice of a strategy is plastic [21] and can be modified by the environment (here the fermentation temperature).
Materials Materials: The HeterosYeast dataset. The genetic material for the experimental study consisted of 7 strains of S. cerevisiae and 4 strains of S. uvarum associated with various food processes (enology, brewing, cider fermentation and distilling) or isolated from the natural environment (oak exudates). The 11 parental lines were selfed and pairwise crossed, which resulted in a half-diallel design with a total of 66 strains: 11 inbred lines, 27 intraspecific hybrids (21 for S. cerevisiae, noted S. c. × S. c. , and 6 for S. uvarum, noted S. u. × S. u. ) and 28 inter-specific (noted S. c. × S. u). The 66 strains were grown in triplicate in fermentors at two temperatures, 26˚C and 18˚, in a medium similar to enological conditions (Sauvignon blanc grape juice, [21]). From a total of 396 alcoholic fermentations (66 strains, 2 temperatures, 3 replicates), 31 failed due to the poor fermenting ability of certain strains.
The design was set up as a block of two sets of 27 fermentations (26 plus a control without yeast to check for contamination), one carried out at 26˚C and the other at 18˚. The distribution of the strains in the block design was randomized to minimize residual variance of the estimators of the strain and temperature effects, as described in  [34]. For each alcoholic fermentation, two types of phenotypic traits were measured or estimated from sophisticated data adjustment models: 35 fermentation/life-history traits and 615 protein abundances. Fermentation/life history traits were classified into four categories [21]: • Kinetics parameters, computed from the CO 2 release curve modeled as a Weibull function fitted on CO 2 release quantification monitored by weight loss of bioreactors: the fermentation lag-phase, t-lag (h); the time to reach the inflection point out of the fermentation lagphase, t-V max (h); the fermentation time at which 45 gL -1 and 75 gL -1 of CO 2 was released, out of the fermentation lag-phase, t-45 (h) and t-75 (h) respectively; the time between t-lag and the time at which the CO 2 emission rate became less than, or equal to, 0.
For the fermentation traits, the model did not converge for most ethyl esters (Ethyl-propanoate, Ethylbutanoate, Ethyl-hexanoate, Ethyl-octanoate and Ethyl-decanoate), and for Acetate Ratio and Acetic acid. These traits were removed from the final analysis, which in the end included 28 traits. Protein functional annotation. Cross-referencing MIPS micro-organism protein classification  [40] , KEGG pathway classification  [41] [42] [43]  and Saccharomyces Genome database  [44]  information, we attributed a single functional category to each protein. The first two hierarchical levels of MIPS functional annotation were taken into account to assign proteins to 34 different categories. All secondary levels were used for the 01. metabolism, 02. energy and 10. cell cycle and DNA processing categories, resulting in 20 different functional categories.
Following these rules for each of the 11 yeast strains and the 55 hybrids at both temperatures, we estimated the protein abundances associated with the reactions of the DynamoYeast model, resulting in 33 reaction weightings out of 70. Among the 70 reactions, only six were associated to protein complexes with an AND Boolean relationship, and four of them (in bold in the right column) matched two proteins of the HeterosYeast proteomic dataset (proteins in bold, Table  2 ). We checked that pairwise correlations between proteins involved in the same complex were either positive or null, suggesting the absence of post-translational regulation that would change the overall enzyme concentration.
The algorithm implemented by  [17]  was translated into R code. Extraction of the stoichiometric matrix from the DynamosYeast model was performed with the sybil package in R  [46]. Prediction of metabolic fluxes from proteomic data. In living systems, most metabolic reactions are catalyzed by enzymes, and quantitative proteomic data retain information about enzyme abundances. Therefore, the metabolism of a cell at a given time is characterized by the set of fluxes, metabolites and protein abundances ðv; m; EÞ t where E = (E 1 , E 2 ,. , E N ), and E i is the abundance of enzyme i associated with the reaction flux v i. Indeed, even though reaction rates are not directly proportional to enzyme abundances, a degree of covariation between protein abundance and flux reaction rate can be expected at the scale of the metabolic network.
• The objective function was calculated only on the subset of observed enzyme abundances: Z ¼ 1 p v X N obs i¼1 ðE i À jv i jÞ 2 Predictions of metabolic fluxes were performed by coupling the DynamoYeast model to our experimental data (protein abundances and CO 2 reaction rates, the only measured flux in our study). We constrained the solution space L by considering the maximum CO 2 release rate, measured at the same time point as the one used in the proteomics analyses  [20]. For each strain × temperature combination a unique solution was obtained by minimizing the objective function, defined in Eq 11, from the set of observations. Testing the prediction algorithm. The prediction algorithm is based on the assumption that fluxes and enzyme abundances covary. Indeed, any reaction rate can be expressed as a more or less complex function of enzyme abundances, kinetic constants and metabolite concentrations  [47] : v i ¼ k cat i E i f ðκ; m; EÞ where k cat is the catalytic constant, κ is a set of other kinetic constants, E is the set of enzymes abundances other than enzyme i.
2. This procedure was performed until the complete set of reactions was selected. In total, simulations were run a thousand times for different values of N s and N obs. Testing the CBM coverage with the observed proteomic dataset. We performed a modular decomposition of the DynamoYeast model via the analysis of the null space (or kernel) of model stoichiometric matrix to check that the 33 reactions that matched with observations of protein abundances were distributed over the main pathways. The stoichiometric matrix S in the DynamoYeast model is of dimension 60 × 70 and of rank 54, meaning that 16 among the 70 reactions are strictly coupled. Thus, the dimension of the null space is Ker(S) = 16, roughly meaning that we need to know 16 independent flux combinations to predict changes in metabolic rates.
Each reaction is associated with a 16 dimensional row vector in the null-space 70 × 16 matrix. We ran a hierarchical clustering using as a metric the symmetric matrix of the angles between the reaction vectors in the null-space  [28]  and we found, as previously shown, that the reactions cluster within metabolic pathways (S4 Fig). We found that the 33 reactions that matched with observations of protein abundances were distributed over the main pathways (S4 Fig). In addition, we used our simulation method to check the correlations between initial and predicted fluxes by choosing as being observed the same 33 fluxes as in the experimental dataset.
Statistical analysis In order to study the main features characterizing fermentation and life-history traits in the HeterosYeast dataset, we analyzed the variation components of data from three different levels of cellular organization: protein abundances E, metabolic fluxes V and fermentation/life-history traits T: D ¼ ðE; V; TÞ The total number of observations was 127 strain × temperature combinations (66 strains, 2 temperatures, minus 5 missing data due to the poor fermenting abilities of some strains). The whole dataset consisted of 615 protein abundances, 70 metabolic fluxes and 28 fermentation and life-history traits. Two types of analyses were performed using different multivariate approaches: analyses at a single phenotypic level and analyses integrating the different levels. We ran Principal Component Analyses (PCA) to identify the most important sources of variation within the datasets and the similarities/differences between the different phenotypic levels.
We included prior knowledge regarding the yeast species in order to perform a supervised sparse Partial Least Squares Discriminant Analysis (sPLS-DA) to extract and combine discriminating features that best separate the different groups. The number of selected features was tuned using 3-fold cross-validation repeated 1,000 times. In addition, the three levels of cellular organization were integrated in an unsupervised framework. First, we performed a regularized Canonical Correlation Analysis (rCCA) between fluxes V and integrated traits T, using the mixomics package in R  [33, 51]  to search for the key features that maximized the correlation between metabolic fluxes and fermentation traits. Second, since the correlation matrix between traits and fluxes was clearly structured, we computed the matrix of Euclidean distance between traits, based on their correlation with metabolic fluxes, and clustered traits using the hclust package in R.
This procedure allowed us to define five groups of traits that showed similar correlation patterns with fluxes of the central carbon metabolism. Finally, we stored the linear correlation coefficients between proteins (P = 615) and traits (T = 28) in a (T × P) matrix and ran a Linear Discriminant Analysis to search for proteins that best discriminate between groups of traits, considering traits as individuals. Pearson's chi-square test of enrichment was computed on protein functional category frequencies taking as prior probability the expected category frequency found in the MIPS database. The same procedure was applied ignoring the flux level. rCCA was run between fluxes V and integrated traits T, and a hierarchical clustering was performed to group traits according to their correlations with protein abundances.
negative) correlation threshold of 0. 85 (resp. -0. 85). The bars represent the proportion of proteins positively/negatively correlated to the LD axis belonging to a functional category divided by the proportion of proteins from the same category found in the MIPS database. Only the categories for which this value is higher than 1 are shown. https://doi. org/10. 1371/journal. pcbi. 1009157. g005.