submodular toolbox がちゃんと使えるようになった

先週、うまく行っていなかった submodular toolbox は単純に入力ミスだった。
やはり、１日置いてから考え直すと分かることも多い。

これによって、最密部分グラフを視覚的レベルではうまく取り出せることが分かった。
数学的に最適解かどうかは確認していないが、今回は厳密な最適解までは必要ないので、視覚的レベルのみで確認している。
ところで、submodular toolbox だと sfo_min_norm_point の関数で xhat を output としていない。
このままだと、

Nagano, K., Kawahara, Y., & Aihara, K. (2011). Size-constrained Submodular Minimization through Minimum Norm Base. Proceedings of the 28th International Conference on Machine Learning,. Bellevue, WA, USA: ICML. Retrieved from http://www.icml-2011.org/papers/506_icmlpaper.pdf

の方法がつかえないので、sfo_min_norm_point の関数冒頭の function 定義で xhat も output にするように変更した。（もともと xhat は内部計算に使われているので、output に追加しただけ。）


アルゴリズムとしてはこれでいいのだが、実用的なレベルではうまく行っていないことに気がついた。
最密部分グラフとしては連結グラフがほしいのだが、密となる部分が数箇所あると、それぞれが密として抽出されるのである。
それぞれの箇所のうち、最大連結成分のものをまずは最密部分グラフとして出力するように、これから改良して行こうと考えている。

今日の作業内容： submodular 3h + 本読み 2h
今日のランチ：つかさ　生サーモン照り焼き
明日の予測作業時間：5h

submodular toolbox がうまく使えていない様子

最密部分グラフを求める計算で、submodular 関数として計算することができるようなので、やってみているのであるが、現状ではうまくいっていない。

今回は submodular function optimization toolbox を利用している。
http://www.mathworks.com/matlabcentral/fileexchange/20504-submodular-function-optimization/content/sfo/sfo_min_norm_point.m
入力として作ったグラフは、前回と同じく

[0,1]x[0,1] の２次元正方形にランダムに100個のノードをばら撒いて、[0.1,0.5]x[0.2,0.6] に 20 個をばら撒く。
このとき、距離が 0.15 以下となるノード同士を枝でつなぐ。

これでできた距離行列を distanceMatrix としたときに、


    [m,n] = size(distanceMatrix);
    [findI,findJ,findV] = find(distanceMatrix);
    findV2 = ones(length(findV),1);
    distanceMatrix1 = sparse(findI, findJ, findV2, m, n);

    fn = @(nodeSet) -full(sum(sum(distanceMatrix1(nodeSet,nodeSet))));
    F = sfo_fn_wrapper(fn);

として関数ハンドラを生成している。
つまり、すべての枝の重さを１として、nodeSet で与えられたノードの集合に両端が含まれる枝の数をマイナスしたものを最小化する。
これで、
[bestSet] = sfo_min_norm_point(F,1:n);
とすると、最密部分グラフを構成するノード集合が得られる予定であったのだが、実際にやってみるとノードが２個か３個になってしまう。

Toolbox を正しく利用できていないのか、toolbox のアルゴリズムに間違いがあるのか、それをまずは切り分ける必要がある。

今日の作業内容： toolbox 2h + 論文読み 3h
今日のランチ：らく 鶏の照り焼き定食
明日の予測作業時間：4h

mendeley で bibtex 管理をしてみる

