Conformational dynamics of α-1 acid glycoprotein (AGP) in cancer: A comparative study of glycosylated and unglycosylated AGP.

α-1 acid glycoprotein (AGP) is one of the most abundant plasma proteins. It fulfills two important functions: immunomodulation, and binding to various drugs and receptors. These different functions are closely associated and modulated via changes in glycosylation and cancer missense mutations. From a structural point of view, glycans alter the local biophysical properties of the protein leading to a diverse ligand-binding spectrum. However, glycans can typically not be observed in the resolved X-ray crystallography structure of AGP due to their high flexibility and microheterogeneity, so limiting our understanding of AGP's conformational dynamics 70 years after its discovery. We here investigate how mutations and glycosylation interfere with AGP's conformational dynamics changing its biophysical behavior, by using molecular dynamics (MD) simulations and sequence-based dynamics predictions. The MD trajectories show that glycosylation decreases the local backbone flexibility of AGP and increases the flexibility of distant regions through allosteric effects. We observe that mutations near the glycosylation site affect glycan's conformational preferences. Thus, we conclude that mutations control glycan dynamics which modulates the protein's backbone flexibility directly affecting its accessibility. These findings may assist in the drug design targeting AGP's glycosylation and mutations in cancer.

PyRETIS 3: Conquering rare and slow events without boundaries.

We present and discuss the advancements made in PyRETIS 3, the third instalment of our Python library for an efficient and user-friendly rare event simulation, focused to execute molecular simulations with replica exchange transition interface sampling (RETIS) and its variations. Apart from a general rewiring of the internal code towards a more modular structure, several recently developed sampling strategies have been implemented. These include recently developed Monte Carlo moves to increase path decorrelation and convergence rate, and new ensemble definitions to handle the challenges of long-lived metastable states and transitions with unbounded reactant and product states. Additionally, the post-analysis software PyVisa is now embedded in the main code, allowing fast use of machine-learning algorithms for clustering and visualising collective variables in the simulation data.

Defining permeability of curved membranes in molecular dynamics simulations.

Many phospholipid membranes in the cell have a high curvature; for instance, in caveolae, mitochondrial crystae, nanotubes, membrane pearls, small liposomes, or exosomes. Molecular dynamics (MD) simulations are a computational tool to gain insight in the transport behavior at the atomic scale. Membrane permeability is a key kinetic property that might be affected in these highly curved membranes. Unfortunately, the geometry of highly curved membranes creates ambiguity in the permeability value, even with an arbitrarily large factor purely based on geometry, caused by the radial flux not being a constant value in steady state. In this contribution, the ambiguity in permeability for liposomes is countered by providing a new permeability definition. First, the inhomogeneous solubility diffusion model based on the Smoluchowski equation is solved analytically under radial symmetry, from which the entrance and escape permeabilities are defined. Next, the liposome permeability is defined guided by the criterion that a flat and curved membrane should have equal permeability, in case these were to be carved out from an imaginary homogeneous medium. With this criterion, our new definition allows for a fair comparison of flat and curved membranes. The definition is then transferred to the counting method, which is a practical computational approach to derive permeability by counting complete membrane crossings. Finally, the usability of the approach is illustrated with MD simulations of diphosphatidylcholine (DPPC) bilayers, without or with some cholesterol content. Our new liposome permeability definition allows us to compare a spherically shaped membrane with its flat counterpart, thus showcasing how the curvature effect on membrane transport may be assessed.

Path sampling with memory reduction and replica exchange to reach long permeation timescales.

