Visium profile of human breast tumour
Introduction
This vignette walks through a typical workflow on a Visium breast
cancer dataset.
Load Package
library(blisa)
library(SpatialExperiment)
library(ggplot2)
Load Example Data
The example dataset is a small SpatialExperiment object
(one Visium breast cancer section) hosted as a GitHub Release asset.
Download it once and cache locally:
data_url <- "https://github.com/ChenLaboratory/example_data/releases/download/v1.0.0/spe_visium_bc_s1.rds"
cache_dir <- tools::R_user_dir("blisa", "cache")
data_file <- file.path(cache_dir, "spe_visium_bc_s1.rds")
if (!file.exists(data_file)) {
dir.create(cache_dir, recursive = TRUE, showWarnings = FALSE)
download.file(data_url, data_file, mode = "wb")
}
spe <- readRDS(data_file)
spe
#> class: SpatialExperiment
#> dim: 18085 4992
#> metadata(2): um_per_pixel coord_unit
#> assays(1): counts
#> rownames(18085): SAMD11 NOC2L ... MT-ND6 MT-CYB
#> rowData names(3): ID Symbol Type
#> colnames(4992): AACACCTACTATCGAA-1 AACACGTGCATCGCAC-1 ... TGTTGGCCAGACCTAC-1 TGTTGGCCTACACGTG-1
#> colData names(6): in_tissue array_row ... n_genes sample_id
#> reducedDimNames(0):
#> mainExpName: NULL
#> altExpNames(0):
#> spatialCoords names(2) : pxl_col_in_fullres pxl_row_in_fullres
#> imgData names(4): sample_id image_id data scaleFactor
Quality control
For a Visium data, there is no need to run hexBinCells()
to bin the data. Here we remove spots with less than 3000 total
counts.
kp <- Matrix::colSums(counts(spe)) > 3000
table(kp)
#> kp
#> FALSE TRUE
#> 580 4412
spe <- spe[,kp]
Run BLISA at the LR level
blisa() does everything in one call: spatial weights, LR
pair filtering against CellChatDB, bivariate Moran’s I per LR pair, and
hotspot identification. However, as Visium data is spot-based, there is
no cell type annotation by default. Therefore, blisa()
doesn’t calculate cell-cell interaction scoring as for imaging-based
data such as Xenium.
res <- blisa(spe)
#> Testing 462 LR pairs...
#> |========================================| 100%
res
#> A blisa object
#> Level : lr
#> LR pairs tested : 462
#> Significant : 462
#> Bins : 4412
#> CCI computed : FALSE
The result is a blisa object with five slots:
Level - whether it is LR level or pathway levelLR_results — one row per LR pair, sorted by number of
hotspot binsbins — the hexagonal grid as an sf
objectspatial_weights — queen and distance-decay weights from
computeSpatialWeights()CCI_scores — wide data frame of sender-receiver
interaction scores per LR pair
Run BLISA at the pathway level
To perform pathway level cell-cell communication analyses, we run
blisaPathway().
res_pathway <- blisaPathway(res)
res_pathway
#> A blisa object
#> Level : pathway
#> Pathways tested : 100
#> Significant : 100
#> Bins : 4412
#> CCI computed : FALSE
By default, the pathway level result is the union of each LR pair’s
significant hotspot spots, with the minimum (best) p-value per spot. Two
alternative methods, Simes’ and Fisher’s methods, are supported to
combine p-values of all LR pair’s across hotspot spots.
Visualization
The plotLRrank function visualizes the top interacting
pathways across the tissue, ranked by the number of significant hotspot
spots.
plotLRrank(res_pathway, flip=TRUE)
Spatial Map of Hotspot Bins
For a chosen pathway (here, the top-ranked pathway by default),
plotHotspots shows which bins are significant hotspots,
coloured by p-value. The H&E image can be shown in the backgound
using plotImage from scider.
p_img <- scider::plotImage(spe)
plotHotspots(res_pathway, index = 1, background = p_img, size = 1.2)
