跳轉到

7.3 應用：影像壓縮與去噪#

本節目標#

  • 建立「影像就是矩陣」的觀念：灰階圖是 \(m \times n\) 矩陣、彩色圖是三層矩陣
  • 用截斷 SVD 壓縮真實圖片，計算壓縮比與重建品質
  • 處理彩色圖：三個通道各自做 SVD
  • 用「保留大奇異值、丟棄雜訊方向」的原理做 SVD 去噪
  • 完成實驗：找出你的圖片需要幾個奇異值才能辨識

影像就是矩陣#

一張灰階圖，就是一個元素為亮度的矩陣：每個像素一個數字。 用 plt.imread 或本書提供的工具箱讀入，就能直接「看見矩陣」：

import numpy as np
import matplotlib.pyplot as plt

# 真實圖片一：scipy 內建的浣熊照片（512 x 512 灰階）
from scipy import datasets
img = datasets.face(gray=True).astype(float) / 255.0   # 標準化到 [0, 1]

# 真實圖片二：自己的圖檔（任何 jpg/png 都行）
# img = plt.imread("my_photo.jpg")
# if img.ndim == 3: img = img[..., 0]        # 先取單通道
# img = img.astype(float) / img.max()

print(img.shape, img.min(), img.max())

fig, axes = plt.subplots(1, 2, figsize=(9, 4))
axes[0].imshow(img, cmap="gray")
axes[0].set_title("原圖")
axes[1].imshow(img[:8, :8], cmap="gray")     # 左上角 8x8 像素放大
axes[1].set_title("放大看：就是一格格數字")
plt.show()

print(np.round(img[:8, :8], 2))              # 這 8x8 就是矩陣本體

矩陣的「秩」在影像語言中就是「影像的結構複雜度」。 平滑的照片秩低、雜亂的雜訊秩高——這正是壓縮與去噪的切入點。

截斷 SVD 影像壓縮#

對影像矩陣 \(A\) 做 SVD，保留前 \(k\) 個奇異值，就得到壓縮版 \(A_k\) 。 儲存 \(A_k\) 只需 \(k(m + n + 1)\) 個數字，壓縮比為

\[ \text{壓縮比} = \frac{k(m + n + 1)}{mn} \]

完整可執行的壓縮實驗：

import numpy as np
import matplotlib.pyplot as plt
from scipy import datasets

img = datasets.face(gray=True).astype(float) / 255.0
m, n = img.shape

U, s, Vt = np.linalg.svd(img, full_matrices=False)

fig, axes = plt.subplots(2, 3, figsize=(13, 8))
for ax, k in zip(axes.ravel(), [5, 10, 20, 50, 100, 200]):
    Ak = (U[:, :k] * s[:k]) @ Vt[:k, :]
    ratio = k * (m + n + 1) / (m * n)
    mse = np.mean((img - Ak)**2)
    ax.imshow(Ak, cmap="gray")
    ax.set_title(f"$k={k}$  壓縮比 {ratio:.2%}  MSE {mse:.2e}")
    ax.axis("off")
fig.suptitle("截斷 SVD 影像壓縮（浣熊照片 512×512）")
plt.tight_layout(); plt.show()

再看奇異值本身——影像的衰減曲線非常陡：

plt.plot(s / s[0], "o-", ms=3)
plt.xlabel("索引 $i$"); plt.ylabel(r"$\sigma_i / \sigma_1$")
plt.title("影像的奇異值衰減：少數方向承載大部分能量")
plt.grid(alpha=0.3); plt.show()

cum = np.cumsum(s**2) / (s**2).sum()
for p in [0.90, 0.95, 0.99]:
    k = np.searchsorted(cum, p) + 1
    print(f"保留 {p:.0%} 能量需要 k = {k}")

典型的結果： \(k \approx 50\) （壓縮比約 20%）時，照片已清晰可辨； 這就是奇異值衰減的直接兌現。

與 JPEG 的關係：JPEG 用的是類似思想的另一套基底（離散餘弦）， 而 SVD 的基底是從這張圖自己算出來的——資料自適應，代價是計算較貴。

彩色圖像：三通道各自壓縮#

彩色圖是 \(m \times n \times 3\) 的陣列。SVD 只吃二維矩陣， 所以標準做法是對 R、G、B 三個通道分別做截斷 SVD，再合回來：

import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import load_sample_image

img = load_sample_image("china.jpg")           # 真實照片：(427, 640, 3)
img = img.astype(float) / 255.0
m, n, _ = img.shape

k = 50
channels = []
for c in range(3):                             # R, G, B 各做一次
    U, s, Vt = np.linalg.svd(img[:, :, c], full_matrices=False)
    channels.append((U[:, :k] * s[:k]) @ Vt[:k, :])
img_k = np.clip(np.stack(channels, axis=-1), 0, 1)   # 合回三通道並夾到 [0,1]

fig, axes = plt.subplots(1, 2, figsize=(11, 4))
axes[0].imshow(img);            axes[0].set_title("原圖")
axes[1].imshow(img_k)
axes[1].set_title(f"三通道各取 $k={k}$（壓縮比 {3*k*(m+n+1)/(m*n*3):.2%}）")
for ax in axes: ax.axis("off")
plt.tight_layout(); plt.show()

注意 np.clip：截斷重建可能讓像素值略微超出 \([0, 1]\) ，顯示前要夾回來。

