我对校正表型一直有些疑惑，以及对育种值准确性的计算有疑惑，这里写一下我困惑的地方和思考的结果。

整理目前的情况

动物的校正表型

校正表型的由来

课本书上写，估计育种值的准确性是估计育种值和真实育种值的相关系数。由于实际数据的真实育种值未知，因此我们需要使用别的基准值来代替真实育种值，比如奶牛上可能会使用 daughter yield deviation 和 Deregressed EBV 之类的值，但是一般来说其他动物上用的最多的是校正表型。

校正表型的另一个作用就是作用GS分析的因变量。

校正表型的计算

育种计算实际数据计算准确性的时候很多时候都会用到校正表型（ yc，有叫 corrected phenotypes ，也有 adjusted phenotypes ），准确性一般使用估计育种值与校正表型的相关。

校正表型的做法一般是用全部表型数据采用传统的常规评估得到的，但是这里还有一点疑义或者说分歧，如下

• 有的文献写的是 y-Xb（表型剔除固定因子） （文献1-2）
• 有的文献写的是 ebv+e （育种值加残差）

这两种做法当随机因子只有加性效应时没有区别，只有当存在其他随机因子时存在不一样（文献2中就是存在非加性随机因子，但是没有剔除这个随机因子）。

如果是重复度量性状，那就是一个个体的多次校正表型求均值，才是最终使用的校正表型。

植物多年多点的BLUE值

植物多年多点表型校正惯例

参考文献3和4的网页，植物的表型一般是多年多点的情况（还有一个可能的因素就是重复 rep，也就是同一年同一个地点存在多个重复），我看植物上的惯例是将品种作为固定因子，其他因子作为随机因子，如年，点，重复（一般用的是 点:重复，因为重复一般是嵌套在地点里的），还可以加入这些因子之间的互作效应。其校正表型就采用品种的 BLUE 值加截距（一般分析软件中当只有品种这一个固定因子的时候得到的 BLUE 值就是包括解决的）。

为什么用BLUE值

其实我还没遇到过其他的数据的情况，只做过动植物，就是上面这两种情况，但是可以讨论一下。

笼统地说，这种情况就是没有系谱，只有表型，因子除了个体本身的效应，就是环境因素，那么在使用 BLUP 分析的情况，就存在：

• 个体效应是作为固定因子，还是随机因子？
• 各种环境因素是作为固定因子，还是环境因子？

这里我感觉环境因素的分类不是重点（植物上貌似都是将环境因素作为随机因子），作为固定因子和随机因子都可以。甚至有想法的话，可以做一下因子的检验。

重点就在于个体效应是应当作为固定因子和随机因子。

Don’t BLUP twice

文献5（《Don’t BLUP twice》）是由植物遗传统计学大牛 James B. Holland 和 Hans-Peter Piepho 在 2024 年发表的短文，堪称近年来植物 GWAS（全基因组关联分析）和 GS（基因组选择）领域最硬核的“避坑指南”。它直接掐断了一个困扰无数研究生和研究员的纠结：“我在第一步校正表型的时候，到底该用 BLUE 还是 BLUP？”

结论：用 BLUE 值

这篇文献的核心结论如下：

• 第一阶段（表型校正）：只是用表型数据，不用基因型数据。必须把品种效应值当作固定效应（Fixed Effect），其他环境效应为随机因子、固定因子或协变量（视具体情况而定），提取品种效应的 BLUE 值。
• 第二阶段（GWAS/GS分析）：将品种效应的 BLUE 值作为表型，此时再把品种效应当作随机效应（含亲缘关系矩阵），计算 BLUP 值。

这篇文章是分情况来讨论的。

当我们使用两步法进行GWAS或GS分析时，第一步对表型进行校正的时候，我们有下面三种选择：

• 品种效应作为固定因子
• 品种效应作为随机因子（假设样本之间不相关，也就是亲缘关系矩阵是一个单位矩阵）
• 品种效应作为加性随机因子（使用G阵作为亲缘关系矩阵）

