Package 'spatialLIBD'

Description

This function re-uses the SpatialExperiment::scaleFactors() from current images when adding new images. This is useful if you take for example a multi-channel VisiumIF image and break into several single-channel images that all have the same dimensions. So you could have a set of images such as channel_01_lowres and channel_02_lowres that have the same dimensions and viewing area as the lowres image produced by SpaceRanger, each with only one channel. Similarly, you might have done some image manipulation for a given image and generated one or more images with the same dimensions as existing images.

Usage

add_images(
  spe,
  image_dir,
  image_pattern,
  image_id_current = "lowres",
  image_id = image_pattern,
  image_paths = locate_images(spe, image_dir, image_pattern)
)

Arguments

Value

A SpatialExperiment-class object with the additional image data in imgData(spe).

See Also

Other Functions for adding non-standard images: locate_images()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Add an image
    SpatialExperiment::imgData(add_images(
        spe,
        image_id_current = "lowres",
        image_id = "lowres_aws",
        image_paths = c("151507" = "https://spatial-dlpfc.s3.us-east-2.amazonaws.com/images/151507_tissue_lowres_image.png")
    ))
}

Description

This function adds spe$key to a SpatialExperiment-class object which is unique across all spots.

Usage

add_key(spe, overwrite = TRUE)

Arguments

Value

A SpatialExperiment-class object with key added to the colData(spe) that is unique across all spots.

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## This object already has a 'key'
    head(spe$key)

    ## We can clean it
    spe$key_original <- spe$key
    spe$key <- NULL

    ## and then add it back
    spe <- add_key(spe)
    head(spe$key)

    ## Note that the original 'key' order was 'sample_id'_'barcode' and we'
    ## have since changed it to 'barcode'_'sample_id'.

    ## Below we restore the original 'key'
    spe$key <- spe$key_original
    spe$key_original <- NULL
    head(spe$key)
}

Description

This function identify spots in a SpatialExperiment-class (SPE) with outlier quality control values: low sum_umi or sum_gene, or high expr_chrM_ratio, utilizing scuttle::isOutlier. Also identifies in-tissue edge spots and distance to the edge for each spot.

Usage

add_qc_metrics(spe, overwrite = FALSE)

Arguments

Details

The initial version of this function lives at https://github.com/LieberInstitute/Visium_SPG_AD/blob/master/code/07_spot_qc/01_qc_metrics_and_segmentation.R.

Value

A SpatialExperiment with added quality control information added to the colData().

scran_low_lib_size

shows spots that have a low library size.

scran_low_n_features

spots with a low number of expressed genes.

scran_high_Mito_percent

spots with a high percent of mitochondrial gene expression.

scran_discard

spots belonging to either scran_low_lib_size, scran_low_n_feature, or scran_high_Mito_percent.

edge_spot

spots that are automatically detected as the edge spots of the in_tissue section.

edge_distance

closest distance in number of spots to either the vertical or horizontal edge.

scran_low_lib_size_edge

spots that have a low library size and are an edge spot.

Author(s)

Louise A. Huuki-Myers

Examples

## Obtain the necessary data
spe_pre_qc <- fetch_data("spatialDLPFC_Visium_example_subset")

## For now, we fake out tissue spots in example data
spe_qc <- spe_pre_qc
spe_qc$in_tissue[spe_qc$array_col < 10] <- FALSE

## adds QC metrics to colData of the spe
spe_qc <- add_qc_metrics(spe_qc, overwrite = TRUE)
vars <- colnames(colData(spe_qc))
vars[grep("^(scran|edge)", vars)]

## visualize edge spots
vis_clus(spe_qc, sampleid = "Br6432_ant", clustervar = "edge_spot")

## specify your own colors
vis_clus(
    spe_qc,
    sampleid = "Br6432_ant",
    clustervar = "edge_spot",
    colors = c(
        "TRUE" = "lightgreen",
        "FALSE" = "pink",
        "NA" = "red"
    )
)
vis_gene(spe_qc, sampleid = "Br6432_ant", geneid = "edge_distance", minCount = -1)

## Visualize scran QC flags

## Check the spots with low library size as detected by scran::isOutlier()
vis_clus(spe_qc, sample_id = "Br6432_ant", clustervar = "scran_low_lib_size")

## Violin plot of library size with low library size highlighted in a
## different color.
scater::plotColData(spe_qc[, spe_qc$in_tissue], x = "sample_id", y = "sum_umi", colour_by = "scran_low_lib_size")

## Check any spots that scran::isOutlier() flagged
vis_clus(spe_qc, sampleid = "Br6432_ant", clustervar = "scran_discard")

## Low library spots that are on the edge of the tissue
vis_clus(spe_qc, sampleid = "Br6432_ant", clustervar = "scran_low_lib_size_edge")

## Use `low_library_size` (or other variables) and `edge_distance` as you
## please.
spe_qc$our_low_lib_edge <- spe_qc$scran_low_lib_size & spe_qc$edge_distance < 5

vis_clus(spe_qc, sample_id = "Br6432_ant", clustervar = "our_low_lib_edge")

## Clean up
rm(spe_qc, spe_pre_qc, vars)

Description

This function adds to a SPE (SpatialExperiment-class) object the output from read10xVisiumAnalysis().

Usage

add10xVisiumAnalysis(spe, visium_analysis)

Arguments

Details

You might want to use read10xVisiumWrapper() instead of using this function directly.

Value

A SpatialExperiment-class object with the clustering results from SpaceRanger added to colData(spe) and the dimension reduction results added to reducedDims(spe). Added data starts with the ⁠10x_⁠ prefix to make them easy to differentiate.

See Also

Other Utility functions for reading data from SpaceRanger output by 10x Genomics: read10xVisiumAnalysis(), read10xVisiumWrapper()

Examples

## See 'Using spatialLIBD with 10x Genomics public datasets' for
## a full example using this function.
if (interactive()) {
    browseVignettes(package = "spatialLIBD")
}

## Note that ?SpatialExperiment::read10xVisium doesn't include all the files
## we need to illustrate read10xVisiumWrapper().

Description

Once you have computed the enrichment t-statistics for your sc/snRNA-seq data using registration_wrapper() and related functions, you can then use layer_stat_cor() and layer_stat_cor_plot() to perform the spatial registartion of your sc/snRNA-seq data. This function helps interpret that matrix and assign layer labels to your clusters.

Usage

annotate_registered_clusters(
  cor_stats_layer,
  confidence_threshold = 0.25,
  cutoff_merge_ratio = 0.25
)

Arguments

Details

If you change the input modeling_results to layer_stat_cor() then the interpretation of this function could change. For example, maybe you have your own spatially-resolved transcriptomics data that doesn't have to be about DLPFC layers.

Value

A data.frame with 3 columns. Your clusters, the layer_confidence which depends on confidence_threshold, and the layer_label.

See Also

Other Layer correlation functions: layer_stat_cor(), layer_stat_cor_plot()

Examples

## Obtain the necessary data
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}

## Compute the correlations
cor_stats_layer <- layer_stat_cor(
    tstats_Human_DLPFC_snRNAseq_Nguyen_topLayer,
    modeling_results,
    model_type = "enrichment"
)

## Obtain labels
annotate_registered_clusters(cor_stats_layer)

## More relaxed merging threshold
annotate_registered_clusters(cor_stats_layer, cutoff_merge_ratio = 1)

Description

This function checks that the modeling_results object has the appropriate structure. For more details please check the vignette documentation.

Usage

check_modeling_results(modeling_results)

Arguments

Value

The input object if all checks are passed.

See Also

Other Check input functions: check_sce(), check_sce_layer(), check_spe()

Examples

if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}

## Check the object
xx <- check_modeling_results(modeling_results)

Description

This function checks that the sce object has the appropriate structure. This is a legacy function and we highly encourage you to use SpatialExperiment-class objects and check them with check_spe().

Usage

check_sce(
  sce,
  variables = c("GraphBased", "ManualAnnotation", "Maynard", "Martinowich",
    paste0("SNN_k50_k", 4:28), "spatialLIBD", "cell_count", "sum_umi", "sum_gene",
    "expr_chrM", "expr_chrM_ratio", "SpatialDE_PCA", "SpatialDE_pool_PCA", "HVG_PCA",
    "pseudobulk_PCA", "markers_PCA", "SpatialDE_UMAP", "SpatialDE_pool_UMAP", "HVG_UMAP",
    "pseudobulk_UMAP", "markers_UMAP", "SpatialDE_PCA_spatial",
    "SpatialDE_pool_PCA_spatial", "HVG_PCA_spatial", "pseudobulk_PCA_spatial",
    "markers_PCA_spatial", "SpatialDE_UMAP_spatial", "SpatialDE_pool_UMAP_spatial", 
    
    "HVG_UMAP_spatial", "pseudobulk_UMAP_spatial", "markers_UMAP_spatial")
)

Arguments

Value

The input object if all checks are passed.

See Also

Other Check input functions: check_modeling_results(), check_sce_layer(), check_spe()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("sce_example")) sce_example <- fetch_data("sce_example")

    ## Check the object
    check_sce(sce_example)
}

Description

This function checks that the sce_layer object has the appropriate structure. For more details please check the vignette documentation.

Usage

check_sce_layer(sce_layer, variables = "spatialLIBD")

Arguments

Value

The input object if all checks are passed.

See Also

Other Check input functions: check_modeling_results(), check_sce(), check_spe()

Examples

