机器学习实战之Logistic回归与判断乳腺癌实例

一、回归的概念与Logistic回归

1.1 回归

假设有一些数据，能在坐标轴上展示，我们利用一条直线对这些点进行拟合(该线称为最佳拟合直线)，这个拟合过程就称作为回归。

它的方程可以表示为：f ( x )= w1x1+w2x2+…+wdxd + b
而我们要做的就是求出wi与b，使其对应的f(x)的值能尽量接近y

1.2 Logistic回归

Logistic回归一种二分类算法，它利用的是Sigmoid函数阈值在[0,1]这个特性。现在首先我们知道以下理论：

sigmoid函数如下：

而我们的函数是这样：

于是，就可以用sigmo函数

最后的问题，就是求解w矩阵（即合适的向量参数）。

1.3 极大似然法求参数

关于极大似然法，就是根据已知的数据，反过来求参数。这里不过多叙述。这里粘上比较好的介绍极大似然法的文章https://www.zhihu.com/question/24124998

设：

就是在已知x和参数的情况下，求样本是正例还是负例的概率。理想情况是每个点的p都是1。但显然是很难实现的。一般情况下，正样本概率>0.5即视为是正例。
将两个公式合并在一起

转换一下，就是这样：

cost是代价函数，y=0，就是求负样本的概率；y=1，求正样本概率。
取对数得：

cost越大，就说明概率最高，因此要求cost的最大值。
设样本互相独立，整个样本集生成的概率就是样本生成概率的积。取对数得到似然函数：

L对θ求偏导，得到极值。
于是要用到梯度上升方法（求最大值）。而如果是求最小值，就是用梯度下降方法。

1.4 梯度上升法

α为步长，控制更新的速度，即每次上升多少。
求解L对θ的偏导过程如下：

得到权重更新公式：

二、代码实现

2.1 梯度上升法

上面已经介绍了梯度上升法的理论，现在用代码将它实现

def gradAscent(dataMatIn, classLabels):
    dataMatrix = np.mat(dataMatIn)  # 转换成numpy的mat
    labelMat = np.mat(classLabels).transpose()  # 转换成numpy的mat,并进行转置
    m, n = np.shape(dataMatrix)  # 返回dataMatrix的大小。m为行数,n为列数。
    # print(m,n)
    alpha = 0.01  # 移动步长,也就是学习速率,控制更新的幅度。
    maxCycles = 500  # 最大迭代次数
    weights = np.ones((n, 1))
    # print(weights)
    # print(dataMatrix * weights)
    # print(dataMatrix.transpose())

    for k in range(maxCycles):
        h = sigmoid(dataMatrix * weights)  # 梯度上升矢量化公式
        # print(dataMatrix * weights)
        # print(h)
        error = labelMat - h
        # print(error)
        weights = weights + alpha * dataMatrix.transpose() * error
    
    return weights.getA()  # 将矩阵转换为数组，并返回

将更新幅度alpha设为0.01。重点在这里：

weights = weights + alpha * dataMatrix.transpose() * error

根据权重和训练样本计算估计值h

error就是原来数据与用sigmoid算出的值的差

最后就得到参数

这样会有一个缺点：每进行一次迭代更新就会计算所有样本，因此得到的模型正确率比较高，但同时计算复杂度高，算法耗时。

2.2 随机梯度上升算法

根据样本数量进行迭代，每计算一个样本就进行一次更新。每次随即取一个样本，用完将它删除，缩短了时间。

def stocGradAscent1(dataMatrix, classLabels, numIter=150):
    m, n = np.shape(dataMatrix)  # 返回dataMatrix的大小。m为行数,n为列数。
    weights = np.ones(n)  # 参数初始化										#存储每次更新的回归系数
    for j in range(numIter):
        dataIndex = list(range(m))
        for i in range(m):
            alpha = 4 / (1.0 + j + i) + 0.01  # 降低alpha的大小，每次减小1/(j+i)。
            randIndex = int(random.uniform(0, len(dataIndex)))  # 随机选取样本
            h = sigmoid(sum(dataMatrix[randIndex] * weights))  # 选择随机选取的一个样本，计算h
            error = classLabels[randIndex] - h  # 计算误差
            weights = weights + alpha * error * dataMatrix[randIndex]  # 更新回归系数
            del (dataIndex[randIndex])  # 删除已经使用的样本
    # print(weights)
    # print("sto ", type(weights))
    return weights  # 返回

三、实例测试

3.1 数据准备

