广义线性混合模型的稳健估计

广义线性混合模型是处理纵向数据和重复测量数据等聚类数据的有力工具。广义线性混合模型或广义部分线性混合模型的推理大体可以通过极大似然估计(MLE)和广义估计方程(GEE)进行。但是MLE和GEE方法对于离群值是敏感的,在面对污染数据时,需要一些稳健的方法。由于大多数稳健方法都是基于回归模型,对于广义线性混合模型的研究较少。本文基于广义线性混合模型,对一些基于似然的稳健方法进行介绍。

对于广义线性混合模型,响应变量服从指数族分布,i=1,\ldots,n, j=1,\ldots,n_i,

f_{y_{i j} \mid b_i}\left(y_{i j} \mid b_i, \beta, \phi\right)=\exp \left\{\left(y_{i j} \theta_{i j}-c\left(\theta_{i j}\right)\right) / \phi+d\left(y_{i j}, \phi\right)\right\},

其中 b_i 表示随机效应,\theta_{ij}为点则参数,\phi为尺度参数,c(\cdot)b(\cdot,\cdot)为已知函数。响应变量的条件均值\mu_{ij}=E[y_{ij}|b_i]=\dot{c}(\theta_{ij})\eta_{ij}=g(\mu_{ij})=x_{ij}^T\beta+z_{ij}^Tb_ig(\cdot)为一个已知的联系函数。当g(\mu_{ij})=\theta_{ij} 时表示为典则连接。一般假设随机效应 b_i\sim N(0,\Sigma_b)。在很多情况下,尺度参数\phi是已知的或者可以估计的,为简单起见,我们假设\phi=1

给出广义线性混合模型的似然函数为

l(\beta,\Sigma_b)=\sum_{i=1}^n\log \int\prod_{j=1}^{n_i}f_{y_{i j} \mid b_i}\left(y_{i j} \mid b_i, \beta, \phi\right)f_{b_i}(b_i|\Sigma_b)db_i.

对于广义线性混合模型,接下来介绍几种稳健估计方法。

1.基于score函数的稳健估计

为极大化似然函数,通过对似然求导有

\bg_white E\left[\frac{\partial \ln f_{y \mid b}(y \mid b, \boldsymbol{\beta})}{\partial \boldsymbol{\beta}} \mid y\right]=\mathbf{0},\qquad E\left[\frac{\partial \ln f_b(b\mid \mathbf{\Sigma_b})}{\partial \boldsymbol{\Sigma_b}} \mid y\right]=\mathbf{0} .

根据前面介绍的广义线性混合模型有

\frac{\partial \ln f_{y_i \mid b}\left(y_i \mid b_i, \boldsymbol{\beta}\right)}{\partial \boldsymbol{\beta}}=\left\{y_i-\mu_i(\boldsymbol{\beta}, b_i)\right\} x_i,

由于上述score函数对异常值敏感,我们给出它的稳健版本,

E\left[\sum _ { i = 1 } ^ { n } \left\{\psi_c\left(\frac{y_i-\mu_i(\boldsymbol{\beta}, \mathbf{b}_i)}{\sigma_i(\boldsymbol{\beta}, \mathbf{b}_i)}\right) \sigma_i(\boldsymbol{\beta}, \mathbf{U}) w\left(\mathbf{x}_i\right) \mathbf{x}_i-\mathbf{q}(\boldsymbol{\beta}, \mathbf{b}_i)\right\} \mid \mathbf{y}\right]=\mathbf{0}

其中\sigma_i^2(\boldsymbol{\beta}, \mathbf{u})=\operatorname{var}\left[Y_i \mid \mathbf{u}\right]=c^{\prime \prime}\left(\theta_i\right)\omega为权重函数,用来降低设计空间中杠杠点的权重,\psi_c可以考虑为Huber的psi函数:\psi_c(x)=\max (-c, \min (x, c))。由上式可以看出,稳健版本的score函数从响应变量和协变量两个角度出发来降低异常值的影响。可以选择权重函数为马氏距离的函数

\omega(\mathbf{x}) =w\left(\mathbf{x}, \mathbf{m}_x, \mathbf{S}_x\right) =\min \left\{1,\left(\frac{b_0}{\left(\mathbf{x}-\mathbf{m}_x\right)^T \mathbf{S}_x^{-1}\left(\mathbf{x}-\mathbf{m}_x\right)}\right)^{\gamma / 2}\right\},

其中\gamma \geq 1b_0选为自由度为\mathbf{x}维数的卡方分布的0.95分位数,\mathbf{m}_x, \mathbf{S}_x分别为\mathbf{x}的位置参数和尺度参数的稳健估计。\mathbf{q}(\boldsymbol{\beta}, \mathbf{b}_i)是为了保证估计量具有一致性的修正项目,

