Оценивание дискретных распределений при дополнительных ограничениях на вероятности некоторых событий (виртуальный семинар)

Материал из MachineLearning.

Перейти к: навигация, поиск

Этот виртуальный семинар посвящён осбуждению некоторых обобщений классической задачи восстановления плотности распределения по конечной выборке данных.

Содержание

Общие постановки задач

Основные особенности рассматриваемых здесь постановок задачи:

  • имеется точная априорная информация о вероятности некоторых событий; это приводит к появлению дополнительных ограничений типа равенств в задаче максимизации правдоподобия;
  • выборка может быть «немного» неоднородной;
  • рассматривается несколько разновидностей задачи: объектами выборки могут быть как элементарные исходы, так и последовательности (временные ряды) элементарных исходов;
  • рассматриваются только дискретные распределения (множество элементарных исходов конечно);

Стационарный однородный случай

Задано конечное множество элементарных исходов \Omega. Для каждого \omega\in\Omega вероятность исхода p(\omega) неизвестна. Имеется информация двух типов:

  1. эмпирические данные: выборка наблюдений \{\omega_1,\ldots,\omega_m\},\; \omega_i\in\Omega, случайных, независимых из распределения \{p(\omega):\: \omega\in\Omega\};
  2. априорные ограничения: известны точные значения P_j вероятностей событий A_j\subseteq\Omega,\; j=1,\ldots,J:
P(A_j) = \sum_{\omega\in A_j} p(\omega_j) = P_j.

Требуется найти оценки вероятностей исходов \hat p(\omega). Эти оценки должны вычисляться достаточно эффективно — за доли секунды при m\sim 10^4,\; |\Omega|\sim 10^2.

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

Обозначим через \nu(\omega) частоту исхода \omega в выборке:

\nu(\omega) = \frac1m \sum_{i=1}^m \left[ \omega_i=\omega \right].

Непараметрическая оценка максимума правдоподобия

Найти оценку максимума правдоподобия, решив оптимизационную задачу

\sum_{\omega\in\Omega} \nu(\omega) \ln \hat p(\omega) \to \max

при ограничениях нормировки

\sum_{\omega\in\Omega}\hat p(\omega) = 1, \quad \hat p(\omega)\ge 0,

и априорных ограничениях-равенствах

\sum_{\omega\in A_j} \hat p(\omega) = P_j,\quad j=1,\ldots,J.

Вопросы:

  • Решается ли данная задача аналитически? (предположительно, да)
  • Обладают ли эти оценки свойствами несмещённости, состоятельности, эффективности? (предположительно, да)
  • Какие свойста этих оценок «испортятся», и насколько сильно, если априорная информация P(A_j)=P_j будет не согласована с неизвестным истинным распределением, то есть с эмпирическими данными? (предположительно, возникнет смещение)
  • Как число априорных ограничений влияет на дисперсию оценок? (предположительно, дисперсия уменьшается с ростом J)

Параметрическая оценка максимума правдоподобия

Эмпирических данных может оказаться не достаточно для получения надёжных оценок, особенно для маловероятных исходов. Тогда вводится ещё один тип информации — параметрическая модель распределения \hat p(\omega) = \phi(\omega,\theta), где \phi — фиксированная функция, \theta — вектор параметров модели. Постановка задачи остаётся той же, только теперь решением задачи является вектор параметров \theta.

Возможен также полупараметрический подход, когда вероятности часто встречающихся исходов (скажем, при \nu(\omega)>\nu_0) оцениваются непараметрически, а маловероятные исходы оцениваются согласно параметрической модели.

Вопросы:

  • Для каких параметрических моделей возможно получить эффективное численное решение?
  • Как определить порог \nu_0 при полупараметрическом оценивании?
  • Как ввести «размытый» порог, чтобы решение определялось моделью в тем большей степени, чем меньше \nu(\omega), без резкого перехода от параметрического оценивания к непараметрическому?

Двухэтапное решение

Для получения вычислительно эффективного метода оценивания предлагается разделить решение задачи на два этапа.

Этап 1. Оценить вероятности исходов \hat p(\omega), параметрически или непараметрически, не учитывая априорные ограничения P(A_j) = P_j. Эта задача решается стандартными методами. Например, при непараметрическом подходе оценка максимума правдоподобия есть просто

