Package: PathwaySpace 1.0.3.1

Overview

This vignette introduces the SpotSpace package, a PathwaySpace extension for the Seurat package (Hao et al. 2024), providing methods for signal propagation and visualization in spatial transcriptomics. SpotSpace extends existing spatial analysis workflows to explore signal patterns in tissue microenvironments. In what follows, we present three step-by-step tutorials describing how to prepare input data for SpotSpace. The results reproduce and refine examples featured in Seurat’s tutorials, so users are encouraged to see how these packages can be used together.

Before you start

This vignette assumes prior experience with Seurat (Hao et al. 2024), especially for handling spatial transcriptomics data.

Note: If you are new to Seurat’s spatial workflows, we recommend reviewing the spatial analysis tutorials before continuing.

Computational requirement:

  • Hardware: RAM >= 16 GB

  • Software: R (>=4.4) and RStudio

Required packages

# Check required packages for this vignette
if (!require("remotes", quietly = TRUE)){
  install.packages("remotes")
}
if (!require("RGraphSpace", quietly = TRUE)){
  remotes::install_github("sysbiolab/RGraphSpace")
}
if (!require("PathwaySpace", quietly = TRUE)){
  remotes::install_github("sysbiolab/PathwaySpace")
}
if (!require("SpotSpace", quietly = TRUE)){
  remotes::install_github("sysbiolab/SpotSpace")
}
if (!require("SeuratData", quietly = TRUE)){
  remotes::install_github("satijalab/seurat-data")
}
# Check versions
if (packageVersion("RGraphSpace") < "1.1.0"){
  message("Need to update 'RGraphSpace' for this vignette")
  remotes::install_github("sysbiolab/RGraphSpace")
}
if (packageVersion("PathwaySpace") < "1.0.3.1"){
  message("Need to update 'PathwaySpace' for this vignette")
  remotes::install_github("sysbiolab/PathwaySpace")
}
if (packageVersion("SpotSpace") < "0.0.2"){
  message("Need to update 'SpotSpace' for this vignette")
  remotes::install_github("sysbiolab/SpotSpace")
}
if (packageVersion("Seurat") < "5.3.1.9999"){
  message("Need to update 'Seurat' for this vignette")
  remotes::install_github("satijalab/Seurat")
}
# Load packages
library(RGraphSpace)
library(PathwaySpace)
library(SpotSpace)
library(Seurat)
library(SeuratObject)
library(SeuratData)
library(patchwork)

Visium v1 dataset

Setting input data

For this tutorial, we will use the stxBrain dataset from the SeuratData package, consisting of spatial transcriptomics data from sagittal mouse brain sections generated with Visium v1 technology. This dataset is commonly used to demonstrate Seurat spatial workflows (Hao et al. 2024). Here, we will preprocess it with Seurat and then extract the relevant data for PathwaySpace downstream analyses.

## Install a Seurat dataset (this step is required only once)
SeuratData::InstallData("stxBrain")
# Check manifest of installed datasets
# SeuratData::InstalledData()

# Load the 'stxBrain' dataset
seurat_obj <- LoadData("stxBrain", type = "anterior1")

The stxBrain dataset is normalized as suggested in Seurat’s spatial_vignette, either using the SCTransform() and NormalizeData() functions.

# Run vst normalization on counts
# seurat_obj <- SCTransform(seurat_obj, assay = "Spatial", verbose = FALSE)

# NOTE: Seurat recommends using SCTransform() for processing this 
# spatial dataset, which may require more computation time. Here,
# we use log-normalization for demonstration purposes.
seurat_obj <- NormalizeData(seurat_obj)

… and then we extract spot coordinates, tissue image, and vst-normalized data.

# Get spot coordinates
spot_coord <- GetTissueCoordinates(seurat_obj, scale = "lowres")

# Get raster image
raster_image <- GetImage(seurat_obj, "raster")

