Skip to contents

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.6

  • magick: 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