# [Rd] constrained optimiser doesn't obey constraints (PR#411)

thomas@biostat.washington.edu thomas@biostat.washington.edu
Fri, 4 Feb 2000 17:46:05 +0100 (MET)

```
---------- Forwarded message ----------
Date: Fri, 4 Feb 2000 01:58:44 +0000 (GMT)
From: Thomas Lumley <thomas@biostat.washington.edu>
To: r-bugs@biostat.ku.dk
Subject: constrained optimiser doesn't obey constraints

The constrained optimiser is calling the function with parameter values
outside the constrained range

In this example the box is from c(0,0,0,0) to c(1,1,.5,1). The function
prints out  the parameter vector when it is out of range. This gives:

[1] "x= -0.001" "x= 1"      "x= 0"      "x= 1"
[1] "x= 0"     "x= 1.001" "x= 0"     "x= 1"
[1] "x= 0"      "x= 1"      "x= -0.001" "x= 1"
[1] "x= 0"     "x= 1"     "x= 0"     "x= 1.001"
[1] "x= NaN" "x= NaN" "x= NaN" "x= NaN"
[1] "x= NaN" "x= NaN" "x= NaN" "x= NaN"
repeated a number of times.

For these starting value it eventually reports convergence, for
some others it goes forever.

---Example code---

revise.mix.binbb_function(x,loss,no.inform){
if (any(is.nan(x)) || any(x<0) || any(x>c(1,1,.5,1))){
print(paste("x=",x))
return(Inf)
}
eta_x[1]
pi1_x[2]
omega1_x[3]
p0_x[4]
log.term_0
for (i in 1:length(no.inform)){
if ( (loss[i]-1) >= 0) {
first.prod_1
for (r in 0:(loss[i]-1)){
first.prod_first.prod*(pi1+r*omega1)
};} else {first.prod_1;}
if ( (no.inform[i]-loss[i]-1) >=0) {
second.prod_1
for (r in 0:(no.inform[i]-loss[i]-1)) {
second.prod_second.prod*(1-pi1+r*omega1)
};} else {second.prod_1;}

third.prod_1
for (r in 0:(no.inform[i]-1)){
third.prod_third.prod*(1+r*omega1)
}

common.term<-choose(no.inform[i],loss[i])
first.part_eta*common.term*first.prod*second.prod/third.prod

whole.term_eta*common.term*first.prod*second.prod/third.prod+(1-eta)*common.term*p0^loss[i]*(1-p0)^(no.inform[i]-loss[i])

log.term_log.term+log(whole.term)

}
## log.term is the log-likelihood.
## For nlminb we need to minimize -log-likelihood which is
## equivalent to maximizing the log-likelihood
return(-log.term)
}

init<-c(0.001, 0.001, 0.001, 0.220)

"sim.ns" <-
c(17, 17, 15, 15, 16, 18, 18, 19, 17, 12, 18, 19, 19, 15, 20,
17, 7, 15, 18, 15, 19, 19, 19, 21, 16, 17, 18, 19, 21, 12, 16,
16, 21, 14, 16, 8, 18, 15)
"loss.dat" <-
structure(c(17, 13, 8, 9, 10, 4, 1, 7, 3, 4, 2, 6, 3, 4, 7, 3,
1, 2, 5, 4, 3, 4, 2, 1, 1, 7, 4, 5, 3, 0, 4, 6, 4, 1, 6, 3, 2,
3), .Names = c("V1", "V2", "V3", "V4", "V5", "V6", "V7", "V8",
"V9", "V10", "V11", "V12", "V13", "V14", "V15", "V16", "V17",
"V18", "V19", "V20", "V21", "V22", "V23", "V24", "V25", "V26",
"V27", "V28", "V29", "V30", "V31", "V32", "V33", "V34", "V35",
"V36", "V37", "V38"))

optim(init,fn=revise.mix.binbb,
lower=c(0.00,0.00,0.00,0.00),upper=c(1,1,.5,1),
no.inform=sim.ns,loss=loss.dat,method="L-BFGS-B")

--please do not edit the information below--

Version:
platform = i686-unknown-linux
arch = i686
os = linux
system = i686, linux
status = Under development (unstable)
major = 0
minor = 99.0
year = 2000
month = January
day = 29
language = R

Search Path: