| 25 |
|
\bibliographystyle{achemso} |
| 26 |
|
|
| 27 |
|
\title{Dipolar ordering in the ripple phases of molecular-scale models |
| 28 |
< |
of lipid membranes} |
| 28 |
> |
of lipid membranes} |
| 29 |
|
\author{Xiuquan Sun and J. Daniel Gezelter \\ |
| 30 |
|
Department of Chemistry and Biochemistry,\\ |
| 31 |
|
University of Notre Dame, \\ |
| 75 |
|
experimental results provide strong support for a 2-dimensional |
| 76 |
|
hexagonal packing lattice of the lipid molecules within the ripple |
| 77 |
|
phase. This is a notable change from the observed lipid packing |
| 78 |
< |
within the gel phase.~\cite{Cevc87} The X-ray diffraction work by |
| 78 |
> |
within the gel phase,~\cite{Cevc87} although Tenchov {\it et al.} have |
| 79 |
> |
recently observed near-hexagonal packing in some phosphatidylcholine |
| 80 |
> |
(PC) gel phases.\cite{Tenchov2001} The X-ray diffraction work by |
| 81 |
|
Katsaras {\it et al.} showed that a rich phase diagram exhibiting both |
| 82 |
|
{\it asymmetric} and {\it symmetric} ripples is possible for lecithin |
| 83 |
|
bilayers.\cite{Katsaras00} |
| 85 |
|
A number of theoretical models have been presented to explain the |
| 86 |
|
formation of the ripple phase. Marder {\it et al.} used a |
| 87 |
|
curvature-dependent Landau-de~Gennes free-energy functional to predict |
| 88 |
< |
a rippled phase.~\cite{Marder84} This model and other related continuum |
| 89 |
< |
models predict higher fluidity in convex regions and that concave |
| 90 |
< |
portions of the membrane correspond to more solid-like regions. |
| 91 |
< |
Carlson and Sethna used a packing-competition model (in which head |
| 92 |
< |
groups and chains have competing packing energetics) to predict the |
| 93 |
< |
formation of a ripple-like phase. Their model predicted that the |
| 94 |
< |
high-curvature portions have lower-chain packing and correspond to |
| 95 |
< |
more fluid-like regions. Goldstein and Leibler used a mean-field |
| 96 |
< |
approach with a planar model for {\em inter-lamellar} interactions to |
| 97 |
< |
predict rippling in multilamellar phases.~\cite{Goldstein88} McCullough |
| 98 |
< |
and Scott proposed that the {\em anisotropy of the nearest-neighbor |
| 99 |
< |
interactions} coupled to hydrophobic constraining forces which |
| 100 |
< |
restrict height differences between nearest neighbors is the origin of |
| 101 |
< |
the ripple phase.~\cite{McCullough90} Lubensky and MacKintosh |
| 102 |
< |
introduced a Landau theory for tilt order and curvature of a single |
| 103 |
< |
membrane and concluded that {\em coupling of molecular tilt to membrane |
| 104 |
< |
curvature} is responsible for the production of |
| 105 |
< |
ripples.~\cite{Lubensky93} Misbah, Duplat and Houchmandzadeh proposed |
| 106 |
< |
that {\em inter-layer dipolar interactions} can lead to ripple |
| 107 |
< |
instabilities.~\cite{Misbah98} Heimburg presented a {\em coexistence |
| 108 |
< |
model} for ripple formation in which he postulates that fluid-phase |
| 109 |
< |
line defects cause sharp curvature between relatively flat gel-phase |
| 110 |
< |
regions.~\cite{Heimburg00} Kubica has suggested that a lattice model of |
| 111 |
< |
polar head groups could be valuable in trying to understand bilayer |
| 112 |
< |
phase formation.~\cite{Kubica02} Bannerjee used Monte Carlo simulations |
| 113 |
< |
of lamellar stacks of hexagonal lattices to show that large head groups |
| 88 |
> |
a rippled phase.~\cite{Marder84} This model and other related |
| 89 |
> |
continuum models predict higher fluidity in convex regions and that |
| 90 |
> |
concave portions of the membrane correspond to more solid-like |
| 91 |
> |
regions. Carlson and Sethna used a packing-competition model (in |
| 92 |
> |
which head groups and chains have competing packing energetics) to |
| 93 |
> |
predict the formation of a ripple-like phase. Their model predicted |
| 94 |
> |
that the high-curvature portions have lower-chain packing and |
| 95 |
> |
correspond to more fluid-like regions. Goldstein and Leibler used a |
| 96 |
> |
mean-field approach with a planar model for {\em inter-lamellar} |
| 97 |
> |
interactions to predict rippling in multilamellar |
| 98 |
> |
phases.~\cite{Goldstein88} McCullough and Scott proposed that the {\em |
| 99 |
> |
anisotropy of the nearest-neighbor interactions} coupled to |
| 100 |
> |
hydrophobic constraining forces which restrict height differences |
| 101 |
> |
between nearest neighbors is the origin of the ripple |
| 102 |
> |
phase.~\cite{McCullough90} Lubensky and MacKintosh introduced a Landau |
| 103 |
> |
theory for tilt order and curvature of a single membrane and concluded |
| 104 |
> |
that {\em coupling of molecular tilt to membrane curvature} is |
| 105 |
> |
responsible for the production of ripples.~\cite{Lubensky93} Misbah, |
| 106 |
> |
Duplat and Houchmandzadeh proposed that {\em inter-layer dipolar |
| 107 |
> |
interactions} can lead to ripple instabilities.~\cite{Misbah98} |
| 108 |
> |
Heimburg presented a {\em coexistence model} for ripple formation in |
| 109 |
> |
which he postulates that fluid-phase line defects cause sharp |
| 110 |
> |
curvature between relatively flat gel-phase regions.~\cite{Heimburg00} |
| 111 |
> |
Kubica has suggested that a lattice model of polar head groups could |
| 112 |
> |
be valuable in trying to understand bilayer phase |
| 113 |
> |
formation.~\cite{Kubica02} Bannerjee used Monte Carlo simulations of |
| 114 |
> |
lamellar stacks of hexagonal lattices to show that large head groups |
| 115 |
|
and molecular tilt with respect to the membrane normal vector can |
| 116 |
< |
cause bulk rippling.~\cite{Bannerjee02} |
| 116 |
> |
cause bulk rippling.~\cite{Bannerjee02} Recently, Kranenburg and Smit |
| 117 |
> |
described the formation of symmetric ripple-like structures using a |
| 118 |
> |
coarse grained solvent-head-tail bead model.\cite{Kranenburg2005} |
| 119 |
> |
Their lipids consisted of a short chain of head beads tied to the two |
| 120 |
> |
longer ``chains''. |
| 121 |
|
|
| 122 |
|
In contrast, few large-scale molecular modeling studies have been |
| 123 |
|
done due to the large size of the resulting structures and the time |
| 349 |
|
along dipole $i$ and $\mathbf{\hat r}_{ij}$ is the unit vector |
| 350 |
|
pointing along the inter-dipole vector $\mathbf{r}_{ij}$. |
| 351 |
|
|
| 352 |
+ |
Since the charge separation distance is so large in zwitterionic head |
| 353 |
+ |
groups (like the PC head groups), it would also be possible to use |
| 354 |
+ |
either point charges or a ``split dipole'' approximation, |
| 355 |
+ |
\begin{equation} |
| 356 |
+ |
V_{ij}({\bf \hat{u}}_{i},{\bf \hat{u}}_{j},{\bf |
| 357 |
+ |
\hat{r}}_{ij})) = \frac{1}{{4\pi \epsilon_0 }}\left[ {\frac{{\mu _i \cdot \mu _j }}{{R_{ij}^3 }} - |
| 358 |
+ |
\frac{{3\left( {\mu _i \cdot r_{ij} } \right)\left( {\mu _i \cdot |
| 359 |
+ |
r_{ij} } \right)}}{{R_{ij}^5 }}} \right] |
| 360 |
+ |
\end{equation} |
| 361 |
+ |
where $\mu _i$ and $\mu _j$ are the dipole moments of sites $i$ and |
| 362 |
+ |
$j$, $r_{ij}$ is vector between the two sites, and $R_{ij}$ is given |
| 363 |
+ |
by, |
| 364 |
+ |
\begin{equation} |
| 365 |
+ |
R_{ij} = \sqrt {r_{ij}^2 + \frac{{d_i^2 }}{4} + \frac{{d_j^2 |
| 366 |
+ |
}}{4}}. |
| 367 |
+ |
\end{equation} |
| 368 |
+ |
Here, $d_i$ and $d_j$ are effect charge separation distances |
| 369 |
+ |
associated with each of the two dipolar sites. This approximation to |
| 370 |
+ |
the multipole expansion maintains the fast fall-off of the multipole |
| 371 |
+ |
potentials but lacks the normal divergences when two polar groups get |
| 372 |
+ |
close to one another. |
| 373 |
+ |
|
| 374 |
|
For the interaction between nonequivalent uniaxial ellipsoids (in this |
| 375 |
|
case, between spheres and ellipsoids), the spheres are treated as |
| 376 |
|
ellipsoids with an aspect ratio of 1 ($d = l$) and with an well depth |
| 379 |
|
et al.} and is appropriate for dissimilar uniaxial |
| 380 |
|
ellipsoids.\cite{Cleaver96} |
| 381 |
|
|
| 382 |
< |
The solvent model in our simulations is identical to one used by |
| 383 |
< |
Marrink {\it et al.} in their dissipative particle dynamics (DPD) |
| 384 |
< |
simulation of lipid bilayers.\cite{Marrink04} This solvent bead is a |
| 385 |
< |
single site that represents four water molecules (m = 72 amu) and has |
| 386 |
< |
comparable density and diffusive behavior to liquid water. However, |
| 387 |
< |
since there are no electrostatic sites on these beads, this solvent |
| 388 |
< |
model cannot replicate the dielectric properties of water. |
| 382 |
> |
The solvent model in our simulations is similar to the one used by |
| 383 |
> |
Marrink {\it et al.} in their coarse grained simulations of lipid |
| 384 |
> |
bilayers.\cite{Marrink04} The solvent bead is a single site that |
| 385 |
> |
represents four water molecules (m = 72 amu) and has comparable |
| 386 |
> |
density and diffusive behavior to liquid water. However, since there |
| 387 |
> |
are no electrostatic sites on these beads, this solvent model cannot |
| 388 |
> |
replicate the dielectric properties of water. Note that although we |
| 389 |
> |
are using larger cutoff and switching radii than Marrink {\it et al.}, |
| 390 |
> |
our solvent density at 300 K remains at 0.944 g cm$^{-3}$, and the |
| 391 |
> |
solvent diffuses at 0.43 $\AA^2 ps^{-1}$ (roughly twice as fast as |
| 392 |
> |
liquid water). |
| 393 |
|
|
| 394 |
|
\begin{table*} |
| 395 |
|
\begin{minipage}{\linewidth} |
| 447 |
|
solvent beads (24 water molecules) per lipid. These configurations |
| 448 |
|
were then equilibrated for another $30$ ns. All simulations utilizing |
| 449 |
|
the solvent were carried out at constant pressure ($P=1$ atm) with |
| 450 |
< |
$3$D anisotropic coupling, and constant surface tension |
| 450 |
> |
$3$D anisotropic coupling, and small constant surface tension |
| 451 |
|
($\gamma=0.015$ N/m). Given the absence of fast degrees of freedom in |
| 452 |
|
this model, a time step of $50$ fs was utilized with excellent energy |
| 453 |
|
conservation. Data collection for structural properties of the |
| 454 |
|
bilayers was carried out during a final 5 ns run following the solvent |
| 455 |
< |
equilibration. All simulations were performed using the OOPSE |
| 456 |
< |
molecular modeling program.\cite{Meineke05} |
| 455 |
> |
equilibration. Orientational correlation functions and diffusion |
| 456 |
> |
constants were computed from 30 ns simulations in the microcanonical |
| 457 |
> |
(NVE) ensemble using the average volume from the end of the constant |
| 458 |
> |
pressure and surface tension runs. The timestep on these final |
| 459 |
> |
molecular dynamics runs was 25 fs. No appreciable changes in phase |
| 460 |
> |
structure were noticed upon switching to a microcanonical ensemble. |
| 461 |
> |
All simulations were performed using the {\sc oopse} molecular |
| 462 |
> |
modeling program.\cite{Meineke05} |
| 463 |
|
|
| 464 |
|
A switching function was applied to all potentials to smoothly turn |
| 465 |
< |
off the interactions between a range of $22$ and $25$ \AA. |
| 465 |
> |
off the interactions between a range of $22$ and $25$ \AA. The |
| 466 |
> |
switching function was the standard (cubic) function, |
| 467 |
> |
\begin{equation} |
| 468 |
> |
s(r) = |
| 469 |
> |
\begin{cases} |
| 470 |
> |
1 & \text{if $r \le r_{\text{sw}}$},\\ |
| 471 |
> |
\frac{(r_{\text{cut}} + 2r - 3r_{\text{sw}})(r_{\text{cut}} - r)^2} |
| 472 |
> |
{(r_{\text{cut}} - r_{\text{sw}})^3} |
| 473 |
> |
& \text{if $r_{\text{sw}} < r \le r_{\text{cut}}$}, \\ |
| 474 |
> |
0 & \text{if $r > r_{\text{cut}}$.} |
| 475 |
> |
\end{cases} |
| 476 |
> |
\label{eq:dipoleSwitching} |
| 477 |
> |
\end{equation} |
| 478 |
|
|
| 479 |
|
\section{Results} |
| 480 |
|
\label{sec:results} |
| 785 |
|
($\mu$).\label{fig:phaseDiagram}} |
| 786 |
|
\end{figure} |
| 787 |
|
|
| 788 |
< |
|
| 789 |
< |
We have also computed orientational diffusion constants for the head |
| 790 |
< |
groups from the relaxation of the second-order Legendre polynomial |
| 791 |
< |
correlation function, |
| 792 |
< |
\begin{eqnarray} |
| 793 |
< |
C_{\ell}(t) & = & \langle P_{\ell}\left({\bf \mu}_{i}(t) \cdot {\bf |
| 794 |
< |
\mu}_{i}(0) \right) \rangle \\ \\ |
| 795 |
< |
& \approx & e^{-\ell(\ell + 1) \theta t}, |
| 796 |
< |
\end{eqnarray} |
| 797 |
< |
of the head group dipoles. In this last line, we have used a simple |
| 798 |
< |
``Debye''-like model for the relaxation of the correlation function, |
| 799 |
< |
specifically in the case when $\ell = 2$. The computed orientational |
| 800 |
< |
diffusion constants are given in table \ref{tab:relaxation}. The |
| 801 |
< |
notable feature we observe is that the orientational diffusion |
| 802 |
< |
constant for the head group exhibits an order of magnitude decrease |
| 803 |
< |
upon entering the rippled phase. Our orientational correlation times |
| 804 |
< |
are substantially in excess of those provided by... |
| 788 |
> |
We have computed translational diffusion constants for lipid molecules |
| 789 |
> |
from the mean-square displacement, |
| 790 |
> |
\begin{equation} |
| 791 |
> |
D = \lim_{t\rightarrow \infty} \frac{1}{6 t} \langle {|\left({\bf r}_{i}(t) - {\bf r}_{i}(0) \right)|}^2 \rangle, |
| 792 |
> |
\end{equation} |
| 793 |
> |
of the lipid bodies. Translational diffusion constants for the |
| 794 |
> |
different head-to-tail size ratios (all at 300 K) are shown in table |
| 795 |
> |
\ref{tab:relaxation}. We have also computed orientational correlation |
| 796 |
> |
times for the head groups from fits of the second-order Legendre |
| 797 |
> |
polynomial correlation function, |
| 798 |
> |
\begin{equation} |
| 799 |
> |
C_{\ell}(t) = \langle P_{\ell}\left({\bf \mu}_{i}(t) \cdot {\bf |
| 800 |
> |
\mu}_{i}(0) \right) |
| 801 |
> |
\end{equation} |
| 802 |
> |
of the head group dipoles. The orientational correlation functions |
| 803 |
> |
appear to have multiple components in their decay: a fast ($12 \pm 2$ |
| 804 |
> |
ps) decay due to librational motion of the head groups, as well as |
| 805 |
> |
moderate ($\tau^h_{\rm mid}$) and slow ($\tau^h_{\rm slow}$) |
| 806 |
> |
components. The fit values for the moderate and slow correlation |
| 807 |
> |
times are listed in table \ref{tab:relaxation}. Standard deviations |
| 808 |
> |
in the fit time constants are quite large (on the order of the values |
| 809 |
> |
themselves). |
| 810 |
|
|
| 811 |
+ |
Sparrman and Westlund used a multi-relaxation model for NMR lineshapes |
| 812 |
+ |
observed in gel, fluid, and ripple phases of DPPC and obtained |
| 813 |
+ |
estimates of a correlation time for water translational diffusion |
| 814 |
+ |
($\tau_c$) of 20 ns.\cite{Sparrman2003} This correlation time |
| 815 |
+ |
corresponds to water bound to small regions of the lipid membrane. |
| 816 |
+ |
They further assume that the lipids can explore only a single period |
| 817 |
+ |
of the ripple (essentially moving in a nearly one-dimensional path to |
| 818 |
+ |
do so), and the correlation time can therefore be used to estimate a |
| 819 |
+ |
value for the translational diffusion constant of $2.25 \times |
| 820 |
+ |
10^{-11} m^2 s^{-1}$. The translational diffusion constants we obtain |
| 821 |
+ |
are in reasonable agreement with this experimentally determined |
| 822 |
+ |
value. However, the $T_2$ relaxation times obtained by Sparrman and |
| 823 |
+ |
Westlund are consistent with P-N vector reorientation timescales of |
| 824 |
+ |
2-5 ms. This is substantially slower than even the slowest component |
| 825 |
+ |
we observe in the decay of our orientational correlation functions. |
| 826 |
+ |
Other than the dipole-dipole interactions, our head groups have no |
| 827 |
+ |
shape anisotropy which would force them to move as a unit with |
| 828 |
+ |
neighboring molecules. This would naturally lead to P-N reorientation |
| 829 |
+ |
times that are too fast when compared with experimental measurements. |
| 830 |
|
|
| 831 |
|
\begin{table*} |
| 832 |
|
\begin{minipage}{\linewidth} |
| 833 |
|
\begin{center} |
| 834 |
< |
\caption{Rotational diffusion constants for the head groups |
| 835 |
< |
($\theta_h$) and molecular bodies ($\theta_b$) as a function of the |
| 836 |
< |
head-to-body width ratio. The orientational mobility of the head |
| 837 |
< |
groups experiences an {\it order of magnitude decrease} upon entering |
| 838 |
< |
the rippled phase, which suggests that the rippling is tied to a |
| 839 |
< |
freezing out of head group orientational freedom. Uncertainties in |
| 840 |
< |
the last digit are indicated by the values in parentheses.} |
| 841 |
< |
\begin{tabular}{lcc} |
| 834 |
> |
\caption{Fit values for the rotational correlation times for the head |
| 835 |
> |
groups ($\tau^h$) and molecular bodies ($\tau^b$) as well as the |
| 836 |
> |
translational diffusion constants for the molecule as a function of |
| 837 |
> |
the head-to-body width ratio (all at 300 K). In all of the phases, |
| 838 |
> |
the head group correlation functions decay with an fast librational |
| 839 |
> |
contribution ($12 \pm 1$ ps). There are additional moderate |
| 840 |
> |
($\tau^h_{\rm mid}$) and slow $\tau^h_{\rm slow}$ contributions to |
| 841 |
> |
orientational decay that depend strongly on the phase exhibited by the |
| 842 |
> |
lipids. The symmetric ripple phase ($\sigma_h / d = 1.35$) appears to |
| 843 |
> |
exhibit the slowest molecular reorientation.} |
| 844 |
> |
\begin{tabular}{lcccc} |
| 845 |
|
\hline |
| 846 |
< |
$\sigma_h / d$ & $\theta_h (\mu s^{-1})$ & $\theta_b (1/fs)$ \\ |
| 846 |
> |
$\sigma_h / d$ & $\tau^h_{\rm mid} (ns)$ & $\tau^h_{\rm |
| 847 |
> |
slow} (\mu s)$ & $\tau_b (\mu s)$ & $D (\times 10^{-11} m^2 s^{-1})$ \\ |
| 848 |
|
\hline |
| 849 |
< |
1.20 & $0.206(1) $ & $0.0175(5) $ \\ |
| 850 |
< |
1.28 & $0.179(2) $ & $0.055(2) $ \\ |
| 851 |
< |
1.35 & $0.025(1) $ & $0.195(3) $ \\ |
| 852 |
< |
1.41 & $0.023(1) $ & $0.024(3) $ \\ |
| 849 |
> |
1.20 & $0.4$ & $9.6$ & $9.5$ & $0.43(1)$ \\ |
| 850 |
> |
1.28 & $2.0$ & $13.5$ & $3.0$ & $5.91(3)$ \\ |
| 851 |
> |
1.35 & $3.2$ & $4.0$ & $0.9$ & $3.42(1)$ \\ |
| 852 |
> |
1.41 & $0.3$ & $23.8$ & $6.9$ & $7.16(1)$ \\ |
| 853 |
|
\end{tabular} |
| 854 |
|
\label{tab:relaxation} |
| 855 |
|
\end{center} |
| 879 |
|
is convex. These structures are held together by the extremely strong |
| 880 |
|
and directional interactions between the head groups. |
| 881 |
|
|
| 882 |
< |
Dipolar head groups are key for the maintaining the bilayer structures |
| 883 |
< |
exhibited by this model. The dipoles are likely to form head-to-tail |
| 884 |
< |
configurations even in flat configurations, but the temperatures are |
| 885 |
< |
high enough that vortex defects become prevalent in the flat phase. |
| 886 |
< |
The flat phase we observed therefore appears to be substantially above |
| 887 |
< |
the Kosterlitz-Thouless transition temperature for a planar system of |
| 888 |
< |
dipoles with this set of parameters. For this reason, it would be |
| 889 |
< |
interesting to observe the thermal behavior of the flat phase at |
| 890 |
< |
substantially lower temperatures. |
| 882 |
> |
The attractive forces holding the bilayer together could either be |
| 883 |
> |
non-directional (as in the work of Kranenburg and |
| 884 |
> |
Smit),\cite{Kranenburg2005} or directional (as we have utilized in |
| 885 |
> |
these simulations). The dipolar head groups are key for the |
| 886 |
> |
maintaining the bilayer structures exhibited by this particular model; |
| 887 |
> |
reducing the strength of the dipole has the tendency to make the |
| 888 |
> |
rippled phase disappear. The dipoles are likely to form attractive |
| 889 |
> |
head-to-tail configurations even in flat configurations, but the |
| 890 |
> |
temperatures are high enough that vortex defects become prevalent in |
| 891 |
> |
the flat phase. The flat phase we observed therefore appears to be |
| 892 |
> |
substantially above the Kosterlitz-Thouless transition temperature for |
| 893 |
> |
a planar system of dipoles with this set of parameters. For this |
| 894 |
> |
reason, it would be interesting to observe the thermal behavior of the |
| 895 |
> |
flat phase at substantially lower temperatures. |
| 896 |
|
|
| 897 |
|
One feature of this model is that an energetically favorable |
| 898 |
|
orientational ordering of the dipoles can be achieved by forming |
