第二十二章 logistic regression 算法（上）

逻辑回归（LogisticRegression）

Logistic regression（逻辑回归）是当前业界比较常用的机器学习方法，用于估计某种事物的可能性。类似某用户购买某商品的可能性，某病人患有某种疾病的可能性啊等等。这个世界是随机的（当然了，人为的确定性系统除外，但也有可能有噪声或产生错误的结果，只是这个错误发生的可能性太小了，小到千万年不遇，小到忽略不计而已），所以万物的发生都可以用可能性或者几率（Odds）来表达。“几率”指的是某事物发生的可能性与不发生的可能性的比值。

Logistic regression可以用来回归，也可以用来分类，主要是二分类。它给我们提供的就是你的这个样本属于正类的可能性是多少。

假设我们的样本是{x, y}，y是0或者1，表示正类或者负类，x是我们的m维的样本特征向量。那么这个样本x属于正类，也就是y=1的“概率”可以通过下面的逻辑函数来表示：

这里θ是模型参数，也就是回归系数，σ是sigmoid函数。实际上这个函数是由下面的对数几率（也就是x属于正类的可能性和负类的可能性的比值的对数）变换得到的：

换句话说，y也就是我们关系的变量，我们将这些对应x1,x2,…, xm的权值叫做回归系数，表达为θ1, θ2,…, θm。他们的加权和就是你的总得分了。

所以说上面的logistic回归就是一个线性分类模型，它与线性回归的不同点在于：为了将线性回归输出的很大范围的数，例如从负无穷到正无穷，压缩到0和1之间，这样的输出值表达为“可能性”才能说服广大民众。当然了，把大值压缩到这个范围还有个很好的好处，就是可以消除特别冒尖的变量的影响。而实现这个伟大的功能其实就只需要平凡一举，也就是在输出加一个logistic函数。另外，对于二分类来说，可以简单的认为：如果样本x属于正类的概率大于0.5，那么就判定它是正类，否则就是负类。

所以说，LogisticRegression 就是一个被logistic方程归一化后的线性回归，仅此而已。

归入到正统的机器学习框架下，模型选好了，只是模型的参数θ还是未知的，我们需要用我们收集到的数据来训练求解得到它。那我们下一步要做的事情就是建立代价函数了。

LogisticRegression最基本的学习算法是最大似然。

假设我们有n个独立的训练样本{(x1, y1) ,(x2, y2),…, (xn, yn)}，y={0, 1}。那每一个观察到的样本(xi, yi)出现的概率是：

上面为什么是这样呢？当y=1的时候，后面那一项是不是没有了，那就只剩下x属于1类的概率，当y=0的时候，第一项是不是没有了，那就只剩下后面那个x属于0的概率（1减去x属于1的概率）。所以不管y是0还是1，上面得到的数，都是(x,y)出现的概率。那我们的整个样本集，也就是n个独立的样本出现的似然函数为（因为每个样本都是独立的，所以n个样本出现的概率就是他们各自出现的概率相乘）：

那最大似然法就是求模型中使得似然函数最大的系数取值θ*。这个最大似然就是我们的代价函数（cost function）了。

OK，那代价函数有了，我们下一步要做的就是优化求解了。我们先尝试对上面的代价函数求导，看导数为0的时候可不可以解出来，也就是有没有解析解，有这个解的时候，就皆大欢喜了，一步到位。如果没有就需要通过迭代了，耗时耗力。

我们先变换下L(θ)：取自然对数，然后化简（不要看到一堆公式就害怕哦，很简单的哦，只需要耐心一点点，自己动手推推就知道了。注：有xi的时候，表示它是第i个样本，下面没有做区分了，相信你的眼睛是雪亮的），得到：

这时候，用L(θ)对θ求导，得到：

然后我们令该导数为0，你会很失望的发现，它无法解析求解。不信你就去尝试一下。所以没办法了，只能借助高大上的迭代来搞定了。这里选用了经典的梯度下降算法。

优化求解

梯度下降(gradient descent)

Gradient descent 又叫 steepest descent，是利用一阶的梯度信息找到函数局部最优解的一种方法，也是机器学习里面最简单最常用的一种优化方法。它的思想很简单，和我开篇说的那样，要找最小值，我只需要每一步都往下走（也就是每一步都可以让代价函数小一点），然后不断的走，那肯定能走到最小值的地方，例如下图所示：

但，我同时也需要更快的到达最小值啊，怎么办呢？我们需要每一步都找下坡最快的地方，也就是每一步我走某个方向，都比走其他方法，要离最小值更近。而这个下坡最快的方向，就是梯度的负方向了。

