求數據的第一主成分
（在notebook中）
將包加載好，再創建出一個虛擬的測試用例，生成的X有兩個特征，特征一為0到100之間隨機分布，共一百個樣本，對於特征二，其和特征一有一個基本的線性關系（為什么要有一個基本的線性關系？是因為含有一個基本的線性關系，這樣對數據降維的效果會更加的明顯）
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])
圖像如下
具體操作
那么我們就開始進行數據的降維的具體操作
首先我們先進行demean的操作，我們設置一個函數為demean，傳來一個X（矩陣），返回均值歸0以后所得到的矩陣，具體計算可以使用X減去X的均值即可，需要注意的是，X代表一個矩陣，每一行代表一個樣本，那么操作實際上就是每一個樣本中的每一個特征都要減去這個特征所對應的均值，所以X減去的應該是個向量（每一個特征所對應的均值，也就是在行方向上去一個均值），其就是每一列的均值，這樣就可以得到均值歸0以后的樣本了
def demean(X):
return X - np.mean(X,axis=0)
那么我們測試一下這個函數
X_demean = demean(X)
然后繪制一下這個圖
plt.scatter(X_demean[:,0],X_demean[:,1])
圖像如下
可以看出來整體分布是個之前一樣的，不過從坐標軸可以看出來，（0，0）大概在中心位置，現在對於這些數據來說，在兩個方向上的均值都為0
那么我們可以驗證一下是不是這樣的
使用np.mean(X_demean[:,0])看一下對於第0維度來說均值是多少，結果如下
使用np.mean(X_demean[:,1])看一下對於第1維度來說均值是多少，結果如下
我們可以看出來，這兩個值都是很小的，十分接近0，可以說這兩個均值是為0的
使用梯度上升法來求解我們數據的主成分
那么做完均值歸零以后，我們就要開始使用梯度上升法來求解我們數據的主成分
求解過程和梯度下降法非常的像
首先我們先設置一下求解目標函數的函數f，傳上w和x，只要將X點乘上w，再進行平方和以后除以樣本數即可
def f(w,X):
return np.sum((X.dot(w)**2))/len(X)
我們將求解目標函數對應的梯度的方法設為df_mean（其實應該是叫df_math，打快了沒注意），其就是計算的公式的實現，即x的轉置點乘上x點乘w乘2再除以樣本數
def df_mean(w,X):
return X.T.dot(X.dot(w))*2. / len(X)
為了驗證我們的這個是正確的，使用這個df_debug這個函數，和線性下降法一樣，使兩個點之間連成的直線不斷的靠近應得的直線，使其斜率相當，注意的是，這里的epsilon取值比較小，是因為在PCA的梯度上升法中，w是一個方向向量，其模為1，所以w的每一個維度其實都很小，那么為了適應，相應的epsilon也要小一些
def df_debug(w,X,epsilon=0.0001):
res = np.empty(len(w))
for i in range(len(w)):
w_1 = w.copy()
w_1[i] += epsilon
w_2 = w.copy()
w_2[i] -= epsilon
res[i] = (f(w_1,X) - f(w_2,X)) / (2*epsilon)
return res
那么下面就開始梯度上升法的過程
每一次循環中先計算一下梯度，然后記錄一下上一次的w是多少，之后計算一下使用梯度上升法后得到的新的w，之后就對新的w和舊的w進行計算，看一下兩次對應的增加的值有沒有超過給定的值，還沒有超過就繼續循環，否則的話，我們就直接break回去，返回最終的w
需要注意的是，w本來是個方向向量，其模應該為1，從之前推導的公式的條件也可以看出來，需要w的模為1，但是在計算過程中，很可能會讓w的模不為1，雖然也可以使用，但是會使得整個計算的過程很不順暢，且還要將eta的值變得特得特別小，這樣會導致程序的效率非常的低
所以每次我們需要對w進行一下將其模變為1的操作，這里使用direction函數將其化成一個只表示方向的單位向量，只要返回w除以w對應的模（求向量的模可以使用np.linalg.norm(w)）就可以，這樣我們在使用的時候首先就要將initial_w轉換成單位向量，在每求到一次新的w之后相應的也調用一下這個函數來將模變成1，這樣就完成了使用梯度上升法求解PCA的代碼
def direction(w):
return w / np.linalg.norm(w)
def gradient_ascent(df,X,initial_w,eta,n_iters=1e4, epsilon=1e-8):
w = direction(initial_w)
cur_iter = 0
while cur_iter < n_iters:
gradient = df(w, X)
last_w = w
w = w + eta * gradient
w = direction(w)
if (abs(f(w,X) - f(last_w,X)) < epsilon):
break
cur_iter += 1
return w
我們想要調用這個函數，還需要注意，我們調用的時候是需要初始化一個w的值，但是這個初始化的值是不可以為0向量的，因為對於這個目標函數來說，w等於0本身就是一個極值點，只不過這個極值點是最小值點，其梯度值也為0，所以帶入0是不可以的，因此我們開始的時候需要隨機化一個向量，在向量中對應的隨機的元素個數就是X這個矩陣的列數，也就是樣本的特征數
initial_w = np.random.random(X.shape[1])
initial_w
結果如下
我們設置eta為0.001
eta=0.001
還有一個注意點，對於PCA問題，不能使用StandardScaler來標准化數據，這是為什么呢？因為對於PCA問題來說，本身就是求一個軸，使樣本映射以后的方差最大，但是一旦對樣本進行標准化了，將樣本的方差打散以后會使樣本的方差為1，這就使得方差的最大值不存在，就會求不出來PCA真正想要得到的結果了
我們先使用df_debug來求一下對應的軸
gradient_ascent(df_debug,X_demean,initial_w,eta)
結果如下
然后可以使用df_mean來求一下
gradient_ascent(df_mean,X_demean,initial_w,eta)
結果如下
可以看出來結果是一樣的，這就說明了我們的結果是正確的
然后我們可視化一下，看一下這個軸對應的方向是怎么個方向，先將原始的數據繪制出來，然后繪制軸對應的方向，方法是對於x和y，首先繪制（0，0）這個點，然后繪制一下w這個點，即w[0]和w[1]對應的點，但是由於w很小，正常顯示可能看不清，所以對其乘上一個系數，最后是其顏色為紅色
w = gradient_ascent(df_mean,X_demean,initial_w,eta)
plt.scatter(X_demean[:,0],X_demean[:,1])
plt.plot([0,w[0]*30],[0,w[1]*30],color='r')
圖像如下
相應的這個軸就是一個主成分，在這里使我們所求出來的第一個主成分，因此我們也稱其為第一主成分
極端情況
我們使用一下極端的結果來試一下，還是使用原來的虛假樣例，但是不添加噪音
X2 = np.empty((100,2))
X2[:,0] = np.random.uniform(0. ,100. , size=100)
X2[:,1] = 0.75 * X2[:,0] + 3.
繪制出圖像
plt.scatter(X2[:,0],X2[:,1])
圖像如下
可以發現其完全在一條直線上，不過依然是分布在一個二維的平面上
然后我們開始求解問題，首先進行demean操作，然后調用gradient_ascent函數進行計算並存在w2中
X2_demean = demean(X2)
w2 = gradient_ascent(df_mean,X2_demean,initial_w,eta)
繪制圖像
plt.scatter(X2_demean[:,0],X2_demean[:,1])
plt.plot([0,w2[0]*30],[0,w2[1]*30],color='r')
圖像如下
這個軸和這些點是重合的，且我們可以得出gradient_ascent(df_mean,X2_demean,initial_w,eta)為0趨近於0.8和趨近於0.6組成的，對應的我們的系數0.75
經過這兩個方面可知是算法正確的
