Python機器學習基礎篇三《無監督學習與預處理》

虛幻私塾發表於2020-12-26

前言

前期回顧: Python機器學習基礎篇二《為什麼用Python進行機器學習》
上面這篇裡面寫了文字和序列相關。

我們要討論的第二種機器學習演算法是無監督學習演算法。無監督學習包括沒有已知輸出、沒 有老師指導學習演算法的各種機器學習。在無監督學習中,學習演算法只有輸入資料,並需要 從這些資料中提取知識。

3.1 無監督學習的型別

本章將研究兩種型別的無監督學習:資料集變換與聚類。

資料集的無監督變換(unsupervised transformation)是建立資料新的表示的演算法,與資料的 原始表示相比,新的表示可能更容易被人或其他機器學習演算法所理解。無監督變換的一個 常見應用是降維(dimensionality reduction),它接受包含許多特徵的資料的高維表示,並找到表示該資料的一種新方法,用較少的特徵就可以概括其重要特性。降維的一個常見應用是為了視覺化將資料降為二維。

無監督變換的另一個應用是找到“構成”資料的各個組成部分。這方面的一個例子就是對文字文件集合進行主題提取。這裡的任務是找到每個文件中討論的未知主題,並學習每個 文件中出現了哪些主題。這可以用於追蹤社交媒體上的話題討論,比如選舉、槍支管制或流行歌手等話題。

與之相反,聚類演算法(clustering algorithm)將資料劃分成不同的組,每組包含相似的物項。思考向社交媒體網站上傳照片的例子。為了方便你整理照片,網站可能想要將同一個 人的照片分在一組。但網站並不知道每張照片是誰,也不知道你的照片集中出現了多少個 人。明智的做法是提取所有的人臉,並將看起來相似的人臉分在一組。但願這些人臉對應 同一個人,這樣圖片的分組也就完成了。

3.2 無監督學習的挑戰

無監督學習的一個主要挑戰就是評估演算法是否學到了有用的東西。無監督學習演算法一般用於不包含任何標籤資訊的資料,所以我們不知道正確的輸出應該是什麼。因此很難判斷一個模型是否“表現很好”。例如,假設我們的聚類演算法已經將所有的側臉照片和所有的正面照片進行分組。這肯定是人臉照片集合的一種可能的劃分方法,但並不是我們想要的那種方法。然而,我們沒有辦法“告訴”演算法我們要的是什麼,通常來說,評估無監督演算法結果的唯一方法就是人工檢查。

因此,如果資料科學家想要更好地理解資料,那麼無監督演算法通常可用於探索性的目的, 而不是作為大型自動化系統的一部分。無監督演算法的另一個常見應用是作為監督演算法的預 處理步驟。學習資料的一種新表示,有時可以提高監督演算法的精度,或者可以減少記憶體佔 用和時間開銷。

在開始學習“真正的”無監督演算法之前,我們先簡要討論幾種簡單又常用的預處理方法。 雖然預處理和縮放通常與監督學習演算法一起使用,但縮放方法並沒有用到與“監督”有關的資訊,所以它是無監督的。

3.3 預處理與縮放

上一章我們學到,一些演算法(如神經網路和 SVM)對資料縮放非常敏感。因此,通常的做法是對特徵進行調節,使資料表示更適合於這些演算法。通常來說,這是對資料的一種簡單的按特徵的縮放和移動。下面的程式碼(圖 3-1)給出了一個簡單的例子:

In[2]:
mglearn.plots.plot_scaling()

在這裡插入圖片描述

圖 3-1:對資料集縮放和預處理的各種方法

3.3.1 不同型別的預處理

在圖 3-1 中,第一張圖顯示的是一個模擬的有兩個特徵的二分類資料集。第一個特徵(x 軸)位於 10 到 15 之間。第二個特徵(y 軸)大約位於 1 到 9 之間。

接下來的 4 張圖展示了 4 種資料變換方法,都生成了更加標準的範圍。scikit-learn 中 的 StandardScaler 確保每個特徵的平均值為 0、方差為 1,使所有特徵都位於同一量級。但這種縮放不能保證特徵任何特定的最大值和最小值。RobustScaler 的工作原理與 StandardScaler 類似,確保每個特徵的統計屬性都位於同一範圍。但 RobustScaler 使用的 是中位數和四分位數 1 ,而不是平均值和方差。這樣 RobustScaler 會忽略與其他點有很大不 同的資料點(比如測量誤差)。這些與眾不同的資料點也叫異常值(outlier),可能會給其他縮放方法造成麻煩。

與之相反,MinMaxScaler 移動資料,使所有特徵都剛好位於 0 到 1 之間。對於二維資料集來說,所有的資料都包含在 x 軸 0 到 1 與 y 軸 0 到 1 組成的矩形中。

最後,Normalizer 用到一種完全不同的縮放方法。它對每個資料點進行縮放,使得特徵向量的歐式長度等於 1。換句話說,它將一個資料點投射到半徑為 1 的圓上(對於更高維度的情況,是球面)。這意味著每個資料點的縮放比例都不相同(乘以其長度的倒數)。如果只有資料的方向(或角度)是重要的,而特徵向量的長度無關緊要,那麼通常會使用這種 歸一化。

3.3.2 應用資料變換

前面我們已經看到不同型別的變換的作用,下面利用 scikit-learn 來應用這些變換。我們 將使用第 2 章見過的 cancer 資料集。通常在應用監督學習演算法之前使用預處理方法(比 如縮放)。舉個例子,比如我們想要將核 SVM(SVC)應用在 cancer 資料集上,並使用 MinMaxScaler 來預處理資料。首先載入資料集並將其分為訓練集和測試集(我們需要分開 的訓練集和資料集來對預處理後構建的監督模型進行評估):

In[3]:
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target,
 random_state=1)
print(X_train.shape)
print(X_test.shape)
Out[3]:
(426, 30)
(143, 30)

提醒一下,這個資料集包含 569 個資料點,每個資料點由 30 個測量值表示。我們將資料 集分成包含 426 個樣本的訓練集與包含 143 個樣本的測試集。

與之前構建的監督模型一樣,我們首先匯入實現預處理的類,然後將其例項化:

In[4]:
from sklearn.preprocessing import MinMaxScaler
scaler = MinMaxScaler()

然後,使用 fit 方法擬合縮放器(scaler),並將其應用於訓練資料。對於 MinMaxScaler 來 說,fit 方法計算訓練集中每個特徵的最大值和最小值。與第 2 章中的分類器和迴歸器 (regressor)不同,在對縮放器呼叫 fit 時只提供了 X_train,而不用 y_train:

In[5]:
scaler.fit(X_train)
Out[5]:
MinMaxScaler(copy=True, feature_range=(0, 1))

為了應用剛剛學習的變換(即對訓練資料進行實際縮放),我們使用縮放器的 transform 方法。在 scikit-learn 中,每當模型返回資料的一種新表示時,都可以使用 transform 方法:

In[6]:
# 變換資料
X_train_scaled = scaler.transform(X_train)
# 在縮放之前和之後分別列印資料集屬性
print("transformed shape: {}".format(X_train_scaled.shape))
print("per-feature minimum before scaling:\n {}".format(X_train.min(axis=0)))
print("per-feature maximum before scaling:\n {}".format(X_train.max(axis=0)))
print("per-feature minimum after scaling:\n {}".format(
 X_train_scaled.min(axis=0)))
print("per-feature maximum after scaling:\n {}".format(
 X_train_scaled.max(axis=0)))
Out[6]:
transformed shape: (426, 30)
per-feature minimum before scaling:
 [ 6.98 9.71 43.79 143.50 0.05 0.02 0. 0. 0.11
 0.05 0.12 0.36 0.76 6.80 0. 0. 0. 0.
 0.01 0. 7.93 12.02 50.41 185.20 0.07 0.03 0.
 0. 0.16 0.06]
per-feature maximum before scaling:
[ 28.11 39.28 188.5 2501.0 0.16 0.29 0.43 0.2
 0.300 0.100 2.87 4.88 21.98 542.20 0.03 0.14
 0.400 0.050 0.06 0.03 36.04 49.54 251.20 4254.00
 0.220 0.940 1.17 0.29 0.58 0.15]
per-feature minimum after scaling:
[ 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
per-feature maximum after scaling:
[ 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.
 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]

變換後的資料形狀與原始資料相同,特徵只是發生了移動和縮放。你可以看到,現在所有特徵都位於 0 到 1 之間,這也符合我們的預期。

為了將 SVM 應用到縮放後的資料上,還需要對測試集進行變換。這可以通過對 X_test 呼叫 transform 方法來完成:

In[7]:
# 對測試資料進行變換
X_test_scaled = scaler.transform(X_test)
# 在縮放之後列印測試資料的屬性
print("per-feature minimum after scaling:\n{}".format(X_test_scaled.min(axis=0)))
print("per-feature maximum after scaling:\n{}".format(X_test_scaled.max(axis=0)))
Out[7]:
per-feature minimum after scaling:
[ 0.034 0.023 0.031 0.011 0.141 0.044 0. 0. 0.154 -0.006
 -0.001 0.006 0.004 0.001 0.039 0.011 0. 0. -0.032 0.007
 0.027 0.058 0.02 0.009 0.109 0.026 0. 0. -0. -0.002]
per-feature maximum after scaling:
[ 0.958 0.815 0.956 0.894 0.811 1.22 0.88 0.933 0.932 1.037
 0.427 0.498 0.441 0.284 0.487 0.739 0.767 0.629 1.337 0.391
 0.896 0.793 0.849 0.745 0.915 1.132 1.07 0.924 1.205 1.631]

你可以發現,對測試集縮放後的最大值和最小值不是 1 和 0,這或許有些出乎意料。有些特徵甚至在 0~1 的範圍之外!對此的解釋是,MinMaxScaler(以及其他所有縮放器)總是對訓練集和測試集應用完全相同的變換。也就是說,transform 方法總是減去訓練集的最小值,然後除以訓練集的範圍,而這兩個值可能與測試集的最小值和範圍並不相同。

3.3.3 對訓練資料和測試資料進行相同的縮放

為了讓監督模型能夠在測試集上執行,對訓練集和測試集應用完全相同的變換是很重要 的。如果我們使用測試集的最小值和範圍,下面這個例子(圖 3-2)展示了會發生什麼:

In[8]:
from sklearn.datasets import make_blobs
# 構造資料
X, _ = make_blobs(n_samples=50, centers=5, random_state=4, cluster_std=2)
# 將其分為訓練集和測試集
X_train, X_test = train_test_split(X, random_state=5, test_size=.1)
# 繪製訓練集和測試集
fig, axes = plt.subplots(1, 3, figsize=(13, 4))
axes[0].scatter(X_train[:, 0], X_train[:, 1],
 c=mglearn.cm2(0), label="Training set", s=60)
axes[0].scatter(X_test[:, 0], X_test[:, 1], marker='^',
 c=mglearn.cm2(1), label="Test set", s=60)
axes[0].legend(loc='upper left')
axes[0].set_title("Original Data")
# 利用MinMaxScaler縮放資料
scaler = MinMaxScaler()
scaler.fit(X_train)
X_train_scaled = scaler.transform(X_train)
X_test_scaled = scaler.transform(X_test)
# 將正確縮放的資料視覺化
axes[1].scatter(X_train_scaled[:, 0], X_train_scaled[:, 1],
 c=mglearn.cm2(0), label="Training set", s=60)
axes[1].scatter(X_test_scaled[:, 0], X_test_scaled[:, 1], marker='^',
 c=mglearn.cm2(1), label="Test set", s=60)
axes[1].set_title("Scaled Data")
# 單獨對測試集進行縮放
# 使得測試集的最小值為0,最大值為1
# 千萬不要這麼做!這裡只是為了舉例
test_scaler = MinMaxScaler()
test_scaler.fit(X_test)
X_test_scaled_badly = test_scaler.transform(X_test)
# 將錯誤縮放的資料視覺化
axes[2].scatter(X_train_scaled[:, 0], X_train_scaled[:, 1],
 c=mglearn.cm2(0), label="training set", s=60)
axes[2].scatter(X_test_scaled_badly[:, 0], X_test_scaled_badly[:, 1],
 marker='^', c=mglearn.cm2(1), label="test set", s=60)
axes[2].set_title("Improperly Scaled Data")
for ax in axes:
ax.set_xlabel("Feature 0")
ax.set_ylabel("Feature 1")

在這裡插入圖片描述

圖 3-2:對左圖中的訓練資料和測試資料同時縮放的效果(中)和分別縮放的效果(右)

第一張圖是未縮放的二維資料集,其中訓練集用圓形表示,測試集用三角形表示。第二張 圖中是同樣的資料,但使用 MinMaxScaler 縮放。這裡我們呼叫 fit 作用在訓練集上,然後 呼叫 transform 作用在訓練集和測試集上。你可以發現,第二張圖中的資料集看起來與第 一張圖中的完全相同,只是座標軸刻度發生了變化。現在所有特徵都位於 0 到 1 之間。你 還可以發現,測試資料(三角形)的特徵最大值和最小值並不是 1 和 0。

第三張圖展示瞭如果我們對訓練集和測試集分別進行縮放會發生什麼。在這種情況下,對 訓練集和測試集而言,特徵的最大值和最小值都是 1 和 0。但現在資料集看起來不一樣。 測試集相對訓練集的移動不一致,因為它們分別做了不同的縮放。我們隨意改變了資料的排列。這顯然不是我們想要做的事情。

再換一種思考方式,想象你的測試集只有一個點。對於一個點而言,無法將其正確地縮放以滿足 MinMaxScaler 的最大值和最小值的要求。但是,測試集的大小不應該對你的處理方式有影響。

3.3.4 預處理對監督學習的作用

現在我們回到 cancer 資料集,觀察使用 MinMaxScaler 對學習 SVC 的作用(這是一種不同的 方法,實現了與第 2 章中相同的縮放)。首先,為了對比,我們再次在原始資料上擬合 SVC:

In[10]:
from sklearn.svm import SVC
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target,
 random_state=0)
