sietraj
sietraj is an alternative to the
AMBER MM-PBSA software. Virtual alanine mutations are also possible. Solvated interaction energies (SIE) are calculated using parameters that have been fitted to reproduce binding free energies of a data set of 99 protein-ligand complexes (Naim et al., 2007). Note that the underlying physics behind the electrostatic components of SIE and MM-PBSA is the same. Hence, all the caveats of implicit solvation apply. The main differences lie in the choice of parameters, surface generation method (Bhat & Purisima, 2007; Chan & Purisima, 1998) and Poisson solver (Purisima & Nilar, 1995; Purisima, 1998).
Click
here for download and installation instructions.
Sample run 1
Two sample prmtop and short trajectory files can be found in $SIEHOME/examples, where $SIEHOME is the path to the sietraj2 folder that was installed. bcl2.mdcrd is a trajectory file with 11 snapshots with all water molecules stripped away. This is a protein-ligand complex with 2611 protein atoms and 43 ligand atoms.
sietraj is invoked as shown below. The output is displayed on the screen as well as written to a file. The screen output is rounded to 3 decimal places. The output file, /tmp/bcl2.out, retains full precision. The output columns are: snapshot frame number, binding free energy, reaction field energy (binding), reaction field energy (complex), reaction field energy (target), reaction field energy (ligand), intermolecular vdW, intermolecular coulomb, change in surface area. Energies are in kcal/mol and areas are in Å
2.
prompt% cd $SIEHOME/examples
prompt% $SIEHOME/bin/sietraj -pt bcl2.prmtop -trj bcl2.mdcrd -sf 1 -ef 11 -inc 4 \
-tr 1-2611 -lr 2612-2654 -o /tmp/bcl2.out -sie
(c) 2007-2019 National Research Council Canada
sietraj 2.0 GLIBCXX_3.4.15
1 -5.706 10.312 -1288.445 -1293.298 -5.459 -29.364 -2.066 -446.976
5 -5.584 8.545 -1299.622 -1302.995 -5.172 -26.969 -1.747 -429.689
9 -5.924 9.447 -1291.449 -1295.816 -5.080 -30.162 -2.219 -467.440
The arguments given to sietraj are explained below. The order does not matter. The start and end atom numbers for the target and ligand can be found by generating a pdb file for the complex and looking for the TER records.
For carrying out an SIE calculation:
-pt prmtop_file (prmtop file to use)
-trj trajectory_file (trajectory coordinate file)
-sf start_frame (First snapshot to include)
-ef end_frame (Last snapshot ≤ end_frame)
-inc increment (Snapshot increment)
-tr target_range (Range of atom numbers comprising the target, separated
by a "-", e.g., 1-1947. If the numbering is not contiguous,
enter all ranges separated with a comma with no spaces,
e.g., 1-1947,2000-2034,2100-2234.)
-lr ligand_range (Range of atom numbers comprising the ligand.
Same format as target_range.)
-o output_file (Output file containing SIE energy terms.)
-sie (Carry out SIE calculation.)
For carrying out a virtual alanine mutation:
-mut2ala resnum (Mutate residue resnum to Ala on the fly.)
For calculating averages from the SIE output:
-ave filename (Process the output of the SIE calculation.
filename is the output file from the SIE calculation.)
The output file can be loaded into a spreadsheet for processing or you can use the sietraj command as shown below. Of course, in an actual application you would be averaging over more than just 3 snapshots. It is useful to plot a time series of the calculated snapshot ΔG values to make sure that you are averaging over a stable part of the trajectory, i.e., with little drift in the mean value over time.
prompt% $SIEHOME/bin/sietraj -ave /tmp/bcl2.out
Energies in kcal/mol:
Average StdErr Stdev
Inter vdW -28.83 0.78 1.36
Inter Coulomb -2.01 0.11 0.20
Reaction Field 9.43 0.42 0.72
Cavity -5.78 0.11 0.20
Constant -2.89 0.00 0.00
-----------------------
Delta G -5.74 0.08 0.14
Sample size= 3
Coefficients used:
alpha= 0.104758
gamma= 0.012894
const= -2.89
Delta_G = alpha * (vdw + Coul + RF + Cav) + constant
Cav = gamma * Delta_SA
The solvated interaction energy (SIE) is calculated for each snapshot via rigid infinite separation of the target and ligand. SIE is the sum of the intermolecular van der Waals and coulomb interactions plus the change in reaction field energy and nonpolar solvation energy. In this approximation, the internal energies of the target and ligand cancel out and do not contribute. The SIE is then scaled by an empirically determined factor, α, obtained by fitting to a training set of 99 protein-ligand complexes. The scaling can be considered a crude treatment of entropy-enthalpy compensation (Naïm et al., 2007; Chen et al., 2004).
Sample run 2
The second example is a protein-protein complex with 3844 atoms. Each protein has about 1900 atoms. 1VET.mdcrd is a trajectory file with 11 snapshots with all water molecules removed. 1VET.prmtop is the corresponding prmtop file for the complex.
prompt% cd $SIEHOME/examples
prompt% $SIEHOME/bin/sietraj -pt 1VET.prmtop -trj 1VET.mdcrd -sf 1 -ef 11 -inc 4 \
-tr 1-1947 -lr 1948-3844 -o /tmp/1VET.out -sie
(c) 2007-2019 National Research Council Canada
sietraj 2.0 GLIBCXX_3.4.15
1 -21.854 98.751 -1139.468 -611.547 -626.672 -151.209 -99.597 -2246.962
5 -22.592 105.858 -1151.303 -631.204 -625.957 -158.472 -107.433 -2173.438
9 -23.499 103.531 -1185.858 -645.690 -643.700 -161.664 -110.797 -2155.880
prompt% $SIEHOME/bin/sietraj -ave /tmp/1VET.out
Energies in kcal/mol:
Average StdErr Stdev
Inter vdW -157.12 2.53 4.37
Inter Coulomb -105.94 2.71 4.69
Reaction Field 102.71 1.71 2.96
Cavity -28.26 0.29 0.51
Constant -2.89 0.00 0.00
-----------------------
Delta G -22.65 0.39 0.67
Sample size= 3
Coefficients used:
alpha= 0.104758
gamma= 0.012894
const= -2.89
Delta_G = alpha * (vdw + Coul + RF + Cav) + constant
Cav = gamma * Delta_SA
Virtual alanine mutation
Virtual alanine mutations can be performed by adding "-mut2ala resnum" to the command line, where resnum is the residue number of the amino acid to be mutated. The mutation is carried out on the fly with partial charges and atom types assigned according to the FF99SB parameters. Note that disulfide bonded CYS cannot be mutated. Also, at the moment GLY cannot be mutated to ALA.
prompt% cd $SIEHOME/examples
prompt% $SIEHOME/bin/sietraj -pt bcl2.prmtop -trj bcl2.mdcrd -sf 1 -ef 11 -inc 4 \
-tr 1-2611 -lr 2612-2654 -o /tmp/bcl2.out -sie -mut2ala 94
(c) 2007-2019 National Research Council Canada
sietraj 2.0 GLIBCXX_3.4.15
Res # 94 (LEU) will be mutated to Ala
1 -5.592 10.330 -1289.129 -1294.000 -5.459 -28.192 -2.075 -454.455
5 -5.496 8.456 -1300.766 -1304.050 -5.172 -25.911 -1.777 -438.023
9 -5.806 9.470 -1292.949 -1297.338 -5.080 -29.161 -2.210 -460.367
Notes
- The prmtop file format is that of Amber 7 or later.
- If you wish to run the calculation on a single-snapshot restrt/inpcrd file, convert the file to a 10-column trajectory file using the Amber's ptraj utility. You can also use the crd2trj.sh and crd2trj.awk scripts found in the apps/sietraj subdirectory.
- The trajectory file must be generated with the nobox option in ptraj. Otherwise, the snapshot frames will not be read in properly.
How to cite
If you use sietraj in your work, please cite the following papers:
- Cui, Q.; Sulea, T.; Schrag, J. D.; Munger, C.; Hung,M.-N.; Naïm, M., Cygler, M.; Purisima, E. O. Molecular Dynamics and Solvated Interaction Energy Studies of Protein-Protein Interactions: the MP1-p14 Scaffolding Complex. J. Mol. Biol. 2008, 379, 787-802. [link]
- Naïm, M.; Bhat, S.; Rankin, K. N.; Dennis, S.; Chowdhury, S. F.; Siddiqi, I.; Drabik, P.; Sulea, T.; Bayly, C.; Jakalian, A.; Purisima, E. O. Solvated interaction energy (SIE) for scoring protein-ligand binding affinities. 1. Exploring the parameter space. J. Chem. Inf. Model. 2007, 47, 122-133. [link]
References
- Bhat, S.; Purisima, E. O. Molecular surface generation using a variable-radius solvent probe. Proteins: Structure, Function, and Bioinformatics 2006, 62, 244-261.
- Chan, S. L.; Purisima, E. O. Molecular Surface Generation Using Marching Tetrahedra. J. Comput. Chem. 1998, 19, 1268-1277.
- Chen, W.; Chang, C.E.; Gilson, M. K.. Calculation of Cyclodextrin Binding Affinities: Energy, Entropy, and Implications for Drug Design. Biophys.J. 2004, 87, 3035-3049.
- Naïm, M.; Bhat, S.; Rankin, K. N.; Dennis, S.; Chowdhury, S. F.; Siddiqi, I.; Drabik, P.; Sulea, T.; Bayly, C.; Jakalian, A.; Purisima, E. O. Solvated interaction energy (SIE) for scoring protein-ligand binding affinities. 1. Exploring the parameter space. J. Chem. Inf. Model. 2007, 47, 122-133.
- Purisima, E. O. Fast Summation Boundary Element Method for Calculating Solvation Free Energies of Macromolecules. J. Comput. Chem. 1998, 19, 1494-1504.
- Purisima, E. O.; Nilar, S. H. A simple yet accurate boundary element method for continuum dielectric calculations. J. Comput. Chem. 1995, 16, 681-689.
How to contact us
Send feedback, suggestions or queries to
mm-admin@nrc-cnrc.gc.ca.
