[R] Non-linear regression problem: R vs JMP (long)

[Douglas Bates]

Andrew Robinson <andrewr at uidaho.edu> writes:

> Summary: 

> A model/data combination achieves convergence in JMP, and survives a
> reasonably rigorous examination (sensible parameter estimates,
> well-behaved surface, confidence intervals exclude 0).  The same
> combination fails to achieve convergence in R despite starting from
> the estimates reported in JMP.  What can I do?

> Full story:

> I am collaborating on a project which presently requires the fit of a seven 
> parameter non-linear function to 874 observations.  The function is:

> freeze.d <- deriv3(~ dbh^b * u^dbh * a1 * smi^a2 * exp(a3*pbal) *
>   (1 + a4*exp(a5*ba)) * cr^a6 * dh5^a7,
>   c("a1","a2","a3","a4","a5","a6","a7","b","u"),
>   function(dbh, smi, pbal, ba, cr, dh5, a1, a2, a3, a4, a5, a6, a7, b, u){})

Isn't that 9 parameters - counting b and u?

> The dbh, smi, pbal, ba, cr, and dh5, are known.

I would not be surprised at the inability to converge on a 9 parameter
nonlinear model with so many multiplicative terms.  Would it be
reasonable to try to fit the logarithm of the response to the
logarithm of the prediction instead? 

Try evaluating the gradient matrix at JMP's converged parameter
estimates and checking the condition number (using the kappa
function).  I expect it will be very large, indicating very poor
conditioning in the gradient.

> The data do not require the level of intricacy reflected in the function. It 
> is debatable whether they support it. However, my collaborator is anxious to 
> avoid linear models because of their relatively poor extrapolative 
> properties.  So, that is an ongoing discussion.

I don't think that an overparameterized nonlinear model is necessarily
the answer.

You could try to reduce the number of parameters by 1 in R using the
plinear algorithm (remove a1 as it is a simple multiplier but state
alg='plinear' in the call to nls).