Lab 4 Image processing¶

0. Environment setup¶

1. Images¶

$$ \def\bbX{{{\mathbf X}}} $$

Images are mathematically modeled as functions of two variables. These variables represent the vertical coordinate $m$ and the horizontal coordinate $n$. These variables can be discrete or continuous but we will consider them here to be discrete. In this case the coordinate pair $(m, n)$ is called a pixel. We restrict each pixel coordinate to be between 0 and $N-1$. We use $\mathbf{X}$ to denote the image and and $x(m, n)$ to represent the value at pixel $(m, n)$. Pixels of an image are arranged in a matrix, $$ \tag{1} \mathbf{X}=\left(\begin{array}{cccccc} x(0,0) & x(0,1) & \cdots & x(0, n) & \cdots & x(0, N) \\ x(1,0) & x(1,1) & \cdots & x(1, n) & \cdots & x(1, N) \\ \vdots & \vdots & & \vdots & & \vdots \\ x(m, 0) & x(m, 1) & \cdots & x(m, n) & \cdots & x(m, N) \\ \vdots & \vdots & & \vdots & & \vdots \\ x(N, 0) & x(N, 1) & \cdots & x(N, n) & \cdots & x(N, N) \end{array}\right) $$ Thus, an image $\mathbf{X}$ is a matrix in which entries $x(m, n)$ represent pixel values.

The interpretation of this matrix is that pixel values $x(m,n)$ represent the luminance of the pixel. The luminance is how much light is reflected by the pixel. This determines how bright the pixel appears in a black and white image. For that reason we do not represent images as functions. We represent them as luminance colormaps.

To represent color images we just consider several matrices with different color contents. For instance, in a red, blue, and green (RGB) representation we have matrices $\bbX_{\text{R}}$, $\bbX_{\text{G}}$, and $ \bbX_{\text{B}}$ that represent the luminance restricted to the red, blue, and green colors. This information can be represented with a tensor,

$$ \tag{2} \bbX = [\bbX_{\text{R}}, \bbX_{\text{G}}, \bbX_{\text{B}}] $$

In this tensor, each of the color matrices is called a channel.

1.1 Handwritten Digits¶

In this lab we work with black and white images that represent handwritten digits. In Figure fig_images we show an image of a handwritten number 1 and a handwritten number 3. The dataset we are given contains pairs $(\bbX_q, y_q)$ of images $\bbX_q$ and human annotations $y_q$ that identify the correct digits.

Task 1¶

Load the data and visualize three images.

Loading the Data¶

We will use torchvision’s MNIST dataset class to download the data.

Visualize Images¶

2. Shifts in Space¶

Before defining convolutions to process images we have to introduce vertical and horizontal shift operators. We name these operators $\mathcal{S}_{\mathrm{v}}$ and $\mathcal{S}_{\mathrm{H}}$ and define them as the operators whose action on an image $\mathbf{X}$ results in a shifting of vertical and horizontal coordinates, respectively. Thus, if applying the vertical shift $\mathcal{S}_{\mathrm{v}}$ to the image $\mathbf{X}$ yields the image $\mathbf{U}_{10}=\mathcal{S}_{\mathrm{v}} \mathbf{X}$, the entries $u_{10}(m, n)$ of the vertically shifted image $\mathbf{U}_{10}$ are

$$ \tag{4} u_{10}(m, n)=x(m-1, n), \text { when } \mathbf{U}_{10}=\mathcal{S}_{\mathrm{v}} \mathbf{X} $$

Likewise if applying the vertical shift $\mathcal{S}_{\mathrm{H}}$ to the image $\mathbf{X}$ yields the image $\mathbf{U}_{01}=\mathcal{S}_{\mathrm{v}} \mathbf{X}$, the entries $u_{01}(m, n)$ of the horizontally shifted image $\mathbf{U}_{01}$ are

$$ u_{01}(m, n)=x(m, n-1), \quad \text { when } \mathbf{U}_{01}=\mathcal{S}_{\mathrm{H}} \mathbf{X} $$

