混合线性模型参数估计与检验研究

山东理工大学 毕业设计（论文）

题 目：

学 院： 理学院 专 业： 统计学 学生姓名： 刘美倩 指导教师：

毕业设计（论文）时间：二ОО 九年 4 月 20 日～ 6月21日 共 九 周

中文摘要

摘要

混合线性模型是一种既包含固定效应又包含随机效应的一类线性模型，它在生物、医学、经济、金融、环境保护、工业设计等都具有广泛应用. 近年来关于这个模型的理论和应用研究都有了很大的发展，本文将综述它在参数（固定效应和方差分量）估计方面的一些重要理论结果，重点介绍该模型的方差分量的估计方法，并以具体模型为例来说明方差分量估计的思路和过程. 方差分量的估计方法主要有方差分析估计（ANOV AE ），极大似然估计（MLE ），限制极大似然估计(RMLE)，最小范数二次无偏估计(MINQUE)及谱分解估计(SDE) .这几种估计有不同的优缺点，我们将在本文中一一讲解. 此外我们将对所做的方差分量的估计进行检验.

关键词:混合线性模型；固定效应；随机效应；方差分量估计

Abstract

Abstract

The linear mixed model is a subclass of linear models which contains both fixed effects and random effects. The model has extensive applications in biology，medical science ，economy ，finance ，environment protection and design of industry and so on. In the recent years，there is great progress in the research on the theory of the model and its applications .The letter will describe some important theory in parameter(fixed effects and random effects)estimator ，especially introduce the estimator method of variance components and explain the process of variance components estimator taking the specific models for example The estimator methods of variance component mainly contain analysis of variance estimate，maximum likelihood estimate，restricted maximum likelihood estimate ，minimum norm quadratic unbiased estimate and spectral decomposition estimate .These estimators have different strong point and weak point ，which will be explained in the letter. Besides ，we will test variance components estimators which we have made.

Keywords: Linear Mixed Model, Fixed Effects ，Random Effects ，Variance Components Estimator

目 录

中文摘要 ................................................................................................... Ⅰ Abstract ……………………………………………………………………………………Ⅱ 目录 . .......................................................................................................... Ⅲ 第一章 引 言 .......................................................................................... 1 1.1 课题的背景和意义 ..................................................................................... 1 1.2 混合线性模型的发展 ................................................................................. 2 第二章 固定效应的估计 ............................................................................ 3 2.1混合线性模型的形式 . ................................................................................. 3 2.2 固定效应的估计 ......................................................................................... 4 第三章 方差分量的估计 ............................................................................ 6 3.1 方差分析估计 ............................................................................................. 6 3.2极大似然估计 . ............................................................................................. 9 3.3 限制极大似然估计 ................................................................................... 12 3.4 最小范数二次无偏估计 ................................................................... 13 3.5 谱分解估计.............................................................................................. 16 第四章 随机效应值的预测及方差分量估计的检验 ................................. 20 4.1 随机效应值的预测 .......................................................................... 20 4.2 方差分量估计的检验 ............................................................................ 21 第五章 方差分量估计方法的比较 .......................................................... 22

结论 . ........................................................................................................ 24 参考文献 ................................................................................................. 25 致谢 . ........................................................................................................ 26

第一章 引言

1.1 课题的背景和意义

方差分析模型 Y=Xb+e,是大家所熟知的，将模型中各种效应分为二大类，固定效应和随机效应，这便得到了下面将要研究的混合线性模型. 混合线性模型分为平衡数据的混合线性模型（对所有因子的水平组合，重复试验次数相同的模型）和非平衡数据的混合线性模型（对所有因子的水平组合，重复试验次数不相同的模型）.

混合线性模型的统计分析包括二方面，首先是估计各项随机效应的方差分

2量σn ，其次是估计固定效应β，随机效应是不可估计的，但可以进行预测，

另外，可以对方差分量进行检验和区间估计.

混合线性模型中，研究最多的是对于方差分量的估计方法，目前来看，主要方法有方差分析法（Henderson 方法三）、极大似然估计、限制极大似然估计、最小范数二乘估计及谱分解估计，不同方法中还有不同的求解方法，此方面的研究还不尽完善. 对于方差分量的估计方法还有待发展，各种方法都有其不足之处.

对于固定效应的估计方法较少，比如最小二乘估计方法，极大似然估计在得到方差分量估计的同时也有固定效应的估计.

最后Markov 链蒙特卡罗（Markov chain Marte Carlo）分析方法，不但可以获得参数的点估计，还能提供这些参数的相应分布及其特征值.

1.2 混合线性模型的发展

近30年来，关于混合线性模型的参数估计一直是线性模型的最活跃的研究方向之一，在这方面已有一些专著将混合线性模型应用于处理纵向数据分析. 对于固定效应的估计方法为最小二乘估计而. 随机效应是不可估计的但可以预测.

对于随机效应的方差分量估计方法有多种：

1953年Henderson 提出了三种估计非平衡数据方差分量的统计分析方法，这些方法都模仿方差分析的原理. 在这三种方法中，方法三运用最为广泛. 适用于析

因试验设计或套试验设计等正规试验设计的方差分析，但对一些特殊的线性模型，采用此法估计的方差分量可能会产生偏差，70年代以来发展的一系列混合线性模型分析方法，包括最大似然法，限制性最大似然法，最小范数二次无偏估计法等. 可以克服Henderson 方法三的局限性.

自70年代开始，Hartley 和Rao(1967)首先提出应用最大似然法分析混合线性模型的非平衡数据. 但这个方法的估计值受到固定效应的影响，可能导致较严重的有偏估计，为了服这一缺点，Patterson 和Thompson(1971)提出了限制性最大似然法，使似然值不包括固定效应.

Rao （1971）提出最小范数二次无偏估计，比最大似然法和限制性最大似然法更简便和优越，它不需要进行迭代运算，对线性模型也没有正态分布的限定.

最近，王松桂和尹素菊提出了同时估计固定效应和方差分量的一种新方法，称为谱分解估计.

第二章 固定效应的估计

2.1 混合线性模型的形式

混合线性模型的一般形式是：

y =Xβ+U 1ξ1+U 2ξ2+„+U k ξk ，

其中y 是n ⨯1观测向量，X 为n ⨯已知设计阵，β 是p ⨯1非随机的参数向量，称为固定效应. U i 为n ⨯t i 已知设计矩阵，ξi 为t i ⨯1随机向量，称为随机效应. 且E(ξi )=0，i=1，„，k .

通常假设

Cov(ξi )=σi 2I t i ，Cov (ξi ，ξj )=0，i ≠j. 于是我们有

E(y )=Xβ， Cov （y ）=∑σi 2U i U I / ∑(σ2) ，

i =1k

这里σ2=(σ12, , σk 2). σi 2称为方差分量，相应的模型也称为方差分量模型. 混合线性模型之所以在生物、医学、经济、金融和工业设计等领域具有十分广泛的应用，是因为在实际问题中，参与试验的个体是随机抽取的，且我们研究的目标不是这些个体本身的特征，而是它们所在的总体特征，这时把个体效应看作随机的而引入模型，可以大大提高模型的精度.

对于混合线性模型，我们感兴趣的参数分两类：固定效应β 和方差分量

σ2=(σ12, , σk 2)，它们分别包含在均值E(y ) 和协方差Cov （y ）中，处理问题

的方法与固定效应模型相比较复杂一些.

在本章里我们将研究固定效应估计方法，以后几章里详细讲解方差分量的估计方法及检验.

2.2固定效应的估计

在考虑到固定效应的估计时，将模型写为如下形式 y =Xβ+Uξ+e，

其中β为固定效应，ξ为随机效应，且E(ξ)=0，E(e)=0,Cov(ξ，e)=0，并且假设

ξ和e 的协方差阵具有较一般的形式Cov(ξ)=D ≥0,Cov(e)=R>0,于是我们有∑=Cov(y )=UDU+R>0.

暂时视D,R 已知应用最小二乘估计法得到正则方程 X '∑-1X β*=X '∑-1y ,

