dataesr/scientific-paragraphs-categorization huggingface dataset:
This suggests that the introduction of a dual treatment scheme, including ceftriaxone, succeeded in slowing the transmission of the second wave of  ST-1901. In fact, the few cases of full ceftriaxone resistance we detected were inferred to have evolved in Europe during the period of ESC monotherapy  (2009) (2010) (2011) (2012). Our reconstruction of the global dissemination of  ST-1901  Clade A points to a central role of antibiotic-induced selection in shaping the evolution and spread of this strain, in line with other recent studies of the gonococcus  (Sánchez-Busó et al. 2019; Golparian et al. 2020). We cannot rule out the existence of additional drivers of the expansion of Clade A, such as prevalence in high-risk groups, as the available metadata is lacking in this regard. In terms of genome dynamics, Clade A does not seem to stand out, as we infer nearly identical evolutionary rates and rates of recombination (Fig.
Sample collection Genomes from the following sources were included: (1) Isolates received at the Norwegian Institute of Public Health between 2016-2018 that matched our typing-specific inclusion criteria (see below); (2) Isolates obtained primarily from Scotland, retained within the Scottish Bacterial STI laboratory as ST-1901 from the period 2015-2016; (3) Isolates from a number of published studies  (Grad et al. 2014; Kwong et al. 2016; Kwong et al. 2018; Lee et al. 2018; Yahara et al. 2018; Peng et al. 2019; Sánchez-Busó et al. 2019). When typing information was not available, the genomes were screened, and those that matched our inclusion criteria retained. NCBI BioProject accessions for these are PRJNA431691, PRJDB6496, PRJNA394216, PRJNA266539, PRJEB4024, PRJEB2999 and PRJEB17738; (4) Isolates available at PathogenWatch (https://pathogen.
watch/) as of May 2019 that matched our inclusion criteria. Only isolates with provided collection year and country of isolation were included. In addition to isolates matching the ST-1901 MLST profile, we included isolates belonging to ST-1579, ST-7360, ST-13590 and ST-8137. The first two STs have been described as single-locus variants of  ST-1901 (Yahara et al. 2018 ) and all four STs were found embedded within the ST-1901 genome-based phylogeny in a previous analysis of Norwegian gonococcal genomes  (Alfsnes et al. 2020). Multilocus sequence typing (MLST) was performed in silico (https://github. com/tseemann/mlst) using the pubMLSTdatabase  (Jolley et al. 2018 ) to determine the sequence type. Extended data on all isolates is available in the Figshare repository 10. 6084/m9. figshare. 12423596.
Unicycler v0. 4. 7  (Wick et al. 2017 ) was used to de novo assemble a high-quality reference genome using jointly the long Oxford Nanopore reads and short accurate Illumina reads. The resulting assembly resolved a circular genome of 2,222,926 base pairs. All short read assemblies were aligned to this closed reference genome using Parsnp v1. 2  (Treangen et al. 2014). A whole-genome multifasta was generated by retaining the reference nucleotide position for all sites by filling all non-core regions with reference nucleotides using an in-house script (https://github. com/krepsen/parsnp2fasta). In total, 1,716,165 sites were retained by Parsnp, and the remaining sites filled with reference nucleotides. Assembled genomes were annotated using Prokka v1. 14. 6. Prior to pangenome reconstruction, the annotated genomes were again screened with Mash  (Ondov et al.
Output from the Panaroo (Table  S8 ) and Scoary (Table  S7 ) analyses are available as supplementary tables. Genome clustering was performed on a large in-house global collection of 9535 assembled genomes, using PopPUNK  (Lees et al. 2019). A Neighbor-joining tree was generated using RapidNJ as implemented in PopPUNK  (Simonsen and Pedersen 2011). The phylogeny was annotated with MLST profiles generated in silico (https://github. com/tseemann/mlst).
17) and tip-date-randomization performed within BactDating demonstrated a clear temporal signal in the data (Figs. S6  and S7 ). 200 million Markov chain Monte Carlo (MCMC) steps were performed to generate a time-resolved tree using the default gamma model for clock rate. Post removal of the first 20 million steps, the resulting evolutionary rate was estimated to 6. 54 substitutions/genome/year (95% CI: 5. 76 -7. 34). To verify the temporal inferences from the above analyses, a separate workflow was employed: Gubbins  (Croucher et al. 2015)  was run on the Clade A portion of the dataset, followed by least squares dating using LSD2  (To et al. 2016)  with 1000 bootstrap replicates for temporal inference. The resulting TMRCA estimate for Clade A was 1838 (CI: 1810 to 1861), which was well aligned with our original BactDating-inferred TMRCA in 1887 (CI: 1798 to 1938).
Effective population size estimates To describe the growth and decline of the waves out of Asia we estimated their effective population size over time using the Skygrowth R-package v0. 2. 0  (Volz and Didelot 2018). Skygrowth uses a nonparametric Bayesian model to estimate the effective population size over time given the estimated genealogies. The model places a prior on the change in the growth rate of the effective pathogen population size, not on the change in the logarithm of the effective population size, like the commonly used skyride models  (Minin et al. 2008). For sections of the phylogenies where data are sparse this has the advantage of using the information on the growth rate learned from the previous parts of the phylogeny, whereas the skyride models has a tendency to estimate stabilizing effective population sizes in such situations.
Identification of resistance alleles To identify mutations and alleles associated with reduced susceptibility or resistance to relevant antibiotics, whole genome assemblies were aligned against gene-specific databases downloaded from pubMLST  (Jolley et al. 2018)  using BLAST  (Altschul et al. 1990 ). This was done for the following genes/segments: penA (NEIS1753), ponA (NG_ponA), porB (NG_porB), mtrR ('mtrR -includes promoter region) 23S (NG_23S), parC (NG_parC), gyrA (NG_gyrA), rplD (NEIS0133) and rpsJ (NEIS0129). In addition, to identify mosaic mtrD alleles associated with azithromycin resistance  (Wadsworth et al. 2018; Williamson et al. 2019) , the assemblies were blasted against the mtrD gene of FA1090, and the alignments visually inspected. For penA alleles, we used NG-STAR nomenclature  (Demczuk et al.
Phylogeography and mutational mapping Geographical transition patterns over time were inferred for the time-calibrated phylogeny using stochastic character mapping  (Bollback 2006)  as implemented in the R-package Phytools v0. 6-99  (Revell 2012). We simulated 1000 stochastic character mappings conditional on the observed geographical states on the tips. We used four continents as categories: "North America", "Europe", "Asia" and "Oceania", to capture large distance transitions that occurred in small time intervals. The posterior probability of the geographical state of a node was calculated as the fraction of the inferred character mappings for which that node was in a given character state (continent). The times of geographical transition events were considered to be the times when the transitions happened on the branches leading to the downstream nodes with posterior probabilities higher than 50% for new locations.
Downsampling for a balanced global representation To rule out that key spatio-temporal inferences were the result of unbalanced sampling across regions and phenotypic drug susceptibility, we generated a downsampled dataset and re-run temporal and phylogeographic analyses. The downsampled dataset included 90 genomes from each of the regions North America, Europe and East Asia, as well as all the genomes from Oceania (n=26), the Philippines, India and Chile (n=4 in total). The European and East Asian collections were downsampled in a random fashion, whereas the North American collection was downsampled in a targeted manner to minimize the effect of targeted sampling of resistant isolates from the USA  (Grad et al. 2014 ): We thus included all 46 phenotypically cefixime susceptible isolates in addition to 30 resistant isolates from the USA.
Acknowledgements The authors would like to thank our colleagues responsible for culturing and whole genome sequencing of gonococci at the  Norwegian Institute of Public Health. We would also like to thank colleagues from  The Scottish Bacterial STI reference lab  for facilitating access and provision of the ST1901 isolates and  Cecile Minder  for the propagation and DNA extraction from these strains. Magnus N. Osnes is supported by  The National Graduate School in Infection Biology and Antimicrobials (IBA)  hosted by the  University of Oslo  ( IBA  project  249062 ). L. v. D and F. B. acknowledge financial support from the  Newton Fund UK-China NSFC initiative  (grant  MR/P007597/1 ) and the  BBSRC  (equipment grant  BB/R01356X/1 ). No funding bodies were involved in the design of the study, the collection, analysis, and interpretation of data and in writing the manuscript.
5e ). Thus, the initial position of stem cells along the central-periodontoblastic axis is associated with its transcriptional state, migratory trajectory, and fates of progeny (Fig. 5c ). Comparisons of composition of growing vs. nongrowing mouse teeth. Although the mouse incisor stands as a model for a growing tooth, molecular features that distinguish it from nongrowing teeth remain unexplored. Therefore, we generated singlecell transcriptional snapshots of a nongrowing adult mouse molar using both 10X Chromium and Smart-seq2 platforms. To leverage total scale of multiple datasets, we analysed them jointly and together with self-renewing incisor datasets using Conos data integration strategy (Supplementary Fig. 4 )  40. Coarse-grained cell-type composition appeared similar between molar and incisor, except for the lack of epithelial populations in adult molars (Supplementary Fig.
In agreement with that, we observe that the Sox2 + /Shh + clusters 12 and 2 contain the majority of TAC's and most likely represent less differentiated states as compared to Sox2 -/Shh + clusters 11, 5, 1, and 13. h, i Expression of wellknown markers corresponding to a different ameloblast stage proving the gradual differentiation from secretory ameloblasts stage (Enam + ) through maturation ameloblast stage (Klk4 + , Odam + ) into postmaturation ameloblast stage (Gm17660 + ). significantly more homogeneous as compared to the pulpal populations of the incisor given the resolution of the current measurements. Joint Conos clustering of incisor and molar datasets shows that molar mesenchyme falls into a single cluster shared with the distal mouse-incisor pulp (Supplementary Fig. 4a, b ). Analysis of mesenchyme heterogeneity using separately 10× and Smart-seq2 platforms corroborated the heterogeneous population structure of mouse incisor and homogeneous distallike population of mouse molar (Supplementary Fig.
4d, e ). At the same time, gene expression programs of mouse molar and distal incisor pulp have noticeable expression differences in 379 genes (p value < 10 -2 , t-test group means comparison and at least two fold change in both Smart-seq2 and 10× Chromium datasets) (Supplementary Fig. 4f ). Mouse-incisor apical genes tend to show high expression in a Smoc2 + compared to Smoc2 -human apical papilla. On the other hand, mouse-incisor distal genes tend to show high expression in a Smoc2 -and not Smoc2 + human apical papilla. In adult teeth, mouse-incisor distal genes are uniformly In-depth single cell analysis of dental mesenchyme expressed in all populations, but incisor apical genes show the affinity to the periodontoblastic pulp. The meaning of this heterogeneity is unknown and requires further investigation.
8c, d ). In that regard, the pulp of human nongrowing molars appeared to be much more transcriptionally diverse compared to the mouse nongrowing molars. In particular, human molar contained a pulp subpopulation that was spatially localized in the periodontoblastic layer, previously morphologically described as cell-free and cell-rich zones, which are absent in mouse (Fig. 8i )  42. We detected a group of proliferative cells in a growing human apical papilla, which showed pronounced transcriptional similarity to a Smoc2 -human apical papilla pulp, and dissimilarity with any subpopulations of nongrowing human molars (Figs. 8e-g, 9e, g). To explore similarity of genetic programs in a mesenchymal compartment of human and mouse teeth, we compiled a set of marker genes that are differentially expressed between apical and distal mouse-incisor subtypes in both Chromium 10x and Smart-seq2 datasets (Supplementary Table  4 ).
Used antibodies: ACTA2 (Protein Tech, 23081-1-AP; 1:500), AIF1 (Novus, NB100-1028; 1:500), CALB1 (Swant; CB-38a; 1:500); COL4 (AbD Serotec, 2150-1470; 1:500), CDH1 (Novus, Af748; 1:500), CSF1 (NSJ, R31901; 1:200), CLDN10 (Sigma-Aldrich, HPA042348; 1:200), DLX5 (LSbio, LS-C352119; 1:200), DPP4 (Novus, AF954; 1:200), EGR1 (Cell Signalling, 4154; 1:200), F4/80 (Abcam, ab6640; 1:200), LYVE1 (Novus, AF2125; 1:200), MKI67 (Zytomed, RBK027-05, 1:200), NOTUM (Sigma-Aldrich, HPA023041; 1:200), PIEZO2 (Sigma-Aldrich, HPA040616; 1:200), POSTN (Novus, NBP1-30042; 1:200), RYR2 (ThermoFisher, PA5-36121; 1:200), S100A13 (DAKO; IS504; 1:500), SALL1 (Abcam, ab31526; 1:200), SMOC2 (MyBioSource, MBS2527784; 1:200), SOX9 (Sigma-Aldrich, HPA001758; 1:200), SOX10 (Santa cruz, sc-365692; 1:200), THBD (RnD systems, MAB3894; 1:200).
After centrifugation, supernatant was removed, the pellet was resuspended in 1 mL 2% FBS in HBSS, suspension was filtered using Tubes with Cell Strainer Snap Cap (Corning, 352235) and FACS (BD FACSAria III; software used: BD FACSDiva 8. 0. 1) was performed. All the work (except of enzymatic digestions at 37 °C) was performed on ice. Numbers of used mice for single-cell RNA-seq analyses. For analysis of adult healthy mouse-incisor 78 mandibular incisors were used out of 39 animals in total. For analysis of mouse molar pulps 48 first molars out of 12 adult animals were used in total. For adult human tooth analyses 7 wisdom molars out of 7 healthy randomly selected males and females of age 18-31 were used and 6 apical papillae out of 3 patients/teeth were used. Single-cell sorting and single-cell transcriptomics.
Following the final centrifugation, cells were resuspended in 500 ul of staining buffer and 1. 5 ug of DAPI (D1306, Invitrogen) added to be used as a dead cell exclusion marker. Samples were then analyzed on BD FACSAria III fusion machine. Data analysis was performed on FlowJo v10 software. Cells were gated based on size using standard SSC-A and FSC-A parameters so that debris is excluded. Following gating of cells, we focused on single cells and excluded doublets using SSC-A and SSC-W parameters. Live cells were then selected as cells identified to be dimly fluorescing in DAPI. Appropriate gating strategies were then used to select cells positive for the antibodies being used as deducted from the use of unstained controls. Data processing of mouse incisor. The reads were aligned to the UCSC mm10 genome assembly, and per-gene read counts in each cell were determined using feature Counts software.
The cells were grouped in 17 clusters using unbiased clustering as determined by PAGODA. t-SNE embedding was generated using Rtsne package and PAGODA-based cell-cell distance with perplexity = 25. Expression of a set of genes, where it is shown, was defined as their average expression for each cell. To characterize gene modules controlling cell-type identities, we selected genes that have at least 2 fpm difference between maximum and mean average expression among clusters. For supplementary heatmap of general dataset genes were clustered using hclust() and cutree() R functions in 20 clusters using hierarchical clustering with Euclidean distance between gene expression profiles using Ward's linkage and cells were arranged using PAGODA clustering described above. Similar procedure was used to characterize gene modules separately in epithelial (552 genes) and other compartments using different cutoffs of 1.
Five clusters representing progressive ameloblast differentiation were identified based on known markers of respective stages. Heterogeneity of ameloblasts cells was modelled as a principal trajectory using crestree R package approach  51  , with parameters (lambda = 100, sigma = 0. 03, M = 100) and cosine-based cell-cell distance. Root of the reconstructed trajectory was selected to biologically correspond to progenitor population of ameloblasts and each cell was assigned pseudotime as a distance from the root along the trajectory. Gene expression levels were modeled as a function of pseudotime using splines of the fifth degree with gam function from mgcv R package. Significance of association was calculated as Benjamini-Hochberg adjusted gam estimates of spline p value. Fitted gene expression levels were used to estimate magnitude of expression levels variation along the trajectory and downstream clustering.
Five hundred and fifty-six genes that had more than 100-fold magnitude differences along the trajectory and adjusted p value < 10 -5 were clustered in nine clusters using hierarchical clustering with Ward linkage based on Euclidean distance. Mesenchymal compartment. Four mesenchymal clusters of the general dataset comprising 1111 cells were reanalyzed separately. Only genes that had at least 20 reads in at least 10 cells were considered for downstream processing using PAGODA (k = 20 and plain batch correction across wells). Mitotic signature was regressed out and top aspects of transcriptional heterogeneity were determined as for general dataset. The cells were grouped in five clusters using unbiased clustering as determined by PAGODA and t-SNE embedding was generated using PAGODAbased cell-cell distance with perplexity = 20.
To clean up non-mesenchymal admixture of cells from other populations, only 1042 cells that had mean correlation of more than 0,2 to 200 the most correlated cells were retained for further analysis. Transcriptional states of mesenchymal cells were modelled as a principal tree using our crestree R package based on the SimplePPT approach  51, 52. Briefly, given a set of data points x_1,. ,x_N in M-dimensional space, a set of principal points z_1,. , z_K are arranged and are connected as a tree in the same space. Positions of principal points and tree structure are learned as alternate convex optimization problem that balances overall proximity of prinicipal tree to data points and stringency of the tree. Tree was learned with parameters (lambda = 2000, sig = 0. 03, M = 200) and cosine-based cell-cell distance.
Mesenchymal cells were grouped in 12 clusters using default PAGODA2 multilevel community detection method. Three follicle clusters were merged, while nine dental pulp clusters were coloured to reflect the most similar cluster colour of Smart-seq2 annotation (see Fig. 4a ). For analysis of fate-specific expression programs, 20 genes of each fate that have the largest mean expression difference between a fate cluster (e. g. apical, distal and pre-odontoblastic clusters) and progenitor clusters were considered as fate-specific markers. Intensity of fate-specific expression program in each cell was estimated as mean expression among 20 fate-specific genes. Cell-cycle score was defined as first principal component of cells transcriptional variability based on cell-cycle-annotated genes from. Cells from a cluster of mitotic cells were projected onto t-SNE embedding of non-mitotic mesenchymal landscape as a mean position of 10-nearest neighbours non-mitotic cells, where neighbours were defined using cosine-based distance in dimensionally reduced space of 20 PCs.
5 for glia, 1 for endothelium, 1 for pericytes)). 1-cor(. ) cell-cell distance was used for clustering and t-SNE (perplexity = 20). Assessment of cell quality. Cell quality was additionally probed using metrics reflecting expression complexity, mitochondrial content and doublet probabilities. Toward that goal, we assessed tradeoffs between number of expressed genes and UMIs (or reads for Smart-seq2) per cells, fraction of total reads from mitochondrial genes and cell doublet probabilities estimated using Scrublet with default parameters  53. We estimated and explored these metrics for four datasets and manually excluded one cluster that had low number of genes compared to UMIs (it was a cluster of spike-in cells, see below) and a number of clusters of joint human analysis that were likely doublets (see Extended Data Fig.
2e-h ). Data preprocessing of 10x Chromium samples. CellRanger-10× Chromium software was used to perform alignment to GRCh38 human genome or mm10 mouse genome assemblies, filtering, barcode counting and UMI counting. For Apical papilla 1, Adult molar 3, Adult molar 4, incisor (10×), mouse molar 1 (10×) datasets preprocessing was performed using CellRanger-2. 2. 0 following by filtering of cells having less than 500 UMIs. For other datasets datasets preprocessing was performed using CellRanger 3. 0. 2 following by default CellRanger 3. 0. 2 filtering of cells. Additionally, a protocol of library preparation used by the facility included spike-in of Jurkat and 32D cells of human and mouse species. Spike-in cells were not used for data processing or analysis and were excluded as Hbb + clusters; they are also easily detectable as having low complexity and forming a separate outlier transcriptional cluster.
Joint analysis of mouse datasets. Two mouse-incisor datasets (10x Chromium, 4236 cells, and Smart-seq2,2889 cells) and three mouse molar datasets (two 10× Chromium, 1460 and 384 cells, and Smart-seq2, 195 cells), each composed of multiple teeth (see chapter "Numbers of used mice for single-cell RNA-seq analyses" in materials and methods and Supplementary Table  2 ), were processed independently using PAGODA2 R package 13 routine basicP2proc(), which performs normalization, log transformation, correction for mean-variance trend of expression levels, dimensionality reduction via PCA and clustering. After filtration of spike-in clusters identified through expression of Hbb, processed datasets were then jointly analysed using CONOS R package, which enables integrative analysis of single-cell datasets across samples and conditions  40.
Joint graph of 9164 cells from all datasets was constructed using CONOS routine buildGraph() with nearest neighbour parameters k = 15, k. self = 15, k. self. weight = 0. 1 in space of 10 common principal components (CPCA) estimated using 1000 overdispersed genes for each pair of samples. Joint graph was layout in 2D using UMAP method through CONOS routine embedGraph() with parameters spread = 1 and min. dist = 0. 05. Graph-based leiden community method with resolution = 1. 0 was used to partition cells in 22 clusters. To provide additional statistical support for reproducible structure of mesenchymal incisor populations and homogeneous distal-like molar state, we computationally isolated clusters of mesenchymal cells and separately explored cells of Smart-seq2 and 10x Chromium platforms. Platform-specific mesenchymal cells were processed using routine basicP2proc() with default PAGODA2 batch correction across samples and dimensionality reduction to 10 principal components based on top 1000 overdispersed genes.