\mathbf{q}(\boldsymbol{\beta}, \mathbf{b})=\frac{1}{n} \sum_{i=1}^n E\left\{\psi_c\left(\frac{Y_i-\mu_i(\boldsymbol{\beta}, \mathbf{b}_i)}{\sigma_i(\boldsymbol{\beta}, \mathbf{b}_i)}\right) \mid \mathbf{u}\right\} \sigma_i(\boldsymbol{\beta}, \mathbf{b}_i) w\left(\mathbf{x}_i\right) \mathbf{x}_i

该方法通过对score函数进行修正,来得到一个稳健版本的scoer函数,通过score函数的期望为0的式子对参数进行求解,得到参数的稳健估计。

2.稳健对数似然函数

上面的方法是基于score函数来获得稳健性估计,具有较好的性质。但变量选择是统计建模中时常需要面对的问题,也是统计学研究中的重点和热点问题之一。上述方法在变量具有稀疏性时无法得到想要的结果。接下来介绍一种基于稳健似然函数的稳健方法,在给出稳健似然函数后,可以根据惩罚似然的方法进行稳健的变量选择。

同样从响应变量和协变量两个方面来降低异常值对于似然的贡献。首先考虑响应变量中的异常值,若y_{ij}为异常值,则通过一个合适的值\vartheta_{ij}y_{ij}进行修改 。若y_{ij}偏大,令y_{ij}^*=y_{ij}-\vartheta_{ij},若y_{ij}偏小,令y_{ij}^*=y_{ij}+\vartheta_{ij}。然后考虑协变量中的异常值,若x_{ij}为异常值,则通过权重函数\omega_{ij}来降低该数据点对似然的贡献。根据上述分析,可以得到稳健对数似然函数为

l^R(\beta,\Sigma_b)=\sum_{i=1}^n\log \int\prod_{j=1}^{n_i}[f_{y_{i j}^* \mid b_i}\left(y_{i j}^* \mid b_i, \beta, f, \phi\right)]^{\omega_{ij}}f_{b_i}(b_i|\Sigma_b)db_i

其中\vartheta_{ij}=\vartheta_{ij}(\tau_{ij})=\text{sign} (\tau_{ij})(|\tau_{ij}|-c)v_{ij}^{1/2} I(|\tau_{ij}|>c)\tau_{ij}=(y_{ij}-\mu_{ij})/v_{ij}^{1/2}为Pearson残差,\mu_{ij}=g^{-1}(x_{ij}^T\beta+z_{ij}^Tb_i)为响应变量的条件均值,v_{ij}为响应变量的条件方差,c为调节常数。则可以得到

y_{ij}^*= \begin{cases} \mu_{ij}-cv_{ij}^{1/2}, &y_{i j}\leq \mu_{ij}-cv_{ij}^{1/2},\\ y_{ij}, & \mu_{ij}-cv_{ij}^{1/2}<y_{i j}\leq\mu_{ij}+cv_{ij}^{1/2}, \\ \mu_{ij}+cv_{ij}^{1/2}, & y_{i j}>\mu_{ij}+cv_{ij}^{1/2}. \end{cases}

同样给出权重函数

\omega(\mathbf{x}) =w\left(\mathbf{x}, \mathbf{m}_x, \mathbf{S}_x\right) =\min \left\{1,\left(\frac{b_0}{\left(\mathbf{x}-\mathbf{m}_x\right)^T \mathbf{S}_x^{-1}\left(\mathbf{x}-\mathbf{m}_x\right)}\right)^{\gamma / 2}\right\},

在给出稳健对数似然后,为保证估计量的一致性,给出调整的稳健对数似然

l^R_a\left(\beta,\Sigma_b\right)=l^R\left(\beta,\Sigma_b\right)-a_n(\beta,\Sigma_b),

其中\frac{\partial}{\partial\psi}a_n(\psi)=E\frac{\partial}{\partial\psi}l^R(\psi)

该方法通过降低响应变量和协变量中异常值对似然的贡献来得到稳健的似然函数,通过极大化稳健似然函数可以得到稳健估计。且在给出似然函数后,可以通过惩罚方法进行稳健变量选择。大概做法是,将随机效应方差进行Cholesky分解\Sigma_b=\Gamma\Gamma^T,通过对广义线性混合模型模型重新参数化,得到\eta_{ij}=g(\mu_{ij})=x_{ij}^T\beta+z_{ij}^T\Gamma b_i^*,其中\Gamma是一个下三角矩阵,然后可以通过对模型中\beta\Gamma进行惩罚来同时获得固定效应和随机效应的稳健变量选择结果。其中对\Gamma惩罚是对下三角矩阵的每一行非零元素进行组惩罚。

由于上述稳健似然函数中存在积分函数,无法直接通过极大化方法进行求解,可以通过EM算法来进行迭代求解。

  • 1
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 1
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值