python實戰之線性回歸、局部加權回歸
1.基本概念與思想
回歸：求回歸方程中回歸系數的過程稱為回歸。
局部加權思想：給待預測點附近的每個點賦予一定的權重。
2.線性回歸
回歸方程的解： Θ=(XTX)-1XTY (1)
其中，Θ表示回歸系數矩陣，X表示樣本矩陣，Y表示樣本類標矩陣。
3.局部加權回歸
回歸方程解： Θ=(XTWX)-1XTWY (2)
其中，（2）與（1）不同的是多了表示局部權重的矩陣W。另外，
w(i,i)=exp((x(i)-x)2/-2k2) (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()
轉載請注明出處，謝謝！