Dielectric saturation of water in a membrane protein channel

Abstract

Water molecules in confined geometries like nanopores and biological ion channels exhibit structural and dynamical properties very different from those found in free solution. Protein channels that open aqueous pores through biological membranes provide a complex spatial and electrostatic environment that decreases the translational and rotational mobility of water molecules, thus altering the effective dielectric constant of the pore water. By using the Booth equation, we study the effect of the large electric field created by ionizable residues of an hour-glass shaped channel, the bacterial porin OmpF, on the pore water dielectric constant, ε_w. We find a space-dependent significant reduction (down to 20) of ε_w that may explain some ad hoc assumptions about the dielectric constant of the protein and the water pore made to reconcile model calculations with measurements of permeation properties and pK_a’s of protein residues. The electric potential calculations based on the OmpF protein atomic structure and the Booth field-dependent dielectric constant show that protein dielectric constants ca. 10 yield good agreement with molecular dynamics simulations as well as permeation experiments.

---

Introduction

The structural and dynamical properties of water are essential to account for many observed features in a wide variety of systems in the micrometer and nanometer range. Typical examples involve the physicochemical descriptions of thin films of aqueous electrolytes, polymeric ion exchange membranes, synthetic nanopores or even solid state materials like fractions of rocks and sandstones.^1–5 Particularly important is the role of water in biological membranes and macromolecules. As an essential part of many metabolic processes, water is vital for photosynthesis and respiration and also central to acid–base neutrality and enzyme function. Acting as a solvent, water has many distinct properties that are critical for the proliferation of life that set it apart from other substances.⁶ In particular, we are interested here in the exchange of charged and neutral solutes across the cell membrane envelope. This physiological function is mostly accomplished by highly conductive ion channels regulating the influx of nutrients and the extrusion of waste products.⁷ In the present study we consider OmpF, a general diffusion bacterial porin, that forms large channels in the outer membrane of E. coli.^8–10 Each monomer forms an hourglass-shaped aqueous pore with diameter in the range 1–4 nm.¹¹ Reconstituted in planar lipid bilayers OmpF homotrimeric channels allow multi-ionic transport and exhibit moderate cationic selectivity at neutral pH in solutions of monovalent salts.^12–14 This selectivity has been reported to be highly sensitive to the charge state of the ionizable residues of the channel,¹⁵ particularly of those lying at the channel constriction.¹⁶

The properties of water in an intricate geometry like a protein channel could differ dramatically from its well-known bulk conditions. In particular, simulation studies of several membrane channels have suggested that pore water under the influence of a large number of ionizable residues and permeating ions exhibits decreased translational and rotational mobility.¹⁷ This means that in high fields the water dipoles become oriented along the electric field direction, in order to minimize the dipole-field interaction energy. The decreased polarizability gives rise to a reduction of the dielectric constant known as dielectric saturation of water. Water is a strongly hydrogen-bonded liquid and its behaviour in confined geometries of nanometer dimensions differs from bulk water not only due to saturation effects but also because the effective cluster that responds to the local field is smaller in size that in bulk water. The reduced cluster size might result in a lowering of the relative permittivity of the same order of magnitude as that due to dielectric saturation effects.^18,19 The goal of the present study is then to show how the dielectric properties of water play a significant role in the electrostatic interactions determining the function of membrane protein channels. In a large number of studies aimed to the evaluation of electrostatic energies in proteins the discussion has focused on the protein dielectric environment,^20–22 while little attention has been paid to the properties of the surrounding water.²³ Thus, ionization of acidic and basic residues of proteins (quantified by their corresponding pK_a values), reduction or oxidation processes in redox centers in proteins and binding of charged ligands have been chosen as relevant benchmarks to discriminate between consistent and inconsistent models on the basis of a single parameter, the protein dielectric constant ε_p.²² This explains the increasing number of papers aimed to understand how dielectric properties of water regulate the electrostatic effects playing a major role in protein systems. Previous studies devoted to ion channel function addressed this problem by assigning to the effective dielectric constant of water in the channel a single value lower than the well known ε_w = 78 bulk value.^17,24 Cheng and co-workers performed finite difference Poisson–Boltzmann calculations showing that dielectric saturation affects the solvation energy of permeating ions in a narrow model channel system (a bundle of four α-helical segments).²⁵ The pore water dielectric constant was assumed to be down to 5 as suggested in previous studies on gramicidin-like channels because of their similar radii.²⁴ For porin-like channels, a water dielectric constant of 24 has been proposed.²⁶ More recent Brownian dynamics simulations, taking advantage of the MthK channel atomic structure, used values of ε_w ranging from 20 to 80 to show that the decrease of dielectric constant results in an increase of repulsive image forces that deepen the channel energy well.²⁷ Note that the use of a single, uniform value for the pore water dielectric constant means that the spatial dependence of ε_w is neglected. If this is analyzed in terms of Poisson’s equation, one can anticipate that the electric field will be underestimated since an additional source of field is lost.⁵ We suggest here a slightly different approach: To use an electric-field-dependent dielectric constant for the water following Booth original treatment.²⁸ We solve the three-dimensional Poisson–Boltzmann equation with a spatially dependent dielectric constant following Booth equation. The calculations of the electric field and the field-dependent water dielectric constant are iterated until the desired accuracy is obtained in the self-consistent solution. By comparison with previous studies done in OmpF channel^12,29–31 we discuss the validity of the approach presented here and the physicochemical meaning of the dielectric constants involved. The coexistence of discrepant values of the dielectric constants and electric potential wells trying to explain the same facts can be understood by taking into account that different methodologies and different simplifying assumptions are made.