## Obtain example data from the HumanPilot project
## (Maynard, Collado-Torres, et al, 2021)
if (!exists("sce_layer")) sce_layer <- fetch_data("sce_layer")

## Check the pseudo-bulked data
check_sce_layer(sce_layer)

Description

This function checks that the spe object has the appropriate structure. For more details please check the vignette documentation.

Usage

check_spe(
  spe,
  variables = c("sum_umi", "sum_gene", "expr_chrM", "expr_chrM_ratio")
)

Arguments

Value

The input object if all checks are passed.

Author(s)

Brenda Pardo, Leonardo Collado-Torres

See Also

Other Check input functions: check_modeling_results(), check_sce(), check_sce_layer()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Check the object
    check_spe(spe)
}

Description

This function creates a clusters.csv file similar to the ones created by SpaceRanger at outs/analysis/clustering but with the key column that combines the barcode and the sample_id, which is needed when the spe object contains data from multiple samples given that the barcodes are duplicated.

Usage

cluster_export(
  spe,
  cluster_var,
  cluster_dir = file.path(tempdir(), "exported_clusters"),
  overwrite = TRUE
)

Arguments

Value

The path to the exported clusters.csv file.

See Also

Other cluster export/import utility functions: cluster_import()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Export two cluster variables
    cluster_export(spe, "spatialLIBD")
    cluster_export(spe, "GraphBased")
}

Description

This function imports previously exported clustering results with cluster_export() and adds them to the colData() slot of your SpatialExperiment-class object.

Usage

cluster_import(
  spe,
  cluster_dir = file.path(tempdir(), "exported_clusters"),
  prefix = "imported_",
  overwrite = TRUE
)

Arguments

Value

A SpatialExperiment-class object with the imported clusters appended on the colData().

See Also

Other cluster export/import utility functions: cluster_export()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Export two cluster variables
    cluster_export(spe, "spatialLIBD")
    cluster_export(spe, "GraphBased")

    ## Re-import them
    colData(cluster_import(spe))
}

Description

This function determines if you have enough RAM memory on your system.

Usage

enough_ram(how_much = 4e+09)

Arguments

Details

If benchmarkme::get_ram() fails, this function will return FALSE as a save bet.

Value

A logical(1) indicating whether your system has enough RAM memory.

Examples

## Do you have ~ 4 GB in your system?
enough_ram(4e9)

## Do you have ~ 100 GB in your system
enough_ram(100e9)

Description

This function downloads from ExperimentHub Visium, Visium Spatial Proteogenomics (Visium-SPG), or single nucleus RNA-seq (snRNA-seq) data and results analyzed by LIBD from multiple projects. If ExperimentHub is not available, this function will download the files from Dropbox using BiocFileCache::bfcrpath() unless the files are present already at destdir. Note that ExperimentHub and BiocFileCache will cache the data and automatically detect if you have previously downloaded it, thus making it the preferred way to interact with the data.

Usage

fetch_data(
  type = c("sce", "sce_layer", "modeling_results", "sce_example", "spe",
    "spatialDLPFC_Visium", "spatialDLPFC_Visium_example_subset",
    "spatialDLPFC_Visium_pseudobulk", "spatialDLPFC_Visium_modeling_results",
    "spatialDLPFC_Visium_SPG", "spatialDLPFC_snRNAseq",
    "Visium_SPG_AD_Visium_wholegenome_spe", "Visium_SPG_AD_Visium_targeted_spe",
    "Visium_SPG_AD_Visium_wholegenome_pseudobulk_spe",
    "Visium_SPG_AD_Visium_wholegenome_modeling_results", "visiumStitched_brain_spe",
    "visiumStitched_brain_spaceranger", "visiumStitched_brain_Fiji_out", 
    
    "LFF_spatial_ERC_SRT", "LFF_spatial_ERC_SRT_pseudobulk",
    "LFF_spatial_ERC_SRT_modeling_results", "LFF_spatial_ERC_snRNAseq",
    "LFF_spatial_ERC_snRNAseq_pseudobulk_broad",
    "LFF_spatial_ERC_snRNAseq_pseudobulk_subcluster",
    "LFF_spatial_ERC_snRNAseq_modeling_results_broad",
    "LFF_spatial_ERC_snRNAseq_modeling_results_subcluster"),
  destdir = tempdir(),
  eh = ExperimentHub::ExperimentHub(),
  bfc = BiocFileCache::BiocFileCache()
)

Arguments

Details

The data was initially prepared by scripts at https://github.com/LieberInstitute/HumanPilot and further refined by https://github.com/LieberInstitute/spatialLIBD/blob/master/inst/scripts/make-data_spatialLIBD.R.

Value

The requested object: sce, sce_layer, ve or modeling_results that you have to assign to an object. If you didn't you can still avoid re-loading the object by using .Last.value.

References

Please always cite the spatialLIBD publication whenever you use fetch_data() to download data.

Pardo B, Spangler A, Weber LM, Page SC, Hicks SC, Jaffe AE, Martinowich K, Maynard KR, Collado-Torres L. spatialLIBD: an R/Bioconductor package to visualize spatially-resolved transcriptomics data. BMC Genomics. 2022 Jun 10;23(1):434. doi: 10.1186/s12864-022-08601-w. PubMed PMID: 35689177; PubMed Central PMCID: PMC9188087.

Additionally, please cite the relevant publication describing the data generation and initial data analysis for the dataset you are using.

• For "sce", "sce_layer", "modeling_results", "sce_example", and "spe" which are files from the HumanPilot study https://github.com/LieberInstitute/HumanPilot please cite: Maynard KR, Collado-Torres L, Weber LM, Uytingco C, Barry BK, Williams SR, Catallini JL 2nd, Tran MN, Besich Z, Tippani M, Chew J, Yin Y, Kleinman JE, Hyde TM, Rao N, Hicks SC, Martinowich K, Jaffe AE. Transcriptome-scale spatial gene expression in the human dorsolateral prefrontal cortex. Nat Neurosci. 2021 Mar;24(3):425-436. doi: 10.1038/s41593-020-00787-0. Epub 2021 Feb 8. PubMed PMID: 33558695; PubMed Central PMCID: PMC8095368.

• For spatialDLPFC files http://research.libd.org/spatialDLPFC/ please cite: Huuki-Myers LA, Spangler A, Eagles NJ, Montgomery KD, Kwon SH, Guo B, Grant-Peters M, Divecha HR, Tippani M, Sriworarat C, Nguyen AB, Ravichandran P, Tran MN, Seyedian A, Hyde TM, Kleinman JE, Battle A, Page SC, Ryten M, Hicks SC, Martinowich K, Collado-Torres L, Maynard KR. A data-driven single-cell and spatial transcriptomic map of the human prefrontal cortex. Science. 2024 May 24;384(6698):eadh1938. doi: 10.1126/science.adh1938. Epub 2024 May 24. PubMed PMID: 38781370; PubMed Central PMCID: PMC11398705.

• Note that spatialDLPFC_snRNAseq from https://github.com/LieberInstitute/Human_DLPFC_Deconvolution was also described in the following publication: Huuki-Myers LA, Montgomery KD, Kwon SH, Cinquemani S, Eagles NJ, Gonzalez-Padilla D, Maden SK, Kleinman JE, Hyde TM, Hicks SC, Maynard KR, Collado-Torres L. Benchmark of cellular deconvolution methods using a multi-assay dataset from postmortem human prefrontal cortex. Genome Biol. 2025 Apr 7;26(1):88. doi: 10.1186/s13059-025-03552-3. PubMed PMID: 40197307; PubMed Central PMCID: PMC11978107.

• For Visium_SPG_AD files https://research.libd.org/Visium_SPG_AD please cite: Kwon SH, Parthiban S, Tippani M, Divecha HR, Eagles NJ, Lobana JS, Williams SR, Mak M, Bharadwaj RA, Kleinman JE, Hyde TM, Page SC, Hicks SC, Martinowich K, Maynard KR, Collado-Torres L. Influence of Alzheimer's disease related neuropathology on local microenvironment gene expression in the human inferior temporal cortex. GEN Biotechnol. 2023 Oct;2(5):399-417. doi: 10.1089/genbio.2023.0019. Epub 2023 Oct 16. PubMed PMID: 39329069; PubMed Central PMCID: PMC11426291.

• For visiumStitched_brain files https://research.libd.org/visiumStitched_brain/ please cite: Eagles NJ, Bach SV, Tippani M, Ravichandran P, Du Y, Miller RA, Hyde TM, Page SC, Martinowich K, Collado-Torres L. Integrating gene expression and imaging data across Visium capture areas with visiumStitched. BMC Genomics. 2024 Nov 13;25(1):1077. doi: 10.1186/s12864-024-10991-y. PubMed PMID: 39533203; PubMed Central PMCID: PMC11559125.

• For LFF_spatial_ERC files https://research.libd.org/LFF_spatial_ERC/ please cite: Huuki-Myers LA, Divecha HR, Bach SV, Valentine MR, Eagles NJ, Mulvey B, Bharadwaj RA, Zhang R, Evans JR, Grant-Peters M, Miller RA, Kleinman JE, Han S, Hyde TM, Page SC, Weinberger DR, Martinowich K, Ryten M, Maynard KR, Collado-Torres L. APOE E4 Alzheimer's Risk Converges on an Oligodendrocyte Subtype in the Human Entorhinal Cortex. bioRxiv. 2025 Nov 20;. doi: 10.1101/2025.11.20.689483. PubMed PMID: 41332786; PubMed Central PMCID: PMC12667772.

