Protein Electrostatics: Rapid Multigrid-Based Newton Algorithm for Solution of the Full Nonlinear Poisson-Boltzmann Equation

Abstract

A new method for solving the full nonlinear Poisson-Boltzmann equation is outlined. This method is robust and efficient, and uses a combination of the multigrid and inexact Newton algorithms. The novelty of this approach lies in the appropriate combination of the two methods, neither of which by themselves are capable of solving the nonlinear problem accurately. Features of the Poisson-Boltzmann equation are fully exploited by each component of the hybrid algorithm to provide robustness and speed. The advantages inherent in this method increase with the size of the problem. The efficacy of the method is illustrated by calculations of the electrostatic potential around the enzyme Superoxide Dismutase. The CPU time required to solve the full nonlinear equation is less than half that needed for a conjugate gradient solution of the corresponding linearized Poisson-Boltzmann equation. The solutions reveal that the field around the active sites is significantly reduced as compared to that obtained by solving the corresponding linearized Poisson-Boltzmann equation. This new method for the nonlinear Poisson-Boltzmann equation will enable fast and accurate solutions of large protein electrostatics problems.     

Comparison of Poisson-Boltzmann and condensation model expression for the colligative properties of cylindrical polyions

Closed-form expression have been derived for the polyelectrolyte contribution to the colligative properties of solutions containing rodlike polyions in the presence of excess added salt. The derivations are based on: the conventional Poisson-Boltzmann equation for cylindrical symmetry; the thermodynamics of the cell model developed by Marcus [J. Chem. Phys. 23 , 1057–1068 (1955)]; and an equation derived from the cylindrical Poisson-Boltzmann cell model by Anderson and Record [Biophys. Chem. 11 , 353–360 (1980)]. Subject to the inherent limitations of the Poisson-Boltzmann approximation [Fixman (1979) J. Chem. Phys. 70 , 4995–5005], the resulting expressions are nevertheless applicable outside the “limit of infinite dilution.” They conform over a range of salt concentrations to the limiting laws deduced by Manning from the hypothesis of counterion condensation [J. Chem. Phys. 51 , 924–933 (1969)]. This hypothesis is found to be compatible with the Poisson-Boltzmann cell model but is not required in the derivation of the thermodynamic coefficients presented here. It is demonstrated that the magnitude of the polyion axial charge density plays a critical role in determining the low-salt limiting forms of the colligative properties obtained from the Poisson-Boltzmann equation, in close analogy with Manning's model.     

The δ-function counterion condensate around a line charge derived from the Fuoss,Katchalsky, and Lifson polyelectrolyte theory

It has recently been proven that the counterion condensate around an isolated line charge in an electrolyte, as characterized by nonlinear Poisson-Boltzmann theory, is an encapsulating δ-function. Here the identical result is shown to hold in the framework of the polyelectrolyte theory of Fuoss, Katchalsky, and Lifson. The proof fully exploits analytic solutions to the differential equation which are not available for the nonlinear, cylindrical Poisson-Boltzmann equation.     

The Poisson-Boltzmann equation for biomolecular electrostatics: a tool for structural biology

Electrostatics plays a fundamental role in virtually all processes involving biomolecules in solution. The Poisson-Boltzmann equation constitutes one of the most fundamental approaches to treat electrostatic effects in solution. The theoretical basis of the Poisson-Boltzmann equation is reviewed and a wide range of applications is presented, including the computation of the electrostatic potential at the solvent-accessible molecular surface, the computation of encounter rates between molecules in solution, the computation of the free energy of association and its salt dependence, the study of pKa shifts and the combination with classical molecular mechanics and dynamics. Theoretical results may be used for rationalizing or predicting experimental results, or for suggesting working hypotheses. An ever-increasing body of successful applications proves that the Poisson-Boltzmann equation is a useful tool for structural biology and complementary to other established experimental and theoretical methodologies.     

A fast and simple method to calculate protonation states in proteins.

A simple model for electrostatic interactions in proteins, based on a distance and position dependent screening of the electrostatic potential, is presented. It is applied in conjunction with a Monte Carlo algorithm to calculate pK(alpha) values of ionizable groups in proteins. The purpose is to furnish a simple, fast, and sufficiently accurate model to be incorporated into molecular dynamic simulations. This will allow for dynamic protonation calculations and for coupling between changes in structure and protonation state during the simulation. The best method of calculating protonation states available today is based on solving the linearized Poisson-Boltzmann equation on a finite difference grid. However, this model consumes far too much computer time to be a practical alternative. Tests are reported for fixed structures on bacteriorhodopsin, lysozyme, myoglobin, and calbindin. The studies include comparisons with Poisson-Boltzmann calculations with dielectric constants 4 and 20 inside the protein, a model with uniform dielectric constant 80 and distance-dependent dielectric models. The accuracy is comparable to that of Poisson-Boltzmann calculations with dielectric constant 20, and it is considerably better than that with epsilon = 4. The time to calculate the protonation at one pH value is at least 100 times less than that of a Poisson-Boltzmann calculation. Proteins 1999;36:474-483.     