Formulation of the problem

Extensive literature is available on the dielectric constant of polar liquids under high electric fields, in cases where polarization is not linear with the field and eventually saturates.^1,32,33 The nonlinear character of polarization is reflected in the decrease of dielectric constant with the field. Experiments in water confirmed that the decrease in the dielectric constant is proportional to the square of the electric field.³⁴ However, because of the technical difficulties of those experiments, only few results have been reported.³⁵ This explains why most of the work done in this subject is either theoretical or, more recently, based on computer simulations.³³ The model employed here was developed by Booth for pure water²⁸ and can be considered as extension of the seminal work of Onsager and Kirkwood. The water dielectric constant ε_w is assumed to change with the electric field as:where ε₀ is the vacuum electric permittivity, n is the optical refractive index, N₀ is the number of molecules per unit volume, μ is the water molecule dipole moment, k is the Boltzmann constant, T is the absolute temperature and L is the Langevin function. Eqn (1) does not assume the dielectric constant of water to be in the linear response regime: the non-linearity of the polarization vector is caused here by the functional dependence of the water dielectric constant on the electric field.³⁶ The use of the Booth theory for the dielectric saturation of water in an ion channel might be questionable from several points of view. One serious objection has to do with the complicated behaviour of water under high electric fields. Computer simulations studying polarization density of water have reported slight deviations from Langevin equation at fields around 5 V nm⁻¹.³⁷ The characteristic flexibility of the water molecule, not considered in the Langevin equation, could change the interaction energy by increasing/reducing the O–H distance and/or the H–O–H angle. More striking effects appear in fields around 10 V nm⁻¹ where water undergoes restructuring phase transitions into a proton-ordered cubic ice-like structure.³⁸ However, for electric fields up to 1–5 V nm⁻¹ the Booth equation is in excellent agreement with MD simulations carried out for water enclosed between charged surfaces.³³

Aside from the commented disparity between bulk and pore water due to the confined geometry and the electrostatic interaction with protein residues, another interesting criticism pays attention to the nature of water in aqueous solutions, stressing that Booth’s original model was developed for pure water and not for an electrolyte solution.³⁹ Keeping in mind that electrophysiological experiments with ion channels are typically performed at decimolar salt concentrations, the average number of ions inside the OmpF channel studied here is so small (two–three ions)²⁹ that we can overlook the above criticism. Only in this sense, such a diluted electrolyte should not differ too much from pure water. In fact, the evaluation of bulk dielectric constant of aqueous electrolyte solutions by means of Booth equation together with the Debye–Hückel theory has early successful antecedents.⁴⁰

Finally, note that Booth’s original treatment was developed for a totally homogeneous electric field, that is not the case of OmpF porin. In regions with a strong spatial dependence of the direction of the electric field it would not make sense to talk about a preferential orientation for the water molecules. However, more elaborated theoretical approaches considering inhomogeneous electric fields have given support to the local use of Booth equation in regions when a moderate spatial dependence in the direction of the electric field is found,^2,41,42 as it is the case of the present study. This point will be addressed later in the paper when the electric field in the OmpF porin is discussed.

Solution procedure