Examples

## Download the SingleCellExperiment object
## at the layer-level
if (!exists("sce_layer")) sce_layer <- fetch_data("sce_layer")

## Explore the data
sce_layer

## How to download and load "spatialDLPFC_snRNAseq"
## A similar process is needed for downloading and loading other
## HDF5SummarizedExperiment files:
## * "LFF_spatial_ERC_SRT"
## * LFF_spatial_ERC_snRNAseq"

## Not run: 
sce_path_zip <- fetch_data("spatialDLPFC_snRNAseq")
sce_path <- unzip(sce_path_zip, exdir = tempdir())
sce <- HDF5Array::loadHDF5SummarizedExperiment(
    file.path(tempdir(), "sce_DLPFC_annotated")
)
sce
#> class: SingleCellExperiment
#> dim: 36601 77604
#> metadata(3): Samples cell_type_colors cell_type_colors_broad
#> assays(2): counts logcounts
#> rownames(36601): MIR1302-2HG FAM138A ... AC007325.4 AC007325.2
#> rowData names(7): source type ... gene_type binomial_deviance
#> colnames(77604): 1_AAACCCAAGTTCTCTT-1 1_AAACCCACAAGGTCTT-1 ... 19_TTTGTTGTCTCATTGT-1 19_TTTGTTGTCTTAAGGC-1
#> colData names(32): Sample Barcode ... cellType_layer layer_annotation
#> reducedDimNames(4): GLMPCA_approx TSNE UMAP HARMONY
#> mainExpName: NULL
#> altExpNames(0):
lobstr::obj_size(sce)
#> 172.28 MB

## End(Not run)

Description

This function is useful for automatically cropping the images. It finds the edge points (min and max on both the X and Y axis) in pixels based on a particular image. This function takes advantage of the known design of Visium slides as documented at https://support.10xgenomics.com/spatial-gene-expression/software/pipelines/latest/output/spatial and https://kb.10xgenomics.com/hc/en-us/articles/360041426992. That is, that for a regular Visium slide, the array row has a range from 0 to 77, the array col from 0 to 127, the capture area has a 6.5 mm edge length, the the fiducial frame area has an edge of 8 mm, and spot centers are about 100 um from each other.

Usage

frame_limits(
  spe,
  sampleid,
  image_id = "lowres",
  visium_grid = list(row_min = 0, row_max = 77, col_min = 0, col_max = 127,
    fiducial_vs_capture_edge = (8 - 6.5) * 1000/2/100)
)

Arguments

Value

A named list with y_min, y_max, x_min, and x_max pixels from the selected image that can be used for cropping the image.

Author(s)

Louise Huuki-Myers and Leonardo Collado-Torres

See Also

Other Spatial cluster visualization functions: vis_clus(), vis_clus_p(), vis_grid_clus(), vis_image()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Obtain the frame limits for one sample
    frame_limits(spe, sampleid = "151673")
}

Description

Using the layer-level (group-level) data, this function evaluates whether list of gene sets (Ensembl gene IDs) are enriched among the significant genes (FDR < 0.1 by default) genes for a given model type result. Test the alternative hypothesis that OR > 1, i.e. that gene set is over-represented in the set of enriched genes. If you want to check depleted genes, change reverse to TRUE.

Usage

gene_set_enrichment(
  gene_list,
  fdr_cut = 0.1,
  modeling_results = fetch_data(type = "modeling_results"),
  model_type = names(modeling_results)[1],
  reverse = FALSE
)

Arguments

Details

Check https://github.com/LieberInstitute/HumanPilot/blob/master/Analysis/Layer_Guesses/check_clinical_gene_sets.R to see a full script from where this family of functions is derived from.

Value

A table in long format with the enrichment results using stats::fisher.test().

• OR odds ratio.

• Pval p-value for fisher.test().

• test group or layer in the modeling_results.

• NumSig Number of genes from the gene set present in modeling_results & with fdr < fdr_cut and t_stat > 0 (unless reverse = TRUE) for test in modeling results.

• SetSize Number of genes from modeling_results present in gene_set.

• ID name of gene set.

• model_type record of input model type from ⁠modeling results⁠.

• fdr_cut record of input frd_cut.

• GeneList List of gene names from input set with fdr < fdr_cut and t_stat > 0, possible that some gene_names are missing from modeling data.

Author(s)

Andrew E Jaffe, Leonardo Collado-Torres

See Also

Other Gene set enrichment functions: gene_set_enrichment_plot()

Examples

## Read in the SFARI gene sets included in the package
asd_sfari <- utils::read.csv(
    system.file(
        "extdata",
        "SFARI-Gene_genes_01-03-2020release_02-04-2020export.csv",
        package = "spatialLIBD"
    ),
    as.is = TRUE
)

## Format them appropriately
asd_sfari_geneList <- list(
    Gene_SFARI_all = asd_sfari$ensembl.id,
    Gene_SFARI_high = asd_sfari$ensembl.id[asd_sfari$gene.score < 3],
    Gene_SFARI_syndromic = asd_sfari$ensembl.id[asd_sfari$syndromic == 1]
)

## Obtain the necessary data
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}

## Compute the gene set enrichment results
asd_sfari_enrichment <- gene_set_enrichment(
    gene_list = asd_sfari_geneList,
    modeling_results = modeling_results,
    model_type = "enrichment"
)

## Explore the results
asd_sfari_enrichment

Description

This function takes the output of gene_set_enrichment() and creates a ComplexHeatmap visualization of the results. Fill of the heatmap represents the -log10(p-val), Odds-ratios are printed for test that pass specified significance threshold ORcut.

Usage

gene_set_enrichment_plot(
  enrichment,
  xlabs = unique(enrichment$ID),
  PThresh = 12,
  ORcut = 3,
  enrichOnly = FALSE,
  mypal = c("white", RColorBrewer::brewer.pal(9, "YlOrRd")),
  plot_SetSize_bar = FALSE,
  gene_list_length = NULL,
  model_sig_length = NULL,
  model_colors = NULL,
  ...
)

Arguments

Details

Includes functionality to plot the size of the input gene sets as barplot annotations.

Check https://github.com/LieberInstitute/HumanPilot/blob/master/Analysis/Layer_Guesses/check_clinical_gene_sets.R to see a full script from where this family of functions is derived from.

Value

A (Heatmap-class) visualizing the gene set enrichment odds ratio and p-value results.

Author(s)

Andrew E Jaffe, Leonardo Collado-Torres, Louise Huuki-Myers

See Also

Other Gene set enrichment functions: gene_set_enrichment()

Examples

## Read in the SFARI gene sets included in the package
asd_sfari <- utils::read.csv(
    system.file(
        "extdata",
        "SFARI-Gene_genes_01-03-2020release_02-04-2020export.csv",
        package = "spatialLIBD"
    ),
    as.is = TRUE
)

## Format them appropriately
asd_safari_geneList <- list(
    Gene_SFARI_all = asd_sfari$ensembl.id,
    Gene_SFARI_high = asd_sfari$ensembl.id[asd_sfari$gene.score < 3],
    Gene_SFARI_syndromic = asd_sfari$ensembl.id[asd_sfari$syndromic == 1]
)

## Obtain the necessary data
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}

## Compute the gene set enrichment results
asd_sfari_enrichment <- gene_set_enrichment(
    gene_list = asd_safari_geneList,
    modeling_results = modeling_results,
    model_type = "enrichment"
)

## Visualize the gene set enrichment results

## Default plot
gene_set_enrichment_plot(
    enrichment = asd_sfari_enrichment
)

## Use a custom green color palette & use shorter gene set names
## (x-axis labels)
gene_set_enrichment_plot(
    asd_sfari_enrichment,
    xlabs = gsub(".*_", "", unique(asd_sfari_enrichment$ID)),
    mypal = c("white", RColorBrewer::brewer.pal(9, "BuGn"))
)

## Add bar plot annotations for SetSize of model genes in the gene_lists
gene_set_enrichment_plot(
    asd_sfari_enrichment,
    xlabs = gsub(".*_", "", unique(asd_sfari_enrichment$ID)),
    plot_SetSize_bar = TRUE
)

## Add stacked bar plot annotations showing SetSize and difference from the
## length of the input gene_list
gene_set_enrichment_plot(
    asd_sfari_enrichment,
    xlabs = gsub(".*_", "", unique(asd_sfari_enrichment$ID)),
    plot_SetSize_bar = TRUE,
    gene_list_length = lapply(asd_safari_geneList, length)
)

## add bar plot annotations for number of enriched genes from layers
if (!exists("sce_layer")) sce_layer <- fetch_data(type = "sce_layer")
sig_genes <- sig_genes_extract(
    modeling_results = modeling_results,
    model = "enrichment",
    sce_layer = sce_layer,
    n = nrow(sce_layer)
)

sig_genes <- sig_genes[sig_genes$fdr < 0.1, ]
n_sig_model <- as.list(table(sig_genes$test))

## add barplot with n significant genes from modeling
gene_set_enrichment_plot(
    asd_sfari_enrichment,
    xlabs = gsub(".*_", "", unique(asd_sfari_enrichment$ID)),
    plot_SetSize_bar = TRUE,
    model_sig_length = n_sig_model
)

## add color annotations
gene_set_enrichment_plot(
    asd_sfari_enrichment,
    xlabs = gsub(".*_", "", unique(asd_sfari_enrichment$ID)),
    plot_SetSize_bar = TRUE,
    model_colors = libd_layer_colors
)

