Dictionary Learning(字典學習、稀疏表示以及其他)

第一部分 字典學習以及稀疏表示的概要

字典學習（Dictionary Learning）和稀疏表示（Sparse Representation）在學術界的正式稱謂應該是 稀疏字典學習（Sparse Dictionary Learning）。該算法理論包含兩個階段：字典構建階段（Dictionary Generate）和利用字典（稀疏的）表示樣本階段（Sparse coding with a precomputed dictionary）。這兩個階段（如下圖）的每個階段都有許多不同算法可供選擇，每種算法的誕生時間都不一樣，以至於稀疏字典學習的理論提出者已變得不可考。筆者嘗試找了Wikipedia和Google Scolar都無法找到這一系列理論的最早發起人。

 

這里有兩個問題是必須要預先解釋清楚：

問題1：我們為什么需要字典學習？

回答這個問題實際上就是要回答“稀疏字典學習 ”中的字典是怎么來的。做一個比喻，句子是人類社會最神奇的東西，人類社會的一切知識無論是已經發現的還是沒有發現的都必然要通過句子來表示出來（從某種意義上講，公式也是句子）。這樣說來，人類懂得的知識可要算是極為浩繁的。有人統計過人類每天新產生的知識可以裝滿一個2T（2048G）大小的硬盤。但無論有多少句子需要被書寫，對於一個句子來說它最本質的特征是什么呢？毫無疑問，是一個個構成這個句子的單詞（對英語來說）或字（對漢語來說）。所以我們可以很傲嬌的這樣認為，無論人類的知識有多么浩繁，也無論人類的科技有多么發達，一本長不過20厘米，寬不過15厘米，厚不過4厘米的新華字典或牛津字典足以表達人類從古至今乃至未來的所有知識，那些知識只不過是字典中字的排列組合罷了！直到這里，我相信相當一部分讀者或許在心中已經明白了字典學習的第一個好處—— 它實質上是對於龐大數據集的一種降維表示。第二，正如同字是句子最質朴的特征一樣， 字典學習總是嘗試學習蘊藏在樣本背后最質朴的特征（假如樣本最質朴的特征就是樣本最好的特征），這兩條原因同時也是這兩年深度學習之風日盛的情況下字典學習也開始隨之升溫的原因。題外話：現代神經科學表明，哺乳動物大腦的初級視覺皮層干就事情就是圖像的字典表示。

問題2：我們為什么需要稀疏表示？

回答這個問題毫無疑問就是要回答“稀疏字典學習”中稀疏兩字的來歷。不妨再舉一個例子。相信大部分人都有這樣一種感覺，當我們在解涉及到新的知識點的數學題時總有一種累心（累腦）的感覺。但是當我們通過艱苦卓絕的訓練將新的知識點牢牢掌握時，再解決與這個知識點相關的問題時就不覺得很累了。這是為什么呢？意大利羅馬大學的Fabio Babiloni教授曾經做過一項實驗，他們讓新飛行員駕駛一架飛機並采集了他們駕駛狀態下的腦電，同時又讓老飛行員駕駛飛機並也采集了他們駕駛狀態下的腦電。如下圖所示：

 

 

 

隨后Fabio教授計算出了兩類飛行員的大腦的活躍狀態，如下圖：

 

 

左圖是新飛行員（不熟練的飛行員）的大腦。圖中黃色的部分，是被認為活躍的腦區。右圖是老飛行員（熟練的飛行員）的大腦，黃色區域相比左邊的圖有明顯的減少。換言之，針對某一特定任務（這里是飛行），熟練者的大腦可以 調動盡可能少的腦區消耗盡可能少的能量進行同樣有效的計算（所以熟悉知識點的你，大腦不會再容易覺得累了），並且由於調動的腦區很少，大腦 計算速度也會變快，這就是我們稱熟練者為熟練者的原理所在。站在我們所要理解的稀疏字典學習的角度上來講就是 大腦學會了知識的稀疏表示。

 

稀疏表示的本質：用盡可能少的資源表示盡可能多的知識，這種表示還能帶來一個附加的好處，即計算速度快。

 

在懂得“字典”和“稀疏”各自的那點事兒以后，我們還要再講講稀疏和字典共同的那點兒事。或許在大腦中“字典”和“稀疏”是兩個不怎么想干的階段，畢竟“字典”涉及初級視覺皮層，而“稀疏”涉及前額葉皮層。但是在計算機中，“字典”和“稀疏”卻是一堆孿生兄弟。在學習樣本字典之初的時候，稀疏條件就已經被加入了。我們希望字典里的字可以盡能的少，但是卻可以盡可能的表示最多的句子。這樣的字典最容易滿足稀疏條件。也就是說，這個“字典”是這個“稀疏”私人訂制的。

 

第二部分 稀疏字典學習的Python實現

 

用Python實現稀疏字典學習需要三個前提條件

1.安裝NumPy

2.安裝SciPy

3.安裝Python機器學習工具包sklearn

為了避免過於麻煩的安裝，這里我干脆建議諸位讀者安裝Python的商業發行版 Anaconda，內含python集成開發環境和數百個常用的python支持包。具體安裝過程和使用細節參見我的博客 附錄D Python接口大法。

 

