1.用线性回归找到最佳拟合直线
我们知道回归方程的公式如下所示:
如果输入向量x满足以下形式:
并且我们有如下形式的回归系数向量w:
那么回归方程又可以写成以下形式:
那么我们如何求取回归系数向量w呢?
一个常用的方法就是找出使得误差最小的w,我们通常采用平方误差。平方误差可以写成:
用矩阵表示还可以写做(y-xw)T(y-xw)。我们对w求导并令其等于零,解出w如下:
值得注意的是,上述式子中的(xTx) 的逆矩阵不一定存在,所以在程序中我们在计算w之前必须要判断矩阵(xTx)是否可逆,即np.linalg.det(xTx)==0.0?
代码如下所示:
import numpy as np
import matplotlib.pyplot as plt
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 standRegres(xArr,yArr):
xMat=np.mat(xArr);yMat=np.mat(yArr).T
xTx=xMat.T*xMat
if np.linalg.det(xTx)==0.0:
print("This matrix is singular,cannot do inverse")
return
ws=xTx.I*(xMat.T*yMat)
return ws
if __name__=="__main__":
xArr,yArr=loadDataSet("ex0.txt")
ws=standRegres(xArr,yArr)
xMat=np.mat(xArr);yMat=np.mat(yArr)
yHat=xMat*ws
print("预测值和真实值之间的相关系数:\n",np.corrcoef(yHat.T,yMat))
fig=plt.figure()
ax=fig.add_subplot(111)
ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0],label="the original data")
xValue=xMat[:,1]
xValue.sort()
ax.plot(xValue,yHat,color="r",linewidth=3,label="fitting straight-line")
plt.grid()
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()
预测值和真实值之间的相关系数:
[[1. 0.98647356]
[0.98647356 1. ]]
我们可以看到该相关系数矩阵的主对角线元素都是1.0,因为yMat(真实值)与自身匹配最完美,而yMat与yHat(预测值)的相关系数是0.98647356。
2.局部加权线性回归
对于上述方法我们发现每个维度上的输入向量的回归系数都是一样的,这就会造成欠拟合的现象,无法得到准确的预测结果。 现在我们将输入向量左乘上一个权值矩阵,之后再右乘上回归系数向量,得到最终的预测结果,即现有如下式子成立:
这里权值矩阵W是对角矩阵,通过上述式子我们发现现在每个输入变量前的回归系数以及回归方程的偏置项都不相同,这就很好地避免了欠拟合。对于w回归系数向量的求法依然跟第一种方法相同,通过对平方误差求导取零,得到w的式子如下:
然而权值矩阵W我们是通过高斯核来对每个点及该点邻近的点赋予更高的权值。高斯核对应的权重如下:
其中k值决定高斯函数的衰减速度。当k值一定时,实际上我们是利用高斯核函数从第一个输入变量开始平移直到最后一个输入变量。这一点有点类似于滤波器,我们将高斯核函数当做滤波器,使得该点及其邻近的点通过滤波器,而远离该点的点将通过滤波器滤除掉,即该点及其邻近的点的权值较大接近于1,而远离该点的点的权值接近于零。我们是要通过该点及其邻近的点来确定该点处的回归系数。
代码如下所示:
import numpy as np
import matplotlib.pyplot as plt
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 lwlr(testPoint,xArr,yArr,k=1.0):
xMat=np.mat(xArr);yMat=np.mat(yArr).T
m=np.shape(xMat)[0]
weigths=np.mat(np.eye(m))
for j in range(m):
diffMat=testPoint-xMat[j,:]
weigths[j,j]=np.exp(diffMat*diffMat.T/(-2.0*k**2))
xTx=xMat.T*(weigths*xMat)
if np.linalg.det(xTx)==0.0: #判断矩阵是否可逆
print("This matrix is singular,cannot do inverse")
return
ws=xTx.I*(xMat.T*(weigths*yMat))
return testPoint*ws
def lwlrTest(testArr,xArr,yArr,k=1.0):
m=np.shape(testArr)[0]
yHat=np.zeros(m)
for i in range(m):
yHat[i]=lwlr(testArr[i],xArr,yArr,k)
return yHat
if __name__=="__main__":
xArr,yArr=loadDataSet("ex0.txt")
yHat=lwlrTest(xArr,xArr,yArr,0.01)
xMat = np.mat(xArr)
yMat = np.mat(yArr)
print("预测值和真实值之间的相关系数:\n", np.corrcoef(yHat.T, yMat))
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(xMat[:, 1].flatten().A[0], yMat.T[:, 0].flatten().A[0], label="the original data")
strInd=xMat[:,1].argsort(0)
xValue=xMat[strInd][:,0,:]
ax.plot(xValue[:,1], yHat[strInd], color="r", linewidth=3, label="weighted linear regression")
plt.grid()
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()
通过上述运行结果我们发现,k=1.0的时候如同将数据视为等权重,得到的拟合直线与第一种方法的结果一致。当k=0.01时,得到了非常好的拟合效果,能够体现处数据的潜在规律。而k=0.003则有些过拟合,拟合结果与数据点过于贴近,不能体现出数据的潜在规律。