前言
最近因为要做一个预测性的工程从而学到了逻辑回归（Logistic Regression），看了很多的资料以及一些论坛上的博文，发现了很多不错的资料，也发现了一些不足甚至是错误的地方。因此，在这里我重新从自己的视角以及理解写一下我对逻辑回归的认识，纠正一些科普资料中的错误，再加上自己的学习过程。
如果有错误或者不足的地方欢迎看到的人批评指正，相互交流。
原理及用途
回归从表现出形式上讲就是给定一些点集，找出一条曲线去拟合。从数学公式的角度讲，就是在已经大概知道函数形式的时候，用给出的训练样本确定函数中未知的参数。所以在进行回归时，根据已有经验找到一个合适形式的预测函数非常重要。假如说预测函数的分布是非线性的，而我们选择了一个线性函数来预测，就会造成较大的偏差。其实逻辑回归就是在线性回归的基础上加了一个逻辑函数，但是这个逻辑函数非常重要，造就了逻辑回归在分类问题中的重要地位。
逻辑回归经常应用于分类问题，特别是“二分类”。就是预测结果只有两种“0”或者“1”。具体的用途有：
- 分析行为：根据某人的行为习惯分析用户某个行为会不会发生。
- 判断预测：根据建立的模型，分析在某些自变量集合下，是否会发生某种疾病。
实现步骤
- 首先是找到一个合适形式的预测函数，称为 hθ 函数，即hypothesis。
- 构造损失函数，主要分为(1)0-1损失函数，称为Cost函数，(2)平方损失函数，称为 j(θ) 。
- j(θ) 函数的值越小表示预测的函数越准确，所以用logistic regression的方法求出 j(θ) 取最小值时函数的参数 θ 。
具体实现
构造预测函数
由于是“二分类”问题，所以这里要用到Logistic函数，也叫做Sigmoid函数，函数的形式为
g(z)=11+e−z
Sigmoid函数的推导过程是：令 y=logp1−p
y=log(p/(1−p))⇔ey=p/(1−p)⇔ey=p+eyp⇔p=ey1+ey⇔p=11+e−y
该函数对应的图像是：
从图像可以看出p在[0,1]之间，在解决二分类问题的时候可以设定一个阈值，当p>阈值的时候分为类别A，当P<阈值的时候分为类别B。
推理过程
下面分析一下可能会遇到的问题，然后分析预测函数h的形式：
线性边界：
首先以线性边界的情况来讨论，线性边界的边界形式为：
θ0+θ1x1+⋯+θnxn=∑i=0nθixi=θTx
构造预测函数 hθ :
hθ=g(θTx)=11+e−θTx
hθ 表示预测结果取1的概率，那么预测结果取0的概率可以表示为 1−hθ ：
P(y=1|x;θ)=hθ(x)P(y=0|x;θ)=1−hθ(x)
该式还可以写成 P(y|x;θ)=(hθ(x))y(1−hθ(x))1−y 。
使用最大似然估计：
L(θ)=∏i=1m(hθ(x(i)))y(i)(1−hθ(x(i)))1−y(i)
两边同时取对数：
l(θ)=logL(θ)=∑i=1m(y(i)hθ(x(i))+(1−y(i))(1−hθ(x(i))))
当 l(θ) 取得最大值时的 θ 就是最佳的 θ 的取值。
根据上式构造0-1损失函数和平方损失函数：
Cost(hθ(x),y)={−loghθ(x)−log(1−hθ(x))y=1y=0
J(θ)=1m∑i=1m(Cost(hθ(x(i))),y(i))=−1m∑i=1m(y(i)loghθ(x(i))+(1−y(i))(1−log(1−hθ(x(i))))
J(θ)=−1ml(θ) ，所以求 l(θ) 的最大值就是求 J(θ) 的最小值。
在这里使用梯度下降的方法求 J(θ) 的最小值，关于梯度下降的介绍见，梯度下降是一个不断地迭代求最小值的方法。简单的形容就是：
WantminθJ(θ)Repeat{θj:=θj−α∂∂θjj(θ)}
其中 θj:=θj−α∂∂θjJ(θ) 就是梯度下降法的公式， θj: 是 θj 经过迭代后的下一次状态。整个迭代过程终止的条件是 ∂∂θjJ(θ) 的结果趋近于0。
化简 ∂∂θjJ(θ) :
∂∂θjJ(θ)=−1m∑i=1m[y(i)1hθ(x(i))∂∂θjhθ(x(i))−(1−y(i))11−hθ(x(i))∂∂θjhθ(x(i))]=−1m∑i=1m(y(i)1g(θTx(i))−(1−y(i))11−g(θTx(i)))∂∂θjg(θTx(i))=−1m∑i=1m(y(i)1g(θTx(i))−(1−y(i))11−g(θTx(i)))g(θTx(i))(1−g(θTx(i)))∂∂θj(θTx(i))=−1m∑i=1my(i)(1−g(θTx(i)))−(1−y(i))g(θTx(i))g(θTx(i))(1−g(θTx(i)))g(θTx(i))(1−g(θTx(i)))∂∂θjθTx(i)(j)=−1m∑i=1m(y(i)−g(θTx(i)))x(i)(j)=1m∑i=1m(hθ(x(i))−y(i))x(i)(j)
上述过程都属于微积分中求偏导的过程，不熟悉的可以再看下微积分的知识。
求偏导过程中要用到的知识：
h(x)∂∂xh(x)=11+eg(x)=1(1+eg(x))2eg(x)∂∂xg(x)=h(x)(1−h(x))∂∂xg(x)
因为 α 是一个常数，所以 α1m 仍然是一个常数，所以式中的 1m 可以省略。所以 θj:=θj−α∂∂θjJ(θ) 可以化简为：
θj:=θj−α∑i=1m(hθ(x(i))−y(i))x(i)(j);(j=0,1,…,n)
上式就是梯度下降法求解的最终公式：
WantminθJ(θ)Repeat{θj:=θj−α∑i=1m(hθ(x(i))−y(i))x(i)(j))}
对于训练样本：
X=⎡⎣⎢⎢⎢⎢⎢⎢x(1)0x(2)0⋮x(m)0x(1)1x(2)1⋮x(m)1……⋮…x(1)nx(2)n⋮x(m)n⎤⎦⎥⎥⎥⎥⎥⎥,Y=⎡⎣⎢⎢⎢⎢⎢y(1)y(2)⋮y(m)⎤⎦⎥⎥⎥⎥⎥
可以求出 θ 的值：
⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪θ0=θ0−α∑mi=1(hθ(x(i))−y(i))x(i)(0))θ1=θ1−α∑mi=1(hθ(x(i))−y(i))x(i)(1))……………………θm=θm−α∑mi=1(hθ(x(i))−y(i))x(i)(m))
至此逻辑回归过程结束。
对于具体的训练样本数据，逻辑回归所起到的作用如表所示：
| 数据 | 分类 | 预测结果1 | 分类结果1 | 预测结果2 | 分类结果2 |
|---|---|---|---|---|---|
| 10 | 1 | 1 | 1 | 0.1 | |
| 9 | 1 | 0.9 | 1 | 0.9 | 1 |
| 8 | 1 | 0.8 | 1 | 0.8 | 1 |
| 7 | 1 | 0.7 | 1 | 0.7 | 1 |
| 6 | 1 | 0.5 | 0.6 | 1 | |
| 5 | 0.6 | 1 | 0.5 | ||
| 4 | 0.4 | 0.4 | |||
| 3 | 0.3 | 0.3 | |||
| 2 | 0.2 | 0.2 | |||
| 1 | 0.1 | 1 | 1 |
并且可以看出训练结果1的效果要好于训练结果2。
代码实现
import matplotlib.pyplot as plt
from numpy import *
def loadDataSet():
dataMat = []; labelMat = []
fr = open(r'C:\Users\crystal\Desktop\testSet-LR.txt')
for line in fr.readlines():
lineArr = line.strip().split()
dataMat.append([, float(lineArr[]), float(lineArr[])])
labelMat.append(int(lineArr[]))
return dataMat,labelMat
def sigmoid(inX):
return /(+exp(-inX))
def gradAscent(dataMatIn, classLabels):
dataMatrix = mat(dataMatIn) #convert to NumPy matrix
labelMat = mat(classLabels).transpose() #convert to NumPy matrix
m,n = shape(dataMatrix)
alpha =
maxCycles =
weights = ones((n,))
for k in range(maxCycles): #heavy on matrix operations
h = sigmoid(dataMatrix*weights) #matrix mult
error = (labelMat - h) #vector subtraction
weights = weights + alpha * dataMatrix.transpose()* error #matrix mult
return weights
def GetResult():
dataMat,labelMat=loadDataSet()
weights=gradAscent(dataMat,labelMat)
print (weights)
plotBestFit(weights)
def plotBestFit(weights):
dataMat,labelMat=loadDataSet()
dataArr = array(dataMat)
n = shape(dataArr)[]
xcord1 = []; ycord1 = []
xcord2 = []; ycord2 = []
for i in range(n):
if int(labelMat[i])== :
xcord1.append(dataArr[i,]); ycord1.append(dataArr[i,])
else:
xcord2.append(dataArr[i,]); ycord2.append(dataArr[i,])
fig = plt.figure()
ax = fig.add_subplot()
ax.scatter(xcord1, ycord1, s=, c='red', marker='s')
ax.scatter(xcord2, ycord2, s=, c='green')
x = arange(-, , )
y=(*x+)/()
# y = (-weights[0]-weights[1]*x)/weights[2]
ax.plot(x,y)
plt.xlabel('X1'); plt.ylabel('X2');
plt.show()
if __name__=='__main__':
GetResult()
训练结果
参考资料：
1.http://www.360doc.com/content/14/0912/08/14875906_408823529.shtml
2.http://blog.csdn.net/pakko/article/details/37878837
3.http://blog.csdn.net/dongtingzhizi/article/details/15962797
4.http://blog.csdn.net/abcjennifer/article/details/7716281
5.http://blog.csdn.net/woxincd/article/details/7040944