Chemical Reactions part1 〜RDKit 直訳 Day17〜

（12/30追記）試訳をまとめたテストサイトを作成しました。よろしければご参照ください。

こちらはRDKit直訳 Advent Calendar 2018 - Adventar 17日目の記事です。基本的な進め方は1日目の記事をご覧ください。

本日の訳出に困った用語
SMARTS-based language: SMARTSベースの記法
MDL rxn files: MDL rxnファイル
unique products: 単一の生成物
sanitize: サニタイズ
the hierarchy of transformations: 変換の階層構造
dummy atom: ダミーアトム

以下、訳

化学反応(Chemical Reactions)

[Link] https://www.rdkit.org/docs/GettingStartedInPython.html#chemical-reactions

RDKitは分子の集合に対して化学反応を適用することもサポートしています。化学反応を構築する一つの方法は、Daylight社のReaction SMILES [脚注11] に似たSMARTSベースの記法を使うことです。:

[脚注11] A more detailed description of reaction smarts, as defined by the rdkit, is in the The RDKit Book.
[link] The RDKit Book

>>> rxn = AllChem.ReactionFromSmarts('[C:1](=[O:2])-[OD1].[N!H0:3]>>[C:1](=[O:2])[N:3]')
>>> rxn
<rdkit.Chem.rdChemReactions.ChemicalReaction object at 0x...>
>>> rxn.GetNumProductTemplates()
1
>>> ps = rxn.RunReactants((Chem.MolFromSmiles('CC(=O)O'),Chem.MolFromSmiles('NC')))
>>> len(ps) # one entry for each possible set of products
1
>>> len(ps[0]) # each entry contains one molecule for each product
1
>>> Chem.MolToSmiles(ps[0][0])
'CNC(C)=O'
>>> ps = rxn.RunReactants((Chem.MolFromSmiles('C(COC(=O)O)C(=O)O'),Chem.MolFromSmiles('NC')))
>>> len(ps)
2
>>> Chem.MolToSmiles(ps[0][0])
'CNC(=O)OCCC(=O)O'
>>> Chem.MolToSmiles(ps[1][0])
'CNC(=O)CCOC(=O)O

MDL rxnファイルから反応を構築することもできます:

>>> rxn = AllChem.ReactionFromRxnFile('data/AmideBond.rxn')
>>> rxn.GetNumReactantTemplates()
2
>>> rxn.GetNumProductTemplates()
1
>>> ps = rxn.RunReactants((Chem.MolFromSmiles('CC(=O)O'), Chem.MolFromSmiles('NC')))
>>> len(ps)
1
>>> Chem.MolToSmiles(ps[0][0])
'CNC(C)=O'

もちろんアミド結合の形成よりももっと複雑な反応を行うこともできます:

>>> rxn = AllChem.ReactionFromSmarts('[C:1]=[C:2].[C:3]=[*:4][*:5]=[C:6]>>[C:1]1[C:2][C:3][*:4]=[*:5][C:6]1')
>>> ps = rxn.RunReactants((Chem.MolFromSmiles('OC=C'), Chem.MolFromSmiles('C=CC(N)=C')))
>>> Chem.MolToSmiles(ps[0][0])
'NC1=CCCC(O)C1'

この場合、テンプレートへの反応物質のマッピングが1通りではないので、得られるものが複数の生成物を含むセットであることに注意してください:

>>> len(ps)
4

カノニカルSMILESとpythonの辞書型を使うことで、それぞれのユニークな生成物を得ることができます。

>>> uniqps = {}
>>> for p in ps:
...   smi = Chem.MolToSmiles(p[0])
...   uniqps[smi] = p[0]
...
>>> sorted(uniqps.keys())
['NC1=CCC(O)CC1', 'NC1=CCCC(O)C1']

化学反応を処理するコードから生成される分子はサニタイズされていないことに注意して下さい。あえて不自然な反応例をつかって示します:

>>> rxn = AllChem.ReactionFromSmarts('[C:1]=[C:2][C:3]=[C:4].[C:5]=[C:6]>>[C:1]1=[C:2][C:3]=[C:4][C:5]=[C:6]1')
>>> ps = rxn.RunReactants((Chem.MolFromSmiles('C=CC=C'), Chem.MolFromSmiles('C=C')))
>>> Chem.MolToSmiles(ps[0][0])
'C1=CC=CC=C1'
>>> p0 = ps[0][0]
>>> Chem.SanitizeMol(p0)
rdkit.Chem.rdmolops.SanitizeFlags.SANITIZE_NONE
>>> Chem.MolToSmiles(p0)
'c1ccccc1'

より高度な反応機能(Advanced Reaction Functionality)

