逻辑回归系数的置信区间是如何得出来的_步步为营KDD：逻辑回归

在介绍如何使用逻辑回归进行分类时，我们首先需要大概了解下什么是回归，什么是逻辑回归。

回归分析(Regression Analysis)

在统计学中，回归分析(regression analysis)是一个用于估算变量之间关系的统计学过程。回归分析关注的焦点是在一个因变量(dependent variable)和一个或多个自变量(independent variable)之间的关系。更明确的说法就是，回归分析帮助人们了解当自变量发生变化时，因变量是如何变化的。回归分析的目标是找到最能够代表所有观测资料的函数，这个叫函数做回归函数(regression function)。

回归模型(Regression Models)

回归模型涉及以下变量：

• β：未知参数，可能是一个标量，也可能是一个向量。
• X：自变量
• Y：因变量

Y是关于X与β的函数：

上面的近似式通常被格式化为E(Y|X)=f(X,β)。为了能够执行回归分析，函数f的形式必须被指定。有时候函数f的形式可从我们对Y和X的认知中获得，而不是依赖于数据。如果不存在这样的认知，那么我们应该选择一个灵活或者方便的函数形式作为函数f。

现在假设β是个长度为k的向量。为了执行回归分析，用户必须提供关于因变量Y的信息：

• 如果有N个形式为(Y,X)数据点被观测，且Nβ。
• 如果N=k个数据点被观测，那么函数f是线性的，Y=f(X,β)将能被准确求的，而不是近似值。
• 最常见的情况就是N>k个数据点被观测。在这种情况下，数据中有足够多的信息用于计算出一个唯一的β以最好地拟合数据，此时的回归模型可以被看作是超定系统(overdetermined system)。

在上述的最后一个场景中，回归分析可以提供以下用途：

• 找到能够最小化被度量值与因变量Y的预测值之间距离的β的解决方案。
• 在特定的统计假设情况下，回归分析使用信息差(the surplus of information)来提供关于β与因变量Y的预测值的统计信息。

回归的分类

回归分析按照涉及的自变量的多少，可分为一元回归分析和多元回归分析；按照自变量和因变量之间的关系类型，可分为线性回归分析和非线性回归分析。如果在回归分析中，只包括一个自变量和一个因变量，且二者的关系可用一条直线近似表示，这种回归分析称为一元线性回归分析。如果回归分析中包括两个或两个以上的自变量，且因变量和自变量之间是线性关系，则称为多元线性回归分析。

逻辑回归(Logistic Regression)

逻辑回归就是回归函数为逻辑函数的回归分析。

逻辑函数是一种常见的S形函数(Sigmoid Function)。它用如下公式表示：

