em.hmm<-function(x, L, maxiter=300) { ##################################################################################### ## USAGE # em2.hmm(x, L) ## ARGUMENTS # x: the observed data # L: num. of components in the non-null mixture # maxiter: the maximum number of iterations ## DETAILS # em0.hmm calculates the MLE for a HMM model with hidden states being 0/1. # the distribution of state 0 is assumed to be known as N(0, 1) # the distribution of state 1 is assumed to be a normal mixture with L components ## VALUES # fuction 'em0.hmm' gives the MLE of model parameters and Lfdr estimates for a HMM # pii: the initial state distribution # A=(a00 a01\\ a10 a11): transition matrix # pc, i from 1 to L: # --probability weights of each component in the non-null mixture # f1: an L by 2 matrix # --specifying the dist. of each component in the non-null mixture # lfdr: the lfdr variables # niter: number of iterations #################################################################################### NUM<-length(x) # precision tolerance level ptol<-1e-3 niter<-0 f0<-c(0, 1) ### initializing model parameters pii.new<-c(1, 0) A.new<-matrix(c(0.95, 0.05, 0.5, 0.5), 2, 2, byrow=T) pc.new<-rep(1, L)/L mus<-seq(from=-1, by=1.5, length=L) sds<-rep(1, L) f1.new<-cbind(mus, sds) diff<-1 ### The E-M Algorithm while(diff>ptol && niter