A Coarse-Grained Molecular Dynamics Simulation Using NAMD Package to Reveal Aggregation Profile of Phospholipids Self-Assembly in Water
The energy profile of self-assembly process of DLPE, DLPS, DOPE, DOPS, DLiPE, and DLiPS in water was investigated by a coarse-grained molecular dynamics simulation using NAMD package. The self-assembly process was initiated from random configurations. The simulation was carried out for 160 ns. This study presented proof that there were three major self-assembled arrangements which became visible for a certain duration when the simulation took place, that is, liposome, deformed liposome, and planar bilayer. The energy profile that shows plateau at the time of these structures emerge confirmed their stability therein. Our findings have highlighted the idea that liposomes and deformed liposomes are metastable phases which eventually will turn into planar bilayer, the stable one.
Solution of phospholipid molecules can demonstrate more than one micellar structures, namely, spherical micelles, rod-like structures, liposomes, bilayers, and others due to their surfactant-like features [1–4]. These structures play important role in drug delivery systems as well as in biological systems [5, 6]. Micellar structures rely on the molecular species, composition, and also on the self-assembly pathways affected by the initial configuration [7–10]. A great deal of experimental research has been done to study self-assembly of phospholipid molecules. Nevertheless the dynamics information about the liposome formation is still hard to achieve experimentally.
Molecular dynamics computer simulation has the ability to deliver more detailed information. It is an impressive device to investigate the mechanism of self-assembly [11–14]. Conventional molecular dynamics (MD) uncovers maximum features however they are limited to small time scales. It demands a long time to equilibrate a real physical system. Hence a coarse-grained molecular dynamics (CGMD) method was built as a simplified model to carry out molecular dynamics. The CGMD models have been used to explore a variety of structural and dynamic properties in large molecular systems. The CGMD method has offered significant outcome when exploring time and length scales further than what is viable with conventional MD. While CGMD has brought important findings for understanding the phospholipid self-assembly [12, 14, 15], there is still limited information accessible for a theoretical perceptive of phospholipid self-assembly pathway.
An all atomic simulation on phospholipid aggregation by Marrink and coworkers  has shown a typical pathway for bilayer formation. However, it did not state the formation of liposomes during the course of aggregation. Applying CGMD method on DLPE, DLPS, DOPE, DOPS, DLiPE, and DLiPS we demonstrate for the first time the aggregation profile of phospholipid molecules which clearly show the formation of liposome as the metastable phase. These phospholipids have been reported as the main phospholipid component of coconut, sesame, and candlenut endosperm  which produced liposomes and planar bilayer during aggregation .
The structure of all phospholipids used in the simulation is presented in Figure 1. The model molecule was prepared using the Open Babel package . In this study, 256 phospholipid molecules were placed randomly in a cube-shaped box with a size of 8 nm using Packmol package . Residue-based coarse-graining was applied on the system based on Martini Force Field ver.2.0 [15, 19] using VMD package . The force field was parameterized to reproduce accurate thermodynamic properties . Each phospholipid molecule was represented as 10–14 beads. Water molecules were modeled by hydrophilic beads; each one represented four real water molecules. Figure 2 presents coarse-grained structure of all phospholipids used in the simulation. In Martini Force Field each bead interacts with the pair wise Lennard-Jones potential (LJ). Screen Coulomb interaction is used to model the electrostatic interaction between the zwitter ionic head groups of phospholipids.
(a) 1,2-Dilauroyl-sn-glycero-3-phosphoethanolamine (DLPE)
(b) 1,2-Dilauroyl-sn-glycero-3-phosphoserine (DLPS)
(c) 1,2-Dioleoyl-sn-glycero-3-phosphoethanolamine (DOPE)
(d) 1,2-Dioleoyl-sn-glycero-3-phosphoserine (DOPS)
(e) 1,2-Dilinoleoyl-sn-glycero-3-phosphoethanolamine (DLiPE)
(f) 1,2-Dilinoleoyl-sn-glycero-3-phosphoserine (DLiPS)
(g) 1,2-Dilinoleoyl-sn-glycero-3-phosphocholine (DLiPC)
Molecular dynamics simulations of phospholipid system began with energy minimization using NAMD package . Energy minimization was done to adjust the structure to the force field, the distribution of solvents, and especially to reduce the steric clashes that might occur in the system. This phase provided the system with the lowest energy to do the simulation. It had been marked by the achievement of energy convergency at the end of minimization, 0.6 ns.
After the minimization simulations were performed with 40 fs time step integration during the effective time of 160 ns. Simulations were conducted on periodic boundary conditions (PBC). The duration of liposome formation and the total energy systems were analyzed from the simulation results. Visualization during the simulation process was also done by VMD. The simulation was also undertaken for larger systems, that is, 1500 phospholipid molecules.
To evaluate our system, before running the simulation on the phospholipids we conducted the simulation on DLiPC which has been recognized to form liposome .
3. Results and Discussion
Molecular dynamics simulation of phospholipid molecules was able to provide an overview of the mechanism of the aggregation process and the relationship between aggregate structure and the total energy of the system. Simulation was performed on a system with 256 molecules of phospholipids in aqueous medium with density of 0.00609 atom/A3 using NAMD package. For all molecules used in this report the simulation began with a random position. In general the molecules then started to form various clusters of phospholipids with hydrocarbon tails directed to their interior. After that they began to form liposomes or half-liposomes which were followed by deformed liposomes or planar bilayer formation at the end. We think the formation of half-liposome is due to shortage of phospholipid molecules supplies in the system.
To examine the influence of the number of molecules the simulation was also performed on the system with 1500 molecules with the same density. These simulations have shown the aggregation process better. It was preceded by the formation of small clusters of phospholipids which then was followed by mergers into larger aggregates, in the form of worm-like, cup-like, tube-like, and other structures, leading to formation of liposome or planar bilayer. For some phospholipids liposomes remained stable until the simulation ends. For others the process was then followed by deformation to produce deformed liposomes that lasted until they all became a planar bilayer. Examples of various aggregate structures observed in the aggregation process are presented as snapshots on Figure 3.
|(a) Clusters of 1500 molecules DLiPE at = 14.4 ns. Water molecules are not presented for clarity purpose only|
|(b) Cup-like structure of 1500 molecules DLPE at = 124.8 ns. Water molecules are not presented for clarity purpose only|
|(c) Tube-like structure of 1500 molecules DLPE at = 96.96 ns. Water molecules are not presented for clarity purpose only|
|(d) Liposome of 1500 molecules DOPE at = 21.6–45.6 ns. This picture clearly shows some water molecules inside the liposome|
|(e) Deformed Liposome of 1500 DLiPS molecules at > 101.6 ns. This picture clearly shows some water molecules that still reside in the deformed liposome|
|(f) Planar bilayer of 1500 molecules DLiPE at = 133.6 ns|
Observation of the total energy changes during the aggregation process showed that the process was accompanied by a decrease in the total energy of the system. The energy decrease occurred in stages before reaching a minimum when the simulation was terminated. This means that metastable structure is formed which is subsequently followed by a stable structure with minimum energy state. The simulations show that liposome and deformed liposome are metastable structures while a planar bilayer is a stable structure. This finding supports the views expressed by Luisi, Lasic, and Laughlin [24–26] that liposomes are metastable structures. For that liposomes also have a finite lifetime. The change of aggregation state which was accompanied by decreases of total energy system throughout the 160 ns simulation time is represented by the aggregation of 256 DOPE molecules. The DOPE molecules showed an energy decrease when the system formed liposome, deformed liposome, and finally planar bilayer. These self-assembly structures appeared at the first, the second, and the last energy plateau on the energy profile of the simulation (Figure 4). It also showed clearly that DOPE liposome is a metastable system with lifetime for 24 ns. Simulations with a larger number of molecules showed that liposome formation occurs at a lower energy.
The transformation and the lifetime of the aggregate structures as well as the changes of system total energy throughout the period of simulation were varied for each molecular species. The data for simulation of 256 molecules of each phospholipid molecular species are presented in Table 1.
Based on the results in Table 1 we propose the aggregation mechanism of phospholipid molecules is through the stages as in Figure 5. Starting from a random configuration, there is a rapid decline in total energy and it is accompanied by formation of irregular phospholipid clusters during the decline. After that the liposome is formed at a relatively stable total energy system, liposomes then change shape (deformed, nonspherical structure) followed by a decrease in the total energy of the system. Deformation of liposomal structures occurs for a certain time and at a relatively stable total energy. In the final stage, deformed liposome releases its total energy to form a planar bilayer which is stable until the end of the simulation.
The simulation suggests that liposomal structure is a metastable structure and the aggregation eventually will produce a planar bilayer at its final stage. It is in line with the notion suggested by Zana  that liposome formed from a single species of phospholipid cannot be thermodynamically stable since the bilayer will have high bending energy in the liposomal structure. On the contrary a planar bilayer can be thermodynamically stable due to its zero curvature.
The simulations also reveal how the hydrophilic head group, the number of carbon atoms in each carbon chain, and the number of double bond affect the aggregation stages experienced by each phospholipid species and the lifespan of liposome produced. Observation on the total system energy of liposomes with ethanolamine (PE) head group shows that they have total system energy higher than the serine (PS) and choline (PC). These results are consistent with Martens and McMahon  findings on mechanism of membrane fusion which suggest that phospholipids with ethanolamine head group prefer negative curvature, while serine and choline are positive.
This work simulated the vesiculation pathway of several phospholipid species, namely, DLPE, DLPS, DOPE, DOPS, DLiPE, and DLiPS from random initial configuration. We employ Marrink’s coarse-grain model and NAMD package for this study. We have provided further evidence that there are three main self-assembled structures which materialize for a certain period of time during the simulation, that is, liposome, deformed liposome, and planar bilayer. The energy profile associated with the appearance of these structures indicated their stability. Our results have lent power to the hypothesis that liposomes and deformed liposomes are metastable phases which will then become the stable one, planar bilayer.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Dwi Hudiyanti gratefully acknowledges the financial support from BPPS Dikti.
J. N. Israelachvili, “Thermodynamic principles of self-assembly,” in Intermolecular and Surface Forces, J. N. Israelachvili, Ed., pp. 503–534, Academic Press, San Diego, Calif, USA, 3rd edition, 2011.View at: Google Scholar
D. Hudiyanti, T. J. Raharjo, and S. Noegrohati, “Investigation on the morphology and properties of aggregate structures of natural phospholipids in aqueous system using cryo-tem,” Indonesian Journal of Chemistry, vol. 12, no. 1, pp. 57–61, 2012.View at: Google Scholar
S. Šegota and D. U. I. Težak, “Spontaneous formation of vesicles,” Advances in Colloid and Interface Science, vol. 121, no. 1–3, pp. 51–75, 2006.View at: Google Scholar
R. A. Karjiban, N. S. Shaari, U. V. Gunasakaran, and M. Basri, “A coarse-grained molecular dynamics study of DLPC, DMPC, DPPC, and DSPC mixtures in aqueous solution,” Journal of Chemistry, vol. 2013, Article ID 931051, 6 pages, 2013.View at: Publisher Site | Google Scholar
A. Samad, Y. Sultana, and M. Aqil, “Liposomal drug delivery systems: an update review,” Current Drug Delivery, vol. 4, no. 4, pp. 297–305, 2007.View at: Publisher Site | Google Scholar
G. Sessa and G. Weissmann, “Phospholipid spherules (liposomes) as a model for biological membranes,” Journal of Lipid Research, vol. 9, no. 3, pp. 310–318, 1968.View at: Google Scholar
K. Yang and Y. Ma, “Computer simulations of vesicle fission induced by external amphipathic inclusions,” Journal of Physical Chemistry B, vol. 113, no. 4, pp. 1048–1057, 2009.View at: Publisher Site | Google Scholar
A. S. Rudolph, B. R. Ratna, and B. Kahn, “Self-assembling phospholipid filaments,” Nature, vol. 352, no. 6330, pp. 52–55, 1991.View at: Publisher Site | Google Scholar
X. He and F. Schmid, “Spontaneous formation of complex micelles from a homogeneous solution,” Physical Review Letters, vol. 100, no. 13, Article ID 137802, 2008.View at: Publisher Site | Google Scholar
A. J. Markvoort, R. A. Van Santen, and P. A. J. Hilbers, “Vesicle shapes from molecular dynamics simulations,” Journal of Physical Chemistry B, vol. 110, no. 45, pp. 22780–22785, 2006.View at: Publisher Site | Google Scholar
M. Venturoli, M. M. Sperotto, M. Kranenburg, and B. Smit, “Mesoscopic models of biological membranes,” Physics Reports, vol. 437, no. 1-2, pp. 1–54, 2006.View at: Publisher Site | Google Scholar
S. J. Marrink, A. H. De Vries, and A. E. Mark, “Coarse grained model for semiquantitative lipid simulations,” Journal of Physical Chemistry B, vol. 108, no. 2, pp. 750–760, 2004.View at: Publisher Site | Google Scholar
S. J. Marrink, E. Lindahl, O. Edholm, and A. E. Mark, “Simulation of the spontaneous aggregation of phospholipids into bilayers,” Journal of the American Chemical Society, vol. 123, no. 35, pp. 8638–8639, 2001.View at: Publisher Site | Google Scholar
S. V. Bennun, M. I. Hoopes, C. Xing, and R. Faller, “Coarse-grained modeling of lipids,” Chemistry and Physics of Lipids, vol. 159, no. 2, pp. 59–66, 2009.View at: Publisher Site | Google Scholar
S. J. Marrink, H. J. Risselada, S. Yefimov, D. P. Tieleman, and A. H. De Vries, “The MARTINI force field: coarse grained model for biomolecular simulations,” Journal of Physical Chemistry B, vol. 111, no. 27, pp. 7812–7824, 2007.View at: Publisher Site | Google Scholar
D. Hudiyanti, T. J. Raharjo, Narsito, and S. Noegrohati, “Isolasi dan Karakterisasi Lesitin Kelapa dan Wijen,” Agritech, vol. 32, no. 1, p. 4, 2012.View at: Google Scholar
N. M. O'Boyle, M. Banck, C. A. James, C. Morley, T. Vandermeersch, and G. R. Hutchison, “Open Babel: an open chemical toolbox,” Journal of Cheminformatics, vol. 3, article 33, 2011.View at: Publisher Site | Google Scholar
L. Martinez, R. Andrade, E. G. Birgin, and J. M. Martínez, “PACKMOL: a package for building initial configurations for molecular dynamics simulations,” Journal of Computational Chemistry, vol. 30, no. 13, pp. 2157–2164, 2009.View at: Publisher Site | Google Scholar
H. I. Ingólfsson, C. A. Lopez, J. J. Uusitalo et al., “The power of coarse graining in biomolecular simulations,” Wiley Interdisciplinary Reviews: Computational Molecular Science, vol. 4, no. 3, pp. 225–248, 2014.View at: Publisher Site | Google Scholar
W. Humphrey, A. Dalke, and K. Schulten, “VMD: visual molecular dynamics,” Journal of Molecular Graphics, vol. 14, no. 1, pp. 33–38, 1996.View at: Publisher Site | Google Scholar
S. J. Marrink and D. P. Tieleman, “Perspective on the Martini model,” Chemical Society Reviews, vol. 42, no. 16, pp. 6801–6822, 2013.View at: Publisher Site | Google Scholar
J. C. Phillips, R. Braun, W. Wang et al., “Scalable molecular dynamics with NAMD,” Journal of Computational Chemistry, vol. 26, no. 16, pp. 1781–1802, 2005.View at: Publisher Site | Google Scholar
H. J. Risselada and S. J. Marrink, “Curvature effects on lipid packing and dynamics in liposomes revealed by coarse grained molecular dynamics simulations,” Physical Chemistry Chemical Physics, vol. 11, no. 12, pp. 2056–2067, 2009.View at: Publisher Site | Google Scholar
P. L. Luisi, “Are micelles and vesicles chemical equilibrium systems?” Journal of Chemical Education, vol. 78, no. 3, pp. 380–384, 2001.View at: Publisher Site | Google Scholar
D. D. Lasic, “Spontaneous vesiculation and spontaneous liposomes,” Journal of Liposome Research, vol. 9, no. 1, pp. 43–52, 1999.View at: Publisher Site | Google Scholar
R. G. Laughlin, “Equilibrium vesicles: fact or fiction?” Colloids and Surfaces A: Physicochemical and Engineering Aspects, vol. 128, no. 1–3, pp. 27–38, 1997.View at: Publisher Site | Google Scholar
R. Zana, Dynamics of Surfactant Self-Assemblies: Micelles, Microemulsions, Vesicles, and Lyotropic Phases, Surfactant Science, Taylor & Francis, Boca Raton, Fla, USA, 2005.
S. Martens and H. T. McMahon, “Mechanisms of membrane fusion: disparate players and common principles,” Nature Reviews Molecular Cell Biology, vol. 9, no. 7, pp. 543–556, 2008.View at: Publisher Site | Google Scholar