1、問題引入
在統計學中,線性回歸是利用稱為線性回歸方程的最小二乘函數對一個或多個自變量和因變量之間關系進行建模的一種回歸分析。這種函數是一個或多個稱為回歸系數的模型參數的線性組合。一個帶有一個自變量的線性回歸方程代表一條直線。我們需要對線性回歸結果進行統計分析。
例如,假設我們已知一些學生年紀和游戲時間的數據,可以建立一個回歸方程,輸入一個新的年紀時,預測該學生的游戲時間。自變量為學生年紀,因變量為游戲時間。當只有一個因變量時,我們稱該類問題為簡單線性回歸。當游戲時間與學生年紀和學生性別有關,因變量有多個時,我們稱該類問題為多元線性回歸。
2、常見的統計量
在研究該問題之前,首先了解下編程中用到的常見的統計量:
序號 |
概念 |
公式 |
算法 |
說明 |
1 |
均值 |
|
整體的均值 |
|
2 |
中位數 |
|
排序后取中間值 |
|
3 |
眾數 |
|
出現次數最多的數 |
出現頻率 |
4 |
方差 |
|
數據的離散程度 |
|
5 |
標准差 |
|
s |
方差的開方 |
2、簡單線性回歸實例及編程實現
研究一個自變量(X)和一個因變量(y)的關系
簡單線性回歸模型定義:
簡單線性回歸方程:
其中:
為回歸線的截距
為回歸線的斜率
通過訓練數據,求取出估計參數建立的直線方程:
實際編程時,主要是根據已知訓練數據,估計出和
的值
和
以下面實例為例,第一列表示每月投放廣告的次數,第二列表示汽車向量,通過Python編程求取線性回歸方程:
投放廣告數 |
汽車銷量 |
1 |
14 |
3 |
24 |
2 |
18 |
1 |
17 |
3 |
27 |
編程關鍵在於如何求取b0和b1的值,我們引入一個方程(sum of square):
當上述方程的值最小時,我們認為求取到線程回歸方程參數的值,對該方程求最小值可以進一步轉化為求導和求極值的問題,求導過程省略,最后結論如下:
實際代碼:
import numpy as np from matplotlib import pylab as pl # 定義訓練數據 x = np.array([1,3,2,1,3]) y = np.array([14,24,18,17,27]) # 回歸方程求取函數 def fit(x,y): if len(x) != len(y): return numerator = 0.0 denominator = 0.0 x_mean = np.mean(x) y_mean = np.mean(y) for i in range(len(x)): numerator += (x[i]-x_mean)*(y[i]-y_mean) denominator += np.square((x[i]-x_mean)) print('numerator:',numerator,'denominator:',denominator) b0 = numerator/denominator b1 = y_mean - b0*x_mean return b0,b1 # 定義預測函數 def predit(x,b0,b1): return b0*x + b1 # 求取回歸方程 b0,b1 = fit(x,y) print('Line is:y = %2.0fx + %2.0f'%(b0,b1)) # 預測 x_test = np.array([0.5,1.5,2.5,3,4]) y_test = np.zeros((1,len(x_test))) for i in range(len(x_test)): y_test[0][i] = predit(x_test[i],b0,b1) # 繪制圖像 xx = np.linspace(0, 5) yy = b0*xx + b1 pl.plot(xx,yy,'k-') pl.scatter(x,y,cmap=pl.cm.Paired) pl.scatter(x_test,y_test[0],cmap=pl.cm.Paired) pl.show()
藍色表示測試數據,橙色表示預測數據。
3、多元線性回歸實例及編程實現
多元線性回歸方程和簡單線性回歸方程類似,不同的是由於因變量個數的增加,求取參數的個數也相應增加,推導和求取過程也不一樣。
y=β0+β1x1+β2x2+ ... +βpxp+ε
對於b0、b1、…、bn的推導和求取過程,引用一個第三方庫進行計算。以如下數據為例,對運輸里程、運輸次數與運輸總時間的關系,建立多元線性回歸模型:
運輸里程 |
運輸次數 |
運輸總時間 |
100 |
4 |
9.3 |
50 |
3 |
4.8 |
100 |
4 |
8.9 |
100 |
2 |
6.5 |
50 |
2 |
4.2 |
80 |
2 |
6.2 |
75 |
3 |
7.4 |
65 |
4 |
6.0 |
90 |
3 |
7.6 |
90 |
2 |
6.1 |
代碼如下:
import numpy as np from sklearn import datasets,linear_model # 定義訓練數據 x = np.array([[100,4,9.3],[50,3,4.8],[100,4,8.9], [100,2,6.5],[50,2,4.2],[80,2,6.2], [75,3,7.4],[65,4,6],[90,3,7.6],[90,2,6.1]]) print(x) X = x[:,:-1] Y = x[:,-1] print(X,Y) # 訓練數據 regr = linear_model.LinearRegression() regr.fit(X,Y) print('coefficients(b1,b2...):',regr.coef_) print('intercept(b0):',regr.intercept_) # 預測 x_test = np.array([[102,6],[100,4]]) y_test = regr.predict(x_test) print(y_test)
如果特征向量中存在分類型變量,例如車型,我們需要進行特殊處理:
運輸里程 |
輸出次數 |
車型 |
隱式轉換 |
運輸總時間 |
100 |
4 |
1 |
010 |
9.3 |
50 |
3 |
0 |
100 |
4.8 |
100 |
4 |
1 |
010 |
8.9 |
100 |
2 |
2 |
001 |
6.5 |
50 |
2 |
2 |
001 |
4.2 |
80 |
2 |
1 |
010 |
6.2 |
75 |
3 |
1 |
010 |
7.4 |
65 |
4 |
0 |
100 |
6.0 |
90 |
3 |
0 |
100 |
7.6 |
100 |
4 |
1 |
010 |
9.3 |
50 |
3 |
0 |
100 |
4.8 |
100 |
4 |
1 |
010 |
8.9 |
100 |
2 |
2 |
001 |
6.5 |
import numpy as np from sklearn.feature_extraction import DictVectorizer from sklearn import linear_model # 定義數據集 x = np.array([[100,4,1,9.3],[50,3,0,4.8],[100,4,1,8.9], [100,2,2,6.5],[50,2,2,4.2],[80,2,1,6.2], [75,3,1,7.4],[65,4,0,6],[90,3,0,7.6], [100,4,1,9.3],[50,3,0,4.8],[100,4,1,8.9],[100,2,2,6.5]]) x_trans = [] for i in range(len(x)): x_trans.append({'x1':str(x[i][2])}) vec = DictVectorizer() dummyX = vec.fit_transform(x_trans).toarray() x = np.concatenate((x[:,:-2],dummyX[:,:],x[:,-1].reshape(len(x),1)),axis=1) x = x.astype(float) X = x[:,:-1] Y = x[:,-1] print(x,X,Y) # 訓練數據 regr = linear_model.LinearRegression() regr.fit(X,Y) print('coefficients(b1,b2...):',regr.coef_) print('intercept(b0):',regr.intercept_)