樣例一：圖片的稀疏字典學習

這段代碼來源於Python的Dictionary Learning的官方文獻教材，主要用途是教會用戶通過字典學習對圖片進行濾波處理。

step1：首先是各種工具包的導入和測試樣例的導入

1. from time import time
3. import matplotlib.pyplot as plt
4. import numpy as np
5. import scipy as sp
7. from sklearn.decomposition import MiniBatchDictionaryLearning
8. from sklearn.feature_extraction.image import extract_patches_2d
9. from sklearn.feature_extraction.image import reconstruct_from_patches_2d
10. from sklearn.utils.testing import SkipTest
11. from sklearn.utils.fixes import sp_version
13. if sp_version < (0, 12):
14.     raise SkipTest("Skipping because SciPy version earlier than 0.12.0 and "
15.                    "thus does not include the scipy.misc.face() image.")
16. try:
17.     from scipy import misc
18.     face = misc.face(gray=True)
19. except AttributeError:
20.     # Old versions of scipy have face in the top level package
21.     face = sp.face(gray=True)

 

第1行：導入time模塊，用於測算一些步驟的時間消耗

第3~5行：導入Python科學計算的基本需求模塊，主要包括NumPy（矩陣計算模塊）、SciPy（科學計算模塊）和matplotlib.pyplot模塊（畫圖）。有了這三個模塊，Python儼然已是基礎版的Matlab。

第7~11行：導入稀疏字典學習所需要的函數，下面分行解釋

第7行：導入MiniBatchDictionaryLearning，MiniBatch是字典學習的一種方法，這種方法專門應用於大數據情況下字典學習。當數據量非常大時，嚴格對待每一個樣本就會消耗大量的時間，而MiniBatch通過降低計算精度來換取時間利益，但是仍然能夠通過大量的數據學到合理的詞典。換言之，普通的DictionaryLearning做的是精品店，量少而精，但是價格高。MiniBatchDictionaryLearning做的是批發市場，量大不精，薄利多銷。

第8行：導入碎片提取函數extract_patches_2d。調用該函數將一張圖片切割為一個一個的pitch。如果一張圖片相當於一篇文章的話，那么該函數的目標就是把文章中的每個句子都找到，這樣才方便提取蘊藏在每個句子中的字。圖片和pitch的關系如下圖所示：

整張頭像的照片是個圖片，通過對圖片的分割可以將圖片分割為一個一個的小塊，也就是一個個Pitch。如果對pitch仍然不了解，只好請你看這個： http://blog.csdn.net/zouxy09/article/details/8775488

第9行：導入圖片復原函數reconstruct_from_patches_2d，它可以通過pitch復原一整張圖片。

第10行：導入測試工具nose下的異常拋出函數SkipTest

第11行：導入SciPy版本檢測函數sp_version用於檢測版本高低，版本低於0.12的SciPy沒有我們需要的樣本測試用例

第13~15行：檢測SciPy版本，如果版本太低就拋出一個異常。程序運行結束

第16~21行：嘗試打開樣本測試用例，如果打不開就拋出一個異常。

 

step2：通過測試樣例計算字典V

