混合效果建模使用MATLAB

线性混合效应(LME)模型是收集和分组总结的数据的线性回归模型的推广。线性混合效应模型为分析分组数据提供了一个灵活的框架，同时考虑了此类数据中经常出现的组内相关性。

这个例子说明了如何拟合基本的分层或多层LME模型。

版权所有(c) 2014, MathWorks, Inc。

内容

数据说明

该数据集包括美国大陆48个州1970年至1986年(17年)的年度观测数据。该数据集由Munnell(1990)提供。

  1. 国内生产总值:按国家划分的国内生产总值(原国家生产总值)
  2. 状态:指示状态的分类变量
  3. 年:国内生产总值记录年份

参考:

加载数据

安全，clc，接近所有loadPublicData

数据进行预处理

将STATE, YR和REGION转换为类别

publicdata。STATE = categorical(publicdata.STATE);计算日志(GDP)publicdata.log_GDP = log(publicdata.GDP);publicdata。GDP = [];

拟合一个线性模型，并将国家生产总值按地区可视化

美元日志\ _GDP_i = \ beta_0 + \ beta_1 YR_i + \ epsilon_i $

$\epsilon_i \sim N(0，\sigma^2)$

数据收集于1970年至1986年间的每个州。拟合线性模型时，截距代表第0年的GDP。这可能导致一个异常低的截距，由一个大的斜率来补偿。这导致斜率的估计和相应的截距估计之间存在高度负相关。我们可以通过将数据居中来消除这种相关性。

publicdata。yrcentric = publicdata。YR - mean(unique(publicdata.YR));Lm = fitlm(公共数据，'log_GDP ~ 1 + yrcentric '）;

可视化数据和拟合模型

图中,gscatter (publicdata.YR publicdata.log_GDP publicdata.STATE, [],“。”, 20);持有在情节(publicdata.YR预测(lm),“k -”，“线宽”2)标题(“简单线性模型拟合”，“字形大小”(15)包含“年”，“字形大小”、15)ylabel (“日志(GDP)”，“字形大小”15)

从图中可以清楚地看出，拟合一个忽略州(或地区)特定信息的简单线性回归模型不足以完全指定log_GDP和YR之间的关系。所有未建模的特定于组(上下文)的信息最终汇集到一个单一的误差项中。这意味着我们忽略了违反线性回归假设的组内的任何误差相关性，并导致有偏差的标准误差。

图,包含(“OLS残差”)箱线图(lm.Residuals.Raw publicdata.STATE,“定位”，“水平”）

按状态划分的OLS残差箱线图表明，状态内残差要么全为正，要么全为负，并且不在0处居中。这违反了OLS的基本假设，由此得出的推论可能是无效的。

每组拟合单独的线性模型，并可视化参数的间隔

图中,plotIntervals (publicdata'log_GDP ~ 1 + yrcentric '，“状态”）

单个置信区间清楚地表明，由于置信区间不重叠，需要随机效应来解释截距和斜率的状态间变化。

为每个状态拟合一个带虚拟变量的单一线性模型

只使用固定效应来合并群体特定效应的一种方法是为STATE和STATE与YEAR之间的相互作用引入虚拟变量，如下所示

Lm_group = fitlm(公共数据，'log_GDP ~ 1 + yrcentric + STATE + yrcentric:STATE'）;显示预测数和估计系数数。disp (＇ ＇) disp ([“变量数量:”， num2str(lm_group.NumPredictors)])估计系数的数量:，.．.num2str (lm_group.NumEstimatedCoefficients)))
变量数:2估计系数数:96

从图中可以看出，这个模型比先前的拟合模型好得多。只有在STATE中有较少的能级和大量的观测时，这才是理想的。这种方法有几个缺点。即使固定效应模型解释了STATE效应，它也不能提供数据的有用表示，因为它只对数据中的特定样本建模，而不能推广到总体。其次，模型中自由参数的数量随着STATE层数的增加而线性增加，这导致自由度的损失。当处理较小的观测集时，这可能是一个问题。

拟合随机截距模型

让我们探索一个混合效应模型，我们允许每个组的截距随机变化，在这种情况下，每个状态。我们对截距和以y为中心有一个固定效应

美元日志\ _GDP_ {ij} = \ beta_ {00} + \ beta_1 YR_ {ij} + b_ j {0} + \ epsilon_ {ij} $

$$b_{0j} \sim N(0，\sigma_{b}^2)$$

$$\epsilon_{ij} \sim N(0，\sigma^2)$$

lme_intercept = fitlme(公共数据，'log_GDP ~ 1 + yrcentric + (1|STATE)'）;

状态随机效应项的标准差，不包括0(下值:0.82594，上值:1.2324)，表明随机效应项显著。

拟合随机截距和斜率模型

美元日志\ _GDP_ {ij} = \ beta_ {00} + \ beta_ {10} YR_ {ij} + b_ j {0} + b_ {1 j} YR_ {ij} + \ epsilon_ {ij} $

左($ $ b = \ \{数组}{c}开始b_ j {0} \ \ b_ {1 j} \结束数组{}\)\ sim N(0 \σ)= N \离开(0 \左(\开始{数组}{cc} \ sigma_ {1} ^ 2 & \ sigma_ {12} ^ 2 \ \ \ sigma_ {12} ^ 2 & \ sigma_{2} \结束数组{}\)\右)$ $

$$\epsilon_{ij} \sim N(0，\sigma^2)$$

lme_intercept_slope = fitlme(公共数据，.．.'log_GDP ~ 1+年轮中心+(1+年轮中心|STATE)'）;[~，~，rEffects] = randomEffects(lme_intercept_slope);图，散射(反射。估计(1:2:结束)，反射。估计(2:2:结束))“随机效应”，“字形大小”(15)包含“拦截”，“字形大小”(15) ylabel“坡”，“字形大小”15)随机效应表中的估计列显示了估计的最佳值%线性无偏预测器(BLUP)随机效应的向量。

用似然比检验比较随机斜率模型和随机截距模型

似然比检验(LRT)统计量可用于确定更复杂的模型是否显著优于更简单的模型。这里我们将随机截距模型与随机斜率和截距模型进行比较，以确定引入随机效应对斜率的意义。

比较(lme_intercept lme_intercept_slope,“CheckNesting”,真正的)
ans =理论似然比TEST Model DF AIC BIC LogLik LRStat lme_intercept 4 -1550.5 -1531.7 779.24 lme_intercept_slope 6 -2085.4 -2057.2 1048.7 538.94 deltaDF pValue 20

p值较小，说明第二个模型在解释随机效应之间可能存在的相关性方面明显优于第一个模型。标记“CheckNesting”尝试检查第一个模型是否嵌套在第二个模型中。

使用矩阵符号拟合混合效应模型

$y = X\beta + Zb + \

如果不能使用公式轻松描述模型，则可以创建设计矩阵来定义固定和随机效果，然后调用fitlmematrix；

y = X* + Z*b + e

我们可以重现具有随机截距和斜率的模型的结果，以及它们之间可能的相关性

%响应(因变量)y = publicdata.log_GDP;固定效果设计矩阵X = [ones(height(publicdata)，1)， publicdata. yrcentric];%随机效应设计矩阵Z = [ones(height(publicdata)，1)， publicdata. yrcentric];%分组变量G = publicdata.STATE;随机效应设计矩阵可以选择包括假编码。%变量。如果提供了“G”组，则fitlmematrix%自动为您执行此操作。%拟合模型lme_matrix = fitlmematrix(X,y,Z,G，“CovariancePattern”，“对角线”）;

预测州GDP

州= {“他”，新泽西的，“马”，“CT”，“犹他”，“公司”，“v”，“为什么”};unbalancedStates = {“他”，“CT”，“为什么”};missingYears = 1970:1983;compareForecasts(公开数据，州，不平衡的州，缺失年份)

该图显示了使用这两种方法的GDP预测的比较。实圈表示用于拟合模型的观测值，空圈表示被排除在模型之外但用于预测的视觉验证的观测值。第一个图显示了固定效应模型的预测，该模型使用STATE和STATE:YEAR的虚拟变量合并了群体特定效应。第二幅图显示了随机斜率和截距模型的LME模型的预测。预测的置信区间表明，即使存在缺失数据，LME模型也能提供更有信心的预测。对WY的仔细观察表明，由于可用的观测数量较少，固定效应模型显示，随着时间的推移，GDP预测会下降。另一方面，LME模型捕捉了数据的整体趋势，并预测了GDP的增长，这与缺失数据所显示的长期趋势一致。


发布与MATLAB®R2014b