R:如何从nlme调用中获取Hessian

Question

library(nlme)
fm1 <- nlme(height ~ SSasymp(age, Asym, R0, lrc),
            data = Loblolly,
            fixed = Asym + R0 + lrc ~ 1,
            random = Asym ~ 1,
            start = c(Asym = -10311111, R0 = 8.5^4, lrc = 0.01),
            verbose = TRUE)

**Iteration 1
LME step: Loglik: -312.2787, nlminb iterations: 23
reStruct  parameters:
    Seed 
10.41021 
Error in nlme.formula(height ~ SSasymp(age, Asym, R0, lrc), data = Loblolly,  : 
  Singularity in backsolve at level 0, block 1

我试图nlme通过观察粗麻布来研究为什么有些模型不能成功.有办法以某种方式提取这个矩阵吗？

我也在查看fdHess函数(也来自同一个pacakge),"使用有限差分评估一个近似Hessian和一个标量函数的梯度"这是否等同于函数中当前实现的内容nlme？

Answer 1

我相信你的问题是由于起点选择不当造成的。向量c(Asym = 103, R0 = -8.5, lrc = -3.3)毫无复杂性地收敛：

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))

#> Nonlinear mixed-effects model fit by maximum likelihood
#>   Model: height ~ SSasymp(age, Asym, R0, lrc) 
#>   Data: Loblolly 
#>   Log-likelihood: -114.7428
#>   Fixed: Asym + R0 + lrc ~ 1 
#>       Asym         R0        lrc 
#> 101.449600  -8.627331  -3.233751 
#> 
#> Random effects:
#>  Formula: Asym ~ 1 | Seed
#>             Asym  Residual
#> StdDev: 3.650642 0.7188625
#> 
#> Number of Observations: 84
#> Number of Groups: 14 

归根结底，模型拟合可以理解为一个优化问题。当您的模型是非线性的（例如混合效应模型）时，必须使用迭代优化算法来解决该问题。因此，起始值的选择非常关键。这是一篇讨论该主题的精彩科学文章：

Balsa-Canto, E.、Alonso, AA 和 Banga, JR 生化网络动态建模的迭代识别程序。BMC 系统生物学 4, 11 (2010) doi:10.1186/1752-0509-4-11