Chapter 11 Nuclei-only analysis

11.1 Introduction

Single-nuclei RNA-seq (snRNA-seq) is another strategy for performing single-cell transcriptomics where individual nuclei instead of cells are captured and sequenced. The major advantage of snRNA-seq over scRNA-seq is that the former does not require the preservation of cellular integrity during sample preparation, especially dissociation. We only need to extract nuclei in an intact state, meaning that snRNA-seq can be applied to cell types and tissues that are not amenable to dissociation. The cost of this flexibility is the loss of cytoplasmic transcripts that might reduce the biological signal.

The computational analysis of snRNA-seq data is very much like that of scRNA-seq data. We have a gene-by-cell matrix of (UMI) counts that requires quality control, normalization and so on. (Technically, the columns correspond to nuclei but we will use these two terms interchangeably in this chapter.) In fact, the biggest difference in processing occurs in the construction of the count matrix itself, where intronic regions must be included in the annotation for each gene to account for the increased abundance of unspliced transcripts. The scrapper analysis itself only requires a few minor adjustments to account for the loss of cytoplasmic transcripts. We demonstrate using a dataset from Wu et al. (2019) involving snRNA-seq on healthy and fibrotic mouse kidneys.

library(scRNAseq)
sce.nuc <- WuKidneyData()
sce.nuc <- sce.nuc[,sce.nuc$Technology=="sNuc-10x"]
sce.nuc
## class: SingleCellExperiment 
## dim: 18249 8231 
## metadata(0):
## assays(1): counts
## rownames(18249): mt-Cytb mt-Nd6 ... Gm44613 Gm38304
## rowData names(0):
## colnames: NULL
## colData names(4): CellBarcode CellType Technology Status
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):

11.2 Quality control for stripped nuclei

We re-use the same functions in Chapter 1 to compute QC metrics in our snRNA-seq dataset (Figure 11.1). Our aim is to remove low-quality cells with low total counts, low numbers of detected genes and high mitochondrial proportions. The biggest difference from scRNA-seq is that the mitochondrial proportion now represents the efficacy of the stripping process. The presence of mitochondrial counts in a library indicates that the removal of the cytoplasm was not complete, possibly introducing irrelevant heterogeneity in downstream analyses.

library(scrapper)
is.mito.nuc <- grep("^mt-", rownames(sce.nuc))
sce.qc.nuc <- quickRnaQc.se(sce.nuc, subsets=list(Mt=is.mito.nuc))
qc.thresh.nuc <- metadata(sce.qc.nuc)$qc$thresholds
qc.thresh.nuc
## $sum
## [1] 192.4676
## 
## $detected
## [1] 181.0179
## 
## $subset.proportion
## Mt 
##  0
library(scater)
gridExtra::grid.arrange(
    plotColData(sce.qc.nuc, y="sum", colour_by="keep") +
        geom_hline(yintercept=qc.thresh.nuc$sum, linetype="dashed", color="red") +
        scale_y_log10() +
        ggtitle("Total count"),
    plotColData(sce.qc.nuc, y="detected", colour_by="keep") +
        geom_hline(yintercept=qc.thresh.nuc$detected, linetype="dashed", color="red") +
        scale_y_log10() +
        ggtitle("Detected features"),
    plotColData(sce.qc.nuc, y="subset.proportion.Mt", colour_by="keep") + 
        geom_hline(yintercept=qc.thresh.nuc$subset.proportion["Mt"], linetype="dashed", color="red") +
        ggtitle("Mito prop"),
    ncol=3
)
Distribution of QC metrics in the snRNA-seq dataset. Each point represents a cell and is colored according to whether it was retained after QC filtering. Dashed lines represent thresholds for each metric.

Figure 11.1: Distribution of QC metrics in the snRNA-seq dataset. Each point represents a cell and is colored according to whether it was retained after QC filtering. Dashed lines represent thresholds for each metric.

That said, it may be too conservative to remove all libraries with non-zero mitochondrial proportions. Some of these counts may be caused by contamination from mitochondrial transcripts in the ambient solution. Even if the stripping was incomplete, a library might still contain enough biological information to be useful. So, we might consider relaxing the QC filter threshold on the mitochondrial proportions:

# Cells retained using the original filter thresholds.
summary(sce.qc.nuc$keep)
##    Mode   FALSE    TRUE 
## logical    2264    5967
# Relaxing the mitochondrial filter to allow a bit of contamination.
qc.thresh.nuc.relaxed <- qc.thresh.nuc
qc.thresh.nuc.relaxed$subset.proportion["Mt"] <- 0.002

# Supplying custom thresholds to our QC function.
sce.qc.nuc <- quickRnaQc.se(
    sce.nuc,
    subsets = list(Mt = is.mito.nuc),
    thresholds = qc.thresh.nuc.relaxed
)

