dataesr/scientific-paragraphs-categorization huggingface dataset:
Introduction Recent technological improvements in microfluidics and low volume sample handling  [1] ha v eenabled the emergence of single-cell transcriptomics  [2, 3]  as a popular tool for analyzing biological systems  [4, 5, 6]. This growing popularity along with a continued increase in the scale of the respective assays  [7]  has resulted in massive amounts of publicly available data and motivated large scale community efforts such as the Human Cell Atlas  [8] , Tabula Muris  [9]  and the BRAIN Initiative Cell Census Network  [10]. The next natural step in the evolution of this field is therefore to integrate many available datasets from related tissues or disease models in order to increase statistical robustness  [11] , achieve consistency and reproducibility among studies  [12, 13] , and ultimately converge to a common ontology of cell states and types  [8, 14].
A fundamental step toward the ideal of a common ontology is data harmonization, namely integration of two or more transcriptomics datasets into a single dataset on which any downstream analysis can be applied. We use the term harmonization rather than batch effect correction in order to emphasize that the input datasets may come from very different sources (e. g. ,t e c hnology, laboratory), and from samples with a different composition of cell types. A wide range of methods have already been developed for this fundamental problem, initially for Microarrays and later on for bulk RNA sequencing, such as ComBat  [15]  and limma  [16]  which rely on generalized linear models with empirical Bayes shrinkage to avoid over-correction. More recently, similar methods have been proposed specifically for single-cell RNA sequencing (scRNA-seq), such as ZINB-WaVE  [17] , which explicitly accounts for the overabundance of zero entries in the data.
However, because of their linear assumptions, these approaches may not be appropriate when provided with a heterogeneous sample that includes different cell states, each of which may be associated with a different sample-to-sample bias  [12]. With these limitations in mind, the next generation of methods turned to non-linear strategies. Broadly speaking, each of these methods includes a combination of two components: (i) joint factorization of the input matrices (each corresponding to a different dataset) to learn a joint low-dimensional latent representation. This is usually done with well established numerical methods, such as integrative non-negative matrix factorization (LIGER  [18] ), singular value decomposition (Scanorama  [19] ), or canonical correlation analysis (Seurat Alignment  [13] ); (ii) additional non-linear transformation of the resulting latent representations so as to optimally "align" them onto each other.
This is usually done using heuristics, such as alignment of mutual nearest neighbors (MNN  [12] , Scanorama  [19]  and Seurat Anchors  [20] ), dynamic time warping (Seurat Alignment  [13] ) or quantile normalization (LIGER  [18] ). While this family of methods has been shown to effectively overlay different datasets, it suffers from two important limitations. First, an explicit alignment procedure may be difficult to tune in a principled manner and consequently result in over-normalization, especially in challenging cases where the cell type composition is different between datasets and when technical differences between samples are confounded with biological differences of interest. Second, the alignment is done in an ad hoc manner and lacks probabilistic interpretability. Consequently, the resulting harmonized dataset is of limited use and cannot be directly applied for differential expression and other probabilistic decision-making tasks.
Another recent line of work makes use of neural networks to learn a joint representation of multiple datasets (SAUCIE  [21] ) or project one dataset into another (maximum mean discrepancy [MMD] ResNet  [22] ). These methods rely on an explicit non-parametric measure of discrepancy between probability distributions (MMD) to match either the latent spaces or directly the gene expression values from pairs of datasets. However, using the MMD with a universal kernel explicitly assumes that the cell type proportion is similar in all the datasets, which may be less suitable in the general case of data harmonization. Besides harmonization, another important and highly related problem is that of automated annotation of cell state. In principle, there are two ways to approach this problem. The first is ab initio labeling of cells based on marker genes or gene signatures  [13, 23, 24].
While this approach is intuitive and straightforward, its performance may be affected in the plausible case where marker genes are absent due to limitations in sensitivity. The second approach is to "transfer" annotations between datasets. In the simplest scenario, we have access to one dataset where states have been annotated either ab initio, or using additional experimental measurements (e. g. , protein expression  [3, 25]  or lineage tracing  [26] ) and another, unannotated dataset from a similar condition or tissue. The goal is to use the labeled data to derive similar annotations for the second dataset, whenever applicable. This task is often complicated by factors such as differences in technology (e. g. , using Smart-Seq2 data to annotate 10x Chromium data), partial overlap in cell type composition (i.
e. , not all labels should be transferred and not all unannotated cells should be assigned a label), complex organization of the labels (e. g. , hierarchy of cell types and sub-types  [27] , continuum along phenotypic or temporal gradients), partial labeling (i. e. , only a subset of cells from the "annotated" dataset can be assigned a label confidently), and the need to handle multiple (more than 2) datasets in a principled and scalable manner. One way to address the annotation problem with this approach is learning a classifier  [27, 28]  in order to predict a fixed stratification of cells. However, this approach might be sensitive to batch effects, which could render a classifier based on a reference dataset less generalizable to an unannotated dataset. Another, more flexible approach is to transfer annotations by first harmonizing the annotated and unannotated datasets, thus also gaining from the benefits of having a single dataset that can be subject to additional, joint, downstream analysis.
In this paper, we propose a strategy to address several of the outstanding hurdles in both of the harmonization and annotation problems. We first demonstrate that single-cell Variational Inference (scVI)  [29]  a deep generative model we previously developed for probabilistic representation of scRNA-seq data -performs well in both harmonization and harmonizationbased annotation, going beyond its previously demonstrated capacity to correct batch effects. We then introduce single-cell ANnotation using Variational Inference (scANVI), a new method that extends scVI and provides a principled way to address the annotation problem probabilistically while leveraging any available label information. scANVI uses a semi-supervised generative model, which can be utilized for both approaches to the annotation problem.
In the first scenario, we are concerned with a single dataset in which only a subset of cells can be confidently labeled (e. g. , based on expression of marker genes) and annotations should then be transferred to other cells, when applicable. In the second scenario, annotated datasets are harmonized with unannotated datasets and then used to assign labels to the unannotated cells. The inference procedure for both of the scVI and scANVI models relies on neural networks, stochastic optimization and variational inference  [30, 31]  and scales to large numbers of cells and datasets. Furthermore, both methods provide a complete probabilistic representation of the data, which non-linearly controls not only for sample-to-sample bias but also for other technical factors of variation such as over-dispersion, library size discrepancies and zero-inflation.
Results In the following we demonstrate that our framework compares favorably to the existing state of the art in the harmonization and annotation problems in terms of accuracy, scalability, and adaptability to various settings. The first part of the paper focuses on the harmonization problem and covers a range of scenarios, including harmonization of datasets with varying levels of biological overlap, handling cases where the data is governed by a continuous (e. g. , pseudotime) rather than discrete (cell types) form of variation, and processing multiple (> 20) datasets. While we demonstrate that scVI performs well in these scenarios, we also demonstrate that the latent space leaned by scANVI provides a proper harmonized representation of the input datasetsa property necessary for guaranteeing its performance in the annotation problem.
In the second part of this manuscript we turn to the annotation problem and study its two main settings, namely transferring labels between datasets and ab-inito labeling. In the first setting we consider the cases of datasets with a complete or partial biological overlap and use both experimentallyand computationally-derived labels to evaluate our performance. In the second settings, we demonstrate how scANVI can be used effectively to annotate a single dataset by propagating high confidence seed labels (based on marker genes) and by leveraging a hierarchical structure of cell state annotations. Finally, we demonstrate that the generative models inferred by scANVI and scVI can be directly applied for hypotheses testing, using differential expression as a case study.
Joint modeling of scRNA-seq datasets We consider K different scRNA-seq datasets (Figure  1a ). A dataset indexed by k 2{1, ••• ,K} consists of a N k ⇥ G k matrix where each entry records the number of transcripts observed for each of G k genes in each of N k cells. In this work, we use a standard heuristic to filter the genes and generate a common (possibly large) gene set of size G (Online Methods), and thus format our data as a unique ( P k N k ) ⇥ G matrix whose individual entries are noted x ng -namely the expression of gene g in cell n. W eu s es n 2{ 1, ••• ,K} to represent the dataset from which each cell n was generated. Furthermore, a subset of the cells may be associated with a cell state annotation c n 2 C, which can describe either discrete cell types or hierarchical cell types.
More complex structures over labels such as gradients are left as a future research direction. Since the problem of data harmonization of single-cell transcriptomics is difficult and can potentially lead to over-correction (Supplementary Figure  1 )  [32] , we propose a fully-generative method as a robust and principled approach to address it. In our previous work  [29] , we built single-cell Variational Inference (scVI), a deep generative model where the expression level x ng is zero-inflated negative binomial (ZINB) when conditioned on the dataset identifier (s n ), and two additional latent random variables. The first, which we denote by l n , is a one-dimensional Gaussian accounting for the variation in capture efficiency and sequencing depth. The second is z n is a low dimensional Gaussian vector that represents the remaining variability (Figure  1b ).
This vector is expected to reflect biological differences between cells, and can be effectively used for visualization, clustering, pseudotime inference and other tasks. Since the scVI model explicitly conditions on the dataset identifier, it provides an effective way of controlling for technical sample-to-sample variability. However, scVI is unsupervised and does not make use of the available annotations c n , which can further guide the inference of an informative latent representation z n. To this end, we present a more refined hierarchical structure for z n. W ed r a wz n as a mixture conditioned on the cell annotation c n and another latent variable u n , accounting for further biological variability within a cell type (Online Methods). We name the resulting approach single-cell ANnotation using Variational Inference (scANVI).
Notably, our work has substantial overlap with machine learning research such as semi-supervised learning in the context of variational-autoencoders  [33, 34]  domain adaptation  [35]  and label propagation  [36]. We review these related research topics in Supplementary Note 1. The variables z n , inferred either with scVI or scANVI, provide an embedding of all cells in a single, joint latent space. Since this latent space is inferred while controlling for the dataset of origin (s n ), it inherently provides a way to address the harmonization problem. The annotation of unlabeled cells can therefore be conducted with scVI using their proximity to annotated cells in the joint latent space (e. g. , using majority vote over the k-nearest neighbors). The scANVI model provides a more principled way to annotate cells, namely through a Bayesian semi-supervised approach.
Once fitted, the model is able to provide posterior estimates for the unobserved cell state c n , which can be particularly useful when labels cannot be entirely trusted. Because the marginal distribution p(x ng ,c n | s n )i fc n observed (resp. p(x ng | s n ) otherwise) is not amenable to exact Bayesian computation, we use variational inference parameterized by neural networks to approximate it  [30] ( Online Methods). Notably, scANVI and scVI both have a certain number of hyperparameters. In the following evaluations, conducted on different datasets and different scenarios, we use the exact same set of hyperparameters in order to demonstrate that our methods can be applied with a minimal requirement of hyperparameter tuning (Online Methods). We provide a robustness study for hyperparameters in the context of harmonization in Supplementary Figure  2.
Datasets We apply our method on datasets generated by a range of technologies (10x Chromium  [3, 37] , plate-based Smart-Seq2  [38] , Fluidigm C1  [39] , MARS-Seq  [40] , inDrop  [41]  and CEL-Seq2  [42] ), spanning different numbers of cells (from a few thousand to over a hundred thousand cells), and originating from various tissues (mouse bone marrow, human peripheral mononuclear blood cells [PBMCs], human pancreas, mouse brain). Datasets are listed and referenced in Supplementary Table  1. Harmonizing pairs of datasets with a discrete population structure We conducted a comparative study of harmonization algorithms on four different instances, each consisting of a pair of datasets. The first pair (PBMC-CITE  [25] , PBMC-8K  [37] ) represents the simplest case, in which the two datasets come from very similar biological settings (i.
We explore three such cases, including a bone marrow 10x and Smart-Seq2 pair from the Tabula Muris project (MarrowTM-10x, MarrowTM-ss2  [9] ), a pancreas inDrop and CEL-Seq2 pair (Pancreas-InDrop, Pancreas-CEL-Seq2  [44] ), and a dentate gyrus 10x and Fluidigm C1 pair (DentateGyrus-10x, DentateGyrus-C1  [45] ). Successful harmonization should satisfy two somewhat opposing criteria (Supplementary Figure  1 ). On the one hand, cells from the different datasets should be well mixed; namely, the set of k-nearest neighbors (kNN) around any given cell (computed e. g. , using euclidean distance in the harmonized latent space) should be balanced across the different datasets. For a fixed value of k, this property can be evaluated using the entropy of batch mixing  [12] , which is akin to evaluating a simple k-nearest neighbors classifier for the batch identifier (Online Methods).
While this property is important, it is not sufficient, since it can be achieved by simply randomizing the data. Therefore, in our evaluations we also consider the extent to which the harmonized data retains the original structure observed at each dataset in isolation. Here, we expect that the set of k-nearest neighbors of any given cell in its original dataset should remain sufficiently close to that cell after harmonization. This property can be evaluated using a measure we call k-nearest neighbors purity (Online Methods). Furthermore, because the conservation of k-nearest neighbors might be more indicative of a local stability of the algorithm and not of actually preserving global aspects of the datasets, we introduced an additional metric, evaluating the conservation of cluster assignments before and after harmonization.
Clearly, these two criteria can be simply optimized by not changing the input datasets, which will result in poor performance with respect to our first measure. Our evaluation therefore relies on both types of measures, namely mixing of data sets and retainment of the original structure. Since our results depend on the neighborhood size k, we consider a range of values -from a high resolution (k = 10) to a coarse (k = 500) view of the data. We compare scVI to several methods, including MNN  [12] , Seurat Alignment  [13] , Com-Bat  [15]  and principal component analysis (PCA). In addition, and in order to compare our methods to unpaired data integration approaches based on generative adversarial networks  [46] , we also tested MAGAN  [47]. However, even after manual tuning of the learning rate hyperparameter, the input datasets remain largely unmixed (Supplementary Figure  3 ).
This might be due to the fact that MAGAN was not directly applied to harmonize pairs of scRNA-seq datasets and need more tuning to be applicable in that context. For each algorithm and pair of datasets, we report embeddings computed via a Uniform Manifold Approximation and Projection (UMAP)  [48] ( Figure  2a , Supplementary Figure  4 5 6 7 ) as well as the three evaluation metrics (Figure  2b -c, Supplementary Table  2 ). Overall, we observed that scVI performs well in terms of mixing, while comparing favorably to the other methods in terms of retainment of the original structure, for a wide range of neighborhood sizes and across all dataset pairs. Re-assuringly, our positive results for preservation of the output of a clustering algorithm indicates that scVI and scANVI are also stable with regards to more general aspects of the data.
As an illustrative example, consider the results of applying scVI and Seurat Alignment to the Tabula Muris bone marrow datasets (Figure  2a ). While scANVI was designed for the problem of cell state annotation, we also wanted to evaluate its ability to harmonize datasets, which can be seen as a prerequisite. To evaluate this, we consider each dataset pair twice, each time using labels from one of the datasets (exploiting the semi-supervision framework of scANVI) and leaving the other one unlabeled. Reassuringly, we found that scANVI is capable of effectively harmonizing the datasets, with a similar performance to that of scVI in terms of entropy of batch mixing and retainment of the original structure (Figure  2b-c , and Supplementary Table  2 ). We further explore the performance of scANVI in the annotation problem in the subsequent sections.
Harmonizing datasets with a different composition of cell types One of the primary challenges of the harmonization problem is handling cases in which the cell types present in the input datasets only partially overlap or do no overlap at all. Since this is a plausible scenario in many applications, it is important to account for it and avoid over-normalizing or "forcing" distinct cell populations onto each other. To evaluate this, we performed several stress tests in which we artificially manipulated the composition of cell types in the input datasets prior to harmonization. As our benchmark method we use Seurat Alignment, which performed better than the remaining benchmark methods in our first round of evaluation (Figure  2 ). As a case study, we used a pair of PBMC datasets (PBMC-CITE  [25] , PBMC-8K  [37] ) that initially contained a similar composition of immune cell types (Supplementary Table  3 ).
We were first interested in the case of no biological overlap (Figure  3a-d ). To test this, for a given cell type c 0 (e. g. , natural killer cells), we only keep cells of this type in the PBMC-CITE dataset and remove all cells of this type from the PBMC-8K dataset. In Figure  3a -b, we show an example of UMAP visualization of the harmonized data, with natural killer cells as the left out cell type c 0. Evidently, when harmonizing the two perturbed datasets with scVI, the natural killer cells appear as a separate cluster and are not wrongly mixed with cells of different types from the other dataset. Conversely, we see a larger extent of mixing in the latent space inferred by Seurat Alignment. A more formal evaluation is provided in Figure  3c-d , which presents our harmonization performance metrics for each cell type averaged across all perturbations (in each perturbation, c 0 is set to a different cell type).
We also included scANVI with the true number of cell types (C = 6) in this analysis, using the cell labels from the PBMC-CITE dataset. Under the ideal scenario of a successful harmonization, we expect both a low entropy of batch mixing (since the datasets do not overlap), and retainment of the original structure. Evidently, both scVI and scANVI exhibit a consistently low level of batch mixing that is better or comparable to that of Seurat Alignment, while retaining the original structure more accurately. As an additional scenario, we investigated the case where the input datasets contain a similar set of cell types, with the exception of one cell type that appears in only one of the datasets. To simulate this, for a given cell type c 0 , we removed cells of this type from the PBMC-8K dataset, and then harmonize the remaining cells with the unaltered PBMC-CITE (which still contains c 0 ).
We show an example of UMAP visualization in Figure  3e -f, removing CD4+ T cells from the PBMC-8K dataset. Evidently, in the scVI latent space, the PBMC-CITE "unique" CD4+ T cell population is not wrongly mixed with cells from the perturbed PBMC-8K dataset, but rather appears as a distinct cluster. For a more formal analysis, Figure  3g -i shows the harmonization statistics for perturbing the six major cell types present in the PBMC datasets. As above, we also evaluated scANVI in this context, using the labels from the unperturbed (PBMC-CITE) dataset. Figure  3g  shows that the entropy of batch mixing from the "unique" populations (averaging over all six perturbations) is low in all three methods (scVI, scANVI and Seurat Alignment), with a slight advantage for scVI and scANVI. Figure  3h -i shows the harmonization statistics for each population, averaging over all shared cell types between the two datasets.
Evidently, for the populations that are indeed common to the two datasets, scVI and scANVI are capable of mixing them properly, while preserving the original structure, comparing favorably to Seurat Alignment on both measures. Overall, the results of this analysis demonstrate that scVI and scANVI are capable of harmonizing datasets with very different compositions, while not forcing erroneous mixing. These results are consistent with the design of scVI and scANVI, which aim to maximize the likelihood of a joint generative model, without making ap r i o r iassumptions about the similarity in the composition of the input datasets. In a similar but more complex experiment, we also study the case when the two datasets both have their own unique cell types but also share several common cell types.
Populations unique to each dataset have low mixing (Supplementary Figure  8a ), especially with scVI and scANVI. Conversely, the shared populations have a substantially higher mixing rate (Supplementary Figure  8c ). Specifically, scANVI and scVI both mixes shared populations better than Seurat, with a better overall performance for scANVI. Finally, the preservation of original structure is higher scVI and scANVI when compared to Seurat across all cell types, especially for B cells, NK cells and FCGR3A+ Monocytes (Supplementary Figure  8b ). Overall, these results demonstrate that our methods do not tend to force wrong alignment of non-overlapping parts of the input datasets.