Ad hoc continuum-atomistic thermostat for modeling heat flow in molecular dynamics simulations

An ad hoc thermostating procedure that couples a molecular dynamics (MD) simulation and a numerical solution to the continuum heat flow equation is presented. The method allows experimental thermal transport properties to be modeled without explicitly including electronic degrees of freedom in a MD simulation. The method is demonstrated using two examples, heat flow from a constant temperature silver surface into a single crystal bulk, and a tip sliding along a silver surface. For the former it is shown that frictional forces based on the Hoover thermostat applied locally to grid regions of the simulation are needed for effective feedback between the atomistic and continuum equations. For fast tip sliding the thermostat results in less surface heating, and higher frictional and normal forces compared to the same simulation without the thermostat.     

Molecular dynamics studies in nanoscale liquid structures: geometry and thermal effects on nanojet development

Conventional macroscopic jet theory relies heavily on experimental correlations which cannot be easily extended to the nanoscale regime. Moreover, the fluid dynamic effects at small length scales and their contribution to the development of nanoscale liquid structures are fundamentally different from their macroscopic counterparts. This coupled with the high spatial and temporal resolution requirements at nanoscale domains make molecular dynamics (MD) an excellent tool for studying such structures. In this study, the formation and breakup of nanojets (NJs) developing from high pressure into vacuum is investigated using MD based on non-Hamiltonian formulations. By ejecting the equilibrated argon atoms through various nozzle geometries and diameters, nanoscale jet flows were generated. The dependence of the jet structure on nozzle geometry and diameter is studied. The influence of geometry on NJ formation is also studied along with issues involved in the equilibration and thermostat coupling parameter. Various thermostats are compared to understand the role they play in MD simulations of liquid nanostructures. Tuning of the thermostat coupling parameter has also been discussed. The jet breakup phenomenon is analysed and a comparative study, vis-à-vis, well-established continuum and stochastic models, is attempted.     

Parametric dependence on shear viscosity of SPC/E water from equilibrium and non-equilibrium molecular dynamics simulations

A parametric dependent study is crucial for the accurate determination of transport coefficients such as shear viscosity. In this study, we calculate the shear viscosity of extended simple point charge water using a transverse current auto-correlation function (TCAF) from equilibrium molecular dynamics (EMD) and the periodic perturbation method from non-equilibrium molecular dynamics (NEMD) simulations for varying coupling time and system sizes. Results show that the shear viscosity calculated using EMD simulations with different thermostats varies significantly with coupling times and system size. The use of Berendsen and velocity-rescale thermostats in NEMD simulations generates a significant drift from the target temperature and results in an inconsistent shear viscosity with coupling time and system size. The use of Nosé–Hoover thermostat in NEMD simulations offers thermodynamic stability which results in a consistent shear viscosity for various coupling times and system sizes.     

TRANSLATIONAL MOBILITY OF THE MEMBRANE INTERCALATED PARTICLES OF HUMAN ERYTHROCYTE GHOSTS : pH-Dependent, Reversible Aggregation

This paper demonstrates the translational movement along the plane of the human erythrocyte ghost of the membrane particles exposed by freeze-fracture. The membrane particles can be aggregated by incubation of the ghosts in media with a pH in the vicinity of 5 5 or 3 5. The particles are disaggregated in neutral and alkaline media (pH 9 5) and also at pH 4.5 Aggregation of the particles at pH 5.5 is reversible, prevented by prefixation in glutaraldehyde and by media of high ionic strength. Particle aggregation occurs within 2–4 min. These results are consistent with the concept that the erythrocyte ghost membrane is a planar fluid domain formed by a bilayer membrane continuum which is interrupted by localized, yet mobile, proteic intercalations.     

Effect of the thermostat in the molecular dynamics simulation on the folding of the model protein chignolin