数据集来自sklearn下的breast_cancer，即乳腺癌数据集
一共有569条数据，结果有两个：0是恶性，1是良性
共有30个特征。

再选取前400个作为训练，剩下的作为测试

sigmoid函数

def sigmoid(inX):
    return 1.0 / (1 + np.exp(-inX))

分类

def classifyVector(inX, weights):
    # print(type(inX), type(weights))
    prob = sigmoid(sum(inX * weights))
    # print(type(prob))
    if prob > 0.5:
        return 1.0
    else:
        return 0.0

3.2 计算错误率

使用两种方法计算出w系数再利用测试集算出错误率

    trainWeights1 = gradAscent(trainSets, trainLabels)  # 梯度上升法
    # 用梯度上升法算出的w矩阵如下：
    # [[9.12640384e+02]
    #  [1.34656613e+03]
    #  [5.43171929e+03]
    #  [3.18052912e+03]
    #  [1.01783385e+01]
    #  [3.91572955e-01]
    #  [-9.25792787e+00]
    #  [-3.95325370e+00]
    #  [1.86576269e+01]
    #  [8.00688379e+00]
    #  [3.37150224e+00]
    #  [9.61580776e+01]
    #  [-2.43123873e+01]
    #  [-2.47572621e+03]
    #  [1.59634684e+00]
    #  [1.00845873e+00]
    #  [7.66821129e-01]
    #  [1.09262847e+00]
    #  [2.80390196e+00]
    #  [1.24131951e+00]
    #  [9.45583653e+02]
    #  [1.74894265e+03]
    #  [5.49163632e+03]
    #  [-3.74922964e+03]
    #  [1.29464734e+01]
    #  [-5.17878865e+00]
    #  [-1.65935799e+01]
    #  [-3.57460592e+00]
    #  [2.61661503e+01]
    #  [8.32518606e+00]]
    # 与改进的随即梯度上升法不同：
    # [2.16342623e+02  6.75941056e+00  1.13458053e+03  2.62966594e+02
    #  2.25865811e+00 - 3.71721911e+00 - 7.06558071e+00 - 2.28140576e+00
    #  3.77389014e+00  2.07068231e+00  6.48377395e-01  9.37684449e+00
    #  - 3.63332065e+01 - 7.36624418e+02  1.00372305e+00 - 4.41795969e-01
    #  - 6.45795207e-01  6.22551836e-01  1.00055534e+00  9.30970745e-01
    #  2.37995214e+02  5.92267223e+00  1.05974442e+03 - 5.43164908e+02
    #  2.17278030e+00 - 1.66792125e+01 - 2.01526270e+01 - 4.47522451e+00
    #  2.12025978e+00  8.62117916e-01]
    # 于是要将它处理一下
    l = trainWeights1.tolist()
    modify = []
    for i in range(len(l)):
        modify.extend(l[i])
    mdfWeights = numpy.array((modify))
    # print(mdfTrainSets)
    # print(trainWeights1)

    trainWeights2 = stocGradAscent1(np.array(trainSets), trainLabels, 500)  # 随机梯度上升法
    print("未改进：")
    print(mdfWeights)
    print("改进：")
    print(trainWeights2)
    errCnt1 = 0
    errCnt2 = 0
    list1 = []
    list2 = []
    print("未改进：")
    for i in range(len(testSets)):
        # print(type(trainSets[i]), type(trainWeights1))
        if int(classifyVector(testSets[i], mdfWeights)) != int(testLabels[i]):
            print(i, end=" ")
            errCnt1 += 1
    print("")
    print("错误数：", errCnt1)
    errorRate1 = (float(errCnt1) / len(testSets)) * 100  # 错误率计算
    print("错误率：", errorRate1)

    print("改进：")
    for i in range(len(testSets)):
        # print(type(trainSets[i]), type(trainWeights1))
        if int(classifyVector(testSets[i], trainWeights2)) != int(testLabels[i]):
            print(i, end=" ")
            errCnt2 += 1
    # print(list1)
    # print(list2)
    print("")
    print("错误数：", errCnt2)
    errorRate2 = (float(errCnt2) / len(testSets)) * 100  # 错误率计算
    print("错误率：", errorRate2)

结果如下：

可以看到，用改进的方法虽然时间上会比原来的梯度上升法快，但是错误率却高很多。在这个数据集下，用来训练的不过是400条数据，而所消耗的时间却相差不多。在数据不多的情况下，显然是梯度上升法更好用。

四、总结

这次实验在代码上相比上次决策树还是比较简单的，难点还是在于理解里面的数学公式。想要彻底吃透还是得要多多复习数学的。