tag:blogger.com,1999:blog-57718703012527218212024-02-02T19:44:31.496+09:00TekTekBLOGPythonやC++のメモhiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comBlogger135125tag:blogger.com,1999:blog-5771870301252721821.post-61381777584190780142011-05-01T07:21:00.003+09:002011-05-01T07:54:33.928+09:00pypyで何も考えずにpythonを高速化する<a href="http://python.org/">Python公式</a>のCPythonよりも高速という噂の<a href="http://codespeak.net/pypy/dist/pypy/doc/">pypy</a>が<a href="http://morepypy.blogspot.com/2011/04/pypy-15-released-catching-up.html">Python2.7.1相当の機能を実装した</a>というので、その実行速度を測定しました。<br/><br/><br /><br /><b>■インストール(OSX)</b><br/><br />pypyのインストールは非常に簡単で、<a href="http://pypy.org/download.html#default-with-a-jit-compiler">本家からbzipをダウンロード</a>して解凍したらできあがりです。<br/><br /><pre><br />$ tar xjvf pypy-1.5-osx64.tar.bz2<br />$ cd pypy-c-jit-43780-b590cf6de419-osx64<br />$ ls<br />LICENSE README bin include lib-python lib_pypy site-packages<br /></pre><br />binの中にpypyというコマンドがありますので、それを実行すると起動します。<br/><br /><pre><br />$ bin/pypy<br /><br />Python 2.7.1 (b590cf6de419, Apr 30 2011, 03:30:00)<br />[PyPy 1.5.0-alpha0 with GCC 4.0.1] on darwin<br />Type "help", "copyright", "credits" or "license" for more information.<br />And now for something completely different: ``how to construct the blackhole<br />interpreter: we reuse the tracing one, add lots of ifs and pray''<br />>>>> <br /></pre><br />>が3つではなく4つ出てくるところが、公式CPythonよりも強そうですよね!<br/><br/><br /><br /><b>■実測</b></br><br />実際にPython標準のパフォーマンス測定ツール<a href="http://www.python.jp/doc/2.4/lib/module-timeit.html">timeit</a>で、実行速度を比較してみました。<br/><br />比較のため、OSX上のCPython2.6.1, 2.7,1とpypyを使いました。<br/><br/><br /><br />まずは単純に配列を作る処理。<br/><br /><pre><br />------- CPython 2.6.1<br />>>> import timeit<br />>>> t = timeit.Timer(stmt="[x*x for x in xrange(1000)]")<br />>>> print "%.2f usec/pass" % (1000000 * t.timeit(number=100000)/100000)<br />184.11 usec/pass<br /><br />------- CPython 2.7.1<br />>>> import timeit<br />>>> t = timeit.Timer(stmt="[x*x for x in xrange(1000)]")<br />>>> print "%.2f usec/pass" % (1000000 * t.timeit(number=100000)/100000)<br />153.30 usec/pass<br /><br />------- pypy 1.5<br />>>>> import timeit<br />>>>> t = timeit.Timer(stmt="[x*x for x in xrange(1000)]")<br />>>>> print "%.2f usec/pass" % (1000000 * t.timeit(number=100000)/100000)<br />120.95 usec/pass<br /></pre><br /><br />CPython2.7は2.6よりもだいぶ速いですが、それをさらにpypyが上回る速度です。<br/><br/><br /><br /><br />次に、集合(set)の中の検索<br/><br /><pre><br />------- CPython 2.6.1<br />>>> t = timeit.Timer(setup="a=set([x*x for x in xrange(100)])", stmt="[(x in a) for x in xrange(100)]")<br />>>> print "%.2f usec/pass" % (1000000 * t.timeit(number=100000)/100000)<br />23.41 usec/pass<br /><br />------- CPython 2.7.1<br />>>> t = timeit.Timer(setup="a=set([x*x for x in xrange(100)])", stmt="[(x in a) for x in xrange(100)]")<br />>>> print "%.2f usec/pass" % (1000000 * t.timeit(number=100000)/100000)<br />19.59 usec/pass<br /><br />------- pypy1.5<br />>>>> t = timeit.Timer(setup="a=set([x*x for x in xrange(100)])", stmt="[(x in a) for x in xrange(100)]")<br />>>>> print "%.2f usec/pass" % (1000000 * t.timeit(number=100000)/100000)<br />11.97 usec/pass<br /></pre><br />これも結果は同様に圧倒的にpypyが高速。(2.7の2倍に迫る速度!)<br/><br/><br /><br />こんなに簡単に高速化できるなら、どこかで試しに使ってみたくなりますよね。<br/>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-24425162266203324232011-02-08T20:39:00.007+09:002011-02-08T21:40:24.411+09:00Ajaxを補足するためのPythonによる手軽なマルチスレッド処理Webアプリの開発で、重い処理はAjaxで別途処理して取得するという手段が取れない場合があります。<br /><br />「画面の一部で矢印がぐるぐる回って待たせるのではなく、本当に一瞬で全部表示したい」<br /><br />などの場合です。<br /><br />そんなときにお手軽に使えるテクニックとして、サーバー処理のマルチスレッドによる並列化があります。<br />Pythonのマルチスレッドでは、高度なCPUの並列処理はできませんが(グローバルインタプリタロックと言います)、少なくともI/O待ちの間別の処理を流す程度の事は可能です。[<a href="http://www.python.jp/doc/2.4/api/threads.html">参考</a>]<br />その特徴を生かして、Webアプリのレスポンスを速めることができます。<br /><br /><br />これを使う場面ですが、たとえば、以下のようなものがあります。<br /><br /><a onblur="try {parent.deselectBloggerImageGracefully();} catch(e) {}" href="https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEjbbaVch-OHv69nQiH7H8Ud-rT6W2KBpkemgLVoGdAIWxeEo4QhZAmT_TDRgkZI2Wh13CX0F0Yt5y3UwwN8nRaWHUqPF1qTeeEZ4ov4ya-TVLQ91Uj5pU6QR16Xv1fV2sOKCr9zDnNmtc3y/s1600/biz.png"><img style="display:block; margin:0px auto 10px; text-align:center;cursor:pointer; cursor:hand;width: 320px; height: 227px;" src="https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEjbbaVch-OHv69nQiH7H8Ud-rT6W2KBpkemgLVoGdAIWxeEo4QhZAmT_TDRgkZI2Wh13CX0F0Yt5y3UwwN8nRaWHUqPF1qTeeEZ4ov4ya-TVLQ91Uj5pU6QR16Xv1fV2sOKCr9zDnNmtc3y/s320/biz.png" border="0" alt=""id="BLOGGER_PHOTO_ID_5571296501616839154" /></a><br /><br />この画面はうちの会社で作っている<a href="http://www.patentresult.co.jp/lp-biz.html">特許分析WEBアプリのBiz Cruncher</a>というもので、特許の詳細な内容を表示しているところです。<br /><br /><br />特許の閲覧というのは企業知財部では古来から「手めくり」と呼ばれ、いかに速く大量の特許の内容を把握するかに心血を注ぎます。<br />そのため、図などは一瞬で見れるべきだし、電子的な画面であろうと、一瞬で次々表示されなければストレスを感じるものです。<br /><br /><br />そういう場合は付属情報はAjax化で別リクエストにするにしても、必須の要素だけは泣く泣く1リクエストに押し込んで、極限まで速度を向上させる必要があります。<br /><br />たとえば、上の特許詳細画面では時間のかかる処理として、テキストマイニング、画像の生成、詳細テキスト(長文)の取得、引用情報・経過情報の取得 etc.をそれこそ一瞬で行う必要がありました。<br /><br />改良前はこの表示に2.8秒かかっていたのですが、それでも不満ということで、冒頭のマルチスレッド化を行って0.7秒にまで縮めることができました。<br /><br />具体的には以下の関数 call_async を使います。<br /><br /><pre>import threading<br />def call_async(target, args=None, kwargs=None):<br /> u"""<br /> 関数非同期呼び出しツール。<br /> 関数を別スレッドで呼び出して、その後結果を取得できます。<br /><br /> #使用例<br /> r = call_async(f,args=(5,))<br /> #ここで別の処理を片付けます。<br /> r["thread"].join() #スレッドが終わるのを待ちます。<br /> if "result" in r:<br /> print "RESULT", r["result"]<br /> else:<br /> print "ERROR", str(r["exception"])<br /> """<br /> res = {}<br /> def f(*inargs, **inkwargs):<br /> try:<br /> inargs[-1]["result"] = target(*(inargs[:-1]), **inkwargs)<br /> except Exception, e:<br /> inargs[-1]["exception"] = e<br /> res["thread"] = threading.Thread(target=f, args=(args+(res,)), kwargs=kwargs)<br /> res["thread"].start()<br /> return res<br /></pre><br /><br /><br />この関数を使うと、関数呼び出しを別スレッドにしながら、その返り値を取得することもできます。<br />使用方法は上記関数のコメント内にありますが、以下のように、先に遅い処理を実行開始しておき、平行してその他の処理を済ませたあとで、ゆっくりと結果を取得します。<br /><pre><br /> #使用例 slow1~3という遅い処理の関数がある場合<br /><br /> #まず冒頭で遅い処理を別スレッドで実行開始しておきます。<br /> r1 = call_async(slow1,args=(5,)) #遅い処理1を開始!<br /> r2 = call_async(slow2,args=("ABC","DEF")) #遅い処理2を開始!<br /> r3 = call_async(slow3,args=([1,2,3,4])) #遅い処理3を開始!<br /><br /> #ここでその他の処理を片付けます。<br /><br /> r1["thread"].join() #スレッドが終わるのを待ちます。<br /> r2["thread"].join() #スレッドが終わるのを待ちます。<br /> r3["thread"].join() #スレッドが終わるのを待ちます。<br /><br /> if "result" in r1: #slow1の結果を使う処理<br /> print "RESULT 1", r1["result"]<br /> else:<br /> print "ERROR 1", str(r1["exception"])<br /><br /> if "result" in r2: #slow2の結果を使う処理<br /> print "RESULT 2", r2["result"]<br /> else:<br /> print "ERROR 2", str(r2["exception"])<br /><br /> if "result" in r3: #slow3の結果を使う処理<br /> print "RESULT 3", r3["result"]<br /> else:<br /> print "ERROR 3", str(r3["exception"])<br /></pre><br /><br />こうすることで、I/O待ちの発生する処理をいくつも並列で実行でき、実行時間の効果的な短縮を実現できます。<br /><br />こうした手段としてPython標準のmultiprocessingが検討にあがりますが、あれはより高度な分散処理用ですので、関数の結果を取得しようとすると、とたんに大仰なものになります。<br /><br />また、もちろんですが、この秒数の短縮で効果的なのはJavascriptの最適化もあるわけです。<br />それをした後に、さらにどうしてもサーバー側も高速化しないといけない場合は、このお手軽なマルチスレッド化も考慮するとよいと思います。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-10395298950578073882010-12-04T16:45:00.003+09:002010-12-04T16:50:48.727+09:00PythonのUnicode->UTF-8変換はSJISやEUCへの変換よりもちょっと速い想像していた通りですが、UnicodeをUTF-8に変換するのは、SJISやEUCに変換するよりもちょっと速いことがわかった。<br/><br/><br /><br /><pre><br />>>> import timeit<br />>>> def unicode_henkan(c):<br />... t = timeit.Timer("a.encode('%s')" % c, "a=u'\u3042'*10000")<br />... print t.timeit(number=1000)<br />... <br />>>> unicode_henkan("utf-8")<br />0.075012922287<br />>>> unicode_henkan("sjis")<br />0.19597196579<br />>>> unicode_henkan("eucjp")<br />0.146716833115<br />>>> unicode_henkan("cp932")<br />0.21554517746<br />>>> <br /></pre><br /><br />ちなみに、timeitは与えられた計測対象のプログラムを内部でpythonのプログラムとして解釈する際に、文字コードのことをあまり考えてくれないようなので、<br /><pre><br />t = timeit.Timer("a.encode('%s')" % c, "a=u'あ'*10000")<br /></pre><br />のように書くとUnicodeErrorになります。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-82220296189414168512010-09-14T11:16:00.002+09:002010-09-14T11:26:56.544+09:00関数呼び出しをタプル形式にして遅延させる関数関数呼び出しをタプルの情報にして保存するしくみを作った。<br />これを使えばWEBアプリでリクエストをまたいで関数呼び出しを遅延させることができるようになります。<br />関数呼び出しをタプル形式にして、memcacheなどに保存しておき、次に来たときに本当に呼び出したり。<br /><br /><pre># coding: utf-8<br />import inspect<br />import imp<br /><br />def freeze_func(func, *args, **kwargs):<br /> u"""<br /> 関数呼び出しを保存可能なタプル形式に変換する。<br /> """<br /> m = inspect.getmodule(func)<br /> return m.__name__, func.__name__, args, kwargs<br /><br />def __importer(m, path):<br /> a = imp.find_module(m[0], path)<br /> b = imp.load_module(m[0], *a)<br /> if len(m) > 1:<br /> return __importer(m[1:], b.__path__)<br /> else:<br /> return b<br /><br />def _importer(modname):<br /> modname = modname.split(".")<br /> return __importer(modname, None)<br /><br /><br /><br />def unfreeze_func(info):<br /> u"""<br /> タプル形式で保存された関数呼び出しを関数に復元して実行する。<br /> """<br /> b = _importer(info[0])<br /> f = getattr(b, info[1])<br /> args = info[2]<br /> kwargs = info[3]<br /> return f(*args, **kwargs)<br /><br />if __name__ == "__main__":<br /> import re<br /> info = freeze_func(re.sub, r"[A-Za-z]", "*", "5-3-1 Asakusabashi")<br /> print info # -> ('re', 'sub', ('[A-Za-z]', '*', '5-3-1 Asakusabashi'), {})<br /> print unfreeze_func(info) # -> "5-3-1 ************"<br /></pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-70585617127341527452008-05-24T03:51:00.003+09:002008-05-24T03:55:11.986+09:00Djangoでblogっぽいのを作ったので移動しますDjangoでblogっぽいのを作ったので移動します。RSSも吐いておりますので宜しくお願いします。<br />→<a href="http://www.tektek.in/d/blog/">新TekTekBlog</a><br /><br />移動の理由:<br /> 1. このBloggerでは一切のJavascriptなどが貼れないので不便。アクセス解析の類も一切無理。<br /> 2. データのエクスポート機能もないので、書き溜めても外部に効率的に持ち出せない。<br /> 3. 自分で作って自分で設置するほうが面白い。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-24080402340603529952008-05-23T02:04:00.005+09:002008-05-23T02:18:16.434+09:00はてブのホットエントリをなんとなーく分類するPythonスクリプトサンプル昨日作った、<a href="http://hiroshiykw.blogspot.com/2008/05/wardpython.html">Ward法のPythonバインディング</a>を使って、はてなブックマークのホットエントリをなんとなーく分類するPythonスクリプトのサンプルを書いてみた。<br /><br />サンプルプログラムは例によって、<a href="http://coderepos.org/share/browser/lang/cplusplus/misc/clustercpp/python/wardcluster/hatenasample.py">CodeReposに置きました。</a><br /><br />試しに今の時点(5/23 深夜2:06)のホットエントリを分類した結果は以下のような感じです。<br />なんとなーく、分類されているような雰囲気。<br />もっと各エントリの情報を取得すれば精度は上がるんでしょうが、今の所はホットエントリのRSS一枚から取得できる情報(今の所タイトル中の名詞と、タグの単語)だけを使って分類しています。<br /><br />----- 以下、その結果 -----<br /><div style="line-height:1.2em;">Group [ 50 ] --------------------<br /><a href='http://business.nikkeibp.co.jp/article/tech/20080520/157860/?P=1'>「トレードオフの概念は日本に無いのか」 三菱東京UFJ銀のシステム一本化報道に思う:NBonline(日経ビジネス オンライン)</a><br /><a href='http://business.nikkeibp.co.jp/article/tech/20080520/157860/'>「トレードオフの概念は日本に無いのか」 三菱東京UFJ銀のシステム一本化報道に思う:NBonline(日経ビジネス オンライン)</a><br /><a href='http://plusd.itmedia.co.jp/pcuser/articles/0805/22/news055.html'>なぜか変換できない vs. なぜか変換できる:最近の「MS-IME」は目に余る――よろしい、ならば「ATOK」だ (1/4) - ITmedia +D PC USER</a><br /><a href='http://plusd.itmedia.co.jp/lifestyle/articles/0805/22/news019.html'>今日から始めるデジカメ撮影術:第97回 一眼レフとボケの関係 (1/3) - ITmedia D LifeStyle</a><br /><a href='http://blog.goo.ne.jp/kotoba_mamoru/e/307f39118eb8847c93c6e27b405debd3'>NHK「クローズアップ現代」児ポ法単純所持処罰に絞って特集 - フィギュア萌え族(仮)犯行説問題ブログ版・サブカル叩き報道を追う</a><br /><a href='http://headlines.yahoo.co.jp/hl?a=20080522-00000015-maiall-soci'><ネットカフェ難民>「しんどい」と訴える妊婦まで 100人の実態調査(毎日新聞) - Yahoo!ニュース</a><br /><br />Group [ other ] --------------------<br /><a href='http://builder.japan.zdnet.com/news/story/0,3800079086,20373608,00.htm'>Steve Jobs氏のようにプレゼンテーションをする方法 - builder by ZDNet Japan</a><br /><a href='http://wiredvision.jp/news/200805/2008052121.html'>調査結果に基づいて合成「人々が最も不愉快になる曲/最も聴きたい曲」:試聴可能 | WIRED VISION</a><br /><a href='http://japan.cnet.com/marketing/story/0,3800080523,20373730,00.htm'>グーグル、検索アルゴリズムを少しずつ明らかに:マーケティング - CNET Japan</a><br /><br />Group [ 59 ] --------------------<br /><a href='http://gigazine.net/index.php?/news/comments/20080522_sex_benefit/'>愛し合うってすばらしい、セックスすると得られる10の効能 - GIGAZINE</a><br /><a href='http://www.hatena.ne.jp/info/hatenaclub'>はてなクラブ - はてな</a><br /><a href='http://d.hatena.ne.jp/mojimoji/20080522/p1'>十字軍はバカに勝てるか - モジモジ君の日記。みたいな。</a><br /><a href='http://www.digitaldj.jp/2008/05/22_180000.html'>「宅配寿司 銀のさら」の自虐CM | DIGITAL DJ</a><br /><a href='http://anond.hatelabo.jp/20080522131842'>変な人があまりにも酷すぎる件について</a><br /><a href='http://d.hatena.ne.jp/komachimania/20080521/p2'>高学歴のお嬢様をお持ちの方、教えて下さい - 発狂小町</a><br /><a href='http://med-legend.com/2008/05/%E7%9B%AE%E3%82%92%E9%96%8B%E3%81%84%E3%81%A6%E3%81%84%E3%82%8B%E3%81%AE%E3%81%AB%E9%96%89%E3%81%98%E3%81%A6%E3%81%84%E3%82%8B%E3%81%A8%E6%84%9F%E3%81%98%E3%82%8B%E9%8C%AF%E8%A6%9A%E4%BD%9C%E6%88%90/'>目を開いているのに閉じていると感じる錯覚作成法 | 医学都市伝説</a><br /><br />Group [ 60 ] --------------------<br /><a href='http://d.hatena.ne.jp/fuku33/20080522/1211444127'>ケーキを売ればいいのに - 福耳コラム</a><br /><a href='http://2chart.fc2web.com/ochiaigundam.html'>落合博満ガンダム宣言</a><br /><a href='http://blitznews.blog97.fc2.com/blog-entry-1021.html'>電撃速報!! Amazonを倉庫代わりにしていた転売厨終了のお知らせ</a><br /><a href='http://d.hatena.ne.jp/GilCrows/20080521/p1'>『らき☆すた』を見た非オタの友人が衝撃の一言 - GilCrowsのペネトレイト・トーク</a><br /><br />Group [ 64 ] --------------------<br /><a href='http://gura5.blog120.fc2.com/blog-entry-137.html'>世界のPhotoshopチュートリアル トップ50 - ITクオリティ</a><br /><a href='http://coliss.com/articles/build-websites/operation/design/1112.html'>光源、光線、光跡を表現するPhotoshopのブラシ -2XGU.NET | コリス</a><br /><a href='http://d.hatena.ne.jp/mkusunok/20080521/wocs'>歌舞伎町で遇ったギークなタクシー運転手の編み出した、群衆の英知を活かした馬券購入法 - 雑種路線でいこう</a><br /><br />Group [ 68 ] --------------------<br /><a href='http://gigazine.net/index.php?/news/comments/20080522_miterudake/'>人間関係が苦手で社会に適合できない人向けのDVDソフト「ミテルだけ」が登場 - GIGAZINE</a><br /><a href='http://gigazine.net/index.php?/news/comments/20080522_youtomb/'>YouTubeから削除された動画の墓場「YouTomb」が登場、日本のアニメの削除件数は圧倒的 - GIGAZINE</a><br /><br />Group [ 72 ] --------------------<br /><a href='http://air-users.jp/'>AIR-users.jp - 日本の AIR ユーザのためのハブサイト</a><br /><a href='http://js-users.jp/'>js-users.jp - 日本の JavaScript ユーザのためのハブサイト</a><br /><a href='http://perl-users.jp/'>perl-users.jp - 日本の Perl ユーザのためのハブサイト</a><br /><br />Group [ 75 ] --------------------<br /><a href='http://blog.livedoor.jp/dqnplus/archives/1130403.html'>痛いニュース(ノ∀`):「美少女ゲーム(エロゲ)・アニメは、人間性を破壊し少女殺害事件につながる。販売規制を」…民主党女性議員が請願</a><br /><a href='http://gigazine.net/index.php?/news/comments/20080522_sangiin_game/'>「美少女ゲーム・アニメをする人は心を破壊され、人間性を失っているので規制すべき」と主張するトンデモ請願が参議院に - GIGAZINE</a><br /><a href='http://blog.livedoor.jp/dqnplus/archives/1130723.html'>痛いニュース(ノ∀`):「ジョジョの奇妙な冒険」、中東で「日本が侮辱」「テレビ局爆破しろ」と批判殺到→集英社、原作の第3部など出荷停止</a><br /><br />Group [ 76 ] --------------------<br /><a href='http://wordpress.rauru-block.org/index.php/1656'>人工知能の叛乱</a><br /><a href='http://www.ideaxidea.com/archives/2008/05/gmail_1.html'>Gmailの検索オプションまとめ | IDEA*IDEA</a><br /><a href='http://d.hatena.ne.jp/tanemori/20080520/MostUsedApps'>マックユーザがインストールしているアプリケーション - soundscape out</a><br /><a href='http://okyuu.com/ja/top'>[okyuu.com] ソーシャルITメディア</a><br /><a href='http://www.designwalker.com/2008/05/transparent-inspiration.html'>透過をきれいに使ったウェブデザインいろいろ - DesignWalker</a><br /><a href='http://gihyo.jp/dev/feature/01/firebug'>特集:Firefox 3とFirebugで始めるJavaScript開発|gihyo.jp … 技術評論社</a><br /><a href='http://d.hatena.ne.jp/koto2/20080518/1211070116'>PHP コード最適化 Best Practices 63+ - カタコト日記</a><br /><br />Group [ 77 ] --------------------<br /><a href='http://www.nikkansports.com/general/news/f-gn-tp0-20080522-362998.html'>「ジョジョの奇妙な冒険」が中東で非難 - 社会ニュース : nikkansports.com</a><br /><a href='http://sankei.jp.msn.com/world/mideast/080522/mds0805220813001-n1.htm'>「コーラン読み殺害指示」日本アニメ、中東で非難 - MSN産経ニュース</a><br /><a href='http://www.asahi.com/business/update/0521/NGY200805210017.html'>asahi.com:トヨタ、「カイゼン」に残業代 業務と認定、来月から - ビジネス</a><br /><br />Group [ 78 ] --------------------<br /><a href='http://alfalfa.livedoor.biz/archives/51299093.html'>ヒトラー名言集:アルファルファモザイク</a><br /><a href='http://alfalfa.livedoor.biz/archives/51299206.html'>傷つきやすい人はどうすれば強くなる?:アルファルファモザイク</a><br /><a href='http://alfalfa.livedoor.biz/archives/51296978.html'>ウマーな牛スジカレーの作り方:アルファルファモザイク</a><br /><a href='http://www.itmedia.co.jp/news/articles/0805/21/news101.html'>ねとらぼ:米WIREDに「やる夫」の特大AA掲載 ひろゆき氏インタビューも - ITmedia News</a><br /><a href='http://www.gamenews.ne.jp/archives/2008/05/1500_1.html'>1日500アクセス以上のブログは全ブログの●%:Garbagenews.com</a><br /><a href='http://portal.nifty.com/2008/05/22/a/'>@nifty:デイリーポータルZ:50年前の新聞広告</a><br /><br />Group [ 79 ] --------------------<br /><a href='http://www.moongift.jp/2008/05/usvn/'>MOONGIFT: � ブラウザベースのSubversion管理ツール「USVN」:オープンソースを毎日紹介</a><br /><a href='http://phpspot.org/blog/archives/2008/05/css_58.html'>CSSで背景画像をリサイズするテクニック:phpspot開発日誌</a><br /><a href='http://www.makonako.com/mt/archives/2008/05/post_689.html'>ナムコの成果主義、有能なゲーム開発者が報われない開発現場</a></div>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-29669505239636493352008-05-22T02:00:00.005+09:002008-05-22T02:14:59.782+09:00Ward法のPythonバインディングを作りました。<a href="http://hiroshiykw.blogspot.com/2008/05/cward-g.html">先日から作っていたクラスタリング(Ward法)</a>のPython用バインディングを作りました。Pythonから一応つないで使えるようになりました。<br /><br />ソースコードの取得とインストール。<br /><pre><br />$ svn export http://svn.coderepos.org/share/lang/cplusplus/misc/clustercpp<br />$ cd clustercpp/python/wardcluster<br />$ python setup.py build_ext --swig-cpp<br />$ sudo python setup.py install<br /></pre><br />サンプル実行。<br />引き続きディレクトリclustercpp/python/wardcluster内で以下を実行してください。<br /><pre><br />$ python sample.py<br />1 + 0 => 4 ( 0.244948974278 )<br />3 + 2 => 5 ( 0.412310562562 )<br />5 + 4 => 6 ( 4.80645253132 )<br />$</pre><br />ちなみに、今上で動いたsample.pyは以下のようなスクリプトで、4本のベクトルをクラスタリングしています。<br /><pre>from wardcluster import Matrix, DoubleVec, Ward<br /><br />m = Matrix()<br />m.append(DoubleVec([1.0, 2.0, 3.0]))<br />m.append(DoubleVec([1.1, 2.2, 2.9]))<br />m.append(DoubleVec([0.1, 1.2, 5.0]))<br />m.append(DoubleVec([0.3, 1.0, 5.3]))<br /><br />w = Ward()<br /><br />history = w.cluster(m, 4) #この4ってのは、ベクトルの本数<br /><br />for h in history:<br /> print h.index1, "+" , h.index2, "=>", h.newindex , "(", h.distance ,")"<br /></pre><br /><br />メモ--------<br /> * setuptoolsを使ったオシャレなモジュールにするのは断念。(swigの結果ファイルとの調整が難しかった)<br /> * swigのtypemap?とか、そういうのはむずかしくて分からないので、めんどくさいインターフェース。<br /> * あくまで実験的なものですよ。<br /> * 相変わらずswigを理解できていない。<br /> * 今回ベクトルの集合を表すために内部でstd::vector<std::vector<double> >を使うようにしたので、パフォーマンスが落ちています。(将来もとのstd::vector<double>*を使って渡す方式に戻したいものです。。。)hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-12406072193924768902008-05-20T00:02:00.003+09:002008-05-20T00:10:54.836+09:00C++のWard法->g++のオプションで驚異的な高速化先日から作っているクラスタリングアルゴリズムの一種<a href="http://hiroshiykw.blogspot.com/2008/05/python-c.html">Ward法のプログラム</a>を<a href="http://coderepos.org/share/changeset/12008">今日もいじっていた。</a><br /><br />こまごまとしたことをしたのだが、コンパイルオプションをいじると、恐ろしく速くなった。<br /><br />まず、比較のために、前回の結果。<br /><pre>件数 処理時間(s) そのうち初期行列作成時間<br />--------------------------------------------<br />2352 33.432 6.39</pre><br />そこに、今回以下のコンパイルオプションをつけてみた。(ちなみにOSX10.4&Intelな白MacBookの環境です。)<br /><pre>g++ -O3 -ffast-math -funroll-loops -Wall -fno-common -march=nocona -o ward -I../include ward.cpp</pre><br />その結果、以下のように高速化された。<br /><pre>件数 処理時間(s) そのうち初期行列作成時間<br />--------------------------------------------<br />2352 12.628 1.06</pre><br />シャアではないが、通常の3倍以上の速さである。<br />C++恐るべし。<br />分散処理とか考えるのはやはりC++で肝心の部分を書いて、その後なんだろうなと実感しました。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-66532500686010742182008-05-16T22:52:00.008+09:002008-05-17T17:17:11.865+09:00Python -> C++への移植でどれぐらい速くなったか<a href="http://hiroshiykw.blogspot.com/2008/05/ward_15.html">先日、Pythonで実装したWard法</a>をC++で再実装し、比較してみた。<br /><br />まずは、前回のPython版の処理時間を再掲。<br /><pre>件数 処理時間(s) そのうち初期行列作成時間<br />--------------------------------------------<br />196 0.677 0.334<br />392 2.158 1.356<br />588 4.636 3.073<br />1176 18.071 12.088<br />2352 79.09 49.507</pre><br /><br />これをC++に移植すると、以下の処理時間に短縮される。<br /><pre><br />件数 処理時間(s) そのうち初期行列作成時間<br />--------------------------------------------<br />196 0.216 0.04<br />392 0.842 0.17<br />588 1.627 0.39<br />1176 7.086 1.57<br />2352 33.432 6.39</pre><br />今回は、まずは行列作成部分に、できるだけC++で可能なチューニングを施してある。<br />樹形図作成部分は、あまり最適化していない。<br />標準のstd::mapとかstd::vectorを使用しており、特別な数値計算用ライブラリなどは一切使っていない。<br />また、コンパイル時の様々なコンパイルオプションでの最適化もまだ行っていない。<br /><br />C++(C)のいい所は、メモリの使い方などをダイレクトに想像しながらチューニングが行える所だろう。<br />今回も、3重ループの効率化を行うだけで、ここまで処理時間が短くなった。<br /><br />しかし、作っていて気づいたのは、最適化を意識しないでC++を書くと、Pythonとたいして変わらない処理時間だったことだ。<br />その証拠に、上の表での処理時間の短縮分のほとんどは、今回意識して最適化を図った行列作成部分のものである。<br /><br />次は、実際の樹形図作成部分がネックになっているので、この部分をより速くしていきたいと思う。<br /><br /><a href="http://coderepos.org/share/browser/lang/cplusplus/misc/clustercpp">今回のソースはcodereposに置いておきました。</a>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-77517704862099951852008-05-15T01:00:00.004+09:002008-05-15T01:13:42.225+09:00行列×行列さえ一瞬でできれば、文書のクラスタリングも一瞬でできそうWard法をここ数日試してみて、十分速くなったのだが、ここまでくるといろいろ思う事がある。<br /><br />1. 最終的に今ネックになっているのは、一番最初の行列×行列の計算である。<br />2. その演算であれば、分散処理での高速化が十分研究されている。<br />3. 最初は樹形図をつくるところがネックかと思っていたが、そうでないのは驚きである。<br />4. dict使いまくりの富豪プログラミングはアルゴリズムの試作に便利。<br /><br />ということで、例えば分散処理で行列×行列を一瞬で計算できるようになれば、直前のエントリのクラスタリングアルゴリズムで1000文書であろうとものの数秒で分類できるということになる。<br /><br />このクラスタリングの計算方法をC++で書き直せば、もっともっと速くなるはずで、たとえばweb経由でアクセスがあった瞬間に分類アルゴリズムを走らせることもできるだろう。<br /><br />樹形図をどこでちょんぎるかのさじ加減は別途研究の必要があるが。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-84901524992046054532008-05-15T00:11:00.008+09:002008-05-15T01:20:22.199+09:00Ward法さらにキャッシュを効かせて高速にさらに昨日の続き。<br /><br />行列の各行の最小距離成分をキャッシュするようにしてみた。<br />この改造によって、昨日O(N**2)でネックとなっていた、距離表からの最小値の探索がほぼO(N)になったことになる。<br />さらにパフォーマンスは向上し、以下のようになった。<br /><br /><pre>件数 処理時間(s) そのうち初期行列作成時間<br />--------------------------------------------<br />196 0.677 0.334<br />392 2.158 1.356<br />588 4.636 3.073<br />1176 18.071 12.088<br />2352 79.09 49.507</pre><br />昨日、おとといに比べると、だいぶ速くなった。例えば昨日は588件の処理に19秒かかっていたのが、今日は4.6秒!<br />しかし、計算は速くなったものの、今度は実は最初の行列作成に大半の時間を使われていることがわかる。<br /><br />ここが速くなれば、あるいは、行列の作成時間を気にしなくていい状況とかであれば、十分使い物になるパフォーマンスになってきた。<br />C/C++にはとうていかなわないでしょうけど、Pythonでも工夫次第でいろいろがんばれるものだ。<br /><br />以下、今日の改造版のソース。<br /><pre>#!/usr/local/python25/bin/python<br /># _*_ coding: utf-8 _*_<br /># Copyright (C) 2007 Ayukawa Hiroshi<br />import sys<br />import time<br />from operator import itemgetter<br />from collections import defaultdict<br />from itertools import repeat, groupby<br />from math import sqrt<br /><br />def distmatrix_by_dict(vectors):<br /> """<br /> 上三角部分のみの距離表を作ります(ベクトルがdict型のsparseベクトルの場合)<br /> 返り値 (距離, ベクトル番号1, ベクトル番号2) の列<br /> """<br /> mat = defaultdict(lambda: defaultdict(float))<br /> for i in range(len(vectors)):<br /> dist = []<br /> for j in range(i):<br /> d = 0.0<br /> for k in set(vectors[i].keys()) & set(vectors[j].keys()): #ユークリッド距離<br /> d += (vectors[i].get(k, 0.0) - vectors[j].get(k, 0.0))**2<br /> dist.append((j, sqrt(d)))<br /> mat[i] = dafaultdict(float, dist)<br /> return mat<br /><br />def distmatrix_by_list(vectors):<br /> """<br /> 上三角部分のみの距離表を作ります(ベクトルが通常のリストの場合)<br /> 返り値 (距離, ベクトル番号1, ベクトル番号2) の列<br /> """<br /> mat = defaultdict(lambda: defaultdict(float))<br /> for i in range(1,len(vectors)):<br /> v = vectors[i]<br /> dist = []<br /> for j in range(i):<br /> dist.append((j, sqrt(sum([(x[0]-x[1])**2 for x in zip(v, vectors[j])])))) #ユークリッド距離<br /> mat[i] = defaultdict(float, dist)<br /> return mat<br /><br />def _initcache(mat):<br /> """<br /> 距離表の各行ごとの最小値をキャッシュする配列を初期化します。<br /> """<br /> cache = []<br /> for i in sorted(mat.keys()):<br /> cache.append((i, min(((v, i, j) for j,v in mat[i].iteritems()))))<br /> return dict(cache)<br /><br />class _Rap:<br /> """<br /> パフォーマンス測定用クラス<br /> """<br /> def __init__(self):<br /> self.t = 0<br /> def b(self):<br /> self.start = time.time()<br /> def e(self):<br /> self.t += time.time() - self.start<br /><br />def cluster(vectors, distfunc=distmatrix_by_list):<br /> """<br /> クラスタリング(ウォード法)<br /> 返り値 クラスタ結合履歴配列<br /> """<br /> rap0 = _Rap()<br /> rap00 = _Rap()<br /><br /> #距離表の用意<br /> rap0.b()<br /> mat = distfunc(vectors)<br /> rap0.e()<br /> <br /> #距離表の最小値キャッシュを初期化する<br /> rap00.b()<br /> cache = _initcache(mat)<br /> rap00.e()<br /><br /> size = len(vectors) #全体のサイズ<br /> sizes = list(repeat(1, size)) # 各クラスタのサイズ<br /><br /> #ノードを結合した履歴を保存する配列(この列を最終的に返します)<br /> hist = [dict(<br /> joined=None, #joined: この回で結合したインデックス<br /> dist=0.0, #dist: この回で結合したインデックス間の距離<br /> newid=None, #newid: この回で新設されたインデックス<br /> )]<br /> rap1 = _Rap()<br /> rap2 = _Rap()<br /> nextnodenum = size #次に使うクラスタ番号<br /> kieta = set() #消えた行番号の保存用<br /><br /> #size-1回繰り返して結合していきます。<br /> while len(hist) < size:<br /> #最小距離のものを選ぶ<br /> rap1.b()<br /> dist, q, p = min(cache.values()) # p < q<br /> rap1.e()<br /><br /> #消えるもの集合を更新<br /> kieta.add(p)<br /> kieta.add(q)<br /><br /> #高速化のための事前計算<br /> pairsize = sizes[p]+sizes[q]<br /> r = mat[q][p]<br /> obj = mat[nextnodenum]<br /> sp = sizes[p]<br /> sq = sizes[q]<br /><br /> #距離の計算しなおし<br /> rap2.b()<br /> #キャッシュから、今回消えるものをあらかじめ排除<br /> if p in cache: cache.pop(p)<br /> if q in cache: cache.pop(q)<br /> tmpdist = [] #新ノードからの最小距離を求めるための一時記憶場所<br /> for i in xrange(nextnodenum):<br /> if i in kieta: continue #すでに消えているなら何もしない<br /><br /> #距離の計算 Ward法<br /> if i < p:<br /> obj[i] = ((sp+sizes[i])*mat[p][i]+(sq+sizes[i])*mat[q][i]-sizes[i]*r)/float(pairsize+sizes[i])<br /> elif p < i and i < q:<br /> obj[i] = ((sp+sizes[i])*mat[i][p]+(sq+sizes[i])*mat[q][i]-sizes[i]*r)/float(pairsize+sizes[i])<br /> mat[i].pop(p) #今後使わなくなる成分の消去<br /> else:<br /> obj[i] = ((sp+sizes[i])*mat[i][p]+(sq+sizes[i])*mat[i][q]-sizes[i]*r)/float(pairsize+sizes[i])<br /> mat[i].pop(p) #今後使わなくなる成分の消去<br /> mat[i].pop(q) #今後使わなくなる成分の消去<br /><br /> #キャッシュ更新<br /> if i in cache:<br /> if obj[i] < cache[i][0]: #新ノードとi番目の距離がより近い場合はi行目のキャッシュを置き換えます。<br /> cache[i] = (obj[i], nextnodenum, i)<br /> elif mat[i] and (p in cache[i][1:] or q in cache[i][1:]):<br /> cache[i] = min(((v, i, j) for j,v in mat[i].iteritems()))<br /> if not mat[i]: #i行が空になってたら、iのキャッシュは不要。matからもけずる<br /> cache.pop(i)<br /> mat.pop(i)<br /> elif not mat[i]: mat.pop(i)<br /> tmpdist.append((obj[i], i))<br /> if p in mat: mat.pop(p) #今後使わなくなる成分の消去<br /> if q in mat: mat.pop(q) #今後使わなくなる成分の消去<br /><br /> #新ノードからの最小距離をキャッシュに記録<br /> if tmpdist:<br /> mindist = min(tmpdist)<br /> cache[nextnodenum] = (mindist[0], nextnodenum, mindist[1])<br /><br /> sizes.append(pairsize) #新クラスタのサイズデータを格納<br /> rap2.e()<br /> #履歴情報に今回の結合結果を格納<br /> hist.append(dict(<br /> joined=(p,q), <br /> dist=dist,<br /> newid=nextnodenum,<br /> ))<br /> nextnodenum += 1 #次のクラスタID <br /><br /> print "rap0", rap0.t<br /> print "rap00", rap00.t<br /> print "rap1", rap1.t<br /> print "rap2", rap2.t<br /> <br /> return hist<br /><br />def main(nottest=False):<br /> """<br /> 引数 nottest パフォーマンステスト用フラグ Falseの場合はテストのために、結果出力しません。<br /> """<br /> ids = [] # データラベル格納用<br /> vectors = [] # ベクトル格納用<br /> for line in sys.stdin: #標準入力からベクトルを読む<br /> data = line.strip().split("\t")<br /> ids.append(data[0])<br /> vectors.append(map(float, data[1:]))<br /><br /> #クラスタリング実行<br /> hist = cluster(vectors) <br /><br /> if nottest: #出力処理<br /> for c, h in enumerate(hist):<br /> j = h["joined"]<br /> if not j: continue<br /> i = "(%s+%s)" % (ids[j[0]], ids[j[1]])<br /> ids.append(i)<br /> print c, ":", h["dist"], ":", i<br /><br />if __name__ == "__main__":<br /> import optparse<br /> <br /> parser = optparse.OptionParser(u"""<br /> Ward法によるクラスタリング。<br /> """)<br /> parser.add_option("-t", "--test", dest="test", help=u"パフォーマンステスト実行", default=False, action="store_true")<br /> (options, args) = parser.parse_args()<br /><br /> if options.test:<br /> #パフォーマンステスト<br /> import hotshot, hotshot.stats<br /> prof = hotshot.Profile("cluster.prof")<br /> prof.runcall(main)<br /> prof.close()<br /> stats = hotshot.stats.load("cluster.prof")<br /> stats.strip_dirs()<br /> stats.sort_stats('time', 'calls')<br /> stats.print_stats(20)<br /> else:<br /> #通常実行<br /> main(True)<br /></pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-91486417771944812062008-05-14T02:03:00.004+09:002008-05-14T02:35:00.223+09:00Ward法再実装で少し速くなった<a href="http://hiroshiykw.blogspot.com/2008/05/pythonward.html">PythonでWard法によるクラスタリング</a><br /><a href="http://hiroshiykw.blogspot.com/2008/05/ward.html">Ward法のパフォーマンスと分散処理のための考察</a><br /><br />このシリーズの続き。<br />プログラムを書き直し、真っ当な行列的データを用意したところ、ずいぶん速くなり、588件でも19秒程度となった。(1100件でも2分程度なので、ずいぶんと自分が想定している実用化に近づいたと思う。)<br /><br />しかし、真っ当な行列データにすると、今度は最小距離の探索に9割の時間を取られるようになったため、ここを改造すべきと強く思う。<br /><br />ここまでくれば、常に行列の成分のソート情報をなんらか記憶しておけばいいように思うのだが、そのソート情報の維持にはまたコストがかかるというジレンマ。<br /><br />今は可変長の行列として、pythonのdefaultdict(lambda: defaultdict(float))を使っている。これで定義した簡易行列クラスは非常に便利なのでそんなにシビアでない計算の目的にはよく使っている。可変長&メモリ節約というのがとくにすばらしい。<br /><br />これがC++であったらどうするべきであろうか。<br />こうした細かい計算アルゴリズムの王道としては最初から使用するメモリをまとめて最小限確保しておき、その上でメモリを上書きしながら計算して行くべきか。<br />それとも今回のPythonでの実装にならってSTLを使って map<size_t, map<size_t, double>*>* などとするべきか。ああ、こういう場合shared_ptrにするべきなんだっけ、、、。ここはいつも覚えてないので参考書を読み直さなくては。<br />はたまたboostになんかいいのないかと見ていたら、<a href="http://www.kmonos.net/alang/boost/iterator.html">いつの間にかboostのコンテナ系が増えている!</a>調査しなければ。<br />うーん、前途多難。<br /><br />以下、今日の改造版。<br /><br /><pre>#!/usr/local/python25/bin/python<br /># _*_ coding: utf-8 _*_<br /># Copyright (C) 2007 Ayukawa Hiroshi<br />import sys<br />import time<br />from operator import itemgetter<br />from collections import defaultdict<br />from itertools import repeat, groupby<br />from math import sqrt<br /><br />def distmatrix_by_dict(vectors):<br /> """<br /> 上三角部分のみの距離表を作ります(ベクトルがdict型のsparseベクトルの場合)<br /> 返り値 (距離, ベクトル番号1, ベクトル番号2) の列<br /> """<br /> for i in xrange(len(vectors)):<br /> for j in xrange(i):<br /> dist = 0.0<br /> for k in set(vectors[i].keys()) & set(vectors[j].keys()): #ユークリッド距離<br /> dist += (vectors[i].get(k, 0.0) - vectors[j].get(k, 0.0))**2<br /> yield (sqrt(dist), i, j)<br /><br />def distmatrix_by_list(vectors):<br /> """<br /> 上三角部分のみの距離表を作ります(ベクトルが通常のリストの場合)<br /> 返り値 (距離, ベクトル番号1, ベクトル番号2) の列<br /> """<br /> for i in xrange(len(vectors)):<br /> for j in xrange(i):<br /> dist = sum([(x[0]-x[1])**2 for x in zip(vectors[i], vectors[j])]) #ユークリッド距離<br /> yield (sqrt(dist), i, j)<br /><br />def minmatrix(mat, size): <br /> d = min((min(((v, j) for j,v in mat[i].iteritems())), i) for i in mat.iterkeys() if mat[i])<br /> return (d[0][0], d[1], d[0][1])<br /><br />def cluster(vectors, distfunc=distmatrix_by_list):<br /> """<br /> クラスタリング(ウォード法)<br /> 返り値 クラスタ結合履歴配列<br /> """<br /> #距離表の用意<br /> mat = defaultdict(lambda : defaultdict(float))<br /> for d in distfunc(vectors):<br /> mat[d[1]][d[2]] = d[0] <br /><br /> size = len(vectors) #全体のサイズ<br /> sizes = list(repeat(1, size)) # 各クラスタのサイズ<br /><br /> #ノードを結合した履歴を保存する配列(この列を最終的に返します)<br /> hist = [dict(<br /> joined=None, #joined: この回で結合したインデックス<br /> dist=0.0, #dist: この回で結合したインデックス間の距離<br /> newid=None, #newid: この回で新設されたインデックス<br /> )]<br /><br /> class Rap:<br /> def __init__(self):<br /> self.t = 0<br /> def b(self):<br /> self.start = time.time()<br /> def e(self):<br /> self.t += time.time() - self.start<br /> rap1 = Rap()<br /> rap2 = Rap()<br /> nextnodenum = size #次に使うクラスタ番号<br /> #size-1回繰り返して結合していきます。<br /> kieta = set()<br /> while len(hist) < size:<br /> #最小距離のものを選ぶ<br /> rap1.b()<br /> dist, q, p = minmatrix(mat, nextnodenum) # p < q<br /> rap1.e()<br /><br /> #消えるもの集合を更新<br /> kieta.add(p)<br /> kieta.add(q)<br /><br /> #高速化のための事前計算<br /> pairsize = sizes[p]+sizes[q]<br /> r = mat[q][p]<br /> obj = mat[nextnodenum]<br /> sp = sizes[p]<br /> sq = sizes[q]<br /><br /> #距離の計算しなおし<br /> rap2.b()<br /> for i in xrange(nextnodenum):<br /> if i in kieta: continue #すでに消えているなら何もしない<br /><br /> #距離の計算 Ward法<br /> if i < p:<br /> obj[i] = ((sp+sizes[i])*mat[p][i]+(sq+sizes[i])*mat[q][i]-sizes[i]*r)/float(pairsize+sizes[i])<br /> elif p < i and i < q:<br /> obj[i] = ((sp+sizes[i])*mat[i][p]+(sq+sizes[i])*mat[q][i]-sizes[i]*r)/float(pairsize+sizes[i])<br /> mat[i].pop(p) #今後使わなくなる成分の消去<br /> else:<br /> obj[i] = ((sp+sizes[i])*mat[i][p]+(sq+sizes[i])*mat[i][q]-sizes[i]*r)/float(pairsize+sizes[i])<br /> mat[i].pop(p) #今後使わなくなる成分の消去<br /> mat[i].pop(q) #今後使わなくなる成分の消去<br /> if p in mat: mat.pop(p) #今後使わなくなる成分の消去<br /> if q in mat: mat.pop(q) #今後使わなくなる成分の消去<br /> sizes.append(pairsize) #新クラスタのサイズデータを格納<br /> rap2.e()<br /> #履歴情報に今回の結合結果を格納<br /> hist.append(dict(<br /> joined=(p,q), <br /> dist=dist,<br /> newid=nextnodenum,<br /> ))<br /> nextnodenum += 1 #次のクラスタID <br /> print "rap1", rap1.t<br /> print "rap2", rap2.t<br /> <br /> return hist<br /><br />def main(nottest=False):<br /> """<br /> 引数 nottest パフォーマンステスト用フラグ Falseの場合はテストのために、結果出力しません。<br /> """<br /> ids = [] # データラベル格納用<br /> vectors = [] # ベクトル格納用<br /> for line in sys.stdin: #標準入力からベクトルを読む<br /> data = line.strip().split("\t")<br /> ids.append(data[0])<br /> vectors.append(map(float, data[1:]))<br /><br /> #クラスタリング実行<br /> hist = cluster(vectors) <br /><br /> if nottest: #出力処理<br /> for c, h in enumerate(hist):<br /> j = h["joined"]<br /> if not j: continue<br /> i = "(%s+%s)" % (ids[j[0]], ids[j[1]])<br /> ids.append(i)<br /> print c, ":", h["dist"], ":", i<br /><br />if __name__ == "__main__":<br /> import optparse<br /> <br /> parser = optparse.OptionParser(u"""<br /> Ward法によるクラスタリング。<br /> """)<br /> parser.add_option("-t", "--test", dest="test", help=u"パフォーマンステスト実行", default=False, action="store_true")<br /> (options, args) = parser.parse_args()<br /><br /> if options.test:<br /> #パフォーマンステスト<br /> import hotshot, hotshot.stats<br /> prof = hotshot.Profile("cluster.prof")<br /> prof.runcall(main)<br /> prof.close()<br /> stats = hotshot.stats.load("cluster.prof")<br /> stats.strip_dirs()<br /> stats.sort_stats('time', 'calls')<br /> stats.print_stats(20)<br /> else:<br /> #通常実行<br /> main(True)<br /></pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-38912043196800443442008-05-13T01:23:00.005+09:002008-05-13T01:37:17.577+09:002ちゃんねる Pythonのお勉強 Part26よりParallel Python関係の議論抜粋<a href="http://pc11.2ch.net/test/read.cgi/tech/1209480428/201-300">Pythonのお勉強 Part26</a><br />からの抜粋。なんともタイムリーなことに、Parallel Pythonが話題になっている。<br />この議論で288が指摘している事は、やはり的を得ているのではないだろうか。<br />Python->C/C++という移植で10倍速くなるなんてことは珍しくない。<br />(例えば、C++からMecabを利用した時の速度は驚異的ですらあるが、Pythonから利用した場合はちょっとのんびりとしたものになってしまう。)<br />numpyがうまくハマっているPythonスクリプトとかなら別ですけど。<br /><br />で、以下その議論抜粋です。。。。<br />---------------------<br /><pre>169 :デフォルトの名無しさん:2008/05/06(火) 21:31:30<br /> threadって計算速度の向上に効果ある?<br /> a = b1 + b2の計算を(a,b1,b2はarray)<br /> b1とb2をthreadで計算して,<br /> それぞれ計算終了後に足すというプログラムを作ったんだけど,<br /> 普通にb1とb2をメインスレッドで順番に求めて足した方が早かった...<br /><br /> ちなみにb1とb2の計算はダミーのforループ10000回です. <br /><br />244 :デフォルトの名無しさん:2008/05/10(土) 04:27:34<br /> >>169<br /> Parallel Python というパッケージがあるのを知ったので参考までに紹介しとく。<br /> http://www.parallelpython.com/<br /> スレッドじゃなく複数プロセスで並列実行する仕組みらしい。<br /> Python のみで実装されていて非常にシンプルなパッケージ構成になっている。<br /> マルチコア環境の場合はスレッドモジュール風に使える。<br /> クラスタ環境では各ノードで計算サーバ(ppserver.py)をしておく仕組みになっている。<br /> 応用プログラムはどちらの環境でも同じにできるっぽい。<br /><br /> 手元の共有メモリマシンで付属サンプル(sum_primes.py)を実行したら下のようになった。<br /> 素数の和を計算するプログラムなのでプロセス間通信はほとんどないけど、<br /> 並列化のオーバヘッドはあるわけで、なかなか良好な結果だと思う。<br /><br /> 並列度,実行時間(秒),速度向上率<br /> 1, 17.54, 1.000<br /> 2, 8.781, 1.998<br /> 4, 4.424, 3.965<br /> 8, 2.261, 7.759 <br /><br />245 :169:2008/05/10(土) 13:53:14<br /> >>244<br /> おおお,サンクス!<br /> 早速いじってみまふ。 <br /><br />280 :244:2008/05/10(土) 23:51:28<br /> 引続き Parallel Python をいじってみたんだけども、このシステムはマスタワーカモデルを前提にしているみたいだ。<br /> つまり、マスタが仕事の集合を持っていて、有限個のワーカに1つずつ仕事を割り当てる。仕事を終えたワーカは<br /> マスタから新しい仕事をもらって処理する。これを仕事の集合が空になるまで続ける。<br /><br /> ワーカ間で通信をする機能は提供されていないし、特定のワーカ(たとえば特定のリモートホスト上のppserver.py)に<br /> 特定の仕事を割り当てることもできないっぽい。Parallel Python が有効かどうかは実現したい並列アルゴリズムが<br /> マスタワーカモデルに適合するかどうかに依る希ガス。 <br /><br />283 :デフォルトの名無しさん:2008/05/11(日) 00:40:31<br /> >>280<br /> 具体的には、どういう用途に使えそうなんですか? <br /><br />285 :デフォルトの名無しさん:2008/05/11(日) 01:34:00<br /> >>283<br /> 最も適しているのは並列処理の分野で embarrassingly parallel と呼ばれるカテゴリに属する種々の計算、<br /> すなわち入力データを複数の独立した部分に分割できて、各部分について他から独立して計算が行なえる<br /> タイプの用途に適している。<br /> http://en.wikipedia.org/wiki/Embarrassingly_parallel<br /> Wikipediaに例があがっている。フラクタル計算とか3DCGのレンダリング(例:レイトレーシング)とか力任せの<br /> 暗号解読とか色々ある。<br /><br /> 並列処理って概して複雑で長時間計算し続ける必要のある応用が多い。そういう応用には C とか Fortran を<br /> 使うことがほとんどなんだけど、きちんと動くようになるまでの開発時間が長くてデバッグがたいへんだったりする。<br /> そういう並列プログラムの試作(プロトタイピング)に Python が使えたらいいなーと個人的には思っている。<br /><br />288 :デフォルトの名無しさん:2008/05/11(日) 09:43:58<br /> >>285<br /> PythonでCPU8個使って並列計算するよりCで書いたプログラムの方が速そうだったり。。。<br /><br />289 :デフォルトの名無しさん:2008/05/11(日) 09:57:09<br /> >>288<br /> 何言ってんだ?<br /><br />290 :デフォルトの名無しさん:2008/05/11(日) 10:28:52<br /> >>289<br /> あ?ヤンのかコラ <br /><br /><br />295 :デフォルトの名無しさん:2008/05/11(日) 14:39:25<br /> >>285<br /> 自分も並列化するときのプロトタイプとしてPythonの並列環境を使ってみたいと<br /> 考えているんだけれど、今のところプロトタイプとしての感触はどう?<br /><br /> データ分割の楽な問題ってCとかFORTRANでもOpenMPを注意深く使うだけでも早くなる事が<br /> 多いし、MPIでも慣れれば苦労せずに書けるけれど、そうでない、スレーブ間での通信が必要<br /> とかうまく分割できないとかそういう問題はMPIでやろうとするとデバッグが(;゚д゚)<br /> なことになりがちで、そういう用途にPythonのプロトタイプが役に立つなら素晴らしいと思う<br /> 実はC並列版とそれほど速度とメモリ使用量が変わらないとかなら最高<br /> C並列版より速ければ神<br /><br />296 :デフォルトの名無しさん:2008/05/11(日) 16:53:12<br /> Parallel Python だけど pp.Server.__scheduler() を適当に書き直せばワーカと仕事の対応(割り当て)を<br /> ユーザ側で決められそうな希ガス。問題はワーカ間のプロセス間通信。これが一番面倒なところなんで<br /> 自前で実装となると Parallel Python を使ううまみがほとんどない・・・。<br /><br /> やっぱ MPI の Python バインディングあたりが一番現実的な解かなー。でもなんか気が重いんだなー。<br /> もっと Python らしく lightweight なソリューションがないかなー。<br /><br /> >>288<br /> Python の並列プログラムより C の逐次プログラムの方が速そうってことだよね。<br /> そういうことは十分(多分頻繁に)あると思われ。<br /><br /> >>295<br /> 残念ながらまだ並列プログラムのプロトタイプ用途には使えてなくて感触を得るところまでいってない。<br /> 理由は単純で、Python で手軽に並列プログラミングを実現できる道具を見つけられていないから。<br /> >>209に書いたようにいろいろ試してるんだけどなかなか・・・。<br /><br /> ただ、個人的には Python を並列化のプロトタイピングに使うのは大いに有望だと思っている。<br /> プロトタイピングの場合、欲しいのは並列度やデータ量を上げたときの実行時間の変化であって、<br /> 実行速度が多少遅くて実験に時間がかかるとしても知りたいことは分かるはずだから。<br /> プロセス間通信にソケットを使うとすると、データ量が大きくなれば Python でも C 等と同じぐらいの<br /> 速度が出る(ハズ)。演算量が多い部分に numpy 等の C/Fortran で書かれた数値カーネルを使うことに<br /> すれば、プロトタイピング用途には十分な速度とメモリ使用量が得られるのではと思う。 </pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-20086985321624028252008-05-13T00:01:00.002+09:002008-05-13T00:06:34.724+09:00Ward法のパフォーマンスと分散処理のための考察昨日実装したWard法の処理だが、200件程度なら2秒で終わるが、詳しく調べると、400件で10秒、700件で60秒と、O(n**2)どころか、「もしかしてO(2^n)?!」ぐらいのできであった(T_T)。。。<br /><br />その99%は距離計算直し対象の発見と選り分けの部分であり、本来のWard法の距離計算の部分は微々たるもの。<br /><br />今後はこの選り分け対象の部分でちゃんとインデックスの効いた選り分けができるように変更すべきか。<br /><br />やはり、まずは王道に従って、きちんとした行列データとして距離を持つべきだろう。<br /><br />作り直しだ。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-33316097164967376692008-05-12T04:29:00.007+09:002008-05-12T05:05:26.524+09:00PythonでWard法によるクラスタリングPythonで、Ward法によるクラスタリング(デンドログラム作成まで)を実装してみた。<br />参考にしたのは、岡山大学が公開している<a href="http://case.f7.ems.okayama-u.ac.jp/statedu/hbw2-book/node124.html">この資料</a>。<br />(岡山大学の先生ありがとう!)<br /><br />今回こんな車輪の再発明をしたのは、距離表と呼ばれる行列データを行列のイメージで保存せず、<br /><pre>[(ベクトル1とベクトル2の距離, ベクトル番号1, ベクトル番号2), ...]</pre><br />という形式のリストにして、距離でソートしておき、そのソート結果を維持したままこのリストを分割、更新していくという方針で作ったらどうかと思ったからです。<br /><br />196個の19次元の非スパースなデータのデンドログラム作成に2.5秒だからなかなかのパフォーマンスではないだろうか。<br />CもC++もNumpy系の数値計算ライブラリも使わないでこれだから、いいですねえ。<br /><br /><b>使い方</b><br />読み込むデータファイルをタブ区切りで用意する。<br />1カラム目はデータの識別用ID文字列、2カラム目以降にベクトルの数値データを記入。<br />あとはそのデータを標準入力から流し込むだけ。<br /><b>例</b><br />以下のsample.txtを用意する。<br /><pre>A 1.0 1.1 0.2<br />B 1.1 0.5 0.3<br />C 3.1 0.5 0.2<br />D 4.0 1.1 1.2</pre><br />以下を実行<br /><pre>ayu@~/work/cluster% ./cluster.py < sample.txt <br />1 : (B+A)<br />2 : (D+C)<br />3 : ((D+C)+(B+A))</pre><br />最初にBとAを結合し、次にDとCを結合し、最後に(D+C)と(B+A)を結合していることが分かる。<br /><br /><b>ソースcluster.py</b><br /><pre>#!/usr/local/Python25/bin/python<br /># _*_ coding: utf-8 _*_<br /># Copyright (C) 2007 Ayukawa Hiroshi<br />import sys<br />from collections import defaultdict<br />from itertools import repeat, groupby<br />from math import sqrt<br /><br />def distmatrix_by_dict(vectors):<br /> """<br /> 上三角部分のみの距離表を作ります(ベクトルがdict型のsparseベクトルの場合)<br /> 返り値 (距離, ベクトル番号1, ベクトル番号2) の列<br /> """<br /> for i in xrange(len(vectors)):<br /> for j in xrange(i):<br /> dist = 0.0<br /> for k in set(vectors[i].keys()) & set(vectors[j].keys()): #ユークリッド距離<br /> dist += (vectors[i].get(k, 0.0) - vectors[j].get(k, 0.0))**2<br /> yield (sqrt(dist), i, j)<br /><br />def distmatrix_by_list(vectors):<br /> """<br /> 上三角部分のみの距離表を作ります(ベクトルが通常のリストの場合)<br /> 返り値 (距離, ベクトル番号1, ベクトル番号2) の列<br /> """<br /> for i in xrange(len(vectors)):<br /> for j in xrange(i):<br /> dist = sum([(x[0]-x[1])**2 for x in zip(vectors[i], vectors[j])]) #ユークリッド距離<br /> yield (sqrt(dist), i, j)<br /><br />def cluster(vectors, distfunc=distmatrix_by_list):<br /> """<br /> クラスタリング(ウォード法)<br /> 返り値 クラスタ結合履歴配列<br /> """<br /> #距離表の用意<br /> data = list(distfunc(vectors))<br /> data.sort()<br /><br /> size = len(vectors) #全体のサイズ<br /> sizes = list(repeat(1, size)) # 各クラスタのサイズ<br /><br /> #ノードを結合した履歴を保存する配列(この列を最終的に返します)<br /> hist = [dict(<br /> data=data, #data: まだ結合されていない成分一覧(常に距離でソート済みの状態で格納)<br /> joined=None, #joined: この回で結合したインデックス<br /> dist=0.0 #dist: この回で結合したインデックス間の距離<br /> )]<br /><br /> nextnodenum = size #次に使うクラスタ番号<br /><br /> #size-1回繰り返して結合していきます。<br /> while len(hist) < size:<br /> #最小距離のものを選ぶ<br /> firstpair = min([h["data"][0] for h in hist if h["data"]])<br /> #最小距離のもののインデックスペア<br /> pairindex = set(firstpair[1:])<br /><br /> movedata = [] #新クラスタ形成のために計算し直し対象とするもの格納用<br /> for h in hist: #履歴をさかのぼって、計算しなおしとそうでないものに分ける。<br /> newdata = [] #計算しなおさない居残り組格納用<br /> for d in h["data"]:<br /> if d == firstpair: #今回結合対象のものは除外<br /> continue<br /> intersect = set(d[1:]) & pairindex<br /> if intersect: #今回結合対象とインデックスがかぶるなら、計算しなおし対象へ<br /> movedata.append(d)<br /> else: #そうでないなら居残り組<br /> newdata.append(d)<br /> h["data"] = newdata #居残り組に更新<br /><br /> newdata = [] #計算しなおしデータ格納用<br /><br /> dpq, p, q = firstpair #改めて今回の結合ペアをp,qなどと名付ける<br /> #計算し直しデータを点ごとにグループ化。この点ごとに新距離を計算しなおします。<br /> movedata = groupby(sorted((((set(x[1:])-pairindex).pop(), x) for x in movedata)), key=lambda x: x[0])<br /> for r, d in movedata:<br /> newdist = 0.0 #新距離計算用<br /> for w in d: #ward法計算式にしたがって計算<br /> w = w[1]<br /> newdist += (sizes[w[1]]+sizes[w[2]])*w[0]<br /> newdist = (newdist-sizes[r] * dpq) / float(sizes[r]+sizes[p]+sizes[q])<br /> newdata.append((newdist, nextnodenum, r)) #新距離を格納<br /> newdata.sort() #距離データは常にソートして整理整頓<br /> sizes.append(sizes[p]+sizes[q]) #新クラスタのサイズデータを格納<br /> nextnodenum += 1 #次のクラスタID<br /> #履歴情報に今回の結合結果を格納<br /> hist.append(dict(<br /> data=newdata,<br /> joined=(p,q), <br /> dist=firstpair[0]<br /> ))<br /> return hist<br /><br />def main(nottest=False):<br /> """<br /> 引数 nottest パフォーマンステスト用フラグ Falseの場合はテストのために、結果出力しません。<br /> """<br /> ids = [] # データラベル格納用<br /> vectors = [] # ベクトル格納用<br /> for line in sys.stdin: #標準入力からベクトルを読む<br /> data = line.strip().split("\t")<br /> ids.append(data[0])<br /> vectors.append(map(float, data[1:]))<br /><br /> #クラスタリング実行<br /> hist = cluster(vectors) <br /><br /> if nottest: #出力処理<br /> for c, h in enumerate(hist):<br /> j = h["joined"]<br /> if not j: continue<br /> i = "(%s+%s)" % (ids[j[0]], ids[j[1]])<br /> ids.append(i)<br /> print c, ":", i<br /> <br /><br />if __name__ == "__main__":<br /> import optparse<br /> <br /> parser = optparse.OptionParser(u"""<br /> Ward法によるクラスタリング。<br /> """)<br /> parser.add_option("-t", "--test", dest="test", help=u"パフォーマンステスト実行", default=False, action="store_true")<br /> (options, args) = parser.parse_args()<br /><br /> if options.test:<br /> #パフォーマンステスト<br /> import hotshot, hotshot.stats<br /> prof = hotshot.Profile("cluster.prof")<br /> prof.runcall(main)<br /> prof.close()<br /> stats = hotshot.stats.load("cluster.prof")<br /> stats.strip_dirs()<br /> stats.sort_stats('time', 'calls')<br /> stats.print_stats(20)<br /> else:<br /> #通常実行<br /> main(True)</pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-65930906164101079912008-05-11T15:02:00.007+09:002008-05-11T16:09:48.118+09:00Parallel Pythonで分散処理Erlangでなくて、Pythonで分散処理を書く意義はなんだろうかと考えた。<br />Erlangはたしかに分散処理が得意なんだろうけど、いろいろ調べた感じでは、複雑な数値計算などの分散処理には向いていないというウワサだ。良く知られているTwitterや通信の例のようなシンプルな処理を膨大な量さばくにはいいようだけど。<br /><br />計算を分散で行う場合、本来はGoogleが採用しているように、C++をベースにすべきだろう。<br />(ただし、GoogleではSawzallという独自言語で記述し、それをC++に変換して実行するそうな。)<br /><br />そうなると、「なぜPythonで分散処理?」というのが重要になる訳だが、おそらく以下のようなことだろうか。<br /><ol><li>既存の豊富なモジュール(しかも多くはCで書かれている)を使える。</li><li>C、C++で書かれたルーチンをswig等でPythonに連結し、それを分散させれば、実質C、C++で実行しているようなものだろう。</li><li>簡単で読みやすい。(重要!)</li></ol><br />そこで、まさに自分が今仕事で考えているものの要件でぴったりの問題があった。<br /><pre>問題<br />文の列(数十文程度)の各文どうしの総当たりのdiff(通常の行単位ではなく、<br />単語単位でのdiff)を取り、diffで異なった割合を数値にした総当たり表を<br />リアルタイムで返す。</pre><br />都合のいいことに、Pythonには標準でdiffを取るためのライブラリ<a href="http://www.python.jp/doc/2.4/lib/module-difflib.html">difflib</a>がついているので、こむずかしいdiff計算はそいつに丸投げできる。<br />これが、文の数が10、20程度なら分散させずに普通に解いても速度上問題はない。<br />しかし、これが40、50になってくると、ちょっと「リアルタイム」とは言えなくなってくる。<br />(文の数をnとすると、対称であることと、自分自身とのdiffは計算しないとしても、n*(n-1)/2回diffを取る必要がある。)<br /><br />ということで、これを分散処理させてみよう。<br /><br /><b>準備</b><br /><a href="http://www.parallelpython.com/">Parallel Python</a>を使う。インストールは通常のpython setup.py install 一発で完了。<br />実験に使ったソースdiffdist.pyは本記事末尾に添付した。<br />以下のように実行できる。<br />分散処理で働いてもらいたいマシン群で、Parallel Pythonのサーバーを起動しておく。<br />ppserver.pyコマンドはParallel Pythonをインストールした時に、pythonコマンドと同じディレクトリ(Linux, OSX)又は、c:¥Python25¥Scritpsフォルダ(Windows)に格納されているものです。<br /><pre>$ /PATH/TO/PYTHONBIN/ppserver.py</pre><br />あとは、diffdist.pyのppservers変数(113行あたり)にサーバー群のアドレスを記入しておいて、<br /><pre>$ ./diffdist.py</pre><br />で実行。<br />(※Parallel Python本家ドキュメントに従ってauto discovery機能を使えば、IPアドレスの設定記入を省く事もできる。)<br />今回は、140文(合計サイズ32kb程度)の英文の総当たり表を計算した。(今回このデータは添付しない。)<br />英文データファイルの各行どうしで、単語ごとのdiffを取るというサンプルになっている。<br /><br /><b>結果</b><br />まず普通に分散させないコードだとどれぐらいになるかというと、<br /><pre>ayu@~/work/difftest% time ./diffdist.py <br />Starting pp with 2 workers<br />./diffdist.py 26.64s user 0.36s system 97% cpu 27.705 total</pre><br />ていう感じ。余裕で20秒以上かかります。<br /><br />これを、自分の手元マシン(MacBook白2GHz Intel Core2 Duo & メモリ2G)でppserverを起動して分散処理させると、<br /><pre>ayu@~/work/difftest% time ./diffdist.py <br />Starting pp with 2 workers<br />./diffdist.py 14.22s user 0.62s system 90% cpu 16.412 total</pre><br />ここに、さらに別マシン(DELL SC1420、Intel Xeon 2.8GHz & メモリ1G)を追加すると、、、<br /><pre>ayu@~/work/difftest% time ./diffdist.py <br />Starting pp with 2 workers<br />./diffdist.py 9.62s user 0.59s system 89% cpu 11.394 total</pre><br />順調に速くなっている。。<br />おまけだが、我が家の物置で眠っていたSharpのメビウス君(AMD Duron 897MHz & メモリ256M)を追加すると、、、<br /><pre>ayu@~/work/difftest% time ./diffdist.py <br />Starting pp with 2 workers<br />./diffdist.py 8.03s user 0.58s system 85% cpu 10.100 total</pre><br />ほんのちょっとだけど速くなるね!<br /><br />この調子で、どんどんマシンを追加していけば、この処理はすごく速くなっていくだろう。<br />明日会社で本格的に試してみよう。<br /><br /><b>ソース diffdist.py</b><br /><pre>#!/usr/local/Python25/bin/python<br /># _*_ coding: utf-8 _*_<br /># Copyright (C) 2007 Ayukawa Hiroshi<br />import sys<br />import difflib<br />from itertools import islice<br />from collections import defaultdict<br /><br />import pp<br /><br />def getdiffrate(src1, src2):<br /> """<br /> diffを取って、異なった割合(diff割合)を返す。<br /> """<br /> SPC = " "<br /> same = 0<br /> diff = 0<br /> for x in difflib.ndiff(src1, src2):<br /> if x[0] == SPC:<br /> same += 2 #一致部分カウントはsrc1とsrc2の分で+2<br /> else:<br /> diff += 1<br /> return 100 * same / float(same+diff)<br /><br />def diffrow(src):<br /> """<br /> 一つの文と、対比する複数の文のdiff割合をバッチ的に計算する。<br /> <br /> 引数 src<br /> [(文番号, 原文, [(対比する文1の番号,対比する文1), (対比する文2の番号,対比する文), ..]), ...]<br /> <br /> 返り値<br /> [(文番号, 対比する文の番号, 原文と対比文のdiff割合), ...]<br /> """<br /> ans = []<br /> for i1, s1, osrc in src:<br /> for i2, s2 in osrc:<br /> rate = getdiffrate(s1, s2)<br /> ans.append((i1, i2, rate))<br /> return ans<br /> <br />def d_diffmatrix(job_server, srcs, unit=140):<br /> """<br /> 文一覧を与えて、各文同士の間のdiff割合行列を返す。(分散処理バージョン)<br /><br /> 引数<br /> job_server 分散処理サーバー<br /> srcs 文一覧(文字列の配列)<br /> unit バッチ処理する単位(文の数)<br /><br /> 返り値<br /> 行列データ(dict型データで、res[1][2]のようにdiff割合を格納しています。)<br /> ただし、対角線成分上には、同じ文同士のdiff割合を格納すべきですが、計算を省略するため、<br /> 0.0を格納しています。<br /> """<br /> #各タスク実行に必要なデータを作成する。<br /> task = [] #タスクに与えるデータ格納用<br /> subtask = [] #サブタスクとりわけ用<br /> sublen = 0 #サブタスクの分量計算用<br /> #全タスクを、文の数に従って、サブタスクに分割します。<br /> #サブタスクごとに分散処理させます。<br /> for i1, s1 in enumerate(srcs):<br /> osrc = list(islice(enumerate(srcs), i1))<br /> subtask.append((i1, s1, osrc))<br /> sublen += len(osrc)<br /> if sublen > unit: #文の数がunitを超えたら、一回分としてtaskに格納。<br /> task.append(subtask)<br /> subtask = []<br /> sublen = 0<br /> if subtask:<br /> task.append(subtask)<br /><br /> #分散処理します。<br /> jobs = [] #実行ジョブ格納用<br /> for t in task:<br /> jobs.append(job_server.submit(diffrow,(t, ), (getdiffrate,), ("difflib",))) <br /> <br /> #結果を取得して格納<br /> diffrates = defaultdict(lambda: defaultdict(float))<br /> for job in jobs:<br /> for i1, i2, rate in job():<br /> #行列状に格納<br /> diffrates[i1][i2] = rate<br /> diffrates[i2][i1] = rate<br /> return diffrates<br /><br />def diffmatrix(srcs):<br /> """<br /> 文一覧を与えて、各文同士の間のdiff割合行列を返す。(比較のための通常処理バージョン)<br /><br /> 引数<br /> srcs 文一覧(文字列の配列)<br /><br /> 返り値<br /> 行列データ(dict型データで、res[1][2]のようにdiff割合を格納しています。)<br /> ただし、対角線成分上には、同じ文同士のdiff割合を格納すべきですが、計算を省略するため、<br /> 0.0を格納しています。<br /> """<br /><br /> #全タスクのジェネレーター<br /> def gen():<br /> for i1, s1 in enumerate(srcs):<br /> osrc = list(islice(enumerate(srcs), i1))<br /> yield (i1, s1, osrc)<br /> #処理しながら、結果を行列に格納<br /> diffrates = defaultdict(lambda: defaultdict(float))<br /> for i1, i2, rate in diffrow(gen()):<br /> diffrates[i1][i2] = rate<br /> diffrates[i2][i1] = rate<br /> return diffrates<br /><br />if __name__ == "__main__":<br /> ppservers = ("192.168.1.21","192.168.1.22","192.168.1.23")<br /><br /> if len(sys.argv) > 1:<br /> ncpus = int(sys.argv[1])<br /> job_server = pp.Server(ncpus, ppservers=ppservers)<br /> else:<br /> job_server = pp.Server(ppservers=ppservers)<br /><br /> print "Starting pp with", job_server.get_ncpus(), "workers"<br /><br /> src = file("samplepatent.txt").readlines()<br /> res = d_diffmatrix(job_server, [x.split(" ") for x in src])<br /><br /> #比較用<br /> #res = diffmatrix([x.split(" ") for x in src])<br /></pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-87143534152725551422008-05-10T22:11:00.004+09:002008-05-10T22:27:59.868+09:00Erlangでベクトルの集合をクラスタリング(k-means)<a href="http://www.tektek.in/d/tk/detail/4756150705/">Erlang本</a>を買ったので、k-means法によるクラスタリングをサンプル的に作ってみた。<br />Erlang簡潔でいいですね。<br /><br /><b>使い方</b><br />データファイルをタブ区切りで用意する。<br />1カラム目はデータの識別用ID文字列、2カラム目以降にベクトルの数値データを記入。改行コードは¥nで。<br /><pre>ayu@~/work% /usr/local/erlang_R12B_2/bin/erl<br />Erlang (BEAM) emulator version 5.6.2 [source] [smp:2] [async-threads:0] [kernel-poll:false]<br />Eshell V5.6.2 (abort with ^G)<br />1> c(kmeans).<br />{ok,kmeans}<br />2> Src = kmeans:read_line("sampledata2.txt").<br />[{"1",[1,9,2,4,4,163.5,54,3,1,4,4,2.5,1,1,1,4,3,3,3,2,5]},<br /> {"2",[2,7,1,5,5,154.2,42,2,2,4,3,24.5,1,1,4,3,4,3,5,2,4]},<br /> {"3",[2,8,2,5,3,153.8,44.8,3.5,1,2,5,24,1,1,5,3,2,5,4,3,2]},<br /> .....<br />3> {L, R} = kmeans:splitvectorset(Src).<br />{[{"210",[2,4,2,3,3,156,42,3,2,3,4,2,2,2,3,2,3,5,4,1,3]},<br /> {"205",[1,6,1,3,5,170,54,2.5,2,3,3,1.5,1,1,4,3,4,1,2,3,5]},<br /> {"203",[1,2,1,2,3,163,53,2,2,4,4,24,2,1,3,2,2,1,3,1,3]},<br /> ....<br />4> <br /></pre><br />以上で、LとRに二分割された集合が格納される。<br /><br /><b>ソース</b><br /><pre>-module(kmeans).<br />-import(lists, [map/2, zip/2, sum/1]).<br />-import(math, [sqrt/1]).<br />-export([splitvectorset/1, read_line/1]).<br /><br />%ベクトルのノルム<br />norm(X)-><br /> sqrt(sum(map(fun(Y)->Y*Y end, X))).<br /><br />%ベクトル集合の重心<br />median([],_)-><br /> [];<br />median(X, L)-><br /> First = lists:nth(1,X),<br /> if<br /> First == [] -> [];<br /> true -> [sum(map(fun erlang:hd/1, X))/L | median(map(fun erlang:tl/1, X), L)]<br /> end.<br />median(X)-><br /> median(X, length(X)).<br /><br />%ベクトルの差<br />sub([], _)-><br /> [];<br />sub(X, Y) -><br /> [hd(X)-hd(Y)|sub(tl(X), tl(Y))].<br /><br />%二つのベクトル間の距離<br />distance(X, Y) -><br /> norm(sub(X,Y)).<br /><br />%与えられた二点X,YのどちらにベクトルZは近いか?<br />% Xに近ければ'Right', Yに近ければ'Left'を返す。<br />which(_,_,[])-> <br /> ok;<br />which(X, Y, Z) -><br /> DX = distance(X,Z),<br /> DY = distance(Y,Z),<br /> if<br /> DX > DY -> 'Left';<br /> true -> 'Right'<br /> end.<br /><br />%[{'Right/Left', ベクトル}, ...]というリストを'Right'/'Left'に従って二つの<br />%リストに分けて返す。<br />splitvector([], R, L)-><br /> {R, L};<br />splitvector([{D,Z}|T],R,L) -><br /> case D of<br /> 'Right' -> splitvector(T, [Z|R], L);<br /> 'Left' -> splitvector(T, R, [Z|L])<br /> end.<br /><br />%集合ZSの2元からなる部分集合列を返す<br />pairs([])-><br /> [];<br />pairs([Z|T]) -><br /> lists:append(map(fun(X)->{Z,X} end, T), pairs(T)).<br />%ベクトル集合の距離&ベクトルペアの列を作る<br />pairdistance([])-><br /> [];<br />pairdistance(ZS) -><br /> map(fun({X,Y})->{distance(X,Y), X, Y} end, pairs(ZS)).<br /><br />%ベクトル集合の最大距離を持つ2点を選ぶ<br />maxdistance(ZS)-><br /> Pairs = pairdistance(ZS),<br /> {_, MX, MY} = maxdistance(Pairs, hd(Pairs)),<br /> {MX, MY}. <br />maxdistance([], M)-> M;<br />maxdistance([Z|T], M)-><br /> {D1,_,_} = M,<br /> {D2,_,_} = Z,<br /> if<br /> D2 > D1 -> maxdistance(T, Z);<br /> true -> maxdistance(T, M)<br /> end.<br /> <br /> <br />%二点X,Yのどっちに近いかで、ベクトル集合ZSを二分割する。<br />splitvectorset(X,Y,ZS, Thre)-><br /> {R,L} = splitvector(zip(map(fun({_,W})->which(X,Y,W) end, ZS), ZS), [], []),<br /> IdRmv = fun ({_, V})->V end,<br /> CR = median(map(IdRmv, R)),<br /> CL = median(map(IdRmv, L)),<br /> MinDist = lists:min([distance(X,CR)+distance(Y, CL), distance(Y, CR)+distance(X, CL)]),<br /> if<br /> MinDist > Thre -> <br /> splitvectorset(CR, CL, ZS, Thre);<br /> true -> {R, L}<br /> end.<br /><br />%公開関数その1 <br />%ベクトル集合ZSを二分割する。<br />%集合ZSは、もう一つの公開関数read_lineで読み込んだデータです。<br />splitvectorset(ZS)-><br /> Threshold = 0.1,<br /> case length(ZS) of<br /> 0 -> [[],[]];<br /> 1 -> [ZS|[[]]];<br /> 2 -> [[lists:nth(1, ZS)], [lists:nth(2, ZS)]];<br /> _ -><br /> {X, Y} = maxdistance(map(fun ({_, X})->X end, ZS)),<br /> splitvectorset(X, Y, ZS, Threshold)<br /> end.<br /><br />%分析するデータファイルを読んで{ID, ベクトル}のリストとして返します。<br />% データファイル仕様<br />% データファイルはタブ区切り、先頭カラムはレコードID文字列、次カラム以降が数値データ<br />% 改行コードは\n<br />read_line(File) -><br /> {ok, IoDevice} = file:open(File, read),<br /> ANS = read_line(IoDevice, 1, []),<br /> file:close(IoDevice),<br /> ANS.<br /><br />read_line(IoDevice, LineNumber, Buf) -><br /> case io:get_line(IoDevice, "") of<br /> eof -><br /> Buf;<br /> Line -><br /> Data = string:tokens(lists:delete(10, Line), "\t"), % 10 = "\n"<br /> if<br /> length(Data) == 0 -> ok;<br /> true -> Vec = {hd(Data), map(fun(X)-> parse_num(X) end, tl(Data))},<br /> read_line(IoDevice, LineNumber + 1, lists:append(Buf, [Vec]))<br /> end<br /> end.<br /><br />%数値を表す文字列をinteger又はfloatに変換します。 <br />parse_num(X)-><br /> case string:str(X, ".") of<br /> 0 -> list_to_integer(X);<br /> _ -> list_to_float(X)<br /> end.</pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-41752289225293992532008-05-08T02:05:00.006+09:002008-05-08T02:14:43.065+09:00pythonで英語の同義語などを取得する<a href="http://nltk.org">nltk</a>には<a href="http://ja.wikipedia.org/wiki/WordNet">wordnet</a>が組み込まれているので、nltkひとつをインストールするだけで、同義語、上位語、下位語などを自由に取得できる。<br /><pre>>>> from nltk import wordnet<br />>>> dog = wordnet.N["dog"]<br />>>> dog.synsets()<br />[{noun: dog, domestic_dog, Canis_familiaris}, {noun: frump, dog}, <br />{noun: dog}, {noun: cad, bounder, blackguard, dog, hound, heel}, <br />{noun: frank, frankfurter, hotdog, hot_dog, dog, wiener, <br />wienerwurst, weenie}, {noun: pawl, detent, click, dog}, <br />{noun: andiron, firedog, dog, dog-iron}]</pre><br /><br />英語の検索エンジンの検索クエリ入力画面を作ったときなんかに、このデータでクエリーを膨らませて検索することとかも考えてみると使い手がありそうですね。<br /><br />より詳しい使い方は<a href="http://nltk.org/doc/api/nltk.wordnet.synset-pysrc.html#demo">nltkのデモ</a>でざっと眺める事ができます。<br /><br />nltkのインストール方法は、<a href="http://hiroshiykw.blogspot.com/2008/04/nltk.html">自分の過去記事にあります。</a>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-89929493507886215622008-05-04T03:03:00.010+09:002008-05-04T03:45:32.425+09:00日本語テキストのトピック分割先日から<a href="http://coderepos.org/share">coderepos</a>に置いている自動要約モジュールに、日本語テキストのトピック分割のソースをコミットしました。<br />(-><a href="http://coderepos.org/share/browser/lang/python/yoyaku/yoyaku/engine/TopicSegmentationJP.py">そのソース</a>)<br /><br />このスクリプトでは、与えられた日本語のテキストを、トピックごとに分割する機能を提供しています。<br /><br />基本的には論文"<a href="http://acl.ldc.upenn.edu/A/A00/A00-2004.pdf">Advances in domain independent linear text segmentation</a>"を参考にしています。<br />この論文では<br /> 1. 文ごとにTFベクトルを計算し、<br /> 2. そのベクトル間でcos正規化された内積を計算して文間の類似度を算出、<br /> 3. 近接する文同士の類似度の変化具合を見て、トピックの変わり目を決定。<br />という方式をとっています。<br />ですが、今回の実装では、上記を日本語にも適用するためにさらに以下の改良を加えています。<br /><br /><h1>動詞も使うことにした。</h1><br />日本語のテキスト、特にブログの記事などは、名詞だけでなく、「節約する」とか「怒る」とかの動詞も考慮にいれたほうがよかろうと思い、動詞も使っています。ただし、動詞は活用しますので、mecabでの処理結果の品詞情報の行から動詞の原型を取得して使うことにしました。<br />(-><a href="http://coderepos.org/share/browser/lang/python/yoyaku/yoyaku/engine/TopicSegmentationJP.py?rev=11041#L25">対応するソースの部分</a>)<br /><br /><h1>表記揺れ対策</h1><br />原論文では英語のporter stemmingをしていましたが、日本語では独自の表記揺れ対策をしなければいけません。<br />表記揺れ対策には以下の方針をとりました。<br /><br /><b>1. カタカナ語は地道に辞書を作る。</b><br />カタカナ語の表記揺れ辞書データは本職で半自動生成したことがありますので、近いうちに作ろうと思います。大量の生テキストデータが必要ですので、どこかからテキストを取得せねばなりません。wikipediaデータでもいいかもしれませんが、できればあんなきれいな文章ではないほうがいいので考えます。<br />(->対応するソース<a href="http://coderepos.org/share/browser/lang/python/yoyaku/yoyaku/engine/HyokiyureData.py">HyokiyureData.py</a>)<br /><br /><b>2. 漢字の単語はunigramも生成する</b><br />名詞&未知語の場合のみ、漢字の単語の表記揺れは、その語のunigramもTFカウントすることで対応します。<br />例えば、通常「人」と「人々」は異なる語と認識されてしまいますが、「人々」の方を文字単位で分割し、「人」と「々」それぞれについてもカウントすることにより、「人」と「人々」が互いに多少は類似していると判定できるようになります。<br />(-><a href="http://coderepos.org/share/browser/lang/python/yoyaku/yoyaku/engine/TopicSegmentationJP.py?rev=11041#L58">対応するソース</a>)<br /><br /><h1>語の重みのバリエーション</h1><br />語の重みを以下のように変えて計算しています。<br /> * 名詞、未知語: 1語 = 重み1<br /> * 名詞、未知語unigram: 1語(1文字) = 重み0.3<br /> * 動詞: 1語 = 重み0.5<br />これは、一文でTFベクトルを作っても、たいていはTF=1の語ばかりで、あまりその文の特徴を表しているとも思えませんでしたので、unigramや動詞も考慮に入れた上で、それらの重みの間に差を付けることで、文の間により特色の違いをだせないかという試みです。検証してませんが、まあ、見たかんじ結果がよくなったように思うので、ひとまずこれで行きます。<br />(-><a href="http://coderepos.org/share/browser/lang/python/yoyaku/yoyaku/engine/TopicSegmentationJP.py?rev=11041#L54">対応するソース</a>)<br /><br /><h1>TFベクトルではなく、TFIDF的ベクトルを使用</h1><br />ブログ記事のような中程度の長さの文を扱うことを想定していますので、より文の間の特徴を鋭く認識したほうがいいかと思いました。<br />そこで、文全体でまんべんなく頻出している単語の効果を間引くために語が現れた文の数の逆数をかけることで、ベクトルの重みを調整しています。<br />(-><a href="http://coderepos.org/share/browser/lang/python/yoyaku/yoyaku/engine/TopicSegmentationJP.py?rev=11041#L86">対応するソース</a>)<br /><br /><br />まだまだ先は長いですが、こつこつ作って行きます。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-12751068521139313702008-05-04T02:33:00.002+09:002008-05-04T02:35:28.447+09:00coderepos<a href="http://coderepos.org">coderepos</a>に自分のページを作った。今後ここに開発中のライブラリのマニュアルなどを書いていこうと思います。<br /><br /><a href="http://coderepos.org/share/wiki/Committers/ayu">http://coderepos.org/share/wiki/Committers/ayu</a>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-60550423501091541342008-05-01T00:35:00.002+09:002008-05-01T00:38:59.836+09:00Jython本を注文したJython本を注文しました。<br /><br />まだAmazonのレビューはない見たいだけど、<a href="http://www.tektek.in/d/tk/detail/4839922829/">ブログでは話題になっているようだ。</a><br /><br />随所で最高の評価を受けているようなので、とても楽しみです!<br /><br />今度「Python禁止!JAVAだけで書け!」みたいなプロジェクトに無理矢理入れられそうなので、ちょっとJythonの準備しとくかな。。。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-91410086228134858682008-04-30T23:57:00.002+09:002008-05-01T00:05:48.495+09:00Python用の自動要約モジュールをcodereposにimportしたよけっこう前からのんびり作っていた自動要約モジュールを<a href="http://coderepos.org">Coderepos</a>にimportしました。<br /><br /><a href="http://coderepos.org/share/browser/lang/python/yoyaku">http://coderepos.org/share/browser/lang/python/yoyaku</a><br /><br />この自動要約モジュールはまだあまり精度はよくないです。<br />が、どのぐらい良くないかを知るには、以下のサイトでいちおうサービスとして設置していますので、お試しいただけます。<br /><br /><a href="http://www.civory.com/">自動要約サービス Civory</a><br /><br />今の所は単純なMMRアルゴリズムで文を抜粋するだけのものですので、まあ、あまり、、、。<br /><br />今取り組んでいるのは、英文の自動要約の際に、文をまずトピックごとに分割して、その上で要約をかけようとしています。<br />トピックに分けるアルゴリズムは以下の論文を参考にしています。<br /><a href="http://acl.ldc.upenn.edu/A/A00/A00-2004.pdf">Advances in domain independent linear text segmentation</a><br />いちおうこれの実装はひとまず終わっていて、それに対応するソースは以下になります。<br /><a href="http://coderepos.org/share/browser/lang/python/yoyaku/yoyaku/engine/TopicSegmentationEN.py">TopicSegmentationEN.py</a><br />今後これにバグがないか注意しながらいろんな例で実験を重ねていく予定です。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-59329088250540164682008-04-30T23:15:00.001+09:002008-04-30T23:17:02.090+09:00NLTKでstemming<a href="http://nltk.org">nltk</a>にはステミングの有名なアルゴリズムであるporterアルゴリズムが用意されている。<br /><br /><pre>>>> import nltk<br />>>> nltk.PorterStemmer().stem("application")<br />'applic'<br />>>><br /></pre>hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-5044936037134061512008-04-23T19:45:00.005+09:002008-04-23T19:52:47.696+09:00自分のサイトでのユーザーの視線の流れを調べてみる<a href="http://coliss.com/articles/build-websites/operation/design/1025.html">プレゼンにも最適、ページ内のユーザーの視線をシミュレートする -Feng GUI heatmap</a><br /><br /><blockquote>Feng GUI heatmapは、独自のアルゴリズムに基づいて、ユーザーの視線がどのようにページ内で移動するかとその頻度を表示するヒートマップを作成するサービスです。</blockquote><br />ということなので、自分のサイト<a href="http://www.tektek.in">TekTek</a>でもユーザーの視線をシミュレートしてみました。<br />対象ページ: <a href="http://www.tektek.in/d/tk/detail/487311361X/">書籍「ハイパフォーマンスWebサイト」の詳細ページ</a><br /><a onblur="try {parent.deselectBloggerImageGracefully();} catch(e) {}" href="http://www.tektek.in/d/tk/detail/487311361X/"><img style="cursor:pointer; cursor:hand;" src="https://blogger.googleusercontent.com/img/b/R29vZ2xl/AVvXsEgxJwb70hPNLhujvZepS4YyEKrgf9pc69xf8yCsbDpXrTVI6rofqRsGn3Ut-WLPsxb8W6G5xTl_thVCAlYrO2phVJ_nruxRvjrwt667_SNX2DFZ20znYlPFZWHbBJyNj_MvC9Z1sVAD1cYz/s400/tektekheat.jpg" border="0" alt=""id="BLOGGER_PHOTO_ID_5192390810959172146" /></a><br />ぱっと見、書籍の画像とAmazonへのリンクに視線が集中するようだ。<br />しかし、このページのメインはレビューを眺めてもらうことなので、そちらがもっと読みやすくなるようにすべきかもしれない。<br />なんかよいアイデアはないだろうか。<br /><br />ついでだが、上記書籍がまだ手元に届かない。早くしてよ〜。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.comtag:blogger.com,1999:blog-5771870301252721821.post-72743002129544896352008-04-20T09:20:00.004+09:002008-04-20T09:26:43.562+09:00NLTKで英文の文末判定英文の文末を判定する簡易なルールベースのアルゴリズム。<br /><br /><a href="http://nedbatchelder.com/blog/200804/separating_sentences.html">Separating sentences</a><br /><br />1年前にこの手のアルゴリズムを実装しようとしたが、この問題は非常にやっかいです。<br />たとえばこんな例:<blockquote>CELLULAR COMMUNICATIONS INC. sold 1,550,000 common shares at $21.75 each<br />yesterday, according to lead underwriter L.F. Rothschild & Co.</blockquote><br />"INC."の直後や"$21.75"、"L.F."などのピリオドを文末と認識しては大間違いになるのです。<br /><br />この問題を解決するのに自分が1年前に着目していた論文は以下のもの。<br /><br /><a href="http://www.linguistics.ruhr-uni-bochum.de/%7Estrunk/ks2005FINAL.pdf">Unsupervised Multilingual Sentence Boundary Detection</a><br /><br />この論文では、特に言語を英語だけに限定しない方法を提案しています。<br />大規模な生のテキストデータから得られる統計情報のみで、文末判定を行えます。<br />難点はルールベースなどとは違って、事前の綿密な統計の作成、統計処理後の各種特別処理の実装が面倒くさいこと。<br />昨年途中まで実装したが、Python用の自然言語処理ライブラリ<a href="http://nltk.org/">NLTK</a>でどうも実装予定との情報を見つけて半端でやめていました。<br /><br />それで久しぶりに調べてみたら、このアルゴリズムがすでに<a href="http://nltk.org/">NLTK</a>で実装されて公開されているではないですか!<br /><br />早速インストールして使ってみました。<br /><br /><b>インストール</b><br /><pre>% wget http://prdownloads.sourceforge.net/nltk/nltk-0.9.2.tar.gz<br />% tar zxvf nltk-0.9.2.tar.gz<br />% cd nltk-0.9.2<br />% sudo python setup.py install<br />% cd ..<br />% sudo mkdir /usr/share/nltk<br />% cd /usr/share/nltk<br />% sudo wget http://prdownloads.sourceforge.net/nltk/nltk-data-0.9.2.zip<br />% sudo unzip nltk-data-0.9.2.zip<br />% sudo chmod -R g+r data<br />% export NLTK_DATA=/usr/share/nltk/data<br />% python<br />>>> import nltk<br />>>> nltk.corpus.brown.words()<br />['The', 'Fulton', 'County', 'Grand', 'Jury', 'said', ...]<br />>>><br /></pre><br />ここまでできればインストール完了。<br /><br /><b>生テキストデータを用意</b><br />例えば、Google News(英語版)などからリンクをたどって、ひたすらニュース本文をファイルnews.txtにコピペする。<br />実験でやるんであれば1000行程度で十分でした。<br /><br /><b>実験!!</b><br />まず生テキストを食わせて学習<br /><pre>>>> from nltk.tokenize import PunktSentenceTokenizer<br />>>> p = PunktSentenceTokenizer()<br />>>> fp = file("news.txt")<br />>>> p.train(fp.read())<br /></pre><br /><br />次に実際に判定してみる。<br />判定に使うテキストは次のような少々意地の悪い例<br /><blockquote>The Finland-based company expects a weaker dollar and slower economic growth in the U.S. and parts of Europe to dampen the overall handset market this year. About half of Nokia's (NOK) sales are in dollars or currencies tied to it; a weaker dollar makes imports more expensive.<br /><br />"What spooked us was its outlook for the industry in general," said Rick Franklin, equities analyst at Edward Jones.<br /><br />Nokia reiterated projections that the industry shipments of handsets will grow 10% this year over last. In the first quarter, though, global shipments rose 17%, suggesting a slowdown in the remainder of the year.<br /><br />For the quarter that ended March 31, Nokia earned $1.9 billion (1.2 euros), up 25% from the same quarter last year but short of an expected $2.3 billion. Overall sales rose 28% to $20.1 billion (12.6 billion euros), roughly in line with views.</blockquote><br /><br /><pre>>>> a=p.tokenize("""The Finland-based company expects<br />a weaker dollar and slower economic growth in the U.S. and parts of Europe to<br />dampen the overall handset market this year. About half of Nokia's (NOK) sales<br />are in dollars or currencies tied to it; a weaker dollar makes imports more<br />expensive.<br /><br />"What spooked us was its outlook for the industry in general," said Rick Franklin,<br />equities analyst at Edward Jones.<br /><br />Nokia reiterated projections that the industry shipments of handsets will grow 10%<br />this year over last. In the first quarter, though, global shipments rose 17%,<br />suggesting a slowdown in the remainder of the year.<br /><br />For the quarter that ended March 31, Nokia earned $1.9 billion (1.2 euros), up 25%<br />from the same quarter last year but short of an expected $2.3 billion. Overall<br />sales rose 28% to $20.1 billion (12.6 billion euros), roughly in line with views.""")<br />>>> for x in a:<br />... print x<br />... print "-"*20<br />The Finland-based company expects a weaker dollar and slower economic growth<br />in the U.S. and parts of Europe to dampen the overall handset market this year.<br />--------------------<br />About half of Nokia's (NOK) sales are in dollars or currencies tied to it;<br />a weaker dollar makes imports more expensive.<br />--------------------<br />"What spooked us was its outlook for the industry in general," said Rick Franklin,<br />equities analyst at Edward Jones.<br />--------------------<br />Nokia reiterated projections that the industry shipments of handsets will<br />grow 10% this year over last.<br />--------------------<br />In the first quarter, though, global shipments rose 17%, suggesting a slowdown<br />in the remainder of the year.<br />--------------------<br />For the quarter that ended March 31, Nokia earned $1.9 billion (1.2 euros), up<br />25% from the same quarter last year but short of an expected $2.3 billion.<br />--------------------<br />Overall sales rose 28% to $20.1 billion (12.6 billion euros), roughly in line<br />with views.<br />--------------------<br />>>><br /></pre><br />きちんと、ピリオドで文をわけつつも、"U.S."や、"1.2 euros"などで区切るのは避けていることが分かります。<br /><br />精度を上げるには、もっともっと大量の生テキストを食わせる必要があります。hiroshiykwhttp://www.blogger.com/profile/01527114074913563640noreply@blogger.com