# More cells retained with the relaxed threshold.
summary(sce.qc.nuc$keep)
##    Mode   FALSE    TRUE 
## logical     432    7799

11.3 Comments on downstream analyses

The rest of the analysis is performed using the same approach described for scRNA-seq (Figure 11.2). Despite the loss of cytoplasmic transcripts, there is usually still enough biological signal to characterize population heterogeneity (Bakken et al. 2018; Wu et al. 2019). In fact, snRNA-seq might even have a higher signal-to-noise ratio than scRNA-seq, as sequencing resources are not wasted on highly abundant but typically uninteresting transcripts for mitochondrial and ribosomal protein genes.

res.full.nuc <- analyze.se(
    sce.nuc,
    rna.qc.subsets=list(Mt=is.mito.nuc),
    # Making use of our filter threshold decisions:
    more.rna.qc.args=list(thresholds=qc.thresh.nuc.relaxed)
)

sce.full.nuc <- res.full.nuc$x
gridExtra::grid.arrange(
    plotReducedDim(sce.full.nuc, dimred="TSNE", colour_by="graph.cluster"),
    plotReducedDim(sce.full.nuc, dimred="TSNE", colour_by="Status"),
    ncol=2
)
$t$-SNE plots of the Wu kidney dataset. Each point is a cell and is colored by its cluster assignment (left) or its disease status (right).

Figure 11.2: \(t\)-SNE plots of the Wu kidney dataset. Each point is a cell and is colored by its cluster assignment (left) or its disease status (right).

We can also apply more complex procedures such as MNN correction (Chapter 8), e.g., to identify shared clusters across healthy and disease samples (Figure 11.3).

# Re-creating our custom filter thresholds for a blocked analysis.
sce.qc.block.nuc <- quickRnaQc.se(
    sce.nuc,
    block=sce.nuc$Status,
    subsets=list(Mt=is.mito.nuc)
)
qc.block.thresh.nuc.relaxed <- metadata(sce.qc.block.nuc)$qc$thresholds
qc.block.thresh.nuc.relaxed$subset.proportion[["Mt"]][] <- 0.002

res.full.block.nuc <- analyze.se(
    sce.nuc,
    block=sce.nuc$Status,
    rna.qc.subsets=list(Mt=is.mito.nuc),
    more.rna.qc.args=list(thresholds=qc.block.thresh.nuc.relaxed)
)

sce.full.block.nuc <- res.full.block.nuc$x
gridExtra::grid.arrange(
    plotReducedDim(sce.full.block.nuc, dimred="TSNE", colour_by="graph.cluster"),
    plotReducedDim(sce.full.block.nuc, dimred="TSNE", colour_by="Status"),
    ncol=2
)
$t$-SNE plots of the Wu kidney dataset after applying MNN correction across samples. Each point is a cell and is colored by its cluster assignment (left) or its disease status (right).

Figure 11.3: \(t\)-SNE plots of the Wu kidney dataset after applying MNN correction across samples. Each point is a cell and is colored by its cluster assignment (left) or its disease status (right).

Similarly, we can perform marker detection on the snRNA-seq expression values as discussed in Chapter 7. For the most part, interpretation of the marker results makes the simplifying assumption that nuclear abundances are a good proxy for the overall expression profile. This is generally reasonable but may not always be true, resulting in some discrepancies in the marker sets between snRNA-seq and scRNA-seq datasets. For example, transcripts for strongly expressed genes might localize to the cytoplasm for efficient translation and subsequently be lost upon stripping, while genes with the same overall expression but differences in the rate of nuclear export may appear to be differentially expressed between clusters.

markers <- res.full.block.nuc$markers$rna
previewMarkers(markers[["8"]])
## DataFrame with 10 rows and 3 columns
##             mean  detected       lfc
##        <numeric> <numeric> <numeric>
## Kcnip4  1.805699  0.755703  1.351806
## Il34    1.712953  0.818187  1.064234
## Pakap   1.460446  0.809479  0.932567
## Esrrg   2.717717  0.973853  1.051882
## Bmp6    1.193670  0.690227  0.893589
## Slc4a4  1.201226  0.752628  0.816459
## Sash1   1.291215  0.778939  0.777754
## Nhs     1.325255  0.784652  0.750925
## Dock10  0.872091  0.591974  0.672949
## Sorcs1  0.839991  0.551595  0.638883

Session information