## add barplot with n significant genes from modeling filled with model color
gene_set_enrichment_plot(
    asd_sfari_enrichment,
    xlabs = gsub(".*_", "", unique(asd_sfari_enrichment$ID)),
    plot_SetSize_bar = TRUE,
    model_sig_length = n_sig_model,
    model_colors = libd_layer_colors
)

Description

This function defines a ggplot2::layer() for visualizing the histology image from Visium. It can be combined with other ggplot2 functions for visualizing the clusters as in vis_clus_p() or gene-level information as in vis_gene_p().

Usage

geom_spatial(
  mapping = NULL,
  data = NULL,
  stat = "identity",
  position = "identity",
  na.rm = FALSE,
  show.legend = NA,
  inherit.aes = FALSE,
  ...
)

Arguments

Value

A ggplot2::layer() for the histology information.

Author(s)

10x Genomics

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Select the first sample and extract the data
    sample_id <- unique(spe$sample_id)[1]
    spe_sub <- spe[, spe$sample_id == sample_id]
    sample_df <- as.data.frame(colData(spe_sub), optional = TRUE)

    ## Obtain the histology image
    img <- SpatialExperiment::imgRaster(spe_sub)

    ## Transform to a rasterGrob object
    grob <- grid::rasterGrob(img, width = grid::unit(1, "npc"), height = grid::unit(1, "npc"))

    ## Make a plot using geom_spatial
    p <- ggplot2::ggplot(
        sample_df,
        ggplot2::aes(
            x = pxl_col_in_fullres * SpatialExperiment::scaleFactors(spe_sub),
            y = pxl_row_in_fullres * SpatialExperiment::scaleFactors(spe_sub),
        )
    ) +
        geom_spatial(
            data = tibble::tibble(grob = list(grob)),
            ggplot2::aes(grob = grob),
            x = 0.5,
            y = 0.5
        )

    ## Show the plot
    print(p)

    ## Clean up
    rm(spe_sub)
}

Description

This function returns a vector of colors based on a vector of cluster names. It can be used to automatically assign colors.

Usage

get_colors(colors = NULL, clusters)

Arguments

Value

A named vector where the values are the colors to use for displaying them different clusters. For some use cases, you might have to either change the names or use unname().

Examples

## Obtain the necessary data
if (!exists("sce_layer")) sce_layer <- fetch_data("sce_layer")

## Example layer colors with the corresponding names
get_colors(libd_layer_colors, sce_layer$layer_guess)
get_colors(libd_layer_colors, sce_layer$layer_guess_reordered_short)

## Example where colors are assigned automatically
## based on a pre-defined set of colors
get_colors(clusters = sce_layer$kmeans_k7)

## Example where Polychrome::palette36.colors() gets used
get_colors(clusters = letters[seq_len(13)])

## What happens if you have a logical variable with NAs?
set.seed(20240712)
log_var <- sample(c(TRUE, FALSE, NA),
    1000,
    replace = TRUE,
    prob = c(0.3, 0.15, 0.55)
)
log_var_sorted <- sort_clusters(log_var)
get_colors(colors = NULL, clusters = log_var_sorted)

Description

This function uses the magick package to edit the color and perform other image manipulations on a background image. It can be useful if you want to highlight certain features of these images.

Usage

img_edit(
  spe,
  sampleid,
  image_id = "lowres",
  channel = NA,
  brightness = 100,
  saturation = 100,
  hue = 100,
  enhance = FALSE,
  contrast_sharpen = NA,
  quantize_max = NA,
  quantize_dither = TRUE,
  equalize = FALSE,
  normalize = FALSE,
  transparent_color = NA,
  transparent_fuzz = 0,
  background_color = NA,
  median_radius = NA,
  negate = FALSE
)

Arguments

Details

The magick functions are used in the sequence represented by the arguments to this function. You can alternatively use this function sequentially. Or directly use the magick package.

Value

A magick image object such as the one returned by magick::image_read.

See Also

Other Image editing functions: img_update(), img_update_all()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Reduce brightness to 25%
    x <- img_edit(spe, sampleid = "151507", brightness = 25)
    plot(x)
}

Description

Edit the image with img_edit() then update the imgData().

Usage

img_update(
  spe,
  sampleid,
  image_id = "lowres",
  new_image_id = paste0("edited_", image_id),
  overwrite = FALSE,
  ...
)

Arguments

Value

A SpatialExperiment-class object with an updated imgData() slot.

See Also

Other Image editing functions: img_edit(), img_update_all()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Reduce brightness to 25% and update the imgData()
    imgData(img_update(spe, sampleid = "151507", brightness = 25))
}

Description

This function uses img_update() for all samples. That is, it loops through every sample and edits the image with img_edit() and then updates the imgData().

Usage

img_update_all(
  spe,
  image_id = "lowres",
  new_image_id = paste0("edited_", image_id),
  overwrite = FALSE,
  ...
)

Arguments

Value

A SpatialExperiment-class object with an updated imgData() slot.

See Also

Other Image editing functions: img_edit(), img_update()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Reduce brightness to 25% for the 'lowres' image for all samples and
    ## update the imgData()
    imgData(img_update_all(spe, brightness = 25))
}

Description

This function uses the output of sig_genes_extract_all() as well as the logcounts from the layer-level (group-level) data to visualize the expression of a given gene and display the modeling results for the given gene.

Usage

layer_boxplot(
  i = 1,
  sig_genes = sig_genes_extract(),
  short_title = TRUE,
  sce_layer = fetch_data(type = "sce_layer"),
  col_bkg_box = "grey80",
  col_bkg_point = "grey40",
  col_low_box = "violet",
  col_low_point = "darkviolet",
  col_high_box = "skyblue",
  col_high_point = "dodgerblue4",
  cex = 2,
  group_var = "layer_guess_reordered_short",
  assayname = "logcounts"
)

Arguments

Value

This function creates a boxplot of the layer-level data (group-level) separated by layer and colored based on the model type from row i of sig_genes.

References

Adapted from https://github.com/LieberInstitute/HumanPilot/blob/master/Analysis/Layer_Guesses/layer_specificity.R

See Also

Other Layer modeling functions: sig_genes_extract(), sig_genes_extract_all()

Examples

## Obtain the necessary data
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}
if (!exists("sce_layer")) sce_layer <- fetch_data(type = "sce_layer")

## Top 2 genes from the enrichment model
sig_genes <- sig_genes_extract_all(
    n = 2,
    modeling_results = modeling_results,
    sce_layer = sce_layer
)

## Example default boxplot
set.seed(20200206)
layer_boxplot(sig_genes = sig_genes, sce_layer = sce_layer)

## Now show the long title version
set.seed(20200206)
layer_boxplot(
    sig_genes = sig_genes,
    short_title = FALSE,
    sce_layer = sce_layer
)

set.seed(20200206)
layer_boxplot(
    i = which(sig_genes$model_type == "anova")[1],
    sig_genes = sig_genes,
    sce_layer = sce_layer
)

set.seed(20200206)
layer_boxplot(
    i = which(sig_genes$model_type == "pairwise")[1],
    sig_genes = sig_genes,
    sce_layer = sce_layer
)

## Viridis colors displayed in the shiny app
library("viridisLite")
set.seed(20200206)
layer_boxplot(
    sig_genes = sig_genes,
    sce_layer = sce_layer,
    col_low_box = viridis(4)[2],
    col_low_point = viridis(4)[1],
    col_high_box = viridis(4)[3],
    col_high_point = viridis(4)[4]
)

## Paper colors displayed in the shiny app
set.seed(20200206)
layer_boxplot(
    sig_genes = sig_genes,
    sce_layer = sce_layer,
    col_low_box = "palegreen3",
    col_low_point = "springgreen2",
    col_high_box = "darkorange2",
    col_high_point = "orange1"
)

## Blue/red colors displayed in the shiny app
set.seed(20200206)
layer_boxplot(
    i = which(sig_genes$model_type == "pairwise")[1],
    sig_genes = sig_genes,
    sce_layer = sce_layer,
    col_bkg_box = "grey90",
    col_bkg_point = "grey60",
    col_low_box = "skyblue2",
    col_low_point = "royalblue3",
    col_high_box = "tomato2",
    col_high_point = "firebrick4",
    cex = 3
)

Description

Layer modeling correlation of statistics

Usage

layer_stat_cor(
  stats,
  modeling_results = fetch_data(type = "modeling_results"),
  model_type = names(modeling_results)[1],
  reverse = FALSE,
  top_n = NULL
)

Arguments

Details

Check https://github.com/LieberInstitute/HumanPilot/blob/master/Analysis/Layer_Guesses/dlpfc_snRNAseq_annotation.R for a full analysis from which this family of functions is derived from.

Value

A correlation matrix between the query stats and the reference statistics using only the ENSEMBL gene IDs present in both tables. The columns are sorted using hierarchical clustering.

Author(s)

Andrew E Jaffe, Leonardo Collado-Torres

See Also

Other Layer correlation functions: annotate_registered_clusters(), layer_stat_cor_plot()

Examples

## Obtain the necessary data
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}

## Compute the correlations
cor_stats_layer <- layer_stat_cor(
    tstats_Human_DLPFC_snRNAseq_Nguyen_topLayer,
    modeling_results,
    model_type = "enrichment"
)

## Explore the correlation matrix
head(cor_stats_layer[, seq_len(3)])
summary(cor_stats_layer)

