Modeling Zinc Enzyme System using the 12-6 LJ Nonbonded Model

IMPORTANT INFORMATION:

Pengfei Li note that the following files in AMBER 2016 to AMBER 2019 have errors about the parameters for +3 and +4 ions: (1) frcmod.ions234lm_hfe_tip3p, (2) frcmod.ions234lm_iod_tip3p, (3) frcmod.ions234lm_126_tip3p, and (4) frcmod.ions234lm_hfe_tip4pew.

These errors are:

(1) the file frcmod.ions234lm_hfe_tip3p should have the HFE parameter set for the +3 and +4 ions in the TIP3P water model, but it has the IOD parameter set for the +3 and +4 ions in the TIP3P water model.

(2) the file frcmod.ions234lm_iod_tip3p should have the IOD parameter set for the +3 and +4 ions in the TIP3P water model, but it has the 12-6-4 parameter set for the +3 and +4 ions in the TIP3P water model.

(3) the file frcmod.ions234lm_126_tip3p should have the IOD parameter set for the +3 and +4 ions in the TIP3P water model, but it has the 12-6-4 parameter set for the +3 and +4 ions in the TIP3P water model.

(4) the file frcmod.ions234lm_hfe_tip4pew should have the HFE parameter set for the +3 and +4 ions in the TIP4P/EW water model, but it has the IOD parameter for the +3 and +4 ions in the TIP4P/EW water model.

Note that only the ion parameters for +3 and +4 ions had errors in these four files, while the ion parameters for +2 ions in these files are correct. Moreover, the ion parameters for +3 and +4 ions in the other ion parameter files (except these four files) are also correct.

Also note that if one wants to use the 12-6-4 model, only the Rmin/2 and epsilon parameters are avaliable in these parameter files, while the C4 terms are still need to be added by the ParmEd program.

The corrected files are attached:

(1) frcmod.ions234lm_hfe_tip3p;

(2) frcmod.ions234lm_iod_tip3p;

(3) frcmod.ions234lm_126_tip3p;

(4) frcmod.ions234lm_hfe_tip4pew.

For users of AMBER2016 to AMBER2019, please use these files to replace the four old files in your AMBER software package under the directory: $AMBERHOME/amber/dat/leap/parm/.

Pengfei Li thank Paulius Kantakevicius for he finding these mistakes. Pengfei Li apologizes for these mistakes he made when merging the ion parameters in 2016 and any inconvicences these mistakes brought to the community. The fortunate part out of all the misfortunes is that the Lennard-Jones parameters in the IOD and 12-6-4 parameter sets are not far away from each other (The Rmin/2 values of a certain ion in the two parameter sets are usually within 0.1 Angstrom from each other), and the HFE parameters for the +3 and +4 ions are usually not used in MD simulations because their big errors for simulating the IODs (ion-oxygen distances) of highly charged ions.

Here we describe the zinc enzyme system (PDB ID: 4V2Y) as an example for modeling an atomic ion using the 12-6 LJ nonbonded model:

1. Prepare the PDB file (effort: several minutes):

Here is the original PDB: 4V2Y.pdb. There are no hydrogen atoms in the file. Even tleap could automatically add them in the modeling step, here we use webserver H++ to add hydrogen atoms, which is a more accurate way to determine the protonation state of residues. Another example could be seen in MCPB.py tutorial.

Clean the PDB file:

awk '$1=="ATOM" || $1=="HETATM" || $1=="TER" || $1=="END"' 4V2Y.pdb > 4V2Y_clean.pdb

Get chain A of the PDB file:

awk 'substr($0, 22, 1) == "A"' 4V2Y_clean.pdb > 4V2Y_clean_A.pdb

Delete the ligand in the PDB file:

awk 'substr($0, 18, 3)!="EF2"' 4V2Y_clean_A.pdb > 4V2Y_clean_A_2.pdb

Add H atoms by using H++ websever: we will upload the 4V2Y_clean_A_2.pdb, use the default settings and obtain the topology and coordinate files for AMBER: 0.15_80_10_pH6.5_4V2Y_clean_A_2.top and 0.15_80_10_pH6.5_4V2Y_clean_A_2.crd. Then we can use these files to generate the PDB file by using ambpdb, for this PDB file has residue names consistent with the AMBER style.

ambpdb -p 0.15_80_10_pH6.5_4V2Y_clean_A_2.top -c 0.15_80_10_pH6.5_4V2Y_clean_A_2.crd > 4V2Y_clean_A_Hpp.pdb

Since H++ will delete all the ions and waters, so we will use awk to pull out the zinc ion and combine it back to the new protein PDB file:

awk 'NR==FNR {if ($1!="END") print $0} NR>FNR {if ($4=="ZN") print $0}' 4V2Y_clean_A_Hpp.pdb 4V2Y_clean_A_2.pdb > 4V2Y_A_prep.pdb

Then we can check the protonation states of the metal site residues and make fixation if necessary. We can see the four CYS residues (with residue numbers as 6, 9, 72, and 75) ligating to zinc ions are all protonated in the newly generated PDB file through a visuialization tool such as VMD. We can manually modify these four CYS residues to CYM residues by deleting the HG atoms inside them and rename their residue names from "CYS" to "CYM". Afterwards we can use pdb4amber to renumber the PDB file.
pdb4amber -i 4V2Y_A_prep.pdb -o 4V2Y_A_fixed.pdb

2. Perform the modeling steps using leap (effort: several minutes):

We need to create an input file manually. Here is the input file: Zn_tleap.in.

We are using the SPC/E water model and the compromise (CM) parameter set for Zn2+: we source the leaprc.water.spce file which will load the frcmod.ions234lm_126_spce file that has the CM parameter set for Zn2+ for the SPC/E water model.

Here is the explanation of the tleap input file:

source leaprc.protein.ff14SB #Source the ff14SB force field for protein source leaprc.water.spce #Source the SPC/E water model and related metal ion force fields mol = loadpdb 4V2Y_A_fixed.pdb #Load the PDB file solvatebox mol SPCBOX 8 #Solvate the system using SPC water box addions mol Cl- 0 #Neutralize the system using Cl- ions savepdb mol Zn_4V2Y.pdb #Save the pdb file saveamberparm mol Zn_4V2Y.prmtop Zn_4V2Y.inpcrd #Save the topology and coordiante files quit #Quit tleap

The atomic_ions.lib file, which contains the library of atomic ions, is loaded when we source the leaprc.water.spce file. Herein we perform "addions mol Cl- 0" command in the input file to neutralize the system with Cl- ions. If the total charge of the protein system is negative, we can use "addions mol Na+ 0" command in the input file to neutralize the system with Na+ ions. If you create a lib file for metal ion(s) by yourself, please use a different name from the "atomic_ions.lib" file, preventing conflict with the existing file. The above example is for AMBER 2016. Note that the names of leaprc files for different force fields and frcmod files for ions may change between different AMBER versions, hence some commands in above input file may need to be changed or deleted. Users can check the Specifying which force field you want in LEaP and Ions sections in the Molecular mechanics force fields chapter in the corresponding AMBER manual and $AMBERHOME/dat/leap/cmd/ directory for more details. Here we can perform the command:

tleap -s -f Zn_tleap.in > Zn_tleap.out

You can also use the addions/addions2 command (for example "addions mol ZN 3") if you want to add several Zinc ions to the system before neutralizing it. After performing the command listed above, we will get the topology and coordinate files for doing minimization and molecular dynamics simulation: (Zn_4V2Y.prmtop and Zn_4V2Y.inpcrd)

3. Check the topology files:

Before performing the minimization and molecular dynamics simualtion, we can use the "printLJMatrix" command in ParmEd to check the metal ion related parameters.

Here is an input file for ParmEd: 126_parmed.in. Here we use ":ZN" (ZN is the residue name of zinc ion in the 4V2Y_clean_A_2.pdb file) to represent the zinc ion. Additional information could be found in the Amber Masks section in the chapter of Atom and Residue Selections in the AMBER manual. We can perform the command:

parmed -i 126_parmed.in -p Zn_4V2Y.prmtop

We can see that charge of zinc ion is 2.0, its LJ Radius is bigger than 1.0 Angstrom, and the A and B coefficients between it and other atoms are normal (only A and B coeffients between zinc ion and HO/HW atom type are zero due to the VDW parameters of HO/HW is zero).

Afterwards you can use the topology and coordinate files to perform molecular dynamics simulations.