Skip to content

nls profiling with algorithm="port" may violate bounds (PR#8508)

4 messages · Spencer Graves, Ben Bolker

#
[posted to R-devel, no discussion:
resubmitting it as a bug, just so it gets
logged appropriately]

   Sorry to report further difficulties with
nls and profiling and constraints ... the problem
this time (which I didn't check for in my last
round of testing) is that the nls profiler doesn't
seem to respect constraints that have been
set when using the port algorithm.
    See test code below ...
    If I can I will try to hack the code, but I will
probably start by redefining my function with
some workarounds to make the fit quadratically "bad" (but well-defined)
when the parameters are negative ...
     As always, please don't hesitate to correct me
if I'm being an idiot ...

    cheers
      Ben Bolker

-----------------------
rm(list=ls())

npts=10
set.seed(1001)

a =2
b =0.5
x= runif(npts)
y = a*x/(1+a*b*x)+rnorm(npts,sd=0.2)

gfun <- function(a,b,x) {
    if (a<0 || b<0) stop("bounds violated")
    a*x/(1+a*b*x)
}

m1 = nls(y~gfun(a,b,x),algorithm="port",
    lower=c(0,0),start=c(a=1,b=1))

try(confint(m1))
----------------

for what it's worth, the logic appears to be OK in mle in the stats4
library:
--------------
library(stats4)

mfun <- function(a,b,s) {
    if (a<0 || b<0 || s<0) stop("bounds violated")
    -sum(dnorm(y,a*x/(1+a*b*x),sd=s,log=TRUE))
}

m2 = mle(minuslogl=mfun,
    start=list(a=1,b=1,s=0.1),
    method="L-BFGS-B",lower=c(0.002,0.002,0.002))

confint(m2)

m2b = mle(minuslogl=mfun,
    fixed=list(b=0),start=list(a=1,s=0.1),
    method="L-BFGS-B",lower=c(0.002,0.002,0.002))
## set boundary slightly above zero to avoid
## boundary cases

dev <- 2*(-logLik(m2b)+logLik(m2))
as.numeric(pchisq(dev,lower.tail=FALSE,df=1))
1 day later
#
Hi, Ben, et al.:

	  The issue Ben identified with confint(nls(... )) generates a hard 
failure for me.  Specifically the command "confint(m1)" in his script 
below under Rgui 2.2.1 first says, "Waiting for profiling to be done..." 
then forces a screen to pop up with heading "R for Windows GUI 
front-end" reading, "R for Windows GUI front-end has encountered a 
problem and needs to close.  We are sorry for the inconvenience.  If you 
were in the middle of something, the information you were working on 
might be lost... ."  When I try the same thing running R under XEmacs 
with ESS, I get essentially the same response, exceppt "R for Windows 
GUI" is replaced by "R for Windows terminal".  In both cases, it kills 
R.  In both cases, sessionInfo() before this command is as follows:

 > sessionInfo()
R version 2.2.1, 2005-12-20, i386-pc-mingw32

attached base packages:
[1] "stats4"    "methods"   "stats"     "graphics"  "grDevices" "utils"
[7] "datasets"  "base"

	  I'm running Windows XP professional version 5.1 on an IBM T30 
notebook computer.

	  Thanks to all of the R Core team for all your hard work to make R 
what it is today, with these kinds of unpleasant surprises to rare.

	  Best Wishes,
	  spencer graves
bolker at zoo.ufl.edu wrote:

            
#
Spencer Graves <spencer.graves <at> pdf.com> writes:
"We" (being Brian Ripley and I) know about this already.
  I'm sorry I failed to specify enough info in my bug report,
but I was using R-devel/2.3.0 of (I think?) 11 January, under Linux.
Your problem is actually PR #8428, which is fixed enough to prevent
a crash in 2.2.1 patched and really fixed in 2.3.0, all thanks to
Brian Ripley.

  cheers
    Ben
#
Dear Professors. Bolker & Ripley:

       Thank you both very much for all your creativity and hard work 
both in your generall contributions to human knowledge and specifically 
for helping make R the great thing it is today.  I had not seen a reply 
to that email in several days, so I made time to check it out.  When I 
replicated the error, I thought it my duty to report same.  I know 
that's more appropraite with "r-help" than with "r-devel", but I thought 
such a comment might help someone.  I certainly did NOT want to add to 
someones' workload a requirement to reply to my comment.

       Thanks again,
       spencer graves
Ben Bolker wrote: