dataesr/scientific-paragraphs-categorization huggingface dataset:
Mapping score clearly indicate the reference cells that bear highest similarity to the projected cells. It was found that in the case of projected astrocytes generated 94. 05% of all mapping in the astrocyte cluster of the reference dataset. Similarly, oligodendrocyte projection generated 99. 49% of all mapping scores in the oligodendrocyte containing cluster of the reference database. Taken together, Scarf presents multiple methods to perform reference-based integration in highly memory efficient way. Scarf allows caching and reusage of training indices for KNN to allow rapid mapping of datasets in cold start scenarios. Using label transfer, co-embedding and mapping scores users can quickly identify relationships between cells across atlas scale datasets.
Memory efficient hierarchical clustering and t-stochastic neighbour embedding for single-cell genomics data Visualization of single-cells on two/three-dimensional space is one of the central tenets of single-cell data analysis. Since single-cell datasets contain tens to hundreds of thousands of features, non-linear dimension reduction algorithms have been considered as ideal solutions for 2D/3D data visualization. Historically, the most popular choice for a variety of single-cell datasets has been t-SNE  42  , while UMAP  14  is the most common current method for visualization on account of its computation efficiency (runtimes) and improved preservation of global structure in the data. Alternatively, FIt-SNE  43  is an improved t-SNE algorithm with runtimes on par with UMAP for larger scale datasets.
Recently, a novel t-SNE algorithm, SG-tSNE  15  , that can directly operate on stochastic KNN graphs was introduced. The implementation of SG-tSNE was shown to have improved runtimes over FIt-SNE and, to be scalable to three-dimensional embedding of large-scale datasets. However, it is not known how the embeddings obtained using SG-tSNE compared to those from UMAP and how SG-tSNE performs on a wider selection of datasets. The neighbourhood graph computed by Scarf can be directly fed into either UMAP or SG-tSNE algorithm. Furthermore, Scarf uses the same initial embedding for both UMAP and SG-tSNE (see Methods), thus restricting the differences between the two methods to the underlying algorithm itself and the parameters used. We benchmarked the run-time of UMAP and SG-tSNE on four atlas scale datasets of 600K, 1M, 2M and 4M cells at multiple iterations (see Methods for details).
Both algorithms run in an iterative manner for a user defined number of steps. Figure  4A  shows 2D UMAP (250 iterations) and SG-tSNE (2,500 iterations) plots obtained of the four atlas scale datasets. SG-tSNE consumed less time for 2500 iterations compared to UMAP's 250 (7. 0, 11. 8, 15. 9, and 29. 5 mins vs. 14. 9, 30. 3, 37. 7, and 96. 2 mins, for the 600K, 1M, 2M, and 4M cell datasets respectively). Moreover, we observed that for the corresponding iterations across datasets, SG-tSNE scaled better to the increasing number of cells; for 4M cell datasets the runtime of 10,000 iterations SG-tSNE (96. 45 minutes) was similar to the runtime of UMAP with 250 iterations (figure  4B ). Annotation of UMAP (supplementary figure  S4A ) and SG-tSNE (supplementary figure  S4B ) plots with Leiden clusters resulted in well separated clustering in both cases and the relative distances between clusters could be easily visualized on the UMAP plot.
Moreover, the clusters were positioned in an equivalent vicinity to other clusters in both embeddings. To quantify this, we calculated the Spearman rank correlation between the nearest neighbour across the clusters in the embedded space and the original cell-cell neighbourhood graph (supplementary figure  S4C ). The Spearman correlations coefficient values across the clusters show that UMAP was able to learn cluster relations quickly and showed only a slight improvement from increased number of iterations. The largest improvement was seen in the case of the 4M cell dataset where the mean Spearman's rho value increased from 0. 78 (100 iterations) to 0. 84 (1,000 iterations). In all four datasets, SG-tSNE showed a smaller mean Spearman's rho value than UMAP at lower number of iterations but improved to be larger than UMAP when higher number of iterations were used.
For example, in the case of the 4M cell dataset, the SG-tSNE's value increased from 0. 63 (1,000 iterations) to 0. 93 (10,000 iterations). Next, we compared UMAP and SG-tSNE for the preservation of local data structure by quantifying the number of KNN neighbours found in the vicinity of each cell in the UMAP or SG-tSNE space (supplementary figure  S4D ). Across the four datasets, SG-tSNE led to improved local structure with the increasing number of iterations while no benefit was observed for UMAP embeddings. UMAP, with 1,000 iterations, had 0. 14, 0. 19, 0. 21, and 0. 24 mean fraction of nearest neighbours preserved in the embedding of 600K, 1M, 2M, and 4M cell datasets, respectively. SG-tSNE, on the other hand, at 10,000 iterations, reached the values of 0. 33, 0. 46, 0. 44, and 0. 40, in the respective datasets.
Clustering is one of the most critical aspects of single-cell data analysis wherein cells are partitioned into distinct groups. Hierarchical clustering-based methods  [44] [45] [46] [47] [48]  can provide accurate partitioning of cells and discover rare cell types in the data but are not scalable to atlas size datasets  49. Alternatively, Graph clustering  [50] [51] [52] [53]  methods have become a popular choice for clustering on single-cell data and Louvain  52  and Leiden 50 algorithms are routinely used through Seurat  29  and Scanpy 16 packages. However, unlike hierarchical clustering methods, Louvain and Leiden are unable to clearly reflect the relationship between clusters. Scarf performs clustering using the highly scalable Leiden algorithm. Like other graph-based methods and unlike hierarchical clustering, Leiden is unable to explicitly indicate the relationships between clusters.
PAGA  54  was introduced as a method that can explicitly indicate relationships between clusters once created by any clustering method. Unfortunately, it can become hard to interpret those relationships on larger datasets where multiple clusters with complex connection patterns are present (figure  4C ). Within Scarf, we introduce a new approach, Paris 13 , that has not previously been applied to single-cell genomics data. Paris is a graph clustering approach that rather than creating flat clusters of cells, creates a dendrogram of cells akin to hierarchical clustering methods. Specifically, Paris computes a binary tree representing the various nested clusters of the graph, at different levels of granularity. The top-level clusters give the core structure of the graph while the lower-level clusters provide the local structure.
Like Louvain and Leiden, Paris is a scalable algorithm applicable to large graphs with millions of cells. High performance is achieved through the nearest neighbour chain algorithm, where pairs of clusters to be merged are locally searched, following the edges of the graph. The successive merged clusters are stored and aggregated to compute the final binary tree. We applied Paris to the six datasets used previously (supplementary figure  S4E ). In order to allow comparison with Leiden clustering (supplementary figure  S4F ), we cut the Paris dendrogram to obtain the same number of clusters as Leiden clustering. Coalesced Paris dendrograms presented a clear structure of cluster relationships (figure  4D ) and the mean concordance between the Paris and Leiden clusters was found to be between 88.
05% and 93. 73% across the datasets (supplementary figure  S4G ). These values were significantly larger (p-value <4. 25e-7; Mann Whitney U test) than randomly shuffled cluster identity wherein mean concordance with Leiden ranged between 9. 92% and 26. 47%. In order to ascertain the biological meaningfulness of Paris dendrograms, we used scRNA-Seq data from cells of the murine nervous system  39. The authors manually curated multiple taxonomy of cell types in the cell atlas. Importantly, Paris was able to retrace the taxonomic relationships between the cell types (figure  4E ) presenting a dendrogram with relevant and clearly separated branches for vascular, immune, glial, and neuronal cell types. Together, these results show that Scarf can perform all critical analysis steps of single-cell genomics data including filtering, normalization, feature selection, linear dimension reduction, neighbourhood graph creation, embedding using UMAP and t-SNE, clustering, downsampling, and cell projection on even the largest available data sets on a regular laptop in a time efficient manner.
DISCUSSION The ability to handle and process large-scale single-cell genomics data with readily available hardware is urgently needed. Most atlas-scale projects resort to providing online interfaces with UMAP/t-SNE plots over which users visualize the expression of different genes. But this can be grossly insufficient for many research tasks. Researchers often generate datasets that need co-analysis with other atlas-scale datasets which necessitates the use of highperformance computing machines. Another common-usage scenario where researchers need to process count matrices of large-scale datasets is when performing custom sub-selection of cells and generating new UMAP/t-SNE and clustering. Feature selection on a sub-selection of cells can improve the resolution of UMAP/t-SNE and increase the sensitivity of rare cell-type and cell-state identification.
We designed Scarf to allow researchers the possibility for analysis and re-analysis of atlas-scale datasets on their laptop computers. Coupled with memory efficient methods, like Kallisto Bustools  55  , for generating cell-gene count matrices, Scarf represents an end-to-end solution for the analysis of single-cell RNA-Seq datasets. During the preparation of this manuscript an R-based memory efficient tool, ArchR  56  , was published for analysis of scATAC-Seq data. In the paper describing this feature-rich tool, the authors benchmarked a simulated PBMC dataset with 1. 2 million cells and found the memory usage to be over 20GB in a small infrastructure setting. With Scarf, we were able to perform corresponding steps in scATAC-Seq dataset with 700K cells within 5GB RAM. Additionally, as of the writing of this manuscript, ArchR needs genome aligned fragment files and unlike Scarf, does not support direct input of precalculated cell-peak matrices that are readily available for many atlas scale datasets.
Just like Scarf has an upstream dependency on the cell-peak matrix, ArchR is dependent on external tools for genome alignment, a step that might not be memory efficient or feasible on laptop computers due to long runtimes and space requirements. Hence, neither of the currently available software yet represents an end-to-end solution for scATAC-Seq analysis. Our novel downsampling method is designed to leverage the fact that a cell-cell neighbourhood graph is already calculated in a memory efficient way. The obtained downsampled set can thereby be exported and used with external tools to apply methods that are not currently directly supported in Scarf. We allow seamless import and export of data between Scarf and Scanpy. Though our comparison metrics indicate that Scarf performs better than GeoSketch, we want to highlight that GeoSketch does not require clustering information that Scarf/TopACeDo needs to perform downsampling.
As previously described  57  , the embedding obtained from t-SNE and UMAP can be sensitive to the parameters used in the analysis. Hence, one method is not necessarily better than the other, but rather complementary techniques for visualization. Without extensive parameter tuning, SG-tSNE can generally reveal the heterogeneity of large-scale datasets more clearly than UMAP, while on the other hand UMAP can be more useful to investigate the relationship between clusters and visualize biological processes like differentiation. Our contribution here was to, for the first time, perform a head-to-head comparison of UMAP and a graph-based t-SNE (SG-tSNE) on multiple datasets by providing the exact same KNN graph and coordinate initialization. In Scarf, we introduce a generic graph-based hierarchical clustering method that can clearly highlight cluster similarity.
Single-cell specific hierarchical clustering algorithms have been shown to not be scalable to large-scale datasets  48. Our benchmarks showed that Paris algorithm applied on Scarf computed cell-cell neighbourhood graph creates a full dendrogram of a 4M cell dataset in 58 minutes. (supplementary figure  S4H-I ). This runtime is much slower compared to Leiden clustering which generated a fixed resolution in less than 6 minutes. However, the two runtimes cannot directly be compared as Paris reveals a full profile of cell-cell relationships while Leiden provides only a static partitioning using an arbitrary 'resolution' parameter. In practice, analysts often run the Leiden algorithm with many different values for the 'resolution' parameter to ascertain an optimal number of clusters. In contrast, the creation of a cell dendrogram in Paris is parameter free and can subsequently be cut into any desired number of clusters within a few seconds.
Efficient usage of CPU resources: We observed that in the Python ecosystem many libraries aggressively use all the available CPU cores. Though this is usually not an issue when running analysis on servers, it can seriously impede user's ability to simultaneously use their laptops/desktops if all CPU cores are being used. Hence, we have extensively placed CPU usage restriction throughout the code so that only a user-determined number of CPU cores are engaged into the analysis. We have further enabled support for running non-linear dimension reduction methods, which can often be longest running step in the pipeline, on multiple CPU cores. The Scarf code was written to facilitate an object-oriented programming experience for the users. As a result, the user needs to interact with only one variable for each dataset.
Multi-assay support Multi-omics experiments like CITE-Seq are well supported by Scarf. Scarf can store assays data from each of the modalities under a separate branch of the Zarr hierarchy. The cell level attributes are stored outside the assays while feature level attributes and assay's count matrix is stored under the assay group. Scarf's programmatic interface mimics this hierarchy. Also, Scarf assigns one of the assays as default when the dataset is loaded, and all the functions will act on this default assay. Almost all the functions associated with `DataStore` object usefrom_assay` parameter which can be used to assign the assay for the operation.
Normalization For scRNA-Seq datasets Scarf performs library size normalization of the cells. Hence, the normalized value of a gene/feature in a cell ‫ݕ‬ ி can be calculated as ‫ݕ‬ ൌ ܵ. ‫ݔ‬ ∑ ‫ݔ‬ ᇲ ᇲ ி where F is the feature set, c represents a cell and S is a scaling factor. The normalized value can optionally be transformed into log scale: ‫ݕ‬ ൌ ݈‫݃‬൫1 ‫ݕ‬ ൯ For scATAC-Seq datasets, TF-IDF (term frequency-inverse document frequency) normalization is performed.
‫ݕ‬ ൌ ‫ݔ‬ ݊ ி. ‫݈݃‬ ቆ1 ܰ ݊ ቇ Where ݊ ி represents the total number of accessible peaks (those with non-zero values) in a given cell, ݊ represents total number of cells where a given peak is accessible. ܰ is total number of cells in the dataset.
Feature selection: For scRNA-Seq datasets, Scarf provides highly variable gene (HVG) selection approach as previously reported  59. For this purpose, the mean expression and variance of genes are calculated across all the cells and converted into log scale. The next step is to remove the mean-variance trend in the feature space. For this, genes are divided into equal sized bins based on their mean expression value; from each bin, the gene with lowest expression value is selected. A lowess curve (using Python's statsmodels package  60  ) is fitted to the selected genes. The fitted curve is used to predict the 'expected variance' for each gene using mean expression as independent variable. The 'expected variance' is subtracted from observed variance to obtain corrected variance based on which HVGs are selected.
Primary dimension reduction: For single-cell RNA-Seq datasets, PCA is applied to perform linear dimension reduction. PCA is normally applied to a subset of features prioritized using a feature selection method like highly variable gene selection but can also be applied to all the features or user designated custom subset of features. In Scarf, Sklearn's incremental PCA implementation  61  is used which allows PCA to be trained iteratively using only a subset of cells at a time. Additionally, using Scarf, users can easily use only a subset of cells to train PCA and the rest of cells can still be transformed into this trained PCA space. The data is standard scaled before performing the PCA. For scATAC-Seq datasets, we use Gensim's iterative latent sematic indexing (LSI)  62  method to perform dimension reduction.
Similar to PCA, the LSI method can be trained iteratively and hence is scalable to extremely large datasets. Cell-cell neighbourhood graph creation: Approximate nearest neighbour algorithm, hierarchical small navigable world (HNSW), as implemented in HNSWlib library is utilized  63. The dimension reduced values are used to build the HNSW index. By default, Scarf uses Euclidean metric as the measure of distances between the cells. The index once trained is saved onto the disk for later use. The nearest neighbour search is then performed using the index, querying for 'k' neighbours. During the query, it is noted if cells are indeed the nearest neighbours of themselves and this is summarized and reported as the recall value. The resulting KNN index and the distances to the KNNs is saved onto the disk.
Calculation of initial embedding Before running UMAP and tSNE, cells are usually ascribed an initial embedding which can be either random coordinates (default for many t-SNE algorithms) or informed measures such as first to/three principal component  57  or Laplacian eigenmaps (default for UMAP). In Scarf, the initial embedding for UMAP and tSNE is calculated by fitting MiniBatch Kmeans algorithm from scikit-learn package  61  during the same step as the graph construction. PCA is performed on the Kmeans clusters centroid matrix which of form n x c, where n is the number of cells and c is the number of centroids.
Benchmarking setup All benchmarks were performed on Lunarc computing cluster comprised of nodes with configuration: 16 2. 6GHz CPUs and 190 GB RAM. Each compute node had a local SSD drive and files were moved to local storage to allow exclusive IO operations for each job. A job comprised of either the Scarf (version 0. 5. 9) or Scanpy (version 1. 6. 1) pipeline which was monitored for memory consumption using Linux's 'ps' utility. For benchmarking of cell downsampling, version 1. 2 of GeoSketch was used.
All the seed cells are assigned a constant prize value. Here, we used value of 10. The edge penalty ‫ܧ‬ for each edge is calculated as following: ‫ܧ‬ ൌ ‫ܧ‬. ‫ܧ‬ ௪ ିா ೢ Wherein, ‫ܧ‬ and ‫ܧ‬ ௪ are user provided parameters, edge cost multiplier and edge bandwidth, respectively and ‫ܧ‬ ௪ is the edge weight in the graph. Higher values for ‫ܧ‬ will make reaching remote cells in the graph more difficult but at the same time will discourage inclusion of nonseed cells in the downsampled set. Higher ‫ܧ‬ ௪ accentuates the difference among edge penalties. Here we used ‫ܧ‬ ൌ 1 and ‫ܧ‬ ௪ ൌ 10 KNN mapping and Integrated embedding KNN mapping is performed using the precalculated HNSW index of the reference dataset. By avoiding recalculation of index, Scarf can quickly perform mapping even when cold started on a dataset.
If the two datasets have identical populations, we suggest usage of domain shift correction method, CORAL, that's built into Scarf. CORAL  34  correction is performed as follows: ‫ܣ‬ ൌ ‫‪ሺܵሻ‬ݒܥ‬ ଵ/ଶ ‫ڄ‬ ‫‪ሺܶሻ‬ݒܥ‬ ଵ/ଶ ܶ ൌ ‫ܣ‬ ‫ڄ‬ ܶ Wherein, ܵ and ܶ are scaled and normalized count matrixes of source and target samples respectively and ܶ is the domain corrected matrix of the target dataset. Because the covariance, ‫,ݒܥ‬ of matrices is calculated for features this normalization is easily scalable to large datasets. All reference graphs were constructed with 11 neighbours and using 25 PCA dimensions (calculated on 2500 HVGs). The unified UMAP/SG-tSNE for the reference and target cells are calculated by spiking the reference cell-cell neighbourhood graph with target cells based on their nearest neighbours (top 3 for all the results here) in the reference set.
Comparison of UMAP and SG-tSNE on atlas scale datasets We chose 4 iteration sets for UMAP (100,250,500 and 1000) and ten times for SG-tSNE. We observed that individual iterations of SG-tSNE had lower runtime that UMAP, hence we chose iterations so that they can have comparable runtimes. Both the methods were run in parallel mode using 16 computing cores. Parameters used for UMAP: min_dist=0. 5, spread=2. 0. Parameters used for SG-tSNE: alpha=50, box_h=1. To find the KNN of cells in the UMAP/SG-tSNE space, we used NearestNeighbors function from scikit-learn package with Euclidean distance metric and KD tree algorithm. To calculate the Spearman's correlation between the clusters in neighbourhood graph and UMAP/tSNE space, we took the following approach. We first created a similarity matrix of dimension (P,P) where P is the number of clusters identified in the dataset.
DATA AND CODE AVAILABILITY The source code for Scarf package is available here: github. com/parashardhapola/scarf. The documentation for installation and usage of Scarf can be found here: scarf. readthedocs. io. Download links for all count matrices can be obtained using following command: `from scarf. bio_data import datasets`. The count matrices can be downloaded using command: `scarf. fetch_dataset(sample)`. Where sample is unique id of dataset in datasets. CC-BY-NC 4. 0 International license available under a was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is made The copyright holder for this preprint (which this version posted May 3, 2021. ; https://doi. org/10. 1101/2021. 05. 02. 441899 doi: bioRxiv preprint.
MaNGA is also releasing a number of VACs (Table  2 ). 3. There is no change in the main survey eBOSS data released since DR16, when a total of 1. 4 million eBOSS spectra were released, completing its main survey goals. However, a number of Value Added Catalogs (VACs) useful for cosmological and other applications are released in DR17. The TDSS survey also released its complete dataset in DR16. However, on-going eBOSS-like observations of Xray sources under the SPIDERS program and continued monitoring of quasars under the reverberation mapping program (SDSS-RM) are released in DR17. 4. DR17 also includes data from all previous SDSS data releases. All MaNGA, BOSS, eBOSS, APOGEE and APOGEE-2 spectra that were previously released have all been reprocessed with the latest reduction and analysis pipelines.
sdss. org/ dr17/software/). The SAS also contains the VACs that science team members have contributed to the data releases (see Table  2 ), as well as raw and intermediate data products. All files available through the SAS have a data model that explains their content (https://data. sdss. org/datamodel/). Data products can be downloaded from the SAS either directly through browsing, or by using methods such as wget, rsync and Globus Online (see https://www. sdss. org/dr17/data_access/bulk, for more details and examples). For large data downloads, we recommend the use of Globus Online. Since SDSS data releases are cumulative, in that data products released in earlier data releases are included in DR17, and will have been processed by the latest available pipelines, we reiterate that users should always use the latest data release, as pipelines have often been updated to improve their output and fix previously known bugs.