Molecular dynamics simulations of doxorubicin in sphingomyelin-based lipid membranes

Doxorubicin (DOX) is one of the most efficient antitumor drugs employed in numerous cancer therapies. Its incorporation into lipid-based nanocarriers, such as liposomes, improves the drug targeting into tumor cells and reduces drug side effects. The carriers' lipid composition is expected to affect the interactions of DOX and its partitioning into liposomal membranes. To get a rational insight into this aspect and determine promising lipid compositions, we use numerical simulations, which provide unique information on DOX-membrane interactions at the atomic level of resolution. In particular, we combine classical molecular dynamics simulations and free energy calculations to elucidate the mechanism of penetration of a protonated Doxorubicin molecule (DOX⁺) into potential liposome membranes, here modeled as lipid bilayers based on mixtures of phosphatidylcholine (PC), sphingomyelin (SM) and cholesterol lipid molecules, of different compositions and lipid phases. Moreover, we analyze DOX⁺ partitioning into relevant regions of SM-based lipid bilayer systems using a combination of free energy methods. Our results show that DOX⁺ penetration and partitioning are facilitated into less tightly packed SM-based membranes and are dependent on lipid composition. This work paves the way to further investigations of optimal formulations for lipid-based carriers, such as those associated with pH-responsive membranes.     

The third leg: Molecular dynamics simulations of lipid binding proteins

Molecular dynamics computer simulations can provide a third leg which balances the contributions of both structural biology and binding studies performed on the lipid binding protein family. In this context, these calculations help to establish a dialogue between all three communities, by relating experimental observables with details of structure. Working towards this connection is important, since experience has shown the difficulty of inferring thermodynamic properties from a single static conformation. The challenge is exemplified by ongoing attempts to interpret the impact of mutagenesis on structure and function (i.e. binding). A detailed atomic-level understanding of this system could be achieved with the support of all three legs, paving the way towards rational design of proteins with novel specificities. This paper provides an outline of the connections possible between experiment and theory concerning lipid binding proteins.     

Molecular dynamics simulations reveal structural insights into inhibitor binding modes and mechanism of casein kinase II inhibitors

Casein kinase-II, a member of protein kinase family, plays significant role in different cellular processes such as cell growth, differentiation, proliferation, gene expression, and embryogenesis. Being a potent suppressor of apoptosis, it serves as a significant link for its association with various types of malignancies such as colorectal and breast cancer. To overcome its pathological role in various cancerous diseases, CK-II procures great consideration as a therapeutic target. This study aimed at understanding the binding mechanism and structural properties of benzimidazole derivatives by utilizing various computational tools including docking simulation, three-dimensional quantitative structure activity relationships and molecular dynamic simulation. Structure-based 3D-QSAR techniques such as CoMFA and CoMSIA models, were established from the conformations gained by protein–ligand docking approach. The attained models have showed a good extrapolative power for internal as well as external validation. Moreover, MD simulation was carried out to explain the detailed binding mechanism and the comparison of inhibitor’s binding mode with diverse biological activities. A good correlation was observed among docking studies, MD results, and contour map analysis. Interestingly new molecules were designed using detail structural information from MD simulation, showed higher potency of inhibition (pIC₅₀ 7.6–7.7) compare to the most active compound of the series.     

Molecular dynamics simulations of T-20 HIV fusion inhibitor interacting with model membranes

T-20 (also known as enfuvirtide) is a fusion inhibitor peptide known to have some effectiveness in the control of progression of HIV infection by inhibiting the fusion of the HIV envelope with the target cell membrane. Recent results indicate that T-20 is able to interact with membranes in the liquid disordered state but not with membranes in an ordered state, which could be linked to its effectiveness. A detailed molecular picture of the interaction of these molecules with membranes is still lacking. To this effect, extensive molecular dynamics simulations (100 ns) were carried out to investigate the interaction between T-20 and bilayers of 1-palmitoyl-2-oleoyl-phosphatidylcholine (POPC) and POPC/cholesterol (1:1). Membrane properties such as area/lipid, density profiles, order parameters and membrane thickness were studied. It was observed that T-20 has the ability to interact to different extents with both model membranes in this study and that peptide interaction with the bilayer surface has a local effect on membrane structure. The formation of hydrogen bonding between certain peptide residues and the POPC phosphate group was observed. However, T-20 showed a more limited extent of interaction with model membranes when compared with other, more efficient, peptides (such as T-1249). This effect is most notable in POPC/Chol membranes in which interaction is especially weak, owing to less peptide residues acting as H bond donors to POPC and virtually no H bonds being formed between T-20 and cholesterol. This lower ability to interact with membranes is probably correlated with its smaller inhibitory efficiency.     

