gganimate　ポアソン分布

pois_anime.R

library(gganimate)
library(ggplot2)
#library(tidyr)

sim_pois1 <- function(t_n, p_n, w, seed){
  set.seed(seed)
  l <- 5 - 0.5*w
  r <- 5 + 0.5*w
  counts <- integer(t_n)
  p_list <- vector("list", t_n)
  for(i in 1:t_n){
    x <- runif(p_n, 0, 10)
    y <- runif(p_n, 0, 10)
    flag <- l < x & x < r & l < y & y < r
    p_list[[i]] <- data.frame(x=x,y=y,hit=flag,time=i)
    counts[i] <- sum(flag)
  }
  return(list(points=do.call("rbind",p_list),
              counts=counts))
}

sim_pois2 <- function(t_n, p, g_n, seed){
  set.seed(seed)
  counts <- integer(t_n)
  g_list <- vector("list", t_n)
  g <- expand.grid(x=seq(0,1,length.out=g_n),
                   y=seq(0,1,length.out=g_n))
  N <- nrow(g)
  for(i in 1:t_n){
    hit <- rbinom(N,1,p)
    g$hit <- hit
    g$time <- i
    g_list[[i]] <- g
    counts[i] <- sum(hit)
  }
  return(list(grid=do.call("rbind",g_list),
              counts=counts))
}


out <- sim_pois1(500,100,1,1)
out$points$facet_dummy <- "raw"
dfc <- do.call("rbind",lapply(1:500, function(i){
  tab <- table(out$counts[1:i])
  data.frame(time=i,
             num=as.integer(names(tab)),
             prob=as.integer(tab)/sum(tab),
             facet_dummy = "prob"
  )}))

dfseg <- data.frame(matrix(c(4.5,4.5,4.5,5.5,
                             4.5,4.5,5.5,4.5,
                             5.5,4.5,5.5,5.5,
                             4.5,5.5,5.5,5.5),4,4,byrow = TRUE),
                    facet_dummy = "raw")

rs <- range(out$counts)
dfprob <- data.frame(num=rs[1]:rs[2],
                      prob=dpois(rs[1]:rs[2],1),
                      facet_dummy = "prob")
ggplot()+
  geom_point(data=out$points, aes(x=x,y=y,colour=hit))+
  scale_color_viridis_d()+
  geom_col(data=dfc, aes(x=num, y=prob), width=0.7, position = "identity")+
  geom_point(data=dfprob, aes(x=num, y=prob))+
  geom_line(data=dfprob, aes(x=num, y=prob),linetype=2)+
  geom_segment(data=dfseg, aes(x=X1,y=X2,xend=X3,yend=X4), linetype=2)+
  facet_wrap(~facet_dummy, scales="free")+
  transition_time(time)+
  labs(title = 'Step: {frame_time}') +
  theme_bw(12)

anim_save("pois1.gif")

######

out2 <- sim_pois2(500,0.01,10,2)
out2$grid$facet_dummy <- "raw"

dfc2 <- do.call("rbind",lapply(1:500, function(i){
  tab <- table(out2$counts[1:i])
  data.frame(time=i,
             num=as.integer(names(tab)),
             prob=as.integer(tab)/sum(tab),
             facet_dummy = "prob"
  )}))



rs <- range(out2$counts)
dfprob2 <- data.frame(num=rs[1]:rs[2],
                     prob=dpois(rs[1]:rs[2],1),
                     facet_dummy = "prob")

ggplot()+
  geom_tile(data=out2$grid, aes(x=x,y=y,fill=hit))+
  scale_fill_viridis_c()+
  geom_col(data=dfc2,aes(x=num,y=prob), width=0.7, position = "identity")+
  geom_point(data=dfprob2, aes(x=num,y=prob))+
  geom_line(data=dfprob2, aes(x=num,y=prob), linetype=2)+
  facet_wrap(~facet_dummy, scales="free")+
  transition_time(time)+
  labs(title = 'Step: {frame_time}') +
  theme_bw(12)

anim_save("pois2.gif")