林业科学  2013, Vol. 49 Issue (3): 1-8   PDF    
DOI: 10.11707/j.1001-7488.20130301
0

文章信息

童春发, 杨立伟, 蒋安纳, 施季森
Tong Chunfa, Yang Liwei, Jiang Anna, Shi Jisen
不平衡巢式设计遗传模型分析
Analyses of Genetic Models for Unbalanced Nested Design
林业科学, 2013, 49(3): 1-8
Scientia Silvae Sinicae, 2013, 49(3): 1-8.
DOI: 10.11707/j.1001-7488.20130301

文章历史

收稿日期：2012-10-16
修回日期：2012-12-30

作者相关文章

童春发
杨立伟
蒋安纳
施季森

引用本文
童春发, 杨立伟, 蒋安纳, 施季森. 2013. 不平衡巢式设计遗传模型分析[J]. 林业科学, 49(3): 1-8.
Tong Chunfa, Yang Liwei, Jiang Anna, Shi Jisen. 2013. Analyses of Genetic Models for Unbalanced Nested Design. Scientia Silvae Sinicae, 49(3): 1-8.  DOI: 10.11707/j.1001-7488.20130301
不平衡巢式设计遗传模型分析
童春发, 杨立伟, 蒋安纳, 施季森    
1. 南京林业大学林木遗传与生物技术省部共建教育部重点实验室 南京 210037;
2. 南京林业大学信息科学技术学院 南京 210037
收稿日期：2012-10-16；修回日期：2012-12-30
基金项目：国家自然科学基金面上项目(31070601);江苏省高校优势学科建设工程资助项目(PAPD)。
通讯作者：施季森
摘要：虽然巢式试验设计在动植物育种中得到广泛应用,但是在不平衡数据条件下有关遗传参数估计和假设检验统计量的计算还存在着很多问题。文中针对不平衡数据条件下的巢式试验设计,使用固定效应模型和随机效应模型估计相关的遗传参数并进行相应的统计假设检验。对于固定效应模型,使用约束线性模型方法推导出亲本配合力估计及亲本间配合力假设检验统计量的计算公式。对于随机效应模型,采用混合线性模型中的方差分析法,推导出方差分量估计的计算公式,并给出方差分量估计标准误以及方差分量假设检验统计量的计算方法,进而给出遗传力计算公式及其标准误的近似计算方法。最后,用VC ++编写多种形式巢式设计遗传模型的各种遗传参数估计和假设检验统计量计算的Windows应用软件,供林木遗传育种工作者使用。
关键词： 巢式设计    一般配合力    方差分量    遗传力    
Analyses of Genetic Models for Unbalanced Nested Design
Tong Chunfa, Yang Liwei, Jiang Anna, Shi Jisen    
1. Key Laboratory of Forest Genetics & Biotechnology of Ministry of Education, Nanjing Forestry University Nanjing 210037;
2. College of Information Science & Technology, Nanjing Forestry University Nanjing 210037
Abstract: Nested design has been widely used in animal and plant breeding, but there still exist some problems in calculating the estimates of genetic parameters and the statistics of some hypothesis tests with unbalanced data. In this paper, both fixed and random linear models were applied to estimate the genetic parameters and to perform tests of corresponding hypotheses with unbalanced data from the nested experiment design. For the fixed-effect model, we applied the methods of restricted linear model to estimate the general combining abilities of parents and to calculate the statistics for testing whether or not the general combining abilities of any two parents are different. For the random effect model, we used the theory of variance component model to derive the formulae for computing the variance components, their standard errors, and the statistics for hypothesis tests. Furthermore, formulae for calculating the heritabilities and estimating their standard errors were also given. Based on these results, a Windows-based software was developed to implement the calculations of all the genetic parameter estimates and the statistics for various hypotheses, and it would be a useful tool for forest breeders.
Key words: nested design    general combining ability    variance components    heritability    

Comstock 等(1948; 1952)提出了3种北卡罗林纳(North Carolina)遗传交配设计，分别称为 NCⅠ、NCⅡ和NC Ⅲ。 其中 NCⅠ设计又称为巢式设计(nested design)、A / B 或 B / A 设计，其特点是从群体中随机选择若干个个体作为父本，然后对于每一个父本再从群体中选择若干个不同的母本与之交配。巢式设计由于其设计简单并能产生大量的子代，已在动物和农作物育种中有较多的应用(Smith et al., 1998; Pereira et al., 2001; Goncalves et al., 2009; 于连洋等，2010)。

在林木育种研究中，国际上巢式设计应用最好的例子是火炬松(Pinus taeda)(Balocchi et al., 1993; Gwaze et al., 2002)，国内在杉木(Cunninghamia lanceolata)和马尾松(Pinus massoniana)等树种中也有不少研究报道(黄金桃，1996; 梁一池，1999; 谭小梅等，2011)。虽然巢氏设计在遗传育种中的应用由来已久，但是由于其遗传分析模型本身具有多种形式以及试验数据常具有不平衡性，使得巢式设计遗传模型中一些参数估计和假设检验存在一定的困难。

