[R] solve.QP with box and equality constraints

Selwyn McCracken selwyn.mccracken at gmail.com
Mon Feb 16 11:25:25 CET 2009


[reposting as a plain text - apologies for the double posting]

Dear list,

I am trying to follow an example that estimates a 2x2 markov
transition matrix across several periods from aggregate data using
restricted least squares.

I seem to be making headway using solve.QP(quadprog) as the
unrestricted solution matches the example I am following, and I can
specify simple equality and inequality constraints. However, I cannot
correctly specify a constraint matrix (Amat) with box constraints per
cell and equality constraints that span multiple cells. Namely the
solution matrix I am aiming for needs to respect the following
conditions:
  - each row must sum to 1
  - each cell must >= 0
  - each cell must <= 1

I understand the general principle of creating box constraints by
creating pairs of positive and negative values within the constraint
matrix (http://tolstoy.newcastle.edu.au/R/help/05/10/14251.html), but
when I pass this expanded constraint matrix to solve.QP it complains
that Amat and dvec are incompatible. How should I expand dvec (and
consequently Dmat) to accomodate the larger Amat? Moreover, I am
unclear how to apply the meq equality constraint across more than one
cell (i.e. rows summing to one) although I have attempted a guess
below.

Any help warmly received.
Selwyn.

################
#examples below
################

library(quadprog)
#pairs of population values over time
mat = cbind(
  cal = c(12988,13408,13834,14267,
14707,15155) ,
  rest = c(152082,154201,156352,158536,160754,163006)
  )

(X = kronecker(diag(1, ncol(mat)), mat[-nrow(mat),]))
(y = c(mat[-1,]))

XX = (t(X) %*% X)
Xy = t(X) %*% y

# a working example of simply constrained solution
# i.e. constrain 1st and 4th values to be <1, 2nd and 3rd >0
(a =  diag(c(-1,1,1,-1)))
(b = c(1,0,0,1))

# this is correct according to these constraints, but I need 0 >= a <=
1 for all a and each row to sum to 1
solve.QP(Dmat = XX,dvec=Xy,Amat = a,bvec=b)


#my guess at the constraint matrix and b_vec...
(a2 = rbind(
   c(1,0,1,0), #specify the two cells in the 'row' I want to sum to one(???)
   c(0,1,0,1), # likewise
   kronecker(diag(rep(1,4)),c(1,-1)) #create positive and negative
constraint pairs
   )
)

(b2 = c(1,1,rep(0:1,times=4)))

solve.QP(Dmat = XX,dvec=Xy,Amat = a2,bvec=b2,meq=2) #complains that
Amat and dvec are incompatible




More information about the R-help mailing list