1. 程式人生 > >scikit-learn : 線性迴歸,多元迴歸,多項式迴歸

scikit-learn : 線性迴歸,多元迴歸,多項式迴歸

匹薩的直徑與價格的資料

%matplotlib inline
import matplotlib.pyplot as plt
def runplt():
    plt.figure()
    plt.title(u'diameter-cost curver')
    plt.xlabel(u'diameter')
    plt.ylabel(u'cost')
    plt.axis([0, 25, 0, 25])
    plt.grid(True)
    return plt

plt = runplt()
X = [[6], [8], [10], [14], [18]]
y = [[7
], [9], [13], [17.5], [18]] plt.plot(X, y, 'k.') plt.show()

這裡寫圖片描述

訓練模型

from sklearn.linear_model import LinearRegression
import numpy as np
# 建立並擬合模型
model = LinearRegression()
model.fit(X, y)
print('預測一張12英寸匹薩價格:$%.2f' % model.predict(np.array([12]).reshape(-1, 1))[0])
預測一張12英寸匹薩價格:$13.68

一元線性迴歸假設解釋變數和響應變數之間存線上性關係;這個線性模型所構成的空間是一個超平面(hyperplane)。

超平面是n維歐氏空間中餘維度等於一的線性子空間,如平面中的直線、空間中的平面等,總比包含它的空間少一維。

在一元線性迴歸中,一個維度是響應變數,另一個維度是解釋變數,總共兩維。因此,其超平面只有一維,就是一條線。

上述程式碼中sklearn.linear_model.LinearRegression類是一個估計器(estimator)。估計器依據觀測值來預測結果。在scikit-learn裡面,所有的估計器都帶有:
- fit()
- predict()

fit()用來分析模型引數,predict()是通過fit()算出的模型引數構成的模型,對解釋變數進行預測獲得的值。
因為所有的估計器都有這兩種方法,所有scikit-learn很容易實驗不同的模型。

一元線性迴歸模型:

y=α+βx

一元線性迴歸擬合模型的引數估計常用方法是:
- 普通最小二乘法(ordinary least squares )
- 線性最小二乘法(linear least squares)

首先,我們定義出擬合成本函式,然後對引數進行數理統計。

plt = runplt()
plt.plot(X, y, 'k.')
X2 = [[0], [10], [14], [25]]
model = LinearRegression()
model.fit(X, y)
y2 = model.predict(X2)
plt.plot(X, y, 'k.')
plt.plot(X2, y2, 'g-')
plt.show()

這裡寫圖片描述

plt = runplt()
plt.plot(X, y, 'k.')
y3 = [14.25, 14.25, 14.25, 14.25]
y4 = y2 * 0.5 + 5
model.fit(X[1:-1], y[1:-1])
y5 = model.predict(X2)
plt.plot(X, y, 'k.')
plt.plot(X2, y2, 'g-.')
plt.plot(X2, y3, 'r-.')
plt.plot(X2, y4, 'y-.')
plt.plot(X2, y5, 'o-')
plt.show()

這裡寫圖片描述

成本函式(cost function)也叫損失函式(loss function),用來定義模型與觀測值的誤差。模型預測的價格與訓練集資料的差異稱為殘差(residuals)或訓練誤差(training errors)。後面我們會用模型計算測試集,那時模型預測的價格與測試集資料的差異稱為預測誤差(prediction errors)或訓練誤差(test errors)。

模型的殘差是訓練樣本點與線性迴歸模型的縱向距離,如下圖所示:

plt = runplt()
plt.plot(X, y, 'k.')
X2 = [[0], [10], [14], [25]]
model = LinearRegression()
model.fit(X, y)
y2 = model.predict(X2)
plt.plot(X, y, 'k.')
plt.plot(X2, y2, 'g-')

# 殘差預測值
yr = model.predict(X)
for idx, x in enumerate(X):
    plt.plot([x, x], [y[idx], yr[idx]], 'r-')

plt.show()

這裡寫圖片描述

我們可以通過殘差之和最小化實現最佳擬合,也就是說模型預測的值與訓練集的資料最接近就是最佳擬合。對模型的擬合度進行評估的函式稱為殘差平方和(residual sum of squares)成本函式。就是讓所有訓練資料與模型的殘差的平方之和最小化,如下所示:

SSres=i=1n(yif(xi))2

其中, yi 是觀測值, f(xi)f(xi) 是預測值。

import numpy as np
print('殘差平方和: %.2f' % np.mean((model.predict(X) - y) ** 2))
殘差平方和: 1.75

解一元線性迴歸的最小二乘法

通過成本函式最小化獲得引數,我們先求相關係數 ββ 。按照頻率論的觀點,我們首先需要計算 xx 的方差和 xx 與 yy 的協方差。
方差是用來衡量樣本分散程度的。如果樣本全部相等,那麼方差為0。方差越小,表示樣本越集中,反正則樣本越分散。方差計算公式如下:

var(x)=ni=1(xix¯)2n1

Numpy裡面有var方法可以直接計算方差,ddof引數是貝塞爾(無偏估計)校正係數(Bessel’s correction),設定為1,可得樣本方差無偏估計量。

