机器算法验证 - 用随机斜率和截距拟合 Poisson GLM 混合模型 - 吾爱随笔录

用随机斜率和截距拟合 Poisson GLM 混合模型

机器算法验证混合模式广义线性模型泊松分布随机效应模型

2022-03-25 08:21:07

我目前正在研究一系列泊松时间序列模型，试图估计计数获取方式的变化（从一种诊断测试切换到另一种诊断测试）的影响，同时控制随时间推移的其他趋势（比如在发病率）。我有许多不同网站的数据。

虽然我也一直在修补 GAM，但我已经将一系列非常基本的 GLM 与其中的时间趋势相结合，然后汇总结果。SAS 中的代码如下所示：

PROC GENMOD data=work.data descending;
  model counts = dependent_variable time time*time / link=log dist = poisson;
run;

或者在 R 中：

glm(counts ~ dependent_variable + time + time*time, family="poisson")

然后进行这些估计，并将它们汇集到各个站点。也有人建议我尝试对每个站点使用具有随机斜率和截距的泊松混合模型，而不是池化。所以基本上你会有dependent_variable的固定效果，然后是截距和时间的随机效果（或者理想的时间和时间^ 2，虽然我知道这有点毛茸茸）。

我的问题是我不知道如何拟合这些模型之一，而且似乎混合模型是每个人的文档突然变得非常不透明的地方。任何人都有一个简单的解释（或代码）关于如何适应我想要的东西，以及要注意什么？

2个回答

在 R 中：

library(lme4)
lmer(counts ~ dependent_variable + (1+time|ID), family="poisson")

在这种情况下并且此代码适合模型 $Y_i \sim {\rm Poisson}(\lambda_i)$

\log (λ_{i}) = β_{0} + β_{1} X_{i} + η_{i 1} + η_{i 2} t

$\log(\lambda_i) = \beta_0 + \beta_1 X_i + \eta_{i1} + \eta_{i2}t$

其中是，是并且是。是固定效应是随机效应，其方差由模型估计。 $X_i$ dependent_variable $t$ time $i$ ID $\beta_0, \beta_1$ $\eta_{i1}, \eta_{i2}$

这是一个带有一些快速模拟数据的示例，其中随机效应方差确实为 0，协变量没有影响，每个结果都是，并且每个人在时间 . ${\rm Poisson}(1)$ $t = 1,...,10$

x = rnorm(100)
t = rep(1:10,each=10)
ID = rep(1:10,10)
y = rpois(100,1)
g <- lmer(y ~ x + (1+t|ID), family="poisson")
summary(g)
Generalized linear mixed model fit by the Laplace approximation 
Formula: y ~ x + (1 + t | ID) 
   AIC   BIC logLik deviance
 108.8 121.9 -49.42    98.85
Random effects:
 Groups Name        Variance  Std.Dev. Corr   
 ID     (Intercept) 0.0285038 0.168831        
        t           0.0027741 0.052669 -1.000 
Number of obs: 100, groups: ID, 10

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.09078    0.11808  -0.769    0.442
x            0.13670    0.08845   1.546    0.122

Correlation of Fixed Effects:
  (Intr)
x -0.127

需要注意的一点 - 该列只是该列Std.Dev.的平方根，而不是方差估计的标准误差！Variance

在 SAS 中：

proc glimmix data = yourdata ic = q;
    class id;
    model y = x / dist = poisson solution;
    random intercept t / subject = id;
run;

但是当然有很多选择，或多或少有用，可以玩。

其它你可能感兴趣的问题

上一篇为什么设计矩阵 X 的秩等于 X'X 的秩？下一篇当 PMF 包含 0 时，如何计算 Kullback-Leibler 散度？