svm = SVC(C=100)
svm.fit(X_train, y_train)
print("Test set accuracy: {:.2f}".format(svm.score(X_test, y_test)))
Out[10]:
Test set accuracy: 0.63

下面先用 MinMaxScaler 對資料進行縮放,然後再擬合 SVC:

In[11]:
# 使用0-1縮放進行預處理
scaler = MinMaxScaler()
scaler.fit(X_train)
X_train_scaled = scaler.transform(X_train)
X_test_scaled = scaler.transform(X_test)
# 在縮放後的訓練資料上學習SVM
svm.fit(X_train_scaled, y_train)
# 在縮放後的測試集上計算分數
print("Scaled test set accuracy: {:.2f}".format(
 svm.score(X_test_scaled, y_test)))
Out[11]:
Scaled test set accuracy: 0.97

正如我們上面所見,資料縮放的作用非常顯著。雖然資料縮放不涉及任何複雜的數學,但良好的做法仍然是使用 scikit-learn 提供的縮放機制,而不是自己重新實現它們,因為即使在這些簡單的計算中也容易犯錯。

你也可以通過改變使用的類將一種預處理演算法輕鬆替換成另一種,因為所有的預處理類都 具有相同的介面,都包含 fit 和 transform 方法:

In[12]:
# 利用零均值和單位方差的縮放方法進行預處理
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
scaler.fit(X_train)
X_train_scaled = scaler.transform(X_train)
X_test_scaled = scaler.transform(X_test)
# 在縮放後的訓練資料上學習SVM
svm.fit(X_train_scaled, y_train)
# 在縮放後的測試集上計算分數
print("SVM test accuracy: {:.2f}".format(svm.score(X_test_scaled, y_test)))
Out[12]:
SVM test accuracy: 0.96

前面我們已經看到了用於預處理的簡單資料變換的工作原理,下面繼續學習利用無監督學習進行更有趣的變換。

3.4 降維、特徵提取與流形學習

前面討論過,利用無監督學習進行資料變換可能有很多種目的。最常見的目的就是可視 化、壓縮資料,以及尋找資訊量更大的資料表示以用於進一步的處理。

為了實現這些目的,最簡單也最常用的一種演算法就是主成分分析。我們也將學習另外兩種 演算法:非負矩陣分解(NMF)和 t-SNE,前者通常用於特徵提取,後者通常用於二維散點 圖的視覺化。

3.4.1 主成分分析

主成分分析(principal component analysis,PCA)是一種旋轉資料集的方法,旋轉後的特徵在統計上不相關。在做完這種旋轉之後,通常是根據新特徵對解釋資料的重要性來選擇它的一個子集。下面的例子(圖 3-3)展示了 PCA 對一個模擬二維資料集的作用:

In[13]:
mglearn.plots.plot_pca_illustration()

在這裡插入圖片描述

圖 3-3:用 PCA 做資料變換

第一張圖(左上)顯示的是原始資料點,用不同顏色加以區分。演算法首先找到方差最大的 方向,將其標記為“成分 1”(Component 1)。這是資料中包含最多資訊的方向(或向量), 換句話說,沿著這個方向的特徵之間最為相關。然後,演算法找到與第一個方向正交(成直 角)且包含最多資訊的方向。在二維空間中,只有一個成直角的方向,但在更高維的空間 中會有(無窮)多的正交方向。雖然這兩個成分都畫成箭頭,但其頭尾的位置並不重要。 我們也可以將第一個成分畫成從中心指向左上,而不是指向右下。利用這一過程找到的方 向被稱為主成分(principal component),因為它們是資料方差的主要方向。一般來說,主 成分的個數與原始特徵相同。

第二張圖(右上)顯示的是同樣的資料,但現在將其旋轉,使得第一主成分與 x 軸平行且 第二主成分與 y 軸平行。在旋轉之前,從資料中減去平均值,使得變換後的資料以零為中 心。在 PCA 找到的旋轉表示中,兩個座標軸是不相關的,也就是說,對於這種資料表示, 除了對角線,相關矩陣全部為零。

我們可以通過僅保留一部分主成分來使用 PCA 進行降維。在這個例子中,我們可以僅保 留第一個主成分,正如圖 3-3 中第三張圖所示(左下)。這將資料從二維資料集降為一維數 據集。但要注意,我們沒有保留原始特徵之一,而是找到了最有趣的方向(第一張圖中從 左上到右下)並保留這一方向,即第一主成分。

最後,我們可以反向旋轉並將平均值重新加到資料中。這樣會得到圖 3-3 最後一張圖中的 資料。這些資料點位於原始特徵空間中,但我們僅保留了第一主成分中包含的資訊。這種 變換有時用於去除資料中的噪聲影響,或者將主成分中保留的那部分資訊視覺化。

  1. 將 PCA 應用於 cancer 資料集並視覺化

PCA 最常見的應用之一就是將高維資料集視覺化。正如第 1 章中所說,對於有兩個以上特 徵的資料,很難繪製散點圖。對於 Iris(鳶尾花)資料集,我們可以建立散點圖矩陣(見 第 1 章圖 1-3),通過展示特徵所有可能的兩兩組合來展示資料的區域性影像。但如果我們想 要檢視乳腺癌資料集,即便用散點圖矩陣也很困難。這個資料集包含 30 個特徵,這就導 致需要繪製 30 * 14 = 420 張散點圖!我們永遠不可能仔細觀察所有這些影像,更不用說試 圖理解它們了。

不過我們可以使用一種更簡單的視覺化方法——對每個特徵分別計算兩個類別(良性腫瘤 和惡性腫瘤)的直方圖(見圖 3-4)。

In[14]:
fig, axes = plt.subplots(15, 2, figsize=(10, 20))
malignant = cancer.data[cancer.target == 0]
benign = cancer.data[cancer.target == 1]
ax = axes.ravel()
for i in range(30):
 _, bins = np.histogram(cancer.data[:, i], bins=50)
 ax[i].hist(malignant[:, i], bins=bins, color=mglearn.cm3(0), alpha=.5)
 ax[i].hist(benign[:, i], bins=bins, color=mglearn.cm3(2), alpha=.5)
 ax[i].set_title(cancer.feature_names[i])
 ax[i].set_yticks(())
ax[0].set_xlabel("Feature magnitude")
ax[0].set_ylabel("Frequency")
ax[0].legend(["malignant", "benign"], loc="best")
fig.tight_layout()

在這裡插入圖片描述

圖 3-4:乳腺癌資料集中每個類別的特徵直方圖

這裡我們為每個特徵建立一個直方圖,計算具有某一特徵的資料點在特定範圍內(叫作 bin)的出現頻率。每張圖都包含兩個直方圖,一個是良性類別的所有點(藍色),一個是 惡性類別的所有點(紅色)。這樣我們可以瞭解每個特徵在兩個類別中的分佈情況,也可 以猜測哪些特徵能夠更好地區分良性樣本和惡性樣本。例如,“smoothness error”特徵似乎 沒有什麼資訊量,因為兩個直方圖大部分都重疊在一起,而“worst concave points”特徵 看起來資訊量相當大,因為兩個直方圖的交集很小。

但是,這種圖無法向我們展示變數之間的相互作用以及這種相互作用與類別之間的關係。 利用 PCA,我們可以獲取到主要的相互作用,並得到稍為完整的影像。我們可以找到前兩 個主成分,並在這個新的二維空間中用散點圖將資料視覺化。

在應用 PCA 之前,我們利用 StandardScaler 縮放資料,使每個特徵的方差均為 1:

In[15]:
from sklearn.datasets import load_breast_cancer
cancer = load_breast_cancer()
scaler = StandardScaler()
scaler.fit(cancer.data)
X_scaled = scaler.transform(cancer.data)

學習並應用 PCA 變換與應用預處理變換一樣簡單。我們將 PCA 物件例項化,呼叫 fit 方 法找到主成分,然後呼叫 transform 來旋轉並降維。預設情況下,PCA 僅旋轉(並移動) 資料,但保留所有的主成分。為了降低資料的維度,我們需要在建立 PCA 物件時指定想要 保留的主成分個數:

In[16]:
from sklearn.decomposition import PCA
# 保留資料的前兩個主成分
pca = PCA(n_components=2)
# 對乳腺癌資料擬合PCA模型
pca.fit(X_scaled)
# 將資料變換到前兩個主成分的方向上
X_pca = pca.transform(X_scaled)
print("Original shape: {}".format(str(X_scaled.shape)))
print("Reduced shape: {}".format(str(X_pca.shape)))
Out[16]:
Original shape: (569, 30)
Reduced shape: (569, 2)

現在我們可以對前兩個主成分作圖(圖 3-5):

In[17]:
# 對第一個和第二個主成分作圖,按類別著色
plt.figure(figsize=(8, 8))
mglearn.discrete_scatter(X_pca[:, 0], X_pca[:, 1], cancer.target)
plt.legend(cancer.target_names, loc="best")
plt.gca().set_aspect("equal")
plt.xlabel("First principal component")
plt.ylabel("Second principal component")

在這裡插入圖片描述

圖 3-5:利用前兩個主成分繪製乳腺癌資料集的二維散點圖

重要的是要注意,PCA 是一種無監督方法,在尋找旋轉方向時沒有用到任何類別資訊。它 只是觀察資料中的相關性。對於這裡所示的散點圖,我們繪製了第一主成分與第二主成分 的關係,然後利用類別資訊對資料點進行著色。你可以看到,在這個二維空間中兩個類別 被很好地分離。這讓我們相信,即使是線性分類器(在這個空間中學習一條直線)也可以 在區分這個兩個類別時表現得相當不錯。我們還可以看到,惡性點比良性點更加分散,這 一點也可以在圖 3-4 的直方圖中看出來。

PCA 的一個缺點在於,通常不容易對圖中的兩個軸做出解釋。主成分對應於原始資料中的 方向,所以它們是原始特徵的組合。但這些組合往往非常複雜,這一點我們很快就會看 到。在擬合過程中,主成分被儲存在 PCA 物件的 components_ 屬性中:

In[18]:
print("PCA component shape: {}".format(pca.components_.shape))
Out[18]:
PCA component shape: (2, 30)

components_ 中的每一行對應於一個主成分,它們按重要性排序(第一主成分排在首位, 以此類推)。列對應於 PCA 的原始特徵屬性,在本例中即為“mean radius”“mean texture” 等。我們來看一下 components_ 的內容:

In[19]:
print("PCA components:\n{}".format(pca.components_))
Out[19]:
PCA components:
[[ 0.219 0.104 0.228 0.221 0.143 0.239 0.258 0.261 0.138 0.064
 0.206 0.017 0.211 0.203 0.015 0.17 0.154 0.183 0.042 0.103
 0.228 0.104 0.237 0.225 0.128 0.21 0.229 0.251 0.123 0.132]
[-0.234 -0.06 -0.215 -0.231 0.186 0.152 0.06 -0.035 0.19 0.367
 -0.106 0.09 -0.089 -0.152 0.204 0.233 0.197 0.13 0.184 0.28
 -0.22 -0.045 -0.2 -0.219 0.172 0.144 0.098 -0.008 0.142 0.275]]

我們還可以用熱圖將係數視覺化(圖 3-6),這可能更容易理解:

In[20]:
plt.matshow(pca.components_, cmap='viridis')
plt.yticks([0, 1], ["First component", "Second component"])
plt.colorbar()
plt.xticks(range(len(cancer.feature_names)),
 cancer.feature_names, rotation=60, ha='left')
plt.xlabel("Feature")
plt.ylabel("Principal components")

在這裡插入圖片描述

圖 3-6:乳腺癌資料集前兩個主成分的熱圖

你可以看到,在第一個主成分中,所有特徵的符號相同(均為正,但前面我們提到過,箭 頭指向哪個方向無關緊要)。這意味著在所有特徵之間存在普遍的相關性。如果一個測量 值較大的話,其他的測量值可能也較大。第二個主成分的符號有正有負,而且兩個主成分 都包含所有 30 個特徵。這種所有特徵的混合使得解釋圖 3-6 中的座標軸變得十分困難。

  1. 特徵提取的特徵臉

前面提到過,PCA 的另一個應用是特徵提取。特徵提取背後的思想是,可以找到一種資料 表示,比給定的原始表示更適合於分析。特徵提取很有用,它的一個很好的應用例項就是 影像。影像由畫素組成,通常儲存為紅綠藍(RGB)強度。影像中的物件通常由上千個像 素組成,它們只有放在一起才有意義。

我們將給出用 PCA 對影像做特徵提取的一個簡單應用,即處理 Wild 資料集 Labeled Faces (標記人臉)中的人臉影像。這一資料集包含從網際網路下載的名人臉部影像,它包含從 21 世紀初開始的政治家、歌手、演員和運動員的人臉影像。我們使用這些影像的灰度版本, 並將它們按比例縮小以加快處理速度。你可以在圖 3-7 中看到其中一些影像:

In[21]:
from sklearn.datasets import fetch_lfw_people
people = fetch_lfw_people(min_faces_per_person=20, resize=0.7)
image_shape = people.images[0].shape
fix, axes = plt.subplots(2, 5, figsize=(15, 8),
 subplot_kw={'xticks': (), 'yticks': ()})
for target, image, ax in zip(people.target, people.images, axes.ravel()):
 ax.imshow(image)
 ax.set_title(people.target_names[target])

在這裡插入圖片描述

圖 3-7:來自 Wild 資料集中 Labeled Faces 的一些影像

一共有 3023 張影像,每張大小為 87 畫素 ×65 畫素,分別屬於 62 個不同的人:

In[22]:
print("people.images.shape: {}".format(people.images.shape))
print("Number of classes: {}".format(len(people.target_names)))
Out[22]:
people.images.shape: (3023, 87, 65)
Number of classes: 62

但這個資料集有些偏斜,其中包含 George W. Bush(小布什)和 Colin Powell(科林 • 鮑威 爾)的大量影像,正如你在下面所見:

In[23]:
# 計算每個目標出現的次數
counts = np.bincount(people.target)
# 將次數與目標名稱一起列印出來
for i, (count, name) in enumerate(zip(counts, people.target_names)):
 print("{0:25} {1:3}".format(name, count), end=' ')
 if (i + 1) % 3 == 0:
 print()
Out[23]:
Alejandro Toledo 39 Alvaro Uribe 35
Amelie Mauresmo 21 Andre Agassi 36
Angelina Jolie 20 Ariel Sharon 77
Arnold Schwarzenegger 42 Atal Bihari Vajpayee 24
Bill Clinton 29 Carlos Menem 21
Colin Powell 236 David Beckham 31
Donald Rumsfeld 121 George Robertson 22
George W Bush 530 Gerhard Schroeder 109
Gloria Macapagal Arroyo 44 Gray Davis 26
Guillermo Coria 30 Hamid Karzai 22
Hans Blix 39 Hugo Chavez 71
Igor Ivanov 20 Jack Straw 28
Jacques Chirac 52 Jean Chretien 55
Jennifer Aniston 21 Jennifer Capriati 42
Jennifer Lopez 21 Jeremy Greenstock 24
Jiang Zemin 20 John Ashcroft 53
John Negroponte 31 Jose Maria Aznar 23
Juan Carlos Ferrero 28 Junichiro Koizumi 60
Kofi Annan 32 Laura Bush 41
Lindsay Davenport 22 Lleyton Hewitt 41
Luiz Inacio Lula da Silva 48 Mahmoud Abbas 29
Megawati Sukarnoputri 33 Michael Bloomberg 20
Naomi Watts 22 Nestor Kirchner 37
Paul Bremer 20 Pete Sampras 22
Recep Tayyip Erdogan 30 Ricardo Lagos 27
Roh Moo-hyun 32 Rudolph Giuliani 26
Saddam Hussein 23 Serena Williams 52
Silvio Berlusconi 33 Tiger Woods 23
Tom Daschle 25 Tom Ridge 33
Tony Blair 144 Vicente Fox 32
Vladimir Putin 49 Winona Ryder 24

為了降低資料偏斜,我們對每個人最多隻取 50 張影像(否則,特徵提取將會被 George W. Bush 的可能性大大影響):

In[24]:
mask = np.zeros(people.target.shape, dtype=np.bool)
for target in np.unique(people.target):
 mask[np.where(people.target == target)[0][:50]] = 1
X_people = people.data[mask]
y_people = people.target[mask]
# 將灰度值縮放到0到1之間,而不是在0到255之間
# 以得到更好的資料穩定性
X_people = X_people / 255.

人臉識別的一個常見任務就是看某個前所未見的人臉是否屬於資料庫中的某個已知人物。 這在照片收集、社交媒體和安全應用中都有應用。解決這個問題的方法之一就是構建一個 分類器,每個人都是一個單獨的類別。但人臉資料庫中通常有許多不同的人,而同一個人 的影像很少(也就是說,每個類別的訓練樣例很少)。這使得大多數分類器的訓練都很困 難。另外,通常你還想要能夠輕鬆新增新的人物,不需要重新訓練一個大型模型。

一種簡單的解決方法是使用單一最近鄰分類器,尋找與你要分類的人臉最為相似的 人臉。這個分類器原則上可以處理每個類別只有一個訓練樣例的情況。下面看一下 KNeighborsClassifier 的表現如何:

In[25]:
from sklearn.neighbors import KNeighborsClassifier
# 將資料分為訓練集和測試集
X_train, X_test, y_train, y_test = train_test_split(
 X_people, y_people, stratify=y_people, random_state=0)
# 使用一個鄰居構建KNeighborsClassifier
knn = KNeighborsClassifier(n_neighbors=1)
knn.fit(X_train, y_train)
print("Test set score of 1-nn: {:.2f}".format(knn.score(X_test, y_test)))
Out[25]:
Test set score of 1-nn: 0.27

我們得到的精度為 26.6%。對於包含 62 個類別的分類問題來說,這實際上不算太差(隨機 猜測的精度約為 1/62=1.5%),但也不算好。我們每識別四次僅正確識別了一個人。

這裡就可以用到 PCA。想要度量人臉的相似度,計算原始畫素空間中的距離是一種相當糟 糕的方法。用畫素表示來比較兩張影像時,我們比較的是每個畫素的灰度值與另一張影像 對應位置的畫素灰度值。這種表示與人們對人臉影像的解釋方式有很大不同,使用這種原 始表示很難獲取到面部特徵。例如,如果使用畫素距離,那麼將人臉向右移動一個畫素將 會發生巨大的變化,得到一個完全不同的表示。我們希望,使用沿著主成分方向的距離可 以提高精度。這裡我們啟用 PCA 的白化(whitening)選項,它將主成分縮放到相同的尺 度。變換後的結果與使用 StandardScaler 相同。再次使用圖 3-3 中的資料,白化不僅對應 於旋轉資料,還對應於縮放資料使其形狀是圓形而不是橢圓(參見圖 3-8):

In[26]:
mglearn.plots.plot_pca_whitening()

在這裡插入圖片描述

圖 3-8:利用啟用白化的 PCA 進行資料變換

我們對訓練資料擬合 PCA 物件,並提取前 100 個主成分。然後對訓練資料和測試資料進行 變換:

In[27]:
pca = PCA(n_components=100, whiten=True, random_state=0).fit(X_train)
X_train_pca = pca.transform(X_train)
X_test_pca = pca.transform(X_test)
print("X_train_pca.shape: {}".format(X_train_pca.shape))
Out[27]:
X_train_pca.shape: (1547, 100)

新資料有 100 個特徵,即前 100 個主成分。現在,可以對新表示使用單一最近鄰分類器來 將我們的影像分類:

In[28]:
knn = KNeighborsClassifier(n_neighbors=1)
knn.fit(X_train_pca, y_train)
print("Test set accuracy: {:.2f}".format(knn.score(X_test_pca, y_test)))
Out[28]:
Test set accuracy: 0.36

我們的精度有了相當顯著的提高,從 26.6% 提升到 35.7%,這證實了我們的直覺,即主成 分可能提供了一種更好的資料表示。

對於影像資料,我們還可以很容易地將找到的主成分視覺化。請記住,成分對應於輸入空 間裡的方向。這裡的輸入空間是 87 畫素 ×65 畫素的灰度影像,所以在這個空間中的方向 也是 87 畫素 ×65 畫素的灰度影像。

我們來看一下前幾個主成分(圖 3-9):

In[29]:
print("pca.components_.shape: {}".format(pca.components_.shape))
Out[29]:
pca.components_.shape: (100, 5655)
In[30]:
fix, axes = plt.subplots(3, 5, figsize=(15, 12),
 subplot_kw={'xticks': (), 'yticks': ()})
for i, (component, ax) in enumerate(zip(pca.components_, axes.ravel())):
 ax.imshow(component.reshape(image_shape),
 cmap='viridis')
 ax.set_title("{}. component".format((i + 1)))

在這裡插入圖片描述

圖 3-9:人臉資料集前 15 個主成分的成分向量

雖然我們肯定無法理解這些成分的所有內容,但可以猜測一些主成分捕捉到了人臉影像的哪些方面。第一個主成分似乎主要編碼的是人臉與背景的對比,第二個主成分編碼的是人 臉左半部分和右半部分的明暗程度差異,如此等等。雖然這種表示比原始畫素值的語義稍強,但它仍與人們感知人臉的方式相去甚遠。由於 PCA 模型是基於畫素的,因此人臉的 相對位置(眼睛、下巴和鼻子的位置)和明暗程度都對兩張影像在畫素表示中的相似程度 有很大影響。但人臉的相對位置和明暗程度可能並不是人們首先感知的內容。在要求人們 評價人臉的相似度時,他們更可能會使用年齡、性別、面部表情和髮型等屬性,而這些屬 性很難從畫素強度中推斷出來。重要的是要記住,演算法對資料(特別是視覺資料,比如人 們非常熟悉的影像)的解釋通常與人類的解釋方式大不相同。

不過讓我們回到 PCA 的具體案例。我們對 PCA 變換的介紹是:先旋轉資料,然後刪除方 差較小的成分。另一種有用的解釋是嘗試找到一些數字(PCA 旋轉後的新特徵值),使我 們可以將測試點表示為主成分的加權求和(見圖 3-10)。

在這裡插入圖片描述

圖 3-10:圖解 PCA:將影像分解為成分的加權求和

這裡 x0、x1 等是這個資料點的主成分的係數,換句話說,它們是影像在旋轉後的空間中的 表示。

我們還可以用另一種方法來理解 PCA 模型,就是僅使用一些成分對原始資料進行重建。 在圖 3-3 中,在去掉第二個成分並來到第三張圖之後,我們反向旋轉並重新加上平均值, 這樣就在原始空間中獲得去掉第二個成分的新資料點,正如最後一張圖所示。我們可以對 人臉做類似的變換,將資料降維到只包含一些主成分,然後反向旋轉回到原始空間。回到 原始特徵空間可以通過 inverse_transform 方法來實現。這裡我們分別利用 10 個、50 個、 100 個和 500 個成分對一些人臉進行重建並將其視覺化(圖 3-11):

In[32]:
mglearn.plots.plot_pca_faces(X_train, X_test, image_shape)

在這裡插入圖片描述

圖 3-11:利用越來越多的主成分對三張人臉影像進行重建

可以看到,在僅使用前 10 個主成分時,僅捕捉到了圖片的基本特點,比如人臉方向和明 暗程度。隨著使用的主成分越來越多,影像中也保留了越來越多的細節。這對應於圖 3-10 的求和中包含越來越多的項。如果使用的成分個數與畫素個數相等,意味著我們在旋轉後 不會丟棄任何資訊,可以完美重建影像。

我們還可以嘗試使用 PCA 的前兩個主成分,將資料集中的所有人臉在散點圖中視覺化 (圖 3-12),其類別在圖中給出。這與我們對 cancer 資料集所做的類似:

In[33]:
mglearn.discrete_scatter(X_train_pca[:, 0], X_train_pca[:, 1], y_train)
plt.xlabel("First principal component")
plt.ylabel("Second principal component")

在這裡插入圖片描述

圖 3-12:利用前兩個主成分繪製人臉資料集的散點圖(cancer 資料集的對應影像見圖 3-5)

如你所見,如果我們只使用前兩個主成分,整個資料只是一大團,看不到類別之間的分 界。這並不意外,因為即使有 10 個成分(正如圖 3-11 所示),PCA 也僅捕捉到人臉非常 粗略的特徵。

3.4.2 非負矩陣分解

非負矩陣分解(non-negative matrix factorization,NMF)是另一種無監督學習演算法,其目 的在於提取有用的特徵。它的工作原理類似於 PCA,也可以用於降維。與 PCA 相同,我 們試圖將每個資料點寫成一些分量的加權求和,正如圖 3-10 所示。但在 PCA 中,我們想 要的是正交分量,並且能夠解釋儘可能多的資料方差;而在 NMF 中,我們希望分量和系 數均為非負,也就是說,我們希望分量和係數都大於或等於 0。因此,這種方法只能應用 於每個特徵都是非負的資料,因為非負分量的非負求和不可能變為負值。

將資料分解成非負加權求和的這個過程,對由多個獨立源相加(或疊加)建立而成的資料 特別有用,比如多人說話的音軌或包含多種樂器的音樂。在這種情況下,NMF 可以識別 出組成合成資料的原始分量。總的來說,與 PCA 相比,NMF 得到的分量更容易解釋,因 為負的分量和係數可能會導致難以解釋的抵消效應(cancellation effect)。舉個例子,圖 3-9 中的特徵臉同時包含正數和負數,我們在 PCA 的說明中也提到過,正負號實際上是任 意的。在將 NMF 應用於人臉資料集之前,我們先來簡要回顧一下模擬資料。

  1. 將 NMF 應用於模擬資料

與使用 PCA 不同,我們需要保證資料是正的,NMF 能夠對資料進行操作。這說明資料相 對於原點 (0, 0) 的位置實際上對 NMF 很重要。因此,你可以將提取出來的非負分量看作是 從 (0, 0) 到資料的方向。

下面的例子(圖 3-13)給出了 NMF 在二維玩具資料上的結果:

In[34]:
mglearn.plots.plot_nmf_illustration()

在這裡插入圖片描述

圖 3-13:兩個分量的非負矩陣分解(左)和一個分量的非負矩陣分解(右)找到的分量

對於兩個分量的 NMF(如左圖所示),顯然所有資料點都可以寫成這兩個分量的正陣列 合。如果有足夠多的分量能夠完美地重建資料(分量個數與特徵個數相同),那麼演算法會 選擇指向資料極值的方向。

如果我們僅使用一個分量,那麼 NMF 會建立一個指向平均值的分量,因為指向這裡可以 對資料做出最好的解釋。你可以看到,與 PCA 不同,減少分量個數不僅會刪除一些方向, 而且會建立一組完全不同的分量! NMF 的分量也沒有按任何特定方法排序,所以不存在 “第一非負分量”:所有分量的地位平等。

NMF 使用了隨機初始化,根據隨機種子的不同可能會產生不同的結果。在相對簡單的情 況下(比如兩個分量的模擬資料),所有資料都可以被完美地解釋,那麼隨機性的影響很 小(雖然可能會影響分量的順序或尺度)。在更加複雜的情況下,影響可能會很大。

  1. 將 NMF 應用於人臉影像

現在我們將 NMF 應用於之前用過的 Wild 資料集中的 Labeled Faces。NMF 的主要引數是 我們想要提取的分量個數。通常來說,這個數字要小於輸入特徵的個數(否則的話,將每 個畫素作為單獨的分量就可以對資料進行解釋)。

首先,我們來觀察分量個數如何影響 NMF 重建資料的好壞(圖 3-14):

In[35]:
mglearn.plots.plot_nmf_faces(X_train, X_test, image_shape)

在這裡插入圖片描述

圖 3-14:利用越來越多分量的 NMF 重建三張人臉影像

反向變換的資料質量與使用 PCA 時類似,但要稍差一些。這是符合預期的,因為 PCA 找 到的是重建的最佳方向。NMF 通常並不用於對資料進行重建或編碼,而是用於在資料中 尋找有趣的模式。

我們嘗試僅提取一部分分量(比如 15 個),初步觀察一下資料。其結果見圖 3-15。

In[36]:
from sklearn.decomposition import NMF
nmf = NMF(n_components=15, random_state=0)
nmf.fit(X_train)
X_train_nmf = nmf.transform(X_train)
X_test_nmf = nmf.transform(X_test)
fix, axes = plt.subplots(3, 5, figsize=(15, 12),
 subplot_kw={'xticks': (), 'yticks': ()})
for i, (component, ax) in enumerate(zip(nmf.components_, axes.ravel())):
 ax.imshow(component.reshape(image_shape))
 ax.set_title("{}. component".format(i))

在這裡插入圖片描述

圖 3-15:使用 15 個分量的 NMF 在人臉資料集上找到的分量

這些分量都是正的,因此比圖 3-9 所示的 PCA 分量更像人臉原型。例如,你可以清楚地看 到,分量 3(component 3)顯示了稍微向右轉動的人臉,而分量 7(component 7)則顯示 了稍微向左轉動的人臉。我們來看一下這兩個分量特別大的那些影像,分別如圖 3-16 和圖 3-17 所示。

In[37]:
compn = 3
# 按第3個分量排序,繪製前10張影像
inds = np.argsort(X_train_nmf[:, compn])[::-1]
fig, axes = plt.subplots(2, 5, figsize=(15, 8),
 subplot_kw={'xticks': (), 'yticks': ()})
for i, (ind, ax) in enumerate(zip(inds, axes.ravel())):
 ax.imshow(X_train[ind].reshape(image_shape))
compn = 7
# 按第7個分量排序,繪製前10張影像
inds = np.argsort(X_train_nmf[:, compn])[::-1]
fig, axes = plt.subplots(2, 5, figsize=(15, 8),
 subplot_kw={'xticks': (), 'yticks': ()})
for i, (ind, ax) in enumerate(zip(inds, axes.ravel())):
 ax.imshow(X_train[ind].reshape(image_shape))

在這裡插入圖片描述

圖 3-16:分量 3 係數較大的人臉

在這裡插入圖片描述

圖 3-17:分量 7 係數較大的人臉

正如所料,分量 3 係數較大的人臉都是向右看的人臉(圖 3-16),而分量 7 係數較大的人 臉都向左看(圖 3-17)。如前所述,提取這樣的模式最適合於具有疊加結構的資料,包括 音訊、基因表達和文字資料。我們通過一個模擬資料的例子來看一下這種用法。

假設我們對一個訊號感興趣,它是三個不同訊號源合成的(圖 3-18):

In[38]:
S = mglearn.datasets.make_signals()
plt.figure(figsize=(6, 1))
plt.plot(S, '-')
plt.xlabel("Time")
plt.ylabel("Signal")

在這裡插入圖片描述

圖 3-18:原始訊號源

不幸的是,我們無法觀測到原始訊號,只能觀測到三個訊號的疊加混合。我們想要將混合 訊號分解為原始分量。假設我們有許多種不同的方法來觀測混合訊號(比如有 100 臺測量 裝置),每種方法都為我們提供了一系列測量結果。

In[39]:
# 將資料混合成100維的狀態
A = np.random.RandomState(0).uniform(size=(100, 3))
X = np.dot(S, A.T)
print("Shape of measurements: {}".format(X.shape))
Out[39]:
Shape of measurements: (2000, 100)

我們可以用 NMF 來還原這三個訊號:

In[40]:
nmf = NMF(n_components=3, random_state=42)
S_ = nmf.fit_transform(X)
print("Recovered signal shape: {}".format(S_.shape))
Out[40]:
Recovered signal shape: (2000, 3)

為了對比,我們也應用了 PCA:

In[41]:
pca = PCA(n_components=3)
H = pca.fit_transform(X)

圖 3-19 給出了 NMF 和 PCA 發現的訊號活動:

In[42]:
models = [X, S, S_, H]
names = ['Observations (first three measurements)',
 'True sources',
 'NMF recovered signals',
 'PCA recovered signals']
 fig, axes = plt.subplots(4, figsize=(8, 4), gridspec_kw={'hspace': .5},
 subplot_kw={'xticks': (), 'yticks': ()})
for model, name, ax in zip(models, names, axes):
 ax.set_title(name)
 ax.plot(model[:, :3], '-')

在這裡插入圖片描述

圖元 3-19:利用 NMF 和 PCA 還原混合訊號源

圖中包含來自 X 的 100 次測量中的 3 次,用於參考。可以看到,NMF 在發現原始訊號源 時得到了不錯的結果,而 PCA 則失敗了,僅使用第一個成分來解釋資料中的大部分變化。 要記住,NMF 生成的分量是沒有順序的。在這個例子中,NMF 分量的順序與原始訊號完 全相同(參見三條曲線的顏色),但這純屬偶然。

