記事の目的

Rを使用して、中華料理店過程を実装します。
参考: ノンパラメトリックベイズ点過程と統計的機械学習の数理

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. alphaの推定を含む実装

alphaの推定値は1.012666でした。

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()