機器學習---用python實現樸素貝葉斯算法（Machine Learning Naive Bayes Algorithm Application）

在《機器學習---樸素貝葉斯分類器（Machine Learning Naive Bayes Classifier）》一文中，咱們介紹了樸素貝葉斯分類器的原理。如今，讓咱們來實踐一下。

 

在這裏，咱們使用一份皮馬印第安女性的醫學數據，用來預測其是否會得糖尿病。文件一共有768個樣本，咱們先剔除缺失值，而後選出20%的樣本做爲測試樣本。

 

文件下載地址：https://raw.githubusercontent.com/jbrownlee/Datasets/master/pima-indians-diabetes.data.csv

 

特徵分別是：

懷孕次數 
口服葡萄糖耐量試驗中血漿葡萄糖濃度
舒張壓（mm Hg） 
三頭肌組織褶厚度（mm） 
2小時血清胰島素（μU/ ml） 
體重指數（kg/（身高(m)）^ 2） 
糖尿病系統功能 
年齡（歲）

（注：特徵值爲0表示缺失值）

 

標籤是：

是否患有糖尿病，0表明沒有糖尿病，1表明患有糖尿病

 

---

 

在開始以前，先回顧一下樸素貝葉斯算法：

 

假設訓練集特徵數爲i，記做x_i，目標爲k個類別，記爲C_k，樣本數爲n，新樣本特徵記做x_{new_sample}。

1，計算出先驗機率P(C_k)，即每一個類別在訓練集中的機率；

2，分別計算出訓練集中每一個特徵在每一個類別下的條件機率P(x_i|C_k)，具體分爲如下三種狀況；

      a）若是特徵數據是離散值，那麼咱們假設其符合多項式分佈，代入公式；

      b）若是特徵數據是布爾類型的離散值，那麼咱們假設其符合伯努利分佈，代入公式；

      c）若是特徵數據是連續值，那麼咱們假設其符合高斯分佈。只須要求出這個特徵在每一個類別下的平均值μ_k和方差σ_k²，而後代入公式；

3，對於須要預測類別的新樣本，分別計算每一個類別下的極大後驗機率：（咱們在實際計算中爲防止數據下溢，將連乘運算取對數變成相加運算）；

4，對上述極大後驗機率進行比較，最大的那個對應的類別即爲結果；

 

---

 

首先導入pandas和numpy庫，讀取csv文件，而後去除缺失值：

import pandas as pd import numpy as np data=pd.read_csv(r"C:\Users\ccav\pima-indians-diabetes.data.csv",header=None,\ names=["懷孕次數","口服葡萄糖耐量試驗中血漿葡萄糖濃度",\ "舒張壓（mm Hg）","三頭肌組織褶厚度（mm）",\ "2小時血清胰島素（μU/ ml）","體重指數（kg/（身高(m)）^ 2）",\ "糖尿病系統功能","年齡（歲）","是否患有糖尿病"]) data.iloc[:,0:8]=data.iloc[:,0:8].applymap(lambda x:np.NaN if x==0 else x) #把屬性值爲0的地方轉換成NaN
 data=data.dropna(how="any",axis=0) #去除有缺失值的行

 

通過處理還只剩336個有效樣本：

print(len(data)) 336

 

經過查看，咱們發如今這336個有效樣本中，225個是沒有糖尿病人的樣本，111個是患有糖尿病人的樣本，顯然原數據裏兩個類別的樣本是不均衡的。這個問題能夠經過不一樣的採樣方法解決，可是鑑於簡便起見，咱們隨機選取大約20%的樣本用於測試：

#隨機選取80%的樣本做爲訓練樣本
data_train=data.sample(frac=0.8,random_state=4,axis=0) #剩下的做爲測試樣本
test_idx=[i for i in data.index.values if i not in data_train.index.values] data_test=data.loc[test_idx,:]

 