## Repeat with top_n set to 10
summary(layer_stat_cor(
    tstats_Human_DLPFC_snRNAseq_Nguyen_topLayer,
    modeling_results,
    model_type = "enrichment",
    top_n = 10
))

Description

This function makes a ComplexHeatmap from the correlation matrix between a reference and query modeling statistics from layer_stat_cor(). For example, between the query statistics from a set of cell cluster/types derived from scRNA-seq or snRNA-seq data (among other types) and the reference layer statistics from the Human DLPFC Visium data (when using the default arguments).

Usage

layer_stat_cor_plot(
  cor_stats_layer,
  color_max = max(cor_stats_layer),
  color_min = min(cor_stats_layer),
  color_scale = RColorBrewer::brewer.pal(7, "PRGn"),
  query_colors = NULL,
  reference_colors = NULL,
  annotation = NULL,
  ...
)

Arguments

Details

Includes functionality to add color annotations, (helpful to match to colors in Visium spot plots), and annotations from annotate_registered_clusters().

Value

(Heatmap-class) plot of t-stat correlations

Author(s)

Louise Huuki-Myers

See Also

Other Layer correlation functions: annotate_registered_clusters(), layer_stat_cor()

Examples

## Obtain the necessary data
## reference human pilot modeling results
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}

## query spatialDLPFC modeling results
query_modeling_results <- fetch_data(
    type = "spatialDLPFC_Visium_modeling_results"
)

## Compute the correlations
cor_stats_layer <- layer_stat_cor(
    stats = query_modeling_results$enrichment,
    modeling_results,
    model_type = "enrichment"
)

## Visualize the correlation matrix

## Default plot with no annotations and defaults for ComplexHeatmap()
layer_stat_cor_plot(cor_stats_layer)

## add Annotation colors
## add libd_layer_colors to reference Human Pilot layers
layer_stat_cor_plot(cor_stats_layer, reference_colors = libd_layer_colors)

## obtain colors for the query clusters
cluster_colors <- get_colors(clusters = rownames(cor_stats_layer))
layer_stat_cor_plot(cor_stats_layer,
    query_colors = cluster_colors,
    reference_colors = libd_layer_colors
)

## Apply additional ComplexHeatmap param
layer_stat_cor_plot(cor_stats_layer,
    cluster_rows = FALSE,
    cluster_columns = FALSE
)

## Add annotation
annotation_df <- annotate_registered_clusters(
    cor_stats_layer,
    confidence_threshold = .55
)
layer_stat_cor_plot(cor_stats_layer, annotation = annotation_df)

## change fill color scale
layer_stat_cor_plot(cor_stats_layer,
    color_scale = RColorBrewer::brewer.pal(2, "PiYG")
)

## All together
layer_stat_cor_plot(
    cor_stats_layer,
    color_scale = RColorBrewer::brewer.pal(5, "PiYG"),
    query_colors = cluster_colors,
    reference_colors = libd_layer_colors,
    annotation = annotation_df,
    cluster_rows = FALSE,
    cluster_columns = FALSE
)

Description

A named vector of colors to use for the LIBD layers designed by Lukas M. Weber with feedback from the spatialLIBD collaborators.

Usage

libd_layer_colors

Format

A vector of length 9 with colors for Layers 1 through 9, WM, NA and a special WM2 that is present in some of the unsupervised clustering results.

Description

Creates a named character() vector that can be helpful for locating image files and used with add_images(). This function is not necessary if the image files don't use the spe$sample_id.

Usage

locate_images(spe, image_dir, image_pattern)

Arguments

Value

A named character() vector with the path to images.

See Also

Other Functions for adding non-standard images: add_images()

Examples

## Not run: 
locate_images(spe, tempdir(), "testImage")

## End(Not run)

Description

This function expands SpatialExperiment::read10xVisium() by reading analysis outputs from SpaceRanger by 10x Genomics.

Usage

read10xVisiumAnalysis(
  samples = "",
  sample_id = paste0("sample", sprintf("%02d", seq_along(samples)))
)

Arguments

Details

You might want to use read10xVisiumWrapper() instead of using this function directly.

Value

A named list() with the information about the clustering and the dimension reduction (projections) from the SpaceRanger output by 10x Genomics.

See Also

Other Utility functions for reading data from SpaceRanger output by 10x Genomics: add10xVisiumAnalysis(), read10xVisiumWrapper()

Examples

## See 'Using spatialLIBD with 10x Genomics public datasets' for
## a full example using this function.
if (interactive()) {
    browseVignettes(package = "spatialLIBD")
}

## Note that ?SpatialExperiment::read10xVisium doesn't include all the files
## we need to illustrate read10xVisiumWrapper().

Description

This function expands SpatialExperiment::read10xVisium() to include analysis results from SpaceRanger by 10x Genomics as well as add information needed by run_app() to visualize the data with the spatialLIBD shiny web application.

Usage

read10xVisiumWrapper(
  samples = "",
  sample_id = paste0("sample", sprintf("%02d", seq_along(samples))),
  type = c("HDF5", "sparse"),
  data = c("filtered", "raw"),
  images = c("lowres", "hires", "detected", "aligned"),
  load = TRUE,
  reference_gtf = NULL,
  chrM = "chrM",
  gtf_cols = c("source", "type", "gene_id", "gene_version", "gene_name", "gene_type"),
  verbose = TRUE
)

Arguments

Value

A SpatialExperiment object with the clustering and dimension reduction (projection) results from SpaceRanger by 10x Genomics as well as other information used by run_app() for visualzing the gene expression data.

See Also

Other Utility functions for reading data from SpaceRanger output by 10x Genomics: add10xVisiumAnalysis(), read10xVisiumAnalysis()

Examples

## See 'Using spatialLIBD with 10x Genomics public datasets' for
## a full example using this function.
if (interactive()) {
    browseVignettes(package = "spatialLIBD")
}

## Note that ?SpatialExperiment::read10xVisium doesn't include all the files
## we need to illustrate read10xVisiumWrapper().

Description

This function computes the block correlation using the sample ID as the blocking factor. This takes into account that cells in scRNA-seq data or spots in spatially-resolved transcriptomics data from Visium (or similar) have a sample ID batch effect.

Usage

registration_block_cor(
  sce_pseudo,
  registration_model,
  var_sample_id = "registration_sample_id"
)

Arguments

Value

A numeric(1) with the block correlation at the sample ID level.

See Also

Other spatial registration and statistical modeling functions: registration_model(), registration_pseudobulk(), registration_stats_anova(), registration_stats_enrichment(), registration_stats_pairwise(), registration_wrapper()

Examples

example("registration_model", package = "spatialLIBD")
block_cor <- registration_block_cor(sce_pseudo, registration_mod)

Description

This function defines the statistical model that will be used for computing the block correlation as well as pairwise statistics. It is useful to check it in case your sample-level covariates need to be casted. For example, an integer() variable might have to be casted into a factor() if you wish to model it as a categorical variable and not a continuous one.

Usage

registration_model(
  sce_pseudo,
  covars = NULL,
  var_registration = "registration_variable"
)

Arguments

Value

The output of model.matrix() which you can inspect to verify that your sample-level covariates are being properly modeled.

See Also

Other spatial registration and statistical modeling functions: registration_block_cor(), registration_pseudobulk(), registration_stats_anova(), registration_stats_enrichment(), registration_stats_pairwise(), registration_wrapper()

Examples

example("registration_pseudobulk", package = "spatialLIBD")
registration_mod <- registration_model(sce_pseudo, "age")
head(registration_mod)

Description

Pseudo-bulk the gene expression, filter lowly-expressed genes, and normalize. This is the first step for spatial registration and for statistical modeling.

Usage

registration_pseudobulk(
  sce,
  var_registration,
  var_sample_id,
  covars = NULL,
  min_ncells = 10,
  pseudobulk_rds_file = NULL,
  filter_expr = TRUE,
  mito_gene = NULL
)

Arguments

Value

A pseudo-bulked SingleCellExperiment-class object. The logcounts() assay are log2-CPM values calculated with edgeR::cpm(log = TRUE). See https://github.com/LieberInstitute/spatialLIBD/issues/106 and https://support.bioconductor.org/p/9161754 for more details about the math behind scuttle::logNormFactors(), edgeR::cpm(), and their differences.

See Also

Other spatial registration and statistical modeling functions: registration_block_cor(), registration_model(), registration_stats_anova(), registration_stats_enrichment(), registration_stats_pairwise(), registration_wrapper()

Examples

## Ensure reproducibility of example data
set.seed(20220907)

## Generate example data
sce <- scuttle::mockSCE()

## Add some sample IDs
sce$sample_id <- sample(LETTERS[1:5], ncol(sce), replace = TRUE)

## Add a sample-level covariate: age
ages <- rnorm(5, mean = 20, sd = 4)
names(ages) <- LETTERS[1:5]
sce$age <- ages[sce$sample_id]

## Add gene-level information
rowData(sce)$gene_id <- paste0("ENSG", seq_len(nrow(sce)))
rowData(sce)$gene_name <- paste0("gene", seq_len(nrow(sce)))

## Pseudo-bulk by Cell Cycle
sce_pseudo <- registration_pseudobulk(
    sce,
    var_registration = "Cell_Cycle",
    var_sample_id = "sample_id",
    covars = c("age"),
    min_ncells = NULL
)
colData(sce_pseudo)
rowData(sce_pseudo)

Description