In (4) and (5) we adopt the convention that $x(m,n)=0$ when either of the arguments $m$ or $n$ are negative.

To understand (4) and (5) it is convenient to represent them in matrix form. Applying the vertical shift to the image (1) yields the matrix

$$ \def\Sv{\mathcal{S}_{\mathrm{V}}} \def\Sh{\mathcal{S}_{\mathrm{H}}} $$

$$ \begin{equation}\tag{6} \Sv\bbX\!=\! \left( \begin{array} {} 0 & 0 & \cdots & 0 & \cdots & 0 \\ x(0,0) & x(0,1) & \cdots & x(0,n) & \cdots & x(0,N) \\ \vdots & \vdots & & \vdots & & \vdots \\ x(m,0) & x(m,1) & \cdots & x(m,n) & \cdots & x(m,N) \\ \vdots & \vdots & & \vdots & & \vdots \\ x(N-1,0) & x(N-1,1) & \cdots & x(N-1,n) & \cdots & x(N-1,N) \\ \end{array} \right) . \end{equation} $$

$$ \def\bbU{{{\mathbf U}}} \def\bbX{{{\mathbf X}}} $$

In the vertically shifted image $\bbU_{10}=\Sv\bbX$ all of the rows of the matrix $\bbX$ are shifted down. We fill the first row with zeros because there are no pixels that can be shifted into the first row.

When we apply the horizontal shift $\Sh$ to the image $\bbX$ in (1) we obtain the image $$ \def \Nmo {N\!-\!1} \begin{equation}\tag{7} \Sh\bbX \!=\! \left( \begin{array} {} 0 & x(0,0) & \cdots & x(0,n) & \cdots & x(0,\Nmo) \\ 0 & x(1,0) & \cdots & x(1,n) & \cdots & x(1,\Nmo) \\ \vdots & \vdots & & \vdots & & \vdots \\ 0 & x(m,0) & \cdots & x(m,n) & \cdots & x(m,\Nmo) \\ \vdots & \vdots & & \vdots & & \vdots \\ 0 & x(N,0) & \cdots & x(N,n) & \cdots & x(N,\Nmo) \\ \end{array} \right) . \end{equation} $$

In the horizontally shifted image $\bbU_{01}=\Sh\bbX$ all of the columns of the matrix $\bbX$ are shifted right. We fill the first column with zeros because there are no pixels that can be shifted into the first column.

2.1 Shift Compositions and Spatial Shift Sequences¶

As in the case of time signals, image shifts can be composed. For instance, applying the vertical shift twice yields the signal

$$ \begin{equation}\tag{8} \bbU_{20} = \Sv \,\bbU_{10} = \Sv^2 \,\bbX. \end{equation} $$

This is a signal in which the entries are $u_{20}(m, n) = x(m-2, n)$. That is, the image $\bbU_{20}$ is one in which the pixels are shifted two pixels down.

As another example consider the application of $l$ shifts in the horizontal direction. We can define define this recursively as

$$ \begin{equation}\tag{9} \bbU_{0l} = \Sh \,\bbU_{0(l-1)} = \Sh^2 \,\bbU_{0(l-2)} = \ldots = \Sh^{l-1} \,\bbU_{01} = \Sh^l \,\bbX. \end{equation} $$

This is a signal in which the entries are $u_{0l}(m, n) = x(m, n-l)$. That is, the image $\bbU_{0l}$ is one in which the pixels are shifted $l$ pixels to the right.

In general, we can compose any number of vertical shifts with any number of horizontal shifts. This results in the definition of the spatial shift sequence composed of images $$ \begin{equation}\tag{10} \bbU_{kl} = \Sv^k \Sh^l \bbX. \end{equation} $$ This is an image in which the entries are $u_{kl}(m, n) = x(m-k, n-l)$. That is, the image $\bbU_{kl}$ is one in which the pixels are shifted $k$ pixels down and $l$ pixels to the right.

It is important to note that the spatial shift sequence in (10) can be computed recursively. To do that we need to define the vertical and horizontal recursions

