Paul,
Here is another approach: I wrote an R function that would generate interior
points as starting values for constrOptim. This might work better than the LP
approach, since the LP approach gives you a starting value that is on the
boundary of the feasible region, i.e a vertex of the polyhedron, whereas this
new approach gives you points on the interior. You can generate as many points
as you wish, but the approach is brute-force and is very inefficient - it takes
on the order of a 1000 tries to find one feasible point.
Hope this helps,
Ravi.
----- Original Message -----
From: Paul Smith <[EMAIL PROTECTED]>
Date: Tuesday, July 3, 2007 7:32 pm
Subject: Re: [R] Fine tunning rgenoud
To: R-help <r-help@stat.math.ethz.ch>
> On 7/4/07, Ravi Varadhan <[EMAIL PROTECTED]> wrote:
> > It should be easy enough to check that your solution is valid (i.e.
> a local
> > minimum): first, check to see if the solution satisfies all the
> > constraints; secondly, check to see if it is an interior point
> (i.e. none of
> > the constraints become equality); and finally, if the solution is an
> > interior point, check to see whether the gradient there is close to
> zero.
> > Note that if the solution is one of the vertices of the polyhedron,
> then the
> > gradient may not be zero.
>
> I am having bad luck: all constraints are satisfied, but the solution
> given by constrOptim is not interior; the gradient is not equal to
> zero.
>
> Paul
>
>
> > -----Original Message-----
> > From: [EMAIL PROTECTED]
> > [ On Behalf Of Paul Smith
> > Sent: Tuesday, July 03, 2007 5:10 PM
> > To: R-help
> > Subject: Re: [R] Fine tunning rgenoud
> >
> > On 7/3/07, Ravi Varadhan <[EMAIL PROTECTED]> wrote:
> > > You had indicated in your previous email that you are having trouble
> > finding
> > > a feasible starting value for constrOptim(). So, you basically
> need to
> > > solve a system of linear inequalities to obtain a starting point.
> Have
> > you
> > > considered using linear programming? Either simplex() in the "boot"
> > package
> > > or solveLP() in "linprog" would work. It seems to me that you
> could use
> > any
> > > linear objective function in solveLP to obtain a feasible
> starting point.
> > > This is not the most efficient solution, but it might be worth a
> try.
> > >
> > > I am aware of other methods for generating n-tuples that satisfy
> linear
> > > inequality constraints, but AFAIK those are not available in R.
> >
> > Thanks, Ravi. I had already conceived the solution that you suggest,
> > actually using "lpSolve". I am able to get a solution for my problem
> > with constrOptim, but I am not enough confident that the solution is
> > right. That is why I am trying to get a solution with rgenoud, but
> > unsuccessfully until now.
> >
> > Paul
> >
> >
> >
> > > -----Original Message-----
> > > From: [EMAIL PROTECTED]
> > > [ On Behalf Of Paul Smith
> > > Sent: Tuesday, July 03, 2007 4:10 PM
> > > To: R-help
> > > Subject: [R] Fine tunning rgenoud
> > >
> > > Dear All,
> > >
> > > I am trying to solve the following maximization problem, but I cannot
> > > have rgenoud giving me a reliable solution.
> > >
> > > Any ideas?
> > >
> > > Thanks in advance,
> > >
> > > Paul
> > >
> > > ----------------------------
> > > library(rgenoud)
> > >
> > > v <- 0.90
> > > O1 <- 10
> > > O2 <- 20
> > > O0 <- v*O1+(1-v)*O2
> > >
> > > myfunc <- function(x) {
> > > U0 <- x[1]
> > > U1 <- x[2]
> > > U2 <- x[3]
> > > q0 <- x[4]
> > > q1 <- x[5]
> > > q2 <- x[6]
> > > p <- x[7]
> > >
> > > if (U0 < 0)
> > > return(-1e+200)
> > > else if (U1 < 0)
> > > return(-1e+200)
> > > else if (U2 < 0)
> > > return(-1e+200)
> > > else if ((U0-(U1+(O1-O0)*q1)) < 0)
> > > return(-1e+200)
> > > else if ((U0-(U2+(O2-O0)*q2)) < 0)
> > > return(-1e+200)
> > > else if ((U1-(U0+(O0-O1)*q0)) < 0)
> > > return(-1e+200)
> > > else if ((U1-(U2+(O2-O1)*q2)) < 0)
> > > return(-1e+200)
> > > else if((U2-(U0+(O0-O2)*q0)) < 0)
> > > return(-1e+200)
> > > else if((U2-(U1+(O1-O2)*q1)) < 0)
> > > return(-1e+200)
> > > else if(p < 0)
> > > return(-1e+200)
> > > else if(p > 1)
> > > return(-1e+200)
> > > else if(q0 < 0)
> > > return(-1e+200)
> > > else if(q1 < 0)
> > > return(-1e+200)
> > > else if(q2 < 0)
> > > return(-1e+200)
> > > else
> > >
> >
> return(p*(sqrt(q0)-(O0*q0+U0))+(1-p)*(v*(sqrt(q1)-(O1*q1+U1))+(1-v)*(sqrt(q2
> > > )-(O2*q2+U2))))
> > >
> > > }
> > >
> >
> genoud(myfunc,nvars=7,max=T,pop.size=6000,starting.values=runif(7),wait.gene
> > > rations=150,max.generations=300,boundary.enforcement=2)
> > >
> > > ______________________________________________
> > > R-help@stat.math.ethz.ch mailing list
> > >
> > > PLEASE do read the posting guide
> >
> > > and provide commented, minimal, self-contained, reproducible code.
> > >
> >
> > ______________________________________________
> > R-help@stat.math.ethz.ch mailing list
> >
> > PLEASE do read the posting guide
> > and provide commented, minimal, self-contained, reproducible code.
> >
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
>
> PLEASE do read the posting guide
> and provide commented, minimal, self-contained, reproducible code.
lineq <- function(a){
r <- rep(NA,14)
v <- a[1]
O1 <- a[2]
O2 <- a[3]
O0 <- v*O1 + (1-v)*O2
r <- rep(-1,14)
kount <- 0
while (any(r < 0)) {
kount <- kount + 1
x <- runif(7)
U0 <- x[1]
U1 <- x[2]
U2 <- x[3]
q0 <- x[4]
q1 <- x[5]
q2 <- x[6]
p <- x[7]
r[1] <- U0
r[2] <- U1
r[3] <- U2
r[4] <- U0 - (U1 + (O1 - O0)*q1)
r[5] <- U0 - (U2 + (O2 - O0)*q2)
r[6] <- U1 - (U0 + (O0 - O1)*q0)
r[7] <- U1 - (U2 + (O2 - O1)*q2)
r[8] <- U2 - (U0 + (O0 - O2)*q0)
r[9] <- U2 - (U1 + (O1 - O2)*q1)
r[10] <- p
r[11] <- 1-p
r[12] <- q0
r[13] <- q1
r[14] <- q2
}
list(x0=x, kount=kount, r=sign(r))
}
a <- c(0.9,10,20)
ans <- lineq(a)
ans$x0 # interior point
______________________________________________
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
and provide commented, minimal, self-contained, reproducible code.
