Thanx Dimitris, Patrick and Berwin!
For other people interested in this problem, here are two valid solutions that work.
1) Re-parameterize e.g.,
EM <- c(100,0,0,0,100,0,0,0,100) W <- array(EM, c(3,3)) d <- c(10, 20, 70)
fn <- function(x){
x <- exp(x) / sum(exp(x))
r <- W%*%x - d
crossprod(r, r)[1,1]
}
opt <- optim(rnorm(3), fn)
res <- exp(opt$par) / sum(exp(opt$par))
res"The first line of the `fn()' function is just a re-pameterization of your problem, i.e., if `y' is a vector of real numbers, then it is straightforward to see that `x = exp(y) / sum(exp(y))' will be real numbers in (0, 1) for which `sum(y)=1'. So instead of finding xs that minimize your function under the constraint (which is more difficult) you just find the ys using the above transformation." (I owe you a drink Dimitris !!!)
2) Or minimize it as a quadratic function under a linear constraint:
EM <- c(100,0,0,0,100,0,0,0,100) W <- array(EM, c(3,3)) d <- c(10, 20, 70) library(quadprog) Dmat <- crossprod(W,W) dvec <- crossprod(d,W) A <-matrix(c(1,1,1),ncol=1) bvec <- 1 solve.QP(Dmat, dvec, A, bvec, meq=1)
This is based on the objective function (i.e. the thing you want to minimise) : min x'C'Cx - 2 d'Cx + d'd where sum(x) = 1 (Thanx Berwin!!)
Kind regards, Stef
______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