print(np.var([6, 8, 10, 14, 18], ddof=1))
23.2

協方差表示兩個變數的總體的變化趨勢。如果兩個變數的變化趨勢一致,也就是說如果其中一個大於自身的期望值,另外一個也大於自身的期望值,那麼兩個變數之間的協方差就是正值。 如果兩個變數的變化趨勢相反,即其中一個大於自身的期望值,另外一個卻小於自身的期望值,那麼兩個變數之間的協方差就是負值。如果兩個變數不相關,則協方差為0,變數線性無關不表示一定沒有其他相關性。協方差公式如下:

cov(x,y)=ni=1(xix¯)(yiy¯)n1

其中, x¯是直徑 x的均值, xi的訓練集的第 i個直徑樣本, y¯是價格y的均值, yi的訓練集的第i個價格樣本, n是樣本數量。Numpy裡面有cov方法可以直接計算方差。

import numpy as np
print(np.cov([6, 8, 10, 14, 18], [7, 9, 13, 17.5, 18])[0][1])
22.65

現在有了方差和協方差,就可以計算相關係統 β 了。

β=cov(x,y)var(x)

算出β後,我們就可以計算α了:

α=y¯βx¯

將前面的資料帶入公式就可以求出α了:

α=12.90.9762931×11.2=1.9655

模型評估

前面我們用學習演算法對訓練集進行估計,得出了模型的引數。有些度量方法可以用來評估預測效果,我們用R方(r-squared)評估匹薩價格預測的效果。R方也叫確定係數(coefficient of determination),表示模型對現實資料擬合的程度。計算R方的方法有幾種。一元線性迴歸中R方等於皮爾遜積矩相關係數(Pearson product moment correlation coefficient或Pearson’s r)的平方。種方法計算的R方一定介於0~1之間的正數。其他計算方法,包括scikit-learn中的方法,不是用皮爾遜積矩相關係數的平方計算的,因此當模型擬合效果很差的時候R方會是負值。下面我們用scikit-learn方法來計算R方。

SStot=i=1n(yiy¯)2
SSres=i=1n(yif(xi))2
R2=1SSresSStot

R方是0.6620說明測試集裡面過半數的價格都可以通過模型解釋。現在,讓我們用scikit-learn來驗證一下。LinearRegression的score方法可以計算R方:

# 測試集
X_test = [[8], [9], [11], [16], [12]]
y_test = [[11], [8.5], [15], [18], [11]]
model = LinearRegression()
model.fit(X, y)
model.score(X_test, y_test)
0.66200528638545164

多元迴歸

from sklearn.linear_model import LinearRegression
X = [[6, 2], [8, 1], [10, 0], [14, 2], [18, 0]]
y = [[7], [9], [13], [17.5], [18]]
model = LinearRegression()
model.fit(X, y)
X_test = [[8, 2], [9, 0], [11, 2], [16, 2], [12, 0]]
y_test = [[11], [8.5], [15], [18], [11]]
predictions = model.predict(X_test)
for i, prediction in enumerate(predictions):
    print('Predicted: %s, Target: %s' % (prediction, y_test[i]))
print('R-squared: %.2f' % model.score(X_test, y_test))
Predicted: [ 10.06250019], Target: [11]
Predicted: [ 10.28125019], Target: [8.5]
Predicted: [ 13.09375019], Target: [15]
Predicted: [ 18.14583353], Target: [18]
Predicted: [ 13.31250019], Target: [11]
R-squared: 0.77

多項式迴歸

上例中,我們假設解釋變數和響應變數的關係是線性的。真實情況未必如此。下面我們用多項式迴歸,一種特殊的多元線性迴歸方法,增加了指數項。現實世界中的曲線關係都是通過增加多項式實現的,其實現方式和多元線性迴歸類似。本例還用一個解釋變數,匹薩直徑。讓我們用下面的資料對兩種模型做個比較:

import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import PolynomialFeatures
X_train = [[6], [8], [10], [14], [18]]
y_train = [[7], [9], [13], [17.5], [18]]
X_test = [[6], [8], [11], [16]]
y_test = [[8], [12], [15], [18]]
# 建立線性迴歸,並用訓練的模型繪圖
regressor = LinearRegression()
regressor.fit(X_train, y_train)
xx = np.linspace(0, 26, 100)
yy = regressor.predict(xx.reshape(xx.shape[0], 1))
plt = runplt()
plt.plot(X_train, y_train, 'k.')
plt.plot(xx, yy)

quadratic_featurizer = PolynomialFeatures(degree=2)
X_train_quadratic = quadratic_featurizer.fit_transform(X_train)
X_test_quadratic = quadratic_featurizer.transform(X_test)
regressor_quadratic = LinearRegression()
regressor_quadratic.fit(X_train_quadratic, y_train)
xx_quadratic = quadratic_featurizer.transform(xx.reshape(xx.shape[0], 1))
plt.plot(xx, regressor_quadratic.predict(xx_quadratic), 'r-')
plt.show()
print(X_train)
print(X_train_quadratic)
print(X_test)
p