# Get vst-normalized gene expression
vst_gexp <- GetAssayData(seurat_obj, layer="data")

# If needed, remove seurat_obj to free memory
rm(seurat_obj)

Next, we create a PathwaySpace object from the spot coordinates and plot the resulting graph overlaid on the tissue image.

# Create a PathwaySpace object from 'spot_coord', mapped to the 'raster_image'
pspace_obj <- buildSpotSpace(spot_coord = spot_coord, raster_image = raster_image)
# Check spots on the top of the 'raster_image'
xy_labs <- labs(x="Spot coordinates 1", y="Spot coordinates 2")
plotGraphSpace(pspace_obj, add.image = TRUE) + xy_labs

Running PathwaySpace

Before projection, we need to specify a distance unit for the signal decay function. This distance unit will affect the extent over which the convolution operation projects the signal, scaled to the coordinate space. We will use the center-to-center distance between spots, which represents 100 µm in the Visium v1 technology.

# Get distance to the nearest spot
nspot <- getNearestNode(pspace_obj)
pdist <- mean(nspot$dist) # average distance
# 'pdist' set as the average center-to-center distance between spots
pdist
# [1] 0.01395144

As an optional step, the silhouetteMapping() function generates an image mask that outlines the graph layout, over which the subsequent methods will project a landscape image. The baseline argument controls the level at which a silhouette is sliced to form the mask. Increasing the baseline (in [0,1]) produces a more detailed, granular silhouette.

# Add a graph silhouette to the PathwaySpace object
pspace_obj <- silhouetteMapping(pspace_obj, baseline = 0.1)
plotPathwaySpace(ps=pspace_obj, theme = "th3", 
  add.image = TRUE, si.alpha = 0.5) + xy_labs

Next, we specify the signal to be projected; for this demonstration, we will use expression data from the Camk2n1 gene. The vertexSignal() accessor function is then used to assign the gene expression values to graph vertices.

# Select a gene of interest (e.g., Camk2n1) and assign its 
# expression values to graph vertices
gene <- "Camk2n1"
vertexSignal(pspace_obj)[colnames(vst_gexp)] <- vst_gexp[gene,]

We then perform the signal projection, setting decay = 0.5. The decay parameter controls how the signal attenuates as a function of distance in pathway space. With decay = 0.5, the signal decreases to half of its initial value at a distance equal to pdist (for additional configuration details, see the modeling signal decay tutorial).

# Project gene signal
pspace_obj <- circularProjection(pspace_obj, k = gs_vcount(pspace_obj), 
  decay.fun = weibullDecay(decay=0.5, pdist = pdist),
  aggregate.fun = signalAggregation("wmean"))

Because each spot produces an independent projection, the resulting projections are aggregated into a unified landscape. Here we use a weighted arithmetic mean, with each projection weighted by its own magnitude (for additional configuration details, see the signal aggregation rules tutorial).

Next, we show the results with minor variations to demonstrate some of the available plot settings.

# Plot tissue image and projection separated 
p1 <- plotPathwaySpace(ps=pspace_obj, theme = "th3", title = gene)
p1$image <- p1$image + xy_labs
p1$graph <- p1$graph + xy_labs

p1$image + p1$graph

# Plot projections overlaid on the tissue image, with alpha = 0.25
p2 <- plotPathwaySpace(ps=pspace_obj, theme = "th3", title = gene, 
  add.image = TRUE, si.alpha = 0.25) + xy_labs

# Plot projections overlaid on the tissue image, with zlim truncated at >=0.5
p3 <- plotPathwaySpace(ps=pspace_obj, theme = "th3", title = gene, 
  add.image = TRUE, si.alpha = 0.25, zlim = c(0.5, 1)) + xy_labs

p2 + p3


Slide-seq v2 dataset

Setting input data

For this tutorial, we will use the ssHippo dataset available from the SeuratData package, consisting of spatial transcriptomics data from mouse hippocampus generated with Slide-seq v2 technology. We will follow the same general steps from our previous spatial tutorial, preprocessing with Seurat and then extracting the relevant data for PathwaySpace downstream analyses. For further details on this dataset, see Seurat’s spatial_vignette.

## Install a Seurat dataset (this step is required only once)
SeuratData::InstallData("ssHippo")
# Check manifest of installed datasets
# SeuratData::InstalledData()

# Load the 'kidneyref' dataset
seurat_obj <- LoadData("ssHippo")
# Run vst normalization on counts
# seurat_obj <- SCTransform(seurat_obj, assay = "Spatial", verbose = FALSE)

# NOTE: Seurat recommends using SCTransform() for processing this 
# spatial dataset, which may require more computation time. Here,
# we use log-normalization for demonstration purposes.
seurat_obj <- NormalizeData(seurat_obj)
# Extract spot coordinates and vst-normalized data
# Get spot coordinates
spot_coord <- GetTissueCoordinates(seurat_obj)

#Note: the `ssHippo` dataset does not include a tissue image

# Get vst-normalized gene expression
vst_gexp <- GetAssayData(seurat_obj, layer="data")

# If needed, remove seurat_obj to free memory
rm(seurat_obj)
# Create a PathwaySpace object from 'spot_coord'
pspace_obj <- buildSpotSpace(spot_coord = spot_coord, nrc=700)

Running PathwaySpace

# Get distance to the nearest spot
nspot <- getNearestNode(pspace_obj)
pdist <- mean(nspot$dist) # average distance
# 'pdist' set as the average center-to-center distance between spots
pdist
# [1] 0.002497466
# Add a graph silhouette to the PathwaySpace object
pspace_obj <- silhouetteMapping(pspace_obj, fill.cavity = FALSE, 
  pdist = max(nspot$dist))

# Check silhouette plot
xy_labs <- labs(x="Spot coordinates 1", y="Spot coordinates 2")
plotPathwaySpace(ps=pspace_obj, theme = "th3", 
  add.image = TRUE, si.alpha = 0.5) + xy_labs

# Choose a gene of interest (e.g., DDN) and assign its 
# expression values to graph vertices
gene <- "DDN"
vertexSignal(pspace_obj)[colnames(vst_gexp)] <- vst_gexp[gene,]

# Project gene signal
pspace_obj <- circularProjection(pspace_obj, k = gs_vcount(pspace_obj), 
  decay.fun = weibullDecay(decay=0.5, pdist = pdist))

# Plot projections
#-- as a suggestion, truncate zlim at the upper limit 
#-- to enhance certain patters
p1 <- plotPathwaySpace(ps=pspace_obj, theme = "th3", 
  title = gene, zlim = c(0, 1)) + xy_labs
# ...another gene (e.g. PCP4)
gene <- "PCP4"
vertexSignal(pspace_obj)[colnames(vst_gexp)] <- vst_gexp[gene,]

# Project gene signal
pspace_obj <- circularProjection(pspace_obj, k = gs_vcount(pspace_obj), 
  decay.fun = weibullDecay(decay=0.5, pdist = pdist))

# Plot projections
p2 <- plotPathwaySpace(ps=pspace_obj, theme = "th3", 
  title = gene, zlim = c(0, 1)) + xy_labs
p1 + p2

Visium HD dataset

Setting input data

Here, we will use a higher-resolution spatial dataset from mouse brain generated with Visium HD technology. This platform provides whole-transcriptome gene expression data at a raw 2-µm resolution, with additional binned versions available at 8 and 16 µm. For this tutorial, we will use the 16-µm binned data. We will follow the same general steps from our previous spatial tutorials, preprocessing with Seurat and then extracting the relevant data for PathwaySpace downstream analyses. For additional details on this dataset, refer to Seurat’s visiumhd_analysis_vignette.

The Visium HD dataset can be downloaded from the 10x Genomics repository:

# Extract the tar.gz and set 'localdir' to the dataset folder
# Use 'bin.size' to choose the data resolution to load (2, 8, or 16 µm)
localdir <- "path/to/data/directory"
seurat_obj <- Load10X_Spatial(data.dir = localdir, bin.size = 16)

# Check default assay
Assays(seurat_obj)
# [1] "Spatial.016um"
# Run log-normalization for spatial data
seurat_obj <- NormalizeData(seurat_obj)
# Get spot coordinates
spot_coord <- GetTissueCoordinates(seurat_obj, scale = "lowres")

# Get raster image
raster_image <- GetImage(seurat_obj, "raster")

# Get normalized gene expression data
norm_gexp <- GetAssayData(seurat_obj, layer="data")

# If needed, remove seurat_obj to free memory
rm(seurat_obj)
# Create a PathwaySpace object from 'spot_coord', mapped to the 'raster_image'
pspace_obj <- buildSpotSpace(spot_coord = spot_coord, 
  raster_image = raster_image, nrc = 700)

Running PathwaySpace

# Get distance to the nearest spot
nspot <- getNearestNode(pspace_obj)
pdist <- mean(nspot$dist) # average distance
# 'pdist' set as the average center-to-center distance between spots
pdist
# [1] 0.002516144
# Add a graph silhouette to the PathwaySpace object
pspace_obj <- silhouetteMapping(pspace_obj, fill.cavity = FALSE, 
  pdist = max(nspot$dist))

# Check silhouette plot
xy_labs <- labs(x="Spot coordinates 1", y="Spot coordinates 2")
plotPathwaySpace(ps=pspace_obj, theme = "th3", 
  add.image = TRUE, si.alpha = 0.5) + xy_labs

# Choose a gene of interest (e.g., Rorb) and assign its 
# expression values to graph vertices
gene <- "Rorb"
vertexSignal(pspace_obj)[colnames(norm_gexp)] <- norm_gexp[gene,]

# Project gene signal
pspace_obj <- circularProjection(pspace_obj, k = gs_vcount(pspace_obj), 
  decay.fun = weibullDecay(decay=0.5, pdist = pdist))

# Plot projections
p1 <- plotPathwaySpace(ps=pspace_obj, theme = "th3", 
  title = gene, add.image = TRUE, zlim = c(0, 1)) + xy_labs
# ...another gene (e.g. Hpca)
gene <- "Hpca"
vertexSignal(pspace_obj)[colnames(norm_gexp)] <- norm_gexp[gene,]

# Project gene signal
pspace_obj <- circularProjection(pspace_obj, k = gs_vcount(pspace_obj), 
  decay.fun = weibullDecay(decay=0.5, pdist = pdist))

# Plot projections
p2 <- plotPathwaySpace(ps=pspace_obj, theme = "th3", 
  title = gene, add.image = TRUE, zlim = c(0, 1)) + xy_labs
p1 + p2

Citation

If you use PathwaySpace, please cite:

Session information