Assessing kinetics in biological processes with molecular dynamics simulations remains a computational and conceptual challenge, given the large time and length scales involved. For kinetic transport of biochemical compounds or drug molecules, the permeability through the phospholipid membranes is a key kinetic property, but long timescales are hindering the accurate computation. Technological advances in high-performance computing therefore need to be accompanied by theoretical and methodological developments. In this contribution, the replica exchange transition interface sampling (RETIS) methodology is shown to give perspective toward observing longer permeation pathways. It is first reviewed how RETIS, a path-sampling methodology that gives in principle exact kinetics, can be used to compute membrane permeability. Next, recent and current developments in three RETIS aspects are discussed: several new Monte Carlo moves in the path-sampling algorithm, memory reduction by reducing pathlengths, and exploitation of parallel computing with CPU-imbalanced replicas. Finally, the memory reduction presenting a new replica exchange implementation, coined REPPTIS, is showcased with a permeant needing to pass a membrane with two permeation channels, either representing an entropic or energetic barrier. The REPPTIS results showed clearly that inclusion of some memory and enhancing ergodic sampling via replica exchange moves are both necessary to obtain correct permeability estimates. In an additional example, ibuprofen permeation through a dipalmitoylphosphatidylcholine membrane was modeled. REPPTIS succeeded in estimating the permeability of this amphiphilic drug molecule with metastable states along the permeation pathway. In conclusion, the presented methodological advances allow for deeper insight into membrane biophysics even if the pathways are slow, as RETIS and REPPTIS push the permeability calculations to longer timescales.

Capric Acid and Myristic Acid Permeability Enhancers in Curved Liposome Membranes.

Liposomes are considered as advanced drug delivery systems for cancer treatment. A generation of pH-sensitive liposomes is being developed that use fatty acids (FAs) as a trigger for drug release in tumor tissues. However, FAs are also known to enhance permeability, and it is unclear whether FAs in liposomes may cause drug leakage or premature drug release. The passive permeability of the drug through the membrane of the liposome is thus a crucial factor for timely drug delivery. To investigate how the curvature and lipid composition of liposomes affect their passive permeability, coarse-grained molecular dynamics were performed. The permeability was determined with a counting method. Flat bilayers and three liposomes with varying diameters were studied, which had varying lipid compositions of dipalmitoylphosphatidylcholine, cholesterol, and deprotonated or neutral saturated FAs. The investigated permeants were water and two other small permeants, which have different free energy profiles (solubility) across the membrane. First, for the curvature effect, our results showed that curvature increases the water permeability by reducing the membrane thickness. The permeability increase for water is about a factor of 1.7 for the most curved membranes. However, a high curvature decreases permeability for permeants with free energy profiles that are a mix of wells and barriers in the headgroup region of the membrane. Importantly, the type of experimental setup is expected to play a dominant role in the permeability value, i.e., whether permeants are escaping or entering the liposomes. Second, for the composition effect, FAs decrease both the area per lipid (APL) and the membrane thickness, resulting in permeability increases of up to 55%. Cholesterol has a similar effect on the APL but has the opposite impact on membrane thickness and permeability. Therefore, FAs and cholesterol have opposing effects on permeability, with cholesterol's effect being slightly stronger in our simulated bilayers. As all permeability values were well within a factor of 2, and with liposomes usually being larger and less curved in experimental applications, it can be concluded that the passive drug release from a pH-sensitive liposome does not seem to be significantly affected by the presence of FAs.

Understanding the Role of Caveolae in Oxygen Buffering: The Effect of Membrane Curvature.

The "oxygen paradox" can be explained as two opposing biological processes with oxygen (O2) as a reactant. On the one hand, oxygen is essential to aerobic metabolism, powering oxidative phosphorylation in mitochondria. On the other hand, an excess supply of oxygen will generate reactive species which are harmful for the cell. In healthy tissues, the first process must be maximized relative to the second one. We have hypothesized that curved and cholesterol-enriched membrane invaginations called caveolae help maintain the proper oxygen level by taking up oxygen and attenuating its release to the mitochondria. The mechanism by which caveolae may help to buffer the oxygen level in cells is still unclear. Here, we aim to assess how structural aspects of caveolae, the curvature of the membrane, influence the local oxygen abundance and the membrane partitioning. We have modelled a flat bilayer and a liposome composed of dipalmitoylphosphatidylcholine (DPPC), using molecular dynamics simulation. Associated changes in the membrane-level oxygen partition coefficient and free energy profiles will be presented.