\hat p(\omega) = \nu(\omega).

Этап 2. Согласовать полученное на этапе 1 решение с априорными ограничениями. При параметрическом подходе согласование сводится к поиску таких оценок \hat p(\omega), которые в точности удовлетворяют априорным ограничениям и как можно лучше приближают модель. Например, можно воспользоваться приближением в среднеквадратичном:

\sum_{\omega\in\Omega} \left( \phi(\omega,\theta) - \hat p(\omega) \right)^2 \to \min,

при ограничениях нормировки \textstyle\sum_{\omega\in\Omega}\hat p(\omega) = 1,\; \hat p(\omega)\ge 0 и априорных ограничениях \textstyle\sum_{\omega\in A_j} \hat p(\omega) = P_j,\; j=1,\ldots,J.

Вопросы:

  • Обосновано ли применение метода наименьших квадратов (или какого-либо другого функционала) на втором этапе, если на первом этапе применяется принцип максимума правдоподобия?
  • Эквивалентно ли двухэтапное решение исходной постановке задачи? (предположительно, нет)
  • Хотя бы асимптотически? (предположительно, да)
  • Что нужно сделать, чтобы они стали эквивалентными?

Стационарный неоднородный случай

Предположим, что объекты выборки \omega_i взяты по-прежнему случайно и независимо, но теперь из разных (неизвестных) распределений \{p_i(\omega):\: \omega\in\Omega\}. Для каждого объекта известны априорные ограничения — точные значения P_{ij} вероятностей событий A_j\subseteq\Omega,\; j=1,\ldots,J. Для некоторого нового объекта \omega\in\Omega, взятого из неизвестного распределения \{p(\omega):\: \omega\in\Omega\}, также заданы априорные ограничения — точные значения P_j вероятностей событий A_j\subseteq\Omega,\; j=1,\ldots,J.

Требуется найти оценки вероятностей исходов \hat p(\omega). Эти оценки должны вычисляться достаточно эффективно.

Чтобы учесть неоднородность выборки, предлагается ввести веса объектов. Вес объекта \omega_i тем меньше, чем сильнее отличаются априорные вероятности P_{ij} для объекта \omega_i от априорных вероятностей P_{j} для объекта \omega. Далее вся методика, разработанная для однородного случая, переносится на неоднородный, с тем отличием, что теперь выборка взвешенная.

Функцию веса можно задать, опираясь на идею ядерного сглаживания:

w_i = K \left( \frac{1}{h} \textstyle \sum_{j=1}^J |P_{ij}-P_{j}| \right),

где K — неотрицательная невозростающая функция, называемая ядром; hширина окна сглаживания.

Вопросы:

  • Каким должно быть ядро?
  • Как подобрать ширину окна, иными словами, как быстро должен убывать вес с возростанием различия априорных вероятностей?
  • Какую метрику использовать для оценивания различия априорных вероятностей?
  • Будет ли оценка состоятельной, несмещённой, эффективной? Как эти свойства зависят от ширины окна?
  • Верна ли догадка, что ядерное сглаживание эквивалентно тихоновской регуляризации — введению штрафа за различия между неизвестными распределениями? Например так:
\sum_{i\neq k}\sum_{\omega}\left(p_i(\omega)-p_k(\omega)\right)^2 \to \min

Нестационарный неоднородный случай

Нестационарная (динамическая) задача является дальнейшим обобщением стационарной (статической).

Теперь объектами являются не элементарные исходы, а последовательности элементарных исходов x_i = \{\omega^1_i,\ldots,\omega^t_i,\ldots,\omega^T_i\}\in\Omega^T. Индекс t=1,\ldots,T будем называть временем. Время считается дискретным.



