Use of SpatialDecon in a Spatial Transcriptomics dataset

Installation


if(!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("SpatialDecon")

Overview

This vignette demonstrates the use of the SpatialDecon package to estimate cell abundance in spatial gene expression studies.

The workflow demonstrated here focuses on analysis of Seurat objects, which are commonly used to store Visium and Spatial Transcriptomics data. See our other vignettes for examples in GeoMx data.

We’ll analyze a Spatial Transcriptomics dataset from a HER2+ breast tumor, looking for abundance of different immune cell types.

Data preparation

First, we load the package:

library(SpatialDecon)
library(SeuratObject)

Let’s load the example ST Seurat object and examine it:

# download data:
con <- gzcon(url("https://github.com/almaan/her2st/raw/master/data/ST-cnts/G1.tsv.gz"))
txt <- readLines(con)
temp <- read.table(textConnection(txt), sep = "\t", header = TRUE, row.names = 1)
# parse data
raw = t(as.matrix(temp))
norm = sweep(raw, 2, colSums(raw), "/") * mean(colSums(raw))
x = as.numeric(substr(rownames(temp), 1, unlist(gregexpr("x", rownames(temp))) - 1))
y = -as.numeric(substr(rownames(temp), unlist(gregexpr("x", rownames(temp))) + 1, nchar(rownames(temp))))
# put into a seurat object:
andersson_g1 = CreateSeuratObject(counts = raw, assay="Spatial")
#> Warning: The following arguments are not used: row.names
andersson_g1@meta.data$x = x
andersson_g1@meta.data$y = y

Cell profile matrices

A “cell profile matrix” is a pre-defined matrix that specifies the expected expression profiles of each cell type in the experiment. The SpatialDecon library comes with one such matrix pre-loaded, the “SafeTME” matrix, designed for estimation of immune and stroma cells in the tumor microenvironment. (This matrix was designed to avoid genes commonly expressed by cancer cells; see the SpatialDecon manuscript for details.)

Let’s take a glance at the safeTME matrix:

data("safeTME")
data("safeTME.matches")

signif(safeTME[seq_len(3), seq_len(3)], 2)
#>       macrophages  mast B.naive
#> A2M        0.8500 0.044  0.0043
#> ABCB1      0.0021 0.023  0.0250
#> ABCB4      0.0044 0.000  0.2200

heatmap(sweep(safeTME, 1, apply(safeTME, 1, max), "/"),
        labRow = NA, margins = c(10, 5))
The safeTME cell profile matrix

The safeTME cell profile matrix

For studies of other tissue types, we have provided a library of cell profile matrices, available on Github and downloadable with the “download_profile_matrix” function.

For a complete list of matrices, see CellProfileLibrary GitHub Page.

Below we download a matrix of cell profiles derived from scRNA-seq of a mouse spleen.

mousespleen <- download_profile_matrix(species = "Mouse",
                                       age_group = "Adult", 
                                       matrixname = "Spleen_MCA")
dim(mousespleen)
#> [1] 11125     9

mousespleen[1:4,1:4]
#>               Dendritic.cell.S100a4.high Dendritic.cell.Siglech.high
#> 0610009B22Rik                 0.02985075                   0.0000000
#> 0610010F05Rik                 0.00000000                   0.0000000
#> 0610010K14Rik                 0.02985075                   0.0000000
#> 0610012G03Rik                 0.08955224                   0.1111111
#>               Granulocyte Macrophage
#> 0610009B22Rik  0.00000000 0.00000000
#> 0610010F05Rik  0.00000000 0.00000000
#> 0610010K14Rik  0.00000000 0.03846154
#> 0610012G03Rik  0.08571429 0.03846154

head(cellGroups)
#> $Dendritic
#> [1] "Dendritic.cell.S100a4.high"  "Dendritic.cell.Siglech.high"
#> 
#> $Granulocyte
#> [1] "Granulocyte"
#> 
#> $Macrophage
#> [1] "Macrophage"
#> 
#> $`Marginal zone B`
#> [1] "Marginal.zone.B.cell"
#> 
#> $Monocyte
#> [1] "Monocyte"
#> 
#> $NK
#> [1] "NK.cell"

metadata
#>    Profile Matrix Tissue Species  Strain        Age Age Group
#> 33            MCA Spleen   Mouse C57BL/6 6-10 weeks     Adult
#>                                          URL
#> 33 https://pubmed.ncbi.nlm.nih.gov/29474909/
#>                                                                                        Citation
#> 33 Han, X. et al. Mapping the Mouse Cell Atlas by Microwell-Seq. Cell172, 1091-1107.e17 (2018).

heatmap(sweep(mousespleen, 1, apply(mousespleen, 1, max), "/"),
        labRow = NA, margins = c(10, 5), cexCol = 0.7)
The Mouse Spleen profile matrix

The Mouse Spleen profile matrix

For studies where the provided cell profile matrices aren’t sufficient or if a specific single cell dataset is wanted, we can make a custom profile matrix using the function create_profile_matrix().

This mini single cell dataset is a fraction of the data from Kinchen, J. et al. Structural Remodeling of the Human Colonic Mesenchyme in Inflammatory Bowel Disease. Cell 175, 372-386.e17 (2018).

data("mini_singleCell_dataset")

mini_singleCell_dataset$mtx@Dim # genes x cells
#> [1] 1814  250

as.matrix(mini_singleCell_dataset$mtx)[1:4,1:4]
#>          ACTGCTCGTAAGTTCC.S90 TGAAAGAAGGCGCTCT.S66 AGCTTGAGTTTGGGCC.S66
#> PLEKHN1                     0                    0                    0
#> PERM1                       0                    0                    0
#> C1orf159                    0                    0                    0
#> TTLL10                      0                    0                    0
#>          ACGGCCATCGTCTGAA.S66
#> PLEKHN1                     0
#> PERM1                       0
#> C1orf159                    0
#> TTLL10                      0

head(mini_singleCell_dataset$annots)
#>                    CellID  LabeledCellType
#> 2660 ACTGCTCGTAAGTTCC.S90     stromal cell
#> 2162 TGAAAGAAGGCGCTCT.S66       glial cell
#> 368  AGCTTGAGTTTGGGCC.S66 endothelial cell
#> 238  ACGGCCATCGTCTGAA.S66     stromal cell
#> 4158 TCTCTAACACTGTTAG.S90     stromal cell
#> 2611 ACGATGTGTGTGGTTT.S90     stromal cell

table(mini_singleCell_dataset$annots$LabeledCellType)
#> 
#>            endothelial cell                  glial cell 
#>                          14                          12 
#>               pericyte cell                 plasma cell 
#>                           3                          30 
#> smooth muscle cell of colon                stromal cell 
#>                           1                         190

Pericyte cell and smooth muscle cell of colon will be dropped from this matrix due to low cell count. The average expression across all cells of one type is returned so the more cells of one type, the better reflection of the true gene expression. The confidence in these averages can be changed using the minCellNum filter.

custom_mtx <- create_profile_matrix(mtx = mini_singleCell_dataset$mtx,            # cell x gene count matrix
                                    cellAnnots = mini_singleCell_dataset$annots,  # cell annotations with cell type and cell name as columns 
                                    cellTypeCol = "LabeledCellType",  # column containing cell type
                                    cellNameCol = "CellID",           # column containing cell ID/name
                                    matrixName = "custom_mini_colon", # name of final profile matrix
                                    outDir = NULL,                    # path to desired output directory, set to NULL if matrix should not be written
                                    normalize = FALSE,                # Should data be normalized? 
                                    minCellNum = 5,                   # minimum number of cells of one type needed to create profile, exclusive
                                    minGenes = 10,                    # minimum number of genes expressed in a cell, exclusive
                                    scalingFactor = 5,                # what should all values be multiplied by for final matrix
                                    discardCellTypes = TRUE)          # should cell types be filtered for types like mitotic, doublet, low quality, unknown, etc.
#> [1] "Creating Atlas"
#> [1] "1 / 6 : stromal cell"
#> [1] "2 / 6 : glial cell"
#> [1] "3 / 6 : endothelial cell"
#> [1] "4 / 6 : plasma cell"
#> [1] "5 / 6 : pericyte cell"
#> Warning in create_profile_matrix(mtx = mini_singleCell_dataset$mtx, cellAnnots = mini_singleCell_dataset$annots, : 
#>  pericyte cell was dropped from matrix because it didn't have enough viable cells based on current filtering thresholds. 
#>                                         If this cell type is necessary consider changing minCellNum or minGenes
#> [1] "6 / 6 : smooth muscle cell of colon"
#> Warning in create_profile_matrix(mtx = mini_singleCell_dataset$mtx, cellAnnots = mini_singleCell_dataset$annots, : 
#>  smooth muscle cell of colon was dropped from matrix because it didn't have enough viable cells based on current filtering thresholds. 
#>                                         If this cell type is necessary consider changing minCellNum or minGenes

head(custom_mtx)
#>          stromal cell glial cell endothelial cell plasma cell
#> PLEKHN1     0.0000000          0        0.0000000  0.05787067
#> PERM1       0.1167131          0        0.0000000  0.00000000
#> C1orf159    0.0000000          0        0.0000000  0.07922668
#> TTLL10      0.0000000          0        0.1853931  0.00000000
#> TAS1R3      0.0000000          0        0.0000000  0.07039008
#> ATAD3C      0.1219237          0        0.0000000  0.07922668

heatmap(sweep(custom_mtx, 1, apply(custom_mtx, 1, max), "/"),
        labRow = NA, margins = c(10, 5), cexCol = 0.7)
Custom profile matrix

Custom profile matrix

Custom matrices can be created from all single cell data classes as long as a counts matrix and cell annotations can be passed to the function. Here is an example of creating a matrix using a Seurat object.

library(SeuratObject)

data("mini_singleCell_dataset")

rownames(mini_singleCell_dataset$annots) <- mini_singleCell_dataset$annots$CellID

seuratObject <- CreateSeuratObject(counts = mini_singleCell_dataset$mtx, meta.data = mini_singleCell_dataset$annots)
Idents(seuratObject) <- seuratObject$LabeledCellType

rm(mini_singleCell_dataset)

annots <- data.frame(cbind(cellType=as.character(Idents(seuratObject)), 
                           cellID=names(Idents(seuratObject))))

custom_mtx_seurat <- create_profile_matrix(mtx = seuratObject@assays$RNA@counts, 
                                           cellAnnots = annots, 
                                           cellTypeCol = "cellType", 
                                           cellNameCol = "cellID", 
                                           matrixName = "custom_mini_colon",
                                           outDir = NULL, 
                                           normalize = FALSE, 
                                           minCellNum = 5, 
                                           minGenes = 10)
#> [1] "Creating Atlas"
#> [1] "1 / 6 : stromal cell"
#> [1] "2 / 6 : glial cell"
#> [1] "3 / 6 : endothelial cell"
#> [1] "4 / 6 : plasma cell"
#> [1] "5 / 6 : pericyte cell"
#> Warning in create_profile_matrix(mtx = seuratObject@assays$RNA@counts, cellAnnots = annots, : 
#>  pericyte cell was dropped from matrix because it didn't have enough viable cells based on current filtering thresholds. 
#>                                         If this cell type is necessary consider changing minCellNum or minGenes
#> [1] "6 / 6 : smooth muscle cell of colon"
#> Warning in create_profile_matrix(mtx = seuratObject@assays$RNA@counts, cellAnnots = annots, : 
#>  smooth muscle cell of colon was dropped from matrix because it didn't have enough viable cells based on current filtering thresholds. 
#>                                         If this cell type is necessary consider changing minCellNum or minGenes

head(custom_mtx_seurat)
#>          stromal cell glial cell endothelial cell plasma cell
#> PLEKHN1     0.0000000          0        0.0000000  0.05787067
#> PERM1       0.1167131          0        0.0000000  0.00000000
#> C1orf159    0.0000000          0        0.0000000  0.07922668
#> TTLL10      0.0000000          0        0.1853931  0.00000000
#> TAS1R3      0.0000000          0        0.0000000  0.07039008
#> ATAD3C      0.1219237          0        0.0000000  0.07922668

paste("custom_mtx and custom_mtx_seurat are identical", all(custom_mtx == custom_mtx_seurat))
#> [1] "custom_mtx and custom_mtx_seurat are identical TRUE"

Deconvolving a Seurat object with the runspatialdecon function

Now our data is ready for deconvolution. First we’ll show how to use spatialdecon under the basic settings, omitting optional bells and whistles.

res = runspatialdecon(object = andersson_g1,
                      bg = 0.01,
                      X = safeTME,
                      align_genes = TRUE)
str(res)
#> List of 10
#>  $ beta            : num [1:18, 1:441] 5.92 4.71 1.94 1.7 1.41 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ sigmas          : num [1:18, 1:18, 1:441] 0.617854 -0.16539 -0.011852 -0.00779 -0.000584 ...
#>   ..- attr(*, "dimnames")=List of 3
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ yhat            : num [1:693, 1:441] 0.933 4.288 0.787 0.09 0.884 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:693] "TNFRSF4" "MXRA8" "PLCH2" "CA6" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ resids          : num [1:693, 1:441] -0.9 -3.1 -0.655 0 -0.822 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:693] "TNFRSF4" "MXRA8" "PLCH2" "CA6" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ p               : num [1:18, 1:441] 5.00e-14 7.11e-10 4.53e-02 1.33e-01 4.85e-05 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ t               : num [1:18, 1:441] 7.53 6.16 2 1.5 4.06 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ se              : num [1:18, 1:441] 0.786 0.765 0.971 1.133 0.347 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ prop_of_all     : num [1:18, 1:441] 0.1461 0.1163 0.048 0.042 0.0348 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ prop_of_nontumor: num [1:18, 1:441] 0.1461 0.1163 0.048 0.042 0.0348 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ X               : num [1:693, 1:18] 0.011117 0.000305 0 0.001419 0.002561 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:693] "TNFRSF4" "MXRA8" "PLCH2" "CA6" ...
#>   .. ..$ : chr [1:18] "macrophages" "mast" "B.naive" "B.memory" ...

