1 Introduction

Single-cell ’omics analysis enables high-resolution characterization of heterogeneous populations of cells by quantifying measurements in individual cells and thus provides a fuller, more nuanced picture into the complexity and heterogeneity between cells. However, the data also present new and significant challenges as compared to previous approaches, especially as single-cell data are much larger and sparser than data generated from bulk sequencing methods. Dimensionality reduction is a key step in the single-cell analysis to address the high dimensionality and sparsity of these data, and to enable the application of more complex, computationally expensive downstream pipelines.

Correspondence analysis (CA) is a matrix factorization method, and is similar to principal components analysis (PCA). Whereas PCA is designed for application to continuous, approximately normally distributed data, CA is appropriate for non-negative, count-based data that are in the same additive scale. corral implements CA for dimensionality reduction of a single matrix of single-cell data.

See the vignette for corralm for the multi-table adaptation of CA for single-cell batch alignment/integration.

corral can be used with various types of input. When called on a matrix (or other matrix-like object), it returns a list with the SVD output, principal coordinates, and standard coordinates. When called on a SingleCellExperiment, it returns the SingleCellExperiment with the corral embeddings in the reducedDim slot named corral. To retrieve the full list output from a SingleCellExperiment input, the fullout argument can be set to TRUE.

2 Loading packages and data

We will use the Zhengmix4eq dataset from the DuoClustering2018 package.

library(corral)
library(SingleCellExperiment)
library(ggplot2)
library(DuoClustering2018)
zm4eq.sce <- sce_full_Zhengmix4eq()

This dataset includes approximately 4,000 pre-sorted and annotated cells of 4 types mixed by Duo et al. in approximately equal proportions (Duò, Robinson, and Soneson, n.d.). The cells were sampled from a “Massively parallel digital transcriptional profiling of single cells” (Zheng et al. 2017).

zm4eq.sce
## class: SingleCellExperiment 
## dim: 15568 3994 
## metadata(1): log.exprs.offset
## assays(3): counts logcounts normcounts
## rownames(15568): ENSG00000237683 ENSG00000228327 ... ENSG00000215700
##   ENSG00000215699
## rowData names(10): id symbol ... total_counts log10_total_counts
## colnames(3994): b.cells1147 b.cells6276 ... regulatory.t1084
##   regulatory.t9696
## colData names(14): dataset barcode ... libsize.drop feature.drop
## reducedDimNames(2): PCA TSNE
## mainExpName: NULL
## altExpNames(0):
table(colData(zm4eq.sce)$phenoid)
## 
##         b.cells  cd14.monocytes naive.cytotoxic    regulatory.t 
##             999            1000             998             997

3 corral on SingleCellExperiment

We will run corral directly on the raw count data:

zm4eq.sce <- corral(inp = zm4eq.sce, 
                    whichmat = 'counts')

zm4eq.sce
## class: SingleCellExperiment 
## dim: 15568 3994 
## metadata(1): log.exprs.offset
## assays(3): counts logcounts normcounts
## rownames(15568): ENSG00000237683 ENSG00000228327 ... ENSG00000215700
##   ENSG00000215699
## rowData names(10): id symbol ... total_counts log10_total_counts
## colnames(3994): b.cells1147 b.cells6276 ... regulatory.t1084
##   regulatory.t9696
## colData names(14): dataset barcode ... libsize.drop feature.drop
## reducedDimNames(3): PCA TSNE corral
## mainExpName: NULL
## altExpNames(0):

We can use plot_embedding to visualize the output:

plot_embedding_sce(sce = zm4eq.sce,
                   which_embedding = 'corral',
                   plot_title = 'corral on Zhengmix4eq',
                   color_attr = 'phenoid',
                   color_title = 'cell type',
                   saveplot = FALSE)

Using the scater package, we can also add and visualize umap and tsne embeddings based on the corral output:

library(scater)
## Loading required package: scuttle
library(gridExtra) # so we can arrange the plots side by side

zm4eq.sce <- runUMAP(zm4eq.sce,
                     dimred = 'corral',
                     name = 'corral_UMAP')
zm4eq.sce <- runTSNE(zm4eq.sce,
                     dimred = 'corral',
                     name = 'corral_TSNE')

