MD Trajectories Research Articles

Benchmarking MD Models of Disordered Proteins using NMR Data on Translational Diffusion

MD trajectories can be viewed as “ultimate models” of disordered proteins. However, these models are in need of careful experimental validation. The coefficient of translational diffusion Dtr, measurable by pulsed-field gradient NMR, offers a measure of compactness for disordered proteins (especially useful for smaller proteins, where SAXS data are unavailable). Here we investigate, both experimentally and via MD modeling, the translational diffusion of 25-residue N-terminal fragment from histone H4 (N-H4).

A Machine Learning Algorithm for the Detection of Function-Related Rare Events in MD Trajectories of Biomolecular Systems

The detection of structural changes occurring in the course of an MD trajectory is a key requirement for the analysis of the dynamics of a molecular system in relation to its function. We have developed a protocol based on Machine Learning algorithms to detect rare events in MD simulation trajectories, and used the identified events to analyze functionally significant structural changes occurring in the trajectory. The algorithm utilizes an unsupervised Machine Learning technique named Non-Negative Matrix Factorization (NMF), that learns a sparse, parts-based representation of the data which tends to be more congruent with human intuition than other dimensionality-reduction methods.

How Effective are Retro-Inverso Peptides? Insights from MD Supported by Paramagnetic NMR Data

Early MD simulations were largely limited to observation of fast structural fluctuations in protein molecules. With recent advent of GPU computing, it became possible to model far more relevant dynamic processes: protein folding, ligand binding, allostery, etc. In our study, we use long unbiased MD simulations to model binding of the so-called retro-inverso (ri) peptide to its protein target. Retro-inverso peptides are designed from native peptides by replacing L-amino acids with D-amino acids and inverting the amino-acid sequence. They are believed to retain the binding affinity of their parent peptides, while being resistant to proteolytic degradation. We investigated the binding of retro-inverso peptide ri-Sos1 (rrrppvpppv) to N-terminal SH3 domain from adaptor protein Grb2. The peptide was initially placed in a random position inside a large simulation cell. During the course of the simulation, we observed the binding of ri-Sos1 at the canonical SH3 binding site. The resulting complex features a number of key hydrogen bonds that are similar to the ones found in the native Sos1/Grb2 N-SH3 structure. However, the pattern of these interactions has undergone some noticeable rearrangements. We argue that the concept of retro-inversion lacks any solid theoretical justification and, therefore, ri-peptides are only poor imitators of their parent peptides. To validate our findings, we conducted a number of NMR experiments. In particular, we used the peptides crrrppvpppv and rrrppvpppvc tagged with paramagnetic label MTSL to measure paramagnetic relaxation enhancements (PRE) in Grb2 N-SH3. The results were compared with the PRE rates rigorously calculated from the MD trajectories involving the two respective MTSL-tagged peptides. The 1HN,15N-HSQC-based titration data showed that ri-Sos1 is a much weaker binder than the original Sos1, thus confirming the tenuous nature of ri-peptide ligands. We acknowledge grant support from SPbU 51142660.

In silico exploration of novel protease inhibitors against coronavirus 2019 (COVID-19).

The spread of SARS-CoV-2 has affected human health globally. Hence, it is necessary to rapidly find the drug-candidates that can be used to treat the infection. Since the main protease (Mpro) is the key protein in the virus's life cycle, Mpro is served as one of the critical targets of antiviral treatment. We employed virtual screening tools to search for new inhibitors to accelerate the drug discovery process. The hit compounds were subsequently docked into the active site of SARS-CoV-2 main protease and ranked by their binding energy. Furthermore, in-silico ADME studies were performed to probe for adoption with the standard ranges. Finally, molecular dynamics simulations were applied to study the protein-drug complex's fluctuation over time in an aqueous medium. This study indicates that the interaction energy of the top ten retrieved compounds with COVID-19 main protease is much higher than the interaction energy of some currently in use protease drugs such as ML188, nelfinavir, lopinavir, ritonavir, and α-ketoamide. Among the discovered compounds, Pubchem44326934 showed druglike properties and was further analyzed by MD and MM/PBSA approaches. Besides, the constant binding free energy over MD trajectories suggests a probable drug possessing antiviral properties. MD simulations demonstrate that GLU166 and GLN189 are the most important residues of Mpro, which interact with inhibitors.

Virtual repurposing of ursodeoxycholate and chenodeoxycholate as lead candidates against SARS-Cov2-Envelope protein: A molecular dynamics investigation

