主成分分析一個(gè)非監(jiān)督學(xué)習(xí)算法，主要用于數(shù)據(jù)降維，通過(guò)降維可以發(fā)現(xiàn)數(shù)據(jù)更容易理解的特征，其他作用也有可視化、降噪等。
假設(shè)現(xiàn)有樣本的分布如圖。
樣本有兩個(gè)特征，如果對(duì)樣本進(jìn)行降維，首先可以考慮基于坐標(biāo)軸進(jìn)行降維。會(huì)有如下兩種方式。
分別對(duì)應(yīng)基于特征2降維和基于特征1降維。就降維之后的效果而言，右圖會(huì)更好一些，因?yàn)闃颖局g的間距比較大，樣本之間有更好的可區(qū)分度。但右圖這種方案不一定是最好的方案。
如果將樣本映射到這樣的一根斜線上效果會(huì)更好，因此，降維的目標(biāo)就是找到樣本間距最大的軸。也就是樣本映射到這個(gè)軸之后，使得樣本間距最大。
在數(shù)學(xué)上，樣本間距離反映了樣本的離散程度，而方差可以很好的表示這種程度，因此我們定義樣本間距離——方差。
主成分分析PCA
PCA算法是最經(jīng)典的數(shù)據(jù)降維算法，步驟大致如下：
①將樣本的均值歸零（demean）
②求一個(gè)軸的方向w，使得映射后的方差值最大
均值歸零
原有樣本分布
歸零之后
歸零的作用是使得樣本的均值為0，這樣就可以在計(jì)算方差時(shí)，將原有公式中的均值置為0，即獲得圖中均值歸零后的方差表達(dá)式。
注意這里的Xi表示的是映射之后樣本的數(shù)據(jù)，也就是降維之后的數(shù)據(jù)。
映射軸和方差
映射的軸可以視為一個(gè)向量，w = (w1, w2, ..., wn)，n即為樣本的特征數(shù)。
方差公式
我們的目標(biāo)就是使得這個(gè)公式得到的方差值最大。
如果將樣本的數(shù)據(jù)映射到坐標(biāo)系中，那么每個(gè)樣本就可以看做是個(gè)向量，而樣本的均值也可以看做是一個(gè)向量。這樣我們的方差公式可以變?yōu)?/p>
也就是求每個(gè)向量和均值向量之間距離的平方和再除以向量個(gè)數(shù)，得到的結(jié)果本質(zhì)上也是方差。由于前一步中的均值歸零操作，使得樣本在各個(gè)維度的均值都是0，那么進(jìn)一步化簡(jiǎn)公式為
假設(shè)在數(shù)據(jù)有兩個(gè)維度，則可以得到一個(gè)二維平面。
Xi是樣本點(diǎn)，有Xi,1和Xi,2兩個(gè)特征，假設(shè)我們求得的最大樣本間距的映射軸是w表示的向量，設(shè)w=(w1, w2)，那么Xi映射到w上的點(diǎn)如圖中藍(lán)色線表示。
高中數(shù)學(xué)中，兩個(gè)向量相乘等于兩個(gè)向量的模乘以兩個(gè)向量的夾角余弦值。
通常在運(yùn)算中，我們會(huì)將向量w設(shè)為單位向量（即模長(zhǎng)為1的向量），因此等式右側(cè)的w的模長(zhǎng)可以省略。
又因?yàn)閄i的模長(zhǎng)乘以?shī)A角余弦值正好是藍(lán)色部分的長(zhǎng)度，因此我們得到最后一個(gè)等式。
因此，我們對(duì)之前求方差的算式進(jìn)行修改得到
這里其實(shí)不應(yīng)使用取模的符號(hào)，因?yàn)閄i樣本有n個(gè)屬性，可以理解為一個(gè)n個(gè)元素的向量，w也是有n個(gè)元素的向量，兩個(gè)向量相乘的結(jié)果是一個(gè)實(shí)數(shù)。而放到矩陣來(lái)說(shuō)，Xi是一個(gè)1 * n的矩陣，w理解為一個(gè)n * 1的矩陣，相乘之后也是一個(gè)實(shí)數(shù)，因此，最佳的表達(dá)式應(yīng)該是
展開(kāi)來(lái)看
合并一下w各個(gè)值
得到這個(gè)算式最大值的過(guò)程可以采用梯度上升法求解。
梯度上升法
回憶梯度下降法，通過(guò)求函數(shù)在一個(gè)點(diǎn)的導(dǎo)數(shù)值，得到函數(shù)在該點(diǎn)的斜率，之后使用減去斜率與學(xué)習(xí)率乘積的方式下降，逐漸靠近極小值點(diǎn)。
梯度上升法與之類似，上升不過(guò)就是將減法變?yōu)榧臃?。第一步還是對(duì)方差函數(shù)f在每個(gè)特征維度求導(dǎo)。
對(duì)括號(hào)中的式子化簡(jiǎn)
提出每個(gè)式子的Xi * w
此時(shí)得到一個(gè)向量和X樣本矩陣，由于向量w和X中每個(gè)樣本都有相乘的操作，因此，可以將其寫(xiě)為2/m * (Xw)，左側(cè)算式可以理解為Xi * w需要和每個(gè)樣本的第i列做一次乘法，因此，可以認(rèn)為最終是要乘以整個(gè)樣本X的矩陣。因此得到最右側(cè)的兩個(gè)式子，下面的算式是為了方便運(yùn)算而對(duì)上面的式子進(jìn)行轉(zhuǎn)換得到的。
主成分分析
單個(gè)主成分分析
先模擬一個(gè)數(shù)據(jù)集
import numpy as np
import matplotlib.pyplot as plt
X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0, 10, size=100)
plt.scatter(X[:, 0], X[:, 1])
plt.show()
該數(shù)據(jù)集有兩個(gè)屬性，打印初始圖像
接下來(lái)對(duì)數(shù)據(jù)進(jìn)行主成分分析，第一步是均值歸零，定義相應(yīng)的函數(shù)
'''均值歸零'''
def demean(X):
# axis=0即按列計(jì)算均值，及每個(gè)屬性的均值，1則是計(jì)算行的均值
return (X - np.mean(X, axis=0))
X_demean = demean(X)
plt.scatter(X_demean[:, 0], X_demean[:, 1])
plt.show()
計(jì)算均值調(diào)用的是numpy自帶的函數(shù)，均值歸零之后打印歸零后的圖像
分布沒(méi)有變化，但從坐標(biāo)軸可以看出均值歸零后的效果。
接下來(lái)就是對(duì)方差函數(shù)和其導(dǎo)數(shù)函數(shù)的定義
'''方差函數(shù)'''
def f(w, X):
return np.sum((X.dot(w)**2)) / len(X)
'''方差函數(shù)導(dǎo)數(shù)'''
def df_math(w, X):
return X.T.dot(X.dot(w)) * 2. / len(X)
'''將向量化簡(jiǎn)為單位向量'''
def direction(w):
return w / np.linalg.norm(w)
'''梯度上升法'''
def gradient_ascent(w, X, eta, n_iter=1e4, epsilon=0.0001):
'''
梯度上升法
:param w:
:param X:
:param eta:
:param n_iter:
:param epsilon:
:return:
'''
#先化簡(jiǎn)w為單位向量，方便運(yùn)算
w = direction(w)
i_iter = 0
while i_iter < n_iter:
gradient = df_math(w, X)
last_w = w
w += gradient * eta
#每次更新后將w化簡(jiǎn)為單位向量
w = direction(w)
if abs(f(w, X) - f(last_w, X)) < epsilon:
break
i_iter += 1
return w
和線性回歸梯度下降尋找極小值點(diǎn)類似，只不過(guò)這次的方向是上升。而且有需要注意限定循環(huán)次數(shù)，避免學(xué)習(xí)率過(guò)大造成陷入無(wú)限循環(huán)中。測(cè)試數(shù)據(jù)降維的結(jié)果
w_init = np.random.random(X_demean.shape[1])
eta = 0.01
w = gradient_ascent(w_init, X_demean, eta)
print(w)
plt.scatter(X_demean[:, 0], X_demean[:, 1])
plt.plot([0, w[0] * 30], [0, w[1] * 30], color='r')
plt.show()
注意一點(diǎn)，線性回歸中，通常將特征系數(shù)θ的值設(shè)為全部為0的向量，但在主成分分析中w的初始值不能為0?。?！
得到映射向量如圖中紅線
向量w的值為
[0.75748163 0.65285648]
多個(gè)主成分分析
剛剛模擬了兩個(gè)屬性中選取一個(gè)主成分的過(guò)程，實(shí)際的降維過(guò)程可能會(huì)涉及到數(shù)據(jù)在多個(gè)維度的降維，需要依次求解多個(gè)主成分。
求解第一個(gè)主成分后，假設(shè)得到映射的軸為w所表示的向量，如果此時(shí)需要求解第二個(gè)主成分怎么做。
需要先將數(shù)據(jù)集在第一個(gè)主成分上的分量去掉，然后在沒(méi)有第一個(gè)主成分的基礎(chǔ)上再尋找第二個(gè)主成分。
由之前的推導(dǎo)所知，藍(lán)色部分的模長(zhǎng)是Xi向量和w向量的乘積，又因?yàn)閣是單位向量，向量的模長(zhǎng)乘以方向上的單位向量就可以得到這個(gè)向量，去掉Xi在w方向的分量得到新的數(shù)據(jù)Xi' = Xi - Xipro。
求第二主成分就是在新的數(shù)據(jù)Xi'上尋找第一主成分。以此類推，之后的主成分求法也是這個(gè)套路。
'''均值歸零'''
def demean(X):
# axis=0即按列計(jì)算均值，及每個(gè)屬性的均值，1則是計(jì)算行的均值
return (X - np.mean(X, axis=0))
'''方差函數(shù)'''
def f(w, X):
return np.sum((X.dot(w)**2)) / len(X)
'''方差函數(shù)導(dǎo)數(shù)'''
def df_ascent(w, X):
return X.T.dot(X.dot(w)) / 2 * len(X)
'''將向量化簡(jiǎn)為單位向量'''
def direction(w):
return w / np.linalg.norm(w)
'''尋找第一主成分'''
def first_component(w_init, X, eta, n_iter=1e4, epsilon=0.0001):
w = direction(w_init)
i_iter = 0
while i_iter < n_iter:
last_w = w
gradient = df_ascent(w, X)
w += eta * gradient
w = direction(w)
if abs(f(w, X) - f(last_w, X)) < epsilon:
break
i_iter += 1
return w
'''取前n個(gè)主成分'''
def first_n_component(n, X, eta, n_iter=1e4, epsilon=0.0001):
'''先對(duì)數(shù)據(jù)均值歸零'''
X = demean(X)
#res記錄每個(gè)主成分
res = []
#進(jìn)行n次
for i in range(n):
#每次初始化w向量，注意不能是0向量
w = np.random.random(X.shape[1])
#尋找當(dāng)前數(shù)據(jù)的第一主成分并記錄到res
w = first_component(w, X, eta)
res.append(w)
#每次減去數(shù)據(jù)在主成分方向的分量獲得新數(shù)據(jù)
X = X - X.dot(w).reshape(-1, 1) * w
return res
代碼中的前5個(gè)函數(shù)和第一主成分分析中的一樣，注意最后一個(gè)函數(shù)。每次求當(dāng)前數(shù)據(jù)的第一主成分，并減去數(shù)據(jù)在第一主成分上的分量獲得的新數(shù)據(jù)繼續(xù)求第一主成分，知道完成n次為止。
出去數(shù)據(jù)在第一主成分的分量，如剛才的算式推導(dǎo)的，對(duì)于每個(gè)樣本Xi減去主成分的過(guò)程可以寫(xiě)為
X_new[i] = X[i] - (X[i].dot(w)) * w
X[i]是第i個(gè)樣本，可以看作是1 * n的向量，w是一個(gè)1 * n的向量，通過(guò)dot()方法可以得到一個(gè)常數(shù)，即為X[i]在w方向上的模長(zhǎng)，模長(zhǎng)再乘以w即為在w方向上的分量。
對(duì)于含有m個(gè)樣本的數(shù)據(jù)集X，可以將整個(gè)過(guò)程視為
for i in range(len(X)):
X_new[i] = X[i] - X[i].dot(w) * w
for循環(huán)方便理解但更好的方式是直接使用矩陣乘法的方式解決，上面的for循環(huán)等價(jià)于
X_new = X - X.dot(w).reshape(-1, 1) * w
X是m * n的矩陣，w是1 * n的矩陣，二者dot()之后得到m * 1的矩陣，這里需要使用reshape(-1, 1)使得X和w相乘之后的矩陣真正是一個(gè)m * 1的矩陣，之后再乘以w，得到的m * n的矩陣即為每個(gè)樣本每個(gè)維度在w上的分量。再用原始數(shù)據(jù)X減去即可。
檢驗(yàn)下第一次降維后的效果
'''模擬數(shù)據(jù)'''
X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0, 10, size=100)
'''均值歸零'''
X_demean = demean(X)
'''獲得第一主成分'''
w_init = np.random.random(X_demean.shape[1])
eta = 0.01
w = first_component(w_init, X_demean, eta)
print(w)
'''畫(huà)出第一主成分所在向量'''
plt.scatter(X_demean[:, 0], X_demean[:, 1])
plt.plot([0, w[0] * 30], [0, w[1] * 30], color='r')
plt.show()
'''X減去X在第一主成分上的分量'''
X2 = X_demean - X_demean.dot(w).reshape(-1, 1) * w
plt.scatter(X2[:, 0], X2[:, 1])
plt.show()
下面在新的數(shù)據(jù)上獲取第二主成分
'''獲取第二主成分'''
w2_init = np.random.random(X2.shape[1])
w2 = first_component(w2_init, X2, eta)
plt.scatter(X2[:, 0], X2[:, 1])
plt.plot([0, w2[0] * 30], [0, w2[1] * 30], color='r')
plt.show()
得到圖像
因?yàn)槭嵌S的數(shù)據(jù)，因此前后所得的兩個(gè)向量應(yīng)該是垂直的關(guān)系，如上文圖中所示，兩個(gè)相互垂直的向量相乘值為0，實(shí)驗(yàn)下
X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0, 10, size=100)
eta = 0.01
res = first_n_component(2, X, eta)
print(res)
print(res[0].dot(res[1]))
調(diào)用之前定義的first_n_component函數(shù)，記錄前后兩次的w，打印res和相乘的結(jié)果
[array([0.77854679, 0.62758656]), array([ 0.62758656, -0.77854679])]
-3.885780586188048e-16
從兩個(gè)向量和他們的乘積來(lái)看，乘積近乎為0，可以認(rèn)為是垂直的。
主成分分析和線性回歸的區(qū)別
這張圖給出的是樣本映射到方差最大的軸上的過(guò)程。很像線性回歸中找一條擬合各個(gè)樣本點(diǎn)的直線的圖像。
這里有幾點(diǎn)區(qū)別：
主成分分析中，坐標(biāo)軸表示的是各個(gè)特征，而線性回歸中，以一個(gè)特征的線性回歸為例
y軸對(duì)應(yīng)的是訓(xùn)練樣本的輸出標(biāo)記。樣本是垂直與特征所在的坐標(biāo)軸，測(cè)量的是實(shí)際值和預(yù)測(cè)值的差距，我們的任務(wù)是使得這個(gè)距離盡可能的小。
主成分分析中坐標(biāo)軸全部對(duì)應(yīng)特征，映射的方向是垂直于方差最大的軸而非特征所在的坐標(biāo)軸，且任務(wù)是使得樣本間距盡可能大。
高維數(shù)據(jù)降維
主成分分析的作用就是選出能使樣本方差最大的維度，選擇完維度之后，進(jìn)入對(duì)數(shù)據(jù)降維的操作。將高維數(shù)據(jù)映射為低維數(shù)據(jù)。
假設(shè)經(jīng)過(guò)主成分分析之后，左側(cè)X還是數(shù)據(jù)樣本，一個(gè)m * n的矩陣，右側(cè)是經(jīng)過(guò)k輪分析之后所獲得的k個(gè)主成分向量，形成一個(gè)k * n的矩陣。
自然想到使用線性代數(shù)將兩個(gè)矩陣相乘，得到m * k的矩陣，將原本的n個(gè)特征降為k個(gè)特征，即達(dá)到了降維映射的目標(biāo)。因此得到高維數(shù)據(jù)映射到低維數(shù)據(jù)的公式：
經(jīng)過(guò)向低維映射得到新的m * k的矩陣Xk如圖。低維數(shù)據(jù)也可以通過(guò)與主成分向量w相乘恢復(fù)成高維數(shù)據(jù)。Xk是m * k的矩陣與k * n的矩陣w相乘正好可以得到一個(gè)m * n的矩陣。但這個(gè)矩陣與原矩陣不一樣！??！
低維數(shù)據(jù)映射回高維數(shù)據(jù)的公式：
代碼實(shí)現(xiàn)PCA完整流程
import numpy as np
import matplotlib.pyplot as plt
class PCA:
def __init__(self, n_component):
assert n_component >= 1, 'n_component is invalidate'
self.n_component = n_component
self.components_ = None
def __repr__(self):
return 'PCA(n_component=%d)' % self.n_component
def fit(self,X, eta, n_iter=1e4, epsilon=0.0001):
'''
主成分分析
:param X:
:param eta:
:param n_iter:
:param epsilon:
:return:
'''
assert X.shape[1] >= self.n_component, 'X is invalidate'
'''均值歸零'''
def demean(X):
return X - np.mean(X, axis=0)
'''方差函數(shù)'''
def f(w, X):
return np.sum(X.dot(w)**2) / len(X)
'''方差函數(shù)導(dǎo)數(shù)'''
def df_ascent(w, X):
return X.T.dot(X.dot(w)) * 2 / len(X)
'''將向量化簡(jiǎn)為單位向量'''
def direction(w):
return w / np.linalg.norm(w)
'''尋找第一主成分'''
def first_component(w, X, eta, n_iter=1e4, epsilon=0.0001):
i_iter = 0
while i_iter < n_iter:
last_w = w
gradient = df_ascent(w, X)
w += eta * gradient
w = direction(w)
if abs(f(w, X) - f(last_w, X)) < epsilon:
break
i_iter += 1
return w
self.components_ = np.empty(shape=(self.n_component, X.shape[1]))
X = demean(X)
for i in range(self.n_component):
w = np.random.random(X.shape[1])
w = first_component(w, X, eta, n_iter, epsilon)
X = X - (X.dot(w)).reshape(-1, 1) * w
self.components_[i, :] = w
return self
def transform(self, X):
'''
將X映射到各個(gè)主成分中
:param X:
:return:
'''
assert X.shape[1] == self.components_.shape[1]
return X.dot(self.components_.T)
def inverse_transform(self, X):
'''
將低維數(shù)據(jù)轉(zhuǎn)回高維
:param X:
:return:
'''
assert X.shape[1] == self.components_.shape[0]
return X.dot(self.components_)
定義一個(gè)PCA類，有兩個(gè)屬性n_comnponent記錄主成分個(gè)數(shù)，components_記錄各個(gè)主成分向量。fit()函數(shù)實(shí)際就是通過(guò)梯度上升法選去n_component個(gè)主成分，并通過(guò)transform()函數(shù)從高維向低維映射，inverse_transform()函數(shù)實(shí)現(xiàn)了低維映射回高維。
使用測(cè)試數(shù)據(jù)進(jìn)行測(cè)試
'''測(cè)試數(shù)據(jù)'''
X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0, 10, size=100)
eta = 0.01
pca = PCA(1)
pca.fit(X, eta)
X_new = pca.transform(X)
print(pca.components_)
print(X_new.shape)
X_inverse = pca.inverse_transform(X_new)
plt.scatter(X[:, 0], X[:, 1], color='b', alpha=0.5)
plt.scatter(X_inverse[:, 0], X_inverse[:, 1], color='r', alpha=0.5)
plt.show()
得到主成分向量和降維后的數(shù)據(jù)集的大小
[[0.76233998 0.64717676]]
(100, 1)
同時(shí)繪制一下原數(shù)據(jù)集和低維數(shù)據(jù)映射回高維數(shù)據(jù)后的新數(shù)據(jù)集X_inverse，進(jìn)行對(duì)比
可見(jiàn)，映射回高維的數(shù)據(jù)和原始數(shù)據(jù)不一樣，說(shuō)明PCA降維會(huì)丟失信息，低維的數(shù)據(jù)不能恢復(fù)成原先高維的數(shù)據(jù)。
sklearn中的PCA實(shí)現(xiàn)
sklearn的PCA在sklearn.decomposition的PCA類中
from sklearn.decomposition import PCA
import numpy as np
import matplotlib.pyplot as plt
import sklearn.datasets as dataset
from sklearn.model_selection import train_test_split
from sklearn.neighbors import KNeighborsClassifier
import time
'''模擬數(shù)據(jù)集'''
X = np.empty((100, 2))
X[:, 0] = np.random.uniform(0., 100., size=100)
X[:, 1] = 0.75 * X[:, 0] + 3. + np.random.normal(0, 10, size=100)
pca = PCA(n_components=1)
pca.fit(X)
print(pca.components_)
X_new = pca.transform(X)
print(X_new.shape)
X_inverse = pca.inverse_transform(X_new)
print(X_inverse.shape)
plt.scatter(X[:, 0], X[:, 1], color='b', alpha=0.5)
plt.scatter(X_inverse[:, 0], X_inverse[:, 1], color='r', alpha=0.5)
plt.show()
初始化一個(gè)PCA對(duì)象，主成分?jǐn)?shù)為1，經(jīng)過(guò)主成分分析后的主成分向量打印
[[0.75763785 0.65267518]]
和之前自定義的PCA類得到的結(jié)果近似。也可以通過(guò)transform()函數(shù)得到降維數(shù)據(jù)，也可以打印inverse_transform()函數(shù)得到映射回高維的數(shù)據(jù)。這里不給出了。
接下來(lái)使用真實(shí)數(shù)據(jù)測(cè)試PCA降維對(duì)性能的影響
'''真實(shí)數(shù)據(jù)'''
data = dataset.load_digits()
X = data.data
y = data.target
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=666)
knn_clf = KNeighborsClassifier()
start = time.time()
knn_clf.fit(X_train, y_train)
end = time.time()
print("without PCA the KNN classifier's time cost is %s" % (end - start))
print("without PCA the KNN classifier's test score is %s" % knn_clf.score(X_test, y_test))
以KNN分類器為例，使用digits數(shù)據(jù)集，在不做降維的情況下進(jìn)行分類的耗時(shí)和準(zhǔn)確度為
without PCA the KNN classifier's time cost is 0.00571894645690918
without PCA the KNN classifier's test score is 0.9866666666666667
下面測(cè)試選2個(gè)主成分進(jìn)行降維的情況
knn_clf = KNeighborsClassifier()
pca = PCA(n_components=2)
pca.fit(X_train)
X_train_new = pca.transform(X_train)
X_test_new = pca.transform(X_test)
start = time.time()
knn_clf.fit(X_train_new, y_train)
end = time.time()
print("with PCA the KNN classifier's time cost is %s" % (end - start))
print("with PCA the KNN classifier's test score is %s" % knn_clf.score(X_test_new, y_test))
得到結(jié)果
with PCA the KNN classifier's time cost is 0.0010619163513183594
with PCA the KNN classifier's test score is 0.6066666666666667
可見(jiàn)時(shí)間上降維的數(shù)據(jù)訓(xùn)練明顯快了，但準(zhǔn)確度降低了。這說(shuō)明只選擇2個(gè)主成分效果并不理想。此時(shí)可以用到PCA類中的explained_variance_ratio_屬性。
print(pca.explained_variance_ratio_)
得到n_component=2時(shí)兩個(gè)主成分對(duì)每個(gè)樣本在主成分方向上的解釋度。
[0.14566817 0.13735469]
第一個(gè)主成分能解釋14%樣本的方差，第二個(gè)主成分能解釋13%樣本的方差。即當(dāng)前主成分量能夠維持原數(shù)據(jù)集方差的百分比，經(jīng)過(guò)2次主成分分析之后，新的數(shù)據(jù)集一共可以覆蓋原始數(shù)據(jù)28%的方差。解釋度越高，預(yù)測(cè)越準(zhǔn)確，解釋度之和相加的結(jié)果不會(huì)超過(guò)1，等于1是最理想的狀態(tài)及可以100%預(yù)測(cè)準(zhǔn)確。相當(dāng)于n_component=2時(shí)，72%的原始信息丟失了，因此預(yù)測(cè)準(zhǔn)確度不高。
如果對(duì)原始數(shù)據(jù)有多少個(gè)屬性求多少個(gè)主成分向量會(huì)得到一下結(jié)果。
knn_clf = KNeighborsClassifier()
pca = PCA(n_components=X_train.shape[1])
start = time.time()
pca.fit(X_train)
X_train_new = pca.transform(X_train)
X_test_new = pca.transform(X_test)
end = time.time()
knn_clf.fit(X_train_new, y_train)
print("with 64 component PCA the KNN classifier's time cost is %s" % (end - start))
print("with 64 component PCA the KNN classifier's test score is %s" % knn_clf.score(X_test_new, y_test))
print(pca.explained_variance_ratio_)
輸出
with 64 component PCA the KNN classifier's time cost is 0.007581233978271484
with 64 component PCA the KNN classifier's test score is 0.9866666666666667
[1.45668166e-01 1.37354688e-01 1.17777287e-01 8.49968861e-02
5.86018996e-02 5.11542945e-02 4.26605279e-02 3.60119663e-02
3.41105814e-02 3.05407804e-02 2.42337671e-02 2.28700570e-02
1.80304649e-02 1.79346003e-02 1.45798298e-02 1.42044841e-02
1.29961033e-02 1.26617002e-02 1.01728635e-02 9.09314698e-03
8.85220461e-03 7.73828332e-03 7.60516219e-03 7.11864860e-03
6.85977267e-03 5.76411920e-03 5.71688020e-03 5.08255707e-03
4.89020776e-03 4.34888085e-03 3.72917505e-03 3.57755036e-03
3.26989470e-03 3.14917937e-03 3.09269839e-03 2.87619649e-03
2.50362666e-03 2.25417403e-03 2.20030857e-03 1.98028746e-03
1.88195578e-03 1.52769283e-03 1.42823692e-03 1.38003340e-03
1.17572392e-03 1.07377463e-03 9.55152460e-04 9.00017642e-04
5.79162563e-04 3.82793717e-04 2.38328586e-04 8.40132221e-05
5.60545588e-05 5.48538930e-05 1.08077650e-05 4.01354717e-06
1.23186515e-06 1.05783059e-06 6.06659094e-07 5.86686040e-07
9.18612290e-34 9.18612290e-34 9.18612290e-34 8.82949950e-34]
訓(xùn)練數(shù)據(jù)集是一個(gè)含有64個(gè)特征的數(shù)據(jù)集，經(jīng)過(guò)主成分分析得到的explained_variance_ratio_也是一個(gè)含有64個(gè)元素的數(shù)據(jù)集，表示每個(gè)主成分對(duì)原始數(shù)據(jù)方差的解釋度，是一個(gè)從大到小的排列。最后幾個(gè)的解釋度幾乎為零，也就是說(shuō)對(duì)原始方差基本沒(méi)有作用。但這種方式雖然耗時(shí)增加了，但分類的準(zhǔn)確度達(dá)到98%。在實(shí)際情況下，可能會(huì)忽略對(duì)原始方差影響小的成分，在時(shí)間和準(zhǔn)確度之間做一個(gè)權(quán)衡。
例如要求降維達(dá)到解釋度0.95即可。
knn_clf = KNeighborsClassifier()
pca = PCA(0.95)
start = time.time()
pca.fit(X_train)
X_train_new = pca.transform(X_train)
X_test_new = pca.transform(X_test)
end = time.time()
knn_clf.fit(X_train_new, y_train)
print("with 0.95 PCA the KNN classifier's time cost is %s" % (end - start))
print("with 0.95 component PCA the KNN classifier's test score is %s" % knn_clf.score(X_test_new, y_test))
print(pca.explained_variance_ratio_.shape)
得到
with 0.95 PCA the KNN classifier's time cost is 0.007498979568481445
with 0.95 component PCA the KNN classifier's test score is 0.98
(28,)
可見(jiàn)，當(dāng)解釋度0.95時(shí)，耗時(shí)縮短了一些，分類準(zhǔn)確度還是98%，此時(shí)的主成分向量數(shù)為28，省去了其他36個(gè)特征維度。
n_component低可能造成信息丟失嚴(yán)重，但并非沒(méi)有用。例如還是n_component=2進(jìn)行降維。
pca = PCA(n_components=2)
pca.fit(X_train)
X_train_new = pca.transform(X_train)
for i in range(10):
plt.scatter(X_train_new[y_train == i, 0], X_train_new[y_train == i, 1], alpha=0.8)
plt.show()
對(duì)降維后的數(shù)據(jù)進(jìn)行分類別繪制得到
在二維空間中，類別見(jiàn)的區(qū)分度已經(jīng)相對(duì)明顯，如果僅僅是對(duì)藍(lán)色和紅色兩種類別進(jìn)行分析，則二維可能已經(jīng)足夠了。
PCA的其他用途
用途一：數(shù)據(jù)降噪
文章開(kāi)頭有提到PCA可以降噪，噪聲在生產(chǎn)過(guò)程中可能因?yàn)楦鞣N原因出現(xiàn)，影響數(shù)據(jù)的準(zhǔn)確性，PCA通過(guò)選取主成分將原有數(shù)據(jù)映射到低維數(shù)據(jù)再映射回高維數(shù)據(jù)的方式進(jìn)行一定程度的降噪，以digits數(shù)據(jù)為例。
import numpy as np
import matplotlib.pyplot as plt
import sklearn.datasets as dataset
from sklearn.decomposition import PCA
'''加載數(shù)據(jù)'''
data = dataset.load_digits()
X = data.data
y = data.target
#有噪聲的數(shù)據(jù)
noisy_digits = X + np.random.normal(0, 4, size=X.shape)
#取100個(gè)樣例數(shù)據(jù)
example_digits = noisy_digits[y == 0, :][:10]
for i in range(1, 10):
X_num = noisy_digits[y == i, :][:10]
example_digits = np.vstack([example_digits, X_num])
'''繪圖函數(shù)'''
def plot_digits(data):
fig, axes = plt.subplots(10, 10, figsize=(10, 10),
subplot_kw={'xticks':[], 'yticks':[]},
gridspec_kw=dict(hspace=0.1, wspace=0.1))
for i, ax in enumerate(axes.flat):
ax.imshow(data[i].reshape(8, 8), cmap='binary', interpolation='nearest', clim=(0, 16))
'''繪制帶有噪聲的數(shù)據(jù)'''
plot_digits(example_digits)
plt.show()
在digits數(shù)據(jù)特征集中人工加入噪聲，每類樣本選去10個(gè)進(jìn)行打印，打印的效果如
此時(shí)通過(guò)PCA進(jìn)行主成分分析，獲得映射后的低維數(shù)據(jù)，在映射回高維數(shù)據(jù)并以相同的方式打印
'''PCA數(shù)據(jù)降噪'''
#這里選用n_component=0.5進(jìn)行主成分分析
pca = PCA(0.5)
pca.fit(example_digits)
#映射到低維數(shù)據(jù)example_digits_new
example_digits_new = pca.transform(example_digits)
#再?gòu)牡途S數(shù)據(jù)映射回高維數(shù)據(jù)并繪制
example_digits_no_noise = pca.inverse_transform(example_digits_new)
plot_digits(example_digits_no_noise)
plt.show()
得到的圖像
明顯的好于原始數(shù)據(jù)圖像，降噪有一定效果。
用途二：人臉識(shí)別與特征臉
假設(shè)原始的人像數(shù)據(jù)是一個(gè)m * n的矩陣，在經(jīng)過(guò)主成分分析之后，一個(gè)k * n的主成分矩陣。如果將主成分矩陣也看成是由k個(gè)樣本組成的矩陣，那么可以理解為第一個(gè)樣本是最重要的樣本，第二個(gè)次之，依次往后。原始數(shù)據(jù)矩陣可以視為有m個(gè)人臉樣本的集合，如果將主成分矩陣每一行也看做是一個(gè)樣本的話，每一行相當(dāng)于是一個(gè)由原始人臉數(shù)據(jù)矩陣經(jīng)過(guò)主成分分析得到的特征臉矩陣，這個(gè)矩陣含有k個(gè)特征臉。每個(gè)特征臉表達(dá)了原有樣本中人臉的部分特征。
使用sklearn的fetch_lfw_people數(shù)據(jù)集，選去36張人臉進(jìn)行測(cè)試。
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import fetch_lfw_people
from sklearn.decomposition import PCA
import ssl
ssl._create_default_https_context = ssl._create_unverified_context
'''加載數(shù)據(jù)'''
face = fetch_lfw_people()
#face數(shù)據(jù)集的特征集是13233 * 2914的矩陣，打亂順序，取其中36個(gè)進(jìn)行測(cè)試
random_index = np.random.permutation(len(face.data))
X = face.data[random_index]
example_face = X[:36, :]
#繪制原始的人臉
'''繪圖函數(shù)'''
def plot_faces(data):
fig, axes = plt.subplots(10, 10, figsize=(10, 10),
subplot_kw={'xticks':[], 'yticks':[]},
gridspec_kw=dict(hspace=0.1, wspace=0.1))
for i, ax in enumerate(axes.flat):
ax.imshow(data[i].reshape(62, 47), cmap='binary', interpolation='nearest', clim=(0, 16))
plot_faces(example_face)
plt.show()
#數(shù)據(jù)集大，使用隨機(jī)的方式求解pca，提升效率
pca = PCA(svd_solver='randomized')
pca.fit(example_face)
agent_faces = pca.components_[:36, :]
plot_faces(agent_faces)
plt.show()
選取36張人臉繪制原始圖像，再根據(jù)這些數(shù)據(jù)進(jìn)行PCA操作，再繪制特征臉。
特征臉