ggplot_umap <- plot_embedding_sce(sce = zm4eq.sce,
                                  which_embedding = 'corral_UMAP',
                                  plot_title = 'Zhengmix4eq corral with UMAP',
                                  color_attr = 'phenoid',
                                  color_title = 'cell type',
                                  returngg = TRUE,
                                  showplot = FALSE,
                                  saveplot = FALSE)

ggplot_tsne <- plot_embedding_sce(sce = zm4eq.sce,
                                  which_embedding = 'corral_TSNE',
                                  plot_title = 'Zhengmix4eq corral with tSNE',
                                  color_attr = 'phenoid',
                                  color_title = 'cell type',
                                  returngg = TRUE,
                                  showplot = FALSE,
                                  saveplot = FALSE)

multiplot(ggplot_umap, ggplot_tsne, cols = 2)
## Warning: 'multiplot' is deprecated.
## Use 'gridExtra::grid.arrange' instead.
## See help("Deprecated")

The corral embeddings stored in the reducedDim slot can be used in downstream analysis, such as for clustering or trajectory analysis.

corral can also be run on a SummarizedExperiment object.

4 corral on matrix

corral can also be performed on a matrix (or matrix-like) input.

zm4eq.countmat <- assay(zm4eq.sce,'counts')
zm4eq.countcorral <- corral(zm4eq.countmat)

The output is in a list format, including the SVD output (u,d,v), the standard coordinates (SCu,SCv), and the principal coordinates (PCu,PCv).

zm4eq.countcorral
## corral output summary===========================================
##   Output "list" includes standard coordinates (SCu, SCv),
##   principal coordinates (PCu, PCv), & SVD output (u, d, v)
## Variance explained----------------------------------------------
##                           PC1  PC2  PC3  PC4  PC5  PC6  PC7  PC8
## percent.Var.explained    0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00
## cumulative.Var.explained 0.01 0.02 0.02 0.02 0.03 0.03 0.03 0.03
## 
## Dimensions of output elements-----------------------------------
##   Singular values (d) :: 30
##   Left singular vectors & coordinates (u, SCu, PCu) :: 15568 30
##   Right singular vectors & coordinates (v, SCv, PCv) :: 3994 30
##   See corral help for details on each output element.
##   Use plot_embedding to visualize; see docs for details.
## ================================================================

We can use plot_embedding to visualize the output: (the embeddings are in the v matrix because these data are by genes in the rows and have cells in the columns; if this were reversed, with cells in the rows and genes/features in the column, then the cell embeddings would instead be in the u matrix.)

celltype_vec <- colData(zm4eq.sce)$phenoid
plot_embedding(embedding = zm4eq.countcorral$v,
               plot_title = 'corral on Zhengmix4eq',
               color_vec = celltype_vec,
               color_title = 'cell type',
               saveplot = FALSE)

The output is the same as above with the SingleCellExperiment, and can be passed as the low-dimension embedding for downstream analysis. Similarly, UMAP and tSNE can be computed for visualization. (Note that in performing SVD, the direction of the axes doesn’t matter so they may be flipped between runs, as corral and corralm use irlba to perform fast approximation.)

5 Session information

