線性回歸模型

1、戴之瑤, 魏凌云, 王楊, 周仕君中國(guó)人民大學(xué)統(tǒng)計(jì)學(xué)院2015級(jí)碩士生教材：孟生旺，回歸模型，中國(guó)人民大學(xué)出版社，2015( 2 )主要內(nèi)容主要內(nèi)容l1.1 模型結(jié)構(gòu)和假設(shè)l1.2 解釋變量l1.3 參數(shù)估計(jì)l1.4 異方差與加權(quán)最小二乘估計(jì)l1.5 假設(shè)檢驗(yàn)l1.6 模型診斷和改進(jìn)l1.7 模型的評(píng)價(jià)與比較l1.8 應(yīng)用示例 1 1.1.1 模型結(jié)構(gòu)和假設(shè)模型結(jié)構(gòu)和假設(shè)( 3 )( 4 )1.1模型結(jié)構(gòu)和假設(shè)模型結(jié)構(gòu)和假設(shè) 假設(shè)我們感興趣的變量是 , 希望建立它與其他 個(gè)解釋變量 之間的函數(shù)關(guān)系 。最一般的函數(shù)形式可以表示為：式中 是隨機(jī)誤差。 在線性回歸模型中, 設(shè) 是一個(gè)線性函數(shù), 可得

2、線性回歸模型為： 如果對(duì)因變量和解釋變量有n次觀測(cè), 第i次觀測(cè)值記為 和 則相應(yīng)的線性回歸模型可以表示為：yk12,kx xx12( ,)kyf x xx12( ,)kf x xx01 1kkyxxiy12,iikixxx01 1Tiikkiiiiyxxx為方便起見, 線性回歸模型可以表示為矩陣形式:式中, ( 5 )yX11110121221211(1)(1) 1111,1kknnknknnnkknyxxyxxyxx yX( 6 )基本假設(shè)基本假設(shè)(1)誤差項(xiàng)的均值為零, 且與解釋變量相互獨(dú)立, 即(2)誤差項(xiàng)獨(dú)立同分布, 即每個(gè)誤差項(xiàng)之間相互獨(dú)立且每個(gè)誤差項(xiàng)的方差都相等：(3)解釋變量之

3、間線性無(wú)關(guān)(4)正態(tài)假設(shè), 即假設(shè)誤差項(xiàng)服從正態(tài)分布：在上述假設(shè)下, 可得： ,TEE0X022,0,( ),1,2,ijiiCovij Varin 2(0,)iN01 122()()()()( ),0(,)TTiiiiikkiTiiiiijijTiiE yExxVar yVarVarCov y yCovyN xxxx1 1.2.2 解釋變量解釋變量( 7 )( 8 )1.2解釋變量解釋變量1.2.1分類解釋分類解釋變量變量 如果解釋變量是分類變量, 在建模過(guò)程中需要把分類解釋變量轉(zhuǎn)化為虛擬變量。 為避免解釋變量之間出現(xiàn)完全共線性, 虛擬變量個(gè)數(shù)等于分類變量的水平數(shù)減去1。( 9 )Examp

4、le： 車型是一個(gè)分類解釋變量, 有A, B, C, D四個(gè)水平, 可以轉(zhuǎn)化為x1, x2, x3三個(gè)虛擬變量, 定義如下表所示：車型車型x1x2x3A100B010C001D000( 10 ) 假設(shè)車型是模型中唯一的解釋變量, 則線性回歸模型的擬合值表示為： 根據(jù)模型, 可以求得不同車型條件下對(duì)因變量的擬合值為： 在模型中, 車型D是基準(zhǔn)水平, 也稱參照水平。為了預(yù)測(cè)結(jié)果的穩(wěn)定性, 通常選擇觀測(cè)值較多的水平為基準(zhǔn)水平。在R中的實(shí)現(xiàn)：type = factor(c(A, B, C, 0D)model.matrix(type)01 12233xxx0102030車型 = A車型 = B車型 =

5、C車型 = D( 11 )1.2.2交互效應(yīng)交互效應(yīng) 交互效應(yīng)交互效應(yīng)是指一個(gè)解釋變量對(duì)因變量的影響與另一個(gè)解釋變量有關(guān)。 譬如, 不同性別的駕駛?cè)? 其年齡對(duì)索賠頻率的影響是不同的, 即年齡和性別之間存在交互效應(yīng)。( 12 )1.2.3變量的標(biāo)準(zhǔn)化變量的標(biāo)準(zhǔn)化 為了消除量綱的影響, 可以考慮對(duì)變量進(jìn)行標(biāo)準(zhǔn)化標(biāo)準(zhǔn)化處理, 即： 式中, ,jjjyjxxyyyxss221111()()11,11nnijijnniijijiyjiiyyxxyy xxssnnnn( 13 )Example: 表示汽車保險(xiǎn)的索賠頻率； 表示駕駛?cè)说哪挲g, 是一個(gè)連續(xù)變量； 表示性別, 是一個(gè)虛擬變量, 值為0 表示男

