非負画像の畳み込み演算における負の値のリスクとその対策

## はじめに
畳み込み演算は画像処理においてよく使われる計算過程であり、エッジ検出、ぼかし、フィルタリングなど、さまざまな画像処理手法で広く使用されています。しかし、非負の画像を畳み込み演算する場合でも、結果として負の値が含まれてしまうことがあります。この問題は、特に高速フーリエ変換（FFT）を用いた畳み込みの効率化手法において発生しやすいです。

高速フーリエ変換は、畳み込み演算を効率的に行うための強力なツールですが、計算精度の問題や数値誤差が原因で、非負の入力画像に対しても負の値が含まれる結果を生成することがあります。このような負の値は、数値計算において重大な問題を引き起こす危険性があります。例えば、後続の処理である画像の正規化や閾値処理において予期せぬ動作を引き起こす危険性があります。

本記事では、畳み込み演算によって負の値が生じるケースを示し、その対処法について説明します。

本記事に用いた各種バージョンは、numpy 1.25.2、scipy 1.11.4となります。また、実装コードは、Google Colab[こちら](https://colab.research.google.com/drive/1yvHPyFAHIA6zXgVKRd8Yzh3m6WtIaAAN?usp=sharing)からもご覧いただけますので、よろしければご参照ください。

## 負の値が生じるケースの一例
さっそく、負の値が出るケースの例を見ていきましょう。非負の画像の典型例としてポアソン分布に従う擬似乱数画像を用意します。カーネル側も非負として畳み込み演算を行います。そして、畳み込んだ画像の負の値のピクセルを判定して描画して確認します。


```python
import numpy as np
from scipy.signal import fftconvolve
import matplotlib.pyplot as plt

# ポアソン分布に従う擬似的な乱数の画像を生成
np.random.seed(2024)

lam = 0.5  # ポアソン分布のλパラメータ
image = np.random.poisson(lam, size=(64, 64)).astype('int32')


# カーネルを生成
kernel = np.ones((3, 3), dtype='float64')

# 畳み込み
convolved_image = fftconvolve(image, kernel, mode='same')

# 負の値を持つピクセルの座標を取得
negative_value_coords = np.argwhere(convolved_image < 0)
# xとyの座標の反転（matplotlibに描画のため）
flipped_negative_value_coords = negative_value_coords[:, [1, 0]]

# 画像を描画
plt.figure(figsize=(17, 4))

# 元の画像を描画
plt.subplot(1, 3, 1)
plt.title("Original Image")
plt.imshow(image, cmap='gray')
plt.colorbar()

# 畳み込んだ画像を描画
plt.subplot(1, 3, 2)
plt.title("Convolved Image")
plt.imshow(convolved_image, cmap='gray')
plt.colorbar()

# 負の値を持つピクセルの座標をプロット
plt.subplot(1, 3, 3)
plt.title("Negative Value Coordinates")
plt.imshow(convolved_image, cmap='gray')
plt.colorbar()
if flipped_negative_value_coords.size > 0:
    plt.scatter(flipped_negative_value_coords[:, 0], flipped_negative_value_coords[:, 1], color='red', s=10)

plt.show()
```
結果
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1599258/c1ceea57-80d0-7ee0-1dca-7c65d99c22c8.png)

このコードでは、ポアソン分布に従う擬似乱数画像を生成し、3x3の全て1のカーネルで畳み込み演算を行います。
畳み込み演算には[`scipy.signal.fftconvolve()`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.fftconvolve.html)を使用しています。

結果の画像は、左から元画像、畳み込み後の画像、畳み込み後の画像に負の値を持つ画素を赤でプロットしています。この結果から負の値が数個ほど含まれていることがわかります。

## 対処法
畳み込み演算によって生じる負の値を避けるための基本的な対処法を2つ紹介します。1つ目は、単純に畳み込み後の負の値をゼロに置換する方法です。この方法は実装が簡単であり、負の値を排除する点において間違いが少ないですが、情報の損失が発生することに注意が必要です。2つ目は、高速フーリエ変換を使用せず、直接的な畳み込み計算を行う方法です。この方法は計算コストが増加するものの、負の値が生成されるリスクを低減できます。以下に、それぞれの方法の詳細と実装例を説明します。

### 1. 非負の値をゼロに置換する方法
畳み込み演算後の画像の負の値をゼロに置換する方法があります。この方法は結果に対して作用させるため、負の値を単純に排除できますが、情報を損失することに注意が必要です。

追記する箇所は、畳み込み後の負の値のゼロへの置換のコードのみです。

```python
# 負の値をゼロに置換
convolved_image[convolved_image < 0] = 0
```

<details><summary>この処理のコード全容</summary>

