浏览本文须要的背景知识点：线性回归、最大似然预计、一丢丢编程常识

一、引言

  后面几节咱们学习了规范线性回归，而后介绍了三种正则化的办法 – 岭回归、Lasso 回归、弹性网络回归，这些线性模型解决的都是回归的问题。最开始还介绍了两种简略的算法 -PLA 与口袋算法，他们解决的是分类问题。

  那么咱们能应用回归的形式来解决分类问题么，答案是必定的，这就是上面要介绍的模型 – 对数几率回归算法 ¹（Logistic Regression Algorithm），也有被直译为逻辑回归。

二、模型介绍

对数几率回归的模型函数

  既然要通过回归的形式来解决分类的问题，能够通过先进行回归剖析，而后通过一个函数将间断的后果映射成离散的分类后果，例如上面的函数表达式：

$$
y=\left\{\begin{array}{ll}
0 & \left(f\left(w^{T} x\right) \leq C\right) \\
1 & \left(f\left(w^{T} x\right)>C\right)
\end{array}\right.
$$

  在介绍对数几率回归模型之前，先来看看一类被称为 S 函数 ²（Sigmoid function）的函数，这类函数的图像成 S 型，其中最常见的一种函数为逻辑函数 ³（Logistic function），函数图像如下图所示：

  通过图像能够看到，这个函数当自变量越大时，函数值越趋近于 1，当自变量最小时，函数值越趋近于 0，当自变量为 0 时，函数值为 0.5。当下面表达式中的 C 等于 0.5 时，能够看作将间断的后果映射成离散的后果。

  逻辑函数的函数表达式为：

$$
f(z)=\frac{1}{1+e^{-z}}
$$

  将线性方程带入逻辑函数中，就失去了对数几率回归的函数方程：

$$
f(x)=\frac{1}{1+e^{-w^Tx}}
$$

对数几率回归的代价函数

  咱们须要一个代价函数来示意数据拟合的状况，这时很容易想到的是应用与线性回归一样的均方差 ⁴（mean-square error/MSE）的办法来做为代价函数

$$
\operatorname{Cost}(w)=\sum_{i=1}^{N}(y-\hat{y})^{2}
$$

  但对数几率回归是应用似然函数 ⁵（likelihood function）的对数模式来作为其代价函数，前面会阐明为什么应用这种形式比 MSE 更适合。

最大似然预计：

思考一个抛硬币的例子。假如这个硬币侧面跟背面轻重不同。咱们把这个硬币抛 80 次（侧面记为 H，背面记为 T）。并把抛出一个侧面的概率记为 p，抛出一个背面的概率记为 1 – p。假如咱们抛出了 49 个侧面，31 个背面，即 49 次 H，31 次 T。假如这个硬币是咱们从一个装了三个硬币的盒子外头取出的。这三个硬币抛出侧面的概率别离为 1/3、1/2、2/3。这些硬币没有标记，所以咱们无奈晓得哪个是哪个。应用最大似然预计，基于二项分布中的概率品质函数公式，通过这些试验数据（即采样数据），咱们能够计算出哪个硬币的可能性最大。咱们能够看到当 p = 2/3 时，似然函数获得值最大。

$$
\begin{array}{l}
L\left(p=\frac{1}{3} \mid H=49, T=31\right)=C_{80}^{49}\left(\frac{1}{3}\right)^{49}\left(1-\frac{1}{3}\right)^{31} \approx 0.000 \\
L\left(p=\frac{1}{2} \mid H=49, T=31\right)=C_{80}^{49}\left(\frac{1}{2}\right)^{49}\left(1-\frac{1}{2}\right)^{31} \approx 0.012 \\
L\left(p=\frac{2}{3} \mid H=49, T=31\right)=C_{80}^{49}\left(\frac{2}{3}\right)^{49}\left(1-\frac{2}{3}\right)^{31} \approx 0.054
\end{array}
$$

第一种示意形式：

  在二元分类的状况下，例如 y 取 0 或者 1，将对数几率回归的函数方程的后果看作概率，能够写作下式：

$$
\left\{\begin{array}{ll}
P(y=1 \mid x, w)= & h(x) \\
P(y=0 \mid x, w)= & 1-h(x)
\end{array}\right.
$$

  因为 y 只能取 0 或者 1，能够将下面两个式写成一个式子：

$$
P(y \mid x, w)=h(x)^{y}+[1-h(x)]^{1-y}
$$

  写出似然函数，其中 N 为样本总数量，将每一个概率累乘就失去了似然函数：

$$
L(w)=\prod_{i=1}^{N} h\left(X_{i}\right)^{y_{i}}\left[1-h\left(X_{i}\right)\right]^{1-y_{i}}
$$

  在 w 的所有可能取值中找一个使得似然函数取到最大值，这时求得的解就是 w 的最大似然预计 ⁶（maximum likelihood estimation/MLE）

$$
w=\underset{w}{\operatorname{argmax}}(L(w))=\underset{w}{\operatorname{argmax}}\left(\prod_{i=1}^{N} h\left(X_{i}\right)^{y_{i}}\left[1-h\left(X_{i}\right)\right]^{1-y_{i}}\right)
$$

  因为带连乘运算的代价函数不好优化，这里咱们对似然函数取自然对数并且取反，S 函数的取值为（0，1），似然函数的取值也是（0，1），对其取对数不影响其枯燥性。这样就失去了对数几率回归的代价函数：

$$
\begin{aligned}
\operatorname{Cost}(w) &=-\ln L(w) \\
&=-\sum_{i=1}^{N}\left(y_{i} \ln h\left(X_{i}\right)+\left(1-y_{i}\right) \ln \left(1-h\left(X_{i}\right)\right)\right)
\end{aligned}
$$

  因为加了一步取反的操作，这是就不是求最大值，而是将其改为求最小值：

$$
w=\underset{w}{\operatorname{argmin}}\left(-\sum_{i=1}^{N}\left(y_{i} \ln h\left(X_{i}\right)+\left(1-y_{i}\right) \ln \left(1-h\left(X_{i}\right)\right)\right)\right)
$$

第二种示意形式：

  咱们先来看下 S 函数的一个性质：

$$
\begin{aligned}
f(z) &=\frac{1}{1+e^{-z}}=\frac{e^{z}}{1+e^{z}} \\
f(-z) &=\frac{1}{1+e^{z}}=1-f(z)
\end{aligned}
$$

  在二元分类的状况下，当 y 的值取 -1 或者 1 时，将对数几率回归的函数方程的后果看作概率，能够写作下式

$$
\left\{\begin{array}{l}
P(y=1 \mid x, w)=h(x)=\frac{1}{1+e^{-w^{T} x}} \\
P(y=-1 \mid x, w)=1-h(x)=h(-x)=\frac{1}{1+e^{w^{T} x}}
\end{array}\right.
$$

  因为 y 只能取 -1 或者 1，能够将下面两个式写成一个式子：

$$
P(y \mid x, w)=\frac{1}{1+e^{-y w^{T} x}}
$$

  写出似然函数，其中 N 为样本总数量，将每一个概率累乘就失去了似然函数：

$$
L(w)=\prod_{i=1}^{N} \frac{1}{1+e^{-y_{i} w^{T} X_{i}}}
$$

  还是求最大似然预计：

$$
w=\underset{w}{\operatorname{argmax}}\left(\prod_{i=1}^{N} \frac{1}{1+e^{-y_{i} w^{T} X_{i}}}\right)
$$

  这里与第一种形式一样，咱们对似然函数取自然对数并且取反，就失去了对数几率回归的代价函数：

$$
\begin{aligned}
\operatorname{Cost}(w) &=-\ln L(w) \\
&=-\sum_{i=1}^{N} \ln \frac{1}{1+e^{-y_{i} w^{T} X_{i}}} \\
&=-\sum_{i=1}^{N} \ln 1-\ln \left(1+e^{-y_{i} w^{T} X_{i}}\right) \\
&=\sum_{i=1}^{N} \ln \left(1+e^{-y_{i} w^{T} X_{i}}\right)
\end{aligned}
$$

  也是求代价函数的最小值：

$$
w=\underset{w}{\operatorname{argmin}}\left(\sum_{i=1}^{N} \ln \left(1+e^{-y_{i} w^{T} X_{i}}\right)\right)
$$

  在 sklearn 中应用的就是下面这个代价函数。

对数几率回归的正则化：

  与线性回归一样，对数几率回归的代价函数也能够加上正则化的惩办项，有两种形式来增加正则项，一个是给正则项增加系数来管制大小，另一个是给代价函数增加系数来管制大小，实质作用是一样的，在 sklearn 中应用的是（2）式中的模式。

  L1 正则化：

$$
\begin{aligned}
\operatorname{Cost}(w)_{L 1} &=\operatorname{Cost}(w)+\lambda\|w\|_{1} & (1) \\
&=C \cdot \operatorname{Cost}(w)+\|w\|_{1} & (2)
\end{aligned}
$$

  L2 正则化：

$$
\begin{aligned}
\operatorname{Cost}(w)_{L2} &=\operatorname{Cost}(w)+\lambda\|w\|_{2}^2 & (1) \\
&=C \cdot \operatorname{Cost}(w)+\|w\|_{2}^2 & (2)
\end{aligned}
$$

  弹性网络正则化：

$$
\begin{aligned}
\operatorname{Cost}(w)_{E N} &=\operatorname{Cost}(w)+\lambda \rho\|w\|_{1}+\frac{\lambda(1-\rho)}{2}\|w\|_{2}^{2} \\
&=C \cdot \operatorname{Cost}(w)+\rho\|w\|_{1}+\frac{(1-\rho)}{2}\|w\|_{2}^{2}
\end{aligned}
$$

  对数几率回归的代价函数最小化的优化办法有多种，例如后面几节介绍过的坐标降落法。无正则化或者 L2 正则化的状况下能够用梯度降落法 ⁷（Gradient descent）、牛顿法 ⁸（Newton’s method）来等等。

  在 scikit-learn 中能够看到更多的优化办法，其中大多是后面提到算法的减速版本或是变体，例如随机均匀梯度降落法（Stochastic Average Gradient/SAG）、随机均匀梯度降落减速法（SAGA）、L-BFGS 算法（Limited-memory Broyden–Fletcher–Goldfarb–Shanno/L-BFGS），前面会一一介绍这些算法。

三、算法步骤

梯度降落法

  梯度降落法与坐标降落法的思路一样，都是通过一次一次更新权重系数 w，一步一步的迫近代价函数的最小值。坐标降落法是通过固定某一个坐标轴，求出部分最优解，而后更新 w 的值。梯度降落法令是求出函数的梯度后，沿着梯度的反方向再找到一个适合的步长系数迭代更新 w 的值。

  能够了解为碗中的一个小球从某一点自在滚落，必然会沿着变动量最大的方向挪动，最初达到最低点。

<center> 坐标降落法 VS 梯度降落法 </center>

  具体步骤：

（1）初始化权重系数 w，例如初始化为零向量

（2）计算降落方向，梯度降落法将梯度的反方向作为降落方向

$$
\Delta w = -\frac{\partial \operatorname{Cost}(w) }{\partial w}
$$

（3）应用线搜寻 ⁹ 办法来计算降落的步长 α，以使每一步迭代都满足 Wolfe 条件 ¹⁰

（4）更新权重系数 w

$$
w=w+\alpha \Delta w
$$

（5）反复步骤（2）～（4）直到梯度的大小足够小

牛顿法

  牛顿法与梯度降落法一样，也是一种降落办法，两个算法的步骤大致相同，区别就在于抉择的降落方向不一样。

<center> 绿色为梯度降落，红色为牛顿法，牛顿法的门路更加间接 </center>

  具体步骤：

（1）初始化权重系数 w，例如初始化为零向量

（2）计算降落方向，牛顿法应用二阶泰勒开展，将黑塞矩阵的逆矩阵与梯度的乘积的反方向作为降落方向

$$
\Delta w=-H^{-1} \frac{\partial \operatorname{Cost}(w)}{\partial w}
$$

（3）应用线搜寻 ⁹ 办法来计算降落的步长 α，以使每一步迭代都满足 Wolfe 条件 ¹⁰

（4）更新权重系数 w

$$
w=w+\alpha \Delta w
$$

（5）反复步骤（2）～（4）直到梯度的大小足够小

四、原理证实

对数几率回归的代价函数是凸函数

  凸函数除了通过后面几节中介绍的办法判断外，如果函数是二阶可导的，当函数的黑塞矩阵是半正定的，该函数就为凸函数。

多元二次可微的连续函数在凸集上是凸的，当且仅当它的黑塞矩阵在凸集的外部是半正定的。

• 第一种代价函数：

（1）将代价函数连加运算去掉失去一个新函数

（2）带入 h(x)

（3）将 e 的指数符号变换一下

（4）对数除法开展为对数的减法

（5）开展括号

（6）第二项与第四项一样，化简失去

$$
\begin{aligned}
f(w) &=-(y \ln h(x)+(1-y) \ln (1-h(x))) & (1)\\
&=-\left(y \ln \frac{1}{1+e^{-w^{T} x}}+(1-y) \ln \left(1-\frac{1}{1+e^{-w^{T} x}}\right)\right) & (2)\\
&=-\left(y \ln \frac{e^{w^{T} x}}{1+e^{w^{T} x}}+(1-y) \ln \left(\frac{1}{1+e^{w^{T} x}}\right)\right) & (3)\\
&=-y\left(w^{T} x-\ln \left(1+e^{w^{T} x}\right)\right)-(1-y)\left(0-\ln \left(1+e^{w^{T} x}\right)\right) & (4) \\
&=-y w^{T} x+y\left(1+e^{w^{T} x}\right)+\ln \left(1+e^{w^{T} x}\right)-y\left(1+e^{w^{T} x}\right) & (5)\\
&=\ln \left(1+e^{w^{T} x}\right)-y w^{T} x & (6)
\end{aligned}
$$

（1）下面失去的 f(w)

（2）依据求导公式求一阶导数

（3）利用 S 函数的性质，将分子化简一下并合并同类项

（4）依据求导公式求二阶导数

（5）整顿一下后果，能够看到二阶导数是一个负数乘以一个由 x 向量与 x 向量组成的矩阵

$$
\begin{aligned}
f(w) &=\ln \left(1+e^{w^{T} x}\right)-y w^{T} x & (1)\\
\frac{\partial f(w)}{\partial w} &=\frac{1}{1+e^{w^{T} x}} e^{w^{T} x} x-y x & (2)\\
&=\left(\frac{1}{1+e^{-w^{T} x}}-y\right) x & (3)\\
\frac{\partial^{2} f(w)}{\partial w \partial w^{T}} &=-\frac{1}{\left(1+e^{-w^{T} x}\right)^{2}} e^{-w^{T} x}x(-x^{T}) & (4)\\
&=\frac{e^{-w^{T} x}}{\left(1+e^{-w^{T} x}\right)^{2}} x x^{T} & (5)
\end{aligned}
$$

• 第二种代价函数：

（1）将代价函数连加运算去掉失去一个新函数

（2）依据求导公式求一阶导数

（3）利用 S 函数的性质，将分子化简一下

（4）依据求导公式求二阶导数

（5）y 的平方为 1，整顿一下后果，能够看到二阶导数也是一个负数乘以一个由 x 向量与 x 向量组成的矩阵

$$
\begin{aligned}
f(w) &=\ln \left(1+e^{-y w^{T} x}\right) & (1)\\
\frac{\partial f(w)}{\partial w} &=\frac{1}{1+e^{-y w^{T} x}} e^{-y w^{T} x}\left(-y x\right) & (2)\\
&=-\frac{y}{1+e^{y w^{T} x}} x & (3)\\
\frac{\partial^{2} f(w)}{\partial w \partial w^{T}} &=\frac{y}{\left(1+e^{y w^{T} x}\right)^{2}} e^{y w^{T} x} x \left(yx^{T}\right) & (4)\\
&=\frac{e^{y w^{T} x}}{\left(1+e^{y w^{T} x}\right)^{2}} x x^{T} & (5)
\end{aligned}
$$

（1）令二阶导数为 H 矩阵

（2）矩阵正定的判断办法，带入二阶导数，两种二阶导数的常数局部都用 C 示意，并且 C>0

（3）利用转置的性质改写一下

（4）能够看到最初的后果大于等于 0，阐明 H 为半正定矩阵

$$
\begin{aligned}
H &=\frac{\partial^{2} f(w)}{\partial w \partial w^{T}} & (1) \\
v^{T} H v &=C \cdot v^{T} x x^{T} v \quad(C>0) & (2) \\
&=C \cdot\left(x^{T} v\right)^{T}\left(x^{T} v\right) & (3) \\
&=C \cdot\left|x^{T} v\right|^{2} \geq 0 & (4)
\end{aligned}
$$

  两种新函数的黑塞矩阵都是半正定的，那么新函数是凸函数，多个凸函数相加的函数仍然是凸函数，则阐明代价函数是一个凸函数，证毕。

为什么不必 MSE 作为代价函数 ¹¹

• MSE 作为代价函数对谬误数据的惩办力度有余

  当预测数据齐全谬误时（例如理论值为 0 预测值为 1，或理论值为 1 预测值为 0）

  用 MSE 作为代价函数：

$$
\text {Cost}=(1-0)^{2}=1
$$

  应用对数似然函数作为代价函数（第一种）：

$$
\text {Cost}=-(1 \cdot \ln 0+(1-1) \cdot \ln (1-0))=\infty
$$

  应用对数似然函数作为代价函数（第二种）：
（1）理论值为 1 预测值为 -1
（2）预测值为 -1，其对应的 h(x) = 0
（3）这时线性组合的值为负无穷大
（4）这时的代价为无穷大

$$
\begin{array}{l}
y=1, \hat{y}=-1 & (1)\\
h(x)=0 & (2)\\
w^{T} x=-\infty & (3)\\
\text {Cost}=\ln \left(1+e^{-1 \cdot(-\infty)}\right)=\infty & (4)
\end{array}
$$

  能够看到 MSE 的代价值远远小于对数似然函数的代价值，阐明 MSE 的代价函数不会重重地惩办谬误分类，哪怕是齐全谬误的分类。

• MSE 作为代价函数不是一个凸函数

  先来看下下面 S 函数导函数的一个性质，上面会用到这个性质：

$$
\begin{aligned}
h(z) &=\frac{1}{1+e^{-z}} \\
\frac{\partial h(z)}{\partial z} &=-\frac{1}{\left(1+e^{-z}\right)^{2}} e^{-z}(-1) \\
&=\frac{e^{-z}}{\left(1+e^{-z}\right)^{2}} \\
&=\frac{1}{1+e^{-z}} \frac{e^{-z}}{1+e^{-z}} \\
&=h(z)(1-h(z))
\end{aligned}
$$

（1）将 MSE 作为代价函数

（2）求一阶导数，因为间接对 w 求导较为简单，所以拆解为两步求导

（3）先对 y^ 求导，前面再对 w 求导，用到了下面说的 S 函数的性质，留神最初还有一个 x

（4）开展括号

（5）求二阶导数，也是拆解为两步求导

（6）先对 y^ 求导，前面再对 w 求导

$$
\begin{aligned}
\operatorname{Cost}(w) &=(y-\hat{y})^{2} & (1)\\
\frac{\partial \operatorname{Cost}(w)}{\partial w} &=\frac{\partial \operatorname{Cost}(w)}{\partial \hat{y}} \frac{\partial \hat{y}}{\partial w} & (2)\\
&=2(y-\hat{y})(-1) \hat{y}(1-\hat{y}) x & (3)\\
&=-2\left(y \hat{y}-\hat{y}^{2}-y \hat{y}^{2}+\hat{y}^{3}\right) x & (4)\\
\frac{\partial^{2} \operatorname{Cost}(w)}{\partial w\partial w^T} &=\frac{\partial \frac{\partial \operatorname{Cost}(w)}{\partial w}}{\partial \hat{y}} \frac{\partial \hat{y}}{\partial w^T} & (5)\\
&=-2\left(y-2 \hat{y}-2 y \hat{y}+3 \hat{y}^{2}\right) \hat{y}(1-\hat{y}) xx^{T} & (6)
\end{aligned}
$$

  同下面证实一样，当初只须要看常数项是否是大于等于 0

（1）当 y = 1 时，带入二阶导函数中

（2）因式分解

（3）能够看到，当 y^ 取不同范畴的时候，常数项不是始终都是大于等于 0 的

$$
\begin{aligned}
g(\hat{y}) &=-2\left(1-4 \hat{y}+3 \hat{y}^{2}\right) \quad(y=1) & (1) \\
&=-2(3 \hat{y}-1)(\hat{y}-1) & (2) \\
& \Rightarrow\left\{\begin{array}{ll}
\hat{y} \in\left[0, \frac{1}{3}\right) & g(\hat{y})<0 \\
\hat{y} \in\left[\frac{1}{3}, 1\right] & g(\hat{y}) \geq 0
\end{array}\right. & (3)
\end{aligned}
$$

（1）当 y = 0 时，带入二阶导函数中

（2）因式分解

（3）能够看到，当 y^ 取不同范畴的时候，常数项也不是始终都是大于等于 0 的

$$
\begin{array}{rlr}
g(\hat{y}) & =-2\left(3 \hat{y}^{2}-2 \hat{y}\right) \quad (y=0) & (1)\\
& =-2 \hat{y}(3 \hat{y}-2) & (2)\\
& \Rightarrow\left\{\begin{aligned}
\hat{y} \in\left[0, \frac{2}{3}\right] & g(\hat{y}) \geq 0 \\
\hat{y} \in\left(\frac{2}{3}, 1\right] & g(\hat{y})<0
\end{aligned}\right. & (3)
\end{array}
$$

  因为应用 MSE 作为代价函数的黑塞矩阵不能保障是半正定的，阐明其代价函数不是一个凸函数。

梯度降落法的降落方向

  先来看下一元的泰勒公式 ¹²，其中 o(…) 示意括号内式子的高阶无穷小，多元泰勒公式就是后面的一、二阶导数改写成雅可比矩阵和黑塞矩阵的模式：

$$
f(x)=f(a)+\frac{f^{\prime}(a)}{1 !}(x-a)+\frac{f^{(2)}(a)}{2 !}(x-a)^{2}+\cdots+\frac{f^{(n)}(a)}{n !}(x-a)^{n}+o\left((x-a)^{n}\right)
$$

  我的指标是每次迭代后的代价函数值比上次都小：

$$
\operatorname{Cost}(w_{n+1}) \lt \operatorname{Cost}(w_{n})
$$

  将对数几率回归的代价函数用多元的一阶泰勒公式近似，雅可比矩阵为代价函数的梯度：

$$
\operatorname{Cost}(w) \approx \operatorname{Cost}(A)+\nabla \operatorname{Cost}(A)^{T}(w-A)
$$

（1）将 w_n 带入 A，w_n+1 带入 w

（2）w_n+1 – w_n = Δw，移项能够看到须要保障梯度乘 Δw 小于 0

（3）当 Δw 为负的梯度时，能够保障梯度乘 Δw 小于 0

（4）这样就失去了梯度降落的降落方向

$$
\begin{aligned}
\operatorname{Cost}\left(w_{n+1}\right) & \approx \operatorname{Cost}\left(w_{n}\right)+\nabla \operatorname{Cost}\left(w_{n}\right)^{T}\left(w_{n+1}-w_{n}\right) & (1) \\
\operatorname{Cost}\left(w_{n+1}\right)-\operatorname{Cost}\left(w_{n}\right) & \approx \nabla \operatorname{Cost}\left(w_{n}\right)^{T} \Delta w<0 & (2) \\
\nabla \operatorname{Cost}\left(w_{n}\right)^{T} \Delta w &=-\nabla \operatorname{Cost}\left(w_{n}\right)^{T} \nabla \operatorname{Cost}\left(w_{n}\right)<0 & (3) \\
\Delta w &=-\nabla \operatorname{Cost}\left(w_{n}\right) & (4)
\end{aligned}
$$

牛顿法的降落方向

  将对数几率回归的代价函数用多元的二阶泰勒公式近似，雅可比矩阵为代价函数的梯度，H 为代价函数的黑塞矩阵：

$$
\operatorname{Cost}(w) \approx \operatorname{Cost}(A)+\nabla \operatorname{Cost}(A)^{T}(w-A)+\frac{1}{2}(w-A)^{T} H(w-A)
$$

（1）将 w_n 带入 A，并对近似函数求梯度，另梯度等于 0 向量，这时获得以后的极小值

（2）这样就失去了牛顿法的降落方向

$$
\begin{aligned}
\nabla \operatorname{Cost}(w) &=\nabla \operatorname{Cost}\left(w_{n}\right)+H\left(w-w_{n}\right)=0 & (1) \\
w &=w_{n}-H^{-1} \nabla \operatorname{Cost}\left(w_{n}\right) & (2)
\end{aligned}
$$

五、代码实现

应用 Python 实现对数几率回归算法（梯度降落法）：

import numpy as np

c_1 = 1e-4
c_2 = 0.9

def cost(X, y, w):
    """
    对数几率回归的代价函数
    args:
        X - 训练数据集
        y - 指标标签值
        w - 权重系数
    return:
        代价函数值
    """
    power = -np.multiply(y, X.dot(w))
    p1 = power[power <= 0]
    p2 = -power[-power < 0]
    # 解决 python 计算 e 的指数幂溢出的问题
    return np.sum(np.log(1 + np.exp(p1))) + np.sum(np.log(1 + np.exp(p2)) - p2)

def dcost(X, y, W):
    """
    对数几率回归的代价函数的梯度
    args:
        X - 训练数据集
        y - 指标标签值
        w - 权重系数
    return:
        代价函数的梯度
    """
    return X.T.dot(np.multiply(-y, 1 / (1 + np.exp(np.multiply(y, X.dot(w))))))

def direction(d):
    """
    更新的方向
    args:
        d - 梯度
    return:
        更新的方向
    """
    return -d

def sufficientDecrease(X, y, w, step):
    """
    判断是否满足充沛降落条件（sufficient decrease condition）args:
        X - 训练数据集
        y - 指标标签值
        w - 权重系数
        step - 步长
    return:
        是否满足充沛降落条件
    """
    d = dcost(X, y, w)
    p = direction(d)
    return cost(X, y, w + step * p) <= cost(X, y, w) + c_1 * step * p.T.dot(d)

def curvature(X, y, w, step):
    """
    判断是否满足曲率条件（curvature condition）args:
        X - 训练数据集
        y - 指标标签值
        w - 权重系数
        step - 步长
    return:
        是否满足曲率条件
    """
    d = dcost(X, y, w)
    p = direction(d)
    return -p.T.dot(dcost(X, y, w + step * p)) <= -c_2 * p.T.dot(d)

def select(step_low, step_high):
    """
    在范畴内抉择一个步长，间接取中值
    args:
        step_low - 步长范畴开始值
        step_high - 步长范畴完结值
    return:
        步长
    """
    return (step_low + step_high) / 2

def lineSearch(X, y, w, step_init, step_max):
    """
    线搜寻步长，使其满足 Wolfe 条件
    args:
        X - 训练数据集
        y - 指标标签值
        w - 权重系数
        step_init - 步长初始值
        step_max - 步长最大值
    return:
        步长
    """
    step_i = step_init
    step_low = step_init
    step_high = step_max
    i = 1
    d = dcost(X, y, w)
    p = direction(d)
    while (True):
        # 不满足充沛降落条件或者前面的代价函数值大于前一个代价函数值
        if (not sufficientDecrease(X, y, w, step_i) or (cost(X, y, w + step_i * p) >= cost(X, y, w + step_low * p) and i > 1)):
            # 将以后步长作为搜寻的右边界
            # return search(X, y, w, step_prev, step_i)
            step_high = step_i
        else:
            # 满足充沛降落条件并且满足曲率条件，即曾经满足了 Wolfe 条件
            if (curvature(X, y, w, step_i)):
                # 间接返回以后步长
                return step_i
            step_low = step_i
        # 抉择下一个步长
        step_i = select(step_low, step_high)
        i = i + 1

def logisticRegressionGd(X, y, max_iter=1000, tol=1e-4, step_init=0, step_max=10):
    """
    对数几率回归，应用梯度降落法（gradient descent）args:
        X - 训练数据集
        y - 指标标签值
        max_iter - 最大迭代次数
        tol - 变动量容忍值
        step_init - 步长初始值
        step_max - 步长最大值
    return:
        w - 权重系数
    """
    # 初始化 w 为零向量
    w = np.zeros(X.shape[1])
    # 开始迭代
    for it in range(max_iter):
        # 计算梯度
        d = dcost(X, y, w)
        # 当梯度足够小时，完结迭代
        if np.linalg.norm(x=d, ord=1) <= tol:
            break
        # 应用线搜寻计算步长 
        step = lineSearch(X, y, w, step_init, step_max)
        # 更新权重系数 w
        w = w + step * direction(d)
    return w

应用 Python 实现对数几率回归算法（牛顿法）：

import numpy as np

c_1 = 1e-4
c_2 = 0.9

def cost(X, y, w):
   """
   对数几率回归的代价函数
   args:
       X - 训练数据集
       y - 指标标签值
       w - 权重系数
   return:
       代价函数值
   """
   power = -np.multiply(y, X.dot(w))
   p1 = power[power <= 0]
   p2 = -power[-power < 0]
   # 解决 python 计算 e 的指数幂溢出的问题
   return np.sum(np.log(1 + np.exp(p1))) + np.sum(np.log(1 + np.exp(p2)) - p2)

def dcost(X, y, w):
   """
   对数几率回归的代价函数的梯度
   args:
       X - 训练数据集
       y - 指标标签值
       w - 权重系数
   return:
       代价函数的梯度
   """
   return X.T.dot(np.multiply(-y, 1 / (1 + np.exp(np.multiply(y, X.dot(w))))))

def ddcost(X, y, w):
   """
   对数几率回归的代价函数的黑塞矩阵
   args:
       X - 训练数据集
       y - 指标标签值
       w - 权重系数
   return:
       代价函数的黑塞矩阵
   """
   exp = np.exp(np.multiply(y, X.dot(w)))
   result = np.multiply(exp, 1 / np.square(1 + exp))
   X_r = np.zeros(X.shape)
   for i in range(X.shape[1]):
       X_r[:, i] = np.multiply(result, X[:, i])
   return X_r.T.dot(X)

def direction(d, H):
   """
   更新的方向
   args:
       d - 梯度
       H - 黑塞矩阵
   return:
       更新的方向
   """
   return - np.linalg.inv(H).dot(d)

def sufficientDecrease(X, y, w, step):
   """
   判断是否满足充沛降落条件（sufficient decrease condition）args:
       X - 训练数据集
       y - 指标标签值
       w - 权重系数
       step - 步长
   return:
       是否满足充沛降落条件
   """
   d = dcost(X, y, w)
   H = ddcost(X, y, w)
   p = direction(d, H)
   return cost(X, y, w + step * p) <= cost(X, y, w) + c_1 * step * p.T.dot(d)

def curvature(X, y, w, step):
   """
   判断是否满足曲率条件（curvature condition）args:
       X - 训练数据集
       y - 指标标签值
       w - 权重系数
       step - 步长
   return:
       是否满足曲率条件
   """
   d = dcost(X, y, w)
   H = ddcost(X, y, w)
   p = direction(d, H)
   return -p.T.dot(dcost(X, y, w + step * p)) <= -c_2 * p.T.dot(d)

def select(step_low, step_high):
   """
   在范畴内抉择一个步长
   args:
       step_low - 步长范畴开始值
       step_high - 步长范畴完结值
   return:
       步长
   """
   return (step_low + step_high) / 2

def lineSearch(X, y, w, step_init, step_max):
   """
   线搜寻步长，使其满足 Wolfe 条件
   args:
       X - 训练数据集
       y - 指标标签值
       w - 权重系数
       step_init - 步长初始值
       step_max - 步长最大值
   return:
       步长
   """
   step_i = step_init
   step_low = step_init
   step_high = step_max
   i = 1
   d = dcost(X, y, w)
   H = ddcost(X, y, w)
   p = direction(d, H)
   while (True):
       # 不满足充沛降落条件或者前面的代价函数值大于前一个代价函数值
       if (not sufficientDecrease(X, y, w, step_i) or (cost(X, y, w + step_i * p) >= cost(X, y, w + step_low * p) and i > 1)):
           # 将以后步长作为搜寻的右边界
           # return search(X, y, w, step_prev, step_i)
           step_high = step_i
       else:
           # 满足充沛降落条件并且满足曲率条件，即曾经满足了 Wolfe 条件
           if (curvature(X, y, w, step_i)):
               # 间接返回以后步长
               return step_i
           step_low = step_i
       # 抉择下一个步长
       step_i = select(step_low, step_high)
       i = i + 1

def logisticRegressionNewton(X, y, max_iter=1000, tol=1e-4, step_init=0, step_max=10):
   """对数几率回归，应用牛顿法（newton's method）args:
       X - 训练数据集
       y - 指标标签值
       max_iter - 最大迭代次数
       tol - 变动量容忍值
       step_init - 步长初始值
       step_max - 步长最大值
   return:
       w - 权重系数
   """
   # 初始化 w 为零向量
   w = np.zeros(X.shape[1])
   # 开始迭代
   for it in range(max_iter):
       # 计算梯度
       d = dcost(X, y, w)
       # 计算黑塞矩阵
       H = ddcost(X, y, w)
       # 当梯度足够小时，完结迭代
       if np.linalg.norm(d) <= tol:
           break
       # 应用线搜寻计算步长 
       step = lineSearch(X, y, w, step_init, step_max)
       # 更新权重系数 w
       w = w + step * direction(d, H)
   return w

六、第三方库实现

scikit-learn¹³ 实现无正则化的对数几率回归：

from sklearn.linear_model import LogisticRegression

# 初始化对数几率回归器，无正则化
reg = LogisticRegression(penalty="none")
# 拟合线性模型
reg.fit(X, y)
# 权重系数
w = reg.coef_
# 截距
b = reg.intercept_

scikit-learn¹³ 实现 L1 正则化的对数几率回归：

from sklearn.linear_model import LogisticRegression

# 初始化对数几率回归器，L1 正则化，应用坐标降落法
reg = LogisticRegression(penalty="l1", C=10, solver="liblinear")
# 拟合线性模型
reg.fit(X, y)
# 权重系数
w = reg.coef_
# 截距
b = reg.intercept_

scikit-learn¹³ 实现 L2 正则化的对数几率回归：

from sklearn.linear_model import LogisticRegression

# 初始化对数几率回归器，L2 正则化
reg = LogisticRegression(penalty="l2", C=10)
# 拟合线性模型
reg.fit(X, y)
# 权重系数
w = reg.coef_
# 截距
b = reg.intercept_

scikit-learn¹³ 实现弹性网络正则化的对数几率回归：

from sklearn.linear_model import LogisticRegression

# 初始化对数几率回归器，弹性网络正则化
reg = LogisticRegression(penalty="elasticnet", C=10, l1_ratio=0.5, solver="saga")
# 拟合线性模型
reg.fit(X, y)
# 权重系数
w = reg.coef_
# 截距
b = reg.intercept_

scikit-learn 每种优化办法所反对的正则项：

七、动画演示

  下图展现了演示数据，其中红色示意标签值为 1、蓝色示意标签值为 -1：

  下图为应用梯度降落法拟合数据的后果，其中浅红色示意拟合后依据权重系数计算出预测值为 1 的局部，浅蓝色示意拟合后依据权重系数计算出预测值为 - 1 的局部：

  下图比照了不同的正则化形式对后果的影响，正则化惩办项系数 C = 0.01，弹性网络中 ρ = 0.5。

八、思维导图

九、参考文献

1. https://en.wikipedia.org/wiki…
2. https://en.wikipedia.org/wiki…
3. https://en.wikipedia.org/wiki…
4. https://en.wikipedia.org/wiki…
5. https://en.wikipedia.org/wiki…
6. https://en.wikipedia.org/wiki…
7. https://en.wikipedia.org/wiki…
8. https://en.wikipedia.org/wiki…
9. https://en.wikipedia.org/wiki…
10. https://en.wikipedia.org/wiki…
11. https://towardsdatascience.co…
12. https://en.wikipedia.org/wiki…
13. https://scikit-learn.org/stab…

残缺演示请点击这里

注：本文力求精确并通俗易懂，但因为笔者也是初学者，程度无限，如文中存在谬误或脱漏之处，恳请读者通过留言的形式批评指正

本文首发于——AI 导图 ，欢送关注