参考资料
\textbf{参考资料}
参考资料
(1)DEEP LEARNING (GoodFellow, Bengio, Aaron)
(2)Quantum Machine Learning What Quantum Computing Means to Data Mining (Peter Wittek)
监督学习与无监督学习
机器学习可分为监督学习和五监督学习两大类,监督学习需要依赖于人为标注的数据,通过选择合适地目标函数令模型的预测不断地贴近监督信息,从而得到精准的预测。无监督学习,意在学习到数据中的潜在特征,并不涉及到监督数据。
机器学习框架与深度学习框架
当今有一些框架可以让我们不必再重复地造轮子,我们称这种程序为深度学习框架,在机器学习中也有机器学习框架。
我们在日常生活中可以常常听到,机器学习比较容易掌握的论调,那是因为通过机器学习框架进行工作,基本除了数据处理,工作量并不是很大。而且有一些工作范式,不需要很多创造力。
关于框架,在机器学习中有sklearn,深度学习中有pytorch以及tensorflow等。
下文中的推导和代码是底层方面的内容,对于仅仅想应用的人,机器学习方面话只需学会sklearn,深度学习方面,一般都有source code。
其他需要说明的问题
本文所提供的代码是基于数学公式进行实现,并不是基于机器学习框架进行介绍。在代码实现过程中,仅仅做了简单实现,如果应用,还是使用机器学习框架最佳。
线性回归
define: 给定特征数据集 X X X,给定标签 Y Y Y,基于类似于 y = k x + b y=kx+b y=kx+b这种类似的函数,通过目标函数去减小模型的预测 Y ^ \hat{Y} Y^同标签 Y Y Y之间的差距。
优化指标
(1)梯度下降
梯度下降是用来最小化目标函数的一种数学方法,假设目标函数为
M
S
E
e
r
r
o
r
(
Y
;
θ
,
X
)
=
∣
p
r
e
(
X
;
θ
)
−
Y
∣
=
∣
Y
^
−
Y
∣
MSE_{error}(Y;\theta,X)=|pre(X;\theta)-Y|=|\hat{Y}-Y|
MSEerror(Y;θ,X)=∣pre(X;θ)−Y∣=∣Y^−Y∣,我们试图最小化
M
S
E
e
r
r
o
r
MSE_{error}
MSEerror,那么我们必须对参数
θ
\theta
θ进行更新,所以我们对每一个
θ
\theta
θ处求偏导数,让
θ
\theta
θ在偏导数方向减小,从而可以最快地减小目标函数,而减小
M
S
E
e
r
r
o
r
MSE_{error}
MSEerror就是在最小化误差的过程,所以梯度下降就是在优化模型的过程。
(2)
M
S
E
e
r
r
o
r
MSE_{error}
MSEerror
M
S
E
e
r
r
o
r
=
∣
Y
^
−
Y
∣
MSE_{error}=|\hat Y-Y|
MSEerror=∣Y^−Y∣反应真实空间和预测值之间的距离。在数学上我们将这个公式形容为
L
2
n
o
r
m
L2 norm
L2norm, 或者说欧氏距离.
如果在偏导数上存在疑问可参考我之前一篇文章
https://blog.youkuaiyun.com/captainAAAjohn/article/details/118215413
\url{https://blog.youkuaiyun.com/captainAAAjohn/article/details/118215413}
https://blog.youkuaiyun.com/captainAAAjohn/article/details/118215413
∇
w
1
m
∑
i
=
0
m
∣
Y
^
i
−
Y
i
∣
=
0
\nabla_w\frac{1}{m}\sum_{i=0}^{m}|\hat{Y}_{i}-Y_{i}| =0
∇wm1∑i=0m∣Y^i−Yi∣=0
→
∇
w
1
m
∑
i
=
0
m
∣
w
X
i
−
Y
i
∣
=
∇
w
1
m
∑
i
=
0
m
(
X
i
w
−
Y
i
)
T
(
X
i
w
−
Y
i
)
\rightarrow \nabla_w\frac{1}{m}\sum_{i=0}^{m}|wX_i-Y_i|=\nabla_w\frac{1}{m}\sum_{i=0}^{m}(X_iw-Y_i)^T(X_iw-Y_i)
→∇wm1∑i=0m∣wXi−Yi∣=∇wm1∑i=0m(Xiw−Yi)T(Xiw−Yi)
→
∇
w
1
m
∑
i
=
0
m
(
w
T
X
i
T
−
Y
i
T
)
(
X
i
w
−
Y
i
)
=
∇
w
1
m
∑
i
=
0
m
[
w
T
X
i
T
X
i
w
−
w
T
X
i
T
Y
i
−
Y
i
T
X
i
w
+
Y
i
T
Y
i
]
\rightarrow \nabla_w\frac{1}{m}\sum_{i=0}^{m}(w^TX_i^T-Y_i^T)(X_iw-Y_i)=\nabla_w\frac{1}{m}\sum_{i=0}^{m}[w^TX_i^TX_iw-w^TX_i^TY_i-Y_i^TX_iw+Y_i^TY_i]
→∇wm1∑i=0m(wTXiT−YiT)(Xiw−Yi)=∇wm1∑i=0m[wTXiTXiw−wTXiTYi−YiTXiw+YiTYi]
→
∇
w
1
m
∑
i
=
0
m
[
w
T
X
i
T
X
i
w
−
w
T
X
i
T
Y
i
−
Y
i
T
X
i
w
+
Y
i
T
Y
i
]
=
1
m
∇
w
∑
i
=
0
m
[
w
T
X
i
T
X
i
w
−
2
w
T
X
i
T
Y
i
+
Y
i
T
Y
i
]
\rightarrow \nabla_w\frac{1}{m}\sum_{i=0}^{m}[w^TX_i^TX_iw-w^TX_i^TY_i-Y_i^TX_iw+Y_i^TY_i]=\frac{1}{m} \nabla_w \sum_{i=0}^{m}[w^TX_i^TX_iw-2w^TX_i^TY_i+Y_i^TY_i]
→∇wm1∑i=0m[wTXiTXiw−wTXiTYi−YiTXiw+YiTYi]=m1∇w∑i=0m[wTXiTXiw−2wTXiTYi+YiTYi]
→
1
m
∑
i
=
0
m
[
2
X
i
T
X
i
w
−
2
X
i
T
Y
i
]
\rightarrow\frac{1}{m}\sum_{i=0}^{m}[2X_i^TX_iw-2X_i^TY_i]
→m1∑i=0m[2XiTXiw−2XiTYi]
有些情况下模型还需带有偏置,通过添加偏置进一步推导可得:
∇
b
M
S
E
e
r
r
o
r
=
2
X
w
T
+
b
−
2
Y
\nabla_b MSE_{error}=2Xw^T+b-2Y
∇bMSEerror=2XwT+b−2Y
(3)分析
那么这样做为何是合理的?为何这样做就可以减少预测和真实标签之间的差异。
理解1:
M
S
E
e
r
r
o
r
MSE_{error}
MSEerror减小了真实值和预测值之间的差异,所被减小的这种差异是具体到预测向量中每一个location中的value,由于一个预测向量的空间中全部预测的数值相对于真实值的误差都被减小了,因此模型变得精确了。
理解2: 从概率密度角度去理解
从概率密度角度证明MSE的合理性
线性回归可以被视作通过参数
w
w
w与特征空间
X
X
X去预测
Y
^
\hat{Y}
Y^的过程。
Y
^
\hat{Y}
Y^相当于在
X
i
X_i
Xi处模型预测的平均值,因此高斯正态分布函数可以被描绘为:
N
(
Y
;
Y
^
;
σ
2
)
=
1
σ
2
π
exp
(
−
1
2
∑
i
=
0
m
(
Y
i
−
Y
i
^
σ
)
2
)
\mathcal N(Y; \hat{Y};\sigma^2)=\frac{1}{\sigma \sqrt{2\pi}} \exp(-\frac{1}{2} \sum_{i=0}^{m}(\frac{Y_i-\hat{Y_i}}{\sigma})^2)
N(Y;Y^;σ2)=σ2π1exp(−21∑i=0m(σYi−Yi^)2)
(1)数值下溢
在机器学习中我们常常会遇到数值下溢的情况,即目标值几乎为0,而机器学习的参数更新模式一般为梯度下降,如果遇到数值下溢的情况,则不能得到一个足够大的梯度令模型的参数得到更新,所以我们一般会对目标函数做数学处理从而使得在原本梯度下溢的情况下,所计算到的梯度足够大。通过取对数的方式可以有效避免数值下溢的现象。
由于高斯函数为指数函数恒大于0,所以可直接取对数,通过对数转化可得以下形式:
→
log
(
1
σ
2
π
)
−
1
2
∑
i
=
0
m
(
Y
i
−
Y
i
^
σ
)
2
=
−
1
2
σ
2
∑
i
=
0
m
(
Y
i
−
Y
i
^
)
2
+
log
(
1
σ
2
π
)
\rightarrow \log(\frac{1}{\sigma \sqrt{2\pi}} )-\frac{1}{2} \sum_{i=0}^{m}(\frac{Y_i-\hat{Y_i}}{\sigma})^2=-\frac{1}{2\sigma^2} \sum_{i=0}^{m}(Y_i-\hat{Y_i})^2+\log(\frac{1}{\sigma \sqrt{2\pi}} )
→log(σ2π1)−21∑i=0m(σYi−Yi^)2=−2σ21∑i=0m(Yi−Yi^)2+log(σ2π1)
(2)最大似然估计
我们引入最大似然估计的概念:
θ
=
arg max
θ
P
r
e
d
i
c
t
i
o
n
(
Y
∣
X
;
θ
)
\theta=\argmax_{\theta}Prediction(Y|X;\theta)
θ=θargmaxPrediction(Y∣X;θ)
从这个数学公式中我们可以直观地了解到,最大似然估计即为找到令模型在特征空间上得到最准确预测的一组参数。
(3)最大似然估计和
M
S
E
e
r
r
o
r
MSE_{error}
MSEerror
ok,我们已经得到了模型预测的密度函数形式
−
1
2
σ
2
∑
i
=
0
m
(
Y
i
−
Y
i
^
)
2
+
log
(
1
σ
2
π
)
-\frac{1}{2\sigma^2} \sum_{i=0}^{m}(Y_i-\hat{Y_i})^2+\log(\frac{1}{\sigma \sqrt{2\pi}} )
−2σ21∑i=0m(Yi−Yi^)2+log(σ2π1)
最大化这个密度函数的过程相当于最小化
1
2
σ
2
∑
i
=
0
m
(
Y
i
−
Y
i
^
)
2
\frac{1}{2\sigma^2} \sum_{i=0}^{m}(Y_i-\hat{Y_i})^2
2σ21∑i=0m(Yi−Yi^)2的过程,因此我们从概率密度的角度又一次证明了
M
S
E
e
r
r
o
r
MSE_{error}
MSEerror作为目标函数的合理性。
代码实现:
import numpy as np
from sklearn import datasets
import matplotlib.pyplot as plt
class LinearRegression:
def __init__(self,input_dim,output_dim,weight=None,bias=None,use_bias=False):
self.input_dim=input_dim
self.output_dim=output_dim
self.weight=weight
if self.weight!=None:
w,h=self.weight.shape
assert w==input_dim and h==output_dim,print('dim mistake')
else:
self.weight=np.random.uniform(-1,1,size=(output_dim,input_dim))
self.use_bias=use_bias
if self.use_bias:
if bias!=None:
assert self.bias.shape[1]==output_dim,print('dim mistake')
self.bias=bias
else:
self.bias=np.zeros((1,output_dim))
def MSE_target_weight(self,x,y,w):
return np.mean((x.T@x@w.T)-(x.T@y))
def MSE_target_bias(self,x,y,w,b):
return np.mean(2*x@w.T+b-2*y)
def MSE(self,y_,y):
return np.mean((y_-y)**2)
def train(self,x,y,epoches,learning_rate=0.001):
if not isinstance(x,np.ndarray):
x=np.array(x)
for i in range(epoches):
pre=self.__call__(x)
mse_error=self.MSE(pre,y)
print(f'Epoch: {i}---------training loss: {mse_error}')
target_w=self.MSE_target_weight(x,y,self.weight)
self.weight-=learning_rate*target_w.T
if self.use_bias:
target_b=self.MSE_target_bias(x,y,self.weight,self.bias)
self.bias-=target_b
def __call__(self,x):
if not isinstance(x,np.ndarray):
x=np.array(x)
assert x.shape[1]==self.weight.shape[1],print(f"dim 1 of x: {x.shape[1]} shold be equal to dim 1 of weight :{self.weight.shape[0]}")
if self.use_bias:
return np.matmul(x,self.weight.T)+self.bias
else:
return np.matmul(x,self.weight.T)
if __name__=='__main__':
D=datasets.load_iris()
x=D.data
y=D.target+1
y=y.reshape((y.shape[0],1))
print(x.shape,y.shape)
model = LinearRegression(input_dim=4, output_dim=1, use_bias=True)
model.train(x, y, epoches=1000,learning_rate=0.000001)
plt.scatter(np.array(np.arange(x.shape[0])),y,label='true')
plt.scatter(np.array(np.arange(x.shape[0])),model(x),label='pre')
plt.legend()
plt.show()
ok,关于密度分布这个话题其实还没有结束。
最大后验分布 MAP
我们已经引入了极大似然的概念,同时常常还有一个对应的概念“最大后验估计”,也就是MAP。
按照MAP的概念我们无法直接求得
P
(
θ
;
x
)
P(\theta;x)
P(θ;x),我们需要基于贝叶斯统计根据
P
(
x
;
θ
)
P(x;\theta)
P(x;θ)以及先验
P
(
θ
)
P(\theta)
P(θ)去表示
P
(
θ
;
x
)
P(\theta;x)
P(θ;x),也就是
P
(
θ
;
x
)
=
P
(
x
;
θ
)
P
(
θ
)
P
(
x
)
P(\theta;x)=\frac{P(x;\theta)P(\theta)}{P(x)}
P(θ;x)=P(x)P(x;θ)P(θ).为了方便处理,正常情况下会作对数处理也就是
ln
P
(
θ
;
x
)
+
P
(
θ
)
−
ln
P
(
x
)
\ln P(\theta;x)+P(\theta)-\ln P(x)
lnP(θ;x)+P(θ)−lnP(x).由于需要考虑整体的数据所以
P
(
x
)
P(x)
P(x)为1,所以上公式也可以被描绘为
ln
P
(
θ
;
x
)
+
P
(
θ
)
\ln P(\theta;x)+P(\theta)
lnP(θ;x)+P(θ)。
所以MAP过程可以被描绘为如下:
θ
M
L
=
arg max
θ
p
(
θ
∣
x
)
→
arg max
θ
log
p
(
x
∣
θ
)
+
arg max
t
h
e
t
a
log
p
(
θ
)
\theta_{ML}=\argmax_{\theta}p(\theta|x) \rightarrow \argmax_{\theta} \log p(x|\theta)+ \argmax_{theta} \log p(\theta)
θML=θargmaxp(θ∣x)→θargmaxlogp(x∣θ)+thetaargmaxlogp(θ)
所以MAP和最大似然区别在多了一项先验
P
(
θ
)
P(\theta)
P(θ)。
PS:这个MAP和mean Average Precision不是同一个东西.
牛顿法
上文我们提到了通过
M
S
E
e
r
r
o
r
MSE_{error}
MSEerror去优化函数的梯度下降法,其实还有一种通过微分去求解的方式。
(1)hessian matrix
我们求导时,求一次导可以得悉函数的局部斜率,求两次导可以得悉函数的曲率。当函数有多个变量时,我们求二阶导数可以通过hessian矩阵进行表示。通过结合hessian矩阵我们可判断某点是哪一种局部极值点,正定hessian矩阵对应着局部最小值点,负定hessian矩阵对应着局部最大值点。如果有一部分是正的另一部分是负的那么就会形成一种马鞍的形状。
(2)泰勒公式
ok,我们明白了hessian matrix的一系列原理后,就可以根据其判断局部最小值的性质找到直接跳跃到最小点的表示形式。
思路如下:
假设hessian matrix正定处的点为
x
0
x_0
x0,我们在
x
0
x_0
x0处展开
f
(
x
)
f(x)
f(x),可得:
f
(
x
)
=
f
(
x
0
)
+
0.5
(
x
−
x
0
)
T
∇
x
f
(
x
0
)
+
0.5
(
x
−
x
0
)
T
H
(
f
)
(
x
0
)
(
x
−
x
0
)
f(x)=f(x_0)+0.5(x-x_0)^T\nabla_x f(x_0)+0.5(x-x_0)^TH(f)(x_0)(x-x_0)
f(x)=f(x0)+0.5(x−x0)T∇xf(x0)+0.5(x−x0)TH(f)(x0)(x−x0)
PS:
H
(
f
)
(
x
0
)
H(f)(x_0)
H(f)(x0)为
f
(
x
)
f(x)
f(x)在
x
0
x_0
x0处的hessian矩阵。
∂
f
∂
x
=
0.5
∇
x
f
(
x
0
)
+
0.5
H
(
f
)
(
x
0
)
(
x
−
x
0
)
=
0
\frac{\partial f}{\partial x}=0.5\nabla_x f(x_0)+0.5H(f)(x_0)(x-x_0)=0
∂x∂f=0.5∇xf(x0)+0.5H(f)(x0)(x−x0)=0
→
x
=
x
0
−
H
(
f
)
(
x
0
)
T
∇
x
f
(
x
0
)
\rightarrow x=x_0-H(f)(x_0)^T\nabla_x f(x_0)
→x=x0−H(f)(x0)T∇xf(x0)
由此可以一次性得出局部最小点的位置。
支持向量机(SVM)
KKT条件
给定一个
f
(
x
)
f(x)
f(x)基于等式
g
(
x
)
g(x)
g(x)以及不等式
h
(
x
)
h(x)
h(x)去最大化
f
(
x
)
f(x)
f(x).很明显我们无法直接最大化
f
(
x
)
f(x)
f(x)。因为
f
(
x
)
f(x)
f(x)和
g
(
x
)
g(x)
g(x)是条件,这些条件相当于限制了x的定义域,因此也限制了
f
(
x
)
f(x)
f(x)值域,这时我们就需要引入两个因子找出
f
(
x
)
f(x)
f(x)同
f
(
x
)
f(x)
f(x),
g
(
x
)
g(x)
g(x),
h
(
x
)
h(x)
h(x)这三个函数组成的广义拉格朗日函数之间的等价关系:
g
(
x
)
=
0
g(x)=0
g(x)=0
h
(
x
)
<
=
0
h(x)<=0
h(x)<=0
a
>
=
0
a>=0
a>=0
→
L
(
x
,
λ
,
a
)
=
max
min
λ
max
a
λ
g
(
x
)
+
a
h
(
x
)
+
f
(
x
)
\rightarrow L(x,\lambda,a)= \max \min_{\lambda} \max_{a}\lambda g(x)+ah(x)+f(x)
→L(x,λ,a)=maxminλmaxaλg(x)+ah(x)+f(x)
如此一来我们发现上公式等价于
max
f
(
x
)
\max f(x)
maxf(x)的条件是
a
g
(
x
)
=
0
ag(x)=0
ag(x)=0
这里不难理解,因为
h
(
x
)
<
0
h(x)<0
h(x)<0所以最大化
a
h
(
x
)
ah(x)
ah(x)只能
a
g
(
x
)
=
0
ag(x)=0
ag(x)=0.而
h
(
x
)
=
=
0
h(x)==0
h(x)==0所以
h
(
x
)
h(x)
h(x)无须考虑。
在满足上述条件的情况下对
L
(
x
,
λ
,
a
)
L(x,\lambda,a)
L(x,λ,a)求导数,即可找到极值点。
总结下KKT条件:
1. 梯度为0时为最优解
2. 所有引入的因子都满足既定关系
3.
a
g
(
x
)
=
0
ag(x)=0
ag(x)=0
超平面
可以想象为一个平面,将所有的点分配到各自的位置,从而实现分类。
几何距离
二维平面为点到直线距离,多于二维可以通过线性代数表示。
初中时我们学过
y
=
A
X
+
C
A
2
y=\frac{AX+C}{\sqrt{A^2}}
y=A2AX+C
如果多于2维:
y
=
A
T
X
+
C
A
A
T
y=\frac{A^TX+C}{\sqrt{AA^T}}
y=AATATX+C
在分类问题中我们希望最大化几何距离,因为几何距离越大则分类越明显。
同时,对于二分类样本可以通过正负标签进行分类所以原形式可进一步描绘为: y = Y ^ ( A X + C ) A 2 y=\frac{\hat Y(AX+C)}{\sqrt{A^2}} y=A2Y^(AX+C),其中 Y ^ \hat Y Y^代表着标签,取值一般为-1或者+1.这样无论是正样本还是负样本几何距离都为正。
我们似乎无法直接对上述公式进行优化,接下来我们将引入函数距离。
函数距离
y
′
=
Y
^
(
A
X
+
C
)
y'=\hat Y(AX+C)
y′=Y^(AX+C)
我们设置一个约束条件,即函数距离需要大于1。
同时我们假设此时的函数距离为一个常数,由于最大最小化优化和正常数系数无关,所以我们将几何距离的分子设置为1,那么几何距离就可以被表示为: 1 ∣ A ∣ \frac{1}{|A|} ∣A∣1,同时我们希望最大化几何距离,这个命题的等价形式为最小化 1 2 A 2 \frac{1}{2}A^2 21A2.
加之我们已经设置了约束条件,所以此优化命题为:
min
A
,
B
1
2
A
2
,
s
.
t
.
1
−
Y
^
(
A
X
+
C
)
<
0
\min_{A,B} \frac{1}{2}A^2,s.t. \, 1-\hat Y(AX+C)<0
minA,B21A2,s.t.1−Y^(AX+C)<0
对偶问题
我们之前已经引入了广义拉格朗日的概念,套用模板,引入一个因子
a
>
=
0
a>=0
a>=0
那么拉格朗日函数
L
(
A
,
C
,
a
)
=
1
2
A
2
+
a
(
1
−
Y
^
(
A
X
+
C
)
)
\mathscr{L}(A,C,a)= \frac{1}{2}A^2+a( 1-\hat Y(AX+C))
L(A,C,a)=21A2+a(1−Y^(AX+C))
由于我们是对一组数作处理,所以需要表示为求和形式:
L
(
A
,
C
,
a
)
=
1
2
A
2
+
∑
a
i
−
∑
a
^
i
Y
i
(
A
X
i
+
C
)
\mathscr{L}(A,C,a)= \frac{1}{2}A^2+\sum a_i-\sum\hat a_iY_i(AX_i+C)
L(A,C,a)=21A2+∑ai−∑a^iYi(AXi+C)
原命题是一个最小最大问题:
min
A
,
C
max
a
L
(
A
,
C
,
a
)
\min_{A,C} \max_a \mathscr{L}(A,C,a)
minA,CmaxaL(A,C,a)
原命题的对偶问题为最大最小问题:
max
a
min
A
,
C
L
(
A
,
C
,
a
)
\max_a \min_{A,C}\mathscr{L}(A,C,a)
maxaminA,CL(A,C,a)
对偶问题是原命题的下界,但是我们需要证明一下。
max
a
L
>
=
L
>
=
min
A
,
C
L
\max_a \mathscr{L}>=\mathscr{L}>=\min_{A,C}\mathscr{L}
maxaL>=L>=minA,CL
→
max
a
L
>
=
min
A
,
C
L
\rightarrow \max_a \mathscr{L}>=\min_{A,C}\mathscr{L}
→maxaL>=minA,CL
如果想让上公式成立那么左边的下界应大于右边的上界所以我们可以得到:
→
min
A
,
C
max
a
L
>
=
max
a
min
A
,
C
L
\rightarrow \min_{A,C}\max_a \mathscr{L}>=\max_a\min_{A,C}\mathscr{L}
→minA,CmaxaL>=maxaminA,CL
所以我们只要是能找到上式的界限就可找到原命题的下界。
我们分别对对偶问题的
A
,
C
A,C
A,C求偏导数,可得:
A
=
∑
a
i
Y
i
X
i
A=\sum a_iY_iX_i
A=∑aiYiXi,
∑
a
i
Y
i
=
0
\sum a_iY_i=0
∑aiYi=0
将上述两个公式带入到对偶问题中可得:
−
1
2
∑
a
i
Y
i
X
i
∑
a
j
Y
j
X
j
+
∑
a
i
-\frac{1}{2}\sum a_iY_iX_i\sum a_jY_jX_j +\sum a_i
−21∑aiYiXi∑ajYjXj+∑ai
这里声明一点
(
a
1
+
a
2
+
a
i
)
(
b
1
+
b
2
+
b
i
)
(a_1+a_2+~a_i)(b_1+b_2+~b_i)
(a1+a2+ ai)(b1+b2+ bi)这种类似的计算除了分别计算出连个子项再相乘以外,也可以将两项中的每一项相乘再加起来,参考这种思路可得:
−
1
2
∑
∑
a
i
a
j
Y
i
Y
j
X
i
X
j
T
+
∑
a
i
-\frac{1}{2}\sum \sum a_i a_j Y_i Y_j X_iX_j^T +\sum a_i
−21∑∑aiajYiYjXiXjT+∑ai
最后,原问题最终转化为
max
a
−
1
2
∑
∑
a
i
a
j
Y
i
Y
j
X
i
X
j
T
+
∑
a
i
s
.
t
.
a
>
=
0
,
∑
a
i
Y
i
=
0
\max_a -\frac{1}{2}\sum \sum a_i a_j Y_i Y_j X_iX_j^T +\sum a_i \qquad s.t.\,\,a>=0,\sum a_iY_i=0
maxa−21∑∑aiajYiYjXiXjT+∑ais.t.a>=0,∑aiYi=0
其他技巧
(1)松弛因子
当我们用一个笔直的线将两类点进行分类时,会遇到一些问题,比如少量点错分,且误差不大,那么我们可以引入一个松弛因子
ε
\varepsilon
ε稍微调整一下此点的几何距离,因此几何距离就变成了
Y
i
(
A
X
i
+
C
)
>
1
−
ε
Y_i(AX_i+C)>1-\varepsilon
Yi(AXi+C)>1−ε, 这样就让平面拐弯从而实现了正确分类。
OK,我们已经在约束条件中引入了松弛因子。但是如何描绘松弛因子的重要性,我们需要引入一个重要性因子
λ
\lambda
λ,在原命题中添加正则项:
λ
∑
ε
i
\lambda \sum \varepsilon_i
λ∑εi。同时我们还需给
−
ε
i
-\varepsilon_i
−εi一个乘子
β
\beta
β至此我们得到了广义拉格朗日函数的新形势:
L
(
A
,
C
,
a
,
λ
)
=
−
1
2
A
2
+
∑
(
1
−
ε
i
)
a
i
−
∑
a
^
i
Y
i
(
A
X
i
+
C
)
+
λ
∑
ε
i
−
∑
β
ε
i
\mathscr{L}(A,C,a,\lambda)= -\frac{1}{2}A^2+\sum (1-\varepsilon_i)a_i-\sum\hat a_iY_i(AX_i+C)+\lambda \sum \varepsilon_i-\sum \beta \varepsilon_i
L(A,C,a,λ)=−21A2+∑(1−εi)ai−∑a^iYi(AXi+C)+λ∑εi−∑βεi
转化为对偶问题后分别对
λ
.
A
,
ε
\lambda.A, \varepsilon
λ.A,ε求导可得以下关系:
A
=
∑
a
i
Y
i
X
i
A=\sum a_iY_iX_i
A=∑aiYiXi
∑
a
i
Y
i
=
0
\sum a_iY_i=0
∑aiYi=0
λ
−
β
−
a
=
0
\lambda-\beta-a =0
λ−β−a=0
有意思的是我们仍可得到同未添加松弛变量相同的优化结果:
−
1
2
∑
∑
a
i
a
j
Y
i
Y
j
X
i
X
j
T
+
∑
a
i
-\frac{1}{2}\sum \sum a_i a_j Y_i Y_j X_iX_j^T +\sum a_i
−21∑∑aiajYiYjXiXjT+∑ai
只不过条件不同:
s
.
t
.
A
=
∑
a
i
Y
i
X
i
;
∑
a
i
Y
i
=
0
;
λ
=
>
a
>
=
0
s.t. \qquad A=\sum a_iY_iX_i;\sum a_iY_i=0;\lambda=>a>=0
s.t.A=∑aiYiXi;∑aiYi=0;λ=>a>=0
(2)核技巧
一些数据难以通过平面进行分类,比如由一组数据所组成的封闭的图形,像是圆形,无论你使用什么平面去划分这个封闭的圆形基本都不顶用,所以我们需要对原有的数据做一个非线性变换映射到其他的空间从而为实现分类提供契机。
怎么优化
思路如下:
根据最终求导出的形式:
max
a
−
1
2
∑
∑
a
i
a
j
Y
i
Y
j
X
i
X
j
T
+
∑
a
i
s
.
t
.
a
>
=
0
,
∑
a
i
Y
i
=
0
\max_a -\frac{1}{2}\sum \sum a_i a_j Y_i Y_j X_iX_j^T +\sum a_i \qquad s.t.\,\,a>=0,\sum a_iY_i=0
maxa−21∑∑aiajYiYjXiXjT+∑ais.t.a>=0,∑aiYi=0 找到一组
a
a
a。
然后根据
A
=
∑
a
i
Y
i
X
i
A=\sum a_iY_iX_i
A=∑aiYiXi,
∑
a
i
−
∑
a
^
i
Y
i
(
A
X
i
+
C
)
=
0
\sum a_i-\sum\hat a_iY_i(AX_i+C)=0
∑ai−∑a^iYi(AXi+C)=0计算A,C
最后根据
y
=
A
X
+
C
y=AX+C
y=AX+C得到超平面。
详细的方法一个是梯度下降另一个是SMO。
坐标下降
坐标下降思想和梯度下降类似,通过求偏导数后对原函数进行优化,但是还是有点区别,区别在梯度下降一次性全优化,而坐标下降是一次性优化一个变量。
几何上理解,梯度下降类似于一个物体目前受到多个力,速度改变(函数改变)方向是合力方向。而坐标下降是每个力不同时作用,用一个个力去改变物体的运动方向。
SMO
为何需要SMO
(1)一次仅仅改变一个参数不能在恒等约束中成立,比如
∑
a
i
y
i
=
0
\sum a_iy_i=0
∑aiyi=0一次仅仅调整一个
a
i
a_i
ai上述公式将不再成立。
(2)我们可以直接得到解析解,就避免了梯度下降中那种迭代式优化。
这里先偷个懒,因为需要写的东西太多,所以只给出解题思路:
推荐这个老师的视频,解地太详细了:
https://www.bilibili.com/video/BV1mE411p7HE?from=search&seid=1859996609451338719
\url{https://www.bilibili.com/video/BV1mE411p7HE?from=search&seid=1859996609451338719}
https://www.bilibili.com/video/BV1mE411p7HE?from=search&seid=1859996609451338719
解题思路如下:
\textbf{解题思路如下:}
解题思路如下:
(1)我们每一次仅仅更新两个alpha,所以两个更新的alpha满足:
a
1
y
1
+
a
2
y
2
=
−
∑
i
=
3
N
a
i
y
i
=
M
a_1y_1+a_2y_2=-\sum_{i=3}^N a_iy_i=\mathcal{M}
a1y1+a2y2=−∑i=3Naiyi=M
将这个公式乘以
y
1
y_1
y1可以得到
a
1
+
a
2
y
2
=
M
y
2
a_1+a_2y_2=\mathcal{M}y_2
a1+a2y2=My2
如此一来将上述关系带入到
max
a
−
1
2
∑
∑
a
i
a
j
Y
i
Y
j
X
i
X
j
T
+
∑
a
i
s
.
t
.
a
>
=
0
,
∑
a
i
Y
i
=
0
\max_a -\frac{1}{2}\sum \sum a_i a_j Y_i Y_j X_iX_j^T +\sum a_i \qquad s.t.\,\,a>=0,\sum a_iY_i=0
maxa−21∑∑aiajYiYjXiXjT+∑ais.t.a>=0,∑aiYi=0 中就变成了单变量问题。
将带入后的公式对
a
2
a_2
a2求导数,并且连理
f
(
x
)
=
∑
a
i
y
i
x
i
x
+
b
f(x)=\sum a_iy_ix_ix+b
f(x)=∑aiyixix+b可以得到以下条件:
a n e w = a 2 o l d + y 2 ( E 1 − E 2 ) ϵ a^{new}=a_2^{old}+\frac{y_2(E_1-E_2)}{\epsilon} anew=a2old+ϵy2(E1−E2)
其中:
a
n
e
w
a^{new}
anew和
a
o
l
d
a^{old}
aold分别是更新后和更新前的alpha,而
E
E
E为
f
(
x
)
−
y
f(x)-y
f(x)−y,
ϵ
{\epsilon}
ϵ为
X
1
X
1
+
X
2
X
2
−
2
X
1
X
2
X_1X_1+X_2X_2-2X_1X_2
X1X1+X2X2−2X1X2.
此时我们的推导还没结束,我们目前有两个约束条件:
0
<
=
a
i
<
=
λ
0<=a_i<=\mathcal{\lambda}
0<=ai<=λ
a
1
y
1
+
a
2
y
2
=
∑
3
N
a
i
y
i
=
M
a_1y_1+a_2y_2=\sum_{3}^Na_iy_i=\mathcal{M}
a1y1+a2y2=∑3Naiyi=M
关于
y
1
,
y
2
y_1,y_2
y1,y2有两种情况第一种就是相同第二种是不相同。
我们对
a
1
y
1
+
a
2
y
2
=
∑
3
N
a
i
y
i
=
M
a_1y_1+a_2y_2=\sum_{3}^Na_iy_i=\mathcal{M}
a1y1+a2y2=∑3Naiyi=M
(1)当
y
1
=
=
y
2
y_1==y_2
y1==y2时:
a
1
y
1
y
1
+
a
2
y
2
y
1
=
a
1
+
a
2
=
M
a_1y_1y_1+a_2y_2y_1=a_1+a_2=\mathcal{M}
a1y1y1+a2y2y1=a1+a2=M
→
a
2
=
M
−
a
1
\rightarrow a_2=\mathcal{M}-a_1
→a2=M−a1
因为
a
1
∈
[
0
,
C
]
a_1 \in[0,C]
a1∈[0,C]
所以
a
2
∈
[
M
−
C
,
M
]
a_2 \in[\mathcal{M}-C,\mathcal{M}]
a2∈[M−C,M] s.t.
a
1
+
a
2
=
M
a_1+a_2=\mathcal{M}
a1+a2=M
此时我们可以得到结论:
a
2
∈
[
a
1
+
a
2
−
C
,
a
1
+
a
2
]
a_2 \in[a_1+a_2-C,a_1+a_2]
a2∈[a1+a2−C,a1+a2]
(2)同理可得,当
y
1
!
=
y
2
y_1!=y_2
y1!=y2时:
a
1
−
a
2
=
M
a_1-a_2=\mathcal{M}
a1−a2=M
→
a
2
=
a
1
−
M
\rightarrow a_2=a_1-\mathcal{M}
→a2=a1−M
所以
a
2
∈
[
a
2
−
a
1
,
C
−
(
a
1
−
a
2
)
]
a_2 \in[a_2-a_1,C-(a_1-a_2)]
a2∈[a2−a1,C−(a1−a2)]
有了上述两个条件后,我们可以先算出
a
2
a_2
a2,然后和上界下界比一比,最后根据
a
1
y
1
+
a
2
y
2
=
a
1
N
e
w
y
1
+
a
2
N
e
w
y
2
a_1y_1+a_2y_2=a_{1New}y_1+a_{2New}y_2
a1y1+a2y2=a1Newy1+a2Newy2算出
a
1
a_1
a1.
代码实现
import numpy as np
import matplotlib.pyplot as plt
class SMO_SVM:
def __init__(self, features, labels, C=-1, random_seed=33):
"""
Params
:param features:
:param labels:
:param C:
:param random_seed:
Introduction:
-0.5*SUM(aaxxyy)+SUM(a(1-relaxation))-SUM(lamb*relaxation)+C*SUM(relaxation)
"""
np.random.seed(random_seed)
if not isinstance(features, np.ndarray) or isinstance(features, np.ndarray):
self.x = np.array(features)
self.y = np.array(labels)
self.num_features = self.x.shape[0]
self.dim = self.x.shape[1]
self.alpha = self.init_alpha()
self.relaxation = self.init_relx()
self.lamb=self.init_lamb()
self.w,self.b=self.init_super_plane()
self.C = C
def init_lamb(self):
return np.random.uniform(-1, 1, size=(self.num_features, 1))
def init_alpha(self):
return np.random.uniform(-1, 1, size=(self.num_features, 1))
def init_relx(self):
return np.random.uniform(0, 1, size=(self.num_features, 1))
def dot_sum(self, x, y):
if not isinstance(x, np.ndarray) or isinstance(y, np.ndarray):
x = np.array(x)
y = np.array(y)
try:
c = y @ x
except:
c=y*x
return c
def init_super_plane(self):
weight = np.sum(
np.array([self.alpha[i] * self.dot_sum(self.y[i], self.x[i]) for i in range(self.num_features)]),
axis=0).reshape(self.dim,1)
bias = np.array([self.y[i]*(1-self.relaxation[i]) for i in range(self.num_features)]).reshape(self.num_features,1) - self.x @ weight
return weight, bias
def E(self, index):
return self.pre(self.x[index]) - self.y[index]
def div(self,index):
return self.x[index]@self.x[index].T+self.x[index+1]@self.x[index+1].T-self.x[index]@self.x[index+1]
def compare(self,index):
if self.y[index]==self.y[index+1]:
down=self.alpha[index]+self.alpha[index+1]
uper=self.alpha[index]+self.alpha[index+1]-self.C
else:
down=self.alpha[index+1]-self.alpha[index]
uper=self.C-(self.alpha[index+1]-self.alpha[index])
if down<0:
down=0
if uper<down:
uper=down
return uper,down
def SMO(self):
for i in np.arange(0,self.num_features,2):
a1=self.alpha[i]
a2=self.alpha[i+1]
a2_new=a2+self.y[i+1]*sum((self.E(i)-self.E(i+1)))/self.div(i)
u_bond,l_bond=self.compare(i)
if a2_new>u_bond:
a2_new=u_bond
if a2_new<l_bond:
a2_new=l_bond
a1_new=a1*self.y[i]+a2*self.y[i+1]-a2*a2_new
self.alpha[i]=a1_new
self.lamb[i]=self.C-self.alpha[i]
self.lamb[i+1] = self.C- self.alpha[i+1]
self.alpha[i + 1]=a2_new
self.w, self.b = self.init_super_plane()
def pre(self,x):
return x @ self.w+self.b
if __name__=='__main__':
import random
a=np.random.uniform(-1,1,size=(10,5))
y=np.array([random.choice([-1,1]) for i in range(10)])
m=SMO_SVM(features=a,labels=y)
m.SMO()
b=m.pre(a)
b=[1 if i>0 else -1 for i in np.squeeze(b)]
print(b-y)
K近邻算法(KNN)
KNN(K-nearst-neighboor), 顾名思义这个算法是基于聚类的一种思想去进行分类,但是聚类过程中依赖训练数据中的标注,所以严格意义上该算法还是属于监督学习算法。
欧氏距离
欧氏距离可以用来衡量两个向量之间的差异程度,其数学形式为 1 2 ∣ ∣ v 1 − v 2 ∣ ∣ = 1 2 ∑ ( v 1 i − v 2 i ) 2 \frac{1}{2}||v1-v2||=\frac{1}{2} \sum(v_{1i}-v_{2i})^2 21∣∣v1−v2∣∣=21∑(v1i−v2i)2.
加权距离
如果我们给一个评价标准,A,B,C,D。所有达到A指标的产品都是好产品。但是尽管所有A类产品都是同一等级,但是它们之间也有所差异,仅仅通过等级划分也不能很好地描述这种情况。加权距离可以更详细地反映差距,加权顾名思义需要权重,权重则体现在个体在总体中所占据的比例。因此通过样本的个体距离比上所有样本距离之和可以反映样本的重要性。
KNN
(1) 给定一个样本
x
x
x,给定一个已标注数据集
D
\mathcal{D}
D,依次计算数据集
D
\mathcal{D}
D中所有样本和
x
x
x之间的欧式距离。
(2)做升序排序,取前K个样本,计算不同类别的加权距离。
(3)样本
x
x
x属于加权距离最大的样本类别
怎么选择K
K一般选择3-5,也有人通过数学推导得出30-45之间为最佳(参考资料:量子机器学习中数据挖掘的量子计算方法,page75)。通过iris前70个数据作为预测数据后80个数据作为训练数据在K为5时准确度达到0.85,而K为35的情况下准确度只有0.65
代码实现
import numpy as np
from sklearn import datasets
class KNN:
"""
input:[{class: , feature1: ,feature2: }]
output:class
"""
def __init__(self, train_data, K):
self.train_data = train_data
self.K = K
def Euclid_distance(self, v1, v2):
SUM=0
for i in v1.keys():
if i!='class':
SUM+=0.5 * (float(v1[i]) - float(v2[i])) **2
return SUM
def Distance_compute(self, sample):
distances = []
for unit in self.train_data:
distance_record = {}
distance_record['class'] = unit['class']
distance_record['distance'] = self.Euclid_distance(sample,unit)
distances.append(distance_record)
return distances
def weight_averaged_distance(self,distances):
SUM=0
for item in distances:
SUM+=item['distance']
results={}
for item in distances:
if item['class'] not in results:
results[item['class']]=0
results[item['class']]+=(item['distance']/(SUM+1e-8))*item['distance']
return results
@property
def data_process(self,features,lable):
data = []
for i in range(features.shape[0]):
feature = {
}
feature['class'] = lable[i]
for index, feature_i in enumerate(train[i]):
feature[f'f{index}'] = feature_i
data.append(feature)
return data
def pre(self, sample):
distances = self.Distance_compute(sample)
K_distance = sorted(distances, key=lambda unit: unit['distance'])[:self.K]
Results=self.weight_averaged_distance(K_distance)
predict=sorted(Results,reverse=True)
return predict[0]
if __name__=='__main__':
#laod data
D=datasets.load_iris()
x=D.data
y=D.target+1
y=y.reshape((y.shape[0],1))
#eval dataset
eval=[]
eval_r=y[70:]
for i in range(eval_r.shape[0]):
feature={}
feature['class']=eval_r[i][0]
for index,feature_i in enumerate(x[70:][i]):
feature[f'f_{index}']=feature_i
eval.append(feature)
# train dataset
train=x[80:]
label=y[80:]
data=[]
for i in range(train.shape[0]):
feature={}
feature['class']=label[i][0]
for index,feature_i in enumerate(train[i]):
feature[f'f_{index}']=feature_i
data.append(feature)
# init model
m=KNN(train_data=data,K=5)
# evaluation
correct=0
false=0
for i in range(len(eval)):
predict=m.pre(eval[i])
if predict==eval_r[i][0]:
correct+=1
else:
false+=1
acc=correct/(correct+false)
print(f"acc:{acc}")