Pythonの話というより、グラフ理論アルゴリズムの話になりつつあって、ちょっとマニア度が高くなってきましたが、NetworkXのライブラリに少し手を入れて、すべての最短経路を求められるようなメソッドを作って見ました。
重みつき無向グラフの最短経路(shortest path)を求めるアルゴリズムとして、ダイクストラ法(Dijkstra's algorithm)は有名ですが、これを実装したNetworkXのpathsモジュールのメソッドsingle_source_dijkstra(G,source,target=None): は、sourceから、あるtargetまでの最短経路のうち1つしか返してくれません。同じ最短経路長で、いくつかの経路がある場合、それが全部知りたいと思うのは、自然な発想かと。 また、ダイクストラ法は、そのアルゴリズムの本質的な部分で、すべての可能な最短経路を計算しているので、改変はいたって簡単。ソースを貼り付けておきます。 def single_source_dijkstra_all_paths(G,source,target=None): if source==target: return (0, [source]) dist = {} # dictionary of final distances paths = {source:[[source]]} # dictionary of paths seen = {source:0} fringe=[] # use heapq with (distance,label) tuples heapq.heappush(fringe,(0,source)) if not G.is_directed(): G.successors=G.neighbors # if unweighted graph, set the weights to 1 on edges by # introducing a get_edge method # NB: for the weighted graphs (XGraph,XDiGraph), the data # on the edge (returned by get_edge) must be numerical if not hasattr(G,"get_edge"): G.get_edge=lambda x,y:1 while fringe: (d,v)=heapq.heappop(fringe) if v in dist: continue # already searched this node. dist[v] = seen[v] if v == target: break for w in G.successors(v): vw_dist = dist[v] + G.get_edge(v,w) if w in dist: if vw_dist < dist[w]: raise ValueError,\ "Contradictory paths found: negative weights?" elif w not in seen or vw_dist < seen[w]: seen[w] = vw_dist heapq.heappush(fringe,(vw_dist,w)) for _each_path in paths[v]: paths.setdefault(w,[]).append( _each_path + [w] ) elif vw_dist == seen[w]: for _each_path in paths[v]: paths[w].append( _each_path + [w] ) return (dist,paths) pathsの値に、経路のリストを保持するようにして、ループの際に、これまで見つかっている経路長と同じでも、その経路を保持するようにしただけです。 小規模なグラフで動くことは確認したけど、大規模になったときどうかは、確認する必要があるかも。 スポンサーサイト
|
Python使っている人には当たり前なのかもしれないけど、なんか実行を高速化してくれるツールがあるらしい。Psycoがそれです。とりあえず、ダウンロードでインストール。windows用のバイナリがあったので、バージョン1.5を選択。psyco-1.5.win32-py2.4.exe
仕組みは少しこれから詳しくなる必要がありそうだけど、使い方は簡単。 import psyco psyco.full() と、スクリプトの冒頭に入れるだけ。 Cと同じくらい速いとか、実行速度が10倍になるとか書いてるので、ネットワークの中心性(centrality)を計算するnetworkx.centrality.betweenness_centrality( Graph )を使ってプログラムの実行速度を調べてみました。使ったマシンは、PentiumD 3.2GHz(2G RAM) WindowsXPsp2 XPだと、Microsoftが出しているResource Kitに含まれるtimeitというツールを使える。Windows Server 2003 Resource Kit Tools C:\Program Files\Windows Resource Kits\ToolsにPATHを通すの忘れずに。 node数5897、edge数15950のネットワークのBetweennessによるCentralityを求める。確かに速度が2倍くらいに。 とりあえず、時間がかかる計算をするときは試してみる価値あるね。ループがないと意味は薄いようだけど。たとえば、巨大なリストからのランダムサンプリングとかは、ダメな模様。 以下は参考までに、timeitの出力結果。 psyco適用前 Version Number: Windows NT 5.1 (Build 2600) Exit Time: 8:05 pm, Monday, October 16 2006 Elapsed Time: 0:09:15.562 Process Time: 0:09:12.031 System Calls: 1228887 Context Switches: 256903 Page Faults: 815529 Bytes Read: 105270080 Bytes Written: 11654265 Bytes Other: 3072382 psyco適用後 Version Number: Windows NT 5.1 (Build 2600) Exit Time: 8:15 pm, Monday, October 16 2006 Elapsed Time: 0:05:33.968 Process Time: 0:05:33.062 System Calls: 2064202 Context Switches: 465510 Page Faults: 1811956 Bytes Read: 6464132 Bytes Written: 808205 Bytes Other: 1981950 しかし、Psycoって何の略なんだ? |
複雑ネットワークを解析するプログラムをPythonで書きたいと思っていて、Java言うyFilesのように、グラフとしてのネットワークの特性を扱うライブラリが無かったら、自分で書こうと思っていたけど、ものすごくいいのがありました。それが、networkXで、もしかしたらyFilesより高機能かも。しかも、フリー。yFilesはアカデミアでも20万くらいとるから、Pythonやっぱり万歳か?
利用環境はPython2.4.3で、Windowsのバイナリなら、installはいたって簡単。最新版である0.32を入れました。networkx-0.32.win32.exe これだけで、ネットワークを扱うプログラムはすぐに書けるようになりますが、描画できる方がいいかと思い、matplotlibも導入。matplotlib-0.87.6.win32-py2.4.exe >>> import pylab Traceback (most recent call last): File " File "C:\Python24\Lib\site-packages\pylab.py", line 1, in ? from matplotlib.pylab import * File "C:\Python24\Lib\site-packages\matplotlib\pylab.py", line 197, in ? import cm File "C:\Python24\Lib\site-packages\matplotlib\cm.py", line 5, in ? import colors File "C:\Python24\Lib\site-packages\matplotlib\colors.py", line 33, in ? from numerix import array, arange, take, put, Float, Int, where, \ File "C:\Python24\Lib\site-packages\matplotlib\numerix\__init__.py", line 73, in ? import numpy ImportError: No module named numpy 数学ライブラリがが必要なよう。matplotlibのサイトに、 For standard python installations, you will also need to install either Numeric or numarray in addition to the matplotlib installer. matplotlib provides installers for Numeric and numarray users. matplotlib has a numerix setting in the matplotlib rc file (which by default resides in c:\python23\share\matplotlitb\matplotlibrc) and you should make sure this setting corresponds to your preferred array package. と、書いてるので、先日Rとの連携で入れたNumericをそのまま使うことに。Python界の流れ的には、numarrayが主流になりつつあるようで、matplotlibのデフォルトもそれを使うようになっていますが、 C:\Python24\Lib\site-packages\matplotlib\mpl-data の、matplotlibrcを編集。 numerix : numpy # numpy, Numeric or numarray となっているところをコメントアウトし、 numerix : Numeric に変更することで、無事importできるようになります。 networkxのエライところの1つに、何でもnodeにできるところがあります。 たとえば、以下のようにtupleをnodeにできる。そのまま描画するとわかりやすくていいです。 >>> import networkx as nx >>> import pylab >>> g = nx.Graph() >>> g.add_edge( ('node1',3) , ('node2', 8)) >>> g.add_edge( ('node2',8) , ('node3', 9)) >>> g.add_edge( ('node2',8) , ('node4', 5)) >>> nx.draw( g ) >>> pylab.savefig('net.png') ![]() |
いくら、Pythonに充実した標準ライブラリが備わっているといえども、統計処理をしようと思うと、なかなかすぐにできるものではありません。順位和検定をしたいと思い、Pythonのライブラリを少し探したけど、見つからないのでR(http://www.r-project.org/)との連携を模索したらすぐにいいのが見つかりました。
RPy(http://rpy.sourceforge.net/) 現時点で、RC1が発表されていました。2.4系に対応しています。ちなみに試した環境は、2.4.3 RPyは、rpy-1.0-RC1.win32-py2.4.exeを選択。 インストールは、言われるがままに、 Python for Windows extensions(pywin32-210.win32-py2.4.exe) Numerical Python(Numeric-24.2.win32-py2.4.exe ) をそれぞれ、インストール。R.dllのある場所を、windowsの環境変数PATHに追加。これで、OK!と思いきや、実行すると >>> import rpy RHOME= C:\Program Files\R\R-2.3.1 RVERSION= 2.3.1 RVER= 2031 RUSER= C:\Documents and Settings\tsuji Loading the R DLL C:\Program Files\R\R-2.3.1\bin\R.dll .. Done. Loading Rpy version 2031 .. Traceback (most recent call last): File " import rpy File "C:\Python24\Lib\site-packages\rpy.py", line 112, in -toplevel- exec("import _rpy%s as _rpy" % RVER) File " ImportError: No module named _rpy2031 といったようなメッセージが出る。ここでは、エラーを再現するためにわざとやっているけど、結局、rpyがインストールされたフォルダ(たとえば、C:\Python24\Lib\site-packages)にモジュール(_rpy2031.pyd)がないとダメみたい。今回使ったRC1は、デフォルトでR2.3.1に対応しているので、Rのバージョンをこれにあわせておくと、楽かも。 import rpy wc_result = rpy.r.wilcox_test( data1 , data2, paired = rpy.r.TRUE ); wc_result.get( "p.value" ) で、計算結果のP値を得られる。TRUE,FALSEは、そのまま文字列では書けないので注意。 |
| ホーム |
|