This function computes the gene ANOVA F-statistics (at least one group is different from the rest). These F-statistics can be used for spatial registration with layer_stat_cor() and related functions. Although, they are more typically used for identifying ANOVA-marker genes.

Usage

registration_stats_anova(
  sce_pseudo,
  block_cor,
  covars = NULL,
  var_registration = "registration_variable",
  var_sample_id = "registration_sample_id",
  gene_ensembl = NULL,
  gene_name = NULL,
  suffix = ""
)

Arguments

Value

A data.frame() with the ANOVA statistical results. This is similar to fetch_data("modeling_results")$anova.

See Also

Other spatial registration and statistical modeling functions: registration_block_cor(), registration_model(), registration_pseudobulk(), registration_stats_enrichment(), registration_stats_pairwise(), registration_wrapper()

Examples

example("registration_block_cor", package = "spatialLIBD")
results_anova <- registration_stats_anova(sce_pseudo,
    block_cor, "age",
    gene_ensembl = "ensembl", gene_name = "gene_name", suffix = "example"
)
head(results_anova)

## Specifying `block_cor = NaN` then ignores the correlation structure
results_anova_nan <- registration_stats_anova(sce_pseudo,
    block_cor = NaN, "age",
    gene_ensembl = "ensembl", gene_name = "gene_name", suffix = "example"
)
head(results_anova_nan)

## Note that you can merge multiple of these data.frames if you run this
## function for different sets. For example, maybe you drop one group
## before pseudo-bulking if you know that there are many differences between
## that group and others. For example, we have dropped the white matter (WM)
## prior to computing ANOVA F-statistics.

## no covariates
results_anova_nocovar <- registration_stats_anova(sce_pseudo,
    block_cor,
    covars = NULL,
    gene_ensembl = "ensembl", gene_name = "gene_name", suffix = "nocovar"
)
head(results_anova_nocovar)

## Merge both results into a single data.frame, thanks to having different
## 'suffix' values.
results_anova_merged <- merge(results_anova, results_anova_nocovar)
head(results_anova_merged)

Description

This function computes the gene enrichment t-statistics (one group > the rest). These t-statistics are the ones typically used for spatial registration with layer_stat_cor() and related functions.

Usage

registration_stats_enrichment(
  sce_pseudo,
  block_cor,
  covars = NULL,
  var_registration = "registration_variable",
  var_sample_id = "registration_sample_id",
  gene_ensembl = NULL,
  gene_name = NULL
)

Arguments

Value

A data.frame() with the enrichment statistical results. This is similar to fetch_data("modeling_results")$enrichment.

See Also

Other spatial registration and statistical modeling functions: registration_block_cor(), registration_model(), registration_pseudobulk(), registration_stats_anova(), registration_stats_pairwise(), registration_wrapper()

Examples

example("registration_block_cor", package = "spatialLIBD")
results_enrichment <- registration_stats_enrichment(sce_pseudo,
    block_cor, "age",
    gene_ensembl = "ensembl", gene_name = "gene_name"
)
head(results_enrichment)

## Specifying `block_cor = NaN` then ignores the correlation structure
results_enrichment_nan <- registration_stats_enrichment(sce_pseudo,
    block_cor = NaN, "age",
    gene_ensembl = "ensembl", gene_name = "gene_name"
)
head(results_enrichment_nan)

Description

This function computes the gene pairwise t-statistics (one group > another, for all combinations). These t-statistics can be used for spatial registration with layer_stat_cor() and related functions. Although, they are more typically used for identifying pairwise-marker genes.

Usage

registration_stats_pairwise(
  sce_pseudo,
  registration_model,
  block_cor,
  var_registration = "registration_variable",
  var_sample_id = "registration_sample_id",
  gene_ensembl = NULL,
  gene_name = NULL
)

Arguments

Value

A data.frame() with the pairwise statistical results. This is similar to fetch_data("modeling_results")$pairwise.

See Also

Other spatial registration and statistical modeling functions: registration_block_cor(), registration_model(), registration_pseudobulk(), registration_stats_anova(), registration_stats_enrichment(), registration_wrapper()

Examples

example("registration_block_cor", package = "spatialLIBD")
results_pairwise <- registration_stats_pairwise(sce_pseudo,
    registration_mod, block_cor,
    gene_ensembl = "ensembl", gene_name = "gene_name"
)
head(results_pairwise)

## Specifying `block_cor = NaN` then ignores the correlation structure
results_pairwise_nan <- registration_stats_pairwise(sce_pseudo,
    registration_mod,
    block_cor = NaN,
    gene_ensembl = "ensembl", gene_name = "gene_name"
)
head(results_pairwise_nan)

Description

This function is provided for convenience. It runs all the functions required for computing the modeling_results. This can be useful for finding marker genes on a new spatially-resolved transcriptomics dataset and thus using it for run_app(). The results from this function can also be used for performing spatial registration through layer_stat_cor() and related functions of sc/snRNA-seq datasets.

Usage

registration_wrapper(
  sce,
  var_registration,
  var_sample_id,
  covars = NULL,
  gene_ensembl = NULL,
  gene_name = NULL,
  suffix = "",
  min_ncells = 10,
  pseudobulk_rds_file = NULL
)

Arguments

Details

We chose a default of min_ncells = 10 based on OSCA from section 4.3 at http://bioconductor.org/books/3.15/OSCA.multisample/multi-sample-comparisons.html. They cite https://doi.org/10.1038/s41467-020-19894-4 as the paper where they came up with the definition of "very low" being 10. You might want to use registration_pseudobulk() and manually explore sce_pseudo$ncells to choose the best cutoff.

Value

A list() of data.frame() with the statistical results. This is similar to fetch_data("modeling_results").

See Also

Other spatial registration and statistical modeling functions: registration_block_cor(), registration_model(), registration_pseudobulk(), registration_stats_anova(), registration_stats_enrichment(), registration_stats_pairwise()

Examples

## Ensure reproducibility of example data
set.seed(20220907)

## Generate example data
sce <- scuttle::mockSCE()

## Add some sample IDs
sce$sample_id <- sample(LETTERS[1:5], ncol(sce), replace = TRUE)

## Add a sample-level covariate: age
ages <- rnorm(5, mean = 20, sd = 4)
names(ages) <- LETTERS[1:5]
sce$age <- ages[sce$sample_id]

## Add gene-level information
rowData(sce)$ensembl <- paste0("ENSG", seq_len(nrow(sce)))
rowData(sce)$gene_name <- paste0("gene", seq_len(nrow(sce)))

## Compute all modeling results
example_modeling_results <- registration_wrapper(
    sce,
    var_registration = "Cell_Cycle",
    var_sample_id = "sample_id",
    covars = c("age"),
    gene_ensembl = "ensembl",
    gene_name = "gene_name",
    suffix = "wrapper"
)

## Explore the results from registration_wrapper()
class(example_modeling_results)
length(example_modeling_results)
names(example_modeling_results)
lapply(example_modeling_results, head)

Description

This function runs the shiny application that allows users to interact with the Visium spatial transcriptomics data from LIBD (by default) or any other data that you have shaped according to our object structure.

Usage

run_app(
  spe = fetch_data(type = "spe"),
  sce_layer = fetch_data(type = "sce_layer"),
  modeling_results = fetch_data(type = "modeling_results"),
  sig_genes = sig_genes_extract_all(n = nrow(sce_layer), modeling_results =
    modeling_results, sce_layer = sce_layer),
  docs_path = system.file("app", "www", package = "spatialLIBD"),
  title = "spatialLIBD",
  spe_discrete_vars = c("spatialLIBD", "GraphBased", "ManualAnnotation", "Maynard",
    "Martinowich", paste0("SNN_k50_k", 4:28), "SpatialDE_PCA", "SpatialDE_pool_PCA",
    "HVG_PCA", "pseudobulk_PCA", "markers_PCA", "SpatialDE_UMAP", "SpatialDE_pool_UMAP",
    "HVG_UMAP", "pseudobulk_UMAP", "markers_UMAP", "SpatialDE_PCA_spatial",
    "SpatialDE_pool_PCA_spatial", "HVG_PCA_spatial", "pseudobulk_PCA_spatial",
    "markers_PCA_spatial", "SpatialDE_UMAP_spatial", "SpatialDE_pool_UMAP_spatial",
    "HVG_UMAP_spatial", "pseudobulk_UMAP_spatial", 
     "markers_UMAP_spatial"),
  spe_continuous_vars = c("cell_count", "sum_umi", "sum_gene", "expr_chrM",
    "expr_chrM_ratio"),
  default_cluster = "spatialLIBD",
  auto_crop_default = TRUE,
  is_stitched = FALSE,
  ...
)

Arguments

Details

If you don't have the pseudo-bulked analysis results like we computed them in our project https://doi.org/10.1038/s41593-020-00787-0 you can set sce_layer, modeling_results and sig_genes to NULL. Doing so will disable the pseudo-bulked portion of the web application. See the examples for one such case as well as the vignette that describes how you can use spatialLIBD with public data sets provided by 10x Genomics. That vignette is available at http://research.libd.org/spatialLIBD/articles/TenX_data_download.html.

Value

A shiny.appobj that contains the input data.

Examples