6、性, 1表示女性； 表示年齡和性別的交互效應(yīng) ； 則線性回歸模型的擬合值可以表示為： 在模型中, 年齡每增加一單位, 對(duì)索賠頻率擬合值的影響是一個(gè)跟性別有關(guān)的值, 即：y1x2x12x x01 1223 12xxx x1321yxx( 14 ) 基于標(biāo)準(zhǔn)化以后的數(shù)據(jù)建立的回歸模型為： 回歸系數(shù)間有下述關(guān)系： 標(biāo)準(zhǔn)化回歸系數(shù)的絕對(duì)值大小度量了解釋變量的相對(duì)重要性, 值越大, 表明該解釋變量對(duì)因變量的影響越大。 在R中, 用scale(data)實(shí)現(xiàn)標(biāo)準(zhǔn)化11iikkiyxx01,1,2, ,kyiiiiiisikyxs( 15 )1.2.4變量變換變量變換 解釋變量與因變量之間如果是非線性關(guān)系,

7、 可以考慮對(duì)解釋變量進(jìn)行變換或建立多項(xiàng)式回歸模型。 多項(xiàng)式回歸多項(xiàng)式回歸是把一個(gè)解釋變量的冪變換作為新的解釋變量引入回歸模型。 為簡(jiǎn)化表述, 不妨假設(shè)只有一個(gè)原始解釋變量, 則m次多項(xiàng)式回歸模型的基本形式如下： 所以多項(xiàng)式回歸也屬于線性回歸模型。2012mkyxxx( 16 ) 在普通多項(xiàng)式回歸中, 多項(xiàng)式的階數(shù)不同 , 參數(shù)估計(jì)結(jié)果也不同。 為了克服這種缺陷, 可以使用正交多項(xiàng)式回歸正交多項(xiàng)式回歸模型 , 即把原來(lái)的解釋變量 轉(zhuǎn)化為新的正交解釋變量。 譬如 , 三階正交多項(xiàng)式回歸模型為表示為：01 133yzz( 17 )在R中的實(shí)現(xiàn)：set.seed(10)x = 1:20y = 2 +

