PyLiger: scalable single-cell multi-omic data integration in Python

2.1 Python implementation of LIGER

We translated the complete, established LIGER framework into Python. Key functions include integration of multiple single-cell datasets using integrative non-negative matrix factorization, joint clustering, visualization and differential expression testing (Fig. 1A). We carefully compared outputs to ensure that function outputs from the R and Python versions are identical to within the limits of numerical precision. The only exceptions are cases when external packages called by PyLiger, such as UMAP and Leiden, produce slightly different results between R and Python.

As an additional feature, we embedded new functionality for gene ontology (GO) enrichment analysis within PyLiger. This makes it much easier to formulate hypotheses about the functions of key genes that are differentially expressed across cell types or biological conditions. Specifically, PyLiger incorporates GOATOOLS (Klopfenstein et al., 2018) for GO enrichment testing and GO-Figure! (Reijnders and Waterhouse, 2021) for visualizing enriched GO terms. Given a list of genes of interest, users can easily run a PyLiger function to identify a list of significantly enriched GO terms. The genes of interest may be derived by e.g. finding genes differentially expressed between clusters; genes differentially expressed across datasets within a cluster or genes with high contribution to metagene factors. Users may then further visualize the GO terms by semantic similarity scatterplots (Supplementary Fig. 1D). The plotting functions allow full customization of colormap, labels and other plot elements.

2.2 PyLiger adapts AnnData format to interoperate with existing packages

We designed the structure of the PyLiger class to smoothly interface with the widely used AnnData format. The AnnData package was initially introduced along with Scanpy offering a convenient way to store data matrices and annotations together. We store cell factor loading matrices (Hⁱ), shared metagenes (W) and dataset-specific metagenes (Vⁱ) as annotations of the raw matrix (Fig. 1B). The use of AnnData format also facilitates interoperability with existing single-cell analysis tools such as Scanpy and scVelo (Bergen et al., 2020). We use the naming rules from Scanpy to name our annotations (UMAP coordinates, for instance) so that each individual AnnData object can be plugged into Scanpy easily.

2.3 Python implementation reduces runtimes

To demonstrate the performance of PyLiger, we tested our functions using a dataset of 100 000 cells sampled from the adult mouse cortex (Saunders et al., 2018). We confirmed that the results from PyLiger are identical (to within numerical precision) to those from the LIGER R package (External packages called by PyLiger, such as UMAP and Leiden, produce slightly different results between R and Python in some cases.) (Supplementary Fig. 1A and B). Moreover, PyLiger functions run 1.5–5 times faster than their R counterparts (Supplementary Fig. 2). In particular, the most time-consuming step—matrix factorization (using in-memory mode iNMF)—is ∼3 times faster in Python than our previous R implementation. This is particularly impressive because many of the R functions are implemented using Rcpp, whereas all of PyLiger are simply implemented in native Python.

2.4 PyLiger scales to arbitrarily large single-cell datasets using fixed memory

PyLiger uses the HDF5 file format for on-demand loading of datasets stored on disk. We found that in AnnData objects, only the raw matrix allows HDF5-based backing, but not other processed matrices stored as layers. Therefore, we store data matrices in a separate HDF5 file while matrix annotations are still stored in AnnData objects. We compared the on-disk mode to the in-memory mode using the same dataset of 100 000 cells. By sacrificing a little processing efficiency (about 2 s on a dataset of 100 000 cells), the on-disk mode functions can process arbitrarily large datasets using fixed memory. Note that the function create_liger in the on-disk mode of PyLiger is slightly slower than on-disk mode of LIGER due to new feature implementation.

Moreover, we implemented the online iNMF algorithm (Gao et al., 2021) in combination with HDF5 file format, providing scalable and efficient data integration as well as significant memory savings. The online iNMF algorithm scales to arbitrarily large numbers of cells but still uses fixed memory and can incorporate new data without recalculating from scratch. To benchmark the performance, we did a comparison of online iNMF between PyLiger and LIGER using five datasets of increasing sizes (ranging from 10 000 to 255 353 cells in total) sampled from the same adult mouse cortex dataset. The PyLiger implementation of online iNMF achieves a 2.3× speedup on average in comparison to its R counterpart (Supplementary Fig. 2B). To demonstrate that the online iNMF algorithm can perform data integration using fixed memory, we tested both the in-memory and the on-disk mode of online iNMF on the same five datasets. With the increasing of number of cells, the on-disk mode of online iNMF has constant peak memory consumption (Supplementary Fig. 2C).