AutoDock Vina — 蛋白质-配体分子对接工具¶

# 方法1：conda 安装（推荐）
conda install -c bioconda autodock-vina  # 安装 AutoDock Vina（当前 v1.2.7）
vina --version  # 确认版本

# 安装辅助工具（准备对接文件）
conda install -c conda-forge -c bioconda openbabel  # 格式转换（SDF→PDBQT）
pip install meeko   # 现代配体准备工具（推荐，替代 MGLTools）

# MGLTools（传统工具，含 prepare_receptor4.py）
conda install -c bioconda mgltools  # 或手动下载

# 安装 AutoDock-GPU（GPU 加速版，高通量筛选时用）
conda install -c bioconda autodock-gpu  # GPU 版本

# 安装 Gnina（深度学习增强版 Vina）
conda install -c bioconda gnina  # 含卷积神经网络评分

# === 步骤 1：准备受体蛋白（.pdbqt 格式）===

# 先用 PyMOL/ChimeraX 处理 PDB 文件：去水、去配体、加氢
# 然后用 prepare_receptor 脚本转换为 pdbqt 格式
python prepare_receptor4.py \
  -r protein_clean.pdb \     # 清理后的蛋白质 PDB
  -o receptor.pdbqt \        # 输出 pdbqt 格式
  -A hydrogens \             # 添加氢原子
  -U nphs_lps               # 去掉非极性氢

# === 步骤 2：准备配体（.pdbqt 格式）===

# 方法A：用 meeko（现代推荐）
mk_prepare_ligand.py -i ligand.sdf -o ligand.pdbqt  # SDF → PDBQT

# 方法B：用 openbabel
obabel ligand.sdf -O ligand.pdbqt --partialcharge gasteiger  # 加 Gasteiger 电荷

# === 步骤 3：运行对接 ===
vina \
  --receptor receptor.pdbqt \   # 受体文件
  --ligand ligand.pdbqt \       # 配体文件
  --center_x -12.5 \            # 对接盒子中心 X 坐标（Å）
  --center_y 15.2 \             # 对接盒子中心 Y 坐标（Å）
  --center_z 8.7 \              # 对接盒子中心 Z 坐标（Å）
  --size_x 25 \                 # 对接盒子大小 X（Å）
  --size_y 25 \                 # 对接盒子大小 Y（Å）
  --size_z 25 \                 # 对接盒子大小 Z（Å）
  --out ligand_docked.pdbqt \   # 对接结果输出文件
  --log docking_log.txt \       # 日志文件
  --exhaustiveness 8 \          # 搜索穷举程度（8=默认，值越大越准但越慢）
  --num_modes 9                 # 输出构象数量（默认9）

# 方法1：PyMOL 中查看活性位点坐标
# 在 PyMOL 命令行：get_position   → 获取当前视角中心坐标

# 方法2：用已知配体定义盒子中心
python3 << 'EOF'
from rdkit import Chem            # 导入 RDKit（需安装）
from rdkit.Chem import AllChem   # 化学功能模块
import numpy as np

# 从配体 SDF 文件获取质心作为盒子中心
mol = Chem.SDMolSupplier("reference_ligand.sdf")[0]  # 读取已知配体
conf = mol.GetConformer()         # 获取构象
coords = conf.GetPositions()      # 获取原子坐标 (N, 3)
center = coords.mean(axis=0)      # 计算质心

print(f"盒子中心: X={center[0]:.2f}, Y={center[1]:.2f}, Z={center[2]:.2f}")
# 盒子大小：配体坐标范围 + 各方向各加 8-10 Å
extent = coords.max(axis=0) - coords.min(axis=0) + 10  # 各方向+10Å
print(f"建议盒子大小: X={extent[0]:.1f}, Y={extent[1]:.1f}, Z={extent[2]:.1f}")
EOF

# 准备化合物库（转为单个 pdbqt 文件）
mkdir ligand_pdbqt/
for sdf_file in library/*.sdf; do
    basename=$(basename $sdf_file .sdf)   # 提取文件名（不含后缀）
    mk_prepare_ligand.py \
      -i "$sdf_file" \                    # 输入 SDF
      -o "ligand_pdbqt/${basename}.pdbqt" # 输出 PDBQT
done

# 批量对接脚本
python3 << 'EOF'
import subprocess  # 运行系统命令
import pandas as pd  # 结果整理
import re          # 正则表达式提取分数
import os

results = []  # 存储所有结果

for pdbqt in os.listdir("ligand_pdbqt/"):
    if not pdbqt.endswith(".pdbqt"):  # 只处理 pdbqt 文件
        continue

    name = pdbqt.replace(".pdbqt", "")   # 化合物名称
    output = f"docked/{name}_out.pdbqt"  # 对接输出文件

    # 运行 vina
    cmd = [
        "vina",
        "--receptor", "receptor.pdbqt",
        "--ligand", f"ligand_pdbqt/{pdbqt}",
        "--center_x", "-12.5",           # 盒子中心坐标
        "--center_y", "15.2",
        "--center_z", "8.7",
        "--size_x", "25",                # 盒子大小
        "--size_y", "25",
        "--size_z", "25",
        "--out", output,
        "--exhaustiveness", "8"          # 穷举程度
    ]

    result = subprocess.run(cmd, capture_output=True, text=True)
    output_text = result.stdout + result.stderr

    # 提取最佳对接能量（第一个模型的分数）
    match = re.search(r'1\s+([-\d.]+)\s', output_text)
    if match:
        score = float(match.group(1))    # 对接能量（kcal/mol）
        results.append({"compound": name, "docking_score": score})
        print(f"{name}: {score:.2f} kcal/mol")

# 按对接分数排序（越负越好）
df = pd.DataFrame(results).sort_values("docking_score")
df.to_csv("screening_results.csv", index=False)  # 保存排名
print("\n前10名：")
print(df.head(10))  # 显示前10名
EOF

# 安装
conda install -c bioconda autodock-vina
pip install meeko   # 现代配体准备

# 配体准备
mk_prepare_ligand.py -i mol.sdf -o mol.pdbqt

# 受体准备
python prepare_receptor4.py -r protein.pdb -o receptor.pdbqt -A hydrogens

# 基本对接
vina --receptor receptor.pdbqt --ligand ligand.pdbqt \
     --center_x X --center_y Y --center_z Z \
     --size_x 25 --size_y 25 --size_z 25 \
     --out out.pdbqt --exhaustiveness 8

# 结合能解读
# < -9 kcal/mol → 强结合，值得关注
# -7 ~ -9       → 中等结合
# > -5           → 弱结合，一般不考虑

# 当前版本：AutoDock Vina v1.2.7