下面来说为什么不能用另外两种选择

为什么不能用 blup 值

如果在第一阶段把基因型当随机效应算 BLUP，由于混合线性模型中随机效应的分布假设（均值为0），与 BLUE 值相比，品种的 BLUP 值会向群体均值发生收缩。这种收缩效应取决于其估计可靠性，而估计可靠性又由每个品种的性状观测数、性状遗传力以及根据遗传相似性加权后的亲属信息量决定（如果品种效应加入了亲缘关系矩阵）。

如果品系效应未通过亲缘关系建模（即遗传值被视为独立同分布），且数据完全平衡（无缺失值、区组完整、所有品系在每次重复和每个环境中均有观测），那么 BLUE 和 BLUP 的排序将完全一致。在此情形下，某品系 BLUP 值与总体均值的偏差，将等于其 BLUE 值与总体均值偏差乘以品系均值的广义遗传力。因此，两者唯一的区别仅在于：受收缩效应影响，品系值与总体均值的偏离程度发生了缩放。

然而在现实实验中，数据常因缺失值或不完整区组设计而呈现不平衡状态。试想一个不平衡数据集：某些品系在特定重复区组中缺失。此时，各品系在数据集中的效应信息量存在差异。若将品系视为独立随机效应，其BLUP值将根据各自的重复次数和测试环境数量，以不同程度向群体均值收缩。测试环境越多的品系，其遗传值预测的收缩程度越低，这反映了更多信息带来的更高精度。极端情况下，对于毫无观测数据的品系，其BLUP预测将直接等于群体均值——这源于随机效应的分布假设：它们是从某个更大的参考群体随机抽取而来，且先验期望均值为0。

相反，若将品系视为固定效应，则无数据支持的品系无法估计其BLUE值：我们对这类效应不作任何分布假设，它们无需服从正态分布，也不存在一个更大的参考群体。

就选种目的而言，这种反映品系BLUP预测精度的差异化收缩效应是有利的（Robinson，1991）：它能降低不可靠预测的权重；但对于任何第二阶段分析而言，这却是弊端——因为品系的BLUP值不仅反映遗传信息，还掺杂了数据结构和可能的试验设计特征。

为什么不能用加性blup值

如果第一阶段不仅把基因型当随机，还加入了亲缘关系矩阵，也就是使用 GBLUP 的 BLUP 值，那么会进一步加重上面的问题。此时的 GBLUP 会呈现差异化收缩：收缩程度既取决于该品系的直接测量数据量，也取决于数据集中近缘品系的数量。

此外，某个特定品系的 GBLUP 值，不仅反映该品系自身的实测数据，还会按基因组亲缘关系加权，融合其他品系的观测信息。如果 GBLUP 模型有效，其预测的品系排名通常会与 BLUE 值不同。就选种而言，这正是我们想要的：GBLUP 能最优预测多基因效应（育种值）。然而，GWAS 的核心目的却截然不同——它旨在剔除多基因背景效应和群体结构干扰后，识别仍有显著效应的单个的遗传变异。因此，如果我们输入 GWAS 的是已经过“多基因背景校正”的预测值，就等于把 GWAS 要搜寻的大部分信号提前剔除了！此前，Hadfield 等人（2010）也曾在进化生物学研究中指出过滥用 BLUP 带来的其他问题。

用大白话解释， GBLUP 会根据亲戚的表现来猜你的表现。如果一群有某个突变的亲戚都长得高，系统就会把大家都拔高。但这样一来，第二阶段 GWAS 再去找这个突变时，会发现“哇，大家矫正后都一样高了”，直接把真正的因果变异信号给抹平了。

使用blue值看似矛盾实则不然

在第一阶段将品系效应视为统计上的固定效应，而在第二阶段又将其视为随机效应，这看似自相矛盾。但实际上，为了获得接近单阶段分析结果的解，这种做法必不可少（Smith等，2001；Piepho等，2012；Damesa等，2017；Verbyla，2023）。这种表面的矛盾，可能正是导致人们在两阶段GWAS中如何正确处理遗传效应时感到长期困惑的原因。

