新闻详情

新闻详情

首页 / 资讯中心 / 详情

BLUPs与收缩效应:混合模型中个体预测的智能校准

发布时间:2026/6/17 13:38:44
BLUPs与收缩效应:混合模型中个体预测的智能校准
1. BLUPs与收缩效应混合模型中那些被“拉回均值”的预测值你有没有遇到过这样的情况在分析重复测量数据时比如同一头种公猪在不同月龄的精液量、同一学生在多次考试中的成绩、同一工厂在不同班次的次品率——每个个体都有一条自己的变化轨迹但单独为每个个体拟合一条回归线结果要么太“飘”要么太“散”一放到新数据上就崩盘这时候老手们往往不会急着调参或加正则项而是会说一句“试试混合模型看看BLUPs怎么把极端估计往回收一收。”BLUPBest Linear Unbiased Prediction这个词听起来像统计学黑话但它背后的核心动作其实特别朴素用群体信息给个体估计“托底”和“纠偏”。它不是强行让所有人一样也不是放任各自为政而是在“全盘统一”和“完全独立”之间找到一条有依据、可计算、抗干扰的中间路径。这正是收缩shrinkage的本质——不是削足适履而是基于数据本身的变异结构自动判断哪些个体偏离是真实信号哪些更可能是噪声。本文聚焦的不是如何跑通SAS代码而是带你一层层拆开BLUPs的生成逻辑为什么它必须同时包含固定效应和随机效应为什么“部分池化”partial pooling能天然抑制过拟合为什么同一个体的BLUP估计值其标准误比普通OLS估计小得多我会用129头种公猪在4个时间点的精液量数据作为贯穿始终的实操案例从原始散点图开始一步步展示“无池化模型”如何疯狂贴合个别异常点、“全池化模型”如何抹平所有个体差异再到混合模型如何用方差分量作为“调节旋钮”把每头猪的截距和斜率稳稳地拉向群体中心。所有解释都基于真实建模过程中的数值输出和图形反馈不堆砌公式不空谈理论只讲你在SAS或R里跑出结果后真正该盯住哪几行数字、哪几张图、哪几个诊断指标。如果你正在处理纵向数据、多水平数据、遗传评估、教育测评或任何带嵌套结构的数据这篇内容就是你调试混合模型时最该放在手边的“操作备忘录”。2. 混合模型设计思路为什么必须“固定随机”双轨并行2.1 三种建模哲学的底层逻辑对比理解BLUPs的第一步是彻底厘清“全池化”pooled、“无池化”no-pooled和“部分池化”partially pooled这三种策略的根本分歧。它们不是技术选项而是对数据生成机制的不同信念。以129头种公猪的精液量数据为例我们观测的是每头猪在4个时间点比如月龄12、14、16、18的精液体积单位mL。目标是刻画每头猪的生长轨迹——是线性上升还是先快后慢的二次曲线关键在于我们相信这些轨迹是“完全一致”、“完全无关”还是“大同小异”全池化模型Pooled Model它的信念是“天下猪一般黑”。它假设所有129头猪共享同一个截距、同一个线性斜率、同一个二次斜率。模型形式是 $Y_{ij} \beta_0 \beta_1 \cdot \text{time}{ij} \beta_2 \cdot \text{time}{ij}^2 \varepsilon_{ij}$其中下标 $i$ 表示猪$j$ 表示时间点。这里没有“猪”的标识符所有数据被当作一个大池子搅匀了算。它的优势是参数极少、估计稳定、标准误小劣势是它完全无视个体差异如果某头猪天生精液量就高或者发育节奏明显不同这个模型就会系统性地低估或高估它的表现。它给出的预测是“如果这头猪是平均猪它应该长什么样”而不是“这头特定的猪它自己会怎么长”。无池化模型No-Pooled Model它的信念是“每头猪都是宇宙中心”。它为每头猪单独拟合一条完整的二次曲线$Y_{ij} \beta_{0i} \beta_{1i} \cdot \text{time}{ij} \beta{2i} \cdot \text{time}{ij}^2 \varepsilon{ij}$。这里$\beta_{0i}, \beta_{1i}, \beta_{2i}$ 都是针对第 $i$ 头猪的独立参数。模型自由度爆炸式增长——129头猪光截距就要估计129个它的优势是灵活性极强能完美捕捉每头猪的独特模式劣势是极度脆弱。当某头猪的4个观测点中有一个异常值比如采样误差导致的超高读数这条专属曲线就会被这个点剧烈扭曲导致对其他3个时间点的预测也严重失真。更致命的是它无法泛化你永远无法回答“一头新来的、从未测过的猪它的典型轨迹是什么”因为它压根没学习任何关于“猪群”的共性知识。部分池化模型Partially Pooled Model / Mixed Model它的信念是“猪群有共性个体有特性”。它承认存在一个“猪群基准线”由固定效应 $\beta_0, \beta_1, \beta_2$ 描述但同时也允许每头猪在这个基准线上进行“个性化微调”由随机效应 $u_{0i}, u_{1i}, u_{2i}$ 描述。模型形式是 $Y_{ij} (\beta_0 u_{0i}) (\beta_1 u_{1i}) \cdot \text{time}{ij} (\beta_2 u{2i}) \cdot \text{time}{ij}^2 \varepsilon{ij}$。这里的 $u_{0i}, u_{1i}, u_{2i}$ 不是待估的固定参数而是被假定服从正态分布的随机变量$u_{0i} \sim N(0, \sigma^2_{u0}), u_{1i} \sim N(0, \sigma^2_{u1}), u_{2i} \sim N(0, \sigma^2_{u2})$。这个正态分布的方差 $\sigma^2_{u0}, \sigma^2_{u1}, \sigma^2_{u2}$ 才是真正的待估参数它们量化了“猪群内部在截距、线性斜率、二次斜率上的变异程度有多大”。这才是混合模型的精髓所在它不直接估计129个截距而是估计一个“截距变异度”它不直接估计129个斜率而是估计一个“斜率变异度”。BLUPs即 $u_{0i}, u_{1i}, u_{2i}$ 的最佳线性无偏预测值正是从这个变异度的框架中“推导”出来的而非“拟合”出来的。2.2 收缩效应Shrinkage的数学直觉与现实意义那么“收缩”是如何发生的为什么BLUPs会比无池化模型的估计值更“靠拢”群体均值这并非人为设定的规则而是上述正态分布假设下的必然结果。我们可以用一个最简化的例子来建立直觉假设我们只关心每头猪的截距即基线精液量忽略斜率。无池化模型给出的第 $i$ 头猪的截距估计 $\hat{\beta}{0i}^{NP}$就是它4个时间点观测值的平均值。而混合模型的BLUP估计 $\hat{u}{0i}^{BLUP}$其经典表达式在平衡设计下是 $$\hat{u}{0i}^{BLUP} \frac{\sigma^2{u0}}{\sigma^2_{u0} \sigma^2_\varepsilon / n_i} \cdot (\bar{Y}_i - \hat{\beta}_0)$$ 其中$\bar{Y}i$ 是第 $i$ 头猪的平均观测值$\hat{\beta}0$ 是固定效应截距即所有猪的总体平均值$\sigma^2{u0}$ 是随机截距的方差代表猪群间的真实变异$\sigma^2\varepsilon$ 是残差方差代表测量误差或个体内变异$n_i4$ 是每头猪的观测次数。这个公式揭示了收缩的全部秘密收缩因子Shrinkage Factor$\frac{\sigma^2_{u0}}{\sigma^2_{u0} \sigma^2_\varepsilon / n_i}$ 这个分数永远在0到1之间。当猪群变异大、误差小、数据多时如果 $\sigma^2_{u0}$ 很大猪和猪差别确实很大而 $\sigma^2_\varepsilon$ 很小测量很准$n_i$ 很大观测点多那么分母接近 $\sigma^2_{u0}$整个分数接近1。这意味着 $\hat{u}_{0i}^{BLUP} \approx \bar{Y}_i - \hat{\beta}_0$即BLUP几乎等于无池化估计收缩很弱——因为数据足够好足以支撑一个独特的个体估计。当猪群变异小、误差大、数据少时如果 $\sigma^2_{u0}$ 很小猪和猪其实差不多而 $\sigma^2_\varepsilon$ 很大测量噪音大$n_i$ 很小只有1-2个点那么分母远大于 $\sigma^2_{u0}$整个分数接近0。这意味着 $\hat{u}_{0i}^{BLUP} \approx 0$即BLUP几乎为零个体估计被完全“拉回”到群体均值 $\hat{\beta}_0$——因为数据太弱不足以证明这头猪真的与众不同最稳妥的猜测就是“它大概率就是平均猪”。因此收缩不是一种武断的惩罚而是一种基于证据强度的理性妥协。它自动地、连续地在“完全信任个体数据”和“完全信任群体先验”之间滑动。在129头猪的数据中那些观测值离群、波动剧烈的猪其BLUP截距会被大幅收缩而那些数据稳定、且明显高于或低于平均值的猪其BLUP截距则收缩得很少。这种机制天然地、无需额外正则化项就实现了对过拟合的免疫。2.3 方差分量混合模型的“决策中枢”在混合模型中固定效应 $\beta$ 告诉我们“平均而言发生了什么”而随机效应的方差分量 $\sigma^2_{u0}, \sigma^2_{u1}, \sigma^2_{u2}, \sigma^2_\varepsilon$ 则决定了“个体差异有多重要”以及“收缩力度有多大”。它们是整个模型的“决策中枢”其估计值的合理性直接决定了BLUPs的可信度。在SAS的PROC MIXED输出中Covariance Parameter Estimates表格就是这个中枢的仪表盘。对于我们的精液量数据我们需要重点关注Intercept(orUN(1,1))这是随机截距的方差 $\sigma^2_{u0}$。如果它被估计为0或非常接近0比如 0.001说明猪与猪之间的基线精液量几乎没有差异强行加入随机截距不仅无益反而会增加模型复杂度可能导致收敛失败。此时一个简单的全池化模型可能就足够了。time(orUN(2,2))这是随机线性斜率的方差 $\sigma^2_{u1}$。它告诉我们猪与猪在精液量随月龄增长的“速度”上变异有多大。如果它显著大于0说明不同猪的发育节奏确实不同需要为每头猪估计一个个性化的斜率。time*time(orUN(3,3))这是随机二次斜率的方差 $\sigma^2_{u2}$。它衡量的是“加速度”的变异即每头猪的生长曲线弯曲程度是否不同。在我们的初步探索中作者发现这个值很小暗示大多数猪的曲线形状如先快后慢是相似的不需要为每头猪都估计一个独特的曲率。Residual这是残差方差 $\sigma^2_\varepsilon$代表了在考虑了所有固定和随机效应后仍无法解释的变异主要来自测量误差和个体内的随机波动。一个健康的混合模型其方差分量应该呈现出清晰的层级$\sigma^2_{u0} \sigma^2_{u1} \sigma^2_{u2} \sigma^2_\varepsilon$这符合生物学直觉——猪的起始点差异最大其次是生长速度再次是生长加速度最后才是测量误差。如果某个方差分量被估计为负值SAS会报告为0或其标准误极大这通常是一个危险信号表明该随机效应可能并不必要或者数据不足以支持其估计应果断将其从模型中移除。3. 核心细节解析从数据可视化到BLUPs提取的完整实操链3.1 数据初探用图形诊断“是否值得做混合模型”在敲下第一个PROC MIXED命令之前最不能跳过的一步是用最原始的图形去“感受”数据。这比任何统计检验都更能告诉你混合模型是否是正确的工具。对于129头猪×4个时间点的数据我习惯性地绘制三张图第一张所有猪的轨迹叠加图All Trajectories Plotproc sgplot databoar_data; series xtime yvolume / groupboar_id lineattrs(thickness0.5); xaxis values(12 to 18 by 2); yaxis labelSemen Volume (mL); title All 129 Boars: Individual Trajectories; run;这张图的目的是看“森林”而非“树木”。如果所有线条都紧紧簇拥在一条主干周围几乎没有交叉和发散那说明个体差异极小全池化模型足矣。反之如果线条像一团乱麻有的从低处陡升有的从高处缓降有的平直如线有的弯曲如弓这就强烈暗示了巨大的个体变异混合模型大有可为。在我们的数据中这张图显示了明显的发散趋势尤其是截距起始点的离散度远大于斜率的离散度这已经为“随机截距随机斜率”模型投下了第一张赞成票。第二张无池化估计值的分布图No-Pooled Estimates Distribution/* 先用PROC REG为每头猪拟合二次模型输出截距、斜率、二次斜率 */ proc reg databoar_data outestno_pooled_estimates noprint; by boar_id; model volume time time*time; quit; /* 绘制三个参数的分布 */ proc sgpanel datano_pooled_estimates; panelby _name_ / layoutlattice columns3; histogram estimate; rowaxis labelEstimate Value; colaxis labelParameter; title Distribution of No-Pooled Estimates; run;这张图将无池化模型得到的129个截距、129个线性斜率、129个二次斜率分别画成三个直方图。它的核心价值在于识别异常值和评估变异幅度。如果截距的直方图跨度很大比如从50mL到150mL而斜率的直方图却非常窄比如集中在0.5到1.5 mL/月这就印证了之前的观察猪与猪的起点差异巨大但生长速率相对一致。更重要的是如果某个截距估计值远远甩开其他所有值比如一头猪的估计截距是200mL而第二名是120mL这头猪很可能就是数据录入错误或极端异常值。在混合模型中这样的点会被强力收缩但如果它真实存在我们就需要在建模前决定是保留它并接受其BLUP被大幅拉回还是剔除它如果确认是错误。第三张无池化估计 vs. 全池化估计的散点图Shrinkage Preview/* 将全池化模型的固定效应估计值合并到无池化估计表中 */ data no_pooled_with_fixed; merge no_pooled_estimates (ina) pooled_fixed_estimates (inb); by _name_; if a; diff estimate - fixed_estimate; run; proc sgscatter datano_pooled_with_fixed; plot diff*fixed_estimate / group_name_ markerattrs(size5); xaxis labelFixed Effect Estimate; yaxis labelDeviation from Fixed Effect; title Deviations of No-Pooled Estimates from Pooled Mean; run;这张图是收缩效应的“预演”。横轴是全池化模型给出的固定效应值即群体均值纵轴是每头猪的无池化估计值减去这个均值即个体偏差。如果所有点都紧密地分布在纵轴y0附近说明个体差异小收缩空间小。如果点散布得很开尤其是截距的点_name_Intercept形成了一条长长的水平带而斜率的点_name_time则聚集在短带中这就直观地展示了收缩将主要作用于截距对斜率的作用较小。这张图本身不是最终答案但它为我们选择随机效应的结构比如是否要包含随机二次斜率提供了最直接、最不容忽视的视觉证据。3.2 模型构建与比较五种方案的SAS实现与解读基于上述图形诊断我们构建并比较五种嵌套模型。所有模型都使用相同的固定效应部分volume time | time即包含截距、time、time*time区别仅在于随机效应部分。以下是SAS代码的核心骨架及关键解读模型1全池化模型Pooledproc mixed databoar_data methodml; class boar_id; model volume time | time / solution; /* 无random语句即无随机效应 */ ods output SolutionFpooled_fixed; run;解读这是所有模型的基准线。SolutionF表给出了 $\hat{\beta}_0, \hat{\beta}_1, \hat{\beta}_2$。注意其标准误StdErr列会非常小因为它利用了全部数据。但它的预测区间Pred会非常宽因为它把所有个体差异都归为“未解释的残差”导致对单个猪的预测信心不足。模型2随机截距模型Random Interceptproc mixed databoar_data methodml; class boar_id; model volume time | time / solution; random intercept / subjectboar_id; ods output SolutionRrandom_intercept_blups; run;解读random intercept / subjectboar_id告诉SAS为每个boar_id估计一个随机截距 $u_{0i}$。SolutionR表输出的就是这129个BLUPs $\hat{u}{0i}$。此时Covariance Parameter Estimates表会给出 $\sigma^2{u0}$ 和 $\sigma^2_\varepsilon$。一个关键检查是$\sigma^2_{u0}$ 是否显著大于0可以通过查看其Wald Z值Estimate / StdErr是否大于2来粗略判断。如果Z值很小说明加入随机截距并未显著提升模型拟合可以考虑放弃。模型3随机截距随机斜率模型Random Intercept Slopeproc mixed databoar_data methodml; class boar_id; model volume time | time / solution; random intercept time / subjectboar_id; ods output SolutionRrandom_int_slope_blups; run;解读random intercept time / subjectboar_id表示为每个boar_id估计一个随机截距 $u_{0i}$ 和一个随机线性斜率 $u_{1i}$。SolutionR表现在有两行一行是所有 $\hat{u}{0i}$另一行是所有 $\hat{u}{1i}$。Covariance Parameter Estimates表现在有三行UN(1,1)($\sigma^2_{u0}$),UN(2,2)($\sigma^2_{u1}$), 和Residual($\sigma^2_\varepsilon$)。UN(1,2)行是截距和斜率的协方差如果为负说明起点高的猪其生长速度往往较慢这是一个有趣的生物学发现。模型4随机截距随机斜率随机二次斜率模型Full Randomproc mixed databoar_data methodml; class boar_id; model volume time | time / solution; random intercept time time*time / subjectboar_id; /* 注意此模型极易不收敛 */ run;解读这是最复杂的模型试图为每个猪估计三个随机效应。但在实践中由于数据点有限每猪仅4点它几乎总是失败。SAS会报错WARNING: Stopped because of infinite likelihood或ERROR: Unable to make progress with the optimization.。这并非软件缺陷而是数据在“说话”它告诉你用4个点去估计3个独立的、可能相关的参数信息量严重不足。此时我们必须回到图形诊断如果二次斜率的无池化估计分布非常窄就果断放弃这个随机效应。模型5协方差结构更复杂的模型如UNproc mixed databoar_data methodml; class boar_id; model volume time | time / solution; random intercept time / subjectboar_id typeun; /* typeun 表示估计一个2x2的非结构化协方差矩阵 */ run;解读typeun比默认的typevc方差成分更灵活它会估计 $\sigma^2_{u0}, \sigma^2_{u1}$ 和它们的协方差 $\sigma_{u01}$而不是假设它们独立。这能捕捉到截距和斜率之间的相关性。但它的代价是增加了3个待估参数vs. 2个同样面临收敛风险。只有当typevc模型的UN(1,2)协方差估计值很大且符号明确时才值得尝试typeun。3.3 BLUPs的提取、验证与可视化让“收缩”看得见得到SolutionR表后BLUPs只是第一步。真正的价值在于如何使用和验证它们。以下是我在实际项目中必做的三件事第一步将BLUPs与固定效应相加得到个体特异性预测/* 假设我们有random_int_slope_blups表其中包含u0和u1的BLUPs */ /* 我们需要将其与固定效应beta0, beta1, beta2合并 */ proc sql; create table final_predictions as select a.boar_id, a.time, /* 个体特异性预测 固定效应 随机效应 */ (b.beta0 c.u0) (b.beta1 c.u1)*a.time b.beta2*a.time*a.time as pred_volume, a.volume as actual_volume from boar_data as a left join pooled_fixed as b on 11 /* 获取固定效应 */ left join random_int_slope_blups as c on a.boar_id c.subject; quit;这个计算过程至关重要。它让我们看到混合模型的最终预测并非一个单一的全局曲线而是129条“被收缩过”的、彼此平行如果只含随机截距或不平行如果含随机斜率的个体曲线。这些曲线不再像无池化模型那样“各走各路”也不像全池化模型那样“千篇一律”而是呈现出一种和谐的多样性。第二步绘制“收缩图”Shrinkage Plot/* 创建一个数据集包含每头猪的无池化截距、BLUP截距、以及固定效应截距 */ proc sql; create table shrinkage_data as select np.boar_id, np.estimate as np_intercept, blup.estimate as blup_intercept, pf.estimate as fixed_intercept from no_pooled_estimates as np inner join random_int_slope_blups as blup on np.boar_id blup.subject and np._name_ Intercept and blup.effect Intercept inner join pooled_fixed as pf on np._name_ pf._name_ and pf._name_ Intercept; quit; proc sgplot datashrinkage_data; scatter xnp_intercept yblup_intercept / markerattrs(size3); lineparm x0 y0 slope1 / lineattrs(patterndash colorgray); refline fixed_intercept / axisx lineattrs(colorred); refline fixed_intercept / axisy lineattrs(colorred); xaxis labelNo-Pooled Intercept Estimate; yaxis labelBLUP Intercept Estimate; title Shrinkage of Intercept Estimates; run;这张图是BLUPs的灵魂所在。图中的虚线是yx线代表“无收缩”。所有点都应该落在这条线的下方如果固定效应是均值且BLUPs被收缩向它。红色参考线是固定效应截距 $\hat{\beta}_0$。你会清晰地看到那些远离红线的无池化估计左上角和右下角的点其对应的BLUP点被明显地“拉向”了红线而那些本来就靠近红线的点则几乎不动。这就是收缩效应的可视化证明。第三步计算并比较标准误验证“收缩”的收益BLUPs带来的最大好处之一是其标准误Standard Error of Prediction, SEP远小于无池化估计的标准误。这是因为BLUPs利用了群体信息降低了不确定性。在SAS的SolutionR输出中StdErr Pred列就是SEP。我们可以计算一个简单的指标收缩比率Shrinkage RatioSEP_BLUP / SE_NoPooled。一个健康的模型这个比率应该在0.3到0.7之间。比率越小说明收缩带来的精度提升越大。如果比率接近1说明收缩几乎没有发生模型可能过于保守如果比率接近0说明收缩过度模型可能过于激进。这个比率是评估混合模型是否“恰到好处”的一个黄金指标。4. 实操过程与核心环节实现从SAS输出到专业解读的全流程4.1 关键SAS输出表格的逐行解码运行一个成功的混合模型后SAS会输出大量表格。新手常被淹没在信息中而老手只盯住几张表。以下是我每天都会打开的三张核心表格及其解读要点表格1Covariance Parameter Estimates方差分量估计表| Cov Parm | Subject | Estimate | Standard Error | Z Value | Pr |Z| | | :--- | :--- | :--- | :--- | :--- | :--- | | UN(1,1) | boar_id | 125.3 | 18.7 | 6.70 | .0001 | | UN(2,2) | boar_id | 2.1 | 0.8 | 2.63 | 0.0086 | | Residual | | 8.9 | 0.5 | 17.8 | .0001 |UN(1,1)这是随机截距的方差 $\sigma^2_{u0} 125.3$。其Z值6.70远大于2P值0.0001说明猪与猪在基线精液量上的差异是真实存在的、统计显著的。这个值的大小125.3意味着如果我们随机抽取两头猪它们的基线精液量之差的期望方差约为 $2 \times 125.3 250.6$标准差约为15.8mL。UN(2,2)这是随机线性斜率的方差 $\sigma^2_{u1} 2.1$。Z值2.63P值0.0086也达到显著性说明猪与猪在生长速度上确实存在可测量的差异但其变异幅度标准差约1.45 mL/月远小于截距的变异标准差约11.2 mL。Residual残差方差 $\sigma^2_\varepsilon 8.9$。这是模型无法解释的“噪音”其标准差约为3.0 mL与生物学测量误差的预期相符。一个经验法则是$\sigma^2_\varepsilon$ 应该显著小于 $\sigma^2_{u0}$否则说明模型未能捕捉到主要的变异来源。提示如果UN(1,1)的P值不显著比如Pr|Z| 0.05不要急于删除随机截距。先检查数据质量再考虑是否是模型设定问题如时间尺度不合适。有时一个不显著的方差分量恰恰说明你的研究设计非常成功——个体差异被控制得很好。表格2Solution for Fixed Effects固定效应解| Effect | time | timetime | Estimate | Standard Error | DF | t Value | Pr |t| | | :--- | :--- | :--- | :--- | :--- | :--- | :--- | :--- | | Intercept | . | . | 85.2 | 2.1 | 128 | 40.57 | .0001 | | time | . | . | 3.8 | 0.4 | 382 | 9.50 | .0001 | | timetime | . | . | -0.05 | 0.01 | 382 | -5.00 | .0001 |这张表给出了群体层面的“平均轨迹”平均基线是85.2mL平均每月增长3.8mL但增长速度每月减缓0.05mL即曲线是倒U型先快后慢。DF自由度列很有意思Intercept的DF是128129-1这是基于“猪”这个水平计算的而time和time*time的DF是382129×4 - 3 - 129这是基于“观测点”这个水平计算的。这反映了混合模型的多水平本质。表格3Solution for Random Effects随机效应解即BLUPs| Effect | Subject | Estimate | Std Err Pred | DF | t Value | Pr |t| | | :--- | :--- | :--- | :--- | :--- | :--- | :--- | | Intercept | 101 | 15.3 | 2.8 | 3 | 5.46 | 0.0112 | | Intercept | 102 | -8.7 | 2.8 | 3 | -3.11 | 0.0512 | | ... | ... | ... | ... | ... | ... | ... | | time | 101 | 0.2 | 0.3 | 3 | 0.67 | 0.5521 | | time | 102 | -0.1 | 0.3 | 3 | -0.33 | 0.7621 |这张表就是BLUPs的宝库。Estimate列是 $\hat{u}{0i}$ 和 $\hat{u}{1i}$ 的值。Std Err Pred预测标准误是关键它比普通OLS的标准误小得多。例如猪101的截距BLUP是15.3mL其SEP是2.8mL而如果用无池化模型估计其标准误可能高达5-6mL。这2-3mL的精度提升在育种值评估中可能就意味着选留决策的天壤之别。t Value和Pr |t|列在这里没有传统意义上的假设检验意义。因为BLUPs是预测值不是参数估计值我们不检验“这个猪的截距是否为零”而是直接使用它。所以我通常会忽略这两列只关注Estimate和Std Err Pred。4.2 收敛性诊断与模型简化实战指南混合模型最大的“拦路虎”不是理解而是收敛失败。SAS报错ERROR: Did not converge.或WARNING: The final Hessian matrix is not positive definite.是家常便饭。这不是你的错而是数据在提醒你“这个模型太贪心了。”以下是我总结的、屡试不爽的“四步简化法”第一步检查数据结构与缺失值在运行任何模型前先执行proc freq databoar_data; tables boar_id * time / list missing; run;确保每头猪在4个时间点都有数据。如果有缺失PROC MIXED默认会删除整条记录listwise deletion这可能导致有效样本量锐减。对于少量缺失可以考虑用多重插补对于大量缺失则需重新审视实验设计。第二步从最简模型开始逐步添加复杂度永远不要一上来就跑random intercept time time*time。我的标准流程是跑random intercept。如果收敛记录 $\sigma^2_{u0}$。在此基础上添加
网站建设 高端定制 企业官网