Molecular dynamics simulations of the model protein chignolin with explicit solvent were carried out, in order to analyze the influence of the Berendsen thermostat on the evolution and folding of the peptide. The dependence of the peptide behavior on temperature was tested with the commonly employed thermostat scheme consisting of one thermostat for the protein and another for the solvent. The thermostat coupling time of the protein was increased to infinity, when the protein is not in direct contact with the thermal bath, a situation known as minimally invasive thermostat. In agreement with other works, it was observed that only in the last situation the instantaneous temperature of the model protein obeys a canonical distribution. As for the folding studies, it was shown that, in the applications of the commonly utilized thermostat schemes, the systems are trapped in local minima regions from which it has difficulty escaping. With the minimally invasive thermostat the time that the protein needs to fold was reduced by two to three times. These results show that the obstacles to the evolution of the extended peptide to the folded structure can be overcome when the temperature of the peptide is not directly controlled.     

Analysis of energy and friction coefficient fluctuations of a Lennard-Jones liquid coupled to the Nosé–Hoover thermostat

Thermal fluctuations of a many-body system coupled to a Nosé–Hoover thermostat depend on the strength of the coupling parameter τ*. A wrong choice may bring non-ergodic features and non-canonical fluctuations. Here, we analyse by means of molecular dynamics simulations both the energy fluctuations and the spectrum of the friction coefficient ζ^* of an extended Lennard-Jones system in a wide range of τ* at liquid density. Three ranges of τ* are identified – small, intermediate and large, plus their transitions. For τ* values in the intermediate range, ζ^* shows chaotic behaviour, and the particle system requires reasonable computing time for its thermalisation. As a result, the extended system is ergodic and energy fluctuations are canonical and stable in time. On the contrary, small and large ranges of τ* reveal clear evidence of periodicity in the thermostat variable, for instance by propagating the initial temperature condition. For these two ranges, the extended system shows non-ergodic features and energy fluctuations are non-canonical. For large τ*, micro-canonical fluctuations are occasionally obtained. For small to intermediate and intermediate to large ranges of τ*, the thermostat variable exhibits beat waves and is thus unable to reach equilibrium no matter how extended in time the simulations are. Here, we compare our results with previous work and explain the differences.     

Domain Expansion and Functional Diversification in Vertebrate Reproductive Proteins

The rapid evolution of fertilization proteins has generated remarkable diversity in molecular structure and function. Glycoproteins of vertebrate egg coats contain multiple zona pellucida (ZP)-N domains (1–6 copies) that facilitate multiple reproductive functions, including species-specific sperm recognition. In this report, we integrate phylogenetics and machine learning to investigate how ZP-N domains diversify in structure and function. The most C-terminal ZP-N domain of each paralog is associated with another domain type (ZP-C), which together form a “ZP module.” All modular ZP-N domains are phylogenetically distinct from nonmodular or free ZP-N domains. Machine learning–based classification identifies eight residues that form a stabilizing network in modular ZP-N domains that is absent in free domains. Positive selection is identified in some free ZP-N domains. Our findings support that strong purifying selection has conserved an essential structural core in modular ZP-N domains, with the relaxation of this structural constraint allowing free N-terminal domains to functionally diversify.     

Computer Simulation of Two-Dimensional Continuum Flows by the Direct Simulation Monte Carlo Method

Abstract

Computer simulations using particles are an attractive method to extract microscopic information of flow phenomena [1]. The molecular dynamics (MD) method, in which Newton's equations of motions are integrated, gives the temporal development of the system. In the MD simulation of fluid flows, the computational region is limited to atomistic scales [2]. On the other hand, the direct simulation Monte Carlo (DSMC) method, in which collisions of particles are made on a probabilistic basis, has a potential of treating a realistic system with a macroscopic scale length retaining the atomistic details. The DSMC method provides an efficient way to integrate the Boltzmann equation from the rarefied gas to the near-continuum region. Bird clarified the validity of the DSMC method in the near-continuum flow region [3]. However, the DSMC method has not been applied to the continuum region and compared with the continuum hydrodynamics.     

Coupling between the voltage-sensing and pore domains in a voltage-gated potassium channel

