Overview

This tutorial demonstrates how to coerce GeoMxSet objects into Seurat or SpatialExperiment objects and the subsequent analyses. For more examples of what analyses are available in these objects, look at these Seurat or SpatialExperiment vignettes.

Data Processing

Data Processing should occur in GeomxTools. Due to the unique nature of the regions of interest (ROIs), it is recommended to use the preproccesing steps available in GeomxTools rather than the single-cell made preprocessing available in Seurat.

library(GeomxTools)
library(Seurat)
library(SpatialDecon)
library(patchwork)
datadir <- system.file("extdata", "DSP_NGS_Example_Data",
                       package="GeomxTools")
DCCFiles <- dir(datadir, pattern=".dcc$", full.names=TRUE)
PKCFiles <- unzip(zipfile = file.path(datadir,  "/pkcs.zip"))
SampleAnnotationFile <- file.path(datadir, "annotations.xlsx")

demoData <-
  suppressWarnings(readNanoStringGeoMxSet(dccFiles = DCCFiles,
                                          pkcFiles = PKCFiles,
                                          phenoDataFile = SampleAnnotationFile,
                                          phenoDataSheet = "CW005",
                                          phenoDataDccColName = "Sample_ID",
                                          protocolDataColNames = c("aoi",
                                                                   "cell_line",
                                                                   "roi_rep",
                                                                   "pool_rep",
                                                                   "slide_rep")))

After reading in the object, we will do a couple of QC steps.

  1. Shift all 0 counts by 1
  2. Flag low quality ROIs
  3. Flag low quality probes
  4. Remove low quality ROIs and probes
demoData <- shiftCountsOne(demoData, useDALogic=TRUE)
demoData <- setSegmentQCFlags(demoData, qcCutoffs = list(percentSaturation = 45))
demoData <- setBioProbeQCFlags(demoData)

# low sequenced ROIs
lowSaturation <- which(protocolData(demoData)[["QCFlags"]]$LowSaturation)

# probes that are considered outliers 
lowQCprobes <- which(featureData(demoData)[["QCFlags"]]$LowProbeRatio | 
                       featureData(demoData)[["QCFlags"]]$GlobalGrubbsOutlier)

# remove low quality ROIs and probes
passedQC <- demoData[-lowQCprobes, -lowSaturation]

dim(demoData)
## Features  Samples 
##     8707       88
dim(passedQC)
## Features  Samples 
##     8698       83

Objects must be aggregated to Target level data before coercing. This changes the row (gene) information to be the gene name rather than the probe ID.

featureType(passedQC)
## [1] "Probe"
data.frame(assayData(passedQC)[["exprs"]][seq_len(3), seq_len(3)])
DSP.1001250002642.A02.dcc DSP.1001250002642.A03.dcc DSP.1001250002642.A04.dcc
RTS0039454 294 239 6
RTS0039455 270 281 6
RTS0039456 255 238 3 
target_demoData <- aggregateCounts(passedQC)

featureType(target_demoData)
## [1] "Target"
data.frame(assayData(target_demoData)[["exprs"]][seq_len(3), seq_len(3)])
DSP.1001250002642.A02.dcc DSP.1001250002642.A03.dcc DSP.1001250002642.A04.dcc
ACTA2 328.286182 323.490808 6.081111
FOXA2 4.919019 4.919019 6.942503
NANOG 2.954177 4.128918 8.359554

It is recommended to normalize using a GeoMx specific model before coercing. The normalized data is now in the assayData slot called “q_norm”.

norm_target_demoData <- normalize(target_demoData, norm_method="quant",
                                  desiredQuantile = .75, toElt = "q_norm")

assayDataElementNames(norm_target_demoData)
## [1] "exprs"  "q_norm"
data.frame(assayData(norm_target_demoData)[["q_norm"]][seq_len(3), seq_len(3)])
DSP.1001250002642.A02.dcc DSP.1001250002642.A03.dcc DSP.1001250002642.A04.dcc
ACTA2 349.571598 344.257297 3.968122
FOXA2 5.237958 5.234795 4.530208
NANOG 3.145720 4.393974 5.454880

Seurat

Seurat Coercion

The three errors that can occur when trying to coerce to Seurat are:

  1. object must be on the target level
  2. object should be normalized, if you want raw data you can set forceRaw to TRUE
  3. normalized count matrix name must be valid
as.Seurat(demoData)
## Error in as.Seurat.NanoStringGeoMxSet(demoData): Data must be on Target level before converting to a Seurat Object
as.Seurat(target_demoData, normData = "exprs")
## Error in as.Seurat.NanoStringGeoMxSet(target_demoData, normData = "exprs"): It is NOT recommended to use Seurat's normalization for GeoMx data. 
##              Normalize using GeomxTools::normalize() or set forceRaw to TRUE if you want to continue with Raw data
as.Seurat(norm_target_demoData, normData = "exprs_norm")
## Error in as.Seurat.NanoStringGeoMxSet(norm_target_demoData, normData = "exprs_norm"): The normData name "exprs_norm" is not a valid assay name. Valid names are: exprs, q_norm

