《统计学习方法》阅读笔记

前言

阅读李航《统计学习方法》一书笔记整理

统计学习方法-第一章:统计学习方法概论

统计学习方法是基于数据结构统计模型从而对数据进行预测分析。

统计学习方法由 监督学习非监督学习半监督学习强化学习 等组成。

模型、策略和算法 称为统计学习方法的三要素

联合概率分布于假设空间

监督学习假设输入输出的随机变量X和Y遵循联合概率分布P(X,Y)

P(X,Y)表示 分布函数分布密度函数

输入空间和输出空间的映射集合,就是 假设空间

监督学习的模型可以是 概率模型非概率模型,由 **条件概率分布P(Y|X)**或 **决策函数Y=f(X)**表示

损失函数

常用的损失函数有: 0-1损失函数平方损失函数绝对损失函数对数(似然)损失函数

  • 损失函数的期望:
    • Rexp(f)=Ep[L(Y,f(x))]=xyL(y,f(x))P(x,y)dxdyR_{exp}(f) = E_p[L(Y,f(x))] = \int_{x*y}L(y,f(x))P(x,y){\rm d}x{\rm d}y

    • 称为 风险函数期望损失

学习的目标就是选择期望风险最小的模型

  • 模型f(x)关于训练数据集的平均损失:
    • Remp(f)=1Ni=1NL(yi,f(xi))R_{emp}(f) = \frac{1}{N} \sum_{i=1}^{N}L(y_i,f(x_i))

    • 称为 经验风险经验损失
  • 当N趋于无穷时,经验风险R_emp趋于期望风险R_exp

经验风险最小化与结构风险最小化

  • 当假设空间、损失函数以及训练数据集确定的情况下,经验风险函数就可以确定
  • 经验风险最小化策略认为:
    • 经验风险最小的模型就是最优的模型
    • min1Ni=1NL(yi,f(xi))min \frac{1}{N} \sum_{i=1}^{N} L(y_i,f(x_i))

    • 极大似然估计*就是经验风险最小化的一个例子
    • 当样本容量很小时,容易过拟合
  • 结构风险最小化策略:
    • 是为了防止过拟合提出来,等价于正则化,结果风险最小的模型就是最优模型
    • 在经验风险上加上表示模型复杂度的 正则化项罚项
    • min1Ni=1NL(yi,f(xi))+λJ(f)min \frac{1}{N} \sum_{i=1}^{N} L(y_i,f(x_i)) + λJ(f)

      • J(f)是模型的复杂度,λ用于权衡经验风险和模型复杂度
    • 贝叶斯估计的最大后验概率估计就是结构风险最小化的一个例子

正则化

正则化可以取不同的形式

  • L2范数:
    • L(w)=1Ni=1N(f(xi;w)yi)2+λ/2w2L(w) = \frac{1}{N} \sum_{i=1}^{N}(f(x_i;w)-y_i)^2 + λ/2 ||w||^2

    • ||w||表示参数向量的L2范数
  • L1范数:
    • L(w)=1Ni=1N(f(xi;w)yi)2+λw1L(w) = \frac{1}{N} \sum_{i=1}^{N}(f(x_i;w)-y_i)^2 + λ||w||_1

    • ||w||_1表示参数向量的L1范数

上述公式中,第1项的经验风险较小的模型可能较复杂(有多个非零的参数),这时第2项的模型复杂度会较大。正则化的作用是选择经验风险与模型复杂度同时较小的模型。

交叉验证

交叉验证的基本想法是重复地使用数据,把给定数据进行切分,将切分的数据集组合成 训练集测试集,在此基础上反复地进行训练、测试以及模型选择。

  • 简单交叉验证

    • 数据分成2部分,训练集和测试集
    • 用训练集在各种条件下训练得到不同模型,在测试集上评价各个模型的测试误差,选出测试误差最小的模型
  • S折交叉验证

    • 将数据切分成S个 互不相交大小相同的子集
    • 利用S-1个子集训练模型,剩下的测试模型
    • 将这一过程重复不同的S次,最后选出S次平均测试误差最小的模型
  • 留一交叉验证

    • S交叉验证的特殊情况,S=N,N是给定数据集的容量

泛化能力

指由该方法学习到的模型对位置数据的预测能力

该模型对未知数据预测的误差就是泛化误差

  • fNf_N 的泛化能力

    • R(fN)=E[L(Y,fN(X))]R(f_N) = E[L(Y,f_N(X))]

  • 泛化误差上界

    • R(f)<=R^(f)+ε(d,N,δ)R(f) <= \hat{R}(f) + ε(d,N,δ)

    • R^(f)\hat{R}(f),是经验风险
    • ε(d,N,δ)=12N(logd+log1δ)ε(d,N,δ) = \sqrt{\frac1{2N}(\log{d} + \log{\frac1δ})}

    • 上述不等式左端是R(f)的泛化误差,右端是泛化误差上界。
    • 泛化误差上界中,第1项是训练误差,训练误差越小,泛化误差越小;第2项是N的单调递减函数,假设空间中包含的函数越多,其值越大。

注:PDF P32公式推导

生成模型与判别模型

  • 生成模型

    • 生成方法由数据学习联合概率分布P(X,Y),然后求出条件概率P(Y|X)作为预测的模型,即生成模型
    • 模型表示了给定输入X产生输出Y的生成关系
    • 朴素贝叶斯算法和马尔科夫模型
  • 判别模型

    • 判别方法由数据直接学习决策函数f(x)或者条件概率分布P(Y|X)作为预测模型,即判别模型
    • 模型是对给定输入X,应该预测什么样的输出Y
    • k近邻、感知机、决策树、逻辑斯蒂回归模型、最大熵模型、支持向量机、提升方法和条件随机场

分类问题和标注问题

  • TP:正类预测为正类数
  • FN:正类预测为负类数
  • FP:负类预测为正类数
  • TN:负类预测为负类数
  • 精准率:$$P = \frac{TP}{TP+FP}$$
  • 召回率:$$R = \frac{TP}{TP+FN}$$
  • 精准率和召回率的调和均值 F1F_1:$$F1 = \frac{2TP}{2TP+FP+FN}$$
    • 精准率和召回率都高时, F1F_1 值也会高

回归问题

  • 回归用于预测输入变量(自变量)和输出变量(因变量)之间的关系。
  • 回归模型表示从输入变量到输出变量之间映射的函数。
  • 按照输入变量的个数分为: 一元回归多元回归
  • 按照模型的类型分为: 线性回归非线性回归
  • 回归问题最常用的损失函数是 平方损失函数

统计学习方法-第二章:感知机

感知机模型

f(x) = sign(w·x + b),称为感知机,w和b是感知机模型参数,w为 权重权值向量,b是偏置,w·x是w和x的内积,sign是函数符号。

感知机是一种 线性分类器,属于 判别模型

对于特征空间中的一个超平面S,b是超平面的截距。这个超平面将特征空间划分为两个部分。位于两个部分的点,分别被分成正负两类。 因此超平面S被称为 分离超平面

感知机学习策略

  • 数据集的线性可分性

    • 给定一个数据集,xiRn,yi+1,1x_i∈R^n,y_i∈{+1,-1},如果存在某个超平面S能够将数据集的正实例点和负实例点完全正确地分到S的两侧,则称数据集T为 线性可分数据集,否则为 线性不可分
  • 感知机学习策略

    • 假设训练数据集是线性可分的,感知机学习的目标是找到分离超平面
    • sign(w·x+b)学习的损失函数为:$$L(w,b) = -\sum_{x_i∈M}y_i(w·x_i+b)$$
      • 其中,M为误分类点的集合,这个损失函数就是感知机学习的经验风险函数

感知机学习算法

感知机学习问题转换为求解损失函数式的最优化问题,最优化的方法是随机梯度下降法

  • 感知学习算法的原始形式

    1. 选取初值 w0,b0w_0,b_0
    2. 在训练集中选取数据 (xi,yi)(x_i,y_i)
    3. 如果yi(wxi+b)<=0y_i(w·x_i+b)<=0:$$w \leftarrow w + ηy_ix_i$$$$b \leftarrow b + ηy_i$$
    4. 转至(2),直至训练集中没有误分类点
  • 算法的收敛性

