一切从逻辑回归开始

JSong @2016.06.13

本系列文章不适合入门，是作者综合各方资源和个人理解而得. 另外最好有数学基础, 因为数学人一言不合就会上公式.

简单模型的魅力在于它能从各个角度去欣赏. 逻辑回归是最简单的二分类模型之一，实际应用中二分类最常见，如判定是否是垃圾邮件，是否是人脸，是否值得借贷等, 而概率模型对于这类问题有得天独厚的优势. 本文将从各个角度来理解逻辑回归，并指出它是一个概率模型、对数线性模型、交叉熵模型.

我们先从线性回归谈起。 考察 m 个变量和 y 之间的线性关系：

\[y\sim x_1+x_2+\cdots +x_m \]

根据要求我们需要找到 m+1 个回归系数 θ_i，使得

\[\min ||y-\hat{y}||_{F},\quad \hat{y}=h_{\theta}(x)=\sum_{i=0}^m \theta_i x_i=\theta^T x \]

其中 x_0=1.

1. 线性代数下的回归模型

如果每个变量都有 n 个样本，即 x_i ∈ R^n ，则上述问题等价于求解线性方程组：

\[X\theta =y, \quad X=(x_0,x_1,\ldots,x_m) \]

一般来讲，上述方程有唯一解(最小二乘法)：

\[\theta=X^{+}y=(X^TX)^{-1}X^Ty \]

讲到这，车就可以开始上路了。

2. 二分类问题简介

分类模型比较常见，例如是否是垃圾邮件等等。而最简单的分类就是二分类，且多分类一定程度也可以转化成二分类问题。 给定一组有标签的训练数据

\[(X,y)=\{(x^{i}\in \mathbb{R}^n, y^{i}\in \{0,1\}): i=1,2,\ldots,m\} \]

在下图中我们随机给了一组散点图，红色代表y=1，蓝色代表y=0.

我们的目的就是找到一种模型将这两类样本点分开(当然实际应用中可能没有这么好的情况。 这里只是给一个Demo)。 最直接的想法就是找一个超平面

\[\theta^T x=0 \]

使得红色点和蓝色点分别分布于超平面的两端. 此时对于任意一个样本点 x^{(i)}, 不失一般性

\[\theta^T x^{(i)} \]

可以代表样本点到该超平面的距离。 因为红色样本点在超平面上方，即红点对应距离为正数，同理蓝点对应的距离为负数。 为避免区分红点或者蓝点，我们可以用

\[(2y-1)\theta^T x \]

来代替， 此时只要被正确分类， 该距离则为正数，否则为负数。 于是可得相应的损失函数为

\[\max(0,(2y-1)\theta^T x) \]

综合一下可得最简单的线性分类模型为：

\[\arg_\theta \min \sum_{i=1}^{m} \max(0,(2y_i-1)\theta^T x^{(i)}) \]

这个损失函数有很多不好的地方，而且也不可导。接下来我们介绍更好的逻辑回归模型，其背后有很多的解释。

2.1 逻辑回归模型

令

\[\theta^T x=\theta_0 +\sum_{i=1}^{m} \theta_i x_i \]

且我们假定

\[P(y=1|x; \theta)=h_{\theta}(x)=g(\theta^T x)=\frac{1}{1+e^{-\theta^T x}} \]

\[P(y=0|x; \theta)=1-h_{\theta}(x) \]

其中

\[g(z)=\frac{1}{1+e^{-z}} \]

是 logistic 分布函数，它的图像和密度函数

\[e^{-z}/(1+e^{-z})^2 \]

图像见下图. 原则上任何连续的分布函数都可以，其把实数映射到0到1之间. 至于为啥选择logistic函数, 我们之后再讨论.

这里我们采用极大似然估计. 由于y取值的特殊性, 上面的模型假设等价于

\[P(y|x;\theta)=(h_{\theta}(x))^y \,(1-h_{\theta}(x))^{1-y} \]

这样我们便有似然函数：

\[l(\theta)=\log L(\theta)=\log \prod_{i=1}^{m} p(y^{(i)}|x^{(i)};\theta) \]