Drug repurposing is an apt choice to combat the currently prevailing global threat of COVID-19, caused by SARS-Cov2in absence of any specific medication/vaccine. The present work employs state of art computational methods like homology modelling, molecular docking and molecular dynamics simulations to evaluate the potential of two widely used surfactant drugs namely chenodeoxycholate(CDC) and ursodeoxycholate (UDC), to bind to the envelope protein of SARS-Cov2(SARS-Cov2-E).The monomeric unit of SARS-Cov2-E was modelled from a close homologue (>90% sequence identity) and a pentameric assembly was modelled using symmetric docking, followed by energy minimization in a DPPC membrane environment. The minimized structure was used to generate best scoring SARS-Cov2-E–CDC/UDC complexes through blind docking. These complexes were subjected to 230 ns molecular dynamics simulations in triplicates in a DPPC membrane environment. Comparative analyses of structural properties and molecular interaction profiles from the MD trajectories revealed that, both CDC and UDC could stably bind to SARS-Cov2-E through H-bonds, water-bridges and hydrophobic contacts with the transmembrane-channelresidues.T30 was observed to be a key residue for CDC/UDC binding. CDC/UDC binding affected the H-bonding pattern between adjacent monomeric chains, slackening the compact transmembrane region of SARS-Cov2-E. Additionally, the polar functional groups of CDC/UDC facilitated entry of a large number of water molecules into the channel. These observations suggest CDC/UDC as potential candidates to hinder the survival of SARS-Cov2 by disrupting the structure of SARS-Cov2-E and facilitating the entry of solvents/polar inhibitors inside the viral cell. Communicated by Ramaswamy H. Sarma

Hierarchical Graph Representation of Pharmacophore Models.

For the investigation of protein-ligand interaction patterns, the current accessibility of a wide variety of sampling methods allows quick access to large-scale data. The main example is the intensive use of molecular dynamics simulations applied to crystallographic structures which provide dynamic information on the binding interactions in protein-ligand complexes. Chemical feature interaction based pharmacophore models extracted from these simulations, were recently used with consensus scoring approaches to identify potentially active molecules. While this approach is rapid and can be fully automated for virtual screening, additional relevant information from such simulations is still opaque and so far the full potential has not been entirely exploited. To address these aspects, we developed the hierarchical graph representation of pharmacophore models (HGPM). This single graph representation enables an intuitive observation of numerous pharmacophore models from long MD trajectories and further emphasizes their relationship and feature hierarchy. The resulting interactive depiction provides an easy-to-apprehend tool for the selection of sets of pharmacophores as well as visual support for analysis of pharmacophore feature composition and virtual screening results. Furthermore, the representation can be adapted to include information involving interactions between the same protein and multiple different ligands. Herein, we describe the generation, visualization and use of HGPMs generated from MD simulations of two x-ray crystallographic derived structures of the human glucokinase protein in complex with allosteric activators. The results demonstrate that a large number of pharmacophores and their relationships can be visualized in an interactive, efficient manner, unique binding modes identified and a combination of models derived from long MD simulations can be strategically prioritized for VS campaigns.

Candidate Binding Sites for Allosteric Inhibition of the SARS-CoV-2 Main Protease from the Analysis of Large-Scale Molecular Dynamics Simulations.

We analyzed a 100 μs MD trajectory of the SARS-CoV-2 main protease by a non-parametric data analysis approach which allows characterizing a free energy landscape as a simultaneous function of hundreds of variables. We identified several conformations that, when visited by the dynamics, are stable for several hundred nanoseconds. We explicitly characterize and describe these metastable states. In some of these configurations, the catalytic dyad is less accessible. Stabilizing them by a suitable binder could lead to an inhibition of the enzymatic activity. In our analysis we keep track of relevant contacts between residues which are selectively broken or formed in the states. Some of these contacts are formed by residues which are far from the catalytic dyad and are accessible to the solvent. Based on this analysis we propose some relevant contact patterns and three possible binding sites which could be targeted to achieve allosteric inhibition.

Combination of QSAR, molecular docking, molecular dynamic simulation and MM-PBSA: analogues of lopinavir and favipiravir as potential drug candidates against COVID-19