After coercing to a Seurat object all of the metadata is still accessible.

demoSeurat <- as.Seurat(norm_target_demoData, normData = "q_norm")

demoSeurat # overall data object
## An object of class Seurat 
## 1821 features across 83 samples within 1 assay 
## Active assay: GeoMx (1821 features, 0 variable features)
head(demoSeurat, 3) # most important ROI metadata
orig.ident nCount_GeoMx nFeature_GeoMx slide.name scan.name panel roi segment area NegGeoMean_Six.gene_test_v1_v1.1 NegGeoMean_VnV_GeoMx_Hs_CTA_v1.2 NegGeoSD_Six.gene_test_v1_v1.1 NegGeoSD_VnV_GeoMx_Hs_CTA_v1.2 q_norm_qFactors SampleID aoi cell_line roi_rep pool_rep slide_rep
DSP-1001250002642-A02.dcc 67643.36 1821 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 1 Geometric Segment 31318.73 1.487738 3.722752 1.560397 1.796952 0.9391100 DSP-1001250002642-A02 Geometric Segment-aoi-001 HS578T 1 1 1
DSP-1001250002642-A03.dcc 66360.01 1821 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 2 Geometric Segment 31318.73 2.518775 3.068217 1.820611 1.806070 0.9396774 DSP-1001250002642-A03 Geometric Segment-aoi-001 HS578T 2 1 1
DSP-1001250002642-A04.dcc 53749.39 1821 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 3 Geometric Segment 31318.73 2.847315 3.556275 1.654831 1.762066 1.5324910 DSP-1001250002642-A04 Geometric Segment-aoi-001 HEL 1 1 1 
demoSeurat@misc[1:8] # experiment data
## $PKCFileName
##            VnV_GeoMx_Hs_CTA_v1.2            Six-gene_test_v1_v1.1 
## "VnV Cancer Transcriptome Atlas"           "Six gene test custom" 
## 
## $PKCModule
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##    "VnV_GeoMx_Hs_CTA"    "Six-gene_test_v1" 
## 
## $PKCFileVersion
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                   1.2                   1.1 
## 
## $PKCFileDate
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##              "200518"              "200707" 
## 
## $AnalyteType
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                 "RNA"                 "RNA" 
## 
## $MinArea
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                 16000                 16000 
## 
## $MinNuclei
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                   200                   200 
## 
## $shiftedByOne
## [1] TRUE
head(demoSeurat@misc$sequencingMetrics) # sequencing metrics
FileVersion SoftwareVersion Date Plate_ID Well SeqSetId Raw Trimmed Stitched Aligned umiQ30 rtsQ30 DeduplicatedReads NTC_ID NTC Trimmed (%) Stitched (%) Aligned (%) Saturated (%)
DSP-1001250002642-A02.dcc 0.1 1.0.0 2020-07-14 1001250002642 A02 VH00121:3:AAAG2YWM5 646250 646250 616150 610390 0.9785 0.9804 312060 DSP-1001250002642-A01.dcc 7 100 95.34236 94.45106 48.87531
DSP-1001250002642-A03.dcc 0.1 1.0.0 2020-07-14 1001250002642 A03 VH00121:3:AAAG2YWM5 629241 629241 603243 597280 0.9784 0.9811 305528 DSP-1001250002642-A01.dcc 7 100 95.86836 94.92071 48.84677
DSP-1001250002642-A04.dcc 0.1 1.0.0 2020-07-14 1001250002642 A04 VH00121:3:AAAG2YWM5 831083 831083 798188 791804 0.9785 0.9801 394981 DSP-1001250002642-A01.dcc 7 100 96.04191 95.27376 50.11632
DSP-1001250002642-A05.dcc 0.1 1.0.0 2020-07-14 1001250002642 A05 VH00121:3:AAAG2YWM5 884485 884485 849060 842133 0.9796 0.9814 424162 DSP-1001250002642-A01.dcc 7 100 95.99484 95.21168 49.63242
DSP-1001250002642-A06.dcc 0.1 1.0.0 2020-07-14 1001250002642 A06 VH00121:3:AAAG2YWM5 781936 781936 751930 744669 0.9779 0.9803 355121 DSP-1001250002642-A01.dcc 7 100 96.16260 95.23401 52.31156
DSP-1001250002642-A07.dcc 0.1 1.0.0 2020-07-14 1001250002642 A07 VH00121:3:AAAG2YWM5 703034 703034 674815 668726 0.9776 0.9797 341008 DSP-1001250002642-A01.dcc 7 100 95.98611 95.12001 49.00632 
head(demoSeurat@misc$QCMetrics$QCFlags) # QC metrics
LowReads LowTrimmed LowStitched LowAligned LowSaturation LowNegatives HighNTC LowArea
DSP-1001250002642-A02.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A03.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A04.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A05.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A06.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A07.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE 
head(demoSeurat@assays$GeoMx@meta.features) # gene metadata
TargetName Module CodeClass GeneID SystematicName Negative
ACTA2 ACTA2 VnV_GeoMx_Hs_CTA_v1.2 Endogenous 59 ACTA2 FALSE
FOXA2 FOXA2 VnV_GeoMx_Hs_CTA_v1.2 Endogenous 3170 FOXA2 FALSE
NANOG NANOG VnV_GeoMx_Hs_CTA_v1.2 Endogenous 79923, 388112 NANOG, NANOGP8 FALSE
TRAC TRAC VnV_GeoMx_Hs_CTA_v1.2 Endogenous NA TRAC FALSE
TRBC1/2 TRBC1/2 VnV_GeoMx_Hs_CTA_v1.2 Endogenous NA TRBC1 FALSE
TRDC TRDC VnV_GeoMx_Hs_CTA_v1.2 Endogenous NA TRDC FALSE