8、x + x2 + runif(20)*50mod1 = lm(yx)mod2 = lm(ypoly(x, 2)mod3 = lm(ypoly(x, 19)plot(y x, yaxs = i, pch = 19, ylim = c(0, 500), xlim = c(0, 21), xaxs = i, las = 1)abline(mod1)points(x, fitted(mod2), col = 2, type = l, lty = 4, pch = )points(x, fitted(mod3), col = 4, type = l, lty = 5, pch = )legend(1,

9、450, c(一元線性回歸, 二次多項(xiàng)式回歸, 19次多項(xiàng)式回歸), lty = c(1, 2, 3), col = c(1, 2, 4) ( 18 ) 如果解釋變量取值較大, 多項(xiàng)式模型中高次項(xiàng)可能會(huì)導(dǎo)致計(jì)算溢出, 從而使得對(duì)其參數(shù)的估計(jì)值出現(xiàn)下溢。解決這一問題的常用方法是對(duì)解釋變量進(jìn)行下述變換： 如果模型中包含多個(gè)自變量(譬如兩個(gè)), 則模型可表示為： 在多項(xiàng)式回歸模型中, 如果已經(jīng)包含高次項(xiàng), 則所有的低次項(xiàng)通常也要保留在模型中。( 19 )*maxminmaxmin2xxxxxx2201 12211 12221212yxxxxx x1 1.3.3 參數(shù)估計(jì)參數(shù)估計(jì)( 20 )( 21

10、)1.3參數(shù)估計(jì)參數(shù)估計(jì)1.3.1最小二乘估計(jì)最小二乘估計(jì) 回歸參數(shù)的最小二乘估計(jì)可以通過(guò)最小化殘差平方和求得： 對(duì)S關(guān)于 求偏導(dǎo), 并令其等于零, 即得 得到回歸參數(shù)最小二乘估計(jì)最小二乘估計(jì)：21()() ()nTTiiiSyxyXyX220TTS X yX X1TT X XX y( 22 )1.3.2極大似然估計(jì)極大似然估計(jì) 假設(shè)誤差項(xiàng)服從正態(tài)分布, 則因變量 的密度函數(shù)為： 線性回歸模型的對(duì)數(shù)似然函數(shù)可以表示為： 對(duì)數(shù)似然函數(shù)的最大化等價(jià)于殘差平方和的最小化, 即在正態(tài)分布假設(shè)下, 回歸參數(shù)的極大似然估計(jì)等價(jià)于最小二乘估計(jì)。2222()1(; ;)exp22Tiiiyf y x22211

11、()ln(2)22nTiiinly x( 23 ) 對(duì)上式關(guān)于 求導(dǎo), 并用回歸參數(shù)的極大似然估計(jì)值 代入上式 , 并令上式等于零, 則可以求得方差參數(shù) 的極大似然估計(jì)值為： 這個(gè)方差估計(jì)是有偏的, 在實(shí)際中很少用。22221111()nnTiiiiiynnx( 24 )1.3.3方差參數(shù)的無(wú)偏估計(jì)方差參數(shù)的無(wú)偏估計(jì) 在線性回歸模型中, 對(duì)因變量的預(yù)測(cè)值可以表示為： 線性回歸模型的殘差可以表示為： 殘差平方和表示為： 在正態(tài)性假設(shè)下, 有: 方差的無(wú)偏估計(jì)方差的無(wú)偏估計(jì)：1(),TTyXHX X XXyHy()yyIH y () ()()TTTTSSE yIHIH yyIH y22(1)SSE

12、nk22111niink( 25 )1.3.4最小二乘參數(shù)估計(jì)值的性質(zhì)最小二乘參數(shù)估計(jì)值的性質(zhì)(1)無(wú)偏性無(wú)偏性 在 的所有線性無(wú)偏估計(jì)中, 最小二乘估計(jì)的方差最小。(2)協(xié)方差矩陣協(xié)方差矩陣 式中, 是 對(duì)角線上元素, 是解釋變量 與其他解釋變量之間的復(fù)相關(guān)系數(shù)。( )E2122221( )()()(1)()TjjjnjjijiVarVarcRxxX X21. .()1njjjjjiicsecnkjjc2jR1()TX Xjx( 26 )(3)正態(tài)分布假設(shè)下的性質(zhì)正態(tài)分布假設(shè)下的性質(zhì) 如果進(jìn)一步假設(shè)誤差項(xiàng)服從正態(tài)分布, 則有： 此外, 在誤差項(xiàng)服從正態(tài)分布的假設(shè)下 最小二乘估計(jì)值與其真實(shí)值之

13、間的加權(quán)距離服從自由度為k + 1 的卡方分布。21( ,() )TN X X2(,)N yXI1 1.4.4 異方差與加權(quán)最小二乘估計(jì)異方差與加權(quán)最小二乘估計(jì)( 27 )( 28 )1 1.4.4 異方差與加權(quán)最小二乘估計(jì)異方差與加權(quán)最小二乘估計(jì)異方差的產(chǎn)生原因：1. 誤差學(xué)習(xí)模型（error-learning models）2. 數(shù)據(jù)采集技術(shù)的改進(jìn), 使得 減小3. 產(chǎn)生于異常觀測(cè)（outliers）4. 對(duì)CLRM假定的破壞, 即回歸模型的設(shè)定有誤5. 模型中一個(gè)或多個(gè)回歸元的分布偏態(tài)（skewness）6. 其他, 例如：不正確的數(shù)據(jù)變形（比率、一階差分變化等）, 不正確的函數(shù)形式（線

14、性與對(duì)數(shù)線性的變換）異方差性問題在截面數(shù)據(jù)中比在時(shí)序數(shù)據(jù)中更為常見 i2( 29 )1 1.4.4 異方差與加權(quán)最小二乘估計(jì)異方差與加權(quán)最小二乘估計(jì)在線性回歸模型中, 如果誤差項(xiàng)的方差互不相同, 即可以表示為則式中的最小二乘估計(jì)值雖然是無(wú)偏估計(jì), 但不再是最優(yōu)線性無(wú)偏估計(jì), 即在所有的線性無(wú)偏估計(jì)中, 不能保證上式的方差是最小的（有效性）。 Var(i) i22fi (XTX)1XTy( 30 )1 1.4.4. .1 1 異方差異方差條件下參數(shù)估計(jì)值的標(biāo)準(zhǔn)誤條件下參數(shù)估計(jì)值的標(biāo)準(zhǔn)誤 ( )(|) (|)VarE VarXVar EX11111(|)()|)()()|)()|)()( )()|

15、)TTTTTTTTTE VarXE VarX XX y XE VarX XXXXE VarX XXXE X XX VarX X XX211( )()()TTTVarX XX WX X X Var() 2W異方差不影響最小二乘估計(jì)值的無(wú)偏性 (|)0Var EX在異方差條件下, 假設(shè)誤差項(xiàng)的協(xié)方差矩陣可以表示： 則最小二乘估計(jì)值的方差可以表示為： ( 31 )1 1.4.4. .1 1 異方差異方差條件下參數(shù)估計(jì)值的標(biāo)準(zhǔn)誤條件下參數(shù)估計(jì)值的標(biāo)準(zhǔn)誤 21211( )() ()()nTTTiiiiVarX Xx xX Xl矩陣W往往是未知的, 不過(guò), 在大樣本條件下可以通過(guò)下式進(jìn)行估計(jì)： l其中 是

