Contents

1 Introduction

The SpatialExperiment package provides classes and methods for single cell spatial data handling.

This vignettes shows how to create a SpatialExperiment class by loading a seqFISH spatial dataset.

The seqFISH data are part of the BIRS Biointegration Workshop Hackathon previously published in Zhu et al. 2018.

2 Installation

if(!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("SpatialExperiment")

3 Loading data and libraries

We stored the provided seqFish and scRNA-seq data inside the SpatialExperiment external data seqFISH folder.

library(SpatialExperiment)

3.1 seqFISH data

fishCoordFile <- system.file(file.path("extdata", "seqFISH",
                            "fcortex.coordinates.txt"), 
                            package="SpatialExperiment")
fishCoordinates <- read.table(fishCoordFile, header=FALSE, sep=" ")
colnames(fishCoordinates) <- c("Cell_ID", "Irrelevant", "x", "y")

fishCellLabsFile <- system.file(file.path("extdata", "seqFISH", 
                            "seqfish_cell_labels.tsv"),
                            package="SpatialExperiment")
fishCellLabels <- read.table(file=fishCellLabsFile, header=FALSE, sep="\t")
colnames(fishCellLabels) <- c("Cell_ID", "cluster", "class", "classID", 
                                "Irrelevant", "Prob")
fishFeatCountsFile <- system.file(file.path("extdata", "seqFISH",
                            "seqfish_normalized_cortex_b2_testing.txt"), 
                            package="SpatialExperiment")
fishFeaturesCounts <- read.table(file=fishFeatCountsFile, 
                                header=FALSE, sep="\t", row.names=1)

4 SpatialExperiment package

4.1 Class handling

We can load the data inside the SpatialExperiment class. We then use the show method for looking at the class description.

se <- SpatialExperiment(rowData=rownames(fishFeaturesCounts),
                        colData=fishCellLabels,
                        assays=SimpleList(counts=as.matrix(fishFeaturesCounts)),
                        spatialCoords=fishCoordinates)
show(se)
## class: SpatialExperiment 
## dim: 113 1597 
## metadata(0):
## assays(1): counts
## rownames(113): abca15 abca9 ... zfp715 zfp90
## rowData names(1): X
## colnames(1597): V2 V3 ... V1597 V1598
## colData names(6): Cell_ID cluster ... Irrelevant Prob
## reducedDimNames(0):
## altExpNames(0):
## spatialCoordinates(4): Cell_ID Irrelevant x y

4.2 Spatial Coordinates

With the aid of the spatialCoords methods we can get and set the spatial coordinates.

4.2.1 Getter

spatialCoords(se)
## DataFrame with 1597 rows and 4 columns
##        Cell_ID Irrelevant         x         y
##      <integer>  <integer> <numeric> <numeric>
## 1            1        100    265.76   -231.14
## 2            2        100    290.48   -261.52
## 3            3        100    257.12   -133.35
## 4            4        100    753.46   -261.14
## 5            5        100    700.01   -169.05
## ...        ...        ...       ...       ...
## 1593      1593        100   1129.06  -1669.57
## 1594      1594        100   1044.02  -1872.66
## 1595      1595        100   1388.76  -1880.47
## 1596      1596        100   5172.85  -1340.96
## 1597      1597        100   5220.60  -1523.37

4.2.2 Setter

We create a fake fish coordinates data frame and overwrite the old one, showing that the y coordinates are identical to the x ones.

fakeFishCoords <- cbind(fishCoordinates[,c(1:3)], fishCoordinates[,3])
colnames(fakeFishCoords) <- c("Cell_ID", "Irrelevant", "x", "y")
spatialCoords(se) <- fakeFishCoords
spatialCoords(se)
## DataFrame with 1597 rows and 4 columns
##        Cell_ID Irrelevant         x         y
##      <integer>  <integer> <numeric> <numeric>
## 1            1        100    265.76    265.76
## 2            2        100    290.48    290.48
## 3            3        100    257.12    257.12
## 4            4        100    753.46    753.46
## 5            5        100    700.01    700.01
## ...        ...        ...       ...       ...
## 1593      1593        100   1129.06   1129.06
## 1594      1594        100   1044.02   1044.02
## 1595      1595        100   1388.76   1388.76
## 1596      1596        100   5172.85   5172.85
## 1597      1597        100   5220.60   5220.60
sessionInfo()
## R version 4.0.3 (2020-10-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.5 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.12-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.12-bioc/R/lib/libRlapack.so
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        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] parallel  stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] rjson_0.2.20                Matrix_1.2-18              
##  [3] SpatialExperiment_1.0.0     SingleCellExperiment_1.12.0
##  [5] SummarizedExperiment_1.20.0 Biobase_2.50.0             
##  [7] GenomicRanges_1.42.0        GenomeInfoDb_1.26.0        
##  [9] IRanges_2.24.0              S4Vectors_0.28.0           
## [11] BiocGenerics_0.36.0         MatrixGenerics_1.2.0       
## [13] matrixStats_0.57.0          BiocStyle_2.18.0           
## 
## loaded via a namespace (and not attached):
##  [1] knitr_1.30             XVector_0.30.0         magrittr_1.5          
##  [4] zlibbioc_1.36.0        lattice_0.20-41        rlang_0.4.8           
##  [7] stringr_1.4.0          tools_4.0.3            grid_4.0.3            
## [10] xfun_0.18              htmltools_0.5.0        yaml_2.2.1            
## [13] digest_0.6.27          bookdown_0.21          GenomeInfoDbData_1.2.4
## [16] BiocManager_1.30.10    bitops_1.0-6           RCurl_1.98-1.2        
## [19] evaluate_0.14          rmarkdown_2.5          DelayedArray_0.16.0   
## [22] stringi_1.5.3          compiler_4.0.3