We’re most interested in “beta”, the matrix of estimated cell abundances.

heatmap(res$beta, cexCol = 0.5, cexRow = 0.7, margins = c(10,7))
Cell abundance estimates

Cell abundance estimates

Using the advanced settings of spatialdecon

spatialdecon has several abilities beyond basic deconvolution:

  1. If given the nuclei counts for each region/observation, it returns results on the scale of total cell counts. This option is generally not available for ST/Visium data.
  2. If given the identities of pure tumor regions/observations, it infers a handful of tumor-specific expression profiles and appends them to the cell profile matrix. Doing this accounts for cancer cell-derived expression from any genes in the cell profile matrix, removing contaminating signal from cancer cells. This operation is complicated in ST/Visium studies, since they generally do not contain regions previously marked as pure tumor. As a heuristic for identifying pure tumor regions, we identify a small subset of regions with relatively low counts from genes in the safeTME matrix vs. from the rest of the transcriptome.
  3. If given raw count data, it derives per-data-point weights. For Visium/ST data, we assume Poisson error.
  4. If given a “cellmatches” argument, it sums multiple closely-related cell types into a single score. E.g. if the safeTME matrix is used with the cell-matching data object “safeTME.matches”, it e.g. sums the “T.CD8.naive” and “T.CD8.memory” scores into a single “CD8.T.cells” score.