還有許多其他演算法可用於將每個資料點分解為一系列固定分量的加權求和,正如 PCA 和 NMF 所做的那樣。討論所有這些演算法已超出了本書的範圍,而且描述對分量和係數的約 束通常要涉及概率論。如果你對這種型別的模式提取感興趣,我們推薦你學習 scikitlearn 使用者指南中關於獨立成分分析(ICA)、因子分析(FA)和稀疏編碼(字典學習)等 的內容,所有這些內容都可以在關於分解方法的頁面中找到(http://scikit-learn.org/stable/ modules/decomposition.html)。

3.4.3 用t-SNE進行流形學習

雖然 PCA 通常是用於變換資料的首選方法,使你能夠用散點圖將其視覺化,但這一方法 的性質(先旋轉然後減少方向)限制了其有效性,正如我們在 Wild 資料集 Labeled Faces 的散點圖中所看到的那樣。有一類用於視覺化的演算法叫作流形學習演算法(manifold learning algorithm),它允許進行更復雜的對映,通常也可以給出更好的視覺化。其中特別有用的一 個就是 t-SNE 演算法。

流形學習演算法主要用於視覺化,因此很少用來生成兩個以上的新特徵。其中一些演算法(包 括 t-SNE)計算訓練資料的一種新表示,但不允許變換新資料。這意味著這些演算法不能用 於測試集:更確切地說,它們只能變換用於訓練的資料。流形學習對探索性資料分析是很 有用的,但如果最終目標是監督學習的話,則很少使用。t-SNE 背後的思想是找到資料的 一個二維表示,儘可能地保持資料點之間的距離。t-SNE 首先給出每個資料點的隨機二維 表示,然後嘗試讓在原始特徵空間中距離較近的點更加靠近,原始特徵空間中相距較遠的 點更加遠離。t-SNE 重點關注距離較近的點,而不是保持距離較遠的點之間的距離。換句 話說,它試圖儲存那些表示哪些點比較靠近的資訊。

我們將對 scikit-learn 包含的一個手寫數字資料集 2 應用 t-SNE 流形學習演算法。在這個數 據集中,每個資料點都是 0 到 9 之間手寫數字的一張 8×8 灰度影像。圖 3-20 給出了每個 類別的一個例子。

In[43]:
from sklearn.datasets import load_digits
digits = load_digits()
fig, axes = plt.subplots(2, 5, figsize=(10, 5),
 subplot_kw={'xticks':(), 'yticks': ()})
for ax, img in zip(axes.ravel(), digits.images):
 ax.imshow(img)

在這裡插入圖片描述

圖 3-20:digits 資料集的示例影像

我們用 PCA 將降到二維的資料視覺化。我們對前兩個主成分作圖,並按類別對資料點著 色(圖 3-21):

In[44]:
# 構建一個PCA模型
pca = PCA(n_components=2)
pca.fit(digits.data)
# 將digits資料變換到前兩個主成分的方向上
digits_pca = pca.transform(digits.data)
colors = ["#476A2A", "#7851B8", "#BD3430", "#4A2D4E", "#875525",
 "#A83683", "#4E655E", "#853541", "#3A3120", "#535D8E"]
plt.figure(figsize=(10, 10))
plt.xlim(digits_pca[:, 0].min(), digits_pca[:, 0].max())
plt.ylim(digits_pca[:, 1].min(), digits_pca[:, 1].max())
for i in range(len(digits.data)):
 # 將資料實際繪製成文字,而不是散點
 plt.text(digits_pca[i, 0], digits_pca[i, 1], str(digits.target[i]),
 color = colors[digits.target[i]],
 fontdict={'weight': 'bold', 'size': 9})
plt.xlabel("First principal component")
plt.ylabel("Second principal component")

在這裡插入圖片描述

圖 3-21:利用前兩個主成分繪製 digits 資料集的散點圖

實際上,這裡我們用每個類別對應的數字作為符號來顯示每個類別的位置。利用前兩個主 成分可以將數字 0、6 和 4 相對較好地分開,儘管仍有重疊。大部分其他數字都大量重疊 在一起。

我們將 t-SNE 應用於同一個資料集,並對結果進行比較。由於 t-SNE 不支援變換新資料, 所以 TSNE 類沒有 transform 方法。我們可以呼叫 fit_transform 方法來代替,它會構建模型並立刻返回變換後的資料(見圖 3-22):

In[45]:
from sklearn.manifold import TSNE
tsne = TSNE(random_state=42)
# 使用fit_transform而不是fit,因為TSNE沒有transform方法
digits_tsne = tsne.fit_transform(digits.data)
In[46]:
plt.figure(figsize=(10, 10))
plt.xlim(digits_tsne[:, 0].min(), digits_tsne[:, 0].max() + 1)
plt.ylim(digits_tsne[:, 1].min(), digits_tsne[:, 1].max() + 1)
for i in range(len(digits.data)):
 # 將資料實際繪製成文字,而不是散點
 plt.text(digits_tsne[i, 0], digits_tsne[i, 1], str(digits.target[i]),
 color = colors[digits.target[i]],
 fontdict={'weight': 'bold', 'size': 9})
plt.xlabel("t-SNE feature 0")
plt.xlabel("t-SNE feature 1")

在這裡插入圖片描述

圖 3-22:利用 t-SNE 找到的兩個分量繪製 digits 資料集的散點圖

t-SNE 的結果非常棒。所有類別都被明確分開。數字 1 和 9 被分成幾塊,但大多數類別都 形成一個密集的組。要記住,這種方法並不知道類別標籤:它完全是無監督的。但它能夠 找到資料的一種二維表示,僅根據原始空間中資料點之間的靠近程度就能夠將各個類別明 確分開。

t-SNE 演算法有一些調節引數,雖然預設引數的效果通常就很好。你可以嘗試修改 perplexity 和 early_exaggeration,但作用一般很小。

3.5 聚類

我們前面說過,聚類(clustering)是將資料集劃分成組的任務,這些組叫作簇(cluster)。 其目標是劃分資料,使得一個簇內的資料點非常相似且不同簇內的資料點非常不同。與分 類演算法類似,聚類演算法為每個資料點分配(或預測)一個數字,表示這個點屬於哪個簇。

3.5.1 k均值聚類

k 均值聚類是最簡單也最常用的聚類演算法之一。它試圖找到代表資料特定區域的簇中心 (cluster center)。演算法交替執行以下兩個步驟:將每個資料點分配給最近的簇中心,然後將 每個簇中心設定為所分配的所有資料點的平均值。如果簇的分配不再發生變化,那麼演算法 結束。下面的例子(圖 3-23)在一個模擬資料集上對這一演算法進行說明:

In[47]:
mglearn.plots.plot_kmeans_algorithm()

在這裡插入圖片描述

圖 3-23:輸入資料與 k 均值演算法的三個步驟

簇中心用三角形表示,而資料點用圓形表示。顏色表示簇成員。我們指定要尋找三個簇, 所以通過宣告三個隨機資料點為簇中心來將演算法初始化(見圖中“Initialization”/“初 始化”)。然後開始迭代演算法。首先,每個資料點被分配給距離最近的簇中心(見圖中 “Assign Points (1)”/“分配資料點(1)”)。接下來,將簇中心修改為所分配點的平均值 (見圖中“Recompute Centers (1)”/“重新計算中心(1)”)。然後將這一過程再重複兩次。 在第三次迭代之後,為簇中心分配的資料點保持不變,因此演算法結束。 給定新的資料點,k 均值會將其分配給最近的簇中心。下一個例子(圖 3-24)展示了圖 3-23 學到的簇中心的邊界:

In[48]:
mglearn.plots.plot_kmeans_boundaries()

在這裡插入圖片描述

圖 3-24:k 均值演算法找到的簇中心和簇邊界

用 scikit-learn 應用 k 均值相當簡單。下面我們將其應用於上圖中的模擬資料。我們將 KMeans 類例項化,並設定我們要尋找的簇個數 3 。然後對資料呼叫 fit 方法:

In[49]:
from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans
# 生成模擬的二維資料
X, y = make_blobs(random_state=1)
# 構建聚類模型
kmeans = KMeans(n_clusters=3)
kmeans.fit(X)

演算法執行期間,為 X 中的每個訓練資料點分配一個簇標籤。你可以在 kmeans.labels_ 屬性 中找到這些標籤:

In[50]:
print("Cluster memberships:\n{}".format(kmeans.labels_))
Out[50]:
Cluster memberships:
[1 2 2 2 0 0 0 2 1 1 2 2 0 1 0 0 0 1 2 2 0 2 0 1 2 0 0 1 1 0 1 1 0 1 2 0 2
 2 2 0 0 2 1 2 2 0 1 1 1 1 2 0 0 0 1 0 2 2 1 1 2 0 0 2 2 0 1 0 1 2 2 2 0 1
 1 2 0 0 1 2 1 2 2 0 1 1 1 1 2 1 0 1 1 2 2 0 0 1 0 1]

因為我們要找的是 3 個簇,所以簇的編號是 0 到 2。

你也可以用 predict 方法為新資料點分配簇標籤。預測時會將最近的簇中心分配給每個新 資料點,但現有模型不會改變。對訓練集執行 predict 會返回與 labels_ 相同的結果:

In[51]:
print(kmeans.predict(X))
Out[51]:
[1 2 2 2 0 0 0 2 1 1 2 2 0 1 0 0 0 1 2 2 0 2 0 1 2 0 0 1 1 0 1 1 0 1 2 0 2
 2 2 0 0 2 1 2 2 0 1 1 1 1 2 0 0 0 1 0 2 2 1 1 2 0 0 2 2 0 1 0 1 2 2 2 0 1
 1 2 0 0 1 2 1 2 2 0 1 1 1 1 2 1 0 1 1 2 2 0 0 1 0 1]

可以看到,聚類演算法與分類演算法有些相似,每個元素都有一個標籤。但並不存在真實的標 籤,因此標籤本身並沒有先驗意義。我們回到之前討論過的人臉影像聚類的例子。聚類的 結果可能是,演算法找到的第 3 個簇僅包含你朋友 Bela 的面孔。但只有在檢視圖片之後才能 知道這一點,而且數字 3 是任意的。演算法給你的唯一資訊就是所有標籤為 3 的人臉都是相 似的。

對於我們剛剛在二維玩具資料集上執行的聚類演算法,這意味著我們不應該為其中一組的標 籤是 0、另一組的標籤是 1 這一事實賦予任何意義。再次執行該演算法可能會得到不同的簇 編號,原因在於初始化的隨機性質。

下面又給出了這個資料的影像(圖 3-25)。簇中心被儲存在 cluster_centers_ 屬性中,我 們用三角形表示它們:

In[52]:
mglearn.discrete_scatter(X[:, 0], X[:, 1], kmeans.labels_, markers='o')
mglearn.discrete_scatter(
 kmeans.cluster_centers_[:, 0], kmeans.cluster_centers_[:, 1], [0, 1, 2],
 markers='^', markeredgewidth=2)

在這裡插入圖片描述

圖 3-25:3 個簇的 k 均值演算法找到的簇分配和簇中心

我們也可以使用更多或更少的簇中心(圖 3-26):

In[53]:
fig, axes = plt.subplots(1, 2, figsize=(10, 5))
# 使用2個簇中心:
kmeans = KMeans(n_clusters=2)
kmeans.fit(X)
assignments = kmeans.labels_
mglearn.discrete_scatter(X[:, 0], X[:, 1], assignments, ax=axes[0])
# 使用5個簇中心:
kmeans = KMeans(n_clusters=5)
kmeans.fit(X)
assignments = kmeans.labels_
mglearn.discrete_scatter(X[:, 0], X[:, 1], assignments, ax=axes[1])

在這裡插入圖片描述

圖 3-26:使用 2 個簇(左)和 5 個簇(右)的 k 均值演算法找到的簇分配

  1. k均值的失敗案例

即使你知道給定資料集中簇的“正確”個數,k 均值可能也不是總能找到它們。每個簇僅 由其中心定義,這意味著每個簇都是凸形(convex)。因此,k 均值只能找到相對簡單的形 狀。k 均值還假設所有簇在某種程度上具有相同的“直徑”,它總是將簇之間的邊界剛好畫 在簇中心的中間位置。有時這會導致令人驚訝的結果,如圖 3-27 所示:

In[54]:
X_varied, y_varied = make_blobs(n_samples=200,
 cluster_std=[1.0, 2.5, 0.5],
 random_state=170)
y_pred = KMeans(n_clusters=3, random_state=0).fit_predict(X_varied)
mglearn.discrete_scatter(X_varied[:, 0], X_varied[:, 1], y_pred)
plt.legend(["cluster 0", "cluster 1", "cluster 2"], loc='best')
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")

在這裡插入圖片描述

圖 3-27:簇的密度不同時,k 均值找到的簇分配

你可能會認為,左下方的密集區域是第一個簇,右上方的密集區域是第二個,中間密度較 小的區域是第三個。但事實上,簇 0 和簇 1 都包含一些遠離簇中其他點的點。

k 均值還假設所有方向對每個簇都同等重要。圖 3-28 顯示了一個二維資料集,資料中包含 明確分開的三部分。但是這三部分被沿著對角線方向拉長。由於 k 均值僅考慮到最近簇中 心的距離,所以它無法處理這種型別的資料:

In[55]:
# 生成一些隨機分組資料
X, y = make_blobs(random_state=170, n_samples=600)
rng = np.random.RandomState(74)
# 變換資料使其拉長
transformation = rng.normal(size=(2, 2))
X = np.dot(X, transformation)
# 將資料聚類成3個簇
kmeans = KMeans(n_clusters=3)
kmeans.fit(X)
y_pred = kmeans.predict(X)
# 畫出簇分配和簇中心
plt.scatter(X[:, 0], X[:, 1], c=y_pred, cmap=mglearn.cm3)
plt.scatter(kmeans.cluster_centers_[:, 0], kmeans.cluster_centers_[:, 1],
 marker='^', c=[0, 1, 2], s=100, linewidth=2, cmap=mglearn.cm3)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")

在這裡插入圖片描述

圖 3-28:k 均值無法識別非球形簇

如果簇的形狀更加複雜,比如我們在第 2 章遇到的 two_moons 資料,那麼 k 均值的表現也 很差(見圖 3-29):

In[56]:
# 生成模擬的two_moons資料(這次的噪聲較小)
from sklearn.datasets import make_moons
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)
# 將資料聚類成2個簇
kmeans = KMeans(n_clusters=2)
kmeans.fit(X)
y_pred = kmeans.predict(X)
# 畫出簇分配和簇中心
plt.scatter(X[:, 0], X[:, 1], c=y_pred, cmap=mglearn.cm2, s=60)
plt.scatter(kmeans.cluster_centers_[:, 0], kmeans.cluster_centers_[:, 1],
 marker='^', c=[mglearn.cm2(0), mglearn.cm2(1)], s=100, linewidth=2)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")

在這裡插入圖片描述

圖 3-29:k 均值無法識別具有複雜形狀的簇

這裡我們希望聚類演算法能夠發現兩個半月形。但利用 k 均值演算法是不可能做到這一點的。

  1. 向量量化,或者將k均值看作分解

雖然 k 均值是一種聚類演算法,但在 k 均值和分解方法(比如之前討論過的 PCA 和 NMF) 之間存在一些有趣的相似之處。你可能還記得,PCA 試圖找到資料中方差最大的方向,而 NMF 試圖找到累加的分量,這通常對應於資料的“極值”或“部分”(見圖 3-13)。兩種 方法都試圖將資料點表示為一些分量之和。與之相反,k 均值則嘗試利用簇中心來表示每 個資料點。你可以將其看作僅用一個分量來表示每個資料點,該分量由簇中心給出。這種 觀點將 k 均值看作是一種分解方法,其中每個點用單一分量來表示,這種觀點被稱為向量 量化(vector quantization)。

我們來並排比較 PCA、NMF 和 k 均值,分別顯示提取的分量(圖 3-30),以及利用 100 個 分量對測試集中人臉的重建(圖 3-31)。對於 k 均值,重建就是在訓練集中找到的最近的 簇中心:

In[57]:
X_train, X_test, y_train, y_test = train_test_split(
 X_people, y_people, stratify=y_people, random_state=0)
nmf = NMF(n_components=100, random_state=0)
nmf.fit(X_train)
pca = PCA(n_components=100, random_state=0)
pca.fit(X_train)
kmeans = KMeans(n_clusters=100, random_state=0)
kmeans.fit(X_train)
X_reconstructed_pca = pca.inverse_transform(pca.transform(X_test))
X_reconstructed_kmeans = kmeans.cluster_centers_[kmeans.predict(X_test)]
X_reconstructed_nmf = np.dot(nmf.transform(X_test), nmf.components_)
In[58]:
fig, axes = plt.subplots(3, 5, figsize=(8, 8),
 subplot_kw={'xticks': (), 'yticks': ()})
 fig.suptitle("Extracted Components")
for ax, comp_kmeans, comp_pca, comp_nmf in zip(
 axes.T, kmeans.cluster_centers_, pca.components_, nmf.components_):
 ax[0].imshow(comp_kmeans.reshape(image_shape))
 ax[1].imshow(comp_pca.reshape(image_shape), cmap='viridis')
 ax[2].imshow(comp_nmf.reshape(image_shape))
axes[0, 0].set_ylabel("kmeans")
axes[1, 0].set_ylabel("pca")
axes[2, 0].set_ylabel("nmf")
fig, axes = plt.subplots(4, 5, subplot_kw={'xticks': (), 'yticks': ()},
 figsize=(8, 8))
fig.suptitle("Reconstructions")
for ax, orig, rec_kmeans, rec_pca, rec_nmf in zip(
 axes.T, X_test, X_reconstructed_kmeans, X_reconstructed_pca,
 X_reconstructed_nmf):
 ax[0].imshow(orig.reshape(image_shape))
 ax[1].imshow(rec_kmeans.reshape(image_shape))
 ax[2].imshow(rec_pca.reshape(image_shape))
 ax[3].imshow(rec_nmf.reshape(image_shape))
axes[0, 0].set_ylabel("original")
axes[1, 0].set_ylabel("kmeans")
axes[2, 0].set_ylabel("pca")
axes[3, 0].set_ylabel("nmf")

在這裡插入圖片描述

圖 3-30:對比 k 均值的簇中心與 PCA 和 NMF 找到的分量

在這裡插入圖片描述

圖 3-31:利用 100 個分量(或簇中心)的 k 均值、PCA 和 NMF 的影像重建的對比——k 均值的每 張影像中僅使用了一個簇中心

利用 k 均值做向量量化的一個有趣之處在於,可以用比輸入維度更多的簇來對資料進行編 碼。讓我們回到 two_moons 資料。利用 PCA 或 NMF,我們對這個資料無能為力,因為它 只有兩個維度。使用 PCA 或 NMF 將其降到一維,將會完全破壞資料的結構。但通過使用 更多的簇中心,我們可以用 k 均值找到一種更具表現力的表示(見圖 3-32):

In[59]:
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)
kmeans = KMeans(n_clusters=10, random_state=0)
kmeans.fit(X)
y_pred = kmeans.predict(X)
plt.scatter(X[:, 0], X[:, 1], c=y_pred, s=60, cmap='Paired')
plt.scatter(kmeans.cluster_centers_[:, 0], kmeans.cluster_centers_[:, 1], s=60,
 marker='^', c=range(kmeans.n_clusters), linewidth=2, cmap='Paired')
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
print("Cluster memberships:\n{}".format(y_pred))
Out[59]:
Cluster memberships:
[9 2 5 4 2 7 9 6 9 6 1 0 2 6 1 9 3 0 3 1 7 6 8 6 8 5 2 7 5 8 9 8 6 5 3 7 0
 9 4 5 0 1 3 5 2 8 9 1 5 6 1 0 7 4 6 3 3 6 3 8 0 4 2 9 6 4 8 2 8 4 0 4 0 5
 6 4 5 9 3 0 7 8 0 7 5 8 9 8 0 7 3 9 7 1 7 2 2 0 4 5 6 7 8 9 4 5 4 1 2 3 1 
 8 8 4 9 2 3 7 0 9 9 1 5 8 5 1 9 5 6 7 9 1 4 0 6 2 6 4 7 9 5 5 3 8 1 9 5 6
 3 5 0 2 9 3 0 8 6 0 3 3 5 6 3 2 0 2 3 0 2 6 3 4 4 1 5 6 7 1 1 3 2 4 7 2 7
 3 8 6 4 1 4 3 9 9 5 1 7 5 8 2]

在這裡插入圖片描述

圖 3-32:利用 k 均值的許多簇來表示複雜資料集中的變化

我們使用了 10 個簇中心,也就是說,現在每個點都被分配了 0 到 9 之間的一個數字。我們 可以將其看作 10 個分量表示的資料(我們有 10 個新特徵),只有表示該點對應的簇中心的 那個特徵不為 0,其他特徵均為 0。利用這個 10 維表示,現在可以用線性模型來劃分兩個 半月形,而利用原始的兩個特徵是不可能做到這一點的。將到每個簇中心的距離作為特徵, 還可以得到一種表現力更強的資料表示。可以利用 kmeans 的 transform 方法來完成這一點:

In[60]:
distance_features = kmeans.transform(X)
print("Distance feature shape: {}".format(distance_features.shape))
print("Distance features:\n{}".format(distance_features))
Out[60]:
Distance feature shape: (200, 10)
Distance features:
[[ 0.922 1.466 1.14 ..., 1.166 1.039 0.233]
 [ 1.142 2.517 0.12 ..., 0.707 2.204 0.983]
 [ 0.788 0.774 1.749 ..., 1.971 0.716 0.944]
 ...,
 [ 0.446 1.106 1.49 ..., 1.791 1.032 0.812]
 [ 1.39 0.798 1.981 ..., 1.978 0.239 1.058]
 [ 1.149 2.454 0.045 ..., 0.572 2.113 0.882]]

k 均值是非常流行的聚類演算法,因為它不僅相對容易理解和實現,而且執行速度也相對較 快。k 均值可以輕鬆擴充套件到大型資料集,scikit-learn 甚至在 MiniBatchKMeans 類中包含了 一種更具可擴充套件性的變體,可以處理非常大的資料集。

k 均值的缺點之一在於,它依賴於隨機初始化,也就是說,演算法的輸出依賴於隨機種子。 預設情況下,scikit-learn 用 10 種不同的隨機初始化將演算法執行 10 次,並返回最佳結果。 k 均值還有一個缺點,就是對簇形狀的假設的約束性較強,而且還要求指定所要尋找 的簇的個數(在現實世界的應用中可能並不知道這個數字)。 接下來,我們將學習另外兩種聚類演算法,它們都在某些方面對這些性質做了改進。

3.5.2 凝聚聚類

凝聚聚類(agglomerative clustering)指的是許多基於相同原則構建的聚類演算法,這一原則 是:演算法首先宣告每個點是自己的簇,然後合併兩個最相似的簇,直到滿足某種停止準則 為止。scikit-learn 中實現的停止準則是簇的個數,因此相似的簇被合併,直到僅剩下指 定個數的簇。還有一些連結(linkage)準則,規定如何度量“最相似的簇”。這種度量總 是定義在兩個現有的簇之間。

scikit-learn 中實現了以下三種選項。

ward 預設選項。ward 挑選兩個簇來合併,使得所有簇中的方差增加最小。這通常會得到大 小差不多相等的簇。

average 連結將簇中所有點之間平均距離最小的兩個簇合並。

complete 連結(也稱為最大連結)將簇中點之間最大距離最小的兩個簇合並。

ward 適用於大多數資料集,在我們的例子中將使用它。如果簇中的成員個數非常不同(比 如其中一個比其他所有都大得多),那麼 average 或 complete 可能效果更好。

圖 3-33 給出了在一個二維資料集上的凝聚聚類過程,要尋找三個簇。

In[61]:
mglearn.plots.plot_agglomerative_algorithm()

在這裡插入圖片描述

圖 3-33:凝聚聚類用迭代的方式合併兩個最近的簇

最開始,每個點自成一簇。然後在每一個步驟中,相距最近的兩個簇被合併。在前四個步 驟中,選出兩個單點簇並將其合併成兩點簇。在步驟 5(Step 5)中,其中一個兩點簇被 擴充套件到三個點,以此類推。在步驟 9(Step 9)中,只剩下 3 個簇。由於我們指定尋找 3 個 簇,因此演算法結束。

我們來看一下凝聚聚類對我們這裡使用的簡單三簇資料的效果如何。由於演算法的工作原 理,凝聚演算法不能對新資料點做出預測。因此 AgglomerativeClustering 沒有 predict 方 法。為了構造模型並得到訓練集上簇的成員關係,可以改用 fit_predict 方法。5 結果如圖 3-34 所示。

In[62]:
from sklearn.cluster import AgglomerativeClustering
X, y = make_blobs(random_state=1)
agg = AgglomerativeClustering(n_clusters=3)
assignment = agg.fit_predict(X)
mglearn.discrete_scatter(X[:, 0], X[:, 1], assignment)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")

在這裡插入圖片描述

圖 3-34:使用 3 個簇的凝聚聚類的簇分配

正如所料,演算法完美地完成了聚類。雖然凝聚聚類的 scikit-learn 實現需要你指定希望 演算法找到的簇的個數,但凝聚聚類方法為選擇正確的個數提供了一些幫助,我們將在下 面討論。

  1. 層次聚類與樹狀圖

凝聚聚類生成了所謂的層次聚類(hierarchical clustering)。聚類過程迭代進行,每個點都從一個單點簇變為屬於最終的某個簇。每個中間步驟都提供了資料的一種聚類(簇的個數 也不相同)。有時候,同時檢視所有可能的聚類是有幫助的。下一個例子(圖 3-35)疊加 顯示了圖 3-33 中所有可能的聚類,有助於深入瞭解每個簇如何分解為較小的簇:

In[63]:
mglearn.plots.plot_agglomerative()

在這裡插入圖片描述

圖 3-35:凝聚聚類生成的層次化的簇分配(用線表示)以及帶有編號的資料點(參見圖 3-36)

雖然這種視覺化為層次聚類提供了非常詳細的檢視,但它依賴於資料的二維性質,因此不 能用於具有兩個以上特徵的資料集。但還有另一個將層次聚類視覺化的工具,叫作樹狀圖 (dendrogram),它可以處理多維資料集。

不幸的是,目前 scikit-learn 沒有繪製樹狀圖的功能。但你可以利用 SciPy 輕鬆生成樹狀 圖。SciPy 的聚類演算法介面與 scikit-learn 的聚類演算法稍有不同。SciPy 提供了一個函式, 接受資料陣列 X 並計算出一個連結陣列(linkage array),它對層次聚類的相似度進行編碼。 然後我們可以將這個連結陣列提供給 scipy 的 dendrogram 函式來繪製樹狀圖(圖 3-36)。

In[64]:
# 從SciPy中匯入dendrogram函式和ward聚類函式
from scipy.cluster.hierarchy import dendrogram, ward
X, y = make_blobs(random_state=0, n_samples=12)
# 將ward聚類應用於資料陣列X
# SciPy的ward函式返回一個陣列,指定執行凝聚聚類時跨越的距離
linkage_array = ward(X)
# 現在為包含簇之間距離的linkage_array繪製樹狀圖
dendrogram(linkage_array)
# 在樹中標記劃分成兩個簇或三個簇的位置
ax = plt.gca()
bounds = ax.get_xbound()
ax.plot(bounds, [7.25, 7.25], '--', c='k')
ax.plot(bounds, [4, 4], '--', c='k')
ax.text(bounds[1], 7.25, ' two clusters', va='center', fontdict={'size': 15})
ax.text(bounds[1], 4, ' three clusters', va='center', fontdict={'size': 15})
plt.xlabel("Sample index")
plt.ylabel("Cluster distance")

在這裡插入圖片描述

圖 3-36:圖 3-35 中聚類的樹狀圖(用線表示劃分成兩個簇和三個簇)

樹狀圖在底部顯示資料點(編號從 0 到 11)。然後以這些點(表示單點簇)作為葉節點繪 制一棵樹,每合併兩個簇就新增一個新的父節點。

從下往上看,資料點 1 和 4 首先被合併(正如你在圖 3-33 中所見)。接下來,點 6 和 9 被 合併為一個簇,以此類推。在頂層有兩個分支,一個由點 11、0、5、10、7、6 和 9 組成, 另一個由點 1、4、3、2 和 8 組成。這對應於圖中左側兩個最大的簇。

樹狀圖的 y 軸不僅說明凝聚演算法中兩個簇何時合併,每個分支的長度還表示被合併的簇之 間的距離。在這張樹狀圖中,最長的分支是用標記為“three clusters”(三個簇)的虛線表 示的三條線。它們是最長的分支,這表示從三個簇到兩個簇的過程中合併了一些距離非常 遠的點。我們在影像上方再次看到這一點,將剩下的兩個簇合併為一個簇也需要跨越相對 較大的距離。

不幸的是,凝聚聚類仍然無法分離像 two_moons 資料集這樣複雜的形狀。但我們要學習的 下一個演算法 DBSCAN 可以解決這個問題。

3.5.3 DBSCAN

另一個非常有用的聚類演算法是 DBSCAN(density-based spatial clustering of applications with noise,即“具有噪聲的基於密度的空間聚類應用”)。DBSCAN 的主要優點是它不需要使用者先驗地設定簇的個數,可以劃分具有複雜形狀的簇,還可以找出不屬於任何簇的點。 DBSCAN 比凝聚聚類和 k 均值稍慢,但仍可以擴充套件到相對較大的資料集。

DBSCAN 的原理是識別特徵空間的“擁擠”區域中的點,在這些區域中許多資料點靠近在 一起。這些區域被稱為特徵空間中的密集(dense)區域。DBSCAN 背後的思想是,簇形 成資料的密集區域,並由相對較空的區域分隔開。

在密集區域內的點被稱為核心樣本(core sample,或核心點),它們的定義如下。DBSCAN 有兩個引數:min_samples 和 eps。如果在距一個給定資料點 eps 的距離內至少有 min_ samples 個資料點,那麼這個資料點就是核心樣本。DBSCAN 將彼此距離小於 eps 的核心 樣本放到同一個簇中。

演算法首先任意選取一個點,然後找到到這個點的距離小於等於 eps 的所有的點。如果 距起始點的距離在 eps 之內的資料點個數小於 min_samples,那麼這個點被標記為噪 聲(noise),也就是說它不屬於任何簇。如果距離在 eps 之內的資料點個數大於 min_ samples,則這個點被標記為核心樣本,並被分配一個新的簇標籤。然後訪問該點的所有 鄰居(在距離 eps 以內)。如果它們還沒有被分配一個簇,那麼就將剛剛建立的新的簇標 籤分配給它們。如果它們是核心樣本,那麼就依次訪問其鄰居,以此類推。簇逐漸增大, 直到在簇的 eps 距離內沒有更多的核心樣本為止。然後選取另一個尚未被訪問過的點, 並重復相同的過程。

最後,一共有三種型別的點:核心點、與核心點的距離在 eps 之內的點(叫作邊界點, boundary point)和噪聲。如果 DBSCAN 演算法在特定資料集上多次執行,那麼核心點的聚 類始終相同,同樣的點也始終被標記為噪聲。但邊界點可能與不止一個簇的核心樣本相 鄰。因此,邊界點所屬的簇依賴於資料點的訪問順序。一般來說只有很少的邊界點,這種 對訪問順序的輕度依賴並不重要。

我們將 DBSCAN 應用於演示凝聚聚類的模擬資料集。與凝聚聚類類似,DBSCAN 也不允 許對新的測試資料進行預測,所以我們將使用 fit_predict 方法來執行聚類並返回簇標籤。

In[65]:
from sklearn.cluster import DBSCAN
X, y = make_blobs(random_state=0, n_samples=12)
dbscan = DBSCAN()
clusters = dbscan.fit_predict(X)
print("Cluster memberships:\n{}".format(clusters))
Out[65]:
Cluster memberships:
[-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1]

如你所見,所有資料點都被分配了標籤 -1,這代表噪聲。這是 eps 和 min_samples 預設參 數設定的結果,對於小型的玩具資料集並沒有調節這些引數。min_samples 和 eps 取不同值 時的簇分類如下所示,其視覺化結果見圖 3-37。

In[66]:
mglearn.plots.plot_dbscan()
Out[66]:
min_samples: 2 eps: 1.000000 cluster: [-1 0 0 -1 0 -1 1 1 0 1 -1 -1]
min_samples: 2 eps: 1.500000 cluster: [0 1 1 1 1 0 2 2 1 2 2 0]
min_samples: 2 eps: 2.000000 cluster: [0 1 1 1 1 0 0 0 1 0 0 0]
min_samples: 2 eps: 3.000000 cluster: [0 0 0 0 0 0 0 0 0 0 0 0]
min_samples: 3 eps: 1.000000 cluster: [-1 0 0 -1 0 -1 1 1 0 1 -1 -1]
min_samples: 3 eps: 1.500000 cluster: [0 1 1 1 1 0 2 2 1 2 2 0]
min_samples: 3 eps: 2.000000 cluster: [0 1 1 1 1 0 0 0 1 0 0 0]
min_samples: 3 eps: 3.000000 cluster: [0 0 0 0 0 0 0 0 0 0 0 0]
min_samples: 5 eps: 1.000000 cluster: [-1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1]
min_samples: 5 eps: 1.500000 cluster: [-1 0 0 0 0 -1 -1 -1 0 -1 -1 -1]
min_samples: 5 eps: 2.000000 cluster: [-1 0 0 0 0 -1 -1 -1 0 -1 -1 -1]
min_samples: 5 eps: 3.000000 cluster: [0 0 0 0 0 0 0 0 0 0 0 0]

在這裡插入圖片描述

圖 3-37:在 min_samples 和 eps 引數不同取值的情況下,DBSCAN 找到的簇分配

在這張圖中,屬於簇的點是實心的,而噪聲點則顯示為空心的。核心樣本顯示為較大的標 記,而邊界點則顯示為較小的標記。增大 eps(在圖中從左到右),更多的點會被包含在一 個簇中。這讓簇變大,但可能也會導致多個簇合併成一個。增大 min_samples(在圖中從 上到下),核心點會變得更少,更多的點被標記為噪聲。

