预测数值型数据：回归

regress.py

标准回归

#!/usr/bin/python
# -*- coding: utf-8 -*-
#coding=utf-8

from numpy import *

#导入函数
#前几行为x，最后一行为y
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(lineArr)
labelMat.append(float(curLine[-1]))
return dataMat, labelMat

#标准回归函数
def standRegress(xArr, yArr):
xMat = mat(xArr)
yMat = mat(yArr).T
xTx = xMat.T * xMat
if linalg.det(xTx) == 0.0:  #如果行列式为0，则逆不存在
print "This matrix is singular, cannot do inverse"
return
ws = xTx.I * (xMat.T * yMat)  #回归系数
return ws

#测试标准回归函数
def testStandRetress(xArr, yArr):
xArr, yArr = loadDataSet('ex0.txt')
ws = standRegress(xArr, yArr)
print ws
xMat = mat(xArr)
yMat = mat(yArr)   #真实值
import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(111)
#原始图像
ax.scatter(xMat[:,1].flatten().A[0], yMat.T[:,0].flatten().A[0])
xCopy = xMat.copy()
xCopy.sort(0)
yHat = xCopy * ws  #预测值
#预测图像
ax.plot(xCopy[:,1], yHat)
plt.show()

测试：

>>> import regress
>>> xArr, yArr = loadDataSet('ex0.txt')
>>> testStandRetress(xArr, yArr)
w= [[ 3.00774324]
[ 1.69532264]]

corrcoef= [[ 1.          0.98647356]
[ 0.98647356  1.        ]]   #yHat与yMat的相关系数为0.986

局部加权线性回归

#局部加权线性回归
#k决定了对附近的嗲赋予多大的权重
def lwlr(testPoint, xArr, yArr, k=1.0):
xMat = mat(xArr)
yMat = mat(yArr).T
m = shape(xMat)[0]
weights = mat(eye((m)))  #单位矩阵
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, cannot do inverse"
return
ws = xTx.I * (xMat.T * (weights * yMat))
return testPoint * ws

#测试局部加权线性回归函数
def lwlrTest(xArr, yArr):
print "yArr[0]", yArr[0]
print "yHrr[0], k=0.01", lwlr(xArr[0], xArr, yArr, 0.01)
print "yHrr[0], k=0.5", lwlr(xArr[0], xArr, yArr, 0.5)
print "yHrr[0], k=1", lwlr(xArr[0], xArr, yArr, 1)
m = shape(xArr)[0]
yHat = zeros(m)
for i in range(m):
yHat[i] = lwlr(xArr[i], xArr, yArr, k=0.01)
xMat = mat(xArr)
srtInd= xMat[:,1].argsort(0) #对x排序
xSort = xMat[srtInd][:,0,:]
import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(xSort[:,1], yHat[srtInd])  #预测图像
ax.scatter(xMat[:,1].flatten().A[0], mat(yArr).T.flatten().A[0], s=2, c='red') #真实图像
plt.show()

测试：

>>> import regress
>>> xArr, yArr = loadDataSet('ex0.txt')
>>> lwlrTest(xArr, yArr)
yArr[0] 3.176513
yHrr[0], k=0.01 [[ 3.20366661]]
yHrr[0], k=0.5 [[ 3.12201662]]
yHrr[0], k=1 [[ 3.12204471]]

缩减法：岭回归

缩减法可将系数缩减的很小的值直接缩减为0

#缩减法：岭回归
def ridgeRegres(xMat, yMat, lam=0.2):
xTx = xMat.T * xMat
denom = xTx + eye(shape(xMat[1])[1]) * lam
if linalg.det(denom) == 0.0:
print "This matrix is singular, cannot do inverse"
return
ws = denom.I * (xMat.T * yMat)
return ws

def ridgeTest(xArr, yArr):
xMat = mat(xArr)
yMat = mat(yArr).T
#数据标准化
yMean = mean(yMat, 0) #均值
xMean = mean(xMat, 0)
yMat = yMat - yMean
xVar = var(xMat, 0)  #方差
xMat = (xMat - xMean)/xVar
numTestPts = 30
wMat = zeros((numTestPts, shape(xMat)[1]))
for i in range(numTestPts):  #在30个不同的lambda下调用ridgeRegres
ws = ridgeRegres(xMat, yMat, exp(i-10))  #lambda以指数级变化
wMat[i,:] = ws.T  #将所有的回归系数输出到一个矩阵
print wMat
import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(wMat)
plt.show()

