Зачем нужны модели с латентными переменными
Предположим, что мы делаем анализ данных для банка, и нам предоставили данные о годовых зарплатах клиентов.
В этом графике заметны три моды, которым, наверное, соответствуют три кластера клиентов. Неопытный аналитик мог бы проигнорировать это и попытаться описать график в отчете для руководства с помощью двух чисел — средней зарплаты и стандартного отклонения зарплат.
Однако данные с гистограммы ниже имеют точно такое же среднее и стандартное отклонение, как и мультимодальные данные выше. Распределение совсем другое, правда?
Очевидно, что графики выглядят совершенно по-разному, и правильная интерпретация первого графика может принести бизнесу дополнительные деньги (скажем, если банк научится предлагать клиентам из каждого кластера более кастомизированные предложения). Так что не надо пытаться описывать мультимодальные данные с помощью унимодальных распределений.
Проблема заключается в том, что нам неизвестно, к какому кластеру относится каждый клиент, и неизвестны характеристики кластеров — как же их тогда описать?
Для каждого кластера можно попытаться задать свои параметры (среднее и дисперсию). Но как определить, из какого кластера конкретный клиент в выборке? Более того, один клиент может, например, «на 0.7» относиться к одному кластеру и «на 0.3» к другому.
Как решать такую задачу «мягкой» кластеризации («мягкой», потому что один объект может относиться к нескольким кластерам)? Мы могли бы действовать итерационно. Сначала зададим начальное приближение на параметры распределений.
Например, в нашем случае с клиентами банка из графика можно предположить, что среднее для первой кластера 30, для второго — 40, для третьей — 50, а стандартные отклонения у всех равняются, скажем, 10. Зная эти начальные параметры, мы можем для каждого клиента посчитать степень принадлежности к каждому из трёх кластеров (важно не забыть отнормировать эти числа, чтобы их сумма действительно равнялась единице).
Дальше мы бы могли пересчитать наши средние и дисперсии, «взвешивая» вклад объектов пропорционально степени их принадлежности к каждому кластеру, и таким образом уточнить средние и дисперсии для всех трёх кластеров. Повторяя эти два шага последовательно, мы получили бы средние и дисперсии кластеров, а для каждого объекта — степени принадлежности к кластерам. Это — EM-алгоритм, подробнее о нём мы поговорим ниже.
Кстати, оказывается, что и метод кластеризации K-средних во многом сродни EM-алгоритму (и на самом деле представляет собой его предельный случай). Действительно, мы сначала случайно расставляем центры кластеров.
Затем мы для каждого объекта пересчитываем расстояние до центра каждого кластера, после чего получаем «вес» объекта в каждом кластере («вес» в том смысле, что чем ближе объект к центру кластера, тем больше этот объект учитывается при пересчете центра этого кластера) через нормировку расстояний.
Теперь, если мы применяем настоящий метод K-средних, то приписываем объект к кластеру с самым большим «весом», после чего опять пересчитываем центры кластеров и потом опять измеряем вес для каждого объекта кластера. Строгое же применение EM-алгоритма даёт «мягкую» версию метода K-средних.
Смеси распределений
Говорят, что распределение является смесью распределений, если его плотность имеет вид
где:
- — число компонент;
- — вероятности компонент;
- — функции правдоподобия, то есть функции вероятности компонент (в дискретном случае) или их плотности (в абсолютно непрерывном случае).
Проиллюстрируем это понятие на примере с банком. Будем считать, что распределения компонент смеси принадлежат некоторому параметрическому семейству: (например, гауссовскому с параметром ).
Мы можем говорить, что каждое из распределений задаёт свой кластер, причём каждый кластер имеет некоторую априорную вероятность . Если у нас нет дополнительных данных, разумно положить .
Если же нам, к примеру, известно, что какой-нибудь кластер описывает сравнительно малочисленную группу людей, эти вероятности окажутся различными. Таким образом, мы проинтерпретировали нашу мягкую кластеризацию в терминах смеси распределений.
Как генерировать из смеси распределений
Рассмотрим следующий эксперимент: сначала из дискретного распределения выбирается номер , а затем из распределения выбирается значение . Покажем, что распределение переменной будет представлять собой смесь вида [ссылка на определение смеси]
Введем скрытую переменную , отвечающую за то, к какой компоненте смеси будет относиться очередной . Пусть она представляет собой -мерный бинарный случайный вектор, ровно одна компонента которого равна единице:
Вероятность того, что единице будет равна -я компонента, положим равной :
Запишем распределение сразу всего вектора:
Теперь, когда номер компоненты смеси известен, сгенерируем из распределения :
или, что то же самое,
Проверим, что имеет нужное нам распределение. Запишем совместное распределение переменных и :
Чтобы найти распределение переменной , нужно избавиться от скрытой переменной:
Суммирование здесь ведется по всем возможным значениям , то есть по всем -мерным бинарным векторам с одной единицей:
Мы получили, что распределение сгенерированной переменной в описанном эксперименте представляет собой смесь компонент.
Модели со скрытыми переменными
Рассмотрим вероятностную модель с наблюдаемыми переменными и параметрами , для которой задано правдоподобие .
Предположим, что в модели также существуют скрытые переменные , описывающие её внутреннее состояние и, возможно, недоступные для непосредственного наблюдения (как то, к какому из кластеров относится клиент). Тогда правдоподобие называется неполным, а правдоподобие — полным. Они связаны соотношением
Нашей основной целью будет создать хорошую модель , то есть оценить параметры . И оказывается, что с помощью введения скрытых переменных нередко удаётся существенно упростить правдоподобие и эффективно решить задачу.
Рассмотрим пример со смесями распределений. В качестве наблюдаемых переменных здесь выступает выборка , в качестве скрытых переменных — номера компонент, из которых сгенерированы объекты (здесь каждый из — -мерный вектор), в качестве параметров — априорные вероятности и параметры компонент .
Неполное правдоподобие выглядит так:
Правдоподобие здесь имеет вид логарифма суммы. Если приравнять нулю его градиент, то получатся сложные уравнения, не имеющие аналитического решения. Данное правдоподобие сложно вычислять: оно не является выпуклым (а точнее, вогнутым) и может иметь много локальных экстремумов, поэтому применение обычных итерационных методов для его непосредственной максимизации приводит к медленной сходимости.
Рассмотрим теперь полное правдоподобие:
Оно имеет вид суммы логарифмов, и это позволяет аналитически найти оценки максимального правдоподобия на параметры при известных переменных и . В общем случае также стараются выбирать таким способом, чтобы распределение оказалось «лучше» исходного. В каком именно смысле, мы увидим дальше.
Проблема, впрочем, заключается в том, что нам не известны скрытые переменные , поэтому их необходимо оценивать одновременно с параметрами, что никак не легче максимизации неполного правдоподобия. Осуществить это позволяет EM-алгоритм.
EM-алгоритм
EM-алгоритм решает задачу максимизации полного правдоподобия путём попеременной оптимизации по параметрам и по скрытым переменным.
Опишем сначала наивный способ оптимизации. Зафиксируем некоторое начальное приближение для параметров . При известных наблюдаемых переменных и параметрах мы можем оценить скрытые переменные, найдя их наиболее правдоподобные значения:
Зная скрытые переменные, мы можем теперь найти следующее приближение для параметров:
Повторяя итерации до сходимости, мы получим некоторый итоговый вектор параметров .
Данная процедура, однако, далека от идеальной — и ниже мы предложим решение, которое приводит к более качественным результатам.
Воспользуемся байесовским подходом. Точечные оценки параметров несут меньше информации, чем их распределение; учтём это и будем оптимизировать не , а условное распределение .
Как и прежде, зафиксируем вектор параметров , но вместо точечной оценки вычислим апостериорное распределение на скрытых переменных , которое будет в некотором смысле оптимальным образом описывать распределение при известных и . В этом заключается E-шаг EM-алгоритма.
Отметим, что вычислить аналитически возможно не для всех распределений, и скрытые переменные стоит подбирать так, чтобы это всё-таки получилось.
Теперь мы должны произвести оптимизацию по . Для этого возьмём логарифм полного правдоподобия и усредним его по всем возможным значениям скрытых переменных :
Формально говоря, мы нашли матожидание логарифма полного правдоподобия по апостериорному распределению на скрытых переменных.
На M-шаге новый вектор параметров находится как максимизатор данного матожидания:
EM-алгоритм состоит в чередовании E-шага и M-шага.
Можно показать, что такой итерационной процесс всегда не уменьшает правдоподобие и сходится.
В двух словах: почему это работает
Напомним, что наша цель — найти параметры , максимизирующие . Преобразуем этот логарифм правдоподобия, введя некоторое (пока произвольное) распределение на латентных переменных:
Здесь мы применили то, что в силу равенства . Далее, домножим числитель и знаменатель под логарифмом на :
Второе слагаемое — это расстояние Кульбака-Лейблера между распределениями и . Так как расстояние Кульбака-Лейблера всегда неотрицательно, мы получаем, что . Первое слагаемое, , называется вариационной нижней оценкой на , и, максимизируя его, мы будем также увеличивать и всю сумму — поэтому в дальнейшем мы позволим себе ограничиться работой только с ним.
Вариационную нижнюю оценку мы будем попеременно оптимизировать по (да-да, по распределению) и по .
E-шаг состоит в максимизации по при фиксированных . Так как
а от не зависит, максимизировать по — это то же самое, что минимизировать . Расстояние Кульбака-Лейблера минимально, когда распределения совпадают, то есть мы должны обновить по правилу
На M-шаге мы максимизируем по при фиксированном . Преобразуем вариационную нижнюю оценку:
Заметим, что — это энтропия распределения , и она не зависит от . Таким образом, нам достаточно максимизировать по математическое ожидание
Жёсткий EM-алгоритм
Не всегда получается подобрать латентные переменные так, чтобы можно было выразить аналитически, то есть на E-шаге не удаётся минимизировать по .
В такой ситуации иногда приходится брать оптимум не по всему пространству распределений, а только по некоторому семейству — например, параметрическому, в котором оптимизацию можно проводить градиентными методами. В максимально упрощённой ситуации мы возьмём семейство дельта-функций, то есть вместо распределения на будем брать просто точечную оценку. Такая модификация EM-алгоритма называется жёстким EM-алгоритмом.
В начале параграфа мы упоминали кластеризацию методом K-средних и отмечали, что EM-алгоритм даёт «мягкую» версию алгоритма: на E-шаге мы не приписываем однозначно точку к какому-то из кластеров (то есть не берём точечную оценку скрытой переменной «номер кластера, к которому принадлежит точка»), а сопоставляем ей вероятности принадлежности каждому из кластеров (то есть распределение на скрытые переменные). Настоящий метод K-средних как раз таки соответствует жёсткому EM-алгоритму.
Разделение смеси гауссиан
Пусть теперь нам известно, что точек были семплированы из K разных гауссовских распределений и нам неизвестно, какая точка из какого распределения пришла в выборку. Нам нужно оценить параметры (, ) для первого распределения, (, ) для второго и, соответственно, (, ) для -го распределения.
Если мы знаем, что точка пришла из распределения , то её правдоподобие в равно:
Напомним, что — это латентная переменная, обозначающая номер гауссианы (от 1 до ), из которой была просемплирована точка . Например, если точка просемплирована из распределения с номером 3, то в формулу вместо (, ) нужно подставлять .
Воспользуемся EM-алгоритмом для нахождения параметров (, ). Сначала инициализируем их:
Зная , выполним E-шаг: нужно найти или что то же самое, для каждого объекта найти распределение на вероятности .
Как найти апостериорную вероятность , если мы знаем и у нас есть приближение ?
Ответ — по формуле Байеса:
где — априорная вероятность того, что объекта получен из распределения с номером . На первом шаге априорную вероятность можно положить равной для всех гауссиан.
Введём обозначение
— правдоподобие того, что объект пришел из нормального распределения с параметрами .
Тогда по формуле Байеса:
Вот так для каждого объекта по начальному приближению мы посчитаем распределение — с какими вероятностями объект приходит из той или иной компоненты смеси.
Теперь выведем формулы для М-шага.
Запишем производную и приравняем к , чтобы найти экстремум:
Отсюда
Мы получили конечную формулу для пересчета по и предыдущему значению . Причем у этой формулы есть простая интерпретация — каждый объект мы взвешиваем с его вероятностью принадлежности к этому классу .
Теперь посчитаем производную по (обратите внимание, что именно по квадрату ):
Стало быть,
Мы снова получили интерпретируемый результат: подсчитывая дисперсию для -ой гауссианы, мы учитываем вес каждого объекта при подсчете среднеквадратичноого отклонения. То есть веса — вероятности происхождения из той или иной компоненты смеси. Сравните эту формулу с формулой для подсчета выборочной дисперсии, где каждый из объектов вносит одинаковый вклад в дисперсию с весом :
Вы можете «пощупать» EM-алгоритм в задаче разделения вероятностной смеси с помощью интерактивной визуализации — попробуйте сделать E и M шаги и последить за изменениями параметров: после одной итераций алгоритма можно выбрать точку на графике и наблюдать за вероятностью её принадлежности к разным кластерам.
Вероятностный PCA
Теперь давайте рассмотрим простой пример того, как введение латентных переменных может помочь выделять новые информативные признаки в данных.
Предположим, что мы имеем выборку данных (вектор-строку), где каждый объект имеет признаков (предположим, что число очень большое). Это достаточно типичная ситуация, например, при работе с текстами или изображениями.
Теперь введём следующую вероятностную модель , где — латентный вектор-строка меньшей размерности , а , где — единичная матрица размером , — скаляр больший 0.
Что означает эта модель? Она означает, что наши сложные многоразмерные данные могут иметь более простое малоразмерное представление , а отображение линейно с точностью до нормально распределенного шума.
Заметим, что так как , отсюда следует, что . Зададим априорное распределение на как стандартное нормальное и распишем совместное распределение через условное и априорное:
Чтобы восстановить параметры , и латентные переменные , снова воспользуемся EM-алгоритмом.
На E-шаге мы оцениваем распределение на при фиксированных и :
По формуле Байеса распределение на при условии :
С точностью до констант и слагаемых, которые не зависят от , логарифм правдоподобия равен:
Обозначим , тогда
Если теперь взять от этого экспоненту, увидим, что .
M-шаг.
Теперь мы оптимизировать по и :
Приравняв производные к , можно найти:
Вероятностный PCA хорош тем, что:
-
как любая байесовская модель, может служить промежуточным участком в более сложной вероятностной модели;
-
если в данных есть пропуски, то вероятностный PCA легко обобщается и на этот случай с добавлением дополнительных скрытых переменных;
-
так как параметры и оценки на получаются через итерационный EM-алгоритм, то вероятностный PCA может быть вычислительно эффективнее. Так, в вычислениях и промежуточных формулах нигде не используется матрица , и все рассматриваемые матрицы имеют меньший размер.
Связь с обычным PCA
Как вероятностный PCA связан с обычным, который мы изучили в теме про разложение матриц?
Напомним, что в обычном SVD-разложении мы полагали, что . Давайте опять положим, что разница между и есть гауссовский шум с нулевым средним :
или
Если зададим априорное распределение на как стандартное нормальное , тогда и соответственно .
Теперь сделаем обратную замену и убедимся, что оценка максимального правдоподобия в точности равна .
(напомним, что и это вектор-строки). Заметим, что число есть след матрицы, состоящей из этого числа, поэтому можно преобразовать вторую часть, как
Отсюда следует, что
Приравняв производную по к нулю, найдем:
Оценка максимума правдоподобия на :
Эту оценку можно интерпретировать как среднюю потерю дисперсии по всем проигнорированным сингулярным направлениям. Если же — константа, то при получаем обычный PCA.
Другой способ получить обычный PCA — это вместо обычного EM-алгоритма воспользоваться его жёсткой модификацией.
Теперь предлагаем вам потренировать изученный материал на практике. Предлагаем вам выполнить лабораторную работу, которая покрывает большинство тем главы “Вероятностные модели”. Скачайте ноутбук с лабораторной работой. В нём вы найдете описания заданий и дополнительные материалы. Задания из лабораторной прикреплены к этому параграфу в виде задач в системе Яндекс Контест. Чтобы проверить себя, отправляйте решения по соответствующим задачам в систему. Успехов в практике!