Logistic回归

不念不忘少年蓝@ 2021-10-01 12:20 502阅读 0赞

Logistic回归进行分类的思想：  
根据现有数据对分类边界建立回归公式，以此进行分类。

"""
    函数说明:梯度上升算法测试函数
    求函数f(x) = -x^2 + 4x的极大值
    """
    def Gradient_Ascent_test():
        def f_prime(x_old):                                    #f(x)的导数
            return -2 * x_old + 4
        x_old = -1                                            #初始值，给一个小于x_new的值
        x_new = 0                                            #梯度上升算法初始值，即从(0,0)开始
        alpha = 0.01                                        #步长，也就是学习速率，控制更新的幅度
        presision = 0.00000001                                #精度，也就是更新阈值
        while abs(x_new - x_old) > presision:
            x_old = x_new
            x_new = x_old + alpha * f_prime(x_old)            #上面提到的公式
        print(x_new)                                        #打印最终求解的极值近似值

跃阶函数瞬间跳跃过程难以处理，采用sigmoid函数

![5057999-f0b678f5916a7b93.png][]

import matplotlib.pyplot as plt
    import numpy as np
    
    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]))                                #添加标签
        fr.close()                                                            #关闭文件
        return dataMat, labelMat                                            #返回
    
    def plotDataSet():
        dataMat, labelMat = loadDataSet()                                    #加载数据集
        dataArr = np.array(dataMat)                                            #转换成numpy的array数组
        n = np.shape(dataMat)[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])    #1为正样本
            else:
                xcord2.append(dataArr[i,1]); ycord2.append(dataArr[i,2])    #0为负样本
        fig = plt.figure()
        ax = fig.add_subplot(111)                                            #添加subplot
        ax.scatter(xcord1, ycord1, s = 20, c = 'red', marker = 's',alpha=.5)#绘制正样本
        ax.scatter(xcord2, ycord2, s = 20, c = 'green',alpha=.5)            #绘制负样本
        plt.title('DataSet')                                                #绘制title
        plt.xlabel('x'); plt.ylabel('y')                                    #绘制label
        plt.show()                                                            #显示
    
    def sigmoid(inX):
        return 1.0 / (1 + np.exp(-inX))
    
    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为列数。
        alpha = 0.001  # 移动步长,也就是学习速率,控制更新的幅度。
        maxCycles = 500  # 最大迭代次数
        weights = np.ones((n, 1))
        for k in range(maxCycles):
            h = sigmoid(dataMatrix * weights)  # 梯度上升矢量化公式
            error = labelMat - h
            weights = weights + alpha * dataMatrix.transpose() * error
        return weights.getA()  # 将矩阵转换为数组，返回权重数组
    
    if __name__ == '__main__':
        dataMat, labelMat = loadDataSet()
        print(gradAscent(dataMat, labelMat))

输出结果

[[ 4.12414349]
     [ 0.48007329]
     [-0.6168482 ]]

画出数据集和最佳拟合函数的直线

def plotBestFit(weights):
        dataMat, labelMat = loadDataSet()  # 加载数据集
        dataArr = np.array(dataMat)  # 转换成numpy的array数组
        n = np.shape(dataMat)[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])  # 1为正样本
            else:
                xcord2.append(dataArr[i, 1]);
                ycord2.append(dataArr[i, 2])  # 0为负样本
        fig = plt.figure()
        ax = fig.add_subplot(111)  # 添加subplot
        ax.scatter(xcord1, ycord1, s=20, c='red', marker='s', alpha=.5)  # 绘制正样本
        ax.scatter(xcord2, ycord2, s=20, c='green', alpha=.5)  # 绘制负样本
        x = np.arange(-3.0, 3.0, 0.1)
        y = (-weights[0] - weights[1] * x) / weights[2]
        ax.plot(x, y)
        plt.title('BestFit')  # 绘制title
        plt.xlabel('X1');
        plt.ylabel('X2')  # 绘制label
        plt.show()
    
    if __name__ == '__main__':
        dataMat, labelMat = loadDataSet()
        weights = gradAscent(dataMat, labelMat)
        plotBestFit(weights)