Voltage-dependent potassium (Kv), sodium (Nav), and calcium channels open and close in response to changes in transmembrane (TM) potential, thus regulating cell excitability by controlling ion flow across the membrane. An outstanding question concerning voltage gating is how voltage-induced conformational changes of the channel voltage-sensing domains (VSDs) are coupled through the S4-S5 interfacial linking helices to the opening and closing of the pore domain (PD). To investigate the coupling between the VSDs and the PD, we generated a closed Kv channel configuration from Aeropyrum pernix (KvAP) using atomistic simulations with experiment-based restraints on the VSDs. Full closure of the channel required, in addition to the experimentally determined TM displacement, that the VSDs be displaced both inwardly and laterally around the PD. This twisting motion generates a tight hydrophobic interface between the S4-S5 linkers and the C-terminal ends of the pore domain S6 helices in agreement with available experimental evidence.     

Dissecting Domain-Specific Evolutionary Pressure Profiles of Transient Receptor Potential Vanilloid Subfamily Members 1 to 4

The transient receptor potential vanilloid family includes four ion channels–TRPV1, TRPV2, TRPV3 and TRPV4–that are represented within the vertebrate subphylum and involved in several sensory and physiological processes. These channels are related to adaptation to the environment, and probably under strong evolutionary pressure. Using multiple sequence alignments as source for evolutionary, bioinformatics and statistical analysis, we have analyzed the evolutionary profiles for TRPV1, TRPV2, TRPV3 and TRPV4. The evolutionary pressure exerted over vertebrate TRPV2 sequences compared to the other channels argues for a positive selection profile for TRPV2 compared to TRPV1, TRPV3 and TRPV4. We have analyzed the selective pressure on specific protein domains, observing a common selective pressure trend for the common TRPV scaffold, consisting of the ankyrin repeat domain, the membrane proximal domain, the transmembrane domain, and the TRP domain. Through a more detailed analysis we have identified evolutionary constraints involved in the subunit contact at the transmembrane domain level. Performing evolutionary comparison, we have translated specific channel structural information such as the transmembrane topology, and the interaction between the membrane proximal domain and the TRP box. We have also identified potential common regulatory domains among all TRPV1-4 members, such as protein-protein, lipid-protein and vesicle trafficking domains.     

Collective Dynamics Underlying Allosteric Transitions in Hemoglobin

Hemoglobin is the prototypic allosteric protein. Still, its molecular allosteric mechanism is not fully understood. To elucidate the mechanism of cooperativity on an atomistic level, we developed a novel computational technique to analyse the coupling of tertiary and quaternary motions. From Molecular Dynamics simulations showing spontaneous quaternary transitions, we separated the transition trajectories into two orthogonal sets of motions: one consisting of intra-chain motions only (referred to as tertiary-only) and one consisting of global inter-chain motions only (referred to as quaternary-only). The two underlying subspaces are orthogonal by construction and their direct sum is the space of full motions. Using Functional Mode Analysis, we were able to identify a collective coordinate within the tertiary-only subspace that is correlated to the most dominant motion within the quaternary-only motions, hence providing direct insight into the allosteric coupling mechanism between tertiary and quaternary conformation changes. This coupling-motion is substantially different from tertiary structure changes between the crystallographic structures of the T- and R-state. We found that hemoglobin''s allosteric mechanism of communication between subunits is equally based on hydrogen bonds and steric interactions. In addition, we were able to affect the T-to-R transition rates by choosing different histidine protonation states, thereby providing a possible atomistic explanation for the Bohr effect.     

A Well-Balanced Preexisting Equilibrium Governs Electron Flux Efficiency of a Multidomain Diflavin Reductase

Diflavin reductases are bidomain electron transfer proteins in which structural reorientation is necessary to account for the various intramolecular and intermolecular electron transfer steps. Using small-angle x-ray scattering and nuclear magnetic resonance data, we describe the conformational free-energy landscape of the NADPH-cytochrome P450 reductase (CPR), a typical bidomain redox enzyme composed of two covalently-bound flavin domains, under various experimental conditions. The CPR enzyme exists in a salt- and pH-dependent rapid equilibrium between a previously described rigid, locked state and a newly characterized, highly flexible, unlocked state. We further establish that maximal electron flux through CPR is conditioned by adjustable stability of the locked-state domain interface under resting conditions. This is rationalized by a kinetic scheme coupling rapid conformational sampling and slow chemical reaction rates. Regulated domain interface stability associated with fast stochastic domain contacts during the catalytic cycle thus provides, to our knowledge, a new paradigm for improving our understanding of multidomain enzyme function.     