e是自然对数基(natural logarithm base)，也称为欧拉数字(Euler's number)。标准的逻辑S形函数如下图所示：

如图所示，当x=0时，y=0.5；当x→+∞，y→1；当x→-∞，y→0。

如果横坐标足够大，逻辑函数很像一个阶跃函数(step function)。为了实现逻辑回归分类器，我们可以在每个特征上都乘以一个回归系数(回归系数也可以称为权重，它是一个向量)，然后把所有的结果值相加，将总和代入Sigmoid函数中。如果最终值大于0.5，则将该特征向量划分为1类；如果最终值小于0.5，则将该特征向量划分为0类。

为了区分权重向量与特征向量，将上面Sigmoid函数中的x换为z：

X表示特征向量，W表示权重向量。

新的Sigmoid函数如下所示：

那么如何才能获得最佳的回归系数，能够将测试数据准确地分类？

因为逻辑回归的代价函数是似然函数，需要最大化似然函数，所以需要使用梯度上升算法。

梯度上升法(Gradient Ascent)

梯度上升法是一种最优化方法。梯度上升法的基本思想是：要找到某函数的最大值，最好的方法是沿着该函数的梯度方向探寻。

如果梯度记为∇(nabla)，则函数f(x,y)的梯度由下式表示：

上面的公式表示，梯度要沿着x的方向移动∂f(x,y)/∂x的距离，沿着y的方向移动∂f(x,y)/∂y的距离。其中，函数f(x,y)必须要在待计算的点上有定义并且可微。

梯度算子总是指向函数值增长最快的方向。每次移动的量称为步长(step size)，记做α。用向量来表示的话，梯度算法的迭代公式如下：

该公式将一直被执行，直至达到某个停止条件为止，比如迭代次数达到某个指定值或算法达到某个可以允许的误差范围。

#!/usr/bin/env python# _*_ coding: utf-8 _*_from numpy import *'''加载数据文件中的数据格式如下所示：X           Y           ClassLabel-0.017612   14.053064   0-1.395634   4.662541    1-0.752157   6.538620    0-1.322371   7.152853    00.423363    11.054677   0......'''def loadDataSet():    dataMat = []; labelMat = []    fr = open('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'''S形函数'''def sigmoid(inX):    return 1.0 / (1 + exp(-inX))'''采用梯度上升算法获取权重(回归系数)'''def gradAscent(dataMatIn, classLabels):    dataMatrix = mat(dataMatIn)    '''矩阵转置，矩阵形状为dataMatrix.rowsCount * 1。'''    labelMat = mat(classLabels).transpose()     m, n = shape(dataMatrix)    '''步长'''    alpha = 0.001     '''迭代次数'''    maxCycles = 500     '''初始权重，矩阵形状为dataMatrix.columnsCount * 1，元素为1。此例中，n=3。'''    weights = ones((n,1))     for k in range(maxCycles):        '''        矩阵m*n乘以矩阵n*p，得到矩阵m*p。所以矩阵h的形状为dataMatrix.rowsCount * 1         dataMatrix * weights得到的是每行的z值。        将z值代入sigmoid函数得到每行的类别值。        '''        h = sigmoid(dataMatrix * weights)        '''        得到每行类别的误差值        '''        error = (labelMat - h)        '''        权重 = 权重 + 步长 * 梯度        '''        weights = weights + alpha * dataMatrix.transpose() * error    return weights

随机梯度上升算法(Stochastic Gradient Ascent)

梯度上升算法在每次更新回归系数时都需要遍历整个数据集，该方法在处理大数据量时将变得不可用。一个改进的方法是一次仅用一个样本点来更新回归系数，该方法称为随机梯度上升算法。由于可以在新样本到来时对分类器进行增量式更新，因而该算法是一个在线学习算法。与在线学习算法相对应的，像之前的梯度学习算法这样一次处理所有数据被称作“批处理”。

'''采用随机梯度上升算法获取权重'''def stocGradAscent0(dataMatrix, classLabels):    m,n = shape(dataMatrix)    alpha = 0.01    weights = ones(n)    for i in range(m):        h = sigmoid(sum(dataMatrix[i] * weights))        error = classLabels[i] - h        weights = weights + alpha * error * dataMatrix[i]    return weights

梯度上升算法和随机梯度上升算法的区别有二：一是前者的变量h和error都是向量，后者都是数值；二是前者拥有矩阵的转换过程。

随机梯度上升算法的拟合效果不如梯度上升算法好，但它的迭代次数少，速度快。

由于随机梯度上升算法会出现更多的错误分类的数据点，因此在每次迭代时会引发系数的剧烈改变。我们期望算法能够避免来回波动，从而收敛到某个值，另外，收敛速度也需要加快。

改进的随机梯度上升算法如下：

'''改进后的随机梯度上升算法'''def stocGradAscent1(dataMatrix, classLabels, numIter=150):    m,n = shape(dataMatrix)    weights = ones(n)    for j in range(numIter):        dataIndex = range(m)        for i in range(m):            '''            每次迭代都会调整步长，但不会减少到0.            '''            alpha = 4/(1.0+j+i) + 0.01            '''            随机获取样本来更新回归系数，以减少周期性的波动            '''            randIndex = int(random.uniform(0, len(dataIndex)))            h = sigmoid(sum(dataMatrix[randIndex] * weights))            error = classLabels[randIndex] - h            weights = weights + alpha * error * dataMatrix[randIndex]            del(dataIndex[randIndex])    return weights