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.

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

spe

A SpatialExperiment-class object. See fetch_data() for how to download some example objects or read10xVisiumWrapper() to read in spaceranger --count output files and build your own spe object.

sampleid

A character(1) specifying which sample to plot from colData(spe)$sample_id (formerly colData(spe)$sample_name).

image_id

A character(1) with the name of the image ID you want to use in the background.

visium_grid

A named list with the parameters known about the Visium grid. This can change for Visium HD vs regular Visium for example.

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.

See also

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

Author

Louise Huuki-Myers and Leonardo Collado-Torres

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")
}
#> 2024-10-31 20:26:04.315375 loading file /github/home/.cache/R/BiocFileCache/f6f1918e278_Human_DLPFC_Visium_processedData_sce_scran_spatialLIBD.Rdata%3Fdl%3D1
#> $y_min
#> [1] 64
#> 
#> $y_max
#> [1] 575
#> 
#> $x_min
#> [1] 62
#> 
#> $x_max
#> [1] 588
#>