Pandemic COVID-19 infections have spread throughout the world. There is no effective treatment against this disease. Viral RNA-dependent RNA polymerase (RdRp) catalyzes the replication of RNA from RNA and the main protease (Mpro) has a role in the processing of polyproteins that are translated from the RNA of SARS-CoV-2, and thus these two enzymes are strong candidates for targeting by anti-viral drugs. Small molecules such as lopinavir and favipiravir significantly inhibit the activity of Mpro and RdRp in vitro. Studies have shown that structurally modified lopinavir, favipiravir, and other similar compounds can inhibit COVID-19 main protease (Mpro) and RNA-dependent RNA polymerase (RdRp). In this study, lopinavir and its structurally similar compounds were chosen to bind the main protease, and favipiravir was chosen to target RNA-dependent RNA polymerase. Molecular docking and the quantitative structure-activity relationships (QSAR) study revealed that the selected candidates have favorable binding affinity but less druggable properties. To improve the druggability, four structural analogues of lopinavir and one structural analogue of favipiravir was designed by structural modification. Molecular interaction analyses have displayed that lopinavir and favipiravir analogues interact with the active site residues of Mpro and RdRp, respectively. Absorption, distribution, metabolism, excretion and toxicity (ADMET) properties, medicinal chemistry profile, and physicochemical features were shown that all structurally modified analogues are less toxic and contain high druggable properties than the selected candidates. Subsequently, 50 ns molecular dynamics simulation of the top four docked complexes demonstrated that CID44271905, a lopinavir analogue, forms the most stable complex with the Mpro. Further MMPBSA analyses using the MD trajectories also confirmed the higher binding affinity of CID44271905 towards Mpro. In summary, this study demonstrates a new way to identify leads for novel anti-viral drugs against COVID-19. Communicated by Ramaswamy H. Sarma

Bicyclo[1.1.1]Pentane as Phenyl Substituent in Atorvastatin Drug to improve Physicochemical Properties: Drug-likeness, DFT, Pharmacokinetics, Docking, and Molecular Dynamic Simulation

Hyperlipidemia influences on numerous diseases, such as atherosclerosis, hypertension, and diabetes mellitus, however it still treatable. Atorvastatin is one of the most common statins; yet, it has poor solubility due to the three phenyl rings, which lead to low bioavailability issues. Herein, we used bicyclo [1.1.1] pentane a benzene bioisostere as a replacement for these phenyl groups. In total, seven analogous 1a, 1b, 1c, 2a, 2b, 2c, and 3 were designed and investigated theoretically. Their drug-likeness, physicochemical properties, bioactivity, and pharmacokinetics properties were reported and analyzed. Further DFT calculations at CAM-B3LYP/6-31+g in water were used to report the electronic properties and thermodynamic properties of the suggested analogues. Docking studied were conducted to mimic the interactions of the new derivatives with the active site of HMG-CoA reductase. Since docking results are not reliable, molecular dynamic simulations were carried out on analogous with the best docking score. Finally, binding free energy using MM-GBSA was calculated from the MD trajectories, reveals that compound 1a has promising properties as a drug with a docking score of -8.99 kcal/mol and MM-GBSA of -45.35 kcal/mol.

Examination of Interactions of Hexafluoro-2-propanol with Trp-Cage in Hexafluoro-2-propanol-Water by MD Simulations and Intermolecular Nuclear Overhauser Effects.

All-atom molecular dynamic simulations of the peptide Trp-cage in 30% hexafluoro-2-propanol- water (V/V) at 278 K have been carried out with the goal of exploring peptide hydrogen-solvent fluorine nuclear spin cross relaxation. Force field parameters for HFIP reported by Fioroni et al. along with the fluorine parameters of the TFE5 model reported by this lab were used. Water was represented by the TIP5P-Ew model. Peptide modeling used the AMBER99SB-ILDN force field. Translational diffusion coefficients of solution components at 278 K were predicted to within 35% of experimental values using these parameter sets. The simulations indicate that the solvent mixture is not homogeneous, with HFIP molecules clustered into aggregates as large as 53 fluoroalcohol molecules. The solvent environment of surface atoms of Trp-cage fluctuates between being HFIP-rich and more water-rich about every 10 ns. In accord with previous studies by other groups, the average concentration of HFIP near the surface of the peptide is significantly enhanced over the concentration of HFIP in the bulk solvent. In the simulations, ∼7% of the initial contacts between HFIP molecules and Trp-cage develop into peptide-fluoroalcohol interactions that persist for times as long as 8 ns. Most of the available experimental nuclear spin cross-relaxation rates (ΣHF) for hydrogens of the Trp-cage in 30% HFIP-water are reproduced from the MD trajectories to within uncertainties of the experimental data and the simulations. However, a few calculated ΣHF values for hydrogens of the Trp-cage do not agree with experiment. These tend to be situations where long-lived peptide-HFIP interactions are predicted. The disagreements between observed and calculated ΣHF in these instances signal defects in the modeling parameters and procedures that are presently unrecognized.

Targeting protein tyrosine phosphatase to unravel possible inhibitors for Streptococcus pneumoniae using molecular docking, molecular dynamics simulations coupled with free energy calculations

AimsProtein tyrosine phosphatase (PTP-CPS4B) is a signaling enzyme that is essential for a wide range of cellular processes, like metabolism, proliferation, survival and motility. Studies suggest that PTPs are vital for the production of Wzy-dependent capsule in bacteria, making it a valuable target for the discovery of pneumonia associated anti-virulence antibacterial agents. Present study aims at identifying the potential drug candidates to be exploited in inhibiting the growth of Streptococcus pneumonia targeting PTP-CPS4B. Materials and methodsThe present study exploits the molecular docking potential coupled with molecular dynamic simulation as well as free energy calculations to identify potential inhibitors of PTP-CPS4B. Libraries of known and unknown compounds were docked into the active site of PTP-CPS4B using MOE. The compounds with best binding affinity and orientation were subjected to MD simulations and free energy calculations. FindingsTop three compounds based on their binding energy and well composed interaction pattern obtained from molecular docking study were subjected to MD simulations and were compared to reported antibiotic drugs. MD Simulation studies have shown that the presence of an inhibitor inside the active site reduces protein flexibility as evident from RMSD, RMSF and Principal component analyses. MD simulations identified a transition from extended to bended motional shift in loop α6 of the PTP-CPS4B in ligand bound state. This flexibility was reported in the RMSF analysis and verified by the visual investigation of the loop α6 at different time intervals during the simulation. Free energy of binding affinity (computed using MMPBSA &MMGBSA approach) and the interaction patterns obtained from MD trajectory indicate that compound ZN1 (−31.50 Kcal/mol), ZN2 (−33.14 Kcal/mol) and ZN3 (−26.60 Kcal/mol) are potential drug candidates against PTP-CPS4B. Residue wise decomposition study helped in identifying the role of individual amino acid towards the overall inhibition behavior of the compounds. PCA analysis has led to the conclusion that the behavior of PTP-CPS4B inhibitors causes conformational dynamics that can be used to describe the protein inhibition mechanism. SignificanceThe outcome reveals that this study provide enough evidences for the consideration of ZN1, ZN2, ZN3 as potential PTP-CPS4B inhibitors and further in vitro and in vivo studies may prove their therapeutic potential.

Computational Characterization of the Substrate Activation in the Active Site of SARS-CoV-2 Main Protease

Molecular dynamics simulations with the QM(DFT)/MM potentials are utilized to discriminate between reactive and nonreactive complexes of the SARS-CoV-2 main protease and its substrates. Classification of frames along the molecular dynamic trajectories is utilized by analysis of the 2D maps of the Laplacian of electron density. Those are calculated in the plane formed by the carbonyl group of the substrate and a nucleophilic sulfur atom of the cysteine residue that initiates enzymatic reaction. Utilization of the GPU-based DFT code allows fast and accurate simulations with the hybrid functional PBE0 and double-zeta basis set. Exclusion of the polarization functions accelerates the calculations 2-fold, however this does not describe the substrate activation. Larger basis set with d-functions on heavy atoms and p-functions on hydrogen atoms enables to disclose equilibrium between the reactive and nonreactive species along the MD trajectory. The suggested approach can be utilized to choose covalent inhibitors that will readily interact with the catalytic residue of the selected enzyme.

In-silico prediction of novel drug-target complex of nsp3 of CHIKV through molecular dynamic simulation

Literature reported that nsp3 of CHIKV is an important target for the designing of drug as it involves in the replication, survival etc. Herein, about eighteen million molecules available in the ZINC database are filtered against nsp3 using RASPD. Top five hit drug molecules were then taken from the total screened molecules (6988) from ZINC database. Then, a one pot-three components reaction is designed to get the pyrazolophthalazine and its formation was studied using DFT method. Authors created a library of 200 compounds using the product obtained in the reaction and filtered against nsp3 of CHIKV based on docking using iGEMDOCK, a computational tool. Authors have studied the best molecules after applying the the Lipinski's rule of five and bioactive score. Further, the authors took the best compound i.e. CMPD178 and performed the MD simulations and tdMD simulations with nsp3 protease using AMBER18. MD trajectories were studied to collect the information about the nsp3 of CHIKV with and without screened compound and then, MM-GBSA calculations were performed to calculate change in binding free energies for the formation of complex. The aim of the work is to find the potential candidate as promising inhibitor against nsp3 of CHIKV.