此外，作物遗传学家在二阶段GWAS或GS（基因组选择）分析的第一阶段，可能会因忌惮而不愿将品系视为固定效应。因为他们深知，就选育和育种的目的而言，BLUP值要优于BLUE值（Piepho等，2008），而GBLUP（基因组BLUP）效果更好（Isik等，2017）。因此，研究人员理所当然地希望在第一阶段就获取BLUP或GBLUP值，以此作为“最优”的品系预测值输入到第二阶段的GWAS中。

但是，如果分析的目标是获取用于后续分析的输入“表型”汇总值，这种做法就大错特错了。 相反，在第一阶段分析中，品系应当被拟合为固定效应而非随机效应。随后，在第二阶段分析中，再将品系拟合为随机效应，并将其多基因效应建模为与实现的亲缘系数（realized kinship coefficients）成比例的协方差。这一建议同样适用于二阶段GWAS、GS或QTL分析（Damesa等，2019；Feldmann等，2023；George和Cavanagh，2015）。无论哪种分析，第一阶段模型的作用都是将复杂数据集“提炼”简化，为每条品系生成一个单一的表型汇总值；真正的区别在于第二阶段的分析目标。

在两阶段中都把遗传效应当作随机效应处理是不合适的，因为这会导致对BLUP结果再进行BLUP分析——这就是我们提出“不要两次BLUP（Don’t BLUP Twice）”建议的原因。在动物育种中，人们曾对第一阶段BLUP使用 “deregression”（“解回归”，即“unshrinking”，“反收缩” ）来处理二阶段GS分析（Garrick等，2009）。而Verbyla（2023）证明，本文推荐的第一阶段估算BLUE值的方法，与这种解回归法有着密切的数学联系（Garrick等，2009）。

两步法的改进

当由于数据稀疏导致非遗传效应的方差参数估计准确性至关重要时，通过将基因型视为随机效应来改进模型可能会带来益处（我的理解就是数据量太少，导致无法估计随机环境因子的方差组分，此时可以这么做））。

所建议的精细化方法在第一阶段包含两个具体步骤（George 和 Cavanagh 2015; Smith 等 2001; Verbyla 2023）：

(i) 首先将基因型作为加性效应（加入G阵）拟合混合模型，并保存所有的非遗传方差分量估计值；

(ii) 接着将基因型作为固定效应拟合相同的模型，同时代入第 (i) 步中获得的非遗传方差分量估计值（就是使用第一步非加性效应的方差组分，总感觉这样用不太对），进而计算出品系的 BLUE 值。

第二阶段的 GWAS 分析不变。

总的分析流程示意图如下

这种做法的问题

若数据存在不平衡，将品系 BLUEs 视为精度相同的独立“观测值”的两阶段分析，会造成部分信息流失（就是不同品种的BLUE值可靠程度不一样，这就是核心问题）。针对基因组预测，可利用多种加权分析方案（Piepho 等，2012；Damesa 等，2017；Endelman，2023）来校正由数据不平衡导致的 BLUEs 间精度差异与协方差。

其中最精确的方法是将品系 BLUEs 的完整方差-协方差矩阵带入第二阶段分析，并将其整合为以品系 BLUEs 为响应变量的残差方差-协方差结构。该方法也是网状 Meta 分析（Network Meta-analysis）中的标准流程（Schwarzer 等，2015；Madden 等，2016），但其计算负荷较大，尤其当品系数目庞大时。采用对角近似（diagonal approximations）可节省计算资源，但这会以牺牲部分效率为代价（Damesa 等，2017；Verbyla，2023）。

机器学习上的情况