对logistic Regression来说，梯度下降算法新鲜出炉，如下：

其中，参数α叫学习率，就是每一步走多远，这个参数蛮关键的，该参数前的正负号分别表示梯度上升或梯度下降。

梯度上升、随机梯度上升算法的python2实现：

from numpy import *

def loadDataSet():#打开文本文件并逐行读取

dataMat = []; labelMat = []

fr = open('C:/Users/HZF/Desktop/machinelearninginaction/Ch05/testSet.txt')

for line in fr.readlines():

lineArr = line.strip().split()

dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])

labelMat.append(int(lineArr[2]))

return dataMat,labelMat

def sigmoid(inX):#sigmoid函数

return 1.0/(1+exp(-inX))

def gradAscent(dataMatIn, classLabels):#梯度上升函数，dataMatIn 为特征列

dataMatrix = mat(dataMatIn)       #convert to NumPy matrix，X0假定为1.0，x1,x2为testSet.txt的前两列数据

labelMat = mat(classLabels).transpose() #convert to NumPy matrix，将类标签的矩阵行向量转置成列向量

m,n = shape(dataMatrix)#计算矩阵大小

alpha = 0.001

maxCycles = 500#maxCycles是迭代次数

weights = ones((n,1))

for k in range(maxCycles):       #heavy on matrix operations

h = sigmoid(dataMatrix*weights)   #matrix mult，计算了m*n次，因为计算量大，所以若应用在真实数据时，需要对该方法改进

error = (labelMat - h)       #vector subtraction

weights = weights + alpha * dataMatrix.transpose()* error #matrix mult，梯度上升算法推导

return weights

def plotBestFit(weights):#画出数据集和logistic回归最佳拟合直线的函数,此处的weights=weights(调用函数梯度上升里的weights).getA(),

import matplotlib.pyplot as plt

dataArr = array(dataMat)

n = shape(dataArr)[0]

xcord1 = []; ycord1 = []

xcord2 = []; ycord2 = []

for i in range(n):

if int(labelMat[i])== 1:

xcord1.append(dataArr[i,1]); ycord1.append(dataArr[i,2])

else:

xcord2.append(dataArr[i,1]); ycord2.append(dataArr[i,2])

fig = plt.figure()

ax = fig.add_subplot(111)

ax.scatter(xcord1, ycord1, s=30, c='red', marker='s')

ax.scatter(xcord2, ycord2, s=30, c='green')

x = arange(-3.0, 3.0, 0.1)

y = (-weights[0]-weights[1]*x)/weights[2]#最佳拟合直线，y=0,因为0是两个分类0/1的分界处

ax.plot(x, y)

plt.xlabel('X1'); plt.ylabel('X2');

plt.show()#画出x1与x2特征数据散点图的分界线，x0默认是1.0

#梯度上升算法的改进算法——随机梯度上升算法（一次仅用一个样本点更新回归系数，减少计算复杂度）

def stocGradAscent0(dataMatrix, classLabels):#随机梯度上升算法，该算法是一个在线学习算法。

m,n = shape(dataMatrix)

alpha = 0.01

weights = ones(n)  #initialize to all ones

maxCycles = 200

for k in range(maxCycles):#此行为自己实验添加代码，经实验，随机梯度上升算法迭代了200次后效果才可与梯度上升算法效果相近且都比较好

for i in range(m):

h = sigmoid(sum(dataMatrix[i]*weights))#h为具体值

error = classLabels[i] - h#error为具体值

weights = weights + alpha * error * dataMatrix[i]

return weights

def classifyVector(inX, weights):#分类函数

prob = sigmoid(sum(inX*weights))#weights取自前面几种最优化的weights,InX为测试集上的特征向量

if prob > 0.5: return 1.0

else: return 0.0

if __name__=="__main__":

#dataMat,labelMat=loadDataSet()

#weights=gradAscent(dataMat, labelMat)

#print weights

#wei=weights.getA()

dataMat,labelMat=loadDataSet()

dataArr = array(dataMat)

#weights=gradAscent(dataMat, labelMat)#迭代500次

#plotBestFit(weights.getA())

weights=stocGradAscent0(array(dataMat), labelMat)#迭代200次，读者可以不用迭代看看效果

plotBestFit(weights)

classifyVector(inX, weights)

上篇就写逻辑回归算法的基本理论、梯度上升，随机梯度上升的代码实现这么多！

参考文献：

1、机器学习-逻辑回归理论简介 - xietingcandice的专栏 - 博客频道 - CSDN.NET

2、《机器学习》（书）