| 4 |
|
|
| 5 |
|
\section{\label{lipidSec:Intro}Introduction} |
| 6 |
|
|
| 7 |
< |
In the past 10 years, computer speeds have allowed for the atomistic |
| 8 |
< |
simulation of phospholipid bilayers. These simulations have ranged |
| 9 |
< |
from simulation of the gel phase ($L_{\beta}$) of |
| 10 |
< |
dipalmitoylphosphatidylcholine (DPPC), \cite{Lindahl:2000} to the |
| 7 |
> |
In the past 10 years, increasing computer speeds have allowed for the |
| 8 |
> |
atomistic simulation of phospholipid bilayers for increasingly |
| 9 |
> |
relevant lenghths of time. These simulations have ranged from |
| 10 |
> |
simulation of the gel phase ($L_{\beta}$) of |
| 11 |
> |
dipalmitoylphosphatidylcholine (DPPC),\cite{lindahl00} to the |
| 12 |
|
spontaneous aggregation of DPPC molecules into fluid phase |
| 13 |
< |
($L_{\alpha}$ bilayers. \cite{Marrinck:2001} With the exception of a |
| 14 |
< |
few ambitious |
| 15 |
< |
simulations,\cite{Marrinch:2001b,Marrinck:2002,Lindahl:2000} most |
| 16 |
< |
investigations are limited to 64 to 256 |
| 17 |
< |
phospholipids.\cite{Lindal:2000,Sum:2003,Venable:2000,Gomez:2003,Smondyrev:1999,Marrinck:2001a} |
| 18 |
< |
This is due to the expense of the computer calculations involved when |
| 19 |
< |
performing these simulations. To properly hydrate a bilayer, one |
| 13 |
> |
($L_{\alpha}$) bilayers.\cite{marrink01} With the exception of a few |
| 14 |
> |
ambitious |
| 15 |
> |
simulations,\cite{marrink01:undulation,marrink:2002,lindahl00} most |
| 16 |
> |
investigations are limited to a range of 64 to 256 |
| 17 |
> |
phospholipids.\cite{lindahl00,sum:2003,venable00,gomez:2003,smondyrev:1999,marrink01} |
| 18 |
> |
The expense of the force calculations involved when performing these |
| 19 |
> |
simulations limits the system size. To properly hydrate a bilayer, one |
| 20 |
|
typically needs 25 water molecules for every lipid, bringing the total |
| 21 |
|
number of atoms simulated to roughly 8,000 for a system of 64 DPPC |
| 22 |
< |
molecules. Added to the difficluty is the electrostatic nature of the |
| 23 |
< |
phospholipid head groups and water, requiring the computationally |
| 24 |
< |
expensive Ewald sum or its slightly faster derivative particle mesh |
| 25 |
< |
Ewald sum.\cite{Nina:2002,Norberg:2000,Patra:2003} These factors all |
| 26 |
< |
limit the potential size and time lenghts of bilayer simulations. |
| 22 |
> |
molecules. Added to the difficulty is the electrostatic nature of the |
| 23 |
> |
phospholipid head groups and water, requiring either the |
| 24 |
> |
computationally expensive Ewald sum or the faster, particle mesh Ewald |
| 25 |
> |
sum.\cite{nina:2002,norberg:2000,patra:2003} These factors all limit |
| 26 |
> |
the system size and time scales of bilayer simulations. |
| 27 |
|
|
| 28 |
|
Unfortunately, much of biological interest happens on time and length |
| 29 |
< |
scales unfeasible with current simulation. One such example is the |
| 30 |
< |
observance of a ripple phase ($P_{\beta'}$) between the $L_{\beta}$ |
| 31 |
< |
and $L_{\alpha}$ phases of certain phospholipid |
| 32 |
< |
bilayers.\cite{Katsaras:2000,Sengupta:2000} These ripples are shown to |
| 33 |
< |
have periodicity on the order of 100-200~$\mbox{\AA}$. A simulation on |
| 34 |
< |
this length scale would have approximately 1,300 lipid molecules with |
| 35 |
< |
an additional 25 water molecules per lipid to fully solvate the |
| 36 |
< |
bilayer. A simulation of this size is impractical with current |
| 37 |
< |
atomistic models. |
| 29 |
> |
scales well beyond the range of current simulation technology. One |
| 30 |
> |
such example is the observance of a ripple phase |
| 31 |
> |
($P_{\beta^{\prime}}$) between the $L_{\beta}$ and $L_{\alpha}$ phases |
| 32 |
> |
of certain phospholipid bilayers.\cite{katsaras00,sengupta00} These |
| 33 |
> |
ripples are shown to have periodicity on the order of |
| 34 |
> |
100-200~$\mbox{\AA}$. A simulation on this length scale would have |
| 35 |
> |
approximately 1,300 lipid molecules with an additional 25 water |
| 36 |
> |
molecules per lipid to fully solvate the bilayer. A simulation of this |
| 37 |
> |
size is impractical with current atomistic models. |
| 38 |
|
|
| 39 |
< |
Another class of simulations to consider, are those dealing with the |
| 40 |
< |
diffusion of molecules through a bilayer. Due to the fluid-like |
| 41 |
< |
properties of a lipid membrane, not all diffusion across the membrane |
| 42 |
< |
happens at pores. Some molecules of interest may incorporate |
| 43 |
< |
themselves directly into the membrane. Once here, they may possess an |
| 44 |
< |
appreciable waiting time (on the order of 10's to 100's of |
| 45 |
< |
nanoseconds) within the bilayer. Such long simulation times are |
| 46 |
< |
difficulty to obtain when integrating the system with atomistic |
| 46 |
< |
detail. |
| 39 |
> |
The time and length scale limitations are most striking in transport |
| 40 |
> |
phenomena. Due to the fluid-like properties of a lipid membrane, not |
| 41 |
> |
all diffusion across the membrane happens at pores. Some molecules of |
| 42 |
> |
interest may incorporate themselves directly into the membrane. Once |
| 43 |
> |
here, they may possess an appreciable waiting time (on the order of |
| 44 |
> |
10's to 100's of nanoseconds) within the bilayer. Such long simulation |
| 45 |
> |
times are nearly impossible to obtain when integrating the system with |
| 46 |
> |
atomistic detail. |
| 47 |
|
|
| 48 |
< |
Addressing these issues, several schemes have been proposed. One |
| 49 |
< |
approach by Goetz and Liposky\cite{Goetz:1998} is to model the entire |
| 48 |
> |
To address these issues, several schemes have been proposed. One |
| 49 |
> |
approach by Goetz and Liposky\cite{goetz98} is to model the entire |
| 50 |
|
system as Lennard-Jones spheres. Phospholipids are represented by |
| 51 |
|
chains of beads with the top most beads identified as the head |
| 52 |
|
atoms. Polar and non-polar interactions are mimicked through |
| 53 |
|
attractive and soft-repulsive potentials respectively. A similar |
| 54 |
< |
model proposed by Marrinck \emph{et. al.}\cite{Marrinck:2004}~ uses a |
| 54 |
> |
model proposed by Marrinck \emph{et. al}.\cite{marrink04}~uses a |
| 55 |
|
similar technique for modeling polar and non-polar interactions with |
| 56 |
|
Lennard-Jones spheres. However, they also include charges on the head |
| 57 |
|
group spheres to mimic the electrostatic interactions of the |
| 63 |
|
interactions than the previous two models, while still balancing the |
| 64 |
|
need for simplifications over atomistic detail. The model uses |
| 65 |
|
Lennard-Jones spheres for the head and tail groups of the |
| 66 |
< |
phopholipids, allowing for the ability to scale the parameters to |
| 66 |
> |
phospholipids, allowing for the ability to scale the parameters to |
| 67 |
|
reflect various sized chain configurations while keeping the number of |
| 68 |
|
interactions small. What sets this model apart, however, is the use |
| 69 |
< |
of dipoles to represent the electrosttaic nature of the |
| 69 |
> |
of dipoles to represent the electrostatic nature of the |
| 70 |
|
phospholipids. The dipole electrostatic interaction is shorter range |
| 71 |
< |
than coulombic ($\frac{1}{r^3}$ versus $\frac{1}{r}$), eliminating the |
| 72 |
< |
need for a costly Ewald sum. |
| 71 |
> |
than Coulombic ($\frac{1}{r^3}$ versus $\frac{1}{r}$), and eliminates |
| 72 |
> |
the need for a costly Ewald sum. |
| 73 |
|
|
| 74 |
|
Another key feature of this model, is the use of a dipolar water model |
| 75 |
< |
to represent the solvent. The soft sticky dipole ({\scssd}) |
| 76 |
< |
water \cite{Liu:1996a} relies on the dipole for long range |
| 77 |
< |
electrostatic effects, butalso contains a short range correction for |
| 78 |
< |
hydrogen bonding. In this way the systems in this research mimic the |
| 79 |
< |
entropic contribution to the hydrophobic effect due to hydrogen-bond |
| 80 |
< |
network deformation around a non-polar entity, \emph{i.e.}~ the |
| 81 |
< |
phospholipid. |
| 75 |
> |
to represent the solvent. The soft sticky dipole ({\sc ssd}) water |
| 76 |
> |
\cite{liu96:new_model} relies on the dipole for long range electrostatic |
| 77 |
> |
effects, but also contains a short range correction for hydrogen |
| 78 |
> |
bonding. In this way the systems in this research mimic the entropic |
| 79 |
> |
contribution to the hydrophobic effect due to hydrogen-bond network |
| 80 |
> |
deformation around a non-polar entity, \emph{i.e.}~the phospholipid |
| 81 |
> |
molecules. |
| 82 |
|
|
| 83 |
|
The following is an outline of this chapter. |
| 84 |
< |
Sec.~\ref{lipoidSec:Methods} is an introduction to the lipid model |
| 85 |
< |
used in these simulations. As well as clarification about the water |
| 86 |
< |
model and integration techniques. The various simulation setups |
| 87 |
< |
explored in this research are outlined in |
| 88 |
< |
Sec.~\ref{lipidSec:ExpSetup}. Sec.~\ref{lipidSec:Results} and |
| 89 |
< |
Sec.~\ref{lipidSec:Discussion} give a summary of the results and |
| 90 |
< |
interpretation of those results respectively. Finally, the |
| 91 |
< |
conclusions of this chapter are presented in |
| 92 |
< |
Sec.~\ref{lipidSec:Conclusion}. |
| 84 |
> |
Sec.~\ref{lipidSec:Methods} is an introduction to the lipid model used |
| 85 |
> |
in these simulations, as well as clarification about the water model |
| 86 |
> |
and integration techniques. The various simulations explored in this |
| 87 |
> |
research are outlined in |
| 88 |
> |
Sec.~\ref{lipidSec:ExpSetup}. Sec.~\ref{lipidSec:resultsDis} gives a |
| 89 |
> |
summary and interpretation of the results. Finally, the conclusions |
| 90 |
> |
of this chapter are presented in Sec.~\ref{lipidSec:Conclusion}. |
| 91 |
|
|
| 92 |
|
\section{\label{lipidSec:Methods}Methods} |
| 93 |
|
|
| 96 |
– |
|
| 97 |
– |
|
| 94 |
|
\subsection{\label{lipidSec:lipidModel}The Lipid Model} |
| 95 |
|
|
| 96 |
|
\begin{figure} |
| 97 |
< |
|
| 98 |
< |
\caption{Schematic diagram of the single chain phospholipid model} |
| 99 |
< |
|
| 97 |
> |
\centering |
| 98 |
> |
\includegraphics[width=\linewidth]{twoChainFig.eps} |
| 99 |
> |
\caption[The two chained lipid model]{Schematic diagram of the double chain phospholipid model. The head group (in red) has a point dipole, $\boldsymbol{\mu}$, located at its center of mass. The two chains are eight methylene groups in length.} |
| 100 |
|
\label{lipidFig:lipidModel} |
| 105 |
– |
|
| 101 |
|
\end{figure} |
| 102 |
|
|
| 103 |
|
The phospholipid model used in these simulations is based on the |
| 104 |
|
design illustrated in Fig.~\ref{lipidFig:lipidModel}. The head group |
| 105 |
|
of the phospholipid is replaced by a single Lennard-Jones sphere of |
| 106 |
< |
diameter $fix$, with $fix$ scaling the well depth of its van der Walls |
| 107 |
< |
interaction. This sphere also contains a single dipole of magnitude |
| 108 |
< |
$fix$, where $fix$ can be varied to mimic the charge separation of a |
| 106 |
> |
diameter $\sigma_{\text{head}}$, with $\epsilon_{\text{head}}$ scaling |
| 107 |
> |
the well depth of its van der Walls interaction. This sphere also |
| 108 |
> |
contains a single dipole of magnitude $|\boldsymbol{\mu}|$, where |
| 109 |
> |
$|\boldsymbol{\mu}|$ can be varied to mimic the charge separation of a |
| 110 |
|
given phospholipid head group. The atoms of the tail region are |
| 111 |
< |
modeled by unified atom beads. They are free of partial charges or |
| 112 |
< |
dipoles, containing only Lennard-Jones interaction sites at their |
| 113 |
< |
centers of mass. As with the head groups, their potentials can be |
| 114 |
< |
scaled by $fix$ and $fix$. |
| 111 |
> |
modeled by beads representing multiple methyl groups. They are free |
| 112 |
> |
of partial charges or dipoles, and contain only Lennard-Jones |
| 113 |
> |
interaction sites at their centers of mass. As with the head groups, |
| 114 |
> |
their potentials can be scaled by $\sigma_{\text{tail}}$ and |
| 115 |
> |
$\epsilon_{\text{tail}}$. |
| 116 |
|
|
| 117 |
|
The long range interactions between lipids are given by the following: |
| 118 |
|
\begin{equation} |
| 119 |
< |
EQ Here |
| 119 |
> |
V_{\text{LJ}}(r_{ij}) = |
| 120 |
> |
4\epsilon_{ij} \biggl[ |
| 121 |
> |
\biggl(\frac{\sigma_{ij}}{r_{ij}}\biggr)^{12} |
| 122 |
> |
- \biggl(\frac{\sigma_{ij}}{r_{ij}}\biggr)^{6} |
| 123 |
> |
\biggr] |
| 124 |
|
\label{lipidEq:LJpot} |
| 125 |
|
\end{equation} |
| 126 |
|
and |
| 127 |
|
\begin{equation} |
| 128 |
< |
EQ Here |
| 128 |
> |
V_{\text{dipole}}(\mathbf{r}_{ij},\boldsymbol{\Omega}_{i}, |
| 129 |
> |
\boldsymbol{\Omega}_{j}) = \frac{|\mu_i||\mu_j|}{4\pi\epsilon_{0}r_{ij}^{3}} \biggl[ |
| 130 |
> |
\boldsymbol{\hat{u}}_{i} \cdot \boldsymbol{\hat{u}}_{j} |
| 131 |
> |
- |
| 132 |
> |
3(\boldsymbol{\hat{u}}_i \cdot \mathbf{\hat{r}}_{ij}) % |
| 133 |
> |
(\boldsymbol{\hat{u}}_j \cdot \mathbf{\hat{r}}_{ij} \biggr] |
| 134 |
|
\label{lipidEq:dipolePot} |
| 135 |
|
\end{equation} |
| 136 |
|
Where $V_{\text{LJ}}$ is the Lennard-Jones potential and |
| 139 |
|
parameters which scale the length and depth of the interaction |
| 140 |
|
respectively, and $r_{ij}$ is the distance between beads $i$ and $j$. |
| 141 |
|
In $V_{\text{dipole}}$, $\mathbf{r}_{ij}$ is the vector starting at |
| 142 |
< |
bead$i$ and pointing towards bead $j$. Vectors $\mathbf{\Omega}_i$ |
| 142 |
> |
bead $i$ and pointing towards bead $j$. Vectors $\mathbf{\Omega}_i$ |
| 143 |
|
and $\mathbf{\Omega}_j$ are the orientational degrees of freedom for |
| 144 |
|
beads $i$ and $j$. $|\mu_i|$ is the magnitude of the dipole moment of |
| 145 |
|
$i$, and $\boldsymbol{\hat{u}}_i$ is the standard unit orientation |
| 146 |
< |
vector of $\boldsymbol{\Omega}_i$. |
| 146 |
> |
vector rotated with euler angles: $\boldsymbol{\Omega}_i$. |
| 147 |
|
|
| 148 |
< |
The model also allows for the bonded interactions of bonds, bends, and |
| 149 |
< |
torsions. The bonds between two beads on a chain are of fixed length, |
| 150 |
< |
and are maintained according to the {\sc rattle} algorithm. \cite{fix} |
| 148 |
> |
The model also allows for the bonded interactions bends, and torsions. |
| 149 |
> |
The bond between two beads on a chain is of fixed length, and is |
| 150 |
> |
maintained according to the {\sc rattle} algorithm.\cite{andersen83} |
| 151 |
|
The bends are subject to a harmonic potential: |
| 152 |
|
\begin{equation} |
| 153 |
< |
eq here |
| 153 |
> |
V_{\text{bend}}(\theta) = k_{\theta}( \theta - \theta_0 )^2 |
| 154 |
|
\label{lipidEq:bendPot} |
| 155 |
|
\end{equation} |
| 156 |
< |
where $fix$ scales the strength of the harmonic well, and $fix$ is the |
| 157 |
< |
angle between bond vectors $fix$ and $fix$. The torsion potential is |
| 158 |
< |
given by: |
| 156 |
> |
where $k_{\theta}$ scales the strength of the harmonic well, and |
| 157 |
> |
$\theta$ is the angle between bond vectors |
| 158 |
> |
(Fig.~\ref{lipidFig:lipidModel}). In addition, we have placed a |
| 159 |
> |
``ghost'' bend on the phospholipid head. The ghost bend adds a |
| 160 |
> |
potential to keep the dipole pointed along the bilayer surface, where |
| 161 |
> |
$\theta$ is now the angle the dipole makes with respect to the {\sc |
| 162 |
> |
head}-$\text{{\sc ch}}_2$ bond vector |
| 163 |
> |
(Fig.~\ref{lipidFig:ghostBend}). The torsion potential is given by: |
| 164 |
|
\begin{equation} |
| 165 |
< |
eq here |
| 165 |
> |
V_{\text{torsion}}(\phi) = |
| 166 |
> |
k_3 \cos^3 \phi + k_2 \cos^2 \phi + k_1 \cos \phi + k_0 |
| 167 |
|
\label{lipidEq:torsionPot} |
| 168 |
|
\end{equation} |
| 169 |
|
Here, the parameters $k_0$, $k_1$, $k_2$, and $k_3$ fit the cosine |
| 170 |
|
power series to the desired torsion potential surface, and $\phi$ is |
| 171 |
< |
the angle between bondvectors $fix$ and $fix$ along the vector $fix$ |
| 172 |
< |
(see Fig.:\ref{lipidFig:lipidModel}). Long range interactions such as |
| 173 |
< |
the Lennard-Jones potential are excluded for bead pairs involved in |
| 174 |
< |
the same bond, bend, or torsion. However, internal interactions not |
| 171 |
> |
the angle the two end atoms have rotated about the middle bond |
| 172 |
> |
(Fig.:\ref{lipidFig:lipidModel}). Long range interactions such as the |
| 173 |
> |
Lennard-Jones potential are excluded for atom pairs involved in the |
| 174 |
> |
same bond, bend, or torsion. However, internal interactions not |
| 175 |
|
directly involved in a bonded pair are calculated. |
| 176 |
|
|
| 177 |
+ |
\begin{figure} |
| 178 |
+ |
\centering |
| 179 |
+ |
\includegraphics[width=\linewidth]{ghostBendFig.eps} |
| 180 |
+ |
\caption[Depiction of the ``ghost'' bend]{The ``ghost'' bend is a bend potential added to constrain the motion of the dipole on the {\sc head} group. The potential follows Eq.~\ref{lipidEq:bendPot} where $\theta$ is now the angle that the dipole makes with the {\sc head}-$\text{{\sc ch}}_2$ bond vector.} |
| 181 |
+ |
\label{lipidFig:ghostBend} |
| 182 |
+ |
\end{figure} |
| 183 |
+ |
|
| 184 |
|
All simulations presented here use a two chained lipid as pictured in |
| 185 |
< |
Fig.~\ref{lipidFig:twochain}. The chains are both eight beads long, |
| 185 |
> |
Fig.~\ref{lipidFig:lipidModel}. The chains are both eight beads long, |
| 186 |
|
and their mass and Lennard Jones parameters are summarized in |
| 187 |
|
Table~\ref{lipidTable:tcLJParams}. The magnitude of the dipole moment |
| 188 |
|
for the head bead is 10.6~Debye, and the bend and torsion parameters |
| 189 |
< |
are summarized in Table~\ref{lipidTable:teBTParams}. |
| 189 |
> |
are summarized in Table~\ref{lipidTable:tcBendParams} and |
| 190 |
> |
\ref{lipidTable:tcTorsionParams}. |
| 191 |
|
|
| 192 |
< |
\section{label{lipidSec:furtherMethod}Further Methodology} |
| 192 |
> |
\begin{table} |
| 193 |
> |
\caption{The Lennard Jones Parameters for the two chain phospholipids.} |
| 194 |
> |
\label{lipidTable:tcLJParams} |
| 195 |
> |
\begin{center} |
| 196 |
> |
\begin{tabular}{|l|c|c|c|} |
| 197 |
> |
\hline |
| 198 |
> |
& mass (amu) & $\sigma$($\mbox{\AA}$) & $\epsilon$ (kcal/mol) \\ \hline |
| 199 |
> |
{\sc head} & 72 & 4.0 & 0.185 \\ \hline |
| 200 |
> |
{\sc ch}\cite{Siepmann1998} & 13.02 & 4.0 & 0.0189 \\ \hline |
| 201 |
> |
$\text{{\sc ch}}_2$\cite{Siepmann1998} & 14.03 & 3.95 & 0.18 \\ \hline |
| 202 |
> |
$\text{{\sc ch}}_3$\cite{Siepmann1998} & 15.04 & 3.75 & 0.25 \\ \hline |
| 203 |
> |
{\sc ssd}\cite{liu96:new_model} & 18.03 & 3.051 & 0.152 \\ \hline |
| 204 |
> |
\end{tabular} |
| 205 |
> |
\end{center} |
| 206 |
> |
\end{table} |
| 207 |
|
|
| 208 |
+ |
\begin{table} |
| 209 |
+ |
\caption[Bend Parameters for the two chain phospholipids]{Bend Parameters for the two chain phospholipids. All alkane parameters are based off of those from TraPPE.\cite{Siepmann1998}} |
| 210 |
+ |
\label{lipidTable:tcBendParams} |
| 211 |
+ |
\begin{center} |
| 212 |
+ |
\begin{tabular}{|l|c|c|} |
| 213 |
+ |
\hline |
| 214 |
+ |
& $k_{\theta}$ ( kcal/($\text{mol deg}^2$) ) & $\theta_0$ ( deg ) \\ \hline |
| 215 |
+ |
{\sc ghost}-{\sc head}-$\text{{\sc ch}}_2$ & 0.00177 & 129.78 \\ \hline |
| 216 |
+ |
$x$-{\sc ch}-$y$ & 58.84 & 112.0 \\ \hline |
| 217 |
+ |
$x$-$\text{{\sc ch}}_2$-$y$ & 58.84 & 114.0 \\ \hline |
| 218 |
+ |
\end{tabular} |
| 219 |
+ |
\end{center} |
| 220 |
+ |
\end{table} |
| 221 |
+ |
|
| 222 |
+ |
\begin{table} |
| 223 |
+ |
\caption[Torsion Parameters for the two chain phospholipids]{Torsion Parameters for the two chain phospholipids. Alkane parameters based on TraPPE.\cite{Siepmann1998}} |
| 224 |
+ |
\label{lipidTable:tcTorsionParams} |
| 225 |
+ |
\begin{center} |
| 226 |
+ |
\begin{tabular}{|l|c|c|c|c|} |
| 227 |
+ |
\hline |
| 228 |
+ |
All are in kcal/mol $\rightarrow$ & $k_3$ & $k_2$ & $k_1$ & $k_0$ \\ \hline |
| 229 |
+ |
$x$-{\sc ch}-$y$-$z$ & 3.3254 & -0.4215 & -1.686 & 1.1661 \\ \hline |
| 230 |
+ |
$x$-$\text{{\sc ch}}_2$-$\text{{\sc ch}}_2$-$y$ & 5.9602 & -0.568 & -3.802 & 2.1586 \\ \hline |
| 231 |
+ |
\end{tabular} |
| 232 |
+ |
\end{center} |
| 233 |
+ |
\end{table} |
| 234 |
+ |
|
| 235 |
+ |
|
| 236 |
+ |
\section{\label{lipidSec:furtherMethod}Further Methodology} |
| 237 |
+ |
|
| 238 |
|
As mentioned previously, the water model used throughout these |
| 239 |
< |
simulations was the {\scssd} model of |
| 240 |
< |
Ichiye.\cite{liu:1996a,Liu:1996b,Chandra:1999} A discussion of the |
| 241 |
< |
model can be found in Sec.~\ref{oopseSec:SSD}. As for the integration |
| 242 |
< |
of the equations of motion, all simulations were performed in an |
| 243 |
< |
orthorhombic periodic box with a thermostat on velocities, and an |
| 244 |
< |
independent barostat on each cartesian axis $x$, $y$, and $z$. This |
| 245 |
< |
is the $\text{NPT}_{xyz}$. ensemble described in Sec.~\ref{oopseSec:Ensembles}. |
| 239 |
> |
simulations was the {\sc ssd} model of |
| 240 |
> |
Ichiye.\cite{liu96:new_model,liu96:monte_carlo,chandra99:ssd_md} A |
| 241 |
> |
discussion of the model can be found in Sec.~\ref{oopseSec:SSD}. As |
| 242 |
> |
for the integration of the equations of motion, all simulations were |
| 243 |
> |
performed in an orthorhombic periodic box with a thermostat on |
| 244 |
> |
velocities, and an independent barostat on each Cartesian axis $x$, |
| 245 |
> |
$y$, and $z$. This is the $\text{NPT}_{xyz}$. ensemble described in |
| 246 |
> |
Sec.~\ref{oopseSec:Ensembles}. |
| 247 |
|
|
| 248 |
|
|
| 249 |
|
\subsection{\label{lipidSec:ExpSetup}Experimental Setup} |
| 265 |
|
estimate for the number of waters each lipid would displace in a |
| 266 |
|
simulation box. A target number of waters was then defined which |
| 267 |
|
included the number of waters each lipid would displace, the number of |
| 268 |
< |
waters desired to solvate each lipid, and a fudge factor to pad the |
| 269 |
< |
initialization. |
| 268 |
> |
waters desired to solvate each lipid, and a factor to pad the |
| 269 |
> |
initial box with a few extra water molecules. |
| 270 |
|
|
| 271 |
|
Next, a cubic simulation box was created that contained at least the |
| 272 |
|
target number of waters in an FCC lattice (the lattice was for ease of |
| 273 |
|
placement). What followed was a RSA simulation similar to those of |
| 274 |
|
Chapt.~\ref{chapt:RSA}. The lipids were sequentially given a random |
| 275 |
|
position and orientation within the box. If a lipid's position caused |
| 276 |
< |
atomic overlap with any previously adsorbed lipid, its position and |
| 277 |
< |
orientation were rejected, and a new random adsorption site was |
| 276 |
> |
atomic overlap with any previously placed lipid, its position and |
| 277 |
> |
orientation were rejected, and a new random placement site was |
| 278 |
|
attempted. The RSA simulation proceeded until all phospholipids had |
| 279 |
< |
been adsorbed. After adsorption, all water molecules with locations |
| 280 |
< |
that overlapped with the atomic coordinates of the lipids were |
| 281 |
< |
removed. |
| 279 |
> |
been adsorbed. After placement of all lipid molecules, water |
| 280 |
> |
molecules with locations that overlapped with the atomic coordinates |
| 281 |
> |
of the lipids were removed. |
| 282 |
|
|
| 283 |
< |
Finally, water molecules were removed one by one at random until the |
| 284 |
< |
desired number of waters per lipid was reached. The typical low final |
| 285 |
< |
density for these initial configurations was not a problem, as the box |
| 286 |
< |
would shrink to an appropriate size within the first 50~ps of a |
| 287 |
< |
simulation in the $\text{NPT}_{xyz}$ ensemble. |
| 283 |
> |
Finally, water molecules were removed at random until the desired |
| 284 |
> |
number of waters per lipid was reached. The typical low final density |
| 285 |
> |
for these initial configurations was not a problem, as the box shrinks |
| 286 |
> |
to an appropriate size within the first 50~ps of a simulation in the |
| 287 |
> |
$\text{NPT}_{xyz}$ ensemble. |
| 288 |
|
|
| 289 |
< |
\subsection{\label{lipidSec:Configs}The simulation configurations} |
| 289 |
> |
\subsection{\label{lipidSec:Configs}Configurations} |
| 290 |
|
|
| 291 |
< |
Table ~\ref{lipidTable:simNames} summarizes the names and important |
| 292 |
< |
details of the simulations. The B set of simulations were all started |
| 293 |
< |
in an ordered bilayer and observed over a period of 10~ns. Simulution |
| 294 |
< |
RL was integrated for approximately 20~ns starting from a random |
| 295 |
< |
configuration as an example of spontaneous bilayer aggregation. |
| 296 |
< |
Lastly, simulation RH was also started from a random configuration, |
| 297 |
< |
but with a lesser water content and higher temperature to show the |
| 298 |
< |
spontaneous aggregation of an inverted hexagonal lamellar phase. |
| 291 |
> |
The first class of simulations were started from ordered |
| 292 |
> |
bilayers. They were all configurations consisting of 60 lipid |
| 293 |
> |
molecules with 30 lipids on each leaf, and were hydrated with 1620 |
| 294 |
> |
{\sc ssd} molecules. The original configuration was assembled |
| 295 |
> |
according to Sec.~\ref{lipidSec:ExpSetup} and simulated for a length |
| 296 |
> |
of 10~ns at 300~K. The other temperature runs were started from a |
| 297 |
> |
frame 7~ns into the 300~K simulation. Their temperatures were reset |
| 298 |
> |
with the thermostating algorithm in the $\text{NPT}_{xyz}$ |
| 299 |
> |
integrator. All of the temperature variants were also run for 10~ns, |
| 300 |
> |
with only the last 5~ns being used for accumulation of statistics. |
| 301 |
> |
|
| 302 |
> |
The second class of simulations were two configurations started from |
| 303 |
> |
randomly dispersed lipids in a ``gas'' of water. The first |
| 304 |
> |
($\text{R}_{\text{I}}$) was a simulation containing 72 lipids with |
| 305 |
> |
1800 {\sc ssd} molecules simulated at 300~K. The second |
| 306 |
> |
($\text{R}_{\text{II}}$) was 90 lipids with 1350 {\sc ssd} molecules |
| 307 |
> |
simulated at 350~K. Both simulations were integrated for more than |
| 308 |
> |
20~ns, and illustrate the spontaneous aggregation of the lipid model |
| 309 |
> |
into phospholipid macrostructures: $\text{R}_{\text{I}}$ into a |
| 310 |
> |
bilayer, and $\text{R}_{\text{II}}$ into a inverted rod. |
| 311 |
> |
|
| 312 |
> |
\section{\label{lipidSec:resultsDis}Results and Discussion} |
| 313 |
> |
|
| 314 |
> |
\subsection{\label{lipidSec:scd}$\text{S}_{\text{{\sc cd}}}$ order parameters} |
| 315 |
> |
|
| 316 |
> |
The $\text{S}_{\text{{\sc cd}}}$ order parameter is often reported in |
| 317 |
> |
the experimental charecterizations of phospholipids. It is obtained |
| 318 |
> |
through deuterium NMR, and measures the ordering of the carbon |
| 319 |
> |
deuterium bond in relation to the bilayer normal at various points |
| 320 |
> |
along the chains. In our model, there are no explicit hydrogens, but |
| 321 |
> |
the order parameter can be written in terms of the carbon ordering at |
| 322 |
> |
each point in the chain:\cite{egberts88} |
| 323 |
> |
\begin{equation} |
| 324 |
> |
S_{\text{{\sc cd}}} = \frac{2}{3}S_{xx} + \frac{1}{3}S_{yy} |
| 325 |
> |
\label{lipidEq:scd1} |
| 326 |
> |
\end{equation} |
| 327 |
> |
Where $S_{ij}$ is given by: |
| 328 |
> |
\begin{equation} |
| 329 |
> |
S_{ij} = \frac{1}{2}\Bigl<(3\cos\Theta_i\cos\Theta_j - \delta_{ij})\Bigr> |
| 330 |
> |
\label{lipidEq:scd2} |
| 331 |
> |
\end{equation} |
| 332 |
> |
Here, $\Theta_i$ is the angle the $i$th carbon atom frame axis makes |
| 333 |
> |
with the bilayer normal. The brackets denote an average over time and |
| 334 |
> |
molecules. The carbon atom axes are defined: |
| 335 |
> |
$\mathbf{\hat{z}}\rightarrow$ vector from $C_{n-1}$ to $C_{n+1}$; |
| 336 |
> |
$\mathbf{\hat{y}}\rightarrow$ vector that is perpindicular to $z$ and |
| 337 |
> |
in the plane through $C_{n-1}$, $C_{n}$, and $C_{n+1}$; |
| 338 |
> |
$\mathbf{\hat{x}}\rightarrow$ vector perpindicular to |
| 339 |
> |
$\mathbf{\hat{y}}$ and $\mathbf{\hat{z}}$. |
| 340 |
> |
|
| 341 |
> |
The order parameter has a range of $[1,-\frac{1}{2}]$. A value of 1 |
| 342 |
> |
implies full order aligned to the bilayer axis, 0 implies full |
| 343 |
> |
disorder, and $-\frac{1}{2}$ implies full order perpindicular to the |
| 344 |
> |
bilayer axis. The {\sc cd} bond vector for carbons near the head group |
| 345 |
> |
are usually ordered perpindicular to the bilayer normal, with tails |
| 346 |
> |
farther away tending toward disorder. This makes the order paramter |
| 347 |
> |
negative for most carbons, and as such $|S_{\text{{\sc cd}}}|$ is more |
| 348 |
> |
commonly reported than $S_{\text{{\sc cd}}}$. |
| 349 |
> |
|
| 350 |
> |
|
| 351 |
> |
|
| 352 |
> |
|
| 353 |
> |
\begin{figure} |
| 354 |
> |
\centering |
| 355 |
> |
\includegraphics[width=\linewidth]{scdFig.eps} |
| 356 |
> |
\caption[$\text{S}_{\text{{\sc cd}}}$ order parameter for our model]{A comparison of $|\text{S}_{\text{{\sc cd}}}|$ between our model (blue) and DMPC\cite{petrache00} (black) near 300~K.} |
| 357 |
> |
\label{lipidFig:scdFig} |
| 358 |
> |
\end{figure} |
| 359 |
> |
|
| 360 |
> |
|
| 361 |
> |
\begin{figure} |
| 362 |
> |
\centering |
| 363 |
> |
\includegraphics[width=\linewidth]{densityProfile.eps} |
| 364 |
> |
\caption[The density profile of the lipid bilayers]{The density profile of the lipid bilayers along the bilayer normal. The black lines are the {\sc head} atoms, red lines are the {\sc ch} atoms, green lines are the $\text{{\sc ch}}_2$ atoms, blue lines are the $\text{{\sc ch}}_3$ atoms, and the magenta lines are the {\sc ssd} atoms.} |
| 365 |
> |
\label{lipidFig:densityProfile} |
| 366 |
> |
\end{figure} |
| 367 |
> |
|
| 368 |
> |
|
| 369 |
> |
|
| 370 |
> |
\begin{figure} |
| 371 |
> |
\centering |
| 372 |
> |
\includegraphics[width=\linewidth]{diffusionFig.eps} |
| 373 |
> |
\caption[The lateral difusion constants versus temperature]{The lateral diffusion constants for the bilayers as a function of temperature.} |
| 374 |
> |
\label{lipidFig:diffusionFig} |
| 375 |
> |
\end{figure} |
| 376 |
> |
|
| 377 |
> |
\begin{table} |
| 378 |
> |
\caption[Structural properties of the bilayers]{Bilayer Structural properties as a function of temperature.} |
| 379 |
> |
\begin{center} |
| 380 |
> |
\begin{tabular}{|c|c|c|c|c|} |
| 381 |
> |
\hline |
| 382 |
> |
Temperature (K) & $<L_{\perp}>$ ($\mbox{\AA}$) & % |
| 383 |
> |
$<A_{\parallel}>$ ($\mbox{\AA}^2$) & $<P_2>_{\text{Lipid}}$ & % |
| 384 |
> |
$<P_2>_{\text{{\sc head}}}$ \\ \hline |
| 385 |
> |
270 & 18.1 & 58.1 & 0.253 & 0.494 \\ \hline |
| 386 |
> |
275 & 17.2 & 56.7 & 0.295 & 0.514 \\ \hline |
| 387 |
> |
277 & 16.9 & 58.0 & 0.301 & 0.541 \\ \hline |
| 388 |
> |
280 & 17.4 & 58.0 & 0.274 & 0.488 \\ \hline |
| 389 |
> |
285 & 16.9 & 57.6 & 0.270 & 0.616 \\ \hline |
| 390 |
> |
290 & 17.0 & 57.6 & 0.263 & 0.534 \\ \hline |
| 391 |
> |
293 & 17.5 & 58.0 & 0.227 & 0.643 \\ \hline |
| 392 |
> |
300 & 16.9 & 57.6 & 0.315 & 0.536 \\ \hline |
| 393 |
> |
\end{tabular} |
| 394 |
> |
\end{center} |
| 395 |
> |
\end{table} |
