pygtkとかpycairoを入れるときのメモ

BKChemという分子構造式エディタがあって、pythonで書かれているうえに、スクリーンショットみると構造なんかも綺麗。で、バッチ処理もできるらしいので、ちょいとこれを使ってみるかなとWindowsのバイナリ配布のやつを入れてみた。

だが、ありがちなことに、このバイナリ版はバッチ処理が出来なかったのでソースから入れたけどpycairoを要求してきたのでGTK関連を一式インストール。

GTK+のバイナリはgtk.orgからAll-in-one bundleってのを落としてきて、c:\gtk2とでも名前を変えてパスを切っておく。

で、python bindingはここのを使った。というよりeasy_installではうまく入らなかった。

これでBKChemはすんなりインストールできる。実行ファイルはC:\Python25\Lib\site-packages\bkchem\bkchem.pyだった。

bkchem.py -b batch.py args

という書式でバッチ処理が行われるようになったのでsdf2pngみたいな処理を行わせればcdkの代わりに使えるはず。

pycairo自体も結構面白そうなのでsampleをやってみた。

triangle

他にもチュートリアルがあるのであとで読もう。

pythonでHTTPのgzipデータを読み込む

Fetching PDB files remotely in pure Python codeを見つけて、 数行でpdbのファイルをフェッチできるなんてやるなRCSBとpythonのコンボとか思ったんだけど、コメントに「生でフェッチすんのは環境に悪いでよ」とか書いてあった。

確かにPDB界ではgzip圧縮したpdbファイルをとってきてローカルで展開すんのが昔からのナラワシだよなと、試しにgzとかZをくっつけてwgetしてみるとその通りのファイルがダウンロードできた。をを。

で、おーこれは、数行スクリプトにgzipモジュールかませばいいだけなんじゃなかろうかと思て書き書きしてみたが、これってファイルしかとれないのかな?

zlibならどうだといじってみたけど、

>>> import urllib
>>> import zlib
>>> url = 'http://www.rcsb.org/pdb/files/1ab6.pdb.gz'
>>> text = urllib.urlopen(url).read()
>>> zlib.decompress(text)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
zlib.error: Error -3 while decompressing data: incorrect header check

なんかヘッダがーっておこられた。やっぱgzipでなんとかすべきなのか? ともうちっと追いかけたらGzipFileのほうを使えばよいらしいことに辿りついた。

>>> import gzip, StringIO, urllib
>>> url = 'http://www.rcsb.org/pdb/files/1ab6.pdb.gz'
>>> content = urllib.urlopen(url).read()
>>> sf = StringIO.StringIO(content)
>>> dec = gzip.GzipFile(fileobj=sf)
>>> data = dec.read()

ちと長いがこれで読めた。というわけで、

>>> import urllib, StringIO, gzip
>>> def fetch_pdb(id):
...     url = 'http://www.rcsb.org/pdb/files/%s.pdb.gz' % id
...     content = urllib.urlopen(url).read()
...     sf = StringIO.StringIO(content)
...     return gzip.GzipFile(fileobj=sf).read()

ちょっと長くなったけど帯域にやさしい。

pythonのジェネレータがちょっと分かった

なつたん: Python でSICP4.3 Nondeterministic Computingを読んだら、非決定計算が割りとすんなり理解できた。去年はambわかんねと思ったが、分かってみるとあー非決定計算ねと思えた。

あとジェネレータが便利じゃないかと。それにしても遅延評価みたいなもんかなと思ったら、似たようなシチュエーションでも使えるのね。

再帰とジェネレータ

が面白く読めた。

そういえば非決定計算を使った構造活性相関解析例なんかはあるのだろうか?prologなんかを使った例は見たことあるけど、実用例じゃなかったからな。pythonでかけるなら書いてみたいなぁ。

pylonsにはpaster shell

pylonsで対話的にモデルを操作したい場合にはdevelopment.iniのあるディレクトリで

