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

Я использовал итеративно переоцененные наименьшие квадраты (IRLS), чтобы минимизировать функции следующей формы,

$J(m) = \sum_{i=1}^{N} \rho \left(\left| x_i - m \right|\right)$

где $N$ - количество экземпляров $x_i \in \mathbb{R}$ , $m \in \mathbb{R}$ - надежная оценка, которую я хочу, а $\rho$ - подходящая робастная штрафная функция. Допустим, он выпуклый (хотя не обязательно строго) и на данный момент дифференцируемый. Хорошим примером такого $\rho$ является функция потерь Хьюбера .

То, что я делал, это дифференцирование $J(m)$ отношению к $m$ (и манипулирование) для получения,

$\frac{dJ}{dm}= \sum_{i=1}^{N} \frac{\rho'\left( \left|x_i-m\right|\right) }{\left|x_i-m\right|} \left( x_i-m \right)$

и итеративно решения, установив его равным 0 и фиксации веса при итерации $k$ к $w_i(k) = \frac{\rho'\left( \left|x_i-m{(k)}\right|\right) }{\left|x_i-m{(k)}\right|}$(обратите внимание, что воспринимаемая особенность при$x_i=m{(k)}$действительно является устраняемой сингулярностью во всех$\rho$о которых я мог бы беспокоиться). Тогда я получаю,

$\sum_{i=1}^{N} w_i(k) \left( x_i-m{(k+1)} \right)=0$

и я решаю получить, .$m(k+1) = \frac{\sum_{i=1}^{N} w_i(k) x_i}{ \sum_{i=1}^{N} w_i(k)}$

Я повторяю этот алгоритм с фиксированной точкой до "сходимости". Отмечу, что если вы доберетесь до фиксированной точки, вы оптимальны, поскольку ваша производная равна 0 и это выпуклая функция.

У меня есть два вопроса об этой процедуре:

1. Это стандартный алгоритм IRLS? После прочтения нескольких статей по этой теме (и они были очень разбросаны и расплывчаты в отношении того, что такое IRLS), это наиболее последовательное определение алгоритма, которое я могу найти. Я могу публиковать газеты, если люди хотят, но я на самом деле не хотел никого пристрастить. Конечно, вы можете обобщить эту базовую технику на многие другие типы проблем, связанных с вектором и аргументами, отличными от | x i - m ( k ) | Предоставление аргумента является нормой аффинной функции ваших параметров. Любая помощь или понимание было бы здорово в этом.$x_i$$\left|x_i-m{(k)}\right|$
2. Конвергенция, кажется, работает на практике, но у меня есть несколько опасений по этому поводу. Я еще не видел доказательства этого. После нескольких простых симуляций Matlab я вижу, что одна итерация этого не является отображением сжатия (я сгенерировал два случайных экземпляра и вычислил | m 1 ( k + 1 ) - m 2 ( k + 1 ) $m$и увидел, что это иногда больше, чем 1). Кроме того, отображение, определенное несколькими последовательными итерациями, не является строго сжатым отображением, но вероятность того, что константа Липшица будет больше 1, становится очень низкой. Так есть ли понятиевероятностного картографирования? Какой механизм я бы использовал, чтобы доказать, что это сходится? Это даже сходится?$\frac{\left|m_1(k+1) - m_2(k+1)\right|}{\left|m_1(k)-m_2(k)\right|}$

Любое руководство вообще полезно.

Редактировать: Мне нравится статья о IRLS для разреженного восстановления / обнаружения сжатия, написанная Daubechies et al. 2008 «Итеративно повторно взвешенный метод минимизации наименьших квадратов для разреженного восстановления» на arXiv. Но, похоже, основное внимание уделяется весам невыпуклых задач. Мой случай значительно проще.

Что касается вашего первого вопроса, нужно определить «стандарт» или признать, что «каноническая модель» была постепенно создана. Как указано в комментарии, по крайней мере, кажется, что вы используете IRWLS довольно стандартно.

Что касается вашего второго вопроса, «сжатие по вероятности» может быть связано (хотя и неформально) с конвергенцией «рекурсивных стохастических алгоритмов». Из того, что я прочитал, есть огромная литература по этому вопросу, в основном в области машиностроения. В экономике мы используем его чуть-чуть, особенно основополагающие работы Леннарта Юнга - первой работой была Юнг (1977), в которой показано, что сходимость (или нет) рекурсивного стохастического алгоритма может определяться стабильностью (или нет) связанного обыкновенного дифференциального уравнения.

(то, что следует, было переработано после плодотворного обсуждения с ФП в комментариях)

конвергенция

Я буду использовать в качестве ссылки Сабер Элайди "Введение в разностные уравнения", 2005, 3-е изд. Анализ условно на некоторой выборке данных, так что рассматриваются как фиксированные. $x's$

Условие первого порядка для минимизации целевой функции, рассматривается как рекурсивная функция в , м ( к + 1 ) = N Σ я = 1 V я [ м ( K ) ] х я$m$

$$m(k+1) = \sum_{i=1}^{N} v_i[m(k)] x_i, \;\; v_i[m(k)] \equiv \frac{w_i[m(k)]}{ \sum_{i=1}^{N} w_i[m(k)]} \qquad [1]$$

имеет фиксированную точку (argmin целевой функции). По теореме 1.13 стр. 27-28 из Elaydi, если первая производная по от RHS из [ 1 ] , оцененная в неподвижной точке m ∗ , обозначим ее A ′ ( m ∗ ) , меньше единицы по абсолютной величине , то м * является асимптотически устойчивым (AS). Более того, по теореме 4.3 с.179 мы имеем, что это также означает, что неподвижная точка равномерно AS (UAS). «Асимптотически устойчивый» означает, что для некоторого диапазона значений вокруг неподвижной точки существует окрестность ( m $m$$[1]$$m^*$$A'(m^*)$$m^*$
, не обязательно небольшого размера, неподвижная точкапривлекательна, и поэтому, если алгоритм дает значения в этой окрестности, он будет сходиться. Свойство «равномерный» означает, что граница этой окрестности и, следовательно, ее размер не зависят от начального значения алгоритма. Фиксированная точка становитсяглобальноUAS, если γ = ∞ . Так что в нашем случае, если мы докажем, что$(m^* \pm \gamma)$$\gamma = \infty$

$$|A'(m^*)|\equiv \left|\sum_{i=1}^{N} \frac{\partial v_i(m^*)}{\partial m}x_i\right| <1 \qquad [2]$$

$x$$\bar x$

$$\frac{\partial v_i(m^*)}{\partial m} = \frac {\frac{\partial w_i(m^*)}{\partial m}\sum_{i=1}^{N} w_i(m^*)-w_i(m^*)\sum_{i=1}^{N}\frac{\partial w_i(m^*)}{\partial m}}{\left(\sum_{i=1}^{N} w_i(m^*)\right)^2}$$

$$=\frac 1{\sum_{i=1}^{N} w_i(m^*)}\cdot\left[\frac{\partial w_i(m^*)}{\partial m}-v_i(m^*)\sum_{i=1}^{N}\frac{\partial w_i(m^*)}{\partial m}\right]$$ Then

$$A'(m^*) = \frac 1{\sum_{i=1}^{N} w_i(m^*)}\cdot\left[\sum_{i=1}^{N}\frac{\partial w_i(m^*)}{\partial m}x_i-\left(\sum_{i=1}^{N}\frac{\partial w_i(m^*)}{\partial m}\right)\sum_{i=1}^{N}v_i(m^*)x_i\right]$$

$$=\frac 1{\sum_{i=1}^{N} w_i(m^*)}\cdot\left[\sum_{i=1}^{N}\frac{\partial w_i(m^*)}{\partial m}x_i-\left(\sum_{i=1}^{N}\frac{\partial w_i(m^*)}{\partial m}\right)m^*\right]$$

and

$$|A'(m^*)| <1 \Rightarrow \left|\sum_{i=1}^{N}\frac{\partial w_i(m^*)}{\partial m}(x_i-m^*)\right| < \left|\sum_{i=1}^{N} w_i(m^*)\right| \qquad [3]$$

we have

$$\begin{align}\frac{\partial w_i(m^*)}{\partial m} = &\frac{-\rho''(|x_i-m^*|)\cdot \frac {x_i-m^*}{|x_i-m^*|}|x_i-m^*|+\frac {x_i-m^*}{|x_i-m^*|}\rho'(|x_i-m^*|)}{|x_i-m^*|^2} \\ \\ &=\frac {x_i-m^*}{|x_i-m^*|^3}\rho'(|x_i-m^*|) - \rho''(|x_i-m^*|)\cdot \frac {x_i-m^*}{|x_i-m^*|^2} \\ \\ &=\frac {x_i-m^*}{|x_i-m^*|^2}\cdot \left[\frac {\rho'(|x_i-m^*|)}{|x_i-m^*|}-\rho''(|x_i-m^*|)\right]\\ \\ &=\frac {x_i-m^*}{|x_i-m^*|^2}\cdot \left[w_i(m^*)-\rho''(|x_i-m^*|)\right] \end{align}$$

Inserting this into $[3]$ we have

$$\left|\sum_{i=1}^{N}\frac {x_i-m^*}{|x_i-m^*|^2}\cdot \left[w_i(m^*)-\rho''(|x_i-m^*|)\right](x_i-m^*)\right| < \left|\sum_{i=1}^{N} w_i(m^*)\right|$$

$$\Rightarrow \left|\sum_{i=1}^{N}w_i(m^*)-\sum_{i=1}^{N}\rho''(|x_i-m^*|)\right| < \left|\sum_{i=1}^{N} w_i(m^*)\right| \qquad [4]$$

This is the condition that must be satisfied for the fixed point to be UAS. Since in our case the penalty function is convex, the sums involved are positive. So condition $[4]$ is equivalent to

$$\sum_{i=1}^{N}\rho''(|x_i-m^*|) < 2\sum_{i=1}^{N}w_i(m^*) \qquad [5]$$

If $\rho(|x_i-m|)$ is Hubert's loss function, then we have a quadratic ($q$) and a linear ($l$) branch,

$$\rho(|x_i-m|)=\cases{ (1/2)|x_i- m|^2 \qquad\;\;\;\; |x_i-m|\leq \delta \\ \\ \delta\big(|x_i-m|-\delta/2\big) \qquad |x_i-m|> \delta}$$

and

$$\rho'(|x_i-m|)=\cases{ |x_i- m| \qquad |x_i-m|\leq \delta \\ \\ \delta \qquad \qquad \;\;\;\; |x_i-m|> \delta}$$

$$\rho''(|x_i-m|)=\cases{ 1\qquad |x_i-m|\leq \delta \\ \\ 0 \qquad |x_i-m|> \delta}$$

$$\cases{ w_{i,q}(m) =1\qquad \qquad \qquad |x_i-m|\leq \delta \\ \\ w_{i,l}(m) =\frac {\delta}{|x_i-m|} <1 \qquad |x_i-m|> \delta}$$

Since we do not know how many of the $|x_i-m^*|$'s place us in the quadratic branch and how many in the linear, we decompose condition $[5]$ as ($N_q + N_l = N$)

$$\sum_{i=1}^{N_q}\rho_q''+\sum_{i=1}^{N_l}\rho_l'' < 2\left[\sum_{i=1}^{N_q}w_{i,q} +\sum_{i=1}^{N_l}w_{i,l}\right]$$

$$\Rightarrow N_q + 0 < 2\left[N_q +\sum_{i=1}^{N_l}w_{i,l}\right] \Rightarrow 0 < N_q+2\sum_{i=1}^{N_l}w_{i,l}$$

which holds. So for the Huber loss function the fixed point of the algorithm is uniformly asymptotically stable, irrespective of the $x$'s. We note that the first derivative is smaller than unity in absolute value for any $m$, not just the fixed point.

What we should do now is either prove that the UAS property is also global, or that, if $m(0) = \bar x$ then $m(0)$ belongs to the neighborhood of attraction of $m^*$.