主要内容

Mixed-Effects模型

混合效应模型导论

在统计中,一个效果是在预测变量的特定设置上影响响应变量值的任何东西。效应被转化为模型参数。在线性模型中,效应成为系数,代表模型项的比例贡献。在非线性模型中,效应通常有特定的物理解释,并以更一般的非线性组合出现。

固定的影响表示总体参数,假定每次收集数据时都是相同的。估计固定效应是回归建模的传统领域。随机效应,则为样本相关的随机变量。在建模中,随机效应的作用类似于附加的误差项,并且必须指定它们的分布和协方差。

例如,考虑从血液中消除药物的模型。该模型使用时间T.作为一个预测因子和药物浓度C作为回应。非线性模型术语C0.E.-R.T.结合参数C0.R.,分别表示初始浓度和消除速率。如果数据是在多个个体之间收集的,那么可以合理地假设消除率是一个随机变量R.一世取决于个人一世,围绕人口不同 R. ¯ .这个词C0.E.-R.T.就变成了

C 0. E. [ R. ¯ + R. 一世 R. ¯ ] T. = C 0. E. β + B. 一世 T.

在哪里β= R. ¯ 是一个固定的效果B.一世= R. 一世 R. ¯ 是随机效应。

当数据属于自然分组时,随机效应很有用。在药物消除模型中,组就是被研究的个体。更复杂的模型可以根据个人的年龄、体重、饮食等对数据进行分组。虽然群体不是研究的重点,但在模型中加入随机效应可以扩大推论的可靠性,超出特定的个体样本。

Mixed-effects模型占固定和随机效果。与所有回归模型一样,它们的目的是描述作为预测变量的函数的响应变量。然而,混合效果模型识别样本子组内的相关性。以这种方式,它们在完全忽略数据组之间提供折衷,并用单独的模型拟合每个组。

Mixed-Effects模型层次结构

假设一个非线性回归模型的数据为m独特的群体一世= 1,…,m.为了解释模型中的组,写回应j在集团一世为:

y 一世 j = F φ X 一世 j + ε 一世 j

y一世j是回应,X一世j是预测器的矢量,φ是模型参数的向量,和ε一世j是测量误差或过程误差。该指数j范围从1到N.一世, 在哪里N.一世组内观察的次数是多少一世.功能F指定模型的形式。通常,X一世j仅仅是观察时间吗T.一世j.误差通常被假定是独立的、相同的、正态分布的、具有常数方差。

中参数的估计φ描述所有群体的估计值都是相同的人口。但是,如果估计因组而异,则模型变为

y 一世 j = F φ 一世 X 一世 j + ε 一世 j

在混合效果模型中,φ一世可能是固定效应和随机效应的组合:

φ 一世 = β + B. 一世

随机效应B.一世通常描述为多元正态分布,均值为零,协方差Ψ.估计固定效应β随机效应的协方差Ψ提供对不承担参数的人口的描述φ一世在不同群体中都是一样的。估计随机效应B.一世还给出了数据中特定组的描述。

模型参数不必与个别效应相一致。一般来说,设计矩阵一种B.用于识别具有固定和随机效果的线性组合的参数:

φ 一世 = 一种 β + B. B. 一世

如果设计矩阵在组之间有所不同,则模型变为

φ 一世 = 一种 一世 β + B. 一世 B. 一世

如果设计矩阵在不同的观察中也不同,模型就变成

φ 一世 j = 一种 一世 j β + B. 一世 j B. 一世 y 一世 j = F φ 一世 j X 一世 j + ε 一世 j

一些特定于群体的预测因子X一世j不可能随观察而改变j.打电话给那些V.一世,模型变成

y 一世 j = F φ 一世 j X 一世 j V. 一世 + ε 一世 j

指定Mixed-Effects模型

假设一个非线性回归模型的数据为m独特的群体一世= 1,…,m.(具体地,假设组不嵌套。)为此数据指定一般非线性混合效果模型:

  1. 定义组特定的模型参数φ一世是固定效果的线性组合β和随机效果B.一世

  2. 定义响应值y一世作为一个非线性函数F的参数和群体特定的预测变量X一世