## R version 4.5.1 (2025-06-13)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.3 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so;  LAPACK version 3.12.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
##  [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/Sao_Paulo
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] patchwork_1.3.2                Seurat_5.3.1.9999             
##  [3] SeuratObject_5.2.0             sp_2.2-0                      
##  [5] stxBrain.SeuratData_0.1.2      ssHippo.SeuratData_3.1.4      
##  [7] kidneyref.SeuratData_1.0.2     bonemarrowref.SeuratData_1.0.0
##  [9] SeuratData_0.2.2.9002          SpotSpace_0.0.2               
## [11] PathwaySpace_1.0.3.1           RGraphSpace_1.1.0             
## [13] ggplot2_4.0.0.9000             remotes_2.5.0                 
## [15] fontawesome_0.5.3             
## 
## loaded via a namespace (and not attached):
##   [1] deldir_2.0-4           pbapply_1.7-4          gridExtra_2.3         
##   [4] rlang_1.1.6            magrittr_2.0.4         RcppAnnoy_0.0.22      
##   [7] spatstat.geom_3.6-0    matrixStats_1.5.0      ggridges_0.5.7        
##  [10] compiler_4.5.1         png_0.1-8              vctrs_0.6.5           
##  [13] reshape2_1.4.4         stringr_1.5.2          crayon_1.5.3          
##  [16] pkgconfig_2.0.3        fastmap_1.2.0          promises_1.3.3        
##  [19] rmarkdown_2.30         purrr_1.1.0            xfun_0.53             
##  [22] cachem_1.1.0           jsonlite_2.0.0         goftest_1.2-3         
##  [25] later_1.4.4            spatstat.utils_3.2-0   irlba_2.3.5.1         
##  [28] parallel_4.5.1         cluster_2.1.8.1        R6_2.6.1              
##  [31] ica_1.0-3              spatstat.data_3.1-8    stringi_1.8.7         
##  [34] bslib_0.9.0            RColorBrewer_1.1-3     reticulate_1.43.0     
##  [37] spatstat.univar_3.1-4  parallelly_1.45.1      lmtest_0.9-40         
##  [40] jquerylib_0.1.4        scattermore_1.2        Rcpp_1.1.0            
##  [43] knitr_1.50             tensor_1.5.1           future.apply_1.20.0   
##  [46] zoo_1.8-14             sctransform_0.4.2      httpuv_1.6.16         
##  [49] Matrix_1.7-4           splines_4.5.1          igraph_2.2.0          
##  [52] tidyselect_1.2.1       abind_1.4-8            rstudioapi_0.17.1     
##  [55] yaml_2.3.10            spatstat.random_3.4-2  spatstat.explore_3.5-3
##  [58] codetools_0.2-20       miniUI_0.1.2           listenv_0.9.1         
##  [61] lattice_0.22-5         tibble_3.3.0           plyr_1.8.9            
##  [64] shiny_1.11.1           withr_3.0.2            S7_0.2.0              
##  [67] ROCR_1.0-11            evaluate_1.0.5         Rtsne_0.17            
##  [70] future_1.67.0          fastDummies_1.7.5      survival_3.8-3        
##  [73] polyclip_1.10-7        fitdistrplus_1.2-4     pillar_1.11.1         
##  [76] KernSmooth_2.23-26     plotly_4.11.0          generics_0.1.4        
##  [79] RcppHNSW_0.6.0         scales_1.4.0           globals_0.18.0        
##  [82] xtable_1.8-4           glue_1.8.0             lazyeval_0.2.2        
##  [85] tools_4.5.1            data.table_1.17.8      RSpectra_0.16-2       
##  [88] RANN_2.6.2             dotCall64_1.2          cowplot_1.2.0         
##  [91] grid_4.5.1             tidyr_1.3.1            nlme_3.1-168          
##  [94] cli_3.6.5              rappdirs_0.3.3         spatstat.sparse_3.1-0 
##  [97] spam_2.11-1            viridisLite_0.4.2      dplyr_1.1.4           
## [100] uwot_0.2.3             gtable_0.3.6           sass_0.4.10           
## [103] digest_0.6.37          progressr_0.17.0       ggrepel_0.9.6         
## [106] htmlwidgets_1.6.4      farver_2.1.2           htmltools_0.5.8.1     
## [109] lifecycle_1.0.4        httr_1.4.7             mime_0.13             
## [112] MASS_7.3-65

References

Hao, Yuhan, Tim Stuart, Madeline H Kowalski, Saket Choudhary, Paul Hoffman, Austin Hartman, Avi Srivastava, et al. 2024. “Dictionary Learning for Integrative, Multimodal and Scalable Single-Cell Analysis.” Nature Biotechnology 42 (2): 293–304. https://doi.org/10.1038/s41587-023-01767-y.