Dlaczego algorytm EM musi być iteracyjny?

9

Załóżmy, że masz populację N jednostki, każda z losową zmienną XiPoisson(λ). Ty obserwujeszn=Nn0wartości dla dowolnej jednostki, dla której . Chcemy oszacowania .Xi>0λ

Istnieją metody chwil i warunkowe maksymalne prawdopodobieństwo uzyskania odpowiedzi, ale chciałem wypróbować algorytm EM. Otrzymuję algorytm EM: gdzie indeks dolny wskazuje wartość z poprzedniej iteracji algorytmu, a jest stały w odniesieniu do Parametry. (Właściwie uważam, że we frakcji w nawiasach powinno wynosić , ale to nie wydaje się dokładne; pytanie na inny raz).

Q(λ1,λ)=λ(n+nexp(λ1)1)+log(λ)i=1nxi+K,
1Knn+1

Aby uczynić to konkretnym, załóżmy, że , . Oczywiście, i nie są obserwowane i należy oszacować .n=10xi=20Nn0λ

Gdy iteruję następującą funkcję, podając maksymalną wartość z poprzedniej iteracji, docieram do prawidłowej odpowiedzi (zweryfikowanej przez CML, MOM i prostą symulację):

EmFunc <- function(lambda, lambda0){
  -lambda * (10 + 10 / (exp(lambda0) - 1)) + 20 * log(lambda)
}

lambda0 <- 2
lambda  <- 1

while(abs(lambda - lambda0) > 0.0001){
  lambda0 <- lambda
  iter    <- optimize(EmFunc, lambda0 = lambda0, c(0,4), maximum = TRUE)
  lambda  <- iter$maximum
}

> iter
$maximum
[1] 1.593573

$objective
[1] -10.68045

Ale to jest prosty problem; zmaksymalizujmy bez iteracji:

MaxFunc <- function(lambda){
  -lambda * (10 + 10 / (exp(lambda) - 1)) + 20 * log(lambda)
}

optimize(MaxFunc, c(0,4), maximum = TRUE)
$maximum
[1] 2.393027

$objective
[1] -8.884968

Wartość funkcji jest wyższa niż w procedurze iteracyjnej, a wynik jest niezgodny z innymi metodologiami. Dlaczego druga procedura daje inną i (jak sądzę) nieprawidłową odpowiedź?

Charlie
źródło

Odpowiedzi:

6

Kiedy znalazłeś swoją funkcję celu dla algorytmu EM, zakładam, że potraktowałeś liczbę jednostek , które nazywam , jako parametr utajony. W tym przypadku zakładam (ponownie), że reprezentuje zmniejszoną formę oczekiwanej wartości w stosunku do prawdopodobieństwa, które podano . To nie jest taki sam, jak pełna najprawdopodobniej z powodu , że jest treadted podane.xi=0yQy λ1λ1

Dlatego nie możesz użyć dla pełnego prawdopodobieństwa, ponieważ nie zawiera ono informacji o tym, jak zmiana zmienia rozkład (i chcesz wybrać również najbardziej prawdopodobne wartości gdy zmaksymalizujesz pełne prawdopodobieństwo). Właśnie dlatego pełne maksymalne prawdopodobieństwo zerowego skróconego Poissona różni się od funkcji i dlatego otrzymujesz inną (i niepoprawną) odpowiedź, gdy maksymalizujesz .QλyyQf(λ)=Q(λ,λ)

Liczbowo, maksymalizacja z konieczności spowoduje, że funkcja celu będzie co najmniej tak duża, jak wynik EM, i prawdopodobnie większa, ponieważ nie ma gwarancji, że algorytm EM zbiegnie się do maksimum - powinien on tylko zbiegać się do maksimum funkcji wiarygodności !f(λ)f

jayk
źródło