| Structure and Interactions of Aggrecans: Statistical Thermodynamic Approach Biophysical Journal, Volume 95, Issue 10, 15 November 2008, Pages 4570-4583 Rikkert J. Nap and Igal Szleifer Abstract Weak polyelectrolytes tethered to cylindrical surfaces are investigated using a molecular theory. These polymers form a model system to describe the properties of aggrecan molecules, which is one of the main components of cartilage. We have studied the structural and thermodynamical properties of two interacting aggrecans with a molecular density functional theory that incorporates the acid-base equilibrium as well as the molecular properties: including conformations, size, shape, and charge distribution of all molecular species. The effect of acidity and salt concentration on the behavior is explored in detail. The repulsive interactions between two cylindrical-shaped aggrecans are strongly influenced by both the salt concentration and the pH. With increasing acidity, the polyelectrolytes of the aggrecan acquire charge and with decreasing salt concentration those charges become less screened. Consequently the interactions increase in size and range with increasing acidity and decreasing salt concentration. The size and range of the forces offers a possible explanation to the aggregation behavior of aggrecans and for their ability to resist compressive forces in cartilage. Likewise, the interdigitation of two aggrecan molecules is strongly affected by the salt concentration as well as the pH. With increasing pH, the number of charges increases, causing the repulsions between the polymers to increase, leading to a lower interdigitation of the two cylindrical polymer layers of the aggrecan molecules. The low interdigitation in charged polyelectrolytes layers provides an explanation for the good lubrication properties of polyelectrolyte layers in general and cartilage in particular. Abstract | Full Text | PDF (486 kb) |
| Effect of the Ionic Strength and pH on the Equilibrium Structure of a Neurofilament Brush Biophysical Journal, Volume 93, Issue 5, 1 September 2007, Pages 1452-1463 E.B. Zhulina and F.A.M. Leermakers Abstract Using the numerical model of Scheutjens and Fleer, we investigated, on a self-consistent field level, the equilibrium structure of the neurofilament brush formed by projection domains of the constituent NF-H, NF-M, and NF-L proteins. The phosphorylation of such a brush is a major regulatory process that triggers the relocation of the H tails from the NF core to the brush periphery. We explore how the pH and the ionic strength affect the rearrangements in the NF brush structure upon phosphorylation. We demonstrate that the translocation of H tails in an individual NF occurs as a sharp cooperative transition below and up to the physiological salt concentration. Regularities of this process are reminiscent of the collapse-to-stretching transition in a cylindrical polyelectrolyte brush in a poor solvent. The effect of pH at physiological ionic strength is noticeable only in the acidic range and is more pronounced for a dephosphorylated NF. Abstract | Full Text | PDF (248 kb) |
| Polylysine-Induced H NMR-Observable Domains in Phosphatidylserine/Phosphatidylcholine Lipid Bilayers Biophysical Journal, Volume 81, Issue 6, 1 December 2001, Pages 3346-3362 Carla M. Franzin and Peter M. Macdonald Abstract The interaction of three polylysines, Lys (=5), Lys (=30), and Lys (=100), where is the number of lysine residues per chain, with phosphatidylserine-containing lipid bilayer membranes was investigated using H NMR spectroscopy. Lys and Lys added to multilamellar vesicles composed of (70:30) (mol:mol) mixtures of choline-deuterated 1-palmitoyl-2-oleoyl--glycero-3-phosphocholine (POPC)+1-palmitoyl-2-oleoyl--glycero-3-phosphoserine (POPS) produced two resolvable H NMR spectral components under conditions of low ionic strength and for cases where the global anionic lipid charge was in excess over the global cationic polypeptide charge. The intensities and quadrupolar splittings of the two spectral components were consistent with the existence of polylysine-bound domains enriched in POPS, in coexistence with polylysine-free domains depleted in POPS. Lys, however, yielded no H NMR resolvable domains. Increasing ionic strength caused domains to become diffuse and eventually dissipate entirely. At physiological salt concentrations, only Lys yielded H NMR-resolvable domains. Therefore, under physiological conditions of ionic strength, pH, and anionic lipid bilayer content, and in the absence of other, e.g., hydrophobic, contributions to the binding free energy, the minimum number of lysine residues sufficient to produce spectroscopically resolvable POPS-enriched domains on the H NMR millisecond timescale may be fewer than 100, but is certainly greater than 30. Abstract | Full Text | PDF (216 kb) |
Copyright © 2005 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 89, Issue 5, 2972-2987, 1 November 2005
doi:10.1529/biophysj.105.068387
Biophysical Theory and Modeling
Shelly Tzlil and Avinoam Ben-Shaul
, 
Address reprint requests to A. Ben-Shaul, Tel.: 972-2-658-5271.The lipid bilayer, constituting the central structural element of biological membranes, is a two-dimensional fluid mixture, composed typically of many lipid species. Owing to the lateral mobility of lipids at physiological temperatures, the membrane can respond to interactions with integral and peripheral macromolecules by mobilizing those lipids interacting favorably with the macromolecule into the interaction zone. This process leads to local changes in lipid composition around the guest molecules which, under certain conditions, may evolve into larger-scale reorganization of membrane components, resulting in domain formation. The molecular composition and phase characteristics of the domains, as in lipid rafts, are different from those of the surrounding membrane 1,2,3.
The ability of integral proteins to induce local and global changes in lipid composition has been extensively documented experimentally 4, and amply analyzed theoretically 5,6. Similarly, experiments reveal that when charged macromolecules, such as certain kinds of proteins or DNA, are adsorbed onto a mixed membrane containing a small amount of oppositely charged lipids, the charged species migrate toward the adsorbed macromolecule 7,8,9, tending to achieve local electrical neutrality. Although lowering the electrostatic (free) energy of the system, the segregation of charged lipids induced by the peripheral macromolecule may involve a non-negligible entropic penalty. Several recent theoretical studies have carefully analyzed the energetic-entropic balance associated with electrostatic adsorption of rigid macromolecules (e.g., DNA and globular protein) onto fluid membranes 10,11,12,13,14,15,16. It was shown, for instance, that the extent of lipid segregation, the corresponding entropy loss, and the interaction free energy depend sensitively on the shape, charge, and concentration of the adsorbing macromolecule 17,18.
This work focuses on the energetic and structural characteristics of the interaction between flexible, electrically charged, macromolecules (polyelectrolytes), and mixed, oppositely charged, fluid membranes. We shall consider three-component membranes composed of one electrically neutral species and two differently charged lipids. Electrostatic adsorption on such membranes may involve significant changes in the spatial configuration of the adsorbing polyelectrolyte and, consequently, a substantial loss of conformational entropy. The two kinds of entropy loss, those associated with lipid segregation and those which lower the macromolecule’s conformational freedom, tend to offset the gain in electrostatic interaction energy between the oppositely charged molecules. It should be noted, however, that both degrees of freedom, namely, lipid mobility and macromolecule flexibility, enable the interacting complex to select, with higher probability, those mutual membrane-macromolecule configurations of lowest free energy.
A delicate balance between the energetic and entropic contributions to the adsorption free energy on mixed fluid membranes is exhibited in various biological processes 19,20. One important example is the electrostatic-switch mechanism underlying the operation of the myristoylated alanine-rich C kinase substrate (MARCKS), and several other proteins 21,22. MARCKS is a prominent protein C kinase substrate, implicated in a variety of signaling pathways involving, for instance, the control of lipid second messengers and the regulation of cytoskeletal actin. This natively unfolded (and thus flexible) protein binds electrostatically to anionic lipids in the inner leaflet of the plasma membrane. Of special importance in this binding is the multivalent (here z=−4, but generally varying between −3 and −5) anionic lipid (phosphatidylinositol 4,5 bisphosphate, PIP2). The average membrane concentration of PIP2 is typically ∼1%, yet it tends to localize in viral envelopes and membrane rafts, as well as in the binding zones of various proteins involved in signal transduction pathways. Among these proteins is MARCKS, which binds to the plasma membrane through its relatively small (25-residue) but strongly charged effector domain that comprises 13 basic residues. A 150-residue-long flexible polypeptide chain separates the effector domain from the myristoylated N-terminus, and a comparably long and flexible peptide chain connects the effector domain with the C-terminus (for more details, see e.g., Gambhir et al. 9). The myristoyl chain inserts into the hydrophilic core of the lipid bilayer, serving to anchor MARCKS to the membrane.
Experiments reveal that the effector domain sequesters approximately three PIP2 molecules 9,23,24, suggesting that these multivalent (z=−4) lipids provide most of the negative charge required for neutralizing the 13 basic charges of the MARCKS’ effector domain 25. Considering the small average PIP2 concentration in the membrane (∼1%), versus the 10–30% abundance of monovalent acidic lipids (primarily phosphatidyl-serine, PS), it is clear that the protein must import the multivalent lipids from remote membrane regions. Another significant observation is that upon lowering the net charge of the effector domain from +13 to +7, MARCKS detaches from the membrane, thereby exposing PIP2 to cleavage by phospholipase C and initiating a series of signal transduction events 26,27. The change in charge is generally achieved through phosphorylation of three serine residues in the basic domain by protein kinase C.
In addition to electrostatic interactions, MARCKS interacts with the membrane hydrophobically as well, through the myristoyl anchor at the N-terminus and via five phenylalanine side chains within the effector domain 28. A subtle interplay between the electrostatic and hydrophobic interactions, as well as the entropies associated with the long flexible chains on both sides of the effector domain, governs the intricate electrostatic switching process involving MARCKS and PIP2. Motivated by this notion, our goal in this study is to examine the coupling between electrostatic interactions, membrane composition, lipid mobility, and polymer flexibility, and elucidate its role in the adsorption of charged macromolecules onto oppositely charged membranes. Although inspired by the biological relevance of these interactions and processes, our present calculations do not attempt to mimic in detail the behavior of any particular biological system. In fact, our simulations involve a rather short (20-segment) flexible polyelectrolyte chain, interacting with a three-component fluid membrane containing neutral, monovalent, and tetravalent lipids. More specific simulations, modeling MARCKS-membrane interaction, are in progress.
To demonstrate the important role of lipid mobility, our results for the fluid (i.e., annealed) membrane are compared to those obtained for a frozen (i.e., quenched) membrane of the same average lipid composition. The structural and energetic characteristics of a polyelectrolyte interacting with a fluid membrane are qualitatively and quantitatively different from those pertaining to any specific quenched lipid membrane. However, from the purely formal-computational aspect, as argued in the next section, the statistical aspects of polymer adsorption on a fluid membrane can also be derived using a biased superposition of statistical averages corresponding to polymer adsorption on an ensemble of quenched membranes. We shall also show that, in the limit of vanishing polymer concentration (in the bulk solution and hence also on the membrane surface), the average adsorption probability on a Boltzmann-weighted ensemble of quenched membranes equals the adsorption probability on an annealing, fluid membrane. Differences between the two kinds of membranes appear at nonzero concentrations of macromolecules, and will be accounted for using a simple cell model.
In addition to the fluid and quenched lipid membranes, some of our calculations describe polyelectrolyte adsorption on uniformly charged surfaces, such as those of metal oxides or metal electrodes 29,30. We shall see that, in general, uniformly charged surfaces adsorb more weakly than either a quenched or a fluid lipid membrane. Another limiting case of interest is that of a stiff polyelectrolyte interacting with a charged membrane. Contrasting the adsorption behavior of such a molecule with that of a self-avoiding freely jointed chain will emphasize the role of polyelectrolyte flexibility. As we shall see, weak electrostatic interaction may not suffice to overcome the loss of chain flexibility upon adsorption. This will be demonstrated by considering the limiting cases of a weakly charged polyelectrolyte and a membrane containing only monovalent lipids.
In the next section we first describe the basic statistical-thermodynamic background underlying the adsorption of a flexible macromolecule onto the various types of lipid membranes mentioned above. Later in that section we introduce the model system and describe our extended version of the Rosenbluth Monte Carlo (MC) simulation scheme 31,32, which enables efficient simulation of polymer adsorption on annealing, fluid, membranes. We then present and analyze the results of the simulations and close with a brief summary of the main conclusions.
The lipid molecules comprising a fluid (annealed) membrane are mobile and can thus diffuse into and out of the interaction region with the adsorbing macromolecule. The statistical thermodynamics of macromolecule adsorption onto such membranes involves simultaneous averaging over all the polymer and lipid degrees of freedom. In a frozen membrane, on the other hand, the lipids are immobile, and thermodynamic averaging is obtained by tracing over many quenched arrangements of the membrane lipids. Qualitative differences between these two membrane types are reflected in various structural and thermodynamic properties, such as the local lipid composition around the adsorbed macromolecule, and the average adsorption free energies. Formally, however, it can be shown that the various characteristics of macromolecule-lipid interaction on a fluid membrane, can be derived by a biased averaging of the same property over a Boltzmann-weighted ensemble of quenched membranes. (This suggests an indirect way of simulating adsorption onto the fluid membrane using quenched membrane simulations.) We shall see that the differences between average macromolecule adsorption characteristics on fluid and quenched membranes are substantial at nonzero concentrations of macromolecules, but disappear in the limit of vanishing concentration. We open this section with a detailed discussion of the relevant thermodynamic background, and then describe the model system and our method of simulation.
All our simulations involve one adsorbed macromolecule, as formally appropriate to low surface concentrations. However, using a simple cell model, the simulations can also be used to account for the adsorption behavior at higher surface concentrations. In this model, the membrane area A is divided into an array of A/a noninteracting cells, all with the same area, a, and the same lipid composition. The cell area is large enough to comfortably accommodate one adsorbed macromolecule. The model thus, approximately, accounts for excluded area effects but ignores other intermacromolecule interactions. Note that to contain a flexible macromolecule, e.g., an unfolded protein, a membrane cell should typically consist of several hundreds of lipid molecules. For the three-component membranes of interest here, this implies an enormous number of two-dimensional lipid arrangements.
We treat the quenched membrane as an ensemble of independent cells, each characterized by a specific frozen lipid configuration, m. To compare the quenched membrane with a fluid membrane having the same lipid composition, we equate the fraction, Pq(m), of quenched membranes in configuration m, with the Boltzmann weight, P(m), of m in the bare fluid membrane, i.e.,
![]() | (1) |
is the partition function per cell of a bare fluid membrane,![]() | (2) |
The partition function, per cell, of a fluid membrane occupied by an adsorbed macromolecule is given by
![]() | (3) |
Here U(m,p)=U(m)+U(p|m) is the potential energy corresponding to the membrane-polymer configuration m,p. The term U(p|m) stands for the energy of a polymer in state p, interacting with a membrane in a given configuration m. It includes the self-energy of the polymer (i.e., the sum of its intersegment potentials), and its interaction energy with a membrane in state m. By p=α,r we refer to the polymer chain conformation, α, and the position, r, of the polymer relative to the membrane plane (see below). The sum over p thus (tacitly) involves integration over r, implying that
is a configurational partition function, bearing the dimensions of volume. As above, U(m) is the interlipid interaction energy.
The sum of Boltzmann factors,
![]() | (4) |
i.e., the energy of the occupied m-cell is measured relative to the ground state energy U(m). On this energy scale we obtain
for the empty cell. The quantity
introduced in the last equality of Eq. (3) may be interpreted as the average partition function, per cell, in a Boltzmann-weighted ensemble of quenched membranes.From Eq. (3), it also follows that the average of any property A (e.g., adsorption energy, polymer radius of gyration, etc.) of a polymer adsorbed on a fluid membrane, can, in principle, be evaluated as a biased average of A in the ensemble of quenched membranes. Explicitly, let
![]() | (5) |
![]() | (6) |
the weight, (
), of the quenched membrane configuration m, is the product of the fraction of such membranes (P(m)) with the statistical weight (∝
) of all polymer conformations on this membrane. It will be shown below that this formal relationship between fluid and quenched membrane averaging may be given a physical meaning in the limit of vanishing macromolecule concentration. In this limit the probability of finding a macromolecule adsorbed onto a quenched membrane m is proportional to
and hence
is the average of A(m,p) in the ensemble of quenched membranes. Note, however, that this average differs from the simple average 
To account for the adsorption behavior at nonzero surface concentrations we should consider a many-cell membrane in equilibrium with a solution of macromolecules. Suppose the bulk solution is of volume V, and contains Nb macromolecules of chemical potential μ. For simplicity we assume dilute solution behavior, in which case μ=−ln qb+ln φb, where φb=Nb/V is the bulk density of macromolecules, and
![]() | (7) |
The cells comprising a fluid membrane are identical. Treating the membrane as an open system with respect to macromolecule exchange, the grand-canonical partition function of the membrane is
![]() | (8) |
is the two-state (i.e., empty and occupied) partition function of a membrane cell. Using Ns to denote the number of macromolecules adsorbed on the membrane surface, the fraction, θf=Ns/M=Nsa/A, of the membrane area occupied by macromolecules (or, the surface coverage), is given by
We thus obtain a Langmuir-like adsorption equation![]() | (9) |
where v is a volume per macromolecule defined in more detail below. Thus,
may be regarded as the volume fraction of polymers in solution.The grand partition function of the quenched membrane is given by
![]() | (10) |
so that![]() | (11) |
![]() | (12) |
![]() | (13) |
Both θm/(1−θm) and exp(−ΔFm) are convex functions of their arguments. Using Jensen’s inequality of convex functions 34, it thus follows from Eqs. (12) that
![]() | (14) |
From Eq. (12) it follows that the equality 〈θm〉q=θf is obtained only in the limit of vanishing surface coverage, i.e., when
(〈ΔFm〉q → ΔFf requires that all ΔFm are negligibly small, as can be seen from Eq. (13).) Underlying the limiting behavior 〈θm〉q=θf is the fact that although the lipid molecules of a quenched membrane are, indeed, immobile, the adsorbed macromolecules can nevertheless explore all membrane states, m. This may be achieved via lateral diffusion on the membrane surface and/or by desorption from one local membrane region and adsorption into another. Note also that in the limit
we obtain
(see Eq. (11)). Thus, in a Boltzmann-weighted ensemble of quenched membranes, the number of macromolecules bound to membranes of lipid configuration m, is proportional to P(m)θm and hence to
If measured over an ensemble of quenched membranes, the average of a physical observable A(m,p) would then be given by
which, as noted in Eq. (6), is equal to 〈A〉f. Thus, in the θ → 0 limit, both the average surface concentration, and the θ-weighted average of A(m,p) over the quenched membrane ensemble, approach the corresponding quantities in the fluid membrane.
A physical interpretation of the last conclusion can be given in terms of the cell model, as follows. When θ → 0 (and in the absence of kinetic constraints), each of the adsorbed macromolecule can freely and independently visit all membrane environments m, thereby sampling all possible lipid-polymer configurations m,p. Furthermore, since none of the cells is blocked, all m and hence all m,p are sampled according to their Boltzmann weights, just like the states sampled by a macromolecule on a fluid membrane. The difference is, of course, that a macromolecule adsorbed on a fluid membrane need not migrate from one cell to another to sample the entire configuration space. Consequently, θ=θf is the same for all cells of the fluid membrane, whereas a wide distribution of θm values characterizes the ensemble of quenched membranes. Note however, that in practice, even this formal equivalence between a fluid membrane and an ensemble of quenched membranes may not materialize, even when θ → 0, because of kinetic barriers to macromolecule mobility.
Interestingly, it is possible that for some quenched membrane states m, ΔFm≤ΔFf (and hence, θm≥θf). This may appear surprising in view of the fact that the free energy, Ff (not ΔFf), of a macromolecule interacting with a fluid membrane is invariably lower than the free energy, Fm, of a macromolecule adsorbed onto any quenched membrane, m. (This is because
involves summation over both m and p and is therefore larger than any
) Hence,
After adsorption, however, the distribution of lipid arrangements in the fluid membrane is no longer the Boltzmann distribution before adsorption, implying a loss of lipid-mixing entropy. Of course, no loss of lipid entropy is involved upon adsorption onto a quenched membrane, which explains why certain quenched states can be more attractive to macromolecule adsorption than the fluid membrane.
Upon increasing the concentration of macromolecules in solution, the more strongly adsorbing m-values of the quenched membrane will be occupied first. Once these favorable local environments (or cells) are populated, further adsorption is necessarily suppressed. This implies 〈θm〉q≤θf, because in the fluid membrane every cell can independently anneal its lipid distribution, thereby enhancing adsorption.
Our conclusions regarding the relationship between macromolecule adsorption on quenched versus fluid membranes agree with previous works pertaining to polymer statistics in random media. Cates and Ball 35 have studied the behavior of a single long polymer chain in a random medium and concluded that, as long as the environment is infinite, the quenched and annealed averaging will yield the same statistical chain properties. Our fluid and frozen membranes are analogous to the annealed and quenched random potentials in the treatment above. Inequalities valid for the multichain adsorption, analogous to Eq. (12), have been obtained by Andelman and Joanny 36,37 for neutral chains adsorbing on annealed and quenched flat surfaces. The main conclusion there is that the density of polymers on an annealed surface (membrane) is always higher than in the frozen case.
The MC simulations presented in the following sections enable evaluation of all the partition functions encountered above, as well as a variety of relevant structural properties. First, however, we have to clarify what distinguishes an adsorbed macromolecule from a free macromolecule in solution.
We noted above that the sum over p in Eqs. (2) and (3) involves all possible chain conformations, α, as well as all possible positions, r, of the macromolecule, relative to some arbitrary point on the membrane. Identifying the membrane surface with the (x,y) plane, the only relevant coordinate is the distance, z, of the macromolecule from the membrane surface. This distance can be expressed in terms of the normal displacement of any chain segment (or the center of mass) from the membrane plane. In the simulations, we find it convenient to measure this distance in terms of z1, the normal displacement of the first (more precisely, terminal) chain segment (see Fig. 1). Beyond a certain distance from the membrane surface, comparable to the range of membrane-macromolecule interactions, the macromolecule is not affected by the membrane. This cutoff distance, λ, may be defined in terms of z1, such that for z1≤λ, the macromolecule is considered adsorbed, and otherwise as free in solution. An alternative, yet practically equivalent definition of λ can be given in terms of the average segment density profile (see below).
With p≡z1, α, we find from Eq. (3) that the partition function of a macromolecule adsorbed on a fluid membrane is given by
![]() | (15) |
![]() | (16) |
represents the (average) partition function per unit volume of an adsorbed macromolecule, or, equivalently, the average partition function per unit length along the membrane normal. Note that for large z1 (practically for z1≥λ) we have U(m,α;z1)=U(m)+U(α), and hence 
We may now rewrite Eq. (9) in the form
![]() | (17) |
![]() | (18) |