与其他数量遗传模型一样，由于研究目的不同，同一个巢式设计遗传模型会被当作固定效应模型和随机效应模型来处理。固定效应模型，除了随机误差项是随机效应外其他变量都是非随机效应，可用于估算亲本的配合力等遗传参数; 而随机效应模型，除了随机误差外还有其他变量是随机效应，一般用来估计遗传方差分量进而估算遗传力。对于固定遗传效应模型，可用成熟的线性模型理论来分析(Searle，1997; 王松桂等，2007; Rao et al., 2008)，但是在大多数遗传模型中固定效应要满足一定的约束条件，这给研究者运用通用软件计算带来了不少困难。为此，在双列杂交分析中，Wu和Matheson(2000)提出用消除约束条件后的最小二乘法来计算一般配合力和特殊配合力，并提供了SAS 计算程序，这种方法的主要局限性在于研究者需要根据实际模型和数据来修改SAS 代码。最近，Searle，1997; Tong 等(2012)根据约束线性模型方法编写了Windows 软件GSCA 来计算分析双列杂交模型，该软件界面友好、使用方便，能解决大多数双列杂交配合力计算问题。

至于随机遗传效应模型，经典的分析方法是针对平衡试验数据的所谓方差分析法。该方差分析法不同于广泛应用的单因素或多因素方差分析中的那种方法，而是由Henderson(1953)提出的估计遗传模型中方差分量的一种统计方法。对于特定的遗传模型，这种方差分量分析方法在很多文献中被总结为一张方差分析表(Wright，1976; Namkoong，1979;Williams et al., 2002; 续九如，2006)，使得不少研究者忽略了它的线性随机效应模型属性。方差分析表上列有各个遗传因素的均方差及其数学期望，每个均方差数学期望是待求方差分量的线性函数。计算时，先根据试验数据计算均方差，然后令均方差等于各自的数学期望就可以解出方差分量。根据方差分析表计算方差分量的优点在于计算简单，使用计算器甚至手工计算都可以完成，但其缺点是显而易见的，当试验数据不平衡或模型稍微有所改变时，就不能再使用同样的方差分析表来计算了，这给大多数育种研究者带来了许多困惑。为了克服这个缺陷，不少研究者编写了SAS 程序用于多种形式的双列杂交遗传方差分量计算(Zhang et al., 1997; Johnson et al., 1998; Wu et al., 2000; Wu et al., 2001; Xiang et al., 2001; Murray et al., 2003; Zhang et al., 2005)。针对林木半同胞子代测定试验，童春发等(2010)编写了Windows 软件HalfSibSS 用于计算平衡或不平衡数据条件下的方差分量，同时可以计算遗传力等遗传参数，并提供相关估计的标准误和假设检验统计量。

到目前为止，很少有易于使用而且功能强大的程序或软件用于分析不平衡巢式设计的固定效应模型和随机效应模型。唐守正等(2009)编写的统计和生物数学模型计算软件Forstat 包含了主要的林木遗传模型计算问题，其Windows 界面友好，便于使用。该软件虽然能计算不平衡巢式设计模型的方差分量，但对方差分量估计的标准误、假设检验统计量、遗传相关系数以及家系遗传力等计算存在局限性。本文针对不平衡数据条件下巢式试验设计的统计分析模型，采用现代线性统计模型理论(王松桂等，2007; Rao et al., 2008)，给出遗传参数估计和相关假设检验统计量的计算公式，编写相应的Windows 软件，供广大林木育种工作者自由使用。对于固定效应模型，使用约束线性模型方法推导出亲本配合力估计及亲本间配合力假设检验统计量的计算公式。对于随机效应模型，采用混合线性模型中的方差分析法，推导出方差分量估计的计算公式，并给出方差分量估计标准误以及方差分量假设检验统计量的计算方法，进而给出亲本遗传力计算公式及其标准误的近似计算方法。

1 统计模型 1.1 固定效应模型

考虑b 个随机区组r 次重复m 个父本分别与若干个不同的母本交配的不平衡巢式试验设计统计模型，

yijklt = μ + Bi + Ril + Mij + Fijk + FRijkl + eijklt， (1)
式中，yijklt为第i 区组第l 个重复中第j 个父本与第k个母本交配子代的第t 个个体观察值; μ 为总平均值;Bi为第i 个区组效应(i = 1，…，b); Ril为第i 个区组第l 个重复的效应(l = 1，…，r); Mij为第i 区组第j个父本的效应(j = 1，…，m); Fijk为第i 个区组里第j个父本与第k 个母本交配下的母本效应(k = 1，…，fij); FRijkl为Fijk与第l 个重复的互作效应; eijklt为随机误差效应(t = 1，…，n ijkl)，假定其均值为0 方差为σ2。

上述模型中除了eijklt为随机效应外，其他效应均为固定效应，其目的在于估算父本的一般配合力，并检验父本配合力间的差异是否显著。而且这些固定效应应该满足如下约束条件:

$\left\{ \begin{align} & \sum\limits_{i=1}^{b}{{{B}_{i}}=0} \\ & \sum\limits_{l=1}^{r}{{{R}_{il}}=0,\sum\limits_{j=1}^{m}{{{M}_{ij}}=0\left( i=1,\cdots ,b \right)}} \\ & \sum\limits_{k=1}^{{{f}_{ij}}}{{{F}_{jik}}=0\left( i=1,\cdots ,b;j=1,\cdots ,m \right)} \\ & \sum\limits_{l=1}^{r}{F{{R}_{ijkl}}=0\left( i=1,\cdots ,b;j=1,\cdots ,m;k=1,\cdots ,{{f}_{ij}} \right)} \\ & \sum\limits_{j=1}^{m}{\sum\limits_{k=1}^{{{f}_{ij}}}{F{{R}_{ijkl}}=0\left( i=1,\cdots ,b;l=1,\cdots r \right)}} \\ \end{align} \right.$ (2)

模型的不平衡性主要体现在2 个方面: 一是每个父本交配的母本数可能不同，二是小区内个体数由于多种原因可能有缺失。

记所有p 个固定效应所组成的未知参数向量为β，那么固定效应模型(1)及其约束条件(2)可用矩阵形式表示为

$\left\{ \begin{align} & y=X\beta +\varepsilon \\ & L\beta =0 \\ \end{align} \right.$ (3)
式中，y 为n × 1 的观测向量，X 为n × p 设计矩阵，L 为一个列数为p 的矩阵，ε 为n × 1 的随机误差向量。设S(A)表示由矩阵A 的列向量张成的向量空间，R(A)表示矩阵A 的秩。在(3)式中，矩阵X和L 若满足1)S(X')∩ S(L')= { 0}和2)R(X)+R(L)= p，那么矩阵L 或(2)式称为满足边界条 件。一般情况下，选择矩阵L 使其满足边界条件，其时未知参数向量β 的最小二乘估计和σ2 的无偏估计为
$\left\{ \begin{align} & \hat{\beta }={{\left( {X}'X+{L}'L \right)}^{-1}}{X}'y \\ & {{{\hat{\sigma }}}^{2}}=SSE/\left( n-s \right) \\ \end{align} \right.$ (4)
式中，SSE = y'y -${\hat \beta }$'X'y 为残差平方和，s = R(X)。此外，还可以得到${\hat \beta }$的协方差矩阵为
cov(${\hat \beta }$)= σ2(X'X + L'L)-1X'X(X'X + L'L)-1。 (5)

相关假设检验的零假设可表示为

H0: Hβ = 0， (6)
这里H 是一个h × p 矩阵。例如，若要检验区组效应是否显著为0，则零假设为
H0: B1 = B2 = … = Bb = 0，
其矩阵表达式则为
${{H}_{0}}:\left( \begin{array}{*{35}{l}} 0 & 1 & 0 & \cdots & 0 & 0 & \cdots & 0 \\ 0 & 0 & 1 & \cdots & 0 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots & \cdots & \vdots \\ 0 & 0 & 0 & \cdots & 1 & 0 & \cdots & 0 \\ \end{array} \right)\beta =0.$
在零假设(6)下，模型(3)变为
$\left\{ \begin{align} & y=X\beta +\varepsilon \\ & C\beta =0 \\ \end{align} \right.$ (7)
式中，C =$\left( \begin{matrix} L \\ H \\ \end{matrix} \right)$。此时，矩阵C 一般不满足边界条件，但仍然可以给出β 的最小二乘估计为
(8)
${\hat \beta }$H =(T -1 - T -1C'Q- CT -1)X'y，(8) 这里T = X'X + L'L + H'H，Q = CT -1C'，而Q- 为Q 的广义逆矩阵。进一步可以推导出检验统计量(Wang et al., 1994)为
$\begin{align} & F=\frac{\left( SSH-SSE \right)/\left( d{{f}_{1}}-d{{f}_{2}} \right)}{SSE/\left( n-d{{f}_{1}} \right)}\tilde{\ } \\ & F\left( d{{f}_{1}}-d{{f}_{2}}n-d{{f}_{1}} \right) \\ \end{align}$ (9)
式中，SSH = y'y - ${\hat \beta }$'HX'y 为模型(7)的残差平方和，df1 = R(X)，df2 = p - R$\left( \begin{matrix} L \\ H \\ \end{matrix} \right)$ 。

1.2 随机效应模型

如果利用巢式设计估算亲本的方差分量进而估 算遗传力，那么参与交配的亲本应视为群体中的一 个抽样，模型(1)中与亲本有关的效应，即Mij、Fijk和FRijkl应看作是随机变量，此时模型(1)变成了随机效应模型或混合线性模型，其矩阵表达式可表示为

y = Xβ + U1M + U2F + U3FR + ε， (10)
式中，β =(μ，B1，…，Bb，R11，…，Rbr)'，M =(M11，…，Mbm)'，F =(F111，…，Fbmfbm)'，FR =(FR1111，…，FRbmfbmr)'，X，U1，U2和U3分别为对应效应向量的系数矩阵。

假定随机效应Mij、Fijk、FRijkl和eijklt的均值都为0，方差分别为σM2、σF2、σFR2和σε2，并且两两之间相互独立，可以使用方差分析法(ANOVA)求这些方差分量的估计。设D1 = I - X(X'X)- X'，D2 = I -(X U1)［(X U1)'(X U1)］-(X U1)'，D3 = I -(X U1 U2)［(X U1 U2)'(X U1 U2)］-(X U1 U2)'，D4 = I -(X U1 U2 U3)［(X U1 U2 U3)'(X U1 U2 U3)］-(X U1 U2 U3)'，A0 = I - D1，A1 = D1 - D2，A2 =D2 - D3，A3 = D3 - D4，A4 = D4，那么表型性状的平方和可以分解为

y'y = y' A0 y + ESSM + ESSF + ESSFR +ESSε， (11)
式中，ESSM = y' A1 y 为父本间平方和，ESSF =y' A2 y 为父本内母本间平方和，ESSFR = y' A3 y 为 重复× 母本间平方和，ESSε = y' A4 y 为误差平方 和。这些平方和的数学期望分别为
$\left\{ \begin{align} & E\left( ES{{S}_{M}} \right)=Tr\left( {{A}_{1}}\sum \right)=Tr\left( {{A}_{1}} \right)\sigma _{\varepsilon }^{2}+Tr(U_{3}^{'}{{A}_{1}}{{U}_{3}})\sigma _{FR}^{2}+ \\ & Tr(U_{2}^{'}{{A}_{1}}{{U}_{2}})\sigma _{F}^{2}+Tr(U_{1}^{'}{{A}_{1}}{{U}_{1}})\sigma _{M}^{2} \\ & E\left( ES{{S}_{F}} \right)=Tr\left( {{A}_{2}}\sum \right)=Tr\left( {{A}_{2}} \right)\sigma _{\varepsilon }^{2}+Tr(U_{3}^{'}{{A}_{2}}{{U}_{3}})\sigma _{FR}^{2}+ \\ & Tr(U_{2}^{'}{{A}_{2}}{{U}_{2}})\sigma _{F}^{2} \\ & E\left( ES{{S}_{FR}} \right)=Tr\left( {{A}_{3}}\sum \right)=Tr\left( {{A}_{3}} \right)\sigma _{\varepsilon }^{2}+Tr(U_{3}^{'}{{A}_{3}}{{U}_{3}})\sigma _{FR}^{2} \\ & E\left( ES{{S}_{\varepsilon }} \right)=Tr\left( {{A}_{4}}\sum \right)=Tr\left( {{A}_{4}} \right)\sigma _{\varepsilon }^{2}. \\ \end{align} \right.$ (12)
其中Tr 表示矩阵的迹，
Σ = cov(y)= U1 U'1 σM2 + U2 U'2 σF2 + U3 U'3 σFR2 + Iσε2 。 (13)

令这些平方和的数学期望等于各自平方和，并解线性方程组可得方差分量的估计为

$\begin{align} & \left( \begin{matrix} \hat{\sigma }_{M}^{2} \\ \hat{\sigma }_{F}^{2} \\ \hat{\sigma }_{FR}^{2} \\ \hat{\sigma }_{\eta }^{2} \\ \end{matrix} \right)={{\left( \begin{matrix} Tr\left( U_{1}^{'}{{A}_{1}}{{U}_{1}} \right) & Tr\left( U_{2}^{'}{{A}_{1}}{{U}_{2}} \right) & Tr\left( U_{3}^{'}{{A}_{1}}{{U}_{3}} \right) & Tr\left( {{A}_{1}} \right) \\ 0 & Tr\left( U_{2}^{'}{{A}_{2}}{{U}_{2}} \right) & Tr\left( U_{3}^{'}{{A}_{2}}{{U}_{3}} \right) & Tr\left( {{A}_{2}} \right) \\ 0 & 0 & Tr\left( U_{3}^{'}{{A}_{3}}{{U}_{3}} \right) & Tr\left( {{A}_{3}} \right) \\ 0 & 0 & 0 & Tr\left( {{A}_{4}} \right) \\ \end{matrix} \right)}^{-1}}\left( \begin{matrix} ES{{S}_{M}} \\ ES{{S}_{F}} \\ ES{{S}_{FR}} \\ ES{{S}_{\varepsilon }} \\ \end{matrix} \right)= \\ & \left( \begin{matrix} {{\lambda }_{11}} & {{\lambda }_{12}} & {{\lambda }_{13}} & {{\lambda }_{14}} \\ 0 & {{\lambda }_{22}} & {{\lambda }_{23}} & {{\lambda }_{24}} \\ 0 & 0 & {{\lambda }_{33}} & {{\lambda }_{34}} \\ 0 & 0 & 0 & {{\lambda }_{44}} \\ \end{matrix} \right)\left( \begin{matrix} ES{{S}_{M}} \\ ES{{S}_{F}} \\ ES{{S}_{FR}} \\ ES{{S}_{\varepsilon }} \\ \end{matrix} \right). \\ \end{align}$ (14)
根据二次型理论(Searle，1997)，进一步可以得到这些方差分量估计的抽样方差为
$\left\{ \begin{align} & Var\left( \hat{\sigma }_{M}^{2} \right)=2\lambda _{11}^{2}Tr{{\left( {{A}_{1}}\sum \right)}^{2}}+2\lambda _{12}^{2}Tr{{\left( {{A}_{2}}\sum \right)}^{2}}+2\lambda _{13}^{2}Tr{{\left( {{A}_{3}}\sum \right)}^{2}}+2\lambda _{14}^{2}Tr{{\left( {{A}_{4}}\sum \right)}^{2}} \\ & Var\left( \hat{\sigma }_{F}^{2} \right)=2\lambda _{22}^{2}Tr{{\left( {{A}_{2}}\sum \right)}^{2}}+2\lambda _{23}^{2}Tr{{\left( {{A}_{3}}\sum \right)}^{2}}+2\lambda _{24}^{2}Tr{{\left( {{A}_{4}}\sum \right)}^{2}} \\ & Var\left( \hat{\sigma }_{FR}^{2} \right)=2\lambda _{33}^{2}Tr{{\left( {{A}_{3}}\sum \right)}^{2}}+2\lambda _{34}^{2}Tr{{\left( {{A}_{4}}\sum \right)}^{2}} \\ & Var\left( \hat{\sigma }_{\eta }^{2} \right)=2\lambda _{44}^{2}Tr{{\left( {{A}_{4}}\sum \right)}^{2}} \\ \end{align} \right.$ (15)

下面使用Ofversten(1993)的方法考虑方差分量的假设检验问题。首先存在一个正交矩阵Γ，使得

$\Gamma y=\left( \begin{array}{*{35}{l}} {{t}_{1}} \\ {{t}_{2}} \\ {{t}_{3}} \\ {{t}_{4}} \\ {{t}_{5}} \\ \end{array} \right)=\left( \begin{matrix} G \\ 0 \\ \end{matrix} \right)\left( \begin{matrix} \beta \\ M \\ F \\ FR \\ \end{matrix} \right)+\Gamma \varepsilon =\left( \begin{array}{*{35}{l}} {{G}_{11}} & {{G}_{12}} & {{G}_{13}} & {{G}_{14}} \\ 0 & {{G}_{22}} & {{G}_{23}} & {{G}_{24}} \\ 0 & 0 & {{G}_{33}} & {{G}_{34}} \\ 0 & 0 & 0 & {{G}_{44}} \\ 0 & 0 & 0 & 0 \\ \end{array} \right)\left( \begin{matrix} \beta \\ M \\ F \\ FR \\ \end{matrix} \right)+\Gamma \varepsilon ,$ (16)
式中，G 为行满秩矩阵，并且满足R(G11 G12 G13G14)= d1，R(G22 G23 G24)= d2 - d1，R(G33 G34)=d3 - d2，R(G44)= d4 - d3，这里d1 = R(X)，d2 =R(X U1)，d3 = R(X U1 U2)，d4 = R(X U1 U2 U3)。对于零假设H0: σFR 2= 0，其检验的统计量可表示为
$F=\frac{t_{4}^{'}{{t}_{4}}/\left( {{d}_{4}}-{{d}_{3}} \right)}{t_{5}^{'}{{t}_{5}}\left( n-{{d}_{4}} \right)}\tilde{\ }F\left( {{d}_{4}}-{{d}_{3}},n-{{d}_{4}} \right).$ (17)

其次考虑零假设H0: σF2 = 0 的检验问题。设$\left( {{K}_{1}}{{K}_{2}}{{K}_{3}} \right)=\left( \begin{matrix} {{G}_{22}} & {{G}_{23}} & {{G}_{24}} \\ 0 & {{G}_{33}} & {{G}_{34}} \\ 0 & 0 & {{G}_{44}} \\ \end{matrix} \right)$，可以证明K3为行满秩矩阵，因此存在一个正交矩阵P，使得K3 K'3 =PΛP'，其中Λ 为对角矩阵。设L = Λ-1 /2P'，λ 为Λ中对角线元素倒数的最大值，ts为t5中含有s 个元素的子向量，这里s = d4 - d1，令z = L$L\text{ }\left( \begin{matrix} {{t}_{2}} \\ {{t}_{3}} \\ {{t}_{4}} \\ \end{matrix} \right)$u = z +(λ Is - Λ-1)1 /2 ts。再设N1 = L K1，N2 = L K2，那么存在一个正交矩阵Q 使得

$Q\left( {{N}_{1}}{{N}_{1}} \right)=\left( \begin{array}{*{35}{l}} {{S}_{11}} & {{S}_{12}} \\ 0 & {{S}_{22}} \\ 0 & 0 \\ \end{array} \right),$ (18)
式中，S11和S22皆为行满秩矩阵，并且R(S11 S22)=d2 - d1，R(S22)= d3 - d2。令w = Qu =$\left( \begin{matrix} {{w}_{1}} \\ {{w}_{2}} \\ {{w}_{3}}\text{ } \\ \end{matrix} \right)$其中w2和w3的维数分别为d3 - d2和d4 - d3。可以验证，在零假设H0: σF2 = 0 成立的条件下，有如下的检验统计量
$F=\frac{w_{2}^{'}{{w}_{2}}/\left( {{d}_{3}}-{{d}_{2}} \right)}{w_{3}^{'}{{w}_{3}}/\left( {{d}_{4}}-{{d}_{3}} \right)}~\tilde{\ }F\left( {{d}_{3}}-{{d}_{2}},{{d}_{4}}-{{d}_{3}} \right).$ (19)

最后，对于零假设H0: σ2M = 0 的检验问题，将上述σ2F的检验问题作适当的调整即可。将(18)式中的N1和N2交换位置以得到新的正交矩阵Q，此时w2的维数为d2 - d1，而对应的F 检验统计量为

$F=\frac{w_{2}^{'}{{w}_{2}}/\left( {{d}_{2}}-{{d}_{1}} \right)}{w_{3}^{'}{{w}_{3}}/\left( {{d}_{4}}-{{d}_{3}} \right)}~\tilde{\ }F\left( {{d}_{2}}-{{d}_{1}},{{d}_{4}}-{{d}_{3}} \right).$ (20)

根据模型(10)可以计算出单株遗传力和父本家系遗传力。单株遗传力比较简单，可直接表示为

$h_i^2 = \frac{{4\sigma _M^2}}{{Var\left({{y_{ijklt}}} \right)}} = \frac{{4\sigma _M^2}}{{\sigma _M^2 + \sigma _F^2 + \sigma _{FR}^2 + \sigma _\varepsilon ^2}}.$ (21)
父本j 家系均值的方差可表示为
$Var\left({{{\bar y}_{ \cdot j \cdot \cdot }}} \right)= {k_{1j}}\sigma _M^2 + {k_{2j}}\sigma _F^2 + {k_{3j}}\sigma _{FR}^2 + {k_{4j}}\sigma _\varepsilon ^2，$
式中，k1j =$\frac{1}{{n_{ \cdot j \cdot \cdot }^2}}\sum\limits_{i = 1}^b {n_{ \cdot j \cdot \cdot }^2}，{k_{2j}} = \frac{1}{{n_{ \cdot j \cdot \cdot }^2}}\sum\limits_{i = 1}^b {\sum\limits_{k = 1}^{{f_{ij}}} {n_{ijk}^2}，} {k_{3j}} = \frac{1}{{n_{ \cdot j \cdot \cdot }^2}}\sum\limits_{i = 1}^b {\sum\limits_{k = 1}^{{f_{ij}}} {\sum\limits_{l = 1}^r {n_{ijk}^2} }，} {k_{4j}}\frac{1}{{n_{ \cdot j \cdot \cdot }^{}}}$这里nijkl表示第i 区组第l 个重复中第j 个父本与第k 个母本交配子代的个数。取所有父本家系均值方差的平均值作为家系均值的方差，于是父本家系遗传力可表示为
$h_M^2 = \frac{{\sigma _M^2}}{{\sigma _M^2 + {k_2}/{k_1}\sigma _F^2 + {k_3}/{k_1}\sigma _{FR}^2 + {k_4}/{k_1}\sigma _\varepsilon ^2}}，$ (22)
式中，${k_1} = \sum\limits_{j = 1}^m {{k_{1j}}，} {k_2} = \sum\limits_{j = 1}^m {{k_{2j}}，} {k_3} = \sum\limits_{j = 1}^m {{k_{3j}}，} {k_4} = \sum\limits_{j = 1}^m {{k_{4j}}。} $

至于遗传力的标准误和不同性状遗传相关系数及其标准误，由于计算复杂繁琐，这里不再赘述，可以参考童春发等(2010)相关内容所用方法。

2 结果与讨论

本文针对包含区组和重复的不平衡巢式试验设计，基于约束线性模型方法给出了亲本配合力的计算公式以及相关假设检验的统计量; 基于混合线性模型方法给出了遗传方差分量的方差分析法估计公式及其假设检验统计量，同时还给出了遗传力的计算公式。这些结果对于平衡的巢式试验设计也成立。本文的统计分析方法也适用于简化的巢式试验设计模型，如仅包含区组效应或既不包含区组也不包含重复效应的统计模型。

鉴于巢式设计模型参数估计和统计量计算的复杂性，根据本文的统计分析结果，用VC + +2008 编写了Windows 计算软件WinNC1 1. 0。该软件可以计算3 种巢式设计统计模型，即包含区组和重复效应、包含区组效应但不包含重复效应以及均不包含区组和重复效应的统计模型。不管是平衡的还是不平衡的巢式试验设计数据，WinNC1 1. 0 不但能计算诸如亲本配合力、遗传方差分量、遗传力和不同性状遗传相关系数，而且还能给出这些遗传参数估计的标准误，同时还给出配合力和方差分量等假设检验统计量的p 值。软件及其使用说明可以自由下载使用，网址为http: ∥fgbio． njfu． edu． cn / tong /WinNC1 /WinNC1． htm。

与其他统计分析软件包(如SAS 等)相比，WinNC1 1. 0 有2 个主要特点: 一是使用的方便性，二是能提供更多的遗传参数计算结果。SAS 是功能强大并且是最广泛流行的统计分析软件，其中GLM、VARCOMP和MIXED 等过程可用于巢式设计统计分析。但是SAS 不是针对特定的遗传模型，使用时需要用户具有一定的SAS 代码编写技术，这给林木育种工作者带来诸多的不便。与此相反，WinNC1 1. 0 是针对巢式设计的，用户只需将试验数据按一定的格式存放在文件里并用WinNC1 打开，然后点击分析菜单就可获得所有的计算结果。

下面分别使用SAS和WinNC1 来分析白榆(Ulmus pumila)巢氏杂交试验数据(朱之悌，1990)，以这个简单的例子说明WinNC1 在巢氏设计数据分析上的简易和实用等特点。该巢氏设计有4个父本，各与2 个母本交配，每个全同胞家系选取4个单株。为了使用不平衡数据，本文将第1 个家系的前2 个单株和第3 个家系的第1 个单株删除，具体数据列于图 1。这样的数据若以Excel 格式存储则可以作为WinNC1 的输入数据文件，运行WinNC1后选用模型Ⅲ分析，其输出结果列于图 2。同样的数据用SAS 分析需要编写如图 3 的SAS 程序，其中REG 过程分析固定效应模型，而VARCOM 过程分析随机效应模型，其输出结果如图 4 所示。对于固定效应模型，SAS和WinNC1 都能给出相同的参数估计和假设检验结果，但是使用SAS 时需要在程序中按(2)式编写约束条件代码并提供每种假设检验的程序代码，特别当数据结构比较复杂时，如包含地点、区组等，固定效应数量和约束条件数量会异常多，使得程序代码多而复杂，很容易出现错误，并且对于不同的数据程序代码需要修改，即程序不具有通用性。对于随机效应模型，使用方差分析法，SAS仅能给出方差分量估计，而WinNC1 不但能给出相同的分析结果，而且还能给出方差分量估计的标准差以及假设检验统计量，同时还给出家系和单株的遗传力及其标准差。总之，同SAS 相比，WinNC1 在分析林木巢氏设计遗传模型时易于使用并能给出更多遗传统计指标。

图1 WinNC1 1. 0 软件输入数据的Excel 格式 Fig.1 The Excel format of input data in the software of WinNC1 1. 0

图2 白榆数据WinNC1 1. 0 分析输出结果 Fig.2 Output of the genetic analysis for the White Elm data with WinNC1 1. 0

图3 白榆数据遗传分析的SAS 程序 Fig.3 SAS codes for genetic analysis of the White Elm data

图4 白榆数据SAS 分析输出结果 Fig.4 SAS output of the genetic analysis for the White Elm data

值得一提的是ASReml 软件(Gilmour et al., 2009)近来在林木育种中应用较多，然而同SAS 一样，ASReml 需要用户根据模型本身编写代码，输出结果与使用者对该软件的熟练程度和对模型的统计理论理解程度有很大的关系。在统计估计方法上，ASReml 使用限制性极大似然法估计方差分量，而WinNC1 1. 0 使用常用的方差分析法，这2 种方法孰优孰劣现在还没有定论。

参考文献(References)
[1] 黄金桃.1996.杉木巢式交配设计生长性状的遗传分析.福建林业科技,23(3):1-6.(1)
[2] 梁一池.1999.杉木配合力第一轮回的选择:Ⅰ遗传参数估计及配合力效应分析.福建林学院学报,19(3):193-196.(1)
[3] 谭小梅,金国庆,邵纪清,等.2011.马尾松巢式交配设计子代生长和木材基本密度遗传分析.林业科学,47(6):30-35.(1)
[4] 唐守正,郎奎建,李海奎.2009.统计和生物数学模型计算(ForStat教程).北京:科学出版社.(1)
[5] 童春发,卫巍,尹辉,等.2010.林木半同胞子代测定遗传模型分析.林业科学,46(1):29-35.(2)
[6] 王松桂,史建红,尹素菊,等.2007.线性模型引论.北京:科学出版社.(2)
[7] 续九如.2006.林木数量遗传学.北京:高等教育出版社.(1)
[8] 于连洋,凌去非,刘洋,等.2010.48个泥鳅家系鱼种阶段的生长比较.水生生态杂志,3(6):62-68.(1)
[9] 朱之悌.1990.林木遗传学基础.北京:中国林业出版社,184.(1)
[10] Balocchi C E,Bridgwater F E,Zobel B J,et al. 1993 . Age tren ds in genetic parameters for tree height in a nonselected population of loblolly pine. Forest Science,39 ( 2 ) : 231 - 251 .(1)
[11] Comstock R E， Robinson H F. 1948. The components of genetic variance in populations of biparental progenies and their use in estimating the average degree of dominance. Biometrics，4( 4) : 254 － 266. .(1)
[12] Comstock R E,Robinson H F. 1952 . Estimation of average dominance of genes∥ Gowen J W. Heterosis. Iowa State University Press,Ames, IA,494 - 516 .(1)
[13] Gilmour A R,Gogel B J,Cullis B R,et al. 2009 . ASReml User Guide Release 3. 0 . Hemel Hempstead,UK: VSN International Ltd. .(1)
[14] Goncalves G M， Viana A P， Pereira M G，et al. 2009. Genetic parameter estimates in yellow passion fruit based on design I. Brazilian Archives of Biology and Technology，52( 3) : 523 － 530. .(1)
[15] Gwaze D P,Hardin g K J,Purnell R C,et al. 2002 . Optimum selection age for wood density in loblolly pine. Can J For Res,32 : 1393 -1399 .(1)
[16] Henderson C R. 1953 . Estimation of variance and covariance components. Biometrics,9 ( 2 ) : 226 - 252 .(1)
[17] Johnson G R，King J N. 1998. Analysis of half-diallel mating designs: I - A practical analysis procedure for ANOVA approximation. Silvae Genet，47: 74 － 79.(1)
[18] Murray L W, Ray I M, Dong H,et al. 2003 . Clarification and reevaluation of population-based diallel analyses: Gardner and Eberhart analyses II and III revisited. Crop Sci,43 ( 6 ) : 1930 -1937 .(1)
[19] Namkoong G. 1979. Introduction to quantitative genetics in forestry.U. S. Department of Agriculture Forest Service Technical Bulletin 1588.(1)
[20] Ofversten J. 1993 . Exact tests for variance components in unbalanced mixed linear models. Biometrics,49 ( 1 ) : 45 - 57 .(1)
[21] Pereira M G, Amaral Júnior A T. 2001 . Estimation of genetic components in popcorn based on the nested design. Crop Breeding and Applied Biotechnology,1 ( 1 ) : 3 - 10 .(1)
[22] Searle S R. 1997 . Linear models. New York: John Wiley & Sons,Inc. ,66 .(3)
[23] Smith S E,Kuehl R O,Ray I M,et al. 1998 . Evaluation of simple methods for estimating broad-sense heritability in stands of randomly planted genotypes. Crop Science,38 : 1125 - 1129 .(1)
[24] Tong C F,Liu G X,Yang L W,et al. 2012 . GSCA: new software and algorithms to analyse diallel mating designs based on restricted linear model. Silvae Genet,61 ( 3 ) : 126 - 132 .(1)
[25] Wang S G, Chow S C. 1994 . Advanced linear models: Theory and applications. New York: Marcel Dekker,Inc.(1)
[26] Williams E R,Matheson A C,Harwood C E. 2002 . Experimental design and analysis for tree improvement. Collingwood, VIC,Australia: CSIRO Publishing.(1)
[27] Williams E R,Matheson A C,Harwood C E. 2002 . Experimental design and analysis for tree improvement. Collingwood, VIC,Australia: CSIRO Publishing.(1)
[28] Wright J W. 1976 . Introduction to forest genetics. New York: Academic Press,Inc.(1)
[29] Wu H X,Matheson A C. 2000 . Analysis of half-diallel mating design with missing crosses: Theory and SAS program for testing and estimating GCA and SCA effects. Silvae Genet,49 : 130 - 137 .(2)
[30] Wu H X,Matheson A C. 2001 . Analysis of half-diallel mating design with missing crosses: Theory and SAS program for testing and estimating GCA and SCA variance components. Silvae Genet,50 : 265 - 271 .(1)
[31] Xiang B,Li B. 2001 . A new mixed analytical method for genetic analysis of diallel data. Can J For Res,31 ( 12 ) : 2252 - 2259 .(1)
[32] Zhang Y D，Kang M S. 1997. DIALLEL-SAS: A SAS program for Griffing’s diallel analysis. Agron J，89( 2) : 176 － 182.(1)
[33] Zhang Y D, Kang M S, Lamkey K R. 2005 . DIALLEL-SAS05 : Acomprehensive program for Griffing’s and Gardner-Eberhart analysis. Agron J,97 ( 4 ) : 1097 - 1106 .(1)
[34] Rao C R，Toutenburg H，Shalabh H，et al. 2008. Linear models and generalizations: Least squares and alternatives. Springer， New York.(2)
不平衡巢式设计遗传模型分析
童春发, 杨立伟, 蒋安纳, 施季森