机器学习之Logistic回归与Python实现
logistic回归是一种广义的线性回归,通过构造回归函数,利用机器学习来实现分类或者预测。
一 Logistic回归概述
Logistic回归的主要思想是,根据现有的数据对分类边界建立回归公式,从而实现分类(一般两类)。“回归”的意思就是要找到最佳拟合参数,其中涉及的数学原理和步骤如下:
(1)需要一个合适的分类函数来实现分类【单位阶跃函数、Sigmoid函数】
(2)损失函数(Cost函数)来表示预测值(h(x))与实际值(y)的偏差(h−y),要使得回归最佳拟合,那么偏差要尽可能小(偏差求和或取均值)。
(3)记J(ω)表示回归系数为ω时的偏差,那么求最佳回归参数ω就转换成了求J(ω)的最小值。【梯度下降法】
所以,接下来就围绕这几个步骤进行展开。
1.1 分类函数
假设要实现二分类,那么可以找一个函数,根据不同的特征变量,输出0和1,并且只输出0和1,这种函数在某个点直接从0跳跃到1,如:
但是这种函数处理起来,稍微有点麻烦,我们选择另外一个连续可导的函数,也就是Sigmoid函数,函数的公式如下:
这个函数的特点是,当x=0.5时,h(x)=0.5,而x越大,h(x)越接近1,x越小,h(x)越接近0。函数图如下:
这个函数很像阶跃函数,当x>0.5,就可以将数据分入1类;当x<0.5,就可以将数据分入0类。
确定了分类函数,接下来,我们将Sigmoid函数的输入记为z,那么
向量x是特征变量,是输入数据,向量w是回归系数是特征
之后的事情就是如何确定最佳回归系数ω(w0,w1,w2,...,wn)
1.2 Cost函数
现有
对于任意确定的x和w,有:
这个函数可以写成:
取似然函数:
求对数似然函数:
因此,就构造得到了函数J(w)来表示预测值与实际值的偏差,于是Cost函数也可以写成:
所以,我们可以用J(w)来表示预测值与实际值的偏差,也就是Cost函数,接下里的任务,就是如何让偏差最小,也就是J(w)最大
Question:为什么J(w)可以表示预测值与实际值的大小,为什么J(w)最大表示偏差最小。
我们回到J(w)的推导来源,来自
P(y=1|x,w)=hw(x)和P(y=0|x,w)=1−hw(x),
那么显然有
当x>0,此时y=1,1/2<hw(x)<1,所以P(y=1|x,w)=hw(x)>1/2
当x<0,此时y=0,0<hw(x)<1/2,所以P(y=0|x,w)=1−hw(x)>1/2
所以,无论如何取值,P(y=0|x,w)都大于等于1/2,P(y=0|x,w)越大,越接近1,表示落入某分类概率越大,那么分类越准确,预测值与实际值差异就越小。
所以P(y=0|x,w)可以表示预测值与实际值的差异,且P(y=0|x,w)越大表示差异越小,所以其似然函数J(w)越大,预测越准确。
所以,接下来的任务,是如何求解J(w)最大时的w值,其方法是梯度上升法。
1.3 梯度上升法求J(w)最大值
梯度上升法的核心思想是:要找某个函数的最大值,就沿着这个函数梯度方向探寻,如果梯度记为∇,那么函数f(x,y)的梯度是:
梯度上升法中,梯度算子沿着函数增长最快的方向移动(移动方向),如果移动大小为α(步长),那么梯度上升法的迭代公式是:
问题转化成:
首先,我们对J(w)求偏导:
在第四至第五行的转换,用到的公式是:
将求得的偏导公式代入梯度上升法迭代公示:
可以看到,式子中所有函数和输入的值,都是已知的了。接下来,可以通过Python实现Logistic回归了。
二、Python算法实现
2.1 梯度上升法求最佳回归系数
首先,数据取自《机器学习实战》中的数据,部分数据如下:
-0.017612 14.053064 0
-1.395634 4.662541 1
-0.752157 6.538620 0
-1.322371 7.152853 0
0.423363 11.054677 0
0.406704 7.067335 1
先定义函数来获取数去,然后定义分类函数Sigmoid函数,最后利用梯度上升法求解回归系数w。
建立一个logRegres.py文件,输入如下代码:
from numpy import *
#构造函数来获取数据
def loadDataSet():
dataMat=[];labelMat=[]
fr=open('machinelearninginaction/Ch05/testSet.txt')
for line in fr.readlines():
lineArr=line.strip().split()
dataMat.append([1.0,float(lineArr[0]),float(lineArr[1])])#特征数据集,添加1是构造常数项x0
labelMat.append(int(lineArr[-1]))#分类数据集
return dataMat,labelMat
def sigmoid(inX):
return 1/(1+exp(-inX))
def gradAscent(dataMatIn,classLabels):
dataMatrix=mat(dataMatIn) #(m,n)
labelMat=mat(classLabels).transpose() #转置后(m,1)
m,n=shape(dataMatrix)
weights=ones((n,1)) #初始化回归系数,(n,1)
alpha=0.001 #定义步长
maxCycles=500 #定义最大循环次数
for i in range(maxCycles):
h=sigmoid(dataMatrix * weights) #sigmoid 函数
error=labelMat - h #即y-h,(m,1)
weights=weights + alpha * dataMatrix.transpose() * error #梯度上升法
return weights
在python命令符中输入代码对函数进行测试:
In [8]: import logRegres
...:
In [9]: dataArr,labelMat=logRegres.loadDataSet()
...:
In [10]: logRegres.gradAscent(dataArr,labelMat)
...:
Out[10]:
matrix([[ 4.12414349],
[ 0.48007329],
[-0.6168482 ]])
于是得到了回归系数。接下来根据回归系数画出决策边界wTx=0
定义作图函数:
def plotBestFit(weights):
import matplotlib.pyplot as plt
dataMat,labelMat=loadDataSet()
n=shape(dataMat)[0]
xcord1=[];ycord1=[]
xcord2=[];ycord2=[]
for i in range(n):
if labelMat[i]==1:
xcord1.append(dataMat[i][1])
ycord1.append(dataMat[i][2])
else:
xcord2.append(dataMat[i][1])
ycord2.append(dataMat[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,3,0.1)
y=(-weights[0,0]-weights[1,0]*x)/weights[2,0] #matix
ax.plot(x,y)
plt.xlabel('X1')
plt.ylabel('X2')
plt.show()
在Python的shell中对函数进行测试:
In [11]: weights=logRegres.gradAscent(dataArr,labelMat)
In [12]: logRegres.plotBestFit(weights)
...:
2.2 算法改进
(1) 随机梯度上升
上述算法,要进行maxCycles次循环,每次循环中矩阵会有m*n次乘法计算,所以时间复杂度(开销)是maxCycles*m*n,当数据量较大时,时间复杂度就会很大。因此,可以是用随机梯度上升法来进行算法改进。
随机梯度上升法的思想是,每次只使用一个数据样本点来更新回归系数。这样就大大减小计算开销。
代码如下:
def stocGradAscent(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] #array 和list矩阵乘法不一样
return weights
注意:gradAscent函数和这个stocGradAscent函数中的h和weights的计算形式不一样,因为
前者是的矩阵的计算,类型是numpy的matrix,按照矩阵的运算规则进行计算。
后者是数值计算,其类型是list,按照数值运算规则计算。
对随机梯度上升算法进行测试:
In [37]: dataMat,labelMat=logRegres.loadDataSet()
...:
In [38]: weights=logRegres.stocGradAscent(array(dataMat),labelMat)
...:
In [39]: logRegres.plotBestFit(mat(weights).transpose())
...:
输出的样本数据点和决策边界是:
(2)改进的随机梯度上升法
def stocGradAscent1(dataMatrix,classLabels,numIter=150):
m,n=shape(dataMatrix)
weights=ones(n)
for j in range(numIter):
dataIndex=list(range(m))
for i in range(m):
alpha=4/(1+i+j)+0.01#保证多次迭代后新数据仍然具有一定影响力
randIndex=int(random.uniform(0,len(dataIndex)))#减少周期波动
h=sigmoid(sum(dataMatrix[randIndex] * weights))
error=classLabels[randIndex]-h
weights=weights + alpha*dataMatrix[randIndex]*error
del(dataIndex[randIndex])
return weights
在Python命令符中测试函数并画出分类边界:
In [188]: weights=logRegres.stocGradAscent1(array(dataMat),labelMat)
...:
In [189]: logRegres.plotBestFit(mat(weights).transpose())
...:
(3)三种方式回归系数波动情况
普通的梯度上升法:
随机梯度上升:
改进的随机梯度上升
评价算法优劣势看它是或否收敛,是否达到稳定值,收敛越快,算法越优。
三 实例
3.1 通过logistic回归和氙气病症预测马的死亡率
数据取自《机器学习实战》一书中的氙气病症与马死亡的数据,部分数据如下:
2.000000 1.000000 38.500000 66.000000 28.000000 3.000000 3.000000 0.000000 2.000000 5.000000 4.000000 4.000000 0.000000 0.000000 0.000000 3.000000 5.000000 45.000000 8.400000 0.000000 0.000000 0.000000
1.000000 1.000000 39.200000 88.000000 20.000000 0.000000 0.000000 4.000000 1.000000 3.000000 4.000000 2.000000 0.000000 0.000000 0.000000 4.000000 2.000000 50.000000 85.000000 2.000000 2.000000 0.000000
2.000000 1.000000 38.300000 40.000000 24.000000 1.000000 1.000000 3.000000 1.000000 3.000000 3.000000 1.000000 0.000000 0.000000 0.000000 1.000000 1.000000 33.000000 6.700000 0.000000 0.000000 1.000000
#定义分类函数,prob>0.5,则分入1,否则分类0
def classifyVector(inX,trainWeights):
prob=sigmoid(sum(inX*trainWeights))
if prob>0.5:return 1
else : return 0
def colicTest():
frTrain = open('machinelearninginaction/Ch05/horseColicTraining.txt')#训练数据
frTest = open('machinelearninginaction/Ch05/horseColicTest.txt')#测试数据
trainSet=[];trainLabels=[]
for line in frTrain.readlines():
currLine=line.strip().split('\t')
lineArr=[]
for i in range(21):
lineArr.append(float(currLine[i]))
trainSet.append(lineArr)
trainLabels.append(float(currLine[21]))
trainWeights=stocGradAscent1(array(trainSet),trainLabels,500)#改进的随机梯度上升法
errorCount=0;numTestVec=0
for line in frTest.readlines():
numTestVec+=1
currLine=line.strip().split('\t')
lineArr=[]
for i in range(21):
lineArr.append(float(currLine[i]))
if classifyVector(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
for k in range(numTests):
errorSum+=colicTest()
print('after %d iterations the average error rate is:%f'%(numTests,errorSum/float(numTests)))
在控制台命令符中输入命令来对函数进行测试:
In [3]: logRegres.multiTest()
G:\Workspaces\MachineLearning\logRegres.py:19: RuntimeWarning: overflow encountered in exp
return 1/(1+exp(-inX))
the error rate of this test is :0.313433
the error rate of this test is :0.268657
the error rate of this test is :0.358209
the error rate of this test is :0.447761
the error rate of this test is :0.298507
the error rate of this test is :0.373134
the error rate of this test is :0.358209
the error rate of this test is :0.417910
the error rate of this test is :0.432836
the error rate of this test is :0.417910
after 10 iterations the average error rate is:0.368657
分类的错误率是36.9%。
数据分析咨询请扫描二维码
若不方便扫码,搜微信号:CDAshujufenxi
在数据分析领域,Excel作为一种普及率极高且功能强大的工具,无疑为无数专业人士提供了便捷的解决方案。尽管Excel自带了丰富的功 ...
2025-01-17在这个瞬息万变的时代,许多人都在寻找能让他们脱颖而出的职业。而数据分析师,作为大数据和人工智能时代的热门职业,自然吸引了 ...
2025-01-14Python作为一门功能强大的编程语言,已经成为数据分析和可视化领域的重要工具。无论你是数据分析的新手,还是经验丰富的专业人士 ...
2025-01-10完全靠数据决策,真的靠谱吗? 最近几年,“数据驱动”成了商界最火的关键词之一,但靠数据就能走天下?其实不然!那些真正成功 ...
2025-01-09SparkSQL 结构化数据处理流程及原理是什么?Spark SQL 可以使用现有的Hive元存储、SerDes 和 UDF。它可以使用 JDBC/ODB ...
2025-01-09在如今这个信息爆炸的时代,数据已然成为企业的生命线。无论是科技公司还是传统行业,数据分析正在深刻地影响着商业决策以及未来 ...
2025-01-08“数据为王”相信大家都听说过。当前,数据信息不再仅仅是传递的媒介,它成为了驱动经济发展的新燃料。对于企业而言,数据指标体 ...
2025-01-07在职场中,当你遇到问题的时候,如果感到无从下手,或者抓不到重点,可能是因为你掌握的思维模型不够多。 一个好用的思维模型, ...
2025-01-06在现代企业中,数据分析师扮演着至关重要的角色。每天都有大量数据涌入,从社交媒体到交易平台,数据以空前的速度和规模生成。面 ...
2025-01-06在职场中,许多言辞并非表面意思那么简单,有时需要听懂背后的“潜台词”。尤其在数据分析的领域里,掌握常用术语就像掌握一门新 ...
2025-01-04在当今信息化社会,数据分析已成为各行各业的核心驱动力。它不仅仅是对数字进行整理与计算,而是在数据的海洋中探寻规律,从而指 ...
2025-01-03又到一年年终时,各位打工人也迎来了展示成果的关键时刻 —— 年终述职。一份出色的年终述职报告,不仅能全面呈现你的工作价值, ...
2025-01-03在竞争激烈的商业世界中,竞品分析对于企业的发展至关重要。今天,我们就来详细聊聊数据分析师写竞品分析的那些事儿。 一、明确 ...
2025-01-03在数据分析的江湖里,有两个阵营总是争论不休。一派信奉“大即是美”,认为数据越多越好;另一派坚守“小而精”,力挺质量胜于规 ...
2025-01-02数据分析是一个复杂且多维度的过程,从数据收集到分析结果应用,每一步都是对信息的提炼与升华。可视化分析结果,以图表的形式展 ...
2025-01-02在当今的数字化时代,数据分析师扮演着一个至关重要的角色。他们如同现代企业的“解密专家”,通过解析数据为企业提供决策支持。 ...
2025-01-02数据分析报告至关重要 一份高质量的数据分析报告不仅能够揭示数据背后的真相,还能为企业决策者提供有价值的洞察和建议。 年薪 ...
2024-12-31数据分析,听起来好像是技术大咖的专属技能,但其实是一项人人都能学会的职场硬核能力!今天,我们来聊聊数据分析的核心流程,拆 ...
2024-12-31提到数据分析,你脑海里可能会浮现出一群“数字控”抱着电脑,在海量数据里疯狂敲代码的画面。但事实是,数据分析并没有你想象的 ...
2024-12-31关于数据分析师是否会成为失业高危职业,近年来的讨论层出不穷。在这个快速变化的时代,技术进步让人既兴奋又不安。今天,我们从 ...
2024-12-30