[R-sig-ME] nlme vs nlme4

Cole, Tim tim.cole at ucl.ac.uk
Mon Jan 23 10:46:49 CET 2017


Terry,

I've also recently reported that same error, at
https://stat.ethz.ch/pipermail/r-sig-mixed-models/2016q4/025213.html.

BTW where is your function dfun defined? - should it be dfun2?

You say you had to add z to the start vector init2, but you've defined it
as a fixed effect when it's random. Try defining init2 as a list with
components nlpars for the fixed effects p1, p2, p3 and p4, and theta for
var(z), and it might work better (note the error message refers to theta).

Best wishes,
Tim Cole

--- 
Tim.cole at ucl.ac.uk Phone +44(0)20 7905 2666 Fax +44(0)20 7905 2381
Population, Policy and Practice Programme
UCL Great Ormond Street Institute of Child Health, London WC1N 1EH, UK





>Date: Thu, 19 Jan 2017 12:08:51 -0600
>From: "Therneau, Terry M., Ph.D." <therneau at mayo.edu>
>To: r-sig-mixed-models at r-project.org
>Subject: [R-sig-ME] nlme vs nlme4
>
>I have a fit that works in nlme and fails in nlme4; I'd like to convert
>to the latter 
>since I now need to use case weights and varFixed() + nlme leads to a lot
>of "wrong 
>length" warnings.  (Google search suggests that this has sometimes been
>asked about but I 
>don't find any anwsers).
>
>As always, my first hypothesis is that I have some silly setup mistake.
>
># nlme code
>tfun2 <- function(age, p1, p2,p3,p4, z)
>     p1 + p2*(age+z- 60) + .5*(p4-p2)*((age+z-p3) + sqrt((age+z-p3)^2
>+20))
>
>dfun2 <- deriv(body(tfun2), c("p1", "p2", "p3", "p4", "z"),
>               function.arg=tfun2)
>init <- c(p1=1.2, p2=.03, p3=1, p4=1)
>
>nfit1 <- nlme(y ~ dfun(age, p1, p2, p3, p4, z),
>               fixed= p1 + p2 + p3 + p4  ~ 1,
>               random= z~ 1|clinic, data=mydata, start=init,
>               na.action=na.omit)
>nfit1
>   Fixed: p1 + p2 + p3 + p4 ~ 1
>          p1          p2          p3          p4
>  1.21516330  0.00237005 73.58059382  0.08424922
>
>Random effects:
>  Formula: z ~ 1 | clinic
>                z   Residual
>StdDev: 9.418258 0.06865345
>
>---------------
>
># nlmer code
>init2 <- c(p1=1.2, p2=.03, p3=1, p4=1, z=0)
>mfit1 <- nlmer(y ~ dfun(age, p1, p2, p3, p4, z) ~
>               0+ p1 + p2 + p3 + p4 + (0+z | clinic),
>                data=mydata, start=init2, verbose=TRUE)
>Error in devfun(rho$pp$theta) : Downdated VtV is not positive definite
>
>------------------
>
>I had to add z=0 to the start vector or else it complains that z is not
>found, apparently 
>the derivative component of dfun is not enough to inform.
>
>For those who want to know more:  This is essentially a change point
>model fit with a 
>hyperbola.  The biomarker in question rises slowly with age, but for a
>substantial subset
>there is a fairly abrupt upward change in the slope sometime between the
>ages of 60 and 
>90.  The trend line is 1.2 + .0024*(age-60), the average age of turn is
>73.5, and the new 
>slope is .086, and the std of turning age is 9.4.  The constant "20" in
>the formula is an 
>arbitrary value that controls how sharp a turn the hyperbola makes.
>
>I can provide an anonymized data set to someone (if this isn't just a
>silly error).  There
>are about 3000 data points.
>
>Terry Therneau
>
>>



More information about the R-sig-mixed-models mailing list