Calculating Effective Charges

The electrostatic forces between our proteins are found by calculating the electric field from the electrostatic potential, and the F = q*E. To speed up our calculations, we will use what are called effective charges. These are a set of charges (relatively few in number) that reproduce the electrostatic potential (for more details see: R.R. Gabdoulline and R.C. Wade. Effective charges for macromolecules in solvent J.Phys.Chem. (1996) 100, 3868-3878). Since the number of effective charges is about two orders of magnitude less than the number of atoms is the protein, the force calculations run much faster.

The first step in determining effective charges is to choose the sites where the charges will be. This is completely arbitrary, but since the acidic and basic groups on the protein are concentrated sites of charge, we choose these are our locations. Specifically, we choose the terminal atoms on each of the side chains:

These sites are automatically set up using the program ecm_mksites. This program takes your protein structure from stdin, and directs the output to stdout, i.e.

ecm_mksites < barnase.pdb > barnase.sites

In the above line, I'm doing this in my ECM directory, but change as needed. Look at the resulting file. I think you should get 28 sites for barnase and 26 for barstar. Now you need to set up an input file for the effective charge calculations. This file looks like:

------------------ pdb file name
barnase.pdb
------------------ file with test charges for a molecule
barnase.sites
------------------ grid file name and its form (0-binary, 1 -asci) UHBD format
barnase.50mM.grd
0
------------------ probe, skin: expansion happens in [probe; probe+skin] interval
4.0, 3.0
------------------ ionic strength, solvent dielectric
50., 78.4
------------------ file to write effective charges
barnase.echa
------------------ nothing else

The input parameters are documented, but the most important factors you may need to change are the structure file, the sites file, the potential file (make sure this is you binary grid created using asc2bin), the ionic strength (must be the same as you apbs calculation), and the effective charge file name. Using this file, run the command

ecm_expand150 < barnase.echa.in

This first step finds the Coulomb grids and calculates the charge matrix. Next, you need to regularize the matrix using

ecm_regularize

Finally, you calculate a set of effective charges using

ecm_mkecharges

This program will prompt you for a level of regularization. If you choose 1, you will do the best job of fitting the potential, but the individual charge values may be quite large. If you choose the maximum level (28 in the case of barnase), you will simply recover the original set of charges produced by ecm_mksites. The best option is to choose a middle value, and I typically choose the midpoint (if there are 28 sites, I choose 14). This does a reasonable job of fitting the potential, but still gives reasonable charge values.

Once you have your effective charges for both barnase and barstar, move on to determining reaction sites.