[Link] https://www.rdkit.org/docs/GettingStartedInPython.html#advanced-reaction-functionality

原子の保護(Protecting Atoms)

[Link] https://www.rdkit.org/docs/GettingStartedInPython.html#protecting-atoms

ありえない生成物を生み出してしまわないように十分に正確に反応を表現することは、時折難しく、特にrxnファイルを扱っているときにこのような問題が生じます。RDKitでは原子を"保護する"方法を提供しており、その原子が反応に参加しないようにすることができます。

上のセクジョンで用いたアミド結合形成の反応を再度用いてデモンストレーションしてみましょう。アミンのクエリが十分に明確ではなく、少なくとも一つの水素原子の結合した窒素原子はどれでもマッチしてしまいます。したがって、この反応をすでにアミド結合を有する分子に適用した場合、アミドのNも反応点として扱われてしまいます:

>>> rxn = AllChem.ReactionFromRxnFile('data/AmideBond.rxn')
>>> acid = Chem.MolFromSmiles('CC(=O)O')
>>> base = Chem.MolFromSmiles('CC(=O)NCCN')
>>> ps = rxn.RunReactants((acid,base))
>>> len(ps)
2
>>> Chem.MolToSmiles(ps[0][0])
'CC(=O)N(CCN)C(C)=O'
>>> Chem.MolToSmiles(ps[1][0])
'CC(=O)NCCNC(C)=O'

最初の生成物がアミドのNで反応が起きたものに相当します。

全てのアミドのNを保護することで、このような事態が生じるのを防ぐことができます。ここでは、アミドとチオアミドにマッチする部分構造クエリをつかって保護を行い、マッチする原子に"_protected"プロパティをセットしましょう:

>>> amidep = Chem.MolFromSmarts('[N;$(NC=[O,S])]')
>>> for match in base.GetSubstructMatches(amidep):
...     base.GetAtomWithIdx(match[0]).SetProp('_protected','1')

これで、たった一つだけ生成物を生じるようになります:

>>> ps = rxn.RunReactants((acid,base))
>>> len(ps)
1
>>> Chem.MolToSmiles(ps[0][0])
'CC(=O)NCCNC(C)=O'

Recapの実装(Recap Implementation)

[Link] https://www.rdkit.org/docs/GettingStartedInPython.html#recap-implementation

化学反応機能と関連して、Recapアルゴリズムの実装を取り上げます。[脚注8] Recapは、実験室でよく使われる反応を模倣した化学変換のセットを使って、分子を一連の合理的なフラグメントに分解します。

[脚注8] Lewell, X.Q.; Judd, D.B.; Watson, S.P.; Hann, M.M. “RECAP-Retrosynthetic Combinatorial Analysis Procedure: A Powerful New Technique for Identifying Privileged Molecular Fragments with Useful Applications in Combinatorial Chemistry” J. Chem. Inf. Comp. Sci. 38:511-22 (1998).

RDKitのrdkit.Chem.Recapの実装は、適用される変換の階層構造の流れを保存します。

[link] rdkit.Chem.Recap

>>> from rdkit import Chem
>>> from rdkit.Chem import Recap
>>> m = Chem.MolFromSmiles('c1ccccc1OCCOC(=O)CC')
>>> hierarch = Recap.RecapDecompose(m)
>>> type(hierarch)
<class 'rdkit.Chem.Recap.RecapHierarchyNode'>

階層構造の起点は元々の分子です：

>>> hierarch.smiles
'CCC(=O)OCCOc1ccccc1'

それぞれのノードはSMILESをキーとする辞書型を使って、子ノードの情報を保持します。

>>> ks=hierarch.children.keys()
>>> sorted(ks)
['*C(=O)CC', '*CCOC(=O)CC', '*CCOc1ccccc1', '*OCCOc1ccccc1', '*c1ccccc1']

階層構造の末端のノード（葉ノード）は簡単にアクセスすることができ、これもSMILESをキーとする辞書型となっています。

>>> ks=hierarch.GetLeaves().keys()
>>> ks=sorted(ks)
>>> ks
['*C(=O)CC', '*CCO*', '*CCOc1ccccc1', '*c1ccccc1']

分子がフラグメント化された場所では、ダミーアトムをつかって印がつけられていることに注意して下さい。

ノード自体には関連づけられた分子が入っています:

>>> leaf = hierarch.GetLeaves()[ks[0]]
>>> Chem.MolToSmiles(leaf.mol)
'*C(=O)CC'

12/17/2018

このブログ記事のライセンス

以上はRDKit Documentationの試訳です。
ライセンスはCC BY-SA 4.0で、Greg Landrum氏の功績の元に作成しています。
Creative Commons — Attribution-ShareAlike 4.0 International — CC BY-SA 4.0