Logistic回归
1. 基于Logistic回归和Sigmoid函数的分类
- logistic回归又称logistic回归分析，是一种广义的线性回归分析模型，常用于数据挖掘，疾病自动诊断，经济预测等领域。例如，探讨引发疾病的危险因素，并根据危险因素预测疾病发生的概率等。
- logistic回归的因变量可以是二分类的，也可以是多分类的，但是二分类的更为常用，也更加容易解释。
假设现在有一些数据点，我们利用一条直线对这些点进行拟合(该线称为最佳拟合直线)，这个拟合过程就称作为回归，如下图所示：
利用Logistic回归进行分类的主要思想是：根据现有数据对分类边界线建立回归公式，依次进行分类。“回归”源于最佳拟合，表示要找到最佳拟合参数集。训练分类器时的做法就是寻找最佳拟合参数，使用的是最优化算法，后面将会介绍，其中包括基本的梯度上升算法和一个改进的随机梯度上升法。
我们想要的函数应该是，能接受所有的输入然后预测出类别。例如海维塞德阶跃函数，或直接称为单位阶跃函数，图像如下：
该函数的问题在于不连续、不可微，该函数在跳跃点上从0瞬间跳跃到1，有时很难处理，而Sigmoid函数也有这样的性质，且数学上更易处理。
Logistic回归一种二分类算法，它利用的是Sigmoid函数阈值在[0,1]这个特性。其实，Logistic本质上是一个基于条件概率的判别模型(Discriminative Model)。所以我们需要先了解Sigmoid函数 ，我们也可以称它为Logistic函数。具体计算公式如下：
σ
(
z
)
=
1
1
+
e
?
z
\sigma (z)=\frac{1}{1+e^{-z}}
σ(z)=1+e?z1?
下图给出了Sigmoid函数在不同坐标尺度下的两条曲线图。当x为0时，Sigmoid函数值为0.5，随着x增大，对应的sigmoid值将逼近与1；而随着x减小，Sigmoid值将逼近于0。如果横坐标刻度足够大（第二张图），Sigmoid函数看起来很像一个阶跃函数。
因此，为了实现Logistic回归分类器，我们可以在每个特征上都乘以一个回归系数，然后把所有的结果值相加，将这个总和代入Sigmoid函数中，进而得到一个范围在0~1之间的数值。任何大于0.5的数据被分入1类，小于0.5即被归入0类。所以，Logistic回归也可以被看成是一种概率估计。
Sigmoid函数的输入记为z，
z
=
w
T
x
+
b
z = w^Tx + b
z=wTx+b，其中向量x是分类器的输入数据，向量
w
w
w也就是我们要找到的最佳系数，b为常数。代入Sigmoid函数中得到
y
=
1
1
+
e
?
(
w
T
x
+
b
)
y=\frac{1}{1+e^{-(w^Tx+b)}}
y=1+e?(wTx+b)1?该预测函数的值表示y=1的概率，所以有y=1和y=0的概率分别为
p
(
y
=
1
∣
x
)
=
e
w
T
x
+
b
1
+
e
w
T
x
+
b
p(y=1|x)=\frac{e^{w^Tx+b}}{1+e^{w^Tx+b}}
p(y=1∣x)=1+ewTx+bewTx+b?
p
(
y
=
0
∣
x
)
=
1
1
+
e
w
T
x
+
b
p(y=0|x)=\frac{1}{1+e^{w^Tx+b}}
p(y=0∣x)=1+ewTx+b1?
我们便可以输入数据代入z中再由Sigmoid函数映射到0~1之间，得到一个这样的数值之后就可以进行分类，大于0.5的分入1类，反之，分入0类。
现在的问题就在于最佳回归系数是多少？即 w w w和b。
2. 基于最优化方法的最佳回归系数确定
2.1 极大似然估计
- 确定待求解的未知参数 θ 1 , θ 2 , . . . . . . , \theta_1,\theta_2,......, θ1?,θ2?,......,如均值、方差或特定分布函数的参数等
- 计算每个样本 X 1 , X 2 , . . . . . . , X n X_1,X_2,......,X_n X1?,X2?,......,Xn?的概率密度 f ( X i ; θ 1 , . . . . . . , θ m ) f(X_i;\theta_1,......,\theta_m) f(Xi?;θ1?,......,θm?)
- 根据样本的概率密度累乘构造似然函数： L ( θ 1 , . . . , θ m ) = ∏ f ( x i ; θ 1 , . . . , θ m ) L(\theta_1,...,\theta_m)=\prod f(x_i;\theta_1,...,\theta_m) L(θ1?,...,θm?)=∏f(xi?;θ1?,...,θm?)
- 通过似然函数最大化（求导为0），求解未知参数 θ \theta θ（为了降低计算难度，可采用对数加法替换概率乘法，通过导数为0/极大值来求解未知参数）
根据以上方法，我们可以得到最大化对数似然函数
l
(
w
,
b
)
=
∑
i
=
1
m
l
n
p
(
y
i
∣
x
i
;
w
i
,
b
)
l(w,b)=\sum_{i=1}^{m}lnp(y_i|x_i;w_i,b)
l(w,b)=i=1∑m?lnp(yi?∣xi?;wi?,b)
记
β
=
(
w
;
b
)
,
x
^
=
(
x
;
1
)
\beta=(w;b),\hat{x}=(x;1)
β=(w;b),x^=(x;1)，则
w
T
x
+
b
w^Tx+b
wTx+b可简写为：
β
T
x
^
\beta^T\hat{x}
βTx^
令
p
1
(
x
i
^
;
β
)
=
p
(
y
=
1
∣
x
^
;
β
)
p_1(\hat{x_i};\beta)=p(y=1|\hat{x};\beta)
p1?(xi?^?;β)=p(y=1∣x^;β)
p
0
(
x
i
^
;
β
)
=
p
(
y
=
0
∣
x
^
;
β
)
=
1
?
p
1
(
x
i
^
;
β
)
p_0(\hat{x_i};\beta)=p(y=0|\hat{x};\beta)=1-p_1(\hat{x_i};\beta)
p0?(xi?^?;β)=p(y=0∣x^;β)=1?p1?(xi?^?;β)则合二为一，当
y
i
y_i
yi?=1时，第二项为0；当
y
i
y_i
yi?=0时，第一项为0
p
(
y
i
∣
x
i
;
w
i
,
b
)
=
y
i
p
1
(
x
i
^
;
β
)
+
(
1
?
y
i
)
p
0
(
x
i
^
;
β
)
p(y_i|x_i;w_i,b)=y_ip_1(\hat{x_i};\beta)+(1-y_i)p_0(\hat{x_i};\beta)
p(yi?∣xi?;wi?,b)=yi?p1?(xi?^?;β)+(1?yi?)p0?(xi?^?;β)
根据Sigmoid函数，代入求得，似然函数可重写为：
l
(
β
)
=
∑
i
=
1
m
(
l
n
(
y
i
e
z
1
+
e
z
+
(
1
?
y
i
)
1
1
+
e
z
)
)
l(\beta)=\sum_{i=1}^{m}(ln(y_{i}\frac{e^{z}}{1+e^{z}}+(1-y_{i})\frac{1}{1+e^{z}}))
l(β)=i=1∑m?(ln(yi?1+ezez?+(1?yi?)1+ez1?))
通分得到，
l
(
β
)
=
∑
i
=
1
m
l
n
(
y
i
e
z
+
(
1
?
y
i
)
1
+
e
z
)
l(\beta)=\sum_{i=1}^{m}ln(\frac{y_ie^z+(1-y_i)}{1+e^z})
l(β)=i=1∑m?ln(1+ezyi?ez+(1?yi?)?)
根据对数性质得到，
l
(
β
)
=
∑
i
=
1
m
(
l
n
(
y
i
e
z
+
(
1
?
y
i
)
)
?
l
n
(
1
+
e
z
)
)
l(\beta)=\sum_{i=1}^{m}(ln(y_ie^z+(1-y_i))-ln(1+e^z))
l(β)=i=1∑m?(ln(yi?ez+(1?yi?))?ln(1+ez))
由于
y
i
?
y_i\epsilon
yi??{0,1}，所以当
y
i
y_i
yi?=1时，减号左边等于z；当
y
i
y_i
yi?=0时，左边等于0，所以可以转换为：
y
i
z
y_iz
yi?z
再对括号里的部分取负号并将
z
=
w
T
x
+
b
=
β
T
x
^
z = w^Tx + b=\beta^T\hat{x}
z=wTx+b=βTx^代入，等价于最后要求最大化该方程：
l
(
β
)
=
∑
i
=
1
m
(
y
i
β
T
x
i
^
?
l
n
(
1
+
e
β
T
x
i
^
)
)
l(\beta)=\sum_{i=1}^{m}(y_i\beta^T\hat{x_i}-ln(1+e^{\beta^T\hat{x_i}}))
l(β)=i=1∑m?(yi?βTxi?^??ln(1+eβTxi?^?))
2.2. 梯度下降算法
根据上面的方程，我们应该去求方程的最小值，然后得到此时的最佳参数是多少。求函数的最小值我们可以用梯度下降算法，迭代公式为： w : = w ? α ▽ w f ( w ) w:=w-\alpha \triangledown_wf(w) w:=w?α▽w?f(w)
2.3. 梯度上升算法
这里我们选用《机器学习实战》中的方法，**（梯度下降和梯度上升其实差不多，取决于你要求函数的最大值还是最小值）**梯度上升法基于的思想是：要找到某函数的最大值，最好的方法是沿着该函数的梯度方向探寻。如果梯度记为
▽
\triangledown
▽，则函数
f
(
x
,
y
)
f(x,y)
f(x,y)的梯度由下式表示：
▽
f
(
x
,
y
)
=
[
?
f
(
x
,
y
)
?
x
?
f
(
x
,
y
)
?
y
]
\triangledown f(x,y)=\begin{bmatrix}\frac{\partial f(x,y)}{\partial x}\\ \frac{\partial f(x,y)}{\partial y} \end{bmatrix}
▽f(x,y)=[?x?f(x,y)??y?f(x,y)??]
这个梯度意味着要沿x的方向移动
?
f
(
x
,
y
)
?
x
\frac{\partial f(x,y)}{\partial x}
?x?f(x,y)?，沿y的方向移动
?
f
(
x
,
y
)
?
y
\frac{\partial f(x,y)}{\partial y}
?y?f(x,y)?。其中，函数
f
(
x
,
y
)
f(x,y)
f(x,y)必须要在待计算的点上有定义并且可微。如下图：
上图展示的，梯度上升算法到达每个点后都会重新估计移动的方向。从 P0 开始，计算完该点的梯度，函数就根据梯度移动到下一点 P1。在 P1 点，梯度再次被重新计算，并沿着新的梯度方向移动到 P2 。如此循环迭代，直到满足停止条件。迭代过程中，梯度算子总是保证我们能选取到最佳的移动方向。
上图中的梯度上升算法沿梯度方向移动了一步。可以看到，梯度算子总是指向函数值增长最快的方向。这里所说的是移动方向，而未提到移动量的大小。该量值称为步长，记作 α 。用向量来表示的话，梯度上升算法的迭代公式如下:
w
:
=
w
+
α
▽
w
f
(
w
)
w:=w+\alpha \triangledown_wf(w)
w:=w+α▽w?f(w)
该公式将一直被迭代执行，直至达到某个停止条件为止，比如迭代次数达到某个指定值或算法达到某个可以允许的误差范围。
所以，接下来我们要做的就是对
l
(
β
)
l(\beta)
l(β)求偏导：
?
?
β
i
l
(
β
)
=
∑
i
=
1
m
(
y
i
?
p
1
(
x
i
;
β
^
)
)
x
i
^
\frac{\partial }{\partial \beta _{i}}l(\beta)=\sum_{i=1}^{m}(y_i-p_1(\hat{x_i;\beta}))\hat{x_i}
?βi???l(β)=i=1∑m?(yi??p1?(xi?;β^?))xi?^?
则最终迭代公式写成：
w
i
:
=
w
i
+
α
∑
i
=
1
m
(
y
i
?
p
1
(
x
i
;
β
^
)
)
x
i
^
w_i:=w_i+\alpha \sum_{i=1}^{m}(y_i-p_1(\hat{x_i;\beta}))\hat{x_i}
wi?:=wi?+αi=1∑m?(yi??p1?(xi?;β^?))xi?^?
2.3.1. 数据准备
一个简单的数据集：
可以看出数据具有二维特征，因此可以将数据在一个二维平面上展示出来。我们可以将第一列数据(X1)看作x轴上的值，第二列数据(X2)看作y轴上的值。而最后一列数据即为分类标签。根据标签的不同，对这些点进行分类。
编写代码来看看分类情况：
import numpy as np
import matplotlib.pyplot as plt
"""
函数说明:加载数据
Parameters:
无
Returns:
dataMat - 数据列表
labelMat - 标签列表
"""
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 #返回
"""
函数说明:绘制数据集
Parameters:
无
Returns:
无
"""
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() #显示
if __name__ == '__main__':
plotDataSet()
运行结果如下，从图可以看出数据的分布情况。假设Sigmoid函数的输入记为z，那么z=w0x0 + w1x1 + w2x2，即可将数据分割开。其中，x0为全是1的向量，x1为数据集的第一列数据，x2为数据集的第二列数据。另z=0，则0=w0 + w1x1 + w2x2。横坐标为x1，纵坐标为x2。这个方程未知的参数为w0，w1，w2，也就是我们需要求的回归系数(最优参数)。
2.3.2. 训练算法：使用梯度上升找到最佳参数
接下来我们将采用梯度上升法找到Logistic回归分类器在此数据集上的最佳回归系数，代码如下：
"""
函数说明:加载数据
Parameters:
无
Returns:
dataMat - 数据列表
labelMat - 标签列表
"""
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 #返回
"""
函数说明:sigmoid函数
Parameters:
inX - 数据
Returns:
sigmoid函数
"""
def sigmoid(inX):
return 1.0 / (1 + np.exp(-inX))
"""
函数说明:梯度上升算法
Parameters:
dataMatIn - 数据集
classLabels - 数据标签
Returns:
weights.getA() - 求得的权重数组(最优参数)
"""
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))
运行结果如下，由此我们就得到了回归系数[w0,w1,w2]，我们可以由此参数确定决策边界，画出分割线。
2.3.3. 分析数据：画出决策边界
得到了回归系数，我们就可以开始画出决策边界，代码如下：
"""
函数说明:加载数据
Parameters:
无
Returns:
dataMat - 数据列表
labelMat - 标签列表
"""
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 #返回
"""
函数说明:sigmoid函数
Parameters:
inX - 数据
Returns:
sigmoid函数
"""
def sigmoid(inX):
return 1.0 / (1 + np.exp(-inX))
"""
函数说明:梯度上升算法
Parameters:
dataMatIn - 数据集
classLabels - 数据标签
Returns:
weights.getA() - 求得的权重数组(最优参数)
"""
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() #将矩阵转换为数组，返回权重数组
"""
函数说明:绘制数据集
Parameters:
weights - 权重参数数组
Returns:
无
"""
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)
运行结果如下，从分类结果可以看出，分错的点只有个别几个，效果整体上来说还是很不错的。但是尽管数据集简单，还是需要300次的乘法计算。
2.4. 随机梯度上升算法
假设，我们使用的数据集一共有100个样本。那么，dataMatrix就是一个1003的矩阵。每次计算h的时候，都要计算dataMatrixweights这个矩阵乘法运算，要进行1003次乘法运算和1002次加法运算。同理，更新回归系数(最优参数)weights时，也需要用到整个数据集，要进行矩阵乘法运算。总而言之，该方法处理100个左右的数据集时尚可，但如果有数十亿样本和成千上万的特征，那么该方法的计算复杂度就太高了。因此，需要对算法进行改进，我们每次更新回归系数(最优参数)的时候，能不能不用所有样本呢？一次只用一个样本点去更新回归系数(最优参数)？这样就可以有效减少计算量了，这种方法就叫做随机梯度上升算法。
代码如下：
"""
函数说明:随机梯度上升算法
Parameters:
dataMatrix - 数据数组
classLabels - 数据标签
Returns:
weights - 求得的回归系数数组(最优参数)
"""
def stocGradAscent0(dataMatrix, classLabels):
m,n = shape(dataMatrix)
alpha = 0.01
weights = ones(n) #initialize to all ones
for i in range(m):
h = sigmoid(sum(dataMatrix[i]*weights))
error = classLabels[i] - h
weights = weights + alpha * error * dataMatrix[i]
return weights
可以看到，随机梯度上升算法和梯度上升算法在代码上很相似，但也有一些区别：第一，后者的变量h和误差error都是向量，而前者则全是数值；第二，前者没有矩阵的转换过程，所有变量的数据类型都是Numpy数组。
测试一下：
if __name__ == '__main__':
dataMat, labelMat = loadDataSet()
weights = stocGradAscent0(array(dataMat), labelMat)
plotBestFit(weights)
运行结果如下，效果就比较一般了，错分了将近三分之一。
2.5. 改进的随机梯度上升算法
相比于梯度上升算法，增加了两处代码来进行改进。代码如下：
"""
函数说明:改进的随机梯度上升算法
Parameters:
dataMatrix - 数据数组
classLabels - 数据标签
numIter - 迭代次数
Returns:
weights - 求得的回归系数数组(最优参数)
"""
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]) #删除已经使用的样本
return weights #返回
该算法第一个改进之处在于，alpha在每次迭代的时候都会调整，并且，虽然alpha会随着迭代次数不断减小，但永远不会减小到0，因为这里还存在一个常数项。必须这样做的原因是为了保证在多次迭代之后新数据仍然具有一定的影响。如果需要处理的问题是动态变化的，那么可以适当加大上述常数项，来确保新的值获得更大的回归系数。另一点值得注意的是，在降低alpha的函数中，alpha每次减少1/(j+i)，其中j是迭代次数，i是样本点的下标。第二个改进的地方在于跟新回归系数(最优参数)时，只使用一个样本点，并且选择的样本点是随机的，每次迭代不使用已经用过的样本点。这样的方法，就有效地减少了计算量，并保证了回归效果。
接下来，让我们看看改进后的随机梯度上升算法的分类效果：
代码如下：
if __name__ == '__main__':
dataMat, labelMat = loadDataSet()
weights = stocGradAscent1(array(dataMat), labelMat)
plotBestFit(weights)
运行结果如下，效果和梯度上升算法已经差不多了，但是用的计算量更少。改进后的随机梯度上升算法没有出现周期性波动，收敛得更快。可从下面的回归系数和迭代次数的关系图看出。
3. 从疝气病症状预测病马的死亡率
3.1. 数据准备
horseColicTest.txt和horseColicTraining.txt
horseColicTest.txt
horseColicTraining.txt
3.2. 基于改进后的随机梯度上升算法的Logistic分类器
"""
函数说明:使用Python写的Logistic分类器做预测
Parameters:
无
Returns:
无
"""
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)
"""
函数说明:分类函数
Parameters:
inX - 特征向量
weights - 回归系数
Returns:
分类结果
"""
def classifyVector(inX, weights):
prob = sigmoid(sum(inX*weights))
if prob > 0.5: return 1.0
else: return 0.0
if __name__ == '__main__':
colicTest()
错误率较高的原因是因为数据量较小，只有几百。
接下来我们试试梯度上升算法，看看效果如何？
3.3. 基于梯度上升算法的Logistic分类器
"""
函数说明:使用Python写的Logistic分类器做预测
Parameters:
无
Returns:
无
"""
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 = gradAscent(np.array(trainingSet), trainingLabels) #使用改进的随即上升梯度训练
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[:,0]))!= int(currLine[-1]):
errorCount += 1
errorRate = (float(errorCount)/numTestVec) * 100 #错误率计算
print("测试集错误率为: %.2f%%" % errorRate)
此时错误率减少了。因此我们可以得出结论，
- 当数据集较小时，我们选用梯度上升算法；
- 当数据集较大时，我们选用随机梯度上升算法。
4. 总结
- Logistic回归无需事先假设数据分布
- 可得到“类别”的近似概率预测（概率值还可用于后续应用）
- 可直接应用现有数值优化算法（如牛顿法）求取最优解，具有快速、高效的特点