pLSIを試してみた
これまでにK-means++とfuzzy c-meansを使用したクラスタリングを試してきましたが、今回はpLSI(probabilistic latent semantic indexing, 潜在的意味インデキシング)によるクラスタリングを試してみようと思います。
pLSIは確率・統計的な枠組みで次元縮約を行う枠組みで、なかなか精度がよいらしく色々な論文で見かけます。Google NewsのレコメンドでもpLSIを使用しており、MapReduceで処理を並列化させて高速に実行しているそうです(論文読んでないので間違っているかも)。また入力ベクトルをあらかじめ重み付けしておく必要がなく、文書であれば単語の頻度をそのまま入力として使用できるのもうれしいところです。
より詳しくは以下のWikipediaのエントリか、書籍をご参照下さい。(書籍は処理結果の表8.4が並びがグチャグチャになってるので注意！)
- 作者: 新納浩幸
- 出版社/メーカー: オーム社
- 発売日: 2007/11
- メディア: 単行本
- 購入: 9人 クリック: 207回
- この商品を含むブログ (29件) を見る
また以下の工藤拓さんによるC++実装も参考にさせていただきました。
今回作成したコードは以下の通りです。Perlでしかも愚直に作ってあるので、大きなデータはちょっと無理です。
#!/usr/bin/perl # # pLSI (Probabilistic Latent Semantic Indexing) # http://en.wikipedia.org/wiki/PLSI # use strict; use warnings; use Getopt::Long; use Data::Dumper; use constant { DEFAULT_NUM_ITERATION => 50, DEFAULT_BETA => 0.75, }; # Format: id\tval1,val2,val3,... sub read_vectors { my $fh = shift; my @vectors; while (my $line = <$fh>) { chomp $line; next if !$line; my ($id, $vecstr) = split /\t/, $line; if (!$id || !$vecstr) { warn "Illegal format: $line\n"; } my @vector = split /,/, $vecstr; push @vectors, { id => $id, vector => \@vector, }; } return \@vectors; } sub init_random { my ($nrow, $ncol) = @_; my @results; for (my $i = 0; $i < $nrow; $i++) { my @arr; my $sum = 0; for (my $j = 0; $j < $ncol; $j++) { my $r = rand(); push @arr, $r; $sum += $r; } @arr = map { $_ / $sum } @arr; push @results, \@arr; } return \@results; } sub expect { my ($pdz, $pwz, $pz, $num_doc, $num_word, $num_cluster, $beta) = @_; my @scores; for (my $id = 0; $id < $num_doc; $id++) { for (my $iw = 0; $iw < $num_word; $iw++) { my $denominator; for (my $iz = 0; $iz < $num_cluster; $iz++) { $denominator += $pz->[$iz] * (($pwz->[$iw][$iz] * $pdz->[$id][$iz]) ** $beta); } for (my $iz = 0; $iz < $num_cluster; $iz++) { my $numerator = $pz->[$iz] * (($pwz->[$iw][$iz] * $pdz->[$id][$iz]) ** $beta); $scores[$id][$iw][$iz] = $denominator ? $numerator / $denominator : 0; } } } return \@scores; } sub maximize { my ($scores, $documents, $pdz, $pwz, $pz, $num_doc, $num_word, $num_cluster) = @_; my (@pdz_new, @pwz_new, @pz_new); my @denominators; my $sum_denom; for (my $iz = 0; $iz < $num_cluster; $iz++) { my $denominator; for (my $id = 0; $id < $num_doc; $id++) { for (my $iw = 0; $iw < $num_word; $iw++) { $denominator += $documents->[$id][$iw] * $scores->[$id][$iw][$iz]; } } push @denominators, $denominator; # update pdz for (my $id = 0; $id < $num_doc; $id++) { my $numerator; for (my $iw = 0; $iw < $num_word; $iw++) { $numerator += $documents->[$id][$iw] * $scores->[$id][$iw][$iz]; } $pdz_new[$id][$iz] = $numerator / $denominator; } # update pwz for (my $iw = 0; $iw < $num_word; $iw++) { my $numerator; for (my $id = 0; $id < $num_doc; $id++) { $numerator += $documents->[$id][$iw] * $scores->[$id][$iw][$iz]; } $pwz_new[$iw][$iz] = $numerator / $denominator; } } # update pz my $denom_sum; map { $denom_sum += $_ } @denominators; for (my $iz = 0; $iz < $num_cluster; $iz++) { $pz_new[$iz] = $denominators[$iz] / $denom_sum; } return (\@pdz_new, \@pwz_new, \@pz_new); } sub plsi { my ($documents, $num_cluster, $num_iter, $beta) = @_; my $num_doc = scalar @{ $documents }; my $num_word = scalar @{ $documents->[0] }; my $pdz = init_random($num_doc, $num_cluster); my $pwz = init_random($num_word, $num_cluster); my $pz = [ map { 1 / $num_cluster } (0 .. $num_cluster-1) ]; for (my $i = 0; $i < $num_iter; $i++) { my $scores = expect( $pdz, $pwz, $pz, $num_doc, $num_word, $num_cluster, $beta); # print "pdz: ", Dumper($pdz); # print "pwz: ", Dumper($pwz); # print "pz: ", Dumper($pz); # print Dumper($scores); ($pdz, $pwz, $pz) = maximize($scores, $documents, $pdz, $pwz, $pz, $num_doc, $num_word, $num_cluster); } return ($pdz, $pwz, $pz); } sub main { my %option; GetOptions( \%option, 'number=i', 'beta=f', 'iter=i', ); my $path = shift @ARGV; if (!$path || !$option{number} || $option{number} <= 0) { print "Usage: plsi.pl -n ncluster [-b beta | -i niter] inputfile\n"; exit 1; } my $num_cluster = $option{number}; my $num_iter = exists $option{iter} ? $option{iter} : DEFAULT_NUM_ITERATION; my $beta = exists $option{beta} ? $option{beta} : DEFAULT_BETA; die 'iter must be >0' if $num_iter <= 0; die 'beta must be >0' if $beta <= 0; open my $fh, "<$path" or die "cannot open: $path"; my $vectors = read_vectors($fh); my @documents; foreach my $v (@{ $vectors }) { push @documents, $v->{vector}; } if ($num_cluster <= 0 || $num_cluster > scalar(@documents)) { die 'number of centers must be >0 and <=number_of_vectors'; } my ($pdz, $pwz, $pz) = plsi(\@documents, $num_cluster, $num_iter, $beta); for (my $i = 0; $i < scalar @{ $pdz }; $i++) { print $vectors->[$i]{id}."\t"; print join ' ', map { sprintf "%.4f", $_ } @{ $pdz->[$i] }; print "\n"; } } main();
実際に動かしてみます。実行結果はP(d|z)だけ表示しており、これをクラスタリングの結果と捉えることができます。
% cat input.txt 1 2,1,0,3,2 2 3,3,1,1,2 3 0,1,2,1,1 4 2,2,3,1,2 % ./plsi.pl -n 2 input.txt 1 0.0376 0.3015 2 0.1847 0.3372 3 0.3211 0.1026 4 0.4566 0.2587
入力文書1, 2はクラスタ2に主に属し、文書3, 4はクラスタ1に属していることが分かります。 書籍の説明ではEMの各ステップでP(d|z),P(w|z)の合計値を1.0に正規化してはいなかったのですが、C++実装では正規化されてました。ここは揃えた方が気持ちいい、程度の問題なんでしょうか？とりあえず今回作ったプログラムでは各ステップでの正規化はしていません。
もう少しわかりやすいデータでも試してみると、
% cat input2.txt 1 1,1,1,0,0,0 2 1,0,1,0,0,0 3 1,1,0,0,0,0 4 0,0,0,1,1,1 5 0,0,0,1,0,1 6 0,0,0,1,1,0 ./plsi.pl -n 2 /input2.txt 1 0.0000 0.4286 2 0.0000 0.2857 3 0.0000 0.2857 4 0.4286 0.0000 5 0.2857 0.0000 6 0.2857 0.0000
文書1,2,3がクラスタ2、文書4,5,6がクラスタ1にきっちり分かれました。ちゃんと動いてそうですね。
またbetaパラメータはEMのスムージング用のパラメータで、デフォルトは0.75になっていますが、これを低い値にするとよりなめらかな分布になるそうです(TemperedEMというらしい)。詳しくは調べていないので、パラメータの説明はC++実装のREADMEを参照してください。
% ./plsi.pl -n 2 input.txt 1 0.0376 0.3015 2 0.1847 0.3372 3 0.3211 0.1026 4 0.4566 0.2587 ./plsi.pl -n 2 -b 1.0 input.txt 1 0.0000 0.5428 2 0.2168 0.4099 3 0.2731 0.0008 4 0.5101 0.0465 ./plsi.pl -n 2 -b 0.65 input.txt 1 0.2412 0.2435 2 0.3027 0.3033 3 0.1522 0.1509 4 0.3038 0.3023
詳しく調べないで何回か試しただけなのですが、初期値によるブレが大きくて、そんなにいいものではないのかな感じたのですが、どうなんでしょう？10回ぐらい同じ処理をして、その平均を取るとかにした方がいいのかな…。
次はlda(latent Dirichlet allocation, 潜在的ディリクレ配分法)を試すぞ！といいたいところですが、説明読んでもよく分からない…orz なのでとりあえずpLSIをC++で書き直して、bayonにでも組み込んでみようかと思います。
それか、CPANにpLSIのパッケージがなかったので、ちゃんと作り直してAlgorithm::PLSIとかでパッケージングしたらうれしい人いるかな？結局pure Perlじゃ大きなデータは厳しいので、使い物にならないかもしれないけど^^;
