Facebook加入
LINE加入
2021/12/25
本文用了一個經典的例子，從數據探索，模型假設，模型訓練，模型可視化，step by step 讓讀者體驗機器學習完整的流程。
導語在模式分類和機器學習實踐中，線性判別分析（Linear Discriminant Analysis, LDA）方法常被用於數據預處理中的降維（dimensionality reduction）步驟。LDA在保證良好的類別區分度的前提下，將數據集向更低維空間投影，以求在避免過擬合（「維數災難」）的同時，減小計算消耗。
Ronald A. Fisher 在1936年（The Use of Multiple Measurements in Taxonomic Problems）提出了線性判別（Linear Discriminant）方法，它有時也用來解決分類問題。最初的線性判別適用於二分類問題，後來在1948年，被C. R. Rao（The utilization of multiple measurements in problems of biological classification）推廣為「多類線性判別分析」或「多重判別分析」。
一般意義上的LDA與主成分分析（Principal Component Analysis，PCA）十分相似，但不同於PCA尋找使全部數據方差最大化的坐標軸成分，LDA關心的是能夠最大化類間區分度的坐標軸成分。
更多關於PCA的內容，可參考 Implementing a Principal Component Analysis (PCA) in Python step by step。
一言蔽之，LDA將特徵空間（數據集中的多維樣本）投影到一個維度更小的 k 維子空間中（ k≤n−1），同時保持區分類別的信息。 通常情況下，降維不僅降低了分類任務的計算量，還能減小參數估計的誤差，從而避免過擬合。
主成分分析 vs. 線性判別分析線性判別分析（LDA）與主成分分析（PCA）都是常用的線性變換降維方法。PCA是一種「無監督」演算法，它不關心類別標籤，而是致力於尋找數據集中最大化方差的那些方向（亦即「主成分」）；LDA則是「有監督」的，它計算的是另一類特定的方向（或稱線性判別「器」）——這些方向刻畫了最大化類間區分度的坐標軸。
雖然直覺上聽起來，在已知類別信息時，LDA對於多分類任務要優於PCA，但實際並不一定。
比如，比較經過PCA或LDA處理后圖像識別任務的分類準確率，如果樣本數比較少，PCA是要優於LDA的（PCA vs. LDA, A.M. Martinez et al., 2001）。聯合使用LDA和PCA也並不罕見，例如降維時先用PCA再做LDA。
假定我們的目標是將一個 d 維數據集投影到一個 k （ k
在後面，我們會計算數據集的本徵向量（成分），將其歸總到一個所謂的「散布矩陣」（類間散布矩陣和類內散布矩陣）。每一個本徵向量對應一個本徵值，本徵值會告訴我們相應本徵向量的「長度」/「大小」。
如果所有的本徵值大小都很相近，那麼這就表示我們的數據已經投影到了一個「好」的特徵空間上。
而如果一部分本徵值遠大於其他的本徵值，我們可能就考慮只留下最大的那些，因為它們包含了更多關於數據的分佈的信息。反之，接近 0 的本徵值含有的信息量更少，我們就考慮在構造新特徵空間時捨棄它們。
LDA的五個步驟下面列出了執行線性判別分析的五個基本步驟。我們會在後面做更詳細的講解。
1. 計算數據集中不同類別數據的 d 維均值向量。
2. 計算散布矩陣，包括類間、類內散布矩陣。
3. 計算散布矩陣的本徵向量 e1,e2,...,ed 和對應的本徵值 λ1,λ2,...,λd。
4. 將本徵向量按本徵值大小降序排列，然後選擇前 k。 個最大本徵值對應的本徵向量，組建一個 d×k 維矩陣——即每一列就是一個本徵向量。
5. 用這個 d×k-維本徵向量矩陣將樣本變換到新的子空間。這一步可以寫作矩陣乘法 Y=X×W 。 X 是 n×d 維矩陣，表示 n 個樣本； y 是變換到子空間后的 n×k 維樣本。
準備樣本數據集鳶尾花集
接下來要在著名的「鳶尾花」數據集上折騰一番。鳶尾花集可從UCI的機器學習目錄中下載：
引用 ：Bache, K. & Lichman, M. (2013). UCI Machine Learning Repository. Irvine, CA: University of California, School of Information and Computer Science.
該集對三類、一百五十朵鳶尾花做了測量。
三類包括：
1. 山鳶尾（n=50）
2. 變色鳶尾（n=50）
3. 維吉尼亞鳶尾（n=50）
四種特徵包括：
1. 萼片長度（厘米）
2. 萼片寬度（厘米）
3. 花瓣長度（厘米）
4. 花瓣寬度（厘米）
feature_dict = {i:label for i,label in zip( range(4), ('sepal length in cm', 'sepal width in cm', 'petal length in cm', 'petal width in cm', ))} 讀取數據集import pandas as pd df = pd.io.parsers.read_csv( filepath_or_buffer='https://archive.ics.uci.edu/ml/machine-learning-databases/iris/iris.data', header=None, sep=',', ) df.columns = [l for i,l in sorted(feature_dict.items)] + ['class label'] df.dropna(how="all", inplace=True) # to drop the empty line at file-end df.tail
因為數值數據更容易處理，我們使用 scikit-learn 庫中的 LabelEncode 將類別標籤轉換成數字 1，2，3：
from sklearn.preprocessing import LabelEncoder X = df[[0,1,2,3]].values y = df['class label'].values enc = LabelEncoder label_encoder = enc.fit(y) y = label_encoder.transform(y) + 1 label_dict = {1: 'Setosa', 2: 'Versicolor', 3:'Virginica'}
直方圖和特徵選擇把四種特徵的分佈情況在一維直方圖中分別做可視化呈現，可以得到一個對 w1,w2,w3 三類數據的粗略的觀察：
%matplotlib inline from matplotlib import pyplot as plt import numpy as np import math fig, axes = plt.subplots(nrows=2, ncols=2, figsize=(12,6)) for ax,cnt in zip(axes.ravel, range(4)): # set bin sizes min_b = math.floor(np.min(X[:,cnt])) max_b = math.ceil(np.max(X[:,cnt])) bins = np.linspace(min_b, max_b, 25) # plottling the histograms for lab,col in zip(range(1,4), ('blue', 'red', 'green')): ax.hist(X[y==lab, cnt], color=col, label='class %s' %label_dict[lab], bins=bins, alpha=0.5,) ylims = ax.get_ylim # plot annotation leg = ax.legend(loc='upper right', fancybox=True, fontsize=8) leg.get_frame.set_alpha(0.5) ax.set_ylim([0, max(ylims)+2]) ax.set_xlabel(feature_dict[cnt]) ax.set_title('Iris histogram #%s' %str(cnt+1)) # hide axis ticks ax.tick_params(axis="both", which="both", bottom="off", top="off", labelbottom="on", left="off", right="off", labelleft="on") # remove axis spines ax.spines["top"].set_visible(False) ax.spines["right"].set_visible(False) ax.spines["bottom"].set_visible(False) ax.spines["left"].set_visible(False) axes[0][0].set_ylabel('count') axes[1][0].set_ylabel('count') fig.tight_layout plt.show
僅憑這些簡單的圖形化展示，已經足以讓我們得出結論：在四種特徵裡面，花瓣的長度、寬度更適合用來區分三種鳶尾花類別。
實際應用中，比起通過投影降維（此處即LDA），另一種比較好的辦法是做特徵篩選。像鳶尾花這樣的低維數據集，看一眼直方圖就能得到很多信息了。
另一種簡單但有效的方式是使用特徵選擇演算法——如果你有興趣，我在此處詳述了「序列特徵選擇」； scikit-learn 有另一種漂亮的手段。
我還寫過一篇「特徵選擇中的濾波器，封裝器和嵌入方法」，是在更高的層次上對不同方法的總結。
規範性假設需要指出，LDA假設數據服從正態分佈、不同特徵之間互相統計獨立且各類數據的協方差矩陣相等。不過這僅僅指LDA用作分類器的情況，當LDA用於降維時，哪怕數據不符合這些假設，LDA通常也能取得不錯的效果。
而即便對於分類任務，LDA對數據的分佈似乎也是相當魯棒的：
「儘管實際情況常常不符合『不同類別數據間有相同協方差矩陣』和規範性假設，線性判別分析在人臉和物體識別任務中也通常能夠得到不錯的結果。」。
Tao Li, Shenghuo Zhu, and Mitsunori Ogihara. 「Using Discriminant Analysis for Multi-Class Classification: An Experimental Investigation.」 Knowledge and Information Systems 10, no. 4 (2006): 453–72.)
Duda, Richard O, Peter E Hart, and David G Stork. 2001. Pattern Classification. New York: Wiley.
五步實現LDA完成以上幾項準備工作后，我們就可以實際運行LDA了。
第一步：計算數據的 d 維均值向量
首先做一個簡單的計算：分別求三種鳶尾花數據在不同特徵維度上的均值向量 mi：
np.set_printoptions(precision=4) mean_vectors = for cl in range(1,4): mean_vectors.append(np.mean(X[y==cl], axis=0)) print('Mean Vector class %s: %s\n' %(cl, mean_vectors[cl-1]))得到
Mean Vector class 1: [ 5.006 3.418 1.464 0.244] Mean Vector class 2: [ 5.936 2.77 4.26 1.326] Mean Vector class 3: [ 6.588 2.974 5.552 2.026]第二步：計算散布矩陣
計算兩個 4×4 維矩陣：類內散布矩陣和類間散布矩陣。
2.1a 類內散布矩陣 Sw
類內散布矩陣 Sw 由該式算出：
S_W = np.zeros((4,4)) for cl,mv in zip(range(1,4), mean_vectors): class_sc_mat = np.zeros((4,4)) # scatter matrix for every class for row in X[y == cl]: row, mv = row.reshape(4,1), mv.reshape(4,1) # make column vectors class_sc_mat += (row-mv).dot((row-mv).T) S_W += class_sc_mat # sum class scatter matrices print('within-class Scatter Matrix:\n', S_W)得到
within-class Scatter Matrix: [[ 38.9562 13.683 24.614 5.6556] [ 13.683 17.035 8.12 4.9132] [ 24.614 8.12 27.22 6.2536] [ 5.6556 4.9132 6.2536 6.1756]]2.1b 貝塞爾糾偏項的影響
在計算類-協方差矩陣時，可以在類內散布矩陣上添加尺度因數 1N−1，這樣計算式就變為
Ni 是類樣本大小（此處為50），因為三個類別的樣本數是相同的，這裡可以捨去 (Ni−1)。
因為本徵向量是相同的，只是本徵值有一個常數項的尺度變化，所以即便將其忽略不計，最後得到的特徵空間也不會改變（這一點在文末還有體現）。
2.2 類間散布矩陣 SB
mm 是全局均值，而 mmi 和 Ni 是每類樣本的均值和樣本數。
overall_mean = np.mean(X, axis=0) S_B = np.zeros((4,4)) for i,mean_vec in enumerate(mean_vectors): n = X[y==i+1,:].shape[0] mean_vec = mean_vec.reshape(4,1) # make column vector overall_mean = overall_mean.reshape(4,1) # make column vector S_B += n * (mean_vec - overall_mean).dot((mean_vec - overall_mean).T) print('between-class Scatter Matrix:\n', S_B)得到
between-class Scatter Matrix: [[ 63.2121 -19.534 165.1647 71.3631] [ -19.534 10.9776 -56.0552 -22.4924] [ 165.1647 -56.0552 436.6437 186.9081] [ 71.3631 -22.4924 186.9081 80.6041]]第三步：求解矩陣 S−1WSB 的廣義本徵值
接下來求解矩陣 S−1WSB 的廣義本徵值問題，從而得到線性判別「器」：
eig_vals, eig_vecs = np.linalg.eig(np.linalg.inv(S_W).dot(S_B)) for i in range(len(eig_vals)): eigvec_sc = eig_vecs[:,i].reshape(4,1) print('\nEigenvector {}: \n{}'.format(i+1, eigvec_sc.real)) print('Eigenvalue {:}: {:.2e}'.format(i+1, eig_vals[i].real))得到
Eigenvector 1: [[-0.2049] [-0.3871] [ 0.5465] [ 0.7138]] Eigenvalue 1: 3.23e+01 Eigenvector 2: [[-0.009 ] [-0.589 ] [ 0.2543] [-0.767 ]] Eigenvalue 2: 2.78e-01 Eigenvector 3: [[ 0.179 ] [-0.3178] [-0.3658] [ 0.6011]] Eigenvalue 3: -4.02e-17 Eigenvector 4: [[ 0.179 ] [-0.3178] [-0.3658] [ 0.6011]] Eigenvalue 4: -4.02e-17計算結果應當作何解釋呢？從第一節線性代數課開始我們就知道，本徵向量和本徵值表示了一個線性變換的形變程度——本徵向量是形變的方向；本徵值是本徵向量的尺度因數，刻畫了形變的幅度。
如果將LDA用於降維，本徵向量非常重要，因為它們將會組成新特徵子空間的坐標軸。對應的本徵值表示了這些新坐標軸的信息量的多少。
再檢查一遍計算過程，然後對本徵值做進一步討論。
檢查本徵向量-本徵值的運算
要檢查該運算是否正確，看其是否滿足等式：
for i in range(len(eig_vals)): eigv = eig_vecs[:,i].reshape(4,1) np.testing.assert_array_almost_equal(np.linalg.inv(S_W).dot(S_B).dot(eigv), eig_vals[i] * eigv, decimal=6, err_msg='', verbose=True) print('ok')輸出：
ok第四步：選擇線性判別「器」構成新的特徵子空間
4.1 按本徵值降序排列本徵向量
本徵向量只是定義了新坐標軸的方向，它們的大小都是單位長度1。
想要構成更低維的子空間，就得選擇丟棄哪個（些）本徵向量，所以我們得看看對應本徵向量的那些本徵值。 大體上說，對於數據的分佈情況，本徵值最小的那些本徵向量所承載的信息最少，所以要捨棄的就是這些本徵向量。
將本徵向量根據本徵值的大小從高到低排序，然後選擇前 k 個本徵向量：
eig_pairs = [(np.abs(eig_vals[i]), eig_vecs[:,i]) for i in range(len(eig_vals))] # Sort the (eigenvalue, eigenvector) tuples from high to low eig_pairs = sorted(eig_pairs, key=lambda k: k[0], reverse=True) # Visually confirm that the list is correctly sorted by decreasing eigenvalues print('Eigenvalues in decreasing order:\n') for i in eig_pairs: print(i[0])得到
Eigenvalues in decreasing order: 32.2719577997 0.27756686384 5.71450476746e-15 5.71450476746e-15可以看到有兩個本徵值接近0，不是因為它們不含信息量，而是浮點數精度的關係。其實，這后兩個本徵值應該恰好為0。
在LDA中，線性判別器的數目最多是 c−1，c 是總的類別數，這是因為類內散布矩陣 SB 是 c 個秩為1或0的矩陣的和。
注意到很少有完全共線的情況（所有樣本點分佈在一條直線上），協方差矩陣秩為1，這導致了只有一個非零本徵值和一個對應的本徵向量。
現在，把「可釋方差」表達為百分數的形式：
print('Variance explained:\n') eigv_sum = sum(eig_vals) for i,j in enumerate(eig_pairs): print('eigenvalue {0:}: {1:.2%}'.format(i+1, (j[0]/eigv_sum).real))得到
Variance explained: eigenvalue 1: 99.15% eigenvalue 2: 0.85% eigenvalue 3: 0.00% eigenvalue 4: 0.00%第一個本徵對是信息量最大的一組，如果我們考慮建立一個一維的向量空間，使用該本徵對就不會丟失太多信息。
4.2 選擇 k 個最大本徵值對應的本徵向量
按本徵值的大小得到降序排列的本徵對之後，現在就可以組建我們的 k×d-維本徵向量矩陣 W 了（此時大小為 4×2），這樣就從最初的4維特徵空間降到了2維的特徵空間。
W = np.hstack((eig_pairs[0][1].reshape(4,1), eig_pairs[1][1].reshape(4,1))) print('Matrix W:\n', W.real)得到
Matrix W: [[-0.2049 -0.009 ] [-0.3871 -0.589 ] [ 0.5465 0.2543] [ 0.7138 -0.767 ]]第五步：將樣本變換到新的子空間中
使用上一步剛剛計算出的 4×2-維矩陣 W， 將樣本變換到新的特徵空間上：
X_lda = X.dot(W) assert X_lda.shape == (150,2), "The matrix is not 150x2 dimensional."得到
from matplotlib import pyplot as plt def plot_step_lda: ax = plt.subplot(111) for label,marker,color in zip( range(1,4),('^', 's', 'o'),('blue', 'red', 'green')): plt.scatter(x=X_lda[:,0].real[y == label], y=X_lda[:,1].real[y == label], marker=marker, color=color, alpha=0.5, label=label_dict[label] ) plt.xlabel('LD1') plt.ylabel('LD2') leg = plt.legend(loc='upper right', fancybox=True) leg.get_frame.set_alpha(0.5) plt.title('LDA: Iris projection onto the first 2 linear discriminants') # hide axis ticks plt.tick_params(axis="both", which="both", bottom="off", top="off", labelbottom="on", left="off", right="off", labelleft="on") # remove axis spines ax.spines["top"].set_visible(False) ax.spines["right"].set_visible(False) ax.spines["bottom"].set_visible(False) ax.spines["left"].set_visible(False) plt.grid plt.tight_layout plt.show plot_step_lda
上方散點圖展示了我們通過 LDA 構建的新的特徵子空間。可以看到第一個線性判別器「LD1」把不同類數據區分得不錯，第二個線性判別器就不行了。其原因在上面已經做了簡單解釋。
PCA 和 LDA 的對比為了與使用線性判別分析得到的特徵子空間作比較，我們將使用 scikit-learn 機器學習庫中的 PCA 類。
文檔看這裡：
http://scikit-learn.org/stable/modules/generated/sklearn.decomposition.PCA.html。
方便起見，我們直接在 n_components 參數上設定了希望從輸入數據集中得到的成分數量。
n_components : int, None or string Number of components to keep. if n_components is not set all components are kept: n_components == min(n_samples, n_features) if n_components == 『mle』, Minka』s MLE is used to guess the dimension if 0 < n_components < 1, select the number of components such that the amount of variance that needs to be explained is greater than the percentage specified by n_components但是在觀察兩種線性變換的結果之前，讓我們快速複習一下 PCA 和 LDA 的目標：PCA 在整個數據集中尋找方差最大的坐標軸，而 LDA 則尋找對於類別區分度最佳的坐標軸。
實際的降維步驟中，經常是LDA放在PCA之後使用。
from sklearn.decomposition import PCA as sklearnPCA sklearn_pca = sklearnPCA(n_components=2) X_pca = sklearn_pca.fit_transform(X) def plot_pca: ax = plt.subplot(111) for label,marker,color in zip( range(1,4),('^', 's', 'o'),('blue', 'red', 'green')): plt.scatter(x=X_pca[:,0][y == label], y=X_pca[:,1][y == label], marker=marker, color=color, alpha=0.5, label=label_dict[label] ) plt.xlabel('PC1') plt.ylabel('PC2') leg = plt.legend(loc='upper right', fancybox=True) leg.get_frame.set_alpha(0.5) plt.title('PCA: Iris projection onto the first 2 principal components') # hide axis ticks plt.tick_params(axis="both", which="both", bottom="off", top="off", labelbottom="on", left="off", right="off", labelleft="on") # remove axis spines ax.spines["top"].set_visible(False) ax.spines["right"].set_visible(False) ax.spines["bottom"].set_visible(False) ax.spines["left"].set_visible(False) plt.tight_layout plt.grid plt.show然後執行
plot_pca plot_step_lda
回想我們之前討論的內容：PCA 獲取數據集中方差最大的成分，而LDA給出不同類別間方差最大的坐標軸。
使用 scikit-learn 中的 LDA我們已經看到，線性判別分析是如何一步步實現的了。其實通過使用 scikit-learn
機器學習庫中的 LDA ，我們可以更方便地實現同樣的結果。
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as LDA # LDA sklearn_lda = LDA(n_components=2) X_lda_sklearn = sklearn_lda.fit_transform(X, y)作圖：
def plot_scikit_lda(X, title): ax = plt.subplot(111) for label,marker,color in zip( range(1,4),('^', 's', 'o'),('blue', 'red', 'green')): plt.scatter(x=X[:,0][y == label], y=X[:,1][y == label] * -1, # flip the figure marker=marker, color=color, alpha=0.5, label=label_dict[label]) plt.xlabel('LD1') plt.ylabel('LD2') leg = plt.legend(loc='upper right', fancybox=True) leg.get_frame.set_alpha(0.5) plt.title(title) # hide axis ticks plt.tick_params(axis="both", which="both", bottom="off", top="off", labelbottom="on", left="off", right="off", labelleft="on") # remove axis spines ax.spines["top"].set_visible(False) ax.spines["right"].set_visible(False) ax.spines["bottom"].set_visible(False) ax.spines["left"].set_visible(False) plt.grid plt.tight_layout plt.show執行
plot_step_lda plot_scikit_lda(X_lda_sklearn, title='Default LDA via scikit-learn')
關於規範化最近被問到了這個問題，所以我做一下說明：對特徵做尺度變換，比如【規範化】，不會改變LDA的最終結果，所以這是可選步驟。的確，散布矩陣會因為特徵的縮放而發生變化，而且本徵向量也因此而改變。但是關鍵的部分在於，最終本徵值不會變——唯一可見的不同只是成分軸的尺度大小。這在數學上是可以證明的（將來我可能在這兒插入一些公式）。
%matplotlib inline import pandas as pd import matplotlib.pyplot as plt import pandas as pd df = pd.read_csv('https://archive.ics.uci.edu/ml/machine-learning-databases/iris/iris.data', header=None) df[4] = df[4].map({'Iris-setosa':0, 'Iris-versicolor':1, 'Iris-virginica':2}) df.tail
sepal length in cmsepal width in cmpetal length in cmpedal width in cmclass label1456.73.05.22.3Iris-virginica1466.32.55.01.9Iris-virginaica1476.53.05.22.0Iris-virginaica1486.23.45.42.3Iris-virginaica1495.93.05.11.8Iris-virginaica讀取數據集，對 X 中的列做規範化。規範化就是把數據用均值做中心化、用標準差做單位化：
這樣所有的列就都是 0 均值（μxstd=0）、標準差為 1 的了（\sigma_{x_{std}}=1）。
y, X = df.iloc[:, 4].values, df.iloc[:, 0:4].values X_cent = X - X.mean(axis=0) X_std = X_cent / X.std(axis=0)接下來我簡單複製了一下LDA的各個步驟，這些之前我們都討論過了。為簡便都寫成了Python函數。
import numpy as np def comp_mean_vectors(X, y): class_labels = np.unique(y) n_classes = class_labels.shape[0] mean_vectors = for cl in class_labels: mean_vectors.append(np.mean(X[y==cl], axis=0)) return mean_vectors def scatter_within(X, y): class_labels = np.unique(y) n_classes = class_labels.shape[0] n_features = X.shape[1] mean_vectors = comp_mean_vectors(X, y) S_W = np.zeros((n_features, n_features)) for cl, mv in zip(class_labels, mean_vectors): class_sc_mat = np.zeros((n_features, n_features)) for row in X[y == cl]: row, mv = row.reshape(n_features, 1), mv.reshape(n_features, 1) class_sc_mat += (row-mv).dot((row-mv).T) S_W += class_sc_mat return S_W def scatter_between(X, y): overall_mean = np.mean(X, axis=0) n_features = X.shape[1] mean_vectors = comp_mean_vectors(X, y) S_B = np.zeros((n_features, n_features)) for i, mean_vec in enumerate(mean_vectors): n = X[y==i+1,:].shape[0] mean_vec = mean_vec.reshape(n_features, 1) overall_mean = overall_mean.reshape(n_features, 1) S_B += n * (mean_vec - overall_mean).dot((mean_vec - overall_mean).T) return S_B def get_components(eig_vals, eig_vecs, n_comp=2): n_features = X.shape[1] eig_pairs = [(np.abs(eig_vals[i]), eig_vecs[:,i]) for i in range(len(eig_vals))] eig_pairs = sorted(eig_pairs, key=lambda k: k[0], reverse=True) W = np.hstack([eig_pairs[i][1].reshape(4, 1) for i in range(0, n_comp)]) return W先把還沒有做尺度變換的數據的本徵值、本徵向量和投影矩陣列印出來：
S_W, S_B = scatter_within(X, y), scatter_between(X, y) eig_vals, eig_vecs = np.linalg.eig(np.linalg.inv(S_W).dot(S_B)) W = get_components(eig_vals, eig_vecs, n_comp=2) print('EigVals: %s\n\nEigVecs: %s' % (eig_vals, eig_vecs)) print('\nW: %s' % W)結果：
EigVals: [ 2.0905e+01 +0.0000e+00j 1.4283e-01 +0.0000e+00j -2.8680e-16 +1.9364e-15j -2.8680e-16 -1.9364e-15j] EigVecs: [[ 0.2067+0.j 0.0018+0.j 0.4846-0.4436j 0.4846+0.4436j] [ 0.4159+0.j -0.5626+0.j 0.0599+0.1958j 0.0599-0.1958j] [-0.5616+0.j 0.2232+0.j 0.1194+0.1929j 0.1194-0.1929j] [-0.6848+0.j -0.7960+0.j -0.6892+0.j -0.6892-0.j ]] W: [[ 0.2067+0.j 0.0018+0.j] [ 0.4159+0.j -0.5626+0.j] [-0.5616+0.j 0.2232+0.j] [-0.6848+0.j -0.7960+0.j]]作圖：
X_lda = X.dot(W) for label,marker,color in zip( np.unique(y),('^', 's', 'o'),('blue', 'red', 'green')): plt.scatter(X_lda[y==label, 0], X_lda[y==label, 1], color=color, marker=marker)
再對規範化之後的數據集重複上述操作：
S_W, S_B = scatter_within(X_std, y), scatter_between(X_std, y) eig_vals, eig_vecs = np.linalg.eig(np.linalg.inv(S_W).dot(S_B)) W_std = get_components(eig_vals, eig_vecs, n_comp=2) print('EigVals: %s\n\nEigVecs: %s' % (eig_vals, eig_vecs)) print('\nW: %s' % W_std)得到
EigVals: [ 2.0905e+01 1.4283e-01 -6.7207e-16 1.1082e-15] EigVecs: [[ 0.1492 -0.0019 0.8194 -0.3704] [ 0.1572 0.3193 -0.1382 -0.0884] [-0.8635 -0.5155 -0.5078 -0.5106] [-0.4554 0.7952 -0.2271 0.7709]] W: [[ 0.1492 -0.0019] [ 0.1572 0.3193] [-0.8635 -0.5155] [-0.4554 0.7952]]作圖：
X_std_lda = X_std.dot(W_std) X_std_lda[:, 0] = -X_std_lda[:, 0] X_std_lda[:, 1] = -X_std_lda[:, 1] for label,marker,color in zip( np.unique(y),('^', 's', 'o'),('blue', 'red', 'green')): plt.scatter(X_std_lda[y==label, 0], X_std_lda[y==label, 1], color=color, marker=marker)
可以看到，不管做不做尺度變化，本徵值是一樣的（要注意到，W 的秩是2，4維數據集中最小的兩個本徵值本來就應該是0）。並且，你也看到了，除了成分軸的尺度不太一樣，以及圖形做了中心對稱翻轉，最後的投影結果基本沒有區別。
注1：文中出現了線性代數術語「eigenvalue」「eigenvector」，中文教材對應有「特徵值」「本徵值」兩種常見譯法。為了與「feature」相區分，本文使用「本徵」翻譯。
注2：文中提到 「k×d-維本徵向量矩陣」，原文寫作 「k×d-dimensional eigenvector matrix」，指的是「k個d維的本徵向量組成的矩陣」。
中文與公式混合排版之後可能引發歧義，故在此做單獨解釋。
譯者 / 李宇琛 Asher Li
原文來源
關注 AI 研習社（okweiwu），回復 1 領取
【超過 1000G 神經網路／AI／大數據、教程、論文！】
▼▼▼
熱門推薦
留言回覆