Let’s take a look at an example cell matching object:

str(safeTME.matches)
#> List of 14
#>  $ macrophages      : chr "macrophages"
#>  $ mast             : chr "mast"
#>  $ B                : chr [1:2] "B.naive" "B.memory"
#>  $ plasma           : chr "plasma"
#>  $ CD4.T.cells      : chr [1:2] "T.CD4.naive" "T.CD4.memory"
#>  $ CD8.T.cells      : chr [1:2] "T.CD8.naive" "T.CD8.memory"
#>  $ NK               : chr "NK"
#>  $ pDC              : chr "pDCs"
#>  $ mDCs             : chr "mDCs"
#>  $ monocytes        : chr [1:2] "monocytes.C" "monocytes.NC.I"
#>  $ neutrophils      : chr "neutrophils"
#>  $ Treg             : chr "Treg"
#>  $ endothelial.cells: chr "endothelial.cells"
#>  $ fibroblasts      : chr "fibroblasts"

Now let’s mark some regions as pure tumor:

sharedgenes = intersect(rownames(raw), rownames(safeTME))
plot(colSums(raw), colSums(raw[sharedgenes, ]), log = "xy")

hist(colSums(raw[sharedgenes, ]) / colSums(raw), breaks = 20)

alltumor = colSums(raw[sharedgenes, ]) / colSums(raw) < 0.03  # for alma data