$$ \begin{equation}\tag{11} \bbU_{kl} = \Sv \bbU_{(k-1)l}, \qquad \bbU_{kl} = \Sh \bbU_{k(l-1)}, \end{equation} $$

with the initial condition $\bbU_{00} = \bbX$. Recursive computation of the spatial diffusion sequence is important in practical implementations.

2.2 Negative shifts¶

A negative vertical shift is a shift that moves the pixels up. We denote this shift by $\mathcal{S}_{\mathrm{v}}^{-1}$ and define it as the shift that when acting on the image $\mathbf{X}$ yields the image $\mathbf{U}_{(-1) 0}=\mathcal{S}_v^{-1} \mathbf{X}$ whose entries are given by $$ u_{(-1) 0}(m, n)=x(m+1, n), \quad \text { when } \mathbf{U}_{(-1) 0}=\mathcal{S}_{\mathrm{v}}^{-1} \mathbf{X} . $$ In this definition we adopt the convention that $x(N+1, n)=0$. This means that when applying the negative shift $\mathcal{S}_{\mathrm{v}}^{-1}$ we fill the last row of $\mathbf{U}_{(-1) 0}$ with zeros. This is because there are no entries of $\mathbf{X}$ that can be shifted into this row.

Likewise a negative horizontal shift is a shift that moves the pixels left. We denote this shift by $\mathcal{S}_{\mathrm{H}}^{-1}$ and define it as the shift that when acting on the image $\mathbf{X}$ yields the image $\mathbf{U}_{0(-1)}=\mathcal{S}_{\mathrm{H}}^{-1} \mathbf{X}$ whose entries are given by $$ u_{0(-1)}(m, n)=x(m, n+1) \text {, when } \mathbf{U}_{0(-1) 0}=\mathcal{S}_{\mathrm{H}}^{-1} \mathbf{X} . $$ In this definition we adopt the convention that $x(m, N+1)=0$. This means that when applying the negative shift $\mathcal{S}_{\mathrm{H}}^{-1}$ we fill the last column of $\mathbf{U}_{(-1) 0}$ with zeros. This is because there are no entries of $\mathbf{X}$ that can be shifted into this column.

As is the case of the positive shifts $\Sv$ and $\Sh$, the negative shifts $\Sv^{-1}$ and $\Sh^{-1}$ can be composed. It is also possible to compose the positive vertical shift $\Sv$ with the negative horizontal shift $\Sh^{-1}$. The converse composition of the positive horizontal shift $\Sh$ with the negative vertical shift $\Sv^{-1}$ is also possible. We can therefore generalize (10) into $$ \begin{equation}\tag{14} \bbU_{kl} = \Sv^k \Sh^l \bbX. \end{equation} $$ This is, in fact, the exact same equation as (10). The difference is that we now allow for $k$ and $l$ to be positive or negative numbers. The entries of $\bbU_{kl}$ are $ u_{kl}(m,n)=x(m-k,n-l)$. When $k$ is positive this entails shifting pixels down and when $k$ is negative this entails shifting pixels up. When $l$ is positive this entails shifting pixels to the right and when $k$ is negative this entails shifting pixels to the left.

3. Convolutions in space¶

A spatial convolution is a linear combination of the components of the spatial diffusion sequence in (14). For a formal definition we consider a filter range $K$ and define filter coefficients $h_{k l}$ for $k$ and $l$ ranging from $-K$ to $K$. The outcome of applying the filter with coefficients $h_{k l}$ to the image $\mathbf{X}$ is the image $$ \mathbf{Y}=\sum_{k=-K}^K \sum_{l=-K}^K h_{k l} \mathcal{S}_{\mathrm{v}}^k \mathcal{S}_{\mathrm{H}}^l \mathbf{X} . $$ Observe that in this definition we allow for positive and negative shifts. For simplicity we assume that the maximum number of shifts in either direction is the same. They can be made different, but it is standard practice to keep them equal.

