python--线性回归、局部加权回归

python实战之线性回归、局部加权回归

1.基本概念与思想

回归：求回归方程中回归系数的过程称为回归。

局部加权思想：给待预测点附近的每个点赋予一定的权重。

2.线性回归

　　回归方程的解： Θ=(X^TX)^-1X^TY (1)

　　其中，Θ表示回归系数矩阵，X表示样本矩阵，Y表示样本类标矩阵。

3.局部加权回归

　　回归方程解： Θ=(X^TWX)^-1X^TWY (2)

　　其中，（2）与（1）不同的是多了表示局部权重的矩阵W。另外，

　　w(i,i)=exp((x⁽ⁱ⁾-x)²/-2k²) (3)

4.python代码实现

  1 #!/usr/bin/python
  2 #-*- coding:utf-8 -*-
  3 
  4 from numpy import *
  5 import matplotlib.pyplot as plt
  6 '''
  7 解决python matplotlib画图无法显示中文的问题！
  8 '''
  9 from pylab import *
 10 mpl.rcParams['font.sans-serif']=['SimHei']
 11 mpl.rcParams['axes.unicode_minus']=False
 12 #############################################################
 13 
 14 def loadDataSet(fileName):
 15     numFeat=len(open(fileName).readline().split('\t'))-1
 16     dataMat=[];labelMat=[]                     #空列表
 17     #为了减少内存消耗，一行行读入
 18     #fr=open(fileName)
 19     #for line in fr.readlines():
 20     for line in open(fileName):
 21         lineArr=[]
 22         curLine=line.strip().split('\t')
 23         for i in range(numFeat):
 24             lineArr.append(float(curLine[i]))   #list方法append()
 25         dataMat.append(lineArr)
 26         labelMat.append(float(curLine[-1]))
 27     return dataMat,labelMat
 28 #线性回归(lr)主函数
 29 def standRegression(xArr,yArr):
 30     xMat=mat(xArr);yMat=mat(yArr).T #注意此处需要转置
 31     xTx=xMat.T*xMat
 32     if linalg.det(xTx)==0.0:   #若行列式为0，则不可逆
 33         print 'This matrix is singular,cannot do inverse'
 34         return
 35     ws=xTx.I*(xMat.T*yMat)
 36     return ws
 37 #lr绘图
 38 def lrPlot(xArr,yArr,ws):
 39     xMat=mat(xArr);yMat=mat(yArr)
 40     fig=plt.figure()
 41     ax=fig.add_subplot(111)
 42     #scatter()画散点图，yMat绘制原始图
 43     ax.scatter(xMat[:,1].flatten().A[0],yMat.T[:,0].flatten().A[0],label='yMat') #scatter(x,y),x和y必须转化为1-D
 44     #yHat绘制拟合图
 45     xCopy=xMat.copy()
 46     xCopy.sort(0)      #按列排序，画直线图需要排序
 47     yHat=xCopy*ws
 48     #plot()画直线图
 49     ax.plot(xCopy[:,1],yHat,label='yHat')
 50 
 51     plt.legend(loc='down left')   #指定方框位置
 52     plt.show()
 53 
 54 #局部线性回归(lwlr)主函数
 55 def lwlr(testPoint,xArr,yArr,k=1.0):
 56     xMat=mat(xArr);yMat=mat(yArr).T
 57     m=shape(xMat)[0]
 58     weights=mat(eye(m))
 59     for j in range(m):
 60         diffMat=testPoint-xMat[j,:]
 61         weights[j,j]=exp(diffMat*diffMat.T/(-2*k**2))
 62     xTwx=xMat.T*(weights*xMat)  #
 63     if linalg.det(xTwx)==0.0:
 64         print 'This matrix is singular,cannot do inverse'
 65         return
 66     ws=xTwx.I*(xMat.T*(weights*yMat))
 67     return testPoint*ws
 68 def lwlrTest(testArr,xArr,yArr,k=1.0):
 69     #testArr=mat(testArr)     #转换为矩阵
 70     m=shape(testArr)[0]
 71     yHat=zeros(m)
 72     for i in range(m):
 73         yHat[i]=lwlr(testArr[i],xArr,yArr,k)   #testArr[i]列表索引
 74     return yHat
 75 #lwlr绘图测试版
 76 def lwlrPlot(xArr,yArr,yHat):
 77     xMat=mat(xArr)
 78     srtInd=xMat[:,1].argsort(0)   #等价于argsort(xMat[:,1],0)
 79     xSort=xMat[srtInd][:,0,:]     #等价于xMat[srtInd.flatten().A[0]] 
 80 
 81     fig=plt.figure()
 82     ax=fig.add_subplot(111)
 83 
 84     #直线图plt.plot(),画plot前要排序
 85     #ax.plot(xMat[:,1],yHat[:].T)
 86     ax.plot(xSort[:,1],yHat[srtInd]) 
 87 
 88     #画散点图不需要排序
 89     ax.scatter(xMat[:,1].flatten().A[0],mat(yHat).T.flatten().A[0],s=2,c='k')
 90     ax.scatter(xMat[:,1].flatten().A[0],mat(yArr).T.flatten().A[0],s=2,c='r')  #散点图plt.scatter()
 91 
 92     plt.show()
 93 #lwlr绘图比较3个yHat(k取1,0.01,0.003)
 94 def lwlrPlot3(xArr,yArr):               #输入：xArr是n×d矩阵/数组/列表；yArr是n×1
 95 
 96     xMat=mat(xArr)
 97     srtInd=xMat[:,1].argsort(0)         #等价于argsort(xMat[:,1],0)
 98     xSort=xMat[srtInd][:,0,:]           #等价于xMat[srtInd.flatten().A[0]] 
 99 
100     yHat1=lwlrTest(xArr,xArr,yArr,1)    #调用局部加权回归(lwlr)主函数
101     yHat2=lwlrTest(xArr,xArr,yArr,0.01)
102     yHat3=lwlrTest(xArr,xArr,yArr,0.03)
103 
104     fig=plt.figure()
105     ax1=fig.add_subplot(311)
106     ax2=fig.add_subplot(312)
107     ax3=fig.add_subplot(313)
108 
109     #画直线图需要排序
110     #直线图plt.plot(),plot前要排序
111     #ax1.plot(xMat[:,1],yHat[:].T)
112     ax1.plot(xSort[:,1],yHat1[srtInd]) 
113     ax2.plot(xSort[:,1],yHat2[srtInd])
114     ax3.plot(xSort[:,1],yHat3[srtInd])
115 
116     #画散点图不需要排序
117     ax1.scatter(xMat[:,1].flatten().A[0],mat(yArr).T.flatten().A[0],s=2,c='r',label=u'欠拟合')  #散点图plt.scatter()
118     ax2.scatter(xMat[:,1].flatten().A[0],mat(yArr).T.flatten().A[0],s=2,c='r',label=u'最好')
119     ax3.scatter(xMat[:,1].flatten().A[0],mat(yArr).T.flatten().A[0],s=2,c='r',label=u'过拟合')
120 
121     ax1.legend(loc='upper left')
122     ax2.legend(loc='upper left')
123     ax3.legend(loc='upper left')
124 
125     plt.show()

代码主要函数：loadDataSet(), standRegression(), lwlr(), lwlrTest()

转载请注明出处，谢谢！