From: Gabor Grothendieck <ggrothendieck_at_gmail.com>

Date: Fri, 21 Dec 2007 21:26:37 -0500

R-help_at_r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Sat 22 Dec 2007 - 02:30:11 GMT

Date: Fri, 21 Dec 2007 21:26:37 -0500

Check out the drc package.

Also, please read the last line to every message to r-help regarding
reproducible code

and note that function names are case sensitive in R so NLS is not the
same as nls.

On Dec 21, 2007 7:36 PM, <rsposto_at_yahoo.com> wrote:

> I am writing a program for automated (i.e. no user intervention - for biologists) iteratively reweighted least-square fit to a function of the form "reading ~ exp(lm2)/(1 + (dose/exp(lm3))^exp(lm4)" using case weight proportional to the mean, e.g., E(reading). Because for some datasets the solution is sensitive to starting values, I first use OPTIM() with Nelder-Mead to locate the solution, and then plug the solution into NLS() (default algorithm) using the appropriate weights as a way to retrieve SEs, deviance, etc rather than computing these from OPTIM results.

*>
**> The problem I have discovered is the following. OPTIM() arrives at the correct solution that minimizes the objective function, and I have confirmed this in numerous examples in Excel. When this solution is plugged into NLS(), the first evaluation yields the same value of the objective function as OPTIM(). NLS() then attempts steps, and apparently finds a slightly different solution with a smaller value for the objective function. However, this new solution actually does not give a smaller objective function value, but rather a larger one, when verifed against an independent computation. Hence, it appears that NLS is getting confused and taking a step based on erroneous computation of the objective function.
**>
**> An example is below, with an abridged OPTIM output followed by NLS output:
**>
**> -- ************** START OUTPUT **************************
**> -- [1] "Drug = Melphelan"
**> -- [1] "Cell line = CHLA-266"
**> -- Nelder-Mead direct search function minimizer
**> -- function value for initial parameters = 3.052411
**> -- Scaled convergence tolerance is 3.05241e-12
**> -- Stepsize computed as 1.801225
**> -- BUILD 4 288782670697.238953 3.052411
**> -- LO-REDUCTION 6 40.478286 3.052411
**> -- HI-REDUCTION 8 27.015234 3.052411
**> -- .
**> -- .
**> -- .
**> -- LO-REDUCTION 208 0.763370 0.763370
**> -- LO-REDUCTION 210 0.763370 0.763370
**> -- Exiting from Nelder Mead minimizer
**> -- 212 function evaluations used
**> -- m2 m3 m4
**> -- 1 80846342 1.190286e-05 1.049291
**>
**> [***NOTE - the correct solution parameters values are the line above on the absolute scale, which are the values below on the ln scale, representing the first evaluation by NLS. This gives the correct objective function value of 0.763370, and hence agrees with OPTIM at this point ***]
**>
**> -- 0.7633697 : 18.20806090 -11.33873192 0.04811485
**> -- It. 1, fac= 1, eval (no.,total): ( 1, 1): new dev = 0.7504
**> -- 0.7504 : 18.19909405 -11.34554803 0.05560241
**> -- It. 2, fac= 1, eval (no.,total): ( 1, 2): new dev = 0.750387
**> -- 0.7503866 : 18.19933189 -11.34796180 0.05429375
**> -- It. 3, fac= 1, eval (no.,total): ( 1, 3): new dev = 0.750387
**> -- .
**> -- . multiple step halving attempts
**> -- .
**> -- It. 14, fac= 3.8147e-06, eval (no.,total): ( 4, 39): new dev = 0.750387
**> -- It. 14, fac= 1.90735e-06, eval (no.,total): ( 5, 40): new dev = 0.750387
**> --
**> -- Formula: reading ~ exp(lm2)/(1 + (dose/exp(lm3))^exp(lm4))
**> --
**> -- Parameters:
**> -- Estimate Std. Error t value Pr(>|t|)
**> -- lm2 18.19934 0.02084 873.252 <2e-16 ***
**> -- lm3 -11.34795 0.08316 -136.459 <2e-16 ***
**> -- lm4 0.05435 0.04103 1.325 0.191
**> -- ---
**> -- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
**> --
**> -- Residual standard error: 0.1147 on 57 degrees of freedom
**> --
**> -- Number of iterations till stop: 13
**> -- Achieved convergence tolerance: 2.611e-08
**> -- Reason stopped: step factor 9.53674e-007 reduced below 'minFactor' of 1e-006
**> --
**> -- Warning message:
**> -- In nls(reading ~ exp(lm2)/(1 + (dose/exp(lm3))^exp(lm4)), start = list(lm2 =
**> -- log(beta_nm$m2), :
**> -- step factor 9.53674e-07 reduced below 'minFactor' of 1e-06
**> --
**> -- *** END OUTPUT *********************************************************
**>
**> Note the objective function value of .750387 at parameter values 18.19934, -11.34795, 0.05435. However, objective function value at this solution is actually 0.776426587, by Excel computation.
**>
**> Can anyone suggest what I might be doing wrong, or is this a bug or anomaly of the NLS algorithm?
**>
**> Thanks in advance
**>
**> Richard Sposto (Los Angeles).
**>
**>
**>
**> --
**> This message was sent on behalf of rsposto_at_yahoo.com at openSubscriber.com
**> http://www.opensubscriber.com/message/r-help@stat.math.ethz.ch/7530379.html
**>
**> ______________________________________________
**> R-help_at_r-project.org mailing list
**> https://stat.ethz.ch/mailman/listinfo/r-help
**> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
**> and provide commented, minimal, self-contained, reproducible code.
**>
*

R-help_at_r-project.org mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Sat 22 Dec 2007 - 02:30:11 GMT

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.2.0, at Sat 22 Dec 2007 - 06:30:20 GMT.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*