ケムインフォマティクスに虚空投げ

JOELibのDescriptorsは、"Native values","Atom properties","Fingerprints","Transformations"の4つのカテゴリーに分けられています。"Native values"は、Drug-likenessフィルタやQSARで利用可能な基本的なDescriptorsの集まりですので、今回は、この"Native values"の取得方法についてメモしたいと思います。

"Native values"に属するDescriptor名を知るには、DescriptorHelperクラスのインスタンスを生成し、getNativeDescs()メソッドを実行します。

　Vector nativeDes = DescriptorHelper.instance().getNativeDescs();
　String desName;
　for (int i = 0; i < nativeDes.size(); i++)
　{
　　desName = (String)nativeDes.get(i);
　　System.out.println(desName);
　}

ここで得られるdesNameがDescriptorの値を取得する際に必要となります。

実際に、Descriptorの値を取得するには、DescriptorFactoryクラスのgetDescriptor(String name)メソッドでDescriptorを取得し、calculate(JOEMol mol)メソッドによりDescriptorの値の計算を行います。

　Descriptor　descHBA = DescriptorFactory.getDescriptor("Number_of_HBA_2");
　DescResult　results = descHBA.calculate(mol);

Descriptorの値をint型にするには次のように記述します。

　int hba =((IntResult)results).getInt();

logPのようにdouble型のDescriptorの値をもつものは、次のように記述します。
　
　double logP =((DoubleResult)results).getDouble();

CDKと同様に、JOELibでもDescriptorの値を取得すれば、簡単にrule of fiveフィルタが作成でき便利です。
ちなみにJOELibでのlogPの計算は、以下の論文の手法が使われています。

　J. Chem. Inf. Comput. Sci., 1999, 39, 868-873.



人気ブログランキング(クリックして応援してね)

JOELibを用いて分子をPOVRay形式で出力してみます。POVRayは、レイトレーシング法により高品質な3次元画像を生成するフリーな3Dグラフィックスツールです。JOELibはPOVRayのシーン記述ファイルである「POVファイル」を出力することができるため、センスがあれば(私にはありませんが)、非常にかっこいい分子画像を作成することができます。

まずは、POVRayクラスのインスタンスを生成し、initWriter(OutputStream os)メッソドを実行します。

　joelib.io.types.POVRay ray = new joelib.io.types.POVRay();
　ray.initWriter(output);

次に分子の表示形式のセットを行います。
　
　ray.setOutputType(joelib.io.types.POVRay.SPHERE);

ここでは、SPHERE表示を選択しています。その他の表示形式は、POVRayクラスのField Summaryをご覧ください。
POVファイルの出力はwrite(JOEMol mol)メソッドで行います。
　
　ray.write(mol);

出力ファイルをPOVRayにより描画してみます。

SPHERE形式


Ball and Stick形式

　
POVRayで描画させると真っ先にやりたくなる操作は、縮小・拡大/回転だと思います。POVファイルの一番下にある以下のような部分にrotate<>及びscale<>を書き加えることにより実現できます。

object { Molecule1
　　　rotate<90,0,0>
　　　scale<0.8,0.8,0.8>
}



人気ブログランキング(クリックして応援してね)

現在、Protein Data Bank(PDB)には、約37000の立体構造が登録されているようです。年ごとの登録数をみると指数関数的に増加していますね。PDBの成熟に伴い、多くの二次データベースが派生していると思いますが、私自身は低分子とタンパク質の複合体構造に興味があるので、これに関連しそうなデータベースをメモしておきます。

MSDChem：
　PDB中のLigandを部分構造,Fragments,Fingerprintなどで検索ができる。

The PDBbind Database：
　Protein-ligand複合体の実測binding affinityデータが収められている。

Computed Ligand Binding Energy(CLiBE) Database：
　コンピュータで計算したprotein-ligand複合体の相互作用エネルギーが収められている。

HIV Protease Database:
　HIV-1,HIV-2,SIVの立体構造が整理されて収められている。

Protein Ligand Database (PLD)：
　Protein-ligand複合体のbinding affinityデータが実測/計算値共に収められている。




人気ブログランキング(クリックして応援してね)

JOElibでは、JOEMolクラスのaddHydrogens()メソッドにより、H原子を付加することができます。H原子を付加する際、基本的には、以下のコードにより、H原子の個数を決定しています(実際のJOELib中のコードです)。