\[=\sum_{i=1}^m [y^{(i)} \log h_{\theta}(x^{(i)})+ (1-y^{(i)}) \log (1-h_{\theta}(x^{(i)})] \]

于是得到逻辑回归模型为

\[\arg_{\theta} \max l(\theta) \]

2.2 几何意义（超平面分类）

重新考虑模型假设, 我们可以得到

\[\theta^T x=\log \frac{p}{1-p}, \quad p=P(y=1|x) \]

同样，逻辑回归模型也是用超平面来分类的. 且相应的极大似然估计可以写成

\[\min -l(\theta)= \sum_{i=1}^{m} \log[(2y-1)(y-h_{\theta}(x))]=\sum_{i=1}^{m} \log (1+e^{(2y-1)\theta^T x}) \]

也即相应的损失函数为

\[\log(1+\exp((2y-1)\theta^T x)) \]

图像见下方. 由统计假设可知，给定一个样本点，如果它被正确分类，则

\[(2y-1)(y-h_{\theta}(x)) \]

应该越小越好, 如y=1时,h_θ(x)越接近1越好. 相应的 (2y−1)θ^Tx 也是越小越好. 又根据第一节的讨论我们有, 当某样本点被正确分类时

\[(2y-1)\theta^T x= \mbox{distance of sample point and hyperplane}\,\, \theta^T x=0 \]

此处

\[\log(1+\exp(x)) \quad and \quad \max(0,x) \]

的作用类似，但连续性更好. 另外有些文章会写逻辑回归的损失函数为

\[\log(1+\exp(-y\cdot \theta^T x)) \]

这是因为该损失函数中的 y 取值范围为-1和1.

2.3 线性代数意义

在上一节，我们给出了表达式

\[\theta^T x=\log \frac{p}{1-p}, \quad p=P(y=1|x) \]

可以看出它相当于把 y 映射到了[0,1]区间。而且 p/(1-p) 是事件发生与事件不发生的概率之比，称为事件的发生比 (the odds of experiencing an event), 简称为 odds 。

2.4 信息论解释(交叉熵模型)

令随机变量

\[p\in \{y,1-y\}, \quad q \in\{\hat{y},1-\hat{y}\} \]

其中 \hat{y} 为模型拟合得到的y, 则它们之间的交叉熵 (cross entropy) 为

\[H(p,q)= -\sum_i p_i \log q_i =-[y\log \hat{y}+(1-y)\log(1-\hat{y})] \]

我们知道熵常用于度量一个随机变量所包含的信息量， 而交叉熵可以用来度量两个随机变量之间的相似性, 从上面式子可以看出逻辑回归的极大似然等价于最小化交叉熵.

事实上真正等价的是交叉熵与极大似然估计，有兴趣的同学可以自己证明。

2.5 神经网络模型

这个就不多说了，逻辑回归是一个最简单的神经网络模型

2.6 梯度下降法参数求解

按照极大似然的来求导：

\[\bigtriangledown_{\theta} l(\theta)=\sum_{i=1}^{m} (y^{(i)}-h_{\theta}(x^{(i)}))x^{(i)}= X^T[Y-h_{\theta}(X)] \]

其中 x^{(i)} 在 X 中作为行向量存储. 于是可得参数 θ 的迭代式：

\[\theta:=\theta+\alpha X^T[Y-h_{\theta}(X)] \]

其中 α 为步长, 另外因为是最大化 l(θ) ，所以是沿着梯度的正方向寻找。

1.7 python代码

import numpy as np
import matplotlib.pyplot as plt
import time
%pylab inline

# calculate the sigmoid function
def sigmoid(inX):
    return 1.0 / (1 + np.exp(-inX))



def trainLogRegres(train_x, train_y, opts):
    # train a logistic regression model using some optional optimize algorithm
    # train_x is a mat datatype, each row stands for one sample
    # train_y is mat datatype too, each row is the corresponding label
    # opts is optimize option include step and maximum number of iterations

    # calculate training time
    startTime = time.time()
    train_x=np.asmatrix(train_x)
    train_y=np.asmatrix(train_y)
    numSamples, numFeatures =train_x.shape
    alpha = opts['alpha']; maxIter = opts['maxIter']
    weights = np.ones((numFeatures, 1))
    for k in range(maxIter):
        err = train_y - sigmoid(train_x * weights)
        weights = weights + alpha * train_x.T * err
    print 'Congratulations, training complete! Took %fs!' % (time.time() - startTime)
    return weights


# show your trained logistic regression model only available with 2-D data
def showLogRegres(weights, train_x, train_y):
    # notice: train_x and train_y is mat datatype
    numSamples, numFeatures = train_x.shape
    if numFeatures != 3:
        print "Sorry! I can not draw because the dimension of your data is not 2!"
        return 1

    # draw all samples   
    idx1=np.asarray(train_y)==1
    idx2=np.asarray(train_y)==0
    plt.plot(x[:,1][idx1].T,x[:,2][idx1].T,'or')
    plt.plot(x[:,1][idx2].T,x[:,2][idx2].T,'ob')
    # draw the classify line
    min_x = min(train_x[:, 1])[0, 0]
    max_x = max(train_x[:, 1])[0, 0]
    weights = weights.getA()  # convert mat to array
    y_min_x = float(-weights[0] - weights[1] * min_x) / weights[2]
    y_max_x = float(-weights[0] - weights[1] * max_x) / weights[2]
    plt.plot([min_x, max_x], [y_min_x, y_max_x], '-g')
    plt.xlabel('X1'); plt.ylabel('X2')

# 测试样例
len_samples=500
x=np.random.randn(len_samples,2)
idx1=x[:,0]+x[:,1]>0.5
y=np.zeros((len_samples,1))
y[idx1,0]=1
opt={'alpha':0.1,'maxIter':1000}
x=np.hstack((np.ones((len_samples,1)),x))
x=np.asmatrix(x)
y=np.asmatrix(y)
w=trainLogRegres(x,y,opt)
print w
showLogRegres(w,x,y)

输出如下：

Populating the interactive namespace from numpy and matplotlib
Congratulations, training complete! Took 0.047000s!
[[-17.37216188]
 [ 34.35470115]
 [ 34.98381818]]

这只是最简单的样例，在实际应用中我们还可以添加惩罚项

\[\min_{w,c} \|w\|_2+ C\sum_{i=1}^{m} \log(1+\exp(-y^{(i)}(w\cdot x^{(i)}+c))) \]

或者 L_1 范数

\[\min_{w,c} \|w\|_1+ C\sum_{i=1}^{m} \log(1+\exp(-y^{(i)}(w\cdot x^{(i)}+c))) \]

在样本很稀疏的时候，惩罚项很有用. 另外当样本量很大，我们又极其要求速度的时候，梯度下降法也要改进，换成SGD（stochastic gradient descent ）、拟牛顿法、AGD等等

专注写文二十年，原文首发于公众号，欢迎关注