Effect of computational methodology on the conformational dynamics of the protein photosensor LOV1 from Chlamydomonas reinhardtii

LOV domains are the light-sensitive protein domains of plant phototropins and bacteria. They photochemically form a covalent bond between a flavin mononucleotide (FMN) chromophore and a cysteine, attached to the apo-protein, upon irradiation with blue light, which triggers a signal in the adjacent kinase. Although their signaling state has been well characterized through experimental means, their signal transduction pathway as well as dark-state activity are generally only poorly understood. Here we show results from molecular dynamics simulations where we investigated the effect of thermostating and long-range electrostatics on the solution structure and dynamical behavior of the wild-type LOV1 domain from the green algae Chlamydomonas reinhardtii in the dark. We demonstrate that these computational issues can dramatically affect the conformational fluctuations of such protein domains by suppressing configurations far from equilibrium or destabilizing local configurations, leading to artificial changes of the protein secondary structure as well as the H-bond network formed by the amino acids and the FMN. By comparing our calculation results with recent experimental data, we show that the non-invasive thermostating strategy, where the protein solute is only indirectly coupled to the thermostat via the solvent, in conjunction with the particle-mesh Ewald technique, provides dark-state conformers, which are in consistency with experimental observations. Moreover, our calculations indicate that the LOV1 domains can alter the intersystem crossing rate and rate of adduct formation by adjusting the population distribution of these dark-state conformers. This might permit them to function as a modulator of the signal intensity under low light conditions.     

Binding free energy landscape of domain-peptide interactions

Peptide recognition domains (PRDs) are ubiquitous protein domains which mediate large numbers of protein interactions in the cell. How these PRDs are able to recognize peptide sequences in a rapid and specific manner is incompletely understood. We explore the peptide binding process of PDZ domains, a large PRD family, from an equilibrium perspective using an all-atom Monte Carlo (MC) approach. Our focus is two different PDZ domains representing two major PDZ classes, I and II. For both domains, a binding free energy surface with a strong bias toward the native bound state is found. Moreover, both domains exhibit a binding process in which the peptides are mostly either bound at the PDZ binding pocket or else interact little with the domain surface. Consistent with this, various binding observables show a temperature dependence well described by a simple two-state model. We also find important differences in the details between the two domains. While both domains exhibit well-defined binding free energy barriers, the class I barrier is significantly weaker than the one for class II. To probe this issue further, we apply our method to a PDZ domain with dual specificity for class I and II peptides, and find an analogous difference in their binding free energy barriers. Lastly, we perform a large number of fixed-temperature MC kinetics trajectories under binding conditions. These trajectories reveal significantly slower binding dynamics for the class II domain relative to class I. Our combined results are consistent with a binding mechanism in which the peptide C terminal residue binds in an initial, rate-limiting step.     

N-palmitoyl-sulfatide participates in lateral domain formation in complex lipid bilayers

Sulfatides (galactosylceramidesulfates) are negatively charged glycosphingolipids that are important constituents of brain myelin membranes. These membranes are also highly enriched in galactosylceramide and cholesterol. It has been implicated that sulfatides, together with other sphingolipids, take part in lateral domain formation in biological membranes. This study was conducted to characterize the lateral phase behavior of N-palmitoyl-sulfatide in mixed bilayer membranes. Going from simple lipid mixtures with sulfatide as the only sphingolipid in a fluid matrix of POPC, to more complex membranes including other sphingolipids, we have examined 1) ordered domain formation with sulfatide, 2) sterol enrichment in such domains and 3) stabilization of the domains against temperature by the addition of calcium. Using two distinct phase selective fluorescent probes, trans-parinaric acid and cholestatrienol, together with a quencher in the fluid phase, we were able to distinguish between ordered domains in general and ordered domains enriched in sterol. We found that N-palmitoyl-sulfatide formed ordered domains when present as the only sphingolipid in a fluid phospholipid bilayer, but these domains did not contain sterol and their stability was unaffected by calcium. However, at low, physiologically relevant concentrations, sulfatide partitioned favorably into domains enriched in other sphingolipids and cholesterol. These domains were stabilized against temperature in the presence of divalent cations. We conclude that sulfatides are likely to affect the lateral organization of biomembranes.     

