Logistic 回归算法及Python实现

每个回归系数初始化为1
重复R次:
    计算整个数据集的梯度
    使用alpha下的gradient更新回归系数的向量
返回回归系数

#! /usr/bin/env python
# -*- coding: utf-8 -*-
"""
实现logistic回归分类算法， 数据集为: dataset.csv
"""

import numpy as np
import matplotlib.pyplot as plt

def loadDataSet():
   """
   加载数据集
   return: 数据列表， 标签列表
   """
   dataMat = []
   labelMat = []
   #打开数据集
   fr = open('dataset.csv')
   #遍历每一行
   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函数
    @: param intX: 矩阵计算的结果(100x1)
    @: return: 计算结果
    """
    return 1.0 / (1 + np.exp(-inX))

def gradAscent(dataMat, labelMat):
    """
    梯度上升函数
    @: param dataMat: 数据集
    @: param labelMat: 标签集
    @： return: 权重参数矩阵(最佳回归系数)
    """
    # 将数据转为numpy的数组
    dataMatrix = np.mat(dataMat)
    labelMat = np.mat(labelMat).transpose()
    # 获取矩阵的行列数
    m, n = np.shape(dataMatrix)
    # 初始化参数
    alpha = 0.001
    # 初始化迭代次数
    maxCyc = 500
    # 初始化矩阵的权重参数矩阵， 均为1
    weights = np.ones((n, 1))
    # 开始迭代计算
    for k in range(maxCyc):
        h = sigmoid(dataMatrix * weights)
        # 计算误差
        error = labelMat-h
        # 更新迭代参数
        weights = weights + alpha * dataMatrix.transpose() * error
    return weights

def plotBestFit(weights):
    #导入数据
    dataMat, labelMat = loadDataSet()
    #创建数组
    dataArr = np.array(dataMat)
    #获取数组行数
    n = np.shape(dataArr)[0]
    #初始化坐标
    xcord1 = []; ycord1 = []
    xcord2 = []; ycord2 = []
    #遍历每一行数据
    for i in range(n):
        #如果对应的类别标签对应数值1，就添加到xcord1，ycord1中
        if int(labelMat[i]) == 1:
            xcord1.append(dataArr[i,1]); ycord1.append(dataArr[i,2])
        #如果对应的类别标签对应数值0，就添加到xcord2，ycord2中
        else:
            xcord2.append(dataArr[i,1]); ycord2.append(dataArr[i,2])
    #创建空图
    fig = plt.figure()
    #添加subplot，三种数据都画在一张图上
    ax = fig.add_subplot(111)
    #1类用红色标识，marker='s'形状为正方形
    ax.scatter(xcord1, ycord1, s=30, c='red', marker='s')
    #0类用绿色标识，弄认marker='o'为圆形
    ax.scatter(xcord2, ycord2, s=30, c='green')
    #设置x取值，arange支持浮点型
    x = np.arange(-3.0, 3.0, 0.1)
    #配计算y的值
    y = (-weights[0]-weights[1]*x)/weights[2]
    #画拟合直线
    ax.plot(x, y)
    #贴坐标表头
    plt.xlabel('X1'); plt.ylabel('X2')
    #显示结果
    plt.show()

if __name__ == '__main__':
    dataArr, labelMat = loadDataSet()
    weights = gradAscent(dataArr, labelMat)
    print (weights)
    plotBestFit(weights.getA())

所有回归系数初始化为1
对数据集中每个样本
    计算该样本的梯度
    使用alpha * gradient 更新回归系数值
返回回归系数值

def randomGradAscent(dataMat, labelMat):
    """
    随机梯度上升函数
    @: param dataMat: 数据集
    @: param labelMat: 标签集
    @： return: 权重参数矩阵(最佳回归系数)
    """
    dataMatrix = np.array(dataMat)
    m, n = np.shape(dataMatrix)
    # 设置步长
    alpha = 0.01
    #初始化参数
    weights = np.ones(n)
    for i in range(m):
        h = sigmoid(sum(dataMatrix[i]*weights))
        # 计算误差
        error = labelMat[i]-h
        # 更新权重矩阵
        weights = weights + alpha * error * dataMatrix[i]
    return weights

def randomGradAscent2(dataMat, labelMat):
    """
    改进的随机梯度上升算法
    """
    dataMatrix = np.array(dataMat)
    m, n = np.shape(dataMatrix)
    # 初始化参数
    weights = np.ones(n)
    # 迭代次数
    numIter = 500
    for i in range(numIter):
        # 初始化index列表，这里要注意将range输出转换成list
        dataIndex = list(range(m))
        # 遍历每一行数据，这里要注意将range输出转换成list
        for j in list(range(m)):
            # 更新alpha值，缓解数据高频波动
            alpha = 4/(1.0+i+j)+0.0001
            # 随机生成序列号，从而减少随机性的波动
            randIndex = int(np.random.uniform(0, len(dataIndex)))
            # 序列号对应的元素与权重矩阵相乘，求和后再求sigmoid
            h = sigmoid(sum(dataMatrix[randIndex]*weights))
            # 求误差，和之前一样的操作
            error = labelMat[randIndex] - h
            # 更新权重矩阵
            weights = weights + alpha * error * dataMatrix[randIndex]
            # 删除这次计算的数据
            del(dataIndex[randIndex])
    return weights

#! /usr/bin/env python
# -*- coding: utf-8 -*-

import numpy as np
import logistic

def classifyVector(inX, weights):
    """
    分类
    """
    prob = logistic.sigmoid(sum(inX * weights))
    if prob > 0.5:
        return 1.0
    else:
        return 0.0

def colicTest():
    """
    训练和测试模型
    """
    frTrain = open('horse_colic_train.txt')
    frTest = open('horse_colic_test.txt')
    trainingSet = []
    trainingLabels = []
    # --------------训练------------------
    for line in frTrain.readlines():
        currLine = line.strip().split('t')
        lineArr =[]
        for i in range(21):
            lineArr.append(float(currLine[i]))
        trainingSet.append(lineArr)
        trainingLabels.append(float(currLine[21]))
    trainWeights = logistic.randomGradAscent2(np.array(trainingSet), trainingLabels)

    # --------------测试------------------
    numTestVec = 0.0
    errorCount = 0
    for line in frTest.readlines():
        numTestVec += 1.0
        currLine = line.strip().split('t')
        lineArr =[]
        for i in range(21):
            lineArr.append(float(currLine[i]))
        if int(classifyVector(np.array(lineArr), trainWeights))!= int(currLine[21]):
            errorCount += 1
    errorRate = (float(errorCount)/numTestVec)
    print ("测试的错误率为: %f" % errorRate)
    return errorRate

def multiTest():
    """
    多次测试
    """
    numTests = 10 
    errorSum = 0.0
    for k in range(numTests):
        errorSum += colicTest()
    print ("第 %d 次迭代后，错误率为: %f" % (numTests, errorSum / float(numTests)))

if __name__ == '__main__':
    multiTest()

测试的错误率为: 0.417910
测试的错误率为: 0.328358
测试的错误率为: 0.417910
测试的错误率为: 0.268657
测试的错误率为: 0.313433
测试的错误率为: 0.388060
测试的错误率为: 0.417910
测试的错误率为: 0.358209
测试的错误率为: 0.343284
测试的错误率为: 0.298507
第 10 次迭代后，错误率为: 0.355224