![5057999-664802312f1db95b.png][]

--------------------

#### 训练算法：随机梯度上升 ####

梯度上升在处理大样本时计算复杂度太高。改进方法：一次仅用一个样本点来更新回归系数。由于可以再新样本到来时对分类器进行增量式更新，因而随机梯度上升算法是一个在线学习算法。

#随机梯度上升算法
    def stocGradAscent0(dataMatrix,classLabels):
        m,n=np.shape(dataMatrix)
        alpha=0.01
        weights=np.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

随机梯度上升和梯度上升很相似。区别：  
1.后者h,error都是向量，前者都是数值  
2.前者没有矩阵转换过程，所有数据都是Numpy数组

if __name__ == '__main__':
        dataMat, labelMat = loadDataSet()
        # weights = gradAscent(dataMat, labelMat)
        weights = stocGradAscent0(np.array(dataMat), labelMat)
        plotBestFit(weights)

![5057999-713491f3cd3b8e66.png][]

如图可以发现，分类器错过三分之一样本

--------------------

#### 改进的随机梯度上升算法 ####

def stocGradAscent1(dataMatrix,classLabels,numIter=150):
        m,n=np.shape(dataMatrix)
        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
                randIndex=int(np.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
    
    if __name__ == '__main__':
        dataMat, labelMat = loadDataSet()
        # weights = gradAscent(dataMat, labelMat)
        weights = stocGradAscent1(np.array(dataMat), labelMat)
        plotBestFit(weights)

![5057999-4590b6004fdab855.png][]

--------------------

对病马进行预测

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('\t')
            lineArr = []
            for i in range(len(currLine)-1):
                lineArr.append(float(currLine[i]))
            trainingSet.append(lineArr)
            trainingLabels.append(float(currLine[-1]))
        trainWeights = stocGradAscent1(np.array(trainingSet), trainingLabels, 500)        #使用改进的随即上升梯度训练
        errorCount = 0; numTestVec = 0.0
        for line in frTest.readlines():
            numTestVec += 1.0
            currLine = line.strip().split('\t')
            lineArr =[]
            for i in range(len(currLine)-1):
                lineArr.append(float(currLine[i]))
            if int(classifyVector(np.array(lineArr), trainWeights))!= int(currLine[-1]):
                errorCount += 1
        errorRate = (float(errorCount)/numTestVec) * 100                                 #错误率计算
        print("测试集错误率为: %.2f%%" % errorRate)
    
    
    
    if __name__ == '__main__':
        colicTest()

--------------------

#### sklearn下的Logistic ####

官网：[http://scikit-learn.org/dev/modules/generated/sklearn.linear\_model.LogisticRegression.html\#sklearn.linear\_model.LogisticRegression][http_scikit-learn.org_dev_modules_generated_sklearn.linear_model.LogisticRegression.html_sklearn.linear_model.LogisticRegression]

![5057999-2b5e21f7bb42e35b.png][]

solver：优化算法选择参数，只有五个可选参数，即newton-cg,lbfgs,liblinear,sag,saga。默认为liblinear。  
solver参数决定了我们对逻辑回归损失函数的优化方法，有四种算法可以选择，分别是：

*  liblinear：使用了开源的liblinear库实现，内部使用了坐标轴下降法来迭代优化损失函数。
 *  lbfgs：拟牛顿法的一种，利用损失函数二阶导数矩阵即海森矩阵来迭代优化损失函数。
 *  newton-cg：也是牛顿法家族的一种，利用损失函数二阶导数矩阵即海森矩阵来迭代优化损失函数。
 *  sag：即随机平均梯度下降，是梯度下降法的变种，和普通梯度下降法的区别是每次迭代仅仅用一部分的样本来计算梯度，适合于样本数据多的时候。
 *  saga：线性收敛的随机优化算法的的变重。

总结：

*  liblinear适用于小数据集，而sag和saga适用于大数据集因为速度更快。
 *  对于多分类问题，只有newton-cg,sag,saga和lbfgs能够处理多项损失，而liblinear受限于一对剩余(OvR)。啥意思，就是用liblinear的时候，如果是多分类问题，得先把一种类别作为一个类别，剩余的所有类别作为另外一个类别。一次类推，遍历所有类别，进行分类。
 *  newton-cg,sag和lbfgs这三种优化算法时都需要损失函数的一阶或者二阶连续导数，因此不能用于没有连续导数的L1正则化，只能用于L2正则化。而liblinear和saga通吃L1正则化和L2正则化。
 *  同时，sag每次仅仅使用了部分样本进行梯度迭代，所以当样本量少的时候不要选择它，而如果样本量非常大，比如大于10万，sag是第一选择。但是sag不能用于L1正则化，所以当你有大量的样本，又需要L1正则化的话就要自己做取舍了。要么通过对样本采样来降低样本量，要么回到L2正则化。
 *  从上面的描述，大家可能觉得，既然newton-cg, lbfgs和sag这么多限制，如果不是大样本，我们选择liblinear不就行了嘛！错，因为liblinear也有自己的弱点！我们知道，逻辑回归有二元逻辑回归和多元逻辑回归。对于多元逻辑回归常见的有one-vs-rest(OvR)和many-vs-many(MvM)两种。而MvM一般比OvR分类相对准确一些。郁闷的是liblinear只支持OvR，不支持MvM，这样如果我们需要相对精确的多元逻辑回归时，就不能选择liblinear了。也意味着如果我们需要相对精确的多元逻辑回归不能使用L1正则化了。  
    （摘自博客[http://cuijiahua.com/blog/2017/11/ml\_7\_logistic\_2.html][http_cuijiahua.com_blog_2017_11_ml_7_logistic_2.html]）

max\_iter：算法收敛最大迭代次数，int类型，默认为10。仅在正则化优化算法为newton-cg, sag和lbfgs才有用，算法收敛的最大迭代次数

代码：

from sklearn.linear_model import LogisticRegression
    
    def colicSklearn():
        frTrain = open('horseColicTraining.txt')  # 打开训练集
        frTest = open('horseColicTest.txt')  # 打开测试集
        trainingSet = []
        trainingLabels = []
        testSet = [];
        testLabels = []
        for line in frTrain.readlines():
            currLine = line.strip().split('\t')
            lineArr = []
            for i in range(len(currLine) - 1):
                lineArr.append(float(currLine[i]))
            trainingSet.append(lineArr)
            trainingLabels.append(float(currLine[-1]))
        for line in frTest.readlines():
            currLine = line.strip().split('\t')
            lineArr = []
            for i in range(len(currLine) - 1):
                lineArr.append(float(currLine[i]))
            testSet.append(lineArr)
            testLabels.append(float(currLine[-1]))
        classifier = LogisticRegression(solver='liblinear', max_iter=10).fit(trainingSet, trainingLabels)
        test_accurcy = classifier.score(testSet, testLabels) * 100
        print('正确率:%f%%' % test_accurcy)
    
    
    if __name__ == '__main__':
        colicSklearn()

[5057999-f0b678f5916a7b93.png]: /images/20210811/dc7e80e41085480a94c794ef970b4913.png
[5057999-664802312f1db95b.png]: /images/20210811/ee8569526fdc46b5b141c657c0a2c82e.png
[5057999-713491f3cd3b8e66.png]: /images/20210811/30852367c5f8483f962e6297af76261e.png
[5057999-4590b6004fdab855.png]: /images/20210811/a5c77d18c5d441e3930bed67af76724a.png
[http_scikit-learn.org_dev_modules_generated_sklearn.linear_model.LogisticRegression.html_sklearn.linear_model.LogisticRegression]: https://link.jianshu.com?t=http%3A%2F%2Fscikit-learn.org%2Fdev%2Fmodules%2Fgenerated%2Fsklearn.linear_model.LogisticRegression.html%23sklearn.linear_model.LogisticRegression
[5057999-2b5e21f7bb42e35b.png]: /images/20210811/101cefa9a3fd4330b49fe7a2795426d2.png
[http_cuijiahua.com_blog_2017_11_ml_7_logistic_2.html]: https://link.jianshu.com?t=http%3A%2F%2Fcuijiahua.com%2Fblog%2F2017%2F11%2Fml_7_logistic_2.html