方差分量線性回歸模型

1、第四章方差分量線性回歸模型本章考慮的線性模型不僅有固定效應(yīng)、隨機(jī)誤差，而且有隨機(jī)效應(yīng)。我們先從隨機(jī)效應(yīng)角度理解回歸概念，導(dǎo)出方差分量模型，然后研究模型三種主要解法。最后本章介紹關(guān)于方差分量模型的兩個(gè)前沿研究成果，是作者近期在應(yīng)用數(shù)學(xué)學(xué)報(bào)與國(guó)際數(shù)學(xué)雜志 Communications in Statistics 上發(fā)表的。第一節(jié)隨機(jī)效應(yīng)與方差分量模型一、隨機(jī)效應(yīng)回歸模型前面所介紹的回歸模型不僅都是線性的，而且自變量看作是固定效應(yīng)。我們從資料對(duì)Yi , X 1i , X pi 1n 出發(fā)建立回歸模型，過去一直是把Y 看作隨機(jī)的， X1，Xp 看作非隨機(jī)的。但是實(shí)際上， 自變量也經(jīng)常是隨機(jī)的，而并不是

2、我們可以事先設(shè)計(jì)好的設(shè)計(jì)矩陣。我們把自變量也是隨機(jī)變量的回歸模型稱為隨機(jī)效應(yīng)回歸模型。究竟一個(gè)回歸模型的自變量是隨機(jī)的還是非隨機(jī)的，要視具體情況而定。 比如一般情況下消費(fèi)函數(shù)可寫為C C0 b( XT )（ 4.1.1）這里 X 是居民收入， T 是稅收， C0 是生存基本消費(fèi)， b 是待估系數(shù)。加上隨機(jī)擾動(dòng)項(xiàng)，就是一元線性回歸模型CC0b( XT )（ ）那么自變量到底是固定效應(yīng)還是隨機(jī)效應(yīng)?那要看你采樣情況。如果你是按一定收入的家庭去調(diào)查他的消費(fèi)， 那是取設(shè)計(jì)矩陣，固定效應(yīng)。 如果你是隨機(jī)抽取一些家庭，不管他收入如何都登記他的收入與消費(fèi)，那就是隨機(jī)效應(yīng)。對(duì)于隨機(jī)效應(yīng)的回歸模型，我們可以從條

