tonglin0325的个人主页

机器学习——预测数值型数据:回归

线性回归

优点:结果易于理解,计算上不复杂

缺点:对非线性的数据拟合不好

适用数据类型:数值型和标称型数据

 

回归的目的就预测数值型的目标值。最直接的办法就是依据输入写一个目标值的计算公式。这个计算公式就是所谓的回归方程(regression equation),其中的参数就是回归系数,求这些回归系数的过程就是回归

  说道回归,一般都是指线性回归(linear regression)

给定由d个属性描述的示例 ,其中xi是x在第i个属性上的取值,线性模型试图学得一个通过属性组合来进行预测的函数,即

**  **

 

一般用向量形式写成

**  **

其中w={w1;w2;…;wd},w和b学得之后,模型就得以确定。

 

更加一般的情况是,给定一个数据集** **,其中 ,“线性回归”试图学得一个线性模型以尽可能准确地预测实值输出标记。

向量形式写成

**  ,这称为“多元线性回归”**

 

因为是d维的特性,那么w就是一个由回归系数组成的d×1维向量,Xn×(d+1)的矩阵(第一列元素都是1,其余列都是x1…xn)(d+1)×1的向量

  

此时,平方误差(由于误差有正有负,所以不能直接相加)可以表示成

  ,其中 

需要使得平方误差最小,就需要对其求导,并另求导后的式子等于0,求出

  

而这个使得平方误差最小的算法就称为普通最小二乘法(OLS)

 

 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
from numpy import *

def loadDataSet(fileName): #general function to parse tab -delimited floats
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(lineArr)
labelMat.append(float(curLine[-1]))
return dataMat,labelMat

def standRegres(xArr,yArr): #根据公式解出xArr的2维回归系数
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 = linalg.inv(xTx) * (xMat.T*yMat)
ws = xTx.I * (xMat.T*yMat) #.I求自身逆矩阵,.T返回自身的转置,.H返回自身的共轭转置,.A返回自身数据的2维数组的一个视图
return ws

 

 mian.py

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
# coding:utf-8
# !/usr/bin/env python

import regression
import matplotlib.pyplot as plt
from numpy import *

if __name__ == '__main__':
xArr,yArr = regression.loadDataSet("ex0.txt")
ws = regression.standRegres(xArr,yArr)
print ws #输出回归系数
xMat = mat(xArr) #输入的X矩阵,n×(d+1)维
yMat = mat(yArr) #真实的标记,n×1维
yHat = xMat*ws #预测的结果
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0]) #取得矩阵的第二列并转换成list
xCopy = xMat.copy()
xCopy.sort(0) #参数0表示从小到大,1表示从大到小,排序
yHat = xCopy*ws
ax.plot(xCopy[:,1],yHat)
plt.show()
print corrcoef(yHat.T,yMat) #计算相关系数

 

最佳拟合直线方法将数据视为直线进行建模,具有不错的表现。但是还可以根据数据来进行局部调整预测

局部加权线性回归

  线性回归的一个问题是有可能出现欠拟合现象,因为它求的是具有最小均方误差的无偏估计。显而易见,如果模型欠拟合将不能取得最好的预测效果。所以将有些方法允许在估计中引入一些偏差,从而降低预测的均方误差

   其中的一个方法是局部加权线性回归(Locally Weighted Linear Regression,LWLR)。在该算法中,给待遇测点附件的每一个点赋予一定的权重,然后在这个子集上基于最小均方误差来进行普通的回归。与kNN一样,这种算法每次预测均需要事先选取出对应的数据子集。该算法解出回归系数w的形式如下:

  

其中,w是一个矩阵,用来给每个数据点赋予权重

  LWLR使用的核类型可以自由选择,最常用的核就是高斯核(RBF核),高斯核对应的权重如下:

  

 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
def lwlr(testPoint,xArr,yArr,k=1.0):		#局部加权线性回归函数
xMat = mat(xArr); yMat = mat(yArr).T
m = shape(xMat)[0]
weights = mat(eye((m))) #创建m维对角矩阵
for j in range(m):
diffMat = testPoint - xMat[j,:] #testPoint分别减去每一行xMat,testPoint也将是每一行的xMat
weights[j,j] = exp(diffMat*diffMat.T/(-2.0*k**2)) #权重值大小以指数级衰减,weights仍然只有对角线上有值
xTx = xMat.T * (weights * xMat)
if linalg.det(xTx) == 0.0:
print "This matrix is singular, cannot do inverse"
return
ws = xTx.I * (xMat.T * (weights * yMat))
return testPoint * ws #直接返回拟合出来的值,每一个样本的回归系数都不一样

def lwlrTest(testArr,xArr,yArr,k=1.0): #循环所有的样本,对其使用lwlr函数
m = shape(testArr)[0]
yHat = zeros(m)
for i in range(m):
yHat[i] = lwlr(testArr[i],xArr,yArr,k)
return yHat

main.py

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
# coding:utf-8
# !/usr/bin/env python

import regression
import matplotlib.pyplot as plt
from numpy import *

if __name__ == '__main__':
xArr,yArr = regression.loadDataSet("ex0.txt")

xMat = mat(xArr) #输入的X矩阵,n×(d+1)维
yMat = mat(yArr) #真实的标记,n×1维
yHat = regression.lwlrTest(xArr,xArr,yArr,k=0.01)
srtInd = xMat[:,1].argsort(0)
xSort = xMat[srtInd][:,0,:]
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0],s=5,c='red') #取得矩阵的第二列并转换成list
ax.plot(xSort[:,1],yHat[srtInd])
plt.show()

 

 

缩减系数(shrinking)

如果数据的特征比样本点还多的话,那么就不能使用线性回归和之前的方法来做预测了。因为在计算 的时候会出错。

