图像降噪算法——稀疏表达:K-SVD算法
圖像降噪算法——稀疏表達:K-SVD算法
- 圖像降噪算法——稀疏表達:K-SVD算法
- 1. 基本原理
- 2. python代碼
- 3. 結論
圖像降噪算法——稀疏表達:K-SVD算法
為了完善下自己降噪算法知識版圖的完整性,我打算花一個周末的時間再了解下基于稀疏表達和低秩聚類這兩種原理實現的圖像降噪算法,因為學習的時間并不長,也沒有花太多時間去做實驗,所以對算法理解得可能比較膚淺,還愿讀者見諒。
這里我分享幾篇很優秀的博客,我對稀疏表達概念的認識目前也就是來自于這幾篇博客:
字典學習(Dictionary Learning, KSVD)詳解
從稀疏表示到K-SVD,再到圖像去噪
稀疏表達的意義在于?為什么稀疏表達得到廣泛的應用?
1. 基本原理
這里我們先來回答幾個問題:
(1)什么是稀疏表達?
稀疏表達的定義是用較少的基本信號的線性組合來表達大部分或者全部的原始信號。在稀疏表達算法中有"字典"這樣一個概念,這個“字典”和我們現實生活中使用的字典作用是相同的,現實生活中,我們通過對詞條進行組合可以表達我們想要的信息,同樣地,在稀疏表達算法中,我們對“原子”進行線性組合就可以還原我們想要的信號。稀疏表達在機器學習領域有著廣泛的應用,諸如壓縮感知、特征提取等,圖像去噪只是其應用的一個方面。
(2)稀疏表達為什么可以進行圖像去噪?
噪聲圖像是原始圖像和噪聲合成的圖像,原始圖像認為是可稀疏的,即可以通過有限個"原子"進行表示,而噪聲是隨機而不可稀疏的,因此通過提取圖像的稀疏成分,再利用稀疏成分來重構圖像,在這個過程中,噪聲被當做觀測圖像和重構圖像之間的殘差而被丟棄,從而起到降噪的作用。
在圖像去噪領域,盡管CNN的方法已經占據了半壁江山,但是CNN方法的一個弊端就是需要大量訓練數據,在難以采集大量訓練數據的場景中,例如某些醫學圖像,稀疏表達的方法仍然起到重要作用。
(3)稀疏表達和轉換域濾波方法的區別?
轉換域濾波指的是離散傅里葉變換、小波變換這樣一些方法,這些方法獲得的是一系列正交基,而稀疏表達獲得的“字典”則是一些列非正交基。正交基往往只能表示圖像的某一個特征而不能夠同時表示其他特征,因此正交基的稀疏性不及非正交基。
下面咱來開始討論基于稀疏濾波的K-SVD算法,K-SVD算法中最核心的部分就是字典學習,字典學習的過程就是基于樣本集Y={yi}i=1M\mathbf{Y}=\{\mathbf{y}_i\}^M_{i=1}Y={yi?}i=1M?尋找一組“超完備”的基向量作為字典矩陣D={di}i=1K\mathbf{D}=\{\mathbfze8trgl8bvbq_i\}^{K}_{i=1}D={di?}i=1K?其中MMM為樣本數,yi∈RN\mathbf{y}_{i} \in R^{N}yi?∈RN為單個樣本,是一個NNN維的特征向量,MMM個單個樣本按列組成樣本集矩陣,di∈RN\mathbfze8trgl8bvbq_i \in R^{N}di?∈RN為基向量,維度也是NNN維,KKK個基向量按列組成字典矩陣,樣本集中的任意樣本都可以根據字典求得起對應的稀疏表達形式。當K>NK > NK>N時為超完備字典,即我們字典學習的通常情況。
字典學習的目標是學習一個字典矩陣D\mathbf{D}D,使得Y≈D?X\mathbf{Y} \approx \mathbf{D} * \mathbf{X} Y≈D?X同時要XXX盡可能稀疏。這里我們注意D∈RN×K\mathbf{D} \in \mathbf{R}^{N \times K}D∈RN×K, Y∈RN×M\mathbf{Y} \in \mathbf{R}^{N \times M}Y∈RN×M,因此X∈RK×M\mathbf{X} \in \mathbf{R}^{K \times M}X∈RK×M。如下圖所示
D\mathbf{D}D矩陣中不同顏色代表不同的基向量,X\mathbf{X}X矩陣中不同的顏色代表對應基向量的系數,Y\mathbf{Y}Y矩陣中某一樣本(灰色塊)就是由D\mathbf{D}D矩陣的基向量通過X\mathbf{X}X矩陣中黑框中的系數進行線性組合獲得的,由于是稀疏組合,所以中間有很多系數會是零,我們用白色塊表示。
這個問題用數學問題描述為min?D,X∥Y?DX∥F2,s.t.??i,∥xi∥0≤T0\min _{\mathbf{D}, \mathbf{X}}\|\mathbf{Y}-\mathbf{D} \mathbf{X}\|_{F}^{2}, \quad \text { s.t. } \forall i,\left\|\mathbf{x}_{i}\right\|_{0} \leq T_{0} D,Xmin?∥Y?DX∥F2?,?s.t.??i,∥xi?∥0?≤T0?或者min?D,x∑i∥xi∥0,s.t.?min?D,x∥Y?DX∥F2≤?\min _{\mathbf{D}, \mathbf{x}} \sum_{i}\left\|\mathbf{x}_{i}\right\|_{0}, \quad \text { s.t. } \min _{\mathbf{D}, \mathbf{x}}\|\mathbf{Y}-\mathbf{D} \mathbf{X}\|_{F}^{2} \leq \epsilon D,xmin?i∑?∥xi?∥0?,?s.t.?D,xmin?∥Y?DX∥F2?≤?其中xi\mathbf{x}_ixi?為稀疏矩陣X\mathbf{X}X的行向量,代表字典矩陣的系數,和上面彩圖中是對應的。這里值得注意的是∥xi∥0\left\|\mathbf{x}_{i}\right\|_{0}∥xi?∥0?是零階范數,它表示向量中不為0的數的個數。
上面兩個公式是帶有約束的優化問題,可以利用拉格朗日乘子法將其轉化為無約束優化問題min?D,x∥Y?DX∥F2+λ∥xi∥1\min _{\mathbf{D}, \mathbf{x}}\|\mathbf{Y}-\mathbf{D} \mathbf{X}\|_{F}^{2}+\lambda\left\|\mathbf{x}_{i}\right\|_{1} D,xmin?∥Y?DX∥F2?+λ∥xi?∥1?其中主要是將∥xi∥0\left\|\mathbf{x}_{i}\right\|_{0}∥xi?∥0?用∥xi∥1\left\|\mathbf{x}_{i}\right\|_{1}∥xi?∥1?代替,這樣更加便于求解。這里要優化D\mathbf{D}D,X\mathbf{X}X兩個變量,而樣本集Y\mathbf{Y}Y是已知的,我們做法是固定一個變量,優化另一個,交替進行。下面分開來說:
已知D\mathbf{D}D,優化X\mathbf{X}X:這個問題有許多經典算法可以進行求解,例如Lasso,OMP,我們之后再進行補充,K-SVD中比較有特點是下面這種情況。
已知X\mathbf{X}X,優化D\mathbf{D}D:我們需要逐列來更新字典,記dk\mathbfze8trgl8bvbq_kdk?為字典D\mathbf{D}D的第kkk列,其實就是某個基向量,xTk\mathbf{x}_{T}^{k}xTk?為稀疏矩陣X\mathbf{X}X的第kkk行,按照本文之前的定義,其實就是該基向量不同樣本對應的系數。那么有:
∥Y?DX∥F2=∥Y?∑j=1KdjxTj∥F2=∥(Y?∑j≠kdjxTj)?dkxTk∥F2=∥Ek?dkxTk∥F2\begin{aligned}\|\mathbf{Y}-\mathbf{D} \mathbf{X}\|_{F}^{2} &=\left\|\mathbf{Y}-\sum_{j=1}^{K} \mathbfze8trgl8bvbq_{j} \mathbf{x}_{T}^{j}\right\|_{F}^{2} \\ &=\left\|\left(\mathbf{Y}-\sum_{j \neq k} \mathbfze8trgl8bvbq_{j} \mathbf{x}_{T}^{j}\right)-\mathbfze8trgl8bvbq_{k} \mathbf{x}_{T}^{k}\right\|_{F}^{2} \\ &=\left\|\mathbf{E}_{k}-\mathbfze8trgl8bvbq_{k} \mathbf{x}_{T}^{k}\right\|_{F}^{2} \end{aligned}∥Y?DX∥F2??=∥∥∥∥∥?Y?j=1∑K?dj?xTj?∥∥∥∥∥?F2?=∥∥∥∥∥∥????Y?j?=k∑?dj?xTj?????dk?xTk?∥∥∥∥∥∥?F2?=∥∥?Ek??dk?xTk?∥∥?F2??上式中定義殘差Ek=Y?∑j≠kdjxTj\mathbf{E}_{k}=\mathbf{Y}-\sum_{j \neq k} \mathbfze8trgl8bvbq_{j} \mathbf{x}_{T}^{j}Ek?=Y?∑j?=k?dj?xTj?就是除了當前更新基向量外其它基向量組合與樣本集的誤差,我們的目的就是讓該殘差項最小,因此有min?dk,xTk∥Ek?dkxTk∥F2\min _{\mathbfze8trgl8bvbq_{k}, \mathbf{x}_{T}^{k}}\left\|\mathbf{E}_{k}-\mathbfze8trgl8bvbq_{k} \mathbf{x}_{T}^{k}\right\|_{F}^{2} dk?,xTk?min?∥∥?Ek??dk?xTk?∥∥?F2?這個目標函數如何求呢?這里采用的就是SVD方法,如下:Ek′=UΣVT\mathbf{E}_{k}^{\prime}=U \Sigma V^{T} Ek′?=UΣVT其中Σ\SigmaΣ矩陣中的奇異值是從大到小排列的,我們要取最大奇異值對應的向量作為最優的dk,xTk\mathbfze8trgl8bvbq_{k}, \mathbf{x}_{T}^{k}dk?,xTk?,那么就是取矩陣UUU的第一個列向量u1=U(?,1)\mathbf{u}_{1}=U(\cdot, 1)u1?=U(?,1)作為dk\mathbfze8trgl8bvbq_{k}dk?,取矩陣VVV的第一個行向量與第一個奇異值的乘積作為xT′k\mathbf{x}_{T}^{\prime k}xT′k?,將其更新到xTk\mathbf{x}_{T}^{k}xTk?,這里值得注意的是,這里更新xTk\mathbf{x}_{T}^{k}xTk?并不能保證其稀疏性,要保證稀疏性還是需要回到上一種情況使用OMP算法進行稀疏編碼
這里為什么可以通過SVD算法求得最優解呢?關于這一點之后再補充,那么以上就完成了字典學習的過程,那么來看一段代碼印證下KSVD的流程。
2. python代碼
這段代碼是我從博客K-SVD字典學習及其實現(Python)扒過來的,該了一下opencv的接口,這里要注意的是,如果直接按照我這個代碼把整張圖片輸入的話其實降噪效果不理想的,為什么呢?
#!/usr/bin/python # -*- coding: UTF-8 -*- import numpy as np from sklearn import linear_model import scipy.misc from matplotlib import pyplot as plt import cv2class KSVD(object):def __init__(self, n_components, max_iter=30, tol=5000,n_nonzero_coefs=None):"""稀疏模型Y = DX,Y為樣本矩陣,使用KSVD動態更新字典矩陣D和稀疏矩陣X:param n_components: 字典所含原子個數(字典的列數):param max_iter: 最大迭代次數:param tol: 稀疏表示結果的容差:param n_nonzero_coefs: 稀疏度"""self.dictionary = Noneself.sparsecode = Noneself.max_iter = max_iterself.tol = tolself.n_components = n_componentsself.n_nonzero_coefs = n_nonzero_coefsdef _initialize(self, y):"""初始化字典矩陣"""u, s, v = np.linalg.svd(y)self.dictionary = u[:, :self.n_components]print(self.dictionary.shape)def _update_dict(self, y, d, x):"""使用KSVD更新字典的過程"""for i in range(self.n_components):index = np.nonzero(x[i, :])[0]if len(index) == 0:continued[:, i] = 0r = (y - np.dot(d, x))[:, index]u, s, v = np.linalg.svd(r, full_matrices=False)d[:, i] = u[:, 0].Tx[i, index] = s[0] * v[0, :]return d, xdef fit(self, y):"""KSVD迭代過程"""self._initialize(y)for i in range(self.max_iter):x = linear_model.orthogonal_mp(self.dictionary, y, n_nonzero_coefs=self.n_nonzero_coefs)e = np.linalg.norm(y - np.dot(self.dictionary, x))if e < self.tol:breakself._update_dict(y, self.dictionary, x)self.sparsecode = linear_model.orthogonal_mp(self.dictionary, y, n_nonzero_coefs=self.n_nonzero_coefs)return self.dictionary, self.sparsecodeif __name__ == '__main__':im_ascent = cv2.imread("./input.png", 0).astype(np.float)print(im_ascent.shape)ksvd = KSVD(300)dictionary, sparsecode = ksvd.fit(im_ascent)output = dictionary.dot(sparsecode)output = np.clip(output, 0, 255)cv2.imwrite("./output.png", output.astype(np.uint8))從代碼中可以看到,我們輸入的樣本集Y\mathbf{Y}Y矩陣就是整張圖片,這樣的話,單個樣本就是圖片的每一列的數據,這樣合理嗎?當然不合理,我覺得這就是造成僅僅通過上面這段代碼降噪效果不理想的原因,合理的做法應該是從圖像中提取patch,將patch轉換成一列一列的樣本數據,組成樣本集矩陣,然后再利用K-SVD算法從樣本集矩陣學習字典和稀疏矩陣。感興趣的同學可以根據該思路嘗試改下代碼,或者參考代碼nel215/ksvd
3. 結論
其中,Global Trained Dictionary指的是從一堆無噪聲的圖片中抽取的patch學習的字典,Adaptive Dictionary指的是從原噪聲圖像patch中學習的字典,從結果看,Adaptive Dictionary是要優于Global Trained Dictionary的
這篇博客沒有做太多實驗去驗證實驗效果,有問題歡迎交流
此外,這里我寫一個各種算法的總結目錄圖像降噪算法——圖像降噪算法總結,對圖像降噪算法感興趣的同學歡迎參考
總結
以上是生活随笔為你收集整理的图像降噪算法——稀疏表达:K-SVD算法的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 图像降噪算法——DnCNN / FFDN
- 下一篇: 图像降噪算法——低秩聚类:WNNM算法