All Seurat functionality is available after coercing. Outputs might differ if the ident value is set or not.

VlnPlot(demoSeurat, features = "nCount_GeoMx", pt.size = 0.1)

demoSeurat <- as.Seurat(norm_target_demoData, normData = "q_norm", ident = "cell_line")
VlnPlot(demoSeurat, features = "nCount_GeoMx", pt.size = 0.1)

Simple GeoMx data workflow

Here is an example of a typical dimensional reduction workflow.

demoSeurat <- FindVariableFeatures(demoSeurat)
demoSeurat <- ScaleData(demoSeurat)
demoSeurat <- RunPCA(demoSeurat, assay = "GeoMx", verbose = FALSE)
demoSeurat <- FindNeighbors(demoSeurat, reduction = "pca", dims = seq_len(30))
demoSeurat <- FindClusters(demoSeurat, verbose = FALSE)
demoSeurat <- RunUMAP(demoSeurat, reduction = "pca", dims = seq_len(30))

DimPlot(demoSeurat, reduction = "umap", label = TRUE, group.by = "cell_line")

In-depth GeoMx data workflow

Here is a work through of a more indepth DSP dataset. This is a non-small cell lung cancer (nsclc) tissue sample that has an ROI strategy to simulate a visium dataset (55 um circles evenly spaced apart). It was segmented on tumor and non-tumor.

data("nsclc", package = "SpatialDecon")
nsclc
## NanoStringGeoMxSet (storageMode: lockedEnvironment)
## assayData: 1700 features, 199 samples 
##   element names: exprs, exprs_norm 
## protocolData
##   sampleNames: ROI01Tumor ROI01TME ... ROI100TME (199 total)
##   varLabels: Mask.type Raw ... hkFactors (17 total)
##   varMetadata: labelDescription
## phenoData
##   sampleNames: ROI01Tumor ROI01TME ... ROI100TME (199 total)
##   varLabels: Sample_ID Tissue ... istumor (10 total)
##   varMetadata: labelDescription
## featureData
##   featureNames: ABCF1 ABL1 ... LAG3 (1700 total)
##   fvarLabels: TargetName HUGOSymbol ... Negative (9 total)
##   fvarMetadata: labelDescription
## experimentData: use 'experimentData(object)'
## Annotation: kiloplex with cell type spike-in [legacy panel] 
## signature: none
## feature: Target
## analyte: RNA
dim(nsclc)
## Features  Samples 
##     1700      199
data.frame(exprs(nsclc)[seq_len(5), seq_len(5)])
ROI01Tumor ROI01TME ROI02Tumor ROI02TME ROI03Tumor
ABCF1 55 26 47 30 102
ABL1 21 22 27 18 47
ACVR1B 89 30 57 29 122
ACVR1C 9 7 4 8 14
ACVR2A 14 15 9 12 22 
head(pData(nsclc))
Sample_ID Tissue Slide.name ROI AOI.name AOI.annotation x y nuclei istumor
ROI01Tumor ICP20th-L11-ICPKilo-ROI01-Tumor-A02 L11 ICPKilo ROI01 Tumor PanCK 0 8000 572 TRUE
ROI01TME ICP20th-L11-ICPKilo-ROI01-TME-A03 L11 ICPKilo ROI01 TME TME 0 8000 733 FALSE
ROI02Tumor ICP20th-L11-ICPKilo-ROI02-Tumor-A04 L11 ICPKilo ROI02 Tumor PanCK 600 8000 307 TRUE
ROI02TME ICP20th-L11-ICPKilo-ROI02-TME-A05 L11 ICPKilo ROI02 TME TME 600 8000 697 FALSE
ROI03Tumor ICP20th-L11-ICPKilo-ROI03-Tumor-A06 L11 ICPKilo ROI03 Tumor PanCK 1200 8000 583 TRUE
ROI03TME ICP20th-L11-ICPKilo-ROI03-TME-A07 L11 ICPKilo ROI03 TME TME 1200 8000 484 FALSE