在机器学习上，和计算育种值这种分析比较类似的，那就是监督学习的回归问题。但是，二者的目标是完全不同的。回归问题的目标就是预测真实表型，预测值约接近真实表型就越准；而育种分析的目标是预测加性效应，但是加性效应本身是不可知的（真实育种值）。因此，在常规的机器学习问题上，不存在对表型进行校正的情况。

但是当我们用机器学习进行GS分析的时候，我们不能像常规的机器学习那样做。如果按照常规的机器学习的做法，因变量使用真实表型，自变量是所有环境因子加上所有SNP的分型，那么预测得到的结果就是这个个体在特定环境下的表型预测值。当我们需要预测新样本时，我们不仅要提供这个样本的所有SNP的分型，而且还需要提供它的所有环境因子的水平。这就会产生两个问题：

• 第一，预测目标不对，选种要的是加性效应，当然机器学习用的基本是非线性方法，那我们的目标也是预测基因型值，而不是这种某一套特定环境水平下的表型（单纯因为环境不同造成的差异用来选种压根没有意义）；
• 第二，即便这么做，也没法预测新样本，因为没有表型的新样本的环境因子通常是未知的（比如测定年季）。

因此我们不能这么做，我们也需要在机器学习分析时也需要校正掉全部环境因素，因变量使用校正表型，自变量使用所有SNP的分型。

更何况，如果是重复度量表型，我们就必须要将这多个表型用某种方式得到一个汇总的值。

因此对我们这个专业的数据而言，机器学习使用的校正表型需要进行下面的处理：

1. 从表型中剔除所有环境因子的效应
2. 重复度量表型要得到一个唯一的值

按照机器学习的思路，那么最简单的办法就是线性回归，然后将环境因素都加到线性回归模型中，残差e=y-Xb 就是我们需要的校正表型。如果是多次度量性状，那就采用多次残差的均值。当然，多次度量性状也可以将品种效应加到线性回归模型中，就像植物的做法那样，采用品种效应的 BLUE 值作为校正表型。

汇总整理

根据上面的内容，汇总一下我们得到存在下面的做法

• 有系谱的情况下（动物）：使用 y-Xb 或 ebv+e
• 没系谱的情况下（植物）：使用品种效应的 BLUE 值
• 机器学习：采用线性回归的残差

下面我得尝试从原理上进行解释和解惑。首先由于环境因子到底是固定因子还是随机因子对我们这个问题没有影响，因此从简化角度考虑，假设全部环境因子全设为固定因子。

植物的BLUE值做法不普适

先看植物上的 BLUE 值，当所有环境因素也视为固定因子，此时的分析模型就是一个固定因素模型，更准确地说就是一个线性回归模型。如果表型不是重复度量表型，那么此时最小二乘方程组压根就无法求解（特征数目大于表型数目），或者更准确地说没有唯一解。

实际上我们回到正常的混合线性模型中，广义最小二乘也是脱胎于最小二乘，这种情况也无法求解，因此当表型不是重复度量性状时，我们压根就不能将品种效应视为固定因子。（这其实就是固定因子的一个假设，固定因子的水平数目是有限的，隐含条件就是固定因子的水平数目肯定不能超过表型数目）

也就是说，植物上的这种做法不是通用的，只能用于重复度量性状，尤其是这种多年多点的性状，所以不用再考虑这一做法了（不普适）。

再回到植物的 BLUE 的做法，假设不是重复度量的性状，那么这种情况还能怎么处理呢？一种很自然的引申做法就是我们从线性回归模型中将品种效应去除（因为没法估计所以就得剔除），只保留所有的环境因素，然后使用 残差e=y-xb 作为校正表型，这才是普适的做法。这就和上面机器学习的处理思路合流了，二者完全一致。

动物到底是 y-Xb 还是 ebv+e？

下面这段话是我自己盘的，核心就是 ebv+e 间接使用了 BLUP 值，而 BLUP 值是有偏的，不能使用。