据此可得到β的广义最小二乘解β*= (X '∑-1X )X '∑-1y . 因此，任意的可估

-

函数c /β的最佳线性无偏估计（BLU 估计）为

c 'β=c '(X /∑-1X )X '∑-1y .

*

-

而实际上D,R 未知，若用它们的估计D , R 代替，即用∑=U D U '+R 代替

^^^^^

∑，便得到c /β的两步估计

c 'β（∑）=c '(X '∑X ) X '∑y .

-

~

^

^-^-1

在假设Cov(ξi )=σi 2I t i ，Cov (ξi ，ξj )=0，i ≠j 下c 'β的两步估计又可变形为

c 'β（σ）=(X '∑(σ) X ) X '∑(σ) -1y .

在以上模型中若假设e, ξ的联合分布关于原点对称，设σ=σ2（y ）是σ2

的一个估计，它是y 的偶函数，且具有变换不变性. 对一切可估函数c 'β，若E （c 'β（σ））存在，则两步估计c 'β（σ）必为c 'β无偏估计. 以上定理说明了在一定条件下c 'β（σ）为c 'β无偏估计.

关于e, ξ分布的假设在许多情况下是满足的. 此外，对方差分量估计的方法

~

~

~

^2^2

-1-

^2

^2

^2

~

^2

^2

不同，所得到的固定效应的估计也不相同，方差分析法、极大似然法、限制极大似然法、最小范数二乘法, 所产生的估计σi 都是y 的偶函数，且是变换不变的，因此，对于混合线性模型的固定效应以上定理给出了一大类两步估计的无偏性. 我们可以看到方差分量的不同估计，往往会产生不同两步估计，而且除最小二乘

^2

估计外，模型还可能存在另外一些简单估计，如Panel 数据下的Between 估计和within 估计，简约估计等.

第三章 方差分量的估计

3.1 方差分析估计

1953年Henderson 提出了三种估计非平衡数据方差分量的统计分析方法，这些方法都模仿方差分析的原理. 在这三种方法中，方法三运用最为广泛. Searle （1968）将Henderson 的分析方法改成矩阵形式表示. 方差分析方法渊源于固定效应模型的方差分析.

我们用下面的例子来说明它的原理和方法 对于平衡单向分类模型

i =μ+αi +e ij ， i =1，„，a ， j =1，„，b,

其中μ为总均值，是固定效应，α1, , αa 为随机效应. 假定所有αi , e ij 都不相

2关，且其均值为0，方差为Var(αi )=σα,Var(e ij )=σe 2. 记y '=（y 11, , y ab ）. 暂

时先把αi 看作因子A 的i 水平A i 的固定效应，有

RSS(μ)=y .. /(ab ) SS μ, (3.1.1) 其自由度为1，对应于α1, , αa 的平方和，即因子A 的平方和 SS A =RSS(μ, α)-RSS(μ)=∑

i

-2

∑(y

j

-

i .

-y .. ) 2， (3.1.2)

-

其自由度为n-1，而残差平方和为 SS e =y y -RSS(μ, α)=∑

/

I

∑(y

J

ij

-y i . ) 2, (3.1.3)

-

其自由度为

a （b-1）. 由以上三式可以推导出总平方和的分解式 y /y =SS μ+SS A +SS e 将平方和除以自由度，得到均方： Q 0=SS μ, Q 1=SS A /(a -1), Q 2=SS e /[a (b -1)],

再按照αi 为随机效应的假设，求出各均方的均值：

2

E (Q 0) =E (SS u ) =ab μ2+b σα+σe 2, 2 E (Q 1) =E (SS A /(a -1)) =b σα+σe 2,

E (Q 2) =E (SS e /[a (b -1)])=σe 2.

可以看到，后两式的右端为方差分量的线性函数，令E (Q 1) =Q 1， E (Q 2) =Q 2，

2便得到关于σα，σe 2的线性方程组 22

⎧⎪b σa +σe =Q 1,

⎨2

σ=Q . ⎪2⎩e

解此方程组得

σe =Q 2， σα=(Q 1-Q 2)/b.

2

它们就是方差分量σα，σe 2的方差分析估计（ANOVA 估计）.

^2^2

从以上求解过程我们可清楚的看到方差分析法的具体思路，我们可以把方差分析法归纳如下：

1.对一个方差分量模型，现将其随机效应看作固定效应，按通常方差分析方法算出各效应对应的平方和（或均方）.

2.求这些平方和（或均方）的均值（此时的随机效应不再看作固定效应），他们是方差分量的线性函数.

3.令这些平方和（或均方）等于它们各自的均值，得到关于方差分量的一个线性方程组，解此方程组便得到方差分量的估计.

至此我们了解到方差分析的一般思想，将此方法用于一般的混合线性模型. 先考虑一个简单模型（方差分类模型）

y =Xβ+U 1ξ1+U 2ξ2+e,

2

'+ σe 2I ∑(σ2) . 其中Cov(ξi )=σi 2I t i ，i=1，2，Cov(y)=σ12U 1U 1'+σ2U 2U 2

2

得到方差分量σ12，σ2，σe 2的线性方程组

2

⎧a 1σ12+(a 2-a 3) σ2+r 2σe 2=SS ξ1, ⎪⎪22

⎨a 2σ2+σe =SS ξ2,

⎪2(n -r -r -r ) σ=SS e . ⎪123e ⎩

2

解此方程组，得到σ12，σ2，σe 2的估计. 它们就是这些方差分量的ANOV A

估计. 估计如下

σ1={SS ξ1-r 2（SS e /(n -r 1-r 2-r 3) ）

-(a 2-a 3) [（SS ξ2-SS e /(n -r 1-r 2-r 3) ）/a 2]}/a 1,

^2

σ2=（SS ξ2-SS e /(n -r 1-r 2-r 3) ）/a 2， σe =SS e /(n -r 1-r 2-r 3) .

对于一般的混合线性模型，文献中也称为拟合常数法，之所以称其为拟合常数法就是因为在构造估计方程时，我们把随机效应看成固定效应，即常数. 对于平衡数据模型，该方法的平方和分解是惟一的，且可根据方差分析得到. 下面将方差分析方法应用于带有交互效应的两向分类混合模型.

^2

^2

y ijk =μ+αi +βj +γij +e ijk ，

i=1,„，a ，j=1，„，b ，k=1，„，e ，

这里μ，并满足通常的假设，即所有的βj ，αi 为固定效应，βj , γij 为随机效应，

2

，Var(γij )=σγ2，γij ，e ijk 都不相关，且具有均值为0，方差为Var(βj )=σβ

Var(e ijk )=σe 2.

首先，暂时视βj , γij 为固定效应，得到总平方和有如下分解： y 'y =SS μ+SS α+SS β+SS γ+SS e ，

再用随机效应的平方和除以各自的自由度，得到均方Q 1=SS β/(b -1) ，

Q 2=SS γ/(a -1)(b -1) ，Q 3=SS e /[ab (c -1)]，求出它们的均值，并令这些均值等

2

于对应的均方，得到关于σβ，σγ2，σe 2的线性方程组

2

⎧ac σβ+c σγ2+σe 2=Q 1, ⎪⎪22

⎨c σγ+σe =Q 2,

⎪2⎪⎩σe =Q 3.

解此方程组，得到方差分量的估计：

σβ=（Q 1-Q 2）/(ac ) ，σγ=（Q 2-Q 3）/c ，

^2

^2

σe =Q 3.

由于方差分析法给出的估计σ 作为一个线性方程组的解，他们未必是正的. 这是方差分析法的一个缺陷. 至于如何对待方差分量的负估计，目前尚无一致的看法. （1）观点认为，若某个σi

研究表明，除非非负性可自然满足，不然无偏性和非负性几乎不可兼得，在近期的研究中，人们往往放弃无偏性，而要求具有另外的优良性，如均方误差最小 .

^2

^2

^2

3.2 极大似然估计

方差分析方法有一定的缺陷，采用此法估计的方差分量可能会产生偏差，

所以，人们开始寻找新的方法. 自70年代开始，Hartley 和Rao(1967)首先提出应用极大似然法分析混合线性模型的非平衡数据. 极大似然方法最早应用在统计中，后来被应用在混合线性模型中. 此方法能同时获得固定效应和方差分量的估计.

考虑一般的混合线性模型

y =Xβ+U 1ξ1+U 2ξ2+„+U k ξk ， （3.2.1） 这里假设ξi ~N(0，σi 2I t i ) ，i=1，„k ，所有ξi 都相互独立，记V i =U i U i '，

σ2=(σ12, , σk 2) '，于是

Cov （y ）=∑σU i U I '=∑σi 2V i ∑(σ2) .

2i i =1

i =1

k

k

我们假设∑(σ2) >0，因此y ~N n （0，∑(σ2) ），所以未知参数β，σ12，„，σk 2的似然函数为

L(β, σ|y )=(2π)

2

-

n

2

1

|∑(σ) |exp{-(y -X β) /∑(σ2) -1y -X β}，

2

2

-

12

取对数，略去常数项及常数倍，得

l(β, σ2|y )=-ln |∑(σ2) |-(y -X β) '∑(σ2) -1(y -X β)

=-ln |∑(σ2) |-tr ∑(σ2) -1(y -X β)(y -X β) '. （3.2.2） 根据如下公式

∂Ax ∂∂A (t )

]. =A , ln |A (t ) |=tr [A (t ) -1

∂t ∂t ∂x

对σ2，β求导可得

∂l

=-tr (V i ∑(σ2) -1) +tr [(∑(σ2) -1V i ∑(σ2) -1)(y -X β)(y -X β) ']，i=1，„，k, 2

∂σi

∂l

=-2X '∑(σ2) -1X β+2X '∑(σ2) -1y . ∂β

令这些导数等于零，得到似然方程

2-12-1

⎧⎪X '∑(σ) X β=X '∑(σ) y ,

⎨（3.2.3） 2-12-12-1

⎪⎩tr (V i ∑(σ) ) =(y -X β) '(∑(σ) V i ∑(σ) )(y -X β).

i=I,„，k . 进行简化变形后似然方程变为

⎧X β=X (X '∑(σ2) -1X ) -X '∑(σ2) -1y ,

⎪k

i =1, „，k . （3.2.4） ⎨2-122-12-1

''tr [V ∑(σ) ]σ=y (I -P ) (∑(σ) V ∑(σ) )(I -P ) y , i j σi σ⎪∑j =1⎩

这就是我们要求的似然方程. 由（2.2.4）得第一方程，任意可估函数c 'β的ML 估计为c 'β(σ) =X (X '∑(σ) X ) X '∑(σ) y ，其中σ为σ2的ML 估计.

^

^2

^2

-1

-

^2

-1

^2

在一般情况下似然方程（2.2.4）没有显式解，即便在有显式解的情形，σ2的解未必是非负的，若为负值，它就没有落在参数空间内，所以并不是ML 估计. 这时一般采取截断法，即取max{σi ,0}作为ML 估计，在没有显式解的情形只能用迭代法求解.

下面详细介绍两种迭代方法，并简单提出其他方法. Anderson等提出一种迭代法是 σi 这里σ

^2(m )

^2

^2(m +1)

=H (σi

^2(m )

) h (y , σ

-1

^2(m )

) ，

为σ的第m 次迭代值σ

2

^2(m )

. 当σ得两次相邻迭代值相差不大时，迭代

^2

停止，这就得到了方差分量的估计. 带入（3.2.4）的第一方程，便可得到固定效应的估计.

另外一种迭代法是由Hartley 和Rao 提出的，其推广形式是 σi

^2(m +1)

=σi

^2(m )

h (y , σ) i 2(m ) , i =1, , k .

^

tr (∑(σ) -1V i )

^2(m )

这个迭代的一个好处是，当初始值为非负时，后面的迭代值永远不会取负值.

另外还有一些迭代方法，如Newton-Raphson 方法，得分方法，以及EM 算法，这里不做详细研究.

下面以两向分类混合模型为例来具体介绍参数估计的过程 对两向分类混合模型

y ij =μ+αi +βj +e ij , i =1, , a , j =1, , b ,

2 这里μ, αi 为固定效应，βj 为随机效应，βj ~N (0,σβ), e ij ~N (0,σe 2) ，且所

有βj e ij 都相互独立，该模型的矩阵形式为

y =X 1μ+X 2α+U β+e ,

用Kronecker 乘积表示设计阵X 1, X 2和U ： X 1=1ab =1a ⊗1b , X 2=1ab =I a ⊗1b ,

U =1a ⊗I b .

/

固定效应的设计阵为X =(X 1 X 2) ，P X =P X 2=I a ⊗J b , 这里J b =1b 1b /b .

--

得到似然方程

X 1μ+X 2α=P X y ,

-1b -

a σβ+σe =∑(y . j -y .. ) 2，

b j =1

^2

^2^

^

b a σβ+σe

^2

^2

+

(a -1) b

σe

^2

=

1(a σβ+σe ) 2+(a -1) b

^4

^2^2

∑(y

i =1

k

-

2

-y ) . j ..

-

σe

∑∑(y

i

j

ij

-y i . -y . j +y .. ) 2

---

.

求得方程组的显式解

---1

σe =(y ij -y i . -y . j +y .. ) 2， ∑∑(a -1) b i j -1b -1^22

σβ=∑(y . j -y .. ) -σe ，

ab j =1a

^2

^2

^2

在此处，σβ可能取负值.

3.3 限制极大似然估计

方差分量的极大似然估计的一个缺陷是在导出方差分量的估计的过程中没

有考虑到固定效应β的估计所引起的自由度的减少. 为此，Patterson 和Thompon 提出的一种修正方法，称为限制极大似然法. 该方法的思想是基于极大似然估计残差，利用极大似然法导出方差分量的估计. 限制极大似然估计的特点是此方法所求方程的解与方差分析法所得的估计相同.

下面我们用两向分类随机模型来说明限制极大似然法的参数估计过程

y ijk =μ+αi +βj +γij +e ijk ,

i =1, , a , j =1, , b , k =1, , e ,

2

这里为μ为总平均，是固定效应，αi , βj , γij 都为随机效应，假设αi ~N (0,σa ) ，

22

βj ~N (0,σβ) ，γij ~N (0,σγ2), e ijk ~N (0,σα) 且都相互独立. 该模型的矩阵形式为

y =X μ+U 1α+U 2β+U 3γ+e ,

这里X =1a ⊗1b ⊗1c ，U 1=I a ⊗1b ⊗1c ，U 2=1a ⊗I b ⊗1c , U 3=I a ⊗I b ⊗1c ，方差阵

2222

'''''∑(σ2) =σαI α⊗11b b ⊗11c c +σβ1a 1a ⊗I b ⊗11c c +σγI a ⊗I b ⊗11c c +σe I abc .

得到方程组如下

2

bc σα+c σγ2+σe 2=Q 1,

2ac σβ+c σγ2+σe 2=Q 2,

c σγ+σ=Q 3,

2

2

e

σe 2=Q 4,

这里我们记

Q 1=∑(y i .. -y ... ) 2/(a -1),

i =1b a

-

-

Q 2=∑(y . j . -y ...) 2/ac (b -1),

j =1a

--

Q 3=∑∑(y ij . -y i .. -y . j . +y ... ) 2/(a -1)(b -1),

i =1j =1

b

----

Q 4=∑∑∑(y ijk -y ij . ) 2/ab (c -1).

i =1j =1k =1

a b c

-

这与方差分析法所得到的线性方程组相同，因此限制极大似然方程的解与方差分析估计相同. 对平衡数据的混合效应模型，这种现象通常成立.

3.4 最小范数二次无偏估计

方差分量的最小范数二次无偏估计始于Townsend 和Harville ，但C.R.Rao 的工作奠定了重要的基础，此方法与前面介绍的三种方差截然不同. 前三种方法都是先按已有的一定方程式去求估计，至于所得估计有何性质，事先并不知道. 而最小范数二次无偏估计的基本思想是对模型的误差和随机效应的分布没有要求，而是先提出估计应具有的性质，然后把为满足这些性质所加的条件提成一个极值问题，即所谓最小迹问题. 解所得的最小迹问题，便得到所要的估计. Rao 这个方法的缺点是，对N 个不同的人，由于各自选用的方差分量的先验值不同，就会得到N 个不同的估计，具有一定程度的主观随意性.

考虑最一般形式的方差分量模型

y =X β+U 1ξ1+ +U k ξk ， （3.4.1）

这里X n ⨯p ，U i , n ⨯t 为已知设计矩阵，β为p ⨯1固定效应向量，ξi 为t i ⨯1随机效应向量，满足E （ξi ）=0，Cov(ξi )=σi 2I t i , ξi 都不相关，若记

', , ξk ') ， U =(U 1 U 2 U K ) ，ξ'=(ξ1', ξ2则模型（3.4.1）可改为

y =X β+U ξ, E (y ) =X β, Cov (y ) =∑σi 2V i ∑,

i =1k

2

我们的基本目的是估计方差分量σ12, , σk 及其线性函数ϕ=c 'σ2，这里2

''σ2=(σ12, , σk ) ，c =(c 1, c k ) .

首先，看所求的估计量应有的一些性质，因为要估计的参数是方差，所以考虑二次型估计y /Ay ，A 为对称阵，要求这个估计具有下述性质.

(1)不变性 即估计y /Ay 关于参数β具有不变性. (2)无偏性

(3)最小范数准则 欲使y 'Ay 为一个好的估计，那么对一切ξ，矩阵U 'AU 与

∆在某种意义下相差很小. 用矩阵范数||U 'AU -∆||来度量U /AU 与∆相差大小，

则应该选择A 极小化范数||U 'AU -∆||. 若线性函数ϕ=c 'σ2的估计y 'Ay 满足 AX =0,

tr (AV i ) =c i ， i =1, , k , 且使范数||U 'AU -∆||达到极小，则称为最小范数二次无偏估计.

2222

这里采用加权欧氏范数，令权矩阵 W =diag {σ0,1其中为I t 1, , σ0, I }σσi k t k 0, i

的一个预先指定值（先验值），因此W 也就是Cov (ξ) 的一个预先指定阵（先验阵）.

12

12

定义F =W (U 'AU -∆) W , 则加权欧氏范数

||U 'AU -∆||=tr (F 'F ) =tr [W (U 'AU -∆) W (U 'AU -∆) W ]

12

12

=tr (W 'AUWU 'AUW ) -2tr (W 'AUW ∆W ) +tr (∆W ) 2. 利用无偏性，上式第二项

tr (W 'AUW ∆W ) =tr (U 'AUWU 'A ∆W ) =∑

i =1

2

再记V w =∑σ0, i V i ，于是

i =1k

1

2121212

1212

k

4c i σ0, i

t i

(AV i ) =∑

i =1

k

4

c i 2σ0, i

t i

=tr (∆W ) 2.

||U 'AU -∆||=tr (AV w ) 2-tr (∆W ) 2.

这样，对加权欧氏范数求ϕ=c /σ2的最小范数二次无偏估计的问题，归结为求下述极值的解

min tr (AV w ) 2

⎧AX =0

⎨ （3.4.2）

tr (AV ) =c , i =1, , k . i i ⎩它的目标函数是矩阵的迹，称（3.4.2）为最小迹问题.

下面考虑极值问题的解是否存在，有如下定理1 极值问题（3.4.2）的解为

A =B w (∑λi V i ) B w ,

*

i =1k

其中

-1-1-1-1

B w =V w ， -V w X (X 'V w X ) -X 'V w

且λi , i =1, , k 为方程组

∑tr (B w V i B w V j ) λi =c j ，j =1, , k

i =1

2

的解，这里V i =U i U i ，V w =∑σ0, i V i .

i =1k

k

证明定理1等价于证明如下定理2 极值问题（3.4.2）的解为

A =N (∑λi V i ) N ,

*

i =1k

其中N =I -X (X 'X ) -X ，λi , i =1, , k 为方程组

∑tr (NV i NV j ) λi =c j ，j =1, , k

i =1k

的解.

尽管方差分量的估计值取决于人为选择的先验值，只要这些先验值不依赖于试验数据，MINQUE 估计量仍是无偏的，如果用估计值替代先验值重新估计，便可获得新的估计，重复这一过程，直到新的估计非常接近旧的估计为止，这种迭代的估计方法就是上一节介绍的限制极大似然（REML ）方法，其估计结果就是限制极大似然估计，所以在多数情况下，REML 估计和MINQUE 估计是比较接近的.

最后说一下先验值选取的问题，先验值的选择可以凭经验，或根据以往的分析结果. 最简单的方法是取残差效应的先验值为1，其它所有的先验值为0，这种方法称为MINQUE （0）法，所获得的方法分量的估计量是MINQUE （0）估计量. 这种取值的方法的优点是可以省去大矩阵的求逆计算，算法简单，但是由MINQUE （0）方法估计的方差分量的抽样方差往往较大. 另一种简便易算的方法是设所有的先验值为1，这种方法称为MINQUE （1）法，所获得的方差分量是MINQUE （1）估计量.

3.5 谱分解估计

下面我们介绍一种估计方差分量的新方法，即谱分解估计，首先我们引入谱分解的有关知识.

协方差阵的谱分解有多方面的应用，一方面，方差分量的谱分解估计就是以协方差阵∑(σ2) 的谱分解式为基础构造的；另一方面，从∑(σ2) 的谱分解式，可以直接获得∑(σ2) 的逆矩阵和行列式，而这些恰是利用极大似然方法求解方差分量的极大似然估计所必需的.

目前有两种谱分解算法，一种是由Smith 和Hocking 提出的，他们基于完全设计的随机效应模型，给出了谱分解的公式，并将它推广到一般的混合模型的情

形，即将公式中未出现的随机项的方差分量用零取代. 另一种算法是由Searle 和Henderson 提出的，该算法没有利用模型随机效应设计阵U 1, , U k 之间的某些约束关系，而是基于对∑(σ2) 扩充很多项的方法，将∑(σ2) 改写成

∑θj p , j 1(J n p p ⊗ ⊗J n j 11) ，这里（j p , j 1）取遍从（0，0，„0）到（1，1，„1）

j

的所有2p 个二进制向量，即添加（2p -k-1）个θj p , j 1=0的项到∑(σ2) 中，给出了一般形式的谱分解公式. 最近吴密霞，王松桂就一类模型研究了协方差矩阵的谱分解，互异特征值个数的确定及求解问题，将谱分解的结果应用到构造模型参数估计，提出了谱分解估计. 史建红在他的博士论文中把结果推广到一般平衡混合线性模型.

下面我们详细介绍王松桂和尹素菊提出的同时估计固定效应和方差分量的一种新方法，即谱分解估计.

谱分解估计的基本思想是：首先对协方差阵进行谱分解，然后利用谱分解得到的主幂等阵对原模型进行适当的线性变换，获得若干个新的奇异线性模型，这些新模型的特点是它的固定效应与原模型相同，但新模型的协方差阵除了一个因子（这个因子是原模型协方差阵的一个特征值）外，不含未知的方差分量，利用最小二乘统一理论，对每个新模型可以得到固定效应和特征值的一个估计，由于在常见情形下协方差阵的特征值是方差分量的线性函数，因此通过解线性方程组可以获得方差分量的估计.

新方法的突出特点是能同时给出固定效应和方差分量的估计，前者是线性的，后者是二次的，且相互独立. 对于固定效应可以获得若干个谱分解估计，它们都是具有一些好的统计性质的线性估计. 因此，利用这些估计可以对模型做进一步假设检验，区间估计，以及模型诊断等一系列统计推断.

对于任一混合线性模型

y =X β+∑U i ξi +e , （3.5.1）

i =1k

如果y 的协方差阵∑(σ2) 有如下谱分解

∑(σ) =∑λi M i , (3.5.2)

2

i =1q

这里λi ，i=1,„，q 是∑(σ2) 的所有互异非零特征根，它们是σ2线性函数，M i 是特征根λi 对应的主幂等阵（即M =M i , M i M j =0. i ≠j , ∑M i =I ），且独立于未知

2

i

i =1q

参数. 我们分别用幂等阵M i (i=1, „，q) 左乘模型(2.5.1).于是得到变换后的q 个新模型

M i y =M i X β+εi ， εi ~（0，λi M i ），i=1,„，q. （3.5.3） 这些模型的特点是，模型协方差阵除λi 之外，独立于未知参数. 因为M i 是奇异阵，我们可以应用最小二乘统一理论获得原模型参数（固定效应和方差分量）的估计，称为谱分解估计.

下面我们简略介绍一下最小二乘统一理论， 最小二乘统一理论应用于线性模型的参数估计，由著名统计学家Rao 应用推广. 对于线性模型 y =X β+e ,

E (e ) =0. Cov (e ) =σ2∑，

如果|∑|=0，则称该模型为奇异线性模型. 对于奇异线性模型，因为∑-1不存在，于是一般的最小二乘法中用到的Q (β) 无定义，Rao 成功解决了这个问题. 关键是寻找一个新矩阵T ，能够充当∑-1所负担的作用.

T=∑+XUX ', 其中U ≥0,rk(T)=rk(∑ X ), 然后定义

Q (β) =(y -X β) 'T '(y -X β) ,

用最小化Q (β) 求出最小值点

β*=(X 'T -X ) -X 'T -y ,

还能证明对任一可估函数c 'β, c 'β*为其BLU 估计. 这个结论既适用于设计阵X 列满秩或列降秩的情形，又适用于∑奇异阵或非奇异阵的情形。正是由于这个原因，通常把这个结果称为最小二乘统一理论.

对于变换后的新模型（3.5.3）由最小二乘统一理论，便得可估函数c 'β及其方差参数λi 的一个无偏估

c 'β=c '(X 'M i X ) -X 'M i y ，

λi =y '(M i -M i X (X 'M i X ) -XM i ) y /m i ，

2

'i =1, , k +1，这里m i =r i -rk (M i X ) . 得到σ2=(σ12, , σk +1) 的无偏估计

~2

^

~(i )

^

σ=H λ.

-1

c β和σ是模型参数的谱分解估计.

T

~(i ) ^2

第四章 随机效应值的预测及方差分量估计的检验

4.1 随机效应值的预测

在混合线性模型的分析中，通过方差分量的估计可以推断随机效应的变异性，在有些试验中，试验者不但需要估计方差分量，而且需要推断模型中某些随机效应的值，混合线性模型中的随机效应是不可估计的，但却是可预测的.

对于以下模型

y =X β+U ξ+e ，E (e ) =0, Cov (e ) =R >0, E (ξ) =0, Cov (ξ) =D ≥0，

这里假设D , R 都是已知的，如果它们含有未知参数，在实用中用它们的估计代替，这样在随机因素的方差是已知条件下，可以获得随机效应的最佳线性无偏预测（BLUP ）. 得到的估计公式为

ξ=DU '(UDU '+R ) -1(y -X β*) .

按上式预测的随机效应是相依变量的线性函数，因而是多元正态随机变量，其期望值是零. 随机效应的BLUP 预测值是参数的无偏预测，并且在所有的线性无偏预测中具有最小方差，而实际上在运用混合线性模型分析具体试验数据时，方差分量是需要被估计的，因而并不可能获得随机效应值的最佳线性无偏预测，在实际应用时常采用方差分量的估计代替其参数.

^

4.2 方差分量估计的检验

考虑混合线性模型

y =X β+U 1ξ1+U 2ξ2+e ,

22我们假设ξ1~N (0,σ1I s ), ξ2~N (0,σ2I ), e ~N (0,σ) 2q 0I n ，且它们彼此独立，其中

22

'+σ2s , q 分别为已知阵U 1和U 2的列数. 于是Cov (y ) =σ12U 1U 1U 2U '+σ2e I n . 在前面2

基于拟合常数的思想，我们给出了此模型方差分量σ12，σ2和σe 2的估计，现在我

们用同样的方法来构造随机效应ξ2是否存在的检验问题，即

22 H 0:σ2=0H 1:σ2≠0,

的检验统计量.

将以上模型中随机效应ξ1和ξ2暂视为固定效应，模型拟合之后的残差平方和为

SS e =y /y -RSS(β, ξ1, ξ2)=y '(I n -P (X :U 1:U 2) ) y .

2

若σ2=0，则模型变为y =X β+U 1ξ1+e ，同样暂视其中的随机效应ξ1为固定效应

时，模型拟合后的残差平方和变为

SS e 0=y /y -RSS(β, ξ1)=y '(I n -P (X :U 1) ) y .

2直观上，当σ2=0时，模型拟合后的残差平方和SS e 0和SS e 应很接近，即SS e 0-SS e 2相对于SS e 应很小，不然，我们就认为随机效应ξ2作用显著，即接受σ2≠0. 令

F =

(SS e 0-SS e ) /t 2

，

SS e /t 1

这里t 1=n -rk (X :U 1:U 2), t 2=rk (X :U 1:U 2) -rk (X :U 1) .

此外还有一种方法不但可以估计方差分量，还能有效的检测方差的显著性（即对方差分量估计进行检验）. 这种方法叫做刀切法（Miller,1974; 朱军，1992），现以试验处理重复抽样为例，介绍估计遗传参数刀切法估计的方法. 对于遗传参数（方差分量或协方差分量），分析所有观察值所得到的估计. 当从资料中剔除第i 种试验处理世代均值的r 次观察值以后，重新估计而得到新的估计，采用此法从完整的数据中每次剔除r 次观察值，可以获得g 个不同的估计，刀切法估计出的方差分量具有t 分布，从而可以用t 检验对遗传参数作统计检验.

Zhu(1989)用蒙特卡罗模拟比较了估计方差分量估计的方差及参数检验的效益，采用的两种统计分析方法分别是：MINQUE 估计方差的渐进估计法和刀切法抽样估计法. 在遗传分析中一些遗传参数估计值的抽样方差也可以用刀切法算得，并可用t-检验作统计检验.

第五章 方差分量估计的比较

在前面的研究中，我们讲解了混合线性模型的参数（固定效应和方差分量）的估计，并简单介绍了随机效应的预测及方差分量估计的检验. 重点是方差分量的估计，接下来我们对方差分量的各种估计方法进行比较，来给出各自的优缺点.

目前，最常用的方差分量的估计有方差分析估计（ANOV AE ），极大似然估计（MLE ），限制极大似然估计(RMLE)，最小范数二次无偏估计(MINQUE)及谱分解估计(SDE) .

方差分析估计：用一句话来总结它的思想是，令各种平方和等于各自平方和的期望值，得到一个关于方差分量的线性方程组，求解这个方程组便得到方差分量的估计. 对于平衡数据混合模型，这种估计方法的优点是无偏性和方差最小性，但对于非平衡模型，这种方法具有不惟一性，而且计算上很麻烦，最主要的缺陷是有时能给出方差分量的负估计值.

极大似然估计：理论上这种方法能同时获得固定效应和方差分量的估计，但一般来说，我们不能获得方差分量的估计的表达式，只有在平衡数据混合模型下，才能得到显式解，再没有显式解的情形下，只能进行计算机迭代，这时需要有较好的迭代算法. 能够处理大维数的矩阵非线性方程，并且导致整体最大值，同时最大值还必需满足非负要求.

限制极大似然估计:此方法是对极大似然估计的一种完善，极大似然估计方法没有考虑估计固定效应所引起的自由度的减少，限制极大似然估计就是基于这个缺陷所做的改进，这种方法与方差分析估计结果相同.

最小范数二次无偏估计：这种估计方法与前几种相比有很大的不同之处，它在给定了方差分量的一组先验值后，只需求解一个线性方程组，不需要进行繁琐的迭代. 缺点是对N 个不同的人，由于各自选用的方差分量的先验值不同，会得到N 个不同的估计，具有一定程度的主观随意性.

谱分解估计：此估计方法是最近提出的新方法，能同时给出固定效应和方差分量的估计，对于固定效应可以获得若干个谱分解估计，都是具有良好统计性质的线性无偏估计，方差分量的估计是二次不变无偏估计，并且SD 估计与ANOVA

估计一样具有显式解，但SD 估计也不能保证估计的非负性. 已有学者证明，对于一些模型方差分量SD 估计和ANOVA 估计相等，所以SD 估计具有ANOVA 估计的一些优良性质.

结论

本文主要研究了混合线性模型参数的估计及检验，混合线性模型是既包含固定效应又包含随机效应的一类线性模型，此类模型在实际中有多方面的应用，我们既要对固定效应进行估计，又要对随机效应的方差分量进行估计，而随机效应是不可估计的，但实际中我们需要对其进行预测. 此外，我们也要对所估计的方差分量进行检验，通过假设检验，看它所对应的随机效应是否显著有影响，从而判断随机效应在实际中对我们所要研究的变量是否有显著影响.

在本文中我们仅研究了一种估计固定效应的方法，即最小二乘估计法，此方法简单易行，只是在估计过程中我们先暂时假设了随机效应的协方差阵已知，在得到固定效应的估计后，我们还要对协方差阵进行估计，因为实际上它们是未知的，从而得到固定效应的两步估计.

对于方差分量，我们研究了五种方法，方差分析法、极大似然法、限制极大似然法、最小范数二次估计法及谱分解估计方法，其中方差分析法及谱分解估计法都能得到显式解，但不能保证估计的非负性，极大似然法和限制极大似然法在一般情况下需要进行迭代，无法得到显式解，最小范数二次估计需要事先给定一组方差分量的先验值，有很大的主观随机性.

方差分量估计的检验，我们主要研究了一类较简单混合模型的检验问题，构造了一个F 检验，其它检验未做详细说明，对于检验方面的理论研究还不尽完善，还有许多问题需要解决.

参 考 文 献

[1] 王松桂 尹素菊 史建红 吴密霞 线性模型引论 科学出版社 2000 [2] 朱军 线性模型分析原理 科学出版社 2004

[3] 王松桂 邓永旭 方差分量的改进估计 应用数学学报 1999，22：115~125

[4] 吴密霞 王松桂 线性混合模型中固定效应和方差分量同时最优 中国科学 A 辑，2004，34（3）：373~384

[5] 王松桂 尹素菊 线性混合模型参数的一种新估计 中国科学 A 辑，2002，32（5）435~443

[6]SHI JIANHONG WANG SONGGUI Some properties of Spectral Dcomposition Estimate of Variance Components. 应用概率统计，2004，287~294

[7]马铁丰 王松桂 线性混合模型方差分量的检验 应用数学学报，2007，22（4）：443~440

致谢

衷心感谢导师林鹏老师对本人的精心指导，他的言传身教是我受益匪浅，导

师从阅读文献、论文选题等方面提出了宝贵意见，在论文的写作过程中，给了多方面的指导，在百忙之中对论文进行了详细的批阅和审查.

感谢统计系的全体老师和同学的关心和支持！感谢所有关心和帮助过我的人们！

山东理工大学 毕业设计（论文）

题 目：

学 院： 理学院 专 业： 统计学 学生姓名： 刘美倩 指导教师：

毕业设计（论文）时间：二ОО 九年 4 月 20 日～ 6月21日 共 九 周

中文摘要

摘要

混合线性模型是一种既包含固定效应又包含随机效应的一类线性模型，它在生物、医学、经济、金融、环境保护、工业设计等都具有广泛应用. 近年来关于这个模型的理论和应用研究都有了很大的发展，本文将综述它在参数（固定效应和方差分量）估计方面的一些重要理论结果，重点介绍该模型的方差分量的估计方法，并以具体模型为例来说明方差分量估计的思路和过程. 方差分量的估计方法主要有方差分析估计（ANOV AE ），极大似然估计（MLE ），限制极大似然估计(RMLE)，最小范数二次无偏估计(MINQUE)及谱分解估计(SDE) .这几种估计有不同的优缺点，我们将在本文中一一讲解. 此外我们将对所做的方差分量的估计进行检验.

关键词:混合线性模型；固定效应；随机效应；方差分量估计

Abstract

Abstract

The linear mixed model is a subclass of linear models which contains both fixed effects and random effects. The model has extensive applications in biology，medical science ，economy ，finance ，environment protection and design of industry and so on. In the recent years，there is great progress in the research on the theory of the model and its applications .The letter will describe some important theory in parameter(fixed effects and random effects)estimator ，especially introduce the estimator method of variance components and explain the process of variance components estimator taking the specific models for example The estimator methods of variance component mainly contain analysis of variance estimate，maximum likelihood estimate，restricted maximum likelihood estimate ，minimum norm quadratic unbiased estimate and spectral decomposition estimate .These estimators have different strong point and weak point ，which will be explained in the letter. Besides ，we will test variance components estimators which we have made.

Keywords: Linear Mixed Model, Fixed Effects ，Random Effects ，Variance Components Estimator

目 录

中文摘要 ................................................................................................... Ⅰ Abstract ……………………………………………………………………………………Ⅱ 目录 . .......................................................................................................... Ⅲ 第一章 引 言 .......................................................................................... 1 1.1 课题的背景和意义 ..................................................................................... 1 1.2 混合线性模型的发展 ................................................................................. 2 第二章 固定效应的估计 ............................................................................ 3 2.1混合线性模型的形式 . ................................................................................. 3 2.2 固定效应的估计 ......................................................................................... 4 第三章 方差分量的估计 ............................................................................ 6 3.1 方差分析估计 ............................................................................................. 6 3.2极大似然估计 . ............................................................................................. 9 3.3 限制极大似然估计 ................................................................................... 12 3.4 最小范数二次无偏估计 ................................................................... 13 3.5 谱分解估计.............................................................................................. 16 第四章 随机效应值的预测及方差分量估计的检验 ................................. 20 4.1 随机效应值的预测 .......................................................................... 20 4.2 方差分量估计的检验 ............................................................................ 21 第五章 方差分量估计方法的比较 .......................................................... 22

结论 . ........................................................................................................ 24 参考文献 ................................................................................................. 25 致谢 . ........................................................................................................ 26

第一章 引言

1.1 课题的背景和意义

方差分析模型 Y=Xb+e,是大家所熟知的，将模型中各种效应分为二大类，固定效应和随机效应，这便得到了下面将要研究的混合线性模型. 混合线性模型分为平衡数据的混合线性模型（对所有因子的水平组合，重复试验次数相同的模型）和非平衡数据的混合线性模型（对所有因子的水平组合，重复试验次数不相同的模型）.

混合线性模型的统计分析包括二方面，首先是估计各项随机效应的方差分

2量σn ，其次是估计固定效应β，随机效应是不可估计的，但可以进行预测，

另外，可以对方差分量进行检验和区间估计.

混合线性模型中，研究最多的是对于方差分量的估计方法，目前来看，主要方法有方差分析法（Henderson 方法三）、极大似然估计、限制极大似然估计、最小范数二乘估计及谱分解估计，不同方法中还有不同的求解方法，此方面的研究还不尽完善. 对于方差分量的估计方法还有待发展，各种方法都有其不足之处.

对于固定效应的估计方法较少，比如最小二乘估计方法，极大似然估计在得到方差分量估计的同时也有固定效应的估计.

最后Markov 链蒙特卡罗（Markov chain Marte Carlo）分析方法，不但可以获得参数的点估计，还能提供这些参数的相应分布及其特征值.

1.2 混合线性模型的发展

近30年来，关于混合线性模型的参数估计一直是线性模型的最活跃的研究方向之一，在这方面已有一些专著将混合线性模型应用于处理纵向数据分析. 对于固定效应的估计方法为最小二乘估计而. 随机效应是不可估计的但可以预测.

对于随机效应的方差分量估计方法有多种：

1953年Henderson 提出了三种估计非平衡数据方差分量的统计分析方法，这些方法都模仿方差分析的原理. 在这三种方法中，方法三运用最为广泛. 适用于析

因试验设计或套试验设计等正规试验设计的方差分析，但对一些特殊的线性模型，采用此法估计的方差分量可能会产生偏差，70年代以来发展的一系列混合线性模型分析方法，包括最大似然法，限制性最大似然法，最小范数二次无偏估计法等. 可以克服Henderson 方法三的局限性.

自70年代开始，Hartley 和Rao(1967)首先提出应用最大似然法分析混合线性模型的非平衡数据. 但这个方法的估计值受到固定效应的影响，可能导致较严重的有偏估计，为了服这一缺点，Patterson 和Thompson(1971)提出了限制性最大似然法，使似然值不包括固定效应.

Rao （1971）提出最小范数二次无偏估计，比最大似然法和限制性最大似然法更简便和优越，它不需要进行迭代运算，对线性模型也没有正态分布的限定.

最近，王松桂和尹素菊提出了同时估计固定效应和方差分量的一种新方法，称为谱分解估计.

第二章 固定效应的估计

2.1 混合线性模型的形式

混合线性模型的一般形式是：

y =Xβ+U 1ξ1+U 2ξ2+„+U k ξk ，

其中y 是n ⨯1观测向量，X 为n ⨯已知设计阵，β 是p ⨯1非随机的参数向量，称为固定效应. U i 为n ⨯t i 已知设计矩阵，ξi 为t i ⨯1随机向量，称为随机效应. 且E(ξi )=0，i=1，„，k .

通常假设

Cov(ξi )=σi 2I t i ，Cov (ξi ，ξj )=0，i ≠j. 于是我们有

E(y )=Xβ， Cov （y ）=∑σi 2U i U I / ∑(σ2) ，

i =1k

这里σ2=(σ12, , σk 2). σi 2称为方差分量，相应的模型也称为方差分量模型. 混合线性模型之所以在生物、医学、经济、金融和工业设计等领域具有十分广泛的应用，是因为在实际问题中，参与试验的个体是随机抽取的，且我们研究的目标不是这些个体本身的特征，而是它们所在的总体特征，这时把个体效应看作随机的而引入模型，可以大大提高模型的精度.

对于混合线性模型，我们感兴趣的参数分两类：固定效应β 和方差分量

σ2=(σ12, , σk 2)，它们分别包含在均值E(y ) 和协方差Cov （y ）中，处理问题

的方法与固定效应模型相比较复杂一些.

在本章里我们将研究固定效应估计方法，以后几章里详细讲解方差分量的估计方法及检验.

2.2固定效应的估计

在考虑到固定效应的估计时，将模型写为如下形式 y =Xβ+Uξ+e，

其中β为固定效应，ξ为随机效应，且E(ξ)=0，E(e)=0,Cov(ξ，e)=0，并且假设

ξ和e 的协方差阵具有较一般的形式Cov(ξ)=D ≥0,Cov(e)=R>0,于是我们有∑=Cov(y )=UDU+R>0.

暂时视D,R 已知应用最小二乘估计法得到正则方程 X '∑-1X β*=X '∑-1y ,

据此可得到β的广义最小二乘解β*= (X '∑-1X )X '∑-1y . 因此，任意的可估

-

函数c /β的最佳线性无偏估计（BLU 估计）为

c 'β=c '(X /∑-1X )X '∑-1y .

*

-

而实际上D,R 未知，若用它们的估计D , R 代替，即用∑=U D U '+R 代替

^^^^^

∑，便得到c /β的两步估计

c 'β（∑）=c '(X '∑X ) X '∑y .

-

~

^

^-^-1

在假设Cov(ξi )=σi 2I t i ，Cov (ξi ，ξj )=0，i ≠j 下c 'β的两步估计又可变形为

c 'β（σ）=(X '∑(σ) X ) X '∑(σ) -1y .

在以上模型中若假设e, ξ的联合分布关于原点对称，设σ=σ2（y ）是σ2

的一个估计，它是y 的偶函数，且具有变换不变性. 对一切可估函数c 'β，若E （c 'β（σ））存在，则两步估计c 'β（σ）必为c 'β无偏估计. 以上定理说明了在一定条件下c 'β（σ）为c 'β无偏估计.

关于e, ξ分布的假设在许多情况下是满足的. 此外，对方差分量估计的方法

~

~

~

^2^2

-1-

^2

^2

^2

~

^2

^2

不同，所得到的固定效应的估计也不相同，方差分析法、极大似然法、限制极大似然法、最小范数二乘法, 所产生的估计σi 都是y 的偶函数，且是变换不变的，因此，对于混合线性模型的固定效应以上定理给出了一大类两步估计的无偏性. 我们可以看到方差分量的不同估计，往往会产生不同两步估计，而且除最小二乘

^2

估计外，模型还可能存在另外一些简单估计，如Panel 数据下的Between 估计和within 估计，简约估计等.

第三章 方差分量的估计

3.1 方差分析估计

1953年Henderson 提出了三种估计非平衡数据方差分量的统计分析方法，这些方法都模仿方差分析的原理. 在这三种方法中，方法三运用最为广泛. Searle （1968）将Henderson 的分析方法改成矩阵形式表示. 方差分析方法渊源于固定效应模型的方差分析.

我们用下面的例子来说明它的原理和方法 对于平衡单向分类模型

i =μ+αi +e ij ， i =1，„，a ， j =1，„，b,

其中μ为总均值，是固定效应，α1, , αa 为随机效应. 假定所有αi , e ij 都不相

2关，且其均值为0，方差为Var(αi )=σα,Var(e ij )=σe 2. 记y '=（y 11, , y ab ）. 暂

时先把αi 看作因子A 的i 水平A i 的固定效应，有

RSS(μ)=y .. /(ab ) SS μ, (3.1.1) 其自由度为1，对应于α1, , αa 的平方和，即因子A 的平方和 SS A =RSS(μ, α)-RSS(μ)=∑

i

-2

∑(y

j

-

i .

-y .. ) 2， (3.1.2)

-

其自由度为n-1，而残差平方和为 SS e =y y -RSS(μ, α)=∑

/

I

∑(y

J

ij

-y i . ) 2, (3.1.3)

-

其自由度为

a （b-1）. 由以上三式可以推导出总平方和的分解式 y /y =SS μ+SS A +SS e 将平方和除以自由度，得到均方： Q 0=SS μ, Q 1=SS A /(a -1), Q 2=SS e /[a (b -1)],

再按照αi 为随机效应的假设，求出各均方的均值：

2

E (Q 0) =E (SS u ) =ab μ2+b σα+σe 2, 2 E (Q 1) =E (SS A /(a -1)) =b σα+σe 2,

E (Q 2) =E (SS e /[a (b -1)])=σe 2.

可以看到，后两式的右端为方差分量的线性函数，令E (Q 1) =Q 1， E (Q 2) =Q 2，

2便得到关于σα，σe 2的线性方程组 22

⎧⎪b σa +σe =Q 1,

⎨2

σ=Q . ⎪2⎩e

解此方程组得

σe =Q 2， σα=(Q 1-Q 2)/b.

2

它们就是方差分量σα，σe 2的方差分析估计（ANOVA 估计）.

^2^2

从以上求解过程我们可清楚的看到方差分析法的具体思路，我们可以把方差分析法归纳如下：

1.对一个方差分量模型，现将其随机效应看作固定效应，按通常方差分析方法算出各效应对应的平方和（或均方）.

2.求这些平方和（或均方）的均值（此时的随机效应不再看作固定效应），他们是方差分量的线性函数.

3.令这些平方和（或均方）等于它们各自的均值，得到关于方差分量的一个线性方程组，解此方程组便得到方差分量的估计.

至此我们了解到方差分析的一般思想，将此方法用于一般的混合线性模型. 先考虑一个简单模型（方差分类模型）

y =Xβ+U 1ξ1+U 2ξ2+e,

2

'+ σe 2I ∑(σ2) . 其中Cov(ξi )=σi 2I t i ，i=1，2，Cov(y)=σ12U 1U 1'+σ2U 2U 2

2

得到方差分量σ12，σ2，σe 2的线性方程组

2

⎧a 1σ12+(a 2-a 3) σ2+r 2σe 2=SS ξ1, ⎪⎪22

⎨a 2σ2+σe =SS ξ2,

⎪2(n -r -r -r ) σ=SS e . ⎪123e ⎩

2

解此方程组，得到σ12，σ2，σe 2的估计. 它们就是这些方差分量的ANOV A

估计. 估计如下

σ1={SS ξ1-r 2（SS e /(n -r 1-r 2-r 3) ）

-(a 2-a 3) [（SS ξ2-SS e /(n -r 1-r 2-r 3) ）/a 2]}/a 1,

^2

σ2=（SS ξ2-SS e /(n -r 1-r 2-r 3) ）/a 2， σe =SS e /(n -r 1-r 2-r 3) .

对于一般的混合线性模型，文献中也称为拟合常数法，之所以称其为拟合常数法就是因为在构造估计方程时，我们把随机效应看成固定效应，即常数. 对于平衡数据模型，该方法的平方和分解是惟一的，且可根据方差分析得到. 下面将方差分析方法应用于带有交互效应的两向分类混合模型.

^2

^2

y ijk =μ+αi +βj +γij +e ijk ，

i=1,„，a ，j=1，„，b ，k=1，„，e ，

这里μ，并满足通常的假设，即所有的βj ，αi 为固定效应，βj , γij 为随机效应，

2

，Var(γij )=σγ2，γij ，e ijk 都不相关，且具有均值为0，方差为Var(βj )=σβ

Var(e ijk )=σe 2.

首先，暂时视βj , γij 为固定效应，得到总平方和有如下分解： y 'y =SS μ+SS α+SS β+SS γ+SS e ，

再用随机效应的平方和除以各自的自由度，得到均方Q 1=SS β/(b -1) ，

Q 2=SS γ/(a -1)(b -1) ，Q 3=SS e /[ab (c -1)]，求出它们的均值，并令这些均值等

2

于对应的均方，得到关于σβ，σγ2，σe 2的线性方程组

2

⎧ac σβ+c σγ2+σe 2=Q 1, ⎪⎪22

⎨c σγ+σe =Q 2,

⎪2⎪⎩σe =Q 3.

解此方程组，得到方差分量的估计：

σβ=（Q 1-Q 2）/(ac ) ，σγ=（Q 2-Q 3）/c ，

^2

^2

σe =Q 3.

由于方差分析法给出的估计σ 作为一个线性方程组的解，他们未必是正的. 这是方差分析法的一个缺陷. 至于如何对待方差分量的负估计，目前尚无一致的看法. （1）观点认为，若某个σi

研究表明，除非非负性可自然满足，不然无偏性和非负性几乎不可兼得，在近期的研究中，人们往往放弃无偏性，而要求具有另外的优良性，如均方误差最小 .

^2

^2

^2

3.2 极大似然估计

方差分析方法有一定的缺陷，采用此法估计的方差分量可能会产生偏差，

所以，人们开始寻找新的方法. 自70年代开始，Hartley 和Rao(1967)首先提出应用极大似然法分析混合线性模型的非平衡数据. 极大似然方法最早应用在统计中，后来被应用在混合线性模型中. 此方法能同时获得固定效应和方差分量的估计.

考虑一般的混合线性模型

y =Xβ+U 1ξ1+U 2ξ2+„+U k ξk ， （3.2.1） 这里假设ξi ~N(0，σi 2I t i ) ，i=1，„k ，所有ξi 都相互独立，记V i =U i U i '，

σ2=(σ12, , σk 2) '，于是

Cov （y ）=∑σU i U I '=∑σi 2V i ∑(σ2) .

2i i =1

i =1

k

k

我们假设∑(σ2) >0，因此y ~N n （0，∑(σ2) ），所以未知参数β，σ12，„，σk 2的似然函数为

L(β, σ|y )=(2π)

2

-

n

2

1

|∑(σ) |exp{-(y -X β) /∑(σ2) -1y -X β}，

2

2

-

12

取对数，略去常数项及常数倍，得

l(β, σ2|y )=-ln |∑(σ2) |-(y -X β) '∑(σ2) -1(y -X β)

=-ln |∑(σ2) |-tr ∑(σ2) -1(y -X β)(y -X β) '. （3.2.2） 根据如下公式

∂Ax ∂∂A (t )

]. =A , ln |A (t ) |=tr [A (t ) -1

∂t ∂t ∂x

对σ2，β求导可得

∂l

=-tr (V i ∑(σ2) -1) +tr [(∑(σ2) -1V i ∑(σ2) -1)(y -X β)(y -X β) ']，i=1，„，k, 2

∂σi

∂l

=-2X '∑(σ2) -1X β+2X '∑(σ2) -1y . ∂β

令这些导数等于零，得到似然方程

2-12-1

⎧⎪X '∑(σ) X β=X '∑(σ) y ,

⎨（3.2.3） 2-12-12-1

⎪⎩tr (V i ∑(σ) ) =(y -X β) '(∑(σ) V i ∑(σ) )(y -X β).

i=I,„，k . 进行简化变形后似然方程变为

⎧X β=X (X '∑(σ2) -1X ) -X '∑(σ2) -1y ,

⎪k

i =1, „，k . （3.2.4） ⎨2-122-12-1

''tr [V ∑(σ) ]σ=y (I -P ) (∑(σ) V ∑(σ) )(I -P ) y , i j σi σ⎪∑j =1⎩

这就是我们要求的似然方程. 由（2.2.4）得第一方程，任意可估函数c 'β的ML 估计为c 'β(σ) =X (X '∑(σ) X ) X '∑(σ) y ，其中σ为σ2的ML 估计.

^

^2

^2

-1

-

^2

-1

^2

在一般情况下似然方程（2.2.4）没有显式解，即便在有显式解的情形，σ2的解未必是非负的，若为负值，它就没有落在参数空间内，所以并不是ML 估计. 这时一般采取截断法，即取max{σi ,0}作为ML 估计，在没有显式解的情形只能用迭代法求解.

下面详细介绍两种迭代方法，并简单提出其他方法. Anderson等提出一种迭代法是 σi 这里σ

^2(m )

^2

^2(m +1)

=H (σi

^2(m )

) h (y , σ

-1

^2(m )

) ，

为σ的第m 次迭代值σ

2

^2(m )

. 当σ得两次相邻迭代值相差不大时，迭代

^2

停止，这就得到了方差分量的估计. 带入（3.2.4）的第一方程，便可得到固定效应的估计.

另外一种迭代法是由Hartley 和Rao 提出的，其推广形式是 σi

^2(m +1)

=σi

^2(m )

h (y , σ) i 2(m ) , i =1, , k .

^

tr (∑(σ) -1V i )

^2(m )

这个迭代的一个好处是，当初始值为非负时，后面的迭代值永远不会取负值.

另外还有一些迭代方法，如Newton-Raphson 方法，得分方法，以及EM 算法，这里不做详细研究.

下面以两向分类混合模型为例来具体介绍参数估计的过程 对两向分类混合模型

y ij =μ+αi +βj +e ij , i =1, , a , j =1, , b ,

2 这里μ, αi 为固定效应，βj 为随机效应，βj ~N (0,σβ), e ij ~N (0,σe 2) ，且所

有βj e ij 都相互独立，该模型的矩阵形式为

y =X 1μ+X 2α+U β+e ,

用Kronecker 乘积表示设计阵X 1, X 2和U ： X 1=1ab =1a ⊗1b , X 2=1ab =I a ⊗1b ,

U =1a ⊗I b .

/

固定效应的设计阵为X =(X 1 X 2) ，P X =P X 2=I a ⊗J b , 这里J b =1b 1b /b .

--

得到似然方程

X 1μ+X 2α=P X y ,

-1b -

a σβ+σe =∑(y . j -y .. ) 2，

b j =1

^2

^2^

^

b a σβ+σe

^2

^2

+

(a -1) b

σe

^2

=

1(a σβ+σe ) 2+(a -1) b

^4

^2^2

∑(y

i =1

k

-

2

-y ) . j ..

-

σe

∑∑(y

i

j

ij

-y i . -y . j +y .. ) 2

---

.

求得方程组的显式解

---1

σe =(y ij -y i . -y . j +y .. ) 2， ∑∑(a -1) b i j -1b -1^22

σβ=∑(y . j -y .. ) -σe ，

ab j =1a

^2

^2

^2

在此处，σβ可能取负值.

3.3 限制极大似然估计

方差分量的极大似然估计的一个缺陷是在导出方差分量的估计的过程中没

有考虑到固定效应β的估计所引起的自由度的减少. 为此，Patterson 和Thompon 提出的一种修正方法，称为限制极大似然法. 该方法的思想是基于极大似然估计残差，利用极大似然法导出方差分量的估计. 限制极大似然估计的特点是此方法所求方程的解与方差分析法所得的估计相同.

下面我们用两向分类随机模型来说明限制极大似然法的参数估计过程

y ijk =μ+αi +βj +γij +e ijk ,

i =1, , a , j =1, , b , k =1, , e ,

2

这里为μ为总平均，是固定效应，αi , βj , γij 都为随机效应，假设αi ~N (0,σa ) ，

22

βj ~N (0,σβ) ，γij ~N (0,σγ2), e ijk ~N (0,σα) 且都相互独立. 该模型的矩阵形式为

y =X μ+U 1α+U 2β+U 3γ+e ,

这里X =1a ⊗1b ⊗1c ，U 1=I a ⊗1b ⊗1c ，U 2=1a ⊗I b ⊗1c , U 3=I a ⊗I b ⊗1c ，方差阵

2222

'''''∑(σ2) =σαI α⊗11b b ⊗11c c +σβ1a 1a ⊗I b ⊗11c c +σγI a ⊗I b ⊗11c c +σe I abc .

得到方程组如下

2

bc σα+c σγ2+σe 2=Q 1,

2ac σβ+c σγ2+σe 2=Q 2,

c σγ+σ=Q 3,

2

2

e

σe 2=Q 4,

这里我们记

Q 1=∑(y i .. -y ... ) 2/(a -1),

i =1b a

-

-

Q 2=∑(y . j . -y ...) 2/ac (b -1),

j =1a

--

Q 3=∑∑(y ij . -y i .. -y . j . +y ... ) 2/(a -1)(b -1),

i =1j =1

b

----

Q 4=∑∑∑(y ijk -y ij . ) 2/ab (c -1).

i =1j =1k =1

a b c

-

这与方差分析法所得到的线性方程组相同，因此限制极大似然方程的解与方差分析估计相同. 对平衡数据的混合效应模型，这种现象通常成立.

3.4 最小范数二次无偏估计

方差分量的最小范数二次无偏估计始于Townsend 和Harville ，但C.R.Rao 的工作奠定了重要的基础，此方法与前面介绍的三种方差截然不同. 前三种方法都是先按已有的一定方程式去求估计，至于所得估计有何性质，事先并不知道. 而最小范数二次无偏估计的基本思想是对模型的误差和随机效应的分布没有要求，而是先提出估计应具有的性质，然后把为满足这些性质所加的条件提成一个极值问题，即所谓最小迹问题. 解所得的最小迹问题，便得到所要的估计. Rao 这个方法的缺点是，对N 个不同的人，由于各自选用的方差分量的先验值不同，就会得到N 个不同的估计，具有一定程度的主观随意性.

考虑最一般形式的方差分量模型

y =X β+U 1ξ1+ +U k ξk ， （3.4.1）

这里X n ⨯p ，U i , n ⨯t 为已知设计矩阵，β为p ⨯1固定效应向量，ξi 为t i ⨯1随机效应向量，满足E （ξi ）=0，Cov(ξi )=σi 2I t i , ξi 都不相关，若记

', , ξk ') ， U =(U 1 U 2 U K ) ，ξ'=(ξ1', ξ2则模型（3.4.1）可改为

y =X β+U ξ, E (y ) =X β, Cov (y ) =∑σi 2V i ∑,

i =1k

2

我们的基本目的是估计方差分量σ12, , σk 及其线性函数ϕ=c 'σ2，这里2

''σ2=(σ12, , σk ) ，c =(c 1, c k ) .

首先，看所求的估计量应有的一些性质，因为要估计的参数是方差，所以考虑二次型估计y /Ay ，A 为对称阵，要求这个估计具有下述性质.

(1)不变性 即估计y /Ay 关于参数β具有不变性. (2)无偏性

(3)最小范数准则 欲使y 'Ay 为一个好的估计，那么对一切ξ，矩阵U 'AU 与

∆在某种意义下相差很小. 用矩阵范数||U 'AU -∆||来度量U /AU 与∆相差大小，

则应该选择A 极小化范数||U 'AU -∆||. 若线性函数ϕ=c 'σ2的估计y 'Ay 满足 AX =0,

tr (AV i ) =c i ， i =1, , k , 且使范数||U 'AU -∆||达到极小，则称为最小范数二次无偏估计.

2222

这里采用加权欧氏范数，令权矩阵 W =diag {σ0,1其中为I t 1, , σ0, I }σσi k t k 0, i

的一个预先指定值（先验值），因此W 也就是Cov (ξ) 的一个预先指定阵（先验阵）.

12

12

定义F =W (U 'AU -∆) W , 则加权欧氏范数

||U 'AU -∆||=tr (F 'F ) =tr [W (U 'AU -∆) W (U 'AU -∆) W ]

12

12

=tr (W 'AUWU 'AUW ) -2tr (W 'AUW ∆W ) +tr (∆W ) 2. 利用无偏性，上式第二项

tr (W 'AUW ∆W ) =tr (U 'AUWU 'A ∆W ) =∑

i =1

2

再记V w =∑σ0, i V i ，于是

i =1k

1

2121212

1212

k

4c i σ0, i

t i

(AV i ) =∑

i =1

k

4

c i 2σ0, i

t i

=tr (∆W ) 2.

||U 'AU -∆||=tr (AV w ) 2-tr (∆W ) 2.

这样，对加权欧氏范数求ϕ=c /σ2的最小范数二次无偏估计的问题，归结为求下述极值的解

min tr (AV w ) 2

⎧AX =0

⎨ （3.4.2）

tr (AV ) =c , i =1, , k . i i ⎩它的目标函数是矩阵的迹，称（3.4.2）为最小迹问题.

下面考虑极值问题的解是否存在，有如下定理1 极值问题（3.4.2）的解为

A =B w (∑λi V i ) B w ,

*

i =1k

其中

-1-1-1-1

B w =V w ， -V w X (X 'V w X ) -X 'V w

且λi , i =1, , k 为方程组

∑tr (B w V i B w V j ) λi =c j ，j =1, , k

i =1

2

的解，这里V i =U i U i ，V w =∑σ0, i V i .

i =1k

k

证明定理1等价于证明如下定理2 极值问题（3.4.2）的解为

A =N (∑λi V i ) N ,

*

i =1k

其中N =I -X (X 'X ) -X ，λi , i =1, , k 为方程组

∑tr (NV i NV j ) λi =c j ，j =1, , k

i =1k

的解.

尽管方差分量的估计值取决于人为选择的先验值，只要这些先验值不依赖于试验数据，MINQUE 估计量仍是无偏的，如果用估计值替代先验值重新估计，便可获得新的估计，重复这一过程，直到新的估计非常接近旧的估计为止，这种迭代的估计方法就是上一节介绍的限制极大似然（REML ）方法，其估计结果就是限制极大似然估计，所以在多数情况下，REML 估计和MINQUE 估计是比较接近的.

最后说一下先验值选取的问题，先验值的选择可以凭经验，或根据以往的分析结果. 最简单的方法是取残差效应的先验值为1，其它所有的先验值为0，这种方法称为MINQUE （0）法，所获得的方法分量的估计量是MINQUE （0）估计量. 这种取值的方法的优点是可以省去大矩阵的求逆计算，算法简单，但是由MINQUE （0）方法估计的方差分量的抽样方差往往较大. 另一种简便易算的方法是设所有的先验值为1，这种方法称为MINQUE （1）法，所获得的方差分量是MINQUE （1）估计量.

3.5 谱分解估计

下面我们介绍一种估计方差分量的新方法，即谱分解估计，首先我们引入谱分解的有关知识.

协方差阵的谱分解有多方面的应用，一方面，方差分量的谱分解估计就是以协方差阵∑(σ2) 的谱分解式为基础构造的；另一方面，从∑(σ2) 的谱分解式，可以直接获得∑(σ2) 的逆矩阵和行列式，而这些恰是利用极大似然方法求解方差分量的极大似然估计所必需的.

目前有两种谱分解算法，一种是由Smith 和Hocking 提出的，他们基于完全设计的随机效应模型，给出了谱分解的公式，并将它推广到一般的混合模型的情

形，即将公式中未出现的随机项的方差分量用零取代. 另一种算法是由Searle 和Henderson 提出的，该算法没有利用模型随机效应设计阵U 1, , U k 之间的某些约束关系，而是基于对∑(σ2) 扩充很多项的方法，将∑(σ2) 改写成

∑θj p , j 1(J n p p ⊗ ⊗J n j 11) ，这里（j p , j 1）取遍从（0，0，„0）到（1，1，„1）

j

的所有2p 个二进制向量，即添加（2p -k-1）个θj p , j 1=0的项到∑(σ2) 中，给出了一般形式的谱分解公式. 最近吴密霞，王松桂就一类模型研究了协方差矩阵的谱分解，互异特征值个数的确定及求解问题，将谱分解的结果应用到构造模型参数估计，提出了谱分解估计. 史建红在他的博士论文中把结果推广到一般平衡混合线性模型.

下面我们详细介绍王松桂和尹素菊提出的同时估计固定效应和方差分量的一种新方法，即谱分解估计.

谱分解估计的基本思想是：首先对协方差阵进行谱分解，然后利用谱分解得到的主幂等阵对原模型进行适当的线性变换，获得若干个新的奇异线性模型，这些新模型的特点是它的固定效应与原模型相同，但新模型的协方差阵除了一个因子（这个因子是原模型协方差阵的一个特征值）外，不含未知的方差分量，利用最小二乘统一理论，对每个新模型可以得到固定效应和特征值的一个估计，由于在常见情形下协方差阵的特征值是方差分量的线性函数，因此通过解线性方程组可以获得方差分量的估计.

新方法的突出特点是能同时给出固定效应和方差分量的估计，前者是线性的，后者是二次的，且相互独立. 对于固定效应可以获得若干个谱分解估计，它们都是具有一些好的统计性质的线性估计. 因此，利用这些估计可以对模型做进一步假设检验，区间估计，以及模型诊断等一系列统计推断.

对于任一混合线性模型

y =X β+∑U i ξi +e , （3.5.1）

i =1k

如果y 的协方差阵∑(σ2) 有如下谱分解

∑(σ) =∑λi M i , (3.5.2)

2

i =1q

这里λi ，i=1,„，q 是∑(σ2) 的所有互异非零特征根，它们是σ2线性函数，M i 是特征根λi 对应的主幂等阵（即M =M i , M i M j =0. i ≠j , ∑M i =I ），且独立于未知

2

i

i =1q

参数. 我们分别用幂等阵M i (i=1, „，q) 左乘模型(2.5.1).于是得到变换后的q 个新模型

M i y =M i X β+εi ， εi ~（0，λi M i ），i=1,„，q. （3.5.3） 这些模型的特点是，模型协方差阵除λi 之外，独立于未知参数. 因为M i 是奇异阵，我们可以应用最小二乘统一理论获得原模型参数（固定效应和方差分量）的估计，称为谱分解估计.

下面我们简略介绍一下最小二乘统一理论， 最小二乘统一理论应用于线性模型的参数估计，由著名统计学家Rao 应用推广. 对于线性模型 y =X β+e ,

E (e ) =0. Cov (e ) =σ2∑，

如果|∑|=0，则称该模型为奇异线性模型. 对于奇异线性模型，因为∑-1不存在，于是一般的最小二乘法中用到的Q (β) 无定义，Rao 成功解决了这个问题. 关键是寻找一个新矩阵T ，能够充当∑-1所负担的作用.

T=∑+XUX ', 其中U ≥0,rk(T)=rk(∑ X ), 然后定义

Q (β) =(y -X β) 'T '(y -X β) ,

用最小化Q (β) 求出最小值点

β*=(X 'T -X ) -X 'T -y ,

还能证明对任一可估函数c 'β, c 'β*为其BLU 估计. 这个结论既适用于设计阵X 列满秩或列降秩的情形，又适用于∑奇异阵或非奇异阵的情形。正是由于这个原因，通常把这个结果称为最小二乘统一理论.

对于变换后的新模型（3.5.3）由最小二乘统一理论，便得可估函数c 'β及其方差参数λi 的一个无偏估

c 'β=c '(X 'M i X ) -X 'M i y ，

λi =y '(M i -M i X (X 'M i X ) -XM i ) y /m i ，

2

'i =1, , k +1，这里m i =r i -rk (M i X ) . 得到σ2=(σ12, , σk +1) 的无偏估计

~2

^

~(i )

^

σ=H λ.

-1

c β和σ是模型参数的谱分解估计.

T

~(i ) ^2

第四章 随机效应值的预测及方差分量估计的检验

4.1 随机效应值的预测

在混合线性模型的分析中，通过方差分量的估计可以推断随机效应的变异性，在有些试验中，试验者不但需要估计方差分量，而且需要推断模型中某些随机效应的值，混合线性模型中的随机效应是不可估计的，但却是可预测的.

对于以下模型

y =X β+U ξ+e ，E (e ) =0, Cov (e ) =R >0, E (ξ) =0, Cov (ξ) =D ≥0，

这里假设D , R 都是已知的，如果它们含有未知参数，在实用中用它们的估计代替，这样在随机因素的方差是已知条件下，可以获得随机效应的最佳线性无偏预测（BLUP ）. 得到的估计公式为

ξ=DU '(UDU '+R ) -1(y -X β*) .

按上式预测的随机效应是相依变量的线性函数，因而是多元正态随机变量，其期望值是零. 随机效应的BLUP 预测值是参数的无偏预测，并且在所有的线性无偏预测中具有最小方差，而实际上在运用混合线性模型分析具体试验数据时，方差分量是需要被估计的，因而并不可能获得随机效应值的最佳线性无偏预测，在实际应用时常采用方差分量的估计代替其参数.

^

4.2 方差分量估计的检验

考虑混合线性模型

y =X β+U 1ξ1+U 2ξ2+e ,

22我们假设ξ1~N (0,σ1I s ), ξ2~N (0,σ2I ), e ~N (0,σ) 2q 0I n ，且它们彼此独立，其中

22

'+σ2s , q 分别为已知阵U 1和U 2的列数. 于是Cov (y ) =σ12U 1U 1U 2U '+σ2e I n . 在前面2

基于拟合常数的思想，我们给出了此模型方差分量σ12，σ2和σe 2的估计，现在我

们用同样的方法来构造随机效应ξ2是否存在的检验问题，即

22 H 0:σ2=0H 1:σ2≠0,

的检验统计量.

将以上模型中随机效应ξ1和ξ2暂视为固定效应，模型拟合之后的残差平方和为

SS e =y /y -RSS(β, ξ1, ξ2)=y '(I n -P (X :U 1:U 2) ) y .

2

若σ2=0，则模型变为y =X β+U 1ξ1+e ，同样暂视其中的随机效应ξ1为固定效应

时，模型拟合后的残差平方和变为

SS e 0=y /y -RSS(β, ξ1)=y '(I n -P (X :U 1) ) y .

2直观上，当σ2=0时，模型拟合后的残差平方和SS e 0和SS e 应很接近，即SS e 0-SS e 2相对于SS e 应很小，不然，我们就认为随机效应ξ2作用显著，即接受σ2≠0. 令

F =

(SS e 0-SS e ) /t 2

，

SS e /t 1

这里t 1=n -rk (X :U 1:U 2), t 2=rk (X :U 1:U 2) -rk (X :U 1) .

此外还有一种方法不但可以估计方差分量，还能有效的检测方差的显著性（即对方差分量估计进行检验）. 这种方法叫做刀切法（Miller,1974; 朱军，1992），现以试验处理重复抽样为例，介绍估计遗传参数刀切法估计的方法. 对于遗传参数（方差分量或协方差分量），分析所有观察值所得到的估计. 当从资料中剔除第i 种试验处理世代均值的r 次观察值以后，重新估计而得到新的估计，采用此法从完整的数据中每次剔除r 次观察值，可以获得g 个不同的估计，刀切法估计出的方差分量具有t 分布，从而可以用t 检验对遗传参数作统计检验.

Zhu(1989)用蒙特卡罗模拟比较了估计方差分量估计的方差及参数检验的效益，采用的两种统计分析方法分别是：MINQUE 估计方差的渐进估计法和刀切法抽样估计法. 在遗传分析中一些遗传参数估计值的抽样方差也可以用刀切法算得，并可用t-检验作统计检验.

第五章 方差分量估计的比较

在前面的研究中，我们讲解了混合线性模型的参数（固定效应和方差分量）的估计，并简单介绍了随机效应的预测及方差分量估计的检验. 重点是方差分量的估计，接下来我们对方差分量的各种估计方法进行比较，来给出各自的优缺点.

目前，最常用的方差分量的估计有方差分析估计（ANOV AE ），极大似然估计（MLE ），限制极大似然估计(RMLE)，最小范数二次无偏估计(MINQUE)及谱分解估计(SDE) .

方差分析估计：用一句话来总结它的思想是，令各种平方和等于各自平方和的期望值，得到一个关于方差分量的线性方程组，求解这个方程组便得到方差分量的估计. 对于平衡数据混合模型，这种估计方法的优点是无偏性和方差最小性，但对于非平衡模型，这种方法具有不惟一性，而且计算上很麻烦，最主要的缺陷是有时能给出方差分量的负估计值.

极大似然估计：理论上这种方法能同时获得固定效应和方差分量的估计，但一般来说，我们不能获得方差分量的估计的表达式，只有在平衡数据混合模型下，才能得到显式解，再没有显式解的情形下，只能进行计算机迭代，这时需要有较好的迭代算法. 能够处理大维数的矩阵非线性方程，并且导致整体最大值，同时最大值还必需满足非负要求.

限制极大似然估计:此方法是对极大似然估计的一种完善，极大似然估计方法没有考虑估计固定效应所引起的自由度的减少，限制极大似然估计就是基于这个缺陷所做的改进，这种方法与方差分析估计结果相同.

最小范数二次无偏估计：这种估计方法与前几种相比有很大的不同之处，它在给定了方差分量的一组先验值后，只需求解一个线性方程组，不需要进行繁琐的迭代. 缺点是对N 个不同的人，由于各自选用的方差分量的先验值不同，会得到N 个不同的估计，具有一定程度的主观随意性.

谱分解估计：此估计方法是最近提出的新方法，能同时给出固定效应和方差分量的估计，对于固定效应可以获得若干个谱分解估计，都是具有良好统计性质的线性无偏估计，方差分量的估计是二次不变无偏估计，并且SD 估计与ANOVA

估计一样具有显式解，但SD 估计也不能保证估计的非负性. 已有学者证明，对于一些模型方差分量SD 估计和ANOVA 估计相等，所以SD 估计具有ANOVA 估计的一些优良性质.

结论

本文主要研究了混合线性模型参数的估计及检验，混合线性模型是既包含固定效应又包含随机效应的一类线性模型，此类模型在实际中有多方面的应用，我们既要对固定效应进行估计，又要对随机效应的方差分量进行估计，而随机效应是不可估计的，但实际中我们需要对其进行预测. 此外，我们也要对所估计的方差分量进行检验，通过假设检验，看它所对应的随机效应是否显著有影响，从而判断随机效应在实际中对我们所要研究的变量是否有显著影响.

在本文中我们仅研究了一种估计固定效应的方法，即最小二乘估计法，此方法简单易行，只是在估计过程中我们先暂时假设了随机效应的协方差阵已知，在得到固定效应的估计后，我们还要对协方差阵进行估计，因为实际上它们是未知的，从而得到固定效应的两步估计.

对于方差分量，我们研究了五种方法，方差分析法、极大似然法、限制极大似然法、最小范数二次估计法及谱分解估计方法，其中方差分析法及谱分解估计法都能得到显式解，但不能保证估计的非负性，极大似然法和限制极大似然法在一般情况下需要进行迭代，无法得到显式解，最小范数二次估计需要事先给定一组方差分量的先验值，有很大的主观随机性.

方差分量估计的检验，我们主要研究了一类较简单混合模型的检验问题，构造了一个F 检验，其它检验未做详细说明，对于检验方面的理论研究还不尽完善，还有许多问题需要解决.

参 考 文 献

[1] 王松桂 尹素菊 史建红 吴密霞 线性模型引论 科学出版社 2000 [2] 朱军 线性模型分析原理 科学出版社 2004

[3] 王松桂 邓永旭 方差分量的改进估计 应用数学学报 1999，22：115~125

[4] 吴密霞 王松桂 线性混合模型中固定效应和方差分量同时最优 中国科学 A 辑，2004，34（3）：373~384

[5] 王松桂 尹素菊 线性混合模型参数的一种新估计 中国科学 A 辑，2002，32（5）435~443

[6]SHI JIANHONG WANG SONGGUI Some properties of Spectral Dcomposition Estimate of Variance Components. 应用概率统计，2004，287~294

[7]马铁丰 王松桂 线性混合模型方差分量的检验 应用数学学报，2007，22（4）：443~440

致谢

衷心感谢导师林鹏老师对本人的精心指导，他的言传身教是我受益匪浅，导

师从阅读文献、论文选题等方面提出了宝贵意见，在论文的写作过程中，给了多方面的指导，在百忙之中对论文进行了详细的批阅和审查.

感谢统计系的全体老师和同学的关心和支持！感谢所有关心和帮助过我的人们！