When coercing, we can add the coordinate columns allowing for spatial graphing using Seurat.

nsclcSeurat <- as.Seurat(nsclc, normData = "exprs_norm", ident = "AOI.annotation", 
                         coordinates = c("x", "y"))

nsclcSeurat
## An object of class Seurat 
## 1700 features across 199 samples within 1 assay 
## Active assay: GeoMx (1700 features, 0 variable features)
VlnPlot(nsclcSeurat, features = "nCount_GeoMx", pt.size = 0.1)

nsclcSeurat <- FindVariableFeatures(nsclcSeurat)
nsclcSeurat <- ScaleData(nsclcSeurat)
nsclcSeurat <- RunPCA(nsclcSeurat, assay = "GeoMx", verbose = FALSE)
nsclcSeurat <- FindNeighbors(nsclcSeurat, reduction = "pca", dims = seq_len(30))
nsclcSeurat <- FindClusters(nsclcSeurat, verbose = FALSE)
nsclcSeurat <- RunUMAP(nsclcSeurat, reduction = "pca", dims = seq_len(30))

DimPlot(nsclcSeurat, reduction = "umap", label = TRUE, group.by = "AOI.name")

Spatial Graphing

Because this dataset is segmented, we need to separate the tumor and TME sections before using the spatial graphing. These Seurat functions were created for Visium data, so they can only plot the same sized circles.

Here we are showing the gene counts in each ROI separated by segment.

tumor <- suppressMessages(SpatialFeaturePlot(nsclcSeurat[,nsclcSeurat$AOI.name == "Tumor"], 
                                             features = "nCount_GeoMx", pt.size.factor = 12) + 
  labs(title = "Tumor") + 
  theme(legend.position = "none") + 
  scale_fill_continuous(type = "viridis",
                        limits = c(min(nsclcSeurat$nCount_GeoMx), 
                                   max(nsclcSeurat$nCount_GeoMx))))

TME <- suppressMessages(SpatialFeaturePlot(nsclcSeurat[,nsclcSeurat$AOI.name == "TME"], 
                                           features = "nCount_GeoMx", pt.size.factor = 12) + 
  labs(title = "TME") + 
  theme(legend.position = "right") +
  scale_fill_continuous(type = "viridis", 
                        limits = c(min(nsclcSeurat$nCount_GeoMx),
                                   max(nsclcSeurat$nCount_GeoMx))))

wrap_plots(tumor, TME)

Here we show the count for A2M

tumor <- suppressMessages(SpatialFeaturePlot(nsclcSeurat[,nsclcSeurat$AOI.name == "Tumor"], 
                                             features = "A2M", pt.size.factor = 12) + 
  labs(title = "Tumor") + 
  theme(legend.position = "none") + 
  scale_fill_continuous(type = "viridis",
                        limits = c(min(nsclcSeurat@assays$GeoMx@counts["A2M",]), 
                                   max(nsclcSeurat@assays$GeoMx@counts["A2M",]))))

TME <- suppressMessages(SpatialFeaturePlot(nsclcSeurat[,nsclcSeurat$AOI.name == "TME"], 
                                           features = "A2M", pt.size.factor = 12) + 
  labs(title = "TME") + 
  theme(legend.position = "right") +
  scale_fill_continuous(type = "viridis", 
                        limits = c(min(nsclcSeurat@assays$GeoMx@counts["A2M",]),
                                   max(nsclcSeurat@assays$GeoMx@counts["A2M",]))))

wrap_plots(tumor, TME)

Using the FindMarkers built in function from Seurat, we can determine the most differentially expressed genes in Tumor and TME

Idents(nsclcSeurat) <- nsclcSeurat$AOI.name
de_genes <- FindMarkers(nsclcSeurat, ident.1 = "Tumor", ident.2 = "TME")

de_genes <- de_genes[order(abs(de_genes$avg_log2FC), decreasing = TRUE),]
de_genes <- de_genes[is.finite(de_genes$avg_log2FC) & de_genes$p_val < 1e-25,]



