Редкий ограниченный линейный решатель наименьших квадратов
Это великое ТАК ответ указывает на хороший редкий решатель для Ax=b, но у меня есть ограничения на x таким образом, что каждый элемент в x >=0 <=N.
Кроме того, A огромно (вокруг 2e6x2e6), но очень редок с <=4 элементы на строку.
Какие-либо идеи/рекомендации? Я ищу что-то как MATLAB lsqlin но с огромными разреженными матрицами.
Я по существу пытаюсь решить ограниченную переменную проблему наименьших квадратов крупного масштаба на разреженных матрицах:
Править: В CVX:
cvx_begin
variable x(n)
minimize( norm(A*x-b) );
subject to
x <= N;
x >= 0;
cvx_end
4 ответа
Вы пытаетесь решить методом наименьших квадратов с ограничениями прямоугольника. Стандартные алгоритмы разреженных наименьших квадратов включают LSQR, а с недавних пор - LSMR. Для этого нужно только применить произведение матрица-вектор. Чтобы добавить ограничения, имейте в виду, что если вы находитесь внутри бокса (ни одно из ограничений не является «активным»), то вы продолжаете использовать любой метод внутренней точки, который вы выбрали. Для всех активных ограничений следующая итерация, которую вы выполняете, либо деактивирует ограничение, либо заставит вас двигаться по гиперплоскости ограничения. Эти ограничения можно реализовать с помощью некоторых (концептуально относительно простых) подходящих модификаций выбранного алгоритма.
Как правило, вы можете использовать любой пакет выпуклой оптимизации. Я лично решил именно этот тип проблемы, используя пакет Matlab CVX, который использует SDPT3 / SeDuMi в качестве бэкэнда. CVX - это просто очень удобная оболочка для этих полуопределенных программных решателей.
Ваша проблема похожа на неотрицательную задачу наименьших квадратов (NNLS), которую можно сформулировать как
$$ \ min_x || Ax-b || _2 ^ 2 \ text {subject to} x \ ge 0 $$,
для которого, кажется, существует много алгоритмов.
Фактически, ваша проблема может быть более или менее преобразована в проблему NNLS, если в дополнение к вашим исходным неотрицательным переменным $ x $ вы создадите дополнительные переменные $ x '$ и свяжете их линейными ограничениями $ x_i + x_i' = N $. Проблема с этим подходом состоит в том, что эти дополнительные линейные ограничения могут не быть выполнены точно в решении наименьших квадратов - тогда может быть уместно попытаться взвесить их большим числом.
Ваша матрица A ^ T A положительно полуопределена, поэтому ваша задача выпуклая; обязательно воспользуйтесь этим при настройке решателя.
Большинство популярных решателей QP написаны на Фортране и / или не бесплатны; однако я слышал хорошие отзывы о OOQP ( http://www.mcs.anl.gov/research/projects/otc/Tools/OOQP/OoqpRequestForm.html ), хотя это немного больно получение копии.
Как насчет CVXOPT? Он работает с разреженными матрицами, и кажется, что некоторые решатели конусного программирования могут помочь:
http://abel.ee.ucla.edu/cvxopt/userguide/coneprog.html#quadratic-cone-programs
Это простая модификация кода в документе выше, чтобы решить вашу проблему:
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
CVXOPT поддерживает разреженные матрицы, так что он может быть полезен для вас.