机器学习实战笔记(四)：Logist线性回归算法的Python实现

from numpy import *

def loadDataSet():                          #遍历函数，打开文本文件testSet.txt并进行逐行读取
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])]) #将X0的值设为1，因为公式w0+w1x1+w2x2
labelMat.append(int(lineArr[2]))
return dataMat,labelMat

def sigmoid(inX):                             #sigmoid函数
return 1.0/(1+exp(-inX))

def gradAscent(dataMatIn, classLabels):    #梯度上升算法计算最佳回归系数
dataMatrix = mat(dataMatIn)                   #获得输入数据并将样本数组转换为矩阵
labelMat = mat(classLabels).transpose()    #将类标签数组转换为项链并将其转置
m,n = shape(dataMatrix)                     #得到矩阵的大小
alpha = 0.001                        #步长
maxCycles = 500                      #迭代次数
weights = ones((n,1))                #回归系数初始化为1,n*1的向量
for k in range(maxCycles):                  #遍历数组
h = sigmoid(dataMatrix*weights)     #h是一个列向量，元素个数等于样本个数，矩阵相乘
error = (labelMat - h)                     #误差计算，向量减法运算
weights = weights + alpha * dataMatrix.transpose()* error    #矩阵相乘，dataMatrix.transpose()*error就是梯度f(w)
return weights

def plotBestFit(weights):                    #画出训练点
import matplotlib.pyplot as plt
dataMat,labelMat=loadDataSet()      #画点
dataArr = array(dataMat)
n = 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 = arange(-3.0, 3.0, 0.1)                          #画线，x取值区间为[-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()      #显示

def stocGradAscent0(dataMatrix, classLabels):    #随机梯度上升算法
m,n = shape(dataMatrix)
alpha = 0.01
weights = ones(n)          #初始化单位矩阵，维度为数据集列数
for i in range(m):
h = sigmoid(sum(dataMatrix[i]*weights))
error = classLabels[i] - h
weights = weights + alpha * error * dataMatrix[i]
return weights

def stocGradAscent1(dataMatrix, classLabels, numIter=150):       #升级后的梯度上升算法——随机梯度下降
m,n = shape(dataMatrix)
weights = ones(n)           #i初始化单位矩阵
for j in range(numIter):
dataIndex = range(m)
for i in range(m):
alpha = 4/(1.0+j+i)+0.0001             #alpha的值每次迭代时都会进行调整，会缓解数据波动或者高频波动
randIndex = int(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

def classifyVector(inX, weights):
prob = sigmoid(sum(inX*weights))            #判别算法
if prob > 0.5: return 1.0             #prob>0.5,返回为1
else: return 0.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):
lineArr.append(float(currLine[i]))
trainingSet.append(lineArr)          #构建训练数据集
trainingLabels.append(float(currLine[21]))     #构建分类标签训练集
trainWeights = stocGradAscent1(array(trainingSet), trainingLabels, 1000)
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(array(lineArr), trainWeights))!= int(currLine[21]):对输入向量分类，currLine[21]为分类标签
errorCount += 1  #如果不相等，误差数+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))

1.优点：计算代价不高，易于理解和实现

缺点：容易欠拟合，分类精度可能不高

2.对回归方法的认识(这里我是根据Andrew ng的Machine Learning课程的课件整理的内容)

给出一个房屋销售的例子，我们知道房屋价格收到很多因素的影响。现在探讨房屋面积与房屋价格的关系，设X轴为面积，Y轴为价格，可得到一个关系图

，接着我们做出一条曲线去尽可能的拟合这些数据


1>首先给出一些常用的概念和符号：

房屋销售记录表：训练集（训练数据），是流程中的输入数据，称为X（特征）

房屋销售价格：输出数据，称为Y（目标变量）

拟合的函数：假设（模型），一般写作Y=h(x)

训练数据的条目数：输入数据的维度n（特征的个数）

2>机器学习的过程



3>线性回归

我们用x1,x2...,xn来描述特征里的分量，如x1是房屋面积，x2是房屋朝向等等。

于是得出一个估计函数：




接下来我们对h函数进行评估，做出一个损失函数J：


我们的目标是得到最小J(theta),这个函数是对x(i)的估计值与真实值y(i)差的平方和的错误估计函数，i/2用来在求导时消去系数

如何调整theta以得到J(theta)最小值有很多方法，上面的代码我们使用了梯度上升和梯度下降的方法。另外我们还有最小二乘法等。

4>梯度下降法

在选定线性回归模型后，我们只需要确定参数theta，就可以将模型用来预测，而theta要在J(theta)最小时才能确定，故问题转化为求极值的问题。

而梯度下降法最大的问题在于求出的解可能是局部最小值也可能是全局最小值，这与初始点的选取有关。

梯度下降法的流程：





推导过程：











5.基于Logist回归和Sigmoid函数的分类

对于上面的分类器代码实现，我们想要的函数应该是，能接受所有的输入然后预测出类别，例如在两个类的情况下，上述函数输出0或1。(Heaviside step function具有这种性质，然而该函数的问题在于在跳跃点上从0瞬间到1，这个瞬间跳跃过程有时很难处理),Sigmoid函数具有相似的性质，且在数学上更易处理。具体计算公式为：



为实现Logist回归分类器，我们可以在每个特征上都乘以一个回归系数，然后把所有结果值都相加，将总和代入Sigmoid函数中去，进而得到一个0-1之间的数值，任何大于0.5的数据被分入1类，小于0.5即被分入0类。所以Logist回归也可以被看成一种概率估计。



这就结合到我们之前讲的算法了。通过看代码可能可以更好地理解这些理论知识。