Мотивация: метод моментов
Метод моментов – это ещё один способ, наряду с методом максимального правдоподобия, оценки параметров распределения по данным $x_1,\ldots,x_N$. Суть его в том, что мы выражаем через параметры распределения теоретические значения моментов $\mu_k = \mathbb{E}x^k$ нашей случайной величины, затем считаем их выборочные оценки $\widehat{\mu}_k = \frac1N\sum_ix_i^k$, приравниваем их все друг к другу и, решая полученную систему, находим оценки параметров. Можно доказать, что полученные оценки являются состоятельными, хотя могут быть смещены.
Пример 1. Оценим параметры нормального распределения $\mathcal{N}(\mu, \sigma^2)$ с помощью метода моментов.
Попробуйте сделать сами, прежде чем смотреть решение.
$$\mu_1 = \mu,\quad\mu_2 = \sigma^2 + \mu^2$$
Запишем систему:
$$\begin{cases} \widehat{\mu} = \frac1N\sum_i x_i,\ \widehat{\sigma}^2 + \widehat{\mu}^2 = \frac1N\sum_ix_i^2 \end{cases}$$
Из неё очевидным образом находим
$$\widehat{\mu} = \frac1N\sum_ix_i $$
$$\widehat{\sigma}^2 = \frac1N\sum_ix_i^2 - \left(\frac1N\sum_i x_i\right)^2=$$
$$=\frac1N\sum_i\left(x_i - \widehat{\mu}\right)^2$$
Легко видеть, что полученные оценки совпадают с оценками максимального правдоподобия
Пример 2. Оценим параметр $\mu$ логнормального распределения
$$p(x) = \frac1{x\sqrt{2\pi\sigma^2}}\exp\left(-\frac{(\log{x} - \mu)^2}{2\sigma^2}\right)$$
при известном $\sigma^2$. Будет ли оценка совпадать с оценкой, полученной с помощью метода максимального правдоподобия?
Попробуйте сделать сами, прежде чем смотреть решение.
Теперь запишем логарифм правдоподобия:
$$l(X) = -\sum_i\log{x_i} - \sum_i\frac{(\log{x_i} - \mu)^2}{2\sigma^2} + const$$
Дифференцируя по $\mu$ и приравнивая производную к нулю, получаем
$$\widehat{\mu}_{MLE} = \frac1N\sum_i\log{x_i}$$
что вовсе не совпадает с оценкой выше.
Несколько приукрасив ситуацию, можно сделать вывод, что первые два выборочных момента позволяют если не править миром, то уверенно восстанавливать параметры распределений. А теперь давайте представим, что мы посчитали $\frac1N\sum_ix_i$ и $\frac1N\sum_ix_i^2$, а семейство распределений пока не выбрали. Как же совершить этот судьбоносный выбор? Давайте посмотрим на следующие три семейства и подумаем, в каком из них мы бы стали искать распределение, зная его истинные матожидание и дисперсию?
Почему-то хочется сказать, что в первом. Почему? Второе не симметрично – но что нас может заставить подозревать, что интересующее нас распределение не симметрично? С третьим проблема в том, что, выбирая его, мы добавляем дополнительную информацию как минимум о том, что у распределения конечный носитель. А с чего бы? У нас такой инфомации вроде бы нет.
Общая идея такова: мы будем искать распределение, которое удовлетворяет только явно заданным нами ограничениям и не отражает никакого дополнительного знания о нём. Но чтобы эти нестрогие рассуждения превратить в формулы, придётся немного обогатить наш математический аппарат и научиться измерять количество информации.
Энтропия и дивергенция Кульбака-Лейблера
Измерять «знание» можно с помощью энтропии Шэннона. Она определяется как
$$\color{#348FEA}{H(P) = -\sum_xP(x)\log{P(x)}}$$
для дискретного распределения и
$$\color{#348FEA}{H(p) = -\int p(x)\log{p(x)}dx}$$
для непрерывного. В классическом определении логарифм двоичный, хотя, конечно, варианты с разным основанием отличаются лишь умножением на константу.
Неформально можно представлять, что энтропия показывает, насколько сложно предсказать значение случайной величины. Чуть более строго – сколько в среднем бит нужно потратить, чтобы передать информацию о её значении.
Пример 1. Рассмотрим схему Бернулли с вероятностью успеха $p$. Энтропия её результата равна
$$-(1 - p)\cdot\log_2(1 - p) - p\cdot\log_2{p}$$
Давайте посмотрим на график этой функции:
Минимальное значение (нулевое) энтропия принимает при $p\in{0,1}$. В самом деле, для такого эксперимента мы всегда можем наверняка сказать, каков будет его исход; обращаясь к другой интерпретации – чтобы сообщить кому-то о результате эксперимента, достаточно $0$ бит (ведь получатель сообщения и так понимает, что вышло).
Максимальное значение принимается в точке $\frac12$, что вполне соответствует тому, что при $p=\frac12$ предсказать исход эксперимента сложнее всего.
Дополнение для ценителей математики.
Теперь пусть $p$ произвольно. Рассмотрим $N»1$ независимых испытаний $x_1,\ldots, x_N$; среди них будет $n_0\approx (1-p)N$ неудачных и $n_1\approx pN$ удачных. Посчитаем, сколько бит потребуется, чтобы закодировать последовательность $x_i$ для известных $n_0$ и $n_1$. Общее число таких последовательностей равно $C_N^{n_1} = \frac{N!}{n_0!n_1!}$, а чтобы закодировать каждую достаточно будет $\log_2\left(\frac{N!}{n_0!n_1!}\right)$ бит – это количество информации, содержащееся во всей последовательности. Таким образом, в среднем, чтобы закодировать результат одного испытания необходимо
$$\frac1N\log_2\left(\frac{N!}{n_0!n_1!}\right)$$
бит информации. Перепишем это выражение, использовав формулу Стирлинга $\log{N!}\approx N\log{N} - N$:
$$\frac1N\left(\log_2{N!} - \log_2{n_0!} - \log_2{n_1!}\right) \approx $$
$$\approx const\cdot\frac1N\left(N\log_2{N} - N - n_0\log_2{n_0} + n_0 - n_1\log_2{n_1} + n_1\right) =$$
$$=const\cdot\frac1N\left(N\log_2{N} - (1-p)N\log_2{(1-p)N} - pN\log_2{pN}\right) =$$
$$=const\cdot\left(-(1-p)\cdot\log_2(1-p) - p\log_2{p}\right)$$
Вот мы и вывели формулу энтропии!
Пример 2. Энтропия нормального распределения $\mathcal{N}(\mu, \sigma^2)$ равна $\frac12\log(2\pi\sigma^2) + \frac12$, и чем меньше дисперсия, тем меньше энтропия, что и логично: ведь когда дисперсия мала, значения сосредоточены возле матожидания, и они становятся менее «разнообразными».
Энтропия тесно связана с другим важным понятием из теории информации – дивергенцией Кульбака-Лейблера. Она определяется для плотностей $p(x)$ и $q(x)$ как
$$\color{#348FEA}{KL(p,\vert\vert, q) = \int p(x)\log{\frac{p(x)}{q(x)}}dx}$$
в непрерывном случае и точно так же, но только с суммой вместо интеграла в дискретном.
Теоретико-информационный смысл дивергенции Кульбака-Лейблера.
Дивергенцию можно представить в виде разности:
$$KL(p,\vert\vert, q) = (-\int p(x)\log{q(x)}dx) - (-\int p(x)\log{p(x)}dx)$$
Вычитаемое – это энтропия, которая, как мы уже поняли, показывает, сколько в среднем бит требуется, чтобы закодировать значение случайной величины. Уменьшаемое похоже по виду, и можно показать, что оно говорит о том, сколько в среднем бит потребуется на кодирование случайной величины с плотностью $p$ алгоритмом, оптимизированным для кодирования случайной величины $q$. Иными словами, дивергенция Кульбака-Лейблера говорит о том, насколько увеличится средняя длина кодов для значений $p$, если при настройке алгоритма кодирования вместо $p$ использовать $q$. Более подробно вы можете почитать, например, в этом посте.
Дивергенция Кульбака-Лейблера в некотором роде играет роль расстояния между распределениями. В частности, $KL(p,\vert\vert, q)\geqslant0$, причём дивергенция равна нулю, только если распределения совпадают почти всюду. Но при этом она не является симметричной: вообще говоря, $KL(p,\vert\vert, q)\ne KL(q,\vert\vert, p)$.
Вопрос на подумать. Пусть $p(x)$ – распределение, заданное на отрезке $[a, b]$. Выразите его энтропию через дивергенцию Кульбака-Лейблера $p(x)$ с равномерным на отрезке распределением $q_U(x)=\frac1{b-a}\mathbb{I}_{[a,b]}(x)$.
Попробуйте вывести сами, прежде чем смотреть решение.
$$KL(p,\vert\vert, q_U) = -\left(-\int_a^b p(x)\log{p(x)}dx\right) - \int_a^b p(x)\log{\underbrace{q(x)}_{=\frac1{b-a}\text{ на }[a,b]}}dx=$$
$$=\log(b-a) - H(p)$$
Аналогичное соотношение можно выписать и для распределения, заданного на конечном множестве.
Принцип максимальной энтропии
Теперь наконец мы готовы сформулировать, какие распределения мы хотим искать.
Принцип максимальной энтропии. Среди всех распределений на заданном носителе $\mathbb{X}$, удовлетворяющих условиям $\mathbb{E}u_1(x) = \mu_1$, ..., $\mathbb{E}u_k(x) = \mu_k$, где $u_i$ – некоторые функции, мы хотим иметь дело с тем, которое имеет наибольшую энтропию.
В самом деле, энтропия выражает нашу меру незнания о том, как ведёт себя распределение, и чем она больше – тем более «произвольное распределение», по крайней мере, в теории.
Давайте рассмотрим несколько примеров, которые помогут ещё лучше понять, почему некоторые распределения так популярны:
Пример 1. На конечном множестве $1,\ldots,n$ наибольшую энтропию имеет равномерное распределение (носитель – конечное множество из $n$ элементов, других ограничений нет).
Доказательство.
$$KL(p\vert\vert q) = \sum_i p_i\log{p_i} - \sum_i p_i\log{q_i} =$$
$$= -H(p) + \log{n}\underbrace{\sum_ip_i}_{=1}$$
Так как дивергенция Кульбака-Лейблера всегда неотрицательна, получаем, что $H(p)\leqslant\log{n}$. При этом равенство возможно, только если распределения совпадают.
Пример 2. Среди распределений, заданных на всей вещественной прямой и имеющих заданные матожидание $\mu$ и дисперсию $\sigma^2$ наибольшую энтропию имеет нормальное распределение $\mathcal{N}(\mu,\sigma^2)$.
Доказательство.
Пусть $p(x)$ – некоторое распределение, $q(x)\sim\mathcal{N}(\mu, \sigma^2)$. Запишем их дивергенцию Кульбака-Лейблера:
$$KL(p\vert\vert q) = \int p(x)\log{p(x)}dx - \int p(x)\log{q(x)}dx =$$
$$= -H(p) - \int p(x)\left(-\frac12\log(2\pi\sigma^2) - \frac1{2\sigma^2}(x - \mu)^2\right)dx =$$
$$= - H(p) +\frac12\log(2\pi\sigma^2)\cdot\underbrace{\int p(x)dx}_{=1} + \frac1{2\sigma^2}\underbrace{\int(x - \mu)^2p(x)dx}_{=\mathbb{V}p=\sigma^2} =$$
$$= - H(p) + \underbrace{\frac12\log(2\pi\sigma^2) + \frac12}_{=H(q)}$$
Так как дивергенция Кульбака-Лейблера всегда неотрицательна, получаем, что $H(p)\leqslant H(q)$. При этом равенство возможно, только если распределения $p$ и $q$ совпадают почти всюду, а с точки зрения теории вероятностей такие распределения различать не имеет смысла.
Пример 3. Среди распределений, заданных на множестве положительных вещественных чисел и имеющих заданное матожидание $\lambda$ наибольшую энтропию имеет показательное распределение с параметром $\frac1{\lambda}$ (его плотность равна $p(x) = \frac1{\lambda}\exp\left(-\frac1{\lambda}x\right)\mathbb{I}_{(0;+\infty)}(x)$).
Все хорошо знакомые нам распределения, не правда ли? Проблема в том, что они свалились на нас чудесным образом. Возникает вопрос, можно ли их было не угадать, а вывести как-нибудь? И как быть, если даны не эти конкретные, а какие-то другие ограничения? Оказывается, что при некоторых не очень обременительных ограничениях ответ можно записать с помощью распределений экспоненциального класса. Давайте же познакомимся с ними поближе.
Экспоненциальное семейство распределений
Говорят, что семейство распределений относится к экспоненциальному классу, если оно может быть представлено в следующем виде:
$$\color{#348FEA}{p(x\vert\theta) = \frac1{h(\theta)}g(x)\cdot\exp\left(\theta^Tu(x)\right)}$$
где $\theta$ – вектор вещественнозначных параметров (различные значения которых дают те или иные распределения из семейства), $h, g > 0$, $u$ – некоторая вектор-функция, и, разумеется, сумма или интеграл по $x$ равняется единице. Последнее, в частности, означает, что
$$h(\theta) = \int g(x)\exp\left(\theta^Tu(x)\right)dx$$
(или сумма в дискретном случае).
Пример 1. Покажем, что нормальное распределение принадлежит экспоненциальному классу. Для этого мы должны представить привычную нам функцию плотности
$$p(x \vert \mu, \sigma^2) = \frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right)$$
в виде
$$p(x\vert\theta) = \frac{g(x)\cdot\exp\left(\sum_i\text{(параметр)}_i\cdot\text{(функция от x)}_i\right)}{\text{что-то, не зависящее от $x$}}$$
Распишем
$$\frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right) = \frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac1{2\sigma^2}x^2 + \frac{\mu}{\sigma^2}x - \frac{\mu^2}{2\sigma^2}\right)=$$ $$ \frac{\exp\left(-\frac1{2\sigma^2}x^2 + \frac{\mu}{\sigma^2}x\right)}{\sqrt{2\pi}\sigma\exp\left(\frac{\mu^2}{2\sigma^2}\right)}=$$
Определим
$$u_1(x) = x,\qquad u_2(x) = x^2$$
$$\theta_1 = \frac{\mu}{\sigma^2},\quad \theta_2 = -\frac1{2\sigma^2}$$
$$h(\theta) = \sqrt{2\pi}\sigma\exp\left(\frac{\mu^2}{2\sigma^2}\right)$$
Если теперь всё-таки честно выразить $h$ через $\theta$ (это мы оставляем в качестве лёгкого упражнения), то получится
$$p(x \vert \mu, \sigma^2) = \frac1{h(\theta)}\exp\left(\theta^Tu(x)\right)$$
В данном случае функция $g(x)$ просто равна единице.
Пример 2. Покажем, что распределение Бернулли принадлежит экспоненциальному классу. Для этого попробуем преобразовать функцию вероятности (ниже $x$ принимает значения $0$ или $1$):
$$P(x \vert p) = p^x(1 - p)^{1 - x} = \exp\left(x\log{p} + (1 - x)\log(1 - p)\right)$$
Теперь мы можем положить $u(x) = \left(x, 1 - x\right)$, $\theta = \left(\log{p}, \log(1 – p)\right)$, и всё получится. Единственное, что смущает, – это то, что компоненты вектора $u(x)$ линейно зависимы. Хотя это не является формальной проблемой, но всё же хочется с этим что-то сделать. Исправить это можно, если переписать
$$p^x(1 - p)^{1 -x} = (1 - p)\exp\left(x\log{p} + (-x)\log(1 - p)\right) =$$
$$=(1 - p)\exp\left(x\log{\frac{p}{1 - p}}\right)$$
и определить уже минимальное представление с $u(x) = x$, $\theta = \log{\frac{p}{1 - p}}$ (мы ведь уже сталкивались с этим выражением, когда изучали логистическую регрессию, не так ли?).
Вопрос на подумать. Принадлежит ли к экспоненциальному классу семейство равномерных распределений на отрезках $U[a, b]$? Казалось бы, да: так как:
$$p(x) = \frac{1}{b - a}\mathbb{I}_{[a,b]}(x)\exp(0)$$
В чём может быть подвох?
Попробуйте определить сами, прежде чем смотреть ответ.
, и если нас интересует семейство равномерных распределений на отрезках, определяемое параметрами $a$ и $b$, то они не могут быть в функции $g(x)$, они должны быть под экспонентой, а экспонента ни от чего не может быть равна индикатору.
При этом странное и не очень полезное семейство с нулём параметров, состоящее из одинокого распределения $U[0,1]$ можно считать относящимся к экспоненциальному классу: ведь для него формула
$$p(x) = \mathbb{I}_{[0,1]}(x)\exp(0)$$
будет работать.
Как мы увидели, к экспоненциальным семействам относятся как непрерывные, так и дискретные распределения. Вообще, к ним относится большая часть распределений, которыми Вам на практике может захотеться описать $Y \vert X$. В том числе
- нормальное
- распределение Пуассона
- экспоненциальное
- биномиальное, мультиномиальное (с фиксированным числом испытаний)
- геометрическое
- $\chi^2$-распределение
- бета-распределение
- гамма-распределение
- распределение Дирихле
К экспоненциальным семействам не относятся, к примеру: семейство равномерных распределений на отрезке, семейство $t$-распределений Стьюдента, семейство распределений Коши, смеси нормальных распределений.
MLE для семейства из экспоненциального класса
Возможно, вас удивил странный и на первый взгляд не очень естественный вид $p(x\vert\theta)$. Но всё не просто так: оказывается, что оценка максимального правдоподобия параметров распределений из экспоненциального класса устроена очень интригующе.
Запишем функцию правдоподобия выборки $X = (x_1,\ldots,x_N)$:
$$p(X\vert\theta) = h(\theta)^{-N}\cdot\left(\prod_{i=1}^Ng(x_i)\right)\cdot\exp\left(\theta^T\left[\sum_{i=1}^Nu(x_i)\right]\right)$$
Её логарифм равен
$$l(X\vert\theta) = -N\log{h(\theta)} + \sum_{i=1}^N\log{g(x_i)} + \theta^T\left[\sum_{i=1}^Nu(x_i)\right]$$
Дифференцируя по $\theta$, получаем
$$\nabla_{\theta}l(X\vert\theta) = -N\nabla_{\theta}\log{h(\theta)} + \left[\sum_{i=1}^Nu(x_i)\right]$$
Тут нам потребуется следующая
Лемма. $\nabla_{\theta}\log{h(\theta)} = \mathbb{E}u(x)$
Доказательство.
$$h(\theta) = \int g(x)\exp\left(\theta^Tu(x)\right)dx$$
Следовательно,
$$\nabla_{\theta}\log{h(\theta)} = \frac{\nabla_{\theta}\int g(x)\exp\left(\theta^Tu(x)\right)dx}{\int g(x)\exp\left(\theta^Tu(x)\right)dx} =$$
$$= \frac{\int u(x)g(x)\exp\left(\theta^Tu(x)\right)dx}{h(\theta)} =$$
$$=\int u(x)\cdot\frac1{h(\theta)}g(x)\exp\left(\theta^Tu(x)\right)dx = \mathbb{E}u(x)$$
Кстати, можно ещё доказать, что
$$\frac{\partial^2}{\partial \theta_i\partial\theta_j}\log{h(\theta)} = \text{Cov}(u_i(x), u_j(x))$$
Приравнивая $\nabla_{\theta}l(X\vert\theta)$ к нулю и применяя лемму, мы получаем, что
$$\color{#348FEA}{\mathbb{E}u(x) = \frac1N\left[\sum_{i=1}^Nu(x_i)\right]}$$
Таким образом, теоретические матожидания всех компонент $u_i(x)$ должны совпадать с их эмпирическими оценками, а метод максимального правдоподобия совпадает с методом моментов для $\mathbb{E}u_i(x)$ в качестве моментов. И в следующем пункте выяснится, что распределения из семейств, относящихся к экспоненциальному классу, это те самые распределения, которые имеют максимальную энтропию из тех, что имеют заданные моменты $\mathbb{E}u_i(x)$.
Пример.
$$p(x) = \frac1{x\sqrt{2\pi\sigma^2}}\exp\left(-\frac{(\log{x} - \mu)^2}{2\sigma^2}\right) =$$
$$=\frac1{x\sqrt{2\pi\sigma^2}}\exp\left(-\frac1{2\sigma^2}\log^2{x} + \frac{\mu}{\sigma^2}\log{x} - \frac{\mu^2}{2\sigma^2}\right) =$$
$$=\frac1{x\sqrt{2\pi\sigma^2}\exp\left(\frac{\mu^2}{2\sigma^2}\right)}\exp\left(\underbrace{\frac{\mu}{\sigma^2}}_{=\theta_1}\underbrace{\log{x}}_{=u_1(x)} -\underbrace{\frac1{2\sigma^2}}_{=\theta_2}\underbrace{\log^2{x}}_{=u_2(x)} \right) =$$
$$\frac{1}{\sqrt{-\pi\theta_2^{-1}}\cdot\exp{-\frac{\theta_1^2}{4\theta_2}}}\cdot\frac1x\exp\left(\theta_1u_1(x) + \theta_2u_2(x)\right)$$
Как видим, логнормальное распределение тоже из экспоненциального класса. Вас может это удивить: ведь выше мы обсуждали, что для него метод моментов и метод максимального правдоподобия дают разные оценки. Но никакого подвоха тут нет: мы просто брали не те моменты. В данном случае $u_1(x) = \log{x}$, $u_2(x) = \log^2{x}$, их матожидания и надо брать; тогда для параметров, получаемых из MLE, должно выполняться
$$\mathbb{E}\log{x} = \frac1N\sum_i\log{x_i},\quad \mathbb{E}\log^2{x} = \frac1N\sum_i\log^2{x_i}$$
Матожидания в левых частях мы должны выразить через параметры – и нам для этого совершенно не обязательно что-то интегрировать! В самом деле:
$$\mathbb{E}\log{x} = \frac{\partial}{\partial\theta_1}\log{h(\theta)} =$$
$$=\frac{\partial}{\partial\theta_1}\left(-\frac12\log{\pi} + \frac12\log{\theta_2} - \frac{\theta_1^2}{4\theta_2^2}\right) = -\frac{\theta_1}{2\theta_2^2}$$
$$\mathbb{E}\log^2{x} = \frac{\partial}{\partial\theta_2}\log{h(\theta)} = \frac1{2\theta_2} + \frac{\theta_1^2}{2\theta^3}$$
Теорема Купмана-Питмана-Дармуа
Теперь мы наконец готовы сформулировать одно из самых любопытных свойств семейств экспоненциального класса.
В следующей теореме мы опустим некоторые не очень обременительные условия регулярности. Просто считайте, что для хороших дискретных и абсолютно непрерывных распределений, с которыми вы в основном и будете сталкиваться, это так.
Теорема. Пусть $p(x) = \frac1{h(\theta)}\exp\left(\theta^Tu(x)\right)$ – распределение, причём $\theta$ – вектор длины $n$ и $\mathbb{E}u_i(x) = \alpha_i$ для некоторых фиксированных $\alpha_i$, $i=1,\ldots,n$. Тогда распределение $p(x)$ обладает наибольшей энтропией среди распределений с тем же носителем, для которых $\mathbb{E}u_i(x) = \alpha_i$, $i=1,\ldots,n$. При этом оно является единственным с таким свойством (в том смысле, что любое другое распределение, обладающее этим свойством, совпадает с ним почти всюду).
Идея обоснования через оптимизацию.
В дискретном случае у нас есть счётное семейство точек $x_1, x_2,\ldots$, и распределение определяется счётным набором вероятностей $p_i$ принимать значение $x_i$. Мы будем решать задачу
$$\begin{cases} -\sum_j p_j\log{p_j}\longrightarrow\max,\ \sum_jp_ju_i(x_j) = \alpha_i, i = 1,\ldots,n,\ \sum_jp_j = 1,\ p_j\geqslant0 \end{cases}$$
Запишем лагранжиан:
$$\mathcal{L} = \sum_j p_j\log{p_j} + \sum_i\theta_i\left(\alpha_i - \sum_jp_ju_i(x_j)\right)+$$
$$+\theta_0\left(\sum_jp_j - 1\right) - \sum_j\lambda_jp_j$$
Продифференцируем его по $p_j$:
$$\frac{\partial\mathcal{L}}{\partial p_j} = \log{p_j} + 1 - \sum_i\theta_iu_i(x_j) + \theta_0 - \lambda_j$$
Приравнивая это к нулю, получаем
$$p_j = \frac{\exp\left(\langle \theta, u(x_j)\rangle\right)}{\exp\left(\lambda_j - \theta_0 - 1\right)}$$
Числитель уже ровно такой, как и должен быть у распределения из экспоненциального класса; разберёмся со знаменателем. Во-первых, легко видеть, что условие $p_j\geqslant0$ заведомо выполнено (ведь тут сплошные экспоненты), так что его можно было выкинуть из постановки задачи оптимизации или, что то же самое, положить $\lambda_j = 0$. Параметр $\theta_0$ находится из условия $\sum_jp_j = 1$, а точнее, выражается через остальные $\theta_i$, что позволяет записать знаменатель в виде $h(\theta)$.
Идея доказательства «в лоб».
$$\int u_i(x)q(x)dx = \int u_i(x)p(x)dx$$
для всех $i = 1,\ldots,n$. Тогда
$$0\leqslant KL(q\vert\vert p) = \int q(x)\log\left(\frac{q(x)}{p(x)}\right)dx = $$
$$=\underbrace{\int q(x)\log{q(x)}dx}_{-H(q)} - \int q(x)\log{p(x)}dx=$$
$$=-H(q) - \int q(x)\left(-\log{h(\theta)} + \sum_i\theta_iu_i(x)\right)dx =$$
$$=-H(q) - \log{h(\theta)}\underbrace{\int q(x)dx}_{=1=\int p(x)dx} - \sum_i\theta_i\underbrace{\int q(x)u_i(x)dx}_{=\int p(x)u_i(x)dx} =$$
$$=-H(q) - \int p(x)\left(-\log{h(\theta)} + \sum_i\theta_iu_i(x)\right)dx =$$
$$=-H(q) + \int p(x)\log{p(x)}dx = -H(q) + H(p)$$
Таким образом, $H(p)\geqslant H(q)$, причём по уже не раз использованному нами свойству дивергенции Кульбака-Лейблера из $H(p) = H(q)$ будет следовать то, что $p$ и $q$ совпадают почти всюду.
Рассмотрим несколько примеров:
Пример 1. Среди распределений на множестве ${1,2,3,\ldots}$ неотрицательных целых чисел с заданным математическим ожиданием $\mu$ найдём распределение с максимальной энтропией.
В данном случае у нас лишь одна функция $u_1(x) = x$, которая соответствует фиксации матожидания $\mathbb{E}x$. Плотность будет вычисляться только в точках $x=k$, $k=1,2,\ldots$ и будет иметь вид
$$p_k = p(k) = \frac1{h(\theta)}\exp\left(\theta k\right)$$
В этой формуле уже безошибочно угадывается геометрическое распределение с $p = 1 - e^{\theta}$. Параметр $p$ можно подобрать из соображений того, что математическое ожидание равно $\mu$. Матожидание геометрического распределения равно $\frac1p$, так что $p = \frac1{\mu}$. Окончательно,
$$p_k = \frac1{\mu}\left(1 - \frac1{\mu}\right)^{k-1}$$
Пример 2. Среди распределений на всей вещественной прямой с заданным математическим ожиданием $\mu$ найдём распределение с максимальной энтропией.
А сможете ли вы его найти? Решение под катом.
$$p(x) = \frac1{h(\theta)}\exp\left(\theta x\right)$$
но интеграла экспоненты не существует, то есть применение «в лоб» теоремы провалилось. И неспроста: если даже рассмотреть все нормально распределённые случайные величины со средним $\mu$, их энтропии, равные $\frac12 + \frac12\log(2\pi\sigma^2)$ не ограничены сверху, то есть величины с наибольшей энтропией не существует.