```python
## 畳み込み後の画像の負の値をゼロに置換して負の値を回避
import numpy as np
from scipy.signal import fftconvolve
import matplotlib.pyplot as plt

# ポアソン分布に従う擬似的な乱数の画像を生成
np.random.seed(2024)

lam = 0.5  # ポアソン分布のλパラメータ
image = np.random.poisson(lam, size=(64, 64)).astype('int32')

# カーネルを生成
kernel = np.ones((3, 3), dtype='float64')

# 畳み込み
convolved_image = fftconvolve(image, kernel, mode='same')
convolved_image = fftconvolve(image, kernel, mode='same')

# 負の値をゼロに置換
convolved_image[convolved_image < 0] = 0

# 負の値を持つピクセルの座標を取得
negative_value_coords = np.argwhere(convolved_image < 0)
# xとyの座標の反転（matplotlibに描画のため）
flipped_negative_value_coords = negative_value_coords[:, [1, 0]]

# 画像を描画
plt.figure(figsize=(17, 4))

# 元の画像を描画
plt.subplot(1, 3, 1)
plt.title("Original Image")
plt.imshow(image, cmap='gray')
plt.colorbar()

# 畳み込んだ画像を描画
plt.subplot(1, 3, 2)
plt.title("Convolved Image")
plt.imshow(convolved_image, cmap='gray')
plt.colorbar()

# 負の値を持つピクセルの座標をプロット
plt.subplot(1, 3, 3)
plt.title("Negative Value Coordinates")
plt.imshow(convolved_image, cmap='gray')
plt.colorbar()
if flipped_negative_value_coords.size > 0:
    plt.scatter(flipped_negative_value_coords[:, 0], flipped_negative_value_coords[:, 1], color='red', s=10)

plt.show()
```

</details>

### 2. 高速フーリエ変換を使わずに畳み込み演算を実直に行う方法
高速フーリエ変換（FFT）を使用しないで畳み込み演算を行う方法です。これにより、負の値が生じるリスクを低減できますが、画像サイズが大きい場合に計算コストの増加が懸念されます。高速フーリエ変換の数学的な説明や計算コストについては、[こちら](https://qiita.com/51_24_11_/items/11c48395603670ea98d8)を参考してください。

ここでは、実直な計算に対応している[`scipy.signal.convolve()`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.convolve.html)の`method=direct`による畳み込み演算するコードを示し、結果を確認します。

```python
## scipy.signal.convolve(method='direct')を使用して負の値を回避
import numpy as np
from scipy.signal import convolve
import matplotlib.pyplot as plt

# ポアソン分布に従う擬似的な乱数の画像を生成
np.random.seed(2024)

lam = 0.5  # ポアソン分布のλパラメータ
image = np.random.poisson(lam, size=(64, 64)).astype('int32')

# カーネルを生成
kernel = np.ones((3, 3), dtype='float64')

# 畳み込み
convolved_image = convolve(image, kernel, mode='same', method='direct')

# 負の値を持つピクセルの座標を取得
negative_value_coords = np.argwhere(convolved_image < 0)
# xとyの座標の反転（matplotlibに描画のため）
flipped_negative_value_coords = negative_value_coords[:, [1, 0]]

# 画像を描画
plt.figure(figsize=(17, 4))

# 元の画像を描画
plt.subplot(1, 3, 1)
plt.title("Original Image")
plt.imshow(image, cmap='gray')
plt.colorbar()

# 畳み込んだ画像を描画
plt.subplot(1, 3, 2)
plt.title("Convolved Image")
plt.imshow(convolved_image, cmap='gray')
plt.colorbar()

# 負の値を持つピクセルの座標をプロット
plt.subplot(1, 3, 3)
plt.title("Negative Value Coordinates")
plt.imshow(convolved_image, cmap='gray')
plt.colorbar()
if flipped_negative_value_coords.size > 0:
    plt.scatter(flipped_negative_value_coords[:, 0], flipped_negative_value_coords[:, 1], color='red', s=10)

plt.show()
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1599258/413cb68c-9619-b094-9333-dc23e96edd6d.png)

このように、負の値（赤色のプロット）が含まれていないことがわかります。なお、`scipy.signal.convolve()`の`method`引数はデフォルトでは`method='auto'`となっており、画像のサイズなどにより高速フーリエ変換（`method='fft'`）が使われることになります。本コードでは`method='direct'`と指定し、実直な計算が明示的に行われるようにしています。

## まとめ
畳み込み演算は強力なツールであり、高速化などの取り組みが行われていますが、非負の画像を扱う際には負の値が生じるリスクに注意が必要です。本記事では、そのリスクを示し、簡単な対処法について紹介しました。課題に応じた適切な方法を選択し、正確な画像処理を行っていきましょう。

## 参考

https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.fftconvolve.html

https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.convolve.html

https://qiita.com/51_24_11_/items/11c48395603670ea98d8