《机器学习实战》——logistic回归

说明：对书中代码错误部分做了修正，可运行于python3.4

基本原理：现在有一些数据点，用一条直线对这些数据进行拟合，将它们分为两类。这条直线叫做最佳拟合直线，这个拟合过程叫做回归。logistic回归的思想是，利用一个阶跃函数（在某一点突然由0变1），实现分类器。Sigmoid函数近似于阶跃函数：

现在将每个特征乘以一个回归系数，再全部相加，总和带入函数作为输入自变量z，进而得到一个0-1之间输出，四舍五入之后划分为0和1两类。这就是分类器的思想。它的代价函数为：

关键词：数据点，直线，最佳拟合，Sigmoid阶跃函数，回归系数

算法实施：

1.数据。假设文本有m行3列，第一列是第零维系数w=1.0，第2、3列是特征x1,x2——loadDataSet()

2. sigmoid函数。套用公式——sigmoid(inX)

3. 梯度算法。梯度上升算法，求出最佳的w参数(矩阵)，每次计算梯度需要所有数据点信息.——gradAscent(dataMatIn, classLabels)

随机梯度上升算法，改进之后，每次计算梯度只需要带入一个数据点信息——stocGradAscent0(dataMatrix, classLabels)

改进后的随机梯度上升算法——stocGradAscent1(dataMatrix, classLabels,numIter = 150)

4. 计算单词出现的条件概率。首先计算侮辱性句子在所有句子中的比例（即为侮辱性文档的概率）。类别1为侮辱性，将所有次类别的矩阵行向量加和，并统计所有侮辱性类别中出现的单词数；同理，如果是类别0也将行向量加和，统计非侮辱性行中出现的所有单词数；用向量中的各个值除以该类别的总单词数，就是各个词属于侮辱性和非侮辱文档的概率（即各个类别出现该单词的条件概率）；其中侮辱性单词向量中概率最大的词即为最能表征该类别的单词。

测试：为了避免各个单词条件概率联乘中出现的零概率，所有词出现次数（分子）初始化为1，该类别总词数（分母）出现次数初始化为2.0——trainNB0(trainMatrix, trainCategory)

5. 画出数据散点图和拟合直线。——plotBestFit(weights)

算法代码：

# coding utf-8
from numpy import *
def loadDataSet():
	#####数据导入#######
    dataMat = []						# 数据矩阵(第一列是w=1.0，第2、3列是特征)
    labelMat = []  						# 标签矩阵（标识每行数据的类别）
    # fd = open('testSet.txt')			# 格式：2.2 3.0 1
    myArr = [[-3.5, -3, 0],[-2.3, 0, 0],
		[-1.0, -0.1, 0],[-1.3, -1.0, 0],
		[-2.5, 5, 0],[-3.5, 7, 1],
		[-1.5, 16, 1],[1, 10, 1],
		[1, 5, 1],[1, 3, 0]]
    # for item in fd.readlines(): 
    for itemArr in myArr:  
        # itemArr = item.strip().split()  
        dataMat.append([1.0, float(itemArr[0]), float(itemArr[1])])  
        labelMat.append(int(itemArr[2])) 
    return dataMat, labelMat  

def sigmoid(inX):  
	#####计算sigmoid函数（即logistic函数）#####  
    return 1.0/(1+exp(-inX))

def gradAscent(dataMatIn, classLabels):  
	#####梯度上升算法，求出最佳的w参数矩阵#######
    dataMatrix = mat(dataMatIn)          # dataMatIn格式：(1.0，第一特征值，第二特征值)
    labelMat = mat(classLabels).transpose()  # 标签向量转置为列矩阵
    m,n = shape(dataMatrix)  				 # n*3
    alpha = 0.001  							 # 梯度表示移动方向，而alpha表示移动量的大小（步长）
    maxCycles = 500  						 # 最多移动的步数（步数越多越精确）
    weights = ones((n,1))  					 
    for k in range(maxCycles):  
        h = sigmoid(dataMatrix*weights)  
        error = (labelMat - h)  
        weights = weights + alpha * dataMatrix.transpose() * error  
    return weights

def stocGradAscent0(dataMatrix, classLabels):
    #####随机梯度上升算法，一次仅用一个样本点来更新回归系数#######
    dataMatrix = mat(dataMatrix)          # dataMatIn格式：(1.0，第一特征值，第二特征值)
    labelMat = mat(classLabels).transpose()  # 标签向量转置为列矩阵
    m,n = shape(dataMatrix)
    alpha = 0.001
    weights = ones((n,1))
    for i in range(m):
        h = sigmoid(sum(dataMatrix[i]*weights))
        error = classLabels[i] - h
        weights = weights + alpha*error*dataMatrix[i].transpose()
    return weights

def stocGradAscent1(dataMatrix, classLabels,numIter = 150):
    #####随机梯度上升算法，一次仅用一个样本点来更新回归系数#######
    dataMatrix = mat(dataMatrix)          # dataMatIn格式：(1.0，第一特征值，第二特征值)
    labelMat = mat(classLabels).transpose()  # 标签向量转置为列矩阵
    m,n = shape(dataMatrix)
    alpha = 0.001
    weights = ones((n,1))
    for j in range(numIter):
        dataIndex = range(m)
        for i in range(m):
            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].transpose()
    return weights

def plotBestFit(weights):
	#####画出最佳拟合直线#######
    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 = array((-weights[0]-weights[1]*x)/weights[2])[0]
    ax.plot(x,y)
    plt.xlabel('X1');plt.ylabel('X2')
    plt.show()

dataMat, labelMat = loadDataSet()   
# weights = gradAscent(dataMat, labelMat)          # 梯度上升算法。此时weights类型为矩阵
# weights = stocGradAscent0(dataMat, labelMat)     # 随机梯度上升算法。此时weights类型为矩阵
weights = stocGradAscent1(dataMat, labelMat)       # 改进的随机梯度上升算法。此时weights类型为矩阵
plotBestFit(weights) 

结果：

1.梯度上升算法

2.随机梯度上升算法

3.改进的随机梯度上升算法