CRPを試す

確率 プログラミング

はじめに

NowいProfessionalなYoungにバカうけなLMを作ってみたいなぁと思ったので、Chinese Restaurant Process(CRP)についてちょっと調べてみた。

CRPとは

  • 中華料理店過程(Chinese Restaurant Process)
  • 以下のルールによって、集合{1,...,n}の分割Bnを決める離散時間確率過程(時系列)
    • t=1のとき、B1 = {1}
    • t=n+1のとき、それまで生成された分割Bnに対し、次の確率で追加する
      • |ブロックiの要素数|/(n+1)の確率で、そのブロックiに追加する
      • 1/(n+1)の確率で、新たなブロックを作りそれに追加する
  • Btの例(「|」はブロックの切れ目)
B1 = {1}
B2 = {1,2}
B3 = {1,2|3}
B4 = {1,2,4|3}
B5 = {1,2,4|3|5}
...
  • で、これは「中華料理店」での"テーブル"と"お客さん"にして考えられるので「中華料理店過程」と呼ばれているみたい(インドのブッフェ版とかもある)
    • ブロック = テーブル
    • 要素 = お客さん
    • ブロックに追加する = お客さんがテーブルを選んで席につく
    • 「お店にテーブルがいくつもあり、次々とやってくるお客さんは、テーブルに座っているお客さんの数によって座るテーブルを決める」というイメージ

CRPの一般化

  • 上記の考え方はとても有用で、いろんな派生が存在しているみたい
    • 確率の式の変形
    • 料理の概念の追加(Labeled CRP)
    • など？
確率の式の変形
  • 2つのパラメータαとθで一般化する
  • n+1番目のお客さんが座る席を決める
    • 新しいテーブルに座る確率 : \frac{\theta + |B|\alpha}{n+\theta}
    • すでにあるテーブルに座る確率 : \frac{|b|-\alpha}{n+\theta}
      • |B| : 人が座っているテーブルの数
      • |b| : そのテーブルに座っているお客さんの数
      • α,θ : パラメータ(ただし、確率測度であるために制約がある)
料理の概念の追加
  • 各テーブルには一つの料理が置かれ、お客さんは座ったテーブルにおいてある料理を食べる
  • この料理が何になるか？は、基底測度G0によって選ばれる
    • G0は何らかの確率分布(一様分布とか)

2パラメータなCRPを試してみる

パラメータによる変化を見てみるために下記のコードを試してみた。

コード
#include <fstream>
#include <sstream>
#include <iostream>
#include <vector>
#include <algorithm>

class ChineseRestaurantProcess {
  std::vector<int> tables; //テーブルに座っている人数
  int N; //人数
  //ハイパーパラメータ
  double alpha; //discount
  double theta; //strength

  //乱数
  // 注意: longではなくint(32bit)にすべき
  unsigned long xor128(){
    static unsigned long x=123456789,y=362436069,z=521288629,w=88675123;
    unsigned long t;
    t=(x^(x<<11));x=y;y=z;z=w; return( w=(w^(w>>19))^(t^(t>>8)) );
  }
  //現在のテーブルの状態から座るテーブルを決定
  size_t selectTable(){
    double z = (xor128()%1000000001)/1000000000.0;
    size_t ret = 0;
    for(size_t i=0; i<tables.size(); i++){
      double pi = (tables[i]-alpha)/static_cast<double>(N+theta);
      z -= pi;
      if(z < 0.0) return ret;
      ret++;
    }
    return ret;
  }
public:
  ChineseRestaurantProcess(double alpha, double theta):
    alpha(alpha),theta(theta)
  {
    tables.push_back(1); //最初の客は必ず最初のテーブルに座る
    N = 1; //最初はぼっち飯
  }
  //お客さんを一人追加
  void addCustomer(){
    size_t idx = selectTable();
    if(idx>=tables.size()){ //new table
      tables.push_back(1);
    }else{ //otherwise
      tables[idx]++;
    }
    N++;
  }
  //gnuplot用output
  void gnuplot_output(std::string filename){
    std::ofstream ofs(std::string(filename + ".plt").c_str());

    std::vector<int> copy_tables(tables);
    std::sort(copy_tables.begin(), copy_tables.end(), std::greater<int>());

    ofs << "set terminal png size 300,300" << std::endl;
    ofs << "set output \"" << filename << ".png\"" << std::endl;
    ofs << "set title \"" << filename << "\"" << std::endl;
    ofs << "set size ratio 1" << std::endl;
    ofs << "set nokey" << std::endl;
    ofs << "set xrange [1:100000]" << std::endl;
    ofs << "set yrange [1:100000]" << std::endl;
    ofs << "set logscale xy" << std::endl;

    ofs << "plot \"-\" w lp" << std::endl;
    for(size_t i=0; i<copy_tables.size(); i++){
      ofs << copy_tables[i] << std::endl;
    }
    ofs << "e" << std::endl;
  }
};

int main(){
  const double alpha = 1.0;
  const double theta = 10.0;

  ChineseRestaurantProcess crp(alpha, theta);
  //1000000人のお客さんに来てもらう
  for(int i=1; i<1000000; i++){
    crp.addCustomer();
  }

  std::stringstream ss;
  ss << alpha << "_" << theta;
  crp.gnuplot_output(ss.str());
  return 0;
}
結果
  • theta = 10.0
  • 縦軸 : テーブルに座っているお客さんの数
  • 横軸 : 各テーブル(座ってる人の多い順)
  • alpha = 0.0

  • alpha = 0.1

  • alpha = 0.6

  • alpha = 1.0

alpha>0のとき両対数で直線(power-law)になってる

応用

  • CRPの何がいいかというと、「Dirichlet過程やPitman-Yor過程と等価で扱いやすい」
    • この2つの過程はノンパラベイズでとても有用なもの
    • 特に、後者のPitman-YorはPower-lawな性質を持っていて、言語モデルとして利用するととても優れた性能を示すことがわかっている
  • 料理の概念を考えることで階層化を考えることができる

関連モデル

  • Stick Breaking Process
  • Polya's urn model

参考ページ