SpatialExperiment-class {SpatialExperiment} | R Documentation |
The SpatialExperiment
class is designed to represent spatially
resolved transcriptomics (ST) data. It inherits from the
SingleCellExperiment
class and is used in the same manner. In
addition, the class supports storage of spatial information via
spatialData
and spatialCoords
, and storage of
images via imgData
.
... |
Arguments passed to the |
sample_id |
A |
spatialDataNames |
A |
spatialCoordsNames |
A |
spatialData |
A |
spatialCoords |
A numeric matrix containing columns of spatial
coordinates, which will be accessible with |
scaleFactors |
Optional scale factors associated with the image(s). This
can be provided as a numeric value, numeric vector, list, or file path to a
JSON file for the 10x Genomics Visium platform. For 10x Genomics Visium,
the correct scale factor will automatically be selected depending on the
resolution of the image from |
imgData |
Optional |
imageSources |
Optional file path(s) or URL(s) for one or more image sources. |
image_id |
Optional character vector (same length as
|
loadImage |
Logical indicating whether to load image into memory.
Default = |
In this class, rows represent genes, and columns represent spots (for
spot-based ST platforms) or cells (for molecule-based ST platforms). As for
SingleCellExperiment
, counts
and logcounts
can be
stored in the assays
slot, and row and column metadata in
rowData
and colData
. For molecule-based ST data,
the additional measurements per molecule per cell can be stored in a
BumpyMatrix
-formatted assay
named molecules
.
The additional arguments in the constructor documented above (e.g.
spatialData
, spatialCoords
, imgData
, and others)
represent the main extensions to the SingleCellExperiment
class to store associated spatial and imaging information for ST data.
The constructor expects colData
to contain a column named
sample_id
. If this is not present, it will assign the value from the
sample_id
argument. If the imgData
argument is provided, the
constructor expects the imgData
DataFrame
to
already be built. Otherwise, it will build it from the imageSources
and (optional) image_id
arguments. If image_id
is not provided,
this will be assumed from sample_id
and imageSources
instead.
To combine multiple samples within a single object, see
combine
.
For 10x Genomics Visium datasets, the function read10xVisium
can be used to load data and create a SpatialExperiment
object
directly from Space Ranger output files.
A SpatialExperiment
object.
Dario Righelli and Helena L. Crowell
?"SpatialExperiment-methods"
, which includes:
spatialData
, spatialDataNames
,
spatialCoords
, spatialCoordsNames
,
imgData
, scaleFactors
?"SpatialExperiment-assays"
, which includes:
molecules
######################################################### # Example 1: Spot-based ST (10x Visium) using constructor ######################################################### dir <- system.file( file.path("extdata", "10xVisium", "section1"), package = "SpatialExperiment") # read in counts fnm <- file.path(dir, "raw_feature_bc_matrix") sce <- DropletUtils::read10xCounts(fnm) # read in image data img <- readImgData( path = file.path(dir, "spatial"), sample_id="foo") # read in spatial coordinates fnm <- file.path(dir, "spatial", "tissue_positions_list.csv") xyz <- read.csv(fnm, header = FALSE, col.names = c( "barcode", "in_tissue", "array_row", "array_col", "pxl_row_in_fullres", "pxl_col_in_fullres")) # construct observation & feature metadata rd <- S4Vectors::DataFrame( symbol = rowData(sce)$Symbol) # construct 'SpatialExperiment' (spe <- SpatialExperiment( assays = list(counts = assay(sce)), colData = colData(sce), rowData = rd, imgData = img, spatialData=DataFrame(xyz), spatialCoordsNames=c("pxl_col_in_fullres", "pxl_row_in_fullres"), sample_id="foo")) ############################################################# # Example 2: Spot-based ST (10x Visium) using 'read10xVisium' ############################################################# # see ?read10xVisium for details example(read10xVisium) ############################## # Example 3: Molecule-based ST ############################## # create simulated data n <- 1000; ng <- 50; nc <- 20 # sample xy-coordinates in [0,1] x <- runif(n) y <- runif(n) # assign each molecule to some gene-cell pair gs <- paste0("gene", seq(ng)) cs <- paste0("cell", seq(nc)) gene <- sample(gs, n, TRUE) cell <- sample(cs, n, TRUE) # construct data.frame of molecule coodinates df <- data.frame(gene, cell, x, y) # (assure gene & cell are factor so that # missing observations aren't dropped) df$gene <- factor(df$gene, gs) df$cell <- factor(df$cell, cs) # construct BumpyMatrix mol <- BumpyMatrix::splitAsBumpyMatrix( df[, c("x", "y")], row = df$gene, column = df$cell) # get count matrix y <- with(df, table(gene, cell)) y <- as.matrix(unclass(y)) # construct SpatialExperiment spe_mol <- SpatialExperiment( assays = list( counts = y, molecules = mol))