Chapter 3: Search Engines & Scoring Functions: Monte Carlo, Genetic Algorithms & Empirical Potentials
Chapter 3: Search Engines & Scoring Functions: Monte Carlo, Genetic Algorithms & Empirical Potentials
Conformational search algorithms: Systematic, Monte Carlo simulated annealing, and Lamarckian Genetic Algorithm (LGA); scoring function classes: physics-based force fields, empirical AutoDock Vina potentials, knowledge-based statistical potentials, and deep learning scoring.
1. Conformational Search Algorithms
A typical drug molecule possesses 3 translational degrees of freedom ($x, y, z$), 3 rotational degrees of freedom ($\theta, \phi, \psi$), and $N$ rotatable torsional angles ($\chi_1, \chi_2, \dots, \chi_N$). This creates a combinatorial hyper-dimensional conformational search space with $6 + N$ dimensions:
- Systematic Search: Rotates every torsion angle through discrete increments (e.g., 30°). Highly accurate for small, rigid fragments ($N \le 3$), but suffers from combinatorial explosion for drug-like molecules ($N \ge 6$).
- Monte Carlo (MC) Simulated Annealing: Applies random spatial perturbations. Moves that decrease binding energy are unconditionally accepted; moves that increase energy are accepted with Boltzmann probability $P = \exp(-\Delta E / k_B T)$. Temperature $T$ is gradually lowered to escape local energy minima.
- Lamarckian Genetic Algorithm (LGA): Utilized in AutoDock 4. Treats ligand parameters as a chromosome (genes = coordinates and torsion angles). Population undergoes crossover, mutation, and selection. Unlike Darwinian GA, LGA applies local phenotypic energy minimization (Solis-Wets search) and writes the optimized traits directly back into the genome.
- Iterated Local Search (AutoDock Vina): Employs Broyden-Fletcher-Goldfarb-Shanno (BFGS) quasi-Newton gradient optimization combined with successive stochastic mutations. Significantly faster than traditional LGA while delivering superior convergence.
2. The Taxonomy of Scoring Functions
1. Physics-Based / Force Fields
Sum of molecular mechanics potential energy terms:
Examples: AMBER, CHARMM, DOCK. High physical rigor, but computationally demanding due to solvent dielectric screening estimations.
2. Empirical Scoring Functions
Weighted sum of interaction terms calibrated by multiple linear regression against experimental $K_d / K_i$ complexes:
Examples: AutoDock Vina, ChemScore, Glide SP. Balances speed with empirical accuracy; standard choice for high-throughput screening.
3. Knowledge-Based Potentials
Inverse Boltzmann statistical analysis of interatomic distance distributions across tens of thousands of PDB crystal complexes:
Examples: DrugScore, PMF, ITScore. Implicitly accounts for entropic and solvent effects without explicitly modeling individual waters.
4. Machine Learning & 3D-CNNs
Deep learning architectures (3D convolutional neural networks, graph neural networks) trained on spatial atom density voxel grids:
Examples: gnina, RF-Score, OnionNet. Outperforms classical empirical functions in ranking affinity and pose classification.
3. AutoDock Vina Scoring Function Deconstructed
AutoDock Vina calculates total free energy as a conformation-dependent interatomic potential plus a conformation-independent rotatable bond penalty:
c = ∑i < j fti, tj(rij) / (1 + wrot · Nrot)
The interaction functions $f_{t_i, t_j}$ comprise Gaussian steric terms (Gauss 1 and Gauss 2), a quadratic repulsion penalty for steric clashes ($r < 0$), a hydrophobic step function, and a linear directional hydrogen bonding term. The denominator $(1 + w_{\text{rot}} \cdot N_{\text{rot}})$ penalizes flexible ligands with high rotatable bond counts due to loss of conformational entropy upon binding.
Common Pitfalls & Troubleshooting Checklist
- Java Heap Exhaustion (GATK / Picard): If processes terminate with
java.lang.OutOfMemoryError, set your environment variable explicitly:export _JAVA_OPTIONS="-Xmx16g -Djava.io.tmpdir=/tmp". - Index Timestamp Mismatch: If tools complain that a BAM index is missing or out of date, regenerate with
samtools index -@ 4 sample.bam. - Conda Channel Conflicts: Avoid ABI incompatibilities between bioconda and conda-forge by setting strict priority:
conda config --set channel_priority strict. - Chromosome Naming Contamination: Always verify whether reference FASTA headers use
chr1or1; divergent naming causes silent alignment drops.