定义
假设满足所有约束条件的模型集合为:
C
=
{
P
∈
P
∣
E
p
(
f
i
)
=
E
p
~
(
f
i
)
,
i
=
1
,
2
,
⋯
,
n
}
{\tt C}=\left\{ P \in {\tt P} | E_p(f_i)=E_{\tilde{p}(f_i)},i=1,2,\cdots,n\right\}
C={P∈P∣Ep(fi)=Ep~(fi),i=1,2,⋯,n}
E
p
(
f
)
=
∑
x
,
y
P
~
(
x
)
P
(
y
∣
x
)
f
(
x
,
y
)
E_p(f)=\sum_{x,y}\tilde{P}(x)P(y|x)f(x,y)
Ep(f)=∑x,yP~(x)P(y∣x)f(x,y)
E
p
~
(
f
)
=
∑
x
,
y
P
~
(
x
,
y
)
f
(
x
,
y
)
E_{\tilde{p}}(f)=\sum_{x,y}\tilde{P}(x,y)f(x,y)
Ep~(f)=∑x,yP~(x,y)f(x,y)
P
~
(
X
=
x
)
=
v
(
X
=
x
,
Y
=
y
)
N
\tilde{P}(X=x)=\dfrac{{\tt v}(X=x,Y=y)}{N}
P~(X=x)=Nv(X=x,Y=y),
P
~
(
X
=
x
,
Y
=
y
)
=
v
(
X
=
x
)
N
\tilde{P}(X=x,Y=y)=\dfrac{{\tt v}(X=x)}{N}
P~(X=x,Y=y)=Nv(X=x),其中
v
(
X
=
x
,
Y
=
y
)
{\tt v}(X=x,Y=y)
v(X=x,Y=y):训练数据中样本
(
x
,
y
)
(x,y)
(x,y)出现的频次
,
v
(
X
=
x
)
{\tt v}(X=x)
v(X=x):训练数据中输入
x
x
x出现的频次,
N
N
N:训练样本容量。
f
(
x
,
y
)
:
{
1
,
x
与
y
满足某一事实
0
,
否则
f(x,y):\begin{cases} 1,x与y满足某一事实\\ 0,否则 \end{cases}
f(x,y):{1,x与y满足某一事实0,否则
定义在条件概率分布
P
(
Y
∣
X
)
P(Y|X)
P(Y∣X)上的条件熵为:
H
(
P
)
=
−
∑
x
,
y
P
~
(
x
)
P
(
y
∣
x
)
l
o
g
P
(
y
∣
x
)
H(P)=-\sum_{x,y}\tilde{P}(x)P(y|x)log P(y|x)
H(P)=−∑x,yP~(x)P(y∣x)logP(y∣x),则模型集合
C
{\tt C}
C中条件熵
H
(
P
)
H(P)
H(P)最大的模型称为最大熵模型。式中的对数为自然对数。
输入空间
T= { ( x 1 , y 1 ) , ( x 2 , y 2 ) , … , ( x N , y N ) } \left\{(x_1,y_1),(x_2,y_2),\dots,(x_N,y_N)\right\} {(x1,y1),(x2,y2),…,(xN,yN)}
import time
import numpy as np
from collections import defaultdict
def loadData(fileName):
'''
加载Mnist数据集 下载地址:https://download.csdn.net/download/nanxiaotao/89720991)
:param fileName:要加载的数据集路径
:return: 数据集
'''
# 存放数据及标记的list
dataList = []; labelList = []
# 打开文件
fr = open(fileName, 'r')
# 将文件按行读取
for line in fr.readlines():
curLine = line.strip().split(',')
if int(curLine[0]) == 0:
labelList.append(1)
else:
labelList.append(0)
#二值化
dataList.append([int(int(num) > 128) for num in curLine[1:]])
#返回data和label
return dataList, labelList
trainData, trainLabel = loadData('../Mnist/mnist_train.csv')
np.shape(trainData)
特征空间(Feature Space)
trainData[0][0:784]
统计学习方法
模型
P ω ( y ∣ x ) = 1 Z ω ( x ) e x p ( ∑ i = 1 n ω i f i ( x , y ) ) P_\omega(y|x)=\dfrac{1}{Z_\omega(x)}exp(\sum_{i=1}^n\omega_if_i(x,y)) Pω(y∣x)=Zω(x)1exp(∑i=1nωifi(x,y)),其中 Z ω ( x ) = ∑ i = 1 n e x p ( ω i f i ( x , y ) ) Z_\omega(x)= \sum_{i=1}^nexp(\omega_if_i(x,y)) Zω(x)=∑i=1nexp(ωifi(x,y))
策略
L
(
P
,
ω
)
=
∑
x
,
y
P
~
(
x
)
P
(
y
∣
x
)
l
o
g
P
(
y
∣
x
)
+
ω
0
(
1
−
∑
y
P
(
y
∣
x
)
)
+
L(P,\omega)=\sum_{x,y}\tilde{P}(x)P(y|x)logP(y|x)+\omega_0\bigg(1-\sum_{y}P(y|x)\bigg)+
L(P,ω)=∑x,yP~(x)P(y∣x)logP(y∣x)+ω0(1−∑yP(y∣x))+
∑
i
=
1
n
ω
i
(
∑
x
,
y
P
~
(
x
,
y
)
f
i
(
x
,
y
)
−
∑
x
,
y
P
~
(
x
)
P
(
y
∣
x
)
f
i
(
x
,
y
)
)
\sum_{i=1}^n \omega_i\bigg(\sum_{x,y}\tilde{P}(x,y)f_i(x,y)-\sum_{x,y}\tilde{P}(x)P(y|x)f_i(x,y)\bigg)
∑i=1nωi(∑x,yP~(x,y)fi(x,y)−∑x,yP~(x)P(y∣x)fi(x,y))
m
i
n
P
∈
C
m
a
x
ω
L
(
P
,
ω
)
\mathop{min}\limits_{ P \in {\tt C} }\mathop{max}\limits_{ \omega }L(P,\omega)
P∈CminωmaxL(P,ω)
m
a
x
ω
m
i
n
P
∈
C
L
(
P
,
ω
)
\mathop{max}\limits_{ \omega }\mathop{min}\limits_{ P \in {\tt C} }L(P,\omega)
ωmaxP∈CminL(P,ω)
算法
trainDataList = trainData[:20000]
trainLabelList = trainLabel[:20000]
M = 10000
N = len(trainDataList) #训练样本容量
featureNum = len(trainDataList[0]) #特征数量
n = 0 #去重后(x,y)的数量,初始值为0
v ( X = x , Y = y ) {\tt v}(X=x,Y=y) v(X=x,Y=y)
def calc_fixy():
'''
计算(x, y)在训练集中出现过的次数
:return:
'''
global n
#建立特征数目个字典,属于不同特征的(x, y)对存入不同的字典中,保证不被混淆
fixyDict = [defaultdict(int) for i in range(featureNum)]
#遍历训练集中所有样本
for i in range(len(trainDataList)):
#遍历样本中所有特征
for j in range(featureNum):
#将出现过的(x, y)对放入字典中并计数值加1
fixyDict[j][(trainDataList[i][j], trainLabelList[i])] += 1
#对整个大字典进行计数,判断去重后还有多少(x, y)对,写入n
for i in fixyDict:
n = n + len(i)
#返回大字典
return fixyDict
fixy = calc_fixy() #每对(x,y)出现的频数
P ω ( y ∣ x ) = 1 Z ω ( x ) e x p ( ∑ i = 1 n ω i f i ( x , y ) ) P_\omega(y|x)=\dfrac{1}{Z_\omega(x)}exp(\sum_{i=1}^n\omega_if_i(x,y)) Pω(y∣x)=Zω(x)1exp(∑i=1nωifi(x,y))
w = [0] * n # Pw(y|x)中的w
def createSearchDict():
'''
创建查询字典
xy2idDict:通过(x,y)对找到其id,所有出现过的xy对都有一个id
id2xyDict:通过id找到对应的(x,y)对
'''
xy2idDict = [{} for i in range(featureNum)]
id2xyDict = {}
index = 0
#对特征进行遍历
for feature in range(featureNum):
for (x, y) in fixy[feature]:
xy2idDict[feature][(x, y)] = index
id2xyDict[index] = (x, y)
#id加一
index += 1
#返回创建的两个字典
return xy2idDict, id2xyDict
xy2idDict, id2xyDict = createSearchDict() #(x, y)->id和id->(x, y)的搜索字典
E p ~ ( f ) = ∑ x , y P ~ ( x , y ) f ( x , y ) E_{\tilde{p}}(f)=\sum_{x,y}\tilde{P}(x,y)f(x,y) Ep~(f)=∑x,yP~(x,y)f(x,y)
def calcEp_xy():
'''
计算特征函数f(x, y)关于经验分布P_(x, y)的期望值(下划线表示P上方的横线,
同理Ep_xy中的“_”也表示p上方的横线)
:return: 计算得到的Ep_xy
'''
#初始化Ep_xy列表,长度为n
Ep_xy = [0] * n
#遍历每一个特征
for feature in range(featureNum):
#遍历每个特征中的(x, y)对
for (x, y) in fixy[feature]:
#获得其id
id = xy2idDict[feature][(x, y)]
#将计算得到的Ep_xy写入对应的位置中
#fixy中存放所有对在训练集中出现过的次数,处于训练集总长度N就是概率了
Ep_xy[id] = fixy[feature][(x, y)] / N
#返回期望
return Ep_xy
Ep_xy = calcEp_xy() #Ep_xy期望值
E p ( f ) = ∑ x , y P ~ ( x ) P ( y ∣ x ) f ( x , y ) E_p(f)=\sum_{x,y}\tilde{P}(x)P(y|x)f(x,y) Ep(f)=∑x,yP~(x)P(y∣x)f(x,y)
def calcEpxy():
'''
计算特征函数f(x, y)关于模型P(Y|X)与经验分布P_(X, Y)的期望值(P后带下划线“_”表示P上方的横线
程序中部分下划线表示“|”,部分表示上方横线)
:return:
'''
Epxy = [0] * n
#对于每一个样本进行遍历
for i in range(N):
#初始化公式中的P(y|x)列表
Pwxy = [0] * 2
#计算P(y = 0 } X)
#注:程序中X表示是一个样本的全部特征,x表示单个特征,这里是全部特征的一个样本
Pwxy[0] = calcPwy_x(trainDataList[i], 0)
#计算P(y = 1 } X)
Pwxy[1] = calcPwy_x(trainDataList[i], 1)
for feature in range(featureNum):
for y in range(2):
if (trainDataList[i][feature], y) in fixy[feature]:
id = xy2idDict[feature][(trainDataList[i][feature], y)]
Epxy[id] += (1 / N) * Pwxy[y]
return Epxy
P ω ( y ∣ x ) = 1 Z ω ( x ) e x p ( ∑ i = 1 n ω i f i ( x , y ) ) P_\omega(y|x)=\dfrac{1}{Z_\omega(x)}exp(\sum_{i=1}^n\omega_if_i(x,y)) Pω(y∣x)=Zω(x)1exp(∑i=1nωifi(x,y)),其中 Z ω ( x ) = ∑ i = 1 n e x p ( ω i f i ( x , y ) ) Z_\omega(x)= \sum_{i=1}^nexp(\omega_if_i(x,y)) Zω(x)=∑i=1nexp(ωifi(x,y))
def calcPwy_x(X, y):
'''
最大熵模型
:param X: 要计算的样本X(一个包含全部特征的样本)
:param y: 该样本的标签
:return: 计算得到的Pw(Y|X)
'''
#分子
numerator = 0
#分母
Z = 0
#对每个特征进行遍历
for i in range(featureNum):
#如果该(xi,y)对在训练集中出现过
if (X[i], y) in xy2idDict[i]:
#在xy->id字典中指定当前特征i,以及(x, y)对:(X[i], y),读取其id
index = xy2idDict[i][(X[i], y)]
numerator += w[index]
if (X[i], 1-y) in xy2idDict[i]:
#原理与上式相同
index = xy2idDict[i][(X[i], 1-y)]
Z += w[index]
#计算分子的指数
numerator = np.exp(numerator)
#计算分母的z
Z = np.exp(Z) + numerator
#返回Pw(y|x)
return numerator / Z
模型训练
def maxEntropyTrain(iter = 500):
global w
#设置迭代次数寻找最优解
for i in range(iter):
#单次迭代起始时间点
iterStart = time.time()
Epxy = calcEpxy()
#使用的是IIS,所以设置sigma列表
sigmaList = [0] * n
#对于所有的n进行一次遍历
for j in range(n):
#依据“6.3.1 改进的迭代尺度法” 式6.34计算
sigmaList[j] = (1 / M) * np.log(Ep_xy[j] / Epxy[j])
#按照算法6.1步骤二中的(b)更新w
w = [w[i] + sigmaList[i] for i in range(n)]
#单次迭代结束
iterEnd = time.time()
#打印运行时长信息
print('iter:%d:%d'%(i, iter, iterEnd-iterStart ))
maxEntropyTrain()
假设空间(Hypothesis Space)
{ f ∣ f ( x ) = m a x ( m a x ω m i n P ∈ C L ( P , ω ) ) } \left\{f|f(x) = max \big( \mathop{max}\limits_{ \omega }\mathop{min}\limits_{ P \in {\tt C} }L(P,\omega) \big) \right\} {f∣f(x)=max(ωmaxP∈CminL(P,ω))}
输出空间
Y ∈ { 0 , 1 } Y \in \{ 0,1 \} Y∈{0,1}
模型评估
训练误差
testDataList, testLabelList = loadData('../Mnist/mnist_test.csv')
def predict(X):
'''
预测标签
:param X:要预测的样本
:return: 预测值
'''
#因为y只有0和1,所有建立两个长度的概率列表
result = [0] * 2
#循环计算两个概率
for i in range(2):
#计算样本x的标签为i的概率
result[i] = calcPwy_x(X, i)
#返回标签
#max(result):找到result中最大的那个概率值
#result.index(max(result)):通过最大的那个概率值再找到其索引,索引是0就返回0,1就返回1
return result.index(max(result))
def model_test():
'''
对测试集进行测试
:return:
'''
#错误值计数
errorCnt = 0
#对测试集中所有样本进行遍历
for i in range(len(testDataList)):
#预测该样本对应的标签
result = predict(testDataList[i])
#如果错误,计数值加1
if result != testLabelList[i]: errorCnt += 1
#返回准确率
return 1 - errorCnt / len(testDataList)
accuracy = model_test()
print('accuracy:', accuracy)