Optimisation convexe CVX-esque en R ?

Question

J'ai besoin de résoudre (plusieurs fois, pour de nombreuses données, en plus d'un tas d'autres choses) ce qui, selon moi, se résume à une Programme de cônes du second ordre . Elle peut être exprimée succinctement en CVX quelque chose comme ça :

cvx_begin
    variable X(2000);
    expression MX(2000);
    MX = M * X;
    minimize( norm(A * X - b) + gamma * norm(MX, 1) )
  subject to
    X >= 0
    MX((1:500) * 4 - 3) == MX((1:500) * 4 - 2)
    MX((1:500) * 4 - 1) == MX((1:500) * 4)
cvx_end

Les longueurs de données et les modèles de contraintes d'égalité présentés ne sont que des valeurs arbitraires tirées de certaines données de test, mais la forme générale sera sensiblement la même, avec deux termes objectifs - l'un minimisant l'erreur, l'autre encourageant la sparsité - et un grand nombre de contraintes d'égalité sur les éléments d'une version transformée de la variable d'optimisation (elle-même contrainte à être non négative).

Cela semble fonctionner assez bien, bien mieux que mon approche précédente, qui truque les contraintes de manière tout à fait pourrie. Le problème, c'est que tout le reste se passe en R, et qu'il serait assez ennuyeux de devoir le porter à Matlab. Alors, est-il possible de faire cela en R, et si oui, comment ?

Cela se résume en fait à deux questions distinctes :

1) Existe-t-il de bonnes ressources R pour cela ? D'après ce que je peux dire de la Page des tâches du CRAN les options du paquet SOCP sont CLSCOP y DWD qui comprend un solveur SOCP en complément de son classificateur. Les deux ont des interfaces similaires mais assez opaques et sont un peu pauvres en documentation et en exemples, ce qui nous amène à.. :

2) Quelle est la meilleure façon de représenter le problème ci-dessus dans le format de bloc de contraintes utilisé par ces paquets ? La syntaxe CVX ci-dessus cache beaucoup de manipulations fastidieuses avec des variables supplémentaires et d'autres éléments, et je me vois bien passer des heures à faire des calculs. semaines J'essaie de bien faire les choses, donc tous les conseils ou les pointeurs pour me mettre dans la bonne direction seraient les bienvenus...

Answer 1

Vous pouvez trouver le paquet R CVXdeR utile. Cela vous permet de passer un problème d'optimisation à CVX depuis R et de retourner la solution à R.

Answer 2

OK, donc la réponse courte à cette question est : il n'y a pas vraiment de manière satisfaisante de gérer cela dans R. J'ai fini par faire les parties pertinentes dans Matlab avec quelques modifications maladroites entre les deux systèmes, et je vais probablement tout migrer vers Matlab par la suite. (Mon approche actuelle est antérieure à la réponse postée par l'utilisateur 2439686. En pratique, mon problème serait tout aussi gênant en utilisant CVXfromR, mais cela semble être un paquetage utile en général, donc je vais accepter cette réponse).

Les ressources R pour cela sont plutôt minces, mais la article de blog de Vincent Zoonekynd qu'il a mentionné dans les commentaires vaut vraiment la peine d'être lu.

Le solveur SOCP contenu dans le package R DWD est porté depuis le solveur Matlab. SDPT3 (sans les parties SDP), donc l'interface programmatique est fondamentalement la même. Cependant, au moins dans mes tests, il fonctionne beaucoup plus lentement et s'effondre sur les problèmes avec quelques milliers de variables + contraintes, alors que SDPT3 les résout en quelques secondes. (Je n'ai pas fait une comparaison tout à fait juste sur ce point, parce que CVX fait quelques transformations astucieuses sur le problème pour le rendre plus efficace, alors que dans R j'utilise une définition assez naïve, mais quand même).

Une autre alternative possible, surtout si vous pouvez bénéficier d'une licence académique, est d'utiliser la licence commerciale Mosek qui possède un paquet d'interface R Rmosek . Je n'ai pas encore essayé, mais je le ferai peut-être à un moment donné.

(Pour l'anecdote, l'autre solveur fourni avec CVX, SeDuMi, échoue complètement sur le même problème ; les auteurs de CVX ne plaisantent pas lorsqu'ils suggèrent d'essayer plusieurs solveurs. De plus, dans un sous-ensemble significatif de cas, SDTP3 doit passer de la décomposition de Cholesky à la décomposition LU, ce qui rend le traitement plus lent de plusieurs ordres de grandeur, avec une amélioration très marginale de l'objectif par rapport aux étapes pré-LU. J'ai trouvé qu'il valait la peine de réduire la précision demandée pour éviter cela, mais YMMV).

Answer 3

Il existe une nouvelle alternative : CVXR qui vient des mêmes personnes. Il y a un site web , a papier et un projet github .

L'observation de la programmation convexe disciplinée semble gagner en popularité. cvxpy (Python) et Convex.jl (Julia), encore une fois, soutenu par les mêmes personnes.