安徽农业大学：《高级生物统计》课程教学资源（讲义）第六讲 方差分析的最优线性无偏估计法

第六讲方差分析的最优线性无偏估计法 BLUP OF ANOVA (本讲不讲) 1974年美国康奈大学C·R· Henderson博士及其协作者首先提出应用线性模型来评价种公牛育种值的 方法。他们将育种值作为观察值的线性函数采用稍做改良的最小二乘分析法,使其估计值的误差为最小。 并定名为最优线性无偏估计( Best liner unbiased predication简称BLUP法)。由于该法具有对估计效应值 的误差最小、精度最优及无偏的特点,并可由计算机进行大量的繁琐计算,故BLUP法在八十年代就被 些先进国家的学者应用于生产实践中,并在方差分析中得到具体的应用 第一节参数的最小二乘估计 、线性模型 假定观察的随机变量Y与若干个因素X1,X2,…,Xp间存在线性关系,其线性模型为: Y=B1x1+B2X2+…+BpXp+ (5-1) 其中β1,β2…βp为待估的参数,ε为随机误差。若对Ⅺ,X2…Xp,Y作n次观察(实验),则(5 1)式可表示为: y1=x1B1+x12B2 XIpBp 1B1+x22B2 2pB. +a2 y=xmB.+x,B, B 将(5-2)式改写成矩阵形式为: Y=BX+E 式中Y为观察值向量,β为参数向量,Ⅹ为设计矩阵或结构矩阵,ε为误差向量。通常假定ε;~N(o, 02)i=1,2,…,n且相互独立,因此,对于(5-3)式可进一步记为(Y,XB,o21),意为观察值向 量Y有 Er=XB Vary=o/ 对线性模型(5-4),主要对β或β的函数及σ2作出估计并进行假设检验 、效应的最小二乘估计 对于β的估计,根据(5-4)式有: B=X′Y 称(5-5)式为正规方程,它的解为β的最小二乘估计,记作B 1、当(X′ⅹ)为可逆(满秩)矩阵时

55 第六讲 方差分析的最优线性无偏估计法 BLUP OF ANOVA （本讲不讲） 1974 年美国康奈大学 C·R·Henderson 博士及其协作者首先提出应用线性模型来评价种公牛育种值的 方法。他们将育种值作为观察值的线性函数采用稍做改良的最小二乘分析法，使其估计值的误差为最小。 并定名为最优线性无偏估计（Best Liner Unbiased Predication 简称 BLUP 法）。由于该法具有对估计效应值 的误差最小、精度最优及无偏的特点，并可由计算机进行大量的繁琐计算，故 BLUP 法在八十年代就被一 些先进国家的学者应用于生产实践中，并在方差分析中得到具体的应用。 第一节 参数的最小二乘估计 一、线性模型 假定观察的随机变量 Y 与若干个因素 X1，X2，…，XP间存在线性关系，其线性模型为： Y=β1X1+β2X2+…+βpXp+ε （5—1） 其中β1，β2…βp 为待估的参数，ε为随机误差。若对 X1，X2…Xp，Y 作 n 次观察（实验），则（5 —1）式可表示为：        = + + + + = + + + + = + + + + n n n np p n p p p p y x x x y x x x y x x x                 1 1 2 2 2 21 1 22 2 2 2 1 11 1 12 2 1 1 （5—2） 将（5—2）式改写成矩阵形式为： Y=βX+ε （5—3） 式中 Y 为观察值向量，β为参数向量，X 为设计矩阵或结构矩阵，ε为误差向量。通常假定εi～N（o， σ2） i=1，2，…，n 且相互独立，因此，对于（5—3）式可进一步记为（Y，Xβ，σ2 I），意为观察值向 量 Y 有    = = VarY I EY X 2   (5—4) 对线性模型（5—4），主要对β或β的函数及σ2 作出估计并进行假设检验。 二、效应的最小二乘估计 对于β的估计，根据（5—4）式有： (X＇X)β=X＇Y （5—5） 称（5—5））式为正规方程，它的解为β的最小二乘估计，记作  ˆ 。 1、当（X＇X）为可逆（满秩）矩阵时