paster shell

とすればよいらしい。modelっていう変数に色々入っているので

session = model.Session
mol = model.Mol
for m in session.query(mol).filter(mol.MW > 199):
    print m

などとやるべしナ感じ。

jythonプログラミング読了

次の二つの章が特に面白かった。

  • 5 Pythonから見たオブジェクト指向
  • 7 Jythonの応用例

ProductName Jythonプログラミング
西尾 泰和
毎日コミュニケーションズ / ¥ 3,150 ()
在庫あり。

単なるjythonの使い方の本というよりは、jythonのような複数の言語を混ぜ合わせた処理系を通して見えるものを感じ取りましょう的な内容だったように思う。あとその場合に固さやわらかさの長所をうまくのばしてやるようなプログラミングのお作法はこんな感じですよみたいな指針。

perlでもInline系のモジュールを使えば直接他の言語を埋め込んだりできて他言語のライブラリを使えるけど、言語自体が両方の世界に干渉できるともっと楽ができる。

Inline::みたいなやり方でも、jythonみたいなシームレスな言語でも結局複数の言語を理解してないとうまいことやれないので、その分習得コストはかかると思うんだけど、シームレスな言語使うと楽しいのはスケッチするような感覚が得られることなのかなと思ってる。

そういう意味ではGainerも似たような感じかも。

ProductName Make: Technology on Your Time Volume 04
オライリー・ジャパン
オライリージャパン / ¥ 1,575 ()
在庫あり。

あとシームレスな言語を使うっていうことはそれにあわせてライブラリのほうも対応していかないといくのが望ましいのかなとPybel as a generic API for cheminformatics libraries - proof of concept using CDKというエントリを見てふと思った。

chemoinformaticsだとopenbabelっていうc++の大きいライブラリとCDKっていうjavaのライブラリがあって、その上にpython(とruby)がまたがっている感じなのだけど、さらにこれらでいじくったデータを解析にもっていくためにRが必要だったりとかするのでRpyが役に立ったりする。

jythonでopsinを使う

jythonがmacbookにインストールできなかったのでとりあえずlinuxで。

といってもgcjだとエラーを吐くので、sunのjavaをインストール。 ここみて設定。alternativeコマンドを使うとjavaの共存ができるのね。いままで、シンボリックリンクを上書きしてたのでめんどいなーと思ってたけど、これだとらくちん。

jythonをインストールしたら、あとはopsinのjarを落としてきてクラスパスに通す。

>>> import uk.ac.cam.ch.wwmm.opsin as opsin
>>> opsin.NameToStructure().parseToCML("4-iodobenzoic acid").toXML() 

二行でIUPAC名がCMLに。 すばらしい。

jrubyの例もある。

ちなみにjrubyはmacbookにさくっと入って、この例の通りにやれば動いた。

perlで同じ事をやる場合にはInline::Javaを使ってやればいいけど、Javaのライブラリを有効に利用するのはJavaで実装された言語処理系がやっぱ楽だ。

matplotlibでレーダーチャート

元ネタはRadar / Spider Chars

五角形にしたかったのでrule of fiveにPolar Surface Areaを加えておいた。

#!/usr/bin/env python
from matplotlib.projections.polar import PolarAxes
from matplotlib.projections import register_projection
from pylab import *

class RadarAxes(PolarAxes):
    """Class for creating a radar chart (a.k.a. a spider or star chart)        
    http://en.wikipedia.org/wiki/Radar_chart
    """
    name = 'radar'
    # use 1 line segment to connect specified points
    RESOLUTION = 1

    def draw_frame(self, x0, y0, r):
        verts = [(r*cos(t) + x0, r*sin(t) + y0) for t in theta]
        return Polygon(verts, closed=True)

    def set_varlabels(self, labels):
        self.set_thetagrids(theta * 180/pi, labels)

    def get_axes_patch(self):
        x0, y0 = (0.5, 0.5)
        r = 0.5
        return self.draw_frame(x0, y0, r)

