Biophysical Journal 87:1544-1557 (2004)
© 2004 The Biophysical Society
Influence of the Solvent Structure on the Electrostatic Interactions in Proteins
Alexander Rubinstein and
Simon Sherman
Eppley Institute for Research in Cancer and Allied Diseases, University of Nebraska Medical Center, Nebraska Medical Center, Omaha, Nebraska 68198-6805
Correspondence: Address reprint requests to Simon Sherman, E-mail: ssherm{at}unmc.edu.
 |
ABSTRACT
|
|---|
The proper estimation of the influence of the many-body dynamic solvent microstructure on a pairwise electrostatic interaction (PEI) at the protein-solvent interface is very important for solving many biophysical problems. In this work, the PEI energy was calculated for a system that models the interface between a protein and an aqueous solvent. The concept of nonlocal electrostatics for interfacial electrochemical systems was used to evaluate the contribution of a solvent orientational polarization, correlated by the network of hydrogen bonds, into the PEI energy in proteins. The analytical expression for this energy was obtained in the form of Coulomb's law with an effective distance-dependent dielectric function. The asymptotic and numerical analysis carried out for this function revealed several features of dielectric heterogeneity at the protein-solvent interface. For charges located in close proximity to this interface, the values of the dielectric function for the short-distance electrostatic interactions were found to be remarkably smaller than those determined by the classical model, in which the solvent was considered as the uniform dielectric medium of high dielectric constant. Our results have shown that taking into consideration the dynamic solvent microstructure remarkably increases the value of the PEI energy at the protein-solvent interface.
|
INTRODUCTION
|
|---|
Electrostatic interactions play an important role in protein-protein complex formations, molecular recognitions, conformational adaptabilities, stabilities, and the function of proteins (Perutz, 1978
; Warshel, 1981; Warshel and Russell, 1984; Russell et al., 1987; Hendsch and Tidor, 1994; Honig and Nicholls, 1995; Sheinerman et al., 2000; Drozdov-Tikhomirov et al., 2001; Heifetz et al., 2002; Sinha and Smith-Gill, 2002; Kumar and Nussinov, 2002a). Despite a large number of studies focused on the evaluation of electrostatic interactions in proteins (Tanford and Kirkwood, 1957; Warshel and Levitt, 1976; Gilson et al., 1985; Rogers, 1986; Nakamura et al., 1988; Gilson and Honig, 1988; Schaefer and Froemmel, 1990; Rashin and Bukatin, 1994; Smith and Pettitt, 1994; Gilson, 1995; Luo et al., 1997; Papazyan and Warshel, 1997; Levy and Gallicchio, 1998; Sham et al., 1998; Schutz and Warshel, 2001; Simonson, 2001; Wisz and Hellinga, 2003), there are still many questions that have been poorly explored. Levy and Gallicchio formulated some of the questions, the answers of which are only beginning to emerge (Levy and Gallicchio, 1998). The major two questions are following: i), "how can the dielectric properties of boundary layer solvent that behaves differently from the bulk be most accurately captured in a continuum model?"; and ii), "what is the spatial variation of the dielectric response in a protein and is there a simple way to capture that variation within the framework of continuum models?" We believe that developing some treatment of the dielectric response of the protein-solvent interface, which would consider the solvent and protein as condensed polar media with inherent many-body dynamic properties, may help to answer the aforementioned questions and perhaps some other questions associated with the electrostatic estimations in biomolecules.
One of the methods that is most often used to calculate a pairwise electrostatic interaction (PEI) energy in protein (or other macromolecule) is the macroscopic continuum model (Gilson et al., 1985; Papazyan and Warshel, 1997). In the framework of this model, the protein is frequently considered as a uniform low-dielectric medium with a value of the dielectric constant between two and four. These values of the dielectric constant were estimated by experimental data on the dried protein films and powders, and were also predicted by a variety of theoretical estimations based on the theory of dielectrics (see Gilson and Honig, 1986; Nakamura et al., 1988; Simonson and Perahia, 1995; Simonson, 2001; and references cited there). The protein's low-dielectric response could in principle be larger, depending upon the extent to which the protein dipoles are free for reorientation (Gilson et al., 1985). One of the molecular mechanisms rationalizing conclusions on the high dielectric constant at the protein exterior was proposed in the framework of the Frölich-Kirkwood theory of dielectrics and takes into account the fluctuation of charged side chains of the residues concentrated at the protein-solvent interface (Simonson and Perahia, 1995; Simonson and Brooks, 1996). The local high dielectric constant at the protein interior was estimated by experimentally obtained pKa values of buried groups (see Dwyer et al., 2000; and references cited there). The microscopic and semimicroscopic models for calculating the effective dielectric constant in protein, which considers the reorientation effects of protein and solvent dipoles (as the Langevin dipoles) as well as effects of the solvent penetration into protein, was developed by Warshel and his associates (Warshel and Levitt, 1976; Warshel and Russell, 1984; Cutler et al., 1989; Papazyan and Warshel, 1997; Sham et al., 1998; Schutz and Warshel, 2001). The results obtained in these works suggest that proteins cannot be treated as a low-dielectric medium (Warshel et al., 1984). It was also suggested that the fluctuations and reorganization of the protein dipoles, but not the solvent effect, are the major factor that determines the dielectric response in proteins. This concept, however, was discussed critically and in great detail in several other works (see Gilson et al., 1985; Rogers, 1986; Nakamura et al., 1988; Levy and Gallicchio, 1998; Simonson, 2001; and references cited there). The strong argument for the low-dielectric protein model is that the protein dipoles, which are in a favorable but more or less rigid orientation, actually cannot be significantly reoriented in response to an externally applied electric field, and as a result, the medium cannot be viewed as having a high dielectric constant (Gilson et al., 1985). The same point of view that "the highly organized spatial structure of the proteins' polar groups results in the existence of a permanent intraprotein electric field and in proteins' weak dielectric response, i.e., its low dielectric constant," is one of the basic concepts of the insightful model for studying the kinetics of charge transfer in enzymatic reactions (Cannon and Benkovic, 1998; Krishtalik and Topolev, 2000; Mertz and Krishtalik, 2000). In another microscopic model, proposed by Simonson and his co-workers (Simonson et al., 1991), the local dielectric properties of a protein were investigated by the calculation of a generalized susceptibility, which determines the macroscopic susceptibility of continuum electrostatics. The obtained susceptibilities were consistent with values of the macroscopic dielectric constants of 
2–4. Thus, different theoretical and experimental approaches reveal the heterogeneous nature of the protein dielectric properties. The results obtained by these approaches are not simply reconciled with each other and some of them are inconsistent with the results of continuum electrostatic models for proteins.
The uniform low-dielectric model for the whole protein suggests that screening by the surrounding high dielectric solvent is negligible. This model is reasonably accurate for determining the electrostatic interactions between the charges on short distances (short-distance interactions) within the interior of proteins. The use of a low-dielectric constant to calculate the PEI energy between distantly remote charges (long-distance interactions), however, can result in an inaccurate estimation of these interactions within any region of a protein. This is especially obvious in the case of arbitrary distance for interactions within the exterior of a protein.
A geometrical description of the electrostatic field by the field lines gives a simple physical explanation for the inaccuracy of the PEI energy estimation, when the uniform dielectric model of a protein is used. In fact, when the force lines between charges of interest in a protein cross over the solvent, the corresponding electrostatic field undergoes an additional screening by high solvent permittivity. The more the force lines cross over the solvent, the bigger electrostatic screening takes place. It suggests that the effective protein permittivity should depend on the distance between the interacting charges of interest, their submergence in the protein interior, and orientation relative to the protein surface. It also suggests that the corresponding value of the effective dielectric constant for two interacting charges at the interface should depend on heterogeneity of the solvent dielectric properties that differ significantly on the protein surface and in the bulk phase (see Edsall and McKenzie, 1983 and references cited therein).
To take into account the solvent effects and evaluate the electrostatic interactions in biopolymers, the linear and nonlinear distance-dependent dielectric functions were used in the framework of the concept of the effective dielectric constant (Warshel et al., 1984; Hingerty et al., 1985; Rogers, 1986; Buravtsev et al., 1989; Smith and Pettitt, 1994; Mehler, 1996). In particular, the nonlinear sigmoidal functions that approximate the low-dielectric constant (4) characteristic for biopolymers at short distances and approach bulk dielectric constant of water (80) at large distances between charges were used (Hingerty et al., 1985; Ramstein and Lavery, 1988; Mehler, 1996; Hassan and Mehler, 2002; Wang et al., 2002). The other sigmoidal function approximates the value of the dielectric constant in proteins between 20 (at short distances) and 80 (at long distances) (Warshel et al., 1984; Schutz and Warshel, 2001). In principal, the aforementioned functions, however, can be incorrect because they ignore a dielectric boundary between a biopolymer and solvent.
The Tanford-Kirkwood theory (Tanford and Kirkwood, 1957) and its modifications (see, for instance, Gilson et al., 1985; Rogers, 1986; and references there) provide means to calculate the effect of the dielectric boundary on the PEI within proteins. It was shown that this effect is significant (Gilson et al., 1985). These approaches consider a protein medium as a spherical region of low-dielectric constant surrounded by a solvent of high-dielectric constant. Numerical methods were used to extend this approach for the native-like geometries of proteins and to solve the linearized Poisson-Boltzmann equation (Warwicker and Watson, 1982; Klapper et al., 1986; Sharp, 1991; Bharadwaj et al., 1995; Gilson, 1995; Baker and McCammon, 2003; Bordner and Huber, 2003). These calculations take into account the ionic strength effects of the solution, as well as the effects of the reaction field appearing due to the induced surface charges at the interface between two uniform dielectric media with high (80) and low (2–4) dielectric constants.
Dielectric models, utilized in the most aforementioned approaches, represent the solvent and solute as homogeneous dielectric continuums. In these models, in each point, r, of the space filled out by a medium, the electric induction vector D(r) and electric field vector E(r) are related by a simple linear equation:
 | (1) |