这里我们知道当没有环境随机因子时，二者是相同的。所以我们就必须得考虑有环境随机因子的情况，假设存在一个固定环境因子 fix1 ，存在一个随机环境因子 rand1。此时 ebv+e 就等于 y - fix1 - rand1 ，而上面文献告诉我们不能用 BLUP 值，因为 BLUP 值会向群体均值收缩，收缩程度会受到数据结构的影响（数据量越多，收缩程度越小），因此 rand1的 BLUP 值也是这样，同样是受到数据结构的影响，数据不平衡程度越高，这个影响就越严重，因此可以说 rand1的 BLUP 值是有偏的，因此就不能使用，因此 ebv+e = y - fix1 - rand1 也是有偏的，不能使用（如果数据完全平衡，那么可以用，但这在实际数据是不可能的）。

因此，我们要使用 y-Xb ，此时我们就完全没有用到 BLUP 值结果，只用到了 BLUE 值。而 BLUE值 b 是可以通过广义最小二乘方程组估计得到的，因此压根不用构建整个MME，也就不用计算随机因子的BLUP值。

最后再分析一下使用 y-Xb 时随机因子的作用。在广义最小二乘方程组中，对应的模型公式是 y=Xb+e ，所有的随机因子其实都归到了残差里面，加入随机因子，其实就是对 var(y) (表型的协方差矩阵) 做了进一步优化，更符合真实情况，这样算出来的 b 会更准一点（不然正常线性回归就是残差乘以单位矩阵）。

最后统一于y-Xb

因此，这三种方式的最终均统一于 y-Xb ，其模型为
$$
\boldsymbol{y=Xb+e}
$$
又可以分为两种情况：

1. 线性回归：不用讲，采用最小二乘方程组求解
2. 混合线性模型：随机因子仅仅是完善了var(y)，采用广义最小二乘方程组求解

用校正表型作为两步法第二步的因变量这一点没有问题了，作用就用两点：

1. 字面意思，剔除固定因子的影响（我的理解：固定因子水平的值是一个固定的值，因此这个影响是恒定的，比如A厂+5，B厂-3，这是必须要剔除了，不然表型是有偏的；而随机因子是满足均值是0的正态分布的随机数，这就像是随机的噪音，在算不准的情况可以保留，区别就是噪音大一点小一点而已）
2. 多次度量的表型得到一个唯一的值（当数据不平衡的情况下，也会存在上面植物BLUE值的问题，就是不同校正表型精度不一样，重复次数多的精度就高）

准确性为什么与校正表型求相关

这也是我一直最困惑的一点，现在来看的话就是没有更好的指标了，无奈之选。

按照我上面整理的来看，假设没有固定因子（只有截距），那么 y-Xb 其实还是 y ，此时校正表型就是表型。其实校正表型就是剔除了固定因子的表型，它的本质还是表型。这有点车轱辘话，我的意思就是校正表型和理论上的真实育种值没有半毛钱关系。我之前就是想不通这一点。

那么按照理论来说，如果我没有记错的话，育种值和表型值的相关系数就是 h （遗传力的开方），这也就是为什么交叉验证的结果受到遗传力的影响，遗传力越高，交叉验证的准确性越高。

我估计这也就是奶牛为什么一直用什么逆回归育种值，DRP 之类的指标，因为奶牛收集的表型数据比别的物种多得多，用这些指标应该是比校正表型更好。

参考文献

1. Liu T, Qu H, Luo C, et al. Accuracy of genomic prediction for growth and carcass traits in Chinese triple-yellow chickens[J]. BMC genetics, 2014, 15(1): 110.

2. Liu J, Yang G, Kong J, et al. Using single-step genomic best linear unbiased prediction to improve the efficiency of genetic evaluation on body weight in Macrobrachium rosenbergii[J]. Aquaculture, 2020, 528: 735577.

3. https://blog.csdn.net/yijiaobani/article/details/103213440

4. https://www.jianshu.com/p/e6f36bb430ec

5. Holland J B, Piepho H P. Don’t BLUP twice[J]. G3: Genes, Genomes, Genetics, 2024, 14(12): jkae250.