机器学习（六）— logistic回归

　　最近一直在看机器学习相关的算法，今天学习logistic回归，在对算法进行了简单分析编程实现之后，通过实例进行验证。

一 logistic概述

　　个人理解的回归就是发现变量之间的关系，也就是求回归系数，经常用回归来预测目标值。回归和分类同属于监督学习，所不同的是回归的目标变量必须是连续数值型。

　　今天要学习的logistic回归的主要思想是根据现有的数据对分类边界线建立回归公式，以此进行分类。主要在流行病学中应用较多，比较常用的情形是探索某疾病的危险因素，根据危险因素预测某疾病发生的概率等等。logistic回归的因变量可以是二分类的，也可以是多分类的，但是二分类的更为常用，也更加容易解释，所以实际中最为常用的就是二分类的logistic回归。

　　今天我们就二分类进行分析，我们在回归分析中需要一个函数可以接受所有的输入然后预测出类别，假定用0和1分别表示两个类别，logistic函数曲线很像S型，故此我们可以联系sigmoid函数：σ = 1/(1/(1+e^-z))。为了实现logistic回归分类器，我们可以在每个特征上乘以一个回归系数，将所有的乘积相加，将和值代入sigmoid函数中，得到一个范围为0-1之间的数，如果该数值大于0.5则被归入1类，否则被归为0类。

　　基于之前的分析，需要找到回归系数，首先我们可以将sigmoid函数的输入形式记为：z = w₀x₀ + w₁x₁ +...+w_nx_n,其中x为输入数据，相应的w就是我们要求的系数，为了求得最佳系数，结合最优化理论，我们可以选取梯度上升法优化算法。梯度上升法的基本思想是:要找到函数的最大值，最好的方法是沿着该函数的梯度方向寻找。要想更进一步的了解这个方法，建议去看Andrew Ng的机器学习课程，记得在第二节主要讲述的就是梯度下降法，与梯度上升所不同的是它求得的是函数的最小值，不过思想是一致的。

二 python实现

　　基于之前的分析，在本节我们对logistic回归一步一步采用python编程实现，今天我用的是2.7版本的，代码如下：

#coding:utf-8
from numpy import *
import math
import matplotlib.pyplot as plt

#导入数据
def loadDataSet():
  dataMat = []
  labelMat = []
  fr = open('testSet.txt')
  for line in fr.readlines():
    lineArr = line.strip().split()#将文本中的每行中的字符一个个分开，变成list
    dataMat.append([1.0,float(lineArr[0]),float(lineArr[1])])
    labelMat.append(int(lineArr[2]))
  return dataMat,labelMat

#定义sigmoid函数
def sigmoid(inX):
    return 1.0/(1+exp(-inX))

#梯度上升方法求出回归系数
def gradAscent(data,label):
  dataMat = mat(data)
  labelMat = mat(label).transpose()
  m,n = shape(dataMat)
  alpha = 0.001
  maxCycles = 500
  weights = ones((n,1))
  for item in range(maxCycles):
    h = sigmoid(dataMat * weights)
    error = (labelMat - h)#注意labelMat中的元素的数据类型应为int
    weights = weights + alpha * dataMat.transpose() * error
  return weights

'''
#测试
data,label = loadDataSet()
print gradAscent(data,label)
'''

##求出回归系数之后，就确定了不同数据类别之间的分隔线，为了便于理解，可以画出那条线
def plotBestFit(weights):
  dataMat,labelMat = loadDataSet()
  dataArr = array(dataMat)
  n = shape(dataArr)[0]
  xcode1 = []
  ycode1 = []
  xcode2 = []
  ycode2 = []
  for i in range(n):
    if int(labelMat[i]) == 1:
      xcode1.append(dataArr[i,1])
      ycode1.append(dataArr[i,2])
    else:
      xcode2.append(dataArr[i,1])
      ycode2.append(dataArr[i,2])
  fig = plt.figure()
  ax = fig.add_subplot(111)
  ax.scatter(xcode1,ycode1,s = 30,c = 'red',marker = 's')
  ax.scatter(xcode2,ycode2,s = 30,c = 'green')
  x = arange(-3.0,3.0,0.1)
  y = (-weights[0] - weights[1] * x) / weights[2]
  ax.plot(x,y)
  plt.xlabel('x1')
  plt.ylabel('y1')
  plt.show()
  
'''
#测试
data,label = loadDataSet()
weights = gradAscent(data,label)
plotBestFit(weights.getA())
'''
##改进的梯度上升法
def stocGradAscent1(dataMatrix, classLabels, numIter=150):
    m,n = shape(dataMatrix)
    weights = ones(n)   #initialize to all ones
    for j in range(numIter):
        dataIndex = range(m)
        for i in range(m):
            alpha = 4/(1.0+j+i)+0.0001    
            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
''' 
#测试
data,label = loadDataSet()
weights = stocGradAscent1(array(data),label)
plotBestFit(weights)
'''

三实例分析

基于之前的分析，本节采用Logistic回归来预测患有疝病的马的存活问题，代码如下：

def classifyVector(inX,weights):
  prob = sigmoid(sum(inX * weights))
  if prob > 0.5:
    return 1.0
  else:
    return 0.0
def colicTest():
    frTrain = open('horseColicTraining.txt'); frTest = open('horseColicTest.txt')
    trainingSet = []; trainingLabels = []
    for line in frTrain.readlines():
        currLine = line.strip().split('	')
        lineArr =[]
        for i in range(21):
            lineArr.append(float(currLine[i]))
        trainingSet.append(lineArr)
        trainingLabels.append(float(currLine[21]))
    trainWeights = stocGradAscent1(array(trainingSet), trainingLabels, 1000)
    errorCount = 0; numTestVec = 0.0
    for line in frTest.readlines():
        numTestVec += 1.0
        currLine = line.strip().split('	')
        lineArr =[]
        for i in range(21):
            lineArr.append(float(currLine[i]))
        if int(classifyVector(array(lineArr), trainWeights))!= int(currLine[21]):
            errorCount += 1
    errorRate = (float(errorCount)/numTestVec)
    print "错误率是：",errorRate
    return errorRate
def multiTest():
    numTests = 10; errorSum=0.0
    for k in range(numTests):
        errorSum += colicTest()
    print "平均错误率是",(numTests, errorSum/float(numTests))
        
multiTest()

最后可以看出错误率在35%左右，通过调节步长还是可以进一步减小错误率。

　　Logistic回归的目的是寻找到一个非线性sigmoid函数的最佳拟合参数，可以采用梯度上升法优化，而在这个过程中，为了减少时间复杂度，又可以使用随机梯度上升法来简化梯度上升法。