关闭

logistic回归和Python实现

869人阅读 评论(0) 收藏 举报
分类:

一、Logistic Regression的基本内容

通过学习了台湾的林教授和Stanford的课程后发现,他们两个人的基本思路虽然一致,但是具体做法有所差异,下面简单介绍一下两种实现方式。
1、台湾的林教授方法 
     使用Logistic回归的思路是,我们不希望只给输出一个确定的{-1,+1}的结果,而是想给出一个概率值,这时就可以用到Logistic回归,得到概率的输出值在0~1之间。 
   但是实际情况是,我们不能知道准确的概率值,只能知道输出的结果。比如以心脏病复发为例。绿色和红色分别表示理想情况下和实际的结果。 
实际情况理想的情况
 可以将实际的情况看成是含有噪声的理想训练数据。和以前一样,同样用输入加权的形式得到一个分数,加权函数并将此分数通过称为Logistic function或Sigmoid function的函数转换为一个0~1之间的值,在这里用θ(s)表示。s越大患病风险越高。假设函数h为函数h。 
     和平时一样,方法是最小化Logistic回归的错误Ein。而Logistic回归的目标函数为
     假设存在一个数据集D:,则产生这个数据集的概率为各个输入样本对应概率的连乘。通过把目标函数带入可以得到带有f(x)的形式,而得到f(x)也就是根本目的。
    但实际上f(x)是未知的,已知的只有假设函数h。是否可以用假设函数代替f呢?意味着假设函数h产生同样的数据集样本D的可能性有多大。这在数学上称为似然(likelihood)。(注意:这里和Stanford后面讲的就不是太一样了,个人感觉这种思路比较好理解。替代公式变成了
      如果假设函数h和未知函数f很接近(err很小),那么h产生数据样本D的可能性与f产生同样数据D的可能性很接近。由于f(x)产生了样本D,那么可以认为f产生该数据样本的几率很大。因此推断出最好的假设函数g,就是似然性最大的假设函数h。即
         通过h的性质,可以将似然函数写为:,P(x)对似然的大小没有影响,只与函数h对每个样本的连乘有关。并将连乘通过取ln变成加和形式,为了得到Ein通过加负号,变成求最小值。
         这里,定义了交叉熵误差(cross entropy error)来衡量训练误差

       我们得到了训练误差的具体形式,由于这个训练误差函数是可微并且是凸函数,所以依照之前的思路,对这个函数求梯度。对每一个权值向量wi求偏微分,最后合并。
        这时分为两种情况,一种是线性可分时权值θ(.)为0,即y与wTX同号。排除这种特殊情况,当加权求和为0时,不能用类似解线性回归时的闭合求解方式,这时的计算就要用迭代求解的方式了。类似梯度上升、梯度下降和牛顿方法都是用于这种加权和的形式。
        按照类似的思路,一步一步的去修正w,使得Ein的结果越来越小。在logistic回归中,Ein的式子是平滑的。现在我们可以将这个Ein想象成一个山谷,我们要到达山谷的最低点,就是要沿着当前的梯度最大的方向每次迈出一小步,直到到达谷底,使得Ein最小,得到最佳解w。迭代函数为(怎么选择的?)。增加一个η表示步伐,v表示方向,则迭代函数可以变为
         最终目标是最小化Ein(Wt+η'v)。要确定一个w使得Ein最小,不可能一步到位,指导思路还是由繁化简,用线性近似(linear approximation)的方式来解决问题,我们用多维度的泰勒展开公式来近似Ein,只要给定一个小的η,就可以近似这个Ein。
这里Ein(wt)和η都是已知的,Ein的梯度表示了下降的方向,也可以求出来,唯一要考虑的就是最好的向量v该如何选择。于是,对于v和Ein的梯度这两个向量,使得其值是最小的方法就是让v向量的方向和Ein的梯度的方向相反,这样使得两个向量的内积是最小的,另外由于v的模是1,还需要有个归一化的步骤。

为了找到一个适合的η,防止η太小时下降慢,太大时不稳定。让η随着梯度的变化变化(与梯度大小成正比)。所以公式变为这里的η为一个新的值,被称为固定的学习速率。

总结:logistic回归求解最优化解的步骤为:
2、Stanford机器学习方法
       和上面的方法主要的区别从likelihood开始。上面方法主要是想让h能够最好的拟合f(x),使产生相同数据集的几率最大。Stanford的选择为对应参数θ(上面的方法参数名为w)的似然函数。但实际上是一样的,因为上面的方法在求解最大值得过程中略去了固定值P(x),得到的似然函数其实和下面的一样。但是上面的方法用的方式不是像下面这种方法将似然函数分开写。
下面是具体方法:

接着为了求解方便,取了对数。而目的就是最大化似然函数。
最大化似然函数有两种方法,第一种是梯度上升法,一种是牛顿方法。
    (1)梯度上升法
              ,不断更新参数。而梯度求解方法如下:

这里面用了一个公式,对于sigmoid function g(x),
所以,最终的梯度上升规则为
然后repeat所有的j。但是其实不需要for所有的j,可以用矩阵来解决。(下面的做个参考就好,实际代码中一般用矩阵虽然全面,但是相应的速度也会很慢,代码中会用随机梯度上升法来解决这一问题)。
      (2)其他方法
              优化算法除了梯度下降算法外,还包括:
  • Conjugate gradient method(共轭梯度法)
  • Quasi-Newton method(拟牛顿法)
  • BFGS method
  • L-BFGS(Limited-memory BFGS)
             这些方法的有点是不用计算步伐就可以。

        下面简单介绍一下牛顿方法。牛顿法的目的是找到f(x)= 0的点,方法用的是泰勒展开不断的找到△x不断的迭代更新。公式为:。下图为不断迭代的过程。

对应于我们的问题,为了求解最大化的似然函数,我们可以令,得到迭代函数:
但在我们的logistic回归中,θ是向量,要把公式对应到高维的情况下,这时我们得到了:
其中H叫做Hessian。缺点是当特征数量很多(几千)计算的代价是很大的。


3、总结
       两种方法有各自的优缺点,第一种方法计算相对麻烦一点,但是得到的迭代函数的学习率却是随着梯度自适应变化的。第二种的计算方法比较直观,但是学习率不好处理,只能靠经验或测试观察。
      其实还有一种方法是直接用代价函数(cost function)做的。
还有更多的内容,如防止过拟合和欠拟合,多分类等。

4、python实现
      代码主要参考了《机器学习实战》这本书,还有http://blog.csdn.net/zouxy09/article/details/20319673所写的程序,他们的程序多少都有一点错误,这些我做了适当的修改,而且即使没有修改的地方也做了备注。程序主要用的是梯度上升法。
# -*- coding: utf-8 -*-

import random
from numpy import *
from math import *
import matplotlib.pyplot as plt

# 逻辑函数
def sigmoid(inX):
	return 1.0/(1+exp(-inX))
#加载数据
def loadDataSet(filename):
	numFeat = len(open(filename).readline().split('\t'))-1
	dataMat = []; labelMat = []
	fr = open(filename)
	for line in fr.readlines():
		lineArr = []
		curLine = line.strip().split('\t')
		for i in range(numFeat):
			lineArr.append(float(curLine[i]))
		dataMat.append([1]+lineArr)  #这里是为了使 x0 等于 1
		labelMat.append(float(curLine[-1]))
	return mat(dataMat),mat(labelMat).T

#Logistic 回归
def logisticRegression(dataMat,labelMat,opts):
	numSamples,numFeature = shape(dataMat)
	alpha = opts['alpha'] ; maxIter = opts['maxIter']
	weight = ones((numFeature,1))
	for j in range(maxIter):
		#选择方法共三种:梯度下降,随机梯度下降,优化的随机梯度下降
		if opts['optimizeType'] == 'gradDescent':# gradient descent algorilthm
			h = sigmoid(dataMat * weight)  #注意这里会报错,因为exp()的自变量只能是数,不能是矩阵或数组
			error = labelMat - h
			weight = weight + alpha * dataMat.T * error
		elif opts['optimizeType'] == 'stocGradDescent': # stochastic gradient descent
			for i in range(numSamples):
				h = sigmoid(sum(dataMat[i]*weight))
				error = labelMat[i] - h
				weight = weight + alpha * dataMat[i].T* error
		elif opts['optimizeType'] == 'smoothStocGradDescent':# smooth stochastic gradient descent
			dataIndex = range(numSamples)
			for i in range(numSamples):
				randIndex = int(random.uniform(0,len(dataIndex)))#产生0到dataIndex之间的随机数
				alpha = 4.0/(1.0+j+i)+0.01 #j是迭代次数,i是迭代时,第i个选出的样本
				h = sigmoid(sum(dataMat[randIndex]*weight))
				error = labelMat[randIndex] - h
				weight = weight + alpha * dataMat[randIndex].T* error
				del(dataIndex[randIndex])
		else:
			raise NameError('Not support optimize method type!')
	return weight

#准确度计算
def testLogRegres(weights, test_x, test_y):
	numSamples, numFeatures = shape(test_x)
	matchCount = 0
	for i in xrange(numSamples):
		predict = sigmoid(test_x[i, :] * weights)[0, 0] > 0.5
		if predict == bool(test_y[i, 0]):
			matchCount += 1
	accuracy = float(matchCount) / numSamples
	return accuracy


#画出决策边界
def showLogRegres(weight,dataMat,labelMat):
	#输入的x和y为矩阵
	numSamples, numFeatures = shape(dataMat)
	weights = weight.getA() #将矩阵转为数组,经测试应该是对应于[[0],[1],[2]]这种的才有效
	if numFeatures != 3:
		print "Sorry! I can not draw because the dimension of your data is not 2!"
		return 1
		#xrange 与range相比返回的是一个xrange对象,而不是一个数列。xrange的对象是单个元素,可以节省空间。本质为一个生成器
		#在python3中,将xrange命名为range,而取消了range
	for i in xrange(numSamples):
		if int(labelMat[i]) == 0:
			plt.plot(dataMat[i,1],dataMat[i,2],'or')
		elif int(labelMat[i]) == 1:
			plt.plot(dataMat[i,1],dataMat[i,2],'ob')

	#画分类线
	min_x = min(dataMat[:,1])[0,0] # 只有用[0,0]才能得到一个数值
	max_x = max(dataMat[:, 1])[0,0]
	y_min_x = float(-weights[0] - weights[1] * min_x) / weights[2] # 这是令sigmoid函数的输入为零,得到分界处的值
	y_max_x = float(-weights[0] - weights[1] * max_x) / weights[2]
	plt.plot([min_x, max_x], [y_min_x, y_max_x], '-g')
	plt.xlabel = ('X1'); plt.ylabel = ('X2')
	plt.show()

其他的方法会在后面介绍。



参考资料:http://blog.csdn.net/abcjennifer/article/details/7716281(讲的代价函数的方法)
               http://blog.chinaunix.net/uid-20761674-id-4336428.html
                http://www.open-open.com/lib/view/open1421575796468.html
0
0

查看评论
* 以上用户言论只代表其个人观点,不代表CSDN网站的观点或立场
    个人资料
    • 访问:37672次
    • 积分:746
    • 等级:
    • 排名:千里之外
    • 原创:35篇
    • 转载:9篇
    • 译文:0篇
    • 评论:27条
    最新评论