【数学建模】Bayes判别分析——Python

目录

一、何为Bayes判别分析

二、案例引入

三、代码实现

四、总结

---

一、何为Bayes判别分析

        不想说的很复杂，给一堆公式对于初学者来说真的很不友好的。

概率论学过吧？就是在已知存在哪些模式（类别）的情况下判断一个样本属于哪个模式的可能性最大的分析方式。核心是通过概率计算来判断大小的。

二、案例引入

        已知类别：

        待识别对象：试判断下列国家应该被分类为哪一类

三、代码实现

        ① 第一步是导入相关库：

import xlrd as xr
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as LDA
import numpy as np
from sklearn.neighbors import KNeighborsClassifier

        ②第二步是读取excel中的数据

##导入数据
#类别数据
file_location="/Users/lifangjian/Desktop/本科课程资料/数学建模/判别分析案例.xls"
data=xr.open_workbook(file_location)
sheet=data.sheet_by_index(0)
lbstats = [[sheet.cell_value(r,c) for c in range(2,sheet.ncols)] for r in range(1,sheet.nrows)]
lbstats_Frame = pd.DataFrame(lbstats)
print(lbstats_Frame)
#样本数据
sheet1=data.sheet_by_index(1)
ybstats = [[sheet1.cell_value(r,c) for c in range(1,sheet1.ncols)] for r in range(1,sheet1.nrows)]
ybstats_Frame = pd.DataFrame(ybstats)
print(ybstats_Frame)

#数据堆砌
X0=np.array(lbstats_Frame)
x=np.array(ybstats_Frame)
g=np.hstack([np.ones(5),2*np.ones(5)])

        数据结果为:

0     1       2
0  76.0  99.0  5374.0
1  79.5  99.0  5359.0
2  78.0  99.0  5372.0
3  72.1  95.9  5242.0
4  73.8  77.7  5370.0
5  71.2  93.0  4250.0
6  75.3  94.9  3412.0
7  70.0  91.2  3390.0
8  72.8  99.0  2300.0
9  62.9  80.6  3799.0

      0     1       2
0  68.5  79.3  1950.0
1  69.9  96.9  2840.0
2  77.6  93.8  5233.0
3  69.3  90.3  5158.0

[5.7910721  0.2649815  0.03407283]

#其中g=[1. 1. 1. 1. 1. 2. 2. 2. 2. 2.]

         ③ 第三步是计算均值、类内协方差矩阵以及逆矩阵、总体协方差矩阵等

##计算类间协方差及其转置矩阵
l1=np.vstack(X0[0:5])
l2=np.vstack(X0[5:])
vl1=np.cov(l1.T)
vl2=np.cov(l2.T)

##计算总体协方差矩阵以及转置矩阵
v1=np.cov(X0)
v=np.cov(X0.T)

#计算类内均值
l1_mean=l1.mean(axis=0)
l2_mean=l2.mean(axis=0)

#判断行数(一类中的个数)
n1 = l1.shape[0]
n2 = l2.shape[0]
#计算类内协方差,rowvar=False代表计算方式为列，每一列称之为一个变量
s1 = np.cov(l1,rowvar=False)
s2 = np.cov(l2,rowvar=False)
#计算联合协方差
s = ((n1-1)*s1 + (n2-1)*s2)/(n1+n2-2)

        ④第四步是计算bayes回归系数和截距，怎么算不用摆出来莫名其妙看不懂的公式，上代码就好：

##计算系数和常数项
#矩阵求逆np.linalg.inv(s)
#np.dot代表矩阵相乘
B1=np.array(l1_mean).reshape(1,3).dot(np.linalg.inv(s))[0]
print(B1)
#就是yhat=0，xhat=li.mean，d=(yhat-xhat*B)/2
d1=0
for i in range(len(B1)):
    d1-=l1_mean[i]*B1[i]
d1=d1/2
print("y1的Baye判别函数为: y1={}+{}x1+{}x2+{}x3".format(d1,B1[0],B1[1],B1[2]))

B2=np.array(l2_mean).reshape(1,3).dot(np.linalg.inv(s))[0]
print(B2)
d2=0
for i in range(len(B2)):
    d2-=l2_mean[i]*B2[i]
d2=d2/2
print("y2的Baye判别函数为: y2={}+{}x1+{}x2+{}x3".format(d2,B2[0],B2[1],B2[2]))

        输出结果（自行取浮点数就好）：

y1的Baye判别函数为: 

y1=-323.2156752155551+5.791072103680565x1+0.2649815029883881x2+0.034072826317056895x3

y2的Baye判别函数为: y2=-236.0382324407707+5.140342643898282x1+0.25166601634249214x2+0.025334641917697583x3

        ⑤第五步就是做判别分析，这里就暂时先不放代码了。

四、总结

        说实话，R语言真的会很方便，什么时候有没有大佬开个源让python也能便捷实现R、Stata的数分功能。