Introduction

BLISA identifies spatially enriched ligand–receptor interactions from spatial transcriptomics data using local bivariate Moran’s I statistics.

This vignette demonstrates a typical BLISA workflow using Xenium data.

Load Package

library(blisa)
library(scider)
library(CellChat)
library(patchwork)

Load Example Data

Example Xenium dataset distributed with the package:

spe <- readRDS(
  system.file("extdata", "xenium380Cell_spe.rds", package = "blisa")
)

dim(spe)

Visualize Spatial Cell Types

cell_type_colors <- c(
  RColorBrewer::brewer.pal(7, "Set1"),
  RColorBrewer::brewer.pal(6, "Set2")
)

names(cell_type_colors) <- unique(spe$cell_type)

scider::plotSpatial(
  spe,
  group.by = "cell_type",
  pt.size = 0.2,
  pt.alpha = 0.9,
  cols = cell_type_colors
)

Filter Ligand–Receptor Pairs

Restrict CellChat ligand–receptor database to genes present in the panel:

counts_matrix <- SummarizedExperiment::assay(spe, "counts")

LR_pairs_filtered <- filterLRpairs(
  counts = counts_matrix,
  min_ligand = 1,
  min_receptor = 1,
  LR_df = CellChatDB.human$interaction
)

head(LR_pairs_filtered)

Run BLISA

system.time(BLISAres<- runBLISA.spe.isolates.removed(
    spe, 
    LR_df = LR_pairs_filtered))

BLISAres$LR_out

Summarize Ligand–Receptor Interaction Enrichment

plotLRIsum(BLISAres$LR_out)

Visualize Spatial Distribution of Top LR Pairs

plots <- lapply(seq_len(nrow(BLISAres$LR_results)), function(i) {
  plotHotspots(BLISAres, index = i)
})

combined_plot <- patchwork::wrap_plots(plots, ncol = 3)

combined_plot

pdf("xenium_LRI_spatial_plots.pdf", width = 15, height = 12)
print(combined_plot)
dev.off()

Run Cell–Cell Interaction Analysis

When blisa() is called with a group argument (e.g. group = "cell_type"), CCI scores are computed automatically and stored in BLISAres$CCI_scores. To compute or recompute them after the fact, call runCCI() on the blisa object and supply counts_by_group from hexBinCells():

# CCI already computed inside blisa() — retrieve directly
CCIres <- BLISAres$CCI_scores

# Alternatively, compute on an existing blisa object that lacks CCI_scores:
# binned  <- hexBinCells(coords, counts_matrix, bin_size = 50, group = cell_type_vec)
# BLISAres <- runCCI(BLISAres, counts_by_group = binned$counts_by_group)
# CCIres  <- BLISAres$CCI_scores

head(CCIres)

Visualize CCI Heatmap

plotCCI(BLISAres)

plotCCI(BLISAres, sender = "Invasive_Tumor", receiver = "Invasive_Tumor")

Visualize One Ligand–Receptor Pair

p <- plotCCIpair(BLISAres, ligand = "CXCL12", receptor = "CXCR4")

draw(
  p,
  heatmap_legend_side = "right",
  padding = unit(c(5, 5, 10, 10), "mm")
)

CCIspatial(
  spe,
  BLISAres,
  index = 2
)

Session Information

sessionInfo()

BLISA Workflow