Задача состоит в восстановлении дискретной функции плотности вероятности f(\omega_t) (где \omega_t - элементарные исходы, зависящие от времени t \in [0, T], T < \infty, \omega_t = ( \delta_D(t-t^{(1)}_1)) * 1 + \delta_D(t-t^{(1)}_2)) * 1 + ... , \delta_D(t-t^{(2)}_1)) * 1 + \delta_D(t-t^{(2)}_2)) * 1 + ... , ), где \delta_D(.) - дельта-функция Дирака. То есть, проще говоря, события разного вида j происходят в случайные моменты времени t^{(j)}_k) ) при условии, что заданы условия на P(\omega_{A_i}) = X_{A_i} (где \omega_{A_i} - суперпозиция финальных исходов (интегрированных по времени: \omega = \int_{0}^{T} {w_t dt} = (i_1, ...,i_D);\; i_k \in Z_{0,+} \; ( k=1,D ))), P(.) - функция распределения вероятностей, X_{A_i} - заданные вероятности, i = 1,...,K).

Эмпирические частоты для \omega_t заданы.

Для несмещенных оценок вероятностей Pr'\{ X \} в качестве функционала качества предлагается использовать: q(Pr')=(1/n  \sum_ {X \in \Omega_X} {Pr\{ X \} / Pr'\{ X \} } - 1)^2, где Pr' - оценки на вероятности исходов, которые строятся из элементарных исходов интегрированием по времени и суперпозицией получившихся исходов; сумма берется по полному набору исходов (n - полное число исходов в \Omega_X), Pr\{ X \} - истинные значения вероятностей.

Частная постановка задачи

В частном случае: D=2, P(\omega_{A_1}) = \sum_{i>j; i,j \in Z_{0,+}} {Pr\{(i,j)\}} = Q_1, \; P(\omega_{A_2}) = \sum_{i<j; i,j \in Z_{0,+}} {Pr\{(i,j)\}} = Q_2, \; P(\omega_{A_3}) = \sum_{i+j \le T; i,j \in Z_{0,+}} {Pr\{(i,j)\}} = Q_3

В качестве функционала качества можно принять среднее среди функционалов качества для интегральных по времени исходов для деления всего времени на M одинаковых интервалов: q(Pr')= 1/M \sum_{l=1,M}(1/n_l  \sum_ {X_l \in \Omega_{X_l}} {Pr_l\{ X \} / Pr_l'\{ X_l \} } - 1)^2, где X_l = \int_{T/M*(l-1) + \delta_+}^{T/M * l} { (  \omega^{(1)}_t , \omega^{(2)}_t ) dt} ( \delta_+ - положительное бесконечно малое число введено, чтобы не учитывать два раза события на границе интервала). Для M=1 и D=2 множество X_l превращается в множество типа (i_1,j_1), а множество функции плотности вероятности для двух интервалов (M=2) есть ((i_1,j_1),(i_2,j_2)), где (i_1,j_1) - количества событий типа i и j, соответственно, которые произошли в интервале [0,T/2].

Известны результаты реализации этого случайного процесса, из которых можно построить эмпирическую плотность распределения f*(\omega_t).

Проблема выбора функции и параметризации для маргинальных плотностей (декомпозиция) и подгонка совместной плотности для удовлетворения связям (обратная композиция)

Проблемы:

  • допустимость перехода к маргинальных плотностям;
  • выбор класса функции плотности вероятностей;
  • выбор метода оценки параметров. Хорош ли метод максимального правдоподобия для оценки параметров функции плотности распределения, если при оценивании используются только интегральные по времени величины (интегральные эмпирическая функция и интегральная функция распределения вероятностей)).
  • метод подгонки совместной плотности для удовлетворения связям. (1) Для случая одного интервала разбиения: \min_{[ Pr**\{(i,j)\}; \lambda_k ; \sum_{i,j=0,N} {Pr**\{(i,j)\} < 1]}}{\; \; \sum_{i,j = 0,N} {( (Pr_i*\{(i)\} \, Pr_j*\{(j)\}) - Pr**\{(i,j)}\} ) ^2 + \sum_{k=1,3} { \lambda_k \, ( Pr*\{\omega_{A_k}\} - Q_k)} } , Pr**\{(i,j)\} = \int_{0}^{T}{f**(\omega_t) dt}. (2) Для случая двух интервалов разбиения: \min_{[ Pr**\{((i_1,j_1),(i_2,j_2))\}; \lambda_k ; \sum_{i1,j1,i2,j2=0,N} {Pr**\{((i1,j1),(i2,j2))\}} < 1 ]}{\; \; \sum_{i1,i2,j1,j2 = 0,N} {( (Pr_{i_1}*\{(i_1)\} \, Pr_{j_1}*\{(j_1)\} \, Pr_{i_2}*\{(i_2)\} \, Pr_{j_2}*\{(j_2)\}) - Pr**\{((i_1,j_1),(i_2,j_2))}\} ) ^2 + \sum_{k=1,3} { \lambda_k \, ( Pr*\{\omega_{A_k}\} - Q_k)} } .
    1. статистические свойства метода подгонки;
    2. возможность введения в квадратичное выражение констант C_{i,j}, зависящих от эмпирических данных для улучшения статистических свойств оценок.
  • статистические свойства и свойства функционала качества итоговой оценки плотности.

