1、理解偏差(bias)和方差(variance)
首先,我们来理解什么是偏差(bias)和方差(variance):
机器学习目的是找到一个最好函数 f ^ \hat f f^,能够非常准确根据输入的数据得到所期望的结果。为了得到这个最好的函数,我们需要一堆训练数据,找到函数 f ∗ f^* f∗, f ∗ f^* f∗是最优函数 f ^ \hat f f^的估测。衡量两者之间的误差,主要来自两个方面:偏差(bias)和方差(variance)。
从统计学来看,偏差和方差的来源:
-
目的:估计一个变量 x x x的平均值
- 假设 x x x的平均值是 μ \mu μ
- 假设 x x x的方差是 σ 2 \sigma^2 σ2
-
估测 μ \mu μ的方法:
- 取样本,样本为 N N N个点: { x 1 , x 2 , … , x N } \{x^1, x^2,…, x^N\} {x1,x2,…,xN}
- 将 N N N个点求平均值: m = 1 N ∑ N x n ≠ μ m = \frac1N\sum_Nx^n \neq \mu m=N1∑Nxn̸=μ
- 再取多次样本,得到 m 1 , m 2 , m 3 . . . m_1,m_2,m_3... m1,m2,m3...
- 然后求 m m m的期望值: E [ m ] = E [ 1 N ∑ n x n ] = 1 N ∑ n E [ x n ] = μ E[m] = E[\frac1N \sum_nx^n] = \frac1N\sum_nE[x^n]=\mu E[m]=E[N1∑nxn]=N1∑nE[xn]=μ
- 上述的 m n m_n mn都是无偏差的(unbiased),但是散落在 μ \mu μ周围,存在方差 V a r [ m ] = σ 2 N Var[m] = \frac{\sigma^2}{N} Var[m]=Nσ2
- 方差的值取决于样本大小 N N N,如果 N N N比较大,方差就比较小;反之,方差比较大。
-
估测 σ 2 \sigma^2 σ2的方法:
- σ 2 \sigma^2 σ2的估测值是 s 2 = 1 N ∑ n ( x n − m ) 2 s^2 = \frac1N\sum_n(x^n-m)^2 s2=N1∑n(xn−m)2
- 此时多个不同的 s n 2 s_n^2 sn2会散落在 μ \mu μ周围,存在偏差: E [ s 2 ] = N − 1 N σ 2 ≠ σ 2 E[s^2] = \frac{N-1}{N}\sigma^2 \neq \sigma^2 E[s2]=NN−1σ2̸=σ2.这个意思应该就是样本数据相比总体会更加紧密聚集在均值周围。【注:从中心极限定理中的样本估计整体可以理解这一问题。】
- 如果把样本大小 N N N加大,方差估测值 s 2 s^2 s2与方差真实值 σ 2 \sigma^2 σ2的差距会变小。
从机器学习寻找最优函数 f ^ \hat f f^来看,偏差和方差来源:
- 为了找最优函数 f ^ \hat f f^,我们通过训练数据找到了一些函数 f ∗ f^* f∗,这些函数 f ∗ f^* f∗与最优函数在预测数据的时候会产生偏差,这就是误差的来源之一偏差(bias)了。
- 假设找到的所有函数 f ∗ f^* f∗的期望值: E [ f ∗ ] = f ‾ E[f^*]=\overline f E[f∗]=f,这些函数 f ∗ f^* f∗之间又会存在偏离,这些偏离就是方差(variance)。
- 以打靶为例,因为准星没对准靶心而产生的偏离就是偏差,准星对准了,但是因为其它因素导致偏离就是方差。讲义中下面图就是一个很好的例子:
2、偏差(bias)、方差(variance)与Model的关系
2.1 重温模型(Model)和函数(Function)
在上次学习中,没弄清楚Model和function的关系。直到在这次学习中才明白两者的差异:
-
Model:我们在面对一个机器学习问题时设计的一个函数形式,例如: y = b + w . x y = b + w.x y=b+w.x
其中的 b b b和 w w w都是未知的,我们需要求的也就是它们的值。
-
Function:针对设计好的模型,我们使用不同的训练数据,会求到不同的 w w w和 b b b的值,每一组值就对应一个函数(Function)。
所以能够使用同一个Model,得到不同的Function。
2.2 过拟合、欠拟合与偏差(bias)、方差(variance)
首先说明三类函数: 理想的最优函数: f ^ \hat f f^,一个Model通过训练数据获取的函数: f ∗ f^* f∗,函数 f ∗ f^* f∗的期望: f ‾ \overline f f。
- **欠拟合:**根据上一节学习中可知,Model设计得太简单,会造成欠拟合,此时Model误差(error)比较大,偏差(
f
‾
\overline f
f与
f
^
\hat f
f^)也比较大,但是Model的各个函数(
f
∗
f^*
f∗)散布比较紧密,即方差(variance)比较小。
- 原因很简单,简单Model参数比较少,受不同样本数据的影响也就比较小。
- 以一个极端例子为例,设计Model为: f ( x ) = c f(x) = c f(x)=c,( c c c是常数)。由此可见,不管训练数据是什么,Model的函数是不变的,方差为0,但是离最优的function差距很大。
- **过拟合:**Model设计太复杂,造成过拟合。此时误差(error)也很大,但是偏差(
f
‾
\overline f
f与
f
^
\hat f
f^)比较小,此时Model的各个函数(
f
∗
f^*
f∗)散布比较松散,方差比较大。
- 原因类似,复杂Model参数比较多,受不同样本数据的影响也就比较大。
- 方差(variance)与Model复杂程度关系图:
- 偏差((bias)与Model复杂程度关系图:
总结,老师在课程中用一张图说明了它们之间的关系:
最优的函数,应该就是Error最小的那点。
由此可见,我们在训练模型的时候,应该了解到Error的来源是Bias还是Variance。
2.3 诊断Error以及优化Model
- 在训练的时候,如果训练数据误差比较大,此时就是偏差比较大了,处于欠拟合,需要重新设计你的模型:
- 输入更多的特征;
- 设计一个更复杂的模型。
- 如果是方差比较大,训练Model时体现在训练数据偏差小,但是在测试数据上偏差比较大,方差很大。解决办法是:
- 非常有效的方法:收集更多的数据;(可以创造新数据,调角度、上下、左右颠倒等)
- 正则化。
3 、回顾:梯度下降
设计一个损失函数 L L L,其中输入是Model产生的函数 f ∗ f^* f∗(或者参数 w w w和 b b b),通过梯度下降,能够在损失函数上找到一个最优点,使得损失最小,从而找到最佳函数。但是,使用梯度下降时需要注意的几点:
3.1 调整学习率
需要小心设置学习率
η
\eta
η,如果学习率太大,可能造成错过全局最优点在之间来回摆动;如果太小,梯度调整速度太慢。因此,需要画图体现Loss的变化:
因此,动态调整learning rate是比较重要的:
-
一个简单且受欢迎的想法是:每几轮训练之后,通过一些因素减少学习率。
- 训练开始时,距离目标点比较远,所以用比较大的学习率。
- 在几轮训练之后,接近目标点,所以减少学习率。
- 例子:学习率 η \eta η随更新次数 t t t的变化: η t = η / t + 1 \eta^t=\eta/\sqrt{t+1} ηt=η/t+1
-
学习率不应该只设置一个。
-
给不同参数不同的学习率:
-
参数 w w w更新数学推导如下:
首先,定义 t t t为更新次数, η t = η t + 1 \eta^t =\frac{\eta}{\sqrt{t+1}} ηt=t+1η是第 t t t次的学习率, g t = d L ( θ t ) d w g^t = \frac{dL(\theta^t)}{dw} gt=dwdL(θt)是第 t t t次时损失函数对参数 w w w的偏微分,也就是梯度, σ t \sigma^t σt是前 t t t项 g t g^t gt的均方根。
-
w 1 = w 0 − η 0 σ 0 g 0 w^1=w^0-\frac{\eta^0}{\sigma^0}g^0 w1=w0−σ0η0g0, σ 0 = ( g 0 ) 2 \sigma^0 = \sqrt{(g^0)^2} σ0=(g0)2
-
w 2 = w 1 − η 1 σ 1 g 1 w^2 = w^1-\frac{\eta^1}{\sigma^1}g^1 w2=w1−σ1η1g1, σ 1 = 1 2 [ ( g 0 ) 2 + ( g 1 ) 2 ] \sigma^1 = \sqrt{\frac12[(g^0)^2+(g^1)^2]} σ1=21[(g0)2+(g1)2]
-
w 3 = w 2 − η 2 σ 2 g 2 w^3 = w^2-\frac{\eta^2}{\sigma^2}g^2 w3=w2−σ2η2g2, σ 2 = 1 2 [ ( g 0 ) 2 + ( g 1 ) 2 + ( g 2 ) 2 ] \sigma^2 = \sqrt{\frac12[(g^0)^2+(g^1)^2+(g^2)^2]} σ2=21[(g0)2+(g1)2+(g2)2]
…
-
w t + 1 = w t − η t σ t g t w^{t+1} = w^t-\frac{\eta^t}{\sigma^t}g^t wt+1=wt−σtηtgt, σ t = 1 t + 1 ∑ i = 0 t ( g t ) 2 \sigma^t = \sqrt{\frac{1}{t+1}\sum_{i=0}^t(g^t)^2} σt=t+11∑i=0t(gt)2
-
-
最后,推导出式子: w t + 1 = w t − η ∑ t = 0 t ( g i ) 2 g t w^{t+1}=w^t-\frac{\eta}{\sqrt{\sum^t_{t=0}(g^i)^2}}g^t wt+1=wt−∑t=0t(gi)2ηgt
-
3.2 随机梯度下降
损失函数的形式: L = ∑ n ( y ^ n − ( b + ∑ w i x i n ) ) 2 L = \sum_n(\hat y^n-(b+\sum w_ix_i^n))^2 L=∑n(y^n−(b+∑wixin))2。损失函数确实应该考虑所有训练样本数据的误差和。下面说明梯度下降和随机梯度下降的区别:
- 梯度下降: θ i = θ i − 1 − η ∇ L ( θ i − 1 ) \theta^i = \theta^{i-1}-\eta \nabla L(\theta^{i-1}) θi=θi−1−η∇L(θi−1),训练完所有样本之后,再更新参数;
- 随机梯度下降:随机或按照顺序取一个样本
x
n
x^n
xn,只考虑这个样本的损失:
L n = ( y ^ n − ( b + ∑ w i x i n ) ) 2 L^n = (\hat y^n-(b+\sum w_ix_i^n))^2 Ln=(y^n−(b+∑wixin))2 , θ i = θ i − 1 − η ∇ L ( θ i − 1 ) \theta^i = \theta^{i-1}-\eta \nabla L(\theta^{i-1}) θi=θi−1−η∇L(θi−1)
一个样本训练之后,就更新依次参数。
3.3 数据归一化(Feature Scaling)
- **描述:**为了说明这个知识点,按照课程的例子,假设有一个Model:
y
=
b
+
w
1
x
1
+
w
2
x
2
y = b+w_1x_1+w_2x_2
y=b+w1x1+w2x2。如果其中的两个特征值
x
1
x_1
x1和
x
2
x_2
x2的取值范围相差比较大,那就把其中一个数大小缩放到与另一个相同。这样做的原因是什么呢?看下图:
- 为了说明这个知识点,按照课程的例子,假设有一个Model: y = b + w 1 x 1 + w 2 x 2 y = b+w_1x_1+w_2x_2 y=b+w1x1+w2x2。如果其中的两个特征值 x 1 x_1 x1和 x 2 x_2 x2的取值范围相差比较大,那就把其中一个数大小缩放到与另一个相同。这样做的原因是什么呢?看下图:
可以看到如果 x 1 x_1 x1和 x 2 x_2 x2取值范围相差比较大, w 1 w_1 w1和 w 2 w_2 w2对损失函数的影响不一样,会导致梯度下降时参数更新效率变慢,在更新时还需要考虑学习率的改变。反之,两者规模如果差不多,情况会好很多。
怎么做?
Model的输入特征值是 x 1 x 2 x 3 … x r … x R x^1x^2x^3…x^r…x^R x1x2x3…xr…xR,每个特征值有一堆样本数据 x i r x^r_i xir。
m i m_i mi表示全部特征值第 i i i组数据的平均值; σ i \sigma_i σi表示全部特征值第 i i i组数据的标准差。
数据归一化的方式: x i r = x i r − m i σ i x^r_i = \frac{x^r_i-m_i}{\sigma_i} xir=σixir−mi。
4、数学理论
4.1 梯度下降说明
梯度下降的原理其实就是在某个点附近找到一个相对最小的值,如下图所示:
这个的理论基础就是我们所学的泰勒公式了。
4.2 泰勒展开式
-
单变量泰勒展开式
定义:函数 h ( x ) h(x) h(x)在 x 0 x_0 x0处的展开式如下所示
h ( x ) = ∑ k = 0 ∞ h k ( x 0 ) k ! ( x − x 0 ) k = h ( x 0 ) + h ′ ( x 0 ) ( x − x 0 ) + h ′ ′ ( x 0 ) 2 ! ( x − x 0 ) 2 + . . . h(x)= \sum_ {k = 0}^{\infty} \frac {h ^ k(x_0)} {k!}(x-x_0)^ k = h (x_0)+ {h}'(x_0)(x-x_0)+ \frac {h''(x_0)} {2!}(x-x_0)^ 2+... h(x)=k=0∑∞k!hk(x0)(x−x0)k=h(x0)+h′(x0)(x−x0)+2!h′′(x0)(x−x0)2+...
当 x x x趋近于 x 0 x_0 x0的时候,展开式形式换成: h ( x ) ≈ h ( x 0 ) + h ′ ( x 0 ) ( x − x 0 ) h(x)≈ h(x_0)+ {h}'(x_0)(x-x_0) h(x)≈h(x0)+h′(x0)(x−x0) -
多变量泰勒展开式
课件中给出了展开式形式:
-
将泰勒展开式与损失函数相联系
图中损失函数随机取一点 ( a , b ) (a, b) (a,b)
-
损失函数的多变量泰勒展开形式: L ( θ ) ≈ L ( a , b ) = d L ( a , b ) d θ 1 ( θ 1 − a ) + d L ( a , b ) d θ 2 ( θ 2 − b ) L(\theta)≈L(a,b)=\frac{dL(a,b)}{d\theta_1}(\theta_1-a)+\frac{dL(a,b)}{d\theta_2}(\theta_2-b) L(θ)≈L(a,b)=dθ1dL(a,b)(θ1−a)+dθ2dL(a,b)(θ2−b)
-
为了简化,令 s = L ( a , b ) s=L(a,b) s=L(a,b), u = d L ( a , b ) d θ 1 u=\frac{dL(a,b)}{d\theta_1} u=dθ1dL(a,b), v = d L ( a , b ) d θ 2 v=\frac{dL(a,b)}{d\theta_2} v=dθ2dL(a,b)
-
化简后,形式是: L ( θ ) ≈ s + u ( θ 1 − a ) + v ( θ 2 − b ) L(\theta)≈s+u(\theta_1-a)+v(\theta_2-b) L(θ)≈s+u(θ1−a)+v(θ2−b)
-
此时,需要找的是小红圈内最小的损失函数,因此需要满足条件:
( θ 1 − a ) 2 (\theta_1-a)^2 (θ1−a)2+ ( θ 2 − b ) 2 ≤ d 2 (\theta_2-b)^2 \leq d^2 (θ2−b)2≤d2
-
如果要求最小的损失函数,也就是下一个点 ( θ 1 , θ 2 ) (\theta_1,\theta_2) (θ1,θ2),就变成了求两个向量 ( u , v ) (u, v) (u,v)和 ( θ 1 − a , θ 2 − b ) (\theta_1-a,\theta_2-b) (θ1−a,θ2−b)的内积最小值,如图所示,向量 ( u , v ) (u, v) (u,v)已知, ( θ 1 − a , θ 2 − b ) (\theta_1-a,\theta_2-b) (θ1−a,θ2−b)在反方向时,内积最小:
-
整理之后,求最小的损失函数的式子如下:
这就是梯度下降了。
-
**注意:**只有红色圈圈的半径无穷小,式子才能保证成立。即学习率 η \eta η需要足够小,才能满足式子的条件。
5、梯度下降的限制
陷入局部最优,在局部最优点处,参数就停止更新了,此时虽然梯度为0,但这只是极值而不是最小值。如下图,参数处在一个比较高的Loss的时候,但是梯度已经很平缓了,会造成已经到了最优值的地方的假象。
6、学习总结
总结主要包含对此次学习过程中,对一些问题的理解。
6.1 batch、mini-batch与SGD
- Batch:对整个训练集进行梯度下降时,整个数据集就是一个batch。这时候就是batch梯度下降。
- SGD:当每次只对一个样本进行梯度下降的时候,是 随机梯度下降(SGD)。
- mini-batch:当每次处理样本的个数在上面二者之间,就是 mini batch 梯度下降。
- 优缺点分析:当数据集很大时,训练算法是非常慢的,和 batch 梯度下降相比,使用 mini batch 梯度下降更新参数更快,有利于更鲁棒地收敛,避免局部最优。和 stochastic 梯度下降相比,使用 mini batch 梯度下降的计算效率更高,可以帮助快速训练模型。
- 代码Mini-Batch SGD代码:(参考代码)
import numpy as np
import os
def get_data(obj_path_name):
pro_path = os.path.abspath('.')
data_path = str(pro_path + obj_path_name)
print data_path
data = np.loadtxt(data_path)
x = np.asarray(np.column_stack((np.ones((data.shape[0], 1)), data[:, 0:-1])), dtype='double')
y = data[:, -1]
y = np.reshape(y, (data.shape[0], 1))
print x.shape, y.shape
return x, y
def inverse_hessian_mat(x):
# 计算hessian矩阵,并求其逆矩阵
x_hessian = np.dot(x.T, x)
inverse_hessian = np.linalg.inv(x_hessian)
return inverse_hessian
def get_e(x, theta, y):
y_pre = np.dot(x, theta)
e = y_pre - y
return e
def compute_grad(x, e, stand_flag=0):
# batch法必须做梯度归一化
print e.T.shape, x.shape
grad = np.dot(e.T, x)
# grad = np.dot(e.T, x)/x.shape[0]
if stand_flag == 1:
grad = grad/(np.dot(grad, grad.T)**0.5)
return np.reshape(grad, (x.shape[1], 1))
def get_cost(e):
# print e.shape
# 计算当前theta组合点下的各个样本的预测值 y_pre
cost = np.dot(e.T, e) / 2
# cost = np.dot(e.T, e) / (2*e.shape[0])
return cost
def get_cost_l12(e, theta, m, l=1, lmd=10e10):
# print e.shape
if l == 1:
cost = (np.dot(e.T, e) + lmd*sum(abs(theta))) / (2*m)
elif l == 2:
cost = (np.dot(e.T, e) + lmd*np.dot(theta.T*theta)) / (2*m)
else:
cost = (np.dot(e.T, e)) / (2*m)
return cost
def update_batch_theta(theta, grad, alpha):
theta = theta - alpha*grad
return theta
def update_batch_theta_l12(theta, grad, alpha, m, l=1, lmd=10e10):
if l == 1:
theta = theta - alpha * (grad + (lmd/m)*theta)
elif l == 2:
theta = theta - alpha * (grad + (lmd/m))
else:
theta = theta - alpha * grad
return theta
def iter_batch(x, theta, y, out_n, out_e_reduce_rate, alpha):
step = 1
while step < out_n:
"""计算初始的损失值"""
if step == 1:
e = get_e(x, theta, y)
cost_0 = get_cost(e)
"""计算当前theta组合下的grad值"""
grad = compute_grad(x, e, stand_flag=1)
"""依据grad更新theta"""
theta = update_batch_theta(theta, grad, alpha)
"""计算新的损失值"""
e = get_e(x, theta, y)
cost_1 = get_cost(e)
e_reduce_rate = abs(cost_1 - cost_0)/cost_0
print 'Step: %-6s, cost: %s' % (step, cost_1[0, 0])
if e_reduce_rate < out_e_reduce_rate:
flag = 'Finish!\n'
print flag
break
else:
cost_0 = cost_1
step += 1
return theta
def iter_random_batch(x, theta, y, out_n, out_e_reduce_rate, alpha):
step = 0
while step < out_n:
step += 1
for i in range(x.shape[0]):
x_i = np.reshape(x[i, ], (1, x.shape[1]))
y_i = np.reshape(y[i, ], (1, 1))
"""计算初始的损失值"""
e_0 = get_e(x, theta, y)
cost_0 = get_cost(e_0)
"""用一个样本,计算当前theta组合下的grad值"""
e_i = get_e(x_i, theta, y_i)
grad = compute_grad(x_i, e_i, stand_flag=1)
"""依据grad更新theta"""
theta = update_batch_theta(theta, grad, alpha)
"""计算新的损失值"""
e_1 = get_e(x, theta, y)
cost_1 = get_cost(e_1)
e_reduce_rate = abs(cost_1 - cost_0)/cost_0
if e_reduce_rate < out_e_reduce_rate:
flag = 'Finish!\n'
print flag
step = out_n + 1
break
print 'Step: %-6s, cost: %s' % (step, cost_1[0,0])
return theta
def iter_mini_batch(x, theta, y, out_n, out_e_reduce_rate, alpha, batch):
batch_n = x.shape[0]//batch
batch_left_n = x.shape[0] % batch
step = 0
while step < out_n:
step += 1
for i in range(batch_n+1):
"""计算初始的损失值"""
e_0 = get_e(x, theta, y)
cost_0 = get_cost(e_0)
"""选取更新梯度用的batch个样本"""
if i <= batch_n-1:
start = i*batch
x_i = np.reshape(x[start:start+batch, ], (batch, x.shape[1]))
y_i = np.reshape(y[start:start+batch, ], (batch, 1))
else:
if batch_left_n != 0:
x_i = np.reshape(x[-1*batch_left_n:, ], (batch_left_n, x.shape[1]))
y_i = np.reshape(y[-1*batch_left_n:, ], (batch_left_n, 1))
"""用batch个样本,计算当前theta组合下的grad值"""
e_i = get_e(x_i, theta, y_i)
grad = compute_grad(x_i, e_i, stand_flag=1)
"""依据grad更新theta"""
theta = update_batch_theta(theta, grad, alpha)
"""计算新的损失值"""
e_1 = get_e(x, theta, y)
cost_1 = get_cost(e_1)
e_reduce_rate = abs(cost_1 - cost_0)/cost_0
if e_reduce_rate < out_e_reduce_rate:
flag = 'Finish!\n'
print flag
step = out_n
break
print 'Step: %-6s, cost: %s' % (step, cost_1[0, 0])
return theta
def update_newton_theta(x_hessian_inv, theta, grad, alpha):
theta = theta - alpha*np.dot(x_hessian_inv, grad)
# print 'New Theta --> ', theta1
return theta
def iter_newton(x, theta, y, out_n, out_e_reduce_rate, alpha):
e = get_e(x, theta, y)
cost_0 = get_cost(e)
x_hessian_inv = inverse_hessian_mat(x)
step = 1
while step < out_n:
grad = compute_grad(x, e)
theta = update_newton_theta(x_hessian_inv, theta, grad, alpha)
e = get_e(x, theta, y)
cost_1 = get_cost(e)
print 'Step: %-6s, cost: %s' % (step, cost_1[0,0])
e_reduce_rate = abs(cost_1 - cost_0)/cost_0
if e_reduce_rate < out_e_reduce_rate:
flag = 'Finish!\n'
print flag
break
else:
cost_0 = cost_1
step += 1
return theta
def iter_batch_linesearch(x, theta_0, y, out_n, out_e_reduce_rate, alpha):
e_0 = get_e(x, theta_0, y)
cost_0 = get_cost(e_0)
step = 1
while step < out_n:
grad = compute_grad(x, e_0, stand_flag=1)
theta_1, cost_1, e_1, count = line_search(x, y, alpha, theta_0, grad)
e_reduce_rate = abs(cost_1 - cost_0)/cost_0
print 'Step: %-6s, count: %-4s, cost: %s' % (step, count, cost_1[0, 0])
if e_reduce_rate < out_e_reduce_rate:
flag = 'Finish!\n'
print flag
break
else:
cost_0 = cost_1
theta_0 = theta_1
e_0 = e_1
step += 1
return theta_1
def line_search(x, y, alpha, theta_0, grad):
# 不更新梯度,在当前梯度方向上,迭代100次寻找最佳的下一个theta组合点,
e_0 = get_e(x, theta_0, y)
cost_0 = get_cost(e_0)
max_iter, count, a, b = 1000, 0, 0.8, 0.5
while count < max_iter:
# 随count增大,alpha减小,0.8*0.8*0.8*0.8*...
alpha_count = pow(a, count) * alpha
theta_1 = theta_0 - alpha_count*grad
e_1 = get_e(x, theta_1, y)
cost_1 = get_cost(e_1)
# 当前theta组合下的梯度为grad, grad == tan(w) == 切线y变化量/theta变化量, 推出 切线y变化量 = grad * theta变化量
grad_dy_chage = abs(np.dot(grad.T, theta_1-theta_0))
# 实际y变化量为cost0-cost1, 不能加绝对值,如果加了就不收敛了。只有cost_y_chage>0且大于b*grad_dy_chage时才符合条件
cost_y_change = cost_0-cost_1
# 当前梯度方向上,实际y变化量 占 切线y变化量的比率为b,b越大越好,至少大于0.5才行,实际损失减小量要占 dy的一半以上。
if cost_y_change > b*grad_dy_chage:
break
else:
count += 1
return theta_1, cost_1, e_1, count
def iter_bfgs(x, theta_0, y, out_n, out_e_reduce_rate, dfp):
e_0 = get_e(x, theta_0, y)
cost_0 = get_cost(e_0)
grad_0 = np.reshape(compute_grad(x, e_0), (x.shape[1], 1))
hk = np.eye(x.shape[1])
step = 1
while step < out_n:
dk = -1*np.dot(hk, grad_0)
theta_1, cost_1, e_1, l_count = armijo_line_search(x, y, theta_0, grad_0, dk, cost_0)
grad_1 = compute_grad(x, e_1)
# print theta_1
# print grad_1
yk = grad_1 - grad_0
sk = theta_1 - theta_0
condition = np.dot(sk.T, yk)
# 更新以此theta就更新一次hk
if dfp == 1:
if condition > 0:
hk = hk - X2(X3(hk, yk, yk.T), hk)/X3(yk.T, hk, yk)+(X2(sk, sk.T))/condition
else:
if condition > 0:
hk = hk + (1+X3(yk.T, hk, yk)/condition)*(X2(sk, sk.T)/condition)-(X3(sk, yk.T, hk)+X3(hk, yk, sk.T))/condition
e_reduce_rate = abs(cost_1 - cost_0) / cost_0
print 'Step: %-6s, l_count: %-6s, cost: %s' % (step, l_count, cost_1[0,0])
if e_reduce_rate < out_e_reduce_rate:
flag = 'Finish!\n'
print flag
break
cost_0, grad_0, theta_0 = cost_1, grad_1, theta_1
step += 1
return theta_1
def iter_lbfgs(x, theta_0, y, out_n, out_e_reduce_rate):
limit_n = 5
ss, yy = [], []
step = 1
while step < out_n:
if step == 1:
e_0 = get_e(x, theta_0, y)
cost_0 = get_cost(e_0)
grad_0 = compute_grad(x, e_0)
dk = -1*grad_0
theta_1, cost_1, e_1, l_count = armijo_line_search(x, y, theta_0, grad_0, dk, cost_0)
grad_1 = compute_grad(x, e_1)
if len(ss) > limit_n and len(yy) > limit_n:
del ss[0]
del yy[0]
yk = grad_1 - grad_0
sk = theta_1 - theta_0
ss.append(sk)
yy.append(yk)
qk = grad_1
k = len(ss)
condition = X2(yk.T, sk)
alpha = []
for i in range(k):
# t 4->0 倒着计算,倒着存alpha
t = k-i-1
pt = 1 / X2(yy[t].T, ss[t])
alpha.append(pt*X2(ss[t].T, qk))
qk = qk - alpha[i] * yy[t]
z = qk
for i in range(k):
t = k - i - 1
# i 0->4 正着计算,正着存beta
pi = 1 / (X2(yy[i].T, ss[i])[0, 0])
# pi数
beta = pi*X2(yy[i].T, z)
# beta[i]数
z = z+ss[i]*(alpha[t] - beta)
if condition > 0:
dk = -z
e_reduce_rate = abs(cost_1 - cost_0) / cost_0
print 'Step: %-6s, l_count: %-6s, cost: %s' % (step, l_count, cost_1[0, 0])
if e_reduce_rate < out_e_reduce_rate:
flag = 'Finish!\n'
print flag
break
cost_0, grad_0, theta_0 = cost_1, grad_1, theta_1
step += 1
return theta_1
def armijo_line_search(x, y, theta_0, grad_0, dk_0, cost_0):
# print theta_0.shape, grad_0.shape, dk_0.shape, cost_0.shape
# 不更新梯度,在当前梯度方向上,迭代100次寻找最佳的下一个theta组合点,
max_iter, count, countk, a, b = 100, 0, 0, 0.55, 0.4
while count < max_iter:
"""
batch方法使用梯度方向grad更新theta,newton等使用牛顿方向dk来更新theta
newton:hessian逆*grad; dfp:当前步dfpHk; bfgs:当前步bfgsHk
当前梯度方向上,实际y变化量 占 切线y变化量的比率为b,b越大越好,至少大于0.5才行,实际损失减小量要占 dy的一半以上。
"""
"""更新theta"""
theta_1 = theta_0 + pow(a, count)*dk_0
"""计算损失"""
e_1 = get_e(x, theta_1, y)
cost_1 = get_cost(e_1)
cost_y_change = cost_1 - cost_0
dy_change = b * pow(a, count) * np.dot(grad_0.T, dk_0)
# if cost_1 < cost_0 + b * pow(a, count) * np.dot(grad_0.T, dk_0):
if cost_y_change < dy_change:
"""如果一直不满足条件,那么一直没有将count赋值给countk,countk仍为0"""
countk = count
break
count += 1
theta_1 = theta_0 + pow(a, countk) * dk_0
e_1 = get_e(x, theta_1, y)
cost_1 = get_cost(e_1)
return theta_1, cost_1, e_1, count
if __name__ == '__main__':
x1, y1 = get_data('/optdata/line_sample_data.txt')
# x1, y1 = get_data('/optdata/airfoil_self_noise.txt')
theta1 = np.zeros(x1.shape[1]).reshape(x1.shape[1], 1)
res_theta = iter_batch(x1, theta1, y1, out_n=1e5, out_e_reduce_rate=1e-6, alpha=0.02)
# res_theta = iter_random_batch(x1, theta1, y1, out_n=1e5, out_e_reduce_rate=1e-6, alpha=0.02)
# res_theta = iter_mini_batch(x1, theta1, y1, out_n=1e5, out_e_reduce_rate=1e-6, alpha=0.01, batch=10)
# res_theta = iter_batch_linesearch(x1, theta1, y1, out_n=1e5, out_e_reduce_rate=1e-6, alpha=1)
# res_theta = iter_newton(x1, theta1, y1, out_n=1e5, out_e_reduce_rate=1e-6, alpha=1)
# res_theta = iter_bfgs(x1, theta1, y1, out_n=100, out_e_reduce_rate=1e-6, dfp=1)
# res_theta = iter_bfgs(x1, theta1, y1, out_n=100, out_e_reduce_rate=1e-6, dfp=0)
# res_theta = iter_lbfgs(x1, theta1, y1, out_n=100, out_e_reduce_rate=1e-6)
print 'Res_theta:', np.reshape(res_theta, (1, res_theta.shape[0]))[0]
# def iter_bfgs(x, theta_0, y, out_n, out_e_reduce_rate):
6.2 交叉验证
交叉验证的基本思想是把在某种意义下将原始数据(dataset)进行分组,一部分做为训练集(train set),另一部分做为验证集(validation set or test set),首先用训练集对分类器进行训练,再利用验证集来测试训练得到的模型(model),以此来做为评价分类器的性能指标。代码如下所示:
from sklearn.cross_validation import cross_val_score # K折交叉验证模块
from sklearn.datasets import load_iris # iris数据集
from sklearn.model_selection import train_test_split # 分割数据模块
from sklearn.neighbors import KNeighborsClassifier # K最近邻(kNN,k-NearestNeighbor)分类算法
#加载iris数据集
iris = load_iris()
X = iris.data
y = iris.target
#分割数据并
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=4)
#建立模型
knn = KNeighborsClassifier()
#训练模型
knn.fit(X_train, y_train)
#将准确率打印出
print(knn.score(X_test, y_test))
# 0.973684210526 基础验证的准确率
#使用K折交叉验证模块
scores = cross_val_score(knn, X, y, cv=5, scoring='accuracy')
#将5次的预测准确率打印出
print(scores)
#将5次的预测准确平均率打印出
print(scores.mean())
6.3 数据归一化
如3.3节所示,如果某个特征的方差比其他特征大几个数量级,那么它就会在学习算法中占据主导位置,导致学习器并不能像我们说期望的那样,从其他特征中学习。数据归一化是让不同维度之间的特征在数值上有一定比较性,可以大大提高分类器的准确性。好处如下:
- 提升模型精度
- 提升收敛速度
6.4 回归模型评价指标
对于回归模型效果的判断指标经过了几个过程,从SSE到R-square再到Ajusted R-square, 是一个完善的过程:
-
SSE(误差平方和)
计算公式: S S E = ∑ ( y ^ − y ∗ ) 2 SSE=\sum(\hat y-y^*)^2 SSE=∑(y^−y∗)2,其中的 y ^ \hat y y^是真实值, y ∗ y^* y∗是预测值。
-
R-square(决定系数)
计算公式: R 2 = 1 − ( y ^ − y ∗ ) 2 ( y ^ − y ‾ ) 2 R^2 = 1-\frac{(\hat y-y^*)^2}{(\hat y-\overline y)^2} R2=1−(y^−y)2(y^−y∗)2,其中 y ‾ \overline y y是 y y y的平均值。
-
数学理解: 分母理解为原始数据的离散程度,分子为预测数据和原始数据的误差,二者相除可以消除原始数据离散程度的影响。
-
其实“决定系数”是通过数据的变化来表征一个拟合的好坏。
-
理论上取值范围(-∞,1], 正常取值范围为[0 1] 。实际操作中通常会选择拟合较好的曲线计算R²,因此很少出现 − ∞ -\infty −∞.
-
越接近1,表明方程的变量对y的解释能力越强,这个模型对数据拟合的也较好
越接近0,表明模型拟合的越差
经验值:>0.4, 拟合效果好
-
缺点:数据集的样本越大,R²越大,因此,不同数据集的模型结果比较会有一定的误差。
-
-
Adjusted R-Square (校正决定系数)
计算公式: R _ a d j u s t e d = 1 − ( 1 − R ) 2 ( n − 1 ) n − p − 1 R\_adjusted=1-\frac{(1-R)^2(n-1)}{n-p-1} R_adjusted=1−n−p−1(1−R)2(n−1),其中 n n n为样本数量, p p p为特征数量。
- 为了消除了样本数量和特征数量的影响。