Integrating NMR,SAXS, and Atomistic Simulations: Structure and Dynamics of a Two-Domain Protein

Multidomain proteins with two or more independently folded functional domains are prevalent in nature. Whereas most multidomain proteins are linked linearly in sequence, roughly one-tenth possess domain insertions where a guest domain is implanted into a loop of a host domain, such that the two domains are connected by a pair of interdomain linkers. Here, we characterized the influence of the interdomain linkers on the structure and dynamics of a domain-insertion protein in which the guest LysM domain is inserted into a central loop of the host CVNH domain. Expanding upon our previous crystallographic and NMR studies, we applied SAXS in combination with NMR paramagnetic relaxation enhancement to construct a structural model of the overall two-domain system. Although the two domains have no fixed relative orientation, certain orientations were found to be preferred over others. We also assessed the accuracies of molecular mechanics force fields in modeling the structure and dynamics of tethered multidomain proteins by integrating our experimental results with microsecond-scale atomistic molecular dynamics simulations. In particular, our evaluation of two different combinations of the latest force fields and water models revealed that both combinations accurately reproduce certain structural and dynamical properties, but are inaccurate for others. Overall, our study illustrates the value of integrating experimental NMR and SAXS studies with long timescale atomistic simulations for characterizing structural ensembles of flexibly linked multidomain systems.     

Hybrid and Rogue Kinases Encoded in the Genomes of Model Eukaryotes

The highly modular nature of protein kinases generates diverse functional roles mediated by evolutionary events such as domain recombination, insertion and deletion of domains. Usually domain architecture of a kinase is related to the subfamily to which the kinase catalytic domain belongs. However outlier kinases with unusual domain architectures serve in the expansion of the functional space of the protein kinase family. For example, Src kinases are made-up of SH2 and SH3 domains in addition to the kinase catalytic domain. A kinase which lacks these two domains but retains sequence characteristics within the kinase catalytic domain is an outlier that is likely to have modes of regulation different from classical src kinases. This study defines two types of outlier kinases: hybrids and rogues depending on the nature of domain recombination. Hybrid kinases are those where the catalytic kinase domain belongs to a kinase subfamily but the domain architecture is typical of another kinase subfamily. Rogue kinases are those with kinase catalytic domain characteristic of a kinase subfamily but the domain architecture is typical of neither that subfamily nor any other kinase subfamily. This report provides a consolidated set of such hybrid and rogue kinases gleaned from six eukaryotic genomes–S.cerevisiae, D. melanogaster, C.elegans, M.musculus, T.rubripes and H.sapiens–and discusses their functions. The presence of such kinases necessitates a revisiting of the classification scheme of the protein kinase family using full length sequences apart from classical classification using solely the sequences of kinase catalytic domains. The study of these kinases provides a good insight in engineering signalling pathways for a desired output. Lastly, identification of hybrids and rogues in pathogenic protozoa such as P.falciparum sheds light on possible strategies in host-pathogen interactions.     

Two photon fluorescence microscopy of coexisting lipid domains in giant unilamellar vesicles of binary phospholipid mixtures

Images of giant unilamellar vesicles (GUVs) formed by different phospholipid mixtures (1,2-dipalmitoyl-sn-glycero-3-phosphocholine/1, 2-dilauroyl-sn-glycero-3-phosphocholine (DPPC/DLPC) 1:1 (mol/mol), and 1,2-dipalmitoyl-sn-glycero-3-phosphoethanolamine/1, 2-dipalmitoyl-sn-glycero-3-phosphocholine (DPPE/DPPC), 7:3 and 3:7 (mol/mol) at different temperatures were obtained by exploiting the sectioning capability of a two-photon excitation fluorescence microscope. 6-Dodecanoyl-2-dimethylamino-naphthalene (LAURDAN), 6-propionyl-2-dimethylamino-naphthalene (PRODAN), and Lissamine rhodamine B 1,2-dihexadecanoyl-sn-glycero-3-phosphoethanolamine (N-Rh-DPPE) were used as fluorescent probes to reveal domain coexistence in the GUVs. We report the first characterization of the morphology of lipid domains in unsupported lipid bilayers. From the LAURDAN intensity images the excitation generalized polarization function (GP) was calculated at different temperatures to characterize the phase state of the lipid domain. On the basis of the phase diagram of each lipid mixture, we found a homogeneous fluorescence distribution in the GUV images at temperatures corresponding to the fluid region in all lipid mixtures. At temperatures corresponding to the phase coexistence region we observed lipid domains of different sizes and shapes, depending on the lipid sample composition. In the case of GUVs formed by DPPE/DPPC mixture, the gel DPPE domains present different shapes, such as hexagonal, rhombic, six-cornered star, dumbbell, or dendritic. At the phase coexistence region, the gel DPPE domains are moving and growing as the temperature decreases. Separated domains remain in the GUVs at temperatures corresponding to the solid region, showing solid-solid immiscibility. A different morphology was found in GUVs composed of DLPC/DPPC 1:1 (mol/mol) mixtures. At temperatures corresponding to the phase coexistence, we observed the gel domains as line defects in the GUV surface. These lines move and become thicker as the temperature decreases. As judged by the LAURDAN GP histogram, we concluded that the lipid phase characteristics at the phase coexistence region are different between the DPPE/DPPC and DLPC/DPPC mixtures. In the DPPE/DPPC mixture the coexistence is between pure gel and pure liquid domains, while in the DLPC/DPPC 1:1 (mol/mol) mixture we observed a strong influence of one phase on the other. In all cases the domains span the inner and outer leaflets of the membrane, suggesting a strong coupling between the inner and outer monolayers of the lipid membrane. This observation is also novel for unsupported lipid bilayers.     

Inherent Regulation of EAL Domain-catalyzed Hydrolysis of Second Messenger Cyclic di-GMP

The universal second messenger cyclic di-GMP (cdG) is involved in the regulation of a diverse range of cellular processes in bacteria. The intracellular concentration of the dinucleotide is determined by the opposing actions of diguanylate cyclases and cdG-specific phosphodiesterases (PDEs). Whereas most PDEs have accessory domains that are involved in the regulation of their activity, the regulatory mechanism of this class of enzymes has remained unclear. Here, we use biophysical and functional analyses to show that the isolated EAL domain of a PDE from Escherichia coli (YahA) is in a fast thermodynamic monomer-dimer equilibrium, and that the domain is active only in its dimeric state. Furthermore, our data indicate thermodynamic coupling between substrate binding and EAL dimerization with the dimerization affinity being increased about 100-fold upon substrate binding. Crystal structures of the YahA-EAL domain determined under various conditions (apo, Mg²⁺, cdG·Ca²⁺ complex) confirm structural coupling between the dimer interface and the catalytic center. The built-in regulatory properties of the EAL domain probably facilitate its modular, functional combination with the diverse repertoire of accessory domains.     

Domain mapping of the Rad51 paralog protein complexes

The five human Rad51 paralogs are suggested to play an important role in the maintenance of genome stability through their function in DNA double-strand break repair. These proteins have been found to form two distinct complexes in vivo, Rad51B–Rad51C–Rad51D–Xrcc2 (BCDX2) and Rad51C–Xrcc3 (CX3). Based on the recent Pyrococcus furiosus Rad51 structure, we have used homology modeling to design deletion mutants of the Rad51 paralogs. The models of the human Rad51B, Rad51C, Xrcc3 and murine Rad51D (mRad51D) proteins reveal distinct N-terminal and C-terminal domains connected by a linker region. Using yeast two-hybrid and co-immunoprecipitation techniques, we have demonstrated that a fragment of Rad51B containing amino acid residues 1–75 interacts with the C-terminus and linker of Rad51C, residues 79–376, and this region of Rad51C also interacts with mRad51D and Xrcc3. We have also determined that the N-terminal domain of mRad51D, residues 4–77, binds to Xrcc2 while the C-terminal domain of mRad51D, residues 77–328, binds Rad51C. By this, we have identified the binding domains of the BCDX2 and CX3 complexes to further characterize the interaction of these proteins and propose a scheme for the three-dimensional architecture of the BCDX2 and CX3 paralog complexes.