学习《机器学习实战》(Logistic回归)
这一章讲了很多优化算法，记不太清了，也是借着写博客来回顾一下学到的知识
这会是激动人心的一章，因为我们将首次接触到最优化算法。 ---- 《机器学习实战》
先介绍几个基本概念
什么是回归
假设现在有一些数据点，我们用一条直线对这些点进行拟合（该线称为最佳拟合直线），这个拟合过程就称作回归。
Sigmoid函数
为了使函数，接受所有的输入都能预测出0或1的类别，而相比海维赛德阶跃函数，sigmoid在临界处有平滑的过渡。
Sogmoid公式
Logistic回归分类器
为了实现Logistic回归分类器，我们把每个特征都乘以一个回归系数，并把结果相加，将这个总和代入Sigmoid函数中，进而得到一个范围在0~1之间的数值。将大于0.5的归为1类，小于0.5的归为0类。
我们可以简单理解为:
Logistic回归的目的是寻找一个非线性Sigmoid的最佳拟合参数
Sigmoid函数的输入记为z,由下面的公式得出:
简化为向量计算 z =
w
T
w^T
wTx, w的转置与x直接相乘得到的就是上述的结果，这是矩阵乘法运算。
第一个优化算法-梯度上升法
梯度上升法基于的思想是: 要找到某函数的最大值，最好的方法就是沿着该函数的梯度方向探寻。如果梯度记为
∇
\nabla
∇,则函数f(x,y)的梯度由下式表示:
\
这个梯度意味着往x、y综合的一个变化方向，类似于f(x)=ax的导数f’(x)，表示函数变化的趋势，书上的图表示一下,这个图很生动了 :
梯度总是指向函数值变化最快的方向，这只是一个方向，还需要关心移动量，也称为步长，在深度学习中称为学习速率，记作 α \alpha α。梯度算法的迭代公式就可以这样表示:
这个公式最为重要，等下代码实现中都会用。
梯度下降算法
学深度学习时，第一个知道的就是这个梯度下降算法，和梯度上升算法基本思想一样，只需要把公式中的加号变为减号即可
接下来就可以开始实验了， 训练算法: 使用梯度上升找到最佳参数
以下是原始数据
其中每个点包含两个数值型特征: X1和X2。在此数据集上，我们将通过使用梯度上升法找到最佳回归系数,也就是拟合出Logistic回归模型的最佳参数。
代码如下:
# Logistic回归梯度上升优化算法
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.0 + np.exp(-inX))
def gradAscent(dataMatIn, classLabels):
dataMatrix = np.mat(dataMatIn) # 转换成矩阵
labelMat = np.mat(classLabels).transpose() # np.transpose矩阵转置
m, n = np.shape(dataMatrix)
alpha = 0.001 # 步长
maxCycles = 500 # 次数
weights = np.ones((n, 1)) # 参数初始化为1
for k in range(maxCycles):
h = sigmoid(dataMatrix*weights) # N行1列
error = (labelMat - h)
weights = weights + alpha * dataMatrix.transpose() * error # 梯度上升迭代公式 N行M列 * M行1列 = N行1列
return weights
这里面涉及了梯度上升法的迭代更新公式，如下
与代码中相对应,通过似然函数推导而来，这个不是重点就先略过。
画出决策边界
根据上面的算法，可以求出参数向量就可以唯一确定一条直线,这条直线是确定不同类别数据之间的分割线。
直接上代码
# 画出最佳拟合直线
def plotBestFit(wei):
import matplotlib.pyplot as plt
weights = wei
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] # 拟合曲线为0 = w0*x0+w1*x1+w2*x2, 故x2 = (-w0*x0-w1*x1)/w2,其中x0 = 1
ax.plot(x, y)
plt.xlabel('X1'); plt.ylabel('X2')
plt.show()
都挺好理解的，唯一难理解的那句代码也加上了注释，其中weights[0]、weights[1]和weights[2]都是一个数值。
看下用梯度上升法的到的边界
不得不说这效果已经很棒了，基本区分开了不同的类别
随机梯度上升
在梯度上升的基础上，想要减少运算次数，梯度上升每次更新梯度都要整个矩阵参与运算，也就是整个矩阵大小的运算量，很费时，所以诞生了随机梯度。随机梯度上升法仅用一个样本点来更新回归系数，这样有一个好处，就是当新的样本点到来时，只需要用新样本去更新一下分类器即可，不同于梯度上升算法，每次新样本到来都要整个矩阵参与更新，因而随机梯度上升算法是一个在线学习算法。
代码
# 随机梯度上升算法
def stocGradAscent(dataMatrix, classLabels):
m, n = np.shape(dataMatrix)
alpha = 0.01
weights = np.ones(n) # 初始为1
for i in range(m):
h = sigmoid(sum(dataMatrix[i]*weights)) # 都是1行N列 ndarray类型 相乘求和
error = classLabels[i] - h # 一个数值
weights = weights + alpha * error * dataMatrix[i] # alpha error都是数值
return weights
代码和梯度上升算法基本一致，不同的在于更新公式的差异，h和error都变成了数值(之前是向量),而且没有矩阵转换过程，所有数据类型都说Nunpy数组。
看下效果
效果差强人意，错分了约三分之一的样本，但原因不在于这个算法，而是迭代次数，梯度上升我们迭代了500次，而这次只是完成了行遍历。
书中开始尝试测试不同迭代次数，并记录下图所示:
总之就是经过超级多次迭代后才趋于稳定值，并且会有很明显的周期波动。书上解释产生这种周期波动的原因是存在一些不能正确分类的样本点(数据集并非线性可分),所以每次迭代都会引发系数的剧烈改变。
想要让算法避免波动，收敛到某个稳定值，并提高收敛速度，也就是减少迭代次数。
改进的随机梯度上升算法
# 改进的随机梯度上升算法
def stocGradAscent(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
这个算法改进了三个地方，一个是增加了迭代次数numIter来控制外层循环，二是学习率alpha初始很大，并随着迭代次数的增加逐渐减小，但不会到0，提高了收敛速度，三是内层循环不再按样本顺序更新参数，而是随机顺序。
书中讲alpha每次减少1/(j+i)，这样当j<<max(i)时，alpha不是严格下降的。不是很理解什么是严格下降???, 以及通过随机选取样本来更新回归系数，能减少周期性波动。
得到下方效果图:
注意横坐标的迭代次数，减少的不是一点点,说明这样改进是非常优秀的。
试一下
也是非常完美的效果。
接下来看案例
病马的生死预测问题
1、准备数据: 处理数据中的缺失值
有很多可选的做法，需要根据实际情况采取不同的措施
- 使用可用特征的均值来填补缺失值
- 使用特殊值来填补缺失值，如-1
- 忽略有缺失值的样本
- 使用相似样本的均值添补缺失值
- 使用另外的机器学习算法预测缺失值
本例中采用缺失值取0的方式
2、测试算法: 用Logistic回归进行分类
所需要做的就是把测试集上每个特征向量乘以最优化得来的回归系数，再将乘积求和，最后输入Sigmoid函数即可。如果对应的Sigmoid值大于0.5就预测类别标签为1，否则为0。
# 分类
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(21): # 前21列为特征
lineArr.append(float(currLine[i]))
trainingSet.append(lineArr)
trainingLabels.append(float(currLine[21])) # 最后一列为标签
trainWeights = stocGradAscent(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(21):
lineArr.append(float(currLine[i]))
if int(classifyVector(np.array(lineArr), trainWeights)) != int(currLine[21]):
errorCount += 1
errorRate = (float(errorCount)/numTestVec)
print('the error rate of this test is: %f'% errorRate)
return errorRate
def multiTest():
numTests = 10; errorSum = 0.0
for k in range(numTests):
errorSum += colicTest()
print("after %d iterations the average error rate is: %f" % (numTests, errorSum/float(numTests)))
代码挺简单的，都看的很明白, 看下结果
错误率38%，毕竟有1/3的数据缺失，如果增加colicTest中的迭代次数还可以缩小错误率，应该能降到20%左右。
总结，随机梯度上升算法和梯度上升算法的效果相当，但占用更少的计算资源。此外，随机梯度上升是一个在线算法，它可以在新数据到来时就完成参数更新，而不需要重新读取整个数据集来进行批处理运算。
下一章介绍与Logistic回归类似的另一种分类算法: 支持向量机，它被认为是目前最好的现成算法之一。
扫一扫
209
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