引數 eps 在某種程度上更加重要,因為它決定了點與點之間“接近”的含義。將 eps 設定 得非常小,意味著沒有點是核心樣本,可能會導致所有點都被標記為噪聲。將 eps 設定得 非常大,可能會導致所有點形成單個簇。

設定 min_samples 主要是為了判斷稀疏區域內的點被標記為異常值還是形成自己的簇。如 果增大 min_samples,任何一個包含少於 min_samples 個樣本的簇現在將被標記為噪聲。因 此,min_samples 決定簇的最小尺寸。在圖 3-37 中 eps=1.5 時,從 min_samples=3 到 min_ samples=5,你可以清楚地看到這一點。min_samples=3 時有三個簇:一個包含 4 個點,一 個包含 5 個點,一個包含 3 個點。min_samples=5 時,兩個較小的簇(分別包含 3 個點和 4 個點)現在被標記為噪聲,只保留包含 5 個樣本的簇。

雖然 DBSCAN 不需要顯式地設定簇的個數,但設定 eps 可以隱式地控制找到的簇的個數。 使用 StandardScaler 或 MinMaxScaler 對資料進行縮放之後,有時會更容易找到 eps 的較好 取值,因為使用這些縮放技術將確保所有特徵具有相似的範圍。

圖 3-38 展示了在 two_moons 資料集上執行 DBSCAN 的結果。利用預設設定,演算法找到了 兩個半圓形並將其分開:

In[67]:
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)
# 將資料縮放成平均值為0、方差為1
scaler = StandardScaler()
scaler.fit(X)
X_scaled = scaler.transform(X)
dbscan = DBSCAN()
clusters = dbscan.fit_predict(X_scaled)
# 繪製簇分配
plt.scatter(X_scaled[:, 0], X_scaled[:, 1], c=clusters, cmap=mglearn.cm2, s=60)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")

在這裡插入圖片描述

圖 3-38:利用預設值 eps=0.5 的 DBSCAN 找到的簇分配

由於演算法找到了我們想要的簇的個數(2 個),因此引數設定的效果似乎很好。如果將 eps 減小到 0.2(預設值為 0.5),我們將會得到 8 個簇,這顯然太多了。將 eps 增大到 0.7 則 會導致只有一個簇。

在使用 DBSCAN 時,你需要謹慎處理返回的簇分配。如果使用簇標籤對另一個資料進行 索引,那麼使用 -1 表示噪聲可能會產生意料之外的結果。

3.5.4 聚類演算法的對比與評估

在應用聚類演算法時,其挑戰之一就是很難評估一個演算法的效果好壞,也很難比較不同演算法 的結果。在討論完 k 均值、凝聚聚類和 DBSCAN 背後的演算法之後,下面我們將在一些現 實世界的資料集上比較它們。

  1. 用真實值評估聚類

有一些指標可用於評估聚類演算法相對於真實聚類的結果,其中最重要的是調整 rand 指數 (adjusted rand index,ARI)和歸一化互資訊(normalized mutual information,NMI),二者 都給出了定量的度量,其最佳值為 1,0 表示不相關的聚類(雖然 ARI 可以取負值)。

下面我們使用 ARI 來比較 k 均值、凝聚聚類和 DBSCAN 演算法。為了對比,我們還新增了 將點隨機分配到兩個簇中的影像(見圖 3-39)。

In[68]:
from sklearn.metrics.cluster import adjusted_rand_score
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)
# 將資料縮放成平均值為0、方差為1
scaler = StandardScaler()
scaler.fit(X)
X_scaled = scaler.transform(X)
fig, axes = plt.subplots(1, 4, figsize=(15, 3),
 subplot_kw={'xticks': (), 'yticks': ()})
# 列出要使用的演算法
algorithms = [KMeans(n_clusters=2), AgglomerativeClustering(n_clusters=2),
 DBSCAN()]
# 建立一個隨機的簇分配,作為參考
random_state = np.random.RandomState(seed=0)
random_clusters = random_state.randint(low=0, high=2, size=len(X))
# 繪製隨機分配
axes[0].scatter(X_scaled[:, 0], X_scaled[:, 1], c=random_clusters,
 cmap=mglearn.cm3, s=60)
axes[0].set_title("Random assignment - ARI: {:.2f}".format(
 adjusted_rand_score(y, random_clusters)))
for ax, algorithm in zip(axes[1:], algorithms):
 # 繪製簇分配和簇中心
 clusters = algorithm.fit_predict(X_scaled)
 ax.scatter(X_scaled[:, 0], X_scaled[:, 1], c=clusters,
  cmap=mglearn.cm3, s=60)
 ax.set_title("{} - ARI: {:.2f}".format(algorithm.__class__.__name__,
 adjusted_rand_score(y, clusters)))

在這裡插入圖片描述

圖 3-39:利用監督 ARI 分數在 two_moons 資料集上比較隨機分配、k 均值、凝聚聚類和 DBSCAN

調整 rand 指數給出了符合直覺的結果,隨機簇分配的分數為 0,而 DBSCAN(完美地找到 了期望中的聚類)的分數為 1。

用這種方式評估聚類時,一個常見的錯誤是使用 accuracy_score 而不是 adjusted_rand_ score、normalized_mutual_info_score 或其他聚類指標。使用精度的問題在於,它要求分 配的簇標籤與真實值完全匹配。但簇標籤本身毫無意義——唯一重要的是哪些點位於同一 個簇中。

In[69]:
from sklearn.metrics import accuracy_score
# 這兩種點標籤對應於相同的聚類
clusters1 = [0, 0, 1, 1, 0]
clusters2 = [1, 1, 0, 0, 1]
# 精度為0,因為二者標籤完全不同
print("Accuracy: {:.2f}".format(accuracy_score(clusters1, clusters2)))
# 調整rand分數為1,因為二者聚類完全相同
print("ARI: {:.2f}".format(adjusted_rand_score(clusters1, clusters2)))
Out[69]:
Accuracy: 0.00
ARI: 1.00
  1. 在沒有真實值的情況下評估聚類

我們剛剛展示了一種評估聚類演算法的方法,但在實踐中,使用諸如 ARI 之類的指標有一個 很大的問題。在應用聚類演算法時,通常沒有真實值來比較結果。如果我們知道了資料的正 確聚類,那麼可以使用這一資訊構建一個監督模型(比如分類器)。因此,使用類似 ARI 和 NMI 的指標通常僅有助於開發演算法,但對評估應用是否成功沒有幫助。

有一些聚類的評分指標不需要真實值,比如輪廓係數(silhouette coeffcient)。但它們在實 踐中的效果並不好。輪廓分數計算一個簇的緊緻度,其值越大越好,最高分數為 1。雖然 緊緻的簇很好,但緊緻度不允許複雜的形狀。

下面是一個例子,利用輪廓分數在 two_moons 資料集上比較 k 均值、凝聚聚類和 DBSCAN (圖 3-40):

In[70]:
from sklearn.metrics.cluster import silhouette_score
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)
# 將資料縮放成平均值為0、方差為1
scaler = StandardScaler()
scaler.fit(X)
X_scaled = scaler.transform(X)
fig, axes = plt.subplots(1, 4, figsize=(15, 3),
 subplot_kw={'xticks': (), 'yticks': ()})
# 建立一個隨機的簇分配,作為參考
random_state = np.random.RandomState(seed=0)
random_clusters = random_state.randint(low=0, high=2, size=len(X))
# 繪製隨機分配
axes[0].scatter(X_scaled[:, 0], X_scaled[:, 1], c=random_clusters,
 cmap=mglearn.cm3, s=60)
axes[0].set_title("Random assignment: {:.2f}".format(
 silhouette_score(X_scaled, random_clusters)))
algorithms = [KMeans(n_clusters=2), AgglomerativeClustering(n_clusters=2),
 DBSCAN()]
for ax, algorithm in zip(axes[1:], algorithms):
 clusters = algorithm.fit_predict(X_scaled)
 # 繪製簇分配和簇中心
 ax.scatter(X_scaled[:, 0], X_scaled[:, 1], c=clusters, cmap=mglearn.cm3,
 s=60)
 ax.set_title("{} : {:.2f}".format(algorithm.__class__.__name__,
 silhouette_score(X_scaled, clusters)))

在這裡插入圖片描述

圖 3-40:利用無監督的輪廓分數在 two_moons 資料集上比較隨機分配、k 均值、凝聚聚類和 DBSCAN(更符合直覺的 DBSCAN 的輪廓分數低於 k 均值找到的分配)

如你所見,k 均值的輪廓分數最高,儘管我們可能更喜歡 DBSCAN 的結果。對於評估聚 類,稍好的策略是使用基於魯棒性的(robustness-based)聚類指標。這種指標先向資料中 新增一些噪聲,或者使用不同的引數設定,然後執行演算法,並對結果進行比較。其思想 是,如果許多演算法引數和許多資料擾動返回相同的結果,那麼它很可能是可信的。不幸的 是,在寫作本書時,scikit-learn 還沒有實現這一策略。

即使我們得到一個魯棒性很好的聚類或者非常高的輪廓分數,但仍然不知道聚類中是否有任何語義含義,或者聚類是否反映了資料中我們感興趣的某個方面。我們回到人臉影像的 例子。我們希望找到類似人臉的分組,比如男人和女人、老人和年輕人,或者有鬍子的人 和沒鬍子的人。假設我們將資料分為兩個簇,關於哪些點應該被聚類在一起,所有演算法的 結果一致。我們仍不知道找到的簇是否以某種方式對應於我們感興趣的概念。演算法找到的 可能是側檢視和正面檢視、夜間拍攝的照片和白天拍攝的照片,或者 iPhone 拍攝的照片和 安卓手機拍攝的照片。要想知道聚類是否對應於我們感興趣的內容,唯一的辦法就是對簇 進行人工分析。

  1. 在人臉資料集上比較演算法

我們將 k 均值、DBSCAN 和凝聚聚類演算法應用於 Wild 資料集中的 Labeled Faces,並查 看它們是否找到了有趣的結構。我們將使用資料的特徵臉表示,它由包含 100 個成分的 PCA(whiten=True) 生成:

In[71]:
# 從lfw資料中提取特徵臉,並對資料進行變換
from sklearn.decomposition import PCA
pca = PCA(n_components=100, whiten=True, random_state=0)
pca.fit_transform(X_people)
X_pca = pca.transform(X_people)

我們之前見到,與原始畫素相比,這是對人臉影像的一種語義更強的表示。它的計算速度 也更快。這裡有一個很好的練習,就是在原始資料上執行下列實驗,不要用 PCA,並觀察 你是否能找到類似的簇。

用 DBSCAN 分析人臉資料集。我們首先應用剛剛討論過的 DBSCAN:

In[72]:
# 應用預設引數的DBSCAN
dbscan = DBSCAN()
labels = dbscan.fit_predict(X_pca)
print("Unique labels: {}".format(np.unique(labels)))
Out[72]:
Unique labels: [-1]

我們看到,所有返回的標籤都是 -1,因此所有資料都被 DBSCAN 標記為“噪聲”。我們可 以改變兩個引數來改進這一點:第一,我們可以增大 eps,從而擴充套件每個點的鄰域;第二, 我們可以減小 min_samples,從而將更小的點組視為簇。我們首先嚐試改變 min_samples:

In[73]:
dbscan = DBSCAN(min_samples=3)
labels = dbscan.fit_predict(X_pca)
print("Unique labels: {}".format(np.unique(labels)))
Out[73]:
Unique labels: [-1]

即使僅考慮由三個點構成的組,所有點也都被標記為噪聲。因此我們需要增大 eps:

In[74]:
dbscan = DBSCAN(min_samples=3, eps=15)
labels = dbscan.fit_predict(X_pca)
print("Unique labels: {}".format(np.unique(labels)))
Out[74]:
Unique labels: [-1 0]

使用更大的 eps(其值為 15),我們只得到了單一簇和噪聲點。我們可以利用這一結果找出 “噪聲”相對於其他資料的形狀。為了進一步理解發生的事情,我們檢視有多少點是噪聲, 有多少點在簇內:

In[75]:
# 計算所有簇中的點數和噪聲中的點數。
# bincount不允許負值,所以我們需要加1。
# 結果中的第一個數字對應於噪聲點。
print("Number of points per cluster: {}".format(np.bincount(labels + 1)))
Out[75]:
Number of points per cluster: [ 27 2036]

噪聲點非常少——只有 27 個,因此我們可以檢視所有的噪聲點(見圖 3-41):

In[76]:
noise = X_people[labels==-1]
fig, axes = plt.subplots(3, 9, subplot_kw={'xticks': (), 'yticks': ()},
 figsize=(12, 4))
for image, ax in zip(noise, axes.ravel()):
 ax.imshow(image.reshape(image_shape), vmin=0, vmax=1)

在這裡插入圖片描述

圖 3-41:人臉資料集中被 DBSCAN 標記為噪聲的樣本

將這些影像與圖 3-7 中隨機選擇的人臉影像樣本進行比較,我們可以猜測它們被標記為 噪聲的原因:第 1 行第 5 張影像顯示一個人正在用玻璃杯喝水,還有人戴帽子的影像, 在最後一張影像中,人臉前面有一隻手。其他影像都包含奇怪的角度,或者太近或太寬 的剪下。

這 種 類 型 的 分 析 —— 嘗 試 找 出“ 奇 怪 的 那 一 個 ” —— 被 稱 為 異常值檢測(outlier detection)。如果這是一個真實的應用,那麼我們可能會嘗試更好地裁切影像,以得到更加 均勻的資料。對於照片中的人有時戴著帽子、喝水或在面前舉著某物,我們能做的事情很 少。但需要知道它們是資料中存在的問題,我們應用任何演算法都需要解決這些問題。