Simple Model of Protein Energetics To Identify Ab Initio Folding Transitions from All-Atom MD Simulations of Proteins.

A fundamental requirement to predict the native conformation, address questions of sequence design and optimization, and gain insights into the folding mechanisms of proteins lies in the definition of an unbiased reaction coordinate that reports on the folding state without the need to compare it to reference values, which might be unavailable for new (designed) sequences. Here, we introduce such a reaction coordinate, which does not depend on previous structural knowledge of the native state but relies solely on the energy partition within the protein: the spectral gap of the pair nonbonded energy matrix (ENergy Gap, ENG). This quantity can be simply calculated along unbiased MD trajectories. We show that upon folding the gap increases significantly, while its fluctuations are reduced to a minimum. This is consistently observed for a diverse set of systems and trajectories. Our approach allows one to promptly identify residues that belong to the folding core as well as residues involved in non-native contacts that need to be disrupted to guide polypeptides to the folded state. The energy gap and fluctuations criteria are then used to develop an automatic detection system which allows us to extract and analyze folding transitions from a generic MD trajectory. We speculate that our method can be used to detect conformational ensembles in dynamic and intrinsically disordered proteins, revealing potential preorganization for binding.

Computational study\xa0of the water-driven graphene wrinkle life-cycle\xa0towards\xa0applications in flexible\xa0electronics

The ubiquitous presence of wrinkles in two-dimensional materials alters their properties significantly. It is observed that during the growth process of graphene, water molecules, sourced from ambient humidity or transferred method used, can get diffused in between graphene and the substrate. The water diffusion causes/assists wrinkle formation in graphene, which influences its properties. The diffused water eventually dries, altering the geometrical parameters and properties of wrinkled graphene nanoribbons. Our study reveals that the initially distributed wrinkles tend to coalesce to form a localized wrinkle whose configuration depends on the initial wrinkle geometry and the quantity of the diffused water. The movement of the localized wrinkle is categorized into three modes—bending, buckling, and sliding. The sliding mode is characterized in terms of velocity as a function of diffused water quantity. Direct bandgap increases linearly with the initial angle except the highest angle considered (21°), which can be attributed to the electron tunneling effect observed in the orbital analysis. The system becomes stable with an increase in the initial angle of wrinkle as observed from the potential energy plots extracted from MD trajectories and confirmed with the DOS plot. The maximum stress generated is less than the plastic limit of the graphene.

The temperature of maximum density for amino acid aqueous solutions. An experimental and molecular dynamics study

The temperature of maximum density (TMD) for aqueous solutions of seven amino acids has been experimentally determined by means of density measurements versus temperature. The selected amino acids have been arginine, cysteine, glutamic acid, glutamine, lysine, methionine, and threonine. The TMD dependence against composition has been obtained from the experimental data and characterized through the Despretz constant. All amino acids induce a depression in the TMD as compared with that of pure water. If the mole fraction is selected as composition variable, a clear dependence against amino acid molar mass is observed, which disappears when TMDs are represented versus mass fraction; almost all data shrink onto a single straight line. It must be pointed out that the TMD depressions for all studied amino acids are quite larger than those previously observed for proteins. This suggests that TMDs for proteins cannot be explained as a simple, additive result of its constituents, and, therefore complex cooperative phenomena seem to take place. The partial molar volume at infinite dilution has been obtained from density data, and a consistency test between its temperature dependence and that of temperature of maximum density versus composition has been performed, obtaining satisfactory results. A molecular dynamics study for all the studied systems has been also carried out. Amino acids have been modeled through the OPLS-AA force field, whereas the TIP4P/2005 model was used for water. The temperature of maximum density and partial molar volume have been calculated from the simulated density, and the results are compared with experimental data. Although the agreement is only fair, similar qualitative trends were obtained. The simulated Despretz constants are smaller than the experimental ones, a result that was already previously observed for methanol aqueous solutions. A structural analysis of water molecules in solution along the MD trajectories showed no enhancement of ice-like structures in complete agreement with the TMD decrease with concentration.

Quantum-Chemical and Molecular Dynamics Investigations of Magnesium Chloride Complexes in Dimethoxyethane Solutions.

