Contents

library(stJoincount)
library(pheatmap)
library(ggplot2)

v0.99.16

stJoincount: Quantification tool for spatial correlation between clusters in spatial transcriptomics preprocessed data using join count statistic approach.

Introduction

Spatial dependency is the relationship between location and attribute similarity. The measure reflects whether an attribute of a variable observed at one location is independent of values observed at neighboring locations. Positive spatial dependency exists when neighboring attributes are more similar than what could be explained by chance. Likewise, a negative spatial dependency is reflected by a dissimilarity of neighboring attributes. Join count analysis allows for quantification of the spatial dependencies of nominal data in an arrangement of spatially adjacent polygons.

This tool requires data produced with the 10X Genomics Visium Spatial Gene Expression platform with customized clusters. The purpose of this R package is to perform join count analysis for spatial correlation between clusters.

Installation

Users can install stJoincount with:

if (!requireNamespace("BiocManager", quietly = TRUE)) {
      install.packages("BiocManager")
  }
BiocManager::install("stJoincount")

Examples of how to run this tool are below:

Preprocessing

In this vignette, we are going to use an human breast cancer spatial transcriptomics sample.

fpath <- system.file("extdata", "dataframe.rda", package="stJoincount")
load(fpath)
head(humanBC)
#>                    imagecol imagerow Cluster
#> AATTGCAGCAATCGAC-1 431.2129 476.8069       4
#> ACCAGGAGTGTGATCT-1 273.0446 117.8218       9
#> ACCTCCGCCCTCGCTG-1 448.2178 423.9109       7
#> AGGTGTATCGCCATGA-1 144.2822 317.5000       1
#> ATAGTTCCACCCACTC-1 431.5099 323.9109       7
#> CCGTATTAGCGCAGTT-1 462.1535 200.4950       2

Within the ‘extdata’ user can find a dataframe “humanBC.rda”. This example data is a data.frame that comes from a Seurat object of a human breast cancer sample. It contains the following information that is essential to this algorithm - barcode (index), cluster (they could either be categorical or numerical labels), central pixel location (imagerow and imagecol). This dataframe is simplified after combining metadata with spatial coordinates. The index contains barcodes, and at least three other columns that have these information are required and the column names should be the same as following: imagerow: The row pixel coordinate of the center of the spot imagecol: The column pixel coordinate of the center of the spot Cluster: The label that corresponding to this barcode

The following codes demonstrate how to generate the described data.frame from Seurat/spatialExperiment Objects.

An example data preparation from Seurat:

fpath <- system.file("extdata", "SeuratBC.rda", package="stJoincount")
load(fpath)
df <- dataPrepFromSeurat(seuratBC, "label")

An example data preparation from SpatialExperiment object:

fpath <- system.file("extdata", "SpeBC.rda", package="stJoincount")
load(fpath)
df2 <- dataPrepFromSpE(SpeObjBC, "label")

Raster processing

This tool first converts a labeled spatial tissue map into a raster object, in which each spatial feature is represented by a pixel coded by label assignment. This process includes automatic calculation of optimal raster resolution and extent for the sample.

resolutionList <- resolutionCalc(humanBC)
resolutionList
#> [1] 152.89604  64.20792
mosaicIntegration <- rasterizeEachCluster(humanBC)
#> No optimal number found, using n = 110 instead.
#> In this case, there may be minor deviations in the subsequent calculation process.
#>         The results are for reference only

Visualization

After the labeled spatial sample being converted, the raster map can be visualized by:

mosaicIntPlot(humanBC, mosaicIntegration)

Join count analysis

A neighbors list is then created from the rasterized sample, in which adjacent and diagonal neighbors for each pixel are identified. After adding binary spatial weights to the neighbors list, a multi-categorical join count analysis is performed to tabulate “joins” between all possible combinations of label pairs. The function returns the observed join counts, the expected count under conditions of spatial randomness, and the variance calculated under non-free sampling.

joincount.result <- joincountAnalysis(mosaicIntegration)

The z-score is then calculated as the difference between observed and expected counts, divided by the square root of the variance. A heatmap of z-scores represents the result from the join count analysis for all possible label pairs.

matrix <- zscoreMatrix(humanBC, joincount.result)
zscorePlot(matrix)

