逻辑回归原理

逻辑回归1.基本原理LogisticRegression和LinearRegression的原理是相似的，可以简单的描述为这样的过程：（1）找一个合适的预测函数，一般表示为h函数，该函数就是我们需要找的分类函数，它用来预测输入数据的判断结果。这个过程时非常关键的，需要对数据有一定的了解或分析，知道或者猜测预测函数的“大概”形式，比如是线性函数还是非线性函数。（2）构造一个Cost函数（损失函数），该函数表示预测的输出（h）与训练数据类别（y）之间的偏差，可以是二者之间的差（h-y）或者是其他的形式。综合考虑所有训练数据的“损失”，将Cost求和或者求平均，记为J(θ)函数，表示所有训练数据预测值与实际类别的偏差。（3）显然，J(θ)函数的值越小表示预测函数越准确（即h函数越准确），所以这一步需要做的是找到J(θ)函数的最小值。找函数的最小值有不同的方法，LogisticRegression实现时有的是梯度下降法（GradientDescent）。2.具体过程1.1构造预测函数LogisticRegression虽然名字里带“回归”，但是它实际上是一种分类方法，用于两分类问题（即输出只有两种）。根据步骤，需要先找到一个预测函数（h），显然，该函数的输出必须是两个值（分别代表两个类别），所以利用了Logistic函数（或称为Sigmoid函数），函数形式为：zez11)(g对应的函数图像是一个取值在0和1之间的S型曲线（图1）。图1接下来需要确定数据划分的边界类型，对于图2和图3中的两种数据分布，显然图2需要一个线性的边界，而图3需要一个非线性的边界。接下来我们只讨论线性边界的情况。图2图3对于线性边界的情况，边界形式如下：构造预测函数为：hθ(x)函数的值有特殊的含义，它表示结果取1的概率，因此对于输入x分类结果为类别1和类别0的概率分别为：3.2构造Cost函数AndrewNg在课程中直接给出了Cost函数及J(θ)函数如式（5）和（6），但是并没有给出具体的解释，只是说明了这个函数来衡量h函数预测的好坏是合理的。实际上这里的Cost函数和J(θ)函数是基于最大似然估计推导得到的。下面详细说明推导的过程。（4）式综合起来可以写成：取似然函数为：对数似然函数为：最大似然估计就是要求得使l(θ)取最大值时的θ，其实这里可以使用梯度上升法求解，求得的θ就是要求的最佳参数。但是，在AndrewNg的课程中将J(θ)取为（6）式，即：因为乘了一个负的系数-1/m，所以J(θ)取最小值时的θ为要求的最佳参数。3.3梯度下降法求J(θ)的最小值求J(θ)的最小值可以使用梯度下降法，根据梯度下降法可得θ的更新过程：式中为α学习步长，下面来求偏导：上式求解过程中用到如下的公式：因此，（11）式的更新过程可以写成：因为式中α本来为一常量，所以1/m一般将省略，所以最终的θ更新过程为：另外，补充一下，3.2节中提到求得l(θ)取最大值时的θ也是一样的，用梯度上升法求（9）式的最大值，可得：观察上式发现跟（14）是一样的，所以，采用梯度上升发和梯度下降法是完全一样的，这也是《机器学习实战》中采用梯度上升法的原因。3.4梯度下降过程向量化关于θ更新过程的vectorization，AndrewNg的课程中只是一带而过，没有具体的讲解。《机器学习实战》连Cost函数及求梯度等都没有说明，所以更不可能说明vectorization了。但是，其中给出的实现代码确是实现了vectorization的，图4所示代码的32行中weights（也就是θ）的更新只用了一行代码，直接通过矩阵或者向量计算更新，没有用for循环，说明确实实现了vectorization，具体代码下一章分析。文献[3]中也提到了vectorization，但是也是比较粗略，很简单的给出vectorization的结果为：且不论该更新公式正确与否，这里的Σ(...)是一个求和的过程，显然需要一个for语句循环m次，所以根本没有完全的实现vectorization，不像《机器学习实战》的代码中一条语句就可以完成θ的更新。下面说明一下我理解《机器学习实战》中代码实现的vectorization过程。约定训练数据的矩阵形式如下，x的每一行为一条训练样本，而每一列为不同的特称取值：约定待求的参数θ的矩阵形式为：先求x.θ并记为A：求hθ(x)-y并记为E：g(A)的参数A为一列向量，所以实现g函数时要支持列向量作为参数，并返回列向量。由上式可知hθ(x)-y可以由g(A)-y一次计算求得。再来看一下（15）式的θ更新过程，当j=0时：同样的可以写出θj，综合起来就是：综上所述，vectorization后θ更新的步骤如下：（1）求A=x.θ；（2）求E=g(A)-y；（3）求θ:=θ-α.x'.E,x'表示矩阵x的转置。也可以综合起来写成：y)-θ)(()1(-θθxgxmTα：前面已经提到过：1/m是可以省略的。4.代码分析图4中是《机器学习实战》中给出的部分实现代码。图4sigmoid函数就是前文中的g(z)函数，参数inX可以是向量，因为程序中使用了Python的numpy。gradAscent函数是梯度上升的实现函数，参数dataMatin和classLabels为训练数据，23和24行对训练数据做了处理，转换成numpy的矩阵类型，同时将横向量的classlabels转换成列向量labelMat，此时的dataMatrix和labelMat就是（18）式中的x和y。alpha为学习步长，maxCycles为迭代次数。weights为n维（等于x的列数）列向量，就是（19）式中的θ。29行的for循环将更新θ的过程迭代maxCycles次，每循环一次更新一次。对比3.4节最后总结的向量化的θ更新步骤，30行相当于求了A=x.θ和g(A)，31行相当于求了E=g(A)-y，32行相当于求θ:=θ-α.x'.E。所以这三行代码实际上与向量化的θ更新步骤是完全一致的。4.1代码示例使用鸢尾花数据集进行代码运行fromsklearn.datasetsimportload_irisimportmatplotlib.pyplotaspltimportnumpyasnpiris=load_iris()data=iris.datatarget=iris.target#printdata[:10]#printtarget[10:]X=data[0:100,[0,2]].y=target[0:100]printX[:5]printy[-5:]label=np.array(y)index_0=np.where(label==0)plt.scatter(X[index_0,0],X[index_0,1],marker='x',color='b',label='0',s=15)index_1=np.where(label==1)plt.scatter(X[index_1,0],X[index_1,1],marker='o',color='r',label='1',s=15)plt.xlabel('X1')plt.ylabel('X2')plt.legend(loc='upperleft')plt.show()逻辑回归模型的类：importnumpyasnpclasslogistic(object):def__init__(self):self.W=Nonedeftrain(self,X,y,learn_rate=0.01,num_iters=5000):num_train,num_feature=X.shape#inittheweightself.W=0.001*np.random.randn(num_feature,1).reshape((-1,1))loss=[]foriinrange(num_iters):error,dW=self.compute_loss(X,y)self.W+=-learn_rate*dWloss.append(error)ifi%200==0:print'i=%d,error=%f'%(i,error)returnlossdefcompute_loss(self,X,y):num_train=X.shape[0]h=self.output(X)loss=-np.sum((y*np.log(h)+(1-y)*np.log((1-h))))loss=loss/num_traindW=X.T.dot((h-y))/num_trainreturnloss,dWdefoutput(self,X):g=np.dot(X,self.W)returnself.sigmod(g)defsigmod(self,X):return1/(1+np.exp(-X))defpredict(self,X_test):h=self.output(X_test)y_pred=np.where(h=0.5,1,0)returny_pred训练测试，并且可视化跟踪的损失lossimportmatplotlib.pyplotasplty=y.reshape((-1,1))#addthex0=1one=np.ones((X.shape[0],1))X_train=np.hstack((one,X))classify=logistic()loss=classify.train(X_train,y)printclassify.Wplt.plot(loss)plt.xlabel('Iterationnumber')plt.ylabel('Lossvalue')plt.show()可视化决策边界label=np.array(y)index_0=np.where(label==0)plt.scatter(X[index_0,0],X[index_0,1],marker='x',color='b',label='0',s=15)index_1=np.where(label==1)plt.scatter(X[index_1,0],X[index_1,1],marker='o',color='r',label='1',s=15)#showthedecisionboundaryx1=np.arange(4,7.5,0.5)x2=(-classify.W[0]-classify.W[1]*x1)/classify.W[2]plt.plot(x1,x2,color='black')plt.xlabel('X1')plt.ylabel('X2')plt.legend(loc='upperleft')plt.show()4.1代码示例第二部分1.数据准备下面的数据来自《机器学习实战》中的示例：-0.01761214.0530640-1.3956344.6625411-0.7521576.5386200-1.3223717.15285300.42336311.05467700.4067047.06733510.66739412.7414520-2.4601506.86680510.5694119.5487550-0.02663210.4277430上面的数据一共是3列10行，其中前两列为x1和x2的值，第3列表示y的值；10行表示取了10个样本点。我们可以将这些数据当做训练模型参数的训练样本。见到训练样本就可以比较直观的理解算法的输入，以及我们如何利用这些数据来训练逻辑回归分类器，进而用训练好的模型来预测新的样本（检测样本）。从逻辑回归的参数形式，式子（1）我们可以看到逻辑回归模型中有两个待定参数a（x的系数）和b（常数项），我们现在给出来的数据有两个特征x1,x2，因此整个模型就增加了一项：ax1+cx2+b。为了形式上的统一，我们使用带下标的a表示不同的参数（a0表示常数项b并作x0的参数x0=1，a1、a2分别表示x1和x2的参数），就可以得到：$$a_0x_0+a_1x_1+a_2x_2$$这样统一起来后，就可以使用矩