Molecular dynamics simulations on SDF-1alpha: binding with CXCR4 receptor

Insights into the interacting mode of CXCR4 with SDF-1alpha are crucial in understanding the structural and functional characteristics of CXCR4 receptor. In this paper a computational pipeline, integrating protein structure prediction, molecular dynamics simulations, automated molecular docking, and Brownian dynamics simulations were employed to investigate the dynamic and energetic aspects of CXCR4 associating with SDF-1alpha. The entire simulation revealed the surface distribution feature of electrostatic potentials and conformational "open-close" process of the receptor. The possible binding conformation of CXCR4 was identified, and the CXCR4-SDF-1alpha binding complex was generated. Arg188-Glu277 salt bridge plays an important role for both the extracellular domain conformational change and SDF-1alpha binding. Two binding sites were mapped at the extracellular domain (Site 1) and inside the transmembrane domain (Site 2), which are composed of conserved residues. Sites 1 and 2 contribute approximately 60% and 40% to the binding affinity with SDF-1alpha, respectively. The binding model is in agreement with most of the experimental data. Transmembrane VI has more significant motion in the harmonious conformational transition of CXCR4 during SDF-1alpha binding, which may be possibly associated with signal transduction. Based on the modeling and simulation, a binding mechanism hypothesis between CXCR4 and SDF-1alpha and its relationship to the signal transduction has been proposed.     

Molecular dynamics simulations of lipid nanodiscs

A lipid nanodisc is a discoidal lipid bilayer stabilized by proteins, peptides, or polymers on its edge. Nanodiscs have two important connections to structural biology. The first is associated with high-density lipoprotein (HDL), a particle with a variety of functionalities including lipid transport. Nascent HDL (nHDL) is a nanodisc stabilized by Apolipoprotein A-I (APOA1). Determining the structure of APOA1 and its mimetic peptides in nanodiscs is crucial to understanding pathologies related to HDL maturation and designing effective therapies. Secondly, nanodiscs offer non-detergent membrane-mimicking environments and greatly facilitate structural studies of membrane proteins. Although seemingly similar, natural and synthetic nanodiscs are different in that nHDL is heterogeneous in size, due to APOA1 elasticity, and gradually matures to become spherical. Synthetic nanodiscs, in contrast, should be homogenous, stable, and size-tunable. This report reviews previous molecular dynamics (MD) simulation studies of nanodiscs and illustrates convergence and accuracy issues using results from new multi-microsecond atomistic MD simulations. These new simulations reveal that APOA1 helices take 10–20 μs to rearrange on the nanodisc, while peptides take 2 μs to migrate from the disc surfaces to the edge. These systems can also become kinetically trapped depending on the initial conditions. For example, APOA1 was trapped in a biologically irrelevant conformation for the duration of a 10 μs trajectory; the peptides were similarly trapped for 5 μs. It therefore remains essential to validate MD simulations of these systems with experiments due to convergence and accuracy issues. This article is part of a Special Issue entitled: Emergence of Complex Behavior in Biomembranes edited by Marjorie Longo.     

Molecular dynamics simulations reveal structural insights into inhibitor binding modes and functionality in human Group IIA phospholipase A2

Human Group IIA phospholipase A₂ (hGIIA) promotes inflammation in immune‐mediated pathologies by regulating the arachidonic acid pathway through both catalysis‐dependent and ‐independent mechanisms. The hGIIA crystal structure, both alone and inhibitor‐bound, together with structures of closely related snake‐venom‐derived secreted phospholipase enzymes has been well described. However, differentiation of biological and nonbiological contacts and the relevance of structures determined from snake venom enzymes to human enzymes are not clear. We employed molecular dynamics (MD) and docking approaches to understand the binding of inhibitors that selectively or nonselectively block the catalysis‐independent mechanism of hGIIA. Our results indicate that hGIIA behaves as a monomer in the solution environment rather than a dimer arrangement that is in the asymmetric unit of some crystal structures. The binding mode of a nonselective inhibitor, KH064, was validated by a combination of the experimental electron density and MD simulations. The binding mode of the selective pentapeptide inhibitor FLSYK to hGIIA was stipulated to be different to that of the snake venom phospholipases A₂ of Daboia russelli pulchella (svPLA₂). Our data suggest that the application of MD approaches to crystal structure data is beneficial in evaluating the robustness of conclusions drawn based on crystal structure data alone. Proteins 2017; 85:827–842. © 2016 Wiley Periodicals, Inc.     

An alamethicin channel in a lipid bilayer: molecular dynamics simulations

