gmx_MMPBSA

Calculate binding free energies with MM/PBSA and MM/GBSA

Docs

Input

MD Trajectory file

Topology file

5000 credits

Output

Configure input settings on the left, then click "Submit job"

What is gmx_MMPBSA?

gmx_MMPBSA calculates binding free energies from molecular dynamics trajectories using MM/PBSA (Molecular Mechanics/Poisson-Boltzmann Surface Area) and MM/GBSA (Molecular Mechanics/Generalized Born Surface Area) methods. Developed by Valdés-Tresanco and colleagues, it bridges GROMACS simulations with AMBER's MMPBSA.py engine, enabling end-state free energy calculations for protein-ligand, protein-protein, and protein-nucleic acid complexes.

These methods estimate binding affinity by calculating the free energy difference between bound and unbound states. The calculation combines molecular mechanics energies (van der Waals, electrostatic) with implicit solvation models, offering a balance between the speed of empirical scoring functions and the accuracy of rigorous alchemical methods.

How MM/PBSA and MM/GBSA work

Binding free energy ( $\Delta G_{\text{bind}}$ ) is computed from snapshots extracted from an MD trajectory:

$\Delta G_{\text{bind}} = G_{\text{complex}} - G_{\text{receptor}} - G_{\text{ligand}}$

Each term comprises several energy components:

Component	Description
$\Delta E_{\text{vdW}}$	Van der Waals interactions (Lennard-Jones potential)
$\Delta E_{\text{elec}}$	Coulombic electrostatic interactions
$\Delta G_{\text{polar}}$

MM/PBSA solves the Poisson-Boltzmann equation for electrostatic solvation, providing higher accuracy but at greater computational cost. MM/GBSA uses analytical Generalized Born approximations, running faster while maintaining reasonable accuracy for ranking compounds.

Generalized Born models

Model	Description
`GB-HCT (igb=1)`	Hawkins-Cramer-Truhlar model, fast but less accurate
`GB-OBC1 (igb=2)`	Onufriev-Bashford-Case model, variant 1
`GB-OBC2 (igb=5)`	OBC variant 2, recommended for most applications
`GB-Neck (igb=7)`	Improved treatment of interstitial regions
`GB-Neck2 (igb=8)`	Further refinements to neck region calculations

The OBC2 model (igb=5) provides a good balance of accuracy and speed for protein-ligand systems.

How to use gmx_MMPBSA online

ProteinIQ provides cloud-based access to gmx_MMPBSA, eliminating the need to configure AMBER, GROMACS, or compute environments locally.

Inputs

Input	Description
`MD Trajectory file`	GROMACS trajectory (`.xtc`, `.trr`) or multi-model PDB. Contains atomic coordinates across simulation frames.
`Topology file`	GROMACS topology (`.tpr` or `.top`). Defines force field parameters and system topology.
`Index file`	Optional GROMACS index file (`.ndx`). Specifies receptor and ligand atom groups. Auto-detected if omitted.
`Input parameter file`	Optional gmx_MMPBSA input file. Overrides settings with custom parameters.

Settings

Calculation settings

Setting	Description
`Calculation method`	`MM/GBSA` (faster, suitable for screening), `MM/PBSA` (more accurate), or `Both` for comparison
`Energy decomposition`	When enabled, calculates per-residue energy contributions to identify key binding site residues

Trajectory settings

Setting	Description
`Start frame`	First frame to analyze (1-indexed)
`End frame`	Last frame to analyze. Use `-1` to include all frames
`Frame interval`	Analyze every nth frame. Set to 2 or higher to reduce computation for long trajectories

Solvation settings

Setting	Description
`Salt concentration`	Ionic strength in molar (default 0.15 M, physiological conditions)
`GB model`	Generalized Born variant for MM/GBSA. `GB-OBC2 (igb=5)` recommended

Output

Results include total binding free energy with component breakdown:

Column	Description
`DELTA TOTAL`	Net binding free energy (kcal/mol). More negative indicates stronger binding.
`VDWAALS`	Van der Waals contribution
`EEL`	Electrostatic energy (gas phase)
`EGB` or `EPB`	Polar solvation (GB or PB method)
`ESURF` or `ENPOLAR`	Nonpolar solvation (surface area term)

When energy decomposition is enabled, per-residue contributions identify which amino acids drive binding.

Interpreting binding energies

$\Delta G_{\text{bind}}$ (kcal/mol)	Interpretation
< −15	Very strong binding
−10 to −15	Strong binding (nanomolar range)
−5 to −10	Moderate binding
> −5	Weak binding

These values should be interpreted as relative rankings rather than absolute affinities. MM/PBSA and MM/GBSA excel at comparing similar ligands binding to the same target, but absolute values carry significant uncertainty.

Applications

Lead optimization: Compare binding energies across congeneric compound series to guide medicinal chemistry
Binding hotspot identification: Energy decomposition reveals which residues contribute most to binding, informing mutagenesis studies
Selectivity analysis: Compare binding to on-target versus off-target proteins
Protein-protein interface analysis: Identify key residues stabilizing macromolecular complexes

Limitations

The single-trajectory approach used here assumes identical conformations for complex, receptor, and ligand. This reduces noise but misses conformational reorganization effects. Results depend heavily on the quality and length of the MD simulation—short simulations with poor sampling produce unreliable energies.

Entropy contributions are not calculated by default due to computational cost and large uncertainties. For charged ligands, electrostatic and polar solvation terms often exhibit significant cancellation, amplifying small errors. Results should guide compound prioritization rather than predict absolute binding constants.

Related tools

OpenFE

Run alchemical free energy calculations for drug discovery using Open Free Energy. Supports Absolute Hydration Free Energy (AHFE) and Relative Binding Free Energy (RBFE) calculations with GPU-accelerated OpenMM simulations.

AIMNet2

Predict molecular energies, atomic forces, atomic charges, stress tensors, and Hessians from coordinate-bearing molecular structures with AIMNet2 neural network potentials.

OpenMM

Run GPU-accelerated molecular dynamics simulations using OpenMM. Simulate protein and protein-ligand complex dynamics with industry-standard force fields (AMBER, CHARMM) and OpenFF ligand parameterization.

AutoDock-GPU

GPU-accelerated molecular docking using the AutoDock4 force field. Up to 56x faster than serial AutoDock via CUDA parallelization of the Lamarckian Genetic Algorithm.

AutoDock Vina

AutoDock Vina is a widely-used molecular docking tool that predicts protein-ligand binding modes using physics-based force fields. Fast, reliable, and the gold standard for structure-based drug discovery.

MD Trajectory Analysis

Analyze molecular dynamics trajectories using a ProteinIQ tool pinned to MDAnalysis 2.9.0. Calculate RMSD, residue-aggregated RMSF, radius of gyration, distance tracking, and additional trajectory observables from standard topology and trajectory files.

ORB v3

ORB v3 is a universal interatomic potential (machine learning force field) that predicts energies, forces, and stress tensors for atomic systems. Supports both molecular and materials structures with geometry optimization using conservative and direct model variants.

PandaDock

Open-source molecular docking platform using physics-based scoring functions. CPU-optimized algorithms achieve sub-angstrom accuracy (0.014A RMSD) without GPU requirements.

SMINA

SMINA is a fork of AutoDock Vina with enhanced scoring functions, custom scoring support, and 10-20x faster minimization. Ideal for scoring function development, pose refinement, and high-performance docking workflows.

Boltz-2

Boltz-2 is a biomolecular foundation model for structure and binding affinity prediction. Supports proteins, ligands, DNA, and RNA in multi-component complexes. Automatically scales GPU resources for large complexes. Predicts binding affinity with near-FEP accuracy at 1000x faster speed.

What is gmx_MMPBSA?