该模型是:

φ 一世 = 一种 一世 β + B. 一世 B. 一世 y 一世 = F φ 一世 X 一世 + ε 一世 B. 一世 N. 0. Ψ ε 一世 N. 0. σ 2

非线性混合效应模型的这种配方使用以下符号:

φ一世 组特定模型参数的矢量
β 固定效果的载体,建模人口参数
B.一世 多元正态分布的特定群体随机效应的向量
一种一世 组合固定效果的特定群体设计矩阵
B.一世 组合随机效应的特定群体设计矩阵
X一世 一组特定预测值的数据矩阵
y一世 群体特定响应值的数据向量
F 的一般的实值函数φ一世X一世
ε一世 群特定误差的向量,假定为独立的、相同的、正态分布的和独立的B.一世
Ψ 随机效应的协方差矩阵
σ2 误差方差,假设在观测中是恒定的

例如,考虑从血液中消除药物的模型。该模型包含两个重叠阶段:

  • 初始阶段P.在此期间,药物浓度与周围组织达到平衡

  • 第二阶段问:在此期间,药物从血液中被清除

查看多个个人的数据一世,模型是

y 一世 j = C P. 一世 E. R. P. 一世 T. 一世 j + C 问: 一世 E. R. 问: 一世 T. 一世 j + ε 一世 j

在哪里y一世j是个人观察到的浓度一世在时间T.一世j.该模型允许对不同的个体进行不同的采样时间和不同的观察次数。

消除率R.P.一世R.问:一世必须是积极的,才能对身体有意义。通过引入日志速率来实现这一点R.P.一世=日志(R.P.一世) 和R.问:一世=日志(R.问:一世)和重新参数化模型:

y 一世 j = C P. 一世 E. 经验值 R. P. 一世 T. 一世 j + C 问: 一世 E. 经验值 R. 问: 一世 T. 一世 j + ε 一世 j

选择具有随机效果的模型的哪个参数是构建混合效果模型时的重要考虑因素。一种技术是向所有参数添加随机效应,并使用差异的估计来确定其在模型中的重要性。替代方案是将模型分开拟合到每个组,而无需随机效果,并查看参数估计的变化。如果估计跨组跨越广泛变化,或者如果每个组的置信区间具有最小的重叠,则该参数是随机效果的良好候选者。

引入固定效果β和随机效果B.一世对于所有的模型参数,将模型重新表述如下:

y 一世 j = [ C ¯ P. + C P. 一世 C ¯ P. ] E. 经验值 [ R. ¯ P. + R. P. 一世 R. ¯ P. ] T. 一世 j + [ C ¯ 问: + C 问: 一世 C ¯ 问: ] E. 经验值 [ R. ¯ 问: + R. 问: 一世 R. ¯ 问: ] T. 一世 j + ε 一世 j = β 1 + B. 1 一世 E. 经验值 β 2 + B. 2 一世 T. 一世 j + β 3. + B. 3. 一世 E. 经验值 β 4. + B. 4. 一世 T. 一世 j + ε 一世 j

在通用模型的符号中:

β = β 1 β 4. B. 一世 = B. 一世 1 B. 一世 4. y 一世 = y 一世 1 y 一世 N. 一世 X 一世 = T. 一世 1 T. 一世 N. 一世

在哪里N.一世是个人的观察人数一世.在这种情况下,设计矩阵一种一世B.一世至少最初,4×4个标识矩阵。根据需要,可以改变设计矩阵,以引入单个效果的加权或时间依赖性。

拟合模型并估算协方差矩阵Ψ经常导致进一步的改进。对随机效果方差的相对较小的估计表明它可以从模型中删除。同样地,某些随机效应中的考码率的相对较小的估计表明,不需要完整的协方差矩阵。由于随机效应不观察起来,Ψ必须间接估计。指定的对角线或块对角线协方差模式Ψ可以提高拟合算法的收敛性和效率。

