[R] maximization help :

[Spencer Graves]

	  Have you considered something like the following:

multilogit <- function(p){
	k <- length(p)
	z <- log(p)
	(z[-k]-z[k])
}

inv.multilogit <- function(z){
	k1 <- length(z)
	p. <- exp(z)
	p.i <- (1+sum(p.))
	(c(p., 1)/p.i)
}

multilogit(c(.1,.2,.7))
inv.multilogit(multilogit(c(.1,.2,.7)))
inv.multilogit(1:2)
multilogit(inv.multilogit(1:2))

prodSum <- function(x, A, log.=TRUE,
	trace.=FALSE, neg=TRUE){
	p <- inv.multilogit(x)
	if(trace.)cat("p =", p, ";")
	logP <- sum(log(A%*%p))
	{if(log.){
		if(neg) return(-logP)
		else return(logP)
		}
	else
		return(exp(logP))
	}
}

prodSum(1:2, diag(3), trace.=T)
sum(log(inv.multilogit(1:2)))

prodSum(0:1, diag(3), trace.=T)
sum(log(inv.multilogit(0:1)))

optim(c(0,0), fn=prodSum, hessian=TRUE, A=diag(3),
	method="CG", control=list(trace=999))
optim(1:2, fn=prodSum, hessian=TRUE, A=diag(3),
	method="CG", control=list(trace=999))

A <- array(c(1,1,1,0), dim=c(2,2))
optim(1, fn=prodSum, hessian=TRUE, A=A,
	method="CG", control=list(trace=999))

A <- array(c(1,1,0, 0, 1, 1), dim=c(3,2))
optim(1, fn=prodSum, hessian=TRUE, A=A,
	method="CG", control=list(trace=999))

There may be a more elegant solution based on singular values of A, but 
I don't see it.

	  hope this helps.
	  spencer graves

mingan yang wrote:

> 
> 
> Given a vector  : pvec=(p1,p2,.... p J)   with sum(pvec)=1,   all the 
> elements are non-negative, that is, they are probabilities
> 
> a  matrix   A  ( N* J ), with the elements  alpha(ij)  are 0 or 1
> 
> 
>    I want to MAXIMIZE THE RESULT
> 
>      RESULT=   product( i=1, to N   [ sum (  alpha(ij)* pj , j =1,to J ) 
> ]  )
> 
>    thus, I need to get pvec. how should I do ?
> 
>      for example
>               say, A=  0   1   0   0
>                 1    1  0   0
>                 1   0   0   0
>                0   0   1   0
>                1  0    0   1
>               0   0    0   1
> 
>   that is A is a matrix 6* 4    thus pvec=(p1,p2,p3,p4)
> 
>    I want to get  values of pvec such that , they can maximize
> 
>   p2 *  ( p1 + p2 ) * p1 * p3 * (p1+p4) * p4
> 
> 
>  thanks
> 
> ______________________________________________
> R-help at 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