This function visualizes the gene expression stored in assays(spe)
or any
continuous variable stored in colData(spe)
for one given sample at the
spot-level using (by default) the histology information on the background.
To visualize clusters (or any discrete variable) use vis_clus()
.
vis_gene(
spe,
sampleid = unique(spe$sample_id)[1],
geneid = rowData(spe)$gene_search[1],
spatial = TRUE,
assayname = "logcounts",
minCount = 0,
viridis = TRUE,
image_id = "lowres",
alpha = NA,
cont_colors = if (viridis) viridisLite::viridis(21) else c("aquamarine4",
"springgreen", "goldenrod", "red"),
point_size = 2,
auto_crop = TRUE,
na_color = "#CCCCCC40",
...
)
Defaults to the output of
fetch_data(type = 'spe')
. This is a
SpatialExperiment-class
object with the spot-level Visium data and information required for
visualizing the histology. See fetch_data()
for more details.
A character(1)
specifying which sample to plot from
colData(spe)$sample_id
(formerly colData(spe)$sample_name
).
A character(1)
specifying the gene ID stored in
rowData(spe)$gene_search
or a continuous variable stored in colData(spe)
to visualize. If rowData(spe)$gene_search
is missing, then rownames(spe)
is used to search for the gene ID.
A logical(1)
indicating whether to include the histology
layer from geom_spatial()
. If you plan to use
ggplotly() then it's best to set this to FALSE
.
The name of the assays(spe)
to use for extracting the
gene expression data. Defaults to logcounts
.
A numeric(1)
specifying the minimum gene expression (or
value in the continuous variable) to visualize. Values at or below this
threshold will be set to NA
. Defaults to 0
.
A logical(1)
whether to use the color-blind friendly
palette from viridis or the color palette used
in the paper that was chosen for contrast when visualizing the data on
top of the histology image. One issue is being able to differentiate low
values from NA ones due to the purple-ish histology information that is
dependent on cell density.
A character(1)
with the name of the image ID you want to
use in the background.
A numeric(1)
in the [0, 1]
range that specifies the
transparency level of the data on the spots.
A character()
vector of colors that supersedes the
viridis
argument.
A numeric(1)
specifying the size of the points. Defaults
to 1.25
. Some colors look better if you use 2
for instance.
A logical(1)
indicating whether to automatically crop
the image / plotting area, which is useful if the Visium capture area is
not centered on the image and if the image is not a square.
A character(1)
specifying a color for the NA values.
If you set alpha = NA
then it's best to set na_color
to a color that has
alpha blending already, which will make non-NA values pop up more and the NA
values will show with a lighter color. This behavior is lost when alpha
is
set to a non-NA
value.
Passed to paste0() for making the title of the
plot following the sampleid
.
A ggplot2 object.
This function subsets spe
to the given sample and prepares the
data and title for vis_gene_p()
. It also adds a caption to the plot.
Other Spatial gene visualization functions:
vis_gene_p()
,
vis_grid_gene()
if (enough_ram()) {
## Obtain the necessary data
if (!exists("spe")) spe <- fetch_data("spe")
## Valid `geneid` values are those in
head(rowData(spe)$gene_search)
## or continuous variables stored in colData(spe)
## or rownames(spe)
## Visualize a default gene on the non-viridis scale
p1 <- vis_gene(
spe = spe,
sampleid = "151507",
viridis = FALSE
)
print(p1)
## Use a custom set of colors in the reverse order than usual
p2 <- vis_gene(
spe = spe,
sampleid = "151507",
cont_colors = rev(viridisLite::viridis(21, option = "magma"))
)
print(p2)
## Turn the alpha to 1, which makes the NA values have a full alpha
p2b <- vis_gene(
spe = spe,
sampleid = "151507",
cont_colors = rev(viridisLite::viridis(21, option = "magma")),
alpha = 1
)
print(p2b)
## Turn the alpha to NA, and use an alpha-blended "forestgreen" for
## the NA values
# https://gist.githubusercontent.com/mages/5339689/raw/2aaa482dfbbecbfcb726525a3d81661f9d802a8e/add.alpha.R
# add.alpha("forestgreen", 0.5)
p2c <- vis_gene(
spe = spe,
sampleid = "151507",
cont_colors = rev(viridisLite::viridis(21, option = "magma")),
alpha = NA,
na_color = "#228B2280"
)
print(p2c)
## Visualize a continuous variable, in this case, the ratio of chrM
## gene expression compared to the total expression at the spot-level
p3 <- vis_gene(
spe = spe,
sampleid = "151507",
geneid = "expr_chrM_ratio"
)
print(p3)
## Visualize a gene using the rownames(spe)
p4 <- vis_gene(
spe = spe,
sampleid = "151507",
geneid = rownames(spe)[which(rowData(spe)$gene_name == "MOBP")]
)
print(p4)
## Repeat without auto-cropping the image
p5 <- vis_gene(
spe = spe,
sampleid = "151507",
geneid = rownames(spe)[which(rowData(spe)$gene_name == "MOBP")],
auto_crop = FALSE
)
print(p5)
}
#> 2023-09-05 20:35:53.090938 loading file /github/home/.cache/R/BiocFileCache/152d33414640_Human_DLPFC_Visium_processedData_sce_scran_spatialLIBD.Rdata%3Fdl%3D1