Computation of the electrostatic interaction energy between a protein

A Modified Poisson−Boltzmann Model Including Charge Regulation for the ... Ananthakrishnan Sethuraman, Mina Han, Ravi S. Kane, and Georges Belfort...
0 downloads 0 Views 633KB Size
J . Phys. Chem. 1992, 96, 3130-3134

3130

of da/dE made little difference to wave shape (Figure 4C), but it did sigmfkantly affect curves for AE, >SO0 mV. Figure 5 shows experimental and simulation cyclic voltammograms for two scan rates on a low-defect surface, demonstrating a marked effect of nonconstant CY on the wave shape. The fits of theory and experiment shown in Figures 4A and 5 and in Table IV are quite good, supporting the conclusion that a potential-dependent a is important when AE, is large. Theory and experiment agree less well for intermediate AE, (Figure 4B), perhaps because the surface is quite defective Cf, 1-3%) and irregular. Although these observations do not prove that a potential-dependent a is the sole origin of the distortion from conventional shape, they do indicate that the observations are consistent with a nonconstant CY. Regardless of the details of the shape of the cyclic voltammograms on low-defect surfaces, there is no doubt that the observed kinetics are very slow for Fe(CN)63-/4-on the basal plane. We have reported similarly slow kinetics and basal plane for dopamine4 and solution-phase AQDSS3Sufficient data are not yet available to ascertain the generality of slow basal plane kinetics, and until such data are available, it would be risky to speculate an underlying phenomeon causing slow kinetics. The results on AQDS adsorption, capacitance, and Fe(CN)63-le kinetics reported here establish a means to reduce or eliminate contamination of basal plane behavior by edge plane defects. We are currently studying a wide range of redox systems on low-defect basal plane

surfaces, in order to determine the phenomena controlling the electrochemical behavior of the basal plane. The results will be reported separately. In summary, koobs,Coobs,and robs for AQDS are controlled largely by defects on basal plane H O E , and all three observables vary greatly for different surfaces. Although perfect basal plane has a Cobs below 1.0 rF/cm2, a koob for Fe(CN)63-/4-below lod cm/s, and negligible AQDS adsorption, observed values are usually higher due to adventitious or intentional defects. Since all three phenomena depend on defect density, they correlate with each other provided they are measured on the same plane surface. The voltammetry of Fe(CN)63-/4- on partly defective basal plane is consistent with a potential-dependent transfer coefficient, which is a major factor when PEPexceeds about 500 mV at 1 V/s. Finally, the reason for the very slow kinetics of Fe(CN)63-/4-on basal plane is not yet clear, but is currently being addressed by studying a variety of redox systems on low-defect basal plane surfaces.

-

Acknowledgment. This work was supported by the Air Force. Office of Scientific Research. We acknowledge Dr. Arthur Moore at Union Carbide for the generous gift of HOPG and Dennis Evans for the simulation software. We also thank Christie Allred for many useful discussions and software for electrochemical experiments.

Computation of the Electrostatic Interaction Energy between a Protein and a Charged Surface Byung Jun Yoont and Abraham M. Lenhoff* Department of Chemical Engineering, University of Delaware, Newark, Delaware 1971 6 (Received: November 14, 1991)

Protein-surface interactions play a central part in many situations both in vivo and ex vivo, but it is not yet possible to predict the extent of interaction quantitatively as a function of protein structure and surface and solution properties. A method for doing so is presented for systems dominated by electrostatic interactions, based on a continuum model of the electrostatic potential in and around a protein molecule near a charged surface in an electrolyte solution. The governing equations are solved by a boundary element approach which uses the method of images to treat the presence of an infinite planar charged surface exactly without any truncation. The potential distribution is used to obtain the protein-surface interaction energy; the boundary element approach allows this to be done rigorously, without approximating the self-energy contribution. The method is applied to determining the electrostatic interaction between a positively charged ribonuclease A molecule and a negatively charged surface at various orientations and separations. The strongest electrostatic attraction is found for ribonuclease. A with its active site facing the surface, while in certain orientations the interaction can be repulsive.