if __name__ == '__main__':

    register_projection(RadarAxes)
    N = 5

    theta = 2*pi * linspace(0, 1, N+1)[:-1]
    theta += pi/2
    labels = ['HBA', 'HBD', 'cLogP', 'MWT', 'PSA']
    rule_of_five = [10, 5, 5, 500, 140]
    desc = [12, 3, 3.6, 532, 160]
    desc_rate = [100*desc[i]/float(v) for (i,v) in enumerate(rule_of_five)]

    ax = subplot(111, projection='radar')

    ax.fill(theta, [100]*N)
    ax.fill(theta, desc_rate)

    for patch in ax.patches:
        patch.set_alpha(0.5)

    ax.set_varlabels(labels)
    rgrids((20, 40, 60, 80, 100))

    grid(True)
    show()

rule of five

rule of fiveのようにある範囲内に収まっていること(超えるとリスク)というような指標を表すのにレーダーチャートは適しているんだろうか。つまり充足している事を示すような面積の表現はいいのかなぁ。あと、レンジが負になったりするのでそれもどうかと思う。

再考の余地はあるな。

TODO:多変量がある決まったレンジ内に収まっているかどうかを視覚的に捉えやすい表現手段を探す。

openbabelのAddHydrogens

XYZでreadした分子(CC)もAddHydrogensできない

>>> from openbabel import *
>>> conv = OBConversion()
>>> conv.SetInFormat("xyz")
True
>>> mol = OBMol()
>>> conv.ReadFile(mol, "test.xyz")
True
>>> mol.NumAtoms()
2  
>>> mol.AddHydrogens()
True
>>> mol.NumAtoms()
2

でも分子中のそれぞれの原子を指定して

mol.AddHydrogens(atom)

としてやればきちんと水素が付加されることがわかったので、とりあえずスクリプト中でループをまわせばいいかな

pubchemスクレイピング

mixiでpubchemからSMILESを抜き出すのは?みたいなエントリがたってたのだけど消えちゃったみたい。

なんか、ここは宿題まるなげちゃうわーみたいな厳しいコメントついてたからかな。まぁああいったコメント書くヒトはすっきりするのかもしれんけど、見るほうにとっては情報量ゼロのゴミなんだよなー。まだ、課題まるなげっていう情報のほうが情報量的に有意義。

ちょっと考えたのでここに書いておく。

XMLをBeautifulSoupで

import urllib2,sys
from BeautifulSoup import BeautifulSoup

cid = sys.argv[1]
url = 'http://pubchem.ncbi.nlm.nih.gov/summary/summary.cgi?cid=%s&disopt=DisplayXML' % cid

opener = urllib2.build_opener()
html = opener.open(url).read()
soup = BeautifulSoup(html)

print soup.findAll('pc-infodata_value_sval')[-2].string 

xmlがなんかfindしにくいので配列の要素指定に-2とかやって良くない香り。この程度のだったらxmlじゃなくてSDFを正規表現でいじるな。

import urllib2,sys,re

cid = sys.argv[1]
url = 'http://pubchem.ncbi.nlm.nih.gov/summary/summary.cgi?cid=%s&disopt=DisplaySDF' % cid

html = urllib2.build_opener().open(url).read()
p = re.compile('<PUBCHEM_OPENEYE_CAN_SMILES>\n(.+)\n')
m = p.search(html)
print  m.group(1)

biopythonとかでもいけそうな気がするし、urllib2+openbabelの組み合わせでも良いかもしれない。

オープンソースで始めるゲノム・プロテオーム・メタボローム解析

perl,python,Rでオームな解析をするための本。ツールやデータベースの説明が主なので、ハウツーな感じではなくて、インフォマティクス側からみたバイオロジーとかケミストリーのサービスとかツールのレビューに近い感じかも。

個人的には5章のケモインフォのとこと、9章のRの章が面白かった。