pavel-kravchenko Virtual Screening & ADMET Prediction
When to Use
- Screening a compound library (SMILES/SDF) against a protein target by molecular docking
- Preparing a receptor (PDB) and ligands (SDF) for AutoDock Vina
- Filtering docking hits with ADMET rules (Lipinski Ro5, LogS, CYP inhibition, hERG risk)
- Combining docking score + QSAR prediction + ADMET pass/fail into one ranked hit list
- Building a structure-based virtual screening pipeline for hit-to-lead prioritization
Version Compatibility
AutoDock Vina ≥1.2, meeko ≥0.5 (ligand/receptor PDBQT prep), RDKit ≥2023.09, pandas ≥2.0, Python ≥3.10. Vina 1.2+ takes PDBQT for both receptor and ligand (no separate prepare_receptor4.py/prepare_ligand4.py MGLTools scripts needed).
Prerequisites
pip install vina meeko rdkit pandas(orconda install -c conda-forge vina)- A cleaned receptor structure (waters/ligands stripped) and a defined binding-site box
- Familiarity with
bio-chemoinformatics-molecular-io(SMILES/SDF handling) andbio-structural-biology-structure-io(PDB parsing)
Protein Structure Preparation
Goal: get a receptor PDBQT with hydrogens added and no waters/crystallographic ligands.
Approach: download the PDB, strip heteroatoms, add hydrogens with PDBFixer, then convert to PDBQT with meeko/OpenBabel. The classic Vina tutorial target is PDB 1IEP — the Abl kinase (ABL1) tyrosine kinase domain in complex with imatinib — not EGFR.
# Download and fix the Abl kinase-imatinib structure (PDB 1IEP)
wget -q https://files.rcsb.org/download/1IEP.pdb
pdbfixer 1IEP.pdb --output 1IEP_fixed.pdb --add-hydrogens --remove-heterogens --keep-water=False
# Convert receptor PDB to PDBQT (Vina 1.2+ / meeko)
mk_prepare_receptor.py -i 1IEP_fixed.pdb -o 1IEP_fixed --box_size 20 20 20 \
--box_center 22.5 5.0 18.0
from Bio.PDB import PDBParser, PDBIO, Select
class HetatmWaterSelect(Select):
"""Keep only standard amino-acid residues; drop waters and heteroatoms."""
def accept_residue(self, residue):
return residue.id[0] == " " # blank hetfield = standard residue
def strip_to_protein_only(pdb_in: str, pdb_out: str) -> None:
"""Write a protein-only PDB, removing waters/ligands/ions before docking prep."""
parser = PDBParser(QUIET=True)
structure = parser.get_structure("receptor", pdb_in)
io = PDBIO()
io.set_structure(structure)
io.save(pdb_out, select=HetatmWaterSelect())
Molecular Docking with AutoDock Vina
Goal: dock a ligand (or library of ligands) into the prepared binding site and get a binding-affinity score (kcal/mol; more negative = better predicted binding).
Approach: convert each ligand SDF/SMILES to PDBQT with meeko, then run Vina against the receptor with a box centered on the known active site.
# Ligand prep: SDF -> PDBQT (adds Gasteiger charges, sets rotatable bonds)
mk_prepare_ligand.py -i imatinib.sdf -o imatinib.pdbqt
# Run Vina docking against the Abl kinase ATP-binding site
vina --receptor 1IEP_fixed.pdbqt --ligand imatinib.pdbqt \
--center_x 22.5 --center_y 5.0 --center_z 18.0 \
--size_x 20 --size_y 20 --size_z 20 \
--exhaustiveness 8 \
--out imatinib_docked.pdbqt --log docking.log
grep "REMARK VINA RESULT" imatinib_docked.pdbqt | head -1
from vina import Vina
def dock_ligand(receptor_pdbqt: str, ligand_pdbqt: str, center: tuple, box_size=(20, 20, 20)) -> float:
"""Dock one ligand with the Vina Python API; return the best (lowest) binding score in kcal/mol."""
v = Vina(sf_name="vina")
v.set_receptor(receptor_pdbqt)
v.set_ligand_from_file(ligand_pdbqt)
v.compute_vina_maps(center=list(center), box_size=list(box_size))
v.dock(exhaustiveness=8, n_poses=10)
energies = v.energies(n_poses=1) # [[total, inter, intra, torsion, intra_best]]
return float(energies[0][0])
ADMET Property Prediction
Goal: flag compounds with poor drug-likeness or predicted toxicity before/after docking.
Approach: compute cheap RDKit-based Lipinski Ro5 descriptors as a first-pass filter; use SwissADME (web) or a DeepChem Tox21 model for hERG/CYP/toxicity predictions on the survivors.
from rdkit import Chem
from rdkit.Chem import Descriptors, Lipinski
def passes_lipinski(smiles: str, max_violations: int = 1) -> bool:
"""Return True if a SMILES passes Lipinski's Rule of Five (allowing max_violations)."""
mol = Chem.MolFromSmiles(smiles)
if mol is None:
return False
violations = sum([
Descriptors.MolWt(mol) > 500,
Descriptors.MolLogP(mol) > 5,
Lipinski.NumHDonors(mol) > 5,
Lipinski.NumHAcceptors(mol) > 10,
])
return violations <= max_violations
Hit Prioritization
Goal: merge docking score, QSAR-predicted pIC50, and ADMET pass/fail into one ranked shortlist.
Approach: build a pandas DataFrame per compound, filter on admet_pass, then rank by a composite score (flip the sign on vina_score since more negative = better binding, so higher composite = better hit).
import pandas as pd
def prioritize_hits(smiles_library, docking_scores, qsar_predictions, admet_flags, top_n=10) -> pd.DataFrame:
"""Combine Vina score + QSAR pIC50 + ADMET pass/fail into a ranked top-N hit table."""
hits = pd.DataFrame({
"SMILES": smiles_library,
"vina_score": docking_scores, # kcal/mol, more negative = stronger predicted binding
"qsar_pIC50": qsar_predictions, # higher = more potent
"admet_pass": admet_flags, # bool: passed ADMET filter
})
hits_filtered = hits[hits["admet_pass"]].copy()
hits_filtered["composite"] = -hits_filtered["vina_score"] + hits_filtered["qsar_pIC50"]
top_hits = hits_filtered.nlargest(top_n, "composite")
return top_hits[["SMILES", "vina_score", "qsar_pIC50", "composite"]]
Pitfalls
- Docking score ≠ binding affinity: Vina scores correlate weakly with experimental Kd/IC50; always validate top hits experimentally before further investment
- Receptor flexibility: rigid docking misses induced-fit binding; use ensemble docking (multiple receptor conformations) for flexible targets
- PAINS compounds: filter pan-assay interference structures (
rdkit.Chem.FilterCatalog) before scoring — they show up as false-positive hits across unrelated assays - Box placement: an incorrectly centered/sized docking box silently produces poses outside the real binding site with plausible-looking scores
- Protonation state: skipping hydrogen addition/pH-appropriate tautomers on both receptor and ligand biases scores and pose geometry
See Also
bio-chemoinformatics-virtual-screening— ligand-based screening (similarity, pharmacophore, QSAR) as an alternative/complement to dockingbio-chemoinformatics-admet-prediction— deeper ADMET modeling (LogS, CYP, hERG) beyond Lipinski filtersbio-structural-biology-structure-io— general PDB/mmCIF parsing and structure cleanupbio-chemoinformatics-molecular-io— SMILES/SDF reading, writing, and format conversion