跳到主要内容

逻辑回归

逻辑回归 (Logistic Regression) 是一种非常容易实现的分类模型,对于线性可分的情况有很好的表现。与之前文章中的感知器和 Adaline 类似,本文中的逻辑回归也是二元分类的线性模型。

算法思想

为了解释逻辑回归模型背后的思想,先介绍一个概念:比值比 (Odds ratio,OR)。

比值比的数学表示为 p1p\cfrac p {1-p},其中 pp 为某事件可能发生的几率。

在此基础上,定义 logit 函数,即比值比的对数 (计算机科学中 log\log 表示自然对数):

logit(p)=logp1plogit(p) = \log {\cfrac p {1-p}}

logit 的作用是将范围为 (0,1)(0,1) 的变量,映射到整个实数集。

对于给定的 xx

  • y=1y=1 的概率 p(y=1x)p(y=1 | x)(0,1)(0,1) 之间。
  • 基于神经元模型,输出的概率是由净输入 z=i=0mwixi=wTxz = \sum_{i=0}^m w_i x_i = w^T x 得到的。

基于以上逻辑,建立 p(y=1x)zp(y=1 \mid x) \rightarrow zlogit 映射:

logit(p(y=1x))=z\begin{align} logit(p(y=1 \mid x)) = z \end{align}

因此,logit 的反函数即为 zp(y=1x)z \rightarrow p(y=1 \mid x) 的映射:

ϕ(z)=11+ez=p(y=1x)\begin{align} \phi (z) = \cfrac 1 {1+e^{-z}} = p(y=1 \mid x) \end{align}

这里 ϕ(z)=11+ez\phi (z) = \cfrac 1 {1+e^{-z}} 称为 logistic sigmoid 函数,通常简称为 sigmoid 函数。(sigmoid 意为 S-shaped,即 S 形状的曲线。)

import matplotlib.pyplot as plt
import numpy as np


def sigmoid(z):
return 1.0 / (1.0 + np.exp(-z))


z = np.arange(-7, 7, 0.1)
phi_z = sigmoid(z)

plt.plot(z, phi_z)
plt.axvline(0.0, color="k")
plt.ylim(-0.1, 1.1)
plt.xlabel("z")
plt.ylabel("$\phi (z)$")
plt.yticks([0.0, 0.5, 1.0])
ax = plt.gca()
ax.yaxis.grid(True)
plt.show()

算法模型

Adaline 和逻辑回归的差异,如下图所示:

Adaline:

  • 激活函数:线性函数 ϕ(z)=z\phi (z) = z
  • 阈值函数:阶跃函数 y^={1if ϕ(z)00otherwise\hat y = \begin{cases} 1 & if \ \phi (z) \ge 0 \\[2ex] 0 & otherwise \end{cases}

逻辑回归:

  • 激活函数:sigmoid 函数 ϕ(z)=11+ez\phi (z) = \cfrac 1 {1+e^{-z}}
  • 阈值函数:阶跃函数 y^={1if ϕ(z)0.50otherwise\hat y = \begin{cases} 1 & if \ \phi (z) \ge 0.5 \\[2ex] 0 & otherwise \end{cases}

需要注意的是,逻辑回归中的条件 ϕ(z)0.5\phi (z) \ge 0.5,等价于 z0z \ge 0

损失函数和成本函数

定义损失函数 (Loss Function) 和成本函数 (Cost Function),以训练权重参数 ww

损失函数用于单个样本,代表预测输出 y^(i)\hat y^{(i)} 和实际输出 y(i)y^{(i)} 之间的误差。

L(y^(i),y(i))=[y(i)log(y^(i))+(1y(i))log(1y^(i))]={log(y^(i))if y(i)=1log(1y^(i))if y(i)=0\begin{align} L(\hat y^{(i)}, y^{(i)}) & = -[y^{(i)} \log (\hat y^{(i)}) + (1 - y^{(i)}) \log (1 - \hat y^{(i)})] \\ & = \begin{cases} -\log (\hat y^{(i)}) & \text{if} \ y^{(i)} = 1 \\[2ex] -\log (1 - \hat y^{(i)}) & \text{if} \ y^{(i)} = 0\end{cases} \end{align}
  • y(i)=1y^{(i)} = 1 时,最小化 LL 等价于 y^(i)\hat y^{(i)} 接近于 1。
  • y(i)=0y^{(i)} = 0 时,最小化 LL 等价于 y^(i)\hat y^{(i)} 接近于 0。

成本函数是整个训练集的损失函数的平均值。最小化成本函数,可以得到最优的参数 ww

J(w)=1mi=1mL(y^(i),y(i))=1mi=1m[y(i)log(y^(i))+(1y(i))log(1y^(i))]\begin{align} J(w) & = {\frac 1 m} \sum _{i=1}^m L(\hat y^{(i)}, y^{(i)}) \\ & = -{\frac 1 m} \sum _{i=1}^m [y^{(i)} \log (\hat y^{(i)}) + (1 - y^{(i)}) \log (1 - \hat y^{(i)})] \end{align}

Python 实现

逻辑回归的实现如下,使用梯度下降法进行参数训练:

import numpy as np


class LogisticRegressionGD(object):
def __init__(self, eta=0.05, n_iter=100, random_state=1):
self.eta = eta
self.n_iter = n_iter
self.random_state = random_state

def fit(self, X, y):
rgen = np.random.RandomState(self.random_state)
self.w_ = rgen.normal(loc=0.0, scale=0.01, size=1 + X.shape[1])
self.cost_ = []
for i in range(self.n_iter):
net_input = self.net_input(X)
output = self.activation(net_input)
errors = y - output
self.w_[1:] += self.eta * X.T.dot(errors)
self.w_[0] += self.eta * errors.sum()
cost = -(y.dot(np.log(output)) + ((1 - y).dot(np.log(1 - output)))) / X.shape[0]
self.cost_.append(cost)
return self

def net_input(self, X):
return np.dot(X, self.w_[1:]) + self.w_[0]

def activation(self, z):
return 1.0 / (1.0 + np.exp(-np.clip(z, -250, 250)))

def predict(self, X):
return np.where(self.net_input(X) >= 0.0, 1, 0)

绘制图像

Iris-setosa (标签 0) and Iris-versicolor (标签 1) 为例,采用逻辑回归分类,绘制决策边界如下,

import matplotlib.pyplot as plt
from plot_decision_regions import plot_decision_regions
from sklearn import datasets
from sklearn.model_selection import train_test_split

iris = datasets.load_iris()
X = iris.data[:, [2, 3]]
y = iris.target

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=1, stratify=y)

X_train_01_subset = X_train[(y_train == 0) | (y_train == 1)]
y_train_01_subset = y_train[(y_train == 0) | (y_train == 1)]

lrgd = LogisticRegressionGD(eta=0.05,n_iter=1000,random_state=1)
lrgd.fit(X_train_01_subset,y_train_01_subset)

plot_decision_regions(X=X_train_01_subset,y=y_train_01_subset,classifier=lrgd)
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
plt.show()