Relative dielectric constants and selectivity ratios in open ionic channels

We investigate the effects of the relative dielectric coefficient on ionic flows in open ion channels, using mathematical analysis of reasonably general Poisson-Nernst-Planck type models that can include the finite sizes of ions. The value of the relative dielectric coefficient is of course a crucial parameter for ionic behavior in general. Using the powerful theory of singularly perturbed problems in applied mathematics, we show that some properties of open channels are quite insensitive to variation in the relative dielectric coefficient, thereby explaining such effects seen unexpectedly in simulations. The ratio between the total number of one ion species and that of another ion species, and the ratio between the flux of one ion species and that of another ion species do not depend significantly on the relative dielectric coefficient.

To understand these properties, the atomic structure of the channel protein is critically important (see, e.g., [11,12,13,16,28,37,92,104,106,105]). It determines qualitatively mechanisms of selectivity and permeation, and provides basis for high resolutions of integrating correlations in continuum models. This approach with structural features does not though deal quantitatively with the experimental measurements of current at all as can be seen by the absence of predicted graphs of currents or contents vs. voltage or concentration. The structures being discussed are determined from crystals, with somewhat different locations of atoms and forces (or they would not crystallize), in solutions remote from the physiological solutions in which the channels function, and often at temperatures around 100K. Qualitative discussion is of great importance but it must not be confused with attempts at quantitative models.
The quantitative treatment of selectivity and permeation differs from the usual qualitative treatment in several important respects. The quantitative approach based on continuum models treats selectivity and permeation as expressions of the current flow through a channel under a variety of conditions, namely different membrane potentials, different concentrations and compositions of ionic solutions, and different mutations. The quantities from computations/analyses of models can be compared directly with experimental measurements of current. The quantitative models are dramatically reduced and oversimplified to be sure, but they are precise. Such is the nature of most physical models of condensed phases. The utility of the models is determined not by discussion of their oversimplifications but by their ability to fit wide ranges of data with single sets of parameter, and to design new systems and predict their properties. Physicists tend to emphasize the importance of fitting data; structural biologists tend to emphasize the importance of structural detail. Both approaches are needed and such complement and compliment each other, in our view.
The essential point here is that models of the PNP type we have studied capture correlations well enough to explain the main experimental data from three major types of channels cited above (see [3] for more discussions). Ionic flow through biological channels exhibits extremely rich dynamics ( [19,20,37,38,96], etc.). The flow is multi-scale (in both space and time) and the general behavior depends on many parameters such as channel structures (channel shape and spatial distribution of permanent charges), boundary concentrations, electric potential differences, diffusion coefficients, dielectric properties, ionic sizes, and many others ( [3,4,17,18,21,22]). Some of the parameters affect each other, for example, the dielectric coefficients and diffusion coefficients are by no means constants and they vary as the environment changes. It is a great challenge to understand the full behavior of ionic flows. On the other hand, not all the parameters play equal role for every biologically interesting quantity. For example, the sign or direction of the flux of an ion species depends only on the difference of its electrochemical potentials at the boundary points of the channel ( [27]); in particular, the sign of the flux DOES NOT depend on the permanent charge (either its spatial occupancy and its density) and specifics of the dielectric coefficients; of course, the magnitude of the flux DOES depend on the permanent charge and specifics of dielectric coefficients as well as other physical parameters of the biological system. Another example concerns ion channels without permanent charge: the shape does not play a significant role in the magnitude of individual fluxes -an average property of the channel shape is the main determinant in this case ( [48]). This paper is not meant to produce complete models of any channel. Rather, motivated by the work in [3], we will examine a particular property of selectivity in this paper with the quantitative approach based on Poisson-Nernst-Planck (PNP) type models.
In this work, we are interested in the effects of the relative dielectric coefficient on ionic flow properties. The value of the relative dielectric coefficient is of course a crucial parameter for ionic behaviors in general. From the modeling point of view, the traditional Poisson-Boltzmann system with constant dielectric coefficients as in the Debye-Hückel theory has been modified with concentration-dependent dielectric coefficients and ion-water interaction that show significant improvements for the qualitative property of the mean activity coefficient (see, for example, [99,100]).
We will examine some properties of ionic movements and show they are not sensitive to changes in the relative dielectric coefficient. For example, if we consider the ionic flow through a biological channel (for example, sodium channel in particular the DEKA aspartate-glutamate-lysine-alanine channel responsible for the rising phase of nerve action potentials) involving the mixture of Na + Cl − and K + Cl − , the total number #(Na) of Na + ions and the total number #(K) of K + ions in the channel could depend on the value of the relative dielectric coefficient significantly, but the ratio #(Na)/#(K) between these two numbers does not depend significantly on the value of the relative dielectric coefficient under reasonable conditions. The latter has been observed in extensive Monte Carlo simulations (to the surprise of the authors, we are told, see Fig. 8-10 in [3]). Applying the powerful singular perturbation theory to Poisson-Nernst-Planck (PNP) type models that include finite sizes of ions, we provide analytical justifications when the ratio between fluxes of different ion species and the ratio between total numbers of particles of different ion species of the ionic mixture do not depend significantly on the relative dielectric coefficient, and when the ratios might depend significantly on the relative dielectric coefficient.
We emphasize that the dependences of the ratios on other physical parameters such as radius of the channel, protein structure, etc. are important. In the future, we will examine these dependences which require much more detailed analysis of PNP type models than that used in this paper for the present purpose.
We stress that our PNP type models can include finite sizes of ions and allow general distributions of permanent charges of the channel. The key condition for our results is that the standard dimensionless parameter ε defined in (2.4) is small. This dimensionless parameter ε is a version of the Debye length divided by the length l of the channel that reflects a combined effect of the relative dielectric coefficient ε r , the dielectric constant ε 0 of vacuum, and the characteristic concentration C 0 of the ionic mixture. The parameter ε is small mostly because of the length scale of the channel and the characteristic concentration of the ionic mixture. For general electrolytes problems, the value of ε could vary several orders of magnitudes. Our results do not apply when the value of ε is large or moderate. For ion channel problems, the value of ε would be small when the characteristic concentration C 0 and/or the length l of channel are large. We estimate ε to be about 10 −3 or smaller for a few biological situations (see the end of Section 3.1 and Remark 3.1).
We remark that when the value of the relative dielectric coefficient is varied, we assume other properties such as diffusion coefficients, permanent charge distribution and shape of the channel remain more or less unchanged. The flexibility of these parameters would require more extensive analysis of PNP type models with possible modifications of the models. We will examine these important issues in future work.

