Title: | Visualization functions for spatial transcriptomics data |
---|---|
Description: | Visualization functions for spatial transcriptomics data. Includes functions to generate several types of plots, including spot plots, feature (molecule) plots, reduced dimension plots, spot-level quality control (QC) plots, and feature-level QC plots, for datasets from the 10x Genomics Visium and other technological platforms. Datasets are assumed to be in either SpatialExperiment or SingleCellExperiment format. |
Authors: | Lukas M. Weber [aut, cre] , Helena L. Crowell [aut] , Yixing E. Dong [aut] |
Maintainer: | Lukas M. Weber <[email protected]> |
License: | MIT + file LICENSE |
Version: | 1.13.0 |
Built: | 2024-11-19 03:57:57 UTC |
Source: | https://github.com/bioc/ggspavis |
Plotting functions for spatial transcriptomics data.
plotDimRed( spe, plot_type = c("UMAP", "PCA"), annotate = NULL, feature_names = NULL, assay_name = "counts", update_dimnames = TRUE, pal = NULL, point_size = 0.3, legend_point_size = 3, text_by = NULL, text_by_size = 5, text_by_color = "black" )
plotDimRed( spe, plot_type = c("UMAP", "PCA"), annotate = NULL, feature_names = NULL, assay_name = "counts", update_dimnames = TRUE, pal = NULL, point_size = 0.3, legend_point_size = 3, text_by = NULL, text_by_size = 5, text_by_color = "black" )
spe |
Input data, assumed to be a |
plot_type |
Type of reduced dimension plot. Possible options are "UMAP", "PCA", or any other set of reduced dimensions stored in the input object. Default = "UMAP". |
annotate |
Variable to show as annotations. This may be discrete or
continuous. For a discrete variable (e.g. cluster labels), this should be
the name of a column in |
feature_names |
Name of column in |
assay_name |
Name of |
update_dimnames |
Whether to update column names of |
pal |
Color palette for annotations. Options for discrete values are "libd_layer_colors", "Okabe-Ito", or any vector of color names or hex values. For continuous values, provide a vector of length 2 for the low and high range, e.g. c("gray90", "navy"). |
point_size |
Point size. Default = 0.3. |
legend_point_size |
Legend point size for discrete annotations. Default = 3. |
text_by |
Column name of annotation labels to display over each cluster
of points. This will usually be the same as |
text_by_size |
Text size for annotation labels over each cluster. Default = 5. |
text_by_color |
Color name or hex code for annotation labels. Default = "black". |
Function to create reduced dimension plot (e.g. PCA or UMAP) with additional optional annotations such as cluster labels, expression of a gene, or quality control metrics.
Returns a ggplot object, which may be further modified using ggplot functions.
Lukas M. Weber and Yixing E. Dong
library(STexampleData) spe <- Visium_humanDLPFC() # select spots over tissue spe <- spe[, colData(spe)$in_tissue == 1] # use small subset of data for this example n <- 200 set.seed(123) spe <- spe[, sample(seq_len(ncol(spe)), n)] # calculate logcounts library(scran) spe <- logNormCounts(spe) # identify top highly variable genes (HVGs) is_mito <- grepl("(^MT-)|(^mt-)", rowData(spe)$gene_name) spe <- spe[!is_mito, ] dec <- modelGeneVar(spe) top_hvgs <- getTopHVGs(dec, prop = 0.1) # run dimensionality reduction library(scater) set.seed(123) spe <- runPCA(spe, subset_row = top_hvgs) set.seed(123) spe <- runUMAP(spe, dimred = "PCA") colnames(reducedDim(spe, "UMAP")) <- paste0("UMAP", 1:2) # generate plot plotDimRed(spe, plot_type = "UMAP", annotate = "ground_truth")
library(STexampleData) spe <- Visium_humanDLPFC() # select spots over tissue spe <- spe[, colData(spe)$in_tissue == 1] # use small subset of data for this example n <- 200 set.seed(123) spe <- spe[, sample(seq_len(ncol(spe)), n)] # calculate logcounts library(scran) spe <- logNormCounts(spe) # identify top highly variable genes (HVGs) is_mito <- grepl("(^MT-)|(^mt-)", rowData(spe)$gene_name) spe <- spe[!is_mito, ] dec <- modelGeneVar(spe) top_hvgs <- getTopHVGs(dec, prop = 0.1) # run dimensionality reduction library(scater) set.seed(123) spe <- runPCA(spe, subset_row = top_hvgs) set.seed(123) spe <- runUMAP(spe, dimred = "PCA") colnames(reducedDim(spe, "UMAP")) <- paste0("UMAP", 1:2) # generate plot plotDimRed(spe, plot_type = "UMAP", annotate = "ground_truth")
Plotting functions for spatial transcriptomics data.
plotFeatureQC( spe, plot_type = c("histogram", "violin"), x_metric = NULL, annotate = NULL, n_bins = 100, point_size = 0.1, scale_log1p = TRUE )
plotFeatureQC( spe, plot_type = c("histogram", "violin"), x_metric = NULL, annotate = NULL, n_bins = 100, point_size = 0.1, scale_log1p = TRUE )
spe |
Input data, assumed to be a |
plot_type |
Type of QC plot. Options are "histogram" and "violin". See Details for additional details. |
x_metric |
Name of column in |
annotate |
Name of column in |
n_bins |
Number of bins for histograms. Default = 100. Optional argument used for histograms. |
point_size |
Point size. Default = 0.1. Optional argument for violin plots. |
scale_log1p |
Whether to log1p-scale axes. Default = TRUE. |
Function to create quality control (QC) plots for spatial transcriptomics data.
The following types of QC plots are available for feature-level QC (see
plotSpotQC
for spot-level or cell-level QC):
Histogram (plot_type = "histogram"
) for a single QC metric, e.g.
total UMI counts across all spots per feature. The histogram can optionally
highlight selected features, e.g. low abundance features.
Violin (plot_type = "violin"
) for a single QC metric, e.g. total
UMI counts across all spots per feature. The violin plot can optionally
highlight selected features, e.g. low abundance features.
Returns a ggplot object, which may be further modified using ggplot functions.
Yixing E. Dong and Lukas M. Weber
library(STexampleData) spe <- Visium_humanDLPFC() rowData(spe)$feature_sum <- rowSums(counts(spe)) rowData(spe)$low_abundance <- rowSums(counts(spe) > 0) < 20 plotFeatureQC(spe, plot_type = "histogram", x_metric = "feature_sum", annotate = "low_abundance") plotFeatureQC(spe, plot_type = "violin", x_metric = "feature_sum", annotate = "low_abundance")
library(STexampleData) spe <- Visium_humanDLPFC() rowData(spe)$feature_sum <- rowSums(counts(spe)) rowData(spe)$low_abundance <- rowSums(counts(spe) > 0) < 20 plotFeatureQC(spe, plot_type = "histogram", x_metric = "feature_sum", annotate = "low_abundance") plotFeatureQC(spe, plot_type = "violin", x_metric = "feature_sum", annotate = "low_abundance")
Plotting functions for spatial transcriptomics data.
plotMolecules( spe, molecule = NULL, x_coord = NULL, y_coord = NULL, sample_id = "sample_id", pal = c("gray90", "navy"), point_size = 0.3 )
plotMolecules( spe, molecule = NULL, x_coord = NULL, y_coord = NULL, sample_id = "sample_id", pal = c("gray90", "navy"), point_size = 0.3 )
spe |
(SpatialExperiment) Input data, assumed to be a
|
molecule |
Name of mRNA molecule to plot (assumed to match one of the
row names of |
x_coord |
Name of column in |
y_coord |
Name of column in |
sample_id |
Name of column in |
pal |
Color palette, provided as a vector of length 2 for the low and high range. Default = c("gray90", "navy"). |
point_size |
Point size. Default = 0.3. |
Function to create spot plot for molecule-based datasets, showing spatial locations in x-y coordinates with optional annotations such as expression of a gene.
Returns a ggplot object, which may be further modified using ggplot functions.
Lukas M. Weber
library(STexampleData) spe <- seqFISH_mouseEmbryo() plotMolecules(spe, molecule = "Sox2")
library(STexampleData) spe <- seqFISH_mouseEmbryo() plotMolecules(spe, molecule = "Sox2")
Plotting functions for spatial transcriptomics data.
plotSpotQC( spe, plot_type = c("histogram", "scatter", "spot", "violin"), x_coord = NULL, y_coord = NULL, x_metric = NULL, y_metric = NULL, x_threshold = NULL, y_threshold = NULL, trend = TRUE, marginal = TRUE, annotate = NULL, in_tissue = NULL, legend_point_size = 3, n_bins = 100, point_size = 0.3, y_reverse = TRUE ) plotQC(...)
plotSpotQC( spe, plot_type = c("histogram", "scatter", "spot", "violin"), x_coord = NULL, y_coord = NULL, x_metric = NULL, y_metric = NULL, x_threshold = NULL, y_threshold = NULL, trend = TRUE, marginal = TRUE, annotate = NULL, in_tissue = NULL, legend_point_size = 3, n_bins = 100, point_size = 0.3, y_reverse = TRUE ) plotQC(...)
spe |
Input data, assumed to be a |
plot_type |
Type of QC plot. Options are "histogram", "scatter", "spot", and "violin". See Details for additional details. |
x_coord |
Name of column in |
y_coord |
Name of column in |
x_metric |
Name of column in |
y_metric |
Name of column in |
x_threshold |
QC filtering threshold on x-axis metric to highlight with vertical line. Default = NULL. Optional argument used for scatter plots. |
y_threshold |
QC filtering threshold on y-axis metric to highlight with horizontal line. Default = NULL. Optional argument used for scatter plots. |
trend |
Whether to show smoothed trend line (loess). Default = TRUE. Optional argument used for scatter plots. |
marginal |
Whether to show marginal histograms. Default = TRUE. Optional argument used for scatter plots. |
annotate |
Name of column in |
in_tissue |
Name of column in |
legend_point_size |
Legend point size. Default = 3. Optional argument used for spot plots. |
n_bins |
Number of bins for histograms. Default = 100. Optional argument used for histograms. |
point_size |
Point size. Default = 0.3. Optional argument for scatter plots, spot plots, and violin plots. Suggested values: 0.5 for scatter plots, 0.3 for spot plots, 0.1 for violin plots. |
y_reverse |
Whether to reverse y coordinates. This is usually required for 10x Genomics Visium datasets when using the default coordinate values. Default = TRUE. Set to FALSE if not needed, e.g. for other platforms. Optional argument used for spot plots. |
... |
Not used. |
Function to create quality control (QC) plots for spatial transcriptomics data.
The following types of QC plots are available for spot-level or cell-level QC
(see plotFeatureQC
for feature-level QC):
Histogram (plot_type = "histogram"
) for a single QC metric, e.g.
number of UMI counts per spot. For number of counts per spot, the histogram
can optionally highlight selected spots, e.g. spots with low library size.
Scatter plot (plot_type = "scatter"
) comparing two QC metrics,
e.g. number of detected features vs. number of cells per spot, with optional
horizontal and vertical lines highlighting QC filtering thresholds.
Spot plot (plot_type = "spot"
) showing spots in spatial x-y
coordinates, e.g. highlighting selected spots that do not meet filtering
thresholds.
Violin plot (plot_type = "violin"
) for a single QC metric, e.g.
number of UMI counts per spot. For number of counts per spot, the violin plot
can optionally highlight selected spots, e.g. spots with low library size.
Returns a ggplot object, which may be further modified using ggplot functions.
Lukas M. Weber and Yixing E. Dong
library(STexampleData) spe <- Visium_humanDLPFC() colData(spe)$sum <- colSums(counts(spe)) colData(spe)$low_libsize <- colData(spe)$sum < 400 plotSpotQC(spe, plot_type = "histogram", x_metric = "sum", annotate = "low_libsize") plotSpotQC(spe, plot_type = "scatter", x_metric = "sum", y_metric = "cell_count") plotSpotQC(spe, plot_type = "spot", annotate = "low_libsize", in_tissue = "in_tissue") plotSpotQC(spe, plot_type = "violin", x_metric = "sum", annotate = "low_libsize")
library(STexampleData) spe <- Visium_humanDLPFC() colData(spe)$sum <- colSums(counts(spe)) colData(spe)$low_libsize <- colData(spe)$sum < 400 plotSpotQC(spe, plot_type = "histogram", x_metric = "sum", annotate = "low_libsize") plotSpotQC(spe, plot_type = "scatter", x_metric = "sum", y_metric = "cell_count") plotSpotQC(spe, plot_type = "spot", annotate = "low_libsize", in_tissue = "in_tissue") plotSpotQC(spe, plot_type = "violin", x_metric = "sum", annotate = "low_libsize")
Plotting functions for spatial transcriptomics data.
plotSpots( spe, x_coord = NULL, y_coord = NULL, sample_id = NULL, in_tissue = "in_tissue", annotate = NULL, feature_names = NULL, assay_name = "counts", pal = NULL, point_size = 0.3, legend_position = "right", legend_point_size = 3, show_axes = FALSE, y_reverse = TRUE, text_by = NULL, text_by_size = 5, text_by_color = "black" )
plotSpots( spe, x_coord = NULL, y_coord = NULL, sample_id = NULL, in_tissue = "in_tissue", annotate = NULL, feature_names = NULL, assay_name = "counts", pal = NULL, point_size = 0.3, legend_position = "right", legend_point_size = 3, show_axes = FALSE, y_reverse = TRUE, text_by = NULL, text_by_size = 5, text_by_color = "black" )
spe |
Input data, assumed to be a |
x_coord |
Name of column in |
y_coord |
Name of column in |
sample_id |
Name of column in |
in_tissue |
Name of column in |
annotate |
Variable to show as annotations. This may be discrete or
continuous. For a discrete variable (e.g. cluster labels), this should be
the name of a column in |
feature_names |
Name of column in |
assay_name |
Name of |
pal |
Color palette for annotations. Options for discrete values are "libd_layer_colors", "Okabe-Ito", or any vector of color names or hex values. For continuous values, provide a vector of length 2 for the low and high range, e.g. c("gray90", "navy"). |
point_size |
Point size. Default = 0.3. |
legend_position |
Legend position for discrete annotations. Options are "left", "right", "top", "bottom", and "none". Default = "right". |
legend_point_size |
Legend point size for discrete annotations. Default = 3. |
show_axes |
Whether to show axis titles, text, and ticks. Default = FALSE. |
y_reverse |
Whether to reverse y coordinates. This is usually required for 10x Genomics Visium datasets when using the default coordinate values. Default = TRUE. Set to FALSE if not needed, e.g. for other platforms. |
text_by |
Column name of annotation labels to display over each cluster
of points. This will usually be the same as |
text_by_size |
Text size for annotation labels over each cluster. Default = 5. |
text_by_color |
Color name or hex code for annotation labels. Default = "black". |
Function to create spot plot showing spatial locations in x-y coordinates with optional annotations such as cluster labels, expression of a gene, or quality control metrics.
Returns a ggplot object, which may be further modified using ggplot functions.
Lukas M. Weber and Yixing E. Dong
library(STexampleData) # discrete annotations spe <- Visium_humanDLPFC() plotSpots(spe, annotate = "ground_truth") # continuous annotations spe <- Visium_mouseCoronal() plotSpots(spe, annotate = "Gapdh", feature_names = "gene_name")
library(STexampleData) # discrete annotations spe <- Visium_humanDLPFC() plotSpots(spe, annotate = "ground_truth") # continuous annotations spe <- Visium_mouseCoronal() plotSpots(spe, annotate = "Gapdh", feature_names = "gene_name")
Plots for spatially resolved transcriptomics data from the 10x Genomics Visium platform
plotVisium( spe, spots = TRUE, annotate = NULL, highlight = NULL, facets = "sample_id", image = TRUE, zoom = FALSE, show_axes = FALSE, assay = "counts", trans = "identity", point_size = 1, legend_position = "right", x_coord = NULL, y_coord = NULL, y_reverse = TRUE, sample_ids = NULL, image_ids = NULL, pal = NULL )
plotVisium( spe, spots = TRUE, annotate = NULL, highlight = NULL, facets = "sample_id", image = TRUE, zoom = FALSE, show_axes = FALSE, assay = "counts", trans = "identity", point_size = 1, legend_position = "right", x_coord = NULL, y_coord = NULL, y_reverse = TRUE, sample_ids = NULL, image_ids = NULL, pal = NULL )
spe |
(SpatialExperiment) Input data object. |
spots |
(logical) Whether to display spots (spatial barcodes) as points. Default = TRUE. |
annotate |
(character) Column in |
highlight |
(character) Column in |
facets |
(character) Column in |
image |
(logical) Whether to show histology image as background. Default = TRUE. |
zoom |
(logical) Whether to zoom to area of tissue containing spots. Default = FALSE |
show_axes |
(logical) Whether to show axes and coordinates. Default = FALSE |
assay |
(character) Name of assay data to use when |
trans |
Transformation to apply for continuous scales. Ignored unless
|
point_size |
(numeric) Point size. Default = 1. |
legend_position |
Legend position for annotations. Options are "left", "right", "top", "bottom", and "none". Default = "right". |
x_coord |
(character) Column in |
y_coord |
(character) Column in |
y_reverse |
(logical) Whether to reverse y coordinates, which is often required for Visium data, depending on the orientation of the raw data. Default = TRUE. |
sample_ids |
(character) Samples to show, if multiple samples are available. Default = NULL (show all samples). |
image_ids |
(character) Images to show, if multiple images are available. Default = NULL (show all images). |
pal |
(character) Color palette for points. Options for discrete labels are "libd_layer_colors", "Okabe-Ito", or a custom vector of hex color codes. Options for continuous values are "viridis", a single color name (e.g. "red", "navy", etc), or a vector of length two containing color names for each end of the scale. Default = "libd_layer_colors" for discrete data, and "viridis" for continuous data. |
Function to generate plots for spatially resolved transcriptomics datasets from the 10x Genomics Visium spatially platform.
This function generates a plot for spot-based spatially resolved transcriptomics data from the 10x Genomics Visium platform, with several options available to adjust the plot type and style.
Returns a ggplot object. Additional plot elements can be added as ggplot elements (e.g. title, customized formatting, etc).
Helena L. Crowell, with modifications by Lukas M. Weber and Yixing E. Dong
library(STexampleData) spe <- Visium_mouseCoronal() # color by x coordinate, highlight in-tissue spots plotVisium(spe, annotate = "pxl_col_in_fullres", highlight = "in_tissue") # subset in-tissue spots sub <- spe[, as.logical(colData(spe)$in_tissue)] # color by feature counts, don't include image rownames(sub) <- make.names(rowData(sub)$gene_name) plotVisium(sub, annotate = "Gad2", assay = "counts")
library(STexampleData) spe <- Visium_mouseCoronal() # color by x coordinate, highlight in-tissue spots plotVisium(spe, annotate = "pxl_col_in_fullres", highlight = "in_tissue") # subset in-tissue spots sub <- spe[, as.logical(colData(spe)$in_tissue)] # color by feature counts, don't include image rownames(sub) <- make.names(rowData(sub)$gene_name) plotVisium(sub, annotate = "Gad2", assay = "counts")