## Not run: 
## The default arguments will download the data from the web
## using fetch_data(). If this is the first time you have run this,
## the files will need to be cached by ExperimentHub. Otherwise it
## will re-use the files you have previously downloaded.
if (enough_ram(4e9)) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Create the interactive website
    run_app(spe)

    ## You can also run a custom version without the pseudo-bulked
    ## layer information. This is useful if you are only interested
    ## in the spatial transcriptomics features.
    run_app(spe,
        sce_layer = NULL, modeling_results = NULL, sig_genes = NULL,
        title = "spatialLIBD without layer info"
    )

    ## When using shinyapps.io aim for less than 3 GB of RAM with your
    ## objects. Check each input object with:
    ## lobstr::obj_size(x)
    ## Do not create the large input objects on the app.R script before
    ## subsetting them. Do this outside app.R since the app.R script is
    ## run at shinyapps.io, so subsetting on that script to reduce the
    ## memory load is pointless. You have to do it outside of app.R.
}

## How to run locally the spatialDLPFC Sp09 spatialLIBD app. That is,
## from http://research.libd.org/spatialDLPFC/#interactive-websites
## how to run https://libd.shinyapps.io/spatialDLPFC_Visium_Sp09 locally.
if (enough_ram(9e9)) {
    ## Download the 3 main objects needed
    spe <- fetch_data("spatialDLPFC_Visium")
    sce_pseudo <- fetch_data("spatialDLPFC_Visium_pseudobulk")
    modeling_results <- fetch_data("spatialDLPFC_Visium_modeling_results")

    ## These are optional commands to further reduce the memory required.
    #
    ## Keep only the "lowres" images. Reduces the object from 6.97 GB to 4.59 GB
    # imgData(spe) <- imgData(spe)[imgData(spe)$image_id == "lowres", ]
    ## Drop the regular counts (keep only the logcounts). Reduces the object
    ## from 4.59 GB to 2.45 GB.
    # counts(spe) <- NULL

    ## For sig_genes_extract_all() to work
    sce_pseudo$spatialLIBD <- sce_pseudo$BayesSpace
    ## Compute the significant genes
    sig_genes <- sig_genes_extract_all(
        n = nrow(sce_pseudo),
        modeling_results = modeling_results,
        sce_layer = sce_pseudo
    )
    ## Reduce the memory from 423.73 MB to 78.88 MB
    lobstr::obj_size(sig_genes)
    sig_genes$in_rows <- NULL
    sig_genes$in_rows_top20 <- NULL
    lobstr::obj_size(sig_genes)

    ## Specify the default variable
    spe$BayesSpace <- spe$BayesSpace_harmony_09
    ## Get all variables
    vars <- colnames(colData(spe))

    ## Set default cluster colors
    colors_BayesSpace <- Polychrome::palette36.colors(28)
    names(colors_BayesSpace) <- c(1:28)
    m <- match(as.character(spe$BayesSpace_harmony_09), names(colors_BayesSpace))
    stopifnot(all(!is.na(m)))
    spe$BayesSpace_colors <- spe$BayesSpace_harmony_09_colors <- colors_BayesSpace[m]

    ## Download documentation files we use
    temp_www <- file.path(tempdir(), "www")
    dir.create(temp_www)
    download.file(
        "https://raw.githubusercontent.com/LieberInstitute/spatialDLPFC/main/README.md",
        file.path(temp_www, "README.md")
    )
    download.file(
        "https://raw.githubusercontent.com/LieberInstitute/spatialDLPFC/main/code/deploy_app_k09/www/documentation_sce_layer.md",
        file.path(temp_www, "documentation_sce_layer.md")
    )
    download.file(
        "https://raw.githubusercontent.com/LieberInstitute/spatialDLPFC/main/code/deploy_app_k09/www/documentation_spe.md",
        file.path(temp_www, "documentation_spe.md")
    )
    download.file(
        "https://raw.githubusercontent.com/LieberInstitute/spatialDLPFC/main/img/favicon.ico",
        file.path(temp_www, "favicon.ico")
    )
    download.file(
        "https://raw.githubusercontent.com/LieberInstitute/spatialDLPFC/main/code/deploy_app_k09/www/footer.html",
        file.path(temp_www, "footer.html")
    )
    list.files(temp_www)

    ## Run the app locally
    run_app(
        spe,
        sce_layer = sce_pseudo,
        modeling_results = modeling_results,
        sig_genes = sig_genes,
        title = "spatialDLPFC, Visium, Sp09",
        spe_discrete_vars = c( # this is the variables for the spe object not the sce_pseudo object
            "BayesSpace",
            "ManualAnnotation",
            vars[grep("^SpaceRanger_|^scran_", vars)],
            vars[grep("^BayesSpace_harmony", vars)],
            vars[grep("^BayesSpace_pca", vars)],
            "graph_based_PCA_within",
            "PCA_SNN_k10_k7",
            "Harmony_SNN_k10_k7",
            "manual_layer_label",
            "wrinkle_type",
            "BayesSpace_colors"
        ),
        spe_continuous_vars = c(
            "sum_umi",
            "sum_gene",
            "expr_chrM",
            "expr_chrM_ratio",
            vars[grep("^VistoSeg_", vars)],
            vars[grep("^layer_", vars)],
            vars[grep("^broad_", vars)]
        ),
        default_cluster = "BayesSpace",
        docs_path = temp_www
    )
}
## See also:
## * https://github.com/LieberInstitute/spatialDLPFC/tree/main/code/deploy_app_k09
## * https://github.com/LieberInstitute/spatialDLPFC/tree/main/code/deploy_app_k09_position
## * https://github.com/LieberInstitute/spatialDLPFC/tree/main/code/deploy_app_k09_position_noWM
## * https://github.com/LieberInstitute/spatialDLPFC/tree/main/code/deploy_app_k16
## * https://github.com/LieberInstitute/spatialDLPFC/tree/main/code/analysis_IF/03_spatialLIBD_app


## Example for an object with multiple capture areas stitched together with
## <http://research.libd.org/visiumStitched/>.
spe_stitched <- fetch_data("visiumStitched_brain_spe")

## Inspect this object
spe_stitched

## Notice the use of "exclude_overlapping"
table(spe_stitched$exclude_overlapping, useNA = "ifany")

## Run the app with this stitched data
run_app(
    spe = spe_stitched,
    sce_layer = NULL, modeling_results = NULL, sig_genes = NULL,
    title = "visiumStitched example data",
    spe_discrete_vars = c("capture_area", "scran_quick_cluster", "ManualAnnotation"),
    spe_continuous_vars = c("sum_umi", "sum_gene", "expr_chrM", "expr_chrM_ratio"),
    default_cluster = "scran_quick_cluster",
    is_stitched = TRUE
)

## End(Not run)

Description

This function converts a spot-level SingleCellExperiment-class (SCE) object as generated by fetch_data() to a SpatialExperiment-class (SPE) object.

Usage

sce_to_spe(sce = fetch_data("sce"), imageData = NULL)

Arguments

Details

Note that the resulting object is a bit more complex than a regular SPE because it contains the data from the spatialLIBD project which you might otherwise have to generate for your own data.

Value

A a SpatialExperiment-class object.

Author(s)

Brenda Pardo, Leonardo Collado-Torres

Examples

if (enough_ram()) {
    ## Download the sce data
    sce <- fetch_data("sce")
    ## Transform it to a SpatialExperiment object
    spe <- sce_to_spe(sce)
}

Description

From the layer-level modeling results, this function extracts the top n significant genes. This is the workhorse function used by sig_genes_extract_all() through which we obtain the information that can then be used by functions such as layer_boxplot() for constructing informative titles.

Usage

sig_genes_extract(
  n = 10,
  modeling_results = fetch_data(type = "modeling_results"),
  model_type = names(modeling_results)[1],
  reverse = FALSE,
  sce_layer = fetch_data(type = "sce_layer"),
  gene_name = "gene_name"
)

Arguments

Value

A data.frame() with the top n significant genes (as ordered by their statistics in decreasing order) in long format. The specific columns are described further in the vignette.

References

Adapted from https://github.com/LieberInstitute/HumanPilot/blob/master/Analysis/Layer_Guesses/layer_specificity_functions.R

See Also

Other Layer modeling functions: layer_boxplot(), sig_genes_extract_all()

Examples

## Obtain the necessary data
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}
if (!exists("sce_layer")) sce_layer <- fetch_data(type = "sce_layer")

## anova top 10 genes
sig_genes_extract(
    modeling_results = modeling_results,
    sce_layer = sce_layer
)

## Extract all genes
sig_genes_extract(
    modeling_results = modeling_results,
    sce_layer = sce_layer,
    n = nrow(sce_layer)
)

Description

This function combines the output of sig_genes_extract() from all the layer-level (group-level) modeling results and builds the data required for functions such as layer_boxplot().

Usage

sig_genes_extract_all(
  n = 10,
  modeling_results = fetch_data(type = "modeling_results"),
  sce_layer = fetch_data(type = "sce_layer"),
  gene_name = "gene_name"
)

Arguments

Value

A DataFrame-class with the extracted statistics in long format. The specific columns are described further in the vignette.

See Also

Other Layer modeling functions: layer_boxplot(), sig_genes_extract()

Examples

## Obtain the necessary data
if (!exists("modeling_results")) {
    modeling_results <- fetch_data(type = "modeling_results")
}
if (!exists("sce_layer")) sce_layer <- fetch_data(type = "sce_layer")

## top 10 genes for all models
sig_genes_extract_all(
    modeling_results = modeling_results,
    sce_layer = sce_layer
)

Description

This function takes a vector with cluster labels, recasts it as a factor(), and sorts the factor() levels by frequency such that the most frequent cluster is the first level and so on.

Usage

sort_clusters(clusters, map_subset = NULL)

Arguments

Value

