预测数值型数据——回归
通过前面博文的介绍,对于分类有了一定的了解,由于在之前使用过sigmod函数进行梯度下降,所以对于这类算法也做了一个学习,今天正式学习了《机器学习实战》中的回归一章。分类针对离散数据,而回归针对连续型数据。
回归的目的是预测数值型的目标值。h=w0x0+w1x1.....+wn*xn
这个方程就是regression equation,其中的w0,w1,,,wn被称为regression weights,求这些回归系数的过程就是回归。得到回归系数后,再给定输入,就可以用这个方程进行预测。回归一般指线性回归,文中简化线性回归的说法为回归。如何找到这些regression weights呢?前面的转载博文中提到了最小二乘法等等优化办法。利用普通最小二乘法(OLS,ordinary least squares)可以得到w=(XtX)-1Xty 注意这里面t表示转置,-1表示逆,X为矩阵(一般而言),y,w表示向量,具体推导见之前博文http://www.cnblogs.com/litian0605/p/5239774.html
注意:这里需要使用逆阵,而一个矩阵是否可逆,是有要求的,在上面这篇博文中有提到(可逆阵等价非奇异阵等价行列式不为0)。代码如下:
首先对读文件做一下讲解:
fr=open('ex0.txt')
a=fr.readline() b=a.split('\t') c=a.strip().split('\t')
open函数打开一个指定文件名的文件,a为读取一行,打印a得到:
b为对a直接进行以\t分隔操作,打印b:可以看到得到了一个list,此时直接求长度也是正确的,最后一个字符串中有'\n'
c为对a做strip操作后再以'\t'分隔,对a做strip()操作后得到:,再以‘\t’分隔,也就是c:
现在贴出用最小二乘法得到的公式计算回归系数的代码:
# -*- coding: utf-8 -*- from numpy import * def loadDataSet(fileName): numfeat=len(open(fileName).readline().split('\t'))-1#获取列数减1(最后一列为真实y值,与x值分开存储) 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(lineArr) labelMat.append(float(curLine[-1])) return dataMat,labelMat def standRegres(xArr,yArr): xMat=mat(xArr);yMat=mat(yArr).T xTx=xMat.T*xMat if linalg.det(xTx)==0.0: print "this matrix is singular,cannot do inverse"#奇异阵,无逆阵 return ws=xTx.I*(xMat.T*yMat) return ws
注意:数据集中的x0永远为1,w0为常数项
如此,利用最小二乘法得到了回归系数ws:
(w0=3.007743,w1=1.69532) h=3.0077443+1.69532*x1为我们得到的回归方程
现在我们将数据集的点进行绘制,并绘制预测直线,同时计算真实值与预测值的相关性:可以通过调用numpy中的corrcoef(yEstimat,yActual)来计算预测值与实际值的相关性,保证两个都是行向量
1 def plot(xArr,yArr,ws): 2 xMat=mat(xArr);yMat=mat(yArr) 3 yHat=xMat*ws#yHat为只有一列的矩阵,yMat为只有一行的矩阵 4 fig=plt.figure() 5 ax=fig.add_subplot(111) 6 #绘制数据集中的散点 7 ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0]) 8 #计算相关性 9 corr=corrcoef(yHat.T,yMat) 10 print corr; 11 #绘制预测直线 12 xCopy=xMat.copy() 13 xCopy.sort(0)#绘制直线时如果点的次序混乱,绘图将出现问题,首先要将点按照升序排列 14 #yHat=xCopy*ws#计算预测值 15 yHat=xCopy*ws 16 ax.plot(xCopy[:,1],yHat)#将xCopy矩阵中的下标为1的列(第二列)作为横坐标,yHat为纵坐标 17 plt.show() 18
可以看到对角线数据为1.0,因为yMat与自身匹配是完美的,而yhat与ymat的相关系数为0.98。
线性回归的一个问题在于可能出现欠拟合现象,所求为最小均方误差的无偏估计。所以在估计中引入一些偏差,降低预测的均方误差。给出一个方法:
局部加权线性回归(Locally weighted linear regression LWLR)为待测点附近的每个点赋予一定的权重,然后与之前类似基于最小均方差(利用普通最小二乘法)来进行普通的回归。有点类似与knn,当给定一个样本中的x0,x1,x2~xn后,要求给出预测值时,我们首先利用给定的核函数,计算训练样本中其他样例与待测样例的相近度,得到一个矩阵W,W为对角阵,对角线上元素即为刚才计算出的相近度,然后基于最小均方差进行普通的回归。该算法解出的回归系数w的形式如下:
w=(XtWX)-1XtWy
其中W为对角阵(权重阵),对角线上的元素表示其他点对于当前待测点的权重。LWLR使用“核”来对附近的点赋予更高的权重,核的类型可以自由选择,一般选择高斯核:
W(i,i)=exp(|x(i)-x|/-2k^2)该式表示W中第i行i列的值,可以理解为i第i个点与当前待测点的相近度。
贴上代码:
def lwlr(testPoint,xArr,yArr,k=1.0): xMat=mat(xArr);yMat=mat(yArr).T m=shape(xMat)[0]#得出训练样本数 weights=mat(eye((m)));#得到阶数为m的单位方阵(对角线全部为1) for j in range(m): diffMat=testPoint-xMat[j,:] weights[j,j]=exp(diffMat*diffMat.T/(-2.0*k**2))#这里使用的分子是平方后的值 xTx=xMat.T*weights*xMat if linalg.det(xTx)==0.0: print "this matrix is singular" return ws=xTx.I*xMat.T*weights*yMat return testPoint*ws#输出预测值 def lwlrTest(testArr,xArr,yArr,k=1.0): m=shape(testArr)[0]#获取待预测阵中样例个数 yHat=zeros(m) for i in range(m): yHat=lwlr(testArr[i],xArr,yArr,k) return yHat
xArr,yArr=regression.loadDataSet('ex0.txt')#其中xArr为数据集中不包含真实y值的数据,yArr为数据集中的真实y值
如果想要得到对某一条样例的预测值:regression.lwlr(xArr[0],xArr,yArr,1.0)
如果想要得到对于所有点的估计,那么如下调用:regression.lwlr(xArr,xArr,yArr,1.0)#这里xArr既为训练集,又为测试集。
同样可以绘制图像。如下为k=0.03时:
可以通过调节k的取值得到不同的回归曲线。
但是如果出现一个数据集中特征的数目大于样本点的数目(n>m),那么XtX就会是非满秩阵(因为XtX是n*n的方阵,又因为矩阵乘积的秩等于矩阵秩的成绩,所以乘积矩阵的行列式小于乘数矩阵的秩,所以乘积阵的秩<=m<n,所以非满秩)提出岭回归应对这种因为非满秩而无法求逆的问题。w=(XtX+lambda*I)-1Xty。其中I为单位阵,lambda是一个不断调节的系数。这样可以保证逆阵存在。为了使用岭回归,需要对数据集进行标准化。