HistoImagePlot 0.99.10
library(imageFeatureTCGA)
library(HistoImagePlot)
library(dplyr)
library(SpatialExperiment)
library(ggplot2)
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("waldronlab/HistoImagePlot")
HoverNet is a deep learning model for simultaneous segmentation and classification of nuclei in multi-tissue histology images. This vignette demonstrates how to import HoverNet output files into Bioconductor’s spatial data structures and create visualizations of the segmentation results.
The package supports two main file formats:
.json or .json.gz): Contains cell coordinates, types, and
optional contours.h5ad): AnnData format with additional
computed features like mean intensity and nearest neighbor distanceSpatialExperimentobject from URLs (with Automatic Caching)
The file will be automatically cached using BiocFileCache:
hov_file <- paste0(
"https://store.cancerdatasci.org/hovernet/h5ad/",
"TCGA-23-1021-01Z-00-DX1.F07C221B-D401-47A5-9519-10DE59CA1E9D.h5ad.gz")
thumb_path <- paste0(
"https://store.cancerdatasci.org/hovernet/thumb/",
"TCGA-23-1021-01Z-00-DX1.F07C221B-D401-47A5-9519-10DE59CA1E9D.png")
hn_spe <- HoverNet(hov_file, outClass = "SpatialExperiment") |>
import()
## Registered S3 method overwritten by 'bit64':
## method from
## print.bitstring tools
## Warning in file.rename(file, dest): cannot rename file
## '/tmp/Rtmp08FLe6/file59e76628e7b5c/https:/store.cancerdatasci.org/hovernet/h5ad/TCGA-23-1021-01Z-00-DX1.F07C221B-D401-47A5-9519-10DE59CA1E9D.h5ad.gz'
## to
## '/home/biocbuild/.cache/R/BiocFileCache/https://store.cancerdatasci.org/hovernet/h5ad/TCGA-23-1021-01Z-00-DX1.F07C221B-D401-47A5-9519-10DE59CA1E9D.h5ad.gz',
## reason 'Invalid cross-device link'
The package provides a convenient function to overlay the segmentation on the original tissue thumbnail image.
thumb_path <- paste0(
"https://store.cancerdatasci.org/hovernet/thumb/",
"TCGA-23-1021-01Z-00-DX1.F07C221B-D401-47A5-9519-10DE59CA1E9D.png")
plotHoverNetH5ADOverlay(hn_spe, thumb_path)
## Warning in file.rename(file, dest): cannot rename file
## '/tmp/Rtmp08FLe6/file59e7634bf84ce/https:/store.cancerdatasci.org/hovernet/thumb/TCGA-23-1021-01Z-00-DX1.F07C221B-D401-47A5-9519-10DE59CA1E9D.png'
## to
## '/home/biocbuild/.cache/R/BiocFileCache/https://store.cancerdatasci.org/hovernet/thumb/TCGA-23-1021-01Z-00-DX1.F07C221B-D401-47A5-9519-10DE59CA1E9D.png',
## reason 'Invalid cross-device link'
plotHoverNetH5ADOverlay(
hn_spe,
thumb_path,
title = "Ovarian Cancer Tissue - Cell Segmentation",
point_size = 0.02,
legend_point_size = 3
)
custom_colors <- c(
"no label" = "#808080",
"neoplastic" = "#E31A1C",
"inflammatory" = "#1F78B4",
"stromal" = "#33A02C",
"necrotic" = "#FF7F00",
"benign epithelial" = "#6A3D9A"
)
plotHoverNetH5ADOverlay(
hn_spe,
thumb_path,
color_palette = custom_colors,
title = "Custom Color Scheme"
)
H5AD files contain additional computed features like mean intensity and nearest neighbor distance.
h5ad_coords <- data.frame(spatialCoords(hn_spe), colData(hn_spe))
p1 <- ggplot(h5ad_coords, aes(x = x_centroid, y = y_centroid,
color = mean_intensity)) +
geom_point(size = 0.5) +
scale_color_viridis_c() +
coord_fixed() +
theme_minimal() +
labs(title = "Mean Intensity", color = "Intensity")
p2 <- ggplot(h5ad_coords, aes(x = x_centroid, y = y_centroid,
color = nearest_neighbor_distance)) +
geom_point(size = 0.5) +
scale_color_viridis_c(option = "plasma") +
coord_fixed() +
theme_minimal() +
labs(title = "Nearest Neighbor Distance", color = "Distance")
cowplot::plot_grid(p1, p2, ncol = 2)
Click here for Session Info
sessionInfo()
## R version 4.6.0 alpha (2026-04-05 r89794)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.4 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.23-bioc/R/lib/libRblas.so
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0 LAPACK version 3.12.0
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB LC_COLLATE=C
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## time zone: America/New_York
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] ggplot2_4.0.2 SpatialExperiment_1.21.0
## [3] SingleCellExperiment_1.33.2 SummarizedExperiment_1.41.1
## [5] Biobase_2.71.0 GenomicRanges_1.63.2
## [7] Seqinfo_1.1.0 IRanges_2.45.0
## [9] S4Vectors_0.49.2 BiocGenerics_0.57.1
## [11] generics_0.1.4 MatrixGenerics_1.23.0
## [13] matrixStats_1.5.0 dplyr_1.2.1
## [15] HistoImagePlot_0.99.10 imageFeatureTCGA_0.99.69
## [17] BiocStyle_2.39.0
##
## loaded via a namespace (and not attached):
## [1] tidyselect_1.2.1 viridisLite_0.4.3 farver_2.1.2
## [4] blob_1.3.0 filelock_1.0.3 S7_0.2.1-1
## [7] fastmap_1.2.0 BiocFileCache_3.1.0 digest_0.6.39
## [10] lifecycle_1.0.5 RSQLite_2.4.6 magrittr_2.0.5
## [13] compiler_4.6.0 rlang_1.2.0 sass_0.4.10
## [16] tools_4.6.0 yaml_2.3.12 knitr_1.51
## [19] labeling_0.4.3 S4Arrays_1.11.1 bit_4.6.0
## [22] curl_7.0.0 reticulate_1.46.0 DelayedArray_0.37.1
## [25] RColorBrewer_1.1-3 abind_1.4-8 withr_3.0.2
## [28] purrr_1.2.2 grid_4.6.0 Rhdf5lib_1.33.6
## [31] scales_1.4.0 tinytex_0.59 dichromat_2.0-0.1
## [34] cli_3.6.6 rmarkdown_2.31 otel_0.2.0
## [37] tzdb_0.5.0 rjson_0.2.23 BiocBaseUtils_1.13.0
## [40] rhdf5_2.55.16 DBI_1.3.0 cachem_1.1.0
## [43] BiocManager_1.30.27 XVector_0.51.0 vctrs_0.7.3
## [46] Matrix_1.7-5 jsonlite_2.0.0 bookdown_0.46
## [49] hms_1.1.4 bit64_4.6.0-1 TENxIO_1.13.4
## [52] magick_2.9.1 jquerylib_0.1.4 glue_1.8.1
## [55] codetools_0.2-20 cowplot_1.2.0 gtable_0.3.6
## [58] BiocIO_1.21.0 tibble_3.3.1 pillar_1.11.1
## [61] rhdf5filters_1.23.3 rappdirs_0.3.4 htmltools_0.5.9
## [64] R6_2.6.1 dbplyr_2.5.2 httr2_1.2.2
## [67] evaluate_1.0.5 lattice_0.22-9 readr_2.2.0
## [70] png_0.1-9 memoise_2.0.1 bslib_0.10.0
## [73] rjsoncons_1.3.2 Rcpp_1.1.1-1 SparseArray_1.11.13
## [76] anndataR_1.1.2 xfun_0.57 pkgconfig_2.0.3