🏝️

PyMOLを通してPubChemから化学構造情報を取得する

2024/12/14に公開

Python

pymol

pubchem

tech

PyMOLでSMILESを扱うことで、タンパクや核酸のような高分子だけでなく、低分子や天然物のような化学構造を描画する方法を紹介しました。

例えば会社とかで立体構造をサクッと確認したい時に便利だったりしますが、承認薬だったり有名どころの化合物はPublicのデータベースからAPIで取ってきた方が早くて楽で正確です。　　
SciFinderやReaxysのような商用のデータベースには基本的には敵わないんですが、有名どころをサクッと扱う分にはOpen Souceのデータベースで十分です。
シンプルに化学構造を扱うだけであればPubChemが使いやすいかなと思います。活性情報が欲しい場合はChEMBL辺りが触りやすいです。
https://pubchem.ncbi.nlm.nih.gov/
https://www.ebi.ac.uk/chembl/
今回は一般名などで検索性が高いPubChemを使ってPyMOlに化学構造をLoadしていきます。

PubChem.py

いつも通り、下記のスクリプトを保存したディレクトリでPyMOLを起動してください。

PubChem.py

from rdkit import Chem
from rdkit.Chem import AllChem
import pubchempy as pcp
from openbabel import openbabel


# https://iwatobipen.wordpress.com/2024/05/14/add-hydrogen-with-user-defined-ph-from-python-openbabel-cheminformatics/
obc = openbabel.OBConversion()
obc.SetInAndOutFormats('smi', 'smi')

# adjusting the pH of SMILES molecules
def get_smi_with_pH(smi, pH:float=7.4):
    obmol = openbabel.OBMol()
    obc.ReadString(obmol, smi)
    obmol.CorrectForPH(float(pH))
    return obc.WriteString(obmol)

# calling molecule from PubChem database (with pH)
def pc(name: str, pH:float=None):
    try:
        sdf = pcp.get_sdf(name, 'name')
        if sdf is None:
            print(f"Error: Could not retrieve SDF for {name}")
            return
        mol = Chem.MolFromMolBlock(sdf)
        smile = Chem.MolToSmiles(mol)
        if not pH is None:
            smile = get_smi_with_pH(smile, pH)
        mol = Chem.MolFromSmiles(smile)
        mol = Chem.AddHs(mol)
        AllChem.EmbedMolecule(mol)
        AllChem.UFFOptimizeMolecule(mol)
        pdb_block = Chem.MolToPDBBlock(mol)
        cmd.read_pdbstr(pdb_block, name)
    except Exception as e:
        print(f"An error occurred: {e}")
    cmd.hide(f'({name} and hydro and (elem C extend 1))')
cmd.extend('pc', pc)

# Calling smiles from PubChem database
def pubchem2smiles(compound_name: str):
    try:
        # PubChemから化合物情報を取得
        compound = pcp.get_compounds(compound_name, 'name')
        if not compound:
            print(f"No compounds found for {compound_name}")
            return

        # SMILESを取得
        smiles = compound[0].canonical_smiles
        if not smiles:
            print(f"SMILES not found for {compound_name}")
            return

        # ファイルに保存
        file_name = f"{compound_name.replace(' ', '_')}.smi"
        with open(file_name, 'w') as f:
            f.write(smiles)
        print(f"SMILES for {compound_name} saved in {file_name} as")
        print(f"{smiles}")

    except Exception as e:
        print(f"An error occurred: {e}")
cmd.extend('pubchem2smiles', pubchem2smiles)

PyMOL

run PubChem.py

今更ですが、どこかの別のディレクトリにスクリプトを保存してもらっても、run以降を絶対パスにしてもらえれば問題なくコマンドをインストールできます。　　
例えばDocuments/PyMOLというディレクトリにPubChem.pyのスクリプトを保存したとして、どこのディレクトリからでも下記のコマンドで新しいコマンドをインストールできます。

PyMOL

run ~/Documents/PyMOL/PubChem.py

使い方

使い方は自体はそんなに難しくありません。

PyMOL console

pc [molecule name]
# ex)pc aspirin

molecule name の引数をPubChemで検索してtopの結果から分子を読み込んでいます。

PyMOL console

pubchem2smiles [molecule name]
# ex) pubchem2smiles aspirin

上記と同じく引数をPubChemで検索してtopの結果のSMILESを.smiで保存します。

終わりに

残念ながら上記のコードは時々エラーが出ます。　　
例えばhalichondrin BはPubChemのBest matchで正しい構造がトップにきますが、pubchempyでは情報をちゃんと引っ張ってこれないみたいです。　　
スペースを含まない化合物名でも同様のエラーが出る分子があったので、スペースの処理の問題では内容でした。　　
この辺が解決出来たら追記しようと思います。

PubChem.py

使い方

終わりに

Discussion