　hcount = atom.getImplicitValence() - atom.getValence();

JOELibの便利な機能として、pHを考慮したH原子の付加ができます。具体的には、joelib/data/plain/phmodel.txtに従っています。例えば、carboxylic acidの場合、以下のように記述されています。

TRANSFORM O=CO[#1:1] >> O=CO
TRANSFORM O=C[OQ1-0:1] >> O=C[O-:1]

この記述によれば、COOHはCOO-に変換されるため、負に荷電しており、H原子は付加されません。したがって、phmodel.txtをカスタマイズすることにより、H原子の付加を自由にコントロールできます。記述方法は、Reaction SMARTSだと思います。

もう１つ便利な機能として、polar hydrogensのみの付加があります。例えば、水素結合のパターンを目視で認識する際に便利ですし、Amberにもこれに対応した力場パラメータが存在するので、分子モデリングでも利用することができます。

それでは、JOEMolクラスのaddHydrogens(boolean polaronly, boolean correctForPH)メソッドを用いて、上述した内容をアラニンを用いて確認したいと思います。第一引数が、polar hydrogensのみのH原子付加を行うか、第二引数は、phmodel.txtを用いたH原子付加のコントロールを行うかを示しております。


上図中のtrue,falseは第一引数、第二引数の順番になっております。ちなみに、addHydrogens()メソッドでは、false,trueとなっております。



人気ブログランキング(クリックして応援してね)

昨日本屋でなんとなく本を見ていると多湖輝さんの「頭の体操」シリーズが平積みで置いてあった。小中学生の頃、このシリーズが大好きで、難解な問題をクリアーできたときの爽快感がたまらなかった。それまで知っていた“なぞなぞ”や“パズル”とは、明らかに一線をきす問題に、これを考えつく多湖さんは、当時、私の尊敬する大人のNo.1であった。もう1つこのシリーズの印象に残っているのは、本の表紙で、具体的には思い出せなかったのでネットで探してみると、ココにあった。このなにか特異な絵柄も惹きつけられる要素だったのだと思う。第10集あたりまでは読んだ記憶があるが、だんだん興味が他に移っていき読まなくなった。

あれから20年程たったが、多湖さんにもう一度挑戦してみようと思う。DS脳トレで鍛えているおじさんは、20年たっても脳年齢は21才なはず(笑)なので8割以上はクリアーできるだろう(といいんなあ)と信じている。






人気ブログランキング(クリックして応援してね)

JOELibを使って簡易babelを作成してみます。元々babel自体がOELibで作られていますので、JOELibでbabelっぽいプログラムを作るのは容易だろうと想像できます。

JOELibではSimpleReaderクラスを用いてファイルの入力を行います。
SimpleReaderクラスのインスタンスの生成時に、ファイル形式の情報が必要となるため、IOTypeHolderクラスを用いて、その情報を取得します。

　IOType inType = IOTypeHolder.instance().filenameToType(inputf);

ここでは、インスタンスの数を1つに保つために、GoFのSingletonパターンが使われていると思います。

　sreader = new SimpleReader(inputs,inType);

また、ファイルの出力は、SimpleWriterクラスを用いて行います。

　IOType outType = IOTypeHolder.instance().filenameToType(outputf);
　swriter = new SimpleWriter(outputs,outType);

そして、JOEMolクラスのインスタンスの生成を次のように行います。
　
　JOEMol mol = new JOEMol(inType,outType);

分子の入力は、readNextメソッドで行い、出力はwriteNextメソッドで行っています。

　while(true){
　　if(!sreader.readNext(mol))
　　　　break;
　　swriter.writeNext(mol);
　}

実行
# java jbabel inputfile outputfile

拡張子で、ファイル形式の判定を行っていますので、とてもシンプルになりました。ファイル変換だけをみればbabelと同等の機能をもっていると思います。

以下に今回のソースコードを示しております。変換可能なファイル形式は、ココに記載されています。

簡単にbabelもどきができましたね(当たり前ですか..)。個人的にJOELibのクラス構成をとても気に入っていますので、ちょくちょくメモする予定です。

Jbabel.java




人気ブログランキング(クリックして応援してね)

RuleOfFiveDescriptorクラスには、以下６つのメソッドがあります。

　calculate(AtomContainer mol)
　getParameterNames()
　getParameters()
　getParameterType(java.lang.String name)
　getSpecification()
　setParameters(java.lang.Object[] params)

APIの説明を読むと、getParameterNames()及びTypeあたりでクラス内で使われているdescriptorを知り、setParameters(java.lang.Object[] params)でそのdescriptorの値を任意に指定できるのかなと思いました。デフォルトでは、Lipinski's Ruleだろうけど、Mw<350などに設定したい場合があるので結構親切じゃんっと思い込んだわけです。ただ、calculate(AtomContainer mol)の説明に以下のような不可解な文書が記載されています。

　the method take a boolean checkAromaticity: if the boolean is 　
　true, it means that aromaticity has to be checked

checkAromaticityがtrueならaromaticityのチェックをするんだなということは分かりますが、どうやって設定するの？と思ってしまいます。aromaticityのチェックはしておかないと特にXLogPの計算結果に影響しそうなので...もしかしてデフォルトでFALSEではないよな...など心配になってきたので、ソースコードを見ることにしました。

まずは、setParameters(java.lang.Object[] params)を見ると次のようなコードが書かれていました。

　if (params.length > 1) {
　　　throw new CDKException("RuleOfFiveDescriptor only expects one
　　　parameter");
　 }
if (!(params[0] instanceof Boolean)) {
throw new CDKException("The first parameter must be of type
Boolean");
}
checkAromaticity = ((Boolean) params[0]).booleanValue();

えっ！このメソッドからcheckAromaticityがセットできることは分かりましたが、いろんなパラメータのセットができるという私の推測は見事にはずれてしまいました。ということで、以下のコードでcheckAromaticityをTRUEにできます。

　// checkAromaticityをTRUEにする
　Boolean param[] = new Boolean[1];
　param[0] = Boolean.TRUE;
　try {
　　rule5.setParameters(param);
　} catch (CDKException e1) {
　　System.err.println(e1.toString());
　}

本当は、コード作成者は、複数のパラメータをセットできるようにしたかったのではないかな.....
ちなみにcheckAromaticityはデフォルトでは、falseですので、注意してください。あと、当然ですが、getParameterNames()でも返ってくるのは、"checkAromaticity"だけです。

ということで、RuleOfFiveDescriptorクラスは、checkAromaticityをTRUEかFALSEに設定して、calculate()すればよい(それしかない)と思われます。

前回も記載しましたが、"rule of five"フィルタのコアになるコードは以下のようになります。
　
　//RuleOf5の計算
　dv = rule5.calculate((AtomContainer)mol);
　IntegerResult ir = (IntegerResult)dv.getValue();

　//RuleOf5を満たすものをファイル出力
　if(ir.intValue() == 0){
　　mw.writeMolecule(mol);
　　count++;
　}

これで、"rule of five"フィルタの完成です。
注意点としては、CDKにおける水素結合のドナー、アクセプターの定義を確認しないとせっかくの有望な化合物がこぼれてしまう可能性があること。また、CDKのlogPの計算方法(精度)も把握しておかないと同様の心配がありますので、このあたりもそのうち取り上げてみます。

CDKex2.java




人気ブログランキング(クリックして応援してね)

CDKには、"rule of five"フィルタを簡単に構築できるRuleOfFiveDescriptorクラスがあります。"rule of five"とは、Lipinskiにより提案された化合物のDrug-likeness(薬らしさ)を評価する方法の1つです。

RuleOfFiveDescriptorクラスでは、以下のように"rule of five"の各項目を満たさなければlipinskifailuresに1を加えています。すなわち、"rule of five"を全て満たした化合物は、lipinskifailuresが0となるわけです。

　if(xlogPvalue > 5.0) { lipinskifailures += 1; }
　if(acceptors > 10) { lipinskifailures += 1; }
　if(donors > 5) { lipinskifailures += 1; }
　if(mwvalue > 500.0) { lipinskifailures += 1; }
　if(rotatablebonds > 10.0) { lipinskifailures += 1; }

今回は、sdf形式のファイルを入力し、"rule of five"を満たす化合物のファイルを出力するプログラムを作成したいと思います。

　Usage: CDKex2 data.sdf data_passed.sdf

まずは、概略ですが、data.sdfに含まれる化合物を順次取り込むために、前回用いたMDLReaderクラスではなく、IteratingMDLReaderクラスを用います。

　IteratingMDLReader itr = new IteratingMDLReader(fr);
　while(itr.hasNext())
　　　mol =(Molecule)itr.next();

sdf形式のファイル中にH原子が付加されていない場合は、HydrogenAdderクラスを用いて付加できます。

　HydrogenAdder adder = new HydrogenAdder();
　adder.addExplicitHydrogensToSatisfyValency(mol);

"rule of five"の評価はRuleOfFiveDescriptorクラスのcalculate()メソッドにより実施し、getValue()メソッドを用いて前述したlipinskifailuresの値を取得します。

　RuleOfFiveDescriptor rule5 = new RuleOfFiveDescriptor();
　dv = rule5.calculate((AtomContainer)mol);
　IntegerResult ir = (IntegerResult)dv.getValue();
　if(ir.intValue() == 0)
　　　mw.writeMolecule(mol);

RuleOfFiveDescriptorクラスはちょっとクセ(問題)のあるクラスですので、これも含め、次回、実際のコードを用いてプログラムの解説をしたいと思います。




人気ブログランキング(クリックして応援してね)

普段、Amber/CHARMM/NAMDなどを研究に利用していると、一度は、分子動力学(MD)を自力でプログラミングしたい。もしくは、MDの中身をより詳しく理解したい。このような話をよく聞きます。
私は、学生のときMDを自力で開発していましたが、結局、私の実力ではAmberには勝てないなっと心底悟らされ、今ではMDを与えられた課題を解くツールとして、いかに使いこなすかで満足しています。パラメータを検討した(いじりまわした)末に、実測値と相関が取れたりするとマジっか?と思いながらもうれしいものです。

MD関連の書籍は、多くでていると思いますが、MDの理論体系を見通せる書籍としては、以下の2冊がお勧めです。




実際に、MDをプログラミングするとっかかりとして、非常に参考になったのは、以下の書籍でした。FORTRANで書かれたソースプログラムが掲載されていますので、Javaとかで書き直しただけで、何となくMDを自作した気分になれます(笑)。ただし、実際にプログラミングするといろいろな問題点が把握できると思いますので、MDを利用する上でも役に立つと思います。




人気ブログランキング(クリックして応援してね)

FrownsでのFingerprintの使い方をメモしたいと思います。
　
　import frowns.Fingerprint

によりモジュールを読み込み、generateFingerprint関数によりFingerprintを生成します。

　def generateFingerprint(molecule, numInts=32, pathLength=7):

デフォルトでは、ビット長は32、部分構造のパス長が7となっています。
ここでは、生成したFingerprintを部分構造検索の前工程に利用するプログラムを作成致します。

プログラムの概要ですが、

1.
クエリー構造のFingerprintを作成する。ただし、ここでは、Smiles.smilinで変換可能なSMARTSパターンを入力すること。

2. #Fingerprint generation
標的となる化合物郡のFingerprintを全て生成し、dbリストに収める。

3. #Smarts Search
第3引数であるflagが1ならSMARTS検索の前工程としてFingerprint同士の一致を評価し、0なら評価しない。
その後、SMARTS検索。

実行例1：Fingerprintを利用しない。
$ python smartsearch2.py data.sdf "c1ccccc1OCCCCC" 0

出力：
C(=O)(O)c1c(OC(=O)CCCCC)cccc1
O(CCCCCCCCCCCC)c1ccc(C(=O)O)cc1
Time 4.43799996376
2 Hit compounds

実行例2：Fingerprintを利用する。
$ python smartsearch2.py data.sdf "c1ccccc1OCCCCC" 1

出力：
C(=O)(O)c1c(OC(=O)CCCCC)cccc1
O(CCCCCCCCCCCC)c1ccc(C(=O)O)cc1
Time 1.84400010109
2 Hit compounds

この例では、Fingerprintの利用により約2.4倍、検索速度のUPが確認できました。実際の検索システムでは、Fingerprintは化合物のデータベースへの登録時に、ビット文字列などで格納されることが多いと思います。

smartsearch2.py




人気ブログランキング(クリックして応援してね)