この前の SDPARA の原稿のときに bibtex の便利さに今更ながら気がついたが、今更過ぎてひとつひとつを bibtex で入力するのは面倒である。
できれば、論文の PDF を登録すると自動的に bibtex にしてくれるようなものが便利である。
Mac OS X なら bibdesk でできるのであるが、これを Windows か Linux でできる方法をさがしていたところ、mendeley という Web サービスを発見した。
(JabRef も試したが、PDF から bibtex に自動抽出する方法が良く分からず断念した。）

これは、500MB までならフリーで PDF をアップロードできて、それぞれの PDF から自動的に bibtex の情報を抽出できるというものである。
モノはためし、ということで SDPARA の論文で早速やってみたところ、問題なく bibtex に抽出できた。
ただし、リサーチレポートなど学術雑誌になっていないものや、古い学術雑誌でスキャナで取り込んで PDF になっているものなどは、うまく抽出できない部分があったため、手入力で情報を修正する必要があった。

いずれにしても、これから論文 PDF をダウンロードした場合には、ここにストックすることにしようと考えている。
登録した時間で時系列表示もできるので、いつごろにどんなことを勉強していたかの記録の代わりにもなりそうである。

今日の作業内容： mendeley チェック 3h
今日のランチ：味庵　麻婆豆腐
明日の予測作業時間：5h

最密グラフのヒューリスティクス

金曜日に聞いた話で、submodular を使うと与えられたグラフの最密部分グラフがわかる、ということで、これはSNLの計算につかえるのでは？と思って調べている。

まず、いきなり submodular だとハードルが高すぎるので、比較的わかりやすそうなヒューリスティクスで、雰囲気をつかむことにした。

今回実装したヒューリスティクスは、

http://www.cs.umd.edu/~barna/icalp-final.pdf

である。

アルゴリズムとしては、degree の低い node からはずしていって、外すごとに密度を計算する。

この過程で最大の密度を達成したところが、ヒューリスティクスの算出する最密部分グラフである。

このヒューリスティクスでは２倍の精度で最適解を出力する。

つまり、最大密度が 0.9 なら、このヒューリスティクスでは 0.45 から 0.9 の密度となる部分グラフが求まる。

実際にやってみるとわかったが、２倍の精度というのは自分の感覚と相当ずれている。

たとえば、

[0,1]x[0,1] の２次元正方形にランダムに100個のノードをばら撒いて、[0.1,0.5]x[0.2,0.6] に 20 個をばら撒く。

このとき、距離が 0.15 以下となるノード同士を枝でつなぐ。

想定としては [0.1,0.5]x[0.2,0.6]に集中配備だから、ここが最密部分グラフとなるかと思うのだが、乱数によって、ここが出力されたり、グラフ全体を出力したり、ということが、ままある。

あと、乱数の関係で２ノードだけが他から孤立していて、その２ノードは互いにつながっている、という孤島ができることがあるのだが、このとき何故か片方だけが最密部分グラフに取り込まれたりする。

やはり、２倍の精度というのは、SNL 計算には甘いのかもしれない。

いずれにしても、だいたいの雰囲気はわかったので、これで submodular についても検討してみようと思う。

今日の作業内容：最密部分グラフ 3h

今日のランチ：ちゅらさん

明日の予測作業時間：3h

sparsePOP の C++ 版の Windows バイナリを作成

sparsePOP には C++ 版があって、Matlab がなくても実行できるが、Windows バイナリは準備されていなかった。
sparsePOP-C++ で利用しているライブラリのライセンス関係を調べてみたところ、ライセンスについてまとめたものを作成すれば大丈夫ということが分かったので、Windows バイナリが作れるかどうかを試してみた。

結果的に、SDPA 同様に Debian のクロスコンパイルを使えば、sparsePOP-C++-Windows を作ることができる。
ただし、以下の2点には注意が必要と考えている。


(1) metis は 5.0 以降になると idxtype がなくなっていて、
sparsePOP からは利用できない。
(2) metis 4.0.3 は srand48, drand48 という、 Linux にはあるけど
Windows にはない関数をつかっているので、それぞれ srand, rand で
代用できるように libsdrand48.c を作成


最初のコンパイル成功までにはエラーを削除したりなどで５時間程度かかったが、一度できてしまうとコンパイル自体は３分しかかかっていなかった。


現状としては上記のライセンスをまとめたものを準備していないので配布はできないが、これを準備すれば公開できる手はずである。


今日の作業内容：sparsePOP コンパイル 5h
今日のランチ：らく　焼き魚定食
明日の予測作業時間：4h

psファイルのトリミング(bouding box)

TeX で図を入れるとなると eps ファイルとなるが、eps ファイルを Matlab などで生成した場合に軸をTeXの数式で表現したいといったことがある。

このときには、eps ファイルを取り込んだ tex ソースファイルを別途作り、ここで pstricks で数式なり文字列なりを任意の位置に挟み込み、普通にコンパイルした後、dvips で１ページ分だけの ps ファイルにする。
ただ、この方法だと ps ファイルで余白が大きくなり、余白のトリミングが必要である。

余白のトリミング ps ファイルをテキストエディタで開いて bouding box の数字を変更すればできるのだが、できるだけ自動で行いたい。
Ubuntu の場合には、
$ ps2epsi a.ps a.epsi
とすれば、トリミングがされた eps ファイルに変更される。
このあとで
$ mv a.epsi a.eps
とでもすればよい。

ただし、
$ ps2eps a.ps a.eps
とすると、トリミングは行われないようである。

今日の作業内容：SDPA 問い合わせ対応 3h + 論文校正 2h
今日のランチ：しなの 玉子丼セット
明日の予測作業時間：4h

SDPA関係問い合わせが３件

冬休み中に SDPA 関係で３件の問い合わせが来ていた。
ここでは、Q&Aっぽく２件をまとめてみる。

[Q1] SDPA-M を実行すると Matlab ごとクラッシュする。特に Schur が SPARSE になったときにクラッシュするようである。

[A1] この場合、原因は MUMPS の可能性が高い。
MUMPS は一部がスレッドセーフになっていないようで、マルチスレッド計算をするとクラッシュすることがある。
この場合は、
>>> maxNumCompThreads(1);
としておけば、シングルスレッドになるため、（計算時間はかかるが）計算は続行できる。

[Q2] SDPA の callable library を呼んだときと SDPA のバイナリファイルを読んだときだと、数値結果が若干異なることがある。

[A2] リンクしているライブラリ (BLASや LAPACK, GLIBC など）が異なると、数値誤差の影響が出やすい（入力行列が一次従属など）場合では数値結果にも影響が出くることがある。
たとえば、OpenBLAS でも、最新版を git で取ってきて、それまでの版と比較すると、入力行列が一次従属な問題では数値誤差の影響が表面化しやすい。

今日の作業内容：SDPA 対応 3h
今日のランチ：シッダルータ チキンカレー
明日の予測作業時間：4h