Oxygen Storage in Stacked Phospholipid Membranes Under an Oxygen Gradient as a Model for Myelin Sheaths.

Axons in the brain and peripheral nervous system are enveloped by myelin sheaths, which are composed of stacked membrane bilayers containing large fractions of cholesterol, phospholipids, and glycolipids. The oxygen availability to the nearby oxygen consuming cytochrome c oxidase in the mitochondria is essential for the well-functioning of a cell. By constructing a rate network model based on molecular dynamics simulations, and solving it for steady-state conditions, this work calculates the oxygen storage in stacked membranes under an oxygen gradient. It is found that stacking membranes increases the oxygen storage capacity, indicating that myelin can function as an oxygen reservoir. However, it is found that the storage enhancement levels out for stacks with a large number of bilayers, suggesting why myelin sheaths consist of only 10-300 membranes rather than thousands. The presence of additional water between the stacked bilayers, as seen in cancer cells, is shown to diminish myelin oxygen storage enhancement.

Sampling efficiency of the counting method for permeability calculations estimated with the inhomogeneous solubility-diffusion model.

Permeability is a key property in various fields such as membrane technology for chemical separation and transport of substances through cell membranes. At the molecular scale, the counting method uses the number of membrane crossings in a conventional unbiased molecular dynamics simulation to predict the permeability. This contribution investigates under which conditions the counting method has insufficient statistics. An equation is derived for a compartmental model based on the inhomogeneous solubility-diffusion (Smoluchowski) model, giving insight into how the flux correlates with the solubility of permeants. This equation shows that a membrane crossing is a rare event not only when the membrane forms a large free energy barrier but also when the membrane forms a deep free energy well that traps permeants. Such a permeant trap has a high permeability; yet, the counting method suffers from poor statistics. To illustrate this, coarse-grained MD was run for 16 systems of dipalmitoylphosphatidylcholine bilayer membranes with different permeant types. The composition rule for permeability is shown to also hold for fluxes, and it is highlighted that the considered thickness of the membrane causes uncertainty in the permeability calculation of highly permeable membranes. In conclusion, a high permeability in itself is not an effective indicator of the sampling efficiency of the counting method, and caution should be taken for permeants whose solubility varies greatly over the simulation box. A practical consequence relevant in, e.g., drug design is that a drug with high membrane permeability might get trapped by membranes thus reducing its efficacy.

Membrane permeability of small molecules from unbiased molecular dynamics simulations.

Permeation of many small molecules through lipid bilayers can be directly observed in molecular dynamics simulations on the nano- and microsecond timescale. While unbiased simulations provide an unobstructed view of the permeation process, their feasibility for computing permeability coefficients depends on various factors that differ for each permeant. The present work studies three small molecules for which unbiased simulations of permeation are feasible within less than a microsecond, one hydrophobic (oxygen), one hydrophilic (water), and one amphiphilic (ethanol). Permeabilities are computed using two approaches: counting methods and a maximum-likelihood estimation for the inhomogeneous solubility diffusion (ISD) model. Counting methods yield nearly model-free estimates of the permeability for all three permeants. While the ISD-based approach is reasonable for oxygen, it lacks precision for water due to insufficient sampling and results in misleading estimates for ethanol due to invalid model assumptions. It is also demonstrated that simulations using a Langevin thermostat with collision frequencies of 1/ps and 5/ps yield oxygen permeabilities and diffusion constants that are lower than those using Nosé-Hoover by statistically significant margins. In contrast, permeabilities from trajectories generated with Nosé-Hoover and the microcanonical ensemble do not show statistically significant differences. As molecular simulations become more affordable and accurate, calculation of permeability for an expanding range of molecules will be feasible using unbiased simulations. The present work summarizes theoretical underpinnings, identifies pitfalls, and develops best practices for such simulations.

