Наименьшие квадраты (подбор кривой модели) алгоритмы

Определение наименьших квадратов

Наименьшие квадраты, в целом, являются задачей нахождения векторный x, который является локальным минимизатором к функции, которая является суммой квадратов, возможно подвергните некоторым ограничениям:

$\min_{x} {‖ F (x) ‖}_{2}^{2} = \min_{x} \sum_{i} F_{i}^{2} (x)$

таким образом, что A·x ≤ b, Aeq·x = beq, lb ≤ x ≤ ub.

Существует несколько решателей Optimization Toolbox™, доступных для различных типов F (x) и различные типы ограничений:

Решатель	F (x)	Ограничения
`mldivide`	C· – d	'none'
`lsqnonneg`	C· – d	x ≥ 0
`lsqlin`	C· – d	Связанный, линейный
`lsqnonlin`	Общий F (x)	Связанный
`lsqcurvefit`	F (x, xdata) – ydata	Связанный

Существуют алгоритмы на пять наименьших квадратов в решателях Optimization Toolbox, в дополнение к алгоритмам, используемым в mldivide:

lsqlin внутренняя точка
lsqlin активный набор
"Доверительная область, отражающая" (нелинейный или линейный метод наименьших квадратов)
Levenberg-Marquardt (нелинейный метод наименьших квадратов)
Алгоритм используется lsqnonneg

Все алгоритмы кроме lsqlin активный набор является крупномасштабным; смотрите Крупномасштабный по сравнению с Алгоритмами Средней шкалы. Для общего обзора методов нелинейного метода наименьших квадратов смотрите Денниса [8]. Определенные детали о методе Levenberg-Marquardt могут быть найдены в Moré [28].

Линейный метод наименьших квадратов: внутренняя точка или Активный Набор

lsqlin 'interior-point' алгоритм использует выпуклый внутренней точкой quadprog Алгоритм и lsqlin 'active-set' алгоритм использует активный набор quadprog алгоритм. quadprog описание задачи должно минимизировать квадратичную функцию

$\min_{x} \frac{1}{2} x^{T} H x + c^{T} x$

подвергните линейным ограничениям и связанным ограничениям. lsqlin функция минимизирует квадратичную норму векторного Cx – d, удовлетворяющего линейным ограничениям и связанным ограничениям. Другими словами, lsqlin минимизирует

$\begin{matrix} {‖ C x - d ‖}_{2}^{2} = {(C x - d)}^{T} (C x - d) \\ = (x^{T} C^{T} - d^{T}) (C x - d) \\ = (x^{T} C^{T} C x) - (x^{T} C^{T} d - d^{T} C x) + d^{T} d \\ = \frac{1}{2} x^{T} (2 C^{T} C) x + {(- 2 C^{T} d)}^{T} x + d^{T} d . \end{matrix}$

Это помещается в quadprog среда путем установки матрицы H на 2^CTC и вектора c к ^(–2CTd). (Вектор свободных членов ^dTd не оказывает влияния на местоположение минимума.) После этой переформулировки lsqlin проблема, quadprog вычисляет решение.

Примечание

quadprog 'interior-point-convex' алгоритм имеет два пути выполнения кода. Это берет тот, когда матрица Гессиана, из которой H является обычной (полной) матрицей, удваивается, и это берет другой, когда H является разреженной матрицей. Для получения дополнительной информации разреженного типа данных, смотрите Разреженные матрицы (MATLAB). Обычно алгоритм быстрее для больших проблем, которые имеют относительно немного ненулевых условий, когда вы задаете H как sparse. Точно так же алгоритм быстрее для небольших или относительно плотных проблем, когда вы задаете H как full.

Доверительная область отражающие наименьшие квадраты

Доверительная область отражающий алгоритм наименьших квадратов

Многие методы, используемые в решателях Optimization Toolbox, основаны на доверительных областях, простой все же мощной концепции в оптимизации.

Чтобы изучить подход доверительной области к оптимизации, рассмотрите задачу безусловной минимизации, минимизируйте f (x), где функция берет аргументы вектора и возвращает скаляры. Предположим, что вы - в точке x в n - пробел, и вы хотите улучшиться, т.е. переместиться в точку с более низким значением функции. Основная идея состоит в том, чтобы аппроксимировать f более простым функциональным q, который обоснованно отражает поведение функционального f в окружении N вокруг точки x. Это окружение является доверительной областью. Испытательный шаг s вычисляется путем минимизации (или приблизительно минимизации) по N. Это - подпроблема доверительной области,

\min_{s} {q (s), s \in N} .

(1)

Текущая точка обновляется, чтобы быть x + s если f (x + s) < f (x); в противном случае текущая точка остается неизменной, и N, область доверия, уменьшается, и испытательный расчет шага повторяется.

Ключевые вопросы в определении определенного подхода доверительной области к минимизации f (x) состоят в том, как выбрать и вычислить приближение q (заданный в текущей точке x), как выбрать и изменить доверительную область N, и как точно решить подпроблему доверительной области. Этот раздел фокусируется на неограниченной проблеме. Более поздние разделы обсуждают дополнительные осложнения из-за присутствия ограничений на переменные.

В стандартном методе доверительной области ([48]), квадратичное приближение q задан первыми двумя сроками приближения Тейлора к F в x; окружение N является обычно сферическим или эллипсоидальным в форме. Математически подпроблема доверительной области обычно утверждается

\min {\frac{1}{2} s^{T} H s + s^{T} g таким образом , что ‖ D s ‖ \leq Δ},

(2)

где g является градиентом f в текущей точке x, H является матрицей Гессиана (симметрическая матрица вторых производных), D является диагональной матрицей масштабирования, Δ является положительной скалярной величиной и ∥. ∥ является 2-нормой. Хорошие алгоритмы существуют для решения уравнения 2 (см. [48]); такие алгоритмы обычно включают расчет всех собственных значений H, и процесс Ньютона применился к характеристическому уравнению

$\frac{1}{Δ} - \frac{1}{‖ s ‖} = 0.$

Такие алгоритмы предоставляют точное решение уравнения 2. Однако они требуют времени, пропорционального нескольким факторизациям H. Поэтому для проблем доверительной области другой подход необходим. Несколько приближений и эвристических стратегий, на основе уравнения 2, были предложены в литературе ([42] и [50]). Подход приближения, сопровождаемый в решателях Optimization Toolbox, должен ограничить подпроблему доверительной области двумерным подпространством S ([39] и [42]). Однажды подпространство S был вычислен, работа, чтобы решить уравнение 2 тривиальна, даже если полная информация о собственном значении/собственном векторе необходима (поскольку в подпространстве, проблема только двумерна). Доминирующая работа теперь переключила к определению подпространства.

Двумерное подпространство S определяется при помощи предобусловленного процесса метода сопряженных градиентов, описанного ниже. Решатель задает S как линейный пробел, заполненный s ₁ и s ₂, где s ₁ является в направлении градиента g, и s ₂ является любой аппроксимированным направлением Ньютона, т.е. решением

H \cdot s_{2} = - g,

(3)

или направление отрицательной кривизны,

s_{2}^{T} \cdot H \cdot s_{2} < 0.

(4)

Философия позади этого выбора S должна обеспечить глобальную сходимость (через направление наискорейшего спуска или направление отрицательной кривизны) и достигнуть быстро локальной сходимости (через шаг Ньютона, когда это существует).

Эскиз безусловной минимизации с помощью идей доверительной области теперь легко дать:

Сформулируйте двумерную подпроблему доверительной области.
Решите уравнение 2, чтобы определить испытательный шаг s.
Если f (x + s) <f (x), то x = x + s.
Настройте Δ.

Эти четыре шага повторяются до сходимости. Размерность доверительной области Δ настроена согласно стандартным правилам. В частности, это уменьшено, если испытательный шаг не принят, т.е. f (x + s) ≥ f (x). См. [46] и [49] для обсуждения этого аспекта.

Решатели Optimization Toolbox обрабатывают несколько важных особых случаев f со специализированными функциями: нелинейный метод наименьших квадратов, квадратичные функции и линейный метод наименьших квадратов. Однако базовые алгоритмические идеи эквивалентны для общего случая. Эти особые случаи обсуждены в более поздних разделах.

Крупномасштабный нелинейный метод наименьших квадратов

Важный особый случай для f (x) является проблемой нелинейного метода наименьших квадратов

\min_{x} \sum_{i} f_{i}^{2} (x) = \min_{x} {‖ F (x) ‖}_{2}^{2},

(5)

где F (x) является функцией с векторным знаком с i компонента F (x), равный _fi (x). Основной метод, используемый, чтобы решить эту задачу, эквивалентен в общем случае, описанном в Методах Доверительной области для Нелинейной Минимизации. Однако структура проблемы нелинейного метода наименьших квадратов используется, чтобы улучшить КПД. В частности, аппроксимированное направление Ньютона Гаусса, т.е. решение s к

