Энтропия
Информативность наблюдений
Представим, что вы пронаблюдали значение $x$ некоторой случайной величины $\xi$. Как измерить количество информации $h(x)$, которое вы при этом получили? Следующие соображения кажутся в этом плане вполне естественными:
- чем выше вероятность $\mathbb P(\xi = x)$, тем более ожидаемо появление значения $x$ и, соответственно, менее информативно;
- и наоборот, наблюдение маловероятного значения $x$ обычно даёт обильную пищу для размышлений и повышает $h(x)$;
- при наблюдении двух независимых реализаций $x$ и $y$ случайной величины $\xi$ логично складывать полученную информацию: $h((x, y)) = h(x) + h(y)$.
Указанные соображения наводят на мысль, что информацию следует считать убывающей функцией от вероятности: $h(x) = g\big(\mathbb P(\xi = x)\big)$. Кроме того, функция $g$ должна превращать произведение в сумму, поскольку для независимых случайных величин $\xi$ и $\eta$ равенство $h((x, y)) = h(x) + h(y)$ влечёт
$$ g\big(\mathbb P(\xi = x, \eta=y)\big) = g\big(\mathbb P(\xi = x) \mathbb P(\eta = y)\big) = g\big(\mathbb P(\xi = x)\big) + g\big(\mathbb P(\eta = y)\big). $$
На самом деле выбор тут небогат. Единственная непрерывная функция, обладающая такими свойствами, — это логарифм: $g(p) = -\log p$. Основание логарифма может быть любым числом больше единицы. Поскольку информацию измеряют в битах и байтах, в теории информации обычно предпочитают двоичные логарифмы. Однако для вычислений удобнее использовать натуральный логарифм, и по умолчанию мы будем подразумевать именно его (кстати, соответствующую единицу информации называют «нат»).
Энтропия Шеннона
Среднее количество информации, которое несёт в себе значение дискретной случайной $\xi$ с распределением вероятностей $p_k= \mathbb P(\xi = k)$, вычисляется по формуле
$$ \mathbb H\xi = \mathbb E \big(g(p(\xi))\big) = -\mathbb E \big(\log p(\xi)\big)=-\sum\limits_k p_k \log p_k. $$
Это так называемая энтропия (Шэннона).
Небольшое математическое замечание
Пример. Рассмотрим схему Бернулли с вероятностью «успеха» $p$. Энтропия её результата равна
$$ \mathbb H\xi = -(1 - p)\log(1 - p) - p\log p, \quad \xi \sim \mathrm{Bern}(p). $$
Давайте посмотрим на график этой функции:
Минимальное значение (нулевое) энтропия принимает при $p = 0$ или $p=1$. Исход такого вырожденного эксперимента заранее известен, и чтобы сообщить кому-то о его результате, достаточно $0$ бит информации. Иначе говоря, можно вообще ничего не передавать, и так всё предельно ясно.
Максимальное значение энтропии достигается в точке $\frac12$, что вполне соответствует тому, что при $p=\frac12$ предсказать исход эксперимента сложнее всего.
Упражнение. Найдите энтропию геометрического распределения с вероятостью «успеха» $p$: $\xi \sim \mathrm{Geom}(p)$,
$$ \mathbb P(\xi = k) = p(1-p)^{k-1}, \quad k \in \mathbb N, \quad 0 < p \leqslant 1. $$
Решение (не открывайте сразу, попробуйте сначала решить самостоятельно)
$$ \mathbb H\xi = -\sum\limits_{k=1}^\infty p(1-p)^{k-1} (\log p + (k-1)\log(1-p)) = $$ $$ = -\log p -p\log(1-p)\sum\limits_{k=1}^\infty k(1-p)^k. $$
Дифференцированием геометрической прогрессии находим
$$ \sum\limits_{k=1}^\infty kx^{k-1} = \sum\limits_{k=1}^\infty (x^k)' = \frac d{dx}\sum\limits_{k=0}^\infty x^k = \Big(\frac 1{1-x}\Big)' = \frac 1{(1-x)^2}; $$
подставляя сюда $x = 1-p$, окончательно получаем, что
$$ \mathbb H\xi = -\log p -p(1-p)\log(1-p)\cdot \frac 1{p^2} = -\log p - \frac{(1-p)\log(1-p)}p. $$
Вспомним, что случайная величина $\xi\sim\mathrm{Geom}(p)$ равна количеству независимых испытаний с вероятностью «успеха» $p$ до появления первого «успеха». При $p=1$ энтропия $\mathbb H\xi$ минимальна и равна нулю, ведь в этом случае геометрическое распределение становится вырожденным: «успех» наступает сразу же с вероятностью $1$. А вот при $p\to 0+$ серия неудачных испытаний может быть сколь угодно длинной; распределение становится всё более неопределённым и «размазанным» по своему носителю, и его энтропия
стремится к $+\infty$. График энтропии как функции от $p$ выглядит так:
Следующие свойства энтропии дискретной случайной величины $\xi$ вытекают прямо из определения:
- неотрицательность: $\mathbb H\xi \geqslant 0$;
- $\mathbb H\xi = 0 \iff \mathbb P(\xi = a) = 1$ при некотором $a\in\mathbb R$ (нулевую энтропию имеют вырожденные распределения и только они);
- $\mathbb H\xi \leqslant \log n$, если случайная величина имеет конечный носитель мощности $n$.
Последнее свойство выводится из неравенства Йенсена. Применяя его к выпуклой вверх логарифмической функции, с учётом нормировки условия $\sum\limits_{k=1}^n p_k = 1$ получаем
$$ -\sum\limits_{k=1}^n p_k\log p_k = \sum\limits_{k=1}^n p_k\log \frac 1p_k \leqslant \log\bigg(\sum\limits_{k=1}^n p_k\cdot \frac 1{p_k}\bigg) = \log n. $$
Ответ (не открывайте сразу, сначала подумайте самостоятельно)
Дифференциальная энтропия
Чтобы вычислить энтропию непрерывной случайной величины $\xi$, надо, как водится, сумму заменить на интеграл, и получится формула дифференциальной энтропии:
$$ \mathbb H\xi = -\int p_{\xi}(x) \log p_{\xi}(x), dx. $$
Замечание. В дальнейшем мы будем использовать одинаковый термин энтропия как для дискретных, так и для непрерывных случайных величин, для краткости опуская слово дифференциальная в последнем случае. Кроме того, энтропию распределения $p$, заданного через pmf или pdf, будем обозначать $\mathbb H[p]$. Такое обозначение позволяет избежать привязки к случайной величине там, где это излишне. Если $\xi \sim p(x)$, то обозначения $\mathbb H\xi$ и $\mathbb H[p]$ эквивалентны. Также отметим, что энтропию можно записать в виде математического ожидания:
$$ \mathbb H[p] = \mathbb E_{\xi \sim p(x)} \log\frac 1{p(\xi)}. $$
Пример. Найдём энтропию нормального распределения $\mathcal N(\mu, \sigma^2)$. Его плотность равна $p(x) = \frac 1{\sqrt{2\pi} \sigma} e^{-\frac{(x-\mu)^2}{2\sigma^2}}$, следовательно,
Делая в последнем интеграле замену $t = -\frac{(x-\mu)^2}{2\sigma^2}$, получаем, что
По свойству гамма-функции $\Gamma\big(\frac 32\big) = \frac 12 \Gamma\big(\frac 12\big) = \frac{\sqrt \pi}2$. Таким образом,
$$ \mathbb H[p] = \frac 12\ln\big(2\pi\sigma^2\big) + \frac 12. $$
Как видно, энтропия гауссовского распределения $\mathcal N(\mu, \sigma^2)$ неограничена ни сверху, ни снизу:
$$ \lim\limits_{\sigma\to +\infty} \frac 12\ln\big(2\pi\sigma^2\big) = +\infty, \quad \lim\limits_{\sigma\to 0+} \frac 12\ln\big(2\pi\sigma^2\big) = -\infty. $$
И да, в отличие от энтропии дискретного распределения дифференциальная энтропия может быть отрицательной. Это связано с тем, что плотность может принимать значения больше единицы, и поэтому математическое ожидание её логарифма с обратным знаком может оказаться меньше нуля. В частности, с нормальным распределением так происходит, если $\sigma^2 < \frac 1{2\pi e}$.
Упражнение. Найдите энтропию показательного распределения $\mathrm{Exp}(\lambda)$.
Решение (не открывайте сразу, попробуйте сначала решить самостоятельно)
С помощью замены $t = \lambda x$ находим, что
KL-дивергенция
В задачах машинного обучения истинное распределение $p(x)$, из которого приходят наблюдения, обычно неизвестно, и его пытаются приблизить распределением $q(x)$ из некоторого класса модельных распределений. Дивергенция Кульбака-Лейблера (KL-дивергенция, относительная энтропия) позволяет оценить расстояние между распределениями $p$ и $q$:
$$\mathbb{KL}(p\vert\vert q) = \sum\limits_k p_k\log\frac{p_k}{q_k}$$
в дискретном случае и
$$\mathbb{KL}(p\vert\vert q) = \int p(x)\log \frac{p(x)}{q(x)}dx$$
в непрерывном. KL-дивергенцию можно представить в виде разности:
$$ \mathbb{KL}(p\vert\vert q) = \int p(x)\log \frac 1{q(x)}dx - \int p(x)\log\frac 1{p(x)}dx = $$
Здесь вычитаемое – это уже знакомая нам энтропия распределения $p(x)$, которая показывает, сколько в среднем бит требуется, чтобы закодировать значение случайной величины $\xi \sim p(x)$. Уменьшаемое носит название кросс-энтропии распределений $p(x)$ и $q(x)$. Кросс-энтропию можно интерпретировать как среднее число бит для кодирования значения случайной величины $\xi \sim p(x)$ алгоритмом, оптимизированным для кодирования случайной величины $\eta \sim q(x)$. Иными словами, дивергенция Кульбака-Лейблера говорит о том, насколько увеличится средняя длина кодов для значений $p$, если при настройке алгоритма кодирования вместо $p$ использовать $q$. Подробнее об этом вы можете почитать, например, в данном посте.
Дивергенция Кульбака-Лейблера в некотором роде играет роль расстояния между распределениями. В частности, $\mathbb{KL}(p\vert\vert q)\geqslant 0$, причём дивергенция равна нулю, только если распределения совпадают почти всюду (для дискретных и непрерывных распределений это означает, что они просто тождественны). Но при этом она не является симметричной: вообще говоря, $\mathbb{KL}(p\vert\vert q)\ne \mathbb{KL}(q\vert\vert p)$.
Упражнение. Пользуясь неравенством $\ln (1+t) \leqslant t$, $ t > -1$, докажите неотрицательность KL-дивергенции.
Попробуйте доказать самостоятельно, прежде чем смотреть решение
$$ \mathbb{KL}(p\vert\vert q) = -\int p(x)\ln{\frac{q(x)}{p(x)}}dx \geqslant $$ $$ \geqslant \int p(x)\Big(\frac{q(x)}{p(x)} - 1\Big)dx = \int q(x)dx - \int p(x)dx = 0. $$
Пример. С помощью KL-дивергенции измерим расстояние между двумя гауссианами $p(x) = \mathcal N(x \vert \mu_1, \sigma_1^2)$ и $q(x) = \mathcal N(x \vert \mu_2, \sigma_2^2)$.
Подставляя явные выражения для плотностей
$$ p(x) = \frac 1{\sqrt{2\pi}\sigma_1} e^{-\frac{(x-\mu_1)^2}{2\sigma_1^2}}\text{ и } q(x) = \frac 1{\sqrt{2\pi}\sigma_2} e^{-\frac{(x-\mu_2)^2}{2\sigma_2^2}}, $$
находим
$$ \ln \frac{p(x)}{q(x)} = \ln \frac{\sigma_2}{\sigma_1} + \frac{(x-\mu_2)^2}{2\sigma_2^2} -\frac{(x-\mu_1)^2}{2\sigma_1^2} = $$ $$ =\ln \frac{\sigma_2}{\sigma_1} + \Big(\frac 1{2\sigma_2^2} - \frac 1{2\sigma_1^2}\Big)(x-\mu_1)^2 + \frac{\mu_1 - \mu_2}{\sigma_2^2}(x-\mu_1) + \frac{(\mu_1 - \mu_2)^2}{2\sigma_2^2}. $$
Из свойств нормального распределения вытекает, что
$$ \int\limits_{-\infty}^{+\infty}p(x) (x-\mu_1),dx = 0, \quad \int\limits_{-\infty}^{+\infty}p(x) (x-\mu_1)^2,dx = \sigma_1^2. $$
Таким образом,
$$ \mathbb{KL}(p\vert\vert q) = \mathbb E_{\xi \sim p(x)} \ln \frac{p(\xi)}{q(\xi)} = \ln \frac{\sigma_2}{\sigma_1} + \frac{\sigma_1^2 + (\mu_1 - \mu_2)^2}{2\sigma_2^2} - \frac 12. $$
Как и должно быть, полученное выражение равно нулю, если гауссианы совпадают. При равных дисперсиях $\sigma_1^2 = \sigma_2^2 = \sigma^2$ получаем, что $\mathbb{KL}(p\vert\vert q) = \frac{(\mu_1 - \mu_2)^2}{2\sigma^2}$. Это выражение остаётся прежним, если поменять местами $\mu_1$ и $\mu_2$, поэтому в этом случае $\mathbb{KL}(p\vert\vert q) = \mathbb{KL}(q\vert\vert p)$. Если же $\sigma_1 \ne \sigma_2$, то выражение для $\mathbb{KL}(q\vert\vert p)$ явно отличается от $\mathbb{KL}(p\vert\vert q)$, что лишний раз показывает несимметричность KL-дивергенции.
Упражнение. Найдите дивергенцию Кульбака-Лейблера двух показательных распределений $p(x) = \mathrm{Exp}(x \vert \lambda)$ и $q(x) = \mathrm{Exp}(x \vert \mu)$.
Попробуйте решить самостоятельно, прежде чем смотреть решение
$$ \ln \frac{p(x)}{q(x)} = \ln \frac \lambda\mu + (\mu-\lambda)x. $$
Следовательно,
$$ \mathbb{KL}(p\vert\vert q) = \ln \frac \lambda\mu + (\mu-\lambda)\int\limits_0^{+\infty} \lambda x e^{-\lambda x},dx = $$ $$ =\ln\frac \lambda\mu + \frac{\mu-\lambda}{\lambda}\int\limits_0^{+\infty} te^{-t}dt = \ln \frac \lambda\mu + \frac \mu \lambda - 1. $$
И здесь KL-дивергенция равна нулю при $\lambda = \mu$.
Кросс-энтропия
При определении KL-дивергенции мы уже встречались с кросс-энтропией
$$ \mathbb H[p, q] = -\mathbb E_{\xi \sim p(x)} \log q(\xi) $$
В зависимости от типа распределений кросс-энтропия вычисляется по формуле
$$ \mathbb H[p, q] = -\sum\limits_k p_k \log q_k \text{ или } \mathbb H[p, q] = -\int\limits_{-\infty}^{+\infty} p(x) \log q(x), dx. $$
Поскольку
$$ \mathbb{KL}(p\vert\vert q) = \mathbb H[p, q] - \mathbb H[p], $$
от $q(x)$ .
В машинном обучении кросс-энтропию часто используют в качестве функции потерь в задаче классификации на $K>1$ классов. Истинное распределение на каждом обучающем объекте задаётся с помощью one hot encoding и является вырожденным:
$$ \boldsymbol y = (y_1, \ldots, y_K), \quad y_k \in {0,1}, \quad \sum\limits_{k=1}^K y_k = 1. $$
Классификатор обычно выдаёт вероятности принадлежности каждому из классов,
$$ \boldsymbol{\widehat y} = (\widehat y_1, \ldots, \widehat y_K), \quad \widehat y_k = \mathbb P(\text{класс }k). $$
Функция потерь на одном объекте полагается равной кросс-энтропии между истинным и предсказанным распределениями:
$$ \mathcal L(\boldsymbol y, \boldsymbol {\widehat y}) = -\sum\limits_{k=1}^K y_k \log \widehat y_k. $$
И это вполне логично: чем ближе модельное распределение к истинному, тем меньше наши потери. В идеале $\mathcal L(\boldsymbol y, \boldsymbol {\widehat y}) = 0$, если $\boldsymbol y = \boldsymbol {\widehat y}$.
Чтобы вычислить функцию потерь по обучающей выборке из $N$ объектов с метками $\boldsymbol y^{(i)}$, обычно берут усреднённную кросс-энтропию
$$ \frac 1N \sum\limits_{i=1}^N\mathcal L(\boldsymbol y^{(i)}, \boldsymbol{\widehat y}^{(i)}) =-\frac 1N\sum\limits_{i=1}^N\sum\limits_{k=1}^K y^{(i)}_k \log \widehat y^{(i)}_k. $$
Принцип максимальной энтропии
А теперь представим, что мы посчитали эти (или какие-то другие) статистики, а семейство распределений пока не выбрали. Как же совершить этот судьбоносный выбор? Давайте посмотрим на следующие три семейства и подумаем, в каком из них мы бы стали искать распределение, зная его истинные матожидание и дисперсию?
Почему-то хочется сказать, что в первом. Почему? Второе не симметрично – но что нас может заставить подозревать, что интересующее нас распределение не симметрично? С третьим проблема в том, что, выбирая его, мы добавляем дополнительную информацию как минимум о том, что у распределения конечный носитель. А с чего бы? У нас такой инфомации вроде бы нет.
Общая идея такова: мы будем искать распределение, которое удовлетворяет только явно заданным нами ограничениям и не отражает никакого дополнительного знания о нём. Таким образом, искомое распределение должно обладать максимальной неопределённостью при заданных ограничениях, или, говоря более научно, иметь максимально возможную энтропию. В самом деле, энтропия выражает нашу меру незнания о том, как ведёт себя распределение, и чем она больше – тем более «произвольное» распределение, по крайней мере, в теории. В этом и заключается принцип максимальной энтропии для выбора модели машинного обучения.
Как мы уже видели выше, среди распределений с конечным носителем максимальную энтропию имеет равномерное распределение. Примеры геометрического и нормального распределения показывают, что энтропия распределений с бесконечным носителем (счётным или континуальным) может быть сколь угодно большой, и среди них нет какого-то одного распределения с максимальной энтропией. Однако в более узком классе распределений с фиксированным средним и/или дисперсией найти распределение с максимальной энтропией, как правило, можно.
Пример. Покажем, что среди распределений на множестве натуральных чисел $\mathbb N$ и математическим ожиданием $\mu > 1$ максимальную энтропию имеет геометрическое распределение.
Для минимизации энтропии $\mathbb H[p] = -\sum\limits_{n=1}^\infty p_n \log p_n$ с учётом ограничений
$$ \sum\limits_{n=1}^\infty p_n = 1, \quad \sum\limits_{n=1}^\infty np_n = \mu $$
воспользуемся методом множителей Лагранжа, согласно которому требуется минимизировать функцию Лагранжа
$$ \mathcal L(p, a, b) = \sum\limits_{n=1}^\infty p_n \log p_n - a\Big(\sum\limits_{n=1}^\infty p_n - 1\Big) - b\Big(\sum\limits_{n=1}^\infty np_n - \mu\Big). $$
Приравняем к нулю частные производные по $p_n$:
$$ \frac{\partial \mathcal L}{\partial p_n} = 1 + \log p_n - a - bn = 0. $$
Отсюда следует, что $p_n = e^{a-1 + bn} = \alpha \beta^n$, так что распределение действительно получается геометрическое. Параметры $\alpha$ и $\beta$ найдём из уравнений
$$ 1 = \sum\limits_{n=1}^\infty p_n = \sum\limits_{n=1}^\infty \alpha \beta^n = \frac{\alpha\beta}{1-\beta}, $$
$$ \mu = \sum\limits_{n=1}^\infty np_n = \alpha \beta\sum\limits_{n=1}^\infty n \beta^{n-1} = \frac{\alpha\beta}{(1-\beta)^2}. $$
Деля первое уравнение на второе, получаем $\frac 1 \mu = 1 - \beta$, или $\beta = 1 - \frac 1\mu$. Далее из первого уравнения находим $\alpha = \frac {1-\beta}{\beta} = \frac 1{\mu -1}$. Итак,
$$ p_n = \frac 1{\mu -1} \Big(1 - \frac 1\mu\Big)^n = \frac 1\mu \Big(1 - \frac 1\mu\Big)^{n-1}, $$
а это и есть геометрическое распределение с параметром $\frac 1\mu$.
У непрерывных распределений возможны более интересные комбинации из ограничений на носитель и параметры. И конечно же, первую скрипку среди распределений с максимальной энтропией играет гауссовское распределение.
Пример. Докажем, что среди распределений на $\mathbb R$ c математическим ожиданием $\mu$ и дисперсией $\sigma^2$ наибольшую энтропию имеет нормальное распределение $\mathcal{N}(\mu,\sigma^2)$.
Пусть $p(x)$ – некоторое распределение со средним $\mu$ и дисперсией $\sigma^2$, $q(x)\sim\mathcal{N}(\mu, \sigma^2)$. Как было показано выше, $\mathbb H[q] = \frac12\log(2\pi\sigma^2) + \frac12$. Запишем дивергенцию Кульбака-Лейблера:
$$ \mathbb KL(p\vert\vert q) = \int\limits_{-\infty}^{+\infty} p(x)\log{p(x)}dx - \int\limits_{-\infty}^{+\infty} p(x)\log{q(x)}dx = $$ $$ = -\mathbb H[p] - \int\limits_{-\infty}^{+\infty} p(x)\left(-\frac12\log(2\pi\sigma^2) - \frac 1{2\sigma^2}(x - \mu)^2\right)dx = $$ $$ = - \mathbb H[p] +\frac12\log(2\pi\sigma^2)\int\limits_{-\infty}^{+\infty} p(x),dx + \frac1{2\sigma^2}\underbrace{\int\limits_{-\infty}^{+\infty}(x - \mu)^2p(x),dx}_{=\mathbb{V}[p]=\sigma^2} = $$ $$ = - \mathbb H[p] + \frac12\log(2\pi\sigma^2) + \frac12 = \mathbb H[q] - \mathbb H[p]. $$
Так как KL-дивергенция всегда неотрицательна, получаем, что $\mathbb H[p]\leqslant \mathbb H[q]$ при любом распределении $p$, удовлетворяющем заданным ограничениям.
Можно показать, что максимальную энтропию среди многомерных распределений с вектором средних $\boldsymbol \mu$ и матрицей ковариаций $\boldsymbol \Sigma$ имеет также гауссовское распределение $\mathcal N(\boldsymbol \mu, \boldsymbol \Sigma)$.
Упражнение. Докажите, что среди распределений на отрезке $[a,b]$ максимальную энтропию имеет равномерное распределение $U[a,b]$.
Попробуйте доказать самостоятельно, прежде чем смотреть решение.
$$ \mathbb H[p,q] = - \int\limits_a^b p(x)\log q(x),dx = \int\limits_a^b p(x) \log(b-a),dx = \log(b-a). $$
В частности, отсюда следует, что $\mathbb H[q] = \log(b-a)$. Расписывая теперь KL-дивергенцию, получаем
$$ \mathbb{KL}(p\vert\vert q) = \mathbb H[p,q] - \mathbb H[p] = \log(b-a) - \mathbb H[p] = \mathbb H[q] - \mathbb H[p] \geqslant 0, $$
что и требовалось доказать.
Упражнение. Докажите, что среди распределений на промежутке $[0, +\infty)$ с математическим ожиданием $\lambda > 0$ максимальную энтропию имеет показательное распределение $\mathrm{Exp\big(\frac 1\lambda)}$.
Попробуйте доказать самостоятельно, прежде чем смотреть решение.
$$ \mathbb{KL}(p\vert\vert q) = \mathbb H[p,q] - \mathbb H[p] = - \int\limits_0^{+\infty} p(x)\Big(-\log\lambda - \frac x\lambda\Big) ,dx - \mathbb H[p]= $$ $$ =\log\lambda + \frac 1\lambda \int\limits_0^{+\infty} xp(x), dx - \mathbb H[p] = \log\lambda + 1 - \mathbb H[p] = \mathbb H[q]- \mathbb H[p] \geqslant 0. $$
Как выяснилось, многие классические распределения имеют максимальную энтропию при весьма естественных ограничениях. Но как быть, если даны не эти конкретные, а какие-то другие ограничения? Есть ли какой-нибудь надёжный алгоритм вывода распределения с максимальной энтропией, позволяющий избежать случайных озарений и гаданий на кофейной гуще? Оказывается, что при некоторых не очень обременительных ограничениях ответ можно записать с помощью распределений экспоненциального класса.
Экспоненциальное семейство распределений
Говорят, что параметрическое семейство распределений относится к экспоненциальному классу, если его pdf (или pmf) может быть представлена в виде
$$ p(x\vert \boldsymbol\theta) = \frac{g(x)}{h(\boldsymbol\theta)}\exp(\boldsymbol \theta^T \boldsymbol u(x)) = g(x)\exp(\boldsymbol \theta^T \boldsymbol u(x) - A(\boldsymbol \theta)), $$
где
- $\boldsymbol\theta \in\mathbb R^m$ – вектор натуральных параметров распределения;
- $g(x)$ — неотрицательная функция (base measure), часто равная единице;
- $h(\boldsymbol\theta) > 0$ — нормализатор (partition), обеспечивающий суммируемость pmf или интегрируемость pdf в единицу:
$$ h(\boldsymbol\theta) = \int g(x)\exp(\boldsymbol \theta^T \boldsymbol u(x))dx; $$
- $A(\boldsymbol \theta) = \ln h(\boldsymbol\theta)$ — log-partition;
- $\boldsymbol u(x)\in\mathbb R^m$ — вектор достаточных статистик распределения.
Пример. Покажем, что нормальное распределение $\mathcal N(x\vert \mu, \sigma^2)$ принадлежит экспоненциальному классу. Оно имеет два параметра, поэтому такую же размерность имеют $\boldsymbol \theta$ и вектор-функция $\boldsymbol u$.
Распишем плотность:
$$ \frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right) = \frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac 1{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)}. $$
Положим $g(x) = \frac 1{\sqrt{2\pi}}$, $\boldsymbol u(x) = (x, x^2)$,
$$\boldsymbol \theta = (\theta_1, \theta_2),\quad \theta_1 = \frac{\mu}{\sigma^2},\quad \theta_2 = -\frac1{2\sigma^2}.$$
Остаётся выразить функцию $h(\boldsymbol\theta) = \sigma\exp\left(-\frac{\mu^2}{2\sigma^2}\right)$ через $\theta_1$ и $\theta_2$.
Упражнение. Выразите partition $h(\boldsymbol\theta)$ и log-partition $A(\boldsymbol\theta)$ через $\boldsymbol\theta$ и запишите плотность нормального распределения в экспоненциальном виде.
Ответ
$$\mathcal N(x \vert \mu, \sigma^2) = \frac1{\sqrt {2\pi}}\cdot \frac 1{h(\boldsymbol\theta)}\exp(\boldsymbol\theta^T \boldsymbol u(x)) = \frac1{\sqrt {2\pi}} \exp(\boldsymbol \theta^T \boldsymbol u(x) - A(\boldsymbol \theta)). $$
Пример. Покажем, что распределение Бернулли $\mathrm{Bern}(p)$ принадлежит экспоненциальному классу. Его pmf $\mathbb P(\xi =x \vert p)$ можно записать как
$$ p^x(1 - p)^{1 - x} = \exp\big(x\log{p} + (1 - x)\log(1 - p)\big)= (1-p) \exp\Big(x\log\frac p{1-p}\Big). $$
Параметр здесь один, поэтому натуральный параметр $\theta$ тоже один: $\theta = \log\frac p{1-p}$. Такая функция от $p$ называется функцией логитов и активно участвует в построении модели логистической регрессии. Остальные функции положим равными $u(x) = x$, $g(x)=1$, $h(\theta) = \frac 1{1-p}$. Остаётся выразить partition через $\theta$:
$$ \log\frac p{1-p} = \log\Big(-1 + \frac 1{1-p}\Big) = \theta \iff \frac 1{1-p} = 1 + e^\theta. $$
Итак, $h(\theta) = 1+e^\theta$, и экспоненциальный вид распределения Бернулли записывается как
$$ \frac 1{1+e^\theta}\exp(\theta x) = \exp\big(\theta x - \log(1 + e^\theta)\big). $$
Вопрос на подумать. Принадлежит ли к экспоненциальному классу семейство равномерных распределений на отрезках $U[a, b]$? Казалось бы, да, ведь
$$p(x) = \frac{1}{b - a}\mathbb{I}_{[a,b]}(x)\exp(0).$$
В чём может быть подвох?
Попробуйте определить сами, прежде чем смотреть ответ.
При этом странное и не очень полезное семейство с нулём параметров, состоящее из одинокого распределения $U[0,1]$ можно считать относящимся к экспоненциальному классу: ведь для него формула
$$p(x) = \mathbb{I}_{[0,1]}(x)\exp(0)$$
будет работать.
К экспоненциальным семействам относятся многие непрерывные и дискретные распределения из часто встречающихся в теории и на практике, в том числе
- нормальное $\mathcal N(\mu, \sigma^2)$;
- распределение Пуассона $\mathrm{Pois}(\lambda)$;
- экспоненциальное $\mathrm{Exp}(\lambda)$;
- биномиальное $\mathrm{Bin}(n, p)$;
- геометрическое $\mathrm{Geom}(p)$;
- бета-распределение;
- гамма-распределение;
- распределение Дирихле.
Как выглядят натуральные параметры, достаточные статистики и нормализаторы этих и других распределений из экспоненциального класса, можно посмотреть на википедии.
К экспоненциальным семействам не относятся, например, равномерное распределение $U[a,b]$, $t$-распределение Стьюдента, распределение Коши, смесь нормальных распределений.
Дифференцирование log-partition
Если распределение $p(x\vert \boldsymbol\theta)$ принадлежит экспоненциальному классу,
$$ p(x\vert \boldsymbol\theta) = \frac{g(x)}{h(\boldsymbol\theta)}\exp(\boldsymbol \theta^T \boldsymbol u(x)) = g(x)\exp(\boldsymbol \theta^T \boldsymbol u(x) - A(\boldsymbol \theta)), $$
то моменты его достаточных статистик $\boldsymbol u(x)$ могут быть получены дифференцированием функции $A(\boldsymbol \theta) = \log h(\boldsymbol\theta)$.
Утверждение. $\nabla_{\boldsymbol\theta}\log h(\boldsymbol \theta) = \mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi)$.
Доказательство. По правилу дифференцирования сложной функции имеем
$$ \nabla_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} = \frac{\nabla_{\boldsymbol\theta}{h(\boldsymbol\theta)}}{h(\boldsymbol\theta)}. $$
Нормализатор $h(\boldsymbol\theta)$ записывается в виде интеграла
$$ h(\boldsymbol\theta) = \int g(x)\exp\left(\boldsymbol\theta^T\boldsymbol u(x)\right)dx, $$
который мы продифференцируем внесением градиента внутрь под знак интеграла:
$$
\nabla_{\boldsymbol\theta} h(\boldsymbol\theta) = \nabla_{\boldsymbol\theta}
\int g(x)\exp\left(\boldsymbol\theta^T\boldsymbol u(x)\right)dx =
$$ $$
=\int g(x)\nabla_{\boldsymbol\theta}
\exp\left(\boldsymbol\theta^T\boldsymbol u(x)\right)dx =
\int g(x)\boldsymbol u(x)
\exp\big(\boldsymbol\theta^T\boldsymbol u(x)\big)dx.
$$
Таким образом,
Замечание о дифференцировании под знаком интегралом
Использованный в доказательстве приём внесения градиента под знак интеграла называют также правилом Лейбница. Этот же метод используется для почленного дифференцирования ряда, что может быть полезно в случае дискретного распределения $p(x\vert \boldsymbol\theta)$. В математическом анализе имеется ряд теорем, обеспечивающих применимость правила Лейбница, однако, мы не будем на них останавливаться. Будем считать, что все рассматриваемые экспоненциальные семейства таковы, что применение правила Лейбница законно.
Если $u_i(x) = x^i$, то в соответствии с только что доказанным частная производная $\frac{\partial A(\boldsymbol \theta)}{\partial \theta_i}$ даёт $i$-й момент распределения $p(x\vert \boldsymbol\theta)$.
Упражнение. Вычислите производные по натуральным параметрам от log-partition для распределения Бернулли $\mathrm{Bern}(x\vert p)$ и нормального распределения $\mathcal N(\mu, \sigma^2)$ и проверьте, что они совпадают со значениями соответствующих моментов.
Решение
Для $\xi \sim \mathrm{Bern}(x\vert p)$ имеем
$$ A(\theta) = \log(1 + e^\theta), \quad \frac 1{1-p} = 1 + e^\theta, $$
поэтому
$$ A'(\theta) = \frac{e^\theta}{1+e^\theta} = 1 - \frac 1{1+e^\theta} = 1 -1 + p = p = \mathbb E\xi. $$
Для гауссовской случайной величины $\eta\sim\mathcal N(x \vert\mu, \sigma^2)$ получаем
$$ A(\theta_1, \theta_2) = -\frac{\theta_1^2}{4\theta_2} - \frac 12\ln(-2\theta_2), $$ $$ \quad \theta_1 = \frac{\mu}{\sigma^2}, $$ $$ \quad \theta_2 = -\frac1{2\sigma^2}, $$
и, значит,
$$ \frac{\partial A}{\partial\theta_1} = -\frac{\theta_1}{2\theta_2} = \mu = \mathbb E\eta, $$ $$ \frac{\partial A}{\partial\theta_2} = \frac{\theta_1^2}{4\theta_2^2} - \frac 1{2\theta_2} = \mu^2 + \sigma^2 = \mathbb E\eta^2. $$
Кстати, можно продифференцировать ещё раз и доказать, что
$$ \nabla^2_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} = \mathrm{cov}(\boldsymbol u(\xi), \boldsymbol u(\xi)). $$
Доказательство
В предыдущий раз мы доказали, что
$$ \nabla_{\boldsymbol\theta}A(\boldsymbol\theta) = \int \boldsymbol u(x)g(x)\exp(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta)),dx. $$
Теперь возьмём ещё раз градиент по $\boldsymbol \theta$ от этого выражения:
$$ \nabla^2_{\boldsymbol\theta}A(\boldsymbol\theta) = \int \boldsymbol u(x)g(x)\nabla_{\boldsymbol\theta}\exp\big(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta)\big),dx = $$ $$ = \int \boldsymbol u(x)g(x)\big(\boldsymbol u(x) - \nabla_{\boldsymbol\theta}A(\boldsymbol \theta)\big)^T \exp\big(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta)\big),dx= $$ $$ =\int \boldsymbol u(x)g(x)\big(\boldsymbol u(x) - \mathbb{E} \boldsymbol u\big)^T \exp\big(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta)\big),dx = $$ $$ = \mathbb{E} \boldsymbol u \boldsymbol u^T - \mathbb{E} \boldsymbol u \big(\mathbb{E} \boldsymbol u\big)^T = \mathrm{cov}(\boldsymbol u(\xi), \boldsymbol u(\xi)). $$
В последней выкладке для краткости обозначено $\mathbb E \boldsymbol u = \mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi)$.
MLE для семейства из экспоненциального класса
Возможно, вас удивил странный и на первый взгляд не очень естественный вид $p(x\vert \boldsymbol \theta)$. Но всё не просто так: оказывается, что оценка максимального правдоподобия параметров распределений из экспоненциального класса устроена очень интригующе.
Запишем функцию правдоподобия i.i.d. выборки $x_1,\ldots,x_n$:
$$ p(x_1, \ldots, x_n\vert\boldsymbol \theta) = h(\boldsymbol \theta)^{-n}\cdot\bigg(\prod_{i=1}^ng(x_i)\bigg)\cdot\exp\bigg(\boldsymbol \theta^T \sum_{i=1}^n \boldsymbol u(x_i)\bigg). $$
Её логарифм равен
$$ \log p(x_1, \ldots, x_n\vert\boldsymbol \theta) = -n\log{h(\boldsymbol \theta)} + \sum_{i=1}^n\log{g(x_i)} + \boldsymbol \theta^T\sum_{i=1}^n \boldsymbol u(x_i). $$
Дифференцируя по $\boldsymbol \theta$, получаем
$$\nabla_{\boldsymbol\theta}\log (x_1, \ldots, x_n\vert\boldsymbol\theta) = -n\nabla_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} + \sum_{i=1}^n \boldsymbol u(x_i).$$
Приравнивая $\nabla_{\boldsymbol\theta}\log p(x_1,\ldots, x_n\vert\boldsymbol\theta)$ к нулю и пользуясь равенством $\nabla_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} = \mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi)$, находим
Таким образом, теоретические матожидания всех компонент $u_i(\xi)$ должны совпадать с их эмпирическими оценками, а метод максимального правдоподобия совпадает с методом моментов для $\mathbb{E}u_i(\xi)$ в качестве моментов. И в следующем пункте выяснится, что распределения из экспоненциальных семейств обладают максимальной энтропией среди тех, что имеют заданные моменты $\mathbb{E}u_i(\xi)$.
Теорема Купмана—Питмана—Дармуа
Теперь мы наконец готовы сформулировать одно из самых любопытных свойств семейств экспоненциального класса.
В следующей теореме мы опустим некоторые не очень обременительные условия регулярности. Просто считайте, что для хороших дискретных и абсолютно непрерывных распределений, с которыми вы в основном и будете сталкиваться, это так.
Теорема. Пусть параметр $\boldsymbol\theta\in R^m$ распределения $p_{\boldsymbol\theta}(x) = \frac1{h(\boldsymbol\theta)}\exp\big(\boldsymbol\theta^T\boldsymbol u(x)\big)$ выбран так, что
для некоторого фиксированного $\boldsymbol \alpha \in\mathbb R^m$. Тогда распределение $p_{\boldsymbol\theta}(x)$ обладает наибольшей энтропией среди распределений $q$ с тем же носителем, для которых $\mathbb{E}_{\xi \sim q(x)} \boldsymbol u(\xi) = \boldsymbol\alpha$.
Идея обоснования через оптимизацию.
В дискретном случае у нас есть счётное семейство точек $x_1, x_2,\ldots$, и распределение определяется счётным набором вероятностей $p_i$ принимать значение $x_i$. Мы будем решать задачу
Для решения этой оптимизационной задачи нам понадобится обобщение метода множителей Лагранжа, известное также как теорема Каруша—Куна—Таккера. В данном случае задача сводится к минимизации лагранжиана
$$\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$$
при имеющихся ограничениях и условиях дополняющей нежёсткости $\lambda_jp_j = 0$. Приравняем частные производные по $p_j$ к нулю:
$$\frac{\partial\mathcal{L}}{\partial p_j} = \log{p_j} + 1 - \sum_i\theta_iu_i(x_j) + \theta_0 - \lambda_j = 0.$$
Отсюда получаем, что
$$p_j = \frac{\exp(\boldsymbol \theta^T \boldsymbol u(x_j))}{\exp\left(\lambda_j - \theta_0 - 1\right)}.$$
Числитель уже ровно такой, как и должен быть у распределения из экспоненциального класса; разберёмся со знаменателем. Поскольку $p_j >0$ (ведь тут сплошные экспоненты), $\lambda_j = 0$ при всех $j$. Параметр $\theta_0$ находится из условия $\sum_jp_j = 1$, а точнее, выражается через остальные $\theta_i$, что позволяет записать знаменатель в виде $h(\boldsymbol \theta)$.
Идея доказательства «в лоб».
$$\int u_i(x)q(x)dx = \int u_i(x)p_{\boldsymbol\theta}(x)dx$$
для всех $i = 1,\ldots,m$. Тогда
$$0\leqslant \mathbb{KL}(q\vert\vert p) = \int q(x)\log\left(\frac{q(x)}{p_{\boldsymbol\theta}(x)}\right)dx = $$
$$=-\mathbb H[q] - \int q(x)\left(-\log{h(\boldsymbol \theta)} + \sum_i\theta_iu_i(x)\right)dx =$$
$$=-\mathbb H[q] - \int p(x)\left(-\log{h(\boldsymbol\theta)} + \sum_i\theta_iu_i(x)\right)dx =$$
$$=-\mathbb H[q] + \int p_{\boldsymbol\theta}(x)\log{p_{\boldsymbol\theta}(x)}dx = -\mathbb H[q] + \mathbb H[p_{\boldsymbol\theta}].$$
Таким образом, $\mathbb H[p_{\boldsymbol\theta}]\geqslant \mathbb H[q]$, что и требовалось доказать.
Выше мы уже находили обладающее максимальной энтропией распределение на множестве натуральных чисел с заданным математическим ожиданием $\mu>1$. Таковым оказалось геометрическое распределение $\mathrm{Geom}\big(\frac 1\mu\big)$.
Теорема Купмана—Питмана—Дармуа позволяет сделать это гораздо быстрее. В данном случае у нас лишь одна функция $u_1(x) = x$, которая соответствует фиксации математического ожидания $\mathbb E\xi$. Искомое дискретное распределение имеет вид
$$ p_k =\mathbb P(\xi = k) = \frac1{h(\theta)}\exp(\theta k) = \frac1{h(\theta)} \big(e^{\theta}\big)^k. $$
Это уже похоже на геометрическое распределение с параметром $p = 1 - e^{\theta}$. Его математическое ожидание равно $\frac 1p$, что по условию должно равняться $\mu$. Итак, наше распределение с максимальной этропией выглядит так:
$$ p_k = \frac1{\mu}\left(1 - \frac1{\mu}\right)^{k-1},\quad k\in\mathbb N. $$
Пример. Среди распределений на всей вещественной прямой с заданным математическим ожиданием $\mu$ найдём распределение с максимальной энтропией.
А сможете ли вы его найти? Решение под катом.
$$p(x) = \frac1{h(\theta)}\exp\left(\theta x\right),$$
но интеграла экспоненты не существует, то есть применение «в лоб» теоремы провалилось. И неспроста: если даже рассмотреть все нормально распределённые случайные величины со средним $\mu$, их энтропии, равные $\frac12 + \frac12\log(2\pi\sigma^2)$, не ограничены сверху, то есть величины с наибольшей энтропией не существует.