We note that the total number of filter coefficients of a two dimensional convolution is $(2 K+1)^2$. These coefficients are arranged in a matrix $\mathbf{H}$.

Task 2¶

$$ \def\bbH{{{\mathbf H}}} $$

Write a filter object and endow it with a function that takes an image $\bbX$ as an input and returns as an output the result of convolving the image with the filter $\bbH$ that is stored as an attribute of the class. Make sure that this function can operate on a batch of input images $\blacksquare$

Task 3¶

Define filters $\bbH$ with entries $h_{kl} = 1/((2K+1)^2)$. Use the class in Task (2) to instantiate these filters with $K=2$, $K=4$, and $K=8$. Apply these filters to the images you visualized in Task (1).

What is the effect of applying these filters? $\blacksquare$

Task 4¶

Pytorch comes with its own function for computing two dimensional convolutions. Consider the same filters of Task 3 with $K=8$. Execute your own convolution code and the code that comes built in in Pytorch to compute convolutions for a batch of $B=100$ images. Compare the execution times. $\blacksquare$

In solving Task 4 you must have noticed that the Pytorch code is much faster. This is because the Pytorch convolution function is just a wrapper that calls a compiled subroutine written in C. If this were the 1990’s we would say that “true men code in C,” which was indeed a quip that was common in the 1990’s. Lucky for us all we live in a time when we do not require men to be tough, nor women to stay away from engineering, nor persons to define themselves as either men or women. So let us just remember that industry level data sciences still relies on low level languages. Python and Pytorch are prototyping languages.

We will use the Pytorch convolution function in the remainder of Lab 4.

4. Filterbanks¶

$$ \def\bbZ{{{\mathbf Z}}} \def\bbA{{{\mathbf A}}} \def\bbx{{{\mathbf x}}} $$

To use convolutional filters in classification tasks we will build convolutional neural networks (CNNs). An intermediate step is to use the energy content of the outputs of a convolutional filterbank.

Begin then by defining the energy of an image as the sum of the squares of the values of each pixel, $$ \begin{equation}\tag{16} \|\bbX\|^2 = \sum_{m=0}^{M-1} \sum_{n=0}^{M-1} x^2(m,n). \end{equation} $$

Continue by defining a collection of $G$ filters $\bbH^g$ each of which contains filter coefficients $h^g_{kl}$. Processing the input image $\bbX$ with each of these filters produces the output feature, $$ \begin{equation}\tag{17} \bbZ^g = \sum_{~k={-K}}^{K~} \sum_{~l={-K}~}^{K} h^g_{kl} \, \Sv^k \, \Sh^l \, \bbX . \end{equation} $$

Each of the $\bbZ^g$ features is an image. This set of images has to be converted into a set of class scores as we did in Lab 2C. We will do that in two steps. The first step is to create a vector that groups the energies of the outputs of each filterbank

$$ \begin{equation}\tag{18} \bbx_1 = \Big[~ \big\| \bbZ^1 \big\|^2;\, \big\| \bbZ^2 \big\|^2;\, \ldots;\, \big\| \bbZ^G \big\|^2 ~\Big]. \end{equation} $$ The second step is to process this vector with a readout layer. This is just a multiplication with a matrix $\bbA$ that matches the dimension $G$ of the filterbank with the number of classes $C$,

$$ \begin{equation}\tag{19} \bbx_2 = \ \bbA \bbx_1. \end{equation} $$ Each of the entries of $\bbx_2$ is a class score. We will train this filterbank to correctly predict class labels.

Task 5¶

Create a filterbank class. In this class the filter coefficients $H_g$ the matrix A are attributes. The class has a method that takes an image X as an input and implements the filterbank equations in (17) – (19). The method returns the vector of scores $x_2$ as an output.

Task 6¶

Split the handwritten digit dataset into train and test sets. Instantiate the class of Task 5 to train a filterbank that classifies images to the corresponding digit class. Use the crossentropy loss as the optimization objective.

Train¶

We already had train and test datasets from Pytorch’s MNIST loaders. All we need to do is take our usual training and evaluation code.