3、件期望的角度推導(dǎo)出與最小二乘法則等價(jià)的回歸函數(shù)。我們希望通過X 預(yù)測(cè) Y，也就是要尋找一個(gè)函數(shù)YM ( X )M ( X 1 , X p ) ，當(dāng) X 的觀察值為x 時(shí)，這個(gè)預(yù)測(cè)的誤差平均起來應(yīng)達(dá)到最小，即EYM ( X ) 2min EYL( X ) 2（ ）L1這里 min 是對(duì)一切X 的可測(cè)函數(shù)L(X) 取極小。由于當(dāng)M(X)E(Y|X)（ 4.1.4）時(shí)，容易證明EYM (X)M(X)L(X ) 0（ 4.1.5）故當(dāng) M (X) E(Y | X)時(shí)，EY L(X)2EYM (X )2EM (X ) L(X ) 2（ 4.1.6）要使上式左邊極小，只有取 L ( X )M(X)E(Y

4、| X)。這個(gè)結(jié)果告訴我們， 預(yù)測(cè)函數(shù)取作條件期望E(Y|X)時(shí)，可使預(yù)測(cè)誤差最小。我們還可以證明，此時(shí) M(X)=E(Y|X) 與 Y 具有最大相關(guān)，即(Y,M (X )max(Y, L( X )（ 4.1.7）L這里 表示相關(guān)系數(shù)。這是因?yàn)楫?dāng)Cov (Y, M ( X )M ( X )E(Y | X ) 時(shí) ， 易 證 Cov (Y, L ( X )Cov ( M ( X ), L( X ) , 同 時(shí)Cov ( M ( X ), M ( X ) ,于是2Cov 2 (Y, L ( X ) Cov2 ( M , L( X )(Y, L(X )D(Y)D L(X )D(Y)DL(X )Cov

5、2 (M ( X ), L( X )DM(X)DM(X)D M ( X )D L(X )D(Y)DM(X)2(M(X),L(X)2(Y,M (X)2 (Y,M (X )等號(hào)當(dāng)且僅當(dāng)| (M (X),L(X) | 1（ 4.1.8）時(shí)成立，此時(shí)L(X)是 M(X)的線性函數(shù)。(4.1.3) 與(4.1.7) 表達(dá)了 M ( X )E(Y | X ) 的極好性質(zhì)，我們稱YM(X) E(Y | X)（ 4.1.9）為 Y 關(guān)于 X 的回歸曲線。上面的 L(X)可取一切函數(shù)。如果限定L(X)是 X 的線性函數(shù)，即要限定2E| Y( 01 X1m X m ) |2 min（4.1.10）L這里 min

6、是對(duì) X 的一切線性函數(shù)取極小，則稱滿足上式的線性函數(shù)為Y 關(guān)于 X 的回歸直線。L我們可以求出 0 , 1, ,m 的解。記(1 ,m )，則L(0 ,) E|Y(01 X1m X m ) |2 b2RXX2 RXYD(Y)（4.1.11）這里bE(Y)(01EX 1m EX m )RXXE(XEX)(XEX) DX 1Cov ( X1 , X 2 )Cov ( X 1 , X m )Cov ( X m , X1 )Cov ( X m , X 2 )D (X m )RXY(Cov(Y, X1 ), Cov(Y, X m )對(duì) L( 0， )求微分 (矩陣微商公式( X AX )2AX )得：

7、Xb0RXXRXY解得?0EY?E(X )?RXX1 RXY這里當(dāng)然假定RXX1 存在，否則使用廣義逆。此時(shí)的預(yù)測(cè)誤差方差是L( ?0 , ?)E| Y( ?0?1 X1?m X m ) |2 ?RXY ?2 ?RXYD(Y)21YRXY RXX RXY1(RR1R)2/XYXYXXXYY（ ）（ ）（ ）（ ）（ ）（ ）（ ）3為復(fù)相關(guān)系數(shù)。它指出了Y 與多元變量 X X1, X m 之間的線性相關(guān)程度，是一元相關(guān)系數(shù)Cov ( X ,Y )rXY（4.1.19）DXDY的推廣。從條件期望角度我們導(dǎo)出的隨機(jī)效應(yīng)回歸模型的回歸直線表達(dá)式， 與從最小二乘角度導(dǎo)出的固定效應(yīng)的回歸方程，表達(dá)式是等

8、價(jià)的，所以從計(jì)算角度，我們不怎么區(qū)分。二、方差分量模型概念上段我們建立了隨機(jī)效應(yīng)概念，將自變量也視作隨機(jī)變量，這就可以導(dǎo)出方差分量模型。方差分量模型研究工作的奠基人是我國(guó)最早的統(tǒng)計(jì)學(xué)家許寶馭馬錄先生。還是剛才提到的消費(fèi)函數(shù)回歸模型，我們作隨機(jī)抽樣。考慮居民按職業(yè)的分類，如工人、教師、醫(yī)生、律師、 店員等等， 記為 X i , i1, m ，我們從這些職業(yè)中隨機(jī)抽取了n 個(gè)樣本，則模型可寫為CijC0b( X i Ti )ij ,j1, n, i 1, ,m（4.1.20）這里 Xi 可看作是第 i 種職業(yè)對(duì)收入的效應(yīng)。如果我們事先安排好取哪個(gè)職業(yè)的，當(dāng)然Xi 是固定效應(yīng)。 可是我們現(xiàn)在對(duì)職業(yè)選

9、取是隨機(jī)的，而且我們還想研究職業(yè)效應(yīng)的方差，這就導(dǎo)入了方差分量模型，因?yàn)楝F(xiàn)在Cij 的方差由兩部分組成：Var(Cij )2b222（4.1.21）0X為了數(shù)學(xué)符號(hào)統(tǒng)一，我們將經(jīng)濟(jì)學(xué)中的符號(hào)改過來，剛才建立的模型是YijU 1 1iij ,i 1, m, j 1, , n（4.1.22）它有一項(xiàng)固定效應(yīng) ，一項(xiàng)隨機(jī)效應(yīng) 1，一項(xiàng)隨機(jī)誤差 。如果還要考慮地區(qū)因素對(duì)消費(fèi)的影響，還可以加進(jìn)第二個(gè)隨機(jī)效應(yīng) 2，于是可得模型YU 1 1U 22（4.1.23）這次我們省掉了取值的標(biāo)記，Y 的方差由三項(xiàng)組成。一般地，我們建立方差分量模型如下：Y XU 11U m m（4.1.24）這里有固定效應(yīng)向量 ，隨

10、機(jī)效應(yīng)向量( 1 ,2 ,m )（4.1.25）并且將隨機(jī)誤差項(xiàng) 也并入了隨機(jī)效應(yīng)向量去。設(shè)計(jì)矩陣X 以及4U (U 1 ,U 2 , ,U m )（4.1.26）都是已知的。對(duì)于隨機(jī)效應(yīng)i ,i1, m ,合理的假定是E(i )0, Cov ( i ,j ) 0,ij（4.1.27）D ( i )i2 , i 1, , m當(dāng)然以后有時(shí)還可以考慮 i 是向量的情況，不過這里假定每個(gè) i 是一維變量。記Vi U iU i , i1, m,12V1m2Vm，（4.1.28）則方差分量模型可記為E (Y)X ,Var (Y)（4.1.29）模型的主要任務(wù)是要估計(jì)固定效應(yīng)向量 與方差分量12 ,22

11、,m2 。和一般的多元線性回歸模型相比，就是待估的方差多了。通過這些介紹， 我們就可以方便地將各種經(jīng)濟(jì)方面的普通線性回歸模型改造成方差分量模型，當(dāng)然要根據(jù)實(shí)際。第二節(jié)方差分量模型的解法對(duì)于方差分量模型YXU 11U mmn 1n p p 1n p1 p1 1n pm p m 1E(Y) X, Var (Y)m2U iU i（ 4.2.1）ii 1一般都采用二步估計(jì)法，首先估計(jì)方差分量12 , m2 ，然后再估計(jì)固定效應(yīng) 。按照廣義最小二乘*(X ? 1X)X?1Y（ 4.2.2）其中?m?i2U iU i（ 4.2.3）i 1所以方差分量模型解法的關(guān)鍵是估計(jì)方差分量。 以下介紹的方法， 也都是

12、針對(duì)方差分量估計(jì)方法而言的。5一、方差分析法先從一個(gè)簡(jiǎn)單的模型結(jié)合數(shù)據(jù)結(jié)構(gòu)形象地說明方法?？紤]模型Yij0iij ,i1, m, j1, , n（ 4.2.4） 0 為總平均，是固定效應(yīng)， 1,， m是隨機(jī)效應(yīng)，E i0, Cov( i , j )0,i j ，Var( i )A2 ,i 1, , m。對(duì)于隨機(jī)誤差ij , Var(ij )2 。這個(gè)模型如果記作方差分量模型的標(biāo)準(zhǔn)形式是YX 0U（ 4.2.5）其中設(shè)計(jì)陣 X=(1,1, ， 1)，隨機(jī)效應(yīng)矩陣為10010k01U1（ 4.2.6）0101m km我們手中資料只有Y(Y11 ,Y1k ,Y21 ,Y2k ,Ymk )我們采用記法

13、方便一些，將資料Y 排成表ij12k組內(nèi)平均1YYYY111121k2Y21Y22Y2kY2mYm1Ym2YmkYm6方差分析主要掌握三點(diǎn)， 一是計(jì)算組內(nèi)差、 組間差， 二是作平方和分解，三是計(jì)算各自的自由度。先計(jì)算總平均：Y1mk（ 4.2.7）Yijmk i 1 j 1總變差 (全體資料與總平均的偏差平方和)：mkY )2ST(Yij（ 4.2.8）i 1j1各組平均 (各組資料橫向相加并平均)Yi1 kYij ,i1, m（ 4.2.9）k j1組間差 (各組平均數(shù)與總平均數(shù)的偏差平方和)kmSA(YiY )（4.2.10）j 1i1組內(nèi)差 (各組數(shù)據(jù)與本組平均數(shù)的偏差平方和)mkYi

14、)2S(Yij（4.2.11）i 1j1則必有平方和分解STSAS（4.2.12）將各平方和除以各自的自由度。T 有一個(gè)約束Y.(4.2.7)，自由度為 n1 mk 1 ； ASS有 m 組差， 1 個(gè)約束，自由度為m 1； Se 有 mk 組差， m 個(gè)約束，自由度為mk-m。注意有自由度分解：f Tf Af e ,mk1( m 1)(mk m)（4.2.13）于是算出均方：QT11 ST（4.2.14）mkQ A1SA（4.2.15）m1Q1S（4.2.16）mkm因?yàn)榧俣殡S機(jī)效應(yīng)，可以算出各均方的均值：7E(QA )k22（4.2.17）E(Q )2（4.2.18）以 QA 代者 E(

15、QA ) ， Q 代替 E(Q ) ，得方程組：k A22Q A（4.2.19）2Q解得?2Qe ,?2A(QAQ ) / k（）這樣就作好了方差分量的估計(jì)，然后可以按作出 的估計(jì)。因?yàn)檫@里的方差分量是由方差分析法作出的，故稱為方差分析法。推廣到一般的方差分量模型時(shí)，基本原則是類似的。我們不妨考慮方差分量模型Y XU 1 1U 22（4.2.21）Cov(Y)12U1U122U 2U 22I先對(duì)總平方和YY 作平方和分解Y YSS 1S 2 S（4.2.22）其中 S 是在模型 Y=X + 中， 的回歸平方和：SSES ()YX(X X)X Y（4.2.23）S 1是在模型 Y X U11中，

16、消去 影響后 1 的平方和S 1SES (, 1 )SES( )（4.2.24）類似地， S 2 是在模型 YXU 11U 22中消去 和 1 影響后， 2 的平方和：S 2EES (,1, 2)SES (, 1 )（4.2.25）最后的 S 為殘差平方和SY YSES (, 1 ,2 )（4.2.26）可以驗(yàn)證SY (ID)Y（4.2.27）S 1Y (DD1 )Y（4.2.28）8S 2Y (D1 D12 )Y（4.2.29）SY D12Y（4.2.30）這里D IX(XX) X IPX（4.2.31）D1DDU 1(U1DU 1) U 1DDPDU 1（4.2.32）D12D1D1U2

17、(U2D1U 2) U2D1D1PD1U 2（4.2.33）這里 P* 表示關(guān)于 * 的投影陣。下面計(jì)算各平方和的均值。E(S 1)X (D D1) Xtr (D D1 )222I U1U1 1U2U2 2X (D D1)Xtr (U 1 DU1 )tr (U 1 D1U 1 )12tr(U 2 DU 2 )2122tr(U 2 D1U 2 ) 22tr ( D D1 ) 2（4.2.34）因?yàn)?DX0, D1X0，所以上式第一項(xiàng)為 0。在第三項(xiàng)中，tr (U1 D1U 1 ) trU 1 DU 1 U 1 DU 1 (U 1DU 1 ) U 1DU 1 0（4.2.35）在第六項(xiàng)中tr (D

18、D1 )tr D 1U 1 (U 1 D1U 1 )U 1D 1 tr(U 1 D1U 1 ) U 1D1U 1 rk (U 1 D1U 1 )rk (U 1 D 1 )rk (U 1X )rk (X )rk (U 1X )rk ( X )所以最后有2(c2 c3 )22E(S 1 ) c1 12r2其中c1tr(U 1 DU 1 )（ ）（ ）（ ）9類似還可以求得于是我們得到方程組解此方程組，就可以得到固定效應(yīng)的估計(jì)。c2tr (U 2 DU 2 )（4.2.39）c3tr(U 2 D1U 2 )（4.2.40）r1rk ( X ),r1r2rk (U 1| X )（4.2.41）E(S

19、2)c22r32（4.2.42）2E(S )(nr1r2r3 )2（4.2.43）r3rk( XU 1U 2 )r1 r 2（4.2.44）S2(c2221c1 1c3 ) 2r2S2r32（4.2.45）2c2 2S(n r1r2r3 ) 212 ,22 ,2 的估計(jì)。然后進(jìn)入二步估計(jì)的第二步，就可以得到關(guān)于算例市場(chǎng)收益率與股利和換手率的關(guān)系考慮一個(gè)隨機(jī)效應(yīng)的多元線性模型YXUn 1n pp 1n mm 1U 的形式如同 (4.2.6) 。問題的實(shí)際背景是，觀測(cè)對(duì)象被分成了m 組，可能存在一個(gè)隨機(jī)效應(yīng)向量對(duì)各組資料有不同的作用。模型也可以寫作YijX ijiij , i1, m, j1, ,

20、 k數(shù)據(jù)結(jié)構(gòu)及具體數(shù)值如下表所示，m=6,k=6。這些資料采自'96 上海股票市場(chǎng)資料總匯 。我們研究目的一是看過去一年的股利收入與當(dāng)年換手率對(duì)當(dāng)年市場(chǎng)收益率有何影響，二是想知道是否存在一個(gè)潛在的尚未觀測(cè)到的隨機(jī)效應(yīng)，對(duì)行業(yè)有明顯影響。 當(dāng)然這種情況采用方差分量模型比較合適。要注意本例是兩個(gè)方差量，上一章第二節(jié)模型(3.2.10) 也是兩個(gè)待估的方差量。它們的隨機(jī)效應(yīng)作用范圍不一樣，不是一回事。10表 4.2.1 1996年股市資料類別股號(hào)股名1996 年收益率 %1995 年股利 %1996 年日換手率628新世界64.769203.12商631中百一店46.84511.81.686

21、32華聯(lián)商廈41.95811.31.81業(yè)655豫園商城16.19511.21.10類682新百公司79.9115.23.36694大連商場(chǎng)91.3885.84.26602真空電子33.112103.52651飛樂音響8.10801.95電800天津磁卡271.76353.74子839四川長(zhǎng)虹381.686604.41類850華東計(jì)算14.63813.27機(jī)870廈華電子68.5793.94.20617聯(lián)合化纖-21.8710.51.53化618鹵堿化工22.37022.63672廣東化纖11.86004.65工688上海石化179.8173.54.38類886湖北興化236.32852.34

22、.45889南京化纖-33.1222.44.64664哈醫(yī)藥191.66654.324.11671天目藥業(yè)111.135164.11醫(yī)812華北制藥152.0158.44.11藥842中西藥業(yè)13.8212.61.71類849四藥股份17.8922.51.68779四川制藥-24.744012.28608異型鋼管8.389102.24鋼665滬昌特鋼75.391.874.01674四川峨鐵35.93233.64鐵808馬鋼股份86.5280.54.52類845鋼管股份-25.17001.61894廣鋼股份51.3712.77.08604二紡機(jī)14.4102.31機(jī)605輕工機(jī)械6.12203.

23、50610中紡機(jī)0.70102.27械806昆明機(jī)床41.8521.14.22類841上柴股份66.981202.20862南通股份41.093201.3611首先我們作普通最小二乘回歸，得到? ( X X ) 1X Y ，然后計(jì)算 Yij*YijX ij? 。此時(shí)的 Yij* 已消除固定效應(yīng)影響，我們將它排成平面表，以作方差分析，計(jì)算A2 與2 。計(jì)算過程從 (4.2.7) 至(4.2.20) 。從下面計(jì)算過程可以看到， 平方和分解式是滿足的：SST SSA SSE 即 (ST=SA +S )。對(duì)于本例資料，隨機(jī)誤差 4055.6 遠(yuǎn)大于隨機(jī)效應(yīng)方差，組內(nèi)差遠(yuǎn)大于組間差，可以認(rèn)為隨機(jī)效應(yīng)不明

24、顯，即行業(yè)差別不明顯。對(duì)于選定的方差分量模型，回歸結(jié)果是Y0.13734.7053X 13.3258X 2它的標(biāo)準(zhǔn)差很小， 為 1.0084，這正是采用方差分量模型廣義最小二乘意義所在。擬合效果圖 (圖令人滿意。-方差分量模型方差分析法計(jì)算程序, 例 4.2.1.第一列為 Y, 以后各列為 X例 421.D 數(shù)據(jù)文件中 , m=6, k=6, p=2要顯示原始資料嗎?0=不顯示 ,1= 顯示（0）先作普通最小二乘, 并打印結(jié)果 :現(xiàn)在作線性回歸顯著性檢驗(yàn) ,計(jì)算 t,F,R統(tǒng)計(jì)量請(qǐng)輸入顯著性水平a, 通常取 a=0.01, 0.05, 0.10, a=?（ 0.05 ）-線性回歸分析計(jì)算結(jié)果樣

25、本總數(shù) 36自變量個(gè)數(shù)2-回歸方程 Y = b0+b1*X1+.+b2*X2Y = 5.3188 +4.5656 X1 +6.0128 X2回歸系數(shù) b0, b1, b2, ., b25.31884.56566.0128-殘差平方和 : 149400.60回歸平方和 : 131759.30誤差方差的估計(jì):4150.0170標(biāo)準(zhǔn)差 =64.4206-線性回歸顯著性檢驗(yàn)顯著性水平:.050-回歸方程整體顯著性F 檢驗(yàn) , H0:b0=b1=.=b2=0F統(tǒng)計(jì)量 :14.5517 F臨界值 F(2, 33)3.285全相關(guān)系數(shù)R :.6846-12回歸系數(shù)逐一顯著性t 檢驗(yàn) , H0:bi=0, i

26、=1,.,2t臨界值 t( 33)1.6924回歸系數(shù) b1-b 2的 t 值 :5.63181.1073-打印方差分析資料Y(I,J)-50.6210-39.0275-38.6721137.5441-56.0541 -50.4640-22.4490-8.9357-7.89368.0546 37.4223-20.2416-25.8348221.1285-21.4183 83.6328 -4.9701-18.2668-46.872275.9169132.1827-13.650151.74866.137130.6483-14.9082-34.5268-8.9422-40.1694-42.87723

27、3.974420.2007-77.2975-103.9000-8.8455-81.7529計(jì)算各種平均 :總平均 YYba:.0000各組平均 Yba:-16.216-2.34139.04534.244-18.463-36.270計(jì)算各種變差 : 總變差 SST, 組間差 SSA, 組內(nèi)差 SSESST= 149400.6000 SSA=27731.9200 SSE= 121668.7000打印方差分量的估計(jì)SIGMAA= 248.4600 SIGMAE= 4055.6240下面計(jì)算協(xié)方差陣及其逆的一塊, 并作分解計(jì)算廣義最小二乘估計(jì),模型轉(zhuǎn)換 Y=PY, X=PX, 下面打印矩陣P 的一塊

28、, P P=的逆SIG1=.1239-.0007-.0007 -.0007-.0007-.0007SIG1=.0000.1238-.0008-.0008-.0008-.0008SIG1=.0000.0000.1237-.0008-.0008-.0008SIG1=.0000.0000.0000.1237-.0008-.0008SIG1=.0000.0000.0000.0000.1236-.0009SIG1=.0000.0000.0000.0000.0000.1235下面打印的是廣義最小二乘的統(tǒng)計(jì)結(jié)果:現(xiàn)在作線性回歸顯著性檢驗(yàn), 計(jì)算 t,F,R 統(tǒng)計(jì)量請(qǐng)輸入顯著性水平a,通常取 a=0.01,

29、0.05, 0.10, a=?(0.05)-線性回歸分析計(jì)算結(jié)果樣本總數(shù) 36自變量個(gè)數(shù)2-回歸方程 Y = b0+b1*X1+.+b2*X2Y =.7309 +4.5682 X1 +5.7726 X2回歸系數(shù) b0, b1, b2, ., b2.73094.56825.7726-殘差平方和 :2290.69回歸平方和 : 2014.61誤差方差的估計(jì):63.6303標(biāo)準(zhǔn)差 = 7.976913-線性回歸顯著性檢驗(yàn)顯著性水平:.050-回歸方程整體顯著性F 檢驗(yàn) , H0:b0=b1=.=b2=0F統(tǒng)計(jì)量 :14.5113 F臨界值 F(2, 33)3.285全相關(guān)系數(shù)R :.6841-回歸系

30、數(shù)逐一顯著性t 檢驗(yàn) , H0:bi=0, i=1,.,2t臨界值 t( 33)1.6924回歸系數(shù) b1-b 2的 t 值 :5.6411 1.0683-比較殘差平方和:普通最小二乘的: 149400.6000:廣義最小二乘的:2290.6910下面打印的是利用廣義最小二乘的回歸系數(shù)去計(jì)算原始資料的回歸擬合結(jié)果:要作回歸預(yù)測(cè)嗎? 鍵入 0= 不預(yù)測(cè) , 1= 要預(yù)測(cè)(0)回歸系數(shù) :.73094.56825.7726要打印擬合數(shù)據(jù)嗎? 0= 不打印 , 1= 打印(0)計(jì)算結(jié)束。-圖500400300原始數(shù)據(jù)200擬合數(shù)據(jù)1000147036925814111122233-100計(jì)算中需要一

31、些分析與技巧。因?yàn)?2AUU2 I n需要用它及計(jì)算 * (X?1X)1X?1Y14注意 UU 是一個(gè)分塊對(duì)角陣，共有 m 個(gè)對(duì)角塊，每塊是一 k× k 方陣，元素皆為 1。這樣 的計(jì)算就容易了。恰好 X 也是分成 m 塊，每塊為 k× p 陣。于是X 1X ? 1 X( X 1 , X m ) diag ( ?1, ?m )X mX 1 ?1 X1X m ?m X m?22I m而 136× 36 的矩陣 mA 11，這就大大簡(jiǎn)化了計(jì)算。如果真的要計(jì)算的逆，一般微機(jī)是不可能的。這個(gè)程序開始時(shí)調(diào)用了普通多元線性回歸程序?qū)υ假Y料回歸。得到的回歸方程為Y5.3188

32、4.5656X16.0128X 2擬合效果也很好 (見圖 4.2.1.2)。但是，這個(gè)模型的殘差向量的標(biāo)準(zhǔn)差為64.4206，遠(yuǎn)大于方差分量模型的 1.0084。這可以對(duì)比說明方差分量模型的作用。圖500400300原始數(shù)據(jù)200擬合數(shù)據(jù)1000147036925814111122233-100二、最小范數(shù)二次無(wú)偏估計(jì)方法方差分量模型中方差分量的最小范數(shù)二次無(wú)偏估計(jì)(Minimum Norm Quadratic UnbiasedEstimator, MINQUE) 是 C.R.Rao 提出的，思路類似于他處理奇異廣義線性模型的分塊逆矩陣法，先提出估計(jì)應(yīng)滿足的性質(zhì)，根據(jù)這些性質(zhì)去求解(一般是一個(gè)

33、極小值問題 )，若能解出，當(dāng)然就有那些設(shè)定的性質(zhì)?？紤]一般的方差分量模型 (4.2.1) ，我們要估計(jì)方差分量12 , , m2 及其線性函數(shù)15c2 , c (c1 ,cm ),2( 12, m2 )（4.2.46）首先考慮 的估計(jì)應(yīng)具有的形式，因?yàn)槭枪烙?jì)方差，可考慮采用二次型YAY 的形式，即? c ?2?Y AY（4.2.47）A 為待估對(duì)稱矩陣。再來考慮Y AY應(yīng)具有的性質(zhì)：(1) 關(guān)于參數(shù) 的平移不變性。若參數(shù) 有平移：0（）則方差分量的估計(jì)應(yīng)該不變。此時(shí)原模型變?yōu)閅X0XU（4.2.49）其二次型估計(jì)變?yōu)?YX 0) A(YX 0)，應(yīng)該有(Y X 0) A(Y X 0) YAY（

34、4.2.50）這等價(jià)于2YAX 00XAX 0 0（4.2.51）即平移不變性應(yīng)滿足的充要條件是AX0(2) 估計(jì)量的無(wú)偏性。因?yàn)镋(Y AY)tr( A )X AXmi2 tr ( AU iU i )（4.2.53）i1今要求對(duì)一切2 成立：2mE (Y AY)2（4.2.54）cci ii 1則其充要條件應(yīng)為tr ( AU iU i ) ci ,i 1, , m（4.2.55）（ 3)最小范數(shù)準(zhǔn)則。經(jīng)過研究， 滿足平移不變性與無(wú)偏性， 尚不能唯一確定待估對(duì)稱矩陣 A。于是可以再加一個(gè)優(yōu)良性質(zhì)。如果隨機(jī)效應(yīng)向量i ,i=1, ， m 是已知的，則=c 2 的估計(jì)應(yīng)該為16c1 ( 1) c2