where 
(r) is a dielectric function of the medium. The potential of the electric field, 
(r), produced by the charge density of external field sources, 
(r), over the medium is determined by solving the first Maxwell equation:
 | (2) |
where 
is the divergence operator. After the substitution of Eq. 1 with E(r) = –(r), Eq. 2 has the form of the Poisson equation:
 | (3) |
or, in the case of the uniform dielectric, (r) = = const, transforms into the simplest one:
 | (4) |
Equations 3 and 4 can be modified into linear or nonlinear Poisson-Boltzmann equations (Robinson and Stokes, 1955; Bharadwaj et al., 1995; Baker and McCammon, 2003; Bordner and Huber, 2003) to model the ionic strength effects of the solution.
It should be noted, however, that Eq. 1 is only a particular case of the dielectric response for the homogeneous dielectric medium. In a more general case, the linear dielectric response for infinite homogeneous dielectric medium can be presented by an integral relationship between the electric induction D and the electric field E (Landau and Lifschitz, 1980):
 | (5) |
where the integration is taken over the volume V of the medium and the function 
ß(r – r') is the static permittivity tensor determined by the spatial correlation of the polarization of the medium in the space. Equation 5 allows one to introduce the influence of the microstructure of the dielectric medium, ß(r – r'), on the value of the electric induction D(r) (Dogonadze et al., 1977). It should be emphasized that the integral relationship, Eq. 5, is required to study the electric field, which is changing on the scale of the characteristic dynamic-structure distances of the system. For example, the Debye length in electrolyte solutions, as well as the characteristic length of several hydrogen bonds, at which the bound charge density fluctuations are correlated in solvent within the orientational Debye polarization mode, may determine the scale (Kornyshev, 1981).
According to Eq. 5, the values of the field E at point r and in some region around this point (in different points r') determine the induction D at point r. This is significantly different from Eq. 1, where the electric field E in point r determines the electric induction D at the same point. In this work, we are using a model based on Eq. 5. This model will be referred to as a nonlocal model for nonlocal dielectric to distinguish it from the models utilizing Eq. 1; models utilizing Eq. 1 will be referred to as local models for local dielectrics. Dependence of the dielectric function ß(r – r') upon r – r' is associated with an assumption on the continuous translational invariance of the infinite homogeneous dielectric medium, and is referred to as spatial dispersion of the dielectric function (Dogonadze et al., 1977).
In contrast to Eq. 1, which is determined by differential equations (Eqs. 3 and 4), Eq. 5 (after the substitution in Eq. 2) determines the integro-differential equation for potential of the electric field, (r) (Kornyshev et al., 1978; Kornyshev, 1981):
 | (6) |
