高斯模型
1. 通过一个例子来说明
下面是一组人类身体特征的统计资料。
当特征是连续变量的时候,运用多项式模型就会导致很多(不做平滑的情况下),此时即使做平滑,所得到的条件概率也难以描述真实情况。所以处理连续的特征变量,应该采用高斯模型。
性别 | 身高(英尺) | 体重(磅) | 脚掌(英寸) |
---|---|---|---|
男 | 6 | 180 | 12 |
男 | 5.92 | 190 | 11 |
男 | 5.58 | 170 | 12 |
男 | 5.92 | 165 | 10 |
女 | 5 | 100 | 6 |
女 | 5.5 | 150 | 8 |
女 | 5.42 | 130 | 7 |
女 | 5.75 | 150 | 9 |
已知某人身高6英尺、体重130磅,脚掌8英寸,请问该人是男是女?
根据朴素贝叶斯分类器,计算下面这个式子的值。
这时,可以假设男性和女性的身高、体重、脚掌都是正态分布,通过样本计算出均值和方差,也就是得到正态分布的密度函数。有了密度函数,就可以把值代入,算出某一点的密度函数的值。
s1=[[6,180,12,1],[5.92,190,11,1],[5.58,170,12,1],[5.92,165,10,1],[5,100,6,0],[5.5,150,8,0],[5.42,130,7,0],[5.75,150,9,0]]
s1=np.array(s1)
data_nan=s1[s1[:,3]==1]
data_nv=s1[s1[:,3]==0]
data_nan_mean=data_nan.mean(axis=0)
data_nan_std=data_nan.std(axis=0)
print(data_nan_mean)
结果为:
array([ 5.855, 176.25 , 11.25 , 1. ])
pirnt(data_nan_std)
结果为:
array([0.16209565, 9.60143218, 0.8291562 , 0. ])
这是高斯函数:
def gaussiFunc(sigema,u,x):#这里的sigema是平方的形式
outxishu = 1/(np.sqrt(2*np.pi*sigema))
inxishu = -((x-u)**2)/(2*sigema)
return outxishu*np.exp(inxishu)
则由
gaussiFunc(0.16209565**2,5.855,6)
结果为:1.649603574852705
gaussiFunc(9.60143218**2,176.25,130)
结果为:3.802086207462094e-07
gaussiFunc(0.8291562**2,11.25,8)
结果为:
0.00022187199343690303
知,男性的身高为6英尺的概率的相对值等于1.6496(大于1并没有关系,因为这里是密度函数的值,只用来反映各个值的相对可能性)。
P(身高=6|男)= 1 2 π σ 2 e x p ( − ( 6 − μ ) 2 2 σ 2 ) \frac{1}{\sqrt{2\pi\sigma^2}}exp\big(\frac{-(6-\mu)^2}{2\sigma^2}\big) 2πσ21exp(2σ2−(6−μ)2)
对于脚掌和体重同样可以计算
P(体重=130|男)=3.802086207462094e-07
P(脚掌=8|男)=0.00022187199343690303
所以,P(身高=6|男) x P(体重=130|男) x P(脚掌=8|男) x P(男)=1.391e-10
同理:
P(身高=6|女) x P(体重=130|女) x P(脚掌=8|女) x P(女)=0.14420.01930.3228=0.9e-3
可以看到,女性的概率比男性要高出将近$10^7$
倍,所以判断该人为女性。
2. iris数据集为例,自定义高斯贝叶斯函数
2.1 准备数据
首先我们需要一些训练数据 这里使用鸢尾花数据。
这里x是一个(150, 4)2维数组,总共150条数据,打印其中的5条数据看一下:
[[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],
… …]
可以看到每条数据都有4个特征项分别是: 萼片的长度,萼片的宽度,花瓣的长度,花瓣的宽度
y是x里每条数据对应的分类:
[0, 0, 1, 1, 2, …]
可以看到x里对应的分类总共有3种[0,1,2]。
2.2 训练模型
其实是求出了属于每种分类里的数据在每个特征列上的平均值和方差。
计算每种分类里每个特征列的平均值和方差
{0: [[5.1, 3.5, 1.4, 0.2],
[4.9, 3.0, 1.4, 0.2],
… …],
1: [[4.7, 3.2, 1.3, 0.2],
[4.6, 3.1, 1.5, 0,2],
… …],
2: [[5.0, 3.6, 1.4, 0.2],
… …]}
得到平均值结果集
{0: [5.006, 3.418, 1.464, 0.244],
1: [5.936, 2.770, 4.260, 1.326],
2: [6.588, 2.974, 5.552, 2.026]}
得到方差结果集
{0: [0.1210, 0.1420, 0.029, 0.011],
1: [0.2611, 0.9650, 0.216, 0.038],
2: [0.3960, 0.1019, 0.298, 0.073]}
2.3 预测数据
其实是求出待预测数据属于哪种分类的概率更大。
待预测的数据
[3.1, 4.4, 2.1, 3.1]
代码如下:
def getGaussianCoeff(data_X,data_y):
values=list(set(data_y.tolist()))
box=[]
coeff=[]
data=np.column_stack((data_X,data_y))
for value in values:
zhongjie=[]
for i in range(len(data)):
if data[i][-1]==value:
zhongjie.append(data[i].tolist())
box.append(zhongjie)
box=np.array(box)
for i in range(len(values)):
sigma=[]
uVector=np.average(box[i],axis=0)
uVector=uVector.round(3)
uVectorMinusAvg=np.array([(box[i][j]-uVector) for j in range(len(box[i]))])
uVectorMinusAvgTranspose=uVectorMinusAvg.T
xieFangChaMatrix=np.dot(uVectorMinusAvgTranspose,uVectorMinusAvg)/len(box[i])
for j in range(len(xieFangChaMatrix)):
for k in range(len(xieFangChaMatrix)):
if j==k:
sigma.append(round(xieFangChaMatrix[j][k],3))
coeff.append(list(zip(sigma,uVector)))
#print(len(sigma),len(uVector))
return coeff
def gaussiFunc(sigema,u,x):#这里的sigema是平方的形式
outxishu = 1/(np.sqrt(2*np.pi*sigema))
inxishu = -((x-u)**2)/(2*sigema)
return outxishu*np.exp(inxishu)
def trainBayes(data_X,data_y,preData):
coeff=getGaussianCoeff(data_X,data_y)
if len(data_X)!=len(data_y):
raise TypeError
if not isinstance(data_X,list):
data_X=data_X.tolist()
data_y=data_y.tolist()
numOfPreData=len(data_X)
values=list(set(data_y))
preDataY=[]
preData=data_X
counts=[0]*len(values)
for k in range(len(data_X)):
for value in values:
if data_y[k]==value:
counts[values.index(value)]=counts[values.index(value)]+1
P1s=[counts[l]/sum(counts) for l in range(len(values))]
for i in range(numOfPreData):
probabilities=[]
for j in range(len(values)):
P2=1
for m in range(len(data_X[0])):
myCoeff=coeff[j][m]
probaOfOneAttr=gaussiFunc(myCoeff[0],myCoeff[1],preData[i][m])
P2=P2*probaOfOneAttr
probabilities.append(P2*P1s[j])
probabilities=np.array(probabilities)
index=np.argmax(probabilities)
preDataY.append(index)
return preDataY
if __name__ == '__main__':
import numpy as np
import time
from sklearn import datasets
iris = datasets.load_iris()
start1=time.clock()
preData = trainBayes(iris.data,iris.target,iris.data)
end1=time.clock()
print('运行的时间是:',end1-start1)
print((iris.target != preData).sum())
x = np.linspace(-10,10,100)#高斯函数检测
y = gaussiFunc(1,2,x)
import matplotlib.pyplot as plt
plt.plot(x,y)
plt.show()
- 使用sklearn包
import numpy as np
import time
from sklearn import datasets
iris=datasets.load_iris()
from sklearn.naive_bayes import GaussianNB
gnb=GaussianNB()
start1=time.clock()
y_pred=gnb.fit(iris.data,iris.target).predict(iris.data)
end1=time.clock()
print('运行的时间是:',end1-start1)
print(iris.data.shape[0],(iris.target!=y_pred).sum())
- 手写高斯贝叶斯
#如何导入数据集?
f=open('e:/dataset1.txt','r')
d=[]
for i in f.readlines():
d.append(i.strip('\n').split('\t'))
d
f.close()
import numpy as np
d=np.array(d)
data=d[:,0:2]
t=d[:,2]
data=[list(map(float,i)) for i in data]
data=np.array(data)
#需要自己求出data_nan_mean,data_nan_std,data_nv_mean,data_nv_std,data_nan_lv,data_nv_lv
data_nan=data[t=='m']
data_nv=data[t=='f']
data_nan_mean=np.mean(data_nan,axis=0)
data_nan_std=np.std(data_nan,axis=0)
data_nv_mean=np.mean(data_nv,axis=0)
data_nv_std=np.std(data_nv,axis=0)
data_nan_lv=len(data_nan)/len(data)
data_nv_lv=len(data_nv)/len(data)
#高斯函数
def gaussiFunc(sigema,u,x):#这里的sigema是平方的形式
outxishu = 1/(np.sqrt(2*np.pi*sigema))
inxishu = -((x-u)**2)/(2*sigema)
return outxishu*np.exp(inxishu)
#主程序代码
idx=[]
data_mean=[data_nv_mean,data_nan_mean]
data_std=[data_nv_std,data_nan_std]
lv=[data_nv_lv,data_nan_lv]
for k in data:
re=[1,1]
for i in range(2):
g=gaussiFunc(data_std[i]**2,data_mean[i],k)
for j in g:
re[i]=re[i]*j
re[i]=re[i]*lv[i]
idx.append(np.argmax(re))
idx=np.array(idx)
idx
#计算误差
t[t=='f']=0
t[t=='m']=1
t=np.array([int(i) for i in t])
c=0
for i in t==idx:
if i==False:
c+=1
print('分类错误数为:',c)
print('错误率为:',c/len(d))
结果为:
分类错误数为: 39
错误率为: 0.08628318584070796