The OmpF channel coordinates were obtained from the Protein Data Bank (PDB code 2OMF) and charges and radii were assigned using the charmm27 force field. The electrostatic potential created by a macromolecule on its neighbourhood can be obtained from the solution to the Poisson equation:where V is the electric potential, ε is the dielectric constant of the medium and ρ is the volume charge density. For macromolecules in free solution, the source term (ρ) includes fixed charges anchored to the macromolecule, but also free or “mobile” charges present in the ionic solution. These mobile charges rearrange in the presence of potential or electric field. This dependence is usually represented using the Boltzmann equation. If we combine Poisson with the Boltzmann distribution for mobile (solute) ions, we obtain the Poisson–Boltzmann (PB) equation:where c is the solute concentration (we assume a 1 : 1 electrolyte) and e is the elementary charge. Very often this equation is further simplified by taking the first term in the series expansion of the right hand side, leading to the linear Poisson–Boltzmann (LPB) equation. This approximation is usually acceptable for small potentials (or alternatively small fixed charge density onto the protein) or high concentrations. Thus, eqn (3) is commonly referred to as nonlinear Poisson–Boltzmann (NPB) equation. Several computational packages have been designed to solve numerically the PB equation in or around macromolecules. However, some of them are only efficient for highly restricted dielectric map distributions. For example, UHBD⁴³ or APBS⁴⁴ can only handle two regions of uniform ε values, which are usually the protein interior and the solution region (generally 4 and 80 are used, respectively, although other values can be chosen⁴⁵), with a smooth transition between them. These programs use a finite difference (UHBD, APBS, Delphi) or finite element (APBS) method to get the potential distribution and/or the electrostatic energies in a predefined grid. These continuum approaches have been widely used in the literature and usually provide a fast and accurate method to determine electrostatic energies.⁴⁶

The dielectric constant that appears in the left-hand side of eqn (3) is assumed to change with the local electric field as predicted by Booth, i.e.eqn (1). In order to solve the mathematical problem posed by eqns (1) and (3) we have extended to our purpose the original code of the APBS program. This specifically tailored PB algorithm is able to deal not only with arbitrary space distributions of the dielectric constant (hereafter denoted as dielectric maps), but also with dielectric maps that depend on the local electric field, which will increase the nonlinear character of the PB equation. We have selected APBS because of its open source characteristics and stability to deal with very nonuniform dielectric maps (in our tests, UHBD did not converge even for small departures from uniform dielectric constants). The APBS python wrapper is also more convenient due to the high level tools available in python to deal with grid maps (basically arrays) resulting from computations. We proceed as follows:

1. First we read the molecule atomic coordinates into the APBS code and compute the potential map using a first standard keyword input file, i.e., solving the PB equation (eqn (3)) in a two-region domain. The output is the solute access region grid and the dielectric constant and potential maps of the molecule.

2. A program is developed to read the above maps and to recalculate the dielectric constant map corresponding to the solute region using the Booth equation (eqn (1)). The code compares the current dielectric constant map with the previous one. Iterations continue until the following convergence criterion is met: the change in the dielectric constant between two consecutive steps at any grid point is less than 0.01.

3. The second input file points to the new dielectric constant map generated in step 2 (also the solute access region grid, as required by the APBS code). The electric potential is then recalculated.

Steps 2 and 3 are repeated until required convergence is achieved. Boundary conditions for the APBS solver were defined by choosing the option “Zero”, which means that the potential values at the boundaries of the computation box are set to zero. The size of the box was assessed, so that an increase in its magnitude had no effect on the resulting potential. Temperature was set to 298 K. All computations assumed protonation states corresponding to pH 7 and standard pK_a taken from the literature for the ionizable residues. We found this procedure fast and robust: less than 15 iterations lead to the desired accuracy, resulting in the 3D dielectric and potential maps. In all computations the three monomers that assemble together to form the channel were taken into account. The numerical computation was done using a 129 × 129 × 97 grid with a 0.1 nm cell size. The program used to create the dielectric grid (step 2) also represented the membrane that surrounds the channel by a low dielectric region (ε_m = 2), non accessible to ions, located between z = 1.6 nm and z = 5.2 nm. z denotes the axial coordinate along the OmpF channel with origin in the intracellular side. All computations were performed in a Dell Precision 690 workstation (8 CPU and 16 GB RAM). Another code was developed to average, when needed, the electric potential and the water dielectric constant over each channel cross-section. The averaging procedure used slices on the xy plane, spaced 0.1 nm in the z direction.

Results

Dielectric saturation and Booth equation

The dielectric constant of water as a function of the field calculated by means of eqn (1) is shown in Fig. 1. The optical refractive index of water is n = 1.33 and N₀ is calculated from the molar concentration of pure water (55.55 × 10³ mol m⁻³). The standard value of the water molecule dipole moment μ = 1.84 Debye⁷ is slightly changed to μ = 2 Debye (1 Debye = 3.33 × 10⁻³⁰ C m) in order to obtain the correct zero field result of ε_w = 78. (Note that this adjustment was also done by Booth in the original development of eqn (1)).²⁸ The dielectric constant drops to half its conventional value about 0.5 V nm⁻¹ and to one-fourth around 1 V nm⁻¹.

Fig. 1 Dependence of water dielectric constant with the strength of electric field according to the Booth equation.

