多元線性回歸與模型診斷

一、多元線性回歸

1.多元線性回歸的基本表達式

註：額外的假設條件

①解釋變量之間不能存在太強的線性相關關係（一般ρ<0.7）

②其他條件與一元線性回歸類似。

2.回歸方程的模型擬合度

在進行回歸模型之前，我們可以計算總的波動誤差如下：

    在運用回歸模型後，總誤差可以分解為以下兩種：

    1）不能解釋變量誤差(誤差項的平方)

    2）被回歸方程可以解釋的誤差

根據以上解釋我們可以整理出以下方差分析表：

df自由度

SS

MSS（方差）

Explained可以解釋的誤差

k

ESS

ESS/k

Residual 誤差項

n-k-1

RSS

RSS/(n-k-1)

Total 總誤差

n-1

TSS

       3)擬合判定係數R方

註：一元線性回歸的R方比較特殊，為相關係數的平方：

通常情況下R方越大，說明模型擬合的解釋力度越好，但有時通過不斷增加解釋變量可以相應提高R方，這時候擬合效果並不一定是最好，所以提出以下修正R方來進行判斷：

3.回歸方程的推論

①置信區間（同一元線性回歸類似）

②假設檢驗（同一元線性回歸類似）

③聯合假設檢驗

通常F檢驗被用於檢測多元線性回歸的解釋變量斜率是否同時等於0：

判斷準則：如果F(T檢驗計算值)>F（關鍵值），則拒絕原假設H0.

Python案例分析：

#導入cv模塊import numpy as npimport matplotlib as mplimport matplotlib.pyplot as pltimport pandas as pddata=pd.read_excel('D:/歐晨的金學智庫/線性回歸案例.xlsx',header=0,index_col=0)data.head()
import statsmodels.api as sm model3 = sm.OLS(data['weight'],sm.add_constant(data[['age','height']]))result = model3.fit()result.summary()  print(result.summary())
分析:經過調整的R方=0.639，擬合係數一般，並不是特別好。
F統計量=26.67（檢驗解釋變量係數是否全為0？），P-value=4.05*10-7，非常小，拒絕原假設；
常數項估計-140.762，P值很小，說明截距項顯著，不為0；
age項估計-0.0485，P值較大>0.05，說明age項不顯著，可以嘗試剔除該解釋變量；
height項1.1973，P值很小，說明height項顯著，不為0。
二、模型診斷
1.自變量選擇偏差的權衡
（1）丟失重要變量
（2）加入無關變量
變量係數的估計偏差（大樣本，無關變量會收斂於0）
增加了模型參數估計的不確定性
增加了R方的值，但是使得調整的R方減小
（3）兩種合理估計線性回歸係數的方法
①一般情況模型變量的選擇方法
a.將所有變量加入進行回歸；
b.移除擬合效果最差的一個變量（尤其是不顯著的變量）；
c.移除後繼續採用線性回歸模型進行擬合，再次移除不顯著的變量；
d.重複以上步驟，直至所有變量的擬合結果都顯著；
註：通常選擇顯著性α在1%~0.1%（相應t值至少為2.57或3.29）
Python案例分析：在上述案例中，我們得到age項不顯著，可以剔除該解釋變量，只用height進行線性回歸：
import statsmodels.api as sm model1 = sm.OLS(data['weight'], sm.add_constant(data['height']))result = model1.fit()result.summary()  print(result.summary())
分析：經過調整的R方=0.652>0.639，說明剔除age變量後，擬合效果更好。剩餘的截距項和身高的P值均很小，說明顯著不為0，所以應當保留。 
②K折交叉檢驗
    a.確定模型數量（有n個解釋變量——每個變量選擇有或無，通常有2^n個模型）
    b.將數據分成相等數量的k個集合，其中k-1個集合作為訓練集擬合回歸方程，剩下的1個集合作為驗證集；重複進行交叉擬合驗證（總計有k次）。
    c.每個模型都採用b的方式進行驗證。
    d.計算每個模型的總的殘差項（驗證k次）的和，選擇殘差項和或其均值最小的一組模型最為最優模型。
Python案例分析：依舊使用以上案例，有2個解釋變量，所以應當有2^2=4個模型，我們排除解釋變量均為0的情況，來做3個模型的K折交叉檢驗：
y=data.weightX1=data[['age','height']]from sklearn.linear_model import LinearRegression lm = LinearRegression()from sklearn.model_selection import cross_val_score scores1 = cross_val_score(lm, X1, y, cv=10)  print (scores1)
分析：上面的scores都是負數，為什麼均方誤差會出現負數的情況呢？因為這裡的mean_squared_error是一種損失函數，優化的目標的使其最小化，而分類準確率是一種獎勵函數，優化的目標是使其最大化，因而選擇只用height做變量。
print('用age,height做自變量：',np.sqrt(-scores1.mean()))np.sqrt(-scores1.mean())X2=data[['age']]scores2 = cross_val_score(lm, X2, y, cv=10)print('只用age做自變量：',np.sqrt(-scores2.mean()))X3=data[['height']]scores3 = cross_val_score(lm, X3, y, cv=10)print('只用height做自變量：',np.sqrt(-scores3.mean()))lm.fit(X3,y)
print('intercept_:%.3f' % lm.intercept_)print('coef_:%.3f' % lm.coef_)
2.殘差的異方差性
    如果殘差項的方差恆定不變（即為常數），則通常認定為同方差性，反之如果方差一直在變動並未恆定則認定為有異方差性。
註：如果存在異方差性進行線性回歸，則回歸係數的假設檢驗以及置信區間的估計通常是有偏差的。
 
採用懷特檢驗法來驗證異方差性：
例如檢驗有2個自變量的線性回歸方程：
 
②將自變量和自變量之間的組合與殘差的平方建立回歸模型檢驗
如果數據滿足同方差性，則殘差項的平方無法被任何自變量變量解釋，即註：殘差項的方程的檢驗統計量的解釋力度記為nR2（即第②步中計算），其檢驗分布為卡方分布（自由度為——k(k+3)/2）
Python案例分析：依舊以weight~age+height為例
import statsmodels.api as sm model3 = sm.OLS(data['weight'],sm.add_constant(data[['age','height']]))result3 = model3.fit()sm.stats.diagnostic.het_white(result3.resid, exog = result3.model.exog)
分析：
第一個為：LM統計量值
第二個為：響應的p值，0.53遠大於顯著性水平0.05，因此接受原假設，即殘差方差是一個常數
第三個為：F統計量值，用來檢驗殘差平方項與自變量之間是否獨立，如果獨立則表明殘差方差齊性
第四個為：F統計量對應的p值，也是遠大於0.05的，因此進一步驗證了殘差方差的齊性。
3.多重共線性
    自變量之間存在相關係數=1的情況，即一個自變量可以被另一個自變量完全解釋，完全替代。
4.繪製殘差圖與異常值
殘差圖即自變量與殘差之間的散點圖，異常值即偏離正常中心值較大的奇異點。
我們以上節一元線性回歸的案例的身高與體重的回歸結果殘差圖為例：
很明顯上圖中對應x在170的時候存在異常值。
 
異常值的判斷：庫克距離(Cook’sdistance)
Python案例分析：
import statsmodels.api as sm model3 = sm.OLS(data['weight'],sm.add_constant(data[['age','height']]))result3 = model3.fit()fig, ax = plt.subplots(figsize=(19.2, 14.4))fig = sm.graphics.influence_plot(result3, ax=ax, criterion="cooks")plt.grid()
分析：由上圖可以得到，第7個數據是偏離較遠的，氣泡很大。
真實數據如下：