table(alltumor)
#> alltumor
#> FALSE  TRUE 
#>   388    53

Calculate weights:

sd_from_noise <- runErrorModel(counts = raw, platform = "st") 
wts <- 1 / sd_from_noise

Now let’s run spatialdecon:


# run spatialdecon with all the bells and whistles:
restils = runspatialdecon(object = andersson_g1,           # Seurat object 
                          X = safeTME,                     # safeTME matrix, used by default
                          bg = 0.01,                       # Recommended value of 0.01 in Visium/ST data
                          wts = wts,                       # weight
                          cellmerges = safeTME.matches,    # safeTME.matches object, used by default
                          is_pure_tumor = alltumor,        # identities of the Tumor segments/observations
                          n_tumor_clusters = 5)            # how many distinct tumor profiles to append to safeTME

str(restils)
#> List of 12
#>  $ beta            : num [1:14, 1:441] 0.0144 0 0 0 0 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ yhat            : num [1:693, 1:441] 1.002 0.999 0.999 0.994 0.994 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:693] "TNFRSF4" "MXRA8" "PLCH2" "CA6" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ resids          : num [1:693, 1:441] -1.003 -0.998 -0.998 -0.992 -0.992 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:693] "TNFRSF4" "MXRA8" "PLCH2" "CA6" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ p               : num [1:14, 1:441] 0.811 1 1 1 1 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ t               : num [1:14, 1:441] 0.24 0 0 0 0 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ se              : num [1:14, 1:441] 0.05991 0.00633 0.17845 0.0811 0.5856 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ beta.granular   : num [1:23, 1:441] 0.0144 0 0 0 0 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:23] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ sigma.granular  : num [1:23, 1:23, 1:441] 3.59e-03 -1.01e-06 9.92e-04 -1.21e-03 -2.62e-05 ...
#>   ..- attr(*, "dimnames")=List of 3
#>   .. ..$ : chr [1:23] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:23] "macrophages" "mast" "B.naive" "B.memory" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ sigma           : num [1:14, 1:14, 1:441] 3.59e-03 -1.01e-06 -2.13e-04 -2.62e-05 8.03e-04 ...
#>   ..- attr(*, "dimnames")=List of 3
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ prop_of_all     : num [1:14, 1:441] 0.208 0 0 0 0 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ prop_of_nontumor: num [1:14, 1:441] 0.208 0 0 0 0 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:14] "macrophages" "mast" "B" "plasma" ...
#>   .. ..$ : chr [1:441] "10x10" "10x11" "10x12" "10x13" ...
#>  $ X               : num [1:693, 1:23] 0.011117 0.000305 0 0.001419 0.002561 ...
#>   ..- attr(*, "dimnames")=List of 2
#>   .. ..$ : chr [1:693] "TNFRSF4" "MXRA8" "PLCH2" "CA6" ...
#>   .. ..$ : chr [1:23] "macrophages" "mast" "B.naive" "B.memory" ...

