机器学习笔记（五）—— 逻辑回归

逻辑回归算法是二分类问题中最常用的几种分类算法之一，通过变形，也能够在多分类问题中发挥余热。今天我将从向大家揭开这个简单算法的神秘面纱！

一、Sigmoid函数

在回归问题中，我们曾经提到，对于数据集$X$,我们可以找到合适的系数$W$使其通过$WX$来预测结果$Y$。

逻辑回归的思路与回归算法一致，我们需要找到一个合适的系数$W$，通过$WX$来得到一个结果$Y$，然后对$Y$进行判断，从而得到分类结果。因此，如果来判断$Y$便是逻辑回归最重要的一个点。

Sigmoid函数的$y$值分布在$[0,1]$之间，它的函数表达式为$y=\frac{1}{1+e^{-x}}$，图像如下图所示。

也就是说，对于n维数据集$X$，如果我能够找到合适的系数$W$，并将$WX$代入Sigmoid函数中，使得结果$Y=sigmoid(WX)$在$[0,1]$范围内分布，大小代表了样本为正样本的概率。例如$Y(X)=0.8$,说明该样本有80%的可能性为正样本，有20%的可能性为负样本。通过这种方法，我们就能够对样本集进行二分类操作。

二、逻辑回归的损失函数

当我们确定好逻辑回归的模型$$y=\frac{1}{1+e^{-WX}}$$后，我们接下来的任务就是如何去评价这个模型的好坏，也就是找打一个损失函数，来表示预测的输出$Y_{predict}$和训练数据类别$Y_{true}$之间的区别。

2.1 最大似然估计