For electric fields higher than 1 V nm⁻¹ the dielectric constant of water reaches almost a plateau value. This means that the orientational polarization is nearly saturated, and the only degree of freedom in which the electric field may act is the internal geometry of the molecule, changing O–H distances and H–O–H angles and potentially leading to a dissociation of the molecule.³² MD simulations of water performed in a variety of charged systems justify the validity of the Booth equations for fields up to 3–5 V nm⁻¹,^33,37 which, as we will show later in the paper, are considerably higher than the typical electric fields found for the OmpF porin in the present study.

Fig. 2 shows a 2D map of the dielectric constant in a longitudinal cross-section of the OmpF channel. ε_w values are calculated by means of the NPB equation and the Booth equation self-consistently over the whole domain protein + aqueous pore. Boundary effects due to the protein surface lower the dielectric constant of water, which reaches its minimum value around the center of the pore. This could be expected from both geometric and electrostatic reasons. On the one hand, the channel displays a narrow central constriction. Molecular dynamics simulations performed in model hourglass-shaped cavities suggested that self-diffusion coefficients and rotational reorientation rates of water molecules were consistently lower for the central zone of the pore than for the channel mouths.⁴⁷ Moreover, a recently reported first passage time approach analyzing water diffusion through the OmpF channel also shows that a small fraction of the water molecules appear to be trapped by the channel walls for considerable lengths of time especially in the central part of the channel.⁴⁸ On the other hand, a considerable number of ionizable residues are located around the central constriction originating a strong electric field that exerts a noticeable effect on the channel electrostatic properties.²⁹

Fig. 2 Longitudinal cross-section of the OmpF channel with a contour plot of the water dielectric constant calculated by means of the NPB equation and the Booth equation self-consistently over the whole domain protein + aqueous pore (just one of the three monomers is shown although in all computations the three monomers that assemble together to form the channel were taken into account). Electrolyte concentration was set to 150 mM to mimic physiological conditions. The dielectric constant in the protein region was ε_p = 2.

Electric field-dependent vs. uniform water dielectric constant

Interestingly, when the ε_w values predicted by the Booth equation (displayed in Fig. 2) are averaged over the whole OmpF pore, one gets a mean dielectric constant of ε_w≈ 60. This value is in good agreement with recent Brownian dynamics studies performed in the MthK channel where the effective pore water dielectric constant is found employing a learning-based stochastic optimization algorithm.²⁷ Then, the question that arises naturally is whether the dielectric constant of water calculated locally by means of the Booth equation can be replaced by its average along the channel axis or not. To address this point we compare the electric potential profile across the pore, calculated by using uniform values of ε_w and field-dependent values of ε_w. Fig. 3 displays the electric potential averaged over the pore cross-section for different axial positions in the OmpF channel. The plot is restricted to the central constriction of the channel, the region where ions face the greater energy well or barrier and also where the differences in the averaged potential between the two approaches are noticeable. In accordance with ref. 27, the depth of the energy well decreases when ε_w = 78 is changed to an effective value of ε_w = 60,²⁷ although in the present case the quantitative effect is insignificant (less that 1 kT/e).

Fig. 3 Electric potential averaged over the pore cross-section along the OmpF channel central constriction. The 3-D map of the electric potential over the protein and the aqueous pore was calculated by using the NPB equation and three options for the pore water dielectric constant: a uniform value of ε_w = 78 (dashed line), a uniform value of ε_w = 60 (dot line), and a value obtained by solving self-consistently the Booth equation (solid line). The bulk ion concentration was set to 150 mM to mimic physiological salt concentration. The dielectric constant in the protein region was ε_p = 2.

Note that this effect is small when compared to the outcome of the Booth equation, which yields a potential well exceeding by almost 50% the one calculated with ε_w = 78. This can be an indication that the use of average values for the whole span of the channel could be misleading in cases where the dielectric saturation effect is concentrated in a relatively small part of the channel, as is the case of the channel constriction of the OmpF porin. Interestingly, the electric potential well provided by the Booth equation could be attained as well using a single value of ε_w≈ 20, similar to that proposed by Jordan et al.²⁴ or Chen et al.²⁵ In view of such result which would turn the whole pore water in an almost insulating medium, one could wonder then if there is a real physical feature behind this value or it is a consequence of pushing the model beyond its limits.

Linearity vs. nonlinearity of Poisson–Boltzmann equation

Several studies on the effect of water dielectric constant on channel ion permeation, electrostatic solvation energy and protein pK_a determination are based on calculations that use the LPB equation.^24,25 Having in mind that extremely strong electric fields are needed for getting noticeable saturation effects on the dielectric constant according to the Booth equation, one may ask whether the LPB equation is accurate enough or the use of the NPB equation becomes mandatory. Fig. 4A shows the water dielectric constant profile across the OmpF channel when using the linear and nonlinear form of the PB equation. Apparently both forms behave similarly and only small differences are found around the central constriction. Indeed, the 2D map shown previously in Fig. 2 using the NPB equation does not display essential differences when it is recalculated via the LPB equation. From previous studies⁴⁹ it is known that the LPB equation overestimates the electric field, which means a lower dielectric constant in terms of eqn (1). Looking back to Fig. 1 we can see that in the region of low values of ε_w the dependence on the electric field is weak, and even in the case when the LPB equation overestimated noticeably the electric field, the resultant ε_w would not differ significantly from the outcome of the NPB equation calculation, as is the case in Fig. 4A. However the situation changes dramatically if we consider the electric potential profile across the OmpF channel using again both the linear and nonlinear forms of the Poisson–Boltzmann equation. Fig. 4B shows that the depth of the potential well obtained by using the LPB equation is almost twice that predicted by using the NPB equation. A close inspection of the values in Fig. 4B shows that the departure from linearity is evident even outside the central constriction. This could be expected since LPB formalism assumes that the electric potential (in kT/e units) is much smaller than 1, which is not the case here.

Fig. 4 Effect of using the LPB or the NPB equation on the cross-section averaged water dielectric constant (A) and on the cross-section averaged electric potential along the OmpF channel (B). The Booth equation was used in both computations. The other parameters used in the computations are the same as in Fig. 3.

Comparison of dielectric constants obtained using different methodologies

Once we have discussed the use of the PB equation, we address now the issue of suitable values that should be used for both the dielectric constants of water and protein media. The experimental measurement of the dielectric constant faces several technical difficulties that come from its dependence on the shape and size of the sample, temperature, water content and frequency and intensity of the incident electromagnetic field.³⁵ The lack of experimental data may be the reason why in the field of protein electrostatics the dielectric constant has been regarded as an adjustable parameter more than a universal characteristics of each material.²² In this sense we can categorize two main ways of tackling the problem, one focusing on ε_w, while ε_p remains unchanged and the other way round. Thus, on the one hand, we can find approaches that consider the protein as a low-dielectric medium (ε_p = 2–4) following the theory of dielectrics and experimental data on dried protein films and powders²³ whereas dielectric constants for pore water are assigned values in the range 20–78 due to dielectric saturation effects.^27,50,51 This has been the case of the present study up to the present point.

On the other hand, in continuum electrostatics calculations of the protonation state of titratable surface groups of proteins, pore water is regarded as bulk water whereas the dielectric constant of the protein is assumed to be in the range between 2 and 40.⁵²Water penetration into the protein, coupled to structural rearrangements is invoked to account for such an increase in ε_p.⁵³ In those pK_a calculations, the contribution of each charge to the total electric potential is computed by means of the superposition principle or, in other words, the linearity of the electric potential with charge. In such scenario the use of LPB equation is mandatory for the sake of consistency.⁴⁹ (Note that an alternative to the superposition principle invoking the use of NPB is non-viable here since it involves 2^N calculations, N denoting the number of titratable residues, and in the case of OmpF N = 306). Within the framework of pK_a calculations via the LPB equation, the water dielectric constant remains unchanged and the values of ε_p must be adjusted in order to obtain a satisfactory comparison with the experiments.^12,29Fig. 5 shows the electric potential across OmpF using the standard LPB equation for different values of ε_p but the water dielectric constant as in the bulk (ε_w = 78). In the case of OmpF, a value of ε_p ~ 20 is needed to obtain reasonable pK_a values.¹² This combination of parameters provides an electrostatic potential well around 6 kT/e, which allows to reproduce accurately current–voltage curves, channel conductance and selectivity pH dependence, and even the slight asymmetry found in the reversal potential when the direction of the salt gradient is inverted.²⁹

Fig. 5 Effect on the cross-section averaged electric potential across OmpF of the protein dielectric constant, using the standard LPB equation. The other parameters used in the computations are the same as in Fig. 3.

One can wonder now whether it is possible to reconcile such different views: the approach presented in Fig. 2 and 4 focused on the dielectric saturation of water and the procedure derived from pK_a calculations which concentrates in the dielectric properties of the protein, as shown in Fig. 5. Having in mind that different methodologies are used (NPB in Fig. 2 and 4 and LPB in Fig. 5) the values of the dielectric constants used in both approaches are not directly comparable. However, we can rationalize the situation saying that up to some extent both views suffer from the same problem: using the widely accepted values for the dielectric constant in a biological membrane or a protein and in bulk water (ε_p = 2 and ε_w = 78, respectively) the available models overestimate the electrostatic effects involved and fail to be consistent with experiments. Then, the depth of calculated energy wells can be reduced a posteriori in two ways: by reducing ε_w (see Fig. 3) or by rising ε_p (see Fig. 5), yielding values that in some cases have no clear physical meaning. The next goal of the present study is to show that both views can be combined in a more complete and robust approach.

Fig. 6 shows that the potential well shown in Fig. 5 for ε_w = 78, ε_p = 20 and the LPB formalism (which provided a satisfactory explanation of a series of experiments)¹² can also be reproduced with slightly different parameter values that we consider more coherent. First of all, ε_w = 78 is replaced by a field-dependent dielectric constant calculated via the Booth equation within an NPB formalism. Secondly, we recognize that the interior dielectric response of proteins is far away from the tabulated ε_p = 2–4 typical of hydrocarbon phase. Taking into account that values of ε_p≈ 20 are considered too high by many authors²² we propose here an intermediate value of ε_p = 10 that does not aim to be conclusive, but a compromise between the extreme values found in literature that have unclear physical meaning. This would reflect, in a more moderate way, the effects previously reported in globular proteins like certain water penetration and, consequently, structural reorganization. The fact that NPB gives similar results to the LPB formalism could seem counterintuitive in the light of Fig. 4. Note, however, that in Fig. 6 both formalisms are not directly comparable since they make different assumptions about water and protein dielectric properties and use different values for both ε_w and ε_p.

Fig. 6 Comparison between the average electric potential obtained using the LPB equation and a high dielectric constant for the protein region (dot line) and that obtained using the NPB equation with the Booth equation and a lower protein dielectric constant (solid line). The other parameters used in computations are the same as in Fig. 3.

Despite the coincidence shown in Fig. 6 between the electric potential well obtained from two different approaches, the fine details of the electric potential map at the channel constriction, calculated in each case, may reveal useful information on discriminating between them.

Fig. 7A and B show contour maps of the electric potential in the cross-sectional plane (at z = 3.42 nm) of the OmpF channel central constriction computed using the LPB equation, ε_w = 78 and ε_p = 20 (Fig. 7A) and the NPB equation with the Booth equation and ε_p = 10 (Fig. 7B). In both panels the dot line separates the aqueous pore (enclosed region) and the surrounding protein domain, non accessible to the solvent. A visual comparison states that Fig. 7A shows a coarser electrical map of the channel constriction, although one of the essentials physical features of the pore eyelet, the strong negative electric potential due to the channel acid groups in the left hand side part of the figure is captured. However, such picture would make the central channel region almost inaccessible to the permeating negative ions. The more detailed picture given by Fig. 7B shows clearly a region of positive potential located near the so-called cluster of arginines. This is in good agreement with MD simulations reporting that cations and anions follow well-separated permeation pathways along the OmpF porin.³⁰

Fig. 7 Contour plot of the electric potential in a cross-sectional plane (z = 3.42 nm) of the OmpF channel central constriction. The numbers denote the maximum value of the potential on each region in kT/e units. The dot line depicts the limit between the aqueous pore (inner region) and the surrounding protein domain, non accessible to solvent. Top panel, A: electric potential obtained using the LPB equation, ε_w = 78 and ε_p = 20; Bottom panel, B: electric potential obtained using the NPB equation with the Booth equation and ε_p = 10. The other parameters used in the computations are the same as in Fig. 3.

The electric field and the validity of Booth equation

As commented earlier in the paper, one key point determining the coherence of the approach presented here is the validity of the Booth equation. An analysis of the models discussed in Fig. 6 and 7 in terms of the associated electric fields can be useful in this sense.

Fig. 8 displays the cross-section averaged axial and transversal electric field obtained using the LPB equation and a high dielectric constant for the protein region (dot line) and the NPB equation with the Booth equation and a lower protein dielectric constant (solid line). Both approaches agree in showing a huge transverse electric field that would explain the screw-like trajectories of the ions found in MD simulations³⁰ as well as the high residence time of large polar molecules at the channel constriction.⁵⁴ As could be anticipated from the wide equipotential zones in Fig. 7A, the LPB equation provides a lower electric field than the non-linear form of the same equation which, according to Fig. 7B, shows a steeper gradient of the electric potential. Interestingly, the values of the electric field provided by the NPB equation using Booth are in good accordance with a MD simulation of the OmpF channel performed by Tielemann and co-workers.³¹ The water dipole moment computed at the central constriction in such simulations yields a transversal electric field around 1 V nm⁻¹ and an axial component of 0.3 V nm⁻¹. As can be seen in Fig. 8, these values are in excellent agreement with the results obtained using the Booth equation around the central constriction (z∼ 3 nm).

Fig. 8 Cross-section averaged axial (A) and transversal (B) electric field obtained using the LPB equation and a high dielectric constant for the protein region (dot line) and that obtained using the NPB equation with the Booth equation and a lower protein dielectric constant (solid line). The other parameters used in the computations are the same as in Fig. 3.

Fig. 8 also shows that both axial and transversal components of the electric field exhibit a similar spatial dependence, and their maximum values are attained around the same location. This means that the dielectric saturation effects are concentrated around the central narrow constriction, where a strong transverse electric field of well-defined orientation is found. The moderate spatial dependence of the direction of the electric field gives credit to the approach presented here, originally developed for homogeneous electric fields.^41,42 Note finally that in all cases, the maximum value of the electric field obtained via the Booth equation is far away from the problematic values of about 5 V nm⁻¹ (where water could undergo restructuring phase transitions), thus assuring the physical soundness of this approach.

Conclusions

We have shown that continuum electrostatics calculations of the electric field in the bacterial porin OmpF yield values that are high enough to orient water molecule dipoles and consequently decrease the water dielectric constant. By using the Booth equation to account for the dependence of ε_w on the electric field calculated iteratively by means of the NPB equation we have obtained a reduction of ε_w of more than 50% at the channel constriction. We have evidenced that the use of average values of ε_w for the whole span of the channel could be inadequate in the cases where the dielectric saturation effect is concentrated in a relatively small part of the channel, as is the case of the channel constriction of the OmpF porin. The analysis of the electric potential along the pore axis sheds some light on the physical meaning of ε_p values used in previous pK_a calculations of the OmpF channel residues performed using the LPB equation. Finally, we compare two different assumptions about the protein and water pore dielectric properties which yield roughly the same potential well at the channel constriction. It is seen that the electric map (electric potential and field) following from a field-dependent ε_w is more consistent with reported MD simulations and channel permeation properties.

Acknowledgements

The authors acknowledge support from the Spanish Ministry of Education (project FIS2007-60205).

References

1. S. H. Chen and M. C. Bellissent-Funel, NATO ASI Series C: Mathematical and Physical Science, Kluwer Academic Publishing, Dordrecht, 1994, vol. 435.
2. R. Paul and S. J. Paddison, Solid State Ionics, 2004, 168, 245.
3. D. Bratko, C. D. Daub, K. Leung and A. Luzar, J. Am. Chem. Soc., 2007, 129, 2504.
4. M. A. Ricci, F. Bruni, P. Gallo, M. Rovere and A. K. Soper, J. Phys.: Condens. Matter, 2000, 12, A345.
5. K. Lebedev, S. Mafé, A. Alcaraz and P. Ramírez, Chem. Phys. Lett., 2000, 326, 87.
6. B. Alberts, A. Johnson, J. Lewis, M. Raff, K. Roberts and P. Walter, Molecular biology of the cell, Garland Science, New York, 4th edn, 2002.
7. B. Hille, Ion Channels of Excitable Membranes, Sinauer Associates, Sunderland, MA, 3rd edn, 2001.
8. R. Benz, K. Janko and P. Läuger, Biochim. Biophys. Acta, 1979, 551, 238.
9. A. H. Delcour, Front. Biosci., 2003, 8, d1055.
10. H. Nikaido, Microbiol. Mol. Biol. Rev., 2003, 67, 593.
11. S. W. Cowan, R. M. Garavito, J. N. Jansonius, J. A. Jenkins, R. Karlsson, N. Konig, E. F. Pai, R. A. Pauptit, P. J. Rizkallah, J. P. Rosenbusch, G. Rummel and T. Schirmer, Structure, 1995, 3, 1041.
12. A. Alcaraz, E. M. Nestorovich, M. Aguilella-Arzo, V. M. Aguilella and S. M. Bezrukov, Biophys. J., 2004, 87, 943.
13. A. Alcaraz, P. Ramírez, E. García-Giménez, M. L. López, A. Andrio and V. M. Aguilella, J. Phys. Chem. B, 2006, 110, 21205.
14. P. Ramírez, M. Aguilella-Arzo, A. Alcaraz, J. Cervera and V. M. Aguilella, Cell Biochem. Biophys., 2006, 44, 287.
15. E. M. Nestorovich, T. K. Rostovtseva and S. M. Bezrukov, Biophys. J., 2003, 85, 3718.
16. P. S. Phale, A. Philippsen, C. Widmer, V. P. Phale, J. P. Rosenbusch and T. Schirmer, Biochemistry, 2001, 40, 6319.
17. M. S. P. Sansom, G. R. Smith, C. Adcock and P. C. Biggin, Biophys. J., 1997, 73, 2404.
18. G. Lamm and G. R. Pack, J. Phys. Chem. B, 1997, 101, 959.
19. M. Troll and B. H. Zimm, J. Phys. Chem., 1983, 87, 3197.
20. M. F. Perutz, Science, 1978, 201, 1187.
21. A.-S. Yang, M. R. Gunner, R. Sampogna, K. Sharp and B. Honig, Proteins, 1993, 15, 252.
22. C. N. Schutz and A. Warshel, Protein Struc. Funct. Genet., 2001, 44, 400.
23. A. Rubinstein and S. Sherman, Biophys. J., 2004, 87, 1544.
24. P. C. Jordan, in Ion channel reconstitution, ed. C. Miller, Springer, New York, 1986.
25. W. Cheng, C. X. Wang, W. Z. Chen, Y. W. Xu and Y. Y. Shi, Eur. Biophys. J., 1998, 27, 105.
26. M. Gutman, Y. Tsfadia, A. Masad and E. Nachiel, Biochim. Biophys. Acta, 1992, 1109, 141.
27. J. A. Ng, T. Vora, V. Krishnamurthy and S.-H. Chung, Eur. Biophys. J., 2008, 37, 213.
28. F. Booth, J. Chem. Phys., 1951, 19, 391.
29. M. Aguilella-Arzo, J. J. García-Celma, J. Cervera, A. Alcaraz and V. M. Aguilella, Bioelectrochemistry, 2007, 70, 320.
30. W. Im and B. Roux, J. Mol. Biol., 2002, 322, 851.
31. D. P. Tieleman and H. J. C. Berendsen, Biophys. J., 1998, 74, 2786.
32. A. A. Kornyshev and G. Sutmann, Phys. Rev. Lett., 1997, 79, 3435.
33. I.-C. Yeh and M. L. Berkowitz, J. Chem. Phys., 1999, 110, 7935.
34. H. A. Kotodziej and B. P. Jones, J. Chem. Soc., Faraday Trans. 2, 1975, 71, 269.
35. B. E. Conway, Ionic Hydration in Chemistry and Biophysics, Elsevier, Amsterdam, 1981.
36. C. J. F. Böttcher, Theory of Electric Polarization, Elsevier, Amsterdam, 2nd edn, 1993.
37. G. Sutmann, J. Electroanal. Chem., 1998, 450, 289.
38. M. L. Berkovitz, I.-C. Yeh and E. Spohr, in Interfacial Electrochemistry: Theory, Experiment and Applications, ed. Andrzej Wieckowski, Marcel Dekker, New York, 1999.
39. K. J. Irwin, S. M. Barnett and D. L. Freeman, J. Membr. Sci., 1989, 47, 79.
40. E. Glueckauf, Trans. Faraday Soc., 1964, 60, 1637.
41. R. Paul and S. J. Paddison, J. Phys. Chem. B, 2004, 108, 13231.
42. J. R. Bontha and P. N. Pintauro, Chem. Eng. Sci., 1994, 49, 3835.
43. M. E. Davis, J. D. Madura, B. A. Luty and J. A. McCammon, Comput. Phys. Commun., 1991, 62, 187.
44. N. A. Baker, D. Sept, S. Joseph, M. J. Holst and J. A. McCammon, Proc. Natl. Acad. Sci. U. S. A., 2001, 98, 10037.
45. M. K. Gilson and B. H. Honig, Biopolymers, 1986, 25, 2097.
46. F. Fogolari, A. Brigo and H. Molinari, J. Mol. Recognit., 2002, 15, 377.
47. M. S. P. Sansom, I. D. Kerr, J. Breed and R. Sankaramakrishnan, Biophys. J., 1996, 70, 693.
48. V. J. van Hijkoop, A. J. Dammers, K. Malek and M.-O. Coppens, J. Chem. Phys., 2007, 127, 085101.
49. F. Fogolari, P. Zuccato, G. Esposito and P. Viglino, Biophys. J., 1999, 76, 1.
50. A. B. Mamonov, R. D. Coalson, A. Nitzan and M. G. Kurnikova, Biophys. J., 2003, 84, 3646.
51. T. Bastug and S. Kuyucak, Biophys. J., 2003, 84, 2871.
52. M. Lund, B. Jönsson and C. E. Woodward, J. Chem. Phys., 2007, 126, 225103.
53. D. A. Karp, A. G. Gittis, M. R. Stahley, C. A. Fitch, W. E. Stites and B. García-Moreno E, Biophys. J., 2007, 92, 2041.
54. V. Levadny, V. M. Aguilella, M. Aguilella-Arzo and M. Belaya, Phys. Rev. E, 2004, 70, 041912.

---