Сглаживание и подгонка совместной плотности (без декомпозиции к маргинальным плотностям)

Проблемы:

  • сглаживание совместной плотности в зависимости от эмпирических данных и значения связей при сохранении статистических свойств;
  • подгонка совместной плотности для удовлетворения связям.

Достаточно общая аппроксимация для маргинальной плотности (для рассматриваемой задачи)

Стоит задача построить достаточно общую аппроксимацию для плотности вероятностей для рассматриваемой задачи и построить быстрый алгоритм для оценки параметров. Есть основания считать, что зависимость от времени вероятности того, что событие не произойдет слабое.

Достаточно общей аппроксимацией выглядит следующая. Все время разбивается на достаточно большое число равных интервалов и принимается, что вероятность того, что событие произойдет в одной элементарном интервале один раз (q=1) мала (и всеми следующими вероятностями для q>1 можно пренебречь). Вероятности для интервала (для q=(0, 1, 2)): (p0(n), (1-p0(n)) * (1 - beta), (1-p0(n)) * beta), где величина beta будет характеризовать погрешность (и допустимость) данной модели при заданном числе интервалов N.

Теперь, подбирая последовательности p0(n) и beta можно достаточно хорошо аппроксимировать общую плотность.

Принимаем: p0(n) = p0 * Exp(a_n tau). Тогда общая плотность (для принятой гипотезы о p0(n) и beta(n)) будет выражаться в виде (Q=\sum_{i=1,N}{q_i} - полное число событий во всех интервалах). Для случая a_n = n:

P\{Q=0\} = p0^N  Exp( (N*(N+1)/2) tau )

P\{Q=1\} = p0^{(N-1)}  \sum_{ n = 1,N } { Exp( (N*(N+1)/2) - n ) * tau )  (1-beta) * ( 1 - p0 * Exp( n * tau  )  )  }

P\{Q=2\} = p0^{(N-2)}  ( \sum_{ n1 = 1,N; n2>n1 } { Exp( (N*(N+1)/2) - n1 - n2 ) * tau )  (1-beta) * ( 1 - p0 * Exp( n1 * tau )  )  (1-beta) * ( 1 - p0 * Exp( n2 * tau )  )} +  + \sum_{ n1 = 1,N } { Exp( (N*(N+1)/2) - n1 ) * tau ) * beta * ( 1 - p0 * Exp( n1 * tau ) * p0  )}

P\{Q>q,Q<N\}=...

В идеале хотелось бы построить некоторое достаточное разложение функции правдоподобия: log(L) = \sum_{Q=0;\infty} {\nu_q * log(P\{Q=q\})} = log( P\{Q=0\}) + \sum_{Q=1;\infty} {\nu_q * log(P\{Q=q\} / P\{Q=0\})}

, чтобы было возможно найти ее максимальное значение.

На первый взгляд, начальные члены разложения не должны быть слишком сложными (напрашиваются какие-то рекурсивные последовательности).

Однако, на практике оценка плотности через поиск максимума правдоподобия, как оказалось имеет некоторые особенности...

Ссылки

Литература

  1. У. Гренандер, "Вероятности на алгебраических структурах".
  2. Гнеденко, Колмогоров, "Предельные распределения для сумм независимых случайных величин".
  3. Дуб Дж. Л., "Вероятностные процессы", 1956
Личные инструменты