机器学习笔记（二十）：逻辑回归(2)

凌云时刻 · 技术

导读：在这一篇笔记中我们继续来学习分类算法逻辑回归(Logistic Regression)。

作者 | 计缘

来源 | 凌云时刻（微信号：linuxpk）

定义损失函数

回顾完对数后，我们再来看看逻辑回归的损失函数，没错，我们就是要使用对数的函数来表示：

  如果 ，那么损失函数为： 如果 ，那么损失函数为：  

我们来看看这两个函数为什么符合之前我们定义的损失函数的描述‍‍。首先   的曲‍‍线为：

 

之前说了，概率的值域‍‍在(0,1)(0,1)之间，所以上图曲线的x轴以下的曲线是没有意义的，所以‍‍对于   ，‍‍它的曲线是：‍‍

从上面的曲线图可以很容易的发现，‍‍‍‍‍‍‍‍‍‍当   趋近于0的时候，‍‍   趋‍‍近于正无穷，这个正无穷其实就是一个很大的损失惩罚，因‍‍为当   趋近于0时‍‍，   应该也是趋‍‍近于0，但是这‍‍里   ‍‍定义的是1。‍‍当   在‍‍不断趋近1的过程中‍‍，   的‍‍值逐渐减小，既损失逐渐减小，‍‍当   趋‍‍近于1时‍‍，   应‍‍该也是趋近于0，和这里定义的是一致的，所以   的‍‍值是0，说明没有损失。‍‍

下面再来看‍‍看   的‍‍曲线：‍‍‍‍‍‍‍‍‍‍

同理因为概率的‍‍值域在(0,1)之间，‍‍所以   的‍‍曲线‍‍为：

这条曲线同样可以解释我们之前定义的损失趋势。

此时我们找到的损失函数还是根据不同的分类分成了两个，其实将其合成一个也很简单：

   

如此一来‍‍，当   时，损失函‍‍数为   ，‍‍‍‍当   时，损‍‍失函数为   。‍‍‍‍

上面的公式，是针对一个样本数据的，那么如果有多个样本数据，其实就是将这些样本数据的损失值加起来然后在求一下平均值：

   

   

下面我们要做的就是找‍‍到一组   值，使得上面‍‍的   达到最‍‍‍‍小值。

损失函数的梯度

上面的公式是没法像线性回归那样求出一个正规方程解的，所以我们需要使用梯度下降法来求得‍‍使   最小的一‍‍组   。‍‍下面我们先把公式都列出来：

• 多元线性回归公式：   ，‍‍注意这里的   ‍‍是加上了值全部为1的一列的矩阵，而为了方便推导，这里‍‍的   是一‍‍个列向量，就‍‍不写成   了。

• Sigmoid函数：   。

• 逻辑回归概率公式：   

• 逻辑回归损失函数：

     

在第五篇笔记中我们知道，求损失函数的梯度就是对   这个列向量逐个元素求导：

   

我们从里往外来看，先从Sigmoid函数求导入手。

 Sigmoid函数求导

先变换一下Sigmoid函数：

   

然后对Sigmoid函数求导，这里遵循求导链式法则以及求导基本法则：

• 复合函数   的导数   为：   

• 代数函数导数：   

• 数学常数的指数求导还是它自己：   

所以可得Sigmoid函数的导数为：

   

 Sigmoid函数的对数求导

下面再往外扩展，来看一下   的导数。这里遵循的导数法则为：

• 对以2为底的对数求导：   

• 复合函数    的导数   为：   

所以   ‍‍的导数为：

   

   

   

   

   

log(σ(t))‘=1σ(t)⋅σ(t)‘=1(1+et)−1⋅(1+e−t)−2⋅e−t=(1+e−t)−1⋅e−t=e−t1+e−t=1+e−t−11+e−t=1−11+e−t=1−σ(t)

‍‍   的导数为：

   

   

   

   

 逻辑回归损失函数求导

当我们知道了Sigmoid函数和Sigmoid函数的对数的求导结果后，我们对逻辑回归损失函数求导就很容易了（这里对‍‍第   ‍‍个   ‍‍求‍‍导），先来看前半部分：

   

   

最后的那‍‍个   是第   个   ‍‍‍‍‍‍‍‍前‍‍‍‍面的系数，再来看后半部分：

   

   

‍‍

此时对整个损失函数求导就是上面两部分相加：

    

   

 

   

   

所以逻辑回归损失函数的梯度为：

   

大家再来回顾一下第五篇笔记中线性回归的梯度：

   

可以发现这两个梯度在形态上是非常相似的，因为线性回归的梯度通过向量化可以优化为：

   

所以逻辑回归的梯度最终可以写为：

   

实现逻辑回归算法

因为逻辑回归拟合损失函数使用的是梯度下降法，所以我们封装逻辑回归算法时大部分都可以套用我们之前封装的线性回归梯度下降方法，需要修改的只是损失函数、预测和评分里的一些代码。

 