Three-dimensional model
For an ionic mixture with n ion species, PNP reads where r ∈ Ω. Ω is a three-dimensional cylindrical-like domain representing the channel, Q(r) is the permanent charge density, ε r (r) is the relative dielectric coefficient, ε 0 is the vacuum permittivity, e is the elementary charge, k B is the Boltzmann constant, T is the absolute temperature; Φ is the electric potential, and, for the kth ion species, C k is the concentration, z k is the valence (the number of charges per particle), µ k is the electrochemical potential depending on the concentrations {C j } and finite sizes of all ion species, J k is the flux density, and D k (r) is the diffusion coefficient.

Quasi-one-dimensional PNP models
Reduction of three-dimensional PNP systems (2.1) to quasi-one-dimensional models was first proposed in [82] using the fact that ion channels have narrow cross-sections relative to their lengths and was rigorously justified in [68] for special cases. A quasi-one-dimensional PNP model is where X ∈ [0, l] is the coordinate along the axis of the channel, A(X) is the area of cross-section of the channel over the location X. Equipped with system (2.2), we impose the following boundary conditions (see, [26] for a reasoning), for k = 1, 2, · · · , n,

Rescaling of the quasi-one-dimensional PNP model
The following rescaling (see [29]) or its variations have been widely used for convenience of mathematical analysis. To do so, we first introduce several quantities.
Let C 0 be a characteristic concentration of the ion solution; for example, one can take C 0 to be the maximal value of the boundary concentrations and the absolute value of the permanent charge: Using these quantities, we now make a dimensionless re-scaling of the variables in system (2.2) as follows.
In terms of the new variables, the BVP (2.2) and (2.3) becomes with boundary conditions at x = 0 and x = 1 (2.6) Our conclusion will be based on properties of the solutions of the boundary value problem (BVP) associated to equation (2.5) and the boundary conditions (2.6). The mathematical consequence of the key assumption of smallness of ε is that the BVP can be treated as a singularly perturbed problem. A general geometric framework for analyzing the singularly perturbed BVP of PNP type systems has been developed in [26,48,66,67,69] for classical PNP systems and in [47,57] for PNP systems with finite ion sizes.

Relevant properties of the singularly perturbed BVP and its consequences
The essential result for the singularly perturbed BVP (2.5) and (2.6) is (see, e.g., [26,66,67]), for small ε > 0, if is the solution, then, for each 1 ≤ k ≤ n, one can express c k (x; ε) as where c k0 (x) is independent of ε and is continuous at all x except at those points where the permanent charge Q(x) jumps (discontinuity points of Q(x)).
Similarly, for each 1 ≤ k ≤ n, where J k0 is a constant independent of ε.

The above property for the concentrations in (3.1) is illustrated by a numerical solution in
Concerning C k , it is noted from (2.4) that, when C 0 is large, is small, but εC 0 = ε r ε 0 k B T C 0 e 2 l 2 is NOT small in general (although, relative to C 0 , it is small). Therefore, ignoring the second term εC 0 c k1 (x; ε) in (3.3) for the value of C k (x; ε) is tricky. Similar remark applies to the flux J k (ε).

Ratio of total numbers of different ion species
We now apply the results in (3.3) to the ratio of total numbers of two different ion species. For definiteness, we take Na + and K + ions.
Let c 1 (x; ε) and c 2 (x; ε) be the concentrations of Na + and K + in terms of the scaled variables, respectively. Then, the total numbers of Na + and K + ions in the channel are, respectively, It follows from X = lx and C k (X; ε) = C 0 c k (x; ε) in (2.4) that

(3.4)
From extensive Monte Carlo simulations for DEKA Na channel in [3], it shows (Fig. 8-10 in [3]) that the ratio N 1 /N 2 decreases as the channel diameter increases BUT the ratio does not depend onε r significantly, say, asε r changes from 1 to 80. We now provide an explanation to the latter from formulas (3.1) and (3.3).  and hence, the ratio between the total number of Na + and that of K + is    then the coefficient of the ε order term is of order O(1). In particular,

It follows from (3.1) or (3.3) that
with the right-hand side independent of ε and with an error bound in the order O(ε). Note that, in the derivation (3.6), the large quantity C 0 is dropped when ratio N 1 (ε)/N 2 (ε) is considered.
We now use the following specific values of the quantities involved to estimate the value of ε in (2.4) for an assessment of the error bound.
• The length between the two baths containing the channel is l(nm) = l × 10 −9 (m); • e ≈ 1.60 × 10 −19 (C), ε 0 ≈ 8.85 × 10 −12 (F·m −1 ); Then, from (2.4), If we take l = 2.5 (the length of the channel would be 2.5 (nm)) and C 0 = 10 (the maximum of concentrations including the concentration of the permanent charge would be 10 (M)), then ε ≈ √ε r × 5.84 × 10 −3 (3.9) As the relative dielectric coefficientε r changes from 1 to 80, the value of ε changes by a factor of √ 80 ≈ 9. If the length l of the channel is 2.5 (nm) and the maximum of the concentration C 0 of the ionic solution is 10 (M), then ε is of order 10 −3 (see (3.9)) so a factor of 9 in the variation ofε r makes about 5% difference in the approximation (3.8) and is not significant.
Remark 3.1. Our quasi-one-dimensional model of a channel setup includes access regions and baths, following [68,82]. The variable diameter (equivalently, the variable cross sectional area A(X) in (2.2)) of the 'channel' in this setup includes a region with the diameter of the actual pore of the channel protein, including the taper often seen in structures as they join the extracellular space. In our setup, the diameter of the 'channel' is dramatically increased outside the channel protein to mimic (crudely but adequately for most purposes) the huge increase in the cross sectional area available for current flow, as ions leave the channel protein and enter the extracellular space with its macroscopic cross sectional area. What is important is that the resistance to electric current flow and effective 'resistance' to diffusional flow is millions of times less in the bath region than inside the channel protein. Thus, in this treatment the quantity l is the length between the two electrodes where the boundary conditions (2.3) are imposed. The channel would correspond to a subinterval of [0, l], or in the scaled variable x, to a subinterval, say [a, b], of the interval [0, 1]. For this setting, c 1 (x; ε) and c 2 (x; ε) would be different from those in the previous setting. The total numbers of Na + and K + ions in the channel would be, respectively, Following the lines of derivation in (3.6), the ratio between the total number of Na + and that of K + would be Hence, with the right-hand side independent of ε and with an error bound in the order O(ε).
Since the quantity l is taken to be the length between the two electrodes where the boundary conditions (2.3) are imposed, it would be much larger than the length of the channel. For example, if we take l = 10 (nm), which is four times the length 2.5 (nm) used for l in the estimation in (3.9), then the quantity ε in (3.9) would be reduced by a factor of 1/4. As a result, the factor of 9 caused by the change ofε r from 1 to 80 makes about 1.25% difference in the approximation (3.10).
When conclusion (3.8) might fail. Recall that our conclusion, from (3.6) to (3.8), that the ratio does not depend significantly on the dielectric coefficient generally requires also condition (3.7), in addition to ε being small. It is clearly possible that condition (3.7) does not hold; for example, in [3], when selectivity of Na + vs. Ca 2+ was examined for DEKA Na channel, a SMALL amount (1mM) of CaCl 2 was added to NaCl (100 mM) bulk solution so that Ca-concentration inside the channel is of order O(10 −3 ) relative to Na-concentration (see Fig. 3 in [3]). In this case, the conclusion (3.8) might fail. This is not surprise and an explanation for this possible failure goes as follows. Suppose, on the contrary, that the conclusion (3.8) still holds for this case. Since the numerator of the right-hand side of (3.8) is of order O(1) due to the dimensionless scaling, the denominator would be order O(10 −3 ) which is comparable to ε; that is, the condition (3.7) does not hold. Going back to formula (3.6), the second term, due to its much smaller denominator, might not be small relative to the first term anymore. In this situation, the approximation (3.8) for the ratio fails in general; that is, the ratio would depend on the relative dielectric coefficient in more significant way. In fact, this is indeed the case as shown in Fig. 10 in [3]. This is consistent with the requirement of (3.7) in our result. Looking closer at Fig. 10 in [3], from if |J 20 | is not small relative to ε, (3.12) with the right-hand side independent of ε and with an error bound in the order of O(ε). Hence, when ε is small, the ratio of the two fluxes does not depend on the value of dielectric coefficient significantly.

Conclusion
Some continuum theories of ion conduction allow reasonably general conclusions concerning physiological parameters and their dependence on underlying physical parameters. The PNP theory of open ionic channels can be analyzed by quite general mathematics to understand the dependence of selectivity on the relative dielectric coefficients. In this short paper, based on PNP type models that allow finite sizes of ions, we provide theoretical justifications when the ratio of the total numbers of two ion species does not depend on the relative dielectric coefficient in a significant way and when it might. A similar result is established for the ratio of fluxes of two ion species. In particular, the result explains that the selectivity does not depend on the relative dielectric coefficients in a significant way for ionic mixtures of KCl and Na + Cl − but may depend on the relative dielectric coefficients in a significant way for ionic mixtures of Na + Cl − and Ca ++ Cl − 2 which has been observed from extensive Monte Carlo simulations for DEKA Na channel in [3]. We believe that PNP theory of open ionic channels can provide an alternative way to understand ion channel properties.