\min {‖ J s + F ‖}_{2}^{2},

(6)

(где J является якобианом F (x)), используется, чтобы помочь задать двумерное подпространство S. Вторые производные функционального _fi компонента (x) не используются.

В каждой итерации метод предобусловленных методов сопряженных градиентов используется, чтобы приблизительно решить нормальные уравнения, т.е.

$J^{T} J s = - J^{T} F,$

несмотря на то, что нормальные уравнения явным образом не формируются.

Крупномасштабный линейный метод наименьших квадратов

В этом случае функциональный f (x), который будет решен,

$f (x) = {‖ C x + d ‖}_{2}^{2},$

возможно подвергните линейным ограничениям. Алгоритм генерирует строго выполнимый, выполняет итерации схождения, в пределе, к локальному решению. Каждая итерация включает приближенное решение большой линейной системы (порядка n, где n является длиной x). Матрицы итерации имеют структуру матричного C. В частности, метод предобусловленных методов сопряженных градиентов используется, чтобы приблизительно решить нормальные уравнения, т.е.

$C^{T} C x = - C^{T} d,$

несмотря на то, что нормальные уравнения явным образом не формируются.

Метод доверительной области подпространства используется, чтобы определить поисковое направление. Однако вместо того, чтобы ограничить шаг (возможно) одним отражательным шагом, как в нелинейном случае минимизации, кусочный отражающий поиск линии проводится в каждой итерации, как в квадратичном случае. См. [45] для деталей поиска линии. В конечном счете линейные системы представляют подход Ньютона, получая условия оптимальности первого порядка в решении, приводя к сильным локальным уровням сходимости.

Функция умножения якобиана. lsqlin может решить задачу линейно-метода-наименьших-квадратов-с-ограничениями, не используя матричный C явным образом. Вместо этого это использует функцию умножения якобиана jmfun,

W = jmfun(Jinfo,Y,flag)

то, что вы обеспечиваете. Функция должна вычислить следующие продукты для матричного Y:

Если flag == 0 затем W = C'*(C*Y).
Если flag > 0 затем W = C*Y.
Если flag < 0 затем W = C'*Y.

Это может быть полезно, если C является большим, но содержит достаточно структуры, что можно записать jmfun не формируя C явным образом. Для примера смотрите функцию умножения якобиана с Линейным методом наименьших квадратов.

Метод Levenberg-Marquardt

В задаче наименьших квадратов функциональный f (x) минимизирован, который сумма квадратов.

\min_{x} f (x) = {‖ F (x) ‖}_{2}^{2} = \sum_{i} F_{i}^{2} (x) .

(7)

Проблемы этого типа происходят в большом количестве практических применений, особенно когда подходящие функции модели к данным, т.е. нелинейная оценка параметра. Они также распространены в управлении, где вы хотите выход, y (x,t), чтобы следовать за некоторой непрерывной траекторией модели, φ (t), для векторного x и скалярного t. Эта проблема может быть выражена как

\min_{x \in ℜ^{n}} \int_{t_{1}}^{t_{2}} {(y (x, t) - φ (t))}^{2} d t,

(8)

где y (x,t) и φ (t) является скалярными функциями.

Когда интеграл дискретизируется с помощью подходящей формулы квадратуры, вышеупомянутое может быть сформулировано как задача наименьших квадратов:

\min_{x \in ℜ^{n}} f (x) = \sum_{i = 1}^{m} {(\bar{y} (x, t_{i}) - \bar{φ} (t_{i}))}^{2},

(9)

где $\bar{y}$ и $\bar{φ}$ включайте веса квадратурной схемы. Обратите внимание на то, что в этой проблеме векторный F (x)

$F (x) = [\begin{matrix} \bar{y} (x, t_{1}) - \bar{φ} (t_{1}) \\ \bar{y} (x, t_{2}) - \bar{φ} (t_{2}) \\ ... \\ \bar{y} (x, t_{m}) - \bar{φ} (t_{m}) \end{matrix}] .$

В проблемах этого вида невязка ∥F (x) ∥, вероятно, будет мала в оптимуме, поскольку это - общая практика, чтобы установить реалистично достижимые целевые траектории. Несмотря на то, что функция в LS может быть минимизирована с помощью общего метода безусловной минимизации, как описано в Основах Оптимизации без ограничений, определенные характеристики проблемы могут часто использоваться, чтобы повысить итеративную эффективность процедуры решения. Градиент и матрица Гессиана LS имеют специальную структуру.

