菜单

牛顿法和梯度下降法实现对数几率回归

10月 14, 2019 - 机器学习&大数据, 算法

前言

众所周知,对数几率回归,也就是Logistic回归,对于机器学习初学者来说是重要的基础,也是从小白走向大神的第一步。借着学校机器学习课的机会,我也开始尝试啃了啃西瓜书,下面这个例子就是基于自西瓜书课后习题,watermelon3.0alpha数据集进行的对数几率回归。

本文将使用两种方法,牛顿法和梯度下降法分别实现回归过程,从而形成一定的对比。

算法原理推导

对数几率回归原理推导

首先,为实现二分类任务,使用sigmoid函数对线性回归进行激活,sigmoid函数如下:
y = \frac{1}{1+e^{-z}}
其中,z为线性回归的表示形式,即z = \theta^Tx+b

与线性回归中常用的均方误差(MSE)不同,对数几率回归其损失函数的形式为:
\begin{aligned}
Loss(\theta,b) =& -[y\times\log\hat{y} + (1-y)\times\log(1-\hat{y})]\\
=& -[y\times\log{\frac{1}{1+e^{-z}}} + (1-y)\times\log(1-\frac{1}{1+e^{-z}})]\\
=& -[-y\log(1+e^{-z}) + (1-y)\log(\frac{e^{-z}}{1+e^{-z}})]\\
=& -[-y\log(1+e^{-z})+\log(\frac{e^{-z}}{1+e^{-z}})-y\log(\frac{e^{-z}}{1+e^{-z}})]\\
=& -[-y\log(1+e^{-z})+\log(\frac{1}{e^{z}+1})-y\log(e^{-z})+y\log(1+e^{-z})]\\
=& \quad y\times z + \log(e^z+1)
\end{aligned}

至此,我们的问题转化为优化上述损失函数,使得其值最小,分别使用梯度下降法牛顿法优化上述损失函数:

梯度下降法

在数据格式化处理阶段,将参数b与\theta放在同一个矩阵内进行运算,因此,梯度下降的迭代过程可以表示为:
\theta = \theta – \alpha\frac{\partial L}{\partial \theta}
其中\alpha为人工指定的学习率,\frac{\partial L}{\partial \theta}由下述方式求得:
\begin{align}
\frac{\partial L}{\partial \theta} =& \frac{\partial L}{\partial z}\frac{\partial z}{\partial \theta} \\
=& \frac{\partial L}{\partial \hat{y}}\frac{\partial \hat{y}}{\partial z}\frac{\partial z}{\partial \theta}\\
\end{align}

因此,分别对\frac{\partial L}{\partial \hat{y}},\frac{\partial \hat{y}}{\partial z},\frac{\partial z}{\partial \theta}进行求解:
\frac{\partial L}{\partial \hat{y}} = \frac{\hat{y}-y}{\hat{y}(1-\hat{y})}

\begin{align}
\frac{\partial \hat{y}}{\partial z}=&(\frac{1}{1+e^{-z}})^{‘}\\
=&\frac{e^{-z}}{(1+e^{-z})^2} = \frac{1+e^{-z}-1}{(1+e^{-z})^2}\\
=&\frac{1}{(1+e^{-z})}-\frac{1}{(1+e^{-z})^2}\\
=& \hat{y} – \hat{y}^2\\
=& \hat{y}(1-\hat{y})
\end{align}

这一步其实是很关键的,看似复杂的sigmoid函数,对其求导后的结果竟然变成\hat{y}(1-\hat{y}),这也是很多人在没有手动推导公式的情况下看到这个时感到不理解的原因(比如我)。
\frac{\partial z}{\partial \theta} = x

综上可得,\frac{\partial L}{\partial \theta} = (\hat{y}-y)x,即参数\theta的更新方式为:
\theta = \theta – \alpha (\hat{y}-y)x

牛顿法

在牛顿法中,其标准形式的向量参数\theta的更新迭代过程可表示为:
\theta = \theta – H(L)^{-1}\frac{\partial L}{\partial \theta}
其中H(\theta)为二阶导数\frac{\partial^2 L}{\partial \theta^2}矩阵,即Hession矩阵,其构造方法为:
\frac{\partial^2 L}{\partial \theta^2}=X diag \{(1-\hat{y})\hat y\}_{m,m}X^T
为什么会有这样一个看起来很奇怪的构造方法呢,在我的理解中,与Hession矩阵的形式有关,众所周知,Hession矩阵长这个样子:
H(L)= \begin{Bmatrix} \frac{\partial^2 L}{\partial \theta_1^2}&\frac{\partial^2 L}{\partial \theta_1\theta_2}&\cdots&\frac{\partial^2 L}{\partial \theta_1\theta_n}\\ \frac{\partial^2 L}{\partial \theta_2\theta_1}&\frac{\partial^2 L}{\partial \theta_2^2}&\cdots&\frac{\partial^2 L}{\partial \theta_2\theta_n} \\ \vdots&\vdots&\ddots&\vdots \\ \frac{\partial^2 L}{\partial \theta_n\theta_1}&\frac{\partial^2 L}{\partial \theta_n\theta_2}&\cdots&\frac{\partial^2 L}{\partial \theta_n^2} \end{Bmatrix}
又根据上述过程中求出的\frac{\partial L}{\partial \theta} = (\hat{y}-y)x可知,最基本基础的二阶偏导为:
\frac{\partial^2 L}{\partial \theta^2} = x\times y(1-y)\times x
因此使用diag\{(1-\hat{y})\hat{y}\}_{m,m}构建对角矩阵,再分别左乘X,右乘X^T,利用矩阵运算的方式代替求和过程,求出Hession矩阵。

依据该过程进行迭代,直到两次迭代产生的损失函数差值小于10^{-5}

算法的实现

Logistic_Regression函数

首先利用read_data()函数读入训练集数据,并构造两个theta矩阵,分别代表两种方式求得的\theta

分别调用Gradient_Descent()函数和Newtown()函数求解参数值并输出结果。

Gradient_Descent()函数

该函数接收训练集,参数,学习率,训练次数作为其参数,在每一轮训练中,首先利用np.dot()函数计算h值,即z=\theta x+b的值,再将其带入sigmoid函数,存入p中。依据上述原理推导中的迭代公式,求解一阶导数并不断进行梯度下降过程。

完成梯度下降后,进行数据可视化处理。

Newtown()函数

​ 与梯度下降过程类似,对于每一次迭代过程,求解一阶导数temp1,利用np.diag()函数构建对角矩阵,并将其二阶导数矩阵存入temp2中,依据上述原理推导公式进行迭代。

​ 完成牛顿法后,进行数据可视化处理。

结果可视化

源代码

牛顿法和梯度下降法实现对数几率回归》有2个想法

zxx

bb好厉害

回复

发表评论

您的电子邮箱地址不会被公开。 必填项已用*标注