[R] not positive definite D matrix in quadprog

From: Molins, Jordi <Jordi.Molins_at_drkw.com>
Date: Wed 01 Sep 2004 - 20:53:09 EST

Hello to everybody,

I have a quadratic programming problem that I am trying to solve by various methods. One of them is to use the quadprog package in R.

When I check positive definiteness of the D matrix, I get that one of the eigenvalues is negative of order 10^(-8). All the others are positive. When I set this particular eigenvalue to 0.0 and I recheck the eigenvalues in R, the last eigenvalue is positive of order 10^(-12). I try to use solve.QP, but I get an error message that matrix D in quadratic function is not positive definite. For reference, a fully R session is listed below.

Is 10^(-12) too close to 0? i.e., does R consider that with an eigenvalue of order +10^(-12) the matrix is not positive definite but positive semidefinite?

In general, has somebody know a way (in R or outside R, maybe in c++) to solve quadratic programming with positive semidefinite matrices?

In particular, my problem is not so hard: given y an n x 1 matrix, and beta an n x m matrix, I want to find an m x 1 matrix x s. t. sum(y - beta * x)^2 is minimum. The particularity is that I want to impose restrictions on x: all x components should be between 0 and 1, and there are also constraints of the type A x = b, where A and b have the necessary dimensions to ensure consistency.

I have tried with some other packages, and they do not give a correct solution when the system increases in size (e.g., 24 variables and 9 constraint equations) ... some idea?

thanks!

Jordi

---

The problem:
 library(MASS)
 library(quadprog)

 D <-
matrix(c(439.5883658,438.8445615,438.1007572,2430.285506,2426.162884,2422.04 0262,44.21800696,44.14261394,  

438.8445615,438.1020157,437.3594699,2426.173348,2422.057702,2417.942056,44.1 43188,44.06792255,  

438.1007572,437.3594699,449.6727418,2445.212326,2542.83573,2643.780669,50.19 455336,52.04059805,  

2430.285506,2426.173348,2445.212326,13491.19467,13614.55046,13737.90625,253. 4897678,255.745654,  

2426.162884,2422.057702,2542.83573,13614.55046,14687.86142,15923.99043,313.8 180838,336.4239658,  

2422.040262,2417.942056,2643.780669,13737.90625,15923.99043,19107.7405,410.9 729841,472.5104919,  

44.21800696,44.143188,50.19455336,253.4897678,313.8180838,410.9729841,9.5462 51262,11.57677661,  

44.14261394,44.06792255,52.04059805,255.745654,336.4239658,472.5104919,11.57 677661,14.51245153),8,8)

 D.vectors <- eigen(D,only.values=F)$vectors  D.values <- eigen(D,only.values=F)$values

#the last value is negative

 D.values
[1] 4.609489e+04 2.458166e+03 8.232288e+01 1.961199e+00 5.976441e-01 [6] 2.810968e-01 1.253157e-09 -2.685763e-08

 D.quad <- matrix(0,8,8)
 diag(D.quad) <- D.values

#checking that the eigenvalue decomposition works fine
 D.vectors%*%D.quad%*%ginv(D.vectors)

 D.quad[8,8]
[1] -2.685763e-08
 D.quad[8,8] <- 0.0

#checking; nothing changes too much

D.vectors%*%D.quad%*%ginv(D.vectors)

#now all eigenvalues are positive:

 D.values.new <-
eigen(D.vectors%*%D.quad%*%ginv(D.vectors),only.values=F)$values  D.values.new
[1] 4.609489e+04 2.458166e+03 8.232288e+01 1.961199e+00 5.976441e-01 [6] 2.810968e-01 1.253140e-09 1.428534e-12

Dmat <- D.vectors%*%D.quad%*%ginv(D.vectors)

dvec <-
-c(-2910.533769,-2905.609008,-3012.223863,-16274.97455,-17222.46423,-18380.6 391,-357.8878464,-379.6371849)

#this ensures that coefficients are positive:
 Amat <- matrix(0,8,8)
 diag(Amat) <- 1
 bvec <- rep(0,8)

#it says D is not positive definite ...
 solve.QP(Dmat,dvec,Amat,bvec=bvec)
Error in solve.QP(Dmat, dvec, Amat, bvec = bvec) :

        matrix D in quadratic function is not positive definite!

---

The information contained herein is confidential and is inte...{{dropped}}

---

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 Sep 01 21:04:50 2004