统计和机器学习工具箱™功能nlmefitnlmefitsa.将一般非线性混合效应模型适合数据,估计固定和随机效应。该功能还估计协方差矩阵Ψ随机效应。额外的诊断输出允许您评估模型参数数量和拟合优度之间的权衡。

指定协变量模型

如果模型在指定Mixed-Effects模型假设一个群体相关的协变量,如权重(W.)模型变成:

φ 1 φ 2 φ 3. = 1 0. 0. 0. 1 0. 0. 0. 1 W. 一世 0. 0. β 1 β 2 β 3. β 4. + 1 0. 0. 0. 1 0. 0. 0. 1 B. 1 B. 2 B. 3.

因此,参数φ一世对于任何个人一世组是:

φ 1 一世 φ 2 一世 φ 3. 一世 = β 1 + β 4. * W. 一世 β 2 β 3. + B. 1 一世 B. 2 一世 B. 3. 一世

要指定一个协变量模型,请使用'fegroupdesign'选择。

'fegroupdesign'是A.p-by-q-by-m数组指定不同的p-by-q固定效果设计矩阵为每个m团体。使用前一个示例,数组类似于以下内容:

  1. 创建数组。

    %模型中的参数个数(Phi) num_params = 3;%协变量数量num_cov = 1;%假设数据集中的组数为7 num_groups = 7;covariates = [75;52个;66;55;70;58;62); A = repmat(eye(num_params, num_params+num_cov),... [1,1,num_groups]); A(1,num_params+1,1:num_groups) = covariates(:,1)
  2. 创建一个结构体使用指定的设计矩阵。

    options.fegroupdesign = a;
  3. 指定参数nlmefit(要么nlmefitsa.)如图所示使用nlmefit和nlmefitsa的混合效应模型

选择nlmefitnlmefitsa.

统计学和机器学习工具箱提供两个功能,nlmefitnlmefitsa.用于拟合非线性混合效果模型。每个函数都提供不同的功能,可以帮助您决定使用它。

近似的方法

nlmefit提供了以下四种拟合非线性混合效应模型的近似方法:

  • 'lme'- 在当前条件估计下使用线性混合效应模型的可能性βB..这是默认值。

  • “RELME”-线性混合效应模型在当前条件估计时使用限制性似然βB.

  • '-无随机效应的一阶拉普拉斯近似。

  • 'foce'条件估计的一阶拉普拉斯近似B.

nlmefitsa.提供额外的近似方法,随机近似期望 - 最大化(SAEM)[25]有三个步骤:

  1. 模拟:生成随机效果的模拟值B.从后密度来看P.B.|Σ)给出当前参数估计。

  2. 随机近似:更新对数似然函数的期望值,从上一步中取其值,并向模拟随机效应计算的对数似然函数的平均值移动一部分。

  3. 最大化步骤:根据随机效应的模拟值,选择新的参数估计来最大化对数似然函数。

这两个nlmefitnlmefitsa.尝试寻找参数估计以最大化似然函数,这是很难计算的。nlmefit通过以各种方式近似似然函数和最大化近似函数来处理问题。它使用传统的优化技术,这些技术取决于收敛标准和迭代限制等的东西。

nlmefitsa.另一方面,以这样的方式模拟参数的随机值,即长期运行它们会聚到最大化精确似然函数的值。结果是随机的,传统的收敛试验不适用。所以nlmefitsa.提供在模拟过程中绘制结果的选项,以及多次重新启动模拟的选项。您可以使用这些特性来判断结果是否收敛到您想要的精度。

具体参数nlmefitsa.

以下参数特定于nlmefitsa..大多数控制随机算法。

  • COV0.-协方差矩阵的初始值Psi..必须是一个R.——- - - - - -R.正定矩阵。如果为空,则默认值取决于BETA0

  • ComputeStdErrors-真正的计算系数估计的标准误差并将其存储在输出中统计数据结构,或者(默认)省略此计算。

  • LogLikMethod—指定日志可能性的近似方法。

  • Nburnin.-不重新计算参数估计的初始老化迭代次数。默认是5。

  • n-控制算法的三个阶段中每个阶段执行的迭代次数。

  • NMCMCITIS.-马尔可夫链蒙特卡罗(MCMC)迭代次数。

