朴素貝葉斯分類器及Python、MATLAB實現

一、問題引入：

設有A，B，C三種花，且它們在自然界的數量都相同，即在這三類中任意取一花，P(A)=P(B)=P(C)=1/3。現有一枝花，問它屬於哪一類，則在沒有任何提示的情況下，可以得知，它是A（或B或C）的可能性一樣。

但，若此時用它們花萼的長度，花萼的寬度，花瓣的長度，花瓣的寬度，即4維向量表示各自的特征，並且這些特征對我們已知。那么，這時它屬於哪一類的不確定度減少，則我們能夠得知它屬於哪一類的概率就大了。

已知某樣本的特征，判斷它是哪一類，就是模式識別的任務，而已知某樣本的特征，得出它屬於這些類的概率，最大者為所屬的類，就是貝葉斯分類的方法。以A為例，利用貝葉斯公式：

其中，為這三類花中，花萼的長度，花萼的寬度，花瓣的長度，花瓣的寬度的總體密度分布，對於三類來說，都是一致的；

P(A)=1/3，稱為先驗概率，在實踐中有已知的統計；

為類條件密度，即A類花的花萼的長度，花萼的寬度，花瓣的長度，花瓣的寬度服從的分布，在朴素貝葉斯分類中，假設該分布密度為4元高斯分布；

稱為后驗概率。

所以，在求解時，我們只需求解其展開公式的分子即可。

 

二、最小錯誤率的朴素貝葉斯分類：

對於，假設該分布密度為4元高斯分布，則為朴素貝葉斯分類，即分布為：

 

計算過程如下：

對展開公式的分子取對數（注：取對數不影響單調性，和極大值的取值），得：

同理，可以計算B,C條件下的結果：，將這三者中的最大值對應的類，作為決策的類別。

 

 

Python代碼:

 

#實驗的數據可以從http://aima.eecs.berkeley.edu/data/ 下載iris文檔。

import numpy

from sklearn.model_selection import train_test_split 

from sklearn.naive_bayes import GaussianNB

data = numpy.genfromtxt('iris.csv', delimiter=',', usecols=(0,1,2,3))

target = numpy.genfromtxt('iris.csv', delimiter=',', usecols=(4), dtype=str) 

#data[0]維花萼的長度，data[1]為花萼的寬，data[2]為花瓣的長，data[3]花瓣的寬，data[4]是標簽。

 t = numpy.zeros(len(target)) 

t[target == 'setosa'] = 1

t[target == 'versicolor'] = 2

t[target == 'virginica'] = 3

#訓練

clf = GaussianNB()

train, test, t_train, t_test = train_test_split(data, t, test_size=0.4, random_state=0)
clf.fit(train, t_train)

#識別 

print(clf.predict([data[3]]))

#print(classifier.predict_proba(data[1]))

#print(classifier.predict_log_proba(data[1]))

print(clf.score(test,t_test))        #計算識別正確率

 

 

MATLAB代碼:

%需要對數據進行提取，剔除標簽（我用Python實現的），得到的數據可直接作為特征向量用於MATLAB中

f1 = open('iris.txt', 'r')
f2 = open('iris_2.txt', 'w')
for i in range(150) :
    line = f1.readline()
    f2.write(line[0:15]+'\n')
f2.close()
f1.close()

clear;clc;

A=[5.1,3.5,1.4,0.2
4.9,3.0,1.4,0.2
4.7,3.2,1.3,0.2
4.6,3.1,1.5,0.2
5.0,3.6,1.4,0.2
5.4,3.9,1.7,0.4
4.6,3.4,1.4,0.3
5.0,3.4,1.5,0.2
4.4,2.9,1.4,0.2
4.9,3.1,1.5,0.1
5.4,3.7,1.5,0.2
4.8,3.4,1.6,0.2
4.8,3.0,1.4,0.1
4.3,3.0,1.1,0.1
5.8,4.0,1.2,0.2
5.7,4.4,1.5,0.4
5.4,3.9,1.3,0.4
5.1,3.5,1.4,0.3
5.7,3.8,1.7,0.3
5.1,3.8,1.5,0.3
5.4,3.4,1.7,0.2
5.2,4.1,1.5,0.1
5.5,4.2,1.4,0.2
4.9,3.1,1.5,0.1
5.0,3.2,1.2,0.2
5.5,3.5,1.3,0.2
4.9,3.1,1.5,0.1
4.4,3.0,1.3,0.2
5.1,3.4,1.5,0.2
5.0,3.5,1.3,0.3
4.5,2.3,1.3,0.3
4.4,3.2,1.3,0.2
5.0,3.5,1.6,0.6
5.1,3.8,1.9,0.4
4.8,3.0,1.4,0.3
5.1,3.8,1.6,0.2
4.6,3.2,1.4,0.2
5.3,3.7,1.5,0.2
5.0,3.3,1.4,0.2
7.0,3.2,4.7,1.4];

