如何从总结模型拟合中提取信息

Question

library(nlme)
fm1 <- nlme(height ~ SSasymp(age, Asym, R0, lrc),
            data = Loblolly,
            fixed = Asym + R0 + lrc ~ 1,
            random = Asym ~ 1,
            start = c(Asym = 103, R0 = -8.5, lrc = -3.3))
> summary(fm1)
Nonlinear mixed-effects model fit by maximum likelihood
  Model: height ~ SSasymp(age, Asym, R0, lrc) 
 Data: Loblolly 
       AIC      BIC    logLik
  239.4856 251.6397 -114.7428

Random effects:
 Formula: Asym ~ 1 | Seed
            Asym  Residual
StdDev: 3.650642 0.7188625

Fixed effects: Asym + R0 + lrc ~ 1 
         Value Std.Error DF   t-value p-value
Asym 101.44960 2.4616951 68  41.21128       0
R0    -8.62733 0.3179505 68 -27.13420       0
lrc   -3.23375 0.0342702 68 -94.36052       0
 Correlation: 
    Asym   R0    
R0   0.704       
lrc -0.908 -0.827

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-2.23601930 -0.62380854  0.05917466  0.65727206  1.95794425 

Number of Observations: 84
Number of Groups: 14 

我对从NLME拟合的汇总输出中提取信息感兴趣。

我想要提取：

1. 随机效应的标准偏差（即Asym的标准偏差，其值为3.65）。对于这个问题，我尝试过使用fm1$apVar，但没有成功。
2. 固定效应的参数估计值（即Asym = 101.44960，R0 = -8.62733等），可以通过fixef(fm1)进行提取。
3. 固定效应的标准误差（即2.46，0.317，0.034）。对于这个问题，我尝试过使用sqrt(diag(fm1$varFix))，但这些值与固定效应下标准误差列下面的值不完全匹配？
4. 对数似然（即-114.7428，可以使用fm1$logLik进行提取）
5. 残差（即0.7188625，可以使用fm1$Residuals进行提取）

我的最终目标是拟合多个模型，并将它们各自的汇总估计值存储在一个组织良好的data.frame中。

fm1 <- nlme(height ~ SSasymp(age, Asym, R0, lrc),
            data = Loblolly,
            fixed = Asym + R0 + lrc ~ 1,
            random = Asym ~ 1,
            start = c(Asym = 103, R0 = -8.5, lrc = -3.3))

fm2 <- nlme(height ~ SSasymp(age, Asym, R0, lrc),
            data = Loblolly,
            fixed = Asym + R0 + lrc ~ 1,
            random = Asym ~ 1,
            start = c(Asym = 103, R0 = -5.4, lrc = -3.3))

summary(fm1)
summary(fm2)

mylist = list(NULL, summary(fm1), NULL, summary(fm2), NULL, NULL)

假设我的列表对象看起来像mylist。现在我想创建一个类似于下面的data.frame：

model    FixedAsym    FixedAsymStdError   FixedR0      ...     Residual
 1       101.44960        2.4616951       -8.62733            0.7188625
 2       101.44934        2.4616788       -8.62736     ...    0.7188625

为了创建这个数据框（行数对应于我在mylist中有多少模型摘要），我需要系统地从模型摘要输出中提取那些值（编号为1-5）。

- Adrian

---

这可能很有用：http://r4ds.had.co.nz/many-models.html - Joshua Rosenberg

创建一个函数，返回一个命名向量，其中包含您想要为单个模型返回的所有统计信息，然后使用sapply在模型列表上运行它。 - lmo

@lmo 听起来不错。你有什么想法可以提取固定效应的标准误吗？ - Adrian

你可以使用broom包中的tidy和glance来获取大部分所需的模型估计值。除了使用tidy之外，另一种方法是从摘要输出中提取tTable并从中提取所需的内容（如SE）;例如，summary(fm1)$tTable。 - aosmith

对数似然和残差都有提取函数：logLik和residuals。您可以将固定效应标准误与 sqrt(diag(vcov(fm1)))(或 sqrt(diag(fm1$varFix)，但最好使用提取函数)匹配到摘要输出中，如果您不调整它们，则可以将其与 summary(fm1, adjustSigma = FALSE) 输出进行比较，或者使用 method="REML"拟合您的模型-请参见 ?nlme:::summary.lme。 - user20650

Answer 1

这里还有几个片段...

as.numeric(VarCorr(fm1)[,2])
# [1] 3.6506418 0.7188625

summary(fm1)$tTable[,2]
#       Asym         R0        lrc 
# 2.46169512 0.31795045 0.03427017 

# looks like you don't need this one anymore, but here's a way of getting it
summary(fm1)$corFixed
#            Asym         R0        lrc
# Asym  1.0000000  0.7039498 -0.9077793
# R0    0.7039498  1.0000000 -0.8271022
# lrc  -0.9077793 -0.8271022  1.0000000

非常抱歉这不是一个完整的回答 - 制作类似您描述的摘要表格可能会很困难，因为每个潜在行的结构都会有所不同，并且将取决于哪些变量被包括为固定和随机效应。

- Matt Tyers