模型和数据要求

在能力上有一些不同nlmefitnlmefitsa..因此,有些数据和模型可以使用任意一种功能,但有些可能要求您只选择其中之一。

  • 误差模型-nlmefitsa.金宝app支持各种错误模型。例如,响应的标准偏差可以是常数,与函数值成比例,或两者的组合。nlmefit在响应标准差不变的假设下拟合模型。其中一个错误模型,'指数',指定响应的日志具有恒定的标准偏差。你可以用nlmefit通过将日志响应提供为输入,并通过重写模型函数来产生非线性函数值的日志。

  • 随机效应-两个函数都将数据拟合到一个带参数的非线性函数,参数可以是简单的标量值或协变量的线性函数。nlmefit允许线性函数的任何系数同时具有固定和随机效果。nlmefitsa.金宝app仅对线性函数的常数(截距)系数仅支持随机效果,而不是用于斜率系数。所以在这个例子中指定协变量模型nlmefitsa.可以仅将前三个测试值视为随机效应。

  • 模型形式-nlmefit金宝app支持一个非常一般的模型规范,对设计矩阵有很少的限制,这些设计矩阵将固定系数和模型参数的随机效应相关。nlmefitsa.更严格的:

    • 固定效应设计必须在每个组(对每个个体)是恒定的,所以不支持观察依赖设计。金宝app

    • 随机效应设计对整个数据集必须是恒定的,因此既不支持观察依赖的设计,也不支持群体依赖的设计。金宝app

    • 正如前面提到的下随机效应,随机效应设计不得规定边坡系数的随机效应。这意味着设计必须由0和1组成。

    • 随机效果设计不得对多个系数使用相同的随机效果,并且不能对任何单个系数使用多个随机效果。

    • 固定效果设计不得使用相同的多个参数系数。这意味着它可以在每列中具有大多数非零值。

    如果你想使用nlmefitsa.对于协变量效应是随机的数据,在非线性模型表达式中直接包含协变量。不要在固定或随机效应设计矩阵中包含协变量。

  • 收敛- 如下所述模型形式nlmefitnlmefitsa.有不同的测量趋同方法。nlmefit使用传统的优化措施,nlmefitsa.提供诊断,以帮助您判断随机仿真的融合。

在实践中,nlmefitsa.往往更加强大,而且可能在困难的问题上失败。然而,nlmefit在它完全收敛的问题上可能会更快收敛。有些问题可能会受益于联合策略,例如跑步nlmefitsa.有一段时间来获得合理的参数估计,并使用那些作为其他迭代的起点使用nlmefit

使用具有混合效果模型的输出功能

Outputfcn场的选项结构指定在每次迭代后调用的一个或多个功能。通常,您可以使用输出函数在每次迭代时绘制点或从算法中显示优化数量。要设置输出功能:

  1. 将输出功能写为MATLAB®文件函数或本地功能。

  2. 使用实例化设置的值Outputfcn要成为函数句柄,即,函数的名称@的迹象。例如,如果输出函数是outfun.m ..,命令

    options = statset('OutputFcn', @outfun);

    指定outputfcn.做…的把柄outfun..要指定多个输出函数,请使用以下语法:

    options = statset('OutputFcn',{@outfun, @outfun2});
  3. 用调用优化函数选项作为输入参数。

有关输出功能的示例,请参阅示例输出功能

输出功能的结构

输出函数的函数定义行具有以下形式:

