Images and spot coordinates
Last compiled: 18 September 2024
images_and_coordinates.rmd
Load required libraries:
Image transformation
semla
uses the magick
R package to read and
process images. Please have a look at this
vignette if you want to know more about magick
.
Let’s load an H&E image and visualize it:
he <- file.path("https://data.mendeley.com/public-files/datasets/kj3ntnt6vb/files",
"d97fb9ce-eb7d-4c1f-98e0-c17582024a40/file_downloaded")
im <- image_read(he)
# get information about image
info <- image_info(im)
info
## # A tibble: 1 × 7
## format width height colorspace matte filesize density
## <chr> <int> <int> <chr> <lgl> <int> <chr>
## 1 JPEG 1882 2000 sRGB FALSE 1030876 72x72
# Plot HE image
par(mar = c(0, 0, 1, 0))
im |> as.raster() |> plot()
title(paste0(info$width, "x", info$height))
The H&E image is 2000 pixels high and 1882 pixels wide. If we
want to scale this image, we can use image_scale()
from the
magick
R package:
# Scale down and plot HE image
im_small <- im |> image_scale("400")
info_small <- image_info(im_small)
info_small
## # A tibble: 1 × 7
## format width height colorspace matte filesize density
## <chr> <int> <int> <chr> <lgl> <int> <chr>
## 1 JPEG 400 425 sRGB FALSE 0 72x72
par(mar = c(0, 0, 1, 0))
im_small |> as.raster() |> plot()
title(paste0(info_small$width, "x", info_small$height))
We can also apply various types of transformations:
# Rotate image
im_rot <- im_small |> image_rotate(degrees = 45)
# Mirror image along x axis
im_xrefl <- im_small |> image_flop()
# Mirror image along y axis
im_yrefl <- im_small |> image_flip()
par(mfrow = c(2, 2), mar = c(0, 0, 2, 0))
im_small |> as.raster() |> plot()
title(paste0("original ", info_small$width, "x", info_small$height))
im_rot |> as.raster() |> plot()
title(paste0("rotated 45 degrees (",
image_info(im_rot)$width, "x", image_info(im_rot)$height, ")"))
im_xrefl |> as.raster() |> plot()
title(paste0("reflected along x axis (",
image_info(im_xrefl)$width, "x", image_info(im_xrefl)$height, ")"))
im_yrefl |> as.raster() |> plot()
title(paste0("reflected along y axis (",
image_info(im_yrefl)$width, "x", image_info(im_yrefl)$height, ")"))
When working with multiple SRT data sets, it can be useful to align
the tissue sections to have roughly the same orientation and size. In
other words, we might want to register our H&E images to a reference
image. When doing so, we want to apply rotations, reflections and
translations. The issue with the default magick
rotation
function is that it creates a bounding box to hold the entire image,
which makes the image bigger. In the plots above you cans see that the
rotated image is 586x586 pixels in size compared to the 400x425 of the
other images.
The main issue is that if we want to map spots to our H&E images after rotation, the scale has now changed and the spots will appear smaller. Instead, it is more convenient to keep the same input dimensions.
semla
provide two functions to complement the
transformation functions available in the magick
R package:
ImageTranslate()
and ImageTransform()
.
ImageTranslate
can be used to move H&E image. If you
imagine that you have an art board with the same dimensions as your
image, this function can be used to move the H&E image on the art
board. However, only the part of the image that is inside the art board
will be kept.
ImageTransform
can be used to rotate and translate
H&E images. Again, the transformations occur inside the art board
and only the part of the image still inside the art board will be
returned.
# Move image 100 pixels to the right and 100 pixels down
im_rot <- im_small |> ImageTranslate(xy_offset = c(100, 100))
# Rotate image 45 degrees
im_transf <- im_small |> ImageTransform(angle = 45, xy_offset = c(0, 0))
# Rotate image 45 degrees and move image 100 pixels to the left and 100 pixels up
im_transf2 <- im_small |> ImageTransform(angle = 45, xy_offset = c(-100, -100))
par(mfrow = c(2, 2), mar = c(0, 0, 2, 0))
im_small |> as.raster() |> plot()
title(paste0("original (", info_small$width, "x", info_small$height, ")"))
im_rot |> as.raster() |> plot()
title(paste0("moved towards bottom right corner (",
image_info(im_rot)$width, "x", image_info(im_rot)$height, ")"))
im_transf |> as.raster() |> plot()
title(paste0("rotated 45 degrees clockwise (",
image_info(im_transf)$width, "x", image_info(im_transf)$height, ")"))
im_transf2 |> as.raster() |> plot()
title(paste0("rotated 45 degrees clockwise and moved (",
image_info(im_transf2)$width, "x", image_info(im_transf2)$height, ")"))
As you can see, regardless of what transformation we apply, the image dimensions remain the same.
Spot transformation
Let’s load some spot coordinates for the H&E image that we have.
# get example coordinate file
cordinatefile <- system.file("extdata/mousebrain/spatial",
"tissue_positions_list.csv",
package = "semla")
# Load coordinates
xy <- LoadSpatialCoordinates(coordinatefiles = cordinatefile, verbose = T)
xy
## # A tibble: 2,560 × 7
## barcode selected y x pxl_row_in_fullres pxl_col_in_fullres sampleID
## <chr> <int> <int> <int> <int> <int> <int>
## 1 CATACAAA… 1 13 35 4117 6086 1
## 2 CTGAGCAA… 1 15 25 4472 5062 1
## 3 GGGTACCC… 1 14 26 4294 5164 1
## 4 ACGGAATT… 1 15 27 4472 5266 1
## 5 GGGCGGTC… 1 14 28 4294 5369 1
## 6 ATGTTACG… 1 15 29 4472 5471 1
## 7 AACCATGG… 1 14 30 4294 5574 1
## 8 TCGCATCC… 1 15 31 4473 5676 1
## 9 ACTTAGTA… 1 14 32 4295 5778 1
## 10 GAGCTCTC… 1 15 33 4473 5881 1
## # ℹ 2,550 more rows
Here we have access to the Visium grid coordinates and the H&E
image coordinates for the full resolution image used for
spaceranger count
. We also have a column with spot barcodes
and a column called selected
which holds information about
what sots are located under the tissue.
We can illustrate what this means with a plot:
ggplot(xy, aes(pxl_col_in_fullres, pxl_row_in_fullres, color = factor(selected))) +
geom_point()
The spots with a value of 1 correspond to spots under the tissue. But right now, the tissue is upside down in the plot relative to our H&E image. This is because the origin (0, 0) of the plot is located in the bottom left corner, but for images its in the upper right corner.
We can fix this easily by inverting the y axis. However, to do this properly, we need the dimensions of the H&E image…
Unfortunately, we don’t have access to this information right now so we need to load the scalefactors and an H&E image provided in the spaceranger output folder.
scalefactorfile <- system.file("extdata/mousebrain/spatial",
"scalefactors_json.json",
package = "semla")
# read scalefactors
scalefactors <- jsonlite::read_json(scalefactorfile)
scalefactors
## $spot_diameter_fullres
## [1] 143.3171
##
## $tissue_hires_scalef
## [1] 0.1039393
##
## $fiducial_diameter_fullres
## [1] 214.9757
##
## $tissue_lowres_scalef
## [1] 0.03118179
Now we can see that the scaling factor between the original H&E image and the tissue_lowres image is ~0.03. Let’s load the tissue_lowres and convert our coordinates to fit the image:
lowresimagefile <- system.file("extdata/mousebrain/spatial",
"tissue_lowres_image.jpg",
package = "semla")
# Load image
im <- image_read(lowresimagefile)
image_info(im)
## # A tibble: 1 × 7
## format width height colorspace matte filesize density
## <chr> <int> <int> <chr> <lgl> <int> <chr>
## 1 JPEG 565 600 sRGB FALSE 106233 72x72
# Convert coordinates
xy <- xy |>
mutate(across(pxl_col_in_fullres:pxl_row_in_fullres,
~ .x*scalefactors$tissue_lowres_scalef))
There are a couple of important things to pay attention to here. The
x, y coordinates are now transformed to fit our H&E image and
therefore we can set the limits of the plot to be the same as the
H&E image dimensions (see limits
in
scale_*_continuous()
). We also need to set
expand = c(0, 0)
to make sure that the margins are removed
from the plot area. We also need to invert the y axis which we can do
now that we have the H&E image height
(image_info(im)$height
).
g <- im |>
rasterGrob(width = unit(1, "npc"), height = unit(1, "npc"), interpolate = TRUE)
ggplot(xy, aes(pxl_col_in_fullres,
pxl_row_in_fullres,
color = factor(selected))) +
annotation_custom(g, -Inf, Inf, -Inf, Inf) +
geom_point() +
scale_x_continuous(limits = c(0, image_info(im)$width), expand = c(0, 0)) +
scale_y_reverse(limits = c(image_info(im)$height, 0), expand = c(0, 0))
There are also functions available to apply transformations to spots:
CoordTransform()
and CoordMirror()
.
CoordTransform
is equivalent to
ImageTransform
but for spot coordinates, meaning that you
can apply rotations and translations. The main difference is that
CoordTransform
rotates around a predefined center. If you
want apply the same rotation to an H&E image and its corresponding
spot coordinates, you want to set the center
argument to be
the center of the H&E image.
Let’s demonstrate this with our spot coordinates:
# Select only x, y coordinates
xy_coords <- xy |>
head(n = 1500) |>
select(pxl_col_in_fullres, pxl_row_in_fullres)
# Midpoint
c_xy <- colMeans(xy_coords)
# Apply transformation to apot coordinates
xy_transformed <- CoordTransform(xy_coords, angle = 45, xy_offset = c(0, 0))
# Plot spot coordinates
gg <- rbind(cbind(xy_coords |> setNames(c("x", "y")), type = "original"),
cbind(xy_transformed |> setNames(c("x", "y")), type = "transformed"))
ggplot(gg, aes(x, y, color = type)) +
geom_point() +
geom_point(aes(x = c_xy[1], y = c_xy[2]), color = "red",
size = 5, shape = 4, stroke = 3) +
facet_grid(~type) +
scale_x_continuous(limits = c(0, image_info(im)$width), expand = c(0, 0)) +
scale_y_reverse(limits = c(image_info(im)$height, 0), expand = c(0, 0))
## Warning in geom_point(aes(x = c_xy[1], y = c_xy[2]), color = "red", size = 5, : All aesthetics have length 1, but the data has 3000 rows.
## ℹ Did you mean to use `annotate()`?
The spots are rotated around the center of our spots (red cross). Instead, we want to define a center to rotate the spots around, which will the the center of our H&E image.
# Select only x, y coordinates
xy_coords <- xy |>
head(n = 1500) |>
select(pxl_col_in_fullres, pxl_row_in_fullres)
# Apply transformation to apot coordinates
xy_transformed <- CoordTransform(xy_coords, angle = 45, xy_offset = c(0, 0),
center = c(image_info(im)$width/2,
image_info(im)$height - image_info(im)$height/2))
# Plot spot coordinates
gg <- rbind(cbind(xy_coords |> setNames(c("x", "y")), type = "original"),
cbind(xy_transformed |> setNames(c("x", "y")), type = "transformed"))
ggplot(gg, aes(x, y, color = type)) +
geom_point() +
geom_point(aes(x = image_info(im)$width/2,
y = image_info(im)$height/2),
color = "blue", size = 5, stroke = 3, shape = 4) +
facet_grid(~type) +
scale_x_continuous(limits = c(0, image_info(im)$width), expand = c(0, 0)) +
scale_y_reverse(limits = c(image_info(im)$height, 0), expand = c(0, 0))
## Warning in geom_point(aes(x = image_info(im)$width/2, y = image_info(im)$height/2), : All aesthetics have length 1, but the data has 3000 rows.
## ℹ Did you mean to use `annotate()`?
Now let’s have a look at how we can transform coordinates and images at the same time.
Transform images and spots
Rotation
The CoordAndImageTransform()
makes the transformation
process a bit simpler. You can provide the H&E image and its
corresponding spot coordinates and apply transformations to both objects
simultaneously. note that the spot coordinates are still defined the
same way as before, i.e. the origin is in the top left corner. This
means that we still need to invert the y axis.
Since some of the spots are now outside the “art board”, they will be missed when drawing the plot. In this example, we lose 46 spots!
# Select only x, y coordinates
xy_coords <- xy |>
select(pxl_col_in_fullres, pxl_row_in_fullres)
# Apply transformations
transf_res <- CoordAndImageTransform(im, xy_coords, angle = 45, xy_offset_image = c(0, 0))
# Plot results
g <- transf_res$im_transf |>
rasterGrob(width = unit(1, "npc"), height = unit(1, "npc"), interpolate = TRUE)
# Add selected to transf_res$xy_transf
transf_res$xy_transf$selected <- xy$selected
# Note that the y axis still needs to be reversed
ggplot(transf_res$xy_transf, aes(tr_x, tr_y,
color = factor(selected))) +
annotation_custom(g, -Inf, Inf, -Inf, Inf) +
geom_point() +
scale_x_continuous(limits = c(0, image_info(im)$width), expand = c(0, 0)) +
scale_y_reverse(limits = c(image_info(im)$height, 0), expand = c(0, 0))
Rotation + translation
# Select only x, y coordinates
xy_coords <- xy |>
select(pxl_col_in_fullres, pxl_row_in_fullres)
# Apply transformations
transf.res <- CoordAndImageTransform(im, xy_coords, angle = 45, xy_offset_image = c(100, 100))
# Plot results
g <- transf.res$im_transf |>
rasterGrob(width = unit(1, "npc"), height = unit(1, "npc"), interpolate = TRUE)
# Add selected to transf.res$xy_transf
transf.res$xy_transf$selected <- xy$selected
ggplot(transf.res$xy_transf, aes(tr_x, image_info(im)$height - tr_y,
color = factor(selected))) +
annotation_custom(g, -Inf, Inf, -Inf, Inf) +
geom_point() +
scale_x_continuous(limits = c(0, image_info(im)$width), expand = c(0, 0)) +
scale_y_continuous(limits = c(0, image_info(im)$height), expand = c(0, 0))
# Select only x, y coordinates
xy_coords <- xy |>
select(pxl_col_in_fullres, pxl_row_in_fullres)
# Apply transformations
transf.res <- CoordAndImageTransform(im, xy_coords, angle = 45, mirror_x = TRUE)
# Plot results
g <- transf.res$im_transf |>
rasterGrob(width = unit(1, "npc"), height = unit(1, "npc"), interpolate = TRUE)
# Add selected to transf.res$xy_transf
transf.res$xy_transf$selected <- xy$selected
ggplot(transf.res$xy_transf, aes(tr_x, image_info(im)$height - tr_y,
color = factor(selected))) +
annotation_custom(g, -Inf, Inf, -Inf, Inf) +
geom_point() +
scale_x_continuous(limits = c(0, image_info(im)$width), expand = c(0, 0)) +
scale_y_continuous(limits = c(0, image_info(im)$height), expand = c(0, 0))
Package versions
semla
: 1.1.6magick
: 2.8.4
Session info
## R version 4.4.0 (2024-04-24)
## Platform: aarch64-apple-darwin20
## Running under: macOS Sonoma 14.5
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.12.0
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## time zone: Europe/Stockholm
## tzcode source: internal
##
## attached base packages:
## [1] grid stats graphics grDevices datasets utils methods
## [8] base
##
## other attached packages:
## [1] jsonlite_1.8.8 magick_2.8.4 semla_1.1.6 ggplot2_3.5.0
## [5] dplyr_1.1.4 SeuratObject_4.1.4 Seurat_4.3.0.1
##
## loaded via a namespace (and not attached):
## [1] RColorBrewer_1.1-3 rstudioapi_0.16.0 magrittr_2.0.3
## [4] spatstat.utils_3.1-0 farver_2.1.2 rmarkdown_2.28
## [7] fs_1.6.4 ragg_1.3.3 vctrs_0.6.5
## [10] ROCR_1.0-11 spatstat.explore_3.3-2 forcats_1.0.0
## [13] htmltools_0.5.8.1 curl_5.2.2 sass_0.4.9
## [16] sctransform_0.4.1 parallelly_1.38.0 KernSmooth_2.23-24
## [19] bslib_0.8.0 htmlwidgets_1.6.4 desc_1.4.3
## [22] ica_1.0-3 plyr_1.8.9 plotly_4.10.4
## [25] zoo_1.8-12 cachem_1.1.0 igraph_2.0.3
## [28] mime_0.12 lifecycle_1.0.4 pkgconfig_2.0.3
## [31] Matrix_1.7-0 R6_2.5.1 fastmap_1.2.0
## [34] fitdistrplus_1.2-1 future_1.34.0 shiny_1.9.1
## [37] digest_0.6.37 colorspace_2.1-1 patchwork_1.3.0
## [40] tensor_1.5 irlba_2.3.5.1 textshaping_0.4.0
## [43] labeling_0.4.3 progressr_0.14.0 fansi_1.0.6
## [46] spatstat.sparse_3.1-0 httr_1.4.7 polyclip_1.10-7
## [49] abind_1.4-8 compiler_4.4.0 withr_3.0.1
## [52] highr_0.11 MASS_7.3-60.2 tools_4.4.0
## [55] lmtest_0.9-40 httpuv_1.6.15 future.apply_1.11.2
## [58] goftest_1.2-3 glue_1.7.0 dbscan_1.2-0
## [61] nlme_3.1-164 promises_1.3.0 Rtsne_0.17
## [64] cluster_2.1.6 reshape2_1.4.4 generics_0.1.3
## [67] gtable_0.3.5 spatstat.data_3.1-2 tidyr_1.3.1
## [70] data.table_1.16.0 sp_2.1-4 utf8_1.2.4
## [73] spatstat.geom_3.3-2 RcppAnnoy_0.0.22 ggrepel_0.9.6
## [76] RANN_2.6.2 pillar_1.9.0 stringr_1.5.1
## [79] spam_2.10-0 later_1.3.2 splines_4.4.0
## [82] lattice_0.22-6 renv_1.0.2 survival_3.6-4
## [85] deldir_2.0-4 tidyselect_1.2.1 miniUI_0.1.1.1
## [88] pbapply_1.7-2 knitr_1.48 gridExtra_2.3
## [91] scattermore_1.2 xfun_0.47 matrixStats_1.4.1
## [94] stringi_1.8.4 lazyeval_0.2.2 yaml_2.3.10
## [97] evaluate_0.24.0 codetools_0.2-20 tibble_3.2.1
## [100] BiocManager_1.30.25 cli_3.6.3 uwot_0.2.2
## [103] xtable_1.8-4 reticulate_1.39.0 systemfonts_1.1.0
## [106] munsell_0.5.1 jquerylib_0.1.4 Rcpp_1.0.13
## [109] globals_0.16.3 spatstat.random_3.3-1 zeallot_0.1.0
## [112] png_0.1-8 spatstat.univar_3.0-1 parallel_4.4.0
## [115] pkgdown_2.1.0 dotCall64_1.1-1 listenv_0.9.1
## [118] viridisLite_0.4.2 scales_1.3.0 ggridges_0.5.6
## [121] leiden_0.4.3.1 purrr_1.0.2 rlang_1.1.4
## [124] cowplot_1.1.3 shinyjs_2.1.0