| 6 |
|
|
| 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 |
| 9 |
> |
relevant lengths 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 |
| 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 rotated with euler angles: $\boldsymbol{\Omega}_i$. |
| 146 |
> |
vector rotated with Euler angles: $\boldsymbol{\Omega}_i$. |
| 147 |
|
|
| 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 |
| 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}. |
| 246 |
> |
Sec.~\ref{oopseSec:integrate}. |
| 247 |
|
|
| 248 |
|
|
| 249 |
|
\subsection{\label{lipidSec:ExpSetup}Experimental Setup} |
| 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 |
| 309 |
> |
into phospholipid macro-structures: $\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} |
| 314 |
> |
\subsection{\label{lipidSec:diffusion}Lateral Diffusion Constants} |
| 315 |
> |
|
| 316 |
> |
The lateral diffusion constant, $D_L$, is the constant characterizing |
| 317 |
> |
the diffusive motion of the lipid within the plane of the bilayer. It |
| 318 |
> |
is given by the following Einstein relation valid at long |
| 319 |
> |
times:\cite{allen87:csl} |
| 320 |
> |
\begin{equation} |
| 321 |
> |
2tD_L = \frac{1}{2}\langle |\mathbf{r}(t) - \mathbf{r}(0)|^2\rangle |
| 322 |
> |
\end{equation} |
| 323 |
> |
Where $\mathbf{r}(t)$ is the position of the lipid at time $t$, and is |
| 324 |
> |
constrained to lie within a plane. For the bilayer simulations the |
| 325 |
> |
plane of constrained motion was that perpendicular to the bilayer |
| 326 |
> |
normal, namely the $xy$-plane. |
| 327 |
> |
|
| 328 |
> |
Fig.~\ref{lipidFig:diffusionFig} shows the lateral diffusion constants |
| 329 |
> |
as a function of temperature. There is a definite increase in the |
| 330 |
> |
lateral diffusion with higher temperatures, which is exactly what one |
| 331 |
> |
would expect with greater fluidity of the chains. However, the |
| 332 |
> |
diffusion constants are all two orders of magnitude smaller than those |
| 333 |
> |
typical of DPPC.\cite{Cevc87} This is counter-intuitive as the DPPC |
| 334 |
> |
molecule is sterically larger and heavier than our model. This could |
| 335 |
> |
be an indication that our model's chains are too interwoven and hinder |
| 336 |
> |
the motion of the lipid, or that a simplification in the model's |
| 337 |
> |
forces has led to a slowing of diffusive behavior within the |
| 338 |
> |
bilayer. In contrast, the diffusion constant of the {\sc ssd} water, |
| 339 |
> |
$9.84\times 10^{-6}\,\text{cm}^2/\text{s}$, compares favorably with |
| 340 |
> |
that of bulk water ($2.2999\times |
| 341 |
> |
10^{-5}\,\text{cm}^2/\text{s}$\cite{Holz00}). |
| 342 |
|
|
| 343 |
+ |
\begin{figure} |
| 344 |
+ |
\centering |
| 345 |
+ |
\includegraphics[width=\linewidth]{diffusionFig.eps} |
| 346 |
+ |
\caption[The lateral diffusion constants versus temperature]{The lateral diffusion constants for the bilayers as a function of temperature.} |
| 347 |
+ |
\label{lipidFig:diffusionFig} |
| 348 |
+ |
\end{figure} |
| 349 |
+ |
|
| 350 |
+ |
\subsection{\label{lipidSec:densProf}Density Profile} |
| 351 |
+ |
|
| 352 |
+ |
Fig.~\ref{lipidFig:densityProfile} illustrates the densities of the |
| 353 |
+ |
atoms in the bilayer systems normalized by the bulk density as a |
| 354 |
+ |
function of distance from the center of the box. The profile is taken |
| 355 |
+ |
along the bilayer normal, in this case the $z$ axis. The profile at |
| 356 |
+ |
270~K shows several structural features that are largely smoothed out |
| 357 |
+ |
by 300~K. The left peak for the {\sc head} atoms is split at 270~K, |
| 358 |
+ |
implying that some freezing of the structure might already be occurring |
| 359 |
+ |
at this temperature. From the dynamics, the tails at this temperature |
| 360 |
+ |
are very much fluid, but the profile could indicate that a phase |
| 361 |
+ |
transition may simply be beyond the length scale of the current |
| 362 |
+ |
simulation. In all profiles, the water penetrates almost |
| 363 |
+ |
5~$\mbox{\AA}$ into the bilayer, completely solvating the {\sc head} |
| 364 |
+ |
atoms. The $\text{{\sc ch}}_3$ atoms although mainly centered at the |
| 365 |
+ |
middle of the bilayer, show appreciable penetration into the head |
| 366 |
+ |
group region. This indicates that the chains have enough mobility to |
| 367 |
+ |
bend back upward to allow the ends to explore areas around the {\sc |
| 368 |
+ |
head} atoms. It is unlikely that this is penetration from a lipid of |
| 369 |
+ |
the opposite face, as the lipids are only 12~$\mbox{\AA}$ in length, |
| 370 |
+ |
and the typical leaf spacing as measured from the {\sc head-head} |
| 371 |
+ |
spacing in the profile is 17.5~$\mbox{\AA}$. |
| 372 |
+ |
|
| 373 |
+ |
\begin{figure} |
| 374 |
+ |
\centering |
| 375 |
+ |
\includegraphics[width=\linewidth]{densityProfile.eps} |
| 376 |
+ |
\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.} |
| 377 |
+ |
\label{lipidFig:densityProfile} |
| 378 |
+ |
\end{figure} |
| 379 |
+ |
|
| 380 |
+ |
|
| 381 |
+ |
\subsection{\label{lipidSec:scd}$\text{S}_{\text{{\sc cd}}}$ Order Parameters} |
| 382 |
+ |
|
| 383 |
|
The $\text{S}_{\text{{\sc cd}}}$ order parameter is often reported in |
| 384 |
< |
the experimental charecterizations of phospholipids. It is obtained |
| 384 |
> |
the experimental characterizations of phospholipids. It is obtained |
| 385 |
|
through deuterium NMR, and measures the ordering of the carbon |
| 386 |
|
deuterium bond in relation to the bilayer normal at various points |
| 387 |
|
along the chains. In our model, there are no explicit hydrogens, but |
| 393 |
|
\end{equation} |
| 394 |
|
Where $S_{ij}$ is given by: |
| 395 |
|
\begin{equation} |
| 396 |
< |
S_{ij} = \frac{1}{2}\Bigl<(3\cos\Theta_i\cos\Theta_j - \delta_{ij})\Bigr> |
| 396 |
> |
S_{ij} = \frac{1}{2}\Bigl\langle(3\cos\Theta_i\cos\Theta_j |
| 397 |
> |
- \delta_{ij})\Bigr\rangle |
| 398 |
|
\label{lipidEq:scd2} |
| 399 |
|
\end{equation} |
| 400 |
< |
Here, $\Theta_i$ is the angle the $i$th carbon atom frame axis makes |
| 401 |
< |
with the bilayer normal. The brackets denote an average over time and |
| 402 |
< |
molecules. The carbon atom axes are defined: |
| 400 |
> |
Here, $\Theta_i$ is the angle the $i$th axis in the reference frame of |
| 401 |
> |
the carbon atom makes with the bilayer normal. The brackets denote an |
| 402 |
> |
average over time and molecules. The carbon atom axes are defined: |
| 403 |
|
$\mathbf{\hat{z}}\rightarrow$ vector from $C_{n-1}$ to $C_{n+1}$; |
| 404 |
< |
$\mathbf{\hat{y}}\rightarrow$ vector that is perpindicular to $z$ and |
| 404 |
> |
$\mathbf{\hat{y}}\rightarrow$ vector that is perpendicular to $z$ and |
| 405 |
|
in the plane through $C_{n-1}$, $C_{n}$, and $C_{n+1}$; |
| 406 |
< |
$\mathbf{\hat{x}}\rightarrow$ vector perpindicular to |
| 406 |
> |
$\mathbf{\hat{x}}\rightarrow$ vector perpendicular to |
| 407 |
|
$\mathbf{\hat{y}}$ and $\mathbf{\hat{z}}$. |
| 408 |
|
|
| 409 |
|
The order parameter has a range of $[1,-\frac{1}{2}]$. A value of 1 |
| 410 |
|
implies full order aligned to the bilayer axis, 0 implies full |
| 411 |
< |
disorder, and $-\frac{1}{2}$ implies full order perpindicular to the |
| 411 |
> |
disorder, and $-\frac{1}{2}$ implies full order perpendicular to the |
| 412 |
|
bilayer axis. The {\sc cd} bond vector for carbons near the head group |
| 413 |
< |
are usually ordered perpindicular to the bilayer normal, with tails |
| 414 |
< |
farther away tending toward disorder. This makes the order paramter |
| 413 |
> |
are usually ordered perpendicular to the bilayer normal, with tails |
| 414 |
> |
farther away tending toward disorder. This makes the order parameter |
| 415 |
|
negative for most carbons, and as such $|S_{\text{{\sc cd}}}|$ is more |
| 416 |
|
commonly reported than $S_{\text{{\sc cd}}}$. |
| 417 |
|
|
| 418 |
+ |
Fig.~\ref{lipidFig:scdFig} shows the $S_{\text{{\sc cd}}}$ order |
| 419 |
+ |
parameters for the bilayer system at 300~K. There is no appreciable |
| 420 |
+ |
difference in the plots for the various temperatures, however, there |
| 421 |
+ |
is a larger difference between our models ordering, and that of |
| 422 |
+ |
DMPC. As our values are closer to $-\frac{1}{2}$, this implies more |
| 423 |
+ |
ordering perpendicular to the normal than in a real system. This is |
| 424 |
+ |
due to the model having only one carbon group separating the chains |
| 425 |
+ |
from the top of the lipid. In DMPC, with the flexibility inherent in a |
| 426 |
+ |
multiple atom head group, as well as a glycerol linkage between the |
| 427 |
+ |
head group and the acyl chains, there is more loss of ordering by the |
| 428 |
+ |
point when the chains start. |
| 429 |
|
|
| 351 |
– |
|
| 352 |
– |
|
| 430 |
|
\begin{figure} |
| 431 |
|
\centering |
| 432 |
|
\includegraphics[width=\linewidth]{scdFig.eps} |
| 434 |
|
\label{lipidFig:scdFig} |
| 435 |
|
\end{figure} |
| 436 |
|
|
| 437 |
+ |
\subsection{\label{lipidSec:p2Order}$P_2$ Order Parameter} |
| 438 |
|
|
| 439 |
< |
\begin{figure} |
| 440 |
< |
\centering |
| 441 |
< |
\includegraphics[width=\linewidth]{densityProfile.eps} |
| 442 |
< |
\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.} |
| 443 |
< |
\label{lipidFig:densityProfile} |
| 444 |
< |
\end{figure} |
| 439 |
> |
The $P_2$ order parameter allows us to measure the amount of |
| 440 |
> |
directional ordering that exists in the bilayer. Each lipid molecule |
| 441 |
> |
can be thought of as a cylindrical tube with the head group at the |
| 442 |
> |
top. If all of the cylinders are perfectly aligned, the $P_2$ order |
| 443 |
> |
parameter will be $1.0$. If the cylinders are completely dispersed, |
| 444 |
> |
the $P_2$ order parameter will be 0. For a collection of unit vectors, |
| 445 |
> |
the $P_2$ order parameter can be solved via the following |
| 446 |
> |
method.\cite{zannoni94} |
| 447 |
|
|
| 448 |
+ |
Define an ordering matrix $\mathbf{Q}$, such that, |
| 449 |
+ |
\begin{equation} |
| 450 |
+ |
\mathbf{Q} = \frac{1}{N}\sum_i^N % |
| 451 |
+ |
\begin{pmatrix} % |
| 452 |
+ |
u_{ix}u_{ix}-\frac{1}{3} & u_{ix}u_{iy} & u_{ix}u_{iz} \\ |
| 453 |
+ |
u_{iy}u_{ix} & u_{iy}u_{iy}-\frac{1}{3} & u_{iy}u_{iz} \\ |
| 454 |
+ |
u_{iz}u_{ix} & u_{iz}u_{iy} & u_{iz}u_{iz}-\frac{1}{3} % |
| 455 |
+ |
\end{pmatrix} |
| 456 |
+ |
\label{lipidEq:po1} |
| 457 |
+ |
\end{equation} |
| 458 |
+ |
Where the $u_{i\alpha}$ is the $\alpha$ element of the unit vector |
| 459 |
+ |
$\mathbf{\hat{u}}_i$, and the sum over $i$ averages over the whole |
| 460 |
+ |
collection of unit vectors. This allows the matrix element |
| 461 |
+ |
$Q_{\alpha\beta}$ to be written: |
| 462 |
+ |
\begin{equation} |
| 463 |
+ |
Q_{\alpha\beta} = \langle u_{\alpha}u_{\beta} - |
| 464 |
+ |
\frac{1}{3}\delta_{\alpha\beta} \rangle |
| 465 |
+ |
\label{lipidEq:po2} |
| 466 |
+ |
\end{equation} |
| 467 |
|
|
| 468 |
+ |
Having constructed the matrix, diagonalizing $\mathbf{Q}$ yields three |
| 469 |
+ |
eigenvalues and eigenvectors. The eigenvector associated with the |
| 470 |
+ |
largest eigenvalue, $\lambda_{\text{max}}$, is the director for the |
| 471 |
+ |
system of unit vectors. The director is the average direction all of |
| 472 |
+ |
the unit vectors are pointing. The $P_2$ order parameter is then |
| 473 |
+ |
simply |
| 474 |
+ |
\begin{equation} |
| 475 |
+ |
\langle P_2 \rangle = \frac{3}{2}\lambda_{\text{max}} |
| 476 |
+ |
\label{lipidEq:po3} |
| 477 |
+ |
\end{equation} |
| 478 |
|
|
| 479 |
< |
\begin{figure} |
| 480 |
< |
\centering |
| 481 |
< |
\includegraphics[width=\linewidth]{diffusionFig.eps} |
| 482 |
< |
\caption[The lateral difusion constants versus temperature]{The lateral diffusion constants for the bilayers as a function of temperature.} |
| 483 |
< |
\label{lipidFig:diffusionFig} |
| 484 |
< |
\end{figure} |
| 479 |
> |
Table~\ref{lipidTab:blSummary} summarizes the $P_2$ values for the |
| 480 |
> |
bilayers, as well as the dipole orientations. The unit vector for the |
| 481 |
> |
lipid molecules was defined by finding the moment of inertia for each |
| 482 |
> |
lipid, then setting $\mathbf{\hat{u}}$ to point along the axis of |
| 483 |
> |
minimum inertia. For the {\sc head} atoms, the unit vector simply |
| 484 |
> |
pointed in the same direction as the dipole moment. For the lipid |
| 485 |
> |
molecules, the ordering was consistent across all temperatures, with |
| 486 |
> |
the director pointed along the $z$ axis of the box. More |
| 487 |
> |
interestingly, is the high degree of ordering the dipoles impose on |
| 488 |
> |
the {\sc head} atoms. The directors for the dipoles consistently |
| 489 |
> |
pointed along the plane of the bilayer, with the directors |
| 490 |
> |
anti-aligned on the top and bottom leaf. |
| 491 |
|
|
| 492 |
|
\begin{table} |
| 493 |
|
\caption[Structural properties of the bilayers]{Bilayer Structural properties as a function of temperature.} |
| 494 |
+ |
\label{lipidTab:blSummary} |
| 495 |
|
\begin{center} |
| 496 |
|
\begin{tabular}{|c|c|c|c|c|} |
| 497 |
|
\hline |
| 498 |
< |
Temperature (K) & $<L_{\perp}>$ ($\mbox{\AA}$) & % |
| 499 |
< |
$<A_{\parallel}>$ ($\mbox{\AA}^2$) & $<P_2>_{\text{Lipid}}$ & % |
| 500 |
< |
$<P_2>_{\text{{\sc head}}}$ \\ \hline |
| 498 |
> |
Temperature (K) & $\langle L_{\perp}\rangle$ ($\mbox{\AA}$) & % |
| 499 |
> |
$\langle A_{\parallel}\rangle$ ($\mbox{\AA}^2$) & % |
| 500 |
> |
$\langle P_2\rangle_{\text{Lipid}}$ & % |
| 501 |
> |
$\langle P_2\rangle_{\text{{\sc head}}}$ \\ \hline |
| 502 |
|
270 & 18.1 & 58.1 & 0.253 & 0.494 \\ \hline |
| 503 |
|
275 & 17.2 & 56.7 & 0.295 & 0.514 \\ \hline |
| 504 |
|
277 & 16.9 & 58.0 & 0.301 & 0.541 \\ \hline |
| 510 |
|
\end{tabular} |
| 511 |
|
\end{center} |
| 512 |
|
\end{table} |
| 513 |
+ |
|
| 514 |
+ |
\subsection{\label{lipidSec:miscData}Further Bilayer Data} |
| 515 |
+ |
|
| 516 |
+ |
Also summarized in Table~\ref{lipidTab:blSummary}, are the bilayer |
| 517 |
+ |
thickness and area per lipid. The bilayer thickness was measured from |
| 518 |
+ |
the peak to peak {\sc head} atom distance in the density profiles. The |
| 519 |
+ |
area per lipid data compares favorably with values typically seen for |
| 520 |
+ |
DMPC (60.0~$\mbox{\AA}^2$ at 303~K)\cite{petrache00}. Although are |
| 521 |
+ |
values are lower this is most likely due to the shorter chain length |
| 522 |
+ |
of our model (8 versus 14 for DMPC). |
| 523 |
+ |
|
| 524 |
+ |
\subsection{\label{lipidSec:randBilayer}Bilayer Aggregation} |
| 525 |
+ |
|
| 526 |
+ |
A very important accomplishment for our model is its ability to |
| 527 |
+ |
spontaneously form bilayers from a randomly dispersed starting |
| 528 |
+ |
configuration. Fig.~\ref{lipidFig:blImage} shows an image sequence for |
| 529 |
+ |
the bilayer aggregation. After 3.0~ns, the basic form of the bilayer |
| 530 |
+ |
can already be seen. By 7.0~ns, the bilayer has a lipid bridge |
| 531 |
+ |
stretched across the simulation box to itself that will turn out to be |
| 532 |
+ |
very long lived ($\sim$20~ns), as well as a water pore, that will |
| 533 |
+ |
persist for the length of the current simulation. At 24~ns, the lipid |
| 534 |
+ |
bridge is dispersed, and the bilayer is still integrating the lipid |
| 535 |
+ |
molecules from the bridge into itself, and has still been unable to |
| 536 |
+ |
disperse the water pore. |
| 537 |
+ |
|
| 538 |
+ |
\begin{figure} |
| 539 |
+ |
\centering |
| 540 |
+ |
\includegraphics[width=\linewidth]{bLayerImage.eps} |
| 541 |
+ |
\caption[Image sequence of the bilayer aggregation]{Image sequence of the bilayer aggregation. The blue beads are the {\sc head} atoms the grey beads are the chains, and the red and white bead are the water molecules. A box has been drawn around the periodic image.} |
| 542 |
+ |
\label{lipidFig:blImage} |
| 543 |
+ |
\end{figure} |
| 544 |
+ |
|
| 545 |
+ |
\subsection{\label{lipidSec:randIrod}Inverted Rod Aggregation} |
| 546 |
+ |
|
| 547 |
+ |
Fig.~\ref{lipidFig:iRimage} shows a second aggregation sequence |
| 548 |
+ |
simulated in this research. Here the fraction of water had been |
| 549 |
+ |
significantly decreased to observe how the model would respond. After |
| 550 |
+ |
1.5~ns, The main body of water in the system has already collected |
| 551 |
+ |
into a central water channel. By 10.0~ns, the channel has widened |
| 552 |
+ |
slightly, but there are still many sub channels permeating the lipid |
| 553 |
+ |
macro-structure. At 23.0~ns, the central water channel has stabilized |
| 554 |
+ |
and several smaller water channels have been absorbed to main |
| 555 |
+ |
one. However, there are still several other channels that persist |
| 556 |
+ |
through the lipid structure. |
| 557 |
+ |
|
| 558 |
+ |
\begin{figure} |
| 559 |
+ |
\centering |
| 560 |
+ |
\includegraphics[width=\linewidth]{iRodImage.eps} |
| 561 |
+ |
\caption[Image sequence of the inverted rod aggregation]{Image sequence of the inverted rod aggregation. color scheme is the same as in Fig.~\ref{lipidFig:blImage}.} |
| 562 |
+ |
\label{lipidFig:iRimage} |
| 563 |
+ |
\end{figure} |
| 564 |
+ |
|
| 565 |
+ |
\section{\label{lipidSec:Conclusion}Conclusion} |
| 566 |
+ |
|
| 567 |
+ |
We have presented a phospholipid model capable of spontaneous |
| 568 |
+ |
aggregation into a bilayer and an inverted rod structure. The time |
| 569 |
+ |
scales of the macro-molecular aggregations are in excess of 24~ns. In |
| 570 |
+ |
addition the model's bilayer properties have been explored over a |
| 571 |
+ |
range of temperatures through prefabricated bilayers. No freezing |
| 572 |
+ |
transition is seen in the temperature range of our current |
| 573 |
+ |
simulations. However, structural information from the lowest |
| 574 |
+ |
temperature may imply that a freezing event is on a much longer time |
| 575 |
+ |
scale than that explored in this current research. Further studies of |
| 576 |
+ |
this system could extend the time length of the simulations at the low |
| 577 |
+ |
temperatures to observe whether lipid crystallization occurs within the |
| 578 |
+ |
framework of this model. |