There are quite a few readouts here. Let’s review the important ones:

To illustrate the derivation of tumor profiles, let’s look at the cell profile matrix output by spatialdecon:

heatmap(sweep(restils$X, 1, apply(restils$X, 1, max), "/"),
         labRow = NA, margins = c(10, 5))
safeTME merged with newly-derived tumor profiles

safeTME merged with newly-derived tumor profiles

Note the new tumor-specific columns.

Plotting deconvolution results

The “florets” function plotting cell abundances atop some 2-D projection. Here, we’ll plot cell abundances atop the first 2 principal components of the data:

# PCA of the normalized data:
pc = prcomp(t(log2(pmax(norm, 1))))$x[, c(1, 2)]

# run florets function:
par(mar = c(5,5,1,1))
layout(mat = (matrix(c(1, 2), 1)), widths = c(6, 2))
florets(x = pc[, 1], y = pc[, 2],
        b = restils$beta, cex = .5,
        xlab = "PC1", ylab = "PC2")
par(mar = c(0,0,0,0))
frame()
legend("center", fill = cellcols[rownames(restils$beta)], 
       legend = rownames(restils$beta), cex = 0.7)
TIL abundance plotted on PC space

TIL abundance plotted on PC space

# and plot florets in space:
par(mar = c(5,5,1,1))
layout(mat = (matrix(c(1, 2), 1)), widths = c(6, 2))
florets(x = x, y = y,
        b = restils$beta, cex = .5,
        xlab = "", ylab = "")
par(mar = c(0,0,0,0))
frame()
legend("center", fill = cellcols[rownames(restils$beta)], 
       legend = rownames(restils$beta), cex = 0.7)
TIL abundance plotted on physical space