Introduction Interactions between protein molecules and solid-liquid interfaces are important in many applications,Is2 such as chromatographic separation of proteins, the design of biocompatible materials for medical uses, and food processing. In addition, proteinsurface interactions are important in vivo, e.g., in the approach of proteins to membrane surfaces. The proteinsurface interaction is extremely complex, involving the interplay of many physical and chemical factors. In principle, the evolution of configurational changes of the proteinsurface system can be examined theoretically using the equations of motion for each protein atom, and results based on molecular mechanics and dynamics computations are emergingn3 However, because of the structural complexities of the protein molecule, solution, and solid *To whom correspondence should be addressed. 'Present address: Department of Chemical Engineering, POSTECH, P.O. Box 125, Pohang 790-600, Korea.

0022-365419212096-3130$03.00/0

surface, analysis a t the atomic level remains very difficult. A less extreme degree of reductionism, applicable under conditions where the protein structure does not change appreciably due to interaction with the surface, is to treat the protein molecule as a rigid particle. Here information furnished in terms of the interaction potential energy between the molecule and the surface can be used to determine equilibrium states or to follow the dynamics by integration of the equation of motion for the molecule. The interactions are entirely electromagnetic in nature, but they are usually analyzed using continuum models4 in which the overall interaction is decomposed into parts, e.g., van der Waals inter(1) MacRitchie, F. Adu. Protein Chem. 1978, 32, 283-326. (2) Andrade, J. D.; Hlady, V. Adu. Pofym. Sci. 1986, 79, 1-63. (3) Lim, K.;Herron, J. N. In Biomedical Applications of Polyethylene

Glycol Chemistry; Harris, J . M., Ed.; Plenum, New York, in press. Lu, D. R.; Lee, S . J.; Park, K. J. Biomater. Sci. Polym. Ed. 1991, 3, 127-147. (4) Israelachvili, J. N. Intermolecular and Surface Forces; Academic

Press: New York, 1985.

0 1992 American Chemical Society

Electrostatic Protein-Surface Interaction Energy action, electrostatic double-layer interaction, hydrophobic interaction. These are usually treated as being independent and additive despite the fact that they may not be mutually exclusive. The most successful example of this approach is the DLVO (Derjaguin-Landau-Verwey-overbeek) the or^^,^ central to much of colloid science, in which the van der Waals interaction and electrostatic double-layer interaction are combined. Protein adsorption has been analyzed in terms of the DLVO theory,' but with the protein molecules treated as spheres with uniform properties. In this paper we present a more realistic approach which makes use of the full three-dimensional protein structure to account for the anisotropy of the shape and function of the protein molecule in determining its electrostatic double-layer interaction with a charged surface. Although the van der Waals interaction and solvation interactions exist, there are both mathematical and physical reasons to concentrate on the electrostatic double-layer interaction. First, the mathematical formulation of the problem is less ambiguous and more tractable, in terms of a continuum electrostatics description.*-I0 Second, electrostatic interactions are affected by such parameters as pH and ionic strength, which are easily manipulated experimentally. In addition, protein-surface interactions in many systems are dominated by the electrostatic double-layer interaction, e&, ion-exchange chromatography. In such systems van der Waals interactions are always present, but because they are relatively nonspecific, selective adsorption behavior must be controlled by other contributions such as the double-layer interaction. The nonspecificity of the van der Waals interaction is a consequence of its usual treatment based on the Lifshitz theory4 in which each phase is modeled as a structureless continuum characterized by only its dielectric permittivity. Since the dielectric permittivity of each phase is usually treated as a constant, computation of the van der Waals interaction for a given set of materials is essentially a geometrical problem, with the effccts of salt concentration and pH on the dielectric properties of the system neglected. Our work is based on an integral formulation of the governing equations for the electrostatic potential of the system, developed as an extension of our earlier approach" for the potential within and around a single molecule in an unbounded electrolyte solution. Numerical results are presented for ribonuclease A with and without the presence of a charged surface, and the potential field information is used to determine the electrostatic interaction energy. The outcome of the encounter between the negatively charged surface and ribonuclease A, positively charged as a whole, at various orientations and separations, is discussed in terms of the electrostatic interaction energy.

