Минимизация наименьших квадратов в пределах порога в MATLAB

Question

Пакет cvx для MATLAB может решить (казалось бы, невинную) проблему оптимизации ниже, но он довольно медленный для больших полных матриц, с которыми я работаю. Я надеюсь, что это потому, что использование cvx является излишним, и что проблема на самом деле имеет аналитическое решение, или что умное использование некоторых встроенных функций MATLAB может быстрее выполнить эту работу.

Справочная информация: Хорошо известно, что оба x1=A\b а также x2=pinv(A)*b решить проблему наименьших квадратов:

minimize norm(A*x-b)

с той разницей, что norm(x2)<=norm(x1), По факту, x2 является минимально-нормальным решением проблемы, поэтому norm(x2)<=norm(x) для всех возможных решений x,

определяющий D=norm(A*x2-b), (эквивалентно D=norm(A*x1-b)), затем x2 решает проблему

minimize norm(x)
subject to
norm(A*x-b) == D

Проблема: я хотел бы найти решение для:

minimize norm(x)
subject to
norm(A*x-b) <= D+threshold

На словах мне не нужно norm(A*x-b) быть как можно меньше, просто в пределах определенной терпимости. Я хочу минимально-нормальное решение x что получает A*x в D+threshold из b,

Мне не удалось найти аналитическое решение проблемы (например, использовать псевдообратную в классической задаче наименьших квадратов) в Интернете или вручную. Я искал такие вещи, как "наименьшие квадраты с нелинейным ограничением" и "наименьшие квадраты с порогом".

Любая идея будет принята с благодарностью, но я полагаю, что мой реальный вопрос: каков самый быстрый способ решения этой "пороговой" задачи наименьших квадратов в MATLAB?

Answer 1

Интересный вопрос. Я не знаю ответа на ваш точный вопрос, но рабочее решение представлено ниже.

резюмировать

определять res(x) := norm(Ax - b), Как вы заявляете, x2 сводит к минимуму res(x), В переопределенном случае (как правило, A иметь больше рядов, чем col), x2 это уникальный минимум. В недоопределенном случае к нему присоединяется бесконечно много других *. Тем не менее, среди всех этих x2 это уникальный, который сводит к минимуму norm(x),

Подвести итоги, x2 минимизирует (1) res(x) и (2) norm(x)и это происходит в порядке приоритета. На самом деле это характеризует (полностью определяет) x2,

Предельная характеристика

Но, еще одна характеристика x2 является

x2 := limit_{e-->0} x_e

где

x_e := argmin_{x} J(x;e)

где

J(x;e) := res(x) + e * norm(x)

Можно показать, что

x_e = (A A' + e I)^{-1} A' b      (eqn a)

Следует понимать, что эта характеристика x2 довольно волшебно. Ограничения существуют, даже если (A A')^{-1} не. И предел как-то сохраняет приоритет (2) сверху.

Используя e>0

Конечно, для конечного (но маленького) e, x_e не будет минимизировать res(x)(вместо этого он минимизирует J(x;e)). В вашей терминологии разница - это порог. Я переименую его в

gap := res(x_e) - min_{x} res(x).

Уменьшение стоимости e гарантированно уменьшить значение gap, Достижение определенного gap поэтому значение (т. е. пороговое значение) легко достичь путем настройки e.**

Этот тип модификации (добавление norm(x) к res(x) проблема минимизации) известна как регуляризация в статистике и обычно считается хорошей идеей для стабильности (численно и в отношении значений параметров).

---

*: Обратите внимание, что x1 а также x2 отличаются только в недоопределенном случае

**: Это даже не требует каких-либо тяжелых вычислений, потому что обратное в (eqn a) легко вычисляется для любого (положительного) значения e если SVD для A уже был вычислен.