假设对于给定的数据集$$D=[(X^{(1)},y^{(1)}),.....,(X^{(m)},y^{(m)})]$$
预测函数为：$$h(X)=\frac{1}{1+e^{-WX}}$$
则对对于输入X其分类结果为1和为0的概率分别为：
$$P(y=1|X;W)=h(X)$$
$$P(y=0|X;W)=1-h(X)$$（注：默认类别1为正类别）
即：$$P(y|X;W)=(h(X){)}^y(1-h(X){)}^{1-y}$$
那么基于最大似然估计，得到似然函数：
$$L(W)=\prod _{i=1}^{m}P(y^{(i)}|X^{(i)};W)=(h(X^{(i)}){)}^{y^{(i)}}(1-h(X^{(i)}){)}^{1-y^{(i)}}$$
对数似然函数则为：
$$l(W)=logL(W)=\sum _{i=1}^m(y^{(i)}logh(X^{(i)})+(1-y^{(i)})log(1-h(X^{(i)}))$$
也就是说，如果我们希望$Y_{predict}$与$Y_{true}$之间的差别最小，那么我们需要使最大对数似然函数$l(W)$取最大值，此时我们所求得的$W$就是我们所需要的系数。

所以我们的损失函数为：
$$J(W)=-l(W)=0=-\sum _{i=1}^m(y^{(i)}logh(X^{(i)})+(1-y^{(i)})log(1-h(X^{(i)}))$$
接下来的任务就是如何求解损失函数最小时的$W$值。

三、求解系数W

3.1 梯度下降法

其实我们可以注意到，损失函数与对数似然函数只相差一个负号，因此如果使对损失函数求最小值，就用梯度下降法； 如果是对对数似然函数求最大值，就用梯度上升法，并没有什么本质上的区别！

（1）确定梯度

已知：$$h(X)=g(WX)=\frac{1}{1+e^{-WX}}$$ $$J(W)=-l(W)=0=-\sum _{i=1}^m(y^{(i)}logh(X^{(i)})+(1-y^{(i)})log(1-h(X^{(i)}))$$
现在求解$J(W)$对$W$偏导：
$$\frac{\partial J(W)}{\partial W_j}=\frac{\partial J(W)}{\partial g(WX)}\ast \frac{\partial g(WX)}{\partial WX}\ast \frac{\partial WX}{\partial W_j}$$
其中：
$$\frac{\partial WX}{\partial W_j}=\frac{\partial (W_1{x}_1+....+W_n{x}_n)}{\partial W_j}=x_j$$
$$\frac{\partial J(W)}{\partial g(WX)}=y\ast \frac{1}{g(WX)}+(y-1)\ast \frac{1}{1-g(WX)}$$
$$\frac{\partial g(WX)}{\partial WX}=g(WX)(1-g(WX))$$
请参考：
$$f(x)=\frac{1}{1+e^{-x}}$$
$$f^{\prime }x=\frac{e^{-x}}{(1+e^{-x}{)}^2}=\frac{1}{1+e^{-x}}(1-\frac{1}{1+e^{-x}})=f(x)(1-f(x))$$
综上所述：
$$\frac{\partial J(W)}{\partial W_j}=-(y-h(X))x_j$$
因此，梯度下降法的迭代公式为：
$$W_j=W_j+\alpha \frac{1}{m}\sum _{i=1}^m(y^{(i)}-h(X^{(i)}))x_j^{(i)}=W_j-\frac{\alpha }{m}\sum _{i=1}^m(h(X^{(i)})-y^{(i)})x_j^{(i)}$$

3.2 向量化

在上面的梯度下降法的迭代公式$W_j=W_j+\alpha \frac{1}{m}{\sum }_{i=1}^m(y^{(i)}-h(X^{(i)}))x_j^{(i)}$中，存在一个求和的过程，显然如果通过for循环来处理，将会存在大量的计算，通过向量化来简化运算过程是我们优化算法的一条重要思路。
通过Numpy或者Pandas处理，我们使得
$$X=[x^{(1)},x^{(2)},......,x^{(m)}{]}^T$$
即$X$的每一行为一组训练样本，每一列代表不同特征的取值
令Y为:
$$Y=[y^{(1)},y^{(2)},......,y^{(m)}{]}^T$$
令待求的系数$W$的矩阵形式为:
$$W=[{\omega }_1,{\omega }_2,........,{\omega }_n{]}^T$$
则$$W=W-\frac{\alpha }{m}{X}^T(g(XW)-Y)$$
上式就是我们的梯度下降法的向量化表达。

四、Python实现

import numpy as np   #导入必要的包

class LogisticRegression(object):
	def __init__(self):  
		self.coef_ = None;
		self.inter_ = None;
		self._theta = None;
	
	def _Sigmoid(self,x):   #Sigmoid函数
		return 1./(1.+np.exp(-x));
	
	def fit(self,X_train,y_train,eta = 0.01,n_iters=1e4): #训练过程
		def J(theta,X_b,y): # 损失函数
			y_hat = self._Sigmoid(X_b.dot(theta));
			try:
				return -np.sum(y*np.log(y_hat)+(1-y)*np.log(1-y_hat))/len(y);
			except:
				return float("inf");
		def dJ(theta,X_b,y):   #梯度
			return ((X_b.T).dot(self._Sigmoid(X_b.dot(theta))-y))/len(X_b);
			
		def gradientdescent(initial_theta,X_b,y,eta_1 = 0.01,epsilon = 1e-8,n_iters_1 = 1e4): # 梯度下降法
			i_iters = 0;
			theta = initial_theta;
			while i_iters<n_iters_1:
				gradient = dJ(theta,X_b,y);
				last_theta = theta;
				theta = theta - eta_1*gradient;
				if abs(J(theta,X_b,y)-J(last_theta,X_b,y))<epsilon:
					break;
			return theta;
		X_b = np.hstack([np.ones([len(X_train),1]),X_train]);
		theta = np.zeros([X_b.shape[1],1]);
		y_train = y_train.reshape(len(y_train),1);
		self._theta = gradientdescent(theta,X_b,y_train,eta_1=eta,n_iters_1=n_iters);
		self._theta = self._theta.reshape(len(self._theta));
		self.intercept_ = self._theta[0];
		self.coef_ = self._theta[1:];
		return self;
		
	def _predict_proba(self,X_test):  #预测概率
		X_b = np.hstack([np.ones([len(X_test),1]),X_test]);
		y_predict = self._Sigmoid(X_b.dot(self._theta));
		return y_predict;
	
	def predict(self,X_test): #分类结果
		X_b = np.hstack([np.ones([len(X_test),1]),X_test]);
		proba = self._predict_proba(X_test);
		return np.array(proba>=0.5,dtype = 'int');

五、如何将逻辑回归推广到多分类问题中

逻辑回归在二分类问题中发挥中巨大的作用，如此好用的算法如果不能推广至多分类问题中，那真的是暴殄天物。好在，有很多种方法可以帮助那些二分类算法在多分类问题中继续发挥它们的余热。

5.1 一对一 ( one vs one. OvO）

对于N个类别的多分类问题，OvO的思想是每次从中选择两个类别$Y_i,Y_j$来训练一个分类器，共计需要训练$N(N-1)/2$个分类器；在测试阶段，将新样本提交给所有的分类器，于是我们可以得到$N(N-1)/2$个分类结果，最终结果通过投票选择。OvO的存储开销和测试时间开销通常会比较大，但是它的分类结果相对于OvR而言更加准确。

5.2 一对多 (one vs rest. OvR)

OvR的思想是每次将一个类别作为正类别，其余的类别都被归纳至负类别中，因此对于N个类别的多分类问题，我们只需要训练N个逻辑回归二分类器就可以了。在测试时，若仅有一个分类器预测为正类，则对应的类别标记作为最终的分类结果 ； 如果有多个分类器预测为正类，则需要考虑分类器的预测置信度，选择置信度大的类别标记作为分类结果。
由于OvR只需要训练N个分类器，因此存储开销和测试时间开销通常比OvO要小；但是在训练OvR分类器时每个分类器都会使用全部的训练样本，因此训练时间开销上OvR要大于OvO。