Theory Two different approaches are commonly used for the realistic treatment of electrostatics problems in biomolecular systems. In the microscopic model1*the analysis is performed at an atomic level of detail, whereas in the continuum model the biomolecule and the solvent are treated as two separate continuous media with different dielectric constants. Our method uses the continuum approach, which has been successfully employed for a range of biophysical applications.8-I0 The governing field equations in these models are generally solved numerically owing to the structural complexity of biomolecules, with the finite difference method most widely used. A boundary element method is used here, however, in view of several advantages that it offers, notably its more exact (5) Derjaguin, B. V.; Landau, D. L. Acta Physicochim. U.R.S.S. 1941,14, 633-662. (6) Verwey, E. J. W.; Overbeek, J . Th. G. Theory of The Stability of Lyophobic Colloids; Elsevier: Amsterdam, 1948. (7) Prieve, D. C.;Ruckenstein, E. In Colloid and Interface Science, Kerker, M., Ed.; Academic Press: New York, 1976, Vol. 4, pp 73-89. (8) Rogers, N . K. Prog. Biophys. Mol. Biol. 1986,48, 37-66. (9) Harvey, S.C.Proteins 1989,5, 78-92. (10) Sharp, K.A,; Honig, B. Annu. Reu. Biophys. Biophys. Chem. 1990, 19, 301-332. ( I 1) Yoon, B. J.; Lenhoff, A. M. J . Compur. Chem. 1990,II. 1080-1086. (12) Warshel, A.; Russell, S.T. Q.Reo. Biophys. 1984, 17, 283-422.

The Journal of Physical Chemistry, Vol. 96, No. 7, 1992 3131 geometrical treatment of interior charges and the molecular surface. A key problem in its implementation, viz., the incorporation of electrolyte effects, was recently remedied by utilizing the so-called direct boundary element method." Essentially the same formulation is used here, with a key modification included owing to the presence of a planar charged surface. The continuum formulation assumes the potential q5i in the interior B' of the biomolecule to be governed by the Poisson equation, with the dielectric constant taken to have the uniform value cis A set of n point charges q k are placed at positions xk inside the molecule. The linearized Poisson-Boltzmann equation is assumed to govern the potential dein the exterior region Bc (dielectric constant ce, Debye length l / ~ ) The . potential and the normal component of the electric displacement are continuous at the molecular surface as. These aspects of the formulation are identical with those presented previously,l' but in order to describe proteinsurface interactions in the present work, a planar surface aP is included here, bounding the electrolyte solution surrounding the molecule. The boundary condition at d P is either a constant surface potential or a constant charge density,

The integral representation of the Poisson equation is given by1IJ3J4

where G'(~x- x'I) = 1,