16、設(shè)計(jì)矩陣第i行的元素。l l問題在于上式可以以求出最小二乘估計(jì)值的標(biāo)準(zhǔn)誤, 但無(wú)法解決最小二乘估計(jì)的有效性問題。( 32 )1 1.4.4.2.2 加權(quán)最小二乘估計(jì)加權(quán)最小二乘估計(jì)l矩陣W的逆矩陣可以分解為： l用矩陣 可以對(duì)因變量、設(shè)計(jì)變量和誤差項(xiàng)進(jìn)行變換： W1 LTLL L L L,yy XXwyXl則可以建立下述的線性回歸模型： l上式模型滿足線性回歸所有假設(shè)：l零均值l同方差 L L( )()( )0EEE LL LL L LLL 22112( )()()() TTTTTTVarEEW( 33 )1 1.4.4.2.2 加權(quán)最小二乘估計(jì)加權(quán)最小二乘估計(jì)l所以其最小二乘估計(jì)模型的回歸參

17、數(shù)為：1111TTTTwX XX yX WXX Wyl上述估計(jì)值的協(xié)方差矩陣為：11221()TTwVarX XX WXl可以證明, 是回歸參數(shù) 的無(wú)偏估計(jì), 即：w()wEl所以, 在異方差條件下回歸參數(shù)的估計(jì)值可以表示為：111()TTX WXX Wy( 34 )1 1.4.4.2.2 加權(quán)最小二乘估計(jì)加權(quán)最小二乘估計(jì)l協(xié)方差矩陣W是 的矩陣, 難以通過(guò)n個(gè)樣本觀測(cè)值估計(jì), 所以通常設(shè)其為對(duì)角矩陣, 且常假設(shè)有 , 所以誤差項(xiàng)的方差與其樣本量n成反比, 即：nn wi1/nil 是對(duì)因變量 的觀測(cè)次數(shù)ni 2i2/niyil而且假設(shè)W為對(duì)角矩陣則意味著誤差項(xiàng)之間是相互獨(dú)立的 Cov(i,j

18、)0,ij( 35 )1 1.4.4 補(bǔ)充補(bǔ)充l異方差的檢驗(yàn)方法：1. 殘差圖2. ncvTest生成計(jì)分檢驗(yàn), 原假設(shè)為誤差方差不變, 備擇假設(shè)為誤差方差隨擬合值水平的變化而變化#library(car)3. Goldfeld-Quandt Test4. bartlett.test( 36 )1 1.4.4 補(bǔ)充補(bǔ)充l異方差的解決方法：方法一：NeweyWest()函數(shù)可以進(jìn)行異方差和自相關(guān)穩(wěn)健性NeweyWest估計(jì)（sandwich）library(sandwich)NeweyWest(fit)neweywest AIC(model, k = 2) AIC(model, k = log(

19、n) 1.7.4 交叉驗(yàn)證得分交叉驗(yàn)證得分基本原理：1、把原始數(shù)據(jù)集隨機(jī)分解成r個(gè)大小近似相等的子數(shù)據(jù)集；2、把第一個(gè)子數(shù)據(jù)及作為驗(yàn)證數(shù)據(jù)集, 其余r-1個(gè)子數(shù)據(jù)集合并后進(jìn)行模型的參數(shù)估計(jì), 并基于該模型對(duì)驗(yàn)證數(shù)據(jù)集的因變量進(jìn)行預(yù)測(cè), 計(jì)算誤差平方和；3、將第二個(gè)第r 個(gè)子數(shù)據(jù)集分別進(jìn)行步驟二；4、計(jì)算前述r個(gè)預(yù)測(cè)誤差平方和的平均值。 特例：每個(gè)觀察值作為一個(gè)子數(shù)據(jù)集 表示提出第i個(gè)觀察值后用其他所有觀察值建模對(duì)第i個(gè)觀察值的預(yù)測(cè)值。 是帽子矩陣對(duì)角線上的元素； 是基于全體數(shù)據(jù)建立的模型對(duì)第i個(gè)對(duì)象的預(yù)測(cè)值。 211niiiiCVyyniiy2111niiiiiyyCVnhiihiy在R中, 交叉驗(yàn)證得分可編寫以下函數(shù)求出：cv - function(fit, k) require(bootstrap) theta.fit - function(x, y)lsfit(x, y) theta.predict - function(fit, x)cbind(1, x)%*%fit$