Примеры подбрасывания монеты с максимизацией ожидания

Второй PDF-файл для меня не загрузится, но я также посетил страницу википедии http://en.wikipedia.org/wiki/Expectation%E2%80%93maximization_algorithm с дополнительной информацией. http://melodi.ee.washington.edu/people/bilmes/mypapers/em.pdf (который утверждает, что это мягкое введение) тоже заслуживает внимания.

Весь смысл алгоритма EM состоит в том, чтобы найти параметры, которые максимизируют вероятность наблюдаемых данных. Это единственный пункт на странице 8 первого PDF-файла, уравнение для нижнего индекса ML с заглавной буквы Theta.

Алгоритм EM пригодится там, где есть скрытые данные, которые, если бы вы знали, упростили бы проблему. В примере с тремя монетами это результат подбрасывания монеты 0. Если бы вы знали результат этого, вы могли бы (конечно) произвести оценку вероятности выпадения орла на монете 0. Вы также должны знать, подбрасывалась ли монета 1 или монета 2 три раза на следующем этапе, что позволит вам оценить вероятность того, что монета 1 и монета 2 окажутся орлом. Эти оценки будут оправданы, если заявить, что они максимизируют вероятность наблюдаемых данных, которые будут включать не только результаты, которые вам даны, но также скрытые данные, которыми вы не являетесь - результаты для монеты 0. Для монеты, которая получает A орел и решка B вы обнаружите, что максимальная вероятность для вероятности выпадения орла равна A / (A + B) - возможно, вам стоит проработать это подробно, потому что это строительный блок для шага M.

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

Алгоритм EM просит вас найти параметры, максимизирующие взвешенную сумму логарифмических вероятностей, заданных всеми возможными значениями скрытых данных, где веса задаются вероятностью связанных скрытых данных с учетом наблюдений с использованием параметров в начало шага ЭМ. Это то, что почти все, включая алгоритм Википедии, называют Q-функцией. Доказательство, лежащее в основе алгоритма EM, приведенное в статье в Википедии, гласит, что если вы измените параметры так, чтобы увеличить Q-функцию (что является только средством для достижения цели), вы также измените их, чтобы увеличить вероятность наблюдаемых данных (которые вас действительно волнуют). На практике вы, как правило, обнаруживаете, что вы можете максимизировать Q-функцию, используя вариант того, что вы бы сделали, если бы вы знали скрытые данные, но используя вероятности скрытых данных, учитывая оценки в начале EM- шаг, чтобы как-то взвесить наблюдения.

В вашем примере это означает суммирование количества орлов и решек, произведенных каждой монетой. В PDF получается P (Y = H | X =) = 0,6967. Это означает, что вы используете вес 0,6967 для случая Y = H, что означает, что вы увеличиваете счетчики для Y = H на 0,6967 и увеличиваете счетчики для X = H в монете 1 на 3 * 0,6967, а также увеличиваете счетчики для Y = T на 0,3033 и увеличиваем счетчики для X = H в монете 2 на 3 * 0,3033. Если у вас есть подробное обоснование того, почему A / (A + B) является максимальной вероятностью вероятностей монет в стандартном случае, вы должны быть готовы превратить его в обоснование того, почему эта схема взвешенного обновления максимизирует Q-функцию.

Наконец, логическая вероятность наблюдаемых данных (то, что вы максимизируете) дает вам очень полезную проверку. Он должен увеличиваться с каждым шагом EM, по крайней мере, до тех пор, пока вы не приблизитесь к сходимости настолько, чтобы возникла ошибка округления, и в этом случае у вас может быть очень небольшое уменьшение, сигнализирующее о сходимости. Если оно резко уменьшается, значит, у вас ошибка в программе или в математике.

Как назло, я недавно тоже боролся с этим материалом. Вот как я подумал об этом:

Рассмотрим связанный, но отличный алгоритм, называемый алгоритмом классификации-максимизации, который мы могли бы использовать в качестве метода решения проблемы смешанной модели. Проблема смешанной модели - это проблема, в которой у нас есть последовательность данных, которая может быть получена любым из N различных процессов, общая форма которых нам известна (например, гауссовский), но мы не знаем параметров процессов (например, средства и / или отклонения) и могут даже не знать относительную вероятность процессов. (Обычно мы, по крайней мере, знаем количество процессов. Без этого мы попадаем на так называемую «непараметрическую» территорию.) В некотором смысле процесс, который генерирует все данные, является «отсутствующими» или «скрытыми» данными. проблемы.

Теперь этот связанный алгоритм классификации-максимизации начинает с некоторых произвольных предположений о параметрах процесса. Каждая точка данных оценивается в соответствии с каждым из этих процессов параметров, и генерируется набор вероятностей - вероятность того, что точка данных была сгенерирована первым процессом, вторым процессом и т. Д., Вплоть до последнего N-го процесса. Затем каждая точка данных классифицируется в соответствии с наиболее вероятным процессом.

На данный момент наши данные разделены на N разных классов. Итак, для каждого класса данных мы можем с помощью некоторых относительно простых вычислений оптимизировать параметры этого кластера с помощью метода максимального правдоподобия. (Если мы попытаемся сделать это для всего набора данных перед классификацией, это обычно аналитически невозможно.)

Затем мы обновляем наши предположения о параметрах, переклассифицируем, обновляем наши параметры, переклассифицируем и т. Д. До сходимости.

Алгоритм максимизации ожидания аналогичен, но имеет более общий характер: вместо жесткой классификации точек данных на класс 1, класс 2, ... через класс N мы теперь используем мягкую классификацию, в которой каждая точка данных принадлежит каждый процесс с некоторой вероятностью. (Очевидно, что вероятности для каждой точки должны быть суммированы до единицы, поэтому происходит некоторая нормализация.) Я думаю, мы могли бы также думать об этом как о каждом процессе / предположении, имеющем определенную «объяснительную силу» для каждого из данных. точки.

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

С учетом сказанного, есть некоторые предостережения:

1) Звучит просто. Это не так, по крайней мере, для меня. Литература изобилует мешаниной специальных приемов и приемов - использование выражений правдоподобия вместо выражений вероятности, преобразование в логарифмические вероятности, использование индикаторных переменных, преобразование их в базисную векторную форму и включение их в экспоненты и т. Д.

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

2) Какие бы ограничения у вас ни были, может быть сложно включить в структуру. В частности, если вы знаете вероятности каждого из процессов, вы, вероятно, в хорошей форме. Если нет, вы также оцениваете их, и сумма вероятностей процессов должна быть равна единице; они должны жить на вероятностном симплексе. Не всегда очевидно, как сохранить эти ограничения нетронутыми.

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

4) Доказано, что этот метод сходимости, но сходимость не обязательно к глобальному максимуму; будьте осторожны.

Я счел следующую ссылку полезной при разработке приведенной выше интерпретации: Статистические обучающие слайды

И следующая статья подробно описывает некоторые болезненные математические детали: Описание Майкла Коллинза

Я написал приведенный ниже код на Python, который объясняет пример, приведенный в вашем втором примере статьи от До и Бацоглу.

Я рекомендую вам сначала прочитать эту ссылку, чтобы получить четкое объяснение того, как и почему получены значения 'weightA' и 'weightB' в приведенном ниже коде.

Заявление об ограничении ответственности: код действительно работает, но я уверен, что он не оптимален. Я обычно не программист на Python и начал использовать его две недели назад.

import numpy as np
import math

#### E-M Coin Toss Example as given in the EM tutorial paper by Do and Batzoglou* #### 

def get_mn_log_likelihood(obs,probs):
    """ Return the (log)likelihood of obs, given the probs"""
    # Multinomial Distribution Log PMF
    # ln (pdf)      =             multinomial coeff            *   product of probabilities
    # ln[f(x|n, p)] = [ln(n!) - (ln(x1!)+ln(x2!)+...+ln(xk!))] + [x1*ln(p1)+x2*ln(p2)+...+xk*ln(pk)]     

    multinomial_coeff_denom= 0
    prod_probs = 0
    for x in range(0,len(obs)): # loop through state counts in each observation
        multinomial_coeff_denom = multinomial_coeff_denom + math.log(math.factorial(obs[x]))
        prod_probs = prod_probs + obs[x]*math.log(probs[x])