Electrostatics of membrane adhesion.

We consider electrical double layer interaction under the conditions typically encountered during membrane fusion. Within the physiological concentration range of monovalent electrolytes the interaction is repulsive and the Poisson-Boltzmann calculation may be used to evaluate the force. When divalent counterions are added, strong ion-ion correlations make the Poisson-Boltzmann approximation inapplicable. We use the anisotropic hypernetted chain method to show that in the presence of small amounts of divalent counterions in adsorption equilibrium with the surfaces, the double layer interaction turns into attraction. This attractive electrostatic force may be the balancing contribution controlling membrane adhesion.     

Induced coalescence of cations through low-temperature Poisson-Boltzmann calculations

The computational determination of preferred binding regions of divalent counterions to nucleic acids is either inaccurate (standard Poisson-Boltzmann approaches) or extremely time-consuming (Monte Carlo or molecular dynamics simulations). A novel "selective low-temperature" Poisson-Boltzmann method is introduced that, although approximate in nature, qualitatively accounts for ion correlation and charge-transfer effects and allows for the rapid determination of such regions through an "induced coalescence" of divalent ions. The method is illustrated here for the binding of Mg(2+) to a double-helical sequence of B-form DNA (CGCGAATTCGCG) but the technique is readily applicable to locating divalent cations in other systems such as DNA-endonuclease complexes and ribozymes.     

Tests of continuum theories as models of ion channels. I. Poisson-Boltzmann theory versus Brownian dynamics

Continuum theories of electrolytes are widely used to describe physical processes in various biological systems. Although these are well-established theories in macroscopic situations, it is not clear from the outset that they should work in small systems whose dimensions are comparable to or smaller than the Debye length. Here, we test the validity of the mean-field approximation in Poisson-Boltzmann theory by comparing its predictions with those of Brownian dynamics simulations. For this purpose we use spherical and cylindrical boundaries and a catenary shape similar to that of the acetylcholine receptor channel. The interior region filled with electrolyte is assumed to have a high dielectric constant, and the exterior region representing protein a low one. Comparisons of the force on a test ion obtained with the two methods show that the shielding effect due to counterions is overestimated in Poisson-Boltzmann theory when the ion is within a Debye length of the boundary. As the ion gets closer to the boundary, the discrepancy in force grows rapidly. The implication for membrane channels, whose radii are typically smaller than the Debye length, is that Poisson-Boltzmann theory cannot be used to obtain reliable estimates of the electrostatic potential energy and force on an ion in the channel environment.     

The relationship between the poisson-boltzmann model and the condensation hypothesis: an analysis based on the low salt form of the Donnan coefficient

Two common models for the interaction of counterions with cylindrical polyions are considered in the context of the Donnan membrane equilibrium. General analytic expressions are obtained from the Poisson-Boltzmann equation for the Donnan coefficient in terms of the potential at the surface of the polyion or the local concentration of unbound ions at the surface. Analysis based on these expressions shows that if, and only if, the polyion charge density exceeds a certain critical value a large local concentration of ions will persist near the polyion surface at low ionic strengths. We therefore conclude that this principal hypothesis of the condensation model is consistent with the characteristics of the Poisson-Boltzmann potential at the surface of the polyion.     

Average electrostatic potential between the filaments in striated muscle and its relation to a simple Donnan potential.

An approximate analytical solution to the Poisson-Boltzmann equation for a cylindrical particle was used to calculate the relationship between the charge on the filaments and the average electrostatic potential. Both thick and then filaments were considered in the muscle lattice with a filament charge ratio of 4 to 1. Comparing this with a similar relationship obtained using simple Donnan theory showed a discrepancy at high charge where the Poisson-Boltzmann equation leads to saturation of the average potential. However, using two separate experiments from the literature, we have shown that at pH 7.0 muscle must not be close to saturation and thus is in a region of the curve where the two approaches agree.     

Energetics of charge-charge interactions in proteins

