From bolker at zoo.ufl.edu  Wed Jan 18 00:02:09 2006
From: bolker at zoo.ufl.edu (Ben Bolker)
Date: Tue, 17 Jan 2006 18:02:09 -0500
Subject: [Rd] nls profile with port/constraints
Message-ID: <43CD7771.60506@zoo.ufl.edu>


   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))


-- 
620B Bartram Hall                            bolker at zoo.ufl.edu
Zoology Department, University of Florida    http://www.zoo.ufl.edu/bolker
Box 118525                                   (ph)  352-392-5697
Gainesville, FL 32611-8525                   (fax) 352-392-3704



[Rd] nls profile with port/constraints
Ben Bolker