因爲咱們一共有8個特徵，2個類別，所以，咱們一共須要計算16個條件機率以及2個先驗機率，而後根據這些數據計算出2個極大後驗機率。該如何計算這些數據，又該如何合併，若是不事先想清楚怎麼作，很容易被搞暈。此外，考慮到數據是連續值，咱們假設其符合正態分佈，所以在訓練時須要計算出每一個特徵的平均值和方差。

 

讓咱們先來理清一下思路：

1，按類別分隔數據

2，獲取類別總數和類別名稱

3，訓練數據：計算每一個類別的先驗機率，計算每一個類別每一個特徵的平均值和方差

4，預測：計算每一個類別每一個特徵的條件機率，計算每一個類別的極大後驗機率併合並，得出最大可能的類別名稱

 

按照這個步驟，咱們發現如下幾項都是須要複用的：按類別分隔數據，計算每一個類別的先驗機率，計算每一個類別每一個特徵的平均值和方差，計算每一個類別每一個特徵的條件機率。所以，咱們把這幾項分別作成function。

 

按類別分隔數據：

def SepByClass(X, y): ###按類別分隔數據###
    ###輸入未分類的特徵和目標，輸出分類完成的數據（字典形式）###
    num_of_samples=len(y) #總樣本數
 y=y.reshape(X.shape[0],1) data=np.hstack((X,y)) #把特徵和目標合併成完整數據
 data_byclass={} #初始化分類數據，爲一個空字典
    #提取各種別數據，字典的鍵爲類別名，值爲對應的分類數據
    for i in range(len(data[:,-1])): if i in data[:,-1]: data_byclass[i]=data[data[:,-1]==i] class_name=list(data_byclass.keys()) #類別名
        num_of_class=len(data_byclass.keys()) #類別總數
        
    return data_byclass

 

計算每一個類別的先驗機率：

def CalPriorProb(y_byclass): ###計算y的先驗機率（使用拉普拉斯平滑）###
    ###輸入當前類別下的目標，輸出該目標的先驗機率###
    #計算公式：（當前類別下的樣本數+1）/（總樣本數+類別總數）
    return (len(y_byclass)+1)/(num_of_samples+num_of_class)

 

計算每一個類別每一個特徵的平均值和方差：

def CalXMean(X_byclass): ###計算各種別特徵各維度的平均值###
    ###輸入當前類別下的特徵，輸出該特徵各個維度的平均值###
    X_mean=[] for i in range(X_byclass.shape[1]): X_mean.append(np.mean(X_byclass[:,i])) return X_mean def CalXVar(X_byclass): ###計算各種別特徵各維度的方差###
    ###輸入當前類別下的特徵，輸出該特徵各個維度的方差###
    X_var=[] for i in range(X_byclass.shape[1]): X_var.append(np.var(X_byclass[:,i])) return X_var

 

計算每一個類別每一個特徵的條件機率：

def CalGaussianProb(X_new, mean, var): ###計算訓練集特徵（符合正態分佈）在各種別下的條件機率###
    ###輸入新樣本的特徵，訓練集特徵的平均值和方差，輸出新樣本的特徵在相應訓練集中的分佈機率###
    #計算公式：(np.exp(-(X_new-mean)**2/(2*var)))*(1/np.sqrt(2*np.pi*var))
    gaussian_prob=[] for a,b,c in zip(X_new, mean, var): formula1=np.exp(-(a-b)**2/(2*c)) formula2=1/np.sqrt(2*np.pi*c) gaussian_prob.append(formula2*formula1) return gaussian_prob

 

接下來，咱們開始訓練數據。首先按類別分隔數據，而後遍歷每個類別，分別計算每一個類別的先驗機率和每一個類別每一個特徵的平均值和方差，並儲存在列表中。

