dataesr/scientific-paragraphs-categorization huggingface dataset:
We also used data from 11 telenivometers, which measure sub-hourly SWE and SD using gamma ray attenuation and acoustic sensors. These data were provided by the ERHIN program (Estimación de Recursos Hídricos Procedentes de la Nieve) of the Hydrological Ebro River Basin Authority (Navarro-Serrano and López-Moreno, 2017). Ten telenivometers were located in the Pyrenees, and one in the Cantabrian Mountains. A complete description of the telenivometers and their locations can be found at www. saihhebro. com. We also used an SD sensor in the Central System mountain range  (Durán et al. , 2017) , which is from the National Meteorological Agency of Spain (AEMET). We projected the meteorological variables from the WRF simulation to elevations of the different telenivometers for simulations. Figure  3  shows a comparison of the modeled and observed SD time series at these 10 sites.
It must be noted that it is challenging to validate gridded products from ground-based data  (Snauffer et al. , 2016). Snowpack can have large variability over small distances  (López-Moreno et al. , 2015; Meromy et al. , 2013). This implies that punctual measurements may not be representative of the 10 km resolution data, even when comparing a simulation at the same elevation as the telenivometer. In addition, snow measurements always include biases from the different measuring devices  (Kinar and Pomeroy, 2015). Thus, we focused on the temporal patterns of snowpack during the season. More specifically, we compared the accumulation patterns during the season, assuming that accumulation and melting rates were similar in the simulated and observational data, but that SD and SWE likely differ between the telenivometer and the simulation.
Gridded snow dataset: applications and limitations The final products of the models are daily gridded datasets (resolution: 0. 088 • , ∼ 10 km) of SD and SWE at elevations from 500 to 2900 m a. s. l. (100 m intervals) from 1980 to 2014. The datasets (ncdf4 format) cover the entire Iberian Peninsula, including the north side of the Pyrenees in France. Each dataset contains information of the entire Iberian Peninsula and a mask that covers pixels that do not present areas at the elevations of the simulation estimated from a 250 m resolution DEM. This snow database provides new opportunities for studies of snow in the Iberian Peninsula. In particular, the tem-  poral resolution and the duration of the series show significant improvements over previous observational data. Also, the geographic data generated on SD and SWE provides the opportunity to obtain more snow and hydrologically relevant information than available from remote sensing alone.
It is also possible to develop different snow products at different elevations, allowing for comparison of different elevations and different regions. For example, Fig. 6  shows the longterm average interannual maximum SD and SWE at three different elevations. Figure  7  shows examples of other snow variables that can be derived from the database: average number of snowfall  events and percentage of days with snow cover at three elevations. These analyses are particularly useful for the development of different snow climatologies for the whole Iberian Peninsula, or for specific areas, in studies that rely on ecological data (e. g. , phenology or distribution of plants and animals, forest growth), studies that require hydrological pa-rameters for different catchments, and studies that determine risk maps for snow-related events.
It is also possible to extract daily time series for different areas or elevations at each pixel. For example, Fig. 8  compares SWE series at three elevations in the pixel at the highest peak of the Pyrenees (Aneto Peak, 3404 m a. s. l. ). Thus, these series allows for the study of different annual snow accu- mulation and melting patterns at a specific location and how elevation influences snow evolution. Similarly, it enables one to study the existence of temporal trends or the occurrence of extreme snowfall and melting events. The database contains uncertainties that are not easy to quantify, due to the limited amount of observational data. Biases may be due to uncertainty of the boundary conditions from the ERA-Interim reanalysis  (Chaudhuri et al. , 2013)  since errors from the WRF downscaling model are difficult to quantify in mountain areas  (Gutmann et al.
, 2012) , and uncertainties that typically result from simulations of snow mass and energy balance from meteorological data  (Essery et al. , 1999 (Essery et al. , , 2013;; Magnusson et al. , 2015). The use of the standard air temperature lapse rate could also be a source of uncertainty. Although other studies have observed a decrease in the lapse rate during winter months, this effect is result of thermic inversions that are not considered to be due to the spatial resolution of the simulation. Despite these limitations, we had very satisfactory results when testing the duration and the interannual variability of the snowpack against MODIS and telenivometer data, which provided reliable observations during several snow seasons. This way, the database presents a reliable validation for more than a third of the time period generated.
Data availability The data presented here are freely available for download from Zenodo (https://doi. org/10. 5281/zenodo. 854618). SD and SWE datasets are in ncdf4 format, with one file for each elevation band. The observational information used to validate the main data is also available for download. All telenivometer data are in CSV format. Daily snow cover (derived from MODIS) is provided as five multiband GeoTiff files (one file for each mountain range, each band is a date), and a CSV file indicates the date of each band. The FSM code is freely available from https://github. com/ RichardEssery/FSM  (Essery, 2015).
Conclusions We presented a new daily gridded database of SD and SWE for the Iberian Peninsula from 1980 to 2014 period at a resolution of 0. 088 • (∼ 10 km). The database consists of 50 ncdf4 files for SD and SWE from 500 to 2900 m a. s. l. , and another 2 files of WRF simulation DEMs, summing more than 652 000 maps. A mask label of "no data" is included if the grid is not found at the elevation of the simulated elevation band. The scarcity of snow observations in the Iberian Peninsula made it necessary to couple a dynamic downscaling of ERA-Interim reanalysis using the WRF model by use of a snow energy and mass balance model (FSM). Input data of FSM provided directly, or estimated from WRF outputs, were available for the average elevation of each 10 km × 10 km pixel, and these data were transformed to achieve an elevation offset at 100 m intervals.
Despite some uncertainties, the database is consistent with available observational data. More specifically, validation with MODIS data indicated an error of 6. 07 % and an R 2 of 0. 76 from the analysis of the mean presence of snow. The database also provides good representation of the temporal patterns of the telenivometers, with kappa values generally over 0. 6, and above 0. 4 for all analyzed percentiles. This database will be an important resource for studies of many different hydrological, environmental, and economic processes in Mediterranean areas. Thus, we expect the database presented here will be useful for future snow-related studies at regional scales on the Iberian Peninsula, and for a broad community of researchers and land managers working in areas where snowfall occurs.
Background: Cavefish populations belonging to the Mexican tetra species Astyanax mexicanus are outstanding models to study the tempo and mode of adaptation to a radical environmental change. They are currently assigned to two main groups, the so-called "old" and "new" lineages, which would have populated several caves independently and at different times. However, we do not have yet accurate estimations of the time frames of evolution of these populations. Results: We reanalyzed the geographic distribution of mitochondrial and nuclear DNA polymorphisms and we found that these data do not support the existence of two cavefish lineages. Using IMa2, a program that allows dating population divergence in addition to demographic parameters, we found that microsatellite polymorphism strongly supports a very recent origin of cave populations (< 20,000 years).
We identified a large number of single-nucleotide polymorphisms (SNPs) in transcript sequences of pools of embryos (Pool-seq) belonging to Pachón cave population and a surface population from Texas. Based on summary statistics that can be computed with this SNP data set together with simulations of evolution of SNP polymorphisms in two recently isolated populations, we looked for sets of demographic parameters that allow the computation of summary statistics with simulated populations that are similar to the ones with the sampled populations. In most simulations for which we could find a good fit between the summary statistics of observed and simulated data, the best fit occurred when the divergence between simulated populations was less than 30,000 years. Conclusions: Although it is often assumed that some cave populations have a very ancient origin, a recent origin of these populations is strongly supported by our analyses of independent sets of nuclear DNA polymorphism.
In summary, and except an attempt using an allozymes data set unfortunately too small to give accurate estimations, no dating has ever been performed using nuclear markers to directly test the assumption that some cavefish populations are millions of years old. Several observations led us to doubt about its accuracy. In particular, looking at sequences of Pachón cavefish, we did not find many obvious loss-of-function mutations, such as frameshifts and stop codons, in eyespecific crystallin genes  [26]  and opsin genes  [44] [45] [46]  (unpublished results), an unexpected observation if this population was established for several hundred thousand years, and which becomes very unlikely if it was established more than one million years ago  [47]. Indeed, other fish that could be confined into caves for millions of years have fixed loss-of-function mutations in several opsins and crystallins genes  [48] [49] [50].
Dating with SNPs In order to make an estimation of the age of the Pachón cave population with the SNPs that is independent of the estimation made with microsatellite polymorphism, we compared the observed summary statistics of synonymous polymorphism with the summary statistics of neutral polymorphism in simulated populations using the model of evolution implemented in IMa2, i. e. an ancestral population divided into two populations at some point in the past. Each population had its own effective population size (Fig. 4 ). Migrations between populations were allowed. Whereas IMa2 is based on simulation of coalescence of a sample of alleles (backward in time), we did simulations forward in time of genetic drift of allele frequencies in whole populations, which was more adapted to analyze our SNP dataset.
Such methods have been applied to study the evolution of A. mexicanus cavefish  [32-39, 42, 43, 59 ], but neither one is of much use and often misleading if the goal is to develop a full picture that includes estimates of recent separation time and gene flow  [58]. In such case it is necessary to consider non-equilibrium models and methods allowing the joint estimation of demographic parameters (populations sizes and migration rates) and divergence time  [58]. Accordingly, we used IMa2, a widely used program for "isolation with migration" (IM) model analyses  [53]  to estimate divergence time between surface and cave populations with a dataset of multi-locus microsatellite polymorphism. IMa2 is based on backward simulations of coalescence of samples of alleles. For Pachón cave population, we estimated the divergence time using an alternative approach based on forward simulations of evolution of SNPs.
mexicanus fish, the variance of the number of descendants is probably high (in the range of 10 1 to 10 2 ) and thus N e might be one to two orders of magnitude smaller than N c  [61, 62]. If N c is in the order of magnitude of 10 4 , it is expected that N e is in the order of magnitude of 10 2 to 10 3 as found in the present study with IMa2. Previous studies, using the same data set and the program Migrate  [63] , found similar N e for several caves, including Pachón  [33, 35]. The results obtained with IMa2 suggested that the migration rate from the cave to surface is negligible whereas the migration rate from the surface to the cave is low. This is also expected. If fish could easily exit the cave, no evolution of cavefish would have occurred. Moreover, it is difficult to imagine that blind cavefish who found their way to the surface will have a good fitness there.
Dating with IMa2 With a multilocus microsatellite polymorphism dataset  [33]  and the program IMa2  [53]  divergence times of pairs of populations were estimated. The program is based on an isolation-with-migration (IM) model and uses Metropolis-coupled Markov chain (MCMC) techniques to estimate the posterior densities of the time of divergence, population sizes and gene flow  [58]. The model assumes random population samples, a stepwise mutation model, neutral mutation, freely recombining loci and constant population sizes and gene exchange rates. Although modelling constant population sizes and gene exchange rates might not be ideal in the case of cavefish populations, modelling general patterns requires simplifications and this is presently the only option in programs such as IMa2 and Migrate which is the program previously used to estimate gene flow  [33, 63].
Transcriptome assembly and annotation The Astyanax mexicanus transcriptome was assembled with Newbler ver. 2. 8 (Roche 454) sequence analysis software using 454 sequences (2. 10 6 reads) of both the Pachón cave and surface fish pooled embryos. We obtained 33,400 contigs (mean contig length = 824 bp). We also tried to generate a transcriptome assembly using the Illumina sequences, but whereas this resulted in more contigs  (49, 728)  than the 454 sequences, many of them were concatenations of different transcripts and in some cases the same transcript was found in more than one contig. We therefore mapped the Illumina sequences onto the 454 contigs to identify and annotate SNPs. Putative coding sequences in each contig were identified using the zebrafish (Zv9) proteome available at EnsEMBL 73 as a reference  [81].
There could be migrations from the surface to the cave. The smallest probability of migration was set to 0. 00001 and the highest to 0. 1 per year. The smallest percentage of cavefish that are migrants when a migration occurs was set to 0 and the highest was set to 1%. We also took into account that genetic drift could have occurred in the laboratory stocks (effective population size in the lab was set to 10 and the number of generation was set to 10). All mutations were neutral (frequency changes at each locus were driven by genetic drift only) and each locus were evolving independently. For a given set of parameters, each ten generations after the isolation of two populations, we estimated the frequency of SNPs in each category and we estimated a score of goodness of fit with the frequencies obtained with the real SNP data set.
Availability of data and materials Microsatellite data are available from  [33]. SNPs and their annotations are stored in a MySQL database and are available online at http://ngspipelines. toulouse. inra. fr:9022. Perl and R scripts for the data analyses and graphics are available upon request. Illumina and 454 transcriptomic sequences are available on ArrayExpress: https://www. ebi. ac. uk/arrayexpress/ experiments/E-MTAB-5142/. Script used can be found on GitHub: http://github. com/julienfumey.
S Y N T H E T I C C A S E In order to visualize the methodology previously described, I use a simple synthetic case where   5 ), shown in Fig. 2(b ). We then compare the true posterior pdf with its approximation provided by eq. (  8 ) obtained for σ = α max( b i 2 ), with α = 0. 5, 1. 0 and 8. 0. Figs  2(c-e ) show that the approximated posterior pdf rapidly converges towards the true pdf. We further compare the results with the posterior pdf derived from a Monte Carlo sampling. The Monte Carlo results have been obtained using the PyMC package (https://github. com/pymc-devs/pymc), using 5 × 10 5 samples, with burn and thinning values of 1000 and 10 (Fig. 2f ). Figs 2(g) and (h) show the posterior marginal pdf obtained for each of the two parameters using eq. (  11 ) to be compared with the Monte Carlo approach and the true value.
The Bayesian approach was implemented using the PyMC package (https://github. com/pymc-devs/pymc), using 10 7 samples, with burn and thinning values of 5000 and 10 respectively, leaving 9. 95 × 10 5 posterior samples. Uniform prior bounded between 0 and the plate velocity were specified for the 30 parameters. For the TMVN approach, the algorithm described in Section 3. 3 has been applied. A convergence tolerance of 2 per cent has been used, leading to values for α i between 1. 4 and 8. 7. In the following, we compare the TVMN results to the MCMC results for various quantities usually presented for Bayesian inversions. The aim is to assess the ability of the TMVN method to reproduce Bayesian results obtained using MCMC sampling for independent uniform priors for an underdetermined problem. I chose this case because it is the most difficult: the underdetermined case prevents an exact solution and uniform priors require the priors to be also approximated.
1-and 2-D marginal pdf Fig. 4  shows the comparison between the 1-D marginal pdfs obtained from the two methods. The 1-D pdfs were obtained using the dtmvnorm. marginal function from the tmvtnorm R package  (Wilhelm 2015)  which implements the method proposed by  Cartinhour (1990). Fig. 4  shows that the chosen inverse problem includes a great diversity in terms of parameter resolution. Some parameters are well defined and centred around a value between the bounds and show nearly Gaussian shape  (parameters 12, 13, 14, 17, 22 and 23) , while for other parameters, the maximum of the pdf has its maximum close to either 0 per cent of coupling (e. g. 21) or 100 per cent (e. g. 1-6). Parameters 7, 8 and 9 are badly determined, as indicated by a widely spread pdf. We can note that the resolution indicated by the marginal pdf is not a direct function of the distance between the subfault and the GPS sites as indicated by the very good resolution for subfaults 1-5.
In general, the TMVN and MCMC approaches show a very good agreement in terms of the shape and the values of the marginal pdf. The agreement is Downloaded from https://academic. oup. com/gji/article/214/1/366/4969687 by CNRS -ISTO user on 07 March 2022    3  even better for well resolved parameters. Differences increase as parameter are less resolved (e. g. 7, 9 and 29), but are always less than 20 per cent of the pdf value. Fig. 5  shows the posterior pdf for pair of parameters, that is, the 2-D covariance functions. The posterior 2-D covariance functions were obtained using the dtmvnorm. marginal2 function from the tmvtnorm R package  (Wilhelm 2015)  which implements the method proposed by  Leppard & Tallis (1989) , modified by  Wilhelm (2015). Fig. 5  shows that the TMVN method is also able to provide the covariance information with a similar accuracy as the MCMC approach.
The obtained equivalent problem is more suitable for numerical integration because it has ordered the importance of variables for the integration. In a second step, the numerical multidimensional integration is performed over the unit hypercube. For this step, different schemes have been proposed. Both methods tested in this paper to produce the results shown in Sections 4 and 5 use a so-called randomized lattice rule, where projections of integration points onto the axis produce an equidistant division of the axis, before a Monte Carlo integration  (Genz & Bretz 1999; Hothorn et al. 2001). Several implementations are available to users for scientific programming languages. Matlab and FORTRAN codes are available at Alan Genz web page (last accessed on 2017 December 14) http://www. math. wsu.
edu/f aculty/genz/software/sof tware. html. Python wrappers to Genz FORTRAN routine MVNDST. f are available either within the numerical Scipy library (scipy. stats. mvn) or statsmodels module. The most elaborated library is available for the R language with the original mvtnorm package  (Genz et al. 2006). Aside a wrapper to the mvtnorm package, the tmvtnorm package  (Wilhelm 2015)  bundles routines for 1-and 2-D marginal pdfs calculation. Furthermore,  Manjunath & Wilhelm (2009)  extending previous results from Tallis  (1961)  show that the mean for a TMVN over a hyper-rectangle can be expressed as a linear sum of the values of the individual marginal pdf at their bounds  (Manjunath & Wilhelm 2009, eq. 11). Although more complicated, a similar relationship holds for the covariance  (Manjunath & Wilhelm 2009, eq.
DNA Extraction and Metagenome Sequencing Filters duplicates were pooled, then genomic DNA was extracted as described by  Laplante et al. (2013). Library preparation (TruSeq DNA Illumina) of paired-end reads (2 × 100 bp read length) was performed by the McGill University/Genome Quebec Innovation Center for whole metagenomic shotgun sequencing using a HiSeq TM 2000 Sequencing System. A total of 30 Gbps were obtained and the sequencing data summary is shown in Supplementary File S2. The sequence files are available from the Sequence Read Archive (http://www. ncbi. nlm. nih. gov/sra), BioProject ID: PRJNA449990. FIGURE 1 | Bioinformatics analysis pipeline. Two approaches were developed for this work. With the ORF-based approach, de novo assembly was performed on raw reads data using Ray Meta software.
Then, the predicted ORFs (Open Reading Frames) were annotated using Diamond similarity research tool against SEED, which is a curated database. With the reads-based approach, merged (with FLASH algorithm) and filtered reads (with Nesoni) (length ∼200 pb in average) were annotated using the same tool as the ORF-based approach, Diamond algorithm and SEED database. For both approaches, we used the lowest common ancestor (LCA) algorithm in the taxonomic assignment and the subsystems hierarchy in function classification. Diversity measures were computed using mothur software. using the Nesoni Clip tool (https://github. com/Victorian-Bioinformatics-Consortium) version 0. 133. Overall, the quality of forward reads (R1) was better than reverse reads (R2). This difference is related to sequencing quality decrease over the length of reads, in addition to the loss of enzymatic specificity overtime in the paired-end platform technology.
Base calling quality was selected at a Phred or Q score of 33 (Supplementary File S2). FLASH software v1. 2. 11  (Magoč and Salzberg, 2011)  was used with default parameters (10-65 bp overlapping window) to merge paired-end reads. As a second step, following the selection of good quality reads for all five metagenomic samples, a sequence similarity search was performed against the SEED database (Version: May 2015)  (Overbeek et al. , 2014)  using Diamond v0. 7. 9. 58  (Buchfink et al. , 2015). The taxonomic content of each sample was assigned using the Lowest Common Ancestor (LCA) method  (Huson et al. , 2007). Functional abundance was estimated using the SUPERFOCUS software  (Silva et al. , 2015)  with Diamond (1e -12 as p-value, 60 identity as threshold, 30 base pairs as minimum alignment length).