diffdock

diffdock

概要

DiffDock は、拡散生成モデルを使用して、タンパク質構造とリガンドの SMILES からタンパク質とリガンドの結合ポーズを直接予測します。ドッキングを探索問題ではなく生成問題として扱います。AutoDock Vina や Glide などの従来のドッキングツールとは異なり、DiffDock では結合部位を事前に指定する必要がなく、タンパク質表面全体からポーズをサンプリングします。信頼度スコアを付与した、順位付け済みの結合ポーズを出力します。DiffDock はブラインドドッキングタスクに優れ、多様なポーズ仮説を生成するため、新規の結合部位の発見や、難しい標的の解析に役立ちます。

使用する場面

• ブラインドドッキング（結合部位が不明）：リガンドがタンパク質のどこに結合するか分からず、結合部位の候補を見つけたい場合。
• 従来のドッキングではうまくいかない難しい標的：アロステリック部位、柔軟な領域、または標的の結合部位に共結晶構造がないタンパク質。
• 複数の結合モードの探索：結合時の構造的柔軟性を把握するため、多様なポーズのアンサンブルを生成する場合。
• 構造活性相関（SAR）の探索：一連の類縁体をすばやくドッキングして、予測された結合モードを比較する場合。
• フラグメントスクリーニングの仮説生成：フラグメント分子の結合部位候補を特定する場合。
• 結合部位が既知で、タンパク質を剛体として扱う場合は、AutoDock Vina または GNINA のほうが高速で、精度も同等になることがあります。
• 大規模なバーチャルスクリーニング（>10,000 化合物）には、標準の DiffDock ではなく、GNINA または DiffDock-L（大規模向けバージョン）を検討してください。
• 結合ポケットが明確に定義されており、大規模な化合物ライブラリを高速に処理する必要がある場合は、代わりに AutoDock Vina を使用してください。

前提条件

• Python パッケージ：diffdock（conda によるインストールを推奨）、rdkit、torch、biopython、nglview（可視化用）
• システム：GPU（NVIDIA CUDA を推奨）。CPU での推論は低速です（化合物あたり約 5～10 分）。
• データ要件：タンパク質の PDB ファイル（クリーニングおよびプロトン化済み）、リガンドの SMILES 文字列または SDF ファイル
• 環境：CUDA 対応の PyTorch を備えた conda 環境

# Recommended: clone and install from source
git clone https://github.com/gcorso/DiffDock.git
cd DiffDock
conda create -n diffdock python=3.9
conda activate diffdock
pip install torch torchvision --extra-index-url https://download.pytorch.org/whl/cu118
pip install -r requirements.txt

# Download pretrained model weights
python -c "from utils.download import download_pretrained; download_pretrained()"

ワークフロー

ステップ 1：タンパク質構造を準備する

from Bio import PDB
from Bio.PDB import PDBParser, PDBIO, Select

class NonHetSelect(Select):
    """Remove HETATM records (ligands, water) — keep only protein atoms."""
    def accept_residue(self, residue):
        return residue.id[0] == " "

def clean_pdb(input_pdb: str, output_pdb: str):
    parser = PDBParser(QUIET=True)
    structure = parser.get_structure("protein", input_pdb)
    io = PDBIO()
    io.set_structure(structure)
    io.save(output_pdb, NonHetSelect())
    print(f"Cleaned PDB saved to: {output_pdb}")

clean_pdb("raw_protein.pdb", "protein_clean.pdb")

ステップ 2：リガンド入力を準備する

from rdkit import Chem
from rdkit.Chem import AllChem, SDWriter

def smiles_to_sdf(smiles: str, output_sdf: str, n_confs: int = 1):
    """Convert SMILES to 3D SDF for DiffDock input."""
    mol = Chem.MolFromSmiles(smiles)
    mol = Chem.AddHs(mol)
    AllChem.EmbedMolecule(mol, AllChem.ETKDGv3())
    AllChem.MMFFOptimizeMolecule(mol)
    writer = SDWriter(output_sdf)
    writer.write(mol)
    writer.close()
    print(f"Ligand SDF written to: {output_sdf}")
    return mol

# Example: ibuprofen
smiles = "CC(C)Cc1ccc(cc1)C(C)C(=O)O"
mol    = smiles_to_sdf(smiles, "ligand.sdf")
print(f"Ligand formula: {Chem.rdMolDescriptors.CalcMolFormula(mol)}")

ステップ 3：DiffDock 推論を実行する

# Command-line inference (run from the DiffDock directory)
python inference.py \
    --protein_path protein_clean.pdb \
    --ligand       "CC(C)Cc1ccc(cc1)C(C)C(=O)O" \
    --out_dir      results/ \
    --inference_steps 20 \
    --samples_per_complex 40 \
    --batch_size 10 \
    --no_final_step_noise

