最近の僕の興味は探索フェーズで最も複雑な事象であろうPK/PDとかTK/TDのようなものに移ってしまっている。SBDDを量子化学計算の枠組みで解釈して行くのも楽しいんだけど、どうしても基礎研究寄り過ぎてしまうのと、成果主義的な主張をしづらい(理解されにくい)ので、時間に余裕のある時の知的欲求を満たすための余興的なあつかいになってしまっている。

上司に紹介してもらった本が素晴らしい。

Preclinical Development Handbook: ADME and Biopharmaceutical Properties (Pharmaceutical Development Series)

Wiley-Interscience / ￥ 18,701 ()
通常2～4週間以内に発送

概論多めだけど、1200p超え。

この二冊も良かった。

Pharmacokinetic-Pharmacodynamic Modeling and Simulation
Peter Bonate
Springer / ￥ 8,626 ()
通常1～3週間以内に発送

Pharmacometrics: The Science of Quantitative Pharmacology

Wiley-Interscience / ￥ 18,105 ()
在庫あり。

臨床系のための解析的な本は多いんだけど、それを探索段階のDMPK解析、コンピュテーショナルケミストまたはインフォマティストが使うための本があればいいなぁと。

実際、探索段階だと実験計画とかいい加減だったりとか、データポイントがやたら少なかったりして、ハードル高いうえに、そういう知見はノウハウ的なものなので外に出にくくなっているんだろうなぁと。そこが強いというのは、最終的な強さにつながるのかなぁと、最近強く思う。

ぶっちゃけSBDDなんておまけだよな。QSAR,QSPkR,QSTkRとか分析能力を如何に強くするかが重要なんじゃないかなぁと。

how to makeよりもwhat's to make

そういう気概で来年度は頑張る。

chemoinformaticsな用途でbayonを使ってみる。

データセットはPrimary screen for compounds that inhibit Insulin promoter activity in TRM-6 cells.でアクティブだった1153検体

ダウンロードしてきたsdfをopenbabelでフィンガープリントに変換

babel -imol pc_sample.sdf -ofpt pc_sample.fpt -xh -xfFP2

これをbayonでクラスタリングにするためのTSVに変換

python f2bayon.py pc_sample.fpt > pc_sample.tsv

f2bayon.pyのソース

def hex2bin(fingerprint):
   bf = ""
   h2b = {"0":"0000","1":"0001","2":"0010","3":"0011",
          "4":"0100","5":"0101","6":"0110","7":"0111",
          "8":"1000","9":"1001","a":"1010","b":"1011",
          "c":"1100","d":"1101","e":"1110","f":"1111",
          }

   for l in fingerprint:
       for c in l:
           b = h2b.get(c)
           if b: bf += b
   return bf

def convert(file):
   result = ""
   for data in open(file,"r").read().split("\n>"):
       fp = ""
       for list in data.split("\n")[1:]:
           fp += hex2bin(list)
       result += data.split("\n")[0].split(" ")[0] + " " + fp + "\n"
   return result

if __name__ == "__main__":
   import sys
   file = sys.argv[1]
   c = convert(file)
   for l in c.split("\n")[:-1]:
      id,fp = l.split()
      fps = ""
      for num,bit in enumerate(fp):
         if int(bit) > 0:
            fps += "\t%d\t%s" % (num,bit)
      print id + fps

でbayonで10クラスターに分割

$ time bayon -n 10 pc_sample.tsv > pc_sample.clust

real    0m0.859s
user    0m0.839s
sys 0m0.015s

超速い。

新しい本が出るらしい。

Evaluation of Drug Candidates for Preclinical Development: Pharmacokinetics, Metabolism, Pharmaceutics, and Toxicology (Wiley Series in Drug Discovery and Development)

Wiley / ￥ 9,585 ()
近日発売　予約可

• Cover drug transporters, cytochrome P-450 and drug-drug interactions, plasma protein binding, stability, drug formulation, preclinical safety assessment, toxicology, and toxicokinetics
• Address developability issues that challenge pharma companies, moving beyond isolated experimental results
• Reveal connections between the key scientific areas that are critical for successful drug discovery and development
• Inspire forward-thinking strategies and decision-making processes in preclinical evaluation to maximize the potential of drug candidates to progress through development efficiently and meet the increasing demands of the marketplace

Evaluation of Drug Candidates for Preclinical Development serves as an introductory reference for those new to the pharmaceutical industry and drug discovery in particular. It is especially well suited for scientists and management teams in small- to mid-sized pharmaceutical companies, as well as academic researchers and graduate students concerned with the practical aspects related to the evaluation of drug developability.

python2.5に入れる場合はNumericがsourceforgeになかったり、MxTextToolsの古いバージョンを入れたりとか、いろいろとあれなので、ドキュメント読みましょうってことです。

で、RCSBからgz圧縮されたpdbファイルをとりにいってBioPythonで扱うサンプル

import urllib2, StringIO, gzip
from Bio.PDB.PDBParser import PDBParser

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