B=( 若要估计β1,β2,…βp的线性函数 P=C1B1+C2B2+…+CpBp (5-7) 则P=c'B 称为P的最小二乘估计。 对于02的估计记 R2=(Y-XB)(-B)=(-BX)Y-B)由(5-6)式有 Ro=YY-YXB=YY-BXr (5-8) 称R2为残差平方和,记r为X的秩,则 R6/ 为σ2的无偏估计 最小二乘估计的性质有: E(B)=EICXX-XY=(XX)XEY=(XX)XXB=B (5-10) Var(B)=Var(X X)XY=(XX)Xo IX(XX) (5-11) (r-X-YY(XX-1=0(Xr) 类似地 E(CB)=CB Var(CP=oIC(XX)] (5-13) 即:B、C'B分别为β、Cβ的唯一最小方差线性无偏估计,简称BLUE( Best linear unbiased estimator 2、当(X′X)不可逆时,(5-5)式无唯一解,B有无穷多解,此时常对β加上一个约束条件 HB=O (5-14) 使B的各分量不独立,从(5-14)可解出B与其中独立的那部分分量间的关系: B=T B(O) (5-15) βo)是B中在约束(5-14)下,独立的那些分量构成的向量(维数<P),于是原模型(5-4)变成 「EY=(XTB0 (5-16) vAry=oI 对于β0而言,又成为无约束问题了,所以上述(5-5)至(5-13)仍可适用,只须将其中X改为 XT,B改为B0),于是有正规方程 (XT)′(XT)B0)=(XT)′Y (5-17) 解出Bo,并由(5-15)式得到B,这时B、CB是在HB=0约束条件下的唯一的BLUE,记Bm0。 、对于参数的假设检验,即检验:H1B=0 Ro/ox(n-r (5-18) 其中r=rank(X),即X阵的秩 2、在约束H1B=0下 RH=(Y-XBH(Y-XBH) (5-19) 可以证明 R.a2~x2(d1) (5-20)

56  ˆ =(X＇X)—1X＇Y （5—6） 若要估计β1，β2，…βp 的线性函数 P=C1β1+C2β2+…+Cpβp （5—7） 则 ˆ  ˆ P = c  C＇=（C1，C2，…Cp） 称为 P 的最小二乘估计。 对于σ2 的估计记： ) ˆ )( ˆ ) ( ˆ ) ( ˆ ( 2 R0 = Y − X  Y − X = Y − X  Y − X 由（5—6）式有 R = YY −YX = YY − X Y 2 ˆ ˆ 0 （5—8） 称 2 R0 为残差平方和，记 r 为 X 的秩，则 ˆ = R / n − r 2 0 2  （5—9） 为σ2 的无偏估计。 最小二乘估计的性质有：  =   =   =    =  − − − E E X X X Y X X X EY X X X X 1 1 1 ) [( ) ] ( ) ( ) ˆ ( （5—10） 1 1 2 1 2 1 1 1 2 1 ( ) ( ) ( ) ) [( ) ] ( ) ( ) ˆ ( − − − − − − − =   =  =   =    X X X X X X I X X Var Var X X X Y X X X IX X X     （5—11） 类似地： E(C ˆ ) = C （5—12） ) [ ( ) ] ˆ ( 2 1 Var C C X X C −  =    （5—13） 即：  ˆ 、C＇ ˆ 分别为β、C＇β的唯一最小方差线性无偏估计，简称 BLUE（Best Linear Unbiased Estimator） 2、当（X＇X）不可逆时，（5—5）式无唯一解，  ˆ 有无穷多解，此时常对β加上一个约束条件 Hβ=0 （5—14） 使β的各分量不独立，从（5—14）可解出β与其中独立的那部分分量间的关系： β=Tβ(0) （5—15） β(0)是β中在约束（5—14）下，独立的那些分量构成的向量(维数<P)，于是原模型（5—4）变成    = = VarY I EY XT 2 (0) ( )   （5—16） 对于β(0)而言，又成为无约束问题了，所以上述（5—5）至（5—13）仍可适用，只须将其中 X 改为 XT，β改为β(0)，于是有正规方程 （XT）＇（XT）β(0) =（XT）＇Y （5—17） 解出β(0)，并由（5—15）式得到  ˆ ，这时  ˆ 、C  ˆ  是在 Hβ=0 约束条件下的唯一的 BLUE，记 0 ˆ  H 。 三、对于参数的假设检验，即检验：H1β=0 1、 2 2 RO / ～ ( ) 2  n − r （5—18） 其中 r=rank(X)，即 X 阵的秩。 2、在约束 H1β=0 下 ) ˆ ) ( ˆ ( 1 1 1 2 RH Y X H Y − X H = −  （5—19） 可以证明 2 2 / 1 RH  ～ ( ) 1 2  d （5—20）

d1=n-[r(x)-r(H1 (R2-R)a2~x2(d2) d,=r(X)-[()-r(HDI 3、在假设H1B=0成立时(5-18)与(5-20)式独立 (R.-R)/d2 (d2,n-r) (5-21) R2/ 在方差分析中虽然常有约束,如各效应之和为0,即H=(1……1),这是个理想的约束。它的特点是秩 ramk()=P=rmk(X)+mmk(H)即满秩,如 21 X 201 H 101 这时有R=R2,仍可用(5-21)检验H1B=0是否成立。 第二节单因素方差分析的最小二乘估计 、数学模型 yj=μ+at+eii=1,2,…,Pj=1,2,…,n 其中y为A因素第i个水平的第j次观察值,a为A的第i个水平的效应,e是随机误差,服从N(O, 2,u为总体平均数。(1A式可写成 (2A) 其中 Y"=(y1…y1,nl…yplr…yp.np),B′=( a 2 e′=(e1…el,nl…"ep…"ep,mp) 设计矩阵X 00 100 011 、效应的最小二乘估计 正规方程组(XX)B=X7 nn (3A) 其中N= . Y y 常把(3A)写成表格形式

57 [ ( ) ( )] 1 1 H1 d n r r X = − H − 2 2 0 2 ( )/ 1 RH − R  ～ ( ) 2 2  d ( ) [ ( ) ( )] 2 1 H1 d r X r r X = − H − 3、在假设 H1β=0 成立时（5—18）与（5—20）式独立 R n r R R d F H − − = / ( )/ 2 0 2 2 0 2 1 ～ ( , ) 2 F d n − r （5—21） 在方差分析中虽然常有约束，如各效应之和为 0，即 H=（1……1），这是个理想的约束。它的特点是秩 rank( ) P rank(X) rank(H) X H = = + 即满秩，如 3 1 0 1 2 0 1 1 2 1 (1 0 1) 2 0 1 1 2 1 3 3 =           = =         =         =  P rank P H X X H rank 这时有 2 0 2 RH = R ，仍可用（5—21）检验 H1β=0 是否成立。 第二节 单因素方差分析的最小二乘估计 一、数学模型 yij=μ+αi+eij i=1，2，…，P j=1，2，…，n； （1A） 其中 yij 为 A 因素第 i 个水平的第 j 次观察值，αi 为 A 的第 i 个水平的效应，eij 是随机误差，服从 N（0， 2  e ），μ为总体平均数。(1A)式可写成 Y=Xβ+ε （2A） 其中 Y＇=（y11…y1, n1…yp1…yp, np）， β＇=（μ，α1，α2…αp） ε＇=（e11…e1, n1…ep1…ep, np） 设计矩阵 p N X +                       = ( 1) 0 0 0 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1 0 1 1 1 1                                  二、效应的最小二乘估计 正规方程组 X X  = X Y ˆ ( )                 =                                   . 2 1 2 1 2 2 1 1 1 2 .. 0 0 0 0 0 0 P P P P Y Y Y Y n np n n n n N n n n               （3A） 其中 i p i N n =1 =  ij n j i Y y i =1  =  （i=1，2，…，P） ij n j p i i p i Y y y i 1 =1 =1  =  =  =   常把（3A）写成表格形式 （4A）

u a a2 a RHM H: N Y ar1:nnt0… 00 (RHM: Right hand member意为右边数) (3A)式的系数矩阵(X′X)是不可逆的,其秩为P。这时须对a1,a2,…ap加上一个约束条件 ∑a1=0(i=1,2,…,P) (5A) (6A) 于是参数B的估计可用以下三种方法 第一种逆矩阵法 该法的最大优点在于求解求逆过程中的有关数值可直接用于假设检验。 为使系数矩阵可逆,利用ap=-(a1+ax+…+ap-1)可化(3A)或(4A)式为同解方程式 如(4A)式中第一个方程和最后一个方程可作变化 Nu+n a 1+n2 a 2++np-la p-1+np(-ai-a 即:Np+(n-n)a计+(n-n)a2+…+(n-1-np)ap-1=Y 而 即;np以 从第二个方程起均减去最后一个方程(7A)式,得 (8A) 以表格示之的正规方程组 &I n -np n tnp np n2 tnp np (8A)的系数矩阵S为可逆矩阵,记S1,则(8A)式有唯一解。 P-12 则B=S-(X)为约束(5A)下的BLUE。实际上不管用何种方法,解得B都是B的BLUE。 第二种简捷法 逆矩阵法须将(4A)化为(8A),其过程较为繁琐。故可直接从(4A)解得β。从该式中的a;方程知

58  ˆ 1  ˆ 2  ˆ …………  P ˆ RHM μ： N n1 n2 ……………np Y．． 1： n1 n1 0 ……………0 Y1．  2： n2 0 n2 ……………0 Y2． ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶  p ： np 0 0 ……………np Yp． （RHM：Right Hand Member 意为右边数） （3A）式的系数矩阵（X＇X）是不可逆的，其秩为 P。这时须对α1，α2，…αp 加上一个约束条件 0 ( 1,2, , ) 1 i i p p i  = =  =  （5A） 即 i p i ap  1 1 −  =  （6A） 于是参数β的估计可用以下三种方法。 第一种 逆矩阵法 该法的最大优点在于求解求逆过程中的有关数值可直接用于假设检验。 为使系数矩阵可逆，利用αp=－（α1+α2+…+αp－1）可化（3A）或（4A）式为同解方程式。 如（4A）式中第一个方程和最后一个方程可作变化 Nμ+n1α1+n2α2+…+np－1αp－1+np（－α1－α2―…―αp－1）=Y．． 即：Nμ+（n1－np）α1+（n2－np）α2+…+（np―1―np）αp－1=Y．． 而 npμ+0α1+0α2+…+np（－α1－α2－…－αp－1）=Yp． 即：npμ－npα1－npα2―…―npαp－1=Yp． （7A） 从第二个方程起均减去最后一个方程（7A）式，得 （8A） 以表格示之的正规方程组  ˆ 1  ˆ 2  ˆ … … 1 ˆ  p− RHM 1 1 1 1. . 2 2 2 2. . 1 1 1 1. . 1 2 1 : : : : .. p p p p p p p p p p p p p p p p p p p p p p p p p n n n n n n n Y Y n n n n n n n Y Y n n n n n n Y Y N n n n n n n Y − − − − + − − + − − − − − − − − −                   （8A）的系数矩阵 S 为可逆矩阵，记 S ―1，则（8A）式有唯一解。               = − − − − − − − − 1 0 1 1 1 2 1 1 10 11 12 1 1 00 01 02 0 1 1 p p p p p p p C C C C C C C C C C C C S         （9A） 则 ( ) ˆ 1 = s X Y −  为约束（5A）下的 BLUE。实际上不管用何种方法，解得  ˆ 都是β的 BLUE。 第二种 简捷法 逆矩阵法须将（4A）化为（8A），其过程较为繁琐。故可直接从（4A）解得  ˆ 。从该式中的αi 方程知

a=y (10A) 其中y=Y/n,为平均数,因∑a=0,故相加得 =y+2 点(P-1)-Y2-…-n 回代(10A) -+(p-)2-…-F (12A) an=[-F-F2…+(P-1) 将A,a1,…,an-1式中的各平均数前的系数用矩阵K表示 P-1-1-1 (13A) 则由该简捷法相配合的S-1为 S-l=KDK′ (14A) 其中D为 第三种参数变换法 对(4A)方程组作如下变换 a (16A) 则方程(4A)变为 比原式少最后一列和一行的 (17A) u aa HM

59        + = + = + =    P YP Y Y       ˆ ˆ ˆ ˆ ˆ ˆ 2 2 1 1  （10A） 其中 Yi Yi ni /  =  为平均数，因∑αi=0，故相加得 ˆ 1 ( ) = Y1 + Y2 + + YP p   （11A） 回代（10A）        = − − − + − = − + − − − = − − − −          ˆ [ ( 1) ] ˆ [ ( 1) ] ˆ [( 1) ] 1 2 1 1 2 1 2 1 2 1 1 p p p p p P p Y Y p Y Y p Y Y p Y Y Y        （12A） 将 1 1 ˆ , ˆ , , ˆ     p− 式中的各平均数前的系数用矩阵 K 表示 P P P P P P K                  − − − − − − − − − − − − − = 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1              （13A） 则由该简捷法相配合的 S －1 为 S －1=KD－1K＇ （14A） 其中 D 为 p p np n n D                  =   2 1 （15A） 第三种 参数变换法 对（4A）方程组作如下变换             =  = −  = −  = −  = + − − ˆ 0 ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ 1 1 2 2 1 1 p p p p p p p               （16A） 则方程（4A）变为： 比原式少最后一列和一行的 （17A）  ˆ 1  ˆ 2  ˆ ……… 1 ˆ −   p RHM

(17A)也有唯一解,由(16A)及∑a=0得 ∑a (18A) ap=-yp2a 于是由(16A)可求得a1,a2,…,am1和。 对a2的估计由(5-8)式可知其残差平方和为 Ro=YY-BXY=2yi-(ar+a pyp.) ∑y2-i+a(X-Yn)+…+an(ym-n) (19A) ∑y2-R(u,a,) 其中R(ya1)为A因素对总平方和的贡献(即处理间平方和),由(5-9)式得a2的无偏估计 =1∑2-R以月 (20A) 三、假设检验 根据(5-11),在约束an=-(a1+…+an-)下 Vari a 0 S (21A) 即:Wa(G)=Coa2 Var(a)=c,o (22A) 因为 ap=-(0,1…31)(A.Gn,…,am=) 所以 Farn=(-1)(01…,)an(.a…p=)0…y =(01…,1G2ls-(0,…1) =o[c1+c12+…+c1p-1)+(c21+c22+…+C2p-1)+ 其中记cpp=[(c1+C12+…+c1p)+…+(cp11+cp-12+…+Cp1p=1)(23A) cp)(=1,2,…,P-1) (24A)

60 μ： N n1 n2 ………np-1 Y．． α1： n1 n1 0 ………0 Y1． α2： n2 0 n2 ………0 Y2． ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶ ∶ αp-1： np-1 0 0 ………np-1 Yp-1． （17A）也有唯一解，由（16A）及∑αi=0 得 i p p p i i p i  ˆ  ˆ p ˆ p ˆ 1 1 1   =  − = − = − = （18A） i p i  p = − p   − = ˆ ˆ 1 1 1 于是由（16A）可求得 1 2 1 ˆ , ˆ , , ˆ     p− 和μ。 对 2  e 的估计由（5—8）式可知其残差平方和为： ( , ) [ ˆ ˆ ( ) ˆ ( )] ( ˆ ˆ ˆ ) ˆ 2 1 1 1 1 2 1 1 2 2 0 i j i i j i j p p p p i j i j P P i j y R y Y Y Y Y Y R Y Y X Y y Y y y          =  − =  − + − + + − =  −   =  − + + +    − −        （19A） 其中 ( , ) R   i 为 A 因素对总平方和的贡献（即处理间平方和），由（5—9）式得 2  e 的无偏估计 [ ( , )] 2 1 1 2 1 ij i n j p i e n p y R i  =   −   = = − (20A） 三、假设检验 根据（5—11），在约束 ˆ ( ˆ ˆ )  p = − 1 ++ p−1 下 2 1 1 1 ( ˆ ˆ ˆ ) − Var −  = s     p  e （21A） 即： 2 00 ( ˆ) Var  = C  e 2 ( ˆ ) Var  i = Cii e (i = 1,2,  , p −1) （22A） 因为 ˆ (0,1, ,1)( ˆ , ˆ , , ˆ ) 1 1 = −   p      p− 所以 2 11 1 2 1 1 11 1 2 1 1 2 1 2 2 2 1 2 2 1 1 1 2 ( )] [( ) ( ) (0,1, ,1) (0,1, ,1) ˆ ( 1) (0,1, ,1)[ ( ˆ, ˆ , ˆ ) ](0,1, ,1) p p e p p p p e p p e p p c c c c c c c c c c s Var Var        = + + + + = + + + + + + + + =  = −   − − − − − − − −          其中记 [( ) ( )] p p = 11 + 1 2 + + 1 p + + p−1 1 + p−1 2 + + p−1 p−1 c c c  c  c c  c （23A） ( ) ( 1,2, , 1) ci p = − ci 1 + ci 2 ++ ci p−1 i =  p − （24A）

则a+a)=(c0+cn,+20,)2=O(=12…,P)(25A) 由(21A)式 par(a±a)=(cn1+c±2c,)2(,j=12,…,p) σ2[(0,0,…,1…;1…,0)s-(0,0,…1…1…0) A因素的平方和SA可以用SA=Sr-S,计算时可用 (27A) 因此,为检验H0:a1=a2=…=aP可用 F S/p-1 (28A) a,间的多重比较,用SSR法,即检验a=,可用 2 检验是否≥σSSR a(k, df), 若是,则差异显著,若不是,则差异不显著 四、实例分析 抽测不同品种3头母猪的窝产仔数,其资料如表5—1。试检验3头母猪平均窝产仔数的差异是否显著。 表5—1不同品种的3头母猪产仔数(头窝) 1311131012141112 1079877 48 6 8 数学模型: =+a:+e j=1,2,…,n,i=1,23 其中y为第i个品种的第j窝产仔数,a为第i个品种的效应,c;~N(0,a2) N=∑n1=18X.=∑y;=198,依(4A)最小二乘正规方程如表(5-2) 表5-2最小二乘正规方程 a. aa RHM

61 则 ( ˆ ˆ ) ( 2 ) ( 1,2, , ) 2 2 0 0 0 i p n Var c c c i e + i = + i i + i e = =      （25A） 由（21A）式 [(0,0, ,1, ,1, ,0) (0,0, ,1, ,1, ,0) ] [(0, ,1, ,1, ,0) ( ˆ, , ˆ , , ˆ , ,0) ] ( ˆ ˆ ) ( 2 ) ( , 1,2, , ) 2 1 2 =  =   = +  =    −           s Var Var c c c i j p e i j i j i j i i j j i j e        （26A） A 因素的平方和 SA 可以用 SA=ST－Se，计算时可用 ( ˆ ˆ ) ( ˆ ˆ ) 1 1 1 11 1 1 11 1 1 1 1            = − − − − − − − p p p p p A p c c c c S         （27A） 因此，为检验 H0：α1=α2=…=αP 可用 1 1 2 0 − − = R n S p F A ～ F(p－1，n－p) （28A）  i ˆ 间的多重比较，用 SSR 法，即检验  i  j ˆ = ˆ ，可用 i i j j i j i j c c 2c 2 ( ˆ ˆ ) + −  − （29A） 检验是否≥ ( , ) ˆ dfe e k  SSR ，若是，则差异显著，若不是，则差异不显著。 四、实例分析 抽测不同品种 3 头母猪的窝产仔数，其资料如表 5—1。试检验 3 头母猪平均窝产仔数的差异是否显著。 表 5—1 不同品种的 3 头母猪产仔数（头/窝） 品种 窝产仔数 Yi． ni Yi． A1 A2 A3 12 14 13 15 13 11 13 10 12 14 11 12 10 7 9 8 7 7 54 96 48 4 8 6 13.5 12 8 数学模型： yi j =  + i + ei j j = 1,2,  ,n; i = 1,2,3 其中 i j y 为第 i 个品种的第 j 窝产仔数，αi 为第 i 个品种的效应，ei j～N（ 2 0, e ） N = ni = 18 Y  =  yi j = 198 ，依（4A）最小二乘正规方程如表（5—2） 表 5—2 最小二乘正规方程  ˆ 1  ˆ 2  ˆ 3  ˆ RHM

18486198 a1440054 808096 600648 因第一行(列)为其余3行(列)的和,故需加a1+a,+a3=0后,方有唯一解。 解一,逆矩阵法 依(8A)有 Aa, a2 RHM 2198 206 6 2 1448 于是可用求解求逆紧凑法解之。即在系数矩阵S上加一列向量(XY),进行消元变换。 21066 261448 0.055556-0.1111110.1111ll1l 0.11ll119777778622222228 0.111111622222213.77777826 0.0568190.0113640.18181811.31818 A(2)=00113640.10227306363642863636 0.1818180.63636498181808.181818 0.0601850023148-001851911.1667 A)=0231480143519-00648152333 0.018519-0.0648150.1018520.8333 A0)中前3行、3列为S,第4列为(,a1a2 依(23A)、(24A)有 c13=-(c11+c12)=-(0.143519-0064815)=-0.078704 c23=-(c21+c22)=-(-0.04815+0.101852)=-0.037037 c33=-(c13+c23)=0.078704+0.037037=0.15741 u、a1、a2、a3的BLUE为 =11167,a1=2.333,a2=0.833,a3=-(a1+a2)=-3.1666 解二、简捷法 P=3,依(13A)、(14A)(15A)有 D=080 006

62 6 0 0 6 48 8 0 8 0 96 4 4 0 0 54 : 18 4 8 6 198 3 2 1     因第一行（列）为其余 3 行（列）的和，故需加  ˆ 1 + ˆ 2 + ˆ 3 = 0 后，方有唯一解。 解一 ，逆矩阵法： 依（8A）有  ˆ 1  ˆ 2  ˆ RHM 2 6 14 48 2 10 6 6 18 2 2 198 2 1    − − 于是可用求解求逆紧凑法解之。即在系数矩阵 S 上加一列向量（XY），进行消元变换。           − − − − =           − − =           − − =           − − = 0.018519 0.064815 0.101852 0.8333 0.023148 0.143519 0.064815 2.3333 0.060185 0.023148 0.018519 11.1667 0.181818 0.636364 9.818180 8.181818 0.011364 0.102273 0.636364 2.863636 0.056819 0.011364 0.181818 11.318182 0.111111 6.222222 13.777778 26 0.111111 9.777778 6.222222 28 0.055556 0.111111 0.111111 11 2 6 14 48 2 10 6 6 18 2 2 198 (3) (2) (1) (0) A A A A (3) A 中前 3 行、3 列为 −1 S ，第 4 列为 ( ˆ , ˆ , ˆ ) 1 2    。 依（23A）、（24A）有 ( ) 0.078704 0.037037 0.115741 ( ) ( 0.064815 0.101852 ) 0.037037 ( ) (0.143519 0.064815 ) 0.078704 33 13 23 23 21 22 13 11 12 = − + = + = = − + = − − + = − = − + = − − = − c c c c c c c c c μ、α1、α2、α3 的 BLUE 为  ˆ = 11.1667， ˆ 1 = 2.3333，  ˆ 2 = 0.8333，  ˆ 3 = −( ˆ 1 + ˆ 2 ) = −3.1666 解二、简捷法： P = 3 ，依（13A）、（14A）、（15A）有           = 0 0 6 0 8 0 4 0 0 D           − − = − − 1 2 1 2 1 1 1 1 1 3 1 K

0012-1 S-I=KD-K 180 0 0.5416 0.208333-0.16666 0.20833312916 0.16666 0.583 0.916 0.0601850.023148-0018519 0.0231480.143519-0064815 0.018519-0.0648150.101852 再依(11A)、(12A)有 (,a1,…ap)=K(1…-1) 111)(13.5)(11.16 2.333 0.8333 即:=11167,a1=2.333,3a2=0.8333,a3=-(a1+a2)=-3.1666 解三、参数代换法 由(17A)用求解求逆紧凑法作如下变换 1848198 A0=44054 0.0555560.222222044444411 A=-0.22223121217710 0.444444-177777644444488 0.071429-00714290.57142810285715 0.0714290.321429-0.5714283.214286 0.5714280.571428342857213.714287 0.166667-0.166667-0.1666678 0.16666704166670.1666675.5 0.1666670.1666670.2916674 实际上,根据a1+a2+a3=0,可在表5-2的方程组中,加一行、一列,亦能使系数矩阵可逆, 于是可在该增广矩阵上进行消元变换,其求解求逆如下 184860198 4400154 A0=8080196 6006148

63           − − − −                       − − = = − − − − 1 1 1 1 1 2 1 2 1 6 0 0 1 0 8 0 1 0 0 4 1 1 2 1 2 1 1 1 1 1 9 1 1 1 1 S KD K           − − − − = 0.16666 0.583 0.916 0.208333 1.2916 0.583 0.5416 0.208333 0.16666 9 1           − − − − = 0.018519 0.064815 0.101852 0.023148 0.143519 0.064815 0.060185 0.023148 0.018519 再依（11A）、（12A）有 ( ˆ , ˆ , , ˆ ) ( ) 1 1 1  = −      p− K Y  Yp ，           =                     − − = − −           0.8333 2.3333 11.1667 8 12 13.5 1 2 1 2 1 1 1 1 1 3 1 ˆ ˆ ˆ 2 1    即：  ˆ = 11.1667， ˆ 1 = 2.3333， ˆ 2 = 0.8333， ˆ 3 = −( ˆ 1 + ˆ 2 ) = −3.1666 解三、参数代换法 由（17A）用求解求逆紧凑法作如下变换:           = 8 0 8 96 4 4 0 54 18 4 8 198 (0) A           − − = − − 0.444444 1.777776 4.444448 8 0.222222 3.111112 1.777776 10 0.055556 0.222222 0.444444 11 (1) A           − − − − = 0.571428 0.571428 3.428572 13.714287 0.071429 0.321429 0.571428 3.214286 0.071429 0.071429 0.571428 10.285715 (2) A           − − − − = 0.166667 0.166667 0.291667 4 0.166667 0.416667 0.166667 5.5 0.166667 0.166667 0.166667 8 (3) A 实际上，根据  ˆ 1 + ˆ 2 + ˆ 3 = 0 ，可在表 5—2 的方程组中，加一行、一列，亦能使系数矩阵可逆， 于是可在该增广矩阵上进行消元变换，其求解求逆如下： ⎯⎯⎯ ⎯→                 = 第一主元素 0 1 1 1 0 0 6 0 0 6 1 48 8 0 8 0 1 96 4 4 0 0 1 54 18 4 8 6 0 198 (0) A

0.0555560.22222204444440.333333011 0.2222223.1111121.777776-1.333333110 0.444444-17777764444445-266666618 -0.333333-1.333333-2.66666640001-18 0 00 0.071429-0.0714290.5714290.4285710.07142910.285715 0.0714290.321429-0.571429-0.4285710.3214293.214286 A2)=-05714290571429342872-342857157491371429第三主元素, 0.4285710.42857134285723.4285721.42857213.71429 0.071429-0.3214291.5714291428572-0.321429-3.214286 0.166666-0.166666-0.1666661-0.3333338 0.16666604166660.166666-10.5833335.5 A3)=-0.1601602916-10458334第五丰元素 03.00000 0333333-0.583333-04583333.000-1.041667-9.5 0 02-0.020.040.3211.04 0.02009-0.090.680.560.18 02-0090009032044-0.18轴四主元素 0.04-0.680.328642.88 7.36 0.320.560. 2.88-0.969.12 0.0601850023148-0.018519-0.004630-0.3333331166667 0.0231480.143519-0.064815-0.0787040.33332.33333 0.018519-0.0648150.10185200370370.3333330.833333 0.004630-0.078704-0.0370370.1157410.333333-3.166667 0.3333330.333330.3333330.333333 0 0 解得:=11.1667,a1=2333,a2=0.8333,a3=-3.1667和上述完全相同 由A)中知=8,a1=5.5,a2=4,依(18A)、(16A)有 a=-5=-155+4)=-3167 a=-a=8-(-3.1667)=1167 a1=a1=an=5.5-3.1667=2333 a2=a2=an=4-3.1667=0.8333 由此可见,三种解法结果完全一致 对G2的估计G2:由(19A)得残差平方和 B=∑∑y-[y.+a(1-11)+a2(y2-Y) =2290-11.1667×198-2.3333×6-0.8333×48 2290-2265.0048=25 依(20A)

64 ⎯⎯⎯ ⎯→                 − − − − − − − − − − = （ ） 第二主元素 0 1 1 1 0 0 0.333333 1.333333 2.666666 4.000 1 18 0.444444 1.777776 4.444445 2.666666 1 8 0.222222 3.111112 1.777776 1.333333 1 10 0.055556 0.222222 0.444444 0.333333 0 11 1 A ⎯⎯⎯⎯⎯→                 − − − − − − − − − − − − − = （ ） 第三主元素 0.071429 0.321429 1.571429 1.428572 0.321429 3.214286 0.428571 0.428571 3.428572 3.428572 1.428572 13.71429 0.571429 0.571429 3.428572 3.428572 1.571429 13.71429 0.071429 0.321429 0.571429 0.428571 0.321429 3.214286 0.071429 0.071429 0.571429 0.428571 0.071429 10.285715 2 A ⎯⎯⎯⎯⎯→                 − − − − − − − − − − − − = （ ） 第五主元素 0.333333 0.583333 0.458333 3.000 1.041667 9.5 1 1 1 0 3.0000 0 0.166666 0.166666 0.291666 1 0.458333 4 0.166666 0.416666 0.166666 1 0.583333 5.5 0.166666 0.166666 0.166666 1 0.333333 8 3 A ⎯⎯⎯ ⎯→                 − − − − − − − − − − − − − = （ ） 第四主元素 0.32 0.56 0.44 2.88 0.96 9.12 0.04 0.68 0.32 8.64 2.88 27.36 0.02 0.09 0.09 0.32 0.44 0.18 0.02 0.09 0.09 0.68 0.56 0.18 0.06 0.02 0.02 0.04 0.32 11.04 5 A                 − − − − − − − − − − − − − = 0.333333 0.333333 0.333333 0.333333 0 0 0.004630 0.078704 0.037037 0.115741 0.333333 3.166667 0.018519 0.064815 0.101852 0.037037 0.333333 0.833333 0.023148 0.143519 0.064815 0.078704 0.333333 2.333333 0.060185 0.023148 0.018519 0.004630 0.333333 11.166667 （4） A 解得：  ˆ =11.1667,  ˆ 1 = 2.3333,  ˆ 2 = 0.8333,  ˆ 3 = −3.1667 和上述完全相同。 由 (3) A 中知 ˆ 8, ˆ 5.5, ˆ 4,  = 1  = 2  = 依 （18A）、（16A）有 (5.5 4) 3.1667 3 ˆ 1 ˆ 1 2 1 3 = −   = − + = − = i p i    ˆ =  ˆ − ˆ p = 8 − (−3.1667) = 11.1667  ˆ 1 =  ˆ 1  =  ˆ p = 5.5 − 3.1667 = 2.3333  ˆ 2 =  ˆ 2  =  ˆ p = 4 − 3.1667 = 0.8333 由此可见，三种解法结果完全一致。 对 2  e 的估计 2 ˆ  e ：由（19A）得残差平方和 [ ˆ ˆ ( ) ˆ ( )] 1 1 3 2 2 3 2 1 3 1 2 0      = = R =   yi j − y + Y −Y + y −Y n i j i    =2290－11.1667×198－2.3333×6－0.8333×48 =2290－2265.0048=25 依（20A）