TIL abundance plotted on physical space

So we can see that PC1 roughly tracks many vs. few immune cells, and PC2 tracks the relative abundance of lymphoid/myeloid populations.

Other functions

The SpatialDecon library includes several helpful functions for further analysis/fine-tuning of deconvolution results.

Combining cell types:

When two cell types are too similar, the estimation of their abundances becomes unstable. However, their sum can still be estimated easily. The function “collapseCellTypes” takes a deconvolution results object and collapses any colsely-related cell types you tell it to:

Cell abundance estimates with related cell types collapsed

Cell abundance estimates with related cell types collapsed

Session Info

sessionInfo()
#> R version 4.1.1 (2021-08-10)
#> Platform: x86_64-pc-linux-gnu (64-bit)
#> Running under: Ubuntu 20.04.3 LTS
#> 
#> Matrix products: default
#> BLAS:   /home/biocbuild/bbs-3.14-bioc/R/lib/libRblas.so
#> LAPACK: /home/biocbuild/bbs-3.14-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] stats4    stats     graphics  grDevices utils     datasets  methods  
#> [8] base     
#> 
#> other attached packages:
#> [1] GeomxTools_2.0.0        NanoStringNCTools_1.2.0 ggplot2_3.3.5          
#> [4] S4Vectors_0.32.0        Biobase_2.54.0          BiocGenerics_0.40.0    
#> [7] SeuratObject_4.0.2      SpatialDecon_1.4.3     
#> 
#> loaded via a namespace (and not attached):
#>  [1] httr_1.4.2             sass_0.4.0             splines_4.1.1         
#>  [4] jsonlite_1.7.2         R.utils_2.11.0         bslib_0.3.1           
#>  [7] assertthat_0.2.1       highr_0.9              GenomeInfoDbData_1.2.7
#> [10] vipor_0.4.5            cellranger_1.1.0       yaml_2.2.1            
#> [13] outliers_0.14          numDeriv_2016.8-1.1    pillar_1.6.4          
#> [16] lattice_0.20-45        glue_1.4.2             uuid_1.0-3            
#> [19] digest_0.6.28          RColorBrewer_1.1-2     XVector_0.34.0        
#> [22] minqa_1.2.4            colorspace_2.0-2       htmltools_0.5.2       
#> [25] Matrix_1.3-4           R.oo_1.24.0            plyr_1.8.6            
#> [28] pkgconfig_2.0.3        pheatmap_1.0.12        zlibbioc_1.40.0       
#> [31] purrr_0.3.4            scales_1.1.1           lme4_1.1-27.1         
#> [34] tibble_3.1.5           generics_0.1.1         IRanges_2.28.0        
#> [37] ellipsis_0.3.2         withr_2.4.2            magrittr_2.0.1        
#> [40] crayon_1.4.2           readxl_1.3.1           evaluate_0.14         
#> [43] R.methodsS3_1.8.1      fansi_0.5.0            nlme_3.1-153          
#> [46] R.cache_0.15.0         MASS_7.3-54            beeswarm_0.4.0        
#> [49] ggthemes_4.2.4         repmis_0.5             tools_4.1.1           
#> [52] data.table_1.14.2      lifecycle_1.0.1        stringr_1.4.0         
#> [55] munsell_0.5.0          EnvStats_2.4.0         Biostrings_2.62.0     
#> [58] compiler_4.1.1         jquerylib_0.1.4        GenomeInfoDb_1.30.0   
#> [61] systemfonts_1.0.3      rlang_0.4.12           nloptr_1.2.2.3        
#> [64] grid_4.1.1             RCurl_1.98-1.5         rjson_0.2.20          
#> [67] htmlwidgets_1.5.4      bitops_1.0-7           rmarkdown_2.11        
#> [70] boot_1.3-28            gtable_0.3.0           lmerTest_3.1-3        
#> [73] curl_4.3.2             DBI_1.1.1              reshape2_1.4.4        
#> [76] R6_2.5.1               knitr_1.36             dplyr_1.0.7           
#> [79] ggiraph_0.7.10         fastmap_1.1.0          utf8_1.2.2            
#> [82] stringi_1.7.5          ggbeeswarm_0.6.0       parallel_4.1.1        
#> [85] Rcpp_1.0.7             vctrs_0.3.8            tidyselect_1.1.1      
#> [88] xfun_0.27