Effect of Chain Unsaturation and Temperature on Oxygen Diffusion Through Lipid Membranes from Simulations.

Rafts are nanoscale ordered domains in biological membranes that are rich in saturated phospholipids. In this study, the influence of chain unsaturation and temperature on oxygen diffusion through lipid membranes is examined using advanced computational modeling. The studied phospholipids with increasing unsaturation are: 1,2-dipalmitoyl-sn-glycero-3-phosphocholine (DPPC), 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC), and 1,2-dioleoyl-sn-glycero-3-phosphocholine (DOPC). The unsaturation correlates with the area per lipid and the order parameter. Oxygen diffusion is found to be faster at higher temperature, and the solubility of oxygen in the membrane with respect to water decreases. Diffusion varies over a larger range across the membrane at 323 K in DPPC than in DOPC, whereas POPC has intermediate diffusivity. Oxygen diffusion in saturated lipids is faster at the membrane center and slower near the head group region than in unsaturated lipids. Oxygen solubility in DPPC is higher than in unsaturated lipids.

Complex Reaction Environments and Competing Reaction Mechanisms in Zeolite Catalysis: Insights from Advanced Molecular Dynamics.

The methanol-to-olefin process is a showcase example of complex zeolite-catalyzed chemistry. At real operating conditions, many factors affect the reactivity, such as framework flexibility, adsorption of various guest molecules, and competitive reaction pathways. In this study, the strength of first principle molecular dynamics techniques to capture this complexity is shown by means of two case studies. Firstly, the adsorption behavior of methanol and water in H-SAPO-34 at 350 °C is investigated. Hereby an important degree of framework flexibility and proton mobility was observed. Secondly, the methylation of benzene by methanol through a competitive direct and stepwise pathway in the AFI topology was studied. Both case studies clearly show that a first-principle molecular dynamics approach enables unprecedented insights into zeolite-catalyzed reactions at the nanometer scale to be obtained.

Exploring the vibrational fingerprint of the electronic excitation energy via molecular dynamics.

A Fourier-based method is presented to relate changes of the molecular structure during a molecular dynamics simulation with fluctuations in the electronic excitation energy. The method implies sampling of the ground state potential energy surface. Subsequently, the power spectrum of the velocities is compared with the power spectrum of the excitation energy computed using time-dependent density functional theory. Peaks in both spectra are compared, and motions exhibiting a linear or quadratic behavior can be distinguished. The quadratically active motions are mainly responsible for the changes in the excitation energy and hence cause shifts between the dynamic and static values of the spectral property. Moreover, information about the potential energy surface of various excited states can be obtained. The procedure is illustrated with three case studies. The first electronic excitation is explored in detail and dominant vibrational motions responsible for changes in the excitation energy are identified for ethylene, biphenyl, and hexamethylbenzene. The proposed method is also extended to other low-energy excitations. Finally, the vibrational fingerprint of the excitation energy of a more complex molecule, in particular the azo dye ethyl orange in a water environment, is analyzed.

An exact algorithm to find a maximum weight clique in a weighted undirected graph.

We introduce a new algorithm MaxCliqueWeight for identifying a maximum weight clique in a weighted graph, and its variant MaxCliqueDynWeight with dynamically varying bounds. This algorithm uses an efficient branch-and-bound approach with a new weighted graph coloring algorithm that efficiently determines upper weight bounds for a maximum weighted clique in a graph. We evaluate our algorithm on random weighted graphs with node counts up to 10,000 and on standard DIMACS benchmark graphs used in a variety of research areas. Our findings reveal a remarkable improvement in computational speed when compared to existing algorithms, particularly evident in the case of high-density random graphs and DIMACS graphs, where our newly developed algorithm outperforms existing alternatives by several orders of magnitude. The newly developed algorithm and its variant are freely available to the broader research community at http://insilab.org/maxcliqueweight , paving the way for transformative applications in various research areas, including drug discovery.

