기사의 목적

목차

1. 중국 음식점 과정의 수식

2. 구현

n <- 30
alpha <- 1
z <- 1
set.seed(100)
for(i in 1:n){
  n.k <- tapply(z,z,length)
  prob <- append(n.k/(i-1+alpha), alpha/(i-1+alpha))
  z.tmp <- which.max(rmultinom(1, 1, prob))
  z <- append(z, z.tmp)
}
#可視化
library(ggplot2)
ggplot(NULL, aes(x=z)) + geom_bar()

3. 알파 추정을 포함한 구현

n <- 30
alpha <- 1
z <- 1
c1 <- 10
c2 <- 10
set.seed(100)
for(i in 1:n){
  n.k <- tapply(z,z,length)
  prob <- append(n.k/(i-1+alpha), alpha/(i-1+alpha))
  z.tmp <- which.max(rmultinom(1, 1, prob))
  z <- append(z, z.tmp)
  #alphaの推定
  #piのサンプリング
  pi <- rbeta(1, alpha+1, n)
  #sのサンプリング
  s <- rbinom(1, 1, (n/alpha)/(1+n/alpha))
  #alphaのサンプリング
  alpha <- rgamma(1, c1+length(n.k)-s, c2-log(pi))
}
#可視化
library(ggplot2)
ggplot(NULL, aes(x=z)) + geom_bar()