For charge q located at point r0 in the dielectric medium
 | (7) |
where 
(r – r0) is the Dirac -function, the electrostatic potential (r) can be found by solving the integro-differential Eq. 6 using a Fourier transform. This solution can be presented by an analytical expression in the form of Coulomb's law with an effective distance-dependent dielectric function eff (R) (Kornyshev, 1985):
 | (9) |
 | (10) |
where (k) is the longitudinal wavenumber static dielectric function that is the Fourier component
and k = {Kx, Ky, Kz}.
In the case of the homogeneous isotropic infinite medium with uniform dielectric constant ,
 | (11) |
where ß is the Kronecker : ß = 1 if = ß, and ß = 0 if 
ß, Eq. 9 can be reduced to the classical coulombic form:
 | (12) |
According to this example, analytical properties of the dielectric function (k) determine the electrostatic screening in a condensed media (see Kornyshev and Vorotyntsev, 1980; Kornyshev, 1981; for more complicated dielectric functions). Such an approach demonstrates the so-called "dielectric formalism". This formalism, in particular Eqs. 9–10, was widely used in solid-state physics (Harrison, 1970), for the cases when the dielectric function (k) can be calculated on the basis of rather realistic micromodels. For more complicated systems, such as water, electrolyte solutions and other associated polar liquids, the corresponding calculations of (k) on the basis of microscopic ("first-principle") statistical-mechanic models are impeded. It is due to an absence of a self-consistence theory of polar liquids that would take into account not only electrostatic interdipole interactions but also quantum-mechanical short-range interactions, in particular, hydrogen bonding between water molecules (Dogonadze et al., 1977). That is why the microscopic models for polar liquids are rather simplified. Nevertheless, some common properties of spatial dispersion of the dielectric function (k) for homogeneous isotropic polar media can be revealed on the basis of a phenomenological theory of the polar solvent, which has been developed by R. R. Dogonadze and his associates (Dogonadze and Kornyshev, 1972, 1974; Dogonadze et al., 1977, 1985). In the framework of this theory, the relationship between (k) and the Fourier transform of the correlation function determined for the longitudinal polarization fluctuation in space (r) and time (t)—the dynamic structure factor S(
, k)—was stated by an effective Hamiltonian formalism and the fluctuation-dissipation theorem (see Zubarev, 1974; Kornyshev, 1981; and references cited therein):
 | (13) |
 | (14) |
where 
is the Planck constant; T is the absolute temperature in energetic units, 
...
means quantum-statistical average, and 
is the Heisenberg operator for the dynamic variable P(r, t): a position- and time-dependent medium polarization. In such an approach, the dielectric function (k) is a result of averaging of all possible (long- and short-range) interactions that take place in medium undergoing an external electric field, and that are very difficult to be fully taken into account in microscopic models. For an associated polar liquid, such as water, it is possible to select three major zones of the frequencies () of the electromagnetic absorption separated by transparency zones where the structure factor S (, k) 0. These absorption zones are associated with the specific polarization modes that have inherent distinguished frequencies of motions (Fröhlich, 1958; Dogonadze and Kornyshev, 1974; Kornyshev and Vorotyntsev, 1980), namely: 1), high-frequency electronic transitions changing the dipole moment; 2), fast oscillations due to infrared intramolecular vibration of dipoles; and 3), Debye fluctuations of the dipole orientations—hindered rotation of dipoles. Thus, total polarization of the medium is a result of superposition of these polarization modes. Each type of polarization is correlated in space with some effective correlation radius due to strong nonelectrostatic interactions. In the framework of this approach, the effect of the spatial dispersion was estimated, and the dielectric function, (k), was expressed through a set of the spectral functions of spatial correlation of polarization fluctuations determined for each of the corresponding three modes, f(kLi) (Dogonadze and Kornyshev, 1974; Kornyshev, 1981, 1985):
 | (15) |
where C1 = 1 – 1/1, C2 = 1/1 – 1/2, C3 = 1/2 – 1/3; i denotes the frequency-dependent dielectric constants in the different transparency regions, corresponding to electronic, infrared, and orientational (Debye) polarization modes, respectively; and Li denotes the typical correlation length of each mode. The expression for the f(kLi) can be approximated by some functions (Dogonadze and Kornyshev, 1972; Kornyshev, 1981; Kornyshev et al., 1982; Kornyshev and Sutmann, 1996). In the simplest case, this function is 
(i = 1,2,3). The analysis of the experimental data and theoretical consideration allowed for the estimation of the distinguished values Li (L1 
0.5 Å, L2 1 Å, and L3 3–5 Å) and the dielectric constants i (1 2, 2 5–6, and 3 80) for water solvent (see Dogonadze and Kornyshev, 1974; Kornyshev, 1981, 1985; Kornyshev and Volkov, 1984; Kornyshev and Ulstrup, 1986; and references cited therein). For the "pole" approximation, when the spatial dispersion in the high-frequency modes (electronic and infrared) may be neglected (L1 = L2 = 0), the dielectric function (k) for water has a simpler expression compared to Eq. 15, and substituted into Eq. 10, gives (see Kornyshev, 1985):
 | (16) |
The obtained expression for the effective dielectric function has a simple feature: the function eff(R) increases from the low value 2 (at small distances) to bulk water dielectric constant 3 (at large distances R >> L3). The physical origin of this effect was explained by a feature of the distribution of the induced bound charge in the effective spherical layer (thickness L3) around the charge of interest (Kornyshev, 1985; Kornyshev and Sutmann, 1996).
The nonlocal electrostatics in the bulk of a homogeneous infinite medium (Eqs. 5–10) was extended for interfacial electrochemical systems (Dogonadze et al., 1977; Kornyshev et al., 1977, 1978; Kornyshev and Vorotyntsev, 1980; Rubinshtein, 1981, 1983; Kornyshev, 1985). An interfacial system should be referred to as a heterogeneous one because the translational invariance of the bulk phase (if it exists) is broken down in the direction normal to the interface (Dogonadze et al., 1977; Kornyshev and Vorotyntsev, 1980). In other words, the medium structure on the interface is distorted significantly in comparison with the bulk phase. Information about this breach stretches from the interface to the depth of the medium (at least for the spatial dispersion scales) (Dogonadze et al., 1977). Consideration of this effect is very important in electrostatics, when the electric field is changing on the same scales. In this case, the kernel of the integro-differential Eq. 6, the permittivity tensor, becomes dependent upon r and r' separately (ß(r – r') becomes ß(r, r')). To solve Eq. 6, this tensor has to be determined at the interface. To do this, the tensor ß(r, r') can be expressed through the bulk dielectric function ß(r – r') of the given phase. This can be done in the framework of the so-called "sharp boundary model" that was previously used in the electrodynamics of semiinfinite plasma-like systems (Heinrichs, 1973; Platzman and Wolf, 1973). In the framework of this model, the so-called "specular-diffuse reflection approximation" was adapted to the dielectric medium (Dogonadze et al., 1977; Kornyshev et al., 1977, 1978; Kornyshev and Vorotyntsev, 1980). This approximation allows one to express the electrostatic potential through the bulk dielectric functions of the contacting media and, thereby, reveals the influence of the nonlocal bulk dielectric functions of the component on the distribution of the potential at the interface with different geometries (see Kornyshev and Vorotyntsev, 1980; and references cited therein).
The aforementioned nonlocal electrostatics developed for the bulk phase and a polar solvent interface led to the revision of some basic concepts of electrolyte solutions' properties and allowed one to explain several experimental data in the electrolyte theory (Dogonadze and Kornyshev, 1974; Kornyshev, 1981, 1985; Kornyshev and Sutmann, 1996) and interfacial electrochemistry (Kornyshev and Vorotyntsev, 1980; Rubinshtein, 1985, 1986; Kornyshev et al., 2002). The nonlocal electrostatic approach was effectively applied to solve several problems in the computational biophysics (Buravtsev et al., 1989; Leikin and Kornyshev, 1990; Kornyshev and Leikin, 2001; Rubinstein and Sherman, 2001).
In this work, a concept of nonlocal (NL) electrostatics for interfacial electrochemical systems was used to estimate the influence of the orientational Debye polarization of an aqueous solvent on the PEI energy in proteins. The effects of the electrolyte solution associated with Debye-Hückel screening by the mobile counterions were not considered in this work. The major goal of our work is to capture (in the continuum model) the dielectric properties of boundary layer solvent that behaves differently from the bulk phase and estimates the spatial variation of the dielectric response in a protein. The effects of protein structural reorganization were factored out of the continuum electrostatic model. Protein was considered as the uniform low-dielectric medium, in accordance with the classical theory of dielectrics. The general integral expression for the PEI energy at the interface between protein and solvent was obtained in the form of Coulomb's law with an effective distance-dependent dielectric function. It was demonstrated for short-distance electrostatic interactions between the charges, located in close proximity to the interface, that the values of the effective dielectric function were greater than the bulk dielectric constant in proteins, but were remarkably smaller than those determined by the classical model of the solvent (without consideration of the nonelectrostatic solvent dipole spatial correlations). This suggests that the spatial correlation effects of the solvent dipoles can make significant contributions into the PEI energy on the exterior of proteins. These contributions, however, have been underestimated in classical electrostatic approaches.
|
METHODS
|
|---|
In this work, the concept of NL electrostatics previously developed for arbitrary interfacial electrochemical system (Dogonadze et al., 1977; Kornyshev et al., 1977, 1978, 2002; Kornyshev and Vorotyntsev, 1980; Kornyshev, 1981; Rubinshtein, 1983, 1986) was used. The linear dielectric response for each of the media in contact was presented by the nonlocal relationship, Eq. 5, between the electric induction D and the electric field E:
 | (17) |
where m = 1, 2 is the number of the medium; the function m,ß(r, r') is the static permittivity tensor typical for each medium m.
To describe the interface between protein and water, the "sharp-boundary model" was utilized. In this model, the function m,ß(r, r') is nonzero only for points r and r' belonging to the same medium, m. Each medium was considered as a homogeneous isotropic one; the corresponding function of the medium was determined as:
 | (18) |
The simple case of the planar dielectric boundary that models the local regions of the interface between protein and solvent was considered, and a model system of two semi-infinite media was used (Fig. 1). The cylindrical coordinate system, (R, Z), was introduced for this boundary, where Z is the axis perpendicular to the plane passing through the boundary, and R is the two-dimensional radius vector in this plane. The semi-infinite region Z > 0 was assigned for the protein-like medium with an immersed charge Q1 located at point (0, Z0), and the region Z < 0 was assigned for the solvent.
View larger version (9K):
[in this window]
[in a new window]
|
FIGURE 1 Schematic presentation of the interface between solvent (Z < 0) and protein-like medium (Z > 0) in the cylindrical coordinate system (R, Z). Charge Q1 is located at the point (0, Z0), charge Q2-at the point (R, Z); r12 is the distance between charges Q1 and Q2.
|
|
To solve the electrostatic problem (Eq. 2 subject to the usual boundary conditions and Eqs. 17 and 18) and to find the electrostatic potential (R, Z, Z0) created by the charge Q1 in the point R, Z of the protein-like medium, the approach, based on the "specular reflection approximation" (Heinrichs, 1973; Kornyshev et al., 1977; Kornyshev and Vorotyntsev, 1980), was utilized. In this approximation, the desired potential can be presented in terms of the Fourier-transform 1(k) and 2(k) of the dielectric functions, 1(r – r') and 2(r – r'), that characterize the bulk properties of the media in contact (Kornyshev and Vorotyntsev, 1980; Rubinshtein, 1981). Thus, if the dielectric bulk functions of the protein-like medium and the solvent given as 1(r – r') and 2(r – r'), respectively, the general expression for the potential can be presented as:
 | (19) |
where:
 | (20) |
and J0 (KR) is the Bessel function of the first kind.
In this work, the orientational polarization, determined by the rotations of the dipoles in water, which, in turn, are hindered by the hydrogen-bonding chains, was considered by the phenomenological theory of the polar solvent (Dogonadze et al., 1977, 1985). In the framework of this theory the dielectric function for water, 2(k), was presented by the simple Lorentzian form utilized in several electrochemical and biophysical applications of the NL electrostatic approach (Kornyshev et al., 1978; Kornyshev, 1981; Rubinshtein, 1983, 1986; Buravtsev et al., 1989; Leikin and Kornyshev, 1990; Kornyshev and Sutmann, 1996):
 | (21) |
where * = 6 and s = 78.3 are short- and long-wavelength dielectric constants of the solvent at room temperature; L is the correlation length of the water dipoles, which is proportionate to the characteristic length (3–5 Å) of the hydrogen-bonding network of water molecules (Kornyshev, 1985; Kornyshev and Ulstrup, 1986). This dielectric function, which is used for modeling the aqueous solvent as a system of rigid and strongly correlated dipoles (Kornyshev, 1981), gives the macroscopic (bulk) dielectric constant, s, for the small k-limit (distances much larger than L) and the short-wavelength dielectric constant, *, for the large k-limit (distances much smaller than L). In fact, at large wave vectors, k, the spatial dispersion of the dielectric function 2 (k) can be neglected: 2 (k 
) = 1. In this case, the macroscopic theory is actually inapplicable. Therefore, the utilized model of the dielectric response for water assumed that the function 2(k) in Eq. 21 changes from long wavelengths, s, to some value (2 () = * > 1), corresponding to distances smaller than the parameter L, but still exceeding interatomic distances (Rubinshtein, 1983; Kornyshev and Sutmann, 1996). Such peculiarities of spatial dispersion of the dielectric function 2 (k) with inherent key parameters L and *, Eq. 21, resemble the effective distance-dependent dielectric function that can be estimated from experimental data interpreted as the free energy of interaction between protonated amino groups in dibasic amines (Kornyshev and Ulstrup, 1986). The protein-like medium was presented as a uniform dielectric with a low-dielectric constant: 1(k) = 1 = 4.
In the protein-like medium the PEI energy between the point partial charges Q1 = 
1e and Q2 = 2e (1, 2—parts of the electron charge e) was calculated as the product of the electrostatic potential, Eq. 19, created by the charge Q1 at point (R, Z) and the charge Q2, which is located at this point (Fig. 1). All integrals in Eq. 20 were calculated as in the work (Rubinshtein, 1983), where the method of complex contour integration in the complex plane KZ was utilized. Omitting cumbersome computations, the PEI energy can be written in the form of:
 | (22) |
where: r12 = [R2 + (Z – Z0)2]1/2 is the distance between charges Q1 and Q2 in the cylindrical coordinate system,
 | (23) |
|
RESULTS AND DISCUSSION
|
|---|
As one can see from Eq. 22, the PEI energy, U12(R, Z, Z0), looks similar to the expression for Coulomb interactions with the NL effective (distance-dependent) dielectric function, eff-NL(R, Z, Z0). The second and third terms in Eq. 23 mean that, in addition to the direct Coulomb interaction between charges, there is a contribution associated not only with the distance-dependent density distribution of the surface-induced bound charges (
b(R)), but also with the space-induced bound charges of the solvent (b(R, Z)) in close vicinity to the interface (Rubinshtein, 1983, 1985). The contribution to the direct coulombic interactions in proteins associated with the reorganization effects of the surrounding environment (due to reorientation of the protein and the solvent dipoles) was also considered in the framework of the microdielectric model (Warshel et al., 1984; Sham et al., 1998). When the correlations between the water dipoles are not considered (L = 0 or * = s = 78.3), the NL effective dielectric function, eff-NL(R, Z, Z0), is transforming into the limit expression for the classical dielectric function, eff-CL(R, Z, Z0):
 | (24) |
This function can be easily obtained from the corresponding expression for the PEI energy in the case of the classical model of two uniform dielectrics (2(k) = s and 1(k) = 1) with a flat interface (Jackson, 1975). The effective dielectric functions, given by Eqs. 23 and 24, approach the value 1 at limiting large distances Z, Z0 + and finite R: eff-CL (R, Z, Z0) = eff-NL(R, Z, Z0) = 1. According to Eq. 23, eff-NL(R, Z, Z0) depends upon the distance, r12 = [R2 + (Z – Z0)2]1/2, between the interacted charges as well as upon their orientation relative to the interface and remoteness, Z0, from the interface. The eff-NL(R, Z, Z0) function also depends upon the correlation length of the water dipoles, L, which is the spatial parameter for R and Z.
To simplify the analysis of the eff-NL(R, Z, Z0) function, two extreme orientations of the charged couple were considered: i), parallel, (||), to the interface, (R = r12, Z = Z0), and ii), perpendicular, (
), to the interface, (R = 0, Z = Z0 + r12). An asymptotic expansion of the integral in Eq. 23 was performed for limiting large and small values of Z0/L, Z0/r12, and r12/L. This asymptotic analysis suggests that the values r12, Z0, and L are expended (in some cases outside of the physical scales) to select small parameters for corresponding expansions and to reveal the major features of the eff-NL(R, Z, Z0) as a formal function of two variations, r12 and Z0, and one parameter of the system, L. Simple estimations of the behavior of eff-NL(R, Z, Z0) were carried out considering only the first terms of the expansions. As is demonstrated below by numerical analysis, the behavior of the NL dielectric function in the field of physically applicable distances (r12 and Z0), comparable with the correlation length L (5 Å), is consistent with the behavior of this function revealed by the asymptotic analysis.
In the (||) case, the limiting (asymptotic) values for 
were obtained for six spatial regions (A, B, C, D, E, F), each of which is characterized by a unique relationship between distances Z0, r12, and the correlation length of the water dipoles, L. Parametrical characteristics of these regions and obtained estimations of the 
values in these regions are shown in Table 1 and Fig. 2. It is necessary to note that similar asymptotic solutions in A and B regions (as well as in C, D, and E regions) are slightly different in the higher order terms of the expansions.
View this table:
[in this window]
[in a new window]
|
TABLE 1 Asymptotic values of the NL dielectric function,  in the (||) case for six spatial regions of the variables r12 and Z0
|
|
Asymptotic expressions for the classical dielectric function, eff-CL (R, Z, Z0), described by Eq. 24, were also obtained. Table 2 presents the limiting values for this function obtained for the (||) case, 
when only the first terms of the expansions were taken into consideration.
View this table:
[in this window]
[in a new window]
|
TABLE 2 Asymptotic values of the classical dielectric function, in the (||) case for two spatial regions of the variables r12 and Z0
|
|
In the () case, the asymptotic expressions for 
and 
were also obtained. Restricting the use of the first terms of expansions, these expressions were found to be similar to the expressions presented in Tables 1 and 2, respectively.
As can be seen from Table 2, there is a spatial scale, r12 Z0, at which the estimated value for the classical dielectric function is significantly changing. This scale divides the space of the parameters in two parts, r12 >> Z0 and r12 << Z0, where the corresponding asymptotic solutions are significantly different from each other: eff-CL (r12 >> Z0, Z0) / eff-CL (r12 << Z0, Z0) >> 1. It is necessary to note that the corresponding values (1 + s)/2 for r12 >> Z0 and 1 for r12 << Z0 are the classical limiting values obtained in the continuum electrostatics for the effective dielectric constant at the interface between two uniform dielectric media with dielectric constants 1 and s (Gilson et al., 1985; Rogers, 1986).
In the case of the NL electrostatic approach, an additional scale, the dimensional parameter L, arises (see Table 1 and Fig. 2). This scale determines the behavior of the function eff-NL(r12, Z0) in protein-like medium with changing Z0 and r12 for both extreme orientations of the interacted charges. More specifically, the correlation length L of the solvent determines an additional spatial scale of dielectric heterogeneity at the interface of the protein-like medium. As a result, the specific spatial region F is arising (Table 1; Fig. 2), where the limiting value of eff-NL(r12, Z0) is much smaller compared to the value (1 + s)/2 obtained in the classical (local) case (Table 2, spatial region A) and the values predicted by the Tanford-Kirkwood theory for surface groups (Tanford and Kirkwood, 1957; and relevant calculations carried out by Gilson et al., 1985).
Overall, the asymptotic solutions shown in Table 1 and Fig. 2 allow for the formulation of the following heterogeneity features of the NL dielectric function:
- For the region Z0 << L, an increase of the distance r12 leads to an increase of the function eff-NL (r12, Z0) starting from the values 1 (for r12 << Z0) passing through the values (1 + *)/2 in the intermediate region Z0 << r12 << L and reaching the values (1 + s)/2 (for r12 >> L). Because the value (1 + *)/2 is slightly different from the value 1, it suggests that the NL dielectric function, Eq. 23, in close proximity to the interface has low values of 1 in the more extended region when r12 << L, in contrast to the classical case, when the corresponding dielectric function, Eq. 24, is restricted by the values 1 in the region where r12 << Z0 (Table 2, region B). The "classical" behavior of the NL dielectric function takes place only for Z0 >> L (regions B, C, and D in Table 1; Fig. 2).
- For the region r12 << L, a decrease of the Z0 values leads to an insignificant increase of the function eff-NL (r12, Z0) from the values 1 (for Z0 >> L) to values (1 + *)/2 (for r12 << Z0 << L). It is different from the classical case, when the function eff-CL (r12, Z0) is approximately equal to (1 + s)/2 for r12 >> Z0 (Table 2). The eff-NL (r12, Z0) has "classical" behavior only for r12 >> L (regions C, B, and A in Table 1; Fig. 2).
In the asymptotic analysis, the values r12, Z0, and L were expanded to reveal peculiarities in the behavior of the NL dielectric function described by Eq. 23. Because the correlation length L in the considered aqueous solvent model, Eq. 21, is 3–5 Å, it is obvious that for distances r12, much smaller than L, the function eff-NL (r12 0, Z0) does not have physical meaning and is restricted by some values on the distances that exceed interatomic distances. Therefore, the eff-NL(R, Z, Z0) function was analyzed for values r12 that are 
1.5 Å. For a model of the planar interface between solvent and protein-like media the applicability of Eqs. 22 and 23 for the accurate estimation of the electrostatic interactions in proteins is restricted by a range of distances r12, which are smaller than a radius of curvature for the local region of the interface. Such a restriction can be satisfied in globular proteins and membranes. For the intercharge distances, when this restriction is not implemented, the corresponding estimation is less accurate because it is impossible to neglect the curvature of the interface. The lowest limit for Z0 was considered as the van der Waals radius of the aromatic hydrogen in proteins (1 Å).
Numerical calculations of the eff-NL(r12, Z0) and eff-CL(r12, Z0) functions were performed using Eqs. 23 and 24, respectively. The calculations were made for the charged couple having two extreme orientations relative to the interface, (|| and ). Dependence of the dielectric function upon the distance Z0 was analyzed for several intercharge distances, r12, taken in the interval of 3–6 Å (Figs. 3 and 4). This interval is typical for distances between partial charges within amino acid residues as well as for protein salt bridges (Kumar and Nussinov, 2002a). For charges located in close proximity to the interface (Z0 within range 5 Å), dependence upon the intercharge distance r12 was also analyzed (Figs. 5 and 6). Figs. 3–6 present behavior of the effective dielectric function for the classical local model (* = s = 78.3), eff–CL(r12, Z0), and for the NL model (* = 6.0, s = 78.3, L = 5.0 Å), eff–NL(r12, Z0), of the aqueous solvent for both (|| and ) extreme orientations.
View larger version (15K):
[in this window]
[in a new window]
|
FIGURE 4 Behavior of the effective dielectric function depending on the distances Z0 from the interface in the case of parallel (||) orientation of two charges of interest with fixed intercharge distances r12: 2.5, 3.5, 4.5, 5.5, and 6.5 Å—the curves 1, 2, 3, 4, and 5, respectively.
|
|
View larger version (15K):
[in this window]
[in a new window]
|
FIGURE 6 Behavior of the effective dielectric function depending on distance r12 in the case of parallel (||) orientation for two charges of interest with five fixed values Z0: 1.5, 2.5, 3.5, 4.5, and 5.5 Å—the curves 1, 2, 3, 4, and 5, respectively.
|
|
As can be seen from Figs. 3 and 4, the value of the eff–NL(r12, Z0) function for both extreme orientations is changing from the values 4–5 (for large Z0 > r12) and reaching the values 6–9 (for small Z0 < r12). For both (|| and ) orientations, such an increase in the value of eff–NL(r12, Z0) for the small Z0 is much smaller than the increase suggested by the classical dielectric function (see Fig. 3, A and B). The same tendency was observed for all other intercharge distances (data are not shown). Fig. 4 shows the greater the distance between the charges, the greater the value of the effective dielectric constant at the interface, and the greater the distances Z0 when dielectric function approaches the bulk value 4.
Our analysis has shown that for the given parameter r12, in the cases when the charged couples are located in close proximity to the interface (i.e., Z0 is small) and when the distances between the charges are short, the NL effective dielectric function determined by Eq. 23 varies within a narrow interval of values. This interval is limited by the values within the area between the and 
curves. For example, when r12 is equal to 3.5 Å, for any orientation of two charges, the eff–NL(r12, Z0) function for the small Z0 adopts the values within the interval 5–7 (Fig. 3). Analogously, when r12 is equal to 5.5 or 6.5 Å, the eff–NL(r12, Z0) function for the same Z0 adopts the values within the interval 6–9 (Fig. 4, data for the are not shown).
Figs. 5 and 6 demonstrate the strong dependence of the eff–NL(r12, Z0) function upon the distances r12 between two interacting charges. According to Fig. 5, there is a significant difference between the values of the classical (local), eff–CL(r12, Z0), and the NL, eff–NL(r12, Z0), dielectric functions. This is especially clear in the case of the parallel orientation to the interface (Fig. 5 A), when the value of parameter r12 is increasing. Fig. 6 shows that the greater the remoteness of the charged couple from the interface, Z0, the smaller the value of the effective dielectric constant, and the greater the distance r12 when eff–NL(r12, Z0) approaches rather large values that are comparable with the value (1 + s)/2 >> 1. This limiting value obtained by the asymptotic analysis (Fig. 2) is the maximum of the eff–NL(r12, Z0) function at the interface, and it is smaller than the value of the long-wavelength dielectric constant (s) of the solvent. Such behavior of the eff–NL(r12, Z0) function differs essentially from the behavior of the sigmoidal functions previously used (Warshel, et al., 1984; Ramstein and Lavery, 1988; Buravtsev et al., 1989; Hassan and Mehler, 2002; Wang et al., 2002) that ignore the presence of a dielectric interface between the protein and solvent.
As can be seen in Table 1 and Figs. 2–6
, the data of the numerical analysis are in agreement with the solutions of the asymptotic analysis. For example, in Fig. 4, at the small distances Z0 (Z0 < L) and r12 (r12 < L, curves 1 and 2), the values of the NL dielectric function are 5–6 that are comparable with the values (1 + *)/2 estimated by the asymptotic analysis in the region F (Fig. 2; Table 1). In Fig. 6 the similar behavior of the NL dielectric function for the same Z0 and r12 is exhib
TOP
ABSTRACT
INTRODUCTION
