Оценка параметров нормального распределения

[stiv1995]

Добрый день, возник следующий вопрос:

пусть дана реализация случайной матрицы $M \in \mathbb{R}^{n \times m}$, компоненты которой независимы и имеют распределение

$$M_{ij} \sim \mathcal{N}(\mu_i, \sigma_j)$$

Иными словами, случайные величины, принадлежащие одной строке матрицы, имеют одинаковое математическое ожидание, а случайные величины, принадлежащие одному столбцу, имеют одинаковую дисперсию.

Вопрос в следующем: как оптимально оценить параметры $\{\mu_i\}, \{\sigma_i\}$ вместе или хотя бы матожидания? В качестве метрики рассматривается $L_2$.

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

[alisa-lebovski]

Точный результат здесь вряд ли можно получить, но можно действовать методом последовательных приближений. То есть начать, как вы сказали. Потом каждый раз вычислять оценки среднего исходя из предыдущих оценок дисперсий, а оценки дисперсий - исходя из оценок средних, пока к чему-нибудь не сойдется. Полагаю, этот метод эквивалентен максимизации правдоподобия опять-таки методом последовательных приближений.

[B@R5uk]

А разве не очевидно, что тут надо минимизировать величину

$$\[\sum\limits_{i,j}{{{\left( \frac{{{m}_{ij}}-{{\mu }_{i}}}{{{\sigma }_{j}}} \right)}^{2}}}\]$$

варьируя ${\mu }_{i}$ при фиксированных ${\sigma }_{j}$, а затем, получив явное выражение первого через второе, добиться чтобы величины под квадратом в сумме были наиболее похожи на $\mathcal{N}(0, 1)$ (так же очевидно, что тут надо варьировать ${\sigma }_{j}$)?

-- 02.07.2018, 14:27 --

stiv1995 в сообщении #1323926 писал(а):

Метод максимального правдоподобия приводит к системе, которая вряд ли имеет решение.

А можно взглянуть на уравнение, пожалуйста, а то я эту тему подзабыл уже?

[stiv1995]

alisa-lebovski в сообщении #1323934 писал(а):

Точный результат здесь вряд ли можно получить, но можно действовать методом последовательных приближений

Да, это как раз то, что я имел ввиду под "EM алгоритм".

-- 02.07.2018, 16:14 --

B@R5uk в сообщении #1323950 писал(а):

А разве не очевидно, что тут надо минимизировать величину

$$\[\sum\limits_{i,j}{{{\left( \frac{{{m}_{ij}}-{{\mu }_{i}}}{{{\sigma }_{j}}} \right)}^{2}}}\]$$

Метод максимизации правдоподобия ведет к оптимизационной задаче

$$\[ \sum\lmits_{j} \sigma_j + \sum\limits_{i,j}{{{\left( \frac{{{m}_{ij}}-{{\mu }_{i}}}{{{\sigma }_{j}}} \right)}^{2}}} \to \min \] $$

Дифференцируя по параметрам, получим систему уравнений, которая аналитически не решается.

Вообще не очевидно, почему метод максимального правдоподобия (ОМП) это то, что нужно. Например, James–Stein estimator строго лучше ОМП для оценки средних.

[B@R5uk]

stiv1995 в сообщении #1323964 писал(а):

Метод максимизации правдоподобия ведет к оптимизационной задаче

Что-то с размерностью в выражении проблемы.

-- 02.07.2018, 16:34 --

stiv1995 в сообщении #1323964 писал(а):

Например, James–Stein estimator строго лучше ОМП для оценки средних.

Вот это очень забавная ситуация, про которую я никогда не слышал. Почитал википедию, там пишут, что улучшение в целом достигается за счёт того, что ошибки перераспределяются между компонентами оцениваемого вектора. Причём, это не значит, что в среднем каждый элемент оцениваемого вектора будет оцениваться точнее, всё ровно наоборот. Но общая ошибка в среднем будет меньше. Если это то, что вам нужно, то почему бы и нет? Но меня лично бы такая оценка вряд ли бы заинтересовала.