import subprocess

def run_diffdock(protein_pdb: str, ligand_smiles: str, out_dir: str,
                 n_samples: int = 40, n_steps: int = 20):
    cmd = [
        "python", "inference.py",
        "--protein_path",      protein_pdb,
        "--ligand",            ligand_smiles,
        "--out_dir",           out_dir,
        "--inference_steps",   str(n_steps),
        "--samples_per_complex", str(n_samples),
        "--batch_size",        "10",
        "--no_final_step_noise",
    ]
    result = subprocess.run(cmd, capture_output=True, text=True, cwd="DiffDock/")
    if result.returncode == 0:
        print(f"DiffDock complete. Results in: {out_dir}")
    else:
        print(f"Error: {result.stderr}")
    return result

run_diffdock("protein_clean.pdb", "CC(C)Cc1ccc(cc1)C(C)C(=O)O", "results/")

ステップ 4：信頼度スコアを解析して順位付けする

import re
from pathlib import Path
import pandas as pd

def parse_diffdock_results(out_dir: str) -> pd.DataFrame:
    """Parse DiffDock output SDF files and confidence scores."""
    out_path = Path(out_dir)
    records  = []

    # DiffDock names output files: rank{N}_confidence{score}.sdf
    for sdf_file in sorted(out_path.glob("rank*_confidence*.sdf")):
        name = sdf_file.stem
        # Extract rank and confidence from filename
        rank_match = re.search(r"rank(\d+)", name)
        conf_match = re.search(r"confidence(-?[\d.]+)", name)
        if rank_match and conf_match:
            records.append({
                "rank":       int(rank_match.group(1)),
                "confidence": float(conf_match.group(1)),
                "sdf_file":   str(sdf_file),
            })

    df = pd.DataFrame(records).sort_values("rank")
    print(f"Found {len(df)} poses")
    print(df[["rank", "confidence", "sdf_file"]].head(10))
    return df

df_results = parse_diffdock_results("results/")

ステップ 5：上位のポーズを解析する — 結合部位の残基を抽出する

from rdkit import Chem
from rdkit.Chem import AllChem
from Bio.PDB import PDBParser
import numpy as np

def get_binding_residues(protein_pdb: str, ligand_sdf: str, cutoff_angstrom: float = 4.0):
    """Find protein residues within cutoff distance of the top-ranked ligand pose."""
    parser    = PDBParser(QUIET=True)
    structure = parser.get_structure("prot", protein_pdb)
    prot_atoms = [(atom.get_coord(), residue.resname, residue.id[1])
                  for chain in structure for residue in chain
                  for atom in residue.get_atoms()]

    mol = Chem.SDMolSupplier(ligand_sdf, removeHs=False)[0]
    lig_coords = mol.GetConformer().GetPositions()

    contacts = []
    for prot_coord, resname, resnum in prot_atoms:
        dists = np.linalg.norm(lig_coords - prot_coord, axis=1)
        if dists.min() <= cutoff_angstrom:
            contacts.append((resnum, resname))

    contacts = sorted(set(contacts))
    print(f"Binding site residues within {cutoff_angstrom} A: {contacts[:10]}")
    return contacts

# Use top-ranked pose
top_sdf = df_results.loc[df_results.rank == 1, "sdf_file"].iloc[0]
contacts = get_binding_residues("protein_clean.pdb", top_sdf)

ステップ 6：NGLview でポーズを可視化する

import nglview as nv
from rdkit import Chem

# Load protein + top pose in Jupyter notebook
view = nv.NGLWidget()
view.add_pdbfile("protein_clean.pdb")

top_sdf  = df_results.loc[df_results.rank == 1, "sdf_file"].iloc[0]
view.add_component(top_sdf)
view.representations = [
    {"type": "cartoon", "params": {"color": "chainindex"}},
    {"type": "ball+stick", "params": {"sele": "ligand"}},
]
print(f"Visualizing top pose: confidence={df_results.confidence.iloc[0]:.3f}")
view

主なパラメーター

パラメーター	デフォルト	範囲／オプション	効果
--inference_steps	20	10–40	拡散過程の逆ステップ数。ステップ数を増やすと遅くなりますが、精度が向上します
--samples_per_complex	40	10–100	サンプリングするポーズ数。数を増やすと、結合モードをより幅広く網羅できます
--batch_size	10	1–32	GPU のバッチサイズ。OOM エラーが発生した場合は減らしてください
--no_final_step_noise	off	flag	最後の拡散ステップのノイズを除去し、ポーズの品質を向上させます
--actual_steps	inference_steps と同じ	1–inference_steps	実際に実行するステップ数（合計より少なくできます）
--save_visualisation	off	flag	SDF とともに、可視化用の PDB ファイルも保存します
cutoff_angstrom	4.0	3.0–6.0 Å	結合部位の残基を定義する距離のカットオフ