Comparing normal modes across different models and scales: Hessian reduction versus coarse-graining.

Dimension reduction is often necessary when attempting to reach longer length and time scales in molecular simulations. It is realized by constraining degrees of freedom or by coarse-graining the system. When evaluating the accuracy of a dimensional reduction, there is a practical challenge: the models yield vectors with different lengths, making a comparison by calculating their dot product impossible. This article investigates mapping procedures for normal mode analysis. We first review a horizontal mapping procedure for the reduced Hessian techniques, which projects out degrees of freedom. We then design a vertical mapping procedure for the "implosion" of the all-atom (AA) Hessian to a coarse-grained scale that is based upon vibrational subsystem analysis. This latter method derives both effective force constants and an effective kinetic tensor. Next, a series of metrics is presented for comparison across different scales, where special attention is given to proper mass-weighting. The dimension-dependent metrics, which require prior mapping for proper evaluation, are frequencies, overlap of normal mode vectors, probability similarity, Hessian similarity, collectivity of modes, and thermal fluctuations. The dimension-independent metrics are shape derivatives, elastic modulus, vibrational free energy differences, heat capacity, and projection on a predefined basis set. The power of these metrics to distinguish between reasonable and unreasonable models is tested on a toy alpha helix system and a globular protein; both are represented at several scales: the AA scale, a Go-like model, a canonical elastic network model, and a network model with intentionally unphysical force constants.

Analysis of the basis set superposition error in molecular dynamics of hydrogen-bonded liquids: application to methanol.

An efficient protocol is presented to compensate for the basis set superposition error (BSSE) in DFT molecular dynamics (MD) simulations using localized Gaussian basis sets. We propose a classical correction term that can be added a posteriori to account for BSSE. It is tested to what extension this term will improve radial distribution functions (RDFs). The proposed term is pairwise between certain atoms in different molecules and was calibrated by fitting reference BSSE data points computed with the counterpoise method. It is verified that the proposed exponential decaying functional form of the model is valid. This work focuses on hydrogen-bonded liquids, i.e., methanol, and more specific on the intermolecular hydrogen bond, but in principle the method is generally applicable on any type of interaction where BSSE is significant. We evaluated the relative importance of the Grimme-dispersion versus BSSE and found that they are of the same order of magnitude, but with an opposite sign. Upon introduction of the correction, the relevant RDFs, obtained from MD, have amplitudes equal to experiment.

First principle kinetic studies of zeolite-catalyzed methylation reactions.

Methylations of ethene, propene, and butene by methanol over the acidic microporous H-ZSM-5 catalyst are studied by means of state of the art computational techniques, to derive Arrhenius plots and rate constants from first principles that can directly be compared with the experimental data. For these key elementary reactions in the methanol to hydrocarbons (MTH) process, direct kinetic data became available only recently [J. Catal.2005, 224, 115-123; J. Catal.2005, 234, 385-400]. At 350 °C, apparent activation energies of 103, 69, and 45 kJ/mol and rate constants of 2.6 × 10(-4), 4.5 × 10(-3), and 1.3 × 10(-2) mol/(g h mbar) for ethene, propene, and butene were derived, giving following relative ratios for methylation k(ethene)/k(propene)/k(butene) = 1:17:50. In this work, rate constants including pre-exponential factors are calculated which give very good agreement with the experimental data: apparent activation energies of 94, 62, and 37 kJ/mol for ethene, propene, and butene are found, and relative ratios of methylation k(ethene)/k(propene)/k(butene) = 1:23:763. The entropies of gas phase alkenes are underestimated in the harmonic oscillator approximation due to the occurrence of internal rotations. These low vibrational modes were substituted by manually constructed partition functions. Overall, the absolute reaction rates can be calculated with near chemical accuracy, and qualitative trends are very well reproduced. In addition, the proposed scheme is computationally very efficient and constitutes significant progress in kinetic modeling of reactions in heterogeneous catalysis.