SVD 去噪：雜訊住在小奇異值裡#

為什麼丟掉小奇異值能「去掉雜訊」？直覺是：

  • 影像結構（邊緣、大塊色面）能量集中在前幾個奇異方向
  • 雜訊在所有方向上均勻分布，能量攤在數百個小奇異值裡

所以 \(A_k\) 保留了結構、丟掉了雜訊。做個對照實驗：

import numpy as np
import matplotlib.pyplot as plt
from scipy import datasets

rng = np.random.default_rng(0)
img = datasets.face(gray=True).astype(float) / 255.0
noisy = img + 0.1 * rng.standard_normal(img.shape)    # 加入高斯雜訊
noisy = np.clip(noisy, 0, 1)

U, s, Vt = np.linalg.svd(noisy, full_matrices=False)

fig, axes = plt.subplots(1, 3, figsize=(13, 4.5))
for ax, k, title in zip(axes, [20, 60, 200],
                        ["$k=20$：太狠，細節也沒了", "$k=60$：去噪 sweet spot",
                         "$k=200$：太鬆，雜訊又回來"]):
    Ak = (U[:, :k] * s[:k]) @ Vt[:k, :]
    ax.imshow(Ak, cmap="gray")
    ax.set_title(f"{title}，MSE {np.mean((img - Ak)**2):.2e}")
    ax.axis("off")
fig.suptitle("SVD 去噪：對含雜訊圖做截斷 SVD")
plt.tight_layout(); plt.show()

對照組：乾淨圖的 \(A_k\) 與含雜訊圖的 \(A_k\) （同樣的 \(k\) ）——

U0, s0, Vt0 = np.linalg.svd(img, full_matrices=False)
plt.plot(np.cumsum(s0**2)/ (s0**2).sum(), label="乾淨圖")
plt.plot(np.cumsum(s**2)  / (s**2).sum(),   label="含雜訊圖")
plt.xlabel("$k$"); plt.ylabel("累積能量比例"); plt.legend(); plt.show()

含雜訊圖要累積同樣能量需要更大的 \(k\) ，因為雜訊「墊高」了後段奇異值—— 這張圖讓「雜訊住在小奇異值裡」變成看得見的事實。

實驗：你的圖片需要幾個奇異值？#

把整套流程包成函數，換任何一張你自己的照片都能跑：

import numpy as np
import matplotlib.pyplot as plt

def svd_compress_experiment(path_or_array, ks=(5, 20, 50, 100)):
    img = plt.imread(path_or_array) if isinstance(path_or_array, str) \
          else np.asarray(path_or_array, dtype=float)
    if img.max() > 1.5:
        img = img / 255.0
    if img.ndim == 3:
        img = img.mean(axis=2)                  # 轉灰階
    m, n = img.shape
    U, s, Vt = np.linalg.svd(img, full_matrices=False)
    fig, axes = plt.subplots(1, len(ks) + 1, figsize=(3*(len(ks)+1), 3.5))
    axes[0].imshow(img, cmap="gray"); axes[0].set_title("原圖")
    for ax, k in zip(axes[1:], ks):
        Ak = (U[:, :k] * s[:k]) @ Vt[:k, :]
        ax.imshow(np.clip(Ak, 0, 1), cmap="gray")
        ax.set_title(f"$k={k}$（{k*(m+n+1)/(m*n):.1%}）")
    for ax in axes: ax.axis("off")
    plt.show()

# svd_compress_experiment("my_photo.jpg")

操作方式：放進你最喜歡的一張照片，從 \(k = 5\) 開始往上加， 記下「第一次能認出畫面內容」的 \(k\) 值。不同類型的圖差異很大：

影像類型 可辨識所需的 \(k\) 原因
天空、漸層背景 很小（ \(< 10\) ） 結構極簡單，秩很低
人臉、建築 中等（數十） 有中等的邊緣結構
樹葉、草地、雜訊 很大 高頻細節多，秩高

小結與練習#

本節重點

  • 灰階影像 \(=\) 矩陣、彩色影像 \(=\) 三個矩陣；影像越「平滑」秩越低
  • 截斷 SVD 壓縮：儲存 \(A_k\) 需 \(k(m+n+1)\) 個數字，壓縮比 \(= k(m+n+1)/mn\)
  • 彩色圖對三通道各自做 SVD 後合併，記得 np.clip 回 \([0, 1]\)
  • 去噪原理：結構能量集中在前面奇異值、雜訊攤在後面； \(k\) 是「去噪力度」的旋鈕
  • \(k\) 太小丟結構、太大留雜訊——實驗找出 sweet spot 是唯一正解

練習

  1. 對浣熊照片畫「MSE 隨 \(k\) 」曲線，標出你肉眼認定的 sweet spot，兩者一致嗎？
  2. 修改彩色壓縮程式：讓三個通道使用不同的 \(k\) （例如 R 用 80、G 用 50、B 用 30），觀察品質差異並解釋（提示：人眼對綠色最敏感）。
  3. 把雜訊強度從 \(0.05\) 加到 \(0.30\) ，對每個強度掃 \(k\) ，找出最佳 \(k\) 與雜訊強度的關係。
  4. 比較「SVD 壓縮比 20%」與「把圖存成 JPEG 品質 20」的檔案大小與視覺品質，寫下你的觀察。