【机器学习实战——第5章】：Logistic回归

梯度下降法

• 假设有$m$组样本$((x^1,x^2),y)$
• 线性回归方程为：
  $$Y=w_0\cdot 1+w_1{x}^1+w_2{x}^2$$
• 目标为最小化损失函数$f_m$
  $$\begin{gathered} f_m=\sum _{i=1}^m(Y_i-y_i{)}^2 \\ =\sum _{i=1}^m{Y}_i^2-2Y_i{y}_i+y_i^2 \end{gathered}$$
• 每次迭代下降的方向为：
  $$\Delta w_j=\frac{\delta f_m}{\delta w_j}=2\sum _{i=1}^m(Y_i-y_i)x_i^j$$
  当$j=0$时：
  $$\Delta w_0=2\sum _{i=1}^m(Y_i-y_i)\cdot 1$$
  $$\begin{gathered} \text{梯度上升 }w:=w+\alpha \Delta w \\ \text{梯度下降 }w:=w-\alpha \Delta w \end{gathered}$$
• 对于本例$x$是三维向量的情况：
  $$w=w-\left(\begin{array}{ccccc}1& 1& 1& \cdots & 1\\ x_1^1& x_2^1& x_3^1& \cdots & x_m^1\\ x_1^2& x_2^2& x_3^2& \cdots & x_m^2\end{array}\right)\cdot \left(\begin{array}{c}{Y}_1-y_1\\ Y_2-y_2\\ ⋮\\ Y_m-y_m\end{array}\right)$$

Logistic回归

• Sigmoid函数
  $$\sigma (x)=\frac{1}{1+e^{-x}}$$
  
• 对于二分类问题，上文提到的损失函数并非凸函数，所以需要重新定义损失函数，该函数是凸函数
  $$J(\theta )=\frac{1}{m}\sum _{i=1}^{m}Cost(h_{\theta }(x^i),y^i)$$
  $$Cost(h_{\theta }(x^i),y^i)=\{\begin{array}{rc}-log(h_{\theta }(x))& y=1\\ -log(1-h_{\theta }(x))& y=0\end{array}$$
  $$J(\theta )=-\frac{1}{m}[\sum _{i=1}^m{y}^{i}log(h_{\theta }(x^i))+(1-y^i)log(1-h_{\theta }(x^i))]$$

训练算法：使用梯度下降找到最佳参数

#Logistic回归梯度上升优化算法
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]))
    return dataMat,labelMat

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

def gradAscent(dataMatIn,classLabels):
    dataMatrix=np.mat(dataMatIn)
    labelMat=np.mat(classLabels).transpose()
    m,n=np.shape(dataMatrix)
    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

分析数据：画出决策边界

def plotBestFit(weights):
    import matplotlib.pyplot as plt
    dataMat,labelMat=loadDataSet()
    dataArr=np.array(dataMat)
    n=np.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=np.arange(-3.0,3.0,0.1)
    y=(-weights[0]-weights[1]*x)/weights[2]
    ax.plot(x,y)
    plt.xlabel('X1'); plt.ylabel('X2');
    plt.show()

随机梯度下降算法

梯度上升算法每次更新回归系数时，要遍历整个数据集，时间复杂度较高。一种改进方法是一次只用一个样本点来更新回归系数，该方法称为随机梯度上升算法。

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

def stocGradAscent0(dataMatrix,classLabels):
    m,n=np.shape(dataMatrix)
    alpha=0.1
    weights=np.ones(n)
    for i in range(m):
        h=sigmoid(np.sum(dataMatrix[i]*weights))
        error=classLabels[i]-h
        weights=weights+alpha*error*dataMatrix[i]
    return weights

该拟合直线并非最佳分类线

改进的随机梯度下降算法

def stocGradAscent1(dataMatrix,classLabels,numIter=150):
    m,n=np.shape(dataMatrix)
    weights=np.ones(n)
    for j in range(numIter):
        dataIndex=np.array(range(m))
        for i in range(m):
            alpha=4/(1.0+i+j)+0.01
            randIndex=int(np.random.uniform(0,m))
            h=sigmoid(sum(dataMatrix[randIndex]*weights))
            error=classLabels[randIndex]-h
            weights=weights+alpha*error*dataMatrix[randIndex]
            np.delete(dataIndex,randIndex)
    return weights