Permeation Rates of Oxygen through a Lipid Bilayer Using Replica Exchange Transition Interface Sampling.

Several simulations strategies have emerged to predict the permeability of solutes across membranes, which is important for many biological or industrial processes such as drug design. The widespread inhomogeneous solubility-diffusion (ISD) model is based on the Smoluchowski equation and describes permeation as purely diffusive. The counting method, which counts membrane transitions in a long molecular dynamics (MD) trajectory, is free of this diffusive assumption, but it lacks sufficient statistics when the permeation involves high free energy barriers. Metadynamics and variations thereof can overcome such barriers, but they generally lack the kinetics information. The milestoning framework has been used to describe permeation as a rare event, but it still relies on the Markovian assumption between the milestones. Replica Exchange Transition Interface Sampling (RETIS) has been shown to be an effective method for sampling rare events while simultaneously describing the kinetics without assumptions. This paper is the first permeation application of RETIS on an all-atom lipid bilayer consisting of 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) to compute the entrance, escape and complete transition of molecular oxygen. Conventional MD was performed as a benchmark, and the MD rates from counting were converted to rate constants, giving good agreement with the RETIS values. Moreover, a correction factor was derived to convert the collective order parameter in RETIS, which was aimed to improve efficiency, to a single-particle order parameter. With this work, we showed how the exact kinetics of drug molecules permeation can be assessed with RETIS even if the permeation is truly a rare event or if the permeation is non-Markovian. RETIS will therefore be a valuable tool for future permeation studies.

Comparative study of various normal mode analysis techniques based on partial Hessians.

Standard normal mode analysis becomes problematic for complex molecular systems, as a result of both the high computational cost and the excessive amount of information when the full Hessian matrix is used. Several partial Hessian methods have been proposed in the literature, yielding approximate normal modes. These methods aim at reducing the computational load and/or calculating only the relevant normal modes of interest in a specific application. Each method has its own (dis)advantages and application field but guidelines for the most suitable choice are lacking. We have investigated several partial Hessian methods, including the Partial Hessian Vibrational Analysis (PHVA), the Mobile Block Hessian (MBH), and the Vibrational Subsystem Analysis (VSA). In this article, we focus on the benefits and drawbacks of these methods, in terms of the reproduction of localized modes, collective modes, and the performance in partially optimized structures. We find that the PHVA is suitable for describing localized modes, that the MBH not only reproduces localized and global modes but also serves as an analysis tool of the spectrum, and that the VSA is mostly useful for the reproduction of the low frequency spectrum. These guidelines are illustrated with the reproduction of the localized amine-stretch, the spectrum of quinine and a bis-cinchona derivative, and the low frequency modes of the LAO binding protein.

TAMkin: a versatile package for vibrational analysis and chemical kinetics.

TAMkin is a program for the calculation and analysis of normal modes, thermochemical properties and chemical reaction rates. At present, the output from the frequently applied software programs ADF, CHARMM, CPMD, CP2K, Gaussian, Q-Chem, and VASP can be analyzed. The normal-mode analysis can be performed using a broad variety of advanced models, including the standard full Hessian, the Mobile Block Hessian, the Partial Hessian Vibrational approach, the Vibrational Subsystem Analysis with or without mass matrix correction, the Elastic Network Model, and other combinations. TAMkin is readily extensible because of its modular structure. Chemical kinetics of unimolecular and bimolecular reactions can be analyzed in a straightforward way using conventional transition state theory, including tunneling corrections and internal rotor refinements. A sensitivity analysis can also be performed, providing important insight into the theoretical error margins on the kinetic parameters. Two extensive examples demonstrate the capabilities of TAMkin: the conformational change of the biological system adenylate kinase is studied, as well as the reaction kinetics of the addition of ethene to the ethyl radical. The important feature of batch processing large amounts of data is highlighted by performing an extended level of theory study, which TAMkin can automate significantly.