posterior density of the mixture normal distribution

lpmixnorm.R

library(ggplot2)
library(gganimate)

logsumexp=function (logx1,logx2){
  logx1 + log1p(exp(logx2-logx1))
}

llmixnorm <- function(par, y){
  a0 <- par[1]
  b0 <- par[2]
  theta <- par[3]
  sum(logsumexp(log(1-a0)+dnorm(y,log = TRUE), log(a0)+dnorm(y,b0,log = TRUE))) + dbeta(a0, theta, theta, log = TRUE)
}

N <- 100L
set.seed(1)
y11 <- c(rnorm(N/2),rnorm(N/2,0.5))

ggplot()+
  geom_histogram(data=NULL, aes(x=y11), fill="grey70", bins=25)+
  theme_bw(14)+labs(x="y", y="count")

parms <- expand.grid(a=seq(0.3,0.99,length.out = 200),
                     b=seq(0,1,length.out = 200),
                     theta=seq(0.5,3,by=0.1))

l11 <- apply(parms, 1, llmixnorm, y=y11)
dfL11 <- data.frame(parms, value = exp(l11))


ggplot(dfL11,aes(x=b,y=a))+
  geom_point(aes(colour=value))+
  scale_colour_continuous(type = "viridis")+
  labs(title = 'theta: {round(frame_time,2)}') +
  transition_time(theta)+
  theme_bw(14)+theme(legend.position = "none")

anim_save("post_mixnorm.gif")