首页 文章

R:如何从摘要模型拟合中提取信息

提问于
浏览
5
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拟合的摘要输出中提取信息 .

我想提取一下

  • 随机效应的StdDev(即Asym的StdDev,= 3.65)对于这个我试过 fm1$apVar 但没有运气 .

  • 固定效应的参数估计(即Asym = 101.44960,R0 = -8.62733等),可通过 fixef(fm1) 提取

  • 固定效应的Std.Error(即2.46,0.317,0.034) . 对于这个我尝试了 sqrt(diag(fm1$varFix)) 但这些值与固定效果下的Std.Error列下的值不完全匹配?

  • logLikelihood(即-114.7428,可以使用 fm1$logLik 提取)

  • 残差(即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

并创建此data.frame(行数对应于我在 mylist 中有多少个模型摘要)我需要系统地从模型摘要输出中提取这些值(编号为1-5) .

1 回答

  • 0

    这里还有几件......

    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
    

    抱歉这不是一个完整的答案 - 可能很难创建一个像你所描述的汇总表,因为每个潜在行的结构将是不同的,并且将取决于包含哪些变量作为固定和随机效应 .

相关问题