[R] logistic transformation using nlminb

From: John Pitchard <johnpitchard_at_googlemail.com>
Date: Fri, 16 May 2008 16:47:02 +0100


Dear all,

I want to find the optimal values of a vector, x (with 6 elements) say, satisfying the following conditions:

  1. for all x>=0
  2. sum(x)=1
  3. x[5]<=0.5 and x[6]<=0.5

For the minimisation I'm using nlminb and to satisfy the first 2 conditions the logistic transformation is used with box constraints for condition 3. However, I don't seem to be able to get the values x that i'm expecting from fitting a simpler model. For example,

set.seed(0)
# creating a correlation matrix

corr <- diag(5)
corr[lower.tri(corr)] <- 0.5
corr[upper.tri(corr)] <- 0.5

# Data for the minimisation

mat <- rmvnorm(10000, mean=c(3, 2.75, 2.8, 3, 2.9), sd=c(0.1, 0.1,0.12, 0.15, 0.10), cov=corr)

# here is the simple optimisation function that allows the 5th
element to be potentially negative

obj.funA <- function(opt, mat) {

      opt <- c(opt, 1-sum(opt))
      LinearComb <- mat%*%opt

 obj <- -min(LinearComb)
 obj
}

# I want to put an upper boundary constraint on the first variable -
not being greater than 0.35 and the rest being between 0 and 1

opt <- nlminb(rep(0,4), lower=rep(0,4), upper=c(0.35, 1, 1, 1) ,

obj.funA, mat=mat)
opt.x <- opt$parameters
opt.x <- c(opt.x, 1-sum(opt.x))
opt.x

# Result

[1] 0.34999902 0.06475651 0.00000000 0.16561760 0.41962686

However, I don't get the same result from the logistic transformation

obj.funB <- function(opt, mat) {
 opt <- c(opt, 0)
 opt <- exp(opt)/sum(exp(opt))

      LinearComb <- mat%*%opt
 obj <- -min(LinearComb)
 obj
}

opt <- nlminb(rep(0,4), upper=c(log(0.35), NA, NA, NA), obj.funB, mat=mat)

opt.x <- opt$parameters
opt.x <- c(opt.x, 0)
opt.x <- exp(opt.x)/sum(exp(opt.x))
opt.x

# Result

[1] 2.355325e-001 6.339398e-009 1.202751e-004 9.139718e-002 6.729500e-001

I don't know how to obtain the same answer for both optimisations. In reality, my own optimisation typically gives negative values for the standard minimisation- so I have no choice but to use a more advanced method. Also, there appears to be a dependency between the first and the last, i.e. 2.355325e-001/6.729500e-001 =0.35

Does anyone know why the logistic doesn't give the same answer as the simpler method?

Any help is much appreciated.

Regards,
John



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 Fri 16 May 2008 - 16:52:09 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 Fri 16 May 2008 - 17:30:37 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.

list of date sections of archive