Building a System With CHARMM GUI

To start a molecular dynamics simulation, you must first build your system. You can certainly follow the procedure as outlined in the Lysozyme in Water tutorial, but that can get a little tedious. As I am lazy and have a fondness for GUIS, I tend to prefer to build systems with CHARMM-GUI, which is an online interface that will solvate a particular molecule, choose a force field, make any modifications you desire, and generate the equilibration files for whatever MD engine you are using (which we will modify, as described in MDP options and Energy Minimization and Equilibration). Also, as a general note on good documentation, PLEASE write down all of the options that you choose for this process - it will make your life a lot simpler when writing the methods section of any resultant paper. Check the methods sections of previous MD papers we have published for the types of details that are required. 

• First, navigate to the Solution builder by clicking Input Generator → Solution Builder in the sidebar. You should load a .pdb file, normally retrieved from the RCSB database. In this case, let's load the trp-cage mini-protein, arguably the Dinner group's favorite test system. You can either directly upload the PDB, or simply enter 1L2Y into the text box. Proceed to the next step.
• You should see the following screen: 

Select model 1 (38 NMR structures were uploaded to this entry, let's just worry about the first one). The next field is selecting what parts of the pdb file you want to consider. With more complex systems, you will see many more entries. You might see an entry for each chain of the protein, or for waters/solvent that were resolved in a crystal structure, or for metallic heteroatoms, or for FAB fragments that aid in crystallization, etc. What you select depends on what you are studying. Make sure that you include everything relevent fo your system - for example, phenol and zinc ligands are both important for the insulin hexamer, but the fab fragments are not. 

Here, there's only one thing to select, since this is just a small, 20 amino acid chain. Make sure PROA is checked and move on to the next step. 

•  This page includes all of the manipulations you can make to change your structure from the PDB file. Again, what you do on this screen depends on what you're studying. You could cap modify/mutate specific residues, change protonation states to reflect a specific pH, add/remove cysteine bonds, etc. This is an important screen to fully consider - many a Dinner group member has wasted a month or two because they forgot to insert the correct cysteine bonds. For changing the pH, you can determine the likely local pKa of specific residues by using PROPKA, and change protonation states accordingly. For now, just ensure that Preserve hydrogen corrdinates is checked so that the hydrogens are not remade, and proceed to the next step. 
•  The next page should look as follows:

It is good practice to frequently check your structure by downloading the linked files at the top, to ensure that the process is proceeding as you would like it to. After looking at the structure, this page is primarily focused on choosing the size of the water/ion box. Again, your specific choice here strongly depends on the system you are studying. Basically all MD codes you will run will have periodic boundary conditions - you want to choose a box large enough that your protein does not interact with itself across these periodic boundaries. Electrostatic interactions, in particularly, are quite long-ranged.  However, you don't want too large of a box, or else your simulations will take forever to actually run. Specifically, due to the Particle Mesh Ewald summations, computation time scales as _Nlog(N) _with N being the total number of atoms in your system. As you will soon see, the vast majority of your system will be solvent, the number of which scales with box volume (or cubicly with box length). Thus, this balance can often times be a little bit of guess-and-check.

You want your box to be big enough so that, at the largest extension that you anticipate your molecule achieving, you see no pathological self-interaction. Thus, the box size needed depends on the process your studying, not just the initial size of the molecule. For example, if we were looking at the unfolding of trp-cage, then we would need a much larger box, as the unfolded protein can be quite long - indeed, the box size John Strahan needed to study this unfolding process was approximately (5.0 nm)3. The default unit for CHARMM-GUI is Angstroms. We won't be looking at unfolding, so we can choose a smaller box size. If you aren't expecting to see much expansion in your molecule size, it is normally sufficient to fit your waterbox to your protein size, and set the edge distance to 10 A. This will give 10 angstroms of space between the end of the protein and the box in all dimensions (as a note, 10 A might be a bit of overkill, but since this is a small system anyway, we can get away with wasting a bit of computation time). Choosing this should give you an overall box size of (4.4 nm)3. The other thing to choose on this screen is the ion concentration. In general, the default is to ensure overall charge neutrality and to bring up the total concentration to 150 mM KCl. You can, of course, change this, depending on the process you are studying. Monte Carlo placement is preferred. Proceed to the next step - be warned that this may take some time, especially if your system is large. It's OK if you have to exit or do something else - you should just be able to use your Job ID to check to see if it has finished. This solvation took about 5-10 minutes on my laptop. 

• The next screen should show you how your solvation ended - check your structures to make sure everything is right. Check the automatic PME FFT option, then move to the next step. 
• This final screen allows you to choose the force field and MD engine to use. I recommend CHARMM36m over the vanilla CHARMM36, especially if you think even a limited amount of protein disorder could occur in your process of interest. We tend to use GROMACS, so check that box, and choose NVT equilibration. You should probably choose NPT input generation, although it doesn't really matter, as we will just be changing the equilibration procedure anyway (see Energy Minimization, Equilibration, and MDP options). Choose your desired temperature (force fields are best near 300-303.15 K, so it's probably best to choose something in that range, unless you have a compelling reason not to). Move to the next step. 
• Download the resultant files, and upload them to Midway in order to minimize and equilibrate. Or just sit back and admire the very pretty solvated structure you have just created. 

Again, make sure you have written down all of your choices (box size, atom number, ion types and numbers, force field, any mutations/modifications, etc.).