Plot Free Energy Landscape (FEL) for Binding

Overview and Methodology

What It Does

This script constructs a 2D Free Energy Landscape (FEL) to characterize the conformational states of two molecular groups (e.g., a receptor and ligand) during a molecular dynamics (MD) simulation. It maps the joint distribution of two selected reaction coordinates — such as Distance, RMSD, Radius of Gyration (Rg), or Number of Contacts — onto a free energy surface.

The resulting FEL plot highlights:

• Stable states as energy minima (dark basins).
• Transition states as ridges or saddle points.
• Unstable or rare conformations as high-energy regions.

How It Works

• Reaction Coordinate Sampling:
  • Computes two metrics (X and Y axes) for each trajectory frame.
• 2D Histogram:
  • Bins the X–Y coordinate pairs into a 2D histogram.
• Free Energy Conversion:
  • Converts the probability density $( P(x, y) )$ to free energy using the Boltzmann relation:

$[ \Delta G(x, y) = -k_B T \cdot \ln P(x, y) ]$

• Normalization:
  • Subtracts the minimum energy value to set the global minimum at $( \Delta G = 0 )$ for interpretability.

---

Configuration and Inputs

Prerequisites

• Requires a loaded trajectory.

Key Configuration Options

• Selections:

  • group1_sel, group2_sel: Atom selections for receptor and ligand (or domains).
  • group1_name, group2_name: Used to label plot axes and files.

• Reaction Coordinates:

  • Distance: Center-of-mass (COM) distance between Group 1 and Group 2.
  • Rg: Radius of gyration of combined atoms from both groups.
  • RMSD: RMSD of the combined system relative to the first frame.
  • Contacts: Number of interatomic contacts (<4 Å) between groups.

• FEL Parameters:

  • bins: Resolution of the 2D histogram.
  • temperature: Used for free energy conversion (Kelvin).

• Plot Style:

  • cmap: Color scheme for the energy surface.
  • title_label, cbar_label: Labels for plot and colorbar.

---

Output

• FEL Plot

  • File: FEL_{group1_name}_{group2_name}_Binding.png
  • Displays:
    • Energy basins representing dominant conformational states.
    • Axes corresponding to selected reaction coordinates.
    • Colorbar showing relative energy in kcal/mol.

• Console Output:

  • Plot saved confirmation.
  • Basic trajectory info.

---

Interpreting the Results

• Global Minima:

  • Indicates the most stable conformation or binding mode.
  • Corresponds to the highest population in the MD simulation.

• Multiple Basins:

  • Suggests conformational heterogeneity or multiple binding modes.

• Saddle Points or Ridges:

  • Mark transition pathways between states.

• Isolated Peaks:

  • Rare or transient conformations, possibly unbound or misbound states.

---

Example Scenarios

Binding Pathway Analysis

• Coordinates: Contacts vs. Distance
• Observation: Basin forms at low distance with high contact number.
• Interpretation: Binding process completes into a stable pocket.

Protein Domain Reorientation

• Coordinates: Rg vs. RMSD
• Observation: Broad FEL surface with multiple shallow minima.
• Interpretation: Conformational plasticity or semi-disordered domains.

Docking Refinement

• Coordinates: RMSD vs. Distance
• Observation: Sharp minimum near low RMSD and low distance.
• Interpretation: Final pose converges toward native complex structure.

---

Usage Tips

• Choice of Coordinates:

  • Use Distance + Contacts to probe binding/unbinding.
  • Use Rg + RMSD for flexibility or shape shifts.
  • RMSD requires stable atom ordering — best used with rigid segments.

• Filtering Water/Solvent:

  • Internal filtering excludes resname TIP, WAT to focus on molecular interactions.

• Bin Count:

  • Use bins=100–200 for smooth FELs.
  • Lower bin count for quick preview.

• Temperature:

  • Set to match simulation (typically 300–310 K).

• Color Maps:

  • Use viridis or inferno for high-contrast FELs.
  • coolwarm or RdBu_r for symmetrical energy spreads.

---