停止= outfun(β状态状态
在哪里

  • β是当前的固定效应。

  • 状态是包含来自当前迭代的数据的结构。状态的字段详细描述了结构。

  • 状态是算法的当前状态。算法的状态列出可能的值。

  • 停止是国旗吗真正的取决于优化例程是否应退出或继续。看停止标志了解更多信息。

求解器将输入参数的值传递给outfun.在每次迭代。

状态的字段

的字段如下表所示状态结构:

描述
程序
  • “ALT”交替算法优化线性混合效应或限制线性混合效应近似

  • “一圈”- 针对一阶或一阶条件估算的Laplacian近似优化

迭代 从0开始的整数。
内心的

属性中描述内部迭代状态的结构ALT腿上程序,包含以下字段:

  • 程序——当程序“ALT”

    • “PNLS”(惩罚非线性最小二乘)

    • 'lme'(线性mixed-effects估计)

    • “没有”

    程序“一圈”

    • “PNLS”(惩罚非线性最小二乘)

    • 'plm'(异形可能性最大化)

    • “没有”

  • 状态-下列其中一项:

    • “init”

    • “通路”

    • “完成”

    • “没有”

  • 迭代—从0开始的整数.为nlmefitsa.对于老化迭代,输出函数将在每次迭代后调用,且为负值status.teartation.

fval 当前对数可能性
ψ 当前随机效应协方差矩阵
θ. 的当前参数化ψ
均方误差 当前误差方差

算法的状态

下表列出了可能的值状态

状态 描述

“init”

算法在第一次迭代之前处于初始状态。

“通路”

该算法在迭代结束时。

“完成”

算法在最后一次迭代后处于最终状态。

下面的代码演示了输出函数如何使用的值状态决定在当前迭代中执行哪些任务:

切换状态案例'ITER'%将更新为绘图或GUI的绘图或GUI,用于绘图或GUIS案例的绘图或GUIS案例的DENT'%的绘图,GUI或最终绘图的句柄

停止标志

输出参数停止是国旗吗真正的.该标志告诉求解器应该退出还是继续。下面的示例展示了使用的典型方法停止国旗。

停止基于中间结果的优化输出函数可以基于传递到它的参数的值来停止估计。例如,以下代码基于存储在中的日志似然的值设置为true'fval'职位结构字段:

Stop = outfun(beta,status,state) Stop = false;%检查loglikelihood是否大于132。如果状态。fval> -132 stop = true; end

停止基于GUI输入的迭代如果你设计一个GUI来执行nlmefit迭代,您可以在用户点击时使输出功能停止停止按钮。例如,下面的代码实现了一个取消计算的对话框:

函数retval = stop_outfcn(beta,str,状态)持久性h stop;如果isequal(str.inner.state,'none')交换机(状态)案例'init'%初始化对话框stop = false;h = msgbox('按停止取消计算。',...'nlmefit:迭代0');Button = findobj(h,'类型','uicontrol');设置(按钮,'String','停止','回调',@ stopper)pos = get(h,'位置');POS(3)= 1.1 * POS(3);设置(h,'位置',pos)在对话框标题集(h,'name',sprintf('nlmefit:referation%d',... str.ceration)中)绘制案例'iter'%显示迭代号。绘制;案例'完成'%删除对话框删除(h);结束结束如果停止%停止如果已按下Delete(h)ext retval = stop;功能暂停器(varargin)%set标志停止按下按下时停止= true; disp('Calculation stopped.') end end

示例输出功能

nmlefitoutputfcn.样本统计和机器学习工具箱输出函数为nlmefitnlmefitsa..它使用固定效果初始化或更新一个情节(β)和随机效应的方差(诊断的地位。ψ)).为nlmefit,绘图还包括日志似然性(status.fval.).

nlmefitoutputfcn默认输出函数为nlmefitsa..用在一起nlmefit,在options结构中为它指定一个函数句柄:

/ / / / / / / / / / / / / / / / /
为了防止nlmefitsa.在使用此函数时,为输出函数指定一个空值:
opt = statset('outputfcn',[],...)beta = nlmefitsa(...,'选项',选择,...)
nlmefitoutputfcn停止nlmefitnlmefitsa.如果你接近它产生的数字。