... and getting an answer is no assurance that the answer is meaningful. In cases like this (which arise frequently because of insistence on using the accepted mechanistic model paradigm even in the absence of informative data to estimate it), the confidence intervals for (correlated) parameters will usually be so wide as to be useless. That is, for practical purposes, the model is nonidentifiable. This can easily be seen by making small random perturbations in the data and watching how the parameter estimates change -- i.e. performing a sensitivity analysis. Incidentally, the predictions will typically be fine, so the standard scientific practice of graphing the data with an overlaid smooth curve as a "check" on the validity of the estimated parameters is nonsense.
One should not get so lost among the trees of statistical niceties that one loses sight of the scientific forest: the model can't be adequately fit by the data. Either get more data or choose a more appropriate model. Cheers, Bert Gunter Genentech Nonclinical Biostatistics -----Original Message----- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of Ravi Varadhan Sent: Tuesday, July 14, 2009 4:35 PM To: 'UyenThao Nguyen'; 'spencerg' Cc: r-help@r-project.org Subject: Re: [R] nls, reach limit bounds I took a quick look at "drc"package and the "drm" function. The drm() function uses optim ("BFGS" method). So, that is one diffference. However, without looking at your code on how you used drm(), I cannot tell further. The fact that you got an answer using optim() does not necessarily mean that everything is swell. Did you check the Hessian to see if it is positive-definite? You might also want to try the nls.lm() function in the "minpack.lm" package. Ravi. ---------------------------------------------------------------------------- ------- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: rvarad...@jhmi.edu Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty_personal_pages/Varadhan.h tml ---------------------------------------------------------------------------- -------- -----Original Message----- From: UyenThao Nguyen [mailto:ungu...@tethysbio.com] Sent: Tuesday, July 14, 2009 7:07 PM To: Ravi Varadhan; 'spencerg' Cc: r-help@r-project.org Subject: RE: [R] nls, reach limit bounds Hi Ravi and Spencer, Thank you very much for your help. I did plot the data, and saw that the data didn't seem to have an inflection point. Yes, the data contained 6 points of duplicates, which the 4 P logistic regression is appropriate to use. I tried the dose response model (drm in drc package), this model works without any problem. Do you know if the drm has different tolerance in convergent or something else? Let's say if I choose drm to fit the data, Can I get the parameters in the same way nls gives me? I tested a converged data set on both drm and nls, and I can't see any relationship between their parameters although the fits were similar. Thank you. Uyen -----Original Message----- From: Ravi Varadhan [mailto:rvarad...@jhmi.edu] Sent: Monday, July 13, 2009 3:32 PM To: 'spencerg'; UyenThao Nguyen Cc: r-help@r-project.org Subject: RE: [R] nls, reach limit bounds Hi Uyen, You do not have enough information to estimate 4 parameters in your nonlinear model. Even though you have 12 data points, only 6 are contributing independent information (you essentially have 6 replicate points). If you plot the fittted function you will see that it fits your data really well. However, you will also see that the fitted curve is far from reaching the asymptote. An implication of this is that you cannot reliably estimate b0 and b1 separately. So, you need more data, especially for larger x values in the asymptotic region. Ravi. ---------------------------------------------------------------------------- ------- Ravi Varadhan, Ph.D. Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University Ph: (410) 502-2619 Fax: (410) 614-9625 Email: rvarad...@jhmi.edu Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty_personal_pages/Varadhan.h tml ---------------------------------------------------------------------------- -------- -----Original Message----- From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On Behalf Of spencerg Sent: Saturday, July 11, 2009 9:50 PM To: UyenThao Nguyen Cc: r-help@r-project.org Subject: Re: [R] nls, reach limit bounds Have you plotted the data? There are two standard sources of non-convergence problems like this: First, there may not be enough information in your data to estimate all four parameters. Second, if that's not the case, then your starting values are not appropriate. I routinely use "optim" or "nlminb" to find a sensible solution, because these general purpose optimizers are more likely to converge than "nls". To do this, I write a function with a name like "SSElogistic" to compute the sum of squares of residuals for your model. I like to use "optim" with "hessian=TRUE". Then I compute "eigen(fit$hessian, symmetric=TRUE)", with "fit" = the object returned by "optim". If the smallest eigenvalue is negative, it says that optim found a saddle point. If the smallest eigenvalue is less than, e.g., 1e-8 times the largest, it says that the smallest eigenvector is very poorly estimated. Round those numbers off grossly to 1 significant digit, and they will likely suggest which parameter you can fix and drop from the model. Hope this helps. Spencer Graves UyenThao Nguyen wrote: > Hi, > > I am trying to fit a 4p logistic to this data, using nls function. The function didn't freely converge; however, it converged if I put a lower and an upper bound (in algorithm port). Also, the b1.A parameter always takes value of the upper bound, which is very strange. Has anyone experienced about non-convergent of nls and how to deal with this kind of problem? > > Thank you very much. > > > > ########################################################################3 > y x > 1 0.8924619 -0.31875876 > 2 1.1814749 -0.21467016 > 3 1.6148266 0.06069784 > 4 2.2091363 0.54032947 > 5 2.7019079 1.04921802 > 6 3.0679585 1.60745502 > 9 0.9436973 -0.31875876 > 10 1.2201133 -0.21467016 > 11 1.6470043 0.06069784 > 12 2.2090048 0.54032947 > 13 2.6864857 1.04921802 > 14 3.0673523 1.60745502 > > new.cont=nls.control(maxiter = 10000, tol = 1e-05, minFactor = 1e-08, > printEval = FALSE, warnOnly = FALSE) > > > b0.A=.9*min(DAT$y) > b1.A=1.1*max(DAT$y)-b0.A > b2.A=-1*mean(DAT$x) > b3.A=1 > > > b0.A > b1.A > b2.A > b3.A > > nls.mdl.A=nls(y~b0 + b1/(1+exp(-b2-b3*x)),data=DAT,start = > list(b0=b0.A, b1=b1.A, b2=b2.A, b3=b3.A), lower=-10, upper=10, > algorithm="port",trace=T,control=new.cont) > > ################################## > > > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help@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@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@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@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.