RDKit入門⑧：反応式の取り扱い(後編)の記事では自作関数による化学反応の繰り返しをしつつ、
重複する分子を削除してユニークな分子のみを描画する方法を削除しました。
この記事ではカノニカルSMILESをset型に格納することで重複を削除していましたが、
この方法だとそれぞれの分子にいくつの重複があるのか知ることはできません。
実際の化学反応の視点で考えると、等価な分子が生成しうるパターンが多い分子は
生成比が大きくなる可能性があり、重複数には速度論的意味がある場合が存在します。
そこで今回はただ重複を削除するのではなく、重複数を整理して描画する方法を紹介します。
環境
- Windows 10
- python 3.8.8
- RDKit 2021.03.1
まずはパッケージの読み込みを行います。
今回新たにcollectionsというパッケージを読み込み、使用します。
(collectionsに関しては記事末尾に参考元を記載いたします)
from rdkit import Chem
from rdkit.Chem import AllChem, Draw
from itertools import chain
import collections
#パッケージの読み込み
前回の記事と同様に芳香環上のC-Hブロモ化をトルエンに適用する反応を実施します。
今回は多段階反応定義の際に重複を削除しないかたちで関数を定義します。
toluene = Chem.MolFromSmiles("Cc1ccccc1")
#トルエンを定義した。
def multi_bromination(mol):
rxn = AllChem.ReactionFromSmarts("[ch:1]>>[c:1]Br")
rxn.Initialize()
#反応の定義と念の為の初期化。
bromides = list(chain.from_iterable(rxn.RunReactants([mol])))
#molに対して1回目のrxnを実施し、生成物の格納されたリストをbromidesと命名。
for bromide in bromides:
bromide.UpdatePropertyCache(strict = False)
#価数などの更新作業。
if rxn.IsMoleculeReactant(bromide):
bromides.extend(list(chain.from_iterable(rxn.RunReactants([bromide]))))
#bromides内の分子がrxnの基質だった場合、さらに反応を実施して生成物をbromidesに追加
return bromides
products = multi_bromination(toluene)
#multi_brominationにtolueneを適用し、得られたリストをproductsと命名した。
len(products) #productsの大きさは325であった。
続いてproducts内の重複数を数えます。
collections.Counter関数を使用し、要素と出現回数が紐付いた辞書型オブジェクトを作成します。
なお、数える際はmolオブジェクトをSMILESに戻して取り扱うようにします。
products_count = collections.Counter(Chem.MolToSmiles(mol) for mol in products)
#productsの重複数を数え、辞書型オブジェクトproducts_countに格納した。
products_countの中身を確認すると以下のようになっています。
Counter({'Cc1ccccc1Br': 2,
'Cc1cccc(Br)c1': 2,
'Cc1ccc(Br)cc1': 1,
'Cc1cccc(Br)c1Br': 4,
'Cc1c(Br)cccc1Br': 2,
'Cc1ccc(Br)cc1Br': 4,
'Cc1cc(Br)ccc1Br': 4,
'Cc1ccc(Br)c(Br)c1': 4,
'Cc1cc(Br)cc(Br)c1': 2,
'Cc1ccc(Br)c(Br)c1Br': 12,
'Cc1cc(Br)cc(Br)c1Br': 12,
'Cc1c(Br)ccc(Br)c1Br': 12,
'Cc1c(Br)cc(Br)cc1Br': 6,
'Cc1cc(Br)c(Br)cc1Br': 12,
'Cc1cc(Br)c(Br)c(Br)c1': 6,
'Cc1cc(Br)c(Br)c(Br)c1Br': 48,
'Cc1c(Br)cc(Br)c(Br)c1Br': 48,
'Cc1c(Br)c(Br)cc(Br)c1Br': 24,
'Cc1c(Br)c(Br)c(Br)c(Br)c1Br': 120})
あとはproducts_countのkeyとvalueを利用して分子群描画用のリストを作成します。
products_unique = [Chem.MolFromSmiles(smile) for smile in products_count.keys()]
#products_countのkeyをmolオブジェクトに変換してリストに格納し、products_uniqueと命名した。
products_legends = [str(products_count[x]) for x in products_count.keys()]
#products_count内のkeyに対応するvalueを取り出して文字列に変換してリストに格納し、products_legendsと命名した。
あとは整列のための二次元座標計算などの加工を行い、描画します。
AllChem.Compute2DCoords(toluene)
#tolueneの二次元座標を計算した。
for mol in products_unique:
if mol.HasSubstructMatch(toluene):
AllChem.GenerateDepictionMatching2DStructure(mol, toluene)
#tolueneをテンプレートとしてproducts_unique内のmolオブジェクトを整列させた。
Draw.MolsToGridImage(products_unique, molsPerRow = 5, legends = products_legends, maxMols = len(products_unique))
#products_uniqueを重複数のリストproducts_legendsとともに描画した。
おわりに
今回は繰り返し反応の出力の重複数をcollections.Counter関数で整理する方法を紹介いたしました。
RDKitの化学反応を題材に使用方法を紹介いたしましたが、他にも応用の効く便利な関数です。