Image noise reduction¶

In this tutorial, we estimate the original from a noisy image using Fixstars Amplify.

Formulation Overview¶

We will use black-and-white images here as a simple example. We can make the following assumptions about the relationship between the noisy and denoised images.

The noisy and denoised images have relatively large overlaps.
In the denoised image, neighboring pixels are likely to have the same color.

Noise reduction is performed by optimizing these conditions.

Constants and Variables¶

Since the color of a pixel is black or white, the value of each pixel can be expressed by a variable that takes value $-1$ or $1$. In the following, the black pixel corresponds to $-1$ and the white pixel to $1$.

Let $y_{ij}$ denote the value corresponding to the color of each pixel in the noisy image.

$$ y_{ij} = \begin{cases} -1 \quad \text{(if the corresponding pixel is black)} \\ +1 \quad \text{(if the corresponding pixel is white)} \end{cases} $$

Let $s_{ij}$ denote the Ising variable corresponding to the color of each pixel in the denoised image. The denoised image can be obtained by optimizing these Ising variables.

$$ s_{ij} = \begin{cases} -1 \quad \text{(if the corresponding pixel is black)} \\ +1 \quad \text{(if the corresponding pixel is white)} \end{cases} $$

Objective Function¶

The following two conditions about the denoised image are assumed:

The noisy and denoised images have relatively large overlaps.
In the denoised image, neighboring pixels are likely to have the same color.

Let us write down these conditions as a polynomial minimization problem.

First, consider the first assumption that the the noisy and denoised images have relatively large overlaps. $-y_{ij}s_{ij} = -1$ if the pixels at the same position (i, j) in the two images match and $-y{ij}s_{ij} = 1$ otherwise. Therefore, the following polynomial $f_1$, which is the summation of the expression by each pixel, takes smaller value if the two images have larger overlaps.

$$ f_1 = \sum_{i, j} -y_{ij} s_{ij} $$

Next, consider the second assumption that neighboring pixels tend to be the same color. $-s_{ij}s_{i'j'} = -1$ if the pixels at the position (i, j) and (i', j') in the denoised image have the same color and $-s_{ij}s_{i'j'} = 1$ otherwise. Therefore, the following polynomial $f_2$, which is the summation of the expression by each pair of neighboring pixels, takes smaller value if the larger number of neighboring pixels have the same color.

$$ f_2 = \sum_{s_{i, j} \text{and} s_{i', j'} \text{are adjacent}} -s_{i, j} s_{i', j'} $$

By adding these expressions $f_1$ and $f_2$ together with an appropriate weight, the objective function is created, which takes smaller value if the two conditions are satisfied.

$$ \begin{align} f & = f_1 + \eta f_2\\ f_1 &= \sum_{i, j} -y_{ij} s_{ij}, \\ f_2 &= \sum_{s_{i, j} \text{and} s_{i', j'} \text{are adjacent}} -s_{i, j} s_{i', j'} \end{align} $$

Here, we have introduced the parameter $\eta>0$. This allows us to adjust the relative strength of $f_1$ and $f_2$. The larger the value of $\eta$ is, the stronger the effect of the noise reduction term is.

By minimizing this objective function and interpreting the values of the Ising variables $s$ as the values of the pixels, we can obtain an image with reduced noise.

Reference¶

Load Noisy Image¶

We have prepared a noisy black-and-white image (noisy.png).

Our goal is to remove the noise from this image and estimate the original image.

First, we load the noisy image and convert it to a NumPy array.

In [ ]:

from PIL import Image
import numpy as np
    
img = Image.open("noisy.png")
noisy_img = np.where(np.array(img) >= 128, 1, -1)

The NumPy array noisy_img is a two-dimentional array of $81 \times 196$ , the same as the image size, each element of which is -1 (black) or 1 (white). This NumPy array can be visualized using Matplotlib as follows:

In [ ]:

%matplotlib inline

In [ ]:

import matplotlib.pyplot as plt
    
print(noisy_img.shape)
plt.imshow(noisy_img, cmap="gray")
plt.show()

Create the Ising Variable Array¶

Next, we create an array s of Ising variables. s is a two-dimensional array of the same shape as the noisy image, i.e. the same shape as the NumPy array noisy_img.

To create variables, use VariableGenerator. The array method of VariableGenerator creates variables in a multi-dimensional array form, which is useful for calculating th objective function.

In [ ]:

from amplify import VariableGenerator
    
gen = VariableGenerator()
# Create an Ising variable in the shape of the noisy image
s = gen.array("Ising", shape=noisy_img.shape)

Objective Function¶

We construct the objective function using the array of the noisy image data noisy_img and the Ising variable array s corresponding to the denoised image.

First, we create a function that expresses the condition that the the noisy and denoised images have relatively large overlaps. This function is expressed by $f_1 = \sum_{i, j} -y_{ij} s_{ij}$ and takes smaller value when the two images have larger overlaps.

In [ ]:

f1 = -(noisy_img * s).sum()

Next, create a function that expresses the condition that adjacent pixels tend to have the same color. This function is expressed by the following expression, and takes value when the large number of neighboring pixels have the same color.

$$ f_2 = \sum_{s_{i, j} \text{and} s_{i', j'} \text{are adjacent}} -s_{i, j} s_{i', j'} $$

In [ ]:

f2 = -((s[:, :-1] * s[:, 1:]).sum() + (s[:-1, :] * s[1:, :]).sum())

Here, the first term of the above expression f2 is the total sum of the products of the corresponding Ising variables $s_{i, j} s_{i', j'}$ for all pairs of adjacent pixels to the left and right. Similarly, the second term is the total sum of the products of the corresponding Ising variables for all pairs of upper and lower neighboring pixels.

The objective function is created by adding the two functions together with an appropriate weight. In this tutorial, we set the weight of the second term f2 for the first term f1 to 0.4. When the values of the variable array s are determined so that this objective function is as small as possible, the image represented by s is expected to be close to the image before the noise was added.

In [ ]:

eta = 0.4
objective = f1 + eta * f2

Set Up the Client and Run the Machine¶

We use Amplify AE as our combinatorial optimization solver. A solver client that corresponds to the solver (FixstarsClient) is created, and the parameters of the solver are set.

In [ ]:

from amplify import FixstarsClient
from datetime import timedelta
    
# Setting Fixstars Amplify AE as a client
client = FixstarsClient()
client.parameters.timeout = timedelta(milliseconds=2000)
# client.token = "xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx"  # If you use it in a local environment, please enter the access token for Amplify AE

Combinatorial optimization is executed using the created objective function and the solver client. When the following cell is executed, the Amplify AE will be run to minimize the objective function.

In [ ]:

from amplify import solve
    
result = solve(objective, client)

Obtaining the Solutions and Displaying the Results¶

The value of the variable array s in the best solution obtained can be obtained as follows:

In [ ]:

s_values = s.evaluate(result.best.values)

The s_values defined above is a two-dimensional NumPy array that has the same shape as s and represents the value of each pixel of the denoised image.

Finally, let us show the image represented by s_values. Compared to the first image (noisy.png), we can see that the noise has been reduced.

In [ ]:

plt.imshow(s_values, cmap="gray")
plt.show()