SpatialTime Package
Dimitri Sokolowskei, Alex Trostle, Achira Shah, Robert Tower
SpatialTime
Spatial transcriptomics (ST) enables single cell resolution gene expression analysis assessment across histology sections while preserving native tissue morphology. In this context, SpatialTime package empower users to assess differential gene expression and module score pathways across a distance gradient by establishing spatially defined reference points and calculating ST spots distances against reference spots using Euclidean distances. Thus, a distance gradient can be established allowing the evaluation of gene expression profiling based on proximity of any reference point.
Installation
Development version install:
#(install.packages("devtools")
devtools::install_github("dimitrisokolowskei/SpatialTime")Dependencies Requirements
The igraph package version < 2.1.0 are
required for SpatialTime package. For instance:
library(remotes)
remove.packages("igraph") # remove package if previously installed
install_version("igraph", version = "2.0.3")Prerequisites
SpatialTime package requires a spatial dataset already processed and integrated. For help in generating a Seurat object, please visit Seurat vignette.
Introduction
For this vignette tutorial, we provide a ST example dataset of mice subjected to tibia osteotomy obtained from Rios et al., 2024. Additional files to run the vignette can be found at zenodo repository here.
Analysis Approaches
SpatialTime supports different strategies for delimiting spatial reference spots by: 1) manually defining references spots using the packages built-in Shiny App; 2) selecting single cell types based on cluster identities, or 3) manually defining references lines using imaging softwares such as Fiji/ImageJ. For tutorial purposes, we will first show how to select spots using the three different approaches possible. Then, we’ll present how to investigate and visualize gene expression and module score across these established distance gradients.
1. Individual Spots Selection
1.1 First, import required libraries and the processed Seurat object.
# Libraries
library(Seurat)
library(tidyverse)
library(SpatialTime)
library(shiny)
library(plotly)
library(monocle3)
library(reshape2)
# set.seed for reproducibilty
set.seed(12345)
# Import object
demo <- readRDS("fracture_demo.rds")
View([email protected])
# Visualize spatial data colored by specified identity
ratio <- ScaleRatio(demo) # Adjust image scale (OPTIONAL)
SpatialDimPlot(demo, images = "slice1") +
theme(aspect.ratio = ratio)1.2 Export XY coordinates of the spatial spot clusters of interest.
Make sure the clustering of interest is set as the object active
identity. If not, clusters labels from [email protected] can
be added into Idents(). For instance:
Idents(demo) # Check cluster names in Seurat object identity
Idents(demo) <- demo$Tissue # Change active identity to those defined in the "Tissue" meta data column
# Export cell type coordinates to an output directory. For this example, we subset only the repaired cartilage and bone tissue
subset2Labels(demo, cluster = c("Cartilage","HypCH","Woven"), export.all = F, dir.out = "coord_out")1.3 Once the coordinates from desired clusters have been exported, we next need to select our reference for gradient distance calculation. The package’s build-in Shiny app allows the visualization of individual spatial spots colored by cluster. Additionally, hoovering the mouse cursor over individual spots shows the information of which cluster that spots belongs to as well its XY coordinate.
Run Shiny app and select reference spots:
# Selecting reference spots
app <- ShinySelection(demo, id = "Tissue")
app
To navigate trough the app, use shift+left click to move
around the sample. Use the left click to select spatial
spots to be used as reference and ctrl + left click to
unselect the spots. Users can scroll with the mouse wheel or use the
+ and - signs in the top right corner to zoom
in or zoom out, respectively. In the example below, we have manually
selected spatial spots representing the central fracture line of the
bone.
For each spot selected, information on spot barcode, cluster it belongs to and XY coordinates are displayed in the Shiny App UI as below:
# First spots selected
Selected spots:
Barcode: GACCGTCAGGTCGTGA-1_1, Cluster: Cartilage, Coordinates: (379.704, -135.0951)
Barcode: CTAGTATTCGGAATTA-1_1, Cluster: Cartilage, Coordinates: (376.2911, -141.045)
Barcode: CAACTATATCGAATGC-1_1, Cluster: Cartilage, Coordinates: (372.8783, -147.0251)
...After spots selection, click Export Selected Spots and
close the Shiny app. The id parameter have to include
[email protected] column name that contains cell types
ID’s for better visualization in the shiny app.
By importing the desired cluster coordinates and the above reference
spots selected in the shiny app, SpatialTime distance gradient can be
calculated using CoordMerge function by specifying the path
of the coord files and cell cluster name. Additionally, import the
reference spots coordinates using ShinySpots function. For
example, using the cluster labeled as Woven, it
follows:
#Import coordinates of cluster and reference spotsof interest
target <- CoordMerge("coord_out/", pattern = "Woven")
# Import reference spot coordinates
ref <- ShinySpots(demo, coord_file = "selected_spots_random.csv") # image below can be generated by changing the coord_file="selected_spots.csv"
# Calculate SpatialTime between cluster and reference
```sh
sp <- SpatialShinyTime(target, ref)
spNext, we can visualize distance gradients to confirm accuracy for the following steps.
# Visualize spots gradient
gradient_vis <- SpatialVis(demo, sp, spatial.by = "rel", slice = "slice1", return_obj = F)
gradient_vis +
theme(aspect.ratio = ratio)
# Calculate gradient distance values and add Spatialtime (st) values into metadata
gradient_calc <- SpatialVis(demo, sp, spatial.by = "rel", slice = "slice1", return_obj = T)
View([email protected])
# Remove spots outside regions of interest with a calculated distance of 0 to improve visualization
gradient_calc <- subset(gradient_calc, subset = st != 0)
# Continue to in "Gene Expression Assessment & Visualization" section below2. Cluster selection based approach
As an alternative to manually reference spot selection, entire annotated clusters can be used as the reference instead.
2.1 Import Seurat’s object and ensure that cell clusters are
annotated and present in your object Idents().
# Libraries
library(Seurat)
library(tidyverse)
library(SpatialTime)
library(monocle3)
library(reshape2)
# Import NF1 data ----
demo <- readRDS("fracture_demo.rds")
View([email protected])
# Export cell type coordinates to an output directory
subset2Labels(demo, cluster = c("Cartilage","HypCH","Woven"), export.all = F, dir.out = "coord_out")2.2 Import coordinates to be used as references, as well as clusters
of interest to calculate distance in, and run Spatial2Time
for gradient distance calculation:
# Selecting clusters of interest ----
target <- CoordMerge("coord_out", pattern = "Woven")
# Selecting clusters to be used as reference ----
ref <- CoordMerge("coord_out", pattern = "Cartilage")
# Calculate SpatialTime between clusters of interest ----
sp <- Spatial2Time(target, ref)
sp
2.3 Visualize gradient distance:
# Distance gradient visualization ----
gradient_vis <- SpatialVis(demo, sp, spatial.by = "rel", slice = "slice1", return_obj = F)
ratio <- ScaleRatio(demo) # Adjust image scale
gradient_vis +
theme(aspect.ratio = ratio)
# Add st values into metadata
gradient_calc <- SpatialVis(demo, sp, spatial.by = "rel", slice = "slice1", return_obj = T)
View([email protected])
# Subset out spots not selected as target tissue where spatialtime (st) value were not calculated
gradient_calc <- subset(gradient_calc, subset = st != 0)
# Continue to "Gene Expression Assessment & Visualization" section below3. Fiji/ImageJ based approach
Fiji is a broadly used open source scientific image processing tool. Fiji allows drawing the reference lines of interest and exporting each delineated pixel as an XY coordinates later used for euclidean distance calculation. Fiji can be downloaded here and a brief Fiji tutorial on how draw and export references lines can be found here.
3.1 First, import your Seurat ST object using readRDS()
function. Then, we proceed using subsetLabels() which
allows us to export ST clusters coordinates as excel sheets.
# Required packages
library(Seurat)
library(tidyverse)
library(SpatialTime)
library(monocle3)
library(reshape2)
# Import dataset ----
demo <- readRDS("fracture_demo.rds")
View([email protected])
# Export clusters coordinates
subsetLabels("fracture_demo.rds", cluster = c("Cartilage","Woven"), slice.n = "slice1", dir.out = "coord_out")3.2 Import Fiji XY coordinates previously prepared. Furthermore, we
need to match the reference line colors to their respective
tissue/region names. Be aware that the tissue names must exist in your
Seurat object identities. Notice that each slice image have its own
factor value that can be found in the .json file in each
10x output folder. Then, by running SpatialCalc() each XY
coordinate can be transformed accordingly to the slice factor.
fiji <- read.csv("Fiji_C1.csv")
colors <- c("Red", "Blue")
tissue <- c("Cartilage","Woven")
data <- SpatialCalc("Fiji_C1.csv", factor = 0.10067452, colors = colors, tissue = tissue)3.3 At this step, SpatialTime distance are ready to be calculated. First, import your cluster type coordinates as reference and calculate the spots distances against other cluster.
# Calculate minimum spots distances from reference line of interest
file <- read.csv("coord_out/Woven_coordinates.csv")
sp <- SpatialTime(data, reference = "Cartilage", compare = file)
sp3.4 SpatialTime distance gradient can be assessed using:
gradient_vis <- SpatialVis(demo, sp, spatial.by = "rel", slice = "slice1", return_obj = F)
gradient_vis
3.5 Adding
spatialtime values into metadata and subsetting only cells with st
values:
# Visualization spatialtime gradient ----
gradient_calc <- SpatialVis(demo, sp, spatial.by = "rel", slice = "slice1", return_obj = T)
gradient_calc <- subset(gradient_calc, subset = st != 0)Gene Expression Assessement & Visualization
In this step, independently of spots selection approach use, st
values should be included in [email protected]. Consequently,
pseudotime can be calculated and gene expression analysis visualized as
heatmap across distance gradient. For example:
# Calculate pseudotime and visualize genes differentialy expressed across distance gradient
ps <- PseudoM3Time(gradient_calc, assay = "SCT", values = "st", return_obj = F)
hm <- gradient_heatmap(ps, rownames = F)
Genes expressed
across heatmap clusters could be obtained by:
# Obtain gene names in each heatmap clusters
genes <- GeneGet(hm, n_clusters = 3)
genes
Cluster Gene
Jph1 1 Jph1
Mcm3 2 Mcm3
Col9a1 3 Col9a1
Ptpn18 2 Ptpn18
Ankrd23 1 Ankrd23
...By using return_obj = T Seurat object is returned
allowing for individuals gene expression visualization across distances.
Additionally, module score pathways can be assessed as follows:
# Calculate pseudotime
ps <- PseudoM3Time(gradient_calc, assay = "SCT", values = "st", return_obj = T)
GeneVis(file = ps, column = c("Col2a1", "Ihh","Sox9"), span = 0.5, signal = "gene")
# Gene list pathways
bmp <- c("Acvr1","Amh","Bmp2","Bmp4","Bmp5","Bmp6","Bmp7","Bmp8a","Bmpr1a","Bmpr1b",
"Bmpr2","Gdf5","Gdf6","Gdf7","Hamp","Hamp2","Id1","Id2","Id3","Id4",
"Inhbb","Rgma","Rgmb","Smad1","Smad4","Smad5","Smad9")
wnt <- c("Rspo3","Rspo4","Tcf7","Tcf7l1","Tcf7l2","Wnt1","Wnt10a","Wnt10b","Wnt11",
"Wnt16","Wnt2","Wnt2b","Wnt3","Wnt3a","Wnt4","Wnt5a","Wnt6","Wnt7a","Wnt7b",
"Wnt8a","Wnt8b","Wnt9a","Wnt9b","Wnt5b","Axin1","Axin2")
# Add gene modules in object with pseudotime calculated
gene_modules <- list(Module1 = bmp, Module2 = wnt)
module <- AddModuleScore(gradient_calc, features = gene_modules, assay = "SCT", name = "module", nbin = 16)
GeneVis(file = module, column = c("module1", "module2"), signal = "pathway", span = 0.5)
Spatial Transcriptomics Visium HD
SpatialTime package also supports 10x genomics Visium HD analysis allowing user to calculate distance gradients in a true single cell resolution. The SpatialTime rationale for HD assays remains with minor modifications from standard cluster-based or image-based approaches. For instance:
library(Seurat)
library(tidyverse)
library(SpatialTime)
library(monocle3)
library(reshape2)
mm <- readRDS("path/seurat_object.rds")
subset2LabelsHD(mm, cluster = c("CellType1", "CellType2"))
ref <- CoordMergeHD("coord_out", pattern = "CellType1")
other <- CoordMergeHD("coord_out", pattern = "CellType2")
sp <- Spatial2TimeHD(other, ref)
gradient_vis <- SpatialVis(mm, sp, spatial.by = "rel", slice = "slice1", return_obj = F)
gradient_vis
...For Fiji image-based analysis:
library(Seurat)
library(tidyverse)
library(SpatialTime)
library(monocle3)
library(reshape2)
mm <- readRDS("path/seurat_object.rds")
colors <- c("Red")
tissue <- c("CellType1")
data <- SpatialCalc("fiji_output.csv", factor = 1, colors = colors, tissue = tissue)
file <- read.csv("coord_out/{cluster_name}_slice1_coordinates.csv")
sp <- SpatialTimeHD(data, reference = "CellType1", compare = file)
sp
gradient_vis <- SpatialVis(mm, sp, spatial.by = "rel", slice = "slice1", return_obj = F)
gradient_vis
...Monocle Usage
SpatialTime package was originally conceived using
monocle for spatialvalues calculation mainly to its
reliability and heatmap plotting facility. However, monocle’s newer
versions such as monocle3 brought series of performance and
analysis enhancements. Thus, we additionally provide a monocle wrapper
function that allows equivalent spatialtime calculation in the standard
SpatialTime pipeline. For instance:
library(monocle)
ps <- Pseudo2Time(gradient_calc, assay = "SCT", return_obj = F)
p <- plot_pseudotime_heatmap(ps, num_clusters = 2, cores = 2, show_rownames = F, return_heatmap = T)
p
For downstream analysis,
such as module scores, object seurat_obj$meta.data should
be returned instead of a expression matrix. Thus, set
return_obj = T:
ps <- Pseudo2Time(fracture, assay = "SCT", return_obj = T)Session Info
R version 4.3.2 (2023-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19045)
Matrix products: default
locale:
[1] LC_COLLATE=English_United States.utf8 LC_CTYPE=English_United States.utf8 LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C LC_TIME=English_United States.utf8
time zone: America/Chicago
tzcode source: internal
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] SpatialTime_0.0.0.9000 testthat_3.2.0 devtools_2.4.5 usethis_2.2.2
loaded via a namespace (and not attached):
[1] RcppAnnoy_0.0.21 splines_4.3.2 later_1.3.1 DDRTree_0.1.5 tibble_3.2.1
[6] polyclip_1.10-6 rpart_4.1.21 lifecycle_1.0.4 rprojroot_2.0.4 globals_0.16.2
[11] processx_3.8.2 lattice_0.21-9 MASS_7.3-60 backports_1.4.1 magrittr_2.0.3
[16] limma_3.58.1 Hmisc_5.1-1 plotly_4.10.3 rmarkdown_2.28 yaml_2.3.7
[21] remotes_2.4.2.1 httpuv_1.6.12 Seurat_4.2.1 sctransform_0.4.1 askpass_1.2.0
[26] spam_2.10-0 spatstat.sparse_3.0-3 sp_2.1-1 sessioninfo_1.2.2 pkgbuild_1.4.2
[31] reticulate_1.34.0 cowplot_1.1.1 pbapply_1.7-2 RColorBrewer_1.1-3 abind_1.4-5
[36] pkgload_1.3.3 Rtsne_0.16 purrr_1.0.2 BiocGenerics_0.48.1 nnet_7.3-19
[41] ggrepel_0.9.4 irlba_2.3.5.1 listenv_0.9.0 spatstat.utils_3.1-0 pheatmap_1.0.12
[46] umap_0.2.10.0 RSpectra_0.16-1 goftest_1.2-3 spatstat.random_3.2-1 fitdistrplus_1.1-11
[51] parallelly_1.36.0 commonmark_1.9.0 leiden_0.4.3 codetools_0.2-19 xml2_1.3.5
[56] tidyselect_1.2.0 viridis_0.6.5 matrixStats_1.1.0 stats4_4.3.2 base64enc_0.1-3
[61] spatstat.explore_3.2-5 roxygen2_7.3.2 jsonlite_1.8.7 ellipsis_0.3.2 progressr_0.14.0
[66] Formula_1.2-5 ggridges_0.5.4 survival_3.5-7 HSMMSingleCell_1.22.0 tools_4.3.2
[71] ica_1.0-3 Rcpp_1.0.11 glue_1.6.2 gridExtra_2.3 xfun_0.44
[76] dplyr_1.1.4 withr_2.5.2 combinat_0.0-8 fastmap_1.1.1 xopen_1.0.0
[81] fansi_1.0.5 openssl_2.1.1 callr_3.7.3 digest_0.6.33 rcmdcheck_1.4.0
[86] R6_2.5.1 mime_0.12 colorspace_2.1-0 scattermore_1.2 tensor_1.5
[91] spatstat.data_3.0-3 monocle_2.30.1 utf8_1.2.4 tidyr_1.3.0 generics_0.1.3
[96] data.table_1.14.8 prettyunits_1.2.0 httr_1.4.7 htmlwidgets_1.6.2 uwot_0.1.16
[101] pkgconfig_2.0.3 gtable_0.3.4 lmtest_0.9-40 brio_1.1.3 htmltools_0.5.7
[106] profvis_0.3.8 dotCall64_1.1-0 tidyverse_2.0.0 SeuratObject_5.0.2 scales_1.3.0
[111] Biobase_2.62.0 png_0.1-8 knitr_1.48 rstudioapi_0.15.0 tzdb_0.4.0
[116] reshape2_1.4.4 curl_5.1.0 nlme_3.1-163 checkmate_2.3.0 cachem_1.0.8
[121] zoo_1.8-12 stringr_1.5.0 KernSmooth_2.23-22 parallel_4.3.2 miniUI_0.1.1.1
[126] foreign_0.8-85 desc_1.4.2 pillar_1.9.0 grid_4.3.2 fastICA_1.2-3
[131] vctrs_0.6.4 RANN_2.6.1 slam_0.1-50 urlchecker_1.0.1 VGAM_1.1-9
[136] promises_1.2.1 xtable_1.8-4 cluster_2.1.4 waldo_0.5.2 htmlTable_2.4.2
[141] evaluate_0.23 readr_2.1.4 cli_3.6.1 compiler_4.3.2 rlang_1.1.2
[146] crayon_1.5.2 leidenbase_0.1.25 future.apply_1.11.0 ps_1.7.5 plyr_1.8.9
[151] fs_1.6.3 stringi_1.7.12 deldir_1.0-9 viridisLite_0.4.2 munsell_0.5.0
[156] lazyeval_0.2.2 spatstat.geom_3.2-7 Matrix_1.6-4 hms_1.1.3 patchwork_1.2.0.9000
[161] future_1.33.0 ggplot2_3.5.0 statmod_1.5.0 shiny_1.7.5.1 ROCR_1.0-11
[166] igraph_2.0.3 memoise_2.0.1R version 4.3.2 (2023-10-31 ucrt)