for(i in rownames(de_genes)[1:2]){
  print(data.frame(de_genes[i,]))
  
  tumor <- suppressMessages(SpatialFeaturePlot(nsclcSeurat[,nsclcSeurat$AOI.name == "Tumor"], 
                                               features = i, pt.size.factor = 12) + 
  labs(title = "Tumor") + 
  theme(legend.position = "none") + 
  scale_fill_continuous(type = "viridis",
                        limits = c(min(nsclcSeurat@assays$GeoMx@counts[i,]), 
                                   max(nsclcSeurat@assays$GeoMx@counts[i,]))))

  TME <- suppressMessages(SpatialFeaturePlot(nsclcSeurat[,nsclcSeurat$AOI.name == "TME"], 
                                             features = i, pt.size.factor = 12) + 
    labs(title = "TME") + 
    theme(legend.position = "right") +
    scale_fill_continuous(type = "viridis", 
                          limits = c(min(nsclcSeurat@assays$GeoMx@counts[i,]),
                                     max(nsclcSeurat@assays$GeoMx@counts[i,]))))
  
  print(wrap_plots(tumor, TME))
}
##                p_val avg_log2FC pct.1 pct.2    p_val_adj
## CEACAM6 1.756187e-31    754.903     1     1 2.985517e-28

##            p_val avg_log2FC pct.1 pct.2    p_val_adj
## LYZ 9.276639e-32  -615.2034     1     1 1.577029e-28

SpatialExperiment

SpatialExperiment is an S4 class inheriting from SingleCellExperiment. It is meant as a data storage object rather than an analysis suite like Seurat. Because of this, this section won’t have the fancy analysis outputs like the Seurat section had but will show where in the object all the pieces are stored.

library(SpatialExperiment)

SpatialExperiment Coercion

The three errors that can occur when trying to coerce to SpatialExperiment are:

  1. object must be on the target level
  2. object should be normalized, if you want raw data you can set forceRaw to TRUE
  3. normalized count matrix name must be valid
as.SpatialExperiment(demoData)
## Error in as.SpatialExperiment.NanoStringGeoMxSet(demoData): Data must be on Target level before converting to a SpatialExperiment Object
as.SpatialExperiment(target_demoData, normData = "exprs")
## Error in as.SpatialExperiment.NanoStringGeoMxSet(target_demoData, normData = "exprs"): It is NOT recommended to use Seurat's normalization for GeoMx data. 
##              Normalize using GeomxTools::normalize() or set forceRaw to TRUE if you want to continue with Raw data
as.SpatialExperiment(norm_target_demoData, normData = "exprs_norm")
## Error in as.SpatialExperiment.NanoStringGeoMxSet(norm_target_demoData, : The normData name "exprs_norm" is not a valid assay name. Valid names are: exprs, q_norm

After coercing to a SpatialExperiment object all of the metadata is still accessible.

demoSPE <- as.SpatialExperiment(norm_target_demoData, normData = "q_norm")

