582 DOI: 10.17344/acsi.2015.1561
Acta Chim. Slov. 2015, 62, 582-587
Scientific paper
Electrolyte Solution at Zwitterionic Lipid Layer
Jurij Re{~i~1* and Klemen Bohinc2
1 Faculty of Chemistry and Chemical Technology, University of Ljubljana, Ljubljana, Slovenia 1 Faculty of Health Sciences, University of Ljubljana, Ljubljana, Slovenia * Corresponding author: E-mail: jurij.rescic@fkkt.uni-lj.si
Received: 26-03-2015 Dedicated to prof. Jože Koller on the occasion of his 70th birthday.
Abstract
A coarse-grained model of simple monovalent electrolyte solution in contact with a zwitterionic lipid layer in continuum solvent is studied by canonical Monte Carlo computer simulations and extended Poisson-Boltzmann theory. A structure of zwitterionic layer as well as concentration profiles of positively and negatively charged monovalent ions were obtained from simulations and compared to theoretical predictions. A relatively good agreement between the Monte Carlo computer simulations and theory was observed.
Keywords: Poisson-Boltzmann, zwitterionic lipid layer, Monte Carlo simulation, rod-like ion
1. Introduction
Cationic and zwitterionic lipids have been considered in the past because of their potential use as nonviral carriers of DNA into living cells.1-3 The pure cationic lipids in bilayer interact with negatively charged DNA electrostatically and form polyplexes4. The major weaknesses of using cationic lipids are the relatively low trans-fection efficiency and cytotoxicity.5 On the contrary, naturally occurring zwitterionic lipids are biodegradable and nontoxic, but the attractive interaction between zwitterio-nic lipid layers and DNA can be induced by multivalent
6-9
ions.6 9
The classical mean field approach to describe electrostatic interactions between charged systems in electrolyte solutions is the Poisson-Boltzmann (PB) method. Within the PB theory ions are modelled as point charges, while the solvent (water) is accounted for by a uniform dielectric constant. The charged surfaces are considered as uniformly charged. In addition, correlations between charges are not taken into account.10-12
The charged groups in the solution are neither pointlike nor can the solvent be regarded as passive and featu-reless.13 Many attempts have been made how to improve the PB theory. A formidable improvement is the modified PB theory, capable of treating electrolytes with asymmetry in both size and valency of simple ions.14
However, more complex ions which in biological systems mediate interactions between macroions have an internal structure with spatially separated charges and possibly with additional internal degrees of freedom.15,16 First studies were made with rod-like and spherical ions.17,18 Later on the features of zwitterionic lipids were incorporated into PB theory.19,20 In the theory it was assumed that negative charges of zwitterionic lipids reside within the membrane interface. These negative charges are connected to positive charges that have considerable conformational freedom to move around the negative charges without penetrating into the hydrocarbon core of the lipid layer. Recently the adsorption of macroions onto zwitterionic lipid layers was studied.21-23
Being established as a powerful and exact tool for solving model systems, MC computer simulations were used in the present work to test theoretical prediction of a model zwitterionic systems described below. Various static properties can be obtained by MC simulations, including thermodynamics and structure of model systems. First MC simulations for electric double layer were performed in the eighties of the previous century.24,25 Later on, MC simulations were performed for systems which include molecules with spatial charge distribution.17-19,26 A very good agreement between theoretical predictions and MC data for both rod-like and spherical ions was obtained.
The aim of the present work is to test the extended Poisson-Boltzmann theory with Monte Carlo computer simulations, both utilized to study a coarse-grained model system composed of a zwitterionic lipid layer in contact with a simple monovalent electrolyte solution. Concentration profiles of positively charged moieties of zwitterionic lipid layer as well as concentration profiles of positively and negatively charged monovalent ions will be compared for a range of parameters, with some of them being stretched over real values to test the extended Poisson-Boltz-mann theory more rigorously.
2. 1. Model and Methods
A model system consists of an infinite plane with zwitterionic molecules attached to it. Zwitterions are modelled as rigid rodlike molecules with one point-like positive charge and one point-like negative charge at each end, separated by the bond of length l = 0.5 nm. Examples of zwitterionic lipids are phosphatidylcholine (PC) and phosphatidylethanolamine (PE). Zwitterions are anchored to the plane with its negatively charged part. Zwitterion's cross-section a is a variable parameter and represents the surface area per one zwitterion from which the surface charge density a can be calculated. Such a layer is in contact with a reservoir of a simple monovalent salt whose concentration was expressed in terms of the Debye
length, defined as lD -	J^/fjžf) . In this equa-
tion e0 is the permittivity of the vacuum, e the relative permittivity of water, which is 78.4 at 25 °C, kB is the Boltz-mann's constant, T is the absolute temperature, e0 is the elementary charge, n° is the bulk number density of species i and zi its valency. The positive part of a zwitterion can move on a hemisphere with radius l. For this study we used values of a among 0.65, 2.60, and 5.85 nm2. The corresponding surface charge densities a were -0.2465, -0.0616, and -0.02739 C/m2, respectively. The values of lD were 0.5, 1, and 2 nm.
Since computer simulations cannot be performed on systems with point-like positive and negative charges being simultaneously present in the model system, we had to adjust the model accordingly. First, we assigned finite radii of 0.125 nm to both positive and negative charges. Second, the system size had to be limited so that the number of particles was not too large. Due to restrictions of the simulation software zwitterions were modelled as two charges, connected by a harmonic bond with large bond constant of 20 N/m and equilibrium distance of l = 0.5 nm, therefore not being completely rigid. As in the theoretical approach, negative charges were represented by a uniform surface charge density a, while the positive parts were charged explicitly. All charged specii interact via the Coulomb potential, while the solvent was treated as a dielectric continuum with the dielectric permittivity value of
78.4, typical for aqueous solutions at room temperature 298.15 K. One should note that a coarse-grained model without explicit water molecules cannot predict certain zwitterionic lipid layer features such as membrane's electrostatic potential.
2. 2. Theory
We introduce a Cartesian coordinate system whose x axis is oriented perpendicular to the zwitterionic lipid layer, which is located at x = 0. Due to sufficiently large planar surface and the translational invariance of the system along y and z directions, we can describe the systems with functions depending only on the x coordinate. Figure 1 is an illustration of the model system.
Figure 1. Schematic presentation of zwitterionic lipid layer. The length of head group is l, whereas the cross-section area per lipid
molecule is a.
Each zwitterionic lipid molecule consists of one negatively charged phosphate and one positively charged amino group. The phosphate group is linked through a glycerol backbone to the hydrocarbon tails and is spatially constrained with respect to motion along the normal direction of the membrane. Here in our calculations we assume that the phosphate groups lie in the plane x = 0. The phosphate and amino groups are separated by a fixed distance l. Each amino group can freely rotate on a semi-sphere around the adjacent phosphate group. The solution is composed of negatively and positively charged point-like ions.
The electrostatic free energy of the system, F, measured per unit area A and expressed in units of the thermal energy kBT (here kB is the Boltzmann's constant and T is the absolute temperature) is given by20,21
where lR=e\l 4яее0квТ is the Bjerrum length, which is in aqueous solutions at a room temperature eual to 0.715 nm. The prime denotes the derivation with respect to x.
The first term in eq (1) corresponds to the energy stored in the electrostatic field, here expressed in terms of the commonly used dimensionless electrostatic potential ¥ = е0Ф / kBT instead of the electrostatic potential Ф. The function	in the second
term of the integrand in Eq. 1 describes the mixing free energy of the positive (i = +) and negative ions (i = -) and n0 is their bulk value. The last term of Eq. 1 is related to the orientational entropy of the zwitterionic headgroups. The function W(x) expresses the probability to find the projection of the amino group on the lipid layer normal direction within the interval (x, x + dx). ф is the fraction of charged lipids at x = 0.
In thermal equilibrium the free energy F = F[ni,W] adopts its minimum with respect to the local concentration of ions ni and probability distribution W. Employing the first variation of Poisson's equation and demanding the condition of vanishing first variation 5F[n] = 0 we obtain the number densities of charges. Inserting the equilibrium distributions into Poisson's equation leads to the Poisson-Boltzmann equation20,21
I i
where q = - J e~'vdx is a partition function. The last term in
I «
eq. 2 contributes only in the lipid layer region. The Debye length is denoted by lD, whereas lc = a / 2 nlB.
Equation (2) has two boundary conditions. The first follows from the electro-neutrality of the system (¥'(x = 0) = 2 / lc), whereas the second one is given far from the zwitterionic lipid layer (¥'(*>) = 0).
2. 3. Monte Carlo Computer Simulation
Canonical Monte Carlo simulations were performed using the integrated Monte Carlo / Molecular Dynamic / Brownian dynamic simulation suite Molsim27 following the standard Metropolis scheme. Total of 2000 zwitterions were placed on the two opposite walls of the Monte Carlo simulation box, while small ions were placed randomly in it. The number of small ions was determined form theoretical density profiles, and varied from 200 to over 3000. Size of the MC box was calculated according to the chosen surface charge density.
A trial move for negative part of a zwitterion was limited to the two dimensions only, while for the positive part and for small ions the move was performed in 3D. Displacement parameters were chosen to obtain approximately 50% acceptance rate. 20,000 attempted moves per particle were used for equilibration followed by
100,000 attempted moves during production runs. Inter-particle interactions were calculated as described elsew-here.24,25 Long range corrections due to ionic distribution outside the MC box are found to be small and were therefore not used in present simulations. To calculate single particle distributions, the x-axis was always divided into 200 bins. The standard deviation of values in histograms was less than 0.5% for each separate bin in all cases.
3. Results and Discussion
Primary focus of present work is to investigate structure of an electrolyte solution next to the zwitterionic layer. Several parameters were varied including surface area per zwitterion, salt concentration, and fraction of ionized zwitterions. Salt concentration of bulk solution which is in contact with the zwitterionic layer, is most often subjected to changes. In the present study we used three different bulk salt concentrations of 0.37 mol/L, 0.092 mol/L, and 0.023 mol/L, which for aqueous solutions at 298.15K translate to the Debye lengths of 0.5 nm, 1 nm, and 2 nm, respectively. Figure 2 displays density profiles of both cations and anions as a function of perpendicular distance x from the plane with negative parts of zwitterions attached to it. The value of lD equal to 1 nm. Since the plane is negatively charged, anions are repelled from it, while cations tend to accumulate next to it. Anions' density profiles exhibit a broad maximum at x = 0.5 nm, caused by positive parts of zwitterions, while at the same distance the distribution of cations reaches its weak minimum. The same trend can be observed at all surface charge densities. The same structural features are predicted by both Poisson-Boltzmann theory and Monte Carlo computer simulations.
Similar structural features are shared throughout the whole salt concentration interval studied. Figure 3 and Figure 4 show density profiles at lD values 2 nm and 0.5 nm, respectively. Absolute values of each density profile increases with the increasing salt concentration. Interestingly though, the agreement between Monte Carlo computer simulations and Poisson-Boltzmann theory is best at the largest salt concentration.
We have also considered few cases where certain fraction of zwitterions gets its positive part neutralized while the negative part remains charged. This is introduced through the parameter ф. Value of ф equal to zero means that all zwitterions are electroneutral, while at ф = 0.5 half of the zwitterions are replaced by negative charges residing on the plane and corresponding number of mobile positive counterions. Resulting profiles are shown in Figure 5. The distinct features at ф = 0 are maxima and minima at x = 0.5 nm, which disappear when small fraction of zwitterions become charged. On the other hand, the case with ф = 1 is a system with charged planar slit in contact with reservoir of
a)
Figure 2. Density profiles of small mobile ions next to the zwitterio-nic layer at Debye length equal to 1nm. Corresponding monovalent salt concentration is 0.092 mol/L. Symbols represent Monte Carlo data and lines theory. Legend: + and full line : negative ions at a =
0.65 nm2 (c = -0.2465 C/m2); x and---: negative ions at a = 5.85
nm2 (c = -0.0274 C/m2) О and......: positive ions at a=0.65 nm2; □
and - - - : positive ions at a = 5.85 nm2.
0.14
0.12 0.1 0.03 0.06 0.04 0.02
1	1 1 -anions, a=0.65nm	•
»	anions, a=5.05nm2	X -
	cations, a=0.65niin2	a
	cations, a=5.85nrm2	0
■'"til во, i ; ;		-
t Sui .		
0 5
1
x/nm
1.5
Figure 3. Density profiles of small mobile ions next to the zwitte-rionic layer at Debye length equal to 2 nm. Notation is the same as for Figure 2, except for the corresponding monovalent salt concentration, which is 0.023 mol/L.
Figure 4. Density profiles of small mobile ions next to the zwitte-rionic layer at Debye length equal to 0.5 nm. Notation is the same as for Figure 2, except for the corresponding monovalent salt concentration, which is 0.37 mol/L.
b)
Figure 5. Density profiles of small mobile ions next to the zwitte-rionic layer at Debye length equal to 1 nm and at a = 0.65 nm2 for three values of ф. Symbols represent Monte Carlo data and lines
theory. Legend: □ and full line: ф = 0; О and---: ф =0.25; A and
---: ф = 0.5; V and ...: ф = 1. Note the value of ф = 0 means that
all zwitterions are electroneutral, whereas the value 0 means that there are no zwitterions at all, only negative phosphate groups remain on the plate. Panel A displays cations' density profiles while the panel B density of anions next to the zwitterionic layer.
K/nm
Figure 6. Density profiles of positive part of zwitterions at Debye length equal to 1 nm. Corresponding monovalent salt concentration is 0.092 mol/L. Symbols represent Monte Carlo data and lines theory. Legend: □ and full line: a = 0.65 nm2 (c = -0.2465 C/m2); О and - - -: a = 2.60 nm2 (c = -0.0616 C/m2); A and - - -: a = 5.85 nm2 (c = -0.0274 C/m2).
simple 1:1 salt. Density profiles for this system are also shown in Figure 5 for the same value of a and lD. Please note the logarithmic y-axis.
Finally, distribution of positive parts of zwitterions was also studied for all range of parameters lD and a. Results for lD = 1nm and a = 0.65nm2 are presented in Figure 6. Note that the positive part can only move on a hemisphere with radius of 0.5 nm. The agreement between theory and simulation is excellent in this case.
Being a basic ingredient of biological membranes and therefore an attractive model system, zwitterionic lipid layers have been studied previously by several techniques ranging from experimental to theoretical ones. The latter includes molecular dynamic28 and Monte Carlo computer simu-lations29,30 at various levels of model details. An all-atom Monte Carlo computer simulation study of dimrystolphosp-hatidylcholine bilayer in contact with aqueous solution containing tetramethylammonium dimethylphosphate in isot-hermal-isobaric ensemble30 exhibit qualitatively similar radial distribution function (not shown) of both zwitterionic ends in the bilayer as was observed in the present study.
4. Conclusions
A structure of a simple monovalent electrolyte solution next to a lipid zwitterionic layer was studied by both the extended Poisson-Boltzmann theory capable of treating rigid ions with spatially separated charges and by canonical Monte Carlo computer simulations. Several aspects have been considered, including various surface density of lipid molecules, salt concentration, and ionization of zwitterions. Negatively charged phosphate groups form a charged plane, while positive parts of zwitterions as well as small ions in a solution are spatially distributed due to their mobility. A significant in-
crease in cation concentration next to the zwitterionic layer is observed, while anions tend to accumulate next to the positive parts of zwitterions. At a distance of approximately 1 nm measured perpendicular to the zwitterionic layer, concentrations of all mobile ions reach their average value, regardless of salt concentration. Monte Carlo computer simulations were performed on the model system in order to test theoretical predictions. Agreement between them was found to be ranging from semi-quantitive to very good.
5. Acknowledgment
Authors wish to thank the Slovenian Research Agency for support through grant P1-0201.
6. References
1.	P. P. Karmali and A. Chaudhuri, Med. Res. Rev. 2007, 27, 696-722. http://dx.doi.org/10.1002/med.20090
2.	V. B. Teif, and K. Bohinc, Progr. Biophys. Mol. Biol. 2011, 105, 208-222.
http://dx.doi.org/10.1016/j.pbiomolbio.2010.07.002
3.	P. L. Felgner and G. M. Ringold, Nature 1989, 331, 461- 462.
4.	J. O. Raedler, I. Koltover. T. Salditt, and C. R. Safinya, Science 1997, 275, 810-814.
5.	C. R. Dass, J. Pharm. Pharmacol. 2002, 54, 593-601. http://dx.doi.org/10.1211/0022357021778817
6.	J. J. McManus, J. O. Radler, and K. A. Dawson, J. Phys. Chem. B 2003, 107, 9869-9875. http://dx.doi.org/10.1021/jp034463d
7.	M. Pisani, P Bruni, G. Caracciolo, R. Caminiti, and O. Fran-cescangeli, J. Phys. Chem. B 2006, 110, 13203-13211. http://dx.doi.org/10.1021/jp062713v
8.	D. Uhrikova, A. Lengyel, M. Hanulova, S. S. Funari and P. Balgavy, Eur. Biophys. J. 2007, 36, 363-375. http://dx.doi.org/10.1007/s00249-006-0086-2
9.	D. H. Mengistu, K. Bohinc, and S. May, J. Phys. Chem. B 2009, 113, 12277-12282. http://dx.doi.org/10.1021/jp904986j
10.	S. McLaughlin, Ann. Rev. Biophys. Chem. 1989, 18, 113-136. http://dx.doi.org/10.1146/annurev.bb.18.060189.000553
11.	E. J. Vervey and J. T. G. Overbeek: Theory of the Stability of Lyophobic Colloids; Elsevier: Amsterdam, 1948.
12.	L. Guldbrand, B. Jönsson, H. Wennerström, and P. Linse, J. Chem. Phys. 1984, 80, 2221- 2228. http://dx.doi.org/10.1063/L446912
13.	J. Reščič, K. Bohinc, Acta Chim. Slov.. 2012, 59(3), 601608.
14.	C. W. Outhwaite, M. Molero, and L. B. Bhuiyan, J. Chem. Soc. - Faraday Trans. 1991, 87, 3227-3230. http://dx.doi.org/10.1039/ft9918703227
15.	K. Bohinc, A. Iglič, and S. May, Europhys. Lett. 2004, 68(4), 494-500.
http://dx.doi.org/10.1209/epl/i2004-10250-2
16.	S. Maset and K. Bohinc, J. Phys. A: Math. Theor. 2007, 40, 1 http://dx.doi.org/10.1088/1751-8113/40/39/008
17.	K. Bohinc, J. Re{~i~, S. Maset, and S. May, J. Chem. Phys. 2011, 134, 074111. http://dx.doi.org/10.1063/L3552226
18.	S. May, A. Igli~, J. Re{~i~, S. Maset, and K. Bohinc, J. Phys. Chem. B 2008, 112, 1685-1692. http://dx.doi.org/10.1021/jp073355e
19.	K. Bohinc, J. Re{~i~, J. F. Dufreche, and L. Lue, J. Phys. Chem. B, 2013, 117 (36), 10846-10851. http://dx.doi.org/10.1021/jp404822f
20.	E. C. Mbamala, A. Fahr, and S. May, Langmuir 2006, 22, 5129-36. http://dx.doi.org/10.1021/la060180b
21.	D.H. Mengistu, A. Fahr and S. May, J. Chem. Phys. 2008, 129, 121105-3. http://dx.doi.org/10.1063/L2990746
22.	K. Bohinc, G. Brezesinski and S. May, Phys. Chem. Chem. Phys. 2012, 14, 10613-10621. http://dx.doi.org/10.1039/c2cp40923b
23.	A. Haugen and S. May, J. Chem. Phys. 2007, 127, 2151048. http://dx.doi.org/10.1063/L2803075
24.	G. M. Torrie and J. P. Valleau, Chem. Phys. Lett. 1979, 65, 343-346. http://dx.doi.org/10.1016/0009-2614(79)87078-5
25.	B. Jönsson, H. Wennerström, and B. Halle, J. Phys. Chem. 1980, 84, 2179-2185. http://dx.doi.org/10.1021/j100454a014
26.	J. Re{~i~ and P. Linse, J. Phys. Chem. B 2000,104, 78527857. http://dx.doi.org/10.1021/jp0007585
27.	J. Re{~i~ and P. Linse, J. Comput. Chem. 2015, 36, 12591274. http://dx.doi.org/10.1002/jcc.23919
28.	A. A. Polyansky, P.E. Volynsky, D. E. Nolde, A. S. Arseniev, and R. G. Efremov, J. Phys. Chem. B 2005, 109, 1505215059. http://dx.doi.org/10.1021/jp0510185
29.	U. Nilsson, B. Jönsson, and H. Wennerström, Adv. Coll. Int. Sci. 1992, 41, 179-195.
http://dx.doi.org/10.1016/0001-8686(92)80012-M
30.	P. Jedlovszky, L. Partay, and. M. Mezei, J. Mol. Liq. 2007, 131-132, 225-234.
http://dx.doi.org/10.1016/j.molliq.2006.08.040
Povzetek
S kanoni~no Monte Carlo simulacijo in raz{irjeno Poisson-Boltzmannovo teorijo smo preu~evali model zwitterionske lipidne plasti, ki je v stiku z raztopino preprostega 1:1 elektrolita. Z obema pristopoma smo izra~unali porazdelitev kationov in anionov pravokotno na lipidno plast ter strukturo same zwitterionske lipidne plasti. Ugotovili smo, da se po-razdelitvene funkcije kationov in anionov, dobljene s simulacijami, relativno dobro ujemajo z rezultati teorije. Strukturo lipidne plasti pa obe metodi napovesta skoraj identi~no.