c={

1

4 ~ - 1x'I ~

w h e n x E E'

1/2, when x E aE

Here Gi is the Green function of the Poisson equation and n denotes the unit normal vector toward the solution phase Bc. The corresponding integral representation of the linearlized Poisson-Boltzmann equation is derived in the Appendix. It is given by

where Ge(lx - x'I) =

exp(-nlx - x'l) 4 r J x - x'I

c={

exp(-xlx - y'l) 4 ~ - y'I 1 ~

1, w h e n x E B"

In, w h e n x E &3

Here Ge is the Green function of the linearized Poisson-Boltzmann equation in the presence of the planar surface and y' is the mirror-image point of x'. The positive (negative) sign in the Green function corresponds to use of the constant charge (potential) condition on the surface. The symbol $ denotes the potential field in the absence of the protein, Le., the solution of the linearized Poisson-Boltzmann equation in a half-space with the boundary condition ( l a ) or (Ib). (1 3) Kellogg, 0.D. Foundations of Potential Theory; Dover: New York, 1953. (14) Jaswon, M. A.; Symm, G. T. Integral Equation Methods in Potential Theory and Elastostatics; Academic Press: New York, 1977.

3132 The Journal of Physical Chemistry, Vol. 96, No. 7, 1992

Yoon and Lenhoff

Y

\

io n -0.01

, /

X

\ \

.-,

‘. .0.05

0.01

Figure 1. Surface discretization of ribonuclease A, consisting of 2444 nodal points and 4884 triangles. The orientation of the molecule is represented by angles in a spherical coordinate system, with the origin of the molecule defined as the center of the nodal points, which is almost identical with the center of the atoms in the molecule. The view shown is along the axis 0 = 0, while 4 increases counterclockwise from the positive x-axis in the xy plane.

By using the boundary conditions on dB we obtain a set of two boundary integral equations

with the surface potential @ and its exterior normal derivative as unknowns. After solving (4) and (5) for them as discussed one can compute the potential a t any arbitrary points using (2) and (3). The boundary integral equations (4) and (5) were solved numerically after discretizing the molecular surface into a number of flat triangular elements and assuming a linear functional form of the unknowns on each element.” The discretization of the molecular surface was performed using a surface triangulation program obtained from Dr.Michael L. Conn01ly.l~ The output from this program was further processed to eliminate very slender triangles by removing short edges; the remaining edges are all at least about 1 A long. Figure 1 shows the resulting surface of bovine pancreatic ribonuclease A, for which electrostatic computations were performed in this work. The dielectric constants of the protein and electrolyte solution were chosen as 2 and 80, respectively, values typical of those used in molecular electrostatics computations.* All computations presented here were at pH 7 and a Debye length of 9.6 A, which corresponds to a 0.1 M 1:l electrolyte solution at 298 K. Assuming the intrinsic pK values for the ionizable groups of the amino acids, formal charges were assigned at the N-terminus and C-terminus, and at Lys, Arg, Asp, and Glu residues. All His residues were assumed to be uncharged. Under these conditions, which are generally consistent with the titration curve of ribonuclease A,I6 the molecule has 15 positive charges and 11 negative charges, and thus a net charge of +4. The dipole moment of the charge distribution is 461 D with its positive end pointing toward the region around the active site. (IS) Connolly, M . L . J . Appl. Crystallogr. 1985, 18, 499-505. (16) Tanford, C.; Hauenstein, J. D. J . Am. Chem. SOC.1956, 78, 5287-5291.

r

--

l / /

,

Figure 2. Equipotential contours inside and around ribonuclease A. Positive (broken line) and negative (solid line) potentials are shown in ~ 9.6 A, pH = 7, c’ = 2, ce units of k T / e . System parameters are: 1 / = = 80.

Quantitative information regarding the electrostatic interaction between a molecule and a surface is obtained from the electrostatic interaction energy U of the system, defined as the difference between the free energy F of the system at a given configuration (orientation and separation) and that a t an infinite separation. One may compute the difference in terms of either the Helmholtz free energy or the Gibbs free energy, but the two energies are same when the system is incompressible. When the Debye-Huckel approximation is valid the free energy is given by

The surface integral term has the form of the well-known expression for the free energy of a colloidal particle with surface charges? whereas the summation term has the form of the expression for the free energy of a protein molec~le.~’J*In order to determine F, the potentials on the surface and a t each protein charge must be computed; eqs 2 and 3 are used for this purpose. Since we are interested only in the difference in the values of F, the first term of the right-hand side of (2), which remains unchanged at different orientations and separations, is not included in the computation of +i. Exclusion of this term, which becomes infinite in the computation of the free energy, is an advantage of the boundary element method compred to the finite difference method, in which the contribution cannot be separated explicitly. ReSdtS We first compute the potential field for ribonuclease A in an unbounded solution. Figure 2 shows a cross-section of the equipotential contours in units of kT/e, where kT is the Boltzmann energy and e is the elementary charge. At 298 K 1 k T / e is approximately 25.7 mV. The cross section in Figure 2 corresponds to the plane z = 0, which passes within 1 A of the ionizable atoms of the Lys 41, His 12, and His 119 residues in the active site. The region of high positive potential inside the molecule corresponds to the location of the charge on Lys 41. The potential field around the molecule resembles that of an electric dipole with its positive end near the active site, but since the net charge of ribonuclease A is positive, as a whole the positive potential field is dominant. All computations including charged surfaces used the condition of constant charge density. Figure 3 shows the equipotential contours for ribonuclease A in two orientations in the presence of a charged planar surface. Since the Debye-Hijckel approximation used in deriving the linearized Poisson-Boltzmann equation implies small potentials (14el