def fit(X, y): ###訓練數據###
    ###輸入訓練集特徵和目標，輸出目標的先驗機率，特徵的平均值和方差###
    #將輸入的X,y轉換爲numpy數組
    X, y = np.asarray(X, np.float32), np.asarray(y, np.float32) data_byclass=Gaussian_NB.SepByClass(X,y) #將數據分類
        #計算各種別數據的目標先驗機率，特徵平均值和方差
    for data in data_byclass.values(): X_byclass=data[:,:-1] y_byclass=data[:,-1] prior_prob.append(Gaussian_NB.CalPriorProb(y_byclass)) X_mean.append(Gaussian_NB.CalXMean(X_byclass)) X_var.append(Gaussian_NB.CalXVar(X_byclass)) return prior_prob, X_mean, X_var

 

最後，輸入一個新樣本的特徵，預測其所屬的類別。首先，遍歷以前在訓練數據時計算出的每一個類別的先驗機率和每一個類別每一個特徵的平均值和方差，再以此計算條件機率和極大後驗機率。選出最大的極大後驗機率所對應的索引，最後提取其對應的類別名稱。

def predict(X_new): ###預測數據###
    ###輸入新樣本的特徵，輸出新樣本最有可能的目標###
    #將輸入的x_new轉換爲numpy數組
    X_new=np.asarray(X_new, np.float32) posteriori_prob=[] #初始化極大後驗機率
        
    for i,j,o in zip(prior_prob, X_mean, X_var): gaussian=Gaussian_NB.CalGaussianProb(X_new,j,o) posteriori_prob.append(np.log(i)+sum(np.log(gaussian))) idx=np.argmax(posteriori_prob) return class_name[idx]

 

---

 

整理一下以上代碼，咱們把這個高斯樸素貝葉斯分類器作成一個類，完整代碼以下：

import pandas as pd import numpy as np data=pd.read_csv(r"C:\Users\ccav\pima-indians-diabetes.data.csv",header=None,\ names=["懷孕次數","口服葡萄糖耐量試驗中血漿葡萄糖濃度",\ "舒張壓（mm Hg）","三頭肌組織褶厚度（mm）",\ "2小時血清胰島素（μU/ ml）","體重指數（kg/（身高(m)）^ 2）",\ "糖尿病系統功能","年齡（歲）","是否患有糖尿病"]) data.iloc[:,0:8]=data.iloc[:,0:8].applymap(lambda x:np.NaN if x==0 else x) #把屬性值爲0的地方轉換成NaN
 data=data.dropna(how="any",axis=0) #去除有缺失值的行

