模式识别第四版(希腊)西奥多里蒂斯第二章上机实验(2.1-2.8)Python实现
本实验使用Pthon 3.x版本
先将下面Python函数保存到一个文件中, 例如文件名: pattern_classification,下面的Python函数包括一个高斯分布数据生成函数,一个贝叶斯分类器函数,一个欧氏距离分类器函数,一个马氏距离分类器函数,一个最近邻分类器函数,一个计算误判概率函数。
import numpy as np
from scipy.stats import multivariate_normal
def generate_gauss_classes(m,S,P,N):
"""
函数用于从高斯类别生成训练数据集.
输入参数:
m: 二维(l×c)数组, 第i列为第i个类别分布的均值向量, 这里l是字母,
不是数字; c是类别数.
S: l×l×c(三维)数组, 其第i个二维l×l部分是第i类分布的协方差.
P: c-维数组, 包含类别的先验概率.
N: 样本数.
输出:
X: l×N数组, 每列是l-维数据向量.
y: N维-数组, 其第i个元素表示第i个数据向量的类别
"""
c = len(m[:,0])
X = np.empty((1,2))
y_1 = np.empty(1)
for j in np.linspace(0,c-1,c,dtype=np.int):
a = np.fix(P[j]*N)
b = a.astype(np.int)
t = np.random.multivariate_normal(m[j,:],S[j,:,:],b)
X = np.append(X,t,0)
w = np.ones(b)* j
y_1 = np.append(y_1,w,0)
X = np.delete(X,0,0)
y_1 = np.delete(y_1,0)
y = y_1.astype(np.int)
return X, y
def bayes_classifier(m,S,P,X):
"""
函数用于实现(高斯过程的)贝叶斯分类器.
输入参数:
m: 二维(l×c)数组, 第i列为第i个类别分布的均值向量, 这里l是字母,
不是数字; c是类别数.
S: 三维(l×l×c)数组, 其第i个二维l×l部分是第i类分布的协方差.
P: c-维数组, 包含类别的先验概率.
X: 二维(N×l)数组, 每列是l-维数据向量, 这里N为待分类的样本数.
输出参数:
y: N维数组, 其第i个元素表示第i个数据向量的类别.
"""
c = len(m[:,0])
N = len(X[:,0])
t = np.empty((N,c))
y = np.empty(N)
for j in np.linspace(0,c-1,c,dtype=np.int):
t[:,j] = P[j]*multivariate_normal.pdf(X, m[j,:],S[j,:,:])
# Determining the maximum quantity Pi*p(x|wi)
y = np.argmax(t,axis=1)
return y
def euclidean_classifier(m,X):
"""
函数用于实现欧氏距离分类器.
输入参数:
m: 二维(c×l)数组, 第i列为第i个类别分布的均值向量; c是类别数.
X: 二维(N×l)数组, 每列是l-维数据向量, 这里N为待分类的样本数.
输出参数:
y: N×1数组, 其第i个元素表示第i个数据向量的类别.
"""
c = len(m[:,0])
t = np.empty(c)
N = len(X[:,0])
y_1 = np.empty(N)
for i in np.linspace(0,N-1,N,dtype=np.int):
for j in np.linspace(0,c-1,c,dtype=np.int):
t[j] = np.sqrt(np.dot(X[i,:]-m[j,:],np.transpose(X[i,:]-m[j,:])))
y_1[i] = np.argmin(t)
y = y_1.astype(np.int)
return y
def mahalanobis_classifier(m,S,X):
"""
函数用于实现欧氏距离分类器.
输入参数:
m: 二维(l×c)数组, 第i列为第i个类别分布的均值向量; c是类别数.
S: 三维(l×l×c)数组, 其第i个二维l×l部分是第i类分布的协方差.
X: 二维(N×l)数组, 每列是l-维数据向量, 这里N为待分类的样本数.
输出参数:
y: N维-向量, 其第i个元素表示第i个数据向量的类别.
"""
c = len(m[:,0])
t = np.empty(c)
N = len(X[:,0])
y_1 = np.empty(N)
for i in np.linspace(0,N-1,N,dtype=np.int):
for j in np.linspace(0,c-1,c,dtype=np.int):
t[j] = np.sqrt(np.dot(X[i,:]-m[j,:],np.dot(np.linalg.inv(S[1,:,:]),
np.transpose(X[i,:]-m[j,:]))))
y_1[i] = np.argmin(t)
y = y_1.astype(np.int)
return y
def k_nn_classifier(Z,v,k,X):
"""
函数用于k-最近邻分类器.
输入参数:
Z: 二维(N_1×l)数组, 它的行由N_1个训练向量构成.
v: 一维(N_1×1)数组,, 其分量为Z中对应向量所属类别标记,
c个类别标记为0,2, ..., c - 1.
k: kNN分类器参数.
X: 二维(N×l)数组, 每行是l-维数据向量, 这里N为待分类的样本数.
输出参数:
y: 一维(N×1)数组, 其第i个元素表示第i个数据向量的类别.
"""
N_1 = len(Z[:,0])
N = len(X[:,0])
c = max(v)
dist = np.empty(N_1)
y_1 = np.empty(N)
for i in np.linspace(0,N-1,N,dtype=np.int):
dist = np.sum((X[i,:]- Z) ** 2, axis=1)
nearest = np.argsort(dist)
refe = (-1)*np.ones(c+1)
for q in np.linspace(0,k-1,k,dtype=np.int):
classes = v[nearest[q]]
refe[classes] = refe[classes] + 1
y_1[i] = np.argmax(refe)
y = y_1.astype(np.int)
return y
def compute_error(y,y_est):
"""
函数用于计算分类器误判概率.
输入参数:
y: N-维数组, 其第i个分量为一个数据集X的列相应的真实类别.
c个类别标记为0,2, ..., c-1.
y_est: N-维数组, 其第i个分量为根据分类器, 为数据集X的列指定相应的类别.
输出参数:
clas_error: 误判概率.
"""
N = len(y)
clas_error = np.sum((y!=y_est))/N
return clas_error
一、实验2.1
Python代码:(下面画图的代码格式是为了使得轴刻度向内)
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.axisartist.axislines import SubplotZero
from pattern_classification import *
# 下面两行是为了能够在图形中插入美观的数学符号
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
np.random.seed(123456789)
N = 1000
m = np.array([[1, 1],[7, 7],[15, 1]])
S = np.array([[[12,0],[0,1]],[[8,3],[3,2]],[[2,0],[0,2]]])
P_1 = np.array([1/3,1/3,1/3])
P_2 = np.array([0.6,0.3,0.1])
X_1, y_1 = generate_gauss_classes(m,S,P_1,N)
X_2, y_2 = generate_gauss_classes(m,S,P_2,N)
markers = ['.','*', '1', 's', '+', '1', '2', '3', '4']
pale=["red", "green", "blue", "yellow", "m", "c"]
N = len(y_1)
c = len(m[:,0])
fig = plt.figure(1, figsize=(10,4.5))
ax = SubplotZero(fig, 1, 2, 1)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_1[i,0],X_1[i,1],color=pale[y_1[i]],ls ='', marker=markers[y_1[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("equiprobable classes",fontsize=18)
ax = SubplotZero(fig, 1, 2, 2)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_2[i,0],X_2[i,1],color=pale[y_2[i]],ls ='', marker=markers[y_2[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("nonequiprobable classes",fontsize=18)
plt.savefig("Figure2_1.pdf", dpi=100)
plt.show
二、实验2.2
Python代码为:
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.axisartist.axislines import SubplotZero
from pattern_classification import *
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
np.random.seed(123456789)
N = 1000
m = np.array([[1, 1],[12, 8],[16, 1]])
S = np.array([[[4,0],[0,4]],[[4,0],[0,4]],[[4,0],[0,4]]])
P = np.array([1/3,1/3,1/3])
X, y_1 = generate_gauss_classes(m,S,P,N)
y_2 = bayes_classifier(m,S,P,X)
y_3 = euclidean_classifier(m,X)
y_4 = mahalanobis_classifier(m,S,X)
e_B = compute_error(y_1,y_2)
e_E = compute_error(y_1,y_3)
e_M = compute_error(y_1,y_4)
print('---------------------------------------------')
print("Bayes classifier error: %4.4f" % e_B)
print("euclidean classifier error: %4.4f" % e_E)
print("mahalanobis classifier error: %4.4f" % e_M)
print('---------------------------------------------')
markers = ['.','*', '1', 's', '+', '1', '2', '3', '4']
pale=["red", "green", "blue", "yellow", "m", "c"]
N = len(y_1)
c = len(m[:,0])
fig = plt.figure(1, figsize=(10,9))
ax = SubplotZero(fig, 2, 2, 1)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X[i,0],X[i,1],color=pale[y_1[i]],ls ='', marker=markers[y_1[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Data Set",fontsize=18)
ax = SubplotZero(fig, 2, 2, 2)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X[i,0],X[i,1],color=pale[y_2[i]],ls ='', marker=markers[y_2[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Bayes classes",fontsize=18)
ax = SubplotZero(fig, 2, 2, 3)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X[i,0],X[i,1],color=pale[y_3[i]],ls ='', marker=markers[y_3[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Euclid classes",fontsize=18)
ax = SubplotZero(fig, 2, 2, 4)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X[i,0],X[i,1],color=pale[y_4[i]],ls ='', marker=markers[y_4[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Mahalanobis classes",fontsize=18)
plt.savefig("Figure2_2.pdf", dpi=100)
plt.show
计算的误判概率为:
分类器类型 | 误判概率 |
---|---|
贝叶斯分类器 | 0.008 |
欧氏距离分类器 | 0.008 |
马氏距离分类器 | 0.008 |
分类结果图形:
三、实验2.3
实现的Python代码为将实验2.2中的代码m和S改为
m = np.array([[1, 1],[14, 7],[16, 1]])
S = np.array([[[5,3],[3,4]],[[5,3],[3,4]],[[5,3],[3,4]]])
即可。
计算的误判概率为:
分类器类型 | 误判概率 |
---|---|
贝叶斯分类器 | 0.003 |
欧氏距离分类器 | 0.007 |
马氏距离分类器 | 0.003 |
分类结果图形为
四、实验2.4
实现的Python代码为将实验2.2中的代码m和S改为
m = np.array([[1, 1],[8, 6],[13, 1]])
S = np.array([[[6,0],[0,6]],[[6,0],[0,6]],[[6,0],[0,6]]])
即可
计算的误判概率为:
分类器类型 | 误判概率 |
---|---|
贝叶斯分类器 | 0.0611 |
欧氏距离分类器 | 0.0611 |
马氏距离分类器 | 0.0611 |
分类结果的图形为:
五、实验2.5
实现的Python代码为将实验2.2中的代码m和S改为
m = np.array([[1, 1],[10, 5],[11, 1]])
S = np.array([[[7,4],[4,5]],[[7,4],[4,5]],[[7,4],[4,5]]])
即可。
计算的误判概率为:
分类器类型 | 误判概率 |
---|---|
贝叶斯分类器 | 0.0711 |
欧氏距离分类器 | 0.1251 |
马氏距离分类器 | 0.0711 |
分类结果的图形为:
六、实验2.6
实验2.6为实验2.2–2.5的实验结果总结。
七、实验2.7
实现的Python代码为:
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.axisartist.axislines import SubplotZero
from pattern_classification import *
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
np.random.seed(123456789)
N = 1000
m = np.array([[1, 1],[4, 4],[8, 1]])
S = np.array([[[2,0],[0,2]],[[2,0],[0,2]],[[2,0],[0,2]]])
P_1 = np.array([1/3,1/3,1/3])
P_2 = np.array([0.8,0.1,0.1])
X_1, y_1 = generate_gauss_classes(m,S,P_1,N)
X_2, y_2 = generate_gauss_classes(m,S,P_2,N)
w_1 = bayes_classifier(m,S,P_1,X_1)
w_2 = bayes_classifier(m,S,P_2,X_2)
z_1 = euclidean_classifier(m,X_1)
z_2 = euclidean_classifier(m,X_2)
e_B1 = compute_error(y_1,w_1)
e_B2 = compute_error(y_2,w_2)
e_E1 = compute_error(y_1,z_1)
e_E2 = compute_error(y_2,z_2)
print('----------------------------------------------------------')
print("Bayes classifier error(equiprobable): %4.4f" % e_B1)
print("Bayes classifier error(nonequiprobable): %4.4f" % e_B2)
print("euclidean classifier error(equiprobable): %4.4f" % e_E1)
print("euclidean classifier error(equiprobable): %4.4f" % e_E2)
print('----------------------------------------------------------')
markers = ['.','*', '1', 's', '+', '1', '2', '3', '4']
pale=["red", "green", "blue", "yellow", "m", "c"]
N = len(y_1)
c = len(m[:,0])
fig = plt.figure(1, figsize=(10,13.5))
ax = SubplotZero(fig, 3, 2, 1)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_1[i,0],X_1[i,1],color=pale[y_1[i]],ls ='', marker=markers[y_1[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Data Set(equiprobable)",fontsize=18)
ax = SubplotZero(fig, 3, 2, 2)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_2[i,0],X_2[i,1],color=pale[y_2[i]],ls ='', marker=markers[y_2[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Data Set(nonequiprobable)",fontsize=18)
ax = SubplotZero(fig, 3, 2, 3)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_1[i,0],X_1[i,1],color=pale[w_1[i]],ls ='', marker=markers[w_1[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Bayes classes(equiprobable)",fontsize=18)
ax = SubplotZero(fig, 3, 2, 4)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_2[i,0],X_2[i,1],color=pale[w_2[i]],ls ='', marker=markers[w_2[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Bayes classes(nonequiprobable)",fontsize=18)
ax = SubplotZero(fig, 3, 2, 5)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_1[i,0],X_1[i,1],color=pale[z_1[i]],ls ='', marker=markers[z_1[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Euclid classes(equiprobable)",fontsize=18)
ax = SubplotZero(fig, 3, 2, 6)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X_2[i,0],X_2[i,1],color=pale[z_2[i]],ls ='', marker=markers[z_2[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Euclid classes(equiprobable)",fontsize=18)
plt.savefig("Figure2_7.pdf", dpi=100)
plt.show
计算的误判概率为:
分类器类型 | 误判概率 |
---|
贝叶斯 | 等可能情形 | 0.0661 |
---|---|---|
非等可能情形 | 0.0400 | |
欧氏距离 | 等可能情形 | 0.0661 |
非等可能情形 | 0.0740 |
分类结果的图形为:
八、实验2.8
实现的Python代码为(生成的训练样本数为400,不同的训练样本数可能对分类结果略有影响,但训练样本数太小可能有实质性的影响):
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.axisartist.axislines import SubplotZero
from pattern_classification import *
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
np.random.seed(123456789)
N = 1000
N_1 = 400
m = np.array([[1, 1],[8, 6],[13, 1]])
S = np.array([[[6,0],[0,6]],[[6,0],[0,6]],[[6,0],[0,6]]])
P = np.array([1/3,1/3,1/3])
X, y = generate_gauss_classes(m,S,P,N)
Z, v = generate_gauss_classes(m,S,P,N_1)
y_1 = k_nn_classifier(Z,v,1,X)
y_2 = k_nn_classifier(Z,v,11,X)
e_NN = compute_error(y,y_1)
e_11NN = compute_error(y,y_2)
print('---------------------------------------------')
print("NN classifier error: %4.4f" % e_NN)
print("11NN classifier error: %4.4f" % e_11NN)
print('---------------------------------------------')
markers = ['.','*', '1', 's', '+', '1', '2', '3', '4']
pale=["red", "green", "blue", "yellow", "m", "c"]
N = len(y_1)
c = len(m[:,0])
fig = plt.figure(1, figsize=(13,3.5))
ax = SubplotZero(fig, 1, 3, 1)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X[i,0],X[i,1],color=pale[y[i]],ls ='', marker=markers[y[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("Data Set",fontsize=18)
ax = SubplotZero(fig, 1, 3, 2)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X[i,0],X[i,1],color=pale[y_1[i]],ls ='', marker=markers[y_1[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("NN classes",fontsize=18)
ax = SubplotZero(fig, 1, 3, 3)
fig.add_subplot(ax)
for i in np.linspace(0,N-1,N,dtype=np.int):
ax.plot(X[i,0],X[i,1],color=pale[y_2[i]],ls ='', marker=markers[y_2[i]])
for j in np.linspace(0,c-1,c,dtype=np.int):
ax.plot(m[j,0],m[j,1],color="k", marker='+')
ax.set_xlabel(r"$x_1$")
ax.set_ylabel(r"$x_2$")
ax.set_title("11NN classes",fontsize=18)
plt.savefig("Figure2_8.pdf", dpi=100)
plt.show
计算的误判概率为:
分类器类型 | 误判概率 |
---|---|
NN分类器 | 0.1091 |
11NN分类器 | 0.0751 |