dataesr/scientific-paragraphs-categorization huggingface dataset:
Gene-Based GWAS We used Scoary to associate specific genes with the presence or absence of specific traits (Bt species, Bt subspecies Bta, Btk, Bti and Btm, and Bt clusters a to d) with a GWAS analysis. Three GWAS analyses were run to identify a large number of candidate genes for each trait (Table  2 ). A total of 1304 genes were significantly associated with the traits of interest. Interestingly, the majority of the identified genes were associated with subspecies israelensis and morrisoni in the GWAS L2 analysis, with 809 and 246 genes, respectively. However, this large number of associated genes can be explained by the low number of genomes available for these two subspecies in the dataset reducing the analytical strength. Given the low diversity in the dataset for these subspecies, we did not pursue the search for specific markers for these subspecies.
Validation of GWAS Results We then validated the markers selected using GWAS by comparing them in silico with the VD using TBLASTN. For each selected candidate gene associated with a specific trait in the GWAS analysis, a protein sequence was extracted and compared to all VD genomes that shared or did not share the trait in question. Out of the 249 highly associated genes obtained with Scoary, 128 protein sequences met the criteria of coverage, specificity and sensitivity of the TBLASTN analysis (Table  2 ). All validated marker sequences had a sensitivity range between 80% and 100%, and 100% specificity on both SD and VD datasets. When a single marker was not sufficient, combinations of markers were tested to identify all the associated genomes for each trait.
Functions and Distribution of Specific Markers To predict whether the genes coding for the markers were localized or not on a plasmid, a TBLASTN alignment on a Bt plasmid database was carried out. When the function was known, the closest protein identified in Uniprot  [46]  was indicated (Supplementary Table  S2 ). Table  3  lists 14 markers associated with low p-values, and whose combination could detect all genomes of a given trait (i. e. , 100% sensitivity). The Bt species can thus be identified with a combination of only six markers (cwlA, intQ, group_3916, group_20749, group_20361 and sdpR). For example, the cwlA gene located in an operon along with cwlB, encodes a peptidoglycan hydrolase known for its role in cell lysis during Bt spore release  [47]. Regarding the markers that emerged from the GWAS L2 analysis for Bt subspecies, all Bta genomes (n = 56) could be identified based on the apr gene, which encodes a specific chromosomal alkaline serine protease, also called subtilisin, in Bacillus subtilis  [48].
2. 4. Workflow for Bt Identification Based on the 14 selected markers, a workflow was develo   [34]  core-gene alignment of the study dataset (SD) genomes. Visualization and annotation were performed using iTol  [43]. Filled circles indicate the presence of the 14 Bt markers for the successive GWAS analysis levels (L1, L2 and L3) in the corresponding SD genomes. The 14 Bt-specific markers identified with the three levels of GWAS analysis and validated by TBLASTN were mapped on the phylogenetic tree (Figure  3 ). Despite the great distance between some Bt strains, a combination of six markers can identify all of them. Similarly, with two markers specific to Btk and cluster b, all associated genomes can be retrieved. Only one marker was needed for the identification of Bta, clusters a, c and d.
Workflow for Bt Identification Based on the 14 selected markers, a workflow was developed to distinguish Bt from the genome of a putative Bc isolate (Figure  4 ). We wrote an in-house script in Python 3 to automate the screening of the markers in genome assemblies. Using BLASTX, the script detects the markers for Bt, Bt subspecies and Bt clusters using the same thresholds as the TBLASTN validation (90% identity and coverage) in the user-input genomes. The output indicates the presence (+) or absence (-) of all trait-associated markers for each genome query. Then, according to the workflow pathway, based on the presence/absence of the markers within the genome, a putative Bt identification is given. The script is freely available at https://github. com/afelten-Anses/Bt_typing (accessed on 2 February 2022).
Discussion To differentiate Bt from other members of the Bc group and especially in the context of food poisoning, we analyzed a large database of Bc to identify genetic markers specific to Bt. The construction of the Bc pangenome from 230 Bc genomes showed that it was composed of 39,021 protein-coding genes and that the core genome corresponded to 1854 genes. These results differ somewhat from those previously published for the Bc pangenome. A previous pangenome study of 114 Bc highlighted 59,989 genes; however, only 600 were described as core genes  [23]. Normally, increasing the number of genomes in a dataset can lead to a decrease in the core genome size, along with an increase in the accessory genome size. However, the difference in tools used for pangenome deduction here can explain this large difference in gene numbers, as already highlighted in a genomewide analysis of Clostridium difficile pangenome  [53].
For instance, the previous popular pangenome analysis software, named Roary  [54] , uses a default threshold of 95% identity (versus 98% with Panaroo) to create clustered protein sequences. Nonetheless, Panaroo, unlike Roary, clusters potential close family genes together. Moreover, Roary does not take assembly and annotation errors into account, thereby resulting in an overestimation of the total gene number and thus a reduction of the core genome size, especially for datasets with high diversity as in the Bc group. In addition, many Bt genomes (especially from the Btk and Bta subspecies) that were used for this study comprised a large number of accessory genes, which can be explained   3 ). A marker is considered present when its sequence coverage and its identity are both higher than 90% in the assembly tested.
Discussion To differentiate Bt from other members of the Bc group and especially in the context of food poisoning, we analyzed a large database of Bc to identify genetic markers specific to Bt. The construction of the Bc pangenome from 230 Bc genomes showed that it was composed of 39,021 protein-coding genes and that the core genome corresponded to 1854 genes. These results differ somewhat from those previously published for the Bc pangenome. A previous pangenome study of 114 Bc highlighted 59,989 genes; however, only 600 were described as core genes  [23]. Normally, increasing the number of genomes in a dataset can lead to a decrease in the core genome size, along with an increase in the accessory genome size. However, the difference in tools used for pangenome deduction here can explain this large difference in gene numbers, as already highlighted in a genome-wide analysis of Clostridium difficile pangenome  [53].
For instance, the previous popular pangenome analysis software, named Roary  [54] , uses a default threshold of 95% identity (versus 98% with Panaroo) to create clustered protein sequences. Nonetheless, Panaroo, unlike Roary, clusters potential close family genes together. Moreover, Roary does not take assembly and annotation errors into account, thereby resulting in an overestimation of the total gene number and thus a reduction of the core genome size, especially for datasets with high diversity as in the Bc group. In addition, many Bt genomes (especially from the Btk and Bta subspecies) that were used for this study comprised a large number of accessory genes, which can be explained by the high diversity of the genomes in these subspecies in the dataset, and/or the presence of mobile genetic elements, particularly plasmids.
A total of 249 candidate genes were selected based on their sensitivity and specificity for a given trait. As recommended for GWAS analyses  [60] , the results obtained were validated by testing them on a test dataset (here, the VD). The fact that almost 50% of the selected markers did not pass this step demonstrates the importance of testing results. To extend identification to other traits or Bt genomes, marker selection preferentially targeted genes with a chromosomal location and a known or predicted function. Of the 128 validated genes, only 21 had predicted functions, annotated with Prokka  [33]  (Supplementary Table  S3 ). At least one gene associated with a known function and located on the chromosome was identified for L1-and L2-associated traits (Bt species and subspecies, respectively), with the exception of Btk.
However, we showed that Btk could be identified based on the combination of two plasmid markers, which underlines the importance of plasmid content, especially for the discrimination of closely related Bt strains. The phylogenic visualization of the SD dataset illustrates the challenge-even when using genomic approaches-of differentiating Bt populations from each other or other Bc group members, due to their close genetic proximity  [3, 23]. For example, L3 analysis performed for cluster identification could not be performed on the whole SD. The impossibility of identifying specific genes revealed a limitation of gene-based GWAS when comparing extremely close groups. Fortunately, the sequential use of a combination of markers allows the use of cluster-specific markers after confirmation of Bt membership and then one of two subspecies, Bta or Btk.
Furthermore, the use of a highly diverse dataset in this study allowed the identification of specific and high sensitivity markers. For example, among the Bt markers identified, the cwlA gene was found in 93. 8% of SD Bt genomes (n = 144), making it a very reliable marker compared with the previously proposed detection system [61  -63]. Furthermore, this previous system did not include commercial Bt strains, unlike the analyses performed for this study. Nevertheless, it is noteworthy that the identified markers refer to a specific dataset, meaning they may not allow the identification of the entire Bt species. In addition, we cannot exclude the existence of some non-sequenced strains, in particular divergent strains not included in our dataset that may not possess one or more of the selected markers.
Based on these 14 markers, we developed a new tool, associated with a workflow (Figure  4 ) and a script to predict the identity putative Bt among Bc isolates. Moreover, the gene-based GWAS approach developed here demonstrated that the four proximity clusters associated with commercial Bta or Btk strains likely belong to Bt  [17]. Currently, the workflow cannot be used to distinguish certain Bt strains (particularly the insecticide strains belonging to the kurstaki subspecies) within the same cluster. Gene-based approaches are limited in their specific identification of commercial Bt strains used in agriculture. New approaches, notably based on SNP calling, can be conducted to search for differentiation methods at the strain level. In this perspective, our workflow can be easily updated and optimized, if necessary.
Supplementary Materials: The following supporting information can be downloaded at: https: //www. mdpi. com/article/10. 3390/foods11233924/s1, Figure  S1. Rarefaction curve of gene diversity as a function of the number of samples generated with non-parametric incidence-based estimator jack1  [64]. Table  S1. List of the 286 genomes used in this study. The genome set was divided into two datasets, a study dataset (SD) (n = 230) and a validation dataset (VD) (n = 56). a Clusters were determined by SNP calling, using iVARCall2  [17]. b All the sequencing data used in this study are associated with two BioProjects: PRJNA547495  [17]  and PRJNA781790 (this study). The attribution to panC groups was performed after partial sequencing of panC gene, and according to Guinebretière et al. ,  [24]. Abbreviations: FBO = foodborne outbreak.
Figure1. Flowchart of the pan-GWAS analysis. Pangenome construction was first performed using Panaroo [34]  to obtain an estimation of the core and accessory genomes of the study dataset (SD). Then, the presence/absence matrix output from Panaroo was used for the GWAS analysis performed using Scoary, thereby providing a list of associated genes for each trait. Only genes with at least 80% sensitivity and 100% specificity for the SD were selected. A TBLASTN validation step of the selected genes was performed after recovering the protein sequences from Panaroo. The genes were validated when the respective sequence criteria of coverage and identity were met and when exhibiting at least 80% sensitivity and 100% specificity for the validation dataset (VD). The core-gene alignment obtained from Panaroo was used to construct the phylogeny of the 230 genomes, using IQtree [36].
Funding: This project was funded and supported by Ecophyto II (project "BtImpact", reference number OFB. 21. 0450) and ANSES (PPV project "OSABt", reference number PPV19/190160), via the tax on sales of plant protection products. The proceeds of this tax are allocated to ANSES to finance the development of a system to monitor the adverse effects of plant protection products, called "phytopharmacovigilance" (PPV), established by the French Act on the Future of Agriculture, Food and Forestry of 13 October 2014. A. Gallet has been supported by the French government, through the UCAJEDI Investments in the Future project managed by the National Research Agency (ANR) with the reference number ANR-15-IDEX-01. A. Fichant has been funded by an INRAE-ANSES thesis co-funding.
FOR PEER REVIEW 7 of 16 Figure 2. 162 Figure 2. Phylogenetic tree of 230 Bacillus cereus (Bc) genomes alongside the gene presence-absence variation matrix. The maximum-likelihood phylogeny (GTR + F + R10 substitution model) was performed with IQtree [36] using the Panaroo [34] core-gene alignment of the study dataset (SD) genomes and visualized with Phandango [42]. The presence-absence matrix obtained with Panaroo identified 1854 core genes (present in 99% of SD genomes) and 37,167 accessory genes for the 230 Bc genomes.
Foods 2022 ,Figure 3. 20223 Figure 3. Core genome phylogeny of 230 Bacillus cereus (Bc) genom phylogeny (GTR + F + R10 substitution model) was performed with IQ core-gene alignment of the study dataset (SD) genomes. Visualizati formed using iTol [43]. Filled circles indicate the presence of the 14 B GWAS analysis levels (L1, L2 and L3) in the corresponding SD genom.
of Genome(s) Representative Genome NCBI Accession Number Bacillus anthracis 1 Ames NC_003997. 3 B. cereus (sensu stricto) (Bc s. s. ) 2 ATCC 14579 NC_004722. 1 B. cytotoxicus 1 NVH 391-98 NC_009674. 1 B. mycoides 1 DSM 2048 NZ_CM000742. 1 B. pseudomycoides 1 DSM 12442 NZ_CM000745. 1 B. thuringiensis (Bt) Bt ssp. aizawai 56 Leapi01 AMXS00000000. 2 Bt ssp. kurstaki 57 HD-1 NZ_CP004870. 1 Bt ssp. israelensis 5 AM6552 NZ_CP013275. 1 Bt ssp. morrisoni 1 BGSC 4AA1 NZ_CP010577. 1 other and/or unknown 25 ATCC 10792 NZ_CP021061. 1 B. toyonensis 1 BCT-7112 NC_022781. 1 B. weihenstephanensis 1 WSBC 10204 NZ_CP009746. 1 B. wiedmannii 1 MM3 NZ_CM000718. 1 other B. cereus (sensu lato) * 77.
Study design and population We used data from the Elfe cohort, a prospective nationally representative birth cohort study initiated in 2011. The general objective of the Elfe cohort was to examine the determinants of the child's development, health and socialization from birth to adulthood. The study design and protocol have been described in detail  [35]. In brief, 18,329 newborns were included from a random sample of 349 maternity units. They were born after 33 weeks' gestation to mothers > 18 years old who were not planning to move outside of Metropolitan France in the following 3 years and were able to read French, Arabic, Turkish, or English. Participation rate at inclusion was 51%. To correct for non-representativeness, basic information on non-participating mothers was collected and a weighting procedure was developed to account for the sampling plan and non-participation at both maternity ward and individual levels.
Mothers provided written consent for their own and their child's participation. Fathers provided written consent for the child's participation when present at inclusion or were informed about their rights to oppose it. The Elfe study was approved by the Advisory Committee for Treatment of Health Research Information (Comité Consultatif sur le Traitement des Informations pour la Recherche en Santé), the National Data Protection Authority (CNIL) and the National Statistics Council. After a face-to-face questionnaire administered to mothers by midwives at inclusion, follow-up surveys were conducted by phone or mail. Here, we used data collected at 2 years by phone interviews with children's mothers and fathers, when available. The survey at 2 years included responses for 13,528 toddlers. The present analysis was based on a sample of 13,117 toddlers with valid data on screen exposure at age 2 years.
We investigated the associations between parental leisure activities and adherence to the no-screen guideline by using bivariable and multivariable logistic regression models, estimating odds ratios (ORs) and 95% confidence intervals (CIs). We computed the likelihood of adherence to the no-screen guideline for each additional hour of parental screen time. We conducted interaction tests and stratified the analyses on three key sociodemographic characteristics, consistently identified as factors of screen viewing in children in previous studies: maternal age, educational attainment, and parental migration status. Missing data for factors included in the multivariable models were imputed by multiple imputation and the fully conditional specification method. Five datasets were generated, and the estimates and their 95% CIs were pooled.
Thus, in households with low socioeconomic status, adherence to the guideline may be low regardless of parental screen use, whereas in households with high socioeconomic status, toddlers' exposure to screens may occur when parental screen use increases. The main strength of our study is that we analysed data from a nationally representative birth cohort with a very large sample size. Original and extensive measurements of leisure activities of both mothers and fathers allowed us to explore these understudied characteristics. We were able to investigate a broad range of sociodemographic variables related to screen exposure in young children, and explore associations of parental screen time and guideline adherence in relation to key sociodemographic variables. Our study has limitations. Because of social desirability, parents may have under-reported their child's exposure to screens.
Conclusions In this nationally representative sample, 13. 5% of parents adhered to the no-screen guideline for toddlers. We identified sociodemographic factors, parental leisure activities and screen time as factors of guideline adherence. This knowledge may be important on a national scale to promote parental adherence to the no-screen guideline and reframe screen use as an issue that concerns the whole family, not just the child. Intervention studies may benefit from targeting parental screen time and promoting alternative activities at the household level, promoting parental strategies, including modelling, for reducing opportunities for sedentary behaviours, thus increasing opportunities for physical activity. • fast, convenient online submission • thorough peer review by experienced researchers in your field • rapid publication on acceptance • support for research data, including large and complex data types • gold Open Access which fosters wider collaboration and increased citations maximum visibility for your research: over 100M website views per year.
Acknowledgements We thank the scientific coordinators ( B Geay ,  H Léridon ,  C Bois ,  JL Lanoé ,  X Thierry ,  C Zaros ), IT and data managers, statisticians ( T Simeon ,  A Candea ), administrative and family communication staff, study technicians ( C Guevel ,  M Zoubiri ,  G Meyer ,  I, Milan ,  R Popa ) of the Elfe coordination team, as well as the families that gave their time for the study. We thank  Ms. Laura Smales  for reviewing and editing the English language.
Biomedicine, despite its numerous flavors  [11] , is generally understood as "the investigation of pathological situations using biological theories and methods: the laws of the normal applied to pathological bodies"  [12]. Nowadays this includes modern disciplines like bioinformatics, computational and systems biology  [13]  as well as network medicine  [14, 15] , and we use this term to indicate the most extensively recognized medical approach, in other words, the clinical gold standard. To take this into account we use PubMed (https://pubmed. ncbi. nlm. nih. gov/) as our reference database. Second, as an exemplar ailment to drive our analysis we use rheumatoid arthritis (RA). RA is a recognized model of pathology for chronic, inflammatory autoimmune diseases, and hence all related considerations are, in principle, translatable to other NCDs.
gov, NCT03944096). Identification of the levels of the variables (see below) was inferred, where not available, in the trial record from other FMT studies  [33]. For the sake of readability, from now on all 28 items in our selection are referred to as articles. To further enhance the discussion we introduce three additional categories, representing therapy macro-clusters: PHA, EXP and USTD. In particular, the pharmacological gold standard therapy (from now on PHA) includes drugs aimed at the symptomatic control of inflammation, such as non-steroidal anti-inflammatory drugs (NSAIDs), and paracetamol or morphine(-derived) analgesic drugs; corticosteroids, which operate beyond symptoms and attempt to counteract the degeneration associated with the disease; and finally disease modifying anti-rheumatic drugs (DMARDs) able to interfere and block more (biologic, bDMARDs) or less (conventional, cDMARDs) specific targets, the activity of pro-inflammatory molecules (cytokines) or cells, once the disease has been diagnosed  [16] [17] [18].
In more formal terms, for each biomedical category C i , i = 1. N C , N c = 3, each retained socio-anthropological variable v j , j = 1. N v , N v = 6 and each level (value) of the categorical variable l i z j ; z ¼ 1::N l j , with N l j the number of levels in variable v j , the N c = 3 reference vectors' values V i j are built as follows: v i j ¼ V i j ; if V i j ¼ max z¼1::N l j PN C i i¼1 l i z j N C i 0 @ 1 A > N C i 2 NA; V i j � N C i 2 8 > > > > < > > > > : 9 > > > > = > > > > ; Unsupervised. To search for emergent properties, given the low dimensionality of the space and the dataset, statistical methods were not applicable. Therefore, we chose a graphical approach to visualize up to six dimensions (variables). Intuitively, we positioned the articlesrepresented as points-in a six-dimensional space and looked for areas where these points are denser.
Numerical experiments All numerical computations were performed in the Julia programming language (version 1. 4. 2). The bifurcation diagrams were computed using a pseudo-arclength continuation method implemented in the package BifurcationKit. jl [63]  with version v0. 1. 2. The continuation is based on a Newton-Krylov method to solve  (12)  with GMRES linear solver  [5]. The computation of the eigen-elements, to detect the bifurcation points, is based on the Arnoldi algorithm  [52]. The linear and eigen solvers are both implemented in the package KrylovKit. jl. The nonlinear equations  (12)  were solved at tolerance 10 -11 in the supremum norm. The bifurcation points were located using a bisection algorithm (on the number of unstable eigenvalues) leading to a precision of 10 -4 on the value of the parameter γ, see  (2).