Displaying 1 result from an estimated 1 matches for "m_step".
Did you mean:
pm_step
2007 May 11
0
EM covergence problem
...(D==0,log(p1),log(p2))
}
return (probs)
}
#H [individuals][class]
E_step = function(alpha,Beta){#calc posterior of H
tmpH = matrix(,nrow = 1000,ncol =3)
lprobs = logProbInd(Beta)
for(i in 1:3){#classes
tmpH[,i] = alpha[i]*exp(lprobs[,i])
}
H = tmpH /apply(tmpH,1,sum)
return( H)
}
M_step = function(H,Beta){
#first part use direct estimation
aita = apply(H,2,sum)/1000
opt.c = optim(Beta,Q2,H=H,method="BFGS",control = list(fnscale = -1))
lik = opt.c$value
return(c(aita,opt.c$par,lik))
}
#EM loops
alf = c(0.33,0.33,0.33)
Bt = seq(0.1,0.6,by=0.1)
sc = rep(-8000,5)...