In ventri2020/rraysplot: Make 2D Arrays Plots Less Aggravating

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.path = "man/figures/10_images_to_tensors-",
  out.width = "100%"
)

`%>%` <- purrr::`%>%`

library(tidyverse, warn.conflicts = FALSE)
library(ANTsRCore)

library(devtools)
library(lobstr)
library(fs)

library(rraysplot)
packageVersion("rraysplot")

Tensor data for U-Net

Source 2.2.11 Data representations for neural networks in Deep Learning with R by François Chollet with J. J. Allaire.

Images typically have three dimensions:

height × width × color channels

Although grayscale images have only a single color channel and could thus be stored in 2D tensors, by convention image tensors are always 3D, with a one-dimensional color channel for grayscale images.

A batch of 80 grayscale images of size 768×384 could thus be stored in a tensor of shape

(80, 768, 384, 1)

and a batch of 80 color images could be stored in a tensor of shape

(80, 768, 384, 3)

and a batch of 80 mask --

(80, 768, 384, 2)

Note

We will crop all images to the 768x384 size and prepare tensors for U-Net NN in dimensions:

768x384, 384x192 and 192x96

Note that both dimensions 768x384 are divisble by 128 (768/2^7=6, 384/2^7=3).

Helper functions

crop_image_wxhx1 <- function(img, ll = c(34, 56), wh = c(768, 384)) {
  ll = c(ll, 1)
  ur = ll + c(wh, 0) - c(1, 1, 0)
  ANTsRCore::cropIndices(img, ll, ur)
}

get_channel <- function(img, channel = 1) {
  ANTsRCore::splitChannels(img)[[channel]]
}

imageList2sliceList <- function(
    iList,
    channel = 1,
    z_slice = 1,
    do_cropping = TRUE,
    lower_left = c(34, 56),
    width_x_height = c(768, 384)
) {
  if (do_cropping == TRUE) {
    iList <- purrr::map(
      iList,
      crop_image_wxhx1, ll = lower_left, wh = width_x_height
    )
  }
  iList <- purrr::map(iList, get_channel, ch = channel)
  iList <- purrr::map(iList,
    ANTsRCore::extractSlice, slice = z_slice, direction = 3
  )
}

imageList2MaskList <- function(iList, clean_up = 2) {
  purrr::map(iList, ANTsRCore::getMask, cleanup = clean_up)
}

imageList2arrayList <- function(iList) {
  purrr::map(iList, as.array)
}

arrayList2Tensor <- function(aList, dim_names = info_mri$patient) {
  tensor <- array(
    data = NA,
    dim = c(length(aList), dim(aList[[1]]))
  )

  for (i in seq_along(aList)) {
    tensor[i,,] <- aList[[i]]
  }
  dimnames(tensor)[[1]] <- dim_names

  tensor
}

Read MR Images Data into Tibble

info <- images_info("../120_images", extension = "dcm")

info_mri <- info %>% filter(kind == "MRI")
info_mri

info_scat <- info %>% filter(kind == "SCAT")
info_scat

info_vsat <- info %>% filter(kind == "VSAT")
info_vsat

Create MRI, Mask, SCAT, and VSAT Lists

mri_list  <- imageFileNames2ImageList(info_mri[["file_path"]]) %>%
  imageList2sliceList(channel = 1)
mask_list <- imageList2MaskList(mri_list, clean_up = 2)

scat_list <- imageFileNames2ImageList(info_scat[["file_path"]]) %>%
  imageList2sliceList(channel = 1) %>%
  imageList2MaskList(clean_up = 0)

vsat_list <- imageFileNames2ImageList(info_vsat[["file_path"]]) %>%
  imageList2sliceList(channel = 3) %>%
  imageList2MaskList(clean_up = 0)

We are going to remove the background noise from image. Why? see the figure below.

The background noise and artifacts (forearms here) present in image (left) and removed from image (right).

invisible(plot(mri_list[[1]], doCropping = FALSE))
bhistogram(as.array(mri_list[[1]]))

img1r <- mri_list[[1]] * mask_list[[1]]

invisible(plot(img1r, doCropping = FALSE))
bhistogram(as.array(img1r))

Background noise and artifacts removal

Removing the background noise and artifacts from MR images.

mri_list <- purrr::map2(mri_list, mask_list, `*`)
scat_list <- purrr::map2(scat_list, mask_list, `*`)
vsat_list <- purrr::map2(vsat_list, mask_list, `*`)

