Artifact removal for mammogram images

Objectives & Principles

Analyse a mammogram (monochrome) image and remove any artifacts that might result from labelling, noise or other scanning artifacts.

The working principle is that the largest region (contiguous set of pixels) is the principal subject and that any other regions are artifacts.

The steps these:

Create a mask

Create an image mask that can be used to identify regions.

1. First we convert the image to a format suitable for thresholding. For monochrome mammogram images we convert to LAB colorspace and use the L (lightness) band to form the mask.
2. Then we gaussian blur the mask to remove image noise or scanning defects.
3. Finally, we threshold the mask to get pure white for the mask and pure black for the rest.

Find image regions

1. Using Vix.Vips.Operation.labelregions/1 we identify the regions in the image. This returns an image whose content is the region definitions. We need to extract those into bounding boxes.
2. Bounding boxes are extracted through a libvips technique that relies on being able to generate an image whose pixels are their own coordinates with Vix.Vips.Operation.xyz/2
3. Zipping the results of :min and :max Vix.Vips.Operation.hist_find_indexed/2 calls with the mask and coordinate image we can return a list of {[x_min, y_min], [x_max, y_max]} tuples.
4. We drop the first two entries in the list because the list of regions will include the zero-pixel origin as well as a region that is the entire image.

Removing the artifacts

Assuming artifacts are those regions that are smaller than the largest region, we

1. Sort the regions by desceding size and drop the first entry.
2. For each entry, create a new image that is the same size as the artifact, with a color that is equal to the dominant color of the original image (which for a mammogram would likely be black).
3. Then we compose this created image over the top or the original image in the same place as the artifact - effectively masking it out.

Install Dependencoes

Mix.install([:image, :kino])

Select the image

image_input = Kino.Input.image("Image to segment")

Create a mask

# Threshold beyond which the pixel is
# included in the mask. Typically a value 
# between 0 (for pure B&W images) and 60
# will be reasonable. This is highly dependant
# on a specific image. Basically this is the
# value used to decide what is the image and
# what is the background.
threshold = 60

# See Image.blur/2. Blurring the image
# before masking will cater for image noise.
blur_sigma = 1.5
min_amplitude = 0.01

image = 
  image_input
  |> Kino.Input.read()
  |> Image.from_kino!() 
  |> Image.flatten!()

# For monochrome mammograms we may not 
# get much benefit from lab space, but
# working with just the luminance band
# reflects the monochrome nature of the
# image. We might use other strategies to
# create a mask with color data.
lab_image = Image.to_colorspace!(image, :lab)

# Mask just using the luminance band. Since
# we are using a luminance mask we just start
# with the "L" band of the image. Thresholding
# with a value > 0 means we're using non-black
# pixels as the mask.
mask = 
  lab_image[0]
  |> Image.blur!(sigma: blur_sigma, min_amplitude: min_amplitude)
  |> Image.Math.greater_than!(threshold)

Image.Kino.show(mask)

Find the bounding boxes of artifacts

{labels, _regions} = Vix.Vips.Operation.labelregions!(mask)

# Operation.xyz/2 creates an image where each pixel
# values are coordinates.
xy = Vix.Vips.Operation.xyz!(Image.width(labels), Image.height(labels))

hist_xy_min = 
  xy
  |> Vix.Vips.Operation.hist_find_indexed!(labels, combine: :VIPS_COMBINE_MIN) 
  |> Vix.Vips.Image.to_list!()
  |> hd()

hist_xy_max = 
  xy
  |> Vix.Vips.Operation.hist_find_indexed!(labels, combine: :VIPS_COMBINE_MAX) 
  |> Vix.Vips.Image.to_list!()
  |> hd()

# Drop the first two bounding boxes since they
# are the origin and the full image. The remaining
# regions are what we want.
bounding_boxes = 
  hist_xy_min
  |> Enum.zip(hist_xy_max) 
  |> Enum.drop(2)
  |> IO.inspect(label: "Region bounding boxes {[x_min, y_min], [x_max, y_max]}")

# Draw the bounding boxes over the original
# image.
{:ok, overlay} = 
  image
  |> Vix.Vips.Operation.copy!()
  |> Image.mutate(fn i ->
    Enum.reduce(bounding_boxes, i, fn {[x_min, y_min], [x_max, y_max]}, acc ->
      x_min = trunc(x_min)
      y_min = trunc(y_min)
      width = trunc(x_max - x_min)
      height = trunc(y_max - y_min)

      Image.Draw.rect!(acc, x_min, y_min, width, height, fill: false, color: :red, stroke_width: 3)
    end)
  end)

Image.Kino.show(overlay)

Identify and remove artifacts

For the purposes of this example lets assume that the subject of the image is the largest region and that all the other regions are artifacts.

It may be that other strategies are required. Perhaps regions under a certain area to capture small regions only. It may also be useful to Image.trim/2 or Image.find_trim/2 on the original image too.

# Sort the regions by their area, descending in size
# then drop the first one (the largest area)
area = fn {[x_min, y_min], [x_max, y_max]} ->
  width = trunc(x_max - x_min)
  height = trunc(y_max - y_min)
  width * height
end

artifacts =
  bounding_boxes
  |> Enum.sort(fn a, b -> area.(a) > area.(b) end)
  |> Enum.drop(1)

# Get the dominant color on the basis that is
# most likely the right color for a mammogram or
# other medical scan.
background_color = Image.dominant_color!(image)

# For each artifact, create an image of the same
# size as the artifact, that has the color of the
# dominant color, and compose it in the same place
# as the artifact thereby masking it out of the image.
image_bands = Image.bands(image)

masked_image =
  Enum.reduce artifacts, image, fn {[x_min, y_min], [x_max, y_max]}, i ->
    width = trunc(x_max - x_min)
    height = trunc(y_max - y_min)

    {:ok, mask} = Image.new(width, height, color: background_color, bands: image_bands)
    Image.compose!(i, mask, x: trunc(x_min), y: trunc(y_min))
  end

Image.Kino.show(masked_image)