35、 ( 2 2 )cm ( mm )（4.2.56）p1p2pm這里diag( c1 I p, cmI p )（4.2.57）1pmmp1現(xiàn)在用 YAY去估計(jì) =c 2，在滿足不變性條件AX=0 時(shí)，Y AY( XU )A(XU )U AU（4.2.58）我們自然希望(4.2.56)與 (4.2.58) 之間相差很小，這只要求矩陣與 UAU 之間相差很小。我們選用矩陣范數(shù) UAU- 來度量與 UAU之間的差異，即應(yīng)選A 使U AUmin范數(shù)可選歐氏范數(shù)Btr (B B)（ ）（ ）若記 V=UU，則因是分塊對(duì)角陣及無(wú)偏性要求：U AUtrU AU U AU2U AU2 tr ( AU U AU U)2tr( AUU )tr 2tr ( AV2)mcitr( AU iU i )tr22pii1tr ( AV 2 )2n ci citr 2i1pitr ( AV 2 )2tr2tr2tr ( AV )2 - tr 2（）是已知的，則極小化范數(shù)U AU等價(jià)于極小化tr ( AV ) 2 ?？偨Y(jié)上述三項(xiàng)優(yōu)良性要求，求 =c 2 的最小范數(shù)無(wú)偏估計(jì)的問題，歸結(jié)為求下述極值問題：(L1) mintr ( AV ) 2AX0s ttr (