We present the results of 2-ns molecular dynamics (MD) simulations of a hexameric bundle of Alm helices in a 1-palmitoyl-2-oleoylphosphatidylcholine bilayer. These simulations explore the dynamic properties of a model of a helix bundle channel in a complete phospholipid bilayer in an aqueous environment. We explore the stability and conformational dynamics of the bundle in a phospholipid bilayer. We also investigate the effect on bundle stability of the ionization state of the ring of Glu18 side chains. If all of the Glu18 side chains are ionised, the bundle is unstable; if none of the Glu18 side chains are ionized, the bundle is stable. pKA calculations suggest that either zero or one ionized Glu18 is present at neutral pH, correlating with the stable form of the helix bundle. The structural and dynamic properties of water in this model channel were examined. As in earlier in vacuo simulations (Breed et al., 1996 .Biophys. J. 70:1643-1661), the dipole moments of water molecules within the pore were aligned antiparallel to the helix dipoles. This contributes to the stability of the helix bundle.     

Molecular dynamics simulations of indolicidin association with model lipid bilayers

Identifying the mechanisms responsible for the interaction of peptides with cell membranes is critical to the design of new antimicrobial peptides and membrane transporters. We report here the results of a computational simulation of the interaction of the 13-residue peptide indolicidin with single-phase lipid bilayers of dioleoylphosphatidylcholine, distearoylphosphatidylcholine, dioleoylphosphatidylglycerol, and distearoylphosphatidylglycerol. Ensemble analysis of the membrane-bound peptide revealed that, in contrast to the extended, linear backbone structure reported for indolicidin in sodium dodecyl sulphate detergent micelles, the peptide adopts a boat-shaped conformation in both phosphatidylglycerol and phosphatidylcholine lipid bilayers, similar to that reported for dodecylphosphocholine micelles. In agreement with fluorescence and NMR experiments, simulations confirmed that the peptide localizes in the membrane interface, with the distance between phosphate headgroups of each leaflet being reduced in the presence of indolicidin. These data, along with a concomitant decrease in lipid order parameters for the upper-tail region, suggest that indolicidin binding results in membrane thinning, consistent with recent in situ atomic force microscopy studies.     

Molecular dynamics simulations of the adipocyte lipid binding protein reveal a novel entry site for the ligand

The adipocyte lipid binding protein (ALBP) binds fatty acids (FA) in a cavity that is inaccessible from the bulk. Therefore, the penetration of the FA necessitates conformational changes whose nature is still unknown. It was suggested that the lipid first enters through a "portal region" which consists of the alphaII helix and the adjacent tight turns. The initial event in the ligand binding must be the interaction of the lipid with the protein surface. To analyze this interaction, we have carried out three molecular dynamics simulations of the apo-ALBP, with a palmitate ion initially located at different positions near the protein surface. The simulation indicated that the ligand could adsorb to the protein in more than one location. Yet, in one case, the ligand managed to penetrate the protein by entering a newly formed cavity some 10 A deep. The entry site is located near the N-terminus, at the junction between the loops connecting the beta-strands. Further progression of the penetration seems to be arrested by the need for desolvation of the COOH end of the headgroup. Evolutionary analysis showed that amino acids in this entry site are well conserved. On the basis of these observations, we suggest that the ligand may enter the protein from a site other than the portal region. Furthermore, the rate-limiting step is proposed to be the desolvation of the FA polar headgroup.     

Molecular dynamics simulations of local anesthetic articaine in a lipid bilayer

In order to investigate structural and dynamical properties of local anesthetic articaine in a model lipid bilayer, a series of molecular dynamics simulations have been performed. Simulations were carried out for neutral and charged (protonated) forms of articaine inserted in fully hydrated dimyristoylphosphatidylcholine (DMPC) lipid bilayer. For comparison purpose, a fully hydrated DMPC bilayer without articaine was also simulated. The length of each simulation was 200 ns. Various properties of the lipid bilayer systems in the presence of both charged and uncharged forms of articaine taken at two different concentrations have been examined: membrane area per lipid, mass density distributions, order parameters, radial distribution functions, head group tilt, diffusion coefficients, electrostatic potential, etc, and compared with results of previous simulations of DMPC bilayer in the presence of lidocaine. It was shown that addition of both charged and neutral forms of articaine causes increase of the dipole electrostatic potential in the membrane interior.     

Molecular dynamics simulations of unsaturated lipid bilayers: effects of varying the numbers of double bonds