如果特征比样本点还多(n>m),也就是说输入数据的矩阵x不是满秩矩阵。非满秩矩阵在求逆的时候会出现问题。

为了解决这个问题,引入了岭回归(ridge regression)的概念,这是一种缩减方法

简单地说,岭回归就是在矩阵 上加一个 从而使得矩阵非奇异,进而能对 求逆。其中矩阵 I 是一个m×m的单位矩阵,对角线上元素全为1,其他元素全为0。在这种情况下,回归系数的计算公式将变成:

  

岭回归最先用来处理特征数多于样本数的情况,现在也用于在估计中加入偏差,从而得到更好的估计。这里通过引入 来限制了所有 w 之和,通过引入该惩罚项,能够减少不重要的参数,这个技术在统计学中称为缩减(shrinkage)

缩减方法可以去掉不重要的参数,因此能更好地理解数据。此外,与简单的线性回归相比,缩减法能够取得更好的预测效果。

 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
def ridgeRegres(xMat,yMat,lam=0.2):			#岭回归,计算回归系数,λ默认等于0.2
xTx = xMat.T*xMat
denom = xTx + eye(shape(xMat)[1])*lam
if linalg.det(denom) == 0.0: #防止λ等于0
print "This matrix is singular, cannot do inverse"
return
ws = denom.I * (xMat.T*yMat)
return ws

def ridgeTest(xArr,yArr): #用于在一组不同的λ上分别测试结果,返回λ和8个特征值的关系曲线
xMat = mat(xArr); yMat=mat(yArr).T
yMean = mean(yMat,0) #求yMat的均值
yMat = yMat - yMean
#regularize X's
xMeans = mean(xMat,0) #求xMat的均值
xVar = var(xMat,0) #求xMat的方差
xMat = (xMat - xMeans)/xVar #数据归一化
numTestPts = 30 #30个λ的取值
wMat = zeros((numTestPts,shape(xMat)[1])) #shape(xMat)[1]表示特征的数量
for i in range(numTestPts):
ws = ridgeRegres(xMat,yMat,exp(i-10))
wMat[i,:]=ws.T
return wMat

 mian.py

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
# coding:utf-8
# !/usr/bin/env python

import regression
import matplotlib.pyplot as plt
from numpy import *

if __name__ == '__main__':
abX,abY = regression.loadDataSet("abalone.txt")
ridgeWeights = regression.ridgeTest(abX,abY)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(ridgeWeights)
plt.show()



 

30个不同的log(λ)所对应的回归系数。在最左边,即λ最小的时候,可以得到所有系数的原始值(与线性回归一致);而在右边,系数全部缩减到0,缩减了一些不重要的参数;在中间部分的某值将可以取得最好的预测效果。为了找到最佳参数值,还需要进行交叉验证。

下图是岭回归的回归系数变化图。可以在中间某处找到使得预测的结果最好的λ值。

还有一些其他缩减方法,如lasso、LAR、PCA回归以及子集选择等。与岭回归一样,这样方法不仅可以提高预测精确率,而且可以解释回归系数。

 

lasso方法

在增加下面约束的条件下,普通最小二乘法回归(OLS)将会得到与岭回归一样的公式:

  

在使用普通的最下二乘法回归在发那个两个或者更多的特征相关的时候,可能会得到一个很大的正系数和一个很大的负系数。正是因为上述限制条件的存在,使用岭回归可以便面这个问题。

lasso方法中也对系数进行了限定,对应的约束条件如下:

  

 

前向逐步回归

前向逐步回归算法可以得到和lasso算法差不多的效果,但是更加简单。它属于一种贪心算法,即每一步都尽可能减少误差。一开始,所有的权重都设为1,然后每一步所做的决策是对某个权重增加或者减少一个很小的值。

 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
def regularize(xMat):			#按均值为0方差为1进行标准化处理
inMat = xMat.copy()
inMeans = mean(inMat,0) #calc mean then subtract it off
inVar = var(inMat,0) #calc variance of Xi then divide by it
inMat = (inMat - inMeans)/inVar
return inMat

def stageWise(xArr,yArr,eps=0.01,numIt=100): #前向逐步线性回归,eps表示每次迭代需要调整的步长,numIt表示迭代次数
xMat = mat(xArr); yMat=mat(yArr).T
yMean = mean(yMat,0)
yMat = yMat - yMean #can also regularize ys but will get smaller coef
xMat = regularize(xMat) #把特征进行标准化
m,n=shape(xMat) #m是数据数量,n是特征数量
returnMat = zeros((numIt,n)) #testing code remove
ws = zeros((n,1)); wsTest = ws.copy(); wsMax = ws.copy()
for i in range(numIt): #迭代100次,每一次迭代只以步进的长度修改一个特征值
#print ws.T
lowestError = inf;
for j in range(n): #循环所有特征,最后找到步进后均方误差最小的保存wsMax
for sign in [-1,1]: #步进的正负
wsTest = ws.copy()
wsTest[j] += eps*sign
yTest = xMat*wsTest
rssE = rssError(yMat.A,yTest.A) #求均方误差
if rssE < lowestError: #比较
lowestError = rssE
wsMax = wsTest #保存均方误差最小的wsMax
ws = wsMax.copy()
returnMat[i,:]=ws.T #保存回归系数ws的整个修改过程
return returnMat

 

mian.py

1
2
3
4
5
6
7
8
9
10
11
12
13
# coding:utf-8
# !/usr/bin/env python

import regression
import Old_regression
import matplotlib.pyplot as plt
from numpy import *

if __name__ == '__main__':
abX,abY = regression.loadDataSet("abalone.txt")
returnMat = Old_regression.stageWise(abX,abY)
print returnMat