Solveur de moindres carrés linéaires contraints

Votre problème est similaire à un problème problème non négatif problème (NNLS), qui peut être formulé comme

$$ \ min_x || AX-B || _2 ^ 2 \ Text {Sujet à} x \ GE 0 $$,

pour lequel il semble exister de nombreux algorithmes.

En réalité, votre problème peut être plus ou moins converti en un problème NNLS, si, en plus de vos variables non négatives d'origine $ x $ x, vous créez des variables supplémentaires $ x '$ et liez-les avec des contraintes linéaires $ x_i + x_i' = N $. Le problème avec cette approche est que ces contraintes linéaires supplémentaires pourraient ne pas être satisfaites exactement dans la solution des moindres carrés - il pourrait être approprié d'essayer de les pondérer avec un grand nombre.

Vous essayez de résoudre les moindres carrés avec des contraintes de boîte. Les algorithmes des moindres carrés standard incluent LSQR et plus récemment, LSMR. Celles-ci ne nécessitent que vous souhaitez appliquer des produits de matrice-vectoriels. Pour ajouter les contraintes, réalisez que si vous êtes à l'intérieur de la boîte (aucune des contraintes n'est «active»), vous procédez à la méthode de la page intérieure que vous avez choisie. Pour toutes les contraintes actives, la prochaine itération que vous effectuez vous permettra de désactiver la contrainte, ou de vous contraindre à vous déplacer le long de l'hyperplan de contrainte. Avec des modifications appropriées (conceptuellement relativement simples) adaptées à l'algorithme que vous choisissez, vous pouvez implémenter ces contraintes.

Généralement Cependant, vous pouvez utiliser n'importe quel package d'optimisation convexe. J'ai personnellement résolu ce type de problème exact à l'aide du package MATLAB CVX, qui utilise SDTPT3 / Sedumi pour un backend. CVX est simplement un wrapper très pratique autour de ces solveurs de programme semi-riginite.

Que diriez-vous de CVXOPT? Cela fonctionne avec des matrices clairsemées et il semble que certains des solveurs de programmation des cônes puissent aider:

http://abel.ee.uucle.edu/cvxopt/userguide/coneprog.html#quaDrac-cone-programs p>

C'est un simple Modification du code dans la DOC ci-dessus, pour résoudre votre problème: p>

from cvxopt import matrix, solvers
A = matrix([ [ .3, -.4,  -.2,  -.4,  1.3 ],
             [ .6, 1.2, -1.7,   .3,  -.3 ],
             [-.3,  .0,   .6, -1.2, -2.0 ] ])
b = matrix([ 1.5, .0, -1.2, -.7, .0])
N = 2;
m, n = A.size
I = matrix(0.0, (n,n))
I[::n+1] = 1.0
G = matrix([-I, I])
h = matrix(n*[0.0]  + n*[N])
print G
print h
dims = {'l': n, 'q': [n+1], 's': []}
x = solvers.coneqp(A.T*A, -A.T*b, G, h)['x']#, dims)['x']
print x

Si vous avez MATLAB, c'est quelque chose que vous pouvez faire avec TFOCS . Ceci est la syntaxe que vous utiliseriez pour résoudre le problème:

x = tfocs( smooth_quad, { A, -b }, proj_box( 0, N ) );

Si vous reformulez votre modèle comme:

min
Sous réserve de:

Ensuite, il s'agit d'un problème de programmation quadratique standard. Il s'agit d'un type de modèle commun qui peut être facilement résolu avec un solveur commercial tel que CPLEX ou GURROBI. (Disclaimer: Je travaille actuellement pour l'optimisation de Gurrobi et travaillé auparavant pour ILOG, qui a fourni CPLEX).