Re: [R] Problems using quadprog for solving quadratic programming problem

From: Berwin A Turlach <berwin_at_maths.uwa.edu.au>
Date: Tue 06 Jun 2006 - 21:05:29 EST

G'day Fabian,

>>>>> "FB" == Fabian Barth <Fabian.Barth@web.de> writes:

    FB> I'm using the package quadprog to solve the following     FB> quadratic programming problem.

    FB> I want to minimize the function 
    FB> (b_1-b_2)^2+(b_3-b_4)^2
    FB> by the following constraints b_i, i=1,...,4:


    FB> b_1+b_3=1
    FB> b_2+b_4=1

    FB> 0.1<=b_1<=0.2
    FB> 0.2<=b_2<=0.4
    FB> 0.8<=b_3<=0.9
    FB> 0.6<=b_4<=0.8

    FB> In my opinion the solution should be b_1=b_2=0.2 und b_3=b_4=0.8. This could well be the correct solution. For these values for the b_i the function evaluates to zero and you definitely won't get anything smaller than that. :)

    FB> Unfortunately R doesn't find this solution and what's
    FB> surprising to me, evaluation the solution of solve.QP with my
    FB> function doesn't lead the minimal "value" calculated by
    FB> solve.QP

    FB> I would be very happy, if anyone could help and tell me,     FB> where's my mistake.
Mmh, the first code that you are sending seems to involve only three unknowns. Did you try to eliminate one by hand? If so, which one. Also, in that part of the code you do not tell solve.QP that there are equality constraints. So I am a bit confused what problem the first part is supposed to solve. :-)

The second code, labelled `my "non working" sample' seems to implement the above problem directly. The problem with that code is that the matrix Dmat that you construct is not symmetric, there is the (apparently undocumented and unchecked) assumption that Dmat is symmetric. (Note if the matrix D in the quadratic form is not summetric, then you can always replace it by (D+D^T)/2 without changing the objective function. Thus, without loss of generality, the matrix D in the objective function is always assumed to be symmetric). So adding code like

          Dmat <- (Dmat+t(Dmat))/2
in front of your `print(Dmat)' statement would be necessary.

Two minor observations. First, note that the help page for solve.QP states that the function to be minimised is

        -d^T b + 1/2 b^T D b
where D is passed in Dmat and d is passed in dvec. In your case it might not matter, but strictly speaking, you should multiply your Dmat by 2 to take care of the 1/2 factor in the definition of the objective function. Secondly, the FORTRAN code that is used by solve.QP only uses the upper triangular part of the matrix Dmat that is passed to solve.QP. Thus, essentially you were minimizing

        min 1/2 (b_1^2 + b_2^2 + b_3^2 + b_4^2) under all the conditions stated. And, for that problem solve.QP gave the correct answer.

Once you correct the problem, and a faster way of coding your problem is probabably the following:

> Dmat1 <- matrix(c(2, -2, -2, 2), 2, 2)
> Dmat2 <- matrix(0, 2, 2)
> Dmat <- rbind(cbind(Dmat1, Dmat2), cbind(Dmat2, Dmat1))
> dvec <- rep(0,4)
> Amat <- cbind(c(1,0,1,0), c(0,1,0,1), diag(4), -diag(4))
> bvec <- c(1, 1, 0.1, 0.2, 0.8, 0.6, -0.2, -0.4, -0.9, -0.8)

you will see that there is still a problem:

> solve.QP(Dmat,dvec,Amat,bvec=bvec, meq=2)
Error in solve.QP(Dmat, dvec, Amat, bvec = bvec, meq = 2) :

        matrix D in quadratic function is not positive definite!

The Goldfarb-Idnani algorithm starts off by calculating the unconstrained solution. Thus, it requires that the matrix D in the objective function is positive definite. In your problem the matrix is indefinite.

Hope this helps.

Cheers,

        Berwin


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 Received on Wed Jun 07 00:38:04 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Wed 07 Jun 2006 - 02:10:36 EST.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.