Quantum-chemical calculations and classical and ab initio molecular dynamics simulations have been performed to study the Mg2+-conducting electrolytes based on Mg(TFSI)2/MgCl2 solutions in dimethoxyethane. It has been shown that depending on the TFSI/Cl– ratio, the Mg2Cl22+ or Mg3Cl42+ complexes are preferred as stable ion aggregates. In the initial stages of the ion association process, MgCl+, MgCl2, and Mg2Cl3+ are formed as intermediate species. Calculations of harmonic frequencies and simulations of the IR spectrum of the electrolyte from the ab initio MD trajectories have been used to identify the spectral range of vibrations of ion aggregates found in the modeled electrolyte. The results have been discussed in the context of experimental data.

Efficient calculation of excess free energy of pure and mixed alcohol solutions

The calculation of free energy of pure and mixed solution is important in chemistry and related fields. Existing accurate methods use alchemical approaches which require multiple molecular dynamics simulations for a series of artificial intermediate states and are computationally expensive. Here we present an efficient method to calculate excess free energies (both enthalpy and entropy) of pure and mixed alcohol solutions. The present method calculates both enthalpy and entropy of the solution system directly from a single MD trajectory. We performed free energy calculations for several alcohol solutions using two different force fields for comparison with each other and with experimental results. The method is for general mixed solutions with application to mixed water-methanol and water-ethanol solutions at various degrees of mixtures. The calculated numerical results are generally in good agreement with experimental data.

Insights of inhibition mechanism of sifuvirtide and MT-sifuvirtide against wild and mutant HIV-1 envelope glycoprotein41: a molecular dynamics simulation and binding free energy study

ABSTRACTSifuvirtide (SFT) is the second generation inhibitor of HIV-1 gp41 fusion protein, which is under phase III clinical trial for anti-HIV. However, over a period of time virus develops resistance against SFT. To overcome resistance, methionine and threonine (MT) amino acids were implemented in the pocket binding domain of SFT namely MT-Sifuvirtide (MTSFT), which have more helicity; higher melting temperature and stability against resistance-virus. To date, the binding stability and functional role of MT-hook against wild type and mutant form of HIV-1 gp41 fusion protein are not known. In the present study, SFT and MTSFT were individually docked with the wild and V10A/A19I/Q24R mutant form of HIV-1 gp41. Further MD simulation was carried out to understand the stability of the ligand–protein complexes. The free energy values were calculated from the MD trajectories to understand the effect of mutation involved in the binding process. MD results revealed helix to loop type conformational changes in N-heptad repeat (NHR) of gp41 due to mutation. Intermolecular interactions reveal that MTSFT forms more strong interactions with NHR of gp41 than SFT, which is critical for six-helix bundle stabilisation. Hence, the results explained the detailed functional role of MT-hook against gp41 fusion at the molecular level.

Mg2+ vs Ca2+ bound active site of group II intron– A MD study

Group II introns are enzymes which undergo self-splicing and remove itself from pre-messenger RNA. X-ray structures of group II intron of Oceanobacillus iheyensis at various stages of the self-splicing pathway (pre-hydrolytic, post-hydrolytic, and ligand-free state) revealed intricate atomic interaction network in the active site of the intron. It has been confirmed that a heteronuclear metal ion cluster consisting of four metal ions (K1, K2 sites with K+ and M1, M2 sites with Mg2+) are crucial for function. Substitution of Mg2+ by Ca2+ results in loss of enzymatic activity. The X-ray structures not only opens up the possibility of modelling Mg2+ and Ca2+ bound active site of group II intron and quantitatively estimate the energetics of Mg2+ vs Ca2+ preference but also explore the relative structural and dynamical differences in response to divalent metal ion substitution. Thus, using X-ray structures as a template we performed molecular dynamics simulations to compare structural and dynamical differences between Mg2+ and Ca2+ bound active site of group II intron at various stages of the splicing pathway (i.e., pre-hydrolytic, post-hydrolytic, and ligand-free state). Quantitative estimation of Mg2+ vs Ca2+ selectivity at the M1, M2 sites confirmed Mg2+ preference at intron active sites relative to Ca2+. Ca2+ is relatively more hydrated in the intron active site relative to Mg2+. The local environment (bound nucleophilic water, interaction with scissile phosphate) around M2 is strikingly different between Mg2+ and Ca2+ bound pre-hydrolytic state. In the post-hydrolytic state, the exon part of the hydrolysis product is involved in direct interaction with the M1, whereas the intron part is highly flexible in our MD trajectories. Solvent exposure of M1, M2 sites are least in the pre-hydrolytic state, highest in the ligand-free state, and intermediate in the post-hydrolytic state.