#隨機選取80%的樣本做爲訓練樣本
data_train=data.sample(frac=0.8,random_state=4,axis=0) #剩下的做爲測試樣本
test_idx=[i for i in data.index.values if i not in data_train.index.values] data_test=data.loc[test_idx,:] #提取訓練集和測試集的特徵和目標
X_train=data_train.iloc[:,:-1] y_train=data_train.iloc[:,-1] X_test=data_test.iloc[:,:-1] y_test=data_test.iloc[:,-1] class Gaussian_NB: def __init__(self): self.num_of_samples = None self.num_of_class = None self.class_name = [] self.prior_prob = [] self.X_mean = [] self.X_var = [] def SepByClass(self, X, y): ###按類別分隔數據###
        ###輸入未分類的特徵和目標，輸出分類完成的數據（字典形式）###
        self.num_of_samples=len(y) #總樣本數
 y=y.reshape(X.shape[0],1) data=np.hstack((X,y)) #把特徵和目標合併成完整數據
 data_byclass={} #初始化分類數據，爲一個空字典
        #提取各種別數據，字典的鍵爲類別名，值爲對應的分類數據
        for i in range(len(data[:,-1])): if i in data[:,-1]: data_byclass[i]=data[data[:,-1]==i] self.class_name=list(data_byclass.keys()) #類別名
        self.num_of_class=len(data_byclass.keys()) #類別總數
        
        return data_byclass def CalPriorProb(self, y_byclass): ###計算y的先驗機率（使用拉普拉斯平滑）###
        ###輸入當前類別下的目標，輸出該目標的先驗機率###
        #計算公式：（當前類別下的樣本數+1）/（總樣本數+類別總數）
        return (len(y_byclass)+1)/(self.num_of_samples+self.num_of_class) def CalXMean(self, X_byclass): ###計算各種別特徵各維度的平均值###
        ###輸入當前類別下的特徵，輸出該特徵各個維度的平均值###
        X_mean=[] for i in range(X_byclass.shape[1]): X_mean.append(np.mean(X_byclass[:,i])) return X_mean def CalXVar(self, X_byclass): ###計算各種別特徵各維度的方差###
        ###輸入當前類別下的特徵，輸出該特徵各個維度的方差###
        X_var=[] for i in range(X_byclass.shape[1]): X_var.append(np.var(X_byclass[:,i])) return X_var def CalGaussianProb(self, X_new, mean, var): ###計算訓練集特徵（符合正態分佈）在各種別下的條件機率###
        ###輸入新樣本的特徵，訓練集特徵的平均值和方差，輸出新樣本的特徵在相應訓練集中的分佈機率###
        #計算公式：(np.exp(-(X_new-mean)**2/(2*var)))*(1/np.sqrt(2*np.pi*var))
        gaussian_prob=[] for a,b,c in zip(X_new, mean, var): formula1=np.exp(-(a-b)**2/(2*c)) formula2=1/np.sqrt(2*np.pi*c) gaussian_prob.append(formula2*formula1) return gaussian_prob def fit(self, X, y): ###訓練數據###
        ###輸入訓練集特徵和目標，輸出目標的先驗機率，特徵的平均值和方差###
        #將輸入的X,y轉換爲numpy數組
        X, y = np.asarray(X, np.float32), np.asarray(y, np.float32) data_byclass=Gaussian_NB.SepByClass(X,y) #將數據分類
        #計算各種別數據的目標先驗機率，特徵平均值和方差
        for data in data_byclass.values(): X_byclass=data[:,:-1] y_byclass=data[:,-1] self.prior_prob.append(Gaussian_NB.CalPriorProb(y_byclass)) self.X_mean.append(Gaussian_NB.CalXMean(X_byclass)) self.X_var.append(Gaussian_NB.CalXVar(X_byclass)) return self.prior_prob, self.X_mean, self.X_var def predict(self,X_new): ###預測數據###
        ###輸入新樣本的特徵，輸出新樣本最有可能的目標###
        #將輸入的x_new轉換爲numpy數組
        X_new=np.asarray(X_new, np.float32) posteriori_prob=[] #初始化極大後驗機率
        
        for i,j,o in zip(self.prior_prob, self.X_mean, self.X_var): gaussian=Gaussian_NB.CalGaussianProb(X_new,j,o) posteriori_prob.append(np.log(i)+sum(np.log(gaussian))) idx=np.argmax(posteriori_prob) return self.class_name[idx] if __name__=="__main__": Gaussian_NB=Gaussian_NB() #實例化Gaussian_NB
    Gaussian_NB.fit(X_train,y_train) #使用Gaussian_NB模型訓練數據
    acc=0 TP=0 FP=0 FN=0 for i in range(len(X_test)): predict=Gaussian_NB.predict(X_test.iloc[i,:]) target=np.array(y_test)[i] if predict==1 and target==1: TP+=1
        if predict==0 and target==1: FP+=1
        if predict==target: acc+=1
        if predict==1 and target==0: FN+=1
    print("準確率:",acc/len(X_test)) print("查準率:",TP/(TP+FP)) print("查全率:",TP/(TP+FN)) print("F1:",2*TP/(2*TP+FP+FN))

 

結果以下：

準確率: 0.7611940298507462 查準率: 0.7391304347826086 查全率: 0.6296296296296297 F1: 0.68

因爲樣本不均衡，所以光看準確率是不夠的，在這裏一併計算了查準率和查全率以及F1值。能夠看出這個分類器的效果不是很好，一方面是由於樣本數量較少，另外一方面樸素貝葉斯的預測能力確實比不上覆雜模型，只能提供一個粗略的判斷。