# file "REM2" -- simple EM algorithm for a Poisson model with merged data
# data (simulated from Poisson(2.7):
nstar = 244 # the count for 0 *plus* the count for 1
count = c( 0, 227, 219, 158, 83, 46, 14, 8, 0, 1 )
# Here count[i] = the count for "i", for i >= 2
N = nstar + sum(count)
len = length(count)
# initial value:
lambda = 0.1
cat("i =", 0, " lambda =", lambda, "\n")
numits = 20
for (i in 1:numits) {
# E-step:
z = lambda/(1+lambda) * nstar
# M-step:
lambda = (1/N) * ( z + sum((1:len)*count) )
# output as we go:
cat("i =", i, " z =", z, " lambda =", lambda, "\n")
}