Interactive visualization of dimensional reduction, clustering and cell propoerties for scRNA-Seq data analysis

David Barrios, Angela Villaverde and Carlos Prieto

2025-10-03

Getting started

looking4clusters is an R package distributed as part of CRAN. To install the package, start R and enter:

if (!require("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("looking4clusters")

Once looking4clusters is installed, it can be loaded by the following command.

library("looking4clusters")

Introduction

looking4clusters performs interactive data visualization to facilitate the analysis of scRNA-Seq. Input data are projected in two-dimensional representations by applying dimensionality reduction methods such as: PCA, MDS, t-SNE, UMAP, NMF. The representation of the data in the same interface with interconnected graphs gives different perspectives which enable an adequate cell classification. The package also integrates the application of unsupervised clustering techniques whose results can be viewed interactively in the graphical interface. As a result, the package generates an interactive web page developed with JavaScript that can be easily explored with a web browser. In addition to the visualization, this interface allows the manual selection of groups, labeling of cell entities based on processed meta-information, generation of new graphs that show gene expression values of each cell, sample identification, and visually comparison of samples and clusters. Interactive visualization can be performed from a numeric matrix or a Seurat Object.

Interactive visualization from a numeric matrix

# Load sample data from scRNAseq package
library(scRNAseq)
sce <- ReprocessedAllenData("tophat_counts")
counts <- assay(sce, "tophat_counts")

# Perform dimensional reduction and an automatic cluster identification
obj <- looking4clusters(t(counts), groups=colData(sce)[,'dissection_s'])

# Output interactive visualization
l4chtml(obj)

Alternatively you can iniltialize an object and then add clusters and reductions

# Create a new looking for cluster object
obj <- looking4clusters(t(counts), running_all=FALSE)

# Add a sample clasification from input data
groups <- colData(sce)[,'dissection_s']
obj <- addcluster(obj, groups, myGroups=TRUE)

# Perform a PCA and TSNE and add to the object as a dimensional reduction layout
libsizes <- colSums(counts)
size.factors <- libsizes/mean(libsizes)
logcounts(sce) <- log2(t(t(counts)/size.factors) + 1)

pca_data <- prcomp(t(logcounts(sce)), rank=50)
obj <- addreduction(obj, pca_data$x[,1:2], "PCA")

library(Rtsne)
tsne_data <- Rtsne(pca_data$x[,1:50], pca = FALSE)
obj <- addreduction(obj, tsne_data$Y, "TSNE")

# Output interactive visualization
l4chtml(obj)

Interactive visualization from a SingleCellExperiment Object

# Adding PCA and TSNE to the object
reducedDims(sce) <- list(PCA=pca_data$x, TSNE=tsne_data$Y)

# Create a looking4clusters object from a SingleCellExperiment object
obj <- looking4clusters(sce, groups="dissection_s")

# Output interactive visualization
l4chtml(obj)

Interactive visualization from a Seurat Object

library(Seurat)
library(Matrix)

# Load sample data from ZilionisLungData
lung <- ZilionisLungData()
immune <- lung$Used & lung$used_in_NSCLC_immune
lung <- lung[,immune]
lung <- lung[1:10000,1:1000]

exp_mat <- Matrix(counts(lung),sparse = TRUE)
colnames(exp_mat) <- paste0(colnames(exp_mat), seq(1,ncol(exp_mat)))

# Create a new Seurat object
seurat_object <- CreateSeuratObject(counts = exp_mat, 
                                    project = "Zilionis_immune")

# Seurat data processing steps
seurat_object <- NormalizeData(seurat_object)
seurat_object <- ScaleData(seurat_object, features = rownames(seurat_object))

seurat_object <- FindVariableFeatures(seurat_object)
seurat_object <- RunPCA(seurat_object,
    features = VariableFeatures(object = seurat_object))

# Create a looking4clusters object from a Seurat object
obj <- looking4clusters(seurat_object)

# Output interactive visualization
l4chtml(obj)

SessionInfo

sessionInfo()
## R version 4.5.1 Patched (2025-08-23 r88802)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.3 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.22-bioc/R/lib/libRblas.so 
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0  LAPACK version 3.12.0
## 
## 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       
## 
## time zone: America/New_York
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] Rtsne_0.17                  scRNAseq_2.23.0            
##  [3] SingleCellExperiment_1.31.1 SummarizedExperiment_1.39.2
##  [5] Biobase_2.69.1              GenomicRanges_1.61.5       
##  [7] Seqinfo_0.99.2              IRanges_2.43.5             
##  [9] S4Vectors_0.47.4            BiocGenerics_0.55.1        
## [11] generics_0.1.4              MatrixGenerics_1.21.0      
## [13] matrixStats_1.5.0           looking4clusters_0.99.4    
## 
## loaded via a namespace (and not attached):
##   [1] RColorBrewer_1.1-3       jsonlite_2.0.0           magrittr_2.0.4          
##   [4] modeltools_0.2-24        GenomicFeatures_1.61.6   gypsum_1.5.0            
##   [7] farver_2.1.2             rmarkdown_2.30           BiocIO_1.19.0           
##  [10] vctrs_0.6.5              memoise_2.0.1            Rsamtools_2.25.3        
##  [13] RCurl_1.98-1.17          htmltools_0.5.8.1        S4Arrays_1.9.1          
##  [16] AnnotationHub_3.99.6     curl_7.0.0               Rhdf5lib_1.31.0         
##  [19] SparseArray_1.9.1        rhdf5_2.53.5             sass_0.4.10             
##  [22] alabaster.base_1.9.5     bslib_0.9.0              alabaster.sce_1.9.0     
##  [25] httr2_1.2.1              cachem_1.1.0             GenomicAlignments_1.45.4
##  [28] lifecycle_1.0.4          pkgconfig_2.0.3          Matrix_1.7-4            
##  [31] R6_2.6.1                 fastmap_1.2.0            digest_0.6.37           
##  [34] AnnotationDbi_1.71.1     irlba_2.3.5.1            ExperimentHub_2.99.5    
##  [37] RSQLite_2.4.3            filelock_1.0.3           httr_1.4.7              
##  [40] abind_1.4-8              compiler_4.5.1           bit64_4.6.0-1           
##  [43] S7_0.2.0                 BiocParallel_1.43.4      viridis_0.6.5           
##  [46] DBI_1.2.3                dendextend_1.19.1        HDF5Array_1.37.0        
##  [49] alabaster.ranges_1.9.1   alabaster.schemas_1.9.0  MASS_7.3-65             
##  [52] rappdirs_0.3.3           DelayedArray_0.35.3      rjson_0.2.23            
##  [55] tools_4.5.1              prabclus_2.3-4           nnet_7.3-20             
##  [58] glue_1.8.0               h5mread_1.1.1            restfulr_0.0.16         
##  [61] rhdf5filters_1.21.0      grid_4.5.1               cluster_2.1.8.1         
##  [64] gtable_0.3.6             class_7.3-23             ensembldb_2.33.2        
##  [67] XVector_0.49.1           flexmix_2.3-20           BiocVersion_3.22.0      
##  [70] pillar_1.11.1            robustbase_0.99-6        dplyr_1.1.4             
##  [73] BiocFileCache_2.99.6     lattice_0.22-7           FNN_1.1.4.1             
##  [76] rtracklayer_1.69.1       bit_4.6.0                tidyselect_1.2.1        
##  [79] Biostrings_2.77.2        knitr_1.50               gridExtra_2.3           
##  [82] ProtGenerics_1.41.0      xfun_0.53                diptest_0.77-2          
##  [85] DEoptimR_1.1-4           UCSC.utils_1.5.0         lazyeval_0.2.2          
##  [88] yaml_2.3.10              evaluate_1.0.5           codetools_0.2-20        
##  [91] kernlab_0.9-33           tibble_3.3.0             alabaster.matrix_1.9.0  
##  [94] BiocManager_1.30.26      cli_3.6.5                uwot_0.2.3              
##  [97] RcppParallel_5.1.11-1    jquerylib_0.1.4          dichromat_2.0-0.1       
## [100] Rcpp_1.1.0               GenomeInfoDb_1.45.12     dbplyr_2.5.1            
## [103] png_0.1-8                XML_3.99-0.19            parallel_4.5.1          
## [106] ggplot2_4.0.0            blob_1.2.4               mclust_6.1.1            
## [109] AnnotationFilter_1.33.0  parallelDist_0.2.7       bitops_1.0-9            
## [112] alabaster.se_1.9.0       viridisLite_0.4.2        scales_1.4.0            
## [115] crayon_1.5.3             fpc_2.2-13               rlang_1.1.6             
## [118] KEGGREST_1.49.1