An animation of posterior distribution with mixture normal model

mixnorm.R

library(dplyr)
library(tidyr)
library(ggplot2)
library(gganimate)

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

llmixnorm <- function(par, y){
  a0 <- par[1]
  b0 <- par[2]
  theta <- par[3]
  sum(logsumexp2(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.01,0.99,length.out = 60),
                     b=seq(0,1,length.out = 60),
                     theta=seq(0.5,2,by=0.1))

l11 <- apply(parms, 1, llmixnorm, y=y11)
dfL11 <- data.frame(parms, posterior = exp(l11)) %>% 
  mutate(prior = dbeta(a, theta,theta)) %>% 
  pivot_longer(posterior:prior) %>% 
  group_by(theta, name) %>% 
  mutate(value = value/max(value)) %>% 
  ungroup()
  
head(dfL11)

pp <- ggplot(dfL11,aes(x=a,y=b,colour=value))+
  geom_point(pch=15, size=3.5)+
  scale_colour_continuous(type = "viridis")+
  labs(title = 'theta: {round(frame_time,2)}') +
  transition_time(theta)+
  facet_wrap(~name)+
  theme_bw(14)+theme(legend.position = "none")

animate(
  pp,
  renderer = gifski_renderer("post_mixnorm.gif")
)