demoSPE # overall data object
## class: SpatialExperiment 
## dim: 1821 83 
## metadata(11): PKCFileName PKCModule ... sequencingMetrics QCMetrics
## assays(1): GeoMx
## rownames(1821): ACTA2 FOXA2 ... C1orf43 SNRPD3
## rowData names(6): TargetName Module ... SystematicName Negative
## colnames(83): DSP-1001250002642-A02.dcc DSP-1001250002642-A03.dcc ...
##   DSP-1001250002642-H04.dcc DSP-1001250002642-H05.dcc
## colData names(18): slide name scan name ... slide_rep sample_id
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):
## spatialCoords names(0) :
## imgData names(0):
data.frame(head(colData(demoSPE))) # most important ROI metadata
slide.name scan.name panel roi segment area NegGeoMean_Six.gene_test_v1_v1.1 NegGeoMean_VnV_GeoMx_Hs_CTA_v1.2 NegGeoSD_Six.gene_test_v1_v1.1 NegGeoSD_VnV_GeoMx_Hs_CTA_v1.2 q_norm_qFactors SampleID aoi cell_line roi_rep pool_rep slide_rep sample_id
DSP-1001250002642-A02.dcc 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 1 Geometric Segment 31318.73 1.487738 3.722752 1.560397 1.796952 0.9391100 DSP-1001250002642-A02 Geometric Segment-aoi-001 HS578T 1 1 1 sample01
DSP-1001250002642-A03.dcc 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 2 Geometric Segment 31318.73 2.518775 3.068217 1.820611 1.806070 0.9396774 DSP-1001250002642-A03 Geometric Segment-aoi-001 HS578T 2 1 1 sample01
DSP-1001250002642-A04.dcc 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 3 Geometric Segment 31318.73 2.847315 3.556275 1.654831 1.762066 1.5324910 DSP-1001250002642-A04 Geometric Segment-aoi-001 HEL 1 1 1 sample01
DSP-1001250002642-A05.dcc 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 4 Geometric Segment 31318.73 2.632148 3.785600 2.042222 1.793823 1.6725916 DSP-1001250002642-A05 Geometric Segment-aoi-001 HEL 2 1 1 sample01
DSP-1001250002642-A06.dcc 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 5 Geometric Segment 31318.73 2.275970 4.064107 1.812577 1.839165 1.2351225 DSP-1001250002642-A06 Geometric Segment-aoi-001 U118MG 1 1 1 sample01
DSP-1001250002642-A07.dcc 6panel-old-slide1 (PTL-10891) cw005 (PTL-10891) Slide1 (v1.2) VnV Cancer Transcriptome Atlas, (v1.0) Six gene test custom 6 Geometric Segment 31318.73 2.059767 4.153701 1.952628 1.626391 1.2229991 DSP-1001250002642-A07 Geometric Segment-aoi-001 U118MG 2 1 1 sample01 
demoSPE@metadata[1:8] # experiment data
## $PKCFileName
##            VnV_GeoMx_Hs_CTA_v1.2            Six-gene_test_v1_v1.1 
## "VnV Cancer Transcriptome Atlas"           "Six gene test custom" 
## 
## $PKCModule
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##    "VnV_GeoMx_Hs_CTA"    "Six-gene_test_v1" 
## 
## $PKCFileVersion
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                   1.2                   1.1 
## 
## $PKCFileDate
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##              "200518"              "200707" 
## 
## $AnalyteType
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                 "RNA"                 "RNA" 
## 
## $MinArea
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                 16000                 16000 
## 
## $MinNuclei
## VnV_GeoMx_Hs_CTA_v1.2 Six-gene_test_v1_v1.1 
##                   200                   200 
## 
## $shiftedByOne
## [1] TRUE
head(demoSPE@metadata$sequencingMetrics) # sequencing metrics
FileVersion SoftwareVersion Date Plate_ID Well SeqSetId Raw Trimmed Stitched Aligned umiQ30 rtsQ30 DeduplicatedReads NTC_ID NTC Trimmed (%) Stitched (%) Aligned (%) Saturated (%)
DSP-1001250002642-A02.dcc 0.1 1.0.0 2020-07-14 1001250002642 A02 VH00121:3:AAAG2YWM5 646250 646250 616150 610390 0.9785 0.9804 312060 DSP-1001250002642-A01.dcc 7 100 95.34236 94.45106 48.87531
DSP-1001250002642-A03.dcc 0.1 1.0.0 2020-07-14 1001250002642 A03 VH00121:3:AAAG2YWM5 629241 629241 603243 597280 0.9784 0.9811 305528 DSP-1001250002642-A01.dcc 7 100 95.86836 94.92071 48.84677
DSP-1001250002642-A04.dcc 0.1 1.0.0 2020-07-14 1001250002642 A04 VH00121:3:AAAG2YWM5 831083 831083 798188 791804 0.9785 0.9801 394981 DSP-1001250002642-A01.dcc 7 100 96.04191 95.27376 50.11632
DSP-1001250002642-A05.dcc 0.1 1.0.0 2020-07-14 1001250002642 A05 VH00121:3:AAAG2YWM5 884485 884485 849060 842133 0.9796 0.9814 424162 DSP-1001250002642-A01.dcc 7 100 95.99484 95.21168 49.63242
DSP-1001250002642-A06.dcc 0.1 1.0.0 2020-07-14 1001250002642 A06 VH00121:3:AAAG2YWM5 781936 781936 751930 744669 0.9779 0.9803 355121 DSP-1001250002642-A01.dcc 7 100 96.16260 95.23401 52.31156
DSP-1001250002642-A07.dcc 0.1 1.0.0 2020-07-14 1001250002642 A07 VH00121:3:AAAG2YWM5 703034 703034 674815 668726 0.9776 0.9797 341008 DSP-1001250002642-A01.dcc 7 100 95.98611 95.12001 49.00632 
head(demoSPE@metadata$QCMetrics$QCFlags) # QC metrics
LowReads LowTrimmed LowStitched LowAligned LowSaturation LowNegatives HighNTC LowArea
DSP-1001250002642-A02.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A03.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A04.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A05.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A06.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE
DSP-1001250002642-A07.dcc FALSE FALSE FALSE FALSE FALSE TRUE FALSE FALSE 
data.frame(head(rowData(demoSPE))) # gene metadata
TargetName Module CodeClass GeneID SystematicName Negative
ACTA2 ACTA2 VnV_GeoMx_Hs_CTA_v1.2 Endogenous 59 ACTA2 FALSE
FOXA2 FOXA2 VnV_GeoMx_Hs_CTA_v1.2 Endogenous 3170 FOXA2 FALSE
NANOG NANOG VnV_GeoMx_Hs_CTA_v1.2 Endogenous 79923, 388112 NANOG, NANOGP8 FALSE
TRAC TRAC VnV_GeoMx_Hs_CTA_v1.2 Endogenous NA TRAC FALSE
TRBC1/2 TRBC1/2 VnV_GeoMx_Hs_CTA_v1.2 Endogenous NA TRBC1 FALSE
TRDC TRDC VnV_GeoMx_Hs_CTA_v1.2 Endogenous NA TRDC FALSE