经过有限次迭代可以得到一个将训练数据集完全划分的分离超平面及感知机模型

  • 记号:

    • w^=(wT,b)T\hat{w} = (w^T,b)^T
    • x^=(xT,1)T\hat{x} = (x^T,1)^T
    • w^x^=wx+b\hat{w}·\hat{x} = w·x + b
  • 若训练集T是线性可分的,则

    • 存在 w^_opt=1||\hat{w}\_{opt}|| = 1的超平面 w^_optx^=0||\hat{w}\_{opt}·\hat{x} = 0,且存在γ>0:$$y_i(\hat{w}_{opt}·\hat{x} = y_i(w_{opt}·x_i + b_{opt}) >= γ$$
    • R=max1iNx^R = \max \limits_{1\le i\le N}||\hat{x}||,则感知机算法在训练集上误分类次数k满足:$$k\le{(\frac Rγ)}^2$$
  • 感知机算法的对偶形式

    • 通过感知学习算法的原始形式可知:$$w \leftarrow w + ηy_ix_i$$$$b \leftarrow b + ηy_i$$
    • 逐步修改w,b,设修改n次,则w,b关于(xi,yi)(x_i,y_i)的增量分别是aiyixia_iy_ix_iaiyia_iy_i,ai=niηa_i=n_iη,最后学习到的w,b可以分别表示为:$$w=\sum_{i=1}^N{a_iy_ix_i}$$, $$b=\sum_{i=1}^N{a_iy_i}$$
    • ai0a_i \ge 0,当η=1时,表示第i个实例点由于误分而进行更新的次数
    • 感知机模型f(x)=signj=1Najyjxjx+bf(x) = sign{\sum_{j=1}^N{a_jy_jx_j}·x + b}
      1. a0a\rightarrow 0,b0b\rightarrow 0
      2. 在训练集中选取数据(xi,yi)(x_i,y_i)
      3. 如果yi(j=1Najyjxjxi+b)0y_i(\sum_{j=1}^N{a_jy_jx_j}·x_i + b) \le 0,$$a_i \leftarrow a_i +η$$, $$b \leftarrow b + ηy_i$$
      4. 转至2直到没有错误分类
    • 对偶形式中训练实例仅以内积形式出现。可以预先将训练集中实例间的内积计算出来并以矩阵形式存储——Gram矩阵
      • G=[xixj]_NNG = [x_i·x_j]\_{N*N}
  • 由于异或计算不具有线性可分性,所以不能由感知机表示

感知机算法-实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
import numpy as np
def getDataset():
# 数据来源:《统计学习方法》例2.1
dataset = np.array([[3,3],[4,3],[1,1]])
labels = np.array([[1],[1],[-1]])
return dataset,labels

def train_1(dataset,labels,w=np.array([[0,0]]),b=0):
# 感知学习算法的原始形式
(num_y,num_x) = dataset.shape
while True:
tab = 0
for i in range(num_y):
re = labels[i]*(np.dot(w,dataset[i])+b)
if re <= 0:
tab += 1
w += labels[i]*dataset[i]
b += labels[i]
break
if tab == 0:
break
return w,b

def train_2(dataset,labels,w=np.array([[0,0]]),b=0):
# 感知机算法对偶形式
(num_y,num_x) = dataset.shape
dataMat = np.dot(dataset,dataset.T)
a = np.array([[0],[0],[0]])
while True:
tab = 0
for i in range (num_y):
x = np.sum(a*labels*dataMat[:,i:i+1])
re = labels[i]*(x+b)
if re <= 0:
tab += 1
a[i] += 1
b += labels[i]
break
if tab == 0:
break
w = np.sum(a*labels*dataset,axis=0)
return w,b

if __name__ == '__main__':
dataset,labels = getDataset()
w,b = train_1(dataset,labels)
w,b = train_2(dataset,labels)
  • 注:书例2.2求解表格第4列应为:a1=1,a2=0,a3=3,b=0

统计学习方法-第三章:k近邻法

k-NN是一种基本分类与回归的方法,其输入为实例的特征向量,输出为实例的类别,可以取多类。

k近邻算法

  • 算法

    1. 根据给定的距离度量,在训练集T中找出与x最邻近的k个点,涵盖这k个点的x的邻域记作Nk(x)N_k(x)
    2. Nk(x)N_k(x)中,根据分类的决策规则(如多数表决)决定x的类别y:$$y=\arg \max\limits_{c_j} \sum_{x_i∈N_k(x)} I(y_j=c_j), i=1,2,…,N;j = 1,2,…,K$$,其中,I为指示函数,即当yi=ciy_i=c_i时I为1,否则I为0。
  • 当k=1时,称为 最近邻算法,选取与x最邻近的类作为x的类

  • k近邻模型

    • 模型由三个基本要素——距离度量、k值的选择和分类决策规则决定
    • 当训练集、距离度量、k值即分类决策规则确定时,对于任何一个新的输入实例,它所属的类唯一地确定
  • 距离度量

    • k近邻模型的特征空间一般是n维的实数向量空间RnR^n,使用的距离是欧式距离,也可以是其他距离,如LpL_p距离或Minkowski距离
  • k值的选择

    • 若k值选择较小,则容易产生过拟合;相反,选择k值较大,则容易欠拟合
    • 在应用中,k一般取一个较小的值。通常采用 交叉验证来选取最优的k值
  • 分类决策规则

    • 一般为多数表决法,即输入实例的k个邻近的训练实例中的多数类决定输入实例的类
    • 多数表决等价于经验风险最小化

K近邻算法-实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
import numpy as np
from collections import defaultdict

def get_DataSet():
# 该函数用于设定数据集
data = [
('A', [1.0, 1.1]), ('A', [1.0, 1.0]),
('B', [0.0, 0.0]), ('B', [0.0, 0.1]),
]
y, x = zip(*data)
return np.array(x), y

def knnclassify(inputX, dataset, labels, k):
datasetSize = dataset.shape[0]
# 用矩阵来计算欧式距离、算法参照《机器学习实战》
disMat = np.sqrt(np.sum(np.square(inputX-dataset),axis=1))
classdict = defaultdict(int)
Sort = disMat.argsort()
for i in range(k):
label = labels[Sort[i]]
classdict[label] += 1
classSort = sorted(classdict,key=lambda x:classdict[x],reverse=True)
return classSort[0][0]

if __name__ == '__main__':
k,(group, labels) = 3,get_DataSet()
inputX = [0,0]
result = knnclassify(inputX,group, labels,k)
print (result)

k近邻法的实现:kd树

k近邻法最简单的实现方法是线性扫描,但当训练集很大的时候十分耗时

为提高k近邻搜索的效率,可以使用特殊的结构存储方法,例如kd树

  • 构造kd树

    • 是一个二叉树,表示对k维空间的一个划分
    1. 构造根节点,使根节点对应于k维空间中包含所有实例点的超矩形区域
    2. 在超矩形区域(结点)上选择一个坐标轴(对深度为j的结点,选择坐标轴l=j(modk)+1l=j(mod k)+1)和在此坐标轴上的一个切分点,确定一个超平面,这个超平面通过选定的切分点且垂直于 选定的坐标轴,将当前的超矩形区域分为两个左右的区域(子结点)
    3. 递归2,不断对k维空间进行切分,直到子区域没有实例时终止
    4. 通常,选择的切分点为——训练实例点在坐标轴上的中位数,这样得到的kd树是平衡的,但搜索时效率不一定最优
  • 搜索kd树

    • 在kd树中找出包含目标点x的 叶结点
    • 以此叶结点为“当前最近点”
    • 递归地向上回退,在每个结点进行以下操作:
      • 如果该结点保存的实例点比当前最近点距离目标点更近,则以该实例点为“当前最近点”
      • 当前最近点一定存在于该结点一个子结点对应的区域。检查该子结点的兄弟结点对应区域是否有更近的点。若有,则移动到另一个子结点,递归进行最近邻搜索;若无,则向上回退
    • 当回退到根结点则结束,最后的“当前最近点”即为x的最近邻点
  • 如果实例点随机分布,则kd树的平均计算复杂度是O(logN),N是训练实例数

kd树的创建与搜索

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
import numpy as np
def createDataSet():
# 该函数用于设定数据集
group = np.array([[2,3], [5,4],[9,6],[4,7],[8,1],[7,2]])
return group

class kdNode:
def __init__(self,data,depth=0,left=None,right=None):
self.data = data
self.depth,self.left,self.right = depth,left,right

@classmethod
# 创建kd树
def creatKDtree(cls,group,depth):
# 解析当前 space 尺寸
num_samples, num_dimensons = group.shape
if num_samples == 0:
return None
# 选取坐标轴,对应书上算法3.2 (2)
l = depth % num_dimensons + 1
sortGroup = group[group[:,l-1].argsort()]
midindex = sortGroup.shape[0] // 2
# 划分结点
node = cls(sortGroup[midindex],depth)
node.left = cls(sortGroup[:midindex],depth+1)
node.right = cls(sortGroup[midindex+1:],depth)
return node

# 中序遍历打印kd树
def showKDTree(self):
print(self.data,self.depth)
for child in [self.left,self.right]:
if child is not None:
child.showKDTree()

# 搜索kd树
def searchKDTree(self,inputX,k):
# 算法3.3(1),寻找到叶子结点作为当前最近点
lastNode = self
l = root.depth % k + 1
child = self.left if inputX[l-1] < root.data[l-1] else self.right
if child is not None:
lastNode = child.searchKDTree(inputX,k)
# 计算距离
disNode,lastNodedis = 0,0
disNode = np.sum(np.square(inputX[:k]-self.data[:k]))
lastNodedis = np.sum(np.square([inputX[:k]]-lastNode.data[:k]))
# 更新当前最近结点
if disNode < lastNodedis:
lastNode = self
return lastNode

if __name__ == '__main__':
group = createDataSet()
root = kdNode.creatKDtree(group,0)
k = group.shape[1]
root.showKDTree()
inputX = [7,3]
ln = root.searchKDTree(inputX,k)
print(ln.data)

统计学习方法-第四章:朴素贝叶斯

朴素贝叶斯法是基于贝叶斯定理和特征条件独立假设的分类方法

朴素贝叶斯的学习与分类

  • 基本方法
    • 朴素贝叶斯法通过训练数据集学习联合概率分布P(X,Y)
    • 通过学习先验概率分布P(Y=ck),k=1,2,...,KP(Y=c_k),k=1,2,...,K以及条件概率分布$$P(X=x|Y=c_k)=P(X{(1)}=x{(1)},…X{(n)}=x{(n)}|Y=c_k),k=1,2,…,K$$,学习到联合概率分布P(X,Y)
    • 朴素贝叶斯对条件概率分布做了条件独立性的假设,即:$$P(X=x|Y=c_k)==P(X{(1)}=x{(1)},…X{(n)}=x{(n)}|Y=c_k)=\prod_{j=1}nP(X{(j)}=x^{(j)}|Y=c_k)$$
    • 朴素贝叶斯属于生成模型
    • 后验概率计算根据贝叶斯定理进行:$$P(Y=c_k|X=x)=\frac{P(X=x|Y=c_k)P(Y=c_k)}{\sum_kP(X=x|Y=c_k)P(Y=c_k)}$$,即:$$P(Y=c_k|X=x)=\frac{P(Y=c_k)\prod_{j}P(X{(j)}=x{(j)}|Y=c_k)}{\sum_kP(Y=c_k)\prod_{j}P(X{(j)}=x{(j)}|Y=c_k)}$$
    • 则,朴素贝叶斯分类器为:$$y=\arg\max \limits_{c_k}P(Y=c_k)\prod_{j}P(X{(j)}=x{(j)}|Y=c_k)$$

朴素贝叶斯的参数估计

  • 极大似然估计

    • 先验概率P(Y=ck)P(Y=c_k)的极大似然估计是:$$P(Y=c_k)=\frac{\sum_{i=1}^N I(y_i=c_k)}{N},k=1,2,…,N$$
    • 设第j个特征x(j)x_{(j)}可能取值集合为aj1,aj2,...,ajSj{a_{j1},a_{j2},...,a_{jS_j}},条件概率的极大似然估计是: $$P(X{(j)}=a_{jl}|Y=c_k)=\frac{\sum_{i=1}N I(x_i{(j)}=a_{jl},y_i=c_k)}{\sum_{i=1}N I(y_i=c_k)},j=1,2,…,n;l=1,2,…,S_j;k=1,2,…,K $$,其中,xi(j)x_i^{(j)}是i样本的第j个特征;ajla_jl是第j个特征可能取的第l个值;I为指示函数
  • 学习与分类算法

    1. 计算先验概率和条件概率
    2. 对于给定的实例x计算:$$P(Y=c_k)\prod_{j}P(X{(j)}=x{(j)}|Y=c_k), k=1,2,…,K$$
    3. 确定实例x的类
  • 贝叶斯估计

    • 用极大似然估计可能会出现估计的概率值为0的情况,会影响到后验概率的计算结果
    • 所以采用贝叶斯估计,条件概率的贝叶斯估计是:$$P_λ(X{(j)}=a_(jl)|Y=c_k)=\frac{\sum_{i=l}NI(x_i{(j)}=a_{jl},y_i=c_k)+λ}{\sum_{i=1}NI(y_i=c_k)+S_jλ}$$,其中,λ>=0,当λ=0时,就是极大似然估计。常取λ=1。
    • 先验概率的贝叶斯分布:$$P_λ(Y=c_k)=\frac{\sum_{i=1}^NI(y_i=c_k)+λ}{N+Kλ}$$
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
# 朴素贝叶斯·离散数据
import numpy as np
from collections import Counter
def getDataset():
# 数据来源:《统计学习方法》例4.1
group = [[1,'S',-1],[1,'M',-1],[1,'M',1],[1,'S',1],[1,'S',-1],[2,'S',-1],[2,'M',-1],[2,'M',1],[2,'L',1],[2,'L',1],[3,'L',1],[3,'M',1],[3,'M',1],[3,'L',1],[3,'L',-1]]
return group
def classifyNB(inputx, group):
# 估计先验概率
num_k = len(group[0]) - 1
lables = [i[num_k] for i in group]
num_y = len(set(lables))
num_data = len(group)
P_y = []
for i in set(lables):
P_y.append([i,lables.count(i)])
# 估计条件概率
P_x = []
for y in P_y:
# 统计不同类型的X样本
P_x.append([x for x in group if x[2] == y[0]])
L = []
for x_y in P_x:
p_x_y = Counter()
for x in x_y:
for j in range(num_k):
p_x_y[x[j]] += 1
L.append(p_x_y)
print (L)
result = []
for i in range(num_y):
result.append(P_y[i][1]/num_data * L[i][inputx[0]]/P_y[i][1] * L[i][inputx[1]]/P_y[i][1])

print (result)
if __name__ == '__main__':
group = getDataset()
inputx = [2,'S']
classifyNB(inputx,group)

统计学习方法-第5章:决策树

决策树是一种基本的分类与回归方法,呈树形结构。在分类问题中表示基于特征对实例进行分类的过程

优点:模型具有可读性,分类速度快

决策树模型与学习

  • 决策树模型

    • 分类决策树模型是一种描述对实例进行分类是树形结构
    • 由结点(内部结点、叶结点)、有向边组成
    • 内部结点表示一个特征或属性,叶结点表示一个类
  • 决策树和if-then规则

    • 由决策树的根结点到叶结点的每一条路径构建一条规则
    • 路径内部结点的特征对应规则的条件,叶结点对应规则的结论
    • 决策树的路径或其对应的if-then规则集合具有一个重要的性质:互斥且完备
  • 决策树与条件概率分布

    • 决策树还表示给定特征条件下类的条件概率分布
    • 决策树的每一条路径对应划分中的一个单元
    • 决策树表示的条件概率分布由各个单元给定条件下的类的条件概率分布组成
  • 决策树学习

    • 其本质是从训练数据集中归纳出一组分类规则
    • 需要训练一个与训练数据矛盾较小的决策树,同时具有很好的泛化能力
    • 决策树学习也是从训练数据集估计条件概率模型
    • 决策树学习的损失函数通常是正则化的极大似然函数
    • 决策树学习的策略是以损失函数为目标函数的最小化
    • 算法
      • 构建根结点,将所有训练数据都放在根结点
      • 选择一个最优特征,将训练数据集分割成子集,使得各个子集有一个在当前条件下最好的分类
      • 如果还有子集不能被正确分类,那么就对这些子集选择新的最优特征,继续分割
      • 如此递归下去,直到所有子集都被正确分类,或者没有合适的特征为止
      • 为防止过拟合,需要对已生成的树自下而上进行 剪枝。即,去掉过于细分的叶结点,使其回退到父结点

特征选择:

对训练数据集D,计算其每个特征的信息增益,选择信息增益最大的特征。选取对于训练数据具有分类能力的特征。

  • 信息增益
    • 熵是表示随机变量不确定性的度量
      • 设X是一个取有限个值的离散随机变量,其概率分布P(X=xi)=pi,i=1,2,..,nP(X=x_i)=p_i,i=1,2,..,n
      • 则随机变量X的熵为:$$H(X)=-\sum_{i=1}^n p_i \log{p_i}$$
      • 上述式子中,若pi=0p_i=0,则定义0log0=0,其中的log通常以2或e为底
      • 此时熵的单位被称作比特(bit)或者纳特(nat)
      • 熵越大,随机变量的不确定性就越大:0H(p)logn0≤H(p)≤\log{n}
    • 条件熵
      • 设有随机变量(X,Y),其联合概率分布:$$P(X=x_i,Y=y_i)=p_{ij}$$
      • 条件熵H(Y|X)表示已知随机变量X的条件下,随机变量Y的不确定性
      • 随机变量Y的条件熵,其定义为X给定条件下Y的条件概率分布的熵对X的数学期望:$$H(Y|X)=\sum_{i=1}^np_iH(Y|X=x_i)$$, 其中pi=P(X=xi)p_i=P(X=x_i)
      • 当熵和条件熵中的概率由数据估计得到时,称为 经验熵经验条件熵
    • 信息增益
      • 特征A对训练数据集D的信息增益g(D,A),定义为集合D的经验熵H(D)与特征A给定条件下D的经验条件熵H(D|A)只差,即:$$g(D,A)=H(D)-H(D|A)$$
      • g(D,A)又称为训练数据集中 类与特征互信息
      • 其意义是:由于特征A对数据集D分类不确定性减少的程度
      • 特征选择:对训练数据集D,计算其每个特征的信息增益,选择信息增益最大的特征
    • 信息增益算法
      • 计算数据集D的经验熵H(D):$$H(D)=-\sum_{k=1}^k \frac{|C_k|}{|D|} \log_2{\frac{|C_k|}{|D|}}$$
      • 计算特征A对数据集D的经验条件熵H(D|A):$$H(D|A)=\sum_{i=1}^n \frac{|D_i|}{D}H(D_i) = -\sum_{i=1}^n\frac{|D_i|}{D}\log_2{\frac{|D_{ik}|}{|D_i|}}$$
      • 计算信息增益
    • 信息增益比
      • 是特征选择的另一准则
      • gR(D,A)=g(D,A)H(D)g_R(D,A)=\frac{g(D,A)}{H(D)}

决策树生成

  • ID3算法

    • 在决策树各个结点上应用信息增益准则选择特征
    1. 若训练集D中所有实例属于同一类CKC_K,则T为单结点树,并将类CkC_k作为该结点的类标记
    2. 若特征集A=∅,则T为单结点树,并将D中实例数最大的类CkC_k作为该结点的类标记
    3. 否则,计算A中各特征对D的信息增益,选择最大的特征AgA_g
    4. 如果AgA_g的信息增益小于阈值ε,则置T为单结点树,并将D中实例数最大的类CkC_k作为该结点的类标记
    5. 否则,对A_g的每一可能值aia_i,依Ag=aiA_g=a_i将D分割为若干非空子集DiD_i,将DiD_i中实例数最大的类作为标记,构建子结点
    6. 对第i个子结点,依DiD_i为训练集,以AAgA-{A_g}为特征集,递归调用1-5
  • C4.5算法

    • 与ID3算法不同,C4.5算法用信息增益比来选择特征
    • 其余步骤相同

ID3算法和C4.5算法实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
import numpy as np
def getDataset():
# 数据来源:《统计学习方法》例5.1
# 这里使用0,1,2分表表示青年中年老年,和信贷情况的一般,好,非常好
# 只用0,1分别表示否、是
group = [[0,0,0,0,0],[0,0,0,1,0],[0,1,0,1,1],[0,1,1,0,1],[0,0,0,0,0],[1,0,0,0,0],[1,0,0,1,0],[1,1,1,1,1],[1,0,1,2,1],[1,0,1,2,1],[2,0,1,2,1],[2,0,1,1,1],[2,1,0,1,1],[2,1,0,2,1],[2,0,0,0,0]]
name = ['年龄','有工作','有自己的房子','信贷情况']
return group,name

def getEntropy(category,num_y, totalnum):
C_k = [category.count(i) for i in num_y]
H_D = -np.sum([i/totalnum * np.log2(i/totalnum) for i in C_k if i])
return H_D
def getConditionalEntropy(feature, category, num_y,totalnum):
num_k_y = list(set(feature))
fea_A = []
for i in num_k_y:
a = []
for f in range(len(feature)):
if i == feature[f]:
a.append(category[f])
fea_A.append(a)
# print (fea_A)
D_k = [feature.count(i) for i in num_k_y]
H_D_A = np.sum(D_k[i]/totalnum * getEntropy(fea_A[i],num_y,D_k[i]) for i in range(len(D_k)))
return H_D_A

def InformationGain(dataset,category, k):
# 该函数用于计算信息增益值,这里的k表示选取的特征下标
# 经验熵H(D)
num_y = list(set(category))
num_x = len(dataset)
H_D = getEntropy(category,num_y, num_x)

# 经验条件熵H(D|A)
feature = [x[k] for x in dataset]
H_D_A = getConditionalEntropy(feature, category, num_y,num_x)
return (H_D-H_D_A)

def InformationGainRate(dataset,category,k):
# 该函数用于计算信息增益比
# 经验熵H(D)
num_y = list(set(category))
if len(num_y) == 1:
return
num_x = len(dataset)
H_D = getEntropy(category,num_y, num_x)

# 经验条件熵H(D|A)
feature = [x[k] for x in dataset]
H_D_A = getConditionalEntropy(feature, category, num_y,num_x)
return (H_D-H_D_A)/H_D

class DecisionTreeNode:
def __init__(self, data,depth = 0, label = -1):
self.data,self.depth,self.label, self.child = data,depth,label,[]
@classmethod
def CreateTree(cls,dataset,name,depth=0,label=-1):
num_feature = len(name)
category = [x[-1] for x in dataset]

if len(set(category)) == 1 or num_feature == 0:
node = cls("None",depth)
return node
# 求信息增益
F_G = [InformationGain(dataset,category, i) for i in range(num_feature)]
# 求信息增益比
# F_G = [InformationGainRate(dataset,category, i) for i in range(num_feature)]
if F_G is None:
node = cls("None",depth)
return node
A = np.argmax(F_G)
# 划分子集
node = cls(name[A],depth,label)
D = list(set([i[A] for i in dataset]))
subsetData = [[i for i in dataset if i[A] == d] for d in D]

del name[A]
depth += 1
for s in subsetData:
# 此处的label表示child是特征A=label时候的子树
label = s[0][A]
for x in s:
del x[A]
child_node = cls.CreateTree(s,name,depth,label)
node.child.append(child_node)
# print (node.child)
return node
def show(self):
print('%2d: %s: %s' % (self.depth, self.data, self.label))
# print (self.child)
for child in self.child:
child.show()

if __name__ == '__main__':
dataset,name = getDataset()
root = DecisionTreeNode.CreateTree(dataset,name)
root.show()

决策树剪枝

决策树算法递归生成的决策树容易产生过拟合现象,所以需要对生成的决策树进行简化——剪枝

剪枝往往通过极小化决策树的 损失函数代价函数来实现

  • 损失函数

    • Ca(T)=C(T)+aT=t=1Tk=1kNtklogNtkNt+aTC_a(T) = C(T) + a|T| = -\sum_{t=1}^{|T|}\sum_{k=1}^k N_{tk} \log \frac{N_{tk}}{N_t} + a|T|

    • C(T)表示模型对训练数据的预测误差
    • |T|表示模型的复杂度,参数a≥0控制两者之间的影响
    • 该损失函数的极小化等价于正则化的极大似然估计
  • 剪枝算法

    1. 计算每个结点的经验熵
    2. 递归地从树的叶结点向上回缩——设一叶结点回缩到父结点之前和之后的整体树分别为TBT_BTAT_A,则其损失函数:Ca(TA)Ca(TB)C_a(T_A)≤C_a(T_B),则进行剪枝,即将父结点变为新的叶结点
    3. 返回2,直到不能再继续为止

CART算法

分类与回归树(CART)模型,既可以运用于分类也可以运用于回归

是在给定输入随机变量X条件下输出随机变量Y的条件概率分布的学习方法

  • CART算法

    • 决策树生成:基于训练集生成决策树,生成的决策树要尽量大
    • 决策树剪枝:用验证数据集对已生成的树剪枝并选择 最优子树,用 损失函数最小作为剪枝的标准
  • 最小二乘回归树生成算法

    • 回归树模型:f(x)=m=1McmI(xRm)f(x) = \sum_{m=1}^M c_mI(x∈R_m)
    • 假设已将输入空间划分成M个子单元R1,R2,...,RMR_1,R_2,...,R_M,并在每一个单元上有一个固定的输出值cmc_m
    • 单元RmR_mcmc_m的最优值c^_m\hat{c}\_mRmR_m上所有输入实例xix_i对应的输出yiy_i的均值:c^_m=ave(yixiRm)\hat{c}\_m=ave(y_i|x_i∈R_m)
    1. 选择最优切分变量j和切分点s:$$\min\limits_{j,s}[\min \limits{c_1} \sum_{x_i∈R_1(j,s)}(y_i-c_1)^2 + \min \limits_{c_2}\sum_{x_i∈R_2(j,s)} (y_i-c_2)^2]$$(最小损失函数,即平方误差),遍历变量j,对固定的切分变量j扫描切分点s,选择使上述式子达到最小的(j,s)
    2. 用选定的(j,s)划分区域并决定相应的输出值:$$R_1(j,s)={x|x{(j)≤s}},R_2(j,s)={x|x{(j)}>s}$$ $$\hat{c}_m = \frac{1}{N_m}\sum_{x_i∈R_m(j,s)}y_i, x∈R_m,m=1,2$$
    3. 继续对两个子区域调整,重复1,2,直到满足条件
    4. 生成决策树
  • 分类树的生成

    • 基尼指数:分类问题中,假设有K个类,样本点属于第k类的概念为pkp_k,则概率分布的基尼指数定义为:$$Gini§=\sum_{k=1}^k p_k(1-p_k)=1-\sum_{k=1}^k p_k^2$$
      • 对于二分类问题,若样本点属于第1类的概率是p,则:$$Gini§=2p(1-p)$$
      • 对于给定样本D,则:$$Gini(D)=1-\sum_{k=1}k(\frac{|C_k|}{|D|})2$$,其中CkC_k是属于第k类的样本子集,k是类的个数
      • 如果样本集合D根据特征A可以被分成2部分,则在特征条件A下,集合D的基尼指数:$$Gini(D,A)=\frac{|D_1|}{|D|}Gini(D_1)+\frac{|D_2|}{|D|}Gini(D_2)$$
      • Gini(D)表示集合D的不确定性,Gini(D,A)表示经过A分割后的集合D的不确定性。
      • 基尼指数越大,集合样本的不确定性也越大
  • CART生成算法

    1. 设结点的训练数据集为D,对每一个特征A,都计算Gini(D,A)
    2. 对所有可能的特征A以及它们所有可能的切分点a中,选择基尼指数最小的特征及其对应的切分点作为最优特征与最优切分点。依此从现结点生成两个子结点,将训练数据集依特征分配到两个子结点中去
    3. 对两个子结点递归调用1,2
    4. 生成CART决策树

CART算法实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
import numpy as np

def getDataset():
# 数据来源:《统计学习方法》例5.1
# 这里使用0,1,2分表表示青年中年老年,和信贷情况的一般,好,非常好
# 只用0,1分别表示否、是
group = [[0,0,0,0,0],[0,0,0,1,0],[0,1,0,1,1],[0,1,1,0,1],[0,0,0,0,0],[1,0,0,0,0],[1,0,0,1,0],[1,1,1,1,1],[1,0,1,2,1],[1,0,1,2,1],[2,0,1,2,1],[2,0,1,1,1],[2,1,0,1,1],[2,1,0,2,1],[2,0,0,0,0]]
name = ['年龄','有工作','有自己的房子','信贷情况']
return group,name
def GetGini(dataset,k):
dataset_x = [x[k] for x in dataset]
category = [x[-1] for x in dataset]
A_x = list(set(dataset_x))
A_y = list(set(category))
l_x = []
l_y = []
for x in A_x:
l_x.append([dataset_x[i] for i in range(len(dataset_x)) if dataset_x[i] == x])
l_y.append([category[i] for i in range(len(dataset_x)) if dataset_x[i] == x])
num_y = []
for y in l_y:
num_y.append([y.count(i) for i in A_y])
Gini = []
for x in range(len(l_x)):
G_1 = 1.0 - sum([(i/sum(num_y[x]))**2 for i in num_y[x]])
G_2_s = np.sum([num_y[i] for i in range(len(num_y)) if i != x],axis=0)
G_2 = 1.0 - sum([(i/sum(G_2_s))**2 for i in G_2_s])
G = len(l_x[x])/len(dataset_x)*G_1 + (1-len(l_x[x])/len(dataset_x))*G_2
Gini.append(G)
return Gini

class DecisionTreeNode:
def __init__(self, data,depth = 0, label = -1):
self.data,self.depth,self.label, self.child = data,depth,label,[]
@classmethod
def CreateTree(cls,dataset,name,depth=0,label=-1):
num_feature = len(name)
category = [x[-1] for x in dataset]

if len(set(category)) == 1 or num_feature == 0:
y = list(set(category))
# 此处label表示所分的类
label = category[np.argmax([category.count(i) for i in y])]
node = cls("None",depth)
return node
# Cart算法,求基尼指数
Gini = [GetGini(dataset,i) for i in range(num_feature)]
A = np.argmin([np.min(G) for G in Gini])
label = np.argmin(Gini[A])
node = cls(name[A],depth,label)

D = list(set([i[A] for i in dataset]))
subsetData = [[i for i in dataset if i[A] == d] for d in D]

del name[A]
depth += 1
for s in subsetData:
label = s[0][A]
for x in s:
del x[A]
child_node = cls.CreateTree(s,name,depth,label)
node.child.append(child_node)
# print (node.child)
return node

def show(self):
print('%2d: %s: %s' % (self.depth, self.data, self.label))
# print (self.child)
for child in self.child:
child.show()


if __name__ == '__main__':
dataset,name = getDataset()
root = DecisionTreeNode.CreateTree(dataset,name)
root.show()
  • CART剪枝
    • 子树的损失函数为:Ca(T)=C(T)+aTC_a(T)=C(T)+a|T|,C(T)为训练数据的预测误差(如基尼指数),|T|是子树结点个数
    • 对固定的a,一定存在使损失函数Ca(T)C_a(T)最小的子树,表示为TaT_a
    • 当a=0时,整体树最优;当a->∞时,根结点组成的单结点树最优
  • CART剪枝算法
    1. k=0,T=T0k=0, T=T_0
    2. 设a=+∞
    3. 自下而上地对各内部结点t计算C(Tt),TtC(T_t),|T_t|,以及 $$g(t)=\frac{C(t)-C(T_t)}{|T_t|-1}$$ $$a=\min(a,g(t))$$,其中,TtT_t表示以t为根结点的子树,C(Tt)C(T_t)是预测误差,Tt|T_t|TtT_t的叶结点个数
    4. 自下而上访问内部结点t,若有g(t)=a,进行剪枝,并对叶结点t以多数表决法决定其类
    5. k=k+1,ak=a,Tk=Tk=k+1,a_k=a,T_k=T
    6. 如果T不是由根结点单独构成的树,则回到4
    7. 采用交叉验证法,在子树序列T0,T1,...,TnT_0,T_1,...,T_n中选取最优子树TaT_a

统计学习方法-第六章:逻辑斯蒂回归于最大熵模型

逻辑斯蒂回归模型和最大熵模型都属于对数线性模型。

逻辑斯蒂回归模型

  • 逻辑斯蒂分布
    • 设X是连续随机变量,设X服从逻辑斯蒂分布是指X的分布函数和密度函数是:$$F(x)=P(X≤x)=\frac{1}{1+e^{-(x-μ)/γ}}$$ $$f(x)=F’(x)=\frac{e{-(x-μ)/γ}}{γ(1+e{-(x-μ)/γ})^2}$$,其中μ是位置参数,γ>0为形状参数
    • 分布函数F(x)属于逻辑斯蒂函数,是一条S形曲线,以点(μ,12)(μ,\frac{1}{2})为中心对称,满足:F(x+μ)12=F(xμ)+12F(-x+μ)-\frac{1}{2}=-F(x-μ)+\frac{1}{2}
    • 图片见pdf P93 图6.1
  • 二项逻辑斯蒂回归模型
    • 是一种分类模型,由条件概率分布P(Y|X)表示,随机变量X取值实数,随机变量Y变量取值1或0。
    • 条件概率分布:$$P(Y=1|x)\frac{exp(w·x+b)}{1+exp(w·x+b)}$$ $$P(Y=0|x)\frac{1}{1+exp(w·x+b)}$$
    • 比较两个条件概率大小,将x分往概率值较大的一类
    • 一个时间的几率是指该事件发生的概率和该事件不发生的概率的比值,即/fracp1p/frac{p}{1-p},对逻辑斯蒂回归而言则是:$$\log{\frac{P(Y=1|x)}{1-P(Y=1|x)}=w·x}$$
    • 当线性函数的值越接近正无穷,概率值就越接近1;线性函数的值越接近负无穷,概率值越接近0。
  • 模型参数估计
    • 逻辑斯蒂回归模型学习时,可以用极大似然估计法估计模型参数
    • 其对数似然函数为:$$L(w) = \sum_{i=1}^N[y_i(w·x_i)-\log(1+exp(w·x_i))]$$
  • 多项逻辑斯蒂回归
    • 设离散随机变量Y的取值集合是{1,2,…,K},多项逻辑斯蒂回归模型是$$P(Y=k|x)=\frac{exp(w_k·x)}{1+\sum_{k=1}^{K-1}exp(w_k·x)}, k=1,2,…K-1$$ $$P(Y=K|x)=\frac{1}{1+\sum_{k=1}^{K-1}exp(w_k·x)}$$

逻辑斯蒂回归模型梯度下降算法

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
def sigmoid(inX):
return 1.0/(1+exp(-inX))

def gradDrop(group, label):
dataMatrix = np.mat(group)
labelMat = np.mat(label).transpose()
m,n = np.shape(dataMatrix)
alpha = 0.001
maxCycles = 500
weights = np.ones((n,1))
for k in range(maxCycles):
h = sigmoid(dataMatrix*weights)
error = (labelMat - h)
weights = weights - alpha * dataMatrix.transpose() * error
return weights

最大熵模型

  • 最大熵原理
    • 最大熵原理是概率模型学习的一个准则,即在所有可能的概率模型(分布)中,熵最大的模型是最好的模型。
    • 熵:H(P)=xP(x)logP(x)H(P)=-\sum{x}P(x)\log{P(x)},其满足0≤H§≤log|X|,|X|表示X的取值个数,当且仅当X服从均匀分布时,熵最大
    • 最大熵原理认为要选择的概率模型首先要满足 约束条件
    • 在没有更多信息的情况下,不确定的部分认为是 等可能的
  • 最大熵模型的定义
    • 给定训练数据集,可以确定联合分布P(X,Y)的经验分布和边缘分布P(X)的经验分布:$$\tilde{P}(X=x,Y=y)=\frac{v(X=x,Y=y)}{N}$$ $$\tilde{P}(X=x) = \frac{v(X=x)}{N}$$
    • v(X=x,Y=y)表示训练数据中样本(x,y)出现的频数,v(X=x)表示训练数据中输入x出现的频数,N表示训练样本容量
    • 特征函数f(x,y)描述x与y之间的某一个事实,即$$f(x,y)={1,x与y满足某一事实 0,否则}$$
    • 特征函数期望值$$E_p(f)=\sum_{x,y}\tilde {P}(x)P(y|x)f(x,y)$$
    • 如果模型能够获取训练数据中的信息,那么可以假设这两个期望值相等:$$E_p(f)=E_{(\tilde{P})(f)}$$ 或 $$\sum_{x,y}\tilde{P}(x)P(y|x)f(x,y)=\sum_{x,y}\tilde{P}(x,y)f(x,y)$$
    • 一般将上述两个式子作为模型学习的约束条件
  • 最大熵模型学习
    • 拉格朗日推导(P102)
  • 极大似然估计
    • 对偶函数的极大化等价于最大熵模型的极大似然估计(P103)

模型学习的最优化算法

逻辑斯蒂回归模型、最大熵模型学习归结为以似然函数为目标函数的最优化问题,通常通过迭代算法求解

  • 改进的迭代尺度法(IIS)

    • 已知最大熵模型为$$P_w(y|x)=\frac{1}{Z_w(x)} = exp{\sum_{i=1}^n w_if_i(x,y)}$$
    • 其中,Zw(x)=yexp(i=1nwifi(x,y)Z_w(x)=\sum_{y}exp(\sum_{i=1}^n w_if_i(x,y)
    • 对数似然函数为:$$L(w)=\sum_{x,y}\tilde{P}(x,y)\sum_{i=1}^nw_if_i(x,y)-\sum_{x}\tilde{P}(x)\log{Z_w}(x)$$
    1. 对所有i∈{1,2,…,n},取初值wi=0w_i=0
    2. 对每一i∈{1,2,…,n}:
      1. δiδ_i是方程$$\sum_{x,y}\tilde{P}(x)P(y|x)f_i(x,y)exp(δ_if^{\# }(x,y))=E_{\tilde{P}(f_i)}$$ 的解,这里f^{\\# }(x,y)=\sum_{i=1}^nf_i(x,y)
      2. 更新wiw_i, wiwi+δiw_i \leftarrow w_i+δ_i
    3. 如果不是所有的wiw_i都收敛,重复步骤2
    4. 牛顿法通过迭代求得δ_i^\*,使得g(δ_i^\*)=0,迭代公式为:δik+1=δikg(δik)δ_i^{k+1}=-\frac{δ_i^k}{g'(δ_i^k)}
  • 拟牛顿法

    1. 选定初始点w(0)w^{(0)},取B0B_0为正定对称矩阵,置k=0
    2. 计算gk=g(w(k))g_k=g(w^{(k)}),若gk=ε||g_k||=ε,则停止计算,得w=w(k)w^*=w^{(k)},否认转至3
    3. Bkpk=gkB_kp_k=-g_k求出pkp_k
    4. 一维搜索:求λkλ_k使得:$$f(w^{(k)}+λ_kp_k) = \min\limits_{λ≥-}f(w^{(k)}+λp_k)$$
    5. w(k+1)=w(k)+λkpkw^{(k+1)}=w^{(k)}+λ_kp_k
    6. 计算gk+1=g(wk+1)g_k+1=g(w^{k+1}),若gk+1ε||g_{k+1}||<ε,则停止计算,得w=w(k)w^*=w^{(k)};否则求Bk+1B_{k+1}:$$B_{k+1} = B_k+\frac{y_ky_kT}{y_kTδ_k}-\frac{B_kδ_kδ_kTB_k}{δ_kTB_kδ_k}$$,其中$$y_k=g_{k+1}-g_k,δ_k = w{(k+1)}-w{(k)}$$
    7. 置k=k+1,转3

统计学习方法-第七章:支持向量机

支持向量机(SVM)是一种二分类模型,它的基本模型是定义在特征空间上的间隔最大的线性分类器

线性可分支持向量机和硬间隔最大化

  • 线性可分支持向量机

    • 假设输入空间与特征空间为两个不同的空间。输入空间为 欧式空间离散空间,特征空间为 欧式空间希伯尔特空间
    • 线性可分支持向量机、线性支持向量机假设两空间元素一一对应
    • 一般当训练数据集线性可分时,存在无穷个分离超平面; 感知机利用误分类最小的策略,求得分离超平面,但此时解有无穷多; 线性可分支持向量机利用间隔最大化求最优分离超平面,此时解 唯一
  • 函数间隔和几何间隔

    • 对于给定的训练数据集T和超平面(w,b),定义超平面(w,b)和样本点(xi,yi)(x_i,y_i)的函数间隔为:$$\hat{γ}_i = y_i(w·x_i+b)$$
    • 函数间隔可以表示分类预测的正确性及确定度
    • 需要对分离超平面的法向量w加某些约束(例如,||w||=1),这样使得间隔确定,此时函数间隔称为 几何间隔
    • 点A(xi,yi)(x_i,y_i)与超平面(w,b)的距离,记作γiγ_i:$$γ_i=y_i(\frac{w}{||w||}·x_i + \frac{b}{||w||})$$,其中,||w||是w的L2范数。即超平面关于点A的 几何间隔
    • 超平民(w,b)关于训练数据集T的几何间隔为,超平面关于数据集的所有几何间隔最小值
    • 函数间隔与几何间隔的关系:$$γ_i=\frac{\hat{γ}_i}{||w||}$$ $$γ=\frac{\hat{γ}}{||w||}$$
  • 间隔最大化

    • 支持向量机学习的基本想法是求解能够 正确划分训练数据集并且几何间隔最大的分离超平面
  • 最大间隔法

    1. 构造并求解约束最优化问题: $$\min\limits_{w,b} = \frac{1}{2}||w||^2$$,$$s.t. y_i(w·x_i+b)-1≥0,i=1,2,…,N$$
    2. 由此得到分离超平面:$$w*·x+b* = 0$$
  • 若训练数据集T线性可分,则可将训练数据集中的样本点完全分开的最大间隔分离超平面 存在且唯一

  • 支持向量和间隔边界

    • 在线性可分的情况下,训练数据集样本点中与分离超平面距离最近的样本点的实例称为 支持向量。支持向量满足:yi(wxi+b)1=0y_i(w·x_i+b)-1=0
    • H1和H2(图7.3)与w平行,且没有点落在他们中间,H1和H2的距离称为 间隔,等于2w\frac{2}{||w||},H1和H2称为 间隔边界

学习的对偶算法(P119)

应用拉格朗日的对偶性,通过求解对偶问题,得到原始问题的最优解

支持向量一定在间隔边界上

  • 线性可分支持向量机算法
    1. 构造并求解约束最优化问题$$\min_\limits{a} \frac{1}{2}\sum_{i=1}N\sum_{j=1}N a_ia_jy_iy_j(x_i·x_j)-\sum_{i=1}^N a_i$$ $$s.t.\sum_{i=1}^N a_i y_i =0$$ $$a_i≥0, i=1,2,…,N$$ 求得最优解 a^{\*} = (a_1^{\*}, a_2^\*,...,a_N^{\*})^T
    2. 计算$$w{*}=\sum_{i=1}N a_i^{*}y_ix_i$$ 并选择a^{\*}的一个正分量,计算$$b{*}=y_j-\sum_{i=1}N a_i^{*}y_i(x_i·x_j)$$
    3. 求得分离超平面 w^{\*}·x+b=0,即:$$w^{*} = \sum_(i=1)Na_i{*}y_ix_i$$ $$b^{*} = y_j-\sum_{i=1}Na_i{*}y_i(x_i·x_j)$$

线性支持向量机与软间隔最大化

  • 线性支持向量机

    • 需要修改硬间隔最大化,使他成为软间隔最大化,将其扩展到线性不可分问题
    • 线性不可分说明某些样本点 不能满足 函数间隔大于等于1的条件
    • 对每一个样本点引进松弛变量ζi0ζ_i≥0,约束条件变成:$$y_i(w·x_i+b)≥1-ξ$$,目标函数变成 $$\frac{1}{2}||w||2+C\sum_{j=1}Nξ_i$$,其中C>0称为惩罚函数
    • 线性不可分的线性支持向量机学习问题变成:$$\min\limits_{w,b,ξ}\frac{1}{2}+C\sum_{i=1}^Nξ_i$$ $$s.t. y_i(w·x_i+b)≥1-ξ_i, i=1,2,…,N$$ $$ξ_≥0, i=1,2,…,N$$
    • 关于(w,b,ξ)的解存在,w的解唯一,而b不唯一,其存在于一个区间
  • 学习的对偶算法

    • a^\*=(a_1^\*,a_2^\*,...a_N^\*)是对偶问题的一个解,则:$$w*=\sum_{i=1}Na_i^*y_ix_i$$ $$b*=y_j-\sum_{i=1}Ny_ia_i^*(x_i·x_j)$$
    • 分离超平面:$$\sum_{i=1}^N a_i*y_i(x·x_i)+b*=0$$
  • 线性支持向量机学习算法

    1. 选择惩罚函数C>0,构造并解凸二次规划问题:$$\min\limits_{a} \frac{1}{2}\sum_{i=1}N\sum_{j=1}N a_ia_jy_iy_j(x_i·x_j)-\sum_{i=1}^N a_i$$ $$s.t. \sum_{i=1}^Na_iy_i=0 $$ $$0≤a_i≤C, i=1,2,…,N$$
    2. 计算w^\*=\sum_{i=1}a_i^\*y_ix_i,选择 a^\* 的一个分量 a_j^\*适合条件0<a_j^\*<C,计算$$b*=y_j-\sum_{i=1}Ny_ia_i^*(x_i·x_j)$$
    3. 求分离超平面
  • 支持向量

    • 软间隔的支持向量 xix_i 或者在 间隔边界上,或者在间隔边界与分离超平面 之间,或者在 分离超平面误分一侧
    • a^\*=C,则ξi=0ξ_i=0, 则支持向量机在间隔边界上
    • a^\*<C,则0<ξi0ξ_i<0, 则分类正确,支持向量机在间隔边界与分离超平面之间
    • a^\*=C,ξi=1ξ_i=1, 则支持向量机在分离超平面上
    • a^\*=C,ξi1ξ_i>1, 则支持向量机在分离超平面误分一侧
  • 合页损失函数

    • 合页损失函数:$$L(y(w·x+b))=[1-y(w·x+b)]_+$$, "+"表示$$[z] = {z,z>0 0,z≤0}$$
    • 线性支持向量机原始最优化问题:$$\min\limits_{w,b,ξ}\frac{1}{2} ||w||2+C\sum_{i=1}Nξ_i$$ 等价于最优化问题$$\min\limits_{w,b}\sum_{i=1}^N[1-y_i(w·x_i+b)]_+ + λ||w||^2$$

非线性支持向量机与核函数

  • 非线性分类问题

    • 非线性分类问题是指用非线性模型才能很好分类的问题
    • 对于训练数据集T,如果RnR^n中的一个超曲面将正负例正确地分开,则称这个问题为非线性可分问题
    • 采取非线性变换,将非线性问题转换成线性问题,通过解变换后的线性问题的方法求解原来的非线性问题
    • 设原空间为XR2X\subset R^2, x=(x(1),x(2))TXx=(x^{(1)},x^{(2)})^T∈X,新空间ZR2,z=(z(1),z(2))TZZ\subset R^2,z=(z^{(1)},z{(2)})^T∈Z,定义从原空间到新空间的映射$$z=Φ(x)=((x{(1)})2,(x{(2)})2)^T$$
    • 在新空间里用线性分类学习方法从训练数据中学习分类模型
  • 核函数

    • 设X是输入空间,设H是特征空间,如果存在一个从X到H的映射$$Φ(x):X->H$$,使得所有的x,zXx,z∈X,函数K(x,z)满足:$$K(x,z)=Φ(x)·Φ(z)$$则称K(x,z)为核函数,Φ(x)为映射函数
    • 对偶问题 的目标函数可以转换为:$$W(a)=\frac{1}{2}\sum_{i=1}N\sum_{j=1}Na_ia_jy_iy_jK(x_i,x_j)-\sum_{i=1}^Na_i$$
    • 分类决策函数转换为:$$f(x)=sign(\sum_{i=1}{N_S}a_i*y_iK(x_i,x)+b^*)$$
  • 正定核

    • 定义映射,构成向量空间S
      • 定义映射Φ(x)=K(,x)Φ(x)=K(·,x),根据映射定义线性组合$$f(·)=\sum_{i=1}^Na_iK(·,x)$$
    • 在S上定义内积,使其成为内积空间
      • 在S上定义一个运算*:$$f*g·\sum_{i=1}m\sum{j=1}la_iβ_jK(x_i,z_j) = f·g$$
    • 将内积空间S完备化为希尔伯特空间
      • 内积范数:f=ff||f||=\sqrt{f·f}
      • 一个内积空间,当作为一个赋范向量空间是完备的时候,就是希尔伯特空间
      • 这一希尔伯特空间H称为 再生核希尔伯特空间,满足:$$K(·,x)·f=f(x)$$以及$$K(·,x)·K(·,z)=K(x,z)$$称为 再生核
    • 正定核的充要条件
      • 设K:X×X->R是对称函数,则对任意x∈X, K(x,z)对应的Gram矩阵K=[K(xi,xj)]m×mK=[K(x_i,x_j)]_{m×m}半正定矩阵
    • 正定核的等价定义
      • XRnX \subset R^n,K(x,z)是定义在X×X上的对称函数,如果对任意xiXx_i∈X,K(x,z)对应的Gram矩阵=[K(xi,xj)]m×m=[K(x_i,x_j)]_{m×m}是半正定矩阵,则K(x,z)是正定核
  • 常用核函数

    • 多项式核函数$$K(x,z)=(x·z+1)p$$,此时,分类决策函数为:$$f(x)=sign(\sum_{i=1}{N_s}a_i*y_i(x_i·x+1)p+b^*$$
    • 高斯核函数$$K(x,z)=exp(-\frac{||x-z||2}{2δ2})$$,分类决策函数同理
    • 字符串核函数(P139)
      • 核函数不仅可以定义在欧式空间上,还可以定义在离散数据集合上。比如,字符串核是定义在字符串集合上的核函数。
      • 字符串s长度用|s|表示,其元素为s(1)s(2)…s(|s|)。两个字符串s和t的连接记作st,所有字符串的集合记作\sum^\*=\bigcup\limits_{n=0}^∞\sum^\*n
      • u是s的子串,映射在u维上的取值为$$[Φ_n(s)]u=\sum{i:s(i)=u}λ^{l(i)}$$,这里0<λ≤1是一个衰减参数,l(i)表示字符串i的长度,求和在s中所有与u相同的子串中进行
      • 内积:$$k_n(s,t)=\sum_{u∈\sumn}[Φ_n(s)]_u[Φ_n(t)]_u=\sum_{u∈\sumn(i,j)}\sum_{s(i)=t(j)=u}λ{l(i)}λ{l(j)}$$
  • 非线性支持向量分类机

    • 非线性支持向量:$$f(x)=sign(\sum_{i=1}Na_i*y_iK(x,x_j)+b^*)$$,K(x,z)是正定核函数
    • 非线性支持向量机学习算法
      • 选取适当的核函数K(x,z)和适当参数C,构造并求解最优化问题$$\min\limits_a \frac{1}{2}\sum_{i=1}N\sum_{j=1}Na_ia_jy_iy_jK(x_i,x_j)-\sum_{i=1}^Na_i$$ $$s.t. \sum_{i=1}^Na_iy_j=0 0≤a_i≤C$$
      • 选择a^\*的一个正分量0<a_j^\*<C,计算$$b*=y_j-\sum_{i=1}Na_i^*y_iK(x_i·x_j)$$
      • 构造决策函数

序列最小化优化算法(SMO)

整个SMO算法包括两个部分:求解两个变量二次规划的解析方法和选择变量的启发式方法

  • 两个变量二次规划的求解方法(P141)
    • 假设选择的两个变量为a1,a2a_1,a_2,于是SMO的最优化问题的子问题为:$$\min\limits_{a_1,a_2} W(a_1,a_2)=\frac{1}{2}K_{11}a_12+\frac{1}{2}K_{22}a_22+y_1y_2K_{12}a_1a_2-(a_1+a_2)+y_1a_1\sum_{i=3}Ny_ia_iK_{i1}+y_2a_2\sum_{i=3}Ny_ia_iK_{i2}$$ $$s.t. a_1y_1+a_2y_2=-\sum_{i=3}^Ny_ia_i=ζ$$ $$0≤a_i≤C, i=1,2$$
    • 沿着约束方向未经剪辑时的解是$$a_2{new,unc}=a_2{old}+\frac{y_2(E_1-E_2)}{η}$$,其中$$η=K_{11}+K_{22}-2K_{12}=||Φ(x_1)-Φ(x_2)||^2$$
    • 经剪辑后a2的解为$$a_2^{new} = \begin{cases} H, a_2^{new,unc}>H \\a_2{new,unc},L≤a_2{new,unc}≤H \\L,a_2^{new,unc}<L \end{cases}$$
    • a1new=a1old+y1y2(a2olda2new)a_1^{new}=a_1^{old}+y_1y_2(a_2^{old}-a_2^{new})
  • 变量选择的方法(P144)
    • 第一个变量的选择
      • 称为外层循环
      • 检查样本点是否满足KKT条件:$$a_i=0 <=> y_ig(x_i)≥1$$ $$0<a_i<C <=> y_ig(x_i)=1$$ $$a_i=C <=> y_ig(x_i)≤1$$,其中,g(xi)=j=1NajyjK(xi,xj)+bg(x_i) = \sum_{j=1}^Na_jy_jK(x_i,x_j)+b
      • 选取违反KKT条件最严重的样本点作为第一个变量
      • 首先便利间隔边界上的样本点,检查他们是否满足KKT条件。如果全部满足则遍历整个训练样本集
    • 第二个变量的选择
      • 称为内层循环
      • 选择a2a_2使得对应的E1E2|E_1-E_2|最大
      • 如果选不到合适的a2a_2,则放弃第一个变量a1a_1,在通过外层循环寻找另外的变量
    • 计算阈值b和差值EiE_i
  • SMO算法
    1. 取初值a(0)=0,k=0a^{(0)}=0,令k=0
    2. 选取优化变量a1(k)a2(k)a_1^{(k)},a_2^{(k)},解析两个变量的最优化问题,求得最优解a1(k+1)a2(k+1)a_1^{(k+1)},a_2^{(k+1)}
    3. 若在精度ε单位内满足停机条件$$\sum_{i=1}^Na_iy_i=0$$ $$0≤a_i≤C,i=1,2,…,N$$ $$y_i·g(x_i)={≥1, {x_i|a_i=0} =1, {x_i|0<a_i<C} ≤1, {x_i|a_i=C}}$$则转4,否则令k=k+1,转2
    4. a^=a(k+1)\hat{a}=a^{(k+1)}

统计学习方法-第八章:提升方法

提升方法是一种常用的统计学习方法。在分类问题中,它通过改变训练样本的权重,学习多个分类器,并将这些分类器进行线性组合,提高分类性能

提升方法AdaBoost方法

  • AdaBoost算法
    1. 初始化训练数据的权值分布$$D_1=(w_{l1},…,w_{li},…,w_{lN}),w_{li}=\frac{1}{N},i=1,2,…,N$$,假设训练数据集具有均匀的权值分布。
    2. 对m=1,2,…,M
      • 使用具有权值分布的DmD_m的训练数据集学习,得到基本分类器$$G_m(x):X \rightarrow {-1,+1}$$
      • 计算Gm(x)G_m(x)在训练数据集上的分类误差率$$e_m=P(G_m(x_i)≠y_i)=\sum_{i=1}^Nw_{mi}I(G_m(x_i)≠y_i)$$
      • 计算Gm(x)G_m(x)的系数$$a_m=\frac{1}{2}\ln{\frac{1-e_m}{e_m}}$$
      • 更新训练集的权值分布$$D_{m+1}=(w_{m+1,1},…,w_{m+1,i},…,w_{m+1,N})$$ $$w_{m+1,i}=\frac{w_{mi}}{Z_m}exp(-a_my_iG_m(x_i))$$,其中Z是规范因子Zm=i=1Nwmiexp(amyiGm(xi))Z_m=\sum_{i=1}^Nw_{mi}exp(-a_my_iG_m(x_i))
    3. 构建基本分类器的线性组合$$f(x)=\sum_{m=1}Ma_mG_m(x)$$,得到最终分类器$$G(x)=sign(\sum_{m=1}Ma_mG_m(x))$$

AdaBoost算法实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
# 代码关于G_m选择仍有问题,有待修改
import numpy as np
def getDataset():
# 数据来源习题8.1
dataset = [[0,1,3,-1],[0,3,1,-1],[1,2,2,-1],[1,1,3,-1],[1,2,3,-1],[0,1,2,-1],[1,1,2,1],[1,1,1,1],[1,3,1,-1],[0,2,1,-1]]
name = ['身体','业务','潜力']
return dataset,name

def InitWeight(num):
return [1.0/num for i in range(num)]
def comperror(category, r, w):
error = 0
for i in range(len(category)):
if category[i]-r[i] != 0:
error += w[i]
return error
def getCoefficient(error):
a = 0.5 * np.log((1.0-error)/error)
return a
def G(w,inX):
if inX[2] == 2:
return 1
else:
return -1
def AdaBoost(dataset):
num = len(dataset)
feature = [x[:-1] for x in dataset]
category = [x[-1] for x in dataset]
w = InitWeight(num)
# 使用决策树作为初始分类器,以“潜力==2”为条件作为分类
error = 0.3
a_m = []
n = 10
while n:
a = getCoefficient(error)
a_m.append(a)
# 更新权重
Z = np.sum([w[i] * np.exp(-a*category[i]*G(w[i],feature[i])) for i in range(len(w))])
w = [w[i]*1.0/Z * np.exp(-a*category[i]*G(w[i],feature[i])) for i in range(len(w))]
r = [G(w[i],feature[i]) for i in range(len(feature))]
error = comperror(category,r,w)
print (error)
n -= 1
print (a_m)

if __name__ == '__main__':
dataset,name = getDataset()
AdaBoost(dataset)

AdaBoost算法的训练误差分析

  • AdaBoost的训练误差界
    • AdaBoost算法最终分类器的训练误差界为$$\frac{1}{N}\sum_{i=1}^NI(G(x_i)≠y_i)≤\frac{1}{N}\sum_{i}exp(-y_if(x_i))=\prod_{m}Z_m$$
  • 二类分类问题AdaBoost的训练误差界
    • m=1MZm=m=1M[2em(1em)]=m=1M(14γm2)exp(2m=1Mγm2)\prod_{m=1}^MZ_m=\prod_{m=1}^M[2\sqrt{e_m(1-e_m)}]=\prod_{m=1}^M\sqrt{(1-4γ_m^2)≤exp(-2\sum_{m=1}^Mγ_m^2)},其中γm=12emγ_m=\frac{1}{2}-e_m
    • 如果存在γ>0,对所有m有γmγγ_m≥γ,则\frac{1}{N}\sum_{i=1}NI(G(x_i)≠y_i)≤exp(-2Mγ2),表明在此条件下AdaBoost的训练误差是指数速率下降的

AdaBoost算法的解释

可以认为AdaBoost算法是模型为 加法模型,损失函数为 指数函数,学习

  • 前向分布算法
    • 加法模型:f(x)=m=1Mβmb(x;γm)f(x)=\sum_{m=1}^Mβ_mb(x;γ_m),其中,b(x;γm)b(x;γ_m)为基函数,γmγ_m为基函数参数,βmβ_m为基函数系数。
    • 在给定训练数据及损失函数L(y,f(x))的条件下,学习加法模型成为经验风险极小化,即损失函数极小化问题:$$\min\limits_{β_m,γ_m}\sum_{i=1}^N L(y_i,\sum_{m=1}^Mβ_mb(x_i;γ_m))$$
    • 前向分布算法
      • 初始化f0(x)=0f_0(x)=0
      • 对m=1,2,…,M
        • 极小化损失函数$$(β_m,γ_m)=\arg\min\limits_{β,γ}\sum_{i=1}^NL(y_i,f_{m-1}(x_i)+βb(x_i;γ))$$,得到参数βm,γmβ_m,γ_m
        • 更新$$f_m(x)=f_{m-1}(x)+β_mb(x;γ_m)$$
        • 得到加法模型$$f(x)=f_M(x)=\sum_{m=1}^Mβ_mb(c;γ_m)$$
  • 前向分布算法与AdaBoost
    • AdaBoost是前向分布加法算法的特例,此时模型是由基本分类器组成的加法模型,损失函数是指数函数

提升树

提升树是以分类树和或回归树为基本分类器的提升方法,被认为是统计学习中最好的方法之一

  • 提升树模型
    • 提升方法实际采用加法模型与前向分布算法,以 决策树 为基函数的提升方法称为提升树。对分类问题决策树是二叉分类树,对回归问题决策树是二叉回归树
    • 提升树模型可以表示为决策树的加法模型:$$f_M(x)=\sum_{m=1}^MT(x;Θ_m)$$,其中T(x;Θm)T(x;Θ_m)表示决策树,ΘmΘ_m为决策树参数,M为树的个数
  • 提升树算法
    • 初始化f0(x)=0f_0(x)=0
    • 对m=1,2,…,M
      • 计算残差$$r_{mi}=y_i-f_{m-1}(x_i),i=1,2,…,N$$
      • 拟合残差rmir_mi学习一个回归树,得到T(x;Θm)T(x;Θ_m)
      • 更新fm(x)=fm=1(x)+T(x;Θm)f_m(x)=f_{m=1}(x)+T(x;Θ_m)
    • 得到回归问题提升树:$$f_M(x)=\sum_{m=1}^MT(x;Θ_m)$$
  • 梯度提升
    • 初始化f0(x)=argminci=1NL(yi,c)f_0(x)=\arg\min\limits_{c}\sum_{i=1}^NL(y_i,c),此时是一个只有跟结点的树
    • 对m=1,2,…,N
      • 对i=1,2,…,N,计算:$$r_{mi}=-[\frac{\partial L(y_i,f(x_i))}{\partial f(x_i)}]_{f(x)=f_{m-1}(x)}$$,计算损失函数的负梯度在当前模型的值,作为残差的估计。
      • rmir_{mi}拟合一个回归树,得到第m棵树的叶结点区域Rmy,j=1,2,...,JR_{my},j=1,2,...,J,估计回归树叶结点区域,以拟合残差的近似值
      • 对j=1,2,…,J,计算$$c_{mj}=\arg\min\limits_{c}\sum_{x_i∈R_{mj}}L(y_i,f_{m-l}(x_i)+c)$$,利用线性搜索估计叶结点区域的值,使损失函数极小化
      • 更新得到fm(x)=fm1(x)+j=1JcmjI(xRmj)f_m(x)=f_{m-1}(x)+\sum_{j=1}^Jc_{mj}I(x∈R_{mj}),更新回归树
    • 得到回归树$$\hat{f}(x)=f_M(x)=\sum_{m=1}M\sum_{j=1}Jc_{mj}I(x∈R_{mj})$$

统计学习方法-第九章:EM算法及其推广

EM算法是一种迭代算法,用于含有隐变量的概率模型参数的极大似然估计,或极大后验概率估计

EM算法迭代每次由2步组成:E步,求期望;M步,求极大,所以称为 期望极大算法

EM算法的引入

当概率模型的变量都是观测变量,那么给定数据,可以直接使用极大似然估计或贝叶斯估计模型参数。但当模型含有隐变量时,就不能简单地使用这些方法。EM算法就是含有隐变量的概率模型参数的极大似然估计法,或极大后验概率估计法。

  • EM算法

    1. 选择参数的初值θ(0)θ^{(0)},开始迭代
    2. E步:记θ(i)θ^{(i)}为第i此迭代参数θ的估计值,在第i+1此迭代的E步,计算:$$Q(θ,θ{(i)})=E_Z[\log{P(Y,Z|θ)}|Y,θ{(i)}]=\sum_{Z}\log{P(Y,Z|θ)}P(Z|Y,θ{(i)})$$,其中P(Z|Y,θ{(i)})是在给定观测数据Y和当前的参数估计θ(i)θ^{(i)}下隐变量数据Z的条件概率分布。Q(θ,θ(i))Q(θ,θ^{(i)})中的第1个变元表示要极大化的参数,第2个表示参数的当前估计值。每次迭代实际在求Q函数及其极大
    3. M步:求使Q(θ,θ(i))Q(θ,θ^{(i)})极大化的θ,确定第i+1次迭代的参数估计值θ(i+1)θ^{(i+1)}:$$θ{(i+1)}=\arg\max\limits_{θ}Q(θ,θ{(i)})$$
    4. 重复第2和第3步,直到收敛。停止迭代的条件一般是对较小的正数ε1,ε2ε_1,ε_2,若满足$$||θ{(i+1)}-θ{(i)}||<ε_1 或 ||Q(θ{(i+1)},θ{(i)})-Q(θ{(i)},θ{(i)})||<ε_2$$,则停止迭代
  • Q函数

    • 完全数据对数似然函数 logP(Y,Zθ)\log{P(Y,Z|θ)}关于在给定观测数据Y和当前参数θ(i)θ^{(i)}在对未观测数据Z的条件概率分布P(ZY,θ(i))P(Z|Y,θ^{(i)})的期望称为 Q函数,即:$$Q(θ,θ{(i)})=E_Z[\log{P(Y,Z|θ)}|Y,θ{(i)}]$$
  • EM算法的导出(P174)

  • EM算法在非监督学习中的应用

    • EM算法可以用于生成模型的非监督学习。生成模型由联合概率分布P(X,Y)表示,可以认为非监督学习训练数据是联合概率分布产生的数据。X为观测数据,Y为未观测数据

EM算法的收敛性

  • P(Yθ)P(Y|θ)为观测数据的似然函数,θ(i)(i=1,2,...)θ^{(i)}(i=1,2,...)为EM算法得到参数估计序列,P(Yθ(i))P(Y|θ^{(i)})为对应的似然函数序列,则P(Yθ(i))P(Y|θ^{(i)})是单调递增的,即:$$P(Y|θ{(i+1)})≥P(Y|θ{(i)})$$
  • L(θ)=logP(Yθ)L(θ)=\log{P(Y|θ)}为观测数据的对数似然函数,θ(i)(i=1,2,..)θ^{(i)}(i=1,2,..)为EM算法得到的参数估计序列,L(θ(i))L(θ^{(i)})为对应的似然函数序列,
    1. 如果P(Y|θ)有上界,则L(θ(i))=logP(Yθ(i))L(θ^{(i)})=\log{P(Y|θ^{(i)})}收敛到某一值L^*
    2. 在函数Q(θ,θ)Q(θ,θ')与L(θ)满足一定条件下,由EM算法得到的参数估计序列θ(i)θ^{(i)}是L(θ)的稳定点

EM算法在高斯混合模型学习中的应用

  • 高斯混合模型
    • 高斯混合模型是指具有如下形式的概率分布模型:$$P(y|θ)=\sum_{k=1}Ka_kΦ(y|θ_k)$$,其中,$a_k$是系数,$a_k≥0,\sum_{k=1}^Ka_k=1$,$Φ(y|θ_k)$是高斯分布密度,$θ_k=(μ_k,σ_k^2)$,$$Φ(y|θ_k)=\frac{1}{\sqrt{2π}σ_k}exp(-\frac{(y-μ_k)2}{2σ_k^2})$$,称为 第k个分模型
  • 高斯混合模型参数估计的EM算法(P179)
    1. 取参数的初始值开始迭代
    2. EM算法的E步,确定Q函数:依据当前模型,计算分模型k对观测数据yiy_i的响应度:$$\hat{y}_{jk}=\frac{a_kΦ(y_j|θ_k)}{\sum_{k=1}^Ka_kΦ(y_j|θ_k)},j=1,2,…,N, K=1,2,…,N$$
    3. 确定EM算法的M步:计算新一轮迭代的模型参数:$$\hat{μ}_k=\frac{\sum_{j=1}N\hat{γ}_{jk}y_j}{\sum_{j=1}N\hat{γ}_{jk}}$$ $$\hat{σ}_k2=\frac{\sum_{j=1}N\hat{γ}_{jk}(y_j-μ_k)2}{\sum_{j=1}N\hat{γ}_{jk}}$$ $$\hat{a}_k=\frac{\sum_{j=1}^N\hat{γ}_{jk}}{N}$$
    4. 重复第2,第3步,直到收敛

EM算法的推广

EM算法还可以解释为F函数的 极大-极大算法

  • F函数
    • 假设隐变量数据Z的概率分布为P~(Z)\tilde{P}(Z),定义分布P~\tilde{P}与参数θ的函数F(P~,θ)F(\tilde{P},θ)如下:$$F(\tilde{P},θ)=E_{\tilde{P}}[\log{P(Y,Z|θ)}]+H(\tilde{P})$$,其中,H(P~)=EP~logP~(Z)H(\tilde{P})=-E_{\tilde{P}}\log{\tilde{P}(Z)}
    • 对于固定的θ,存在唯一的分布P~_θ\tilde{P}\_θ极大化F(P~,θ)F(\tilde{P},θ),此时:$$\tilde{P}_θ(Z)=P(Z|Y,θ)$$,且P~_θ\tilde{P}\_θ随θ连续变化
    • P~_θ(Z)=P(ZY,θ)\tilde{P}\_θ(Z)=P(Z|Y,θ),则$$F(\tilde{P},θ)=\log{P(Y|θ)}$$
    • L(θ)=logP(Yθ)L(θ)=\log{P(Y|θ)}为观测数据的对数似然函数,θ(i),i=1,2,...θ^{(i)},i=1,2,...为EM算法得到的参数估计序列。如果函数F(P~,θ)F(\tilde{P},θ)\tilde{P}^\*θ^\*有局部(全局)最大值,则L(θ)在θ^\*也有局部(全局)最大值;
    • EM算法的一次迭代可由F函数的极大-极大算法实现
      • θ(i)θ^{(i)}为第i次迭代参数θ的估计,P~(i)\tilde{P}^{(i)}为第i次迭代函数的P~\tilde{P}的估计,在第i+1次迭代的两步为:
      1. 对固定的θiθ^{i},求P~(i+l)\tilde{P}^{(i+l)}使F(P~,θ(i))F(\tilde{P},θ^{(i)})极大化
      2. 对固定的P~(i+l)\tilde{P}^{(i+l)},求θi+1θ^{i+1}使F(P~(i+1),θ)F(\tilde{P}^{(i+1)},θ)极大化
  • GEM算法-1
    1. 初始化参数θ(0)θ^{(0)},开始迭代
    2. 第i+1次迭代,第1步:记θ(i)θ^{(i)}为参数θ的估计值,P~(i)\tilde{P}^{(i)}为函数P~\tilde{P}的估计,求P~(i+1)\tilde{P}^{(i+1)}使 P~\tilde{P}极大化的F(P~,θ(i))F(\tilde{P},θ^{(i)})
    3. 第2步:求θ(i+1)θ^{(i+1)}使F(P~(i+1),θF(\tilde{P}^{(i+1)},θ极大化
    4. 重复2,3步,直到收敛
  • GEM算法-2
    1. 初始化参数θ(0)θ^{(0)},开始迭代
    2. 第i+1次迭代,第1步: 记θ(i)θ^{(i)}为参数θ的估计值,计算$$Q(θ,θ{(i)})=E_Z[\log{P(Y,Z|θ)}|Y,θ{(i)}]=\sum_{Z}P(Z|Y,θ^{(i)})\log{P(Y,Z|θ)}$$
    3. 第2步:求θ(i+1)θ^{(i+1)}使$$Q(θ{(i+1)},θ{(i)})>Q(θ{(i)},θ{(i)})$$
    4. 重复2,3步,直到收敛
  • GEM算法-3
    1. 初始化参数θ(0)θ^{(0)},开始迭代
    2. 第i+1次迭代,第1步,记θ(i)θ^{(i)}为参数θ的估计值,计算$$Q(θ,θ{(i)})=E_Z[\log{P(Y,Z|θ)}|Y,θ{(i)}]=\sum_{Z}P(Z|Y,θ^{(i)})\log{P(Y,Z|θ)}$$
    3. 第2步,进行d次条件极大化:首先,在θ2(i),...,θk(i)θ_2^{(i)},...,θ_k^{(i)}保持不变的条件下求使Q(θ,θ(i))Q(θ,θ^{(i)})达到极大的θ1(i+1)θ_1^{(i+1)};然后在θ1=θ1(i+1)θ_1=θ_1^{(i+1)},θj=θj(i),j=3,4,...,kθ_j=θ_j^{(i)},j=3,4,...,k的条件下达到极大的θ2(i+1)θ_2^{(i+1)}。经过d次条件最大化,得到θ(i+1)=(θii+1,θ2i+1,...,θdi+1)θ^{(i+1)}=(θ_i^{i+1},θ_2^{i+1},...,θ_d^{i+1})使得$$Q(θ{(i+1)},θ{(i)})>Q(θ{(i)},θ{(i)})$$
    4. 重复2,3步,直到收敛

统计学习方法-第十章:隐马尔科夫模型

隐马尔科夫模型是可用于标注问题的统计学习模型,属于生成模型

隐马尔科夫模型的基本概念

  • 隐马尔科夫模型的定义

    • 隐马尔科夫模型是 关于时序的概率模型,描述由一个隐藏的马尔科夫链随机生成不可观测的状态随机序列,再由各个状态生成一个观测而产生观测随机序列的过程
    • 隐藏的马尔科夫链随机生成的状态的序列,称为 序列状态
    • 每个状态生成一个观测,而由此产生的观测的随机序列,称为 观测序列
    • 序列的每一个位置又可以看做一个时刻
    • 隐马尔科夫模型由 初始概率分布状态转移分布、以及观测概率分布确定
  • 定义

    • 设Q是所有可能的状态的集合,V是所有可能的观测的结合:$$Q=\{q_1,q_2,…,q_N\}$$ $$V=\{v_1,v_2,…,v_M\}$$,其中,N是可能的状态数,M是可能的观测数
    • I是长度为T的状态序列,O是对应的观测序列:$$I=(i_1,i_2,…i_T)$$ $$O=(o_1,o_2,…,o_T)$$
    • A是状态转移概率矩阵:$$A=[a_{ij}]_{N×N}$$,其中$$a_{ij}=P(i_{t+1}=q_j|i_t=q_i),i=1,2,…,N;j=1,2,…,N$$表示时刻t处于状态qtq_t的条件在t+1时刻转移到qjq_j的概率
    • B是观测概率矩阵:$$B=[b_j(k)]_{N×M}$$,其中$$b_j(k)=P(o_t=v_k|i_t=q_j),k=1,2,…,M;j=1,2,…,N$$,表示时刻t处在状态qjq_j的条件下生成观测vkv_k的概率
    • π是初始状态概率向量:$$π=(π_i)$$其中,$$π_i=p(i_1=q_i),i=1,2,…,N$$是时刻t=1处于状态q_1的概率
    • 隐马尔科夫模型由π和A确定了隐藏的马尔科夫链,生成不可观测的状态序列;B确定了如何从状态生成观测,与状态序列综合确定了如何产生观测序列。因此,隐马尔科夫模型λ可以用三元符号表示,即$$λ=(A,B,π)$$
  • 马尔科夫模型的2个基本假设

    1. 齐次马尔 科夫性假设。即假设隐藏的马尔科夫链在任意时刻t的状态只 依赖于其前一时刻的状态,与其他时刻的状态及观测无关,也与时刻t无关
    2. 观测独立性假设。即假设任意时刻的观测只依赖于该时刻的马尔科夫链的状态,与其他观测及状态无关
  • 观测序列的生成过程

    1. 按照初始状态分布π产生状态i1i_1
    2. 令t=1
    3. 按照状态iti_t的观测概率分布bit(k)b_{i_t}(k)生成oto_t
    4. 按照状态iti_t的状态转移概率分布aitit+1\\{a_{i_ti_{t+1}}\\}产生状态it+1,it+1=1,2,...,Ni_{t+1},i_{t+1}=1,2,...,N
    5. 令t=t+1;如果t<T,转第3步;否则,终止
  • 隐马尔科夫模型的3个基本问题

    1. 概率计算问题。给定模型λ=(A,B,π)和观测序列O=(o1,o2,...,oT)O=(o_1,o_2,...,o_T),计算在模型λ下观测序列O出现的概率P(O|λ)
    2. 学习问题。已知观测序列O=(o_1,o_2,…,o_T),估计模型λ=(A,B,π)参数,使得在该模型下观测序列概率P(O|λ)最大。即用极大似然估计的方法估计参数
    3. 预测问题。也称为解码问题。已知模型λ=(A,B,π)和观测序列O=(o_1,o_2,…,o_T),求对给定观测序列条件概率P(I|O)最大的状态序列I=(i1,i2,...,iT)I=(i_1,i_2,...,i_T)。即,给定观测序列,求最有可能的对应的状态序列

概率计算算法

  • 直接计算法

    • 给定模型λ=(A,B,π)和观测序列O=(o1,o2,..,oT)O=(o_1,o_2,..,o_T),计算观测序列O出现的概率P(O|λ)
    • 状态序列I=(i1,i2,...,iT)I=(i_1,i_2,...,i_T)的概率是$$P(I|λ)=π_{i_t}a_{i_1i_2}a_{i_2i_3}…a_{i_{T-1}i_T}$$
    • 对固定的状态序列I=(i1,i2,...,iT)I=(i_1,i_2,...,i_T),观测序列O=(o1,o2,...,oT)O=(o_1,o_2,...,o_T)的概率是P(O|I,λ),$$P(O|I,λ)=b_{i_1}(o_1)b_{i_2}(o_2)…b_{i_T}(o_T)$$
    • O和I同时出现的联合概率为$$P(O,I|γ)=P(O|I,λ)P(I|λ)=π_{i_1}b_{i_1}(o_1)a_{i_1i_2}b_{i_2}(o_2)…a_{i_{T-1}i_T}b_{i_T}(o_T)$$
  • 前向算法

    • 前向概率
      • 给定马尔科夫模型λ,定义到时刻t部分观测序列为o1,o2,...,oto_1,o_2,...,o_t且状态为qiq_i的概率为前向概率,记作$$a_t(i)=P(o_1,o_2,…,o_t,i_t=q_i|λ)$$
    • 前向算法
      1. 初值$$a_1(i)=π_ib_i(o_1),i=1,2,…,N$$,初始化前向概率
      2. 递推 对t=1,2,…,T-1 $$a_{t+1}(i)=[\sum_{j=1}^Na_t(j)a_{ji}]b_i(o_{t+1}),i=1,2,…,N$$
      3. 终止$$P(O|λ)=\sum_{i=1}^Na_T(i)$$
  • 后向算法

    • 后向概率
      • 给定隐马尔科夫模型λ,定义在时刻t状态为qiq_i的条件下,从t+1到T部分观测序列为ot+1,ot+2,...,oTo_{t+1},o_{t+2},...,o_{T}的概率为后向概率,记作$$β_t(i)=P(o_{t+1},o_{t+2},…,o_T|i_t=q_i,λ)$$
    • 后向算法
      1. 初始化后向概率,对最终时刻的所有状态qiq_i规定βT(i)=1β_T(i)=1
      2. 对t=T-1,T-2,…,1 $$β_t(i)=\sum_{j=1}^Na_{ij}b_j(o_{t+1})β_{t+1}(j),i=1,2,…,N$$
      3. P(Oλ)=i=1Nπibi(o1)β1(i)P(O|λ)=\sum_{i=1}^Nπ_ib_i(o_1)β_1(i)
  • 一些概率与期望值的计算

    1. 给定模型λ和观测O,在t时刻状态qiq_i的概率:$$γ_t(i)=P(i_t=q_i|O,λ)=\frac{P(i_t=q_i,O|λ)}{P(0|λ)}$$。由前向概率at(i)a_t(i)和后向概率βt(i)β_t(i)可得:γt(i)=at(i)βt(i)j=1Nat(j)βt(j)γ_t(i)=\frac{a_t(i)β_t(i)}{\sum_{j=1}^Na_t(j)β_t(j)}$

    2. 给定模型λ和观测O,在时刻t处于状态qiq_i且在时刻t+1处于状态qjq_j的概率:$$ζ_t(i,j)=\frac{a_t(i)a_{ij}b_j(o_{t+1})β_{t+1}(j)}{\sum_{i=1}N\sum_{j=1}Na_t(i)a_{ij}b_j(o_{t+1})β_{t+1}(j)}$$

    3. γt(i)γ_t(i)ζt(i,j)ζ_t(i,j)对各个时刻t求和,得到:

    4. 在观测O下状态i出现的期望值$$\sum_{t=1}^Tγ_t(i)$$

    5. 在观测O下由状态i转移的期望值$$\sum_{t=1}^{T-1}γ_t(i)$$

    6. 在观测O下由状态i转移到状态j的期望值$$\sum_{t=1}^{T-1}ζ_t(i,j)$$

学习算法

  • 监督学习方法

    • 转移概率aija_{ij}的估计
      • 设样本中时刻t处于状态i,时刻t+1转移到状态j的频数为AijA_{ij},那么状态转移概率aija_{ij}的估计是$$\hat{a}_{ij}=\frac{A_{ij}}{\sum_{i=1}^NA_{ij}}$$
    • 观测概率bj(k)b_j(k)的估计
      • 设样本中状态为j并观测为k的频数是BjkB_{jk},那么状态为j观测为k的概率bj(k)b_j(k)的估计是$$\hat{b}_j(k)=\frac{B_{jk}}{\sum_{k=1}^MB_{jk}}$$
    • 初始状态概率πiπ_i的估计π^_i\hat{π}\_i为S个样本中初始状态为qiq_i的频率
  • Baum-Welch算法

    • 对于隐马尔科夫模型,我们将观测序列数据看做观测数据O,状态序列数据看做不可观测的隐数据I,那么隐马尔科夫模型事实上是一个含有隐变量的概率模型$$P(O|λ)=\sum_{I}P(O|I,λ)P(I|λ)$$,它的参数学习可以由EM算法实现
    1. 确定完全数据的对数似然函数:logP(O,Iλ)\log{P(O,I|λ)}
    2. EM算法E步:求Q函数Q(λ,λˉ)Q(λ,\bar{λ}):$$Q(λ,\bar{λ})=\sum_{I}\log{P(O,I|λ)}P(O,I|\bar{λ})$$,其中λˉ\bar{λ}是隐马尔科夫模型参数的当前估计值,λ是要极大化的隐马尔科夫模型参数
    3. EM算法M步:极大化Q函数求模型参数A,B,π $$π_i=\frac{P(O,i_1=i|\bar{λ})}{P(O|\bar{λ})}$$ $$a_{ij}=\frac{\sum_{t=1}{T-1}P(O,i_t=i,i_{t+1}=j|\bar{λ})}{\sum_{t=1}{T-1}P(O,i_t=i|\bar{λ})}$$ $$b_j(k)=\frac{\sum_{t=1}^TP(O,i_t=j|\bar{λ})I(o_t=v_k)}{\sum_{t=1}P(O,i_t=j|\bar{λ})}$$
  • Baum-Welch算法

    1. 初始化,对n=0,选取aij(0),bj(k),πi(0)a_{ij}^{(0)},b_j(k),π_i^{(0)},得到模型λ(0)=(A(0),B(0),π(0))λ^{(0)}=(A^{(0)},B^{(0)},π^{(0)})
    2. 递推,对n=1,2,…$$a_{ij}=\frac{\sum_{t=1}{T-1}ζ_t(i,j)}{\sum_{t=1}{T-1}γ_t(i)}$$ $$b_j(k){(n+1)}=\frac{\sum_{t=1,o_t=v_k}Tγ_t(j)}{\sum_{t=1}^Tγ_t(j)}$$ $$π_i^{(n+1)}=γ_1(j)$$
    3. 终止,得到模型参数λ(n+1)λ^{(n+1)}

预测算法

  • 近似算法
    • 近似算法:在每个时刻t选择在该时刻最有可能出现的状态i_t^\*,从而得到一个状态序列I^\*=(i_1^\*,i_2^\*,...,i_T^\*),将它作为预测结果
    • 给定隐马尔科夫模型λ和观测序列O,在时刻t处于状态qiq_i的概率γt(i)γ_t(i)是$$γ_t(i)=\frac{a_t{(i)}β_t{(i)}}{\sum_{j=1}Na_t(j)β_t(j)}$$,在每一个时刻t最有可能的状态i_t^\*=\arg\max\limits_{1≤i≤N}[γ_t(i)],t=1,2,...,T
  • 维特比算法
    • 用动态规划求概率最大路径,此时每条路径对应一个状态序列
    1. 初始化$$δ_1(i)=π_ib_i(o_1),i=1,2,…,N$$ $$Ψ_1(i)=0,i=1,2,…,N$$
    2. 递推。对t=2,3,…,T $$δ_1(i)=\max\limits_{1≤j≤N}[δ_{t-1}(j)a_{ji}]b_i(o_t),i=1,2,…,N$$ $$Ψ_t(i)=\arg\max\limits_{1≤j≤N}[δ_{t-1}(j)a_{ji}],i=1,2,…,N$$
    3. 终止$$P^*=\max\limits_{1≤j≤N}δ_T(i)$$ $$i_T^*=\arg\max\limits_{1≤j≤N}[δ_T(i)]$$
    4. 最优路径回溯

统计学习方法-第十一章:条件随机场

条件随机场是给定一组输入随机变量条件下,另一组输出随机变量的条件概率分布模型。

特点:假设输出随机变量构成马尔科夫随机场

概率无向图模型

概率无向图模型又称为马尔科夫随机场,是一个可以由无向图表示的联合概率分布

  • 模型定义
    • 记图G=(V,E),V和E分别是点和边的集合
    • 设有联合概率分布P(Y),Y∈γ是一组随机变量。由无向图G表示概率分布P(Y),即在图G中,v∈V表示一个随机变量Yv,Y=(Yv)_vVY_v,Y=(Y_v)\_{v∈V};边e∈E表示随机变量之间的概率依赖关系
    • 成对马尔科夫性
      • 设u和v是无向图G中任意两个没有边连接的结点,结点u和v分别对应随机变量YuYvY_u和Y_v。其他所有结点为O,对应的随机变量组是YOY_O。成对马尔科夫性是指,给定随机变量组YOY_O的条件下随机变量Yu,YvY_u,Y_v是条件独立的,即$$P(Y_u,Y_v|Y_o)=P(Y_u|Y_O)P(Y_v|Y_O)$$
    • 局部马尔科夫性
      • 设v∈V是无向图G中任意一个结点,W是与v有边连接的所有结点,O是v,W以外的其他所有结点。局部马尔科夫性是指在给定随机变量组YWY_W的条件下随机变量YvY_v和随机变量YOY_O是独立的,即$$P(Y_v,Y_O|Y_W)=P(Y_v|Y_W)P(Y_O|Y_W)$$
      • P(YOYW)>0P(Y_O|Y_W)>0时,有$$P(Y_v|Y_W)=P(Y_v|Y_W,Y_O)$$
    • 全局马尔科夫性
      • 设结点集合A,B是在无向图G中被结点集合C分开的任意结点集合。全局马尔科夫性是指给定随机变量组YCY_C条件下随机变量组YAY_AYBY_B是条件独立的,即$$P(Y_A,Y_B|Y_C)=P(Y_A|Y_C)P(Y_B|Y_C)$$
    • 概率无向图
      • 设有联合概率分布P(Y),由无向图G=(V,E)表示。在G中,结点表示随机变量,边表示随机变量之间的依赖关系。如果联合概率分布P(Y) 满足成对、局部或全局马尔科夫性,就称次联合概率分布为概率无向图模型,或马尔科夫随机场
  • 概率无向图的因子分解
    • 团与最大团
      • 无向图G中任何两个结点均有边连接的结点子集称为
      • 若C是无向图G的一个团,并且不能再加进任何一个结点使其成为一个更大的团,则称C为最大团
    • 将概率无向图的联合概率分布表示为其最大团上的 随机变量的函数的乘积形式的操作,称为 概率无向图模型的因子分解,即$$P(Y)=\frac{1}{Z}\prod_{C}Ψ_C(Y_C)$$,其中Z是规范化因子$$Z=\sum_{Y}\prod_{C}Ψ_C(Y_C)$$
    • ΨC(YC)Ψ_C(Y_C)是势函数,ΨC(YC)=expE(YC)Ψ_C(Y_C)=exp{-E(Y_C)}
    • Hammersley-Clifford定理
      • 概率分布图模型的联合概率分布P(Y)可以表示为:$$P(Y)=\frac{1}{Z}\prod_CΨ_C(Y_C)$$ $$Z=\sum_{Y}\prod_{C}Ψ_C(Y_C)$$

条件随机场的定义与形式

  • 条件随机场的定义
    • 条件随机场是指给定随机变量X条件下,随机变量Y的马尔科夫随机场。
    • 设X与Y是随机变量,P(Y|X)是在给定X的条件下Y的条件概率分布。所随机变量Y构成一个由无向图G=(V,E)表示的马尔科夫随机场,即$$P(Y_v|X,Y_w,w≠v)=P(Y_v|X,Y_w,wv)$$对任意结点v成立,则称条件概率分布P(Y|X)为条件随机场。wv表示在G中与v有边连接的所有结点w
    • 线性链条件随机场
      • 是定义在线性链上的特殊的条件随机场,可以运用在标注等问题
      • 设X=(X_1,X_2,…,X_n),Y=(Y_1,Y_2,…,Y_n)均为线性链表示的随机变量序列,若在给定随机变量序列X的条件下,随机变量序列Y的条件概率分布P(Y|X)构成条件随机场,即满足马尔科夫性$$P(Y_i|X,Y_1,…,Y_{i-1},Y_{i+1},…,Y_n)=P(Y_i|X,Y_{i-1},Y_{i+1})$$ $$i=1,2,…,n(在i=1和n时只考虑单边)$$
  • 条件随机场的参数化形式
    • 线性链条件随机场的参数化形式
      • 设P(Y|X)为线性链条件随机场,则在随机变量X取值为x的条件下,随机变量Y取值为y的条件概率具有如下形式:$$P(y|x)=\frac{1}{Z(x)}exp(\sum_{i,k}λ_kt_k(y_{i-1},y_i,x,i)+\sum_{i,l}u_ls_l(y_i,x,i))$$,其中$$Z(x)=\sum_exp(\sum_{i,k}λ_kt_k(y_{i-1},y_i,x,y)+\sum_{i,l}u_ls_l(y_i,x,i))$$式中,tk,slt_k,s_l是特征函数,λk,ulλ_k,u_l是对应的权值。Z(x)是规范化因子。
      • tkt_k是定义在边上的特征函数,称为 转移特征,依赖于当前和前一个位置
      • sls_l是定义在结点上的特征函数,称为 状态特征,依赖于当前位置
      • tk,slt_k,s_l都依赖于位置,是局部特征函数,通常取值为1(满足条件时)和0
    • 条件随机场的简化形式
      • 设有K1K_1个转移特征,K2K_2个状态特征,K=K1+K2K=K_1+K_2,记$$f_k(y_{i-1},y_i,x,i)=\begin{cases}t_k(y_{i-1},y_i,x,i),k=1,2,…,K_1 \
        s_l(y_i,x,i),k=K_1+l;l=1,2,…,K_2 \end{cases}$$
      • 然后,对转移与状态特征在各个位置i求和,记作$$f_k(y,x)=\sum_{i=1}^n f_k(y_{i-1},y_i,x,i), k=1,2,…,K$$
      • wkw_k表示特征fk(y,x)f_k(y,x)的权值,即$$w_k=\begin{cases}λ_k,k=1,2,…,K_1 \
        u_l, k=K_1+l;l=1,2,…,K_2 \end{cases}$$
      • 条件随机场表示为$$P(y|x)=\frac{1}{Z(x)}exp\sum_{k=1}^Kw_kf_k(y,x)$$ $$Z(x)=\sum_{y}exp\sum_{k=l}^K w_kf_k(y,x)$$
      • w=(w1,w2,...,wK)Tw=(w_1,w_2,...,w_K)^T,F(y,x)=(f1(y,x),f2(y,x),...,fK(y,x))TF(y,x)=(f_1(y,x),f_2(y,x),...,f_K(y,x))^T,则条件随机场为$$P_w(y|x)=\frac{exp(w·F(y,x))}{Z_w(x)}$$,其中Zw(x)=yexp(wF(y,x))Z_w(x)=\sum_{y}exp(w·F(y,x))
    • 条件随机场的矩阵形式
      • 对观测序列x的每一个位置i=1,2,…,n+1,定义一个m阶矩阵(m是标记y_i取值的个数) $$M_i(x)=[M_i(y_{i-1},y_i|x)]$$ $$M_i(y_{i-1},y_i|x)=exp(W_i(y_{i-1},y_i|x))$$ $$W_i(y_{i-k},y_i|x)=\sum_{i=1}^Kw_kf_k(y_{i-1},y_i,x,i)$$
      • 条件概率为:$$P_w(y|x)=\frac{1}{Z_w(x)}\prod_{i=1}^{n+1}M_i(y_{i-1},y_i|x)$$,其中Zw(x)Z_w(x)为规范因子,$$Z_w(x)=(M_1(x)M_2(x)…M_{n+1}(x))_{start,stop}$$

条件随机场的概率计算问题

  • 前向-后向算法
    • 对每个指标i=0,1,…,n+1,定义前向向量ai(x)a_i(x):$$a_0(y|x) =\begin{cases} 1, y=start\\ 0,否则 \end{cases}$$
    • aiT(x)=ai1T(x)Mi(x)a_i^T(x)=a_{i-1}^T(x)M_i(x),ai(x)a_i(x)是m维列向量
    • 对每个指标i=0,1,…,n+1,定义后向向量βi(x)β_i(x)βn+1(yn+1x)={1,yn+1=stop 0,β_{n+1}(y_{n+1}|x) =\begin{cases} 1, y_{n+1}=stop\\\ 0,否则 \end{cases}$
    • βi(x)=Mi+1(x)βi+1(x)β_i(x)=M_{i+1}(x)β_{i+1}(x)
    • 从前向-后向定义可得$$Z(x)=a_nT(x)·1=1T·β_1(x)$$,此处的1是元素均为1的m维列向量
  • 概率计算
    • 标记序列在位置i是标记yiy_i的条件概率和在位置i-1与i是标记yi1yiy_{i-1}和y_i的条件概率:$$P(Y_i=y_i|x)=\frac{a_i^T(y_i|x)β_i(y_i|x)}{Z(x)}$$ $$P(Y_{i-1}=y_{i-1},Y_i=y_i|x)=\frac{a_{i-1}^T(y_{i-1}|x)M_i(y_{i-1},y_i|x)β_i(y_i|x)}{Z(x)}$$,其中Z(x)=anT(x)1Z(x)=a_n^T(x)·1
  • 期望值计算
    • 特征函数fkf_k关于条件分布P(Y|X)的数学期望是$$E_{P(Y|X)}=[f_k]=\sum_{y}P(y|x)f_k(y,x)=\sum_{i=1}{n+1}\sum_{y_{i-1}y_i}f_k(y_{i-1},y_i,x,i)\frac{a_{i-1}T(y_{i-1}|x)M_i(y_{i-1},y_i|x)β_i(y_i|x)}{Z(x)}$$其中,Z(x)=anT(x)1Z(x)=a_n^T(x)·1
    • 假设经验分布为P~(x)\tilde{P}(x),特征函数fkf_k关于联合分布P(X,Y)的数学期望是$$E_{P(X,Y)}[f_k]=\sum_{x}\tilde{P}(x)\sum_{i=1}{n+1}\sum_{y_{i-1},y_i}f_k(y_{i-1},y_i,x,i)\frac{a_{i-1}T(y_{i-1}|x)M_i(y_{i-1},y_i|x)β_i(y_i|x)}{Z(x)}$$

条件随机场的学习算法

  • 改进的迭代进度法
    • 训练数据的对数似然函数:$$L(w)=L_{\tilde{P}(P_w)}=\log{\prod_{x,y}P_w(y|x)^{\tilde{P}(x,y)}}=\sum_{x,y}\tilde{P}(x,y)$$
    • 假设模型的当前参数向量为w=(w1,w2,...,wk)Tw=(w_1,w_2,...,w_k)^T,向量的增量为δ=(δ1,δ2,...,δK)Tδ=(δ_1,δ_2,...,δ_K)^T,更新参数向量为w+δ=(w1+δ1,w2+δ2,...,wK+δK)Tw+δ=(w_1+δ_1,w_2+δ_2,...,w_K+δ_K)^T
    • 转移特征tkt_k的更新方程:$$E_{\tilde{P}[t_k]}=\sum_{x,y}\tilde{P}(x)P(y|x)\sum_{i=1}^{n+1}t_k(y_{i-1},y_i,x,i)exp(δ_kT(x,y)),k=1,2,…,K_1$$
    • 状态特征sls_l的更新方程:$$E_{\tilde{P}}[s_l]=\sum_{x,y}\tilde{P}(x)P(y|x)\sum_{i=1}^ns_l(y_i,x,i)exp(δ_{K_1+1}T(x,y)), l=1,2,…,K_2$$,此处T(x,y)是在数据(x,y)中出现所有特征数的总和:$$T(x,y)=\sum_{k}f_k(x,y)=\sum_{k=1}K\sum_{i=1}{n+1}f_k(y_{i-1},y_i,x,i)$$
  • 条件随机场模型学习的改进的迭代尺度法
    1. 对所有k∈{1,2,…,K},取初值wk=0w_k=0
    2. 对每一k∈{1,2,…,K}:
      • k=1,2,...,K1k=1,2,...,K_1,令δkδ_k是方程:$$\sum_{x,y}\tilde{P}(x)P(y|x)\sum_{i=1}{n+1}t_k(y_{i-1},y_i,x,i)exp(δ_kT(x,y))=E_{\tilde{P}[t_k]}$$的解;当$k=K_1+l,l=1,2,…,K_2$时,令$δ_{K_1+1}$是方程$$\sum_{x,y}\tilde{P}(x)P(y|x)\sum_{i=1}ns_l(y_i,x,i)exp(δ_{K_1+l}T(x,y))=E_{\tilde{P}}[s_l]$$的解
      • 更新wkw_k值,wkwk+δkw_k\rightarrow w_k+δ_k
    3. 如果不是所有w_k都收敛,重复步骤2
  • 松弛特征 (算法S)
    • s(x,y)=Si=1n+1k=1Kfk(yi1,yi,x,i)s(x,y)=S-\sum_{i=1}^{n+1}\sum_{k=1}^Kf_k(y_{i-1},y_i,x,i),其中S是一个常数。选择足够大的常数S使得对训练数据集的所有数据(x,y),s(x,y)≥0成立,此时特征总数可取S
    • 对于转移特征tk,δkt_k,δ_k的更新方程是:$$\sum_{x,y}\tilde{P}(x)P(y|x)\sum_{i=1}^{n+1}t_k(y_{i-1},y_i,x,i)exp(δ_kS)=E_{\tilde{P}}[t_k]$$ $$δ_k=\frac{1}{S}\log{\frac{E_{\tilde{P}[t_k]}}{E_P[t_k]}}$$,其中$$E_P(t_k)=\sum_{x}\tilde{P}(x)\sum_{i=1}{n+1}\sum_{y_{i-1},y_i}T_k(y_{i-1},y_i,x,i)\frac{a_{i-1}T(y_{i-1}|x)M_i(y_{i-1},y_i|x)β_i(y_i|x)}{Z(x)}$$
    • 对于状态特征sl,δks_l,δ_k的更新方程:$$\sum_{x,y}\tilde{P}(x)P(y|x)\sum_{i=1}^ns_l(y_i,x,i)exp(δ_{K_1+l}S)=E_{\tilde{P}}[s_l]$$ $$δ_{K_1+1}=\frac{1}{S}\log{\frac{E_{\tilde{P}[S_l]}}{E_P[S_l]}}$$,其中$$E_P(S_l)=\sum_{x}\tilde{P}(x)\sum_{i=1}n\sum_{y_i}S_l(y_i,x,i)\frac{a_iT(y_i|x)β_i(y_i|x)}{Z(x)}$$
    • 由于算法S中需要S取足够大,这样一来,每一步迭代的增量会变大,算法收敛慢。于是提出了算法T
  • 算法T
    • 对每一个观测序列x计算其特征总数最大值T(x):$$T(x)=\max\limits_{y}T(x,y)$$
    • 转移特征的更新方程:$$E_{\tilde{P}}[t_k]=\sum_{t=0}{T_{max}}a_{k,t}β_kt$$,这里ak,ta_{k,t}是特征tkt_k的期待值,δk=logβkδ_k=\log{β_k}
    • 状态特征的更新方程:$$E_{\tilde{P}}[S_l]=\sum_{t=0}{T_{max}}b_{l,t}γ_lt$$,这里bl,tb_{l,t}是特征SlS_l的期待值,δk=logγtδ_k=\log{γ_t}
  • 拟牛顿法
    1. 选定初始点w(0)w^{(0)},取B0B_0为正定对称矩阵,置k=0
    2. 计算gk=g(w(k))g_k=g(w^{(k)}).若gk=0g_k=0,则停止计算,否则转3
    3. Bkpk=gkB_kp_k=-g_k求出pkp_k
    4. 一维搜索,求λkλ_k使得$$f(w_{(k)}+λ_kp_k)=\min\limits_{λ≥0}f(w^{(k)}+λp_k)$$
    5. w(k+1)=w(k)+λkpkw^{(k+1)}=w^{(k)}+λ_kp_k
    6. 计算gk+1=g(wk+1)g_{k+1}=g(w^{k+1}),若gk=0g_k=0,则停止计算;否则,按下式求出Bk+1B_{k+1}:$$B_{k+1}=B_k+\frac{y_ky_kT}{y_kTδ_k}-\frac{B_kδ_kδ_kTB_k}{δ_kTB_kδ_k}$$,其中$$y_k=g_{k+1}-g_k,δ_k=w{(k+1)}-w{(k)}$$
    7. 置k=k+1,转3

条件随机场的预测算法

条件随机场的预测问题是给定条件随机场P(Y|X)和输入序列(观测序列)x,求条件概率最大的输出序列(标记序列)y*,即对观测序列进行标注

条件随机场的预测问题,可以转为求非规范化概率最大的最优路径问题:$$\max\limits_y(w·F(y,x))$$,其中,$$w=(w_1,w_2,…,w_k)T$$,$$F(y,x)=(f_1(y,x),f_2(y,x),…,f_k(y,x))T$$, f_k(y,x)=\sum_{i=1}^nf_k(y_{i-1},y_i,x,i),k=1,2,…,K

  • 条件随机场预测的维特比算法
    • 初始化
      • δ1(j)=wF1(y0=start,yi=j,x),j=1,2,...,mδ_1(j)=w·F_1(y_0=start,y_i=j,x),j=1,2,...,m
    • 递推。对于i=2,3,…,n
      • δi(l)=max1jmδi1(j)+wFi(yi1=j,yi=l,x),l=1,2,...,mδ_i(l)=\max\limits_{1≤j≤m}\\{δ_{i-1}(j)+w·F_i(y_{i-1}=j,y_i=l,x)\\},l=1,2,...,m
      • Ψi(l)=argmaxijmδi1(j)+wFi(yi1=j,yi=l,x),l=1,2,...,mΨ_i(l)=\arg\max\limits_{i≤j≤m}\\{δ_{i-1}(j)+w·F_i(y_{i-1}=j,y_i=l,x)\\},l=1,2,...,m
    • 终止
    • 返回路径
      • y_i^\*=Ψ_{i+1}(y_{i+1}^\*),i=n-1,n-2,..,1
    • 求得最优路径y^\*=(y_1^\*,y_2^\*,...,y_n^\*)

附录·梯度下降法

梯度下降法或最速下降法是求解无约束最优化问题的一种最常用的方法。

  • 梯度下降算法
    • 取初始值x(0)Rnx^{(0)}∈R^n,置k=0
    • 计算f(x^{(k)})
    • 计算梯度gk=g(xk)g_k=g(x^{k}),当gk=ε||g_k||=ε时,停止迭代,令x^\*=x^{(k)};否则,令x^\*=x^{(k)};否则,令pk=g(x(k))p_k=-g(x^{(k)}),求λkλ_k,使$$f(x{(k)}+λ_kp_k)=\min\limits_{λ≥0}f(x{(k)}+λp_k)$$

附录·牛顿法和拟牛顿法

牛顿法和拟牛顿法也是求解无约束最优化问题的常用方法,有收敛速度快的优点

  • 牛顿法
    1. 取初始点x(0)x^{(0)},置k=0
    2. 计算gk=g(x(k))g_k=g(x^{(k)})
    3. gkε||g_k||<ε,则停止计算,得近似解x^\*=x^{(k)}
    4. 计算Hk=H(x(k))H_k=H(x^{(k)}),并求pkp_k:$$H_kp_k=-g_k$$
    5. x(k+1)=x(k)+pkx^{(k+1)}=x^{(k)}+p_k
    6. 置k=k+1,转2
  • 拟牛顿法
    • 在牛顿法中,需要计算海塞矩阵的逆矩阵H1H^{-1},该计算比较复杂,所以考虑用一个n阶矩阵Gk=G(x(k))G_k=G(x^{(k)})来近似替代H_K{-1}=H{-1}(x^{k})$
    • 在牛顿法迭代中,海塞矩阵HkH_k满足:$$g_{k+1}-g_k=H_k(x{(k+1)}-x{(k)})$$,记yk=gk+1gk,δk=x(k+1)x(k)y_k=g_{k+1}-g_k,δ_k=x^{(k+1)}-x^{(k)},则$$y_k=H_kδ_k$$或$$H_k^{-1}y_k=δ_k$$,此称为拟牛顿条件
    • 拟牛顿法将GkG_k作为Hk1H_k^{-1}的近似,要求矩阵GkG_k是正定的,且满足$$G_{k+1}y_k=δ_k$$
    • 每一次迭代中,拟牛顿法可以选择更新矩阵Gk+1G_{k+1}:$$G_{k+1}=G_k+△G_k$$
  • DFP算法
    • 假设每一步迭代中矩阵Gk+1G_{k+1}是由GkG_k加上两个附加项构成的,即$$G_{k+1}=G_k+P_k+Q_k$$
    • 其中Pk,QkP_k,Q_k是待定矩阵,此时:$$G_{k+1}y_k=G_ky_k+P_ky_k+Q_ky_k$$
    • 为了使Gk+1G_{k+1}满足拟牛顿条件,可使得Pk,QkP_k,Q_k满足:$$P_ky_k=δ_k$$ $$Q_ky_k=-G_ky_k$$
    1. 选定初始点x(0)x^{(0)},取G0G_0为正定对称矩阵,置k=0
    2. 计算g_k=g(x^{(k)}).若gk<ε||g_k||<ε,则停止计算,得近似解x^\*=x^{(k)};否则转3
    3. pk=Gkgkp_k=-G_kg_k
    4. 一维搜索:求λ_k使得$$f(x{(k)}+λ_kp_k)=\min\limits_{λ≥0}f(x{(k)}+λp_k)$$
    5. 置x{(k+1)}=x{(k)}+λ_kp_k
    6. 计算g_{k+1}=g(x^{k+1}),若gk+1<ε||g_{k+1}||<ε,则停止计算,得近似解x^\*=x^{(k_1)};否则,Qk=GkykykTGkykTGkykQ_k=-\frac{G_ky_ky_k^TG_k}{y_k^TG_ky_k}以此算出Gk+1G_{k+1}
    7. 置k=k+1,转3
  • BFGS算法
    • 这是最流行的拟牛顿算法
    • 可以用BkB_k逼近海塞矩阵H,即拟牛顿条件$$B_{k+1}δ_k=y_k$$
    1. 选定初始点x(0)x^(0),取B0B_0为正定对称矩阵,置k=0
    2. 计算gk=g(x(k))g_k=g(x^{(k)})gk<ε||g_k||<ε,则停止计算,得近似解x^\*=x^{(k)};否则转
    3. Bkpk=gkB_kp_k=-g_k,求得pkp_k
    4. 一维搜索:求λ_k使得$$f(x{(k)}+λ_kp_k)=\min\limits_{λ≥0}f(x{(k)}+λp_k)$$
    5. x(k+1)=x(k)+λkpkx^{(k+1)}=x^{(k)}+λ_kp_k
    6. 计算gk+1=g(xk+1)g_{k+1}=g(x^{k+1}),若gk+1<ε||g_{k+1}||<ε,则停止计算,得近似解x^\*=x^{(k_1)};否则,Bk+1=Bk+ykykTykTδkBkδkδkTBkδkTBkδkB_k+1=B_k+\frac{y_ky_k^T}{y_k^Tδ_k}-\frac{B_kδ_kδ_k^TB_k}{δ_k^TB_kδ_k}以此算出Bk+1B_{k+1}
    7. 置k=k+1,转3