5. Filtering

Filtering images is at the core of image processing. The idea is to use a small image that travels across a larger image and does some operation on each sub-region of this larger image. In the simplest case, we can imagine that a small square travels across the image and computes the local mean of the pixels with which it overlaps. This is then a mean filter:

from IPython.display import HTML, display
HTML(filename='illustrations//convol_mean.html')

Multiple properties of the little travelling square can change:

  • the pixel values of the square can be anything. For example it’s common to have a 2D gaussian for smoothing the image (see below)

  • the size and shape of the square can change

  • the operation done by the square does not have to be a simple multiplication-addition (convolution). It can also be another function like a median for example.

Filters can be applied to regular images, as well as to masks (binary images) and are often used to clean-up an image or to detect object. scikit-image implements most of these filters in skimage.filters and skimage.morphology and we are going to look at a few.

import numpy as np
import matplotlib.pyplot as plt
import skimage
import skimage.io
import skimage.morphology

Linear filters

Blurring effect

Linear filters work by convoluting an image with a smaller template. That small template can be anything, but commonly it is a well-defined small matrix. For example it can be a 2D Gaussian. We can have a look at its actual shape by filtering an image that has a single pixel ON:

single_dot = np.zeros((20,20))
single_dot[10,10] = 1
plt.imshow(single_dot, cmap = 'gray');
_images/05-Filtering_8_0.png

The Gaussian has default values and we can simply apply it on an image:

filtered = skimage.filters.gaussian(single_dot)
plt.imshow(filtered, cmap = 'gray');
_images/05-Filtering_11_0.png

We can just plot a cross section to better understand what happened:

plt.plot(filtered[10,:]);
_images/05-Filtering_13_0.png

The effect of this filtering was to “widen” the original single pixel, which is now spread out across the image. We can now also specify the width of the gaussian:

filtered = skimage.filters.gaussian(single_dot,sigma=2)
plt.imshow(filtered, cmap = 'gray');
_images/05-Filtering_15_0.png

The Gaussian has become even wider. Let’s put two dots on our image and see what happens:

double_dot = np.zeros((20,20))
double_dot[10,7] = 1
double_dot[10,14] = 1
plt.imshow(double_dot, cmap = 'gray');
_images/05-Filtering_17_0.png
filtered = skimage.filters.gaussian(double_dot,sigma=2)
plt.imshow(filtered, cmap = 'gray');
_images/05-Filtering_18_0.png

The end effect is to blur our image. If we apply our filter twice we get an even blurrier image:

filtered2 = skimage.filters.gaussian(filtered,sigma=2)
plt.imshow(filtered2, cmap = 'gray');
_images/05-Filtering_20_0.png

Removing noise

Gaussian filtering does not only blur information, it also helps removing noise. Let’s for example take again our blurred image and add some noise to it. For that we use the numpy random module:

noisy = filtered + 0.02*np.random.rand(20,20)
plt.imshow(noisy, cmap = 'gray');
_images/05-Filtering_24_0.png

If we now apply again our gaussian filter:

noise_filtered = skimage.filters.gaussian(noisy,sigma=2)
plt.imshow(noise_filtered, cmap = 'gray');
_images/05-Filtering_26_0.png

We obtain an image which of course is blurred, but where the strong noise has mostly disappeared.

Filtering an image

We have been toying here with singe pixel images, but of course the goal is to apply such a filter to an actual image. Let’s import one and adjust the width of our gaussian to see it’s effect on a real image:

#image_stack = skimage.io.imread('images/46658_784_B12_1.tif')
image_stack = skimage.io.imread('https://github.com/guiwitz/PyImageCourse_beginner/raw/master/images/46658_784_B12_1.tif')
image_nuclei = image_stack[1250:1750,300:800,2]
plt.subplots(figsize=(10,10))
plt.imshow(image_nuclei, cmap = 'gray');
_images/05-Filtering_31_0.png
filtered_image = skimage.filters.gaussian(image_nuclei)
plt.subplots(figsize=(10,10))
plt.imshow(filtered_image, cmap = 'gray');
_images/05-Filtering_33_0.png
filtered_image = skimage.filters.gaussian(image_nuclei, sigma = 5)
plt.subplots(figsize=(10,10))
plt.imshow(filtered_image, cmap = 'gray');
_images/05-Filtering_35_0.png

Image filtering is of course a vast topic, and you should read more about it if you need it in your work. This example just illustrates how to use filters in Python using scikit-image.

Non-linear filters

Another class of filters does not use convolution to filter the image but does some other type of operation on small neighborhoods of the image. For example it can calculate the median value of each 3x3 sub-region in the image. This is extremely useful if our image has pixels that are “dead” or have extremely large values.

For purpose of illustration, we can artificially create such an image using a skimage utility:

image_noisy = skimage.util.random_noise(image_nuclei, mode='s&p')
plt.subplots(figsize=(10,10))
plt.imshow(image_noisy, cmap = 'gray');
_images/05-Filtering_40_0.png