import numpy as np
from .metrics import accuracy_score

class LogisticRegression:

	def __init__(self):
		# 截距theta0
		self.intercept_ = None
		# 系数，theta1 ... thetaN
		self.coef_ = None
		# theta列向量
		self._theta = None

	# 定义Sigmoid私有函数
	def _sigmoid(self, t):
		return 1. / (1. + np.exp(-t))

	# 使用批量梯度下降法，根据训练数据集X_train，y_train训练LogisticRegression模型
	def fit(self, X_train, y_train, is_debug=False, eta=0.01, n_iters=1e4):
		assert X_train.shape[0] == y_train.shape[0], \
			"特征数据矩阵的行数要等于样本结果数据的行数"

		# 定义逻辑回归损失函数
		def L(theta, X_b, y):
			# 定义逻辑回归概率公式
			y_hat = self._sigmoid(X_train.dot(theta))

			try:
				return -np.sum(y*np.log(y_hat)+(1-y)*np.log(1-y_hat)) / len(X_b)
			except:
				return float('inf')

		# 定义逻辑回归梯度
		def dL(theta, X_b, y):
			return X_b.T.dot(self._sigmoid(X_b.dot(theta)) - y) / len(X_b)

		def dL_debug(theta, X_b, y, epsilon=0.01):
			# 开辟大小与theta向量一致的向量空间
			result = np.empty(len(theta))
			# 便利theta向量中的每一个theta
			for i in range(len(theta)):
				# 复制一份theta向量
				theta_1 = theta.copy()
				# 将第i个theta加上一个距离，既求该theta正方向的theta
				theta_1[i] += epsilon
				# 在复制一份theta向量
				theta_2 = theta.copy()
				# 将第i个theta减去同样的距离，既求该theta负方向的theta
				theta_2[i] -= epsilon
				# 求出这两个点连线的斜率，既模拟该theta的导数
				result[i] = (L(theta_1, X_b, y) - L(theta_2, X_b, y)) / (2 * epsilon)
			return result

		# 实现批量梯度下降法
		def gradient_descent(X_b, y, initial_theta, eta, difference=1e-8):
			theta = initial_theta
			i_iter = 0
			while i_iter < n_iters:
				# 当is_debug为True时走debug的求梯度的方法，反之走梯度公式的方法
				if is_debug:
					gradient = dL_debug(theta, X_b, y)
				else:
					gradient = dL(theta, X_b, y)
				last_theta = theta
				theta = theta - eta * gradient

				if (abs(L(theta, X_b, y) - L(last_theta, X_b, y)) < difference):
					break

				i_iter += 1
			return theta

		# 构建X_b
		X_b = np.hstack([np.ones((len(X_train), 1)), X_train])
		# 初始化theta向量为元素全为0的向量
		initial_theta = np.zeros(X_b.shape[1])

		self._theta = gradient_descent(X_b, y_train, initial_theta, eta)
		self.intercept_ = self._theta[0]
		self.coef_ = self._theta[1:]

		return self

	# 计算概率，给定待预测数据集X_predict，返回表示X_predict的结果概率向量
	def predict_probability(self, X_predict):
		assert self.intercept_ is not None and self.coef_ is not None, \
		"截距和系数都不为空，表示已经经过了fit方法"
		assert X_predict.shape[1] == len(self.coef_), \
		"要预测的特征数据集列数要与theta的系数数量相等"

		X_b = np.hstack([np.ones((len(X_predict), 1)), X_predict])

		# 返回0，1之间的浮点数
		return self._sigmoid(X_b.dot(self._theta))

	# 给定待预测数据集X_predict，返回表示X_predict的结果向量
	def predict(self, X_predict):
		assert self.intercept_ is not None and self.coef_ is not None, \
		"截距和系数都不为空，表示已经经过了fit方法"
		assert X_predict.shape[1] == len(self.coef_), \
		"要预测的特征数据集列数要与theta的系数数量相等"

		probability = self.predict_probability(X_predict)
		# 将概率转换为0和1的向量，True对应1，False对应0
		return np.array(probability >= 0.5, dtype='int')

	# 根据测试数据集X_test和y_test确定当前模型的准确度
	def score(self, X_test, y_test):
		y_predict = self.predict(X_test)
		return accuracy_score(y_test, y_predict)

	def __repr__(self):
		return "LinearRegression()"

下面我们在Jupyter Notebook中使用Scikit Learn提供的鸢尾花数据验证我们封装的逻辑回归的方法：

 

import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets

iris = datasets.load_iris()
X = iris.data
y = iris.target

 

因为鸢尾花数据中有三类鸢尾花，而逻辑回归在一开始就讲过是一个解决二分类问题的算法，所以我们只取前两类的鸢尾花数据来验证，并且只用每类鸢尾花的前两个特征，方便绘图：

 

# 只取前两类的鸢尾花数据
X = X[y