在《多元統計分析——聚類分析——K-均值聚類（K-中值、K-衆數）》當中，我們了解了K-均值聚類的原理，也簡單的介紹了K-均值聚類的兩個應用場景：
- 發現異常情況：如果不對資料進行任何形式的轉換，隻是經過中心标準化或級差标準化就進行快速聚類，會根據資料分布特征得到聚類結果。這種聚類會将極端資料單獨聚為幾類。這種方法适用于統計分析之前的異常值剔除，對異常行為的挖掘，比如監控銀行賬戶是否有洗錢行為、監控POS機是有從事套現、監控某個終端是否是電話卡養卡客戶等等。
注：差別于單個次元的異常值：單個變量的異常值一般是離中心點位置超過3倍或者5倍的标準差。
- 将個案資料做劃分：出于客戶細分目的的聚類分析一般希望聚類結果為大緻平均的幾大類（原始資料盡量服從正态分布，這樣聚類出來的簇的樣本點個數大緻接近），是以需要将資料進行轉換，比如使用原始變量的百分位秩、Turkey正态評分、對數轉換等等。在這類分析中資料的具體數值并沒有太多的意義，重要的是相對位置。這種方法适用場景包括客戶消費行為聚類、客戶積分使用行為聚類等等。
以上兩種場景的大緻步驟如下：
現在我們通過案例來分析這兩個場景。
案例：有一份電信使用者的資料集，字段包括：
1、導入資料
資料樣本大體如下：
%matplotlib inline
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
profile_telecom = pd.read_csv('profile_telecom.csv')
profile_telecom.head()
輸出：
2、變量的相關系數矩陣
去除ID列，輸出後四列的相關系數矩陣。
profile_telecom.loc[: ,'cnt_call':].corr()
輸出：
可以看出有些變量之間的相關系數還是比較高的。
3、檢測變量分布
plt.figure(figsize=(8, 3))
for i in range(4):
plt.subplot(220 + i + 1)
plt.hist(profile_telecom.iloc[:, i + 1], bins=20)
plt.show()
輸出：
可以看出，各變量都是一個右偏的分布，在實際工作當中，遇到最多的一般就是正态分布和右偏分布。
4、場景一：檢測異常情況——發現離群點
4.1、變量歸一化——學生标準化
僅使用學生标準化進行預處理，不使用其它會改變資料分布形态的變換。
from sklearn.preprocessing import scale
from sklearn.cluster import KMeans
from sklearn.decomposition import PCA
tele_scaled = scale(profile_telecom.loc[:, 'cnt_call':]) #學生标準化
tele_scaled
輸出：
array([[-0.296453 , 3.68086843, 0.85513907, 2.08171356],
[-0.1905409 , -0.80292101, -0.59755491, -0.64747691],
[-0.56879838, 0.31802635, -0.39579186, -0.08281681],
...,
[ 0.24823778, 0.01231344, 3.2762957 , 2.45815362],
[-0.11488941, 0.01231344, -0.11332358, -0.08281681],
[ 1.01988304, -0.80292101, -0.59755491, -0.64747691]])
4.2、降維——主成分分析
tele_pca = PCA(n_components=2) #PCA降維（2維）
tele_pca_score = tele_pca.fit_transform(tele_scaled)
print('variance_ratio:', tele_pca.explained_variance_ratio_)
輸出：
variance_ratio: [0.62510442 0.24620209]
前兩個主成分，方差占比之和為87%， 是以降到二維是比較合适的。
4.3、K-means聚類
k = 4
tele_kmeans = KMeans(n_clusters=k, n_init=15).fit(tele_pca_score) #聚成四類
tele_kmeans.cluster_centers_ #輸出質心
輸出：
array([[ 5.10963437, -0.21710494],
[ 1.30272379, -0.10395259],
[-0.8669528 , -0.10911188],
[-0.06943007, 4.0981177 ]])
4.4、聚類特征分析
4.4.1、輪廓系數
from sklearn.metrics import silhouette_score
silhouette_score(tele_pca_score, tele_kmeans.labels_)
輸出：
0.5081846457150366
4.4.2、降維之後聚類效果
plt.figure(figsize=[4, 3])
markers = 'xvo+*^dDhs|_<,.>'
k = 4
for cluster, marker in zip(range(k), markers[:k]):
x_axis = tele_pca_score[:, 0][tele_kmeans.labels_ == cluster]
y_axis = tele_pca_score[:, 1][tele_kmeans.labels_ == cluster]
plt.scatter(x_axis, y_axis, marker=marker)
plt.show()
輸出：
在沒有進行資料分布形态轉換的情況之下，比較明顯的看出紅框中為疑似離群點。注意：k=4是我們随意選擇的數字，k我們可以選大一些，如7、8、9...，這樣異常點檢測的效果越明顯。聚類完成之後，可以對資料樣本按聚類标簽進行分組，檢視每個标簽下都有多少個樣本，如下：
4.4.3、聚類之後各類的數量
pd.DataFrame(tele_pca_score).groupby(tele_kmeans.labels_).count()
輸出：
我們會發現，有某些簇的樣本量較小，即疑似異常值點。
5、場景二：将個案資料做劃分
5.1、變量分布轉換——取對數
上面我們已經知道，原始變量都是呈右偏分布的。
# 對變量取對數（也可取rank）
log_telecom = np.log1p(profile_telecom.iloc[:, 1:]) #有很多0值，0值取對數是無意義的，是以+1之後取對數。
plt.figure(figsize=(8, 3))
for i in range(4):
plt.subplot(220 + i + 1)
plt.hist(log_telecom.iloc[:, i], bins=20)
plt.show()
輸出：
右偏分布的資料集為了讓其盡量靠近正态，處理的方法一般是取對數。結合實際業務，變量中有很多是0值，我們先将其+1，因為對0取對數是沒有意義的。
5.2、變量歸一化——學生标準化
log_telecom_scaled=scale(log_telecom)
log_telecom_scaled
輸出：
array([[ 0.04174954, 1.98292009, 1.26385741, 1.60145475],
[ 0.19291257, -1.24657303, -1.11310713, -0.63398781],
[-0.48397217, 0.76013283, 0.06635716, 0.40350962],
...,
[ 0.6609339 , 0.50039436, 1.89829497, 1.71268561],
[ 0.28932817, 0.50039436, 0.57532607, 0.40350962],
[ 1.18246098, -1.24657303, -1.11310713, -0.63398781]])
注意，是否進行标準化，可以結合實際的業務場景，如本例題中，四個變量的量綱是一樣的（機關都為“次數”） ，但是其各變量之間的方差差異是比較大的，是以需要進行标準化。否則後續的PCA效果可能不太好。
不需要标準化是有兩個先決條件的：1、量綱一緻；2、方差相近。
5.3、降維——主成分分析
log_pca = PCA(n_components=2, whiten=True)
log_pca_score = log_pca.fit_transform(log_telecom_scaled)
print('variance_ratio:', log_pca.explained_variance_ratio_)
輸出：
variance_ratio: [0.71197124 0.22755356]
前兩個主成分方差占比超過90%，取兩個主成分正好。前兩個主成分（各變量的權重）如下：
components=log_pca.components_ #獲得轉換後的所有主成分
components
輸出：
array([[ 0.22077367, 0.54471581, 0.5687719 , 0.57536266],
[ 0.97202048, -0.19050882, -0.08855497, -0.10507422]])
5.4、使用輪廓系數或者樣本到類中心的距離和(離差平方和)确定聚類數量
plt.figure(figsize=[8, 2])
Ks = range(2, 10)
rssds = []; silhs = []
for k in Ks:
model = KMeans(n_clusters=k, n_init=15)
model.fit(log_pca_score)
rssds.append(model.inertia_)
silhs.append(silhouette_score(log_pca_score, model.labels_, sample_size=None))
# 輪廓系數計算複雜度高，使用sample_size
plt.subplot(121); plt.plot(Ks, rssds)
plt.subplot(122); plt.plot(Ks, silhs)
plt.show()
輸出：
左邊為樣本到類中心的距離和(離差平方和)，右邊為輪廓系數，結合兩個名額，我們确定聚成3類是最合适的。
5.5、聚類特征分析
5.5.1、輪廓系數
k = 3
log_pca_kmeans = KMeans(n_clusters=k, n_init=15).fit(log_pca_score)
silhouette_score(log_pca_score, log_pca_kmeans.labels_)
輸出：
0.4454977689225694
輪廓系數0.445，效果還行。
5.5.2、降維之後聚類效果
plt.figure(figsize=[4, 3])
for cluster, marker in zip(range(k), markers[:k]):
x_axis = log_pca_score[:, 0][log_pca_kmeans.labels_ == cluster]
y_axis = log_pca_score[:, 1][log_pca_kmeans.labels_ == cluster]
plt.scatter(x_axis, y_axis, marker=marker)
plt.show()
輸出：
我們發現3個簇之間，樣本點相對是比較均勻的。可以根據聚類标簽檢視各類的樣本數。
5.5.3、聚類之後各類的數量
pd.DataFrame(log_pca_score).groupby(log_pca_kmeans.labels_).count()
輸出：
從上面可以看出，各類是相對均衡的。
5.5.4、解釋模型（特征分析）——使用原始資料
結合聚類結果，根據簇的标簽進行分組，看原始資料各簇的質心（平均值）。
co = profile_telecom.iloc[:, 1:5].groupby(log_pca_kmeans.labels_).mean()
co
輸出：
由上面我們可知：第一類：打電話最多，其他項目很少；第二類：所有項目用的都不多；第三類：所有項目都有用，且處于正常水準。
結合上面主成分分析輸出的主成分：
array([[ 0.22077367, 0.54471581, 0.5687719 , 0.57536266],
[ 0.97202048, -0.19050882, -0.08855497, -0.10507422]])
第一個主成分，各變量的權重相對接近， 第二個主成分，第一個變量的權重較高。
通過以下的可視化更加直覺：
co.T.plot(figsize=[4, 3])
plt.show()
輸出：
結論：在電信行業，打電話為基本業務，其他的為增值業務。第0簇：打電話較多，可進行語音包推薦；第1簇：所有項目都不多，屬于低端使用者，語音為剛性需求，可激活一些流量的需求，如閑時的流量包；第2簇：屬于正常使用者。
注意：以上的根據“标簽分組算平均值”的方法适用于變量少的情況，對于變量次元較多的，我們對聚類結果的特征分析時，可以運用到決策樹，從上往下，探索各類使用者最主要的特征，如下：
關于決策樹，詳見《機器學習——有監督——決策樹（分類樹）相關原理及sklearn實作（資訊熵、基尼系數、資訊增益、特征重要程度的量化）》。