If we apply a median filter of size 3x3, pixels which are completely out of range are going to be replaced by the median value of their surroundings. The others, will be mildly affected:

filtered_median = skimage.filters.median(image_noisy)

plt.subplots(figsize=(10,10))
plt.imshow(filtered_median, cmap = 'gray');
_images/05-Filtering_42_0.png

This is a complete success!

Other selection elements

Until now we have relied on the default regions considered for filtering. For example the median filter operates by default on a 3x3 square regions. However, we can specify any other shape or footprint. For example we can use a Numpy array:

footprint = np.ones((5,5))
filtered_median = skimage.filters.median(image_noisy, footprint=footprint)
plt.imshow(filtered_median, cmap = 'gray');
_images/05-Filtering_45_0.png

Or we can use pre-made shapes that are available directly in scikit-image. For example a disk with a certain radius. This is actually what is done in the Demo notebook for the filtering step:

footprint = skimage.morphology.disk(3)
filtered_median = skimage.filters.median(image_noisy, footprint=footprint)
plt.imshow(filtered_median, cmap = 'gray');
_images/05-Filtering_47_0.png

Morphological filters

When working with masks, one can use another class of filters, that just turn pixels ON and OFF depending on their surroundings. The simplest to understand are the erosion and dilation filters. As their name indicate, these filters either shrink or enlarge ON objects. The size and shape of the filters, specify what region is considered around objects to shrink or enlarge them. We use here again the footprint argument to specify such structuring elements.

Let’s create a simple masks using an Otsu threshold:

mask = image_nuclei > skimage.filters.threshold_otsu(image_nuclei)

plt.imshow(mask, cmap = 'gray');
_images/05-Filtering_50_0.png

We can just specify the influence region as a Numpy array. For example a square of size 5:

influence_region = np.ones((5,5))
eroded = skimage.morphology.binary_erosion(mask, footprint=influence_region)
plt.imshow(eroded, cmap = 'gray');
_images/05-Filtering_54_0.png

We can erode even more:

influence_region = np.ones((10,10))

eroded = skimage.morphology.binary_erosion(mask, footprint=influence_region)

plt.imshow(eroded, cmap = 'gray');
_images/05-Filtering_56_0.png

We can do the reseverse an dilate our image:

influence_region = np.ones((10,10))

dilated = skimage.morphology.binary_dilation(mask, footprint=influence_region)

plt.imshow(dilated, cmap = 'gray');
_images/05-Filtering_58_0.png

How are these operators useful in image processing ? Let’s again imagine that we have a noisy image. Our mask would look like this:

mask = image_noisy > skimage.filters.threshold_otsu(image_noisy)

plt.imshow(mask, cmap = 'gray');
_images/05-Filtering_60_0.png

We can now fill the little holes in the nuclei with a small dilation:

influence_region = np.ones((2,2))

step1 = skimage.morphology.binary_dilation(mask, footprint=influence_region)

plt.imshow(step1, cmap = 'gray');
_images/05-Filtering_62_0.png

And erode again to put the nuclei back to their original size:

influence_region = np.ones((2,2))

step2 = skimage.morphology.binary_erosion(step1, footprint=influence_region)

plt.imshow(step2, cmap = 'gray');
_images/05-Filtering_64_0.png

We still have some noise. So we can do another round where we first erode the image (to remove dust) and dilate it again to recover nuclei of correct size:

influence_region = np.ones((2,2))

step3 = skimage.morphology.binary_erosion(step2, footprint=influence_region)
step4 = skimage.morphology.binary_dilation(step3, footprint=influence_region)

plt.imshow(step4, cmap = 'gray');
_images/05-Filtering_66_0.png

The combination of erosion and dilation is defined as specific operators (opening, closing) and are widely used. For example in the Demo notebook, we use the binary_closing filter which applies firts a dilation and then an erosion:

plt.imshow(skimage.morphology.binary_closing(mask, footprint=influence_region), cmap='gray');
_images/05-Filtering_68_0.png

Other morphological operators exist. For example one can shrink as much as possible all the structures of the image using the skeletonize function:

skeleton = skimage.morphology.skeletonize(step4)

plt.imshow(skeleton, cmap = 'gray');
_images/05-Filtering_70_0.png

Additionally, you can find many other tools for flood filling, removing holes, finding maxima etc. in this module. We leave it up to the reader to explore the documentation.

Other packages

You can of course find useful functions in other packages. An example of this is given in the Demo notebook, where we fill holes that we have in our nuclei with a function from the scipy package:

import scipy.ndimage as ndi
image_closed = skimage.morphology.binary_closing(mask, footprint=influence_region)
plt.imshow(image_closed, cmap='gray');
_images/05-Filtering_75_0.png
image_fill = ndi.binary_fill_holes(image_closed, skimage.morphology.disk(5))
plt.imshow(image_fill, cmap='gray');
_images/05-Filtering_77_0.png