[R-sig-ME] nlme bug ?

[Cole, Tim]

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