Bayesian Setsを試してみた
この前YAPC Asia 2009に参加してきたのですが、そこで「はてなブックマークのシステムについて」の発表の中で、「はてブの関連エントリはBayesian Setsを使って計算されている」という話を聞いてBayesian Setsに俄然興味が湧いてきました。Bayesian Setsは以前論文だけ少し読んで、あまりよく分からないまま放置していたのですが、せっかくなのでPerlで作って試してみました。
Bayesian Setsについて詳しくは、以下のリンク先の資料をご参照下さい。
実際に作成したコードは以下の通りです。上記のMatlabのコードを参考にさせていただいています。
#!/usr/bin/perl # # Bayesian Sets # # Usage: # % bayesian_sets.pl input.tsv query1 query2 .. # # Reference # - Paper: http://www.gatsby.ucl.ac.uk/~heller/bsets.pdf # - Matlab code: http://chasen.org/~daiti-m/dist/bsets/ # use strict; use warnings; use constant { MAX_OUTPUT => 20, }; sub read_vectors { my $fh = shift; my %vectors; while (my $line = <$fh>) { chomp $line; my @arr = split /\t/, $line; my $doc_id = shift @arr; my %vector = @arr; map { $vector{$_} = 1 } keys %vector; $vectors{$doc_id} = \%vector; } return \%vectors; } sub get_parameters { my ($vectors, $c) = @_; my $average_vector = _average_vector($vectors); my (%alpha, %beta); while (my ($key, $val) = each %{ $average_vector }) { $alpha{$key} = $c * $val; $beta{$key} = $c * (1 - $val); } return (\%alpha, \%beta); } sub calc_similarities { my ($vectors, $queries, $alpha, $beta) = @_; my %query_vector; my $num_query; foreach my $query (@{ $queries }) { if (!exists $vectors->{$query}) { warn "[WARN] '$query' doesn't exist in input documents. Skip it\n"; next; } while (my ($key, $val) = each %{ $vectors->{$query} }) { $query_vector{$key} += $val; } $num_query++; } return if !%query_vector; my $length_qvec = scalar(keys %query_vector); my %weight_vector; foreach my $key (keys %{ $alpha }) { my $val = $query_vector{$key} || 0; $weight_vector{$key} = log(1 + $val / $alpha->{$key}) - log(1 + ($num_query - $val) / $beta->{$key}); } my %score; while (my ($doc_id, $vector) = each %{ $vectors }) { $score{$doc_id} = _inner_product(\%weight_vector, $vector); } return \%score; } sub _average_vector { my $vectors = shift; my %average_vector; my $num_vector = 0; while (my ($doc_id, $vector) = each %{ $vectors }) { while (my ($key, $val) = each %{ $vector }) { $average_vector{$key} += $val; } $num_vector++; } map { $average_vector{$_} /= $num_vector } keys %average_vector; return \%average_vector; } sub _inner_product { my ($v1, $v2) = @_; return 0 if !$v1 || !$v2; my @keys = scalar(keys %{ $v1 }) < scalar(keys %{ $v2 }) ? keys %{ $v1 } : keys %{ $v2 }; my $prod = 0; foreach my $key (@keys) { $prod += $v1->{$key} * $v2->{$key} if $v1->{$key} && $v2->{$key}; } return $prod; } sub main { my $path = shift @ARGV; my @queries = @ARGV; if (!$path || !@queries) { warn "Usage $0 file query1 query2 ..\n"; exit 1; } open my $fh, $path or die "cannot open file: $path"; my $vectors = read_vectors($fh); return if !$vectors; my ($alpha, $beta) = get_parameters($vectors, 2); # c=2 my $score = calc_similarities($vectors, \@queries, $alpha, $beta); my $count = 0; foreach my $doc_id (sort { $score->{$b} <=> $score->{$a} } keys %{ $score }) { last if ++$count > MAX_OUTPUT; printf "%s\t%.3f\n", $doc_id, $score->{$doc_id}; } } main();
入力データには以下のようなタブ区切りのテキストファイルを用意します。データ例は以前の記事で作成した、Wikipediaの各キーワードをTFIDFで重み付けしたデータです。実際にはvalue部分はいらないのですが、以前の名残りでそのままのフォーマットにしています。
# document_id \t key1 \t val1 \t key2 \t val2 \t ... 龍淵郡 郡 2365 長淵 1857 半島 1555 淵 1460 龍 1313 朝鮮 1137 (/ 1126 九 987 串 971 長山 788 リョンヨン 777 モングンポ 777 クミポ 777 朝鮮民主主義人民共和国 777 美浦 733 派 720 教会 714 救面 708 ざんかん 708 チャンサンゴッ 708 金浦 605 史上 595 松川 555 プロテスタント 529 分割 527 地 526 岬 513 里 506 建設 505 苔 500 年 496 南岸 490 当時 487 珪砂 483 黄海 481 再編 468 夢 460 ら 460 北朝鮮 444 宣教 404 1884 397 相 397 ソレ 382 甕 379 白島 360 徐 359 灘 344 教会堂 341 景勝 336 崙 334 神さまこんにちは ピョンテ 2332 チュンジャ 2332 慶州 1381 映画 1270 人 929 旅 827 チョン・ムソン 777 キム・ボヨン 777 ハン・ミヌ 777 国際 756 ペ・チャンホ 708 脳性 689 ソンギ 639 ミヌ 639 彼 624 祭 615 詩人 614 3 574 妊婦 548 麻痺 539 回 481 念願 468 英 458 39 430 Hello 425 ロード 419 アン 414 家 413 モスクワ 409 電車 399 God 395 16 377 乗車 366 遠足 350 神さま 342 無賃 339 !) 332 さ 331 大人 329 小学生 309 ムービー 298 ベルリン 296 とき 295 実家 289 一 288 原題 286 途方 272 神様 264 韓国 263 公開 259 藤原雅長 年 5208 雅 1285 藤原 1271 近臣 1263 位 1138 従 1118 任 1110 長 957 元年 903 3 861 7 813 院 798 近衛 782 五位上民部権大輔 777 嘉 770 視 764 4 760 下 747 日 747 左 713 隆信 708 月 689 わら 638 叙 629 保 600 8 599 能 595 参議 572 守 556 蔵人 521 久安 515 五 507 少将 506 26 503 顕能 500 在任 460 2 451 ふじ 448 中将 448 19 445 処分 431 久 425 死去 424 任命 422 兼 416 教 406 後期 403 長男 401 この間 377 辞任 377
実際に動かしてみます。入力ドキュメント集合のファイルと、(複数の)クエリを指定します。クエリは入力ドキュメント集合中に含まれている document_id しか使用できません。
% ./bayesian_sets.pl input.tsv トヨタ自動車 トヨタ自動車 140.748 豊田英二 47.205 トヨタ・トヨエース 45.501 5チャンネル 42.303 ポルシェ 40.428 トヨタ・マークX 37.372 日産・マーチ 36.659 フォード・モーター 36.333 三菱自動車工業 33.872 マツダ 33.625 ホンダ・インサイト 32.407 トヨタ・タンドラ 32.004 トヨタ・クルーガー 29.959 ジオ (自動車) 29.293 シボレー・キャバリエ 29.142 トヨタF1 28.830 三菱・ディオン 27.994 いすゞ・エルフ 27.712 ヒュンダイ・ソナタ 27.452 富士重工業 27.111 % ./bayesian_sets.pl input.tsv 原宿警察署 原宿警察署 203.638 滝野川警察署 80.618 牛込警察署 59.605 万世橋警察署 52.241 三田警察署 (東京都) 47.331 城東警察署 (東京都) 43.736 王子警察署 42.404 丸の内警察署 37.594 中野警察署 (東京都) 36.628 広島東警察署 32.846 田無警察署 31.444 豊田幹部交番 31.430 小金井警察署 31.352 尾鷲警察署 30.760 日下部警察署 30.656 山梨県警察 30.290 戸塚町 30.032 阿東幹部交番 29.785 神宮前 (渋谷区) 28.061 倶知安警察署 27.630
ちゃんと関連しているものが取得できていますね。本当は「apple banana」と「apple iPod」のように複数の意味を持つ言葉で、複数のクエリを入れた場合にそのクラスタ内の単語が出てくるのを確かめられるとよかったのですが。
今回使用した入力データのサイズは10万個のドキュメント(wikipediaのキーワード)ですが、実行には25秒ほどかかっています。今回のコードは入力データからのインデックスの作成と、クエリからの類似アイテムの検索の両方を同時に処理させているため動作が遅いですが、インデックスは初めに作成して保存しておけるようにしておけば、そこそこのスピードで検索できるようになると思います。
まだBayesian Setsの肝をあまり理解できていないのですが、
- 与えられたクエリが含まれていると思われるクラスタをオンデマンドで求め、そのクラスタに含まれている他のアイテムを返す
- クエリから求めた重みベクトルと入力ドキュメント集合から構成した行列の積(重みベクトルと全文書の内積)を求めるだけで計算できる
あたりがうれしいところなんでしょうか?一般的な転置インデックスを使用した類似コンテンツを求める方法でも同様のことはできそうな気がするのですが、「オンデマンドにクラスタを求める」というところでゴミが消されるんですかね。ちなみに単純に全文書とのcosine距離で類似したものを取ってきた場合だと、以下のような結果になります。大抵は良好な結果が得られますが、たまにクエリによってはゴミだらけになっちゃってますね。
% ./cos.pl input.tsv トヨタ自動車 トヨタ自動車 1.000 豊田英二 0.625 松井石根 0.610 島崎藤村 0.609 鉄道の歴史 (日本) 0.598 土佐電気鉄道 0.592 日産・マーチ 0.591 新関八洲太郎 0.588 近衛道嗣 0.588 岡村貢 0.587 高田保馬 0.587 8月7日 0.586 仙台東照宮 0.584 3月25日 0.584 飛鳥時代 0.583 三菱自動車工業 0.582 蒲郡競艇場 0.580 富士重工業 0.580 大林組 0.579 ホンダ・インサイト 0.576 % ./cos.pl input.tsv 原宿警察署 原宿警察署 1.000 阿東幹部交番 0.776 周南警察署 0.775 竹原警察署 0.773 高梁警察署 0.771 日下部警察署 0.769 大牟田警察署 0.764 三原警察署 0.758 音戸警察署 0.756 横須賀警察署 0.754 城東警察署 (東京都) 0.753 山口南警察署 0.752 松浦警察署 0.752 王子警察署 0.751 滝野川警察署 0.750 東近江警察署 0.750 つがる警察署 0.748 万世橋警察署 0.743 あわら警察署 0.743 丸亀警察署 0.743
あと、この前のpLSIの記事では結局そのまま放置になってるのですが、せっかくなので今度こそパッケージングしてAlgorithm::BayesianSetsを作ろうと思います。あまり使い道はなさそうですが、ざっくりどんなものか知るのには多少役に立ちそうですし、何より僕自身がCPAN Authorになってみたいので…^^;
追記:weight_vectorを求めるところが一部間違っていたため、コードと実行結果、実行時間を修正しました。ただ実行結果はポイント以外はほとんど違いはないと思います。
東村アキコ『ママはテンパリスト』1、2巻
- 作者: 東村アキコ
- 出版社/メーカー: 集英社
- 発売日: 2008/10
- メディア: コミック
- 購入: 32人 クリック: 319回
- この商品を含むブログ (255件) を見る
- 作者: 東村アキコ
- 出版社/メーカー: 集英社
- 発売日: 2009/06/19
- メディア: コミック
- 購入: 23人 クリック: 70回
- この商品を含むブログ (139件) を見る
同じ作者がモーニングで連載している「ひまわりっ〜健一レジェンド」が非常に好きなので、何となく買ってみた。作者の日々の子育ての話が書かれているんだけど、これが非常に面白い!子供ってこんなに面白い反応するものなのか。
あとちょこちょこ知らなかった話も出てきて、それもまた面白い。出産のときに会陰切開をする場合があるらしいけど、それは麻酔なしで切るらしい。痛そうだ…と思ったけど、出産の痛みの方がすごすぎて、切ってもほとんど気にならないぐらいらしい。あと本に出てくる子供は3歳ぐらいになって、普通の食事を食べるようになってからも母乳を飲んでいる。ってそんなにながく母乳って出続けるものなのか。
続きが非常に待ち遠しい漫画だ。
最近読んだ本
- 作者: 宮尾岳
- 出版社/メーカー: 少年画報社
- 発売日: 2009/08/07
- メディア: コミック
- 購入: 1人 クリック: 8回
- この商品を含むブログ (14件) を見る
アオバ自転車店の新刊。相変わらず安心して読める。最近ちょっとロードレーサーが欲しくなってきた。
- 作者: 二ノ宮知子
- 出版社/メーカー: 講談社
- 発売日: 2009/08/10
- メディア: コミック
- 購入: 20人 クリック: 257回
- この商品を含むブログ (300件) を見る
今回はのだめとシュトレーゼマンとの共演の話がほとんど。なんか最近あまり面白くないなぁ。ネタ切れなんだろうか。内容的にもあと1、2巻で終わりそうな気配。
雑誌での連載も再開したようだけど、「最後に伝えたいこと。」なんて副題が付いていたので、そろそろ締めの準備をしているのかな。
- 作者: 羽海野チカ
- 出版社/メーカー: 白泉社
- 発売日: 2009/08/12
- メディア: コミック
- 購入: 14人 クリック: 343回
- この商品を含むブログ (479件) を見る
この漫画は本当に大好き。読んでて暖かい気持ちになったり、主人公の孤独に同感したり。将棋も、そんなに深くはつっこんでないけど、緊迫感があって面白い。早く続きが読みたい。
- 作者: きらたかし
- 出版社/メーカー: 講談社
- 発売日: 2006/07/06
- メディア: コミック
- 購入: 2人 クリック: 17回
- この商品を含むブログ (7件) を見る
ヤンマガの「赤灯えれじい」「ケッチン」の人が、バイクに付いての作者自身の実体験を1話2ページで書いた本。これ見てるとバイクって楽しそうだなーと思う。本では大体中古でバイクを買っているんだけど、中古だとかなり安く買えていて、これなら自分もすぐに始められるかな?とか妄想したり。ただやっぱりバイクは事故が怖いので、なかなか手が出せないかな…。
シェルパ斉藤のリッター60kmで行く!日本全国スーパーカブの旅
- 作者: 斉藤政喜
- 出版社/メーカー: 小学館
- 発売日: 2009/07
- メディア: 単行本
- 購入: 3人 クリック: 17回
- この商品を含むブログ (10件) を見る
バイクで九州八十八ヶ所巡礼の旅に出たり、北海道を回ったり、東北を回ったり。カブなのでスローな旅ができて、かつ自転車のようにきつくはなくて、まったりと旅ができるのは本当に羨ましい。宿泊はテントや、ライダーハウスで他の旅行者と同じ部屋になって仲良くなったりして、楽しそうだな−と思う。
- 作者: 沼尾ひろ子
- 出版社/メーカー: エイ出版社
- 発売日: 2009/05/09
- メディア: 文庫
- クリック: 2回
- この商品を含むブログ (1件) を見る
こっちは自転車で輪行していろんなところに旅行に行っている。自転車は疲れたら、乗るの止めて電車とか飛行機に載せて移動しちゃえるのがいい。ただこの本はなんか自分にはフィーリングが合わなくて、途中で文章を読むのをやめて、写真だけさらっと見て終わりにしてしまった。
- 作者: 松田公太
- 出版社/メーカー: 新潮社
- 発売日: 2005/03/27
- メディア: 文庫
- 購入: 19人 クリック: 146回
- この商品を含むブログ (79件) を見る
日本のタリーズの創業者が、小さいころからどのような生活を送ってきて、どのようにして日本でタリーズを開店、展開していったかを創業者自身の言葉で書かれた本。この人の人生は格好いいなぁ。自分の人生における使命を「食を通じて文化の架け橋になる」というように早くから決めていて、それに向かってしっかりと人生を進んでいる。この本読んでから、スタバとか他のシアトル系コーヒーショップに行くよりも、どうせならタリーズに行こうかな、とか思うようになった。歴史を知ると愛着がわいてくるんですかね。うちの会社も、社長がこういう本書けばいいのに。
2006年 - 一部役員の画策により、ファンド各社からバイアウトのターゲットとなる。その危機を回避するために、伊藤園に安定株主として過半数の株式を売却。それに伴い、フードエックスグローブ株式会社社長を退任。
なんてことが書いてある。一体何があったんだろう…?
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じゃ大きなデータは厳しいので、使い物にならないかもしれないけど^^;
鈴木ともこ『山登りはじめました』
- 作者: 鈴木ともこ
- 出版社/メーカー: メディアファクトリー
- 発売日: 2009/06/17
- メディア: 単行本(ソフトカバー)
- 購入: 6人 クリック: 123回
- この商品を含むブログ (41件) を見る
なんとなくアウトドアには憧れがあるんです。ただ大変そうなイメージがあるのと、本気でやろうと思うとシューズやらバッグやらレインウェアやら…とすごく金がかかりそうな気がして、なかなか手がでません。あと近頃は全く運動してないので、体力も心配。
この本は元々は運動が苦手でインドア派だった作者が、いろんな山を登るようになって最終的に富士山に登るまでの話が漫画で書かれています。写真も掲載されていて、綺麗な景色も楽しめる。行く先々での食事もすごくおいしそうで、自分も行きたくなってきます。自分のペースでじっくりのんびり行くのが楽しむコツなんですかね。
この本にも書いてありましたが、日本人に生まれたからには一度は富士山に登ってみたい!。でも富士山って7,8月しか登れないんですね。今年はちょっと無理そうなので、来年登るのを目標にして、この本に書いてある場所に行ってみようかなぁと妄想中。まずは高尾山か鎌倉アルプスですかね。
- 作者: 石塚真一
- 出版社/メーカー: 小学館
- 発売日: 2008/01/30
- メディア: コミック
- 購入: 1人 クリック: 39回
- この商品を含むブログ (71件) を見る
そうそう、山の本といえば「岳」。今回この本を読んでて、岳に出てきたフレーズを思い出した。
遠くで見たときずっとキレイになるんだよ 頂上まで登った山は
岳は本格的な山登りなので、自分も登りたいという気持ちにはなかなかならなかったけど、今回の本はもっと気軽に登れそうな気がしてすごくよかった。
STLのvectorとpriority_queueのソート用比較関数は不等号が逆
この前自分のソースを読んでいたら、両方とも降順にソートするために作った比較関数なのに何故か不等号が逆になっていて、「うわ、ひどいバグ作っちゃった?!」って慌ててテストしたら問題なし。調べてみると、STLのvectorとpriority_queueのソート用比較関数は何故か不等号が逆になっちゃうんですね。
#include <cstdlib> #include <vector> #include <queue> #include <iostream> struct StructGreater { bool operator() (int a, int b) { return a < b; } }; bool funcGreater(int a, int b) { return a > b; } int main() { std::priority_queue<int, std::vector<int>, StructGreater> q; std::vector<int> v; size_t max = 10; for (size_t i = 0; i < max; i++) { int n = rand() % 100; q.push(n); v.push_back(n); } sort(v.begin(), v.end(), funcGreater); for (size_t i = 0; i < max; i++) { std::cout << "v:" << v[i] << "\tq:" << q.top() << std::endl; q.pop(); } return 0; }
ランダムに生成した整数列をvectorとpriority_queueに入れて、ソートして出力してみます。StructGreaterとfuncGreaterで不等号は逆なのですが、結果はちゃんと両方とも降順にソートされています。
% g++ sort.cc % ./a.out vector:93 queue:93 vector:92 queue:92 vector:86 queue:86 vector:86 queue:86 vector:83 queue:83 vector:77 queue:77 vector:49 queue:49 vector:35 queue:35 vector:21 queue:21 vector:15 queue:15
知ってる人にとっては当たり前なのかもしれないけど、念のためメモメモ。
Fuzzy c-meansを試してみた
K-meansは各入力ドキュメントがただ1つのクラスタにのみ属するハードクラスタリング手法ですが、fuzzy c-meansは所属度を持って複数のクラスタへの所属を許すソフトクラスタリング手法です。K-meansは以前に作りましたので、今回はfuzzy c-meansを試したいと思います。
fuzzy c-meansの詳しいアルゴリズムは下のページを参照してください。
Cluster analysis - Wikipedia, the free encyclopedia
とりあえずPerlで実装しました。ソースの半分くらいは、以前作ったK-means++からのコピペです^^;
#!/usr/bin/perl # # Fuzzy c-means clustering # http://en.wikipedia.org/wiki/Cluster_Analysis#Fuzzy_c-means_clustering # use strict; use warnings; use Data::Dumper; use List::Util qw(shuffle); use List::MoreUtils qw(any); sub choose_random_centers { my ($vectors, $num_center) = @_; my @ids = keys %{ $vectors }; @ids = shuffle @ids; my @centers = map { $vectors->{$_} } @ids[0 .. $num_center-1]; return \@centers; } sub choose_smart_centers { my ($vectors, $num_center) = @_; my $cur_potential = 0; my @centers; # choose one random center my $vector = (shuffle values %{ $vectors })[0]; push @centers, $vector; my %closest_dist; foreach my $id (keys %{ $vectors }) { $closest_dist{$id} = squared_euclid_distance($vectors->{$id}, $vector); $cur_potential += $closest_dist{$id}; } # choose each center for (my $i = 1; $i < $num_center; $i++) { my $randval = rand() * $cur_potential; my $center_id; foreach my $id (keys %{ $vectors }) { $center_id = $id; last if $randval <= $closest_dist{$id}; $randval -= $closest_dist{$id}; } my $new_potential = 0; foreach my $id (keys %{ $vectors }) { my $dist = squared_euclid_distance( $vectors->{$id}, $vectors->{$center_id}); $closest_dist{$id} = $dist if $dist < $closest_dist{$id}; $new_potential += $closest_dist{$id}; } push @centers, $vectors->{$center_id}; $cur_potential = $new_potential; } return \@centers; } # 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; $vectors{$id} = \@vector; } return \%vectors; } sub squared_euclid_distance { my ($vec1, $vec2) = @_; my $size = scalar @{ $vec1 }; my $dist = 0; for (my $i = 0; $i < $size; $i++) { $dist += ($vec1->[$i] - $vec2->[$i]) * ($vec1->[$i] - $vec2->[$i]); } return $dist; } sub get_distances { my ($vectors, $centers) = @_; my %distances; foreach my $id (keys %{ $vectors }) { foreach my $center (@{ $centers }) { my $dist = squared_euclid_distance($vectors->{$id}, $center); push @{ $distances{$id} }, $dist; } } return \%distances; } sub calc_weights { my ($vectors, $centers, $m) = @_; my $distances = get_distances($vectors, $centers); return if !$distances; my $num_center = scalar @{ $centers }; my %weights; foreach my $id (keys %{ $distances }) { if (any { $_ == 0 } @{ $distances->{$id} }) { foreach my $dist (@{ $distances->{$id} }) { push @{ $weights{$id} }, $dist == 0 ? 1 : 0; } } else { for (my $i = 0; $i < $num_center; $i++) { my $weight; for (my $j = 0; $j < $num_center; $j++) { my $x = $distances->{$id}[$i] / $distances->{$id}[$j]; $weight += $x * $x; } $weight **= (-1) / ($m - 1); push @{ $weights{$id} }, $weight; } } } return \%weights; } sub calc_centers { my ($vectors, $weights, $num_center, $m) = @_; my @centers; # initialize centers my $veclen = scalar @{ (values %{ $vectors })[0] }; for (my $i = 0; $i < $num_center; $i++) { my @arr; for (my $j = 0; $j < $veclen; $j++) { push @arr, 0; } push @centers, \@arr; } # sum of weights my @weight_sums; for (my $i = 0; $i < $num_center; $i++) { push @weight_sums, 0; } foreach my $id (keys %{ $weights} ) { for (my $i = 0; $i < $num_center; $i++) { $weight_sums[$i] += $weights->{$id}[$i] ** 2; } } # calc center position foreach my $id (keys %{ $vectors }) { for (my $i = 0; $i < $num_center; $i++) { for (my $j = 0; $j < $veclen; $j++) { $centers[$i][$j] += $weights->{$id}[$i] ** 2 * $vectors->{$id}[$j] / $weight_sums[$i]; } } } return \@centers; } sub fuzzy_cmeans { my ($vectors, $num_center, $num_iter, $m) = @_; # choose initial centers; #my $centers = choose_random_centers($vectors, $num_center); my $centers = choose_smart_centers($vectors, $num_center); my $weights; for (my $i = 0; $i < $num_iter; $i++) { #warn "Loop: $i\n"; $weights = calc_weights($vectors, $centers, $m); $centers = calc_centers($vectors, $weights, $num_center, $m); } return $weights; } sub main { my ($path, $num_center) = @ARGV; if (!$path || !$num_center) { print "Usage: fuzzy_cmeans.pl inputfile num_center\n"; exit 1; } open my $fh, "<$path" or die "cannot open: $path"; my $vectors = read_vectors($fh); my $num_vectors = scalar keys %{ $vectors }; if ($num_center <= 0 || $num_center > $num_vectors) { die 'number of centers must be >0 and <=number_of_vectors'; } my $num_iter = 10; my $m = 2; my $weights = fuzzy_cmeans($vectors, $num_center, $num_iter, $m); print Dumper($weights); } main() __END__
通常のfuzzy c-meansは初めの中心の選択はランダムに行うのですが、上記のソースではK-means++と同じ手法を使用して、なるべく似ていないものを選択するようにしています。
d1 1,2 d2 2,2 d3 4,2 d4 1,1 d5 2,1 d6 4,1
上記のデータを使用して、実際にクラスタリングを実行してみます。
% ./fuzzy_cmeans.pl input.txt 2 $VAR1 = { 'd3' => [ '0.00146446087457261', '0.998535539125427' ], 'd1' => [ '0.99714951936669', '0.00285048063330992' ], 'd2' => [ '0.986032261281484', '0.013967738718516' ], 'd5' => [ '0.986032261281475', '0.0139677387185249' ], 'd4' => [ '0.997149519366688', '0.00285048063331179' ], 'd6' => [ '0.00146446087457242', '0.998535539125428' ] };
整形して出力してないので見づらいのですが、ポイントが高いほどそのクラスタに適合していると考えられるので、d1,d2,d4,d5がクラスタ1、d3,d6がクラスタ2と大体分けることができます。
また各ドキュメントに対し、最も適合したクラスタ以外のクラスタへの所属度も出力できています。今回のデータではどちらか一方のクラスタへの所属度が大幅に大きくなってしまいましたが、データによっては同じくらいずつの所属度を複数のクラスタに持つこともあります。データは必ずしも1つの属性のみを持つのではなく、複数の属性を持つこともあるので、fuzzy c-meansのように複数のクラスタへの所属度を表示できるのはやはり便利ですね。
- 作者: 新納浩幸
- 出版社/メーカー: オーム社
- 発売日: 2007/11
- メディア: 単行本
- 購入: 9人 クリック: 207回
- この商品を含むブログ (29件) を見る
ちなみに僕は「Rで学ぶクラスタ解析」でfuzzy c-meansを調べました。この本は説明も読みやすく、数値での実際の計算例も掲載されているため、非常にわかりやすいです。
ただ、今回のfuzzy c-meansは数式で一部誤植があるような?他にもEMアルゴリズムの数式も多分間違っているところがあって、サポートページとか誤植の一覧ページがあると助かるのですが…。