[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