利用 Python學習資料探勘【1】
覆蓋使用Python進行資料探勘查詢和描述資料結構模式的實踐工具。
第一節
介紹
資料探勘是一個隱式提取以前未知的潛在有用的資料資訊提取方式。它使用廣泛，並且是眾多應用的技術基礎。
本文介紹那些使用Python資料探勘實踐用於發現和描述結構模式資料的工具。近些年來，Python在開發以資料為中心的應用中被用的越來越多。感謝大型科學計算社群的支援以及大大豐富的資料分析函式庫。尤其是，我們可以看到如何：
• 匯入和視覺化資料
• 資料分類
• 使用迴歸分析和相關測量法發現數據之間的關係
• 資料降維以壓縮和視覺化資料帶來的資訊
• 分析結構化資料
每個主題都會提供程式碼例項，它們基於四個主要的Python資料分析和處理的類庫：numpy，matplotlib，sklearn和networkx。
第二節
資料匯入和視覺化
通常，資料分析的第一步由獲取資料和匯入資料到我們的工作環境組成。我們可以使用以下的Python程式碼簡單的下載資料：
Python
| 1 2 3 4 5 6 |
import urllib2 url = 'http://aima.cs.berkeley.edu/data/iris.csv' u = urllib2.urlopen(url) localFile = open('iris.csv'', 'w') localFile.write(u.read()) localFile.close() |
在以上的程式碼片段中，我們使用了urllib2類庫以獲取伯克利大學網站的一個檔案，並使用標準類庫提供的File物件把它儲存到本地磁碟。資料包含鳶尾花（iris）資料集，這是一個包含了三種鳶尾花（山鳶尾、維吉尼亞鳶尾和變色鳶尾）的各50個數據樣本的多元資料集，每個樣本都有四個特徵（或者說變數），即花萼（sepal）和花瓣（petal）的長度和寬度。以釐米為單位。
資料集以CSV(逗號分割值）的格式儲存。CSV檔案可以很方便的轉化並把其中的資訊儲存為適合的資料結構。此資料集有5列（譯者注：原文是行，但這裡應該是列的意思），前4列包含著特徵值，最後一列代表著樣本型別。CSV檔案很容易被numpy類庫的genfromtxt方法解析：
Python
| 1 2 3 4 5 |
from numpy import genfromtxt, zeros # read the first 4 columns data = genfromtxt('iris.csv',delimiter=',',usecols=(0,1,2,3)) # read the fifth column target = genfromtxt('iris.csv',delimiter=',',usecols=(4),dtype=str) |
在上面的例子中我們建立了一個包含特徵值的矩陣以及一個包含樣本型別的向量。我們可以通過檢視我們載入的資料結構的shape值來確認資料集的大小：
Python
| 1 2 3 4 |
print data.shape (150, 4) print target.shape (150,) |
我們也可以檢視我們有多少種樣本型別以及它們的名字：
Python
| 1 2 |
print set(target) # build a collection of unique elements set(['setosa', 'versicolor', 'virginica']) |
當我們處理新資料的時候，一項很重要的任務是嘗試去理解資料包含的資訊以及它的組織結構。視覺化可以靈活生動的展示資料，幫助我們深入理解資料。
使用pylab類庫（matplotlib的介面）的plotting方法可以建一個二維散點圖讓我們在兩個維度上分析資料集的兩個特徵值：
Python
| 1 2 3 4 5 |
from pylab import plot, show plot(data[target=='setosa',0],data[target=='setosa',2],'bo') plot(data[target=='versicolor',0],data[target=='versicolor',2],'ro') plot(data[target=='virginica',0],data[target=='virginica',2],'go') show() |
上面那段程式碼使用第一和第三維度（花萼的長和寬），結果如下圖所示： 
在上圖中有150個點，不同的顏色代表不同的型別；藍色點代表山鳶尾，紅色點代表變色鳶尾，綠色點代表維吉尼亞鳶尾。
另一種常用的檢視資料的方法是分特性繪製直方圖。在本例中，既然資料被分為三類，我們就可以比較每一類的分佈特徵。下面這個程式碼可以繪製資料中每一型別的第一個特性（花萼的長度）：
Python
| 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 |
from pylab import figure, subplot, hist, xlim, show xmin = min(data[:,0]) xmax = max(data[:,0]) figure() subplot(411) # distribution of the setosa class (1st, on the top) hist(data[target=='setosa',0],color='b',alpha=.7) xlim(xmin,xmax) subplot(412) # distribution of the versicolor class (2nd) hist(data[target=='versicolor',0],color='r',alpha=.7) xlim(xmin,xmax) subplot(413) # distribution of the virginica class (3rd) hist(data[target=='virginica',0],color='g',alpha=.7) xlim(xmin,xmax) subplot(414) # global histogram (4th, on the bottom) hist(data[:,0],color='y',alpha=.7) xlim(xmin,xmax) show() |
結果如下圖：
根據上圖的直方圖，我們可以根據資料型別區分理解資料的特徵。例如，我們可以觀察到，山鳶尾的平均花萼長度小於維吉尼亞鳶尾。
第三節
分類
分類是一個數據挖掘方法，用於把一個數據集中的樣本資料分配給各個目標類。實現這個方法的模組叫做分類器。使用分類器需要以下兩步：訓練和分類。訓練是指採集已知其特定類歸屬的資料並基於這些資料建立分類器。 分類是指使用通過這些已知資料建立的分類器來處理未知的資料，以判斷未知資料的分類情況。
Sklearn類庫包含很多分類器的實現，本章我們將會使用高斯樸素貝葉斯來分析我們在第一章載入的鳶尾花資料，包含山鳶尾、變色鳶尾和維吉尼亞鳶尾。最後我們把字串陣列轉型成整型資料：
Python
| 1 2 3 4 |
t = zeros(len(target)) t[target == 'setosa'] = 1 t[target == 'versicolor'] = 2 t[target == 'virginica'] = 3 |
現在我們已經做好例項化和訓練分類器的準備了：
Python
| 1 2 3 |
from sklearn.naive_bayes import GaussianNB classifier = GaussianNB() classifier.fit(data,t) # training on the iris dataset |
分類器可以由predict方法完成，並且只要輸出一個樣例就可以很簡單的檢測：
Python
| 1 2 3 4 |
print classifier.predict(data[0]) [ 1.] print t[0] 1 |
上例中predicted類包含了一個正確的樣本（山鳶尾），但是在廣泛的樣本上評估分類器並且使用非訓練環節的資料測試是很重要的。最終我們通過從源資料集中隨機抽取樣本把資料分為訓練集和測試集。我們將會使用訓練集的資料來訓練分類器，並使用測試集的資料來測試分類器。train_test_split方法正是實現此功能的：
Python
| 1 2 3 |
from sklearn import cross_validation train, test, t_train, t_test = cross_validation.train_test_split(data, t, … test_size=0.4, random_state=0) |
資料集被分一分為二，測試集被指定為源資料的40%（命名為test_size），我們用它反覆訓練我們的分類器並輸出精確度：
Python
| 1 2 3 |
classifier.fit(train,t_train) # train print classifier.score(test,t_test) # test 0.93333333333333335 |
在此例中，我們的精確度為93%。一個分類器的精確度是通過正確分類樣本的數量除以總樣本的數量得出的。也就是說，它意味著我們正確預測的比例。
另一個估計分類器表現的工具叫做混淆矩陣。在此矩陣中每列代表一個預測類的例項，每行代表一個實際類的例項。使用它可以很容易的計算和列印矩陣：
Python
| 1 2 3 4 5 |
from sklearn.metrics import confusion_matrix print confusion_matrix(classifier.predict(test),t_test) [[16 0 0] [ 0 23 3] [ 0 0 18]] |
在這個混淆矩陣中我們可以看到所有山鳶尾和維吉尼亞鳶尾都被正確的分類了，但是實際上應該是26個的變色鳶尾，系統卻預測其中三個是維吉尼亞鳶尾。如果我們牢記所有正確的猜測都在表格的對角線上，那麼觀測表格的錯誤就很容易了，即對角線以外的非零值。
可以展示分類器效能的完整報告的方法也是很好用的：
Python
| 1 2 3 4 5 6 7 |
from sklearn.metrics import classification_report print classification_report(classifier.predict(test), t_test, target_names=['setosa', 'versicolor', 'virginica']) precision recall f1-score support setosa 1.00 1.00 1.00 16 versicolor 1.00 0.85 0.92 27 virginica 0.81 1.00 0.89 17 avg / total 0.95 0.93 0.93 60 |
以下是該報告使用到的方法總結：
Precision：正確預測的比例
Recall（或者叫真陽性率）：正確識別的比例
F1-Score：precision和recall的調和平均數
以上僅僅只是給出用於支撐測試分類的資料量。當然，分割資料、減少用於訓練的樣本數以及評估結果等操作都依賴於配對的訓練集和測試集的隨機選擇。如果要切實評估一個分類器並與其它的分類器作比較的話，我們需要使用一個更加精確的評估模型，例如Cross Validation。該模型背後的思想很簡單：多次將資料分為不同的訓練集和測試集，最終分類器評估選取多次預測的平均值。這次，sklearn為我們提供了執行模型的方法：
Python
| 1 2 3 4 5 |
from sklearn.cross_validation import cross_val_score # cross validation with 6 iterations scores = cross_val_score(classifier, data, t, cv=6) print scores [ 0.84 0.96 1. 1. 1. 0.96] |
如上所見，輸出是每次模型迭代產生的精確度的陣列。我們可以很容易計算出平均精確度：
Python
| 1 2 3 |
from numpy import mean print mean(scores) 0.96 |
第四章
聚類
通常我們的資料上不會有標籤告訴我們它的樣本型別；我們需要分析資料，把資料按照它們的相似度標準分成不同的群組，群組（或者群集）指的是相似樣本的集合。這種分析被稱為無監督資料分析。最著名的聚類工具之一叫做k-means演算法，如下所示：
Python
| 1 2 3 |
from sklearn.cluster import KMeans kmeans = KMeans(k=3, init='random') # initialization kmeans.fit(data) # actual execution |
上述片段執行k-measn演算法並把資料分為三個群集（引數k所指定的）。現在我們可以使用模型把每一個樣本分配到三個群集中：
Python
| 1 |
c = kmeans.predict(data) |
我們可以估計群集的結果，與使用完整性得分和同質性得分計算而得的標籤作比較：
Python
| 1 2 3 4 5 |
from sklearn.metrics import completeness_score, homogeneity_score print completeness_score(t,c) 0.7649861514489815 print homogeneity_score(t,c) 0.7514854021988338 |
當大部分資料點屬於一個給定的類並且屬於同一個群集，那麼完整性得分就趨向於1。當所有群集都幾乎只包含某個單一類的資料點時同質性得分就趨向於1.
我們可以把叢集視覺化並和帶有真實標籤的做視覺化比較：
Python
| 1 2 3 4 5 6 7 8 9 10 |
figure() subplot(211) # top figure with the real classes plot(data[t==1,0],data[t==1,2],'bo') plot(data[t==2,0],data[t==2,2],'ro') plot(data[t==3,0],data[t==3,2],'go') subplot(212) # bottom figure with classes assigned automatically plot(data[c==1,0],data[tt==1,2],'bo',alpha=.7) plot(data[c==2,0],data[tt==2,2],'go',alpha=.7) plot(data[c==0,0],data[tt==0,2],'mo',alpha=.7) show() |
結果如下圖所示：
觀察此圖我們可以看到，底部左側的群集可以被k-means完全識別，然而頂部的兩個群集有部分識別錯誤。
第五章
迴歸
迴歸是一個用於預測變數之間函式關係調查的方法。例如，我們有兩個變數，一個被認為是解釋，一個被認為是依賴。我們希望使用模型描述兩者的關係。當這種關係是一條線的時候就稱為線性迴歸。
為了應用線性迴歸我們建立一個由上所述的綜合資料集：
Python
| 1 2 3 |
from numpy.random import rand x = rand(40,1) # explanatory variable y = x*x*x+rand(40,1)/5 # depentend variable |
我們可以使用在sklear.linear_model模組中發現的LinearRegression模型。該模型可以通過計算每個資料點到擬合線的垂直差的平方和，找到平方和最小的最佳擬合線。使用方法和我們之前遇到的實現sklearn的模型類似：
Python
| 1 2 3 |
from sklearn.linear_model import LinearRegression linreg = LinearRegression() linreg.fit(x,y) |
我們可以通過把擬合線和實際資料點畫在同一幅圖上來評估結果：
Python
| 1 2 3 4 |
from numpy import linspace, matrix xx = linspace(0,1,40) plot(x,y,'o',xx,linreg.predict(matrix(xx).T),'--r') show() |
圖見下：
觀察該圖我們可以得出結論：擬合線從資料點中心穿過，並可以確定是增長的趨勢。
我們還可以使用均方誤差來量化模型和原始資料的擬合度：
Python
| 1 2 3 |
from sklearn.metrics import mean_squared_error print mean_squared_error(linreg.predict(x),y) 0.01093512327489268 |
該指標度量了預期的擬合線和真實資料之間的距離平方。當擬合線很完美時該值為0。
第六章
相關
我們通過研究相關性來理解成對的變數之間是否相關，相關性的強弱。此類分析幫助我們精確定位被依賴的重要變數。最好的相關方法是皮爾遜積矩相關係數。它是由兩個變數的協方差除以它們的標準差的乘積計算而來。我們將鳶尾花資料集的變數兩兩組合計算出其係數如下所示：
Python
| 1 2 3 4 5 6 7 |
from numpy import corrcoef corr = corrcoef(data.T) # .T gives the transpose print corr [[ 1. -0.10936925 0.87175416 0.81795363] [-0.10936925 1. -0.4205161 -0.35654409] [ 0.87175416 -0.4205161 1. 0.9627571 ] [ 0.81795363 -0.35654409 0.9627571 1. ]] |
corrcoef方法通過輸入行為變數列為觀察值的矩陣，計算返回相關係數的對稱矩陣。該矩陣的每個元素代表著兩個變數的相關性。
當值一起增長時相關性為正。當一個值減少而另一個只增加時相關性為負。特別說明，1代表完美的正相關，0代表不相關，-1代表完美的負相關。
當變數數增長時我們可以使用偽彩色點很方便的視覺化相關矩陣：
Python
| 1 2 3 4 5 6 7 8 |
from pylab import pcolor, colorbar, xticks, yticks from numpy import arrange pcolor(corr) colorbar() # add # arranging the names of the variables on the axis xticks(arange(0.5,4.5),['sepal length', 'sepal width', 'petal length', 'petal width'],rotation=-20) yticks(arange(0.5,4.5),['sepal length', 'sepal width', 'petal length', 'petal width'],rotation=-20) show() |
結果如下： 
看圖右側的彩條，我們可以把顏色點關聯到數值上。在本例中，紅色被關聯為最高的正相關，我們可以看出我們資料集的最強相關是“花瓣寬度”和“花瓣長度”這兩個變數。
第七章
降維
在第一章中我們瞭解瞭如何將鳶尾花資料集的兩個維度視覺化。單獨使用該方法，我們只能看到資料集的部分資料檢視。既然我們可以同時繪製的最高維度數為3，將整個資料集嵌入一系列維度並建立一個整體視覺化檢視是很有必要的。這個嵌入過程就被稱作降維。最著名的降維技術之一就是主成分分析（PCA）。該技術把資料變數轉換為等量或更少的不相關變數，稱為主成分（PCs）。
這次，sklearn滿足我們本次分析的所有需求：
Python
| 1 2 |
from sklearn.decomposition import PCA pca = PCA(n_components=2) |
上述片段中我們例項化了一個PCA物件，用於計算前兩個主成分。轉換計算如下：
Python
| 1 |
pcad = pca.fit_transform(data) |
然後如往常一樣繪製結果：
Python
| 1 2 3 4 |
plot(pcad[target=='setosa',0],pcad[target=='setosa',1],'bo') plot(pcad[target=='versicolor',0],pcad[target=='versicolor',1],'ro') plot(pcad[target=='virginica',0],pcad[target=='virginica',1],'go') show() |
結果如下：
可以注意到上圖和第一章提到的有些相似，不過這次變色鳶尾（紅色的）和維吉尼亞鳶尾（綠色的）的間隔更清晰了。
PCA將空間資料方差最大化，我們可以通過方差比判斷PCs包含的資訊量：
Python
| 1 2 |
print pca.explained_variance_ratio_ [ 0.92461621 0.05301557] |
現在我們知道第一個PC佔原始資料的92%的資訊量而第二個佔剩下的5%。我們還可以輸出在轉化過程中丟失的資訊量：
Python
| 1 2 |
print 1-sum(pca.explained_variance_ratio_) 0.0223682249752 |
在本例中我們損失了2%的資訊量。
此時，我們可以是應用逆變換還原原始資料：
Python
| 1 |
data_inv = pca.inverse_transform(pcad) |
可以證明的是，由於資訊丟失逆變換不能給出準確的原始資料。我們可以估算逆變換的結果和原始資料的相似度：
Python
| 1 2 |
print abs(sum(sum(data - data_inv))) 2.8421709430404007e-14 |
可以看出原始資料和逆變換計算出的近似值之間的差異接近於零。通過改變主成分的數值來計算我們能夠覆蓋多少資訊量是很有趣的：
Python
| 1 2 3 4 |
for i in range(1,5): pca = PCA(n_components=i) pca.fit(data) print sum(pca.explained_variance_ratio_) * 100,'%' |
上述片段輸出如下：
Python
| 1 2 3 4 |
92.4616207174 % 97.7631775025 % 99.481691455 % 100.0 % |
PCs用得越多，資訊覆蓋就越全，不過這段分析有助於我們理解儲存一段特定的資訊需要哪些元件。例如，從上述片段可以看出，只要使用三個PCs就可以覆蓋鳶尾花資料集的幾乎100%的資訊。
第八章
網路挖掘
通常我們分析的資料是以網路結構儲存的，例如我，我們的資料可以描述一群facebook使用者的朋友關係或者科學家的論文的合作者關係。這些研究物件可以使用點和邊描述之間的關係。
本章中我們將會介紹分析此類資料的基本步驟，稱為圖論，一個幫助我們創造、處理和研究網路的類庫。尤其我們將會介紹如何使用特定方法建立有意義的資料視覺化，以及如何建立一組關聯稠密的點。
使用圖論可以讓我們很容易的匯入用於描述資料結構的最常用結構：
Python
| 1 |
G = nx.read_gml('lesmiserables.gml',relabel=True) |
在上述程式碼我們匯入了《悲慘世界》同時出現的單片語成的網路，可以通過https://gephi.org/datasets/lesmiserables.gml.zip免費下載，資料以GML格式儲存。我們還可以使用下面的命令匯入並可視化網路：
Python
| 1 |
nx.draw(G,node_size=0,edge_color='b',alpha=.2,font_size=7) |
結果如下：
上圖網路中每個點代表小說中的一個單詞，兩單詞間的聯絡代表同一章裡兩單詞同時出現了。很容易就發現這圖不是很有幫助。該網路的大部分細節依然一藏著，並且發現很難發現那些重要點。我們可以研究節點度來獲取一些內部細節。節點度是指一個最簡中心測量並由一個點的多個關聯組成。我們可以通過最大值、最小值、中位數、第一四分位數和第三四分位數來總結一個網路的度分佈：
Python
| 1 2 3 4 5 6 7 8 9 10 11 12 |
deg = nx.degree(G) from numpy import percentile, mean, median print min(deg.values()) print percentile(deg.values(),25) # computes the 1st quartile print median(deg.values()) print percentile(deg.values(),75) # computes the 3rd quartile print max(deg.values())10 1 2.0 6.0 10.0 36 |
經過分析我們決定只考慮節點度大於10的點。我們可以建立一副只包含我們需要的點的新圖來展現這些點：
Python
| 1 2 3 4 5 6 |
Gt = G.copy() dn = nx.degree(Gt) for n in Gt.nodes(): if dn[n] <= 10: Gt.remove_node(n) nx.draw(Gt,node_size=0,edge_color='b',alpha=.2,font_size=12) |
結果見下圖：
現在這幅圖可讀性大大提高了。這樣我們就能觀察到相關度最高的單詞以及它們的關係。
通過識別團來研究網路也是很有趣的。團是指一個點和其他所有點相連的群組，極大團是指網路中不屬於其它任何一個團的子集的團。我們可以通過如下方法發現我們網路中所有的極大團：
Python
| 1 2 |
from networkx import find_cliques cliques = list(find_cliques(G)) |
用下面這行命令輸出極大團：
Python
| 1 2 |
print max(cliques, key=lambda l: len(l)) [u'Joly', u'Gavroche', u'Bahorel', u'Enjolras', u'Courfeyrac', u'Bossuet', u'Combeferre', u'Feuilly', u'Prouvaire', u'Grantaire'] |
我們可以看見列表中的絕大部分名字都和上圖中可以看見的資料群中的相同。
第九章
其他資源
- IPython and IPython notebook，Python互動式指令碼，也提供很棒的瀏覽器介面。[http://ipython.org/]
- NLTK,自然語言工具包，研究和開發自然語言處理的模組、資料和文件的套裝。[http://nltk.org/]
- OpenCV,影象處理和計算機視覺的最重要的類庫之一。[http://opencv.willowgarage.com/]
- Pandas, 為處理“關聯的”或者“標籤化的”資料提供最快、最自由、最富表現力的資料結構，使資料處理簡單直接。[http://pandas.pydata.org/]
- Scipy，建立在numpy基礎之上，提供大量的高等資料、訊號處理、最優化和統計的演算法和工具。[http://www.scipy.org/]
- Statsmodels, 大量的描述統計學、統計學檢驗、繪圖函式，輸出不同資料和不同統計的結果。[http://statsmodels.sourceforge.net/]