# 批量虚拟筛选脚本使用说明

本目录放了一个可直接复用的批量虚拟筛选模板：`vs_batch.py`。
它做的事很简单——给一个受体和一个装满配体的文件夹，用多进程并行调用
AutoDock Vina 命令行，再把每个配体的最优打分会总成一张 CSV。

## 一、安装

推荐用 conda / mamba 建独立环境：

```bash
conda create -n vina python=3.10 -y
conda activate vina
conda install -c conda-forge numpy scipy rdkit vina meeko gemmi autogrid -y
pip install molscrub joblib prody
```

装好后确认命令行能调到：

```bash
vina --version
mk_prepare_ligand.py --help
scrub.py --help
```

`requirements.txt` 里列了完整清单，`pip install -r requirements.txt` 也可，
但 `vina` 用 conda 装通常更省事。

## 二、准备输入

### 1. 受体（PDBQT）

用 Meeko 的 `mk_prepare_receptor.py`。以官方 1IEP 例子为例：

```bash
mk_prepare_receptor.py -i 1iep_receptorH.pdb -o 1iep_receptor -p -v \
    --box_size 20 20 20 --box_center 15.190 53.903 16.917
```

会得到 `1iep_receptor.pdbqt`（受体）和 `1iep_receptor.box.txt`（盒子配置）。
`1iep_receptorH.pdb` 是已加好氢的受体坐标，来自 AutoDock-Vina 官方仓库
`example/basic_docking`。

### 2. 配体（PDBQT）

从 SDF 批量转 PDBQT：

```bash
for f in ligands_sdf/*.sdf; do
    b=$(basename "$f" .sdf)
    mk_prepare_ligand.py -i "$f" -o "ligands/$b.pdbqt"
done
```

如果配体只有 SMILES，先用 `scrub.py` 加氢、生成 3D 构象并枚举质子化状态
（这一步不能省，质子化态选错对接结果会跑偏）：

```bash
# 从 .smi 批量生成 SDF
scrub.py ligands.smi -o mols.sdf
# 再从 SDF 批量转 PDBQT（一次处理多分子，Meeko 支持 --multimol_outdir）
mk_prepare_ligand.py -i mols.sdf --multimol_outdir ligands/
```

### 3. 盒子（config）

`--config` 指向一个文本文件，内容形如：

```
center_x = 15.190
center_y = 53.903
center_z = 16.917
size_x   = 20.0
size_y   = 20.0
size_z   = 20.0
```

没有现成文件时，也可以直接：

```bash
--center 15.190 53.903 16.917 --size 20 20 20
```

脚本会自动生成一个临时 config。

## 三、运行

```bash
python vs_batch.py \
    --receptor 1iep_receptor.pdbqt \
    --ligands  ligands/ \
    --config   1iep_receptor.box.txt \
    --out-dir  poses \
    --csv      results.csv \
    --nproc    8 \
    --cpu-per-job 1 \
    --exhaustiveness 16 \
    --seed 42 \
    --top-n 20
```

关键参数：

| 参数 | 作用 | 建议 |
| --- | --- | --- |
| `--ligands` | 配体来源，目录 / glob / 文件都行 | 目录最方便 |
| `--out-dir` | 输出姿态和日志 | 别和输入混在一起 |
| `--csv` | 汇总结果 | 用于后续排序 |
| `--nproc` | 同时跑几个配体 | 设成 CPU 核数 ÷ `--cpu-per-job` |
| `--cpu-per-job` | 每个 Vina 内部用几核 | 分子多时用 1，分子少时调大 |
| `--exhaustiveness` | 搜索强度 | 试跑用 8–16，正式筛选用 32 |
| `--seed` | 随机种子 | 结果要可复现就固定它 |
| `--scoring` | 打分函数 | 默认 `vina`；`ad4` 需先算 affinity maps |
| `--summarize-only` | 只汇总不重跑 | 断点续跑、重新排序时用 |

产出：

- `poses/<配体名>_out.pdbqt`：每个配体的对接姿态（含打分）
- `poses/<配体名>.log`：Vina 原始日志
- `results.csv`：按最优打分升序排列的汇总表，字段为
  `rank, ligand, best_affinity_kcal_mol, rmsd_lb, rmsd_ub, n_poses, out_file`

## 四、结果怎么读

`results.csv` 第一列 `best_affinity_kcal_mol` 数值越负，表示预测结合越强。
但要注意：

1. 打分只能用来**排序**，不是真实的结合自由能，绝对值别当 Kd。
2. `vina` 与 `ad4` 两套打分函数**不能横向比较**，筛选全程只用一套。
3. 光看打分容易假阳性，建议叠加两条过滤：关键残基相互作用（氢键、
   疏水接触）是否命中，以及配体本身有没有可疑官能团（PAINS）。
4. 想更可信，挑前几十个候选做 MM-GBSA 或短 MD 复算，见主笔记后处理一节。

可视化流程：`mk_export.py poses/<配体>_out.pdbqt -s pose.sdf` 转成 SDF，
再用 PyMOL 打开受体 + pose.sdf 看结合模式。

## 五、断点续跑与集群

- 单机跑到一半中断：直接重跑同一条命令即可，脚本会跳过已有结果的逻辑
  需要你手动确认——更稳的做法是加 `--out-dir` 不变，重跑后用
  `--summarize-only` 只汇总。
- 集群（Slurm/PBS）上更推荐用作业数组：每个配体一个任务，
  用 `--cpu-per-job` 控制核数，最后 `--summarize-only` 汇总。
- 配体规模上万时，换 AutoDock-GPU（https://github.com/ccsb-scripps/AutoDock-GPU）
  或 GNINA（https://github.com/gnina/gnina）更划算。