Path B - Module 42: Molecular Mechanics (Energies)#

How do we know if your new disulfide bridge actually stabilizes the enzyme? Geometry only tells us distances, but Potential Energy tells us the physical stability.

In this module, you will learn to use the molecular mechanics engine of MolSysMT to calculate the energy profile of your industrial biocatalyst.

import molsysmt as msm
from molsysmt import systems

# Load our engineered complex
petase = msm.convert('pdb:6EQE', to_form='molsysmt.MolSys', selection='molecule_type=="protein"')
msm.build.add_missing_hydrogens(petase, pH=8.0)

bhet = msm.convert('C1=CC(=CC=C1C(=O)OCCO)C(=O)OCCO', from_form='string:smiles', to_form='molsysmt.MolSys')
msm.build.add_missing_hydrogens(bhet)

molsys = msm.merge([petase, bhet])

1. Calculating Global Potential Energy#

The function get_potential_energy() uses a force field engine (OpenMM by default) to compute the energy of the current coordinates.

# Calculate the total energy of the complex
energy = msm.molecular_mechanics.get_potential_energy(molsys)

print(f"Total Potential Energy: {energy}")

2. Interaction Energy (Binding Strength)#

For an enzyme engineer, the most important number is how much energy is released when the substrate binds. We call this the Interaction Energy.

# Calculate the non-bonded potential energy between enzyme and substrate
interaction = msm.molecular_mechanics.get_non_bonded_potential_energy(molsys, 
                                                                   selection='molecule_type=="protein"', 
                                                                   selection_2='molecule_type=="small molecule"')

print(f"Enzyme-Substrate Interaction Energy: {interaction}")

3. Stress Mapping (Forces)#

Where is the enzyme most “stressed”? Large forces indicate regions that want to move or change conformation.

# Get the forces on every atom
forces = msm.molecular_mechanics.get_forces(molsys)

print(f"Max force acting on the system: {forces.max()}")

🏆 Path B Challenge: The Energy Auditor#

  1. Calculate the energy of the Wild-type PETase vs your Mutated PETase.

  2. Did the potential energy decrease (become more negative)?

  3. Identify the atom index that has the absolute highest force acting on it.

  4. Use msm.get_label() to see which residue that atom belongs to.

You have moved from geometry to physics! In Module 43, we will use these energies to relax the system: Energy Minimization.