機器學習（Machine Learning）是一門涉及多個領域的交叉學科，它涵蓋了計算機科學、統計學、概率論及最優化理論等多門學科的知識。目前，機器學習在疾病的鑑別診斷，藥物的生產研究，臨床試驗研究，放射影像學等醫學領域已有十分廣泛的應用。
今天推出機器學習系列筆記第8期，本期分享內容為：機器學習之降維算法Part2。（筆者使用的是Windows系統）
SVD有一種驚人的數學性質，即是它可以跳過數學神秘的宇宙，不計算協方差矩陣，直接找出一個新特徵向量組成的n維空間，而這個n維空間就是奇異值分解後的右矩陣V^T。該矩陣有著如下性質：
k就是n_components，是我們降維後希望得到的維度。若X為(m,n)的特徵矩陣，V^T就是結構為(n,n)的矩陣，取這個矩陣的前k行（進行切片），即將V轉換為結構為(k,n)的矩陣。而
SVD在矩陣分解中的過程比PCA簡單快速，但是SVD的信息量衡量指標比較複雜，要理解」奇異值「遠不如理解」方差「來得容易，因此，sklearn將降維流程拆成了兩部分：一部分是計算特徵空間V，由奇異值分解完成，另一部分是映射數據和求解新特徵矩陣，由主成分分析完成，實現了用SVD的性質減少計算量，卻讓信息量的評估指標是方差，具體流程如下圖：
也就是說，可以把SVD當作PCA的一種求解方法，其實指的就是在矩陣分解時不使用PCA本身的特徵值分解，而使用奇異值分解來減少計算量。這種方法確實存在，但在sklearn中，矩陣U和Σ雖然會被計算出來（同樣也是一種比起PCA來說簡化非常多的數學過程，不產生協方差矩陣），但完全不會被用到，也無法調取查看或者使用，因此我們可以認為，U和Σ在fit過後就被遺棄了。奇異值分解追求的僅僅是V，只要有了V，就可以計算出降維後的特徵矩陣。在transform過程之後，fit中奇異值分解的結果除了V(k,n)以外，就會被捨棄，而V(k,n)會被保存在屬性components_ 當中，可以調用查看。
PCA(2).fit(X).components_
PCA(2).fit(X).components_.shape
#返回降維後的新特徵空間V（k，n），2行4列的二維數組
參數svd_solver是在降維過程中，用來控制矩陣分解的一些細節的參數。有四種模式可選："auto", "full", "arpack","randomized"，默認」auto"。
"auto"：基於X.shape和n_components的默認策略來選擇分解器：如果輸入數據的尺寸大於500x500且要提取的特徵數小於數據最小維度min(X.shape)的80％，啟用效率更高的」randomized「方法。否則，精確完整的SVD（參考「full」的情況）將被計算，截斷將會在矩陣被分解完成後有選擇地發生。
"full"：從scipy.linalg.svd中調用標準的LAPACK分解器來生成精確完整的SVD，適合數據量比較適中，計算時間充足的情況，生成的精確完整的SVD的結構為：
"arpack"：從scipy.sparse.linalg.svds調用ARPACK分解器來運行截斷奇異值分解(SVD truncated)，分解時就將特徵數量降到n_components中輸入的數值k，可以加快運算速度，適合特徵矩陣很大的時候，但一般用於特徵矩陣為稀疏矩陣的情況（稀疏矩陣：每一列都是由0和1組成，並且大部分都是0），此過程包含一定的隨機性。截斷後的SVD分解出的結構為：
而參數random_state在參數svd_solver的值為"arpack" or "randomized"的時候生效，可以控制這兩種SVD模式中的隨機模式。通常我們就選用」auto「，不必對這個參數糾結太多。
2.3.3 重要屬性components_通常來說，在新的特徵矩陣生成之前，我們無法知曉PCA都建立了怎樣的新特徵向量，新特徵矩陣生成之後也不具有可讀性，我們無法判斷新特徵矩陣的特徵是從原數據中的什麼特徵組合而來，新特徵雖然帶有原始數據的信息，卻已經不是原數據上代表著的含義了。
但是其實，在矩陣分解時，PCA是有目標的：在原有特徵的基礎上，找出能夠讓信息儘量聚集的新特徵向量。在sklearn使用的PCA和SVD聯合的降維方法中，這些新特徵向量組成的新特徵空間其實就是V(k,n)。當V(k,n)是數字時，我們無法判斷V(k,n)和原有的特徵究竟有著怎樣的數學聯繫。但是，如果原特徵矩陣是圖像，V(k,n)這個空間矩陣也可以被可視化的話，我們就可以通過兩張圖來比較，就可以看出新特徵空間究竟從原始數據裡提取了什麼重要的信息。
讓我們來看一個，人臉識別中屬性components_的運用。
1.導入需要的庫和模塊
from sklearn.datasets import fetch_lfw_people #7個人的人臉數據
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import numpy as np
2.實例化數據集，探索數據
faces = fetch_lfw_people(min_faces_per_person=60) #實例化
faces.data.shape
#行是樣本
#列是樣本相關的所有特徵
faces.images.shape
#返回（1277，62，47）
#1277是矩陣中圖像的個數
#62是每個圖像的特徵矩陣的行
#47是每個圖像的特徵矩陣的列
X = faces.data
3.看看圖像什麼樣？將原特徵矩陣進行可視化
#數據本身是圖像，和數據本身只是數字，使用的可視化方法不同
#創建畫布和子圖對象
fig, axes = plt.subplots(3,8 #3行8列24個子圖
,figsize=(8,4)
,subplot_kw = {"xticks":[],"yticks":[]} #不要顯示坐標軸
)
fig #畫好的畫布
axes
#不難發現，axes中的一個對象對應fig中的一個空格
#我們希望，在每一個子圖對象中填充圖像（共24張圖），因此我們需要寫一個在子圖對象中遍歷的循環
#二維結構，可以有兩種循環方式，一種是使用索引，循環一次同時生成一列上的三個圖
#另一種是把數據拉成一維，循環一次只生成一個圖
#在這裡，究竟使用哪一種循環方式，是要看我們要畫的圖的信息，儲存在一個怎樣的結構裡#我們使用子圖對象.imshow 來將圖像填充到空白畫布上
#而imshow要求的數據格式必須是一個(m,n)格式的矩陣，即每個數據都是一張單獨的圖
#因此我們需要遍歷的是faces.images，其結構是(1277, 62, 47)
#要從一個數據集中取出24個圖，明顯是一次性的循環切片[i,:,:]來得便利
#因此我們要把axes的結構拉成一維來循環
axes.flat #降維成一維數組
[*axes.flat] #使用該方法查看惰性對象裡的內容
enumerate(axes.flat) #也是惰性對象，添加了索引
#填充圖像
for i, ax in enumerate(axes.flat):
ax.imshow(faces.images[i,:,:]
,cmap="gray"#選擇色彩的模式
)
#https://matplotlib.org/tutorials/colors/colormaps.html 查看顏色取值
4.建模降維，提取新特徵空間矩陣
#原本有2900維，我們現在來降到150維
pca = PCA(150).fit(X)
V = pca.components_ #V代表用來映射的新特徵向量空間
V.shape
5.將新特徵空間矩陣可視化
fig, axes = plt.subplots(3,8,figsize=(8,4),subplot_kw = {"xticks":[],"yticks":[]})
for i, ax in enumerate(axes.flat): #遍歷循環
ax.imshow(V[i,:].reshape(62,47),cmap="gray") #reshape還原圖像
可以看出比起降維前的數據，新特徵空間可視化後的人臉非常模糊，這是因為原始數據還沒有被映射到特徵空間中。但是可以看出，整體比較亮的圖片，獲取的信息較多，整體比較暗的圖片，卻只能看見黑漆漆的一塊。在比較亮的圖片中，眼睛，鼻子，嘴巴，都相對清晰，臉的輪廓，頭髮之類的比較模糊。
這說明，新特徵空間裡的特徵向量們，大部分是"五官"和"亮度"相關的向量，所以新特徵向量上的信息肯定大部分是由原數據中和"五官"和"亮度"相關的特徵中提取出來的。
我們通過可視化新特徵空間V，解釋了一部分降維後的特徵：雖然顯示出來的數字看著不知所云，但畫出來的圖表示，這些特徵是和」五官「以及」亮度「有關的。這也再次證明了，PCA能夠將原始數據集中重要的數據進行聚集。
2.4 重要接口inverse_transform在sklearn中，我們通過讓原特徵矩陣X右乘新特徵空間矩陣V(k,n)來生成新特徵矩陣X_dr，那理論上來說。讓新特徵矩陣X_dr右乘V(k,n)的逆矩陣V_(k,n)，就可以將新特徵矩陣X_dr還原為X。那sklearn是否這樣做了呢？讓我們來看看下面的案例。
2.4.1 迷你案例：用人臉識別看PCA降維後的信息保存量1. 導入需要的庫和模塊
from sklearn.datasets import fetch_lfw_people #7個人的人臉數據
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import numpy as np
2.導入數據，探索數據
faces = fetch_lfw_people(min_faces_per_person=60)
faces.images.shape
faces.data.shape
X = faces.data
3. 建模降維，獲取降維後的特徵矩陣X_dr
pca = PCA(150) #實例化
X_dr = pca.fit_transform(X) #擬合+提取結果
X_dr.shape
4. 將降維後矩陣用inverse_transform返回原空間
X_inverse = pca.inverse_transform(X_dr)
X_inverse.shape
#期待X_inverse應該和原數據有相同的結果，如果相同，我們就說inverse_transform實現了降維過程的逆轉
5. 將特徵矩陣X和X_inverse可視化
fig, ax = plt.subplots(2,10,figsize=(10,2.5)
,subplot_kw={"xticks":[],"yticks":[]}
)
#和2.3.3節中的案例一樣，我們需要對子圖對象進行遍歷的循環，來將圖像填入子圖中
#那在這裡，我們使用怎樣的循環？
#現在我們的ax中是2行10列，第一行是原數據，第二行是inverse_transform後返回的數據
#所以我們需要同時循環兩份數據，即一次循環畫一列上的兩張圖，而不是把ax拉平
for i in range(10):
ax[0,i].imshow(faces.image[i,:,:],cmap="binary_r")
ax[1,i].imshow(X_inverse[i].reshape(62,47),cmap="binary_r") #imshow只接受二維數組
可以明顯看出，這兩組數據可視化後，由降維後再通過inverse_transform轉換回原維度的數據畫出的圖像和原數據畫的圖像大致相似，但原數據的圖像明顯更加清晰。這說明，inverse_transform並沒有實現數據的完全逆轉。這是因為，在降維的時候，部分信息已經被捨棄了，X_dr中往往不會包含原數據100%的信息，所以在逆轉的時候，即便維度升高，原數據中已經被捨棄的信息也不可能再回來了。所以，降維不是完全可逆的。
2.4.2 迷你案例：用PCA做噪音過濾降維的目的之一就是希望拋棄掉對模型帶來負面影響的特徵，而帶有效信息的特徵的方差應該是遠大於噪音的，所以相比噪音，有效的特徵所帶的信息應該不會在PCA過程中被大量拋棄。
inverse_transform能夠在不恢復原始數據的情況下，將降維後的數據返回到原本的高維空間，即是說能夠實現」保證維度，但去掉方差很小特徵所帶的信息「。利用inverse_transform的這個性質，我們能夠實現噪音過濾。
1. 導入所需要的庫和模塊
from sklearn.datasets import load_digits
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import numpy as np
2. 導入數據，探索數據
digits = load_digits()
digits.data.shape
#(1797,64)
set(digits.target.tolist()) #去掉重複項查看標籤
3.定義畫圖函數
def plot_digits(data):
#data的結構必須是(m,n)，並且n要能夠被分成(8,8)這樣的結構
fig, axes = plt.subplots(4,10,figsize=(10,4)
,subplot_kw = {"xticks":[],"yticks":[]}
)
for i, ax in enumerate(axes.flat):
ax.imshow(data[i].reshape(8,8),cmap="binary")
plot_digits(digits.data)
4. 為數據加上噪音
import numpy as np
rng = np.random.RandomState(42)
#在指定的數據集中，隨機抽取服從正態分布的數據
#兩個參數，分別是指定的數據集，和抽取出來的正態分布的方差
noisy = rng.normal(digits.data,2)
#normal()從輸入的數據集中隨機抽取一個滿足正態分布的數據，2是方差的大小
plot_digits(noisy)
5.降維
pca = PCA(0.5,svd_solver="full").fit(noisy) #讓降維後的數據帶有50%的原始數據
X_dr = pca.transform(noisy)
X_dr.shape
6. 逆轉降維結果，實現降噪
without_noise = pca.inverse_transform(X_dr)
plot_digits(without_noise)
返回的圖片雖然不能達到像第一張圖清晰，但是已經達到了降噪效果。在現實生活中，收上來的數據就是帶有噪音的，此時我們就可以將它降維後再返回原有的數據結構，就相當於是在不調整維度的狀況下進行了噪音過濾，提取出最重要的特徵信息。
2.5 重要接口，參數和屬性總結 3 案例：PCA對手寫數字數據集的降維1.導入需要的模塊和庫
from sklearn.decomposition import PCA
from sklearn.ensemble import RandomForestClassifier as RFC
from sklearn.model_selection import cross_val_score
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
2.導入數據，探索數據
data = pd.read_csv(r"C:\work\learnbetter\micro-class\week 3 Preprocessing\digit recognizor.csv")
X = data.iloc[:,1:]
y = data.iloc[:,0]
X.shape
#(42000, 784)
3. 畫累計方差貢獻率曲線，找最佳降維後維度的範圍
pca_line = PCA().fit(X)
plt.figure(figsize=[20,5])
plt.plot(np.cumsum(pca_line.explained_variance_ratio_))
plt.xlabel("number of components after dimension reduction")
plt.ylabel("cumulative explained variance ratio")
plt.show()
可以看到從0開始，曲線一開始增長非常快，0-200之間有我們需要的取值。
4. 降維後維度的學習曲線，繼續縮小最佳維度的範圍
#======【TIME WARNING：2mins 30s】======#
score = []
for i in range(1,101,10):
X_dr = PCA(i).fit_transform(X)
once = cross_val_score(RFC(n_estimators=10,random_state=0)
,X_dr,y,cv=5).mean()
score.append(once)
plt.figure(figsize=[20,5])
plt.plot(range(1,101,10),score) #x取值，y取值
plt.show()
可以看到，在特徵數為10左右時，score已經達到90左右，可以大致判斷最佳取值在20附近。
5. 細化學習曲線，找出降維後的最佳維度
#======【TIME WARNING：2mins 30s】======#
score = []
for i in range(10,25):
X_dr = PCA(i).fit_transform(X)
once = cross_val_score(RFC(n_estimators=10,random_state=0),X_dr,y,cv=5).mean()
score.append(once)
plt.figure(figsize=[20,5])
plt.plot(range(10,25),score)
plt.show()
大致判斷20-22之間可找到最高score對應的特徵數。
6. 導入找出的最佳維度進行降維，查看模型效果
X_dr = PCA(21).fit_transform(X)
#======【TIME WARNING:1mins 30s】======#
cross_val_score(RFC(n_estimators=10,random_state=0),X_dr,y,cv=5).mean()
#0.9165476190476192
cross_val_score(RFC(n_estimators=100,random_state=0),X_dr,y,cv=5).mean()
#0.9432380952380953
#對比不同n_estimators的分數差異，進行參數選擇
模型效果還好，但有沒有什麼辦法能夠進一步提高模型的表現呢？
7. 突發奇想，特徵數量已經不足原來的3%，換模型怎麼樣？
在之前的建模過程中，因為計算量太大，所以我們一直使用隨機森林，但事實上，我們知道KNN的效果比隨機森林更好，KNN在未調參的狀況下已經達到96%的準確率，而隨機森林在未調參前只能達到93%，這是模型本身的限制帶來的，這個數據使用KNN效果就是會更好。現在我們的特徵數量已經降到不足原來的3%，可以使用KNN了嗎？
#更換模型
from sklearn.neighbors import KNeighborsClassifier as KNN
cross_val_score(KNN(),X_dr,y,cv=5).mean() #k默認值為5
#0.9699761904761905
8. KNN的k值學習曲線
#======【TIME WARNING: 2mins 30s】======#
score = []
for i in range(10):
X_dr = PCA(21).fit_transform(X)
once = cross_val_score(KNN(i+1),X_dr,y,cv=5).mean()
score.append(once)
plt.figure(figsize=[20,5])
plt.plot(range(10),score)
plt.show()
9.定下超參數後，模型效果如何，模型運行時間如何？
cross_val_score(KNN(4),X_dr,y,cv=5).mean()
#0.9686190476190475
#=======【TIME WARNING: 3mins】======#
#%%timelit
#cross_val_score(KNN(4),X_dr,y,cv=5).mean()
可以發現，原本785列的特徵被縮減到21列之後，用KNN跑出了目前位置這個數據集上最好的結果。再進行更細緻的調整，也許可以將KNN的效果調整到98%以上。PCA為我們提供了無限的可能——不用再因為數據量太龐大而被迫選擇更加複雜的模型。