lobstr::obj_sizes(mri_list, mask_list, scat_list, vsat_list)

Basic checks on generated images

testthat::expect_equal(length(mri_list), 120)
testthat::expect_equal(length(mask_list), 120)
testthat::expect_equal(length(scat_list), 120)
testthat::expect_equal(length(vsat_list), 120)

dims <- purrr::map(mri_list, dim)
testthat::expect_setequal(dims[1], dims)
testthat::expect_setequal(dims[1], purrr::map(mask_list, dim))
testthat::expect_setequal(dims[1], purrr::map(scat_list, dim))
testthat::expect_setequal(dims[1], purrr::map(vsat_list, dim))

n <- 53
is_na <- purrr::compose(is.na, as.numeric)
testthat::expect_setequal(is_na(mri_list[[n]]), FALSE)
testthat::expect_setequal(is_na(mask_list[[n]]), FALSE)
testthat::expect_setequal(is_na(scat_list[[n]]), FALSE)
testthat::expect_setequal(is_na(vsat_list[[n]]), FALSE)

in_range <- function(img, rng = 0:255) all(unique(img) %in% rng)

testthat::expect_setequal(
  purrr::map(mri_list, in_range, rng = 0:255), TRUE
)
testthat::expect_setequal(
  purrr::map(mask_list, in_range, rng = 0:1), TRUE
)
testthat::expect_setequal(
  purrr::map(scat_list, in_range, rng = 0:1), TRUE
)
testthat::expect_setequal(
  purrr::map(vsat_list, in_range, rng = 0:1), TRUE
)

Alpha blending masks with images

$$ {\rm blended\ image}{\rm RGB} = (1 - \alpha)\cdot {\rm image}{\rm RGB} + \alpha\cdot {\rm background\ mask}_{\rm RGB} $$

Set alpha == 0 -- for transparent background mask.
Set alpha == 1 -- for fully opaque mask.

An example of blending antsImages.

i = 53
alpha = 0.25

mri <- antsImageClone(mri_list[[i]])
mask <- antsImageClone(mask_list[[i]])
scat <- antsImageClone(scat_list[[i]])
vsat <- antsImageClone(vsat_list[[i]])

# Neon Rainbow
# https://www.schemecolor.com/neon-rainbow.php

plotBlendedImages(mri, mask, alpha = 0.35, title = i, title_size = 24)
plotBlendedImages(mri, scat, alpha = 0.35, title = i, title_size = 24)
plotBlendedImages(mri, vsat, alpha = 0.35, title = i, title_size = 24)

Image Lists to Tensors

mask_alist <- imageList2arrayList(mask_list)
mri_alist <- imageList2arrayList(mri_list)
scat_alist <- imageList2arrayList(scat_list)
vsat_alist <- imageList2arrayList(vsat_list)

lobstr::obj_sizes(mri_alist, mask_alist, scat_alist, vsat_alist)

mri_tensor <- arrayList2Tensor(mri_alist)
mask_tensor <- arrayList2Tensor(mask_alist)
scat_tensor <- arrayList2Tensor(scat_alist)
vsat_tensor <- arrayList2Tensor(vsat_alist)

str(mri_tensor)
lobstr::obj_sizes(mri_tensor, mask_tensor, scat_tensor, vsat_tensor)

plot2_array2d(
  mri_tensor[29,,], mask_tensor[29,,],
  title = dimnames(mri_tensor)[[1]][[29]],
  title_size = 24
)
plot2_array2d(
  mri_tensor[29,,], scat_tensor[29,,],
  title = dimnames(mri_tensor)[[1]][[29]],
  title_size = 24
)
plot2_array2d(
  mri_tensor[29,,], vsat_tensor[29,,],
  title = dimnames(mri_tensor)[[1]][[29]],
  title_size = 24
)

table(mask_tensor[29,,])
table(scat_tensor[29,,])
table(vsat_tensor[29,,])

Save tensors in `data-raw` directory

str(mri_tensor)
str(mask_tensor)
str(scat_tensor)
str(vsat_tensor)

fat120_768x384 = list(
  image = mri_tensor,
  mask = mask_tensor,
  scat = scat_tensor,
  vsat = vsat_tensor
)

# usethis::use_data_raw(name = "fat120_768x384")
usethis::use_data(fat120_768x384, overwrite = TRUE)
cat("Move the fat120_768x384 from data/ to data-raw/ directory.\n")
fs::file_move(
  "data/fat120_768x384.rda", "data-raw/fat120_768x384.rda"
)