Обозначая m-by-n якобиевская матрица F (x) как J (x), вектор градиента f (x) как G (x), матрица Гессиана f (x) как H (x) и матрица Гессиана каждого _Fi (x) как _Hi (x), вы имеете

\begin{array}{l} G (x) = 2 J (^{x) T} F (x) \\ H (x) = 2 J (^{x) T} J (x) + 2 Q (x), \end{array}

(10)

где

$Q (x) = \sum_{i = 1}^{m} F_{i} (x) \cdot H_{i} (x) .$

Матричный Q (x) имеет свойство, что, когда невязка ∥F (x) ∥ имеет тенденцию обнулять, когда _xk приближается к решению, затем Q (x) также имеет тенденцию обнулять. Таким образом, когда ∥F (x) ∥ мал в решении, очень эффективный метод состоит в том, чтобы использовать направление Ньютона Гаусса в качестве основания для процедуры оптимизации.

В методе Ньютона Гаусса поисковое направление, _dk, получено в каждой главной итерации, k, который является решением проблемы линейного метода наименьших квадратов:

\min_{x \in ℜ^{n}} {‖ J (x_{k}) - F (x_{k}) ‖}_{2}^{2} .

(11)

Направление, выведенное из этого метода, эквивалентно направлению Ньютона, когда условия Q (x) могут быть проигнорированы. Поисковое направление _dk может использоваться в качестве части стратегии поиска линии гарантировать в каждой итерации функциональный f (x) уменьшения.

Метод Ньютона Гаусса часто сталкивается с проблемами, когда термин второго порядка Q (x) является значительным. Метод, который преодолевает эту проблему, является методом Levenberg-Marquardt.

Levenberg-Marquardt [25], и [27] метод использует поисковое направление, которое является решением линейной системы уравнений

(J {(x_{k})}^{T} J (x_{k}) + λ_{k} I) d_{k} = - J {(x_{k})}^{T} F (x_{k}),

(12)

или, опционально, уравнений

(J {(x_{k})}^{T} J (x_{k}) + λ_{k} d i a g (J {(x_{k})}^{T} J (x_{k}))) d_{k} = - J {(x_{k})}^{T} F (x_{k}),

(13)

где скалярный _λk управляет и величиной и направлением _dk. Установите опцию ScaleProblem к 'none' выбрать Equation 12 и установить ScaleProblem к 'Jacobian' выбрать Equation 13.

Вы устанавливаете начальное значение параметра λ ₀ использований InitDamping опция. Иногда, 0.01 значение по умолчанию этой опции может быть неподходящим. Если вы находите, что алгоритм Levenberg-Marquardt делает мало начальных успехов, попробуйте установку InitDamping к различному значению, чем значение по умолчанию, возможно, 1e2.

Когда _λk является нулем, направление_{, dk} идентичен тому из метода Ньютона Гаусса. Когда _λk стремится к бесконечности, _dk стремится к направлению наискорейшего спуска с величиной, имеющей тенденцию обнулять. Это подразумевает, что для некоторого достаточно большого _λk, термин F (_xk + _dk) < F (_xk) сохраняется. Термином _λk можно поэтому управлять, чтобы гарантировать спуск, даже когда с условиями второго порядка, которые ограничивают КПД метода Ньютона Гаусса, сталкиваются. То, когда шаг успешен (дает более низкое значение функции), алгоритм устанавливает λ _{k +1} = _λk/10. Когда шаг неудачен, алгоритм устанавливает λ _{k +1} = _λk *10.

Внутренне, алгоритм Levenberg-Marquardt использует допуск оптимальности (останавливающий критерий) 1e-4 времена функциональный допуск.

Метод Levenberg-Marquardt поэтому использует поисковое направление, которое является пересечением направления Ньютона Гаусса и направления наискорейшего спуска. Это проиллюстрировано в рисунке 12-1, Методе Levenberg-Marquardt на Функции Розенброка. Решение для функции Розенброка сходится после 90 вычислений функции по сравнению с 48 для метода Ньютона Гаусса. Более плохой КПД частично, потому что метод Ньютона Гаусса является обычно более эффективным, когда невязка является нулем в решении. Однако такая информация не всегда доступна заранее, и увеличенная робастность метода Levenberg-Marquardt компенсирует свой случайный более плохой КПД.

Рисунок 12-1, метод Levenberg-Marquardt на функции Розенброка

Для большего количества полного описания этого рисунка, включая скрипты, которые генерируют итеративные точки, смотрите Банановую Минимизацию Функции.

Документация