if __name__ == "__main__":
    p=PDBParser(PERMISSIVE=1)
    s=p.get_structure("1bgw", fetch_pdb("1bgw"))

    for model in s.get_list():
        for chain in model.get_list():
            for residue in chain.get_list():
                if residue.has_id("CA"):
                    ca_atom=residue["CA"]
                    print ca_atom.get_coord()

ところで、製薬系のドラッグデザインっていうかchemoinformaticsなSoftwareはpythonで拡張できるものが多いんだけれども、その割にはpythonでガリガリ書くよっていうヒトはあんまし見かけたことないんですよね。

化合物の最適化の過程というもをグラフ表現にしてとっておけば、あとあと色々と使いまわせていいんじゃないかと思っている。

で、エッヂは化合物でいいとして、ノードは有向グラフで表したいが、部分構造を内包するかどうかで親子関係を決めるのがいいのか、それとも、時系列で古いほうを親にするのがいいのかそれとももっといいやり方があるのか悩ましい。

可視化は、Cytoscapeが丁度よさそうだったんだけど、Python から Graphviz を使うを見て、NetworkXというものがあるのを知った。

matplotlibで綺麗に出力できるし、これはなかなか楽しそうだ。

そういえば多数の化合物からなるsdfからExactMassってどうやって計算したらいいのかね？と尋ねられて、

うーんそれはね、perlでChemistry::OpenBabelモジュールを使ってだなぁ、ごにょごにょとやればいいんだよ、ほら楽勝だろ。

とか言ってみたんだけど、別にpythonでもいいので書いてみた。

import openbabel as ob

obconv = ob.OBConversion()
obconv.SetInFormat("sdf");
obmol = ob.OBMol()

notatend = obconv.ReadFile(obmol, "test.sdf");
while notatend:
    print obmol.GetExactMass()
    obmol = ob.OBMol()
    notatend = obconv.Read(obmol)

たぶんPlutoだったら

from pluto import *

for mol in Molecule.readfile("test.sdf"):
    print mol.ExactMass

と書けるようになるはず。

でも今夜はOBChemTsfmとかいうクラスがねーとかいう謎の深みにはまりかけ、SWIGやらC++にやられっぱなしでそれどころじゃないのであった。

こんな感じで。 titleとかはmol.title("newtitle")にしたほうがいいのかこのままでいいのか決めかねているのでそのうちちゃんと考える。

エタンをエチレンに変えてみた

>>> from pluto import *
>>> mol = Molecule.readstring("CC\tethane")
>>> mol.bond([1,2]).order(2)
>>> mol.title = "ethylene"
>>> mol.writestring()
'C=C\tethylene\n'

あとはmolecule同士のの結合とatomの追加と削除を実装すればよさそう

PythonでもPerlMolみたいに化学反応を手軽に扱えるツールが欲しいなぁという欲求があったので書いてみることにした。PybelみたいなOpenBabelのラッパーだけど。

そういえば、PyMOLやPyChemという、それっぽい名前は既に使われてるんですよな。

で、アトムにちなんでみた。

PLUTO (1) (ビッグコミックス)
浦沢 直樹
小学館 / ￥ 550 ()
在庫あり。

でも、国内しかわからなんだろうな。海外だとastro boyなんでしょ？まぁいいや。

>>> from pluto import *
>>> mol = Molecule.readstring("c1ccccc1O")
>>> mol.title = "phenol"
>>> mol.writestring()
'c1c(cccc1)O\tphenol\n'

とりあえず、readとwriteまではできた。SMARTSのクラスも作ってあるのだけどアトムばっか返ってきてボンドの返し方がよく分からん。

diels-alder反応ぐらいまで書けるようにするとこまではとっととやる。

なかなか面白い内容だった。

Pybel: a Python wrapper for the OpenBabel cheminformatics toolkit

Pybelボンドとか扱えないような感じなんで、ちょっと複雑なことやろうとすると結局Openbabelのインターフェースからアクセスしちゃうんだよなー。OpenBabelのラッパーとしてPerlMolみたいなの書けばいいんだろうけど。

とか思いながら、OpenBabelのAPIなどみてたら、PerlMolと同等なものを書けそうな気がしてきたので少しやってみることにする。

というわけで、分子量とかformulaとか出せるようにしておいた。

CCC       propane
CC        ethane
c1ccccc1  benzene
CO        methanol

というsmilesを

babel -ismi test.smi -osdf testsdf

で変換したファイルを使った。ExactMassを計算する場合は

>>> from pluto import *
>>> mols = Molecule.readfile("test.sdf")
>>> for mol in mols:
...   print mol.exactmass
...
44.062600256
30.046950192
78.046950192
32.026214748

ついでに分子量と化学式も

>>> from pluto import *
>>> mols = Molecule.readfile("test.sdf")
>>> for mol in mols:
...   print mol.formula, mol.molwt
...
C3H8 44.09562
C2H6 30.06904
C6H6 78.11184
CH4O 32.04186

openbabelにはPerlMolでいうところのcombineみたいなメソッドがないので反応させるところで悩み中。