NumPy之:多維數組中的線性代數

目錄

• 簡介
• 圖形加載和說明
• 圖形的灰度
• 灰度圖像的壓縮
• 原始圖像的壓縮
• 總結

簡介

本文將會以圖表的形式爲你們講解怎麼在NumPy中進行多維數據的線性代數運算。

多維數據的線性代數一般被用在圖像處理的圖形變換中，本文將會使用一個圖像的例子進行說明。

圖形加載和說明

熟悉顏色的朋友應該都知道，一個顏色能夠用R，G，B來表示，若是更高級一點，那麼還有一個A表示透明度。一般咱們用一個四個屬性的數組來表示。

對於一個二維的圖像來講，其分辨率能夠看作是一個X*Y的矩陣，矩陣中的每一個點的顏色均可以用（R，G，B）來表示。

有了上面的知識，咱們就能夠對圖像的顏色進行分解了。

首先須要加載一個圖像，咱們使用imageio.imread方法來加載一個本地圖像，以下所示：

import imageio
img=imageio.imread('img.png')
print(type(img))

上面的代碼從本地讀取圖片到img對象中，使用type能夠查看img的類型，從運行結果，咱們能夠看到img的類型是一個數組。

class 'imageio.core.util.Array'

經過img.shape能夠獲得img是一個(80, 170, 4)的三維數組，也就是說這個圖像的分辨率是80*170，每一個像素是一個（R，B，G，A）的數組。

最後將圖像畫出來以下所示：

import matplotlib.pyplot as plt
plt.imshow(img)

圖形的灰度

對於三維數組來講，咱們能夠分別獲得三種顏色的數組以下所示：

red_array = img_array[:, :, 0]
green_array = img_array[:, :, 1]
blue_array = img_array[:, :, 2]

有了三個顏色以後咱們可使用下面的公式對其進行灰度變換：

Y=0.2126R + 0.7152G + 0.0722B

上圖中Y表示的是灰度。

怎麼使用矩陣的乘法呢？使用 @ 就能夠了：

img_gray = img_array @ [0.2126, 0.7152, 0.0722]

如今img是一個80 * 170的矩陣。

如今使用cmap="gray"做圖：

plt.imshow(img_gray, cmap="gray")

能夠獲得下面的灰度圖像：

灰度圖像的壓縮

灰度圖像是對圖像的顏色進行變換，若是要對圖像進行壓縮該怎麼處理呢？

矩陣運算中有一個概念叫作奇異值和特徵值。

設A爲n階矩陣，若存在常數λ及n維非零向量x，使得Ax=λx，則稱λ是矩陣A的特徵值，x是A屬於特徵值λ的特徵向量。

一個矩陣的一組特徵向量是一組正交向量。

即特徵向量被施以線性變換 A 只會使向量伸長或縮短而其方向不被改變。

特徵分解（Eigendecomposition），又稱譜分解（Spectral decomposition）是將矩陣分解爲由其特徵值和特徵向量表示的矩陣之積的方法。

假如A是m * n階矩陣，q=min(m,n)，A*A的q個非負特徵值的算術平方根叫做A的奇異值。

特徵值分解能夠方便的提取矩陣的特徵，可是前提是這個矩陣是一個方陣。若是是非方陣的狀況下，就須要用到奇異值分解了。先看下奇異值分解的定義：

\(A=UΣV^T\)

其中A是目標要分解的m * n的矩陣，U是一個 m * m的方陣，Σ 是一個m * n 的矩陣，其非對角線上的元素都是0。\(V^T\)是V的轉置，也是一個n * n的矩陣。

奇異值跟特徵值相似，在矩陣Σ中也是從大到小排列，並且奇異值的減小特別的快，在不少狀況下，前10%甚至1%的奇異值的和就佔了所有的奇異值之和的99%以上了。也就是說，咱們也能夠用前r大的奇異值來近似描述矩陣。r是一個遠小於m、n的數，這樣就能夠進行壓縮矩陣。

經過奇異值分解，咱們能夠經過更加少許的數據來近似替代原矩陣。

要想使用奇異值分解svd能夠直接調用linalg.svd 以下所示：

U, s, Vt = linalg.svd(img_gray)

其中U是一個m * m矩陣，Vt是一個n * n矩陣。

在上述的圖像中，U是一個(80, 80)的矩陣，而Vt是一個(170, 170) 的矩陣。而s是一個80的數組，s包含了img中的奇異值。

若是將s用圖像來表示，咱們能夠看到大部分的奇異值都集中在前的部分：

這也就意味着，咱們能夠取s中前面的部分值來進行圖像的重構。

使用s對圖像進行重構，須要將s還原成80 * 170 的矩陣：

# 重建
import numpy as np
Sigma = np.zeros((80, 170))
for i in range(80):
    Sigma[i, i] = s[i]

使用 U @ Sigma @ Vt 便可重建原來的矩陣，能夠經過計算linalg.norm來比較一下原矩陣和重建的矩陣之間的差別。

linalg.norm(img_gray - U @ Sigma @ Vt)

或者使用np.allclose來比較兩個矩陣的不一樣：

np.allclose(img_gray, U @ Sigma @ Vt)

或者只取s數組的前10個元素，進行從新繪圖，比較一下和原圖的區別：

k = 10
approx = U @ Sigma[:, :k] @ Vt[:k, :]
plt.imshow(approx, cmap="gray")

能夠看到，差別並非很大：

原始圖像的壓縮

上一節咱們講到了如何進行灰度圖像的壓縮，那麼如何對原始圖像進行壓縮呢？

一樣可使用linalg.svd對矩陣進行分解。

可是在使用前須要進行一些處理，由於原始圖像的img_array 是一個(80, 170, 3)的矩陣--這裏咱們將透明度去掉了，只保留了R，B，G三個屬性。

在進行轉換以前，咱們須要把不須要變換的軸放到最前面，也就是說將index=2，換到index=0的位置，而後進行svd操做：

img_array_transposed = np.transpose(img_array, (2, 0, 1))
print(img_array_transposed.shape)

U, s, Vt = linalg.svd(img_array_transposed)
print(U.shape, s.shape, Vt.shape)

一樣的，如今s是一個(3, 80)的矩陣，仍是少了一維，若是重建圖像，須要將其進行填充和處理，最後將重建的圖像輸出：

Sigma = np.zeros((3, 80, 170))

for j in range(3):
    np.fill_diagonal(Sigma[j, :, :], s[j, :])

reconstructed = U @ Sigma @ Vt
print(reconstructed.shape)

plt.imshow(np.transpose(reconstructed, (1, 2, 0)))

固然，也能夠選擇前面的K個特徵值對圖像進行壓縮：

approx_img = U @ Sigma[..., :k] @ Vt[..., :k, :]
print(approx_img.shape)
plt.imshow(np.transpose(approx_img, (1, 2, 0)))

從新構建的圖像以下：

對比能夠發現，雖然損失了部分精度，可是圖像仍是能夠分辨的。

總結

圖像的變化會涉及到不少線性運算，你們能夠以此文爲例，仔細研究。

本文已收錄於 http://www.flydean.com/08-python-numpy-linear-algebra/

最通俗的解讀，最深入的乾貨，最簡潔的教程，衆多你不知道的小技巧等你來發現！

歡迎關注個人公衆號:「程序那些事」,懂技術，更懂你！