RDKit の基本的な使い方
RDKit はケモインフォマティクスを行う上で、必要不可欠なツールです。
大体のことは RDKit で出来てしまうのですが、公式ドキュメントがめちゃくちゃわかりにくいという致命的な欠点があります。
私が使った範囲で、RDKit の機能についてメモを示します。
RDKit のインストール (pip)
pip install rdkit
RDKit のインストール (conda)
conda install -c rdkit rdkit
RDKit のバージョン確認
import rdkit
rdkit.__version__
SMILES から mol オブジェクトへの変換
from rdkit import Chem mol = Chem.MolFromSmiles("CCO")
Mol オブジェクトから SMILES への変換
smiles = Chem.MolToSmiles(mol)
Canonical SMILES への変換
canonical_smiles = Chem.MolToSmiles(Chem.MolFromSmiles("c1ccnc(Cl)c1"))
同じ要領で SMARTS、Inchi、CXSMILES 等も扱える
smarts = Chem.MolToSmarts(mol) inchi = Chem.MolToInchi(mol) cxsmiles = Chem.MolToCXSmiles(mol)
水素原子を付ける
mol_h = Chem.AddHs(mol)
隣接行列の算出
Chem.GetAdjacencyMatrix(mol)
距離行列の算出
Chem.GetDistanceMatrix(mol)
距離行列を用いて Wiener index も算出できる
np.sum(Chem.GetDistanceMatrix(mol)) / 2
分子量の算出
from rdkit.Chem.Descriptors import ExactMolWt ExactMolWt(mol)
QED の算出
from rdkit.Chem.Descriptors import qed qed(mol)
SA score の算出
import sys import os from rdkit.Chem import RDConfig sys.path.append(os.path.join(RDConfig.RDContribDir, "SA_Score")) import sascorer sascorer.calculateScore(Chem.MolFromSmiles("c1ccccc1OC(=O)O"))
TPSA の算出
from rdkit.Chem.rdMolDescriptors import CalcTPSA mol = Chem.MolFromSmiles("CCO") CalcTPSA(mol)
Labute ASA の算出
from rdkit.Chem.MolSurf import LabuteASA mol = Chem.MolFromSmiles("CCO") # 立体配座は不要だが、水素原子の有無で値が変わる。 mol = Chem.AddHs(mol) LabuteASA(mol)
水素原子の授受に関連する官能基
from rdkit.Chem.Lipinski import NumHDonors, NumHAcceptors mol = Chem.MolFromSmiles("CC(=O)c1c(O)cccc1") # 水素原子受容体の数 NumHAcceptors(mol) # 水素原子供与体の数 NumHDonors(mol)
分子中の全炭素原子のうち SP3 炭素の割合
from rdkit.Chem.Lipinski import FractionCSP3 mol = Chem.MolFromSmiles("FC(F)(F)c1cc(O)ccc1C=CC") FractionCSP3(mol)
分子中に含まれる環構造の数を算出。縮合環は別々に数える（例：ナフタレンは 2 とカウント）。
from rdkit.Chem.Lipinski import RingCount mol = Chem.MolFromSmiles("c1c2ccc(CCC3CC3)cc2cc(C(=O)O)c1") RingCount(mol)
立体異性体や幾何異性体を網羅的に発生させる（不斉中心や cis/trans の考慮）
from rdkit.Chem.EnumerateStereoisomers import EnumerateStereoisomers mol = Chem.MolFromSmiles("NC(C=CC)C(=O)O") [m for m in EnumerateStereoisomers(mol)]
指定した部分構造に対応する index を取得する
mol = Chem.MolFromSmiles("c1c(Cl)c(c2nccc(OCCC)c2)ccc1C(=O)NC1CC1") query = Chem.MolFromSmiles("CCC") mol.GetSubstructMatches(query)
分子中の scaffold を取得
from rdkit.Chem.Scaffolds.MurckoScaffold import GetScaffoldForMol, MakeScaffoldGeneric mol = Chem.MolFromSmiles("c1c(OC)c(Cl)ccc1NC(=O)C1CC1CC") # Bemis-Murcko scaffold の取得。環構造と linker とアミドのカルボニル酸素から Bemis-Murcko scaffold は構成される。 GetScaffoldForMol(mol) # 構造式中の原子を全て飽和炭素に置き換えた構造に変換する（グラフが同じアルカン）。 MakeScaffoldGeneric(mol)
共通する部分構造（Maximum Common Substructure: MCS）を抽出
from rdkit.Chem.rdFMCS import FindMCS mols = [Chem.MolFromSmiles(s) for s in ["c1ccccc1OC", "Oc1ccccc1C(=O)O", "Oc1cc(O)cc(O)c1"]] mcs = FindMCS(mols) Chem.MolFromSmarts(mcs.smartsString)
Mol オブジェクトから atom オブジェクトを取り出す。引数は mol オブジェクト内の原子の index。原子番号ではなく、mol オブジェクトによって index が異なる場合があるので注意。
mol.GetAtomWithIdx(0)
Mol オブジェクト中の atom オブジェクトを全て取り出す。
[atom for atom in mol.GetAtoms()]
Atom オブジェクトの質量（原子量）を取り出す。
atom.GetMass()
不斉点の情報（原子の index と R/S）を取得。
mol = Chem.MolFromSmiles("N[C@@H]([C@H](O)C)C(=O)O")
Chem.FindMolChiralCenters(mol)
Mulliken 電荷とグラフ構造から算出される部分電荷である Gasteiger charge を取得。
from rdkit.Chem.rdPartialCharges import ComputeGasteigerCharges # Gasteiger charge の算出 ComputeGasteigerCharges(mol) # Gasteiger charge は atom オブジェクトの property に str 型として記載されている。計算で用いたい場合は型を変換する。 atom = mol.GetAtomWithIdx(0) float(atom.GetProp("_GasteigerCharge"))
Mol オブジェクトから bond オブジェクトを取り出す。引数は mol オブジェクト内の結合の index。Atom オブジェクトと扱いが似ている。
mol.GetBondWithIdx(0)
Mol オブジェクト中の bond オブジェクトを全て取り出す。
[bond for bond in mol.GetBonds()]
Bond オブジェクトの結合次数を取り出す。共役系は 1.5 となる。
bond.GetBondTypeAsDouble()
Bond オブジェクトの結合を構成する原子の atom オブジェクトを取り出す。結合を構成する原子は 2 つあるのでそれぞれ別の関数で取り出す。
bond.GetBeginAtom() bond.GetEndAtom()
分子中の部分構造を検索。出力は部分構造に対応する原子の index のタプル。複数ある場合、index が最も若い一種類だけが返ってくる。反対に、部分構造がヒットしなければ空のタプルが返ってくる。
mol = Chem.MolFromSmiles("c1c(Cl)c(c2nccc(OCCC)c2)ccc1C(=O)NC1CC1") query = Chem.MolFromSmiles("CCC") mol.GetSubstructMatch(query)
検索した部分構造に対応する原子の index の組み合わせ全てを返す。
mol.GetSubstructMatches(query)
Distance geometry 法による立体配座の発生。
from rdkit.Chem.rdDistGeom import ETKDGv3, EmbedMolecule, EmbedMultipleConfs # Distance geometry で立体配座を生成 param = ETKDGv3() # Seed 値は 0 だと固定されない param.randomSeed = 1 # 立体廃座が生成できていれば 0 を返す。入力した mol オブジェクトは立体配座情報が付与される。入力する mol オブジェクトに水素原子を付加することを忘れずに。 EmbedMolecule(mol_h, param)
生成した立体配座の座標を表示。
mol_h.GetConformer().GetPositions()
複数の立体配座をまとめて生成。
from rdkit.Chem.rdDistGeom import EmbedMultipleConfs param = ETKDGv3() # 10 個の立体配座を生成。mol_h に配座情報が加わるので、mol オブジェクトのリストみたいにはならない。 EmbedMultipleConfs(mol_h, 10, param)
立体配座のリストを表示
list(mol_h.GetConformers())
Mol オブジェクトを sdf として保存する。
mols = [Chem.MolFromSmiles(smiles) for smiles in ["CCO", "c1ccccc1", "N[C@@H](C)C(=O)O"]] writer = Chem.SDWriter("saved_mols.sdf") for mol in mols: writer.write(mol) writer.close()
SDF の読み込み
supplier = Chem.SDMolSupplier("saved_mols.sdf") # Supplier は mol オブジェクトのリストに近い操作ができる。 supplier[0]