A factor() version of clusters where the levels are ordered by frequency.

Examples

## Build an initial set of cluster labels
clus <- letters[unlist(lapply(4:1, function(x) rep(x, x)))]

## In this case, it's a character vector
class(clus)

## We see that we have 10 elements in this vector, which is
## an unnamed character vector
clus

## letter 'd' is the most frequent
table(clus)

## Sort them and obtain a factor. Notice that it's a named
## factor, and the names correspond to the original values
## in the character vector.
sort_clusters(clus)

## Since 'd' was the most frequent, it gets assigned to the first level
## in the factor variable.
table(sort_clusters(clus))

## If we skip the first 3 values of clus (which are all 'd'), we can
## change the most frequent cluster. And thus the ordering of the
## factor levels.
sort_clusters(clus, map_subset = seq_len(length(clus)) > 3)

## Let's try with a factor variable
clus_factor <- factor(clus)
## sort_clusters() returns an identical result in this case
stopifnot(identical(sort_clusters(clus), sort_clusters(clus_factor)))

## What happens if you have a logical variable with NAs?
set.seed(20240712)
log_var <- sample(c(TRUE, FALSE, NA),
    1000,
    replace = TRUE,
    prob = c(0.3, 0.15, 0.55)
)
## Here, the NAs are the most frequent group.
table(log_var, useNA = "ifany")

## The NAs are not used for sorting. Since we have more 'TRUE' than 'FALSE'
## then, 'TRUE' becomes the first level.
table(sort_clusters(log_var), useNA = "ifany")

Description

Using the DLPFC snRNA-seq data from Matthew N Tran et al https://doi.org/10.1016/j.neuron.2021.09.001 we computed enrichment t-statistics for the cell clusters. The Tran et al data has been subset to the top 100 DLPFC layer markers found in Maynard, Collado-Torres, et al 2021. This data is used in examples such as in layer_stat_cor_plot(). The Tran et al data is from the pre-print version of that project.

Usage

tstats_Human_DLPFC_snRNAseq_Nguyen_topLayer

Format

A matrix with 692 rows and 31 variables where each column is a given cell cluster from Tran et al and each row is one gene. The row names are Ensembl gene IDs which are used by layer_stat_cor() to match to our modeling results.

Source

https://github.com/LieberInstitute/HumanPilot/blob/master/Analysis/Layer_Guesses/dlpfc_snRNAseq_annotation.R and https://github.com/LieberInstitute/spatialLIBD/blob/master/dev/02_dev.R#L107-L194.

Description

This function visualizes the clusters for one given sample at the spot-level using (by default) the histology information on the background. To visualize gene-level (or any continuous variable) use vis_gene().

Usage

vis_clus(
  spe,
  sampleid = unique(spe$sample_id)[1],
  clustervar,
  colors = c("#b2df8a", "#e41a1c", "#377eb8", "#4daf4a", "#ff7f00", "gold", "#a65628",
    "#999999", "black", "grey", "white", "purple"),
  spatial = TRUE,
  image_id = "lowres",
  alpha = NA,
  point_size = 2,
  auto_crop = TRUE,
  na_color = "#CCCCCC40",
  is_stitched = FALSE,
  guide_point_size = point_size,
  ...
)

Arguments

Details

This function subsets spe to the given sample and prepares the data and title for vis_clus_p().

Value

A ggplot2 object.

See Also

Other Spatial cluster visualization functions: frame_limits(), vis_clus_p(), vis_grid_clus(), vis_image()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Check the colors defined by Lukas M Weber
    libd_layer_colors

    ## Use the manual color palette by Lukas M Weber
    p1 <- vis_clus(
        spe = spe,
        clustervar = "layer_guess_reordered",
        sampleid = "151673",
        colors = libd_layer_colors,
        ... = " LIBD Layers"
    )
    print(p1)

    ## Without auto-cropping the image
    p2 <- vis_clus(
        spe = spe,
        clustervar = "layer_guess_reordered",
        sampleid = "151673",
        colors = libd_layer_colors,
        auto_crop = FALSE,
        ... = " LIBD Layers"
    )
    print(p2)

    ## Without histology
    p3 <- vis_clus(
        spe = spe,
        clustervar = "layer_guess_reordered",
        sampleid = "151673",
        colors = libd_layer_colors,
        ... = " LIBD Layers",
        spatial = FALSE
    )
    print(p3)

    ## With some NA values
    spe$tmp <- spe$layer_guess_reordered
    spe$tmp[spe$sample_id == "151673"][seq_len(500)] <- NA
    p4 <- vis_clus(
        spe = spe,
        clustervar = "tmp",
        sampleid = "151673",
        colors = libd_layer_colors,
        na_color = "white",
        ... = " LIBD Layers"
    )
    print(p4)
    
    ## edit plot point size but keep guide size larger
    p5 <- vis_clus(
        spe = spe,
        clustervar = "layer_guess_reordered",
        sampleid = "151673",
        colors = libd_layer_colors,
        na_color = "white",
        point_size = 1,
        guide_point_size = 3,
        ... = " LIBD Layers"
    )
    print(p5)
         
}

Description

This function visualizes the clusters for one given sample at the spot-level using (by default) the histology information on the background. This is the function that does all the plotting behind vis_clus(). To visualize gene-level (or any continuous variable) use vis_gene_p().

Usage

vis_clus_p(
  spe,
  d,
  clustervar,
  sampleid = unique(spe$sample_id)[1],
  colors,
  spatial,
  title,
  image_id = "lowres",
  alpha = NA,
  point_size = 2,
  auto_crop = TRUE,
  na_color = "#CCCCCC40"
)

Arguments

Value

A ggplot2 object.

See Also

Other Spatial cluster visualization functions: frame_limits(), vis_clus(), vis_grid_clus(), vis_image()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")
    spe_sub <- spe[, spe$sample_id == "151673"]

    ## Use the manual color palette by Lukas M Weber
    ## Don't plot the histology information
    p <- vis_clus_p(
        spe = spe_sub,
        d = as.data.frame(cbind(colData(spe_sub), SpatialExperiment::spatialCoords(spe_sub)), optional = TRUE),
        clustervar = "layer_guess_reordered",
        sampleid = "151673",
        colors = libd_layer_colors,
        title = "151673 LIBD Layers",
        spatial = FALSE
    )
    print(p)

    ## Clean up
    rm(spe_sub)
}

Description

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().

Usage

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",
  multi_gene_method = c("z_score", "pca", "sparsity"),
  is_stitched = FALSE,
  cap_percentile = 1,
  ...
)

Arguments

Details

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.

Value

A ggplot2 object.

See Also

Other Spatial gene visualization functions: vis_gene_p(), vis_grid_gene()

Examples

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)

    #    Define several markers for white matter
    white_matter_genes <- c(
        "ENSG00000197971", "ENSG00000131095", "ENSG00000123560",
        "ENSG00000171885"
    )

    ## Plot all white matter markers at once using the Z-score combination
    ## method. Flatten this quantity at the top 5% of values for plotting
    p6 <- vis_gene(
        spe = spe,
        sampleid = "151507",
        geneid = white_matter_genes,
        multi_gene_method = "z_score",
        cap_percentile = 0.95
    )
    print(p6)

    ## Plot all white matter markers at once using the sparsity combination
    ## method
    p7 <- vis_gene(
        spe = spe,
        sampleid = "151507",
        geneid = white_matter_genes,
        multi_gene_method = "sparsity"
    )
    print(p7)

    ## Plot all white matter markers at once using the PCA combination
    ## method
    p8 <- vis_gene(
        spe = spe,
        sampleid = "151507",
        geneid = white_matter_genes,
        multi_gene_method = "pca"
    )
    print(p8)
}

Description

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. This is the function that does all the plotting behind vis_gene() To visualize clusters (or any discrete variable) use vis_clus_p().

Usage

vis_gene_p(
  spe,
  d,
  sampleid = unique(spe$sample_id)[1],
  spatial,
  title,
  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",
  legend_title = ""
)

Arguments

Value

A ggplot2 object.

See Also

Other Spatial gene visualization functions: vis_gene(), vis_grid_gene()

Examples

if (enough_ram()) {
    ## Obtain the necessary data
    if (!exists("spe")) spe <- fetch_data("spe")

    ## Prepare the data for the plotting function
    spe_sub <- spe[, spe$sample_id == "151673"]
    df <- as.data.frame(cbind(colData(spe_sub), SpatialExperiment::spatialCoords(spe_sub)), optional = TRUE)
    df$COUNT <- df$expr_chrM_ratio

    ## Don't plot the histology information
    p <- vis_gene_p(
        spe = spe_sub,
        d = df,
        sampleid = "151673",
        title = "151673 chrM expr ratio",
        spatial = FALSE
    )
    print(p)

    ## Clean up
    rm(spe_sub)
}

Description

This function visualizes the clusters for a set of samples at the spot-level using (by default) the histology information on the background. To visualize gene-level (or any continuous variable) use vis_grid_gene().

Usage

vis_grid_clus(
  spe,
  clustervar,
  pdf_file,
  sort_clust = TRUE,
  colors = NULL,
  return_plots = FALSE,
  spatial = TRUE,
  height = 24,
  width = 36,
  image_id = "lowres",
  alpha = NA,
  sample_order = unique(spe$sample_id),
  point_size = 2,
  auto_crop = TRUE,
  na_color = "#CCCCCC40",
  is_stitched = FALSE,
  guide_point_size = point_size,
  ...
)

Arguments