The importance of unsaturated, and especially polyunsaturated phosphatidylcholine molecules for the functional properties of biological membranes is widely accepted. Here, the effects of unsaturation on the nanosecond-scale structural and dynamic properties of the phosphatidylcholine bilayer were elucidated by performance of multinanosecond molecular dynamics simulations of all-atom bilayer models. Bilayers of dipalmitoylphosphatidylcholine and its mono-, di-, and tetraunsaturated counterparts were simulated, containing, respectively, oleoyl, linoleoyl, or arachidonoyl chains in the sn-2 position. Analysis of the simulations focused on comparison of the structural properties, especially the ordering of the chains in the membranes. Although the results suggest some problems in the CHARMM force field of the lipids when applied in a constant pressure ensemble, the features appearing in the ordering of the unsaturated chains are consistent with the behaviour known from ²H NMR experiments. The rigidity of the double bonds is compensated by the flexibility of skew state single bonds juxtaposed with double bonds. The presence of double bonds in the sn-2 chains considerably reduces the order parameters of the CH bonds. Moreover, the double bond region of tetraunsaturated chains is shown to span all the way from the bilayer centre to the head group region.     

Molecular dynamics simulations of N-terminal peptides from a nucleotide binding protein

Molecular dynamics (MD) simulations of N-terminal peptides from lactate dehydrogenase (LDH) with increasing length and individual secondary structure elements were used to study their stability in relation to folding. Ten simulations of 1–2 ns of different peptides in water starting from the coordinates of the crystal structure were performed. The stability of the peptides was compared qualitatively by analyzing the root mean square deviation (RMSD) from the crystal structure, radius of gyration, secondary and tertiary structure, and solvent accessible surface area. In agreement with earlier MD studies, relatively short (< 15 amino acids) peptides containing individual secondary structure elements were generally found to be unstable; the hydrophobic α₁-helix of the nucleotide binding fold displayed a significantly higher stability, however. Our simulations further showed that the first βαβ supersecondary unit of the characteristic dinucleotide binding fold (Rossmann fold) of LDH is somewhat more stable than other units of similar length and that the α₂-helix, which unfolds by itself, is stabilized by binding to this unit. This finding suggests that the first βαβ unit could function as an N-terminal folding nucleus, upon which the remainder of the polypeptide chain can be assembled. Indeed, simulations with longer units (βαβα and βαβαββ) showed that all structural elements of these units are rather stable. The outcome of our studies is in line with suggestions that folding of the N-terminal portion of LDH in vivo can be a cotranslational process that takes place during the ribosomal peptide synthesis.     

Molecular dynamics simulations of the effects of ring-saturated thymine lesions on DNA structure

The effect of thymine lesions produced by radiation or oxidative damage on DNA structure was studied by molecular dynamics simulations of native and damaged DNA. Thymine in position 7 of native dodecamer d(CGCGAATTCGCG)₂ was replaced by one of the four thymine lesions 5-hydroxy-5,6-dihydrothymine, 6-hydroxy-5,6-dihydrothymine (thymine photohydrate), 5,6-dihydmxy-5,6-dihydro-thymine (thymine glycol), and 5,6-dihydmthymine. Simulations were performed with Assisted Model Building with Energy Refinement force field. Solvent was represented by a rectangular box of water with periodic boundary conditions applied. A constant temperature and constant volume protocol was used, the observed level of distortions of DNA structure depends on the specific nature of the lesion. The 5,6-dihydrothymine does not cause distinguishable perturbations to DNA. Other lesions produce a dramatic increase in the rise parameter between the lesion and the 5′ adjacent adenine. These changes are accompanied by weakening of Watson–Crick hydrogen bonds in the A6-T19 base pair on the 5′ side of the lesion. The lesioned bases also show negative values of inclination relative to the helical axis. No changes in the pattern of backbone torsional angles are observed with any of the lesions incorporated into DNA. The structural distortions in DNA correlate well with known biological effects of 5,6-dihydrothymine and thymine glycol on such processes as polymerase action or recognition by repair enzymes. © 1995 John Wiley & Sons, Inc.     

Molecular dynamics simulations of stratum corneum lipid models: fatty acids and cholesterol

We report the results of an investigation on stratum corneum lipids, which present the main barrier of the skin. Molecular dynamics simulations, thermal analysis and FTIR measurements were applied. The primary objective of this work was to study the effect of cholesterol on skin structure and dynamics. Two molecular models were constructed, a free fatty acid bilayer (stearic acid, palmitic acid) and a fatty acid/cholesterol mixture at a 1:1 molar ratio. Our simulations were performed at constant pressure and temperature on a nanosecond time scale. The resulting model structures were characterized by calculating surface areas per headgroup, conformational properties, atom densities and order parameters of the fatty acids. Analysis of the simulations indicates that the free fatty acid fraction of stratum corneum lipids stays in a highly ordered crystalline state at skin temperatures. The phase behavior is strongly influenced when cholesterol is added. Cholesterol smoothes the rigid phases of the fatty acids: the order of the hydrocarbon tails (mainly of the last eight bonds) is reduced, the area per molecule becomes larger, the fraction of trans dihedrals is lower and the hydrophobic thickness is reduced. The simulation results are in good agreement with our experimental data from FTIR analysis and NIR-FT Raman spectroscopy.     