よく使うレシピ

レシピ：複数のリガンドをまとめてドッキングする

使用場面：SAR 解析のため、同じタンパク質に対して類縁体のライブラリをドッキングします。

import pandas as pd
import subprocess

smiles_list = [
    ("compound_1", "CC(C)Cc1ccc(cc1)C(C)C(=O)O"),
    ("compound_2", "CC(C)Cc1ccc(cc1)C(C)C(=O)N"),
    ("compound_3", "CC(C)Cc1ccc(cc1)C(C)C(=O)OC"),
]

results = []
for name, smiles in smiles_list:
    out = f"results/{name}"
    cmd = ["python", "inference.py",
           "--protein_path", "protein_clean.pdb",
           "--ligand", smiles,
           "--out_dir", out,
           "--inference_steps", "20",
           "--samples_per_complex", "20"]
    subprocess.run(cmd, cwd="DiffDock/", capture_output=True)
    # Parse top confidence score
    df_r = parse_diffdock_results(out)
    if not df_r.empty:
        top_conf = df_r.loc[df_r.rank == 1, "confidence"].iloc[0]
        results.append({"name": name, "smiles": smiles, "top_confidence": top_conf})

df_batch = pd.DataFrame(results).sort_values("top_confidence", ascending=False)
df_batch.to_csv("batch_docking_results.csv", index=False)
print(df_batch)

レシピ：信頼度のしきい値でポーズを絞り込む

使用場面：以降の解析や可視化に向けて、信頼度の高いポーズだけを残します。

# Confidence > 0 generally indicates a plausible binding pose
# DiffDock confidence scores: higher = more confident; ~0 is marginal; < -1 is poor
high_conf = df_results[df_results["confidence"] > 0.0]
print(f"High-confidence poses: {len(high_conf)} / {len(df_results)}")
print(high_conf[["rank", "confidence", "sdf_file"]])

レシピ：上位のポーズを PDBQT に変換し、Vina で再スコアリングする

使用場面：AutoDock Vina のエネルギー関数を使って、DiffDock のポーズを再スコアリングします。

# Convert SDF to PDBQT using OpenBabel
obabel rank1_confidence0.75.sdf -O rank1_ligand.pdbqt
obabel protein_clean.pdb -O protein.pdbqt -xr

# Rescore (no docking search, just energy evaluation)
vina --receptor protein.pdbqt --ligand rank1_ligand.pdbqt \
     --score_only --out rank1_rescored.pdbqt

期待される出力

• results/rank{N}_confidence{score}.sdf — 信頼度スコア順にランク付けされた3Dリガンドポーズ
• 各ポーズのランク、信頼度スコア、ファイルパスを含む df_results DataFrame
• 結合部位残基のリスト（カットオフ距離以内の残基番号と名称）
• Jupyter ノートブックでの NGLview によるインタラクティブな3D可視化

トラブルシューティング

問題	原因	解決策
CUDA out of memory	バッチサイズがGPUに対して大きすぎます	--batch_size を 4 または 2 に減らします
結果ディレクトリが空です	Protein PDBの解析に失敗しました	PDBにATOMレコードのみが含まれていることを確認します。clean_pdb() でHETATMを削除します
すべての信頼度スコアが < -2 です	リガンドまたはタンパク質の形式に問題があります	RDKitでSMILESを検証し、タンパク質がプロトン化され、完全であることを確認します
推論が非常に遅い（30分超）	CPUで実行されています	GPUの使用を強く推奨します。CUDA環境が正しく設定されている必要があります
ModuleNotFoundError: e3nn	依存関係がインストールされていません	DiffDockのconda環境で pip install e3nn を実行します
ポーズが1つの部位に集中します	--samples_per_complex の値が小さいです	部位をより適切にカバーするには、40～100 に増やします
タンパク質に残基がありません	結晶構造が不完全です	ドッキング前にMODELLERまたはSwiss-Modelを使って欠損部分を補完します

参考資料

• DiffDock GitHub（gcorso/DiffDock） — ソースコード、インストール方法、事前学習済み重み
• Corso et al.（2023）、ICLR — DiffDock論文 — 拡散ドッキングに関する原著論文
• DiffDockウェブサーバー — ローカルにインストールせず、ブラウザー上で推論を実行できます
• PatentsView APIドキュメント — 知的財産分析の参考情報