#交叉测试岭回归
#numVal为算法中交叉验证的次数
def crossValidation(xArr, yArr, numVal = 10):
m = len(yArr)
indexList = range(m)
errorMat = zeros((numVal, 30))
#创建训练集和测试集容器
for i in range(numVal):
trainX = []
trainY = []
testX = []
testY = []
random.shuffle(indexList)
for j in range(m):
#将数据分为训练集和测试集
if j < m * 0.9:  #90%的数据为训练集，10%为测试集
trainX.append(xArr[indexList[j]])
trainY.append(yArr[indexList[j]])
else:
testX.append(xArr[indexList[j]])
testY.append(yArr[indexList[j]])
wMat = ridgeTest(trainX, trainY)  #保存所有回归系数
#用训练时的参数将测试数据标准化
for k in range(30):
matTestX = mat(testX)
matTrainX = mat(trainX)
meanTrainX = mean(matTrainX, 0)
varTrainX = var(matTrainX, 0)
matTestX = (matTestX - meanTrainX) / varTrainX
yEst = matTestX * mat(wMat[k,:]).T + mean(trainY)
errorMat[i,k] = rssError(yEst.T.A, array(testY))  #误差
meanErrors = mean(errorMat, 0)
minMean = float(min(meanErrors))
bestWeights = wMat[nonzero(meanErrors == minMean)]
xMat = mat(xArr)
yMat = mat(yArr).T
meanX = mean(xMat, 0)
varX = var(xMat, 0)
#因为作了标准化，需进行数据还原
unReg = bestWeights / varX
print "the best model from Ridge Regression is:\n", unReg
print "with constant term: ", -1 * sum(multiply(meanX, unReg)) + mean(yMat)

测试：

>>> import regress
>>> xArr, yArr = loadDataSet('abalone.txt')
>>> ridgeTest(xArr, yArr)
[[  4.30441949e-02  -2.27416346e-02   1.32140875e-01   2.07518171e-02 2.22403745e+00  -9.98952980e-01  -1.17254237e-01   1.66229222e-01]
[  4.30441928e-02  -2.27416370e-02   1.32140878e-01   2.07518175e-02 2.22403626e+00  -9.98952746e-01  -1.17254174e-01   1.66229339e-01]
[  4.30441874e-02  -2.27416435e-02   1.32140885e-01   2.07518187e-02 2.22403305e+00  -9.98952110e-01  -1.17254003e-01   1.66229656e-01]

... ...

得到8个特征值的系数与log(lambda)的关系图像。最左端，lambda最小时，所有系数值与线性回归一致；最右边，都缩减成0；在中间的某值可以取得最好的预测效果

crossValidation(xArr, yArr, numVal = 10)

the best model from Ridge Regression is:

[[ 0.0823653 -3.91472054 16.9506189 10.7538564 8.09628256 -19.35008479 -7.85980015 9.56536757]]

with constant term: 2.94541561804

缩减法：前向逐步线性回归

#缩减法：前向逐步线性回归
def stageWise(xArr, yArr, eps=0.01, numItr=100):
xMat = mat(xArr)
yMat = mat(yArr).T
#数据标准化
yMean = mean(yMat, 0)
yMat = yMat - yMean
xMat = regularize(xMat)
m, n = shape(xMat)
returnMat = zeros((numItr, n))
ws = zeros((n, 1))
wsTest = ws.copy()
wsMax = ws.copy()
for i in range(numItr):
print ws.T
lowestError = inf
for j in range(n):  #对于每个特征值
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
ws = wsMax.copy()
returnMat[i, :] = ws.T
return returnMat

def stageWiseTest():
xArr, yArr = loadDataSet('abalone.txt')
wMat = stageWise(xArr, yArr, eps=0.01, numItr=200)
print "前向逐步线性回归，eps=0.01, numTtr=200:", wMat
wMat = stageWise(xArr, yArr, eps=0.001, numItr=5000)
print  "前向逐步线性回归，eps=0.001, numTtr=5000:", wMat
#与最小二乘法比较
xMat = mat(xArr)
yMat = mat(yArr).T
yMean = mean(yMat, 0)
xMat = regularize(xMat)
yMat = yMat - yMean
wMat = standRegress(xMat, yMat.T)
print "标准回归", wMat.T

测试：

>>> import regress
>>> stageWiseTest()

前向逐步线性回归，eps=0.01, numTtr=200:
[[ 0.    0.    0.   ...,  0.    0.    0.  ]
[ 0.    0.    0.   ...,  0.    0.    0.  ]
[ 0.    0.    0.   ...,  0.    0.    0.  ]
...,
[ 0.05  0.    0.09 ..., -0.64  0.    0.36]
[ 0.04  0.    0.09 ..., -0.64  0.    0.36]
[ 0.05  0.    0.09 ..., -0.64  0.    0.36]]
#第2列和第7列都是0，说明不对目标值造成任何影响

前向逐步线性回归，eps=0.001, numTtr=5000:
[[ 0.     0.     0.    ...,  0.     0.     0.   ]
[ 0.     0.     0.    ...,  0.     0.     0.   ]
[ 0.     0.     0.    ...,  0.     0.     0.   ]
...,
[ 0.043 -0.011  0.12  ..., -0.963 -0.105  0.187]
[ 0.044 -0.011  0.12  ..., -0.963 -0.105  0.187]

标准回归：
[[ 0.0430442  -0.02274163  0.13214087  0.02075182  2.22403814 -0.99895312 -0.11725427  0.16622915]]