[R-sig-ME] nlme bug ?

[Phillip Alday]

This could be a bug, but I would check two related things first:

1. Have you tried plotting a profile plot of the likelihood?

2. Have you tried  different starting values?

In my limited experience, these types of models can be surprisingly
sensitive to starting values.

Phillip

On 29/7/20 7:01 pm, Cole, Tim wrote:
> My sitar package fits the Lindstrom-Beath growth curve model, which varies the slope of the fitted spline curve y by scaling x.
>
> The model generally works well, but it fails in the simplest case, with a 1 df spline curve corresponding to a linear regression. The noddy example below fits a random intercept a and scaling factor c, with regression slope s1. The fitnlme function is called 20 times and the values hardly change, and then it produces nonsense numbers and fails.
>
> Might this be a bug? It works fine with 2+ df for the spline curve.
>
> Thanks for your thoughts.
>
> Tim Cole
>
> ``` r
> # install.packages('sitar')
> # also uses splines and nlme
> suppressMessages(library(dplyr))
>
> dat <- sitar::heights %>%
>   select(x = 'age', y = 'height', id = 'id') %>%
>   mutate(x = scale(x, TRUE, FALSE))
>
> start <- lm(y ~ splines::ns(x), data = dat)$coef
> (start <- c(s1 = start[[2]], a = start[[1]], c = 0))
> #>        s1         a         c
> #>  56.82702 128.76694   0.00000
>
> fitnlme <- function(x, s1, a, c) {
>   print(unique(s1))
>   print(unique(a))
>   print(unique(c))
>   a + drop(s1 * splines::ns(x * exp(c)))
> }
>
> nlme::nlme(y ~ fitnlme(x, s1, a, c), fixed = s1 + a + c ~ 1, random = a + c ~ 1 | id, start = start, data = dat)
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 1.490116e-08
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 1.490116e-08
> #> Warning in nlme.formula(y ~ fitnlme(x, s1, a, c), fixed = s1 + a + c ~ 1, :
> #> Iteration 1, LME step: nlminb() did not converge (code = 1). PORT message: false
> #> convergence (8)
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 1.490116e-08
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 0
> #> [1] 56.82702
> #> [1] 128.7669
> #> [1] 1.490116e-08
> #> [1] 48.22096
> #>  [1] 134.7411 130.4434 125.1878 134.2582 139.3737 132.7397 120.4782 128.5690
> #>  [9] 137.1287 125.7633 136.5816 131.8558
> #>  [1]  168940.88  664404.25 1336620.21  242694.18 -481809.06  256605.51
> #>  [7] 1934638.21  901713.94 -251016.54 1287590.98  -99830.76  482104.58
> #> Error in qr.default(t(const)): NA/NaN/Inf in foreign function call (arg 1)
> ```
>
> <sup>Created on 2020-07-29 by the [reprex package](https://reprex.tidyverse.org) (v0.3.0)</sup>
>
> --
> Population Policy and Practice
> UCL Great Ormond Street Institute of Child Health,
> 30 Guilford Street, London WC1N 1EH, UK
>
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models