如果我們想要找到更有趣的簇,而不是一個非常大的簇,那麼需要將 eps 設定得更小,取 值在 15 和 0.5(預設值)之間。我們來看一下 eps 不同取值對應的結果:

In[77]:
for eps in [1, 3, 5, 7, 9, 11, 13]:
 print("\neps={}".format(eps))
 dbscan = DBSCAN(eps=eps, min_samples=3)
 labels = dbscan.fit_predict(X_pca)
 print("Clusters present: {}".format(np.unique(labels)))
 print("Cluster sizes: {}".format(np.bincount(labels + 1)))
Out[77]:
eps=1
Clusters present: [-1]
Cluster sizes: [2063]
eps=3
Clusters present: [-1]
Cluster sizes: [2063]
eps=5
Clusters present: [-1]
Cluster sizes: [2063]
eps=7
Clusters present: [-1 0 1 2 3 4 5 6 7 8 9 10 11 12]
Cluster sizes: [2006 4 6 6 6 9 3 3 4 3 3 3 3 4]
eps=9
Clusters present: [-1 0 1 2]
Cluster sizes: [1269 788 3 3]
eps=11
Clusters present: [-1 0]
Cluster sizes: [ 430 1633]
eps=13
Clusters present: [-1 0]
Cluster sizes: [ 112 1951]

對於較小的 eps,所有點都被標記為噪聲。eps=7 時,我們得到許多噪聲點和許多較小的 簇。eps=9 時,我們仍得到許多噪聲點,但我們得到了一個較大的簇和一些較小的簇。從 eps=11 開始,我們僅得到一個較大的簇和噪聲。

有趣的是,較大的簇從來沒有超過一個。最多有一個較大的簇包含大多數點,還有一些較 小的簇。這表示資料中沒有兩類或三類非常不同的人臉影像,而是所有影像或多或少地都 與其他影像具有相同的相似度(或不相似度)。

eps=7 的結果看起來最有趣,它有許多較小的簇。我們可以通過將 13 個較小的簇中的點全 部視覺化來深入研究這一聚類(圖 3-42):

In[78]:
dbscan = DBSCAN(min_samples=3, eps=7)
labels = dbscan.fit_predict(X_pca)
for cluster in range(max(labels) + 1):
 mask = labels == cluster
 n_images = np.sum(mask)
 fig, axes = plt.subplots(1, n_images, figsize=(n_images * 1.5, 4),
 subplot_kw={'xticks': (), 'yticks': ()})
 for image, label, ax in zip(X_people[mask], y_people[mask], axes):
 ax.imshow(image.reshape(image_shape), vmin=0, vmax=1)
 ax.set_title(people.target_names[label].split()[-1])

在這裡插入圖片描述

圖 3-42:eps=7 的 DBSCAN 找到的簇

有一些簇對應於(這個資料集中)臉部非常不同的人,比如 Sharon(沙龍)或 Koizumi(小泉)。在每個簇內,人臉方向和麵部表情也是固定的。有些簇中包含多個人的面孔,但 他們的方向和表情都相似。

這就是我們將 DBSCAN 演算法應用於人臉資料集的分析結論。如你所見,我們這裡進行了 人工分析,不同於監督學習中基於 R2 分數或精度的更為自動化的搜尋方法。

下面我們將繼續應用 k 均值和凝聚聚類。

用 k 均值分析人臉資料集。我們看到,利用 DBSCAN 無法建立多於一個較大的簇。凝聚 聚類和 k 均值更可能建立均勻大小的簇,但我們需要設定簇的目標個數。我們可以將簇的 數量設定為資料集中的已知人數,雖然無監督聚類演算法不太可能完全找到它們。相反,我 們可以首先設定一個比較小的簇的數量,比如 10 個,這樣我們可以分析每個簇:

In[79]:
# 用k均值提取簇
km = KMeans(n_clusters=10, random_state=0)
labels_km = km.fit_predict(X_pca)
print("Cluster sizes k-means: {}".format(np.bincount(labels_km)))
Out[79]:
Cluster sizes k-means: [269 128 170 186 386 222 237 64 253 148]

如你所見,k 均值聚類將資料劃分為大小相似的簇,其大小在 64 和 386 之間。這與 DBSCAN 的結果非常不同。

我們可以通過將簇中心視覺化來進一步分析 k 均值的結果(圖 3-43)。由於我們是在 PCA 生成的表示中進行聚類,因此我們需要使用 pca.inverse_transform 將簇中心旋轉回到原 始空間並視覺化:

In[80]:
fig, axes = plt.subplots(2, 5, subplot_kw={'xticks': (), 'yticks': ()},
 figsize=(12, 4))
for center, ax in zip(km.cluster_centers_, axes.ravel()):
 ax.imshow(pca.inverse_transform(center).reshape(image_shape),
 vmin=0, vmax=1)

在這裡插入圖片描述

圖 3-43:將簇的數量設定為 10 時,k 均值找到的簇中心

k 均值找到的簇中心是非常平滑的人臉。這並不奇怪,因為每個簇中心都是 64 到 386 張人 臉影像的平均。使用降維的 PCA 表示,可以增加影像的平滑度(對比圖 3-11 中利用 100個 PCA 維度重建的人臉)。聚類似乎捕捉到人臉的不同方向、不同表情(第 3 個簇中心似 乎顯示的是一張笑臉),以及是否有襯衫領子(見倒數第二個簇中心)。

圖 3-44 給出了更詳細的檢視,我們對每個簇中心給出了簇中 5 張最典型的影像(該簇中與 簇中心距離最近的影像)與 5 張最不典型的影像(該簇中與簇中心距離最遠的影像):

In[81]:
mglearn.plots.plot_kmeans_faces(km, pca, X_pca, X_people,
 y_people, people.target_names)

在這裡插入圖片描述

圖 3-44:k 均值為每個簇找到的樣本影像——簇中心在最左邊,然後是五個距中心最近的點,然後 是五個距該簇距中心最遠的點

圖 3-44 證實了我們認為第 3 個簇是笑臉的直覺,也證實了其他簇中方向的重要性。不過 “非典型的”點與簇中心不太相似,而且它們的分配似乎有些隨意。這可以歸因於以下事 實:k 均值對所有資料點進行劃分,不像 DBSCAN 那樣具有“噪聲”點的概念。利用更多 數量的簇,演算法可以找到更細微的區別。但新增更多的簇會使得人工檢查更加困難。

用凝聚聚類分析人臉資料集。下面我們來看一下凝聚聚類的結果:

In[82]:
# 用ward凝聚聚類提取簇
agglomerative = AgglomerativeClustering(n_clusters=10)
labels_agg = agglomerative.fit_predict(X_pca)
print("Cluster sizes agglomerative clustering: {}".format(
 np.bincount(labels_agg)))
Out[82]:
Cluster sizes agglomerative clustering: [255 623 86 102 122 199 265 26 230 155]

凝聚聚類生成的也是大小相近的簇,其大小在 26 和 623 之間。這比 k 均值生成的簇更不 均勻,但比 DBSCAN 生成的簇要更加均勻。

我們可以通過計算 ARI 來度量凝聚聚類和 k 均值給出的兩種資料劃分是否相似:

In[83]:
print("ARI: {:.2f}".format(adjusted_rand_score(labels_agg, labels_km)))
Out[83]:
ARI: 0.13

ARI 只有 0.13,說明 labels_agg 和 labels_km 這兩種聚類的共同點很少。這並不奇怪,原 因在於以下事實:對於 k 均值,遠離簇中心的點似乎沒有什麼共同點。

下面,我們可能會想要繪製樹狀圖(圖 3-45)。我們將限制圖中樹的深度,因為如果分支 到 2063 個資料點,影像將密密麻麻無法閱讀:

In[84]:
linkage_array = ward(X_pca)
# 現在我們為包含簇之間距離的linkage_array繪製樹狀圖
plt.figure(figsize=(20, 5))
dendrogram(linkage_array, p=7, truncate_mode='level', no_labels=True)
plt.xlabel("Sample index")
plt.ylabel("Cluster distance")

在這裡插入圖片描述

圖 3-45:凝聚聚類在人臉資料集上的樹狀圖

要想建立 10 個簇,我們在頂部有 10 條豎線的位置將樹橫切。在圖 3-36 所示的玩具資料 的樹狀圖中,你可以從分支的長度中看出,兩個或三個簇就可以很好地劃分資料。對於 人臉資料而言,似乎沒有非常自然的切割點。有一些分支代表更為不同的組,但似乎沒 有一個特別合適的簇的數量。這並不奇怪,因為 DBSCAN 的結果是試圖將所有的點都聚類在一起。

我們將 10 個簇視覺化,正如之前對 k 均值所做的那樣(圖 3-46)。請注意,在凝聚聚類中 沒有簇中心的概念(雖然我們計算平均值),我們只是給出了每個簇的前幾個點。我們在 第一張影像的左側給出了每個簇中的點的數量:

In[85]:
n_clusters = 10
for cluster in range(n_clusters):
 mask = labels_agg == cluster
 fig, axes = plt.subplots(1, 10, subplot_kw={'xticks': (), 'yticks': ()},
 figsize=(15, 8))
 axes[0].set_ylabel(np.sum(mask))
 for image, label, asdf, ax in zip(X_people[mask], y_people[mask],
 labels_agg[mask], axes):
 ax.imshow(image.reshape(image_shape), vmin=0, vmax=1)
 ax.set_title(people.target_names[label].split()[-1],
 fontdict={'fontsize': 9})

在這裡插入圖片描述

圖 3-46:In[82] 生成的簇中的隨機影像——每一行對應一個簇,左側的數字表示每個簇中影像的數量

雖然某些簇似乎具有語義上的主題,但許多簇都太大而實際上很難是均勻的。為了得到更加均 勻的簇,我們可以再次執行演算法,這次使用 40 個簇,並挑選出一些特別有趣的簇(圖 3-47):

In[86]:
# 用ward凝聚聚類提取簇
agglomerative = AgglomerativeClustering(n_clusters=40)
labels_agg = agglomerative.fit_predict(X_pca)
print("cluster sizes agglomerative clustering: {}".format(np.bincount(labels_agg)))
n_clusters = 40
for cluster in [10, 13, 19, 22, 36]: # 手動挑選“有趣的”簇
 mask = labels_agg == cluster
 fig, axes = plt.subplots(1, 15, subplot_kw={'xticks': (), 'yticks': ()},
 figsize=(15, 8))
 cluster_size = np.sum(mask)
 axes[0].set_ylabel("#{}: {}".format(cluster, cluster_size))
 for image, label, asdf, ax in zip(X_people[mask], y_people[mask],
 labels_agg[mask], axes):
 ax.imshow(image.reshape(image_shape), vmin=0, vmax=1)
 ax.set_title(people.target_names[label].split()[-1],
 fontdict={'fontsize': 9})
 for i in range(cluster_size, 15):
 axes[i].set_visible(False)
Out[86]:
cluster sizes agglomerative clustering:
 [ 58 80 79 40 222 50 55 78 172 28 26 34 14 11 60 66 152 27
 47 31 54 5 8 56 3 5 8 18 22 82 37 89 28 24 41 40
 21 10 113 69]

在這裡插入圖片描述

圖 3-47:將簇的數量設定為 40 時,從凝聚聚類找到的簇中挑選的影像——左側文字表示簇的編號 和簇中的點的總數

這裡聚類挑選出的似乎是“深色皮膚且微笑”“有領子的襯衫”“微笑的女性”“薩達姆” 和“高額頭”。如果進一步詳細分析,我們還可以利用樹狀圖找到這些高度相似的簇。

3.5.5 聚類方法小結

本節的內容表明,聚類的應用與評估是一個非常定性的過程,通常在資料分析的探索階 段很有幫助。我們學習了三種聚類演算法:k 均值、DBSCAN 和凝聚聚類。這三種演算法都 可以控制聚類的粒度(granularity)。k 均值和凝聚聚類允許你指定想要的簇的數量,而 DBSCAN 允許你用 eps 引數定義接近程度,從而間接影響簇的大小。三種方法都可以用於 大型的現實世界資料集,都相對容易理解,也都可以聚類成多個簇。

每種演算法的優點稍有不同。k 均值可以用簇的平均值來表示簇。它還可以被看作一種分解 方法,每個資料點都由其簇中心表示。DBSCAN 可以檢測到沒有分配任何簇的“噪聲點”, 還可以幫助自動判斷簇的數量。與其他兩種方法不同,它允許簇具有複雜的形狀,正如我 們在 two_moons 的例子中所看到的那樣。DBSCAN 有時會生成大小差別很大的簇,這可能 是它的優點,也可能是缺點。凝聚聚類可以提供資料的可能劃分的整個層次結構,可以通 過樹狀圖輕鬆檢視。

3.6 小結與展望

本章介紹了一系列無監督學習演算法,可用於探索性資料分析和預處理。找到資料的正確表 示對於監督學習和無監督學習的成功通常都至關重要,預處理和分解方法在資料準備中具 有重要作用。

分解、流形學習和聚類都是加深資料理解的重要工具,在沒有監督資訊的情況下,也是理 解資料的僅有的方法。即使是在監督學習中,探索性工具對於更好地理解資料性質也很重 要。通常來說,很難量化無監督演算法的有用性,但這不應該妨礙你使用它們來深入理解數 據。學完這些方法,你就已經掌握了機器學習從業者每天使用的所有必要的學習演算法。

我們建議你在 scikit-learn 中包含的二維玩具資料和現實世界資料集(比如 digits、iris 和 cancer 資料集)上嘗試聚類和分解方法。

相關文章