1. 指数分布族
1.1 定义
指数族分布(exponential family of distributions)亦称指数型分布族,在上世纪30年代中期被提出,在概率学和统计学中,它是一些有着特殊形式的概率分布的集合,是统计中最重要的参数分布族,包含了二项分布、正态分布,指数分布、伯努利分布、泊松分布、gamma分布、beta分布等.指数分布族为很多重要而且常用的概率分布提供了统一的框架,这种一般性有助于表达的方便和从更大的宏观尺度理解这些分布.
指数分布族的概率密度函数PDF或者概率质量函数PMF一般表达式:
分布函数中:
1. η η 是自然参数(natural parameter), 通常是一个实数,也可以是向量
2. b(y) b ( y ) 是底层观测值(underlying measure)
3. T(y) T ( y ) 是充分统计量(sufficient statistic),通常为y
4. a(η) a ( η ) 被称为对数规则化(log normalizer)
当参数a,b,T都固定的时候,就定义了一个以 η η 为参数的函数族.
a(η)
a
(
η
)
的作用是使PDF的积分等于1,称其为对数规则化,原因如下, 将(1.1)变形:
对(1.3)两边同时积分:
根据定义,概率密度函数PDF的积分为1,所以:
只要 a(η) a ( η ) 满足上述(1.6)中等式,就能保证概率密度函数PDF的积分为1.
1.2 sigmoid函数
假设y满足二项分布,如下:
将1.7式转化合并为如下形式:
1.8就是二项分布的概率密度函数,所以,我们可以将其转换为指数分布族的通用框架形式:
和1.1公式比较,可知:
由公式1.9可知:
由公式1.11可知,二项分布的概率参数 θ θ `和指数分布族的自然参数` η η 之间的关系.
公式1.11就是sigmoid函数. 将式1.11带入式1.10,可知:
2. Logstic回归
2.1 推导模型
假设:
1. 样本集和样本:
2. 假设分类y在x给定的情况下,满足二项分布,并且其为1的概率是样本特征x的函数,参数 θ θ :
已知二项分布满足指数分布族通用框架,在其他参数固定的情况下,不同的自然参数 η η 代表着不同概率的二项分布,所以,我们假设自然参数和样本特征之间存在线性的关系:
结合公式2.5,可知:
由于公式2.2的假设, 所以,Logstic回归也属于广义线性模型, 这里并不直接假设y和x之间的关系,否则就是线性模型,而是假设在已知x的情况下,y的取值存在一个概率,x和这个概率之间存在关系,而 η η `又代表了不同概率的同一分布,所以很自然假设x和` η η 之间存在线性关系,并由此来训练模型.
2.2 构建模型
已知y的概率分布,将其应用到所有样本,求取参数向量
θT
θ
T
的值,使其应当满足,所有样本的概率之积最大,即假设,已经收集到的样本是最有可能出现的, 未收集到的样本出现的概率相对较低.
如上,根据需求,应该求解2.4关于 θ θ 的函数,使其结果最大.由于连乘积较难计算最大值,所以,左右两边同时求对数,因为对数函数递增,不影响原函数的升降.
为了求解式2.5的最大值,这里选择使用 梯度上升法, 顾名思义,和梯度下降法相反,不断地叠加梯度,进而不断靠近极值点.第一步需要求 θj θ j 的偏导数:
到这里,就需要计算 hθ(x) h θ ( x ) `关于` θj θ j 的偏导数了:
将式2.7和式2.8带入公式2.6中,可得:
结果很优美,就是分类y减去预测分类 hθ(x) h θ ( x ) 的差值,乘以x.
至此,我们能够知道参数 θj θ j 的更新规则:
其中, α α `是步长参数,用以控制每次更新` θ θ 的增量大小,也叫学习速率,这个需要在实践中进行调整,其值过小会导致多次迭代才能收缩,过大会导致越过最优点发生震荡现象.
梯度下降会导致局部极值点的产生,解决方法是随机初始化,寻找多个最优点结果,在这些最优解中寻找最终结果.但是,对于本例,再次对公式2.9求导可知,二次导数小于0,属于凸函数(看定义,有些教材是二次导数大于0,称为凸函数,凸凹函数一样,只有一个最优解),最优解只有一个.
当我们使用Logistic回归时,只需要使用公式2.3和公式2.10即可.
3. 代码实现
3.1 伪代码
梯度下降法(梯度上升法)有两种实现方式:
1. 批梯度下降法 (batch gradient descent)
更新每一参数时,都使用所有的样本进行更新.
Repeat Until Converge{
}
优点:全局最优解;迭代次数相对较少;易于并行实现;
缺点:当样本数目很多时,训练过程会很慢
2. 随机梯度下降法(stochastic gradient descent) 或者 增量梯度下降法(increment gradient descent)
为了解决样本数据量过大的问题,提出了随机梯度下降法,通过每个样本来迭代更新一次,不需要每次都计算全部样本:
Repeat Until Converge{
for i = 1,2,…m{
}
}
如果样本量很大的情况下,例如十几万,可能只需要其中的几万条或者几千条,就已经迭代到最优解,对比上面的批量梯度下降,迭代一次需要用到十几万训练样本,一次迭代不可能最优,如果迭代10次的话就需要遍历训练样本10次。但是,SGD伴随的一个问题是噪音较BGD要多,使得SGD并不是每次迭代都向着整体最优化方向, 如果连续出现相似的样本,导致迭代结果变化不大,迭代可能就此终止在该局部最优解处.
优点:训练速度快
缺点:准确度下降,并不是全局最优;不易于并行实现.
既然如此,实际上还有另外一种方案,结合1,2方案,每次迭代选取10个样本,一直迭代至收缩,这种方法称为小批量梯度下降法(MBGD),即保证训练速度比较快,又保证最终参数的准确率,理解1,2方案,该方法实现就很简单了,不详述.
3.2 python实现
# -- coding: utf-8 --
import numpy as np
import matplotlib.pyplot as plt
def func_sigmoid(x):
sig = 1.0 / (1 + np.exp(-x))
return sig
def draw_split_line(theta):
"""
hθ(z) = logistic(z)
z = θT * x z:特征空间 x:输入空间
θT * X =0 唯一确定一个超平面L, 使数据线性可分, 此时 L平面上方数据对应 z > 0, L平面下方数据对应z < 0
:param theta:
:return:
"""
x = np.linspace(-10, 10, 50)
vector_x = np.row_stack((np.ones((1, 50)), np.mat(x)))
y = -((np.dot(theta[0, 0:2], vector_x)) / theta[0, 2])
plt.plot(x, y.T)
def draw_train_data(train_x, train_y):
num_futures, num_samples = np.shape(train_x)
for i in range(0, num_samples):
color = 'or' if train_y[i, 0] == 1 else 'ob'
plt.plot(train_x[1, i], train_x[2, i], color, markersize=5)
def get_train_data():
train_x = []
train_y = []
data_file = open("../doc/logistic_train_data.txt")
for line in data_file.readlines():
line_array = line.strip().split()
train_x.append([1.0, float(line_array[0]), float(line_array[1])])
train_y.append(int(line_array[2]))
return np.mat(train_x).transpose(), np.mat(train_y).transpose()
def train(theta, train_x, train_y, times, alpha, precision, train_type):
"""
Parameters
----------
theta: 目标参数
train_x:训练集 x维度
train_y:训练集 y维度
times:最大循环次数
alpha:步长
precision:精度,满足精度时停止循环
train_type: 训练模型, base: 全量梯度下降, incremental:随机梯度/增量梯度下降法
----------
"""
for c in range(0, times):
theta_bak = theta.copy()
if train_type == 'base':
theta = train_mode_base(theta, train_x, train_y, alpha)
elif train_type == 'incremental':
theta = train_mode_incremental(theta, train_x, train_y, alpha)
error = theta_bak - theta
# 判断精度是否满足要求
if (np.fabs(error) < precision).any():
print "precision:", precision, "count:", c, "error:", error
break
return theta
# 全量梯度下降法
def train_mode_base(theta, train_x, train_y, alpha):
num_futures, num_samples = np.shape(train_x)
for j in range(0, num_futures):
step_sum = 0
for i in range(0, num_samples):
train_xi = train_x[:, i]
train_xji = train_x[j, i]
train_yi = train_y[i, 0]
theta_x = np.dot(theta, train_xi)[0, 0]
sig = func_sigmoid(theta_x)
error = train_yi - sig
step = error * train_xji
step_sum += step
theta[0, j] += alpha * step_sum
return theta
# 增量/随机梯度下降法
def train_mode_incremental(theta, train_x, train_y, alpha):
num_futures, num_samples = np.shape(train_x)
for i in range(0, num_samples):
for j in range(0, num_futures):
train_xi = train_x[:, i]
train_xji = train_x[j, i]
train_yi = train_y[i, 0]
theta_x = np.dot(theta, train_xi)[0, 0]
sig = func_sigmoid(theta_x)
error = train_yi - sig
step = error * train_xji
theta[0, j] += alpha * step
return theta
def main():
plt.figure("mode")
train_x, train_y = get_train_data()
# 绘制样本数据
draw_train_data(train_x, train_y)
# 训练
theta = np.mat([1.0, 1.0, 1.0])
theta = train(theta, train_x, train_y, 1000, 0.1, 0.001, 'base')
print theta
draw_split_line(theta)
plt.show()
if __name__ == '__main__':
main()
训练数据下载地址: logistic_train_data.txt
参考:
1. https://www.cnblogs.com/maybe2030/p/5089753.html
2. http://blog.csdn.net/saltriver/article/details/55105285
3. http://blog.csdn.net/stdcoutzyx/article/details/9113681
4. http://blog.csdn.net/stdcoutzyx/article/details/9207047