When coercing, we can add the coordinate columns and they will be appended to the correct location in SpatialExperiment

nsclcSPE <- as.SpatialExperiment(nsclc, normData = "exprs_norm", coordinates = c("x", "y"))

nsclcSPE
## class: SpatialExperiment 
## dim: 1700 199 
## metadata(1): sequencingMetrics
## assays(1): GeoMx
## rownames(1700): ABCF1 ABL1 ... TNFSF4 LAG3
## rowData names(9): TargetName HUGOSymbol ... GlobalOutliers Negative
## colnames(199): ROI01Tumor ROI01TME ... ROI100Tumor ROI100TME
## colData names(20): Sample_ID Tissue ... hkFactors sample_id
## reducedDimNames(0):
## mainExpName: NULL
## altExpNames(0):
## spatialCoords names(2) : x y
## imgData names(0):
data.frame(head(spatialCoords(nsclcSPE)))
x y
ROI01Tumor 0 8000
ROI01TME 0 8000
ROI02Tumor 600 8000
ROI02TME 600 8000
ROI03Tumor 1200 8000
ROI03TME 1200 8000

With the coordinates and the metadata, we can create spatial graphing figures similar to Seurat’s

figureData <- as.data.frame(cbind(colData(nsclcSPE), spatialCoords(nsclcSPE)))

figureData <- cbind(figureData, A2M=as.numeric(nsclcSPE@assays@data$GeoMx["A2M",]))

tumor <- ggplot(figureData[figureData$AOI.name == "Tumor",], aes(x=x, y=y, color = A2M))+
  geom_point(size = 6)+
  scale_color_continuous(type = "viridis",
                        limits = c(min(figureData$A2M), 
                                   max(figureData$A2M)))+
  theme(legend.position = "none", panel.grid = element_blank(),
        panel.background = element_rect(fill = "white"),
        axis.title = element_blank(), axis.text = element_blank(), 
        axis.ticks = element_blank(), axis.line = element_blank())+
  labs(title = "Tumor")


TME <- ggplot(figureData[figureData$AOI.name == "TME",], aes(x=x, y=y, color = A2M))+
  geom_point(size = 6)+
  scale_color_continuous(type = "viridis",
                        limits = c(min(figureData$A2M), 
                                   max(figureData$A2M))) +
  theme(panel.grid = element_blank(), 
        panel.background = element_rect(fill = "white"), axis.title = element_blank(), 
        axis.text = element_blank(), axis.ticks = element_blank(), axis.line = element_blank())+
  labs(title = "TME")

wrap_plots(tumor, TME)

Image Overlays

The free-handed nature of Region of Interest (ROI) selection in a GeoMx experiment makes visualization on top of the image difficult in packages designed for different data. We created SpatialOmicsOverlay specifically to visualize and analyze these types of ROIs in a GeoMx experiment and the immunofluorescent-guided segmentation process.