sessionInfo()
#> R version 4.2.1 (2022-06-23)
#> Platform: x86_64-pc-linux-gnu (64-bit)
#> Running under: Ubuntu 20.04.5 LTS
#> 
#> Matrix products: default
#> BLAS:   /home/biocbuild/bbs-3.16-bioc/R/lib/libRblas.so
#> LAPACK: /home/biocbuild/bbs-3.16-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] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] ggplot2_3.3.6     pheatmap_1.0.12   stJoincount_1.0.0 BiocStyle_2.26.0 
#> 
#> loaded via a namespace (and not attached):
#>   [1] utf8_1.2.2                  reticulate_1.26            
#>   [3] R.utils_2.12.1              tidyselect_1.2.0           
#>   [5] htmlwidgets_1.5.4           grid_4.2.1                 
#>   [7] BiocParallel_1.32.0         Rtsne_0.16                 
#>   [9] DropletUtils_1.18.0         munsell_0.5.0              
#>  [11] units_0.8-0                 codetools_0.2-18           
#>  [13] ica_1.0-3                   future_1.28.0              
#>  [15] miniUI_0.1.1.1              withr_2.5.0                
#>  [17] spatstat.random_3.0-0       colorspace_2.0-3           
#>  [19] progressr_0.11.0            Biobase_2.58.0             
#>  [21] highr_0.9                   knitr_1.40                 
#>  [23] Seurat_4.2.0                stats4_4.2.1               
#>  [25] SingleCellExperiment_1.20.0 ROCR_1.0-11                
#>  [27] wk_0.7.0                    tensor_1.5                 
#>  [29] listenv_0.8.0               labeling_0.4.2             
#>  [31] MatrixGenerics_1.10.0       GenomeInfoDbData_1.2.9     
#>  [33] polyclip_1.10-4             farver_2.1.1               
#>  [35] rhdf5_2.42.0                parallelly_1.32.1          
#>  [37] vctrs_0.5.0                 generics_0.1.3             
#>  [39] xfun_0.34                   R6_2.5.1                   
#>  [41] GenomeInfoDb_1.34.0         locfit_1.5-9.6             
#>  [43] bitops_1.0-7                rhdf5filters_1.10.0        
#>  [45] spatstat.utils_3.0-1        cachem_1.0.6               
#>  [47] DelayedArray_0.24.0         assertthat_0.2.1           
#>  [49] promises_1.2.0.1            scales_1.2.1               
#>  [51] rgeos_0.5-9                 gtable_0.3.1               
#>  [53] beachmat_2.14.0             globals_0.16.1             
#>  [55] goftest_1.2-3               rlang_1.0.6                
#>  [57] splines_4.2.1               lazyeval_0.2.2             
#>  [59] spatstat.geom_3.0-3         s2_1.1.0                   
#>  [61] BiocManager_1.30.19         yaml_2.3.6                 
#>  [63] reshape2_1.4.4              abind_1.4-5                
#>  [65] httpuv_1.6.6                tools_4.2.1                
#>  [67] bookdown_0.29               spData_2.2.0               
#>  [69] SpatialExperiment_1.8.0     ellipsis_0.3.2             
#>  [71] spatstat.core_2.4-4         raster_3.6-3               
#>  [73] jquerylib_0.1.4             RColorBrewer_1.1-3         
#>  [75] proxy_0.4-27                BiocGenerics_0.44.0        
#>  [77] ggridges_0.5.4              Rcpp_1.0.9                 
#>  [79] plyr_1.8.7                  sparseMatrixStats_1.10.0   
#>  [81] zlibbioc_1.44.0             classInt_0.4-8             
#>  [83] purrr_0.3.5                 RCurl_1.98-1.9             
#>  [85] rpart_4.1.19                deldir_1.0-6               
#>  [87] pbapply_1.5-0               cowplot_1.1.1              
#>  [89] S4Vectors_0.36.0            zoo_1.8-11                 
#>  [91] SeuratObject_4.1.2          SummarizedExperiment_1.28.0
#>  [93] ggrepel_0.9.1               cluster_2.1.4              
#>  [95] magrittr_2.0.3              data.table_1.14.4          
#>  [97] magick_2.7.3                scattermore_0.8            
#>  [99] lmtest_0.9-40               RANN_2.6.1                 
#> [101] fitdistrplus_1.1-8          matrixStats_0.62.0         
#> [103] patchwork_1.1.2             mime_0.12                  
#> [105] evaluate_0.17               xtable_1.8-4               
#> [107] IRanges_2.32.0              gridExtra_2.3              
#> [109] compiler_4.2.1              tibble_3.1.8               
#> [111] KernSmooth_2.23-20          R.oo_1.25.0                
#> [113] htmltools_0.5.3             mgcv_1.8-41                
#> [115] later_1.3.0                 spdep_1.2-7                
#> [117] tidyr_1.2.1                 DBI_1.1.3                  
#> [119] MASS_7.3-58.1               boot_1.3-28                
#> [121] sf_1.0-8                    Matrix_1.5-1               
#> [123] cli_3.4.1                   R.methodsS3_1.8.2          
#> [125] parallel_4.2.1              igraph_1.3.5               
#> [127] GenomicRanges_1.50.0        pkgconfig_2.0.3            
#> [129] sp_1.5-0                    terra_1.6-17               
#> [131] plotly_4.10.0               scuttle_1.8.0              
#> [133] spatstat.sparse_3.0-0       bslib_0.4.0                
#> [135] dqrng_0.3.0                 XVector_0.38.0             
#> [137] stringr_1.4.1               digest_0.6.30              
#> [139] sctransform_0.3.5           RcppAnnoy_0.0.20           
#> [141] spatstat.data_3.0-0         rmarkdown_2.17             
#> [143] leiden_0.4.3                uwot_0.1.14                
#> [145] edgeR_3.40.0                DelayedMatrixStats_1.20.0  
#> [147] shiny_1.7.3                 rjson_0.2.21               
#> [149] lifecycle_1.0.3             nlme_3.1-160               
#> [151] jsonlite_1.8.3              Rhdf5lib_1.20.0            
#> [153] viridisLite_0.4.1           limma_3.54.0               
#> [155] fansi_1.0.3                 pillar_1.8.1               
#> [157] lattice_0.20-45             fastmap_1.1.0              
#> [159] httr_1.4.4                  survival_3.4-0             
#> [161] glue_1.6.2                  png_0.1-7                  
#> [163] class_7.3-20                stringi_1.7.8              
#> [165] sass_0.4.2                  HDF5Array_1.26.0           
#> [167] dplyr_1.0.10                irlba_2.3.5.1              
#> [169] e1071_1.7-12                future.apply_1.9.1