sessionInfo()
## R version 4.1.0 (2021-05-18)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.2 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.13-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.13-bioc/R/lib/libRlapack.so
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_GB              LC_COLLATE=C              
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## attached base packages:
## [1] parallel  stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] scater_1.20.0               scuttle_1.2.0              
##  [3] DuoClustering2018_1.9.0     ggplot2_3.3.3              
##  [5] SingleCellExperiment_1.14.0 SummarizedExperiment_1.22.0
##  [7] Biobase_2.52.0              GenomicRanges_1.44.0       
##  [9] GenomeInfoDb_1.28.0         IRanges_2.26.0             
## [11] S4Vectors_0.30.0            BiocGenerics_0.38.0        
## [13] MatrixGenerics_1.4.0        matrixStats_0.58.0         
## [15] corral_1.2.0                gridExtra_2.3              
## [17] BiocStyle_2.20.0           
## 
## loaded via a namespace (and not attached):
##   [1] Rtsne_0.15                    ggbeeswarm_0.6.0             
##   [3] colorspace_2.0-1              ellipsis_0.3.2               
##   [5] mclust_5.4.7                  XVector_0.32.0               
##   [7] BiocNeighbors_1.10.0          dichromat_2.0-0              
##   [9] farver_2.1.0                  MultiAssayExperiment_1.18.0  
##  [11] bit64_4.0.5                   RSpectra_0.16-0              
##  [13] interactiveDisplayBase_1.30.0 AnnotationDbi_1.54.0         
##  [15] fansi_0.4.2                   sparseMatrixStats_1.4.0      
##  [17] cachem_1.0.5                  knitr_1.33                   
##  [19] jsonlite_1.7.2                dbplyr_2.1.1                 
##  [21] png_0.1-7                     uwot_0.1.10                  
##  [23] shiny_1.6.0                   BiocManager_1.30.15          
##  [25] mapproj_1.2.7                 compiler_4.1.0               
##  [27] httr_1.4.2                    assertthat_0.2.1             
##  [29] Matrix_1.3-3                  fastmap_1.1.0                
##  [31] BiocSingular_1.8.0            later_1.2.0                  
##  [33] htmltools_0.5.1.1             tools_4.1.0                  
##  [35] rsvd_1.0.5                    gtable_0.3.0                 
##  [37] glue_1.4.2                    GenomeInfoDbData_1.2.6       
##  [39] reshape2_1.4.4                dplyr_1.0.6                  
##  [41] ggthemes_4.2.4                maps_3.3.0                   
##  [43] rappdirs_0.3.3                Rcpp_1.0.6                   
##  [45] jquerylib_0.1.4               vctrs_0.3.8                  
##  [47] Biostrings_2.60.0             ExperimentHub_2.0.0          
##  [49] DelayedMatrixStats_1.14.0     xfun_0.23                    
##  [51] stringr_1.4.0                 beachmat_2.8.0               
##  [53] mime_0.10                     lifecycle_1.0.0              
##  [55] irlba_2.3.3                   AnnotationHub_3.0.0          
##  [57] zlibbioc_1.38.0               scales_1.1.1                 
##  [59] promises_1.2.0.1              yaml_2.2.1                   
##  [61] curl_4.3.1                    memoise_2.0.0                
##  [63] sass_0.4.0                    stringi_1.6.2                
##  [65] RSQLite_2.2.7                 BiocVersion_3.13.1           
##  [67] highr_0.9                     ScaledMatrix_1.0.0           
##  [69] filelock_1.0.2                BiocParallel_1.26.0          
##  [71] pals_1.7                      rlang_0.4.11                 
##  [73] pkgconfig_2.0.3               bitops_1.0-7                 
##  [75] evaluate_0.14                 lattice_0.20-44              
##  [77] purrr_0.3.4                   labeling_0.4.2               
##  [79] transport_0.12-2              bit_4.0.4                    
##  [81] tidyselect_1.1.1              plyr_1.8.6                   
##  [83] magrittr_2.0.1                bookdown_0.22                
##  [85] R6_2.5.0                      magick_2.7.2                 
##  [87] generics_0.1.0                DelayedArray_0.18.0          
##  [89] DBI_1.1.1                     pillar_1.6.1                 
##  [91] withr_2.4.2                   KEGGREST_1.32.0              
##  [93] RCurl_1.98-1.3                tibble_3.1.2                 
##  [95] crayon_1.4.1                  utf8_1.2.1                   
##  [97] BiocFileCache_2.0.0           rmarkdown_2.8                
##  [99] viridis_0.6.1                 grid_4.1.0                   
## [101] data.table_1.14.0             FNN_1.1.3                    
## [103] blob_1.2.1                    digest_0.6.27                
## [105] xtable_1.8-4                  tidyr_1.1.3                  
## [107] httpuv_1.6.1                  munsell_0.5.0                
## [109] beeswarm_0.3.1                viridisLite_0.4.0            
## [111] vipor_0.4.5                   bslib_0.2.5.1

References

Duò, A, MD Robinson, and C Soneson. n.d. “A Systematic Performance Evaluation of Clustering Methods for Single-Cell Rna-Seq Data [Version 2; Peer Review: 2 Approved], Journal = F1000Research, Volume = 7, Year = 2018, Number = 1141, Doi = 10.12688/f1000research.15666.2.”

Zheng, Grace X. Y., Jessica M. Terry, Phillip Belgrader, Paul Ryvkin, Zachary W. Bent, Ryan Wilson, Solongo B. Ziraldo, et al. 2017. “Massively Parallel Digital Transcriptional Profiling of Single Cells.” Nature Communications 8 (1): 14049. https://doi.org/10.1038/ncomms14049.