Building and running a molecular dynamics (MD) simulation of a protein-ligand system on the PlayMolecule® web platform [TUTORIAL]

In this article we will look at how to prepare, parameterize, build and run a molecular dynamics (MD) simulation of a protein-ligand system using the PlayMolecule platform

Gerard Martínez
Nov 4 · 7 min read

Welcome to another PlayMolecule tutorial. In this occasion we will look at one of the most powerful features of PlayMolecule: our pipeline for building and running MD simulations. As you will see, the process is quite straightforward and thanks to PlayMolecule you should be able to do most of the actions required to run an MD simulation directly from your browser.

PlayMolecule consists of a set of applications that have the capability to “talk to each other” — this is, the output of one application can be used as input for another application. For instance, the SystemBuilder application requires a prepared protein structure that can be obtained through the application ProteinPrepare. Therefore, in order to build and run an MD simulation, we must use a number of applications in succession. Here’s an outline of the pipeline we will use in this tutorial:

Figure 1. Pipeline to build and run an MD simulation

You can either run this tutorial from the public account (no login required) or in your personal account by logging in (top right corner) — in the latter case your structures will not be shown to the public. Either way, you must start and end the tutorial with same setup since the jobs from the public account are not visible from the private account and vice-versa.

As an example, in this tutorial we will build the well-known trypsin-benzamidine system. Without further ado, let’s learn how to do it.


1. ProteinPrepare: preparing your protein structure

Essentially, by using ProteinPrepare, we will input our initial structure and we will obtain a protonated and optimized structure that we will be able to use in SystemBuilder in a later step. To build our trypsin-benzamidine system, we will input the 3PTB RSCB identifier into ProteinPrepare and we will specifically select to use the chain A.


2. Parameterize: generate parameters for your ligand

Parameterize is quite complex (you can read the full publication here) but in essence it takes care of:

  1. Computing the partial charge of the molecule atoms
  2. Assigning general parameters (using the GAFF2 forcefield)
  3. Optimizing dihedral parameters

The third step is particularly important because general force-fields like GAFF2 are not really good at guessing the dihedral angle terms.

Figure 2. Dihedral angle. Obtained from http://cbio.bmt.tue.nl/pumma/index.php/Theory/Potentials

Parameterize can obtain the dihedral terms following one of the modes below:

  1. Simple GAFF2 guessing (quick but inaccurate)
  2. Dihedral scanning and parameter fitting using QM calculations (slow but accurate)
  3. Dihedral scanning and parameter fitting using the ANI-1x neural network potential (a balance between speed and accuracy)

You can select which mode you want to use from the graphical user interface. Beware that ANI1-x only accepts C,H,N,O atom types and both ANI-1x and QM modes have a max. limit of 25 atoms in the (free) public website.

Also bare in mind that while Parameterize takes care of generating parameters, it does not take care of “preparing” your ligand, this is it does not take care of its protonation. Actually, in Parameterize you can draw or input a ligand without hydrogens and they will be guessed using RDKit, but its accuracy is dubious and we recommend using external software to come up with the right protonation state for your ligand.

Another remark: if you are planning to build a protein-ligand system with the ligand bound to the binding site, make sure the coordinates of the ligand are the correct (binding) ones. In this tutorial we assume you have the crystallographic or docking pose for your ligand.

Figure 3. We download the SDF file from RSCB from PDB.

In our example, we will simply download and use the SDF of the benzamidine ligand directly from the 3PTB entry in RSCB: https://www.rcsb.org/structure/3PTB and rely on Parameterize for the protonation state. We will use the ANI-1x dihedral mode. Bare in mind these SDF files usually have multiple copies of the ligands for different chains but Parameterize will only use the first one, which in this case corresponds to the benzamidine that binds the chain A.

Figure 4. The parameters for Benzamidine obtained using Parameterize.

3. SystemBuilder: building your protein-ligand system

Figure 5. We select our protein and ligand from the drop-down menus.

In order to include the protein we previously prepared with ProteinPrepare and the ligand we just parameterized, we must select the jobs we previously executed from the drop down lists in the input form. For instance, in our case we select the protein 3ptb and the ligand we identified as Benzamidine.

We leave the other options by default since we do not want to randomize the position of the ligand in the solvent (we want it to stay in the binding pocket as the initial coordinates provided) and we will leave the ionization to ‘0’, which means that ions will be added simply to neutralize the charge of the whole system to 0.

Figure 6. Output of SystemBuilder: our MD-ready protein-ligand system.

4. SimpleRun: running a short MD simulation

  1. Minimization (forces are minimized, clashes are solved, etc.)
  2. Equilibration run (MD is run for few nanoseconds with constraints so that the system doesn’t blow up)
  3. Production run (MD without constraints)
Figure 7. Input options for SimpleRun.

To run SimpleRun, simply select the SystemBuilder job from the drop-down menu, select the length of the simulation and leave the rest of the options by default. Make sure that the ligand selected is the correct one and the constraints applied are “protein-ligand”. This is important because during equilibration we want the MD engine to apply constraints to our ligand as well (otherwise the ligand would be probably forced out of the pocket due to initial extreme forces).

Once you have run SimpleRun, you will be able to see an output like the one below. You can visualize the simulation by clicking on the “Play button”. You can also see some plots regarding the fluctuation of the ligand as well as a PDF report with all the information condensed, such as the previous plots, plus a protein-ligand interaction map (obtained with the PlexView application) and other useful info. Additionally, three different types of clustering are applied to generate “representative” snapshots of your simulation.

Figure 8. SimpleRun output example for the 3UZA system.

Conclusions

Thank you for reading!

PlayMolecule

PlayMolecule blog. News, tutorials and use cases of PlayMolecule.

Gerard Martínez

Written by

Manager and developer of the PlayMolecule drug discovery platform (https://playmolecule.com/)

PlayMolecule

PlayMolecule blog. News, tutorials and use cases of PlayMolecule.

Welcome to a place where words matter. On Medium, smart voices and original ideas take center stage - with no ads in sight. Watch
Follow all the topics you care about, and we’ll deliver the best stories for you to your homepage and inbox. Explore
Get unlimited access to the best stories on Medium — and support writers while you’re at it. Just $5/month. Upgrade