-- 02.07.2018, 16:47 --

stiv1995 в сообщении #1323964 писал(а):

Метод максимизации правдоподобия ведет к оптимизационной задаче

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

[stiv1995]

B@R5uk в сообщении #1323981 писал(а):

Да. Сейчас сам вспомнил, да пересчитал, у вас серьёзная ошибка в формуле.

Да, логарифм потерял, когда печатал (и коэффициент).

[B@R5uk]

То есть оптимизационная задача задача будет всё-таки вот такой:

$$\[k\sum\limits_{j=1\ldots l}{\ln {{\sigma }_{j}}}+\frac{1}{2}\sum\limits_{\begin{smallmatrix} i=1\ldots k \\ j=1\ldots l \end{smallmatrix}}{\frac{{{\left( {{m}_{ij}}-{{\mu }_{i}} \right)}^{2}}}{\sigma _{j}^{2}}}\to \min \]$$

Её можно решить только численно. При этом ${\mu}_{i}$ явно выражаются через ${\sigma}_{j}$, что позволяет уменьшить размерность численной задачи. Так же наиболее эффективно искать именно минимум, а не решение системы, получающейся при дифференцировании. С точки зрения численного счёта это просто шикарная задача: минимум один и он хорошо обусловлен. В качестве начального приближение можно взять то, что вы уже предложили:

stiv1995 в сообщении #1323926 писал(а):

оценить матожидания как среднее по строкам и потом оценить дисперсии относительно этих матожиданий.

То есть вместо

$$\[{{\mu }_{i}}=\frac{\sum\limits_{j=1\ldots l}{\frac{m_{ij}^{2}}{\sigma _{j}^{2}}}}{\sum\limits_{j=1\ldots l}{\frac{1}{\sigma _{j}^{2}}}}\]$$

взять

$$\[{{\mu }_{i}}=\frac{1}{l}\sum\limits_{j=1\ldots l}{m_{ij}^{2}}\]$$

А затем, считая ${\mu}_{i}$ константами рассчитать оценку для ${\sigma}_{j}$ из уравнений:

$$\[\frac{k}{{{\sigma }_{j}}}-\sum\limits_{i=1\ldots k}{\frac{{{\left( {{m}_{ij}}-{{\mu }_{i}} \right)}^{2}}}{\sigma _{j}^{3}}}=0\]$$

stiv1995 в сообщении #1323926 писал(а):

нужна теоретическая оценка

А для чего? Если для анализа, то можно же анализировать неявные решения, хоть это и сложнее.

[alisa-lebovski]

Только квадраты при $m_{ij}$ убрать. Формула для $\mu_i$ с учетом $\sigma_j$ применяется начиная со второго шага. Получается в точности то, что я предлагала.

[B@R5uk]

alisa-lebovski в сообщении #1324150 писал(а):

Только квадраты при $m_{ij}$ убрать.

Спасибо, что указали на опечатку. Копи-паст иногда зло.

alisa-lebovski в сообщении #1324150 писал(а):

Получается в точности то, что я предлагала.

Действительно. Можно начать с ${\sigma}_{j}=1$ и последовательно применять формулы:

$$\[{{\mu }_{i}}=\frac{\sum\limits_{j=1\ldots l}{\frac{m_{ij}}{\sigma _{j}^{2}}}}{\sum\limits_{j=1\ldots l}{\frac{1}{\sigma _{j}^{2}}}}\]$$

$$\[\sigma _{j}^{2}=\frac{1}{k}\sum\limits_{i=1\ldots k}{{{\left( {{m}_{ij}}-{{\mu }_{i}} \right)}^{2}}}\]$$

Если будет быстрая сходимость, то почему бы и нет? Просто я люблю использовать матлабовскую функцию fminsearch для многомерной оптимизации, но если итерационный алгоритм будет сходится быстрее, то, разумеется, лучше использовать его, так как он ещё и концептуально более прозрачный и простой.