Electrostatic interactions between pairs of atoms in proteins are calculated with a model based on the linearized Poisson-Boltzmann equation. The equation is solved accurately by a method that takes into account the detailed shape of the protein. This paper presents applications to several systems. Experimental data for the interaction of ionized residues with an active site histidine in subtilisin BPN' allow the model to be tested, using various assumptions for the electrical properties of the protein and solvent. The electrostatic stabilization of the active site thiolate of rhodanese is analyzed, with attention to the influence of alpha-helices. Finally, relationships between electrostatic potential and charge-charge distance are reported for large and small globular proteins. The above results are compared with those of simpler electrostatic models, including Coulomb's law with both a distance-dependent dielectric constant (epsilon = R) and a fixed dielectric constant (epsilon = 2), and Tanford-Kirkwood theory. The primary conclusions are as follows: 1) The Poisson-Boltzmann model agrees with the subtilisin data over a range of ionic strengths; 2) two alpha-helices generate a large potential in the active site of rhodanese; 3) epsilon = R overestimates weak electrostatic interactions but yields relatively good results for strong ones; 4) Tanford-Kirkwood theory is a useful approximation to detailed solutions of the linearized Poisson-Boltzmann equation in globular proteins; and 5) the modified Tanford-Kirkwood theory over-screens the measured electrostatic interactions in subtilisin.     

Protocol for MM/PBSA molecular dynamics simulations of proteins

Continuum solvent models have been employed in past years for understanding processes such as protein folding or biomolecular association. In the last decade, several attempts have been made to merge atomic detail molecular dynamics simulations with solvent continuum models. Among continuum models, the Poisson-Boltzmann solvent accessible surface area model is one of the oldest and most fundamental. Notwithstanding its wide usage for simulation of biomolecular electrostatic potential, the Poisson-Boltzmann equation has been very seldom used to obtain solvation forces for molecular dynamics simulation. We propose here a fast and reliable methodology to implement continuum forces in standard molecular mechanics and dynamics algorithms. Results for a totally unrestrained 1 ns molecular dynamics simulation of a small protein are quantitatively similar to results obtained by explicit solvent molecular dynamics simulations.     

The theoretical distributions and diffusivities of small ions in chondroitin sulphate and hyaluronate

The electrostatic interactions between polyionic glycosaminoglycans and small mobile ions are investigated using the Poisson-Boltzmann equation and a rod-in-cell model of the polyelectrolyte. Calculations are made for the range of polyelectrolyte concentrations and buffer compositions for which measurements of ion distributions and diffusivities are reported in a companion paper (Maroudas et al., Biophys. Chem. 32 (1988) 257). We conclude that the distribution of mobile ions is largely determined by the 'far-field' potential and is adequately described by the Poisson-Boltzmann theory and also by more approximate theories such as ideal Donnan or 'condensation' theory. The measured variations in cation diffusivities, particularly the increase in diffusivity with increasing matrix concentration at low ionic strengths, are predicted qualitatively using an approximate diffusion theory together with the calculated potential fields. However, the same theory applied to anion diffusion gives qualitatively wrong results.     

Exact solution of the unidimensional Poisson-Boltzmann equation for a 1:2 (2:1) electrolyte.

The unidimensional Poisson-Boltzmann equation for a 1:2 (2:1) electrolyte has been solved analytically. The results have been compared with those obtained from the linearized equation. It is shown that in physiological conditions the difference may be greater than 10%. The value of the derivative of the potential in x=0, (dpsi/dx)x=0, has been used by many authors in the evaluation of the superficial charges of biological membranes. The value of (dpsi/dx)x-0 have also been compared with the ones derived from the linearized equation. The difference may be greater than 25%. Our results suggest that the linearization of the Poisson-Boltzmann equation for a 1:2(2:1) electrolyte may be greatly misleading.     

Electrostatic Potential of B-DNA: Effect of Interionic Correlations

Modified Poisson-Boltzmann (MPB) equations have been numerically solved to study ionic distributions and mean electrostatic potentials around a macromolecule of arbitrarily complex shape and charge distribution. Results for DNA are compared with those obtained by classical Poisson-Boltzmann (PB) calculations. The comparisons were made for 1:1 and 2:1 electrolytes at ionic strengths up to 1 M. It is found that ion-image charge interactions and interionic correlations, which are neglected by the PB equation, have relatively weak effects on the electrostatic potential at charged groups of the DNA. The PB equation predicts errors in the long-range electrostatic part of the free energy that are only ∼1.5 kJ/mol per nucleotide even in the case of an asymmetrical electrolyte. In contrast, the spatial correlations between ions drastically affect the electrostatic potential at significant separations from the macromolecule leading to a clearly predicted effect of charge overneutralization.     