B=[6.4,3.2,4.5,1.5
6.9,3.1,4.9,1.5
5.5,2.3,4.0,1.3
6.5,2.8,4.6,1.5
5.7,2.8,4.5,1.3
6.3,3.3,4.7,1.6
4.9,2.4,3.3,1.0
6.6,2.9,4.6,1.3
5.2,2.7,3.9,1.4
5.0,2.0,3.5,1.0
5.9,3.0,4.2,1.5
6.0,2.2,4.0,1.0
6.1,2.9,4.7,1.4
5.6,2.9,3.6,1.3
6.7,3.1,4.4,1.4
5.6,3.0,4.5,1.5
5.8,2.7,4.1,1.0
6.2,2.2,4.5,1.5
5.6,2.5,3.9,1.1
5.9,3.2,4.8,1.8
6.1,2.8,4.0,1.3
6.3,2.5,4.9,1.5
6.1,2.8,4.7,1.2
6.4,2.9,4.3,1.3
6.6,3.0,4.4,1.4
6.8,2.8,4.8,1.4
6.7,3.0,5.0,1.7
6.0,2.9,4.5,1.5
5.7,2.6,3.5,1.0
5.5,2.4,3.8,1.1
5.5,2.4,3.7,1.0
5.8,2.7,3.9,1.2
6.0,2.7,5.1,1.6
5.4,3.0,4.5,1.5
6.0,3.4,4.5,1.6
6.7,3.1,4.7,1.5
6.3,2.3,4.4,1.3
5.6,3.0,4.1,1.3
5.5,2.5,4.0,1.3
5.5,2.6,4.4,1.2
6.1,3.0,4.6,1.4
5.8,2.6,4.0,1.2
5.0,2.3,3.3,1.0
5.6,2.7,4.2,1.3
5.7,3.0,4.2,1.2
5.7,2.9,4.2,1.3
6.2,2.9,4.3,1.3
5.1,2.5,3.0,1.1
5.7,2.8,4.1,1.3];

C=[6.3,3.3,6.0,2.5
5.8,2.7,5.1,1.9
7.1,3.0,5.9,2.1
6.3,2.9,5.6,1.8
6.5,3.0,5.8,2.2
7.6,3.0,6.6,2.1
4.9,2.5,4.5,1.7
7.3,2.9,6.3,1.8
6.7,2.5,5.8,1.8
7.2,3.6,6.1,2.5
6.5,3.2,5.1,2.0
6.4,2.7,5.3,1.9
6.8,3.0,5.5,2.1
5.7,2.5,5.0,2.0
5.8,2.8,5.1,2.4
6.4,3.2,5.3,2.3
6.5,3.0,5.5,1.8
7.7,3.8,6.7,2.2
7.7,2.6,6.9,2.3
6.0,2.2,5.0,1.5
6.9,3.2,5.7,2.3
5.6,2.8,4.9,2.0
7.7,2.8,6.7,2.0
6.3,3.4,5.6,2.4
6.4,3.1,5.5,1.8
6.0,3.0,4.8,1.8
6.9,3.1,5.4,2.1
6.7,3.1,5.6,2.4
6.9,3.1,5.1,2.3
5.8,2.7,5.1,1.9
6.8,3.2,5.9,2.3
6.7,3.3,5.7,2.5
6.7,3.0,5.2,2.3
6.3,2.5,5.0,1.9
6.5,3.0,5.2,2.0
6.2,3.4,5.4,2.3
5.9,3.0,5.1,1.8];

NA=size(A,1);NB=size(B,1);NC=size(C,1);

A_train=A(1:floor(NA/2),:);%訓練數據取1/2（或者1/3,3/4,1/4）
B_train=B(1:floor(NB/2),:);
C_train=C(1:floor(NC/2),:);
u1=mean(A_train)';u2=mean(B_train)';u3=mean(C_train)';
S1=cov(A_train);S2=cov(B_train);S3=cov(C_train);
S11=inv(S1);S22=inv(S2);S33=inv(S3);
S1_d=det(S1);S2_d=det(S2);S3_d=det(S3);
PA=1/3;PB=1/3;PC=1/3; %假設各類的先驗概率相等，即都為1/3

A_test=A((floor(NA/2)+1):end,:);
B_test=B((floor(NB/2)+1):end,:);
C_test=C((floor(NC/2)+1):end,:);

%test of Sample_A
right1=0;
error1=0;
for i=1:size(A_test,1)
P1=(-1/2)*(A_test(i,:)'-u1)'*S11*(A_test(i,:)'-u1)-(1/2)*log(S1_d)+log(PA);
P2=(-1/2)*(A_test(i,:)'-u2)'*S22*(A_test(i,:)'-u2)-(1/2)*log(S2_d)+log(PB);
P3=(-1/2)*(A_test(i,:)'-u3)'*S33*(A_test(i,:)'-u3)-(1/2)*log(S3_d)+log(PC);
P=[P1 P2 P3];
[Pm,ind]=max(P);
if ind==1
right1=right1+1;
else
error1=error1+1;
end
end

right_rate=right1/size(A_test,1)   %計算出A中測試數據的准確率，同理可以測試B、C