multinomial_coeff = math.log(math.factorial(sum(obs))) -  multinomial_coeff_denom
likelihood = multinomial_coeff + prod_probs
return likelihood

# 1st:  Coin B, {HTTTHHTHTH}, 5H,5T
# 2nd:  Coin A, {HHHHTHHHHH}, 9H,1T
# 3rd:  Coin A, {HTHHHHHTHH}, 8H,2T
# 4th:  Coin B, {HTHTTTHHTT}, 4H,6T
# 5th:  Coin A, {THHHTHHHTH}, 7H,3T
# so, from MLE: pA(heads) = 0.80 and pB(heads)=0.45

# represent the experiments
head_counts = np.array([5,9,8,4,7])
tail_counts = 10-head_counts
experiments = zip(head_counts,tail_counts)

# initialise the pA(heads) and pB(heads)
pA_heads = np.zeros(100); pA_heads[0] = 0.60
pB_heads = np.zeros(100); pB_heads[0] = 0.50

# E-M begins!
delta = 0.001  
j = 0 # iteration counter
improvement = float('inf')
while (improvement>delta):
    expectation_A = np.zeros((5,2), dtype=float) 
    expectation_B = np.zeros((5,2), dtype=float)
    for i in range(0,len(experiments)):
        e = experiments[i] # i'th experiment
        ll_A = get_mn_log_likelihood(e,np.array([pA_heads[j],1-pA_heads[j]])) # loglikelihood of e given coin A
        ll_B = get_mn_log_likelihood(e,np.array([pB_heads[j],1-pB_heads[j]])) # loglikelihood of e given coin B

        weightA = math.exp(ll_A) / ( math.exp(ll_A) + math.exp(ll_B) ) # corresponding weight of A proportional to likelihood of A 
        weightB = math.exp(ll_B) / ( math.exp(ll_A) + math.exp(ll_B) ) # corresponding weight of B proportional to likelihood of B                            

        expectation_A[i] = np.dot(weightA, e) 
        expectation_B[i] = np.dot(weightB, e)

    pA_heads[j+1] = sum(expectation_A)[0] / sum(sum(expectation_A)); 
    pB_heads[j+1] = sum(expectation_B)[0] / sum(sum(expectation_B)); 

    improvement = max( abs(np.array([pA_heads[j+1],pB_heads[j+1]]) - np.array([pA_heads[j],pB_heads[j]]) ))
    j = j+1

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

Дополните каждую последовательность орла / решки двумя двоичными переменными, которые указывают, была ли использована монета 1 или монета 2. Теперь наши данные выглядят следующим образом:

c_11 c_12 c_21 c_22 c_31 c_32 ...

Для каждого i либо c_i1 = 1, либо c_i2 = 1, а другой равен 0. Если бы мы знали значения, которые эти переменные принимали в нашей выборке, оценка параметров была бы тривиальной: p1 было бы долей голов в выборках, где c_i1 = 1, аналогично для c_i2, а \ lambda будет средним значением c_i1s.

Однако нам неизвестны значения этих двоичных переменных. Итак, что мы в основном делаем, это угадываем их (на самом деле, принимаем их ожидания), а затем обновляем параметры в нашей модели, предполагая, что наши предположения были правильными. Таким образом, шаг E - это ожидание значений c_i1s и c_i2s. Шаг M заключается в получении оценок максимального правдоподобия для p_1, p_2 и \ lambda с учетом этих cs.

В этом есть немного больше смысла? Я могу выписать обновления для шагов E и M, если хотите. Тогда EM просто гарантирует, что, следуя этой процедуре, вероятность никогда не уменьшится по мере увеличения количества итераций.