Molecular dynamics simulations of lipid vesicle fusion in atomic detail

The fusion of a membrane-bounded vesicle with a target membrane is a key step in intracellular trafficking, exocytosis, and drug delivery. Molecular dynamics simulations have been used to study the fusion of small unilamellar vesicles composed of a dipalmitoyl-phosphatidylcholine (DPPC)/palmitic acid 1:2 mixture in atomic detail. The simulations were performed at 350-370 K and mimicked the temperature- and pH-induced fusion of DPPC/palmitic acid vesicles from experiments by others. To make the calculations computationally feasible, a vesicle simulated at periodic boundary conditions was fused with its periodic image. Starting from a preformed stalk between the outer leaflets of the vesicle and its periodic image, a hemifused state formed within 2 ns. In one out of six simulations, a transient pore formed close to the stalk, resulting in the mixing of DPPC lipids between the outer and the inner leaflet. The hemifused state was (meta)stable on a timescale of up to 11 ns. Forcing a single lipid into the interior of the hemifusion diaphragm induced the formation and expansion of a fusion pore on a nanosecond timescale. This work opens the perspective to study a wide variety of mesoscopic biological processes in atomic detail.     

Molecular dynamics simulations of a pulmonary surfactant protein B peptide in a lipid monolayer

Pulmonary surfactant is a complex mixture of lipids and proteins that lines the air/liquid interface of the alveolar hypophase and confers mechanical stability to the alveoli during the breathing process. The desire to formulate synthetic mixtures for low-cost prophylactic and therapeutic applications has motivated the study of the specific roles and interactions of the different components. All-atom molecular dynamics simulations were carried out on a model system composed of a monolayer of palmitic acid (PA) and a surfactant protein B peptide, SP-B(1-25). A detailed structural characterization as a function of the lipid monolayer specific area revealed that the peptide remains inserted in the monolayer up to values of specific area corresponding to an untilted condensed phase of the the pure palmitic acid monolayer. The system remains stable by altering the conformational order of both the anionic lipid monolayer and the peptide secondary structure. Two elements appear to be key for the constitution of this phase: an electrostatic interaction between the cationic peptide residues with the anionic headgroups, and an exclusion of the aromatic residues on the hydrophobic end of the peptide from the hydrophilic and aqueous regions.     

Molecular dynamics simulations of charged and neutral lipid bilayers: treatment of electrostatic interactions

Molecular dynamics (MD) simulations complement experimental methods in studies of the structure and dynamics of lipid bilayers. The choice of algorithms employed in this computational method represents a trade-off between the accuracy and real calculation time. The largest portion of the simulation time is devoted to calculation of long-range electrostatic interactions. To speed-up evaluation of these interactions, various approximations have been used. The most common ones are the truncation of long-range interactions with the use of cut-offs, and the particle-mesh Ewald (PME) method. In this study, several multi-nanosecond cut-off and PME simulations were performed to establish the influence of the simulation protocol on the bilayer properties. Two bilayers were used. One consisted of neutral phosphatidylcholine molecules. The other was a mixed lipid bilayer consisting of neutral phosphatidylethanolamine and negatively charged phosphatidylglycerol molecules. The study shows that the cut-off simulation of a bilayer containing charge molecules generates artefacts; in particular the mobility and order of the charged molecules are vastly different from those determined experimentally. In the PME simulation, the bilayer properties are in general agreement with experimental data. The cut-off simulation of bilayers containing only uncharged molecules does not generate artefacts, nevertheless, the PME simulation gives generally better agreement with experimental data.     

Molecular dynamics study of aquaporin-1 water channel in a lipid bilayer

The aquaporin-1 water channel was modeled in a palmitoyl-oleoyl-phosphatidyl-choline lipid bilayer, by means of molecular dynamics simulations. Interaction of the protein with the membrane and inter-monomer interactions were analyzed. Structural features of the channel important for its biological function, including the Asn-Pro-Ala (NPA) motifs, and the diffusion of water molecules into the channels, were investigated. Simulations revealed the formation of single file water inside the channels for certain relative positions of the NPA motifs.