sessionInfo()
## R version 4.6.1 (2026-06-24)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.4 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.24-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] scater_1.41.2               ggplot2_4.0.3              
##  [3] scuttle_1.23.2              scrapper_1.7.8             
##  [5] scRNAseq_2.27.0             SingleCellExperiment_1.35.2
##  [7] SummarizedExperiment_1.43.0 Biobase_2.73.2             
##  [9] GenomicRanges_1.65.3        Seqinfo_1.3.2              
## [11] IRanges_2.47.5              S4Vectors_0.51.9           
## [13] BiocGenerics_0.59.12        generics_0.1.4             
## [15] MatrixGenerics_1.25.0       matrixStats_1.5.0          
## [17] BiocStyle_2.41.0           
## 
## loaded via a namespace (and not attached):
##   [1] RColorBrewer_1.1-3       jsonlite_2.0.0           magrittr_2.0.5          
##   [4] ggbeeswarm_0.7.3         GenomicFeatures_1.65.0   gypsum_1.9.0            
##   [7] farver_2.1.2             rmarkdown_2.31           BiocIO_1.23.3           
##  [10] vctrs_0.7.3              memoise_2.0.1            Rsamtools_2.29.0        
##  [13] RCurl_1.98-1.20          htmltools_0.5.9          S4Arrays_1.13.0         
##  [16] BiocBaseUtils_1.15.1     AnnotationHub_4.3.2      curl_8.0.0              
##  [19] BiocNeighbors_2.7.3      Rhdf5lib_2.1.0           SparseArray_1.13.2      
##  [22] rhdf5_2.57.12            sass_0.4.10              alabaster.base_1.13.2   
##  [25] bslib_0.12.0             alabaster.sce_1.13.0     httr2_1.3.0             
##  [28] cachem_1.1.0             GenomicAlignments_1.49.1 lifecycle_1.0.5         
##  [31] pkgconfig_2.0.3          rsvd_1.0.5               Matrix_1.7-6            
##  [34] R6_2.6.1                 fastmap_1.2.0            digest_0.6.39           
##  [37] AnnotationDbi_1.75.2     irlba_2.3.7              ExperimentHub_3.3.2     
##  [40] RSQLite_3.53.3           beachmat_2.29.2          labeling_0.4.3          
##  [43] filelock_1.0.3           httr_1.4.8               abind_1.4-8             
##  [46] compiler_4.6.1           bit64_4.8.4              withr_3.0.3             
##  [49] S7_0.2.2                 BiocParallel_1.47.0      viridis_0.6.5           
##  [52] DBI_1.3.0                HDF5Array_1.41.3         alabaster.ranges_1.13.0 
##  [55] alabaster.schemas_1.13.0 rappdirs_0.3.4           DelayedArray_0.39.6     
##  [58] rjson_0.2.23             tools_4.6.1              vipor_0.4.7             
##  [61] otel_0.2.0               beeswarm_0.4.0           glue_1.8.1              
##  [64] h5mread_1.5.2            restfulr_0.0.17          rhdf5filters_1.25.4     
##  [67] grid_4.6.1               gtable_0.3.6             ensembldb_2.37.3        
##  [70] BiocSingular_1.29.1      ScaledMatrix_1.21.0      XVector_0.53.0          
##  [73] ggrepel_0.9.8            BiocVersion_3.24.0       pillar_1.11.1           
##  [76] dplyr_1.2.1              BiocFileCache_3.3.0      lattice_0.23-1          
##  [79] rtracklayer_1.73.0       bit_4.6.0                tidyselect_1.2.1        
##  [82] Biostrings_2.81.7        knitr_1.51               gridExtra_2.3.1         
##  [85] bookdown_0.47            ProtGenerics_1.45.0      xfun_0.60               
##  [88] UCSC.utils_1.9.0         lazyeval_0.2.3           yaml_2.3.12             
##  [91] evaluate_1.0.5           codetools_0.2-20         cigarillo_1.3.1         
##  [94] tibble_3.3.1             alabaster.matrix_1.13.1  BiocManager_1.30.27     
##  [97] cli_3.6.6                jquerylib_0.1.4          dichromat_2.0-1         
## [100] Rcpp_1.1.2               GenomeInfoDb_1.49.1      dbplyr_2.6.0            
## [103] png_0.1-9                XML_3.99-0.24            parallel_4.6.1          
## [106] blob_1.3.0               AnnotationFilter_1.37.0  bitops_1.1-0            
## [109] viridisLite_0.4.3        alabaster.se_1.13.0      scales_1.4.0            
## [112] crayon_1.5.3             rlang_1.3.0              cowplot_1.2.0           
## [115] KEGGREST_1.53.6

References

Bakken, T. E., R. D. Hodge, J. A. Miller, Z. Yao, T. N. Nguyen, B. Aevermann, E. Barkan, et al. 2018. “Single-nucleus and single-cell transcriptomes compared in matched cortical cell types.” PLoS ONE 13 (12): e0209648.

Wu, H., Y. Kirita, E. L. Donnelly, and B. D. Humphreys. 2019. “Advantages of Single-Nucleus over Single-Cell RNA Sequencing of Adult Kidney: Rare Cell Types and Novel Cell States Revealed in Fibrosis.” J. Am. Soc. Nephrol. 30 (1): 23–32.