Modeling salt-mediated electrostatics of macromolecules: the discrete surface charge optimization algorithm and its application to the nucleosome

Much progress has been achieved on quantitative assessment of electrostatic interactions on the all-atom level by molecular mechanics and dynamics, as well as on the macroscopic level by models of continuum solvation. Bridging of the two representations-an area of active research-is necessary for studying integrated functions of large systems of biological importance. Following perspectives of both discrete (N-body) interaction and continuum solvation, we present a new algorithm, DiSCO (Discrete Surface Charge Optimization), for economically describing the electrostatic field predicted by Poisson-Boltzmann theory using a discrete set of Debye-Hückel charges distributed on a virtual surface enclosing the macromolecule. The procedure in DiSCO relies on the linear behavior of the Poisson-Boltzmann equation in the far zone; thus contributions from a number of molecules may be superimposed, and the electrostatic potential, or equivalently the electrostatic field, may be quickly and efficiently approximated by the summation of contributions from the set of charges. The desired accuracy of this approximation is achieved by minimizing the difference between the Poisson-Boltzmann electrostatic field and that produced by the linearized Debye-Hückel approximation using our truncated Newton optimization package. DiSCO is applied here to describe the salt-dependent electrostatic environment of the nucleosome core particle in terms of several hundred surface charges. This representation forms the basis for modeling-by dynamic simulations (or Monte Carlo)-the folding of chromatin. DiSCO can be applied more generally to many macromolecular systems whose size and complexity warrant a model resolution between the all-atom and macroscopic levels.     

Composite cylinder models of DNA: application to the electrostatics of the B-Z transition.

We develop and test a Poisson-Boltzmann model of the electrostatics of the B-Z transition of DNA. Starting from the detailed geometries of the two forms, we compute at each radius the fractions of DNA matter, of volume forbidden (for nonpoint-like ions), and of volume accessible to the center of ions. These radial distributions are incorporated in a composite cylinder model; availability to ions (porosity) and the dielectric constant at each radial distance are then obtained. The phosphate charge is distributed with cylindrical symmetry on two layers at the appropriate radial distances. The porous sheath, between the axis and the charge distribution, provides much more room for ions in B-DNA than in Z-DNA. By using previously developed methods, the Poisson-Boltzmann problem of such cylinders is easily solved. The computational load is small, so that results can be obtained for a large set of salt concentrations and for a number of ionic radii. The variation of the electrostatic free energy difference with salt concentration compares favorably with the experimental value (it is half as large). There is also qualitative agreement with experiments on supercoiled DNA, including a maximum of the free energy difference at submolar salt concentrations. The results for this cylinder with porous sheath are in line with those of the earlier simple planar model and of a plain cylinder with sheath, which is also presented here. They are thus insensitive to details of the model. They support the proposition that the main electrostatic feature of the B-Z transition is the better immersion of the B-DNA phosphates into the solution. They also give confidence in the validity of the Poisson-Boltzmann approach, despite the large salt concentrations involved. Prior studies using an approach based on the potential of mean force are discussed.     

Electrostatic potentials in membrane systems

A membrane with an arbitrary distribution of fixed charges inside and on its surfaces is considered. A procedure for calculating the local electrostatic potential at an arbitrary point of the system is described and its validity discussed. This procedure is based on the linearization of the 3-dimensional Poisson-Boltzmann equation around an exact 1-dimensional solution.     

Calculation of the Electrophoretic Mobility of a Particle Bearing Bound Polyelectrolyte Using the Nonlinear Poisson-Boltzmann Equation

A numerical method for determining the electrophoretic mobility of a polyelectrolyte-coated particle is presented. The particle surface is modeled as having a permeable layer of polyelectrolyte molecules anchored to its surface. Fluid flow within the polyelectrolyte layer is subject to Stokes drag arising from the polyelectrolyte segments. The method allows arbitrary distribution of polymer segments and charge density normal to the surface to be used. The hydrodynamic plane of shear may also be varied. The potential profile is determined by a numerical solution to the nonlinearized Poisson-Boltzmann equation. The potential profile is then used in a numerical solution to the Navier-Stokes equation to give the required mobility. The use of the nonlinearized Poisson-Boltzmann equation extends the results to higher charge density/lower ionic strength conditions than previous treatments. The surface potentials and mobilities for three limiting charge distributions are compared for both the linear and nonlinear treatments to delimit the range of validity of the linear treatment. The utility of the numerical, nonlinear treatment is demonstrated by an improved fit to the electrophoretic mobility of human erythrocytes as a function of ionic strength in the range 10 to 150 mM.