如何解决如何在R中存储每个混合成分的观测值
我想手动在R中拟合混合模型。然后,我想分别存储对混合模型每个组成部分的观察。也就是说,我希望我的代码保留从每个组件淹没的观察结果。这是一个使用here中的EM算法的混合模型示例。
e_step <- function(x,mu.vector,sd.vector,alpha.vector) {
comp1.prod <- dnorm(x,mu.vector[1],sd.vector[1]) * alpha.vector[1]
comp2.prod <- dnorm(x,mu.vector[2],sd.vector[2]) * alpha.vector[2]
sum.of.comps <- comp1.prod + comp2.prod
comp1.post <- comp1.prod / sum.of.comps
comp2.post <- comp2.prod / sum.of.comps
sum.of.comps.ln <- log(sum.of.comps,base = exp(1))
sum.of.comps.ln.sum <- sum(sum.of.comps.ln)
list("loglik" = sum.of.comps.ln.sum,"posterior.df" = cbind(comp1.post,comp2.post))
}
m_step <- function(x,posterior.df) {
comp1.n <- sum(posterior.df[,1])
comp2.n <- sum(posterior.df[,2])
comp1.mu <- 1/comp1.n * sum(posterior.df[,1] * x)
comp2.mu <- 1/comp2.n * sum(posterior.df[,2] * x)
comp1.var <- sum(posterior.df[,1] * (x - comp1.mu)^2) * 1/comp1.n
comp2.var <- sum(posterior.df[,2] * (x - comp2.mu)^2) * 1/comp2.n
comp1.alpha <- comp1.n / length(x)
comp2.alpha <- comp2.n / length(x)
list("mu" = c(comp1.mu,comp2.mu),"var" = c(comp1.var,comp2.var),"alpha" = c(comp1.alpha,comp2.alpha))
}
for (i in 1:50) {
if (i == 1) {
# Initialization
e.step <- e_step(wait,wait.summary.df[["mu"]],wait.summary.df[["std"]],wait.summary.df[["alpha"]])
m.step <- m_step(wait,e.step[["posterior.df"]])
cur.loglik <- e.step[["loglik"]]
loglik.vector <- e.step[["loglik"]]
} else {
# Repeat E and M steps till convergence
e.step <- e_step(wait,m.step[["mu"]],sqrt(m.step[["var"]]),m.step[["alpha"]])
m.step <- m_step(wait,e.step[["posterior.df"]])
loglik.vector <- c(loglik.vector,e.step[["loglik"]])
loglik.diff <- abs((cur.loglik - e.step[["loglik"]]))
if(loglik.diff < 1e-6) {
break
} else {
cur.loglik <- e.step[["loglik"]]
}
}
}
loglik.vector
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。