1. # Convert from uint8 representation with values between 0 and 255 to
2. # a floating point representation with values between 0 and 1.
3. face = face / 255.0
5. # downsample for higher speed
6. face = face[::2, ::2] + face[1::2, ::2] + face[::2, 1::2] + face[1::2, 1::2]
7. face = face / 4.0
8. height, width = face.shape
10. # Distort the right half of the image
11. print('Distorting image...')
12. distorted = face.copy()
13. distorted[:, width // 2:] += 0.075 * np.random.randn(height, width // 2)
15. # Extract all reference patches from the left half of the image
16. print('Extracting reference patches...')
17. t0 = time()
18. patch_size = (7, 7)
19. data = extract_patches_2d(distorted[:, :width // 2], patch_size)
20. data = data.reshape(data.shape[0], -1)
21. data -= np.mean(data, axis=0)
22. data /= np.std(data, axis=0)
23. print('done in %.2fs.' % (time() - t0))
24. print('Learning the dictionary...')
25. t0 = time()
26. dico = MiniBatchDictionaryLearning(n_components=100, alpha=1, n_iter=500)
27. V = dico.fit(data).components_
28. dt = time() - t0
29. print('done in %.2fs.' % dt)
31. plt.figure(figsize=(4.2, 4))
32. for i, comp in enumerate(V[:100]):
33.     plt.subplot(10, 10, i + 1)
34.     plt.imshow(comp.reshape(patch_size), cmap=plt.cm.gray_r,
35.                interpolation='nearest')
36.     plt.xticks(())
37.     plt.yticks(())
38. plt.suptitle('Dictionary learned from face patches\n' +
39.              'Train time %.1fs on %d patches' % (dt, len(data)),
40.              fontsize=16)
41. plt.subplots_adjust(0.08, 0.02, 0.92, 0.85, 0.08, 0.23)#left, right, bottom, top, wspace, hspace

 

第3行：讀入的face大小在0~255之間，所以通過除以255將face的大小映射到0~1上去

第6~7行：對圖形進行采樣，把圖片的長和寬各縮小一般。記住array矩陣的訪問方式      array[起始點：終結點（不包括）：步長]

第8行：圖片的長寬大小

第12行：將face的內容復制給distorted，這里不用等號因為等號在python中其實是地址的引用。

第13行：對照片的右半部分加上噪聲，之所以左半部分不加是因為教材想要產生一個對比的效果

第17行：開始計時，並保存在t0中

第18行：tuple格式的pitch大小

第19行：對圖片的左半部分（未加噪聲的部分）提取pitch

第20行：用reshape函數對data(94500,7,7)進行整形，reshape中如果某一位是-1，則這一維會根據（元素個數/已指明的維度）來計算這里經過整形后data變成（94500，49）

第21~22行：每一行的data減去均值除以方差，這是zscore標准化的方法

第26行：初始化MiniBatchDictionaryLearning類，並按照初始參數初始化類的屬性

第27行：調用fit方法對傳入的樣本集data進行字典提取，components_返回該類fit方法的運算結果，也就是我們想要的字典V

第31~41行：畫出V中的字典，下面逐行解釋

第31行：figsize方法指明圖片的大小，4.2英寸寬，4英寸高。其中一英寸的定義是80個像素點

第32行：循環畫出100個字典V中的字

第41行：6個參數與注釋后的6個屬性對應

運行程序，查看輸出結果：

step3：畫出標准圖像和真正的噪聲，方便同之后字典學習學到的噪聲相比較

1. def show_with_diff(image, reference, title):
2.     """Helper function to display denoising"""
3.     plt.figure(figsize=(5, 3.3))
4.     plt.subplot(1, 2, 1)
5.     plt.title('Image')
6.     plt.imshow(image, vmin=0, vmax=1, cmap=plt.cm.gray,
7.                interpolation='nearest')
8.     plt.xticks(())
9.     plt.yticks(())
10.     plt.subplot(1, 2, 2)
11.     difference = image - reference
13.     plt.title('Difference (norm: %.2f)' % np.sqrt(np.sum(difference ** 2)))
14.     plt.imshow(difference, vmin=-0.5, vmax=0.5, cmap=plt.cm.PuOr,
15.                interpolation='nearest')
16.     plt.xticks(())
17.     plt.yticks(())
18.     plt.suptitle(title, size=16)
19.     plt.subplots_adjust(0.02, 0.02, 0.98, 0.79, 0.02, 0.2)
21. show_with_diff(distorted, face, 'Distorted image')

程序輸出如下圖所示：

step4：測試不同的字典學習方法和參數對字典學習的影響

1. print('Extracting noisy patches... ')
2. t0 = time()
3. data = extract_patches_2d(distorted[:, width // 2:], patch_size)
4. data = data.reshape(data.shape[0], -1)
5. intercept = np.mean(data, axis=0)
6. data -= intercept
7. print('done in %.2fs.' % (time() - t0))
9. transform_algorithms = [
10.     ('Orthogonal Matching Pursuit\n1 atom', 'omp',
11.      {'transform_n_nonzero_coefs': 1}),
12.     ('Orthogonal Matching Pursuit\n2 atoms', 'omp',
13.      {'transform_n_nonzero_coefs': 2}),
14.     ('Least-angle regression\n5 atoms', 'lars',
15.      {'transform_n_nonzero_coefs': 5}),
16.     ('Thresholding\n alpha=0.1', 'threshold', {'transform_alpha': .1})]
18. reconstructions = {}
19. for title, transform_algorithm, kwargs in transform_algorithms:
20.     print(title + '...')
21.     reconstructions[title] = face.copy()
22.     t0 = time()
23.     dico.set_params(transform_algorithm=transform_algorithm, **kwargs)
24.     code = dico.transform(data)
25.     patches = np.dot(code, V)
27.     patches += intercept
28.     patches = patches.reshape(len(data), *patch_size)
29.     if transform_algorithm == 'threshold':
30.         patches -= patches.min()
31.         patches /= patches.max()
32.     reconstructions[title][:, width // 2:] = reconstruct_from_patches_2d(
33.         patches, (height, width // 2))
34.     dt = time() - t0
35.     print('done in %.2fs.' % dt)
36.     show_with_diff(reconstructions[title], face,
37.                    title + ' (time: %.1fs)' % dt)
39. plt.show()

 

第3行：提取照片中被污染過的右半部進行字典學習。

第10~16行：四中不同的字典表示策略

第23行：通過set_params對第二階段的參數進行設置

第24行：transform根據set_params對設完參數的模型進行字典表示，表示結果放在code中。code總共有100列，每一列對應着V中的一個字典元素，所謂稀疏性就是code中每一行的大部分元素都是0，這樣就可以用盡可能少的字典元素表示回去。

第25行：code矩陣乘V得到復原后的矩陣patches

第28行：將patches從（94500，49）變回（94500，7，7）

第32行：通過reconstruct_from_patches_2d函數將patches重新拼接回圖片

 

該程序輸出為四中不同轉換算法下的降噪效果：