dataesr/scientific-paragraphs-categorization huggingface dataset:
This leads to a total of 400 simulated datasets. With the Negative Binomial, the overdispersion value has been chosen as VarðY ij Þ=EðY ij Þ ¼ 5:5, in line with real data. The ratio between the lowest and highest risks is defined as 16 for Large differences in Risk levels (LR), and 4 for Slight differences in Risk levels (SR). The large differences are expected to be quite easy to detect, while the slight differences may be more challenging and explore the limits of the method. Similarly, the local overdispersion may complicate outbreak detection and the accuracy of risk estimation. S3 Fig shows an example of a dataset obtained for the scenario (COLR). We remark that in this scenario, cases are very scattered. Congruently with the overdispersion implied by the negative binomial distribution, there is a high proportion of null values (� 92%) and many high values occur (over 10% of the non-null values are over 5).
The Bayesian framework with MCMC estimation implies long calculation durations in the context of intensive simulation studies and requires parallel computing to perform the analyses, thus we used the HTCondor system to allocate jobs to the different cores (about 150). We realized scripts encoded in the C++ language to treat the resulting files. We used the software R to create all the preliminary files needed by BUGS (data, scripts and models) and HTCondor (submission files), to analyze the results and to realize maps. We also used the R libraries CODA, RPostgreSQL, rgdal, spdep and scales to import the results of BUGS, to interpret the results, and to build the maps of the risk.
Results: Best models to fit contagious simulated data We consider a wide range of models (see S1 Table ) to deal with overdispersion and chose the DIC as model comparison method. We compute the parameters for all the 400 datasets (100 replicates associated to each of the 4 scenarios) and all the 60 models. Thus these analyses provide lots of risk values, of risk maps and of DIC values whose interpretation can be delicate. Figs  4 5 6 7 show, for each model, the repartition of the DIC values based on the 100 replicates for each considered scenario. Table  2  describes the best models for each scenario according to the average value of the DIC over all the replicates. We included all the models for which the DIC does not exceed the best DIC by more than 10. We also compute average values of the DIC to determine the most relevant strategy for each scenario (Table  3 ).
Handling local overdispersion We deal with three situations, depending on the type of dataset: overdispersed data (scenarios (COLR) and (COSR)), non-overdispersed data with a weakly contrasted risk (scenario (NOLR)) and non-overdispersed data with a weakly contrasted risk (scenario (NOSR)). The best models are quite different according to this classification. Overdispersed data. The best models for the scenarios (COLR) and (COSR) have similar characteristics (Figs  4  and 5 ); namely all are negative binomial models without the Gaussian noise term (Table  2 ). These models have very similar mean values for the DIC and show the lowest DIC dispersion among all the models; however, the repartitions of the associated DIC values show slight differences (Fig  4 ). The Poisson models without the Gaussian noise � ij have also quite similar distributions despite higher DIC mean values, and they also have quite low DIC variance values.
Moreover, the use of the Gaussian white noise seems redundant when a negative binomial distribution is tested. For the Poisson models, contrary to what one could expect, the Gaussian noise � ij does not seem interesting to handle the local overdipersion, regardless of the values of the risk. Non-overdispersed highly contrasted data. We notice that the two best models for the scenario (NOLR) use a negative binomial distribution (Table  2 ). However, 16 models can be considered as relevant (distance to the best model of average DIC values smaller than 10): 9 of them are negative binomial models and the 7 others are Poisson models (without noise). Nevertheless, if we consider all the models and all the datasets, we can notice the negative binomial models generally better fit the outcome of cases when risk values are highly contrasted (scenarios (COLR) and (NOLR)) (Table  3 ).
Nevertheless, the gain is not really significant in the case of the scenario (NOSR), contrary to the scenarios (COLR), (COSR) and (NOLR), associated to local overdispersion and/or highly contrasted risk values. The best models to describe datasets for which the risk is highly contrasted (scenarios (COLR) and (NOLR)) include the spatiotemporal CAR component V ij. For the scenario (COSR), four of the five best models also include this term. The spatiotemporal CAR component V ij appears as particulary relevant when a Gaussian white noise � ij is included. The only model considered as the most relevant and which contains the Gaussian noise term � ij is nebinnopa_x_x_ST_gaus, which includes the spatiotemporal CAR process V ij. For scenarios (COLR), (COSR) and (NOLR), the best Poisson models with the Gaussian noise term � ij also integrate the spatiotemporal CAR process V ij.
Contrary to what is commonly said in the literature, the Gaussian noise term is not relevant in our case, even to model overdispersed data. However, both analyses (simulations and real data) have slight differences. In particular models with weight coefficients are irrelevant for simulated datasets but they are helpful to quantify the contribution of each CAR component to model the bovine tuberculosis data. Thus in the context of our simulations, the three CAR components provide a similar smoothing of the risk.
For the analyses of simulated data we performed, the weight coefficients were identified as irrelevant. It may be a consequence of the equilibrium between spatial and temporal correlations in our simulations. However, these weights are highly relevant for the real dataset. Thus it appears that such weights can be useful in case of desequilibrium, for instance when the spatial correlations are stronger than the temporal ones, or conversely. However it would be interesting to test them again in other contexts. This study shows the relevance of disease mapping to address rare infectious diseases. Most of the simulated hotspots have been identified, even when the risk values are weakly contrasted. The simulated risk was quite difficult to model and the delineation of the estimated risk is not very smooth.
It may be due to the very low concentration of the simulated cases which are more scattered than the real datasets. We also remarked that our method overestimates the real relative risk in regions which have a very low cattle population, as previously observed in studies based on purely spatial (and non-contagious) data. Moreover, our methodology provided very relevant maps to represent the estimated bovine tuberculosis risk. These maps are very congruent with the literature concerning the repartition of the risk in France and its evolution between 2001 and 2010. In particular, the detection of simultaneous clusters of cases in three regions is very promising and shows the good capabilities of this methodology to analyse contagious data. The resulting maps also correctly describe the evolution of the risk and its concentration over time around very small areas, congruently with the known contamination of wildlife in these regions.
Thus disease mapping appears as a very relevant way to investigate infectious diseases whose incidence is low. In our study, the model selection criterion is a crucial point as we compared 60 models for each dataset. The DIC (Deviance Information Criterion) is systematically used in the context of spatiotemporal disease mapping and largely used in the purely spatial context. For the simulated data, it provided an interesting ranking of the different models. The best models, according to the DIC, appeared as relevant, because the obtained risk was congruent with the distribution of the simulated risk and the resulting maps were interpretable. Thus this indicator provided interesting evidences concerning the structure and the repartition of the risk. However, in the bovine tuberculosis context, this criterion seemed to choose overfitted models, for which the estimated risk appeared as too temporally structured, to the detriment of the spatial smoothing.
| Genome sequencing approach Because all sequencing technologies are known to be biased in both known and unknown ways, we utilized a variety of sequencing technologies to assemble this complex and large genome. At the time of genome construction, both sequencing technologies and assembly algorithms were undergoing rapid change, and as such the generation of an optimal assembly required substantial trial and error. For instance, several of the authors of this paper were involved in the strawberry poison frog genome assembly  (Rogers et al. , 2018)  and the difficulties encountered in the attempt to assemble that genome made it clear that short read data alone are insufficient if the goal is to assemble a highly contiguous and complete genome. To overcome the limitations associated with short read data, we collected linked and long read data from technologies thought to be complementary to one another (Illumina 10X, Oxford Nanopore and Pacific Biosystems).
To do this we first mapped the Illumina short reads to our assembly using bwa version 0. 7. 17-r1188  (Li & Durbin, 2009). Polishing was then conducted using the software pilon version 1. 22  (Walker et al. , 2014). We used a flank size of five bases (thus ignoring the last five well-aligned bases on either end of the alignment) and fixed both bases and gaps with a minimum size of one base. We then used arcs version 3. 82  (Coombe et al. , 2018)  to scaffold our existing assembly using Illumina 10X data. Prior to doing this, we used the "basic" function within the program longranger  (Marks et al. , 2019)  to trim, error correct and identify barcodes in the 10X data. We then used the "arks. mk" makefile provided in the arcs GitHub page (https://github. com/bcgsc/ arcs/blob/maste r/Examp les/arcsmake), to run the arcs software.
We then attempted to fill in gaps within our assembly using lr gapcloser  (Xu et al. , 2019). We ran this program iteratively, first using the Nanopore data and specifying the "-s n" flag. We then filled gaps using lr gapcloser with PacBio reads and the "-s p" flag. Following this, we again polished our assembly using bwa and pilon. We examined genome quality in two main ways. First, we examined the presence of genic content in our genome using Benchmarking Universal Single-Copy Orthologs version 3. 0. 0 (busco;  Simão et al. , 2015)  using the tetrapod database (tetrapoda_odb9; 2016-02-13). busco version 4. 0 was released during this project, so we also examined our final assembly with busco version 4. 0. 6 and the tetrapod database  (tetrapoda_odb10; 2019-11-20). Our second method of genome examination was genome contiguity.
| Repeats and genome annotation We modelled genomic repeats with repeat modeler version 1. 0. 8 (Smit & Hubley, 2008 -2015)  using RepBase database 20170127. We used F I G U R E 2 Flowchart of our genome assembly approach. Colours of boxes represent the type of data used in that step (see internal key) and italicized font indicates the program(s) used the classified consensus output from repeat modeler as input to repeat masker version 4. 1. 0 (Smit et al. , 2013 (Smit et al. , -2015) )  with the Homo sapiens database specified. We used a combined database of Dfam_3. 1 and RepBase-20181026  (Bao et al. , 2015). We annotated our genome using maker version 3. 01. 02 (Campbell et al. , 2014). We used transcript evidence from R. imitator to aid in assembly ("est2genome=1"). We produced transcript evidence from gene expression data of various R.
imitator tissues and populations; for additional details on transcript evidence used to annotate this genome please see the Supplemental Methods. We annotated the genome using both the draft genome as well as a repeat masked genome from repeat modeler and repeat masker in order to identify the best annotation of our genome. We used busco version 4. 0. 6  (Simão et al. , 2015)  on the final version of the R. imitator genome hto locate putative duplicated orthologues within the genome. busco searches a database of universal single copy orthologues to measure genic completeness of the genome. busco identifies whether these single copy orthologues are present in the genome and if they are complete or fragmented. Additionally, it identifies whether these orthologues are present in a single copy or duplicated in the genome.
Given the high proportion of duplicated orthologues present in our genome, we proceeded to investigate the read support for each of these duplicated regions across our sequencing technologies. We aligned PacBio data using minimap2 version 2. 10  (Li, 2018). After we aligned the data, we used samtools version 1. 10  (Li et al. , 2009)  to generate the depth at each base in the genome for data from each sequencing platform. We wrote a custom Python script (see code availability statement for link to code repository) to extract the depth per base within the duplicated regions of the genome. Additionally, we extracted the genic regions identified as single copy genes by busco in the same fashion. This allowed us to examine if there were any discrepancies between the duplicates and the rest of the genomic sequence.
| Differential gene expression We indexed our new R. imitator genome using star version 2. 5. 4a (Dobin et al. , 2013). We removed adaptor sequences from reads using trimmomatic version 0. 39  (Bolger et al. , 2014). We then aligned our trimmed reads to our genome using star version 2. 5. 4a (Dobin et al. , 2013) , allowing 10 base mismatches (--outFilterMismatchNmax 10), a maximum of 20 multiple alignments per read (--outFilterMultimapNmax 20) and discarding reads that mapped at <50% of the read length (--outFilterScoreMinOverLread 0. 5). We then counted aligned reads using htseq-count version 0. 11. 3  (Anders et al. , 2015). Differential expression analyses were conducted in R version 3. 6. 0 (Team, 2019)  using the package DESeq2  (Love et al. , 2014). Some genes in our annotated genome are represented multiple times, and thus the alignment is nearly to the gene level with some exceptions.
As a result, when we imported data into R we corrected for this by merging counts from htseq-count into a gene-level count. We filtered out low-expression genes by removing any gene with a total experiment-wide expression level ≤50 total counts. cDNA libraries for R. imitator were sequenced at a different core facility than those of R. fantastica and R. variabilis, so in order to statistically account for batch effects we analysed the data from each species independently (combining all species and batch effects in our data set produces rank deficient models). For each species we compared two models using likelihood ratio tests, one which tested the effect of colour morph and the other which tested the effect of developmental stage. Both models included sequencing lane, developmental stage and colour morph as fixed effects.
We used a Benjamini and Hochberg (1995) correction for multiple comparisons and used an alpha value of. 01 for significance. We then extracted data from our models for particular a priori colour genes that play a role in colour or pattern production in other taxa. This a priori list was originally used in  Stuckert et al. (2019) , but was updated by searching for genes that have been implicated in coloration in genomics studies from the last 3 years. Plots in this paper were produced using ggplot2  (Wickham, 2011). Finally, we ran an analysis with the specific intent of identifying genes involved in the production of different colour morphs that are convergent between model (R. fantastica or R. variabilis) and mimic (R. imitator). To do this we conducted a Walds test within species between the spotted and striped morph of R.
| Weighted gene correlation network analysis To identify networks of genes that interact in response to differences in colour morphs, species or developmental stage, we used weighted gene co-expression network analysis (WGCNA), an approach used for finding clusters of transcripts with highly correlated expression levels and for relating them to phenotypic traits, using the package WGCNA  (Langfelder & Horvath, 2008). Although WGCNA names modules after colours, these name designations are unrelated to any phenotype in our data. For the WGCNA analysis, we used the variance stabilizing transformed data produced by deseq2 at this stage. WGCNA requires filtering out genes with low expression, which we did prior to all differential expression analyses. We estimated a soft threshold power (β) that fits our data, by plotting this value against mean connectivity to determine the minimum value at which mean connectivity asymptotes, which represents scale-free topology.
| Repetitive elements Our analyses indicate a high proportion of repeats in the R. imitator genome (see Figure  4 ). repeat modeler masked 54. 14% of total bases in the genome, of which 50. 3% consisted of repeat elements (see    S2. A large proportion of tetrapod orthologues (i. e. , those genes in busco's tetrapoda database) are present in duplicate copies that are spread throughout scaffolds in the genome (23. 2%). On average, scaffolds had 1. 92 ± 2. 93 duplicated orthologues, with a median and mode of 1. When normalized per 10,000 bp, the average was 0. 060 ± 0. 108 SD duplicated orthologues per 10,000 bp of a scaffold, with a median of 0. 032. Intriguingly, the majority of orthologues identified as single copies by busco had average coverage greater than the average genome-wide coverage.
| Gene expression We aligned an average of 20. 97 million reads (±6. 7 SD) per sample. On average, 73. 27% of reads were uniquely mapped (±4. 75% SD). Mapping rates were slightly higher in R. imitator than in R. fantastica or R. variabilis, because the libraries were of slightly higher quality in R. imitator. To further test if this was an artefact derived from mapping reads from other species to the R. imitator genome, we also mapped these reads to species-specific transcriptome assemblies. We assembled these using data for each species in this study using the Oyster River Protocol (MacManes, 2018); additional details on this protocol can be found in the Appendix. After mapping our read data to these species-specific transcriptome assemblies, we found similar results to that of our genome-guided mapping.
Thus, our mapping rates are driven primarily by the slightly lower quality of the R. fantastica and R. variabilis cDNA libraries (as exhibited by slightly lower Phred scores towards the 3' end of reads) and not by species-specific differences in coding regions. For data on the number of reads and mapping rates in each sample Table  S3. All gene expression count data can be found in the GitHub repository (https://github. com/AdamS tucke rt/Ranit omeya_imita tor_genom e/tree/maste r/GeneE xpres sion/data). Patterns of gene expression are largely driven by developmental stage (principal component 1; 43% of the variation) and species (principal component 2; 13% of the variance; see Figure  6 ). We note that differences in sample preparation and sequencing localities between R. imitator and R. fantastica/variabilis may be driving a portion of the pattern in principal component 2.
| Weighted gene correlation network analysis We identified 20 individual gene modules in our expression data using WGCNA (see Figure  10  WGCNA modules). Of these, we (r = -. 345; p =. 00017), and salmon (r =. 431; p <. 0001). We identified 10 modules that were significantly correlated with species; green (r = -. 280; p =. 0025), pink (r = -. 566; p <. 0001), purple (r =. 308; p =. 000857), greenyellow (r =. 320; p =. 000529), light cyan (r =. 226; p =. 0156), grey 60 (r =. 196; p =. 0366), cyan (r = -. 209; p =. 0257), magenta (r =. 304; p =. 00103), salmon (r = -. 308; p =. 000869) and grey (r =. 656; p <. 0001). Finally, we identified three modules that were significantly correlated with colour morph, the grey (r =. 403; p <. 0001), pink (r = -. 372; p <. 0001), and purple modules (r =.
F I G U R E 1 0 The relationship between WGCNA modules and treatment groups. Within each cell the top number represents the R value (correlation between the grouping and module expression) and the bottom number represents the adjusted p value. The strength of the colour (i. e. , darker red or darker blue) in the heatmap represents the strength of the association Somewhat surprisingly, there was minimal overlap in genes that were differentially expressed between convergent colour morphs in both the model and mimic. Of those genes that were shared, only three were in our a priori colour gene list (edaradd, pkc-3 and qdpr). This seems to indicate one of three things: (i) the pattern of differential gene expression affecting convergent colour morph development is largely species-specific, (ii) the expression patterns of a small number of genes have a very large effect on colour morph divergence, or (iii) we have insufficient power to identify these convergent genes.
, 2020; Salis et al. , 2019). Hooper et al. (2019)  found that this gene is associated with bill colour in the long-tailed finch, Poephila acuticauda, and they hypothesized that ttc39b is being used to transport hydrophobic carotenoids to their deposition site. We found that this gene was most highly expressed in the yellow-green spotted morph of R. variabilis and the redheaded morph of R. fantastica, consistent with findings from other gene expression studies. Given the repeat occurrence of this gene in pigmentation studies outlined above, functional validation of ttc39b is likely to be informative. While we found very few putative carotenoid genes identified in other taxa that were differentially expressed between our colour morphs, we nevertheless identified the oxidoreductase activity gene ontology as the primary gene ontology term in the purple module associated with colour morphs.
ing 7. 45% GYPSY/DIRS1 repeats. Given the quality of repeat databases and the scarcity of amphibian genomic resources in these databases, our results probably represent an underrepresentation of repeats in the genome as a whole. A fairly large proportion of the genome's repeat elements are likely to be unassembled and missing in the genome, contributing to our underestimation of repeat F I G U R E 3 Summary figures of genome assembly. (a) Distribution of scaffolds by length. (b) The cumulative proportion of the genome represented by scaffolds. Scaffolds were ordered from longest to shortest prior to plotting, and therefore in this panel the longest scaffolds are represented in the left portion of the figure. (c) Violin plot of the distribution of scaffold lengths in the genome. (d) Violin plot of the length of transcripts in the genome as predicted by maker.
Funding information Funding for this project was provided by an  East Carolina University Thomas Harriot College of Arts and Sciences Advancement Council Distinguished Professorship  and  NSF   DEB 165536  to K. S. ,  NSF   DEB 1655585  to M. D. M. , XSEDE XRAC  MCB110134  to M. D. M. and A. M. M. S. , and M. C. received support from  Agence Nationale de la Recherche (ANR)  grants RANAPOSA (ref. ANR-20-CE02-0003 ), from an  "Investissement d'Avenir" grant CEBA (ref. ANR-10-LABX-25-01 ) and from a  Marie Sklodowska-Curie fellowship  ( FITINV, N 655857 ).
multocida serotypes. The aim of the present study was to investigate the potential epitopes of the PlpE protein from P. multocida A:1, one of the causative agents of pasteurellosis. The protein was mapped in terms of immunogenicity and also evaluated in view of hydrophilicity, antigenicity, half time, stability and so on. Our data indicated that three different regions, each consisting of one or more epitopes, are responsible for protection. These regions are suitable to express as a possible recombinant subunit vaccine.