sessionInfo()
## R version 4.2.0 RC (2022-04-19 r82224)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.4 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.15-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.15-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] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] SpatialExperiment_1.6.0     SingleCellExperiment_1.18.0
##  [3] SummarizedExperiment_1.26.1 GenomicRanges_1.48.0       
##  [5] GenomeInfoDb_1.32.1         IRanges_2.30.0             
##  [7] MatrixGenerics_1.8.0        matrixStats_0.62.0         
##  [9] patchwork_1.1.1             SpatialDecon_1.6.0         
## [11] SeuratObject_4.0.4          Seurat_4.1.0               
## [13] ggiraph_0.8.2               EnvStats_2.7.0             
## [15] GeomxTools_3.0.1            NanoStringNCTools_1.4.0    
## [17] ggplot2_3.3.5               S4Vectors_0.34.0           
## [19] Biobase_2.56.0              BiocGenerics_0.42.0        
## 
## loaded via a namespace (and not attached):
##   [1] utf8_1.2.2                reticulate_1.24          
##   [3] R.utils_2.11.0            tidyselect_1.1.2         
##   [5] lme4_1.1-29               htmlwidgets_1.5.4        
##   [7] grid_4.2.0                BiocParallel_1.30.0      
##   [9] Rtsne_0.16                DropletUtils_1.16.0      
##  [11] munsell_0.5.0             codetools_0.2-18         
##  [13] ica_1.0-2                 future_1.25.0            
##  [15] miniUI_0.1.1.1            withr_2.5.0              
##  [17] spatstat.random_2.2-0     colorspace_2.0-3         
##  [19] highr_0.9                 knitr_1.39               
##  [21] uuid_1.1-0                ROCR_1.0-11              
##  [23] tensor_1.5                listenv_0.8.0            
##  [25] labeling_0.4.2            GenomeInfoDbData_1.2.8   
##  [27] polyclip_1.10-0           farver_2.1.0             
##  [29] pheatmap_1.0.12           rhdf5_2.40.0             
##  [31] repmis_0.5                parallelly_1.31.1        
##  [33] vctrs_0.4.1               generics_0.1.2           
##  [35] xfun_0.30                 ggthemes_4.2.4           
##  [37] R6_2.5.1                  ggbeeswarm_0.6.0         
##  [39] locfit_1.5-9.5            rhdf5filters_1.8.0       
##  [41] bitops_1.0-7              spatstat.utils_2.3-0     
##  [43] reshape_0.8.9             DelayedArray_0.22.0      
##  [45] assertthat_0.2.1          promises_1.2.0.1         
##  [47] scales_1.2.0              beeswarm_0.4.0           
##  [49] gtable_0.3.0              beachmat_2.12.0          
##  [51] Cairo_1.5-15              globals_0.14.0           
##  [53] goftest_1.2-3             rlang_1.0.2              
##  [55] systemfonts_1.0.4         splines_4.2.0            
##  [57] lazyeval_0.2.2            spatstat.geom_2.4-0      
##  [59] yaml_2.3.5                reshape2_1.4.4           
##  [61] abind_1.4-5               httpuv_1.6.5             
##  [63] tools_4.2.0               ellipsis_0.3.2           
##  [65] spatstat.core_2.4-2       jquerylib_0.1.4          
##  [67] RColorBrewer_1.1-3        ggridges_0.5.3           
##  [69] Rcpp_1.0.8.3              plyr_1.8.7               
##  [71] sparseMatrixStats_1.8.0   zlibbioc_1.42.0          
##  [73] purrr_0.3.4               RCurl_1.98-1.6           
##  [75] rpart_4.1.16              deldir_1.0-6             
##  [77] pbapply_1.5-0             cowplot_1.1.1            
##  [79] zoo_1.8-10                ggrepel_0.9.1            
##  [81] cluster_2.1.3             magrittr_2.0.3           
##  [83] magick_2.7.3              data.table_1.14.2        
##  [85] RSpectra_0.16-1           scattermore_0.8          
##  [87] lmerTest_3.1-3            lmtest_0.9-40            
##  [89] RANN_2.6.1                fitdistrplus_1.1-8       
##  [91] R.cache_0.15.0            mime_0.12                
##  [93] evaluate_0.15             xtable_1.8-4             
##  [95] readxl_1.4.0              gridExtra_2.3            
##  [97] compiler_4.2.0            tibble_3.1.6             
##  [99] KernSmooth_2.23-20        crayon_1.5.1             
## [101] minqa_1.2.4               R.oo_1.24.0              
## [103] htmltools_0.5.2           mgcv_1.8-40              
## [105] later_1.3.0               tidyr_1.2.0              
## [107] DBI_1.1.2                 MASS_7.3-57              
## [109] boot_1.3-28               Matrix_1.4-1             
## [111] cli_3.3.0                 R.methodsS3_1.8.1        
## [113] parallel_4.2.0            igraph_1.3.1             
## [115] pkgconfig_2.0.3           numDeriv_2016.8-1.1      
## [117] scuttle_1.6.0             plotly_4.10.0            
## [119] spatstat.sparse_2.1-1     vipor_0.4.5              
## [121] bslib_0.3.1               dqrng_0.3.0              
## [123] XVector_0.36.0            stringr_1.4.0            
## [125] digest_0.6.29             sctransform_0.3.3        
## [127] RcppAnnoy_0.0.19          spatstat.data_2.2-0      
## [129] Biostrings_2.64.0         rmarkdown_2.14           
## [131] cellranger_1.1.0          leiden_0.3.10            
## [133] edgeR_3.38.0              uwot_0.1.11              
## [135] DelayedMatrixStats_1.18.0 shiny_1.7.1              
## [137] rjson_0.2.21              nloptr_2.0.0             
## [139] lifecycle_1.0.1           nlme_3.1-157             
## [141] jsonlite_1.8.0            Rhdf5lib_1.18.0          
## [143] viridisLite_0.4.0         limma_3.52.0             
## [145] fansi_1.0.3               pillar_1.7.0             
## [147] lattice_0.20-45           GGally_2.1.2             
## [149] ggrastr_1.0.1             fastmap_1.1.0            
## [151] httr_1.4.2                survival_3.3-1           
## [153] glue_1.6.2                png_0.1-7                
## [155] HDF5Array_1.24.0          stringi_1.7.6            
## [157] sass_0.4.1                dplyr_1.0.9              
## [159] irlba_2.3.5               future.apply_1.9.0