For chain molecules composed of L segments, the segment density in the bulk solution is ρ(z=∞)=ρb=Lφb. Near the membrane surface the segment density, ρ(z), is different from ρb and is given by
![]() | (19) |
The surface excess of adsorbed macromolecule is, by definition,
![]() | (20) |
is the average (three-dimensional) density of chain segments within the surface layer. Note that the upper limit in the integral defining Γ can be replaced by λ (or any larger value). The second equality may also be regarded as the definition of the surface layer thickness λ. Note that in the limit of vanishing surface density, λaρs/L=θ.A slightly different definition of λ, useful in our numerical calculations, can be given in terms of the ratio q(z1)/q(∞)=φ(z1)/φb. Namely, we can choose λ as the smallest value of z1, beyond which this ratio is practically 1. In practice, the two definitions are indistinguishable, because
so that the integral over ρ(z) in Eq. (20) can be replaced by the integral over φ(z). The equality of the two integrals follows from the fact that for practically all z1 within λ, all chain segments will be found inside the surface layer. For chains originating near λ, say at z1=λ−δ (δ ≪ λ), some conformations will cross the z=λ surface, contributing less than L segments to the surface layer density. By symmetry, however, chains originating at z1=λ+δ will compensate for the loss of segments from the z1=λ−δ chains. The near-equivalence of chains originating at z1=λ±δ follows from the fact that these chains are hardly affected by the membrane.
Equations (19) are applicable to the fluid membrane, as well as any quenched membrane state m. For small values of θ we can use φ(z1)=φbq(z1)/q(∞) with q(z1)/q(∞) derived from our single-chain simulations. Approximate density profiles for nonzero surface concentrations can be derived by expressing ρ(z) as the product of the probability (θ) to find the cell occupied and the normalized density profile corresponding to one adsorbed molecule. With the aid of Eqs. (17) we then find that for z≤λ,
![]() | (21) |
Our model system consists of a single polyelectrolyte interacting with a finite size membrane, large compared to the size of the polymer and the range of intermolecular potentials. With the exception of several limiting test cases, in all simulations we consider polyelectrolyte chains composed of L=20 spherical segments of diameter d, interacting with a considerably larger two-dimensional membrane cell consisting of M=2500 lipid headgroups. The lipid membrane is modeled as a perfectly flat and impenetrable two-dimensional hexagonal lattice, with lipid headgroups occupying all its lattice sites. The lattice constant is set equal to d. The membrane may thus be regarded as an hexagonal array of closely packed disks of diameter d, as illustrated schematically in Fig. 1. Using a typical lipid headgroup area of 65Å2 we find d=8.66Å. We simulate three-component membranes, composed of electrically neutral (z=0), monovalent (z=−1), and tetravalent (z=−4) headgroups. These may be regarded as representing, respectively, the phosphatidyl-choline (PC), phosphatidyl-serine (PS), and phosphatidylinositol 4,5 bisphosphate (PIP2) lipids mentioned in the previous section. The lipid charges are treated as point charges residing at the grid points of the hexagonal lattice, and the electrostatic repulsion between them is modeled in the Debye-Hückel (DH) approximation. Explicitly, the interaction potential between lipids of valences z1 and z2 at distance r apart, and in units of kBT, is
![]() | (22) |
In the simulations each polymer bead carries a unit positive charge (z=+1), localized at its center. Although the polymer bond length d is fixed, there are no other restrictions on bond angles, except for those implied by electrostatic and spatial (excluded volume) repulsion between nonbonded segments. For the electrostatic interaction between polymer charges we again use DH potentials. The spatial repulsion is modeled using the shifted and truncated Lennard-Jones potential:
![]() | (23) |
ensures the onset of steep repulsion as soon as r falls below d38.The electrostatic attraction between the oppositely charged polymer and membrane is also modeled using screened DH potentials. In addition, the membrane surface is treated as an impenetrable wall to the polymer, implying a minimal distance of d/2 between polymer and lipid charges. At this distance the electrostatic attraction between a polymer (z=+1) segment and a monovalent (z=−1) lipid headgroup is 1.07 kBT. For comparison, the electrostatic repulsion between neighboring monovalent lipids or adjacent polymer beads, taking the distance of closest approach to be r=d, is 0.35 kBT.
Since the distances between charges in the system are either comparable to or larger than the Debye length, i.e., r≥d ≈ κ−1, screening by counterions is expected to be effective. Under physiological conditions, when κ−1 is small (of the order of few Ångströms), the long-range character of the electrostatic interactions is screened and DH potentials offer a reasonable approximation. These potentials are commonly employed in simulation and theoretical studies of polyelectrolyte-surface interactions (see, e.g., 39,40).
Henceforth, we shall measure all distances in units of d. Recall also that energies are measured in units of kBT.
The Rosenbluth MC method 31, or its configurational-bias variant, provides an efficient means for simulating polymer statistics 32. In this approach, chain conformations are generated, segment after segment, with preference for conformations of large statistical weight. Based on these ideas we present below our extension of the Rosenbluth scheme for modeling polyelectrolyte adsorption on fluid, as well as frozen and uniform membranes.
Consider first a polymer interacting with a membrane of quenched lipid configuration m. The simulation begins by placing the first chain segment at distance z1 above the center of the membrane cell, where its interaction energy with membrane lipids is u(z1,m) (see Fig. 1). We then sample k random directions (and hence positions, r2) for segment 2 and select one, say
with probability
where
is the interaction energy of segment 2 with segment 1 and the membrane, and
is a local partition function. This procedure is continued until all segments of the chain are generated. Repeated applications of this scheme (for the given (m, z1)) yield an ensemble of conformations {α=r2,…,rL;z1,m} with probabilities
![]() | (24) |
is the total interaction energy of polymer segments with each other and with the membrane. The partition function,![]() | (25) |
Since every possible conformation α is sampled with probability proportional to exp[−U(α)]/W(α), proper Boltzmann averaging requires weighting each α by its Rosenbluth factor W(α); i.e., the average (over α, for the given z1,m) of any structural or energetic polymer property A is given by
![]() | (26) |
![]() | (27) |
For z1>λ we have qm(z1)=qm(∞)≡qb. Averaging 〈A(z1,m)〉 over all z1≤λ we obtain the average of A (over all conformations) for molecules adsorbed on a frozen membrane of lipid configuration m,
![]() | (28) |
![]() | (29) |
From Eqs. (3) we know that the thermodynamic and structural properties of a fluid membrane can be modeled based on simulating an ensemble of quenched membranes. However, this procedure is rather indirect and often impractical. Alternatively, adsorption on the fluid membrane could be simulated by combining the Rosenbluth and Metropolis methods. That is, after generating a polymer in conformation p=(α;z1) for a given lipid configuration m, the membrane is allowed to relax to a new configuration m′ through a series of Metropolis moves. Another polymer conformation p′ can then be generated for m′, letting the membrane relax to m″, and so on. The problem here is that the relaxed membrane is no longer the one which served to generate the last polymer conformation. A retracing procedure 32 can be used to improve this scheme, but not fully eliminate its inconsistencies. We have adopted, therefore, an alternative simulation method for the fluid membrane whereby, in the spirit of the Rosenbluth sampling scheme, we generate simultaneously both polymer conformations p and membrane configurations m, as follows.
Any joint polymer-membrane configuration p,m is fully specified by the coordinates of K=L+M(−1)+M(−4) particles; that is, L polymer segments, M(−1) monovalent lipids, and M(−4) tetravalent lipids (M(0)=M−M(−1)−M(−4) neutral lipids occupy all other membrane sites). We now generate a joint (p,m) configuration by randomly adding either a polymer segment or a charged lipid, until all particles have been placed. More explicitly, suppose the new configuration is already partly grown, consisting of a polymer chain of length l, and a partially charged membrane containing m(−1) and m(−4) anionic lipids. One of the remaining (K−l−m(−1)−m(−4)) particles is now randomly selected and added to the system. If this is a polymer segment it is added as the (l+1)th segment of the chain. As before, this segment is placed in one of k possible positions, with probability exp[−u(l+1;l,m(−1),m(−4))]wl+1, and u(l+1;l,m(−1),m(−4)) is the interaction potential of the added particle with all those already placed, and wl+1 is defined as usual. If the new particle is, say, a monovalent lipid, it is placed with probability
in one of n randomly chosen membrane sites, where u(m(−1)+1;l,m(−1),m(−4)) is the interaction energy of this lipid with the rest of the system, and
is the sum of the Boltzmann factors corresponding to the n membrane sites. (In the simulations we usually sample n=1000 sites, some of which are possibly occupied already and thus do not contribute to w.) This procedure is repeated until all chain segments and all charged lipids are placed, resulting in a statistical distribution of p,m configurations, whose probabilities are
![]() | (30) |
![]() | (31) |
As for the quenched membrane, we generally sample many polymer-membrane configurations corresponding to various z1 values and only then average over this variable. The averaging procedure is analogous, e.g., the average of A for a given z1 is
![]() | (32) |
![]() | (33) |
![]() | (34) |
From the simulations we have derived the basic thermodynamic characteristics of macromolecules interacting with fluid, quenched, and uniformly charged membranes. In parallel, for every system considered we have calculated a variety of structural properties, such as the two-dimensional distribution of charged lipids in the membrane plane, or the density profile of chain segments along, as well as perpendicular to, the membrane normal. Two membrane compositions were analyzed in detail:
Note that the average charge, per lipid, corresponding to these membranes (hereafter also referred to as the weakly-charged and the strongly-charged membranes) is
and
respectively.
For both compositions, simulations were performed for fluid, quenched, and uniformly charged membranes. In the uniformly charged membrane all lipids carry the same partial charge, 
We repeat the numerical values of the various parameters in our model: The polyelectrolyte is a freely jointed homopolymer chain composed of L=20 spherical segments, each carrying a z=+1 charge (see Fig. 1). The membrane cell is an hexagonal array of M=50×50 lipid molecules. The headgroup diameter,
is equal to the polymer’s bond length. The distance d also marks the onset of steep excluded volume repulsion between nonbonded chain segments (see Eq. (23) where σ=d/21/6≈7.72Å and
). The Bjerrum and Debye lengths are lB=7.14Å and κ−1=10Å, respectively.
For the sake of comparison, we have also performed a limited number of simulations for a stiff (rodlike) polymer, as well as for a weakly charged
polymer. Recall that simulations are performed for varying values of the first segment position, z1, and that
corresponds to a free polymer in solution. For the three-dimensional case of a polymer in solution we have also carried out, for comparative reasons, one set of simulations for an electrically neutral polymer.
The number of chain-membrane conformations generated for each z1 value of a polymer adsorbed on a fluid membrane is ∼106. The number of chain conformations generated for each z1 value of a given quenched membrane m is ∼103, and the number of membrane configurations is 104. The increments in chain origin positions are Δz1=1. (Recall that distances are measured in units of d.) The number of possible bond directions when generating polymer conformations is k=50. The number of possible positions for lipid addition in our simulation scheme of the fluid membrane is n=1000.
A pictorial illustration of the polymer-membrane configurations generated by our simulations is given in Fig. 2. The figure shows top and side views of two (rather arbitrary) simulation snapshots of a polyelectrolyte interacting with a fluid membrane of composition PC:PS:PIP2=98:1:1. Only part of the membrane is shown, yet it is apparent that the local concentration of charged lipids in the vicinity of the polymer significantly exceeds the membrane average.
Fig. 3 shows how ΔF(z1), the differential adsorption free energy, and ΔE(z1), the differential adsorption energy, vary with the distance (z1) of the chain origin from the surface of the weakly charged (PC:PS:PIP2=98:1:1) membrane. Fig. 4 shows the same quantities for the strongly charged (PC:PS:PIP2=89:10:1) membrane. The value ΔF(z1) is the free energy change, or, the potential of mean force, associated with bringing the first segment of the macromolecule from the bulk solution to distance z1 from the membrane. Then ΔE(z1) and TΔS(z1)=ΔE(z1)−ΔF(z1) are the energetic and entropic components of this free energy difference. More explicitly, for the fluid and uniformly charged membranes ΔEf(z1)=〈U(α,m;z1)〉f−〈U(α)〉b−〈U(m)〉f and ΔEu(z1)=〈U(α,z1)〉−〈U(α)〉b, respectively. For the quenched membrane we show here the average energy change corresponding to the Boltzmann-weighted ensemble of quenched membranes,
![]() | (35) |
with a similar definition of ΔFu(z1). The corresponding free energy change for the quenched membrane is defined here as
It should be noted that the net (or integral) adsorption energy of the quenched membrane is not a simple integral of
Similarly, the net free energy change of all membranes is not a direct integral of ΔF(z1). These issues will be clarified after analyzing Figure 3 and Figure 4.Figure 3 and Figure 4 reveal, as expected, that the interaction (potential of mean force) between the polyelectrolyte and all three types of oppositely charged membranes is attractive, i.e., ΔE(z1)<0. In Fig. 3 we also show the results for a weakly charged polymer (z=+1/2) interacting with a fluid membrane. This figure reveals that although in all cases ΔE(z1)<0, this attractive interaction may not suffice to ensure adsorption. More explicitly, we note that in the case of a z=+1 polymer interacting with the uniformly (
) charged membrane, as well as in the case of a weakly charged (z=+1/2) polymer interacting with a fluid membrane, the free energy change is positive, ΔF(z1)>0. This is because the electrostatic attraction cannot counterbalance the repulsive depletion interaction resulting from the loss of conformational entropy experienced by any flexible molecule near a rigid wall. The weakly charged quenched membrane is, on average, nonadsorbing as well. Only the fluid membrane appears attractive to the peripheral macromolecule, owing to its ability to recruit charged lipids into the interaction zone. However, even this membrane is repulsive when the polymer charge is reduced to z=+1/2. Fig. 4 reveals that, upon increasing the membrane charge (to
per lipid), all membranes become attractive. The strongest binding is to the fluid membrane and the weakest corresponds to the uniformly charged one.
The entropy losses, TΔS(z1)=ΔE(z1)−ΔF(z1), associated with polyelectrolyte adsorption are quite substantial. In the quenched and uniform membrane cases these entropy losses reflect the lower conformational entropy of the adsorbed molecule, compared to that of a polymer in solution. The entropy loss is even higher, reaching ∼70% in the case of the fluid membrane, see Figure 3 and Figure 4. The origin of the enhanced entropy deficit experienced by this membrane is the additional loss of lipid mixing entropy.
From Eq. (13) we know that
An analogous equality is also valid for the differential partition functions, q(z1). That is,
![]() | (36) |
explaining why
in Figure 3 and Figure 4.In Fig. 5 we show, for our three model membranes, how the partition function (equivalently, the first segment density) ratio defined in Eq. (36) varies with z1. It should be emphasized that the partition functions corresponding to the quenched and fluid membranes have been obtained using the two different MC simulation schemes described in the previous section. Apart from the small numerical noise, we indeed find that the partition functions corresponding to the fluid and the ensemble of quenched membranes are essentially identical, reassuring that the different simulation methods indeed yield identical results.
for a macromolecule interacting with weakly charged membranes of composition PC:PS:PIP2=98:1:1 (a), and strongly charged membranes where PC:PS:PIP2=89:10:1 (b). Solid, dashed, and dotted curves correspond to the fluid, quenched, and uniformly charged membranes, respectively.The ratio q(z1)/qb=φ(z1)/φb reveals, as expected, the stronger attraction of the polyelectrolyte to the strongly charged membrane (Figure 5b). Similar behavior is shown by the average segment density profiles, ρ(z), as defined in Eq. (19) and shown in Fig. 6. Again we see that for θ → 0, the density profiles corresponding to the fluid and quenched membranes are the same.
(and hence θ → 0), and the lower one is for 
Figure 5 and Figure 6 convey similar information. Fig. 5 displays the density profile of chain termini, whereas Fig. 6 shows the average density due to all chain segments (see Eq. (19)). Indeed, apart from small differences at the very small (i.e., near the membrane) and large (
) values of z1, the two profiles are quite similar. Unlike φ(z1) (∝ q(z1)), which decreases monotonically with z1, the maximum in ρ(z) occurs slightly away from the membrane surface. This is probably due to the fact that terminal segments can more easily attach and detach from the surface. Comparing Figure 5 and Figure 6 we also note that φ(z1) (∝ q(z1)) decays slightly more slowly than ρ(z), reflecting the fact that although one chain-end may reside relatively far from the membrane, other segments are attracted to the membrane (see Fig. 2). Of course, all end effects become negligible for very long chains and we then expect similar density distributions for all chain segments.
In Fig. 6, we also show representative density profiles corresponding to high surface concentrations of macromolecules (see Eq. (21)). Under these conditions we expect different adsorption probabilities on the fluid and quenched membrane. Indeed, for a macromolecule bulk density of
Eq. (17) yields θf=0.5 for the strongly charged fluid membrane, whereas Eq. (18) implies a much smaller surface density for the quenched membrane,
Additional values are given in Table 1. Note that the average free energy of adsorption in the ensemble of quenched membranes is zero, indicating that some membrane environments must be repulsive (see below). Also repulsive is the weakly charged uniform membrane, as clearly seen in Figure 3b. Indeed, the ratio
), which may be interpreted as the three-dimensional density of macromolecules very near the membrane, is 
| Table 1 Adsorption properties |
PC:PS:PIP2=89:10:1 ( = −0.14) | PC:PS:PIP2=98:1:1 ( ) | |||||||
|---|---|---|---|---|---|---|---|---|
| Fluid | Quenched | Uniform | Fluid | Quenched | Uniform | |||
| ΔE | −12.5 | −7.4 | −3.1 | −5.0 | −2.4 | −0.7 | ||
| ΔF | −3.5 | −1.4 | −1.0 | −0.7 | 0 | 1.3 | ||
| θ | 0.5 | 0.28 | 0.09 | 0.06 | 0.05 | 0.01 | ||
Adsorption energies and free energies for the ) and membranes. For the quenched lipid membrane we list 〈ΔEm〉q and 〈ΔFm〉q. The surface concentrations, θ, in the bottom row (〈θm〉q for the quenched membrane) are for a bulk concentration of macromolecules ![]() |
The adsorption free energy and related thermodynamic functions are calculated using the partition functions appearing in Eqs. (15), whose values depend on the cutoff distance λ. We have determined
as the distance beyond which ρ(z)/ρb≤1.1 for attractive membranes (ΔF<0), or >0.9 for repulsive ones. (This criterion closely satisfies the second equality in Eq. (20).)
Given the λ-values we have calculated, the integral adsorption energies, free energies, and surface concentrations θ for the fluid, quenched, and uniform membranes. Fig. 7 shows the distributions, P(ΔF) and P(ΔE), of adsorption free energies, ΔFm, and energies, ΔEm, for the ensemble of quenched membranes. The adsorption energies are defined here by
and their average is
The integral adsorption energies for the fluid and uniformly charged membranes are
where ΔE(z1) are the differential adsorption energies shown in Figure 3 and Figure 4. The adsorption free energies are given by
with
Their average is
In Figure 3 and Figure 4 we have shown
Thus, as noted above,
is not simply the integral of
The relationship between ΔF(z1) and ΔF of the fluid membrane is different, namely,
with a similar relationship for the uniform membrane.
The numerical values of ΔF and ΔE for the fluid, quenched, and uniform membranes are listed in Table 1. Again we note that the adsorption energy is largest for the fluid membrane and smallest for the uniform membrane. In fact, no adsorption takes place on the weakly charged uniform membrane, (ΔFu>0), and the weakly charged quenched membrane is, on average, nonadsorbing as well. Note, however, that in this case P(ΔF) is bimodal. From Figure 7a, we note that although ΔEm<0 for all m values, the bimodal distribution of ΔFm reflects two distinct classes of quenched environments, corresponding to attractive (ΔFm<0) and repulsive (ΔFm>0) membranes.
The distributions, P(θ), of surface concentrations for the ensembles of weakly and strongly charged quenched membrane are shown in Fig. 8. Also mentioned there (and in Table 1) are the average values of θ for the fluid and uniform membranes, confirming that adsorption onto the fluid membrane is, indeed, the strongest of all. In accordance with the results in Figure 7a, we note in Figure 8a that a large fraction of the local environments comprising a weakly charged quenched membrane are repulsive. On the other hand, a weakly charged fluid membrane is everywhere attractive. This is of course due to the ability of its charged lipids to diffuse and localize at the macromolecule adsorption site. In biological systems, where the interactions are often weak, such subtle differences could be of crucial importance.

The structural and thermodynamic properties of the adsorbed macromolecules are intimately related to each other. For instance, the density profile of chain termini, φ(z1), enters the calculation of partition functions and free energies. In this subsection we present additional information, pertaining to the configurational stat