Twenty-five years of nucleic acid simulations

Unwinding DNA strands by single-walled carbon nanotubes: Molecular docking and MD simulation approach

Despite the growing research into the use of carbon nano-tubes (CNTs) in science and medicine, concerns about their potential toxicity remain insufficiently studied. This study utilizes molecular docking calculations combined by molecular dynamics simulations to investigate the dynamic intricacies of the interaction between single-walled carbon nanotubes (swCNTs) and double-stranded DNA (dsDNA). By examining the influence of swCNT characteristics such as length, radius, and chirality, our findings shed light on the complex interplay that shapes the binding affinity and stability of the dsDNA-swCNT complex. Molecular docking results identify a zigzag swCNT, with a radius of 0.16 Å and a length of 38 Å, as exhibiting the highest binding affinity with dsDNA (-23.9 kcal/mol). Comprehensive analyses, spanning docking results, binding energies, RMSD, radius of gyration, and potential of mean force (PMF) profiles, provide a detailed understanding of the denaturation dynamics. The PMF profiles reveal the thermodynamic feasibility of the DNA-CNT interaction, outlining distinct energy landscapes and barriers: when the selected swCNT binds within the dsDNA groove, the system becomes trapped at the first and second local energy minima, occurring at 1.48 nm and 1.00 nm, respectively. Intramolecular hydrogen bond calculations show a significant reduction, affirming the denaturing effect of swCNTs on DNA. Furthermore, the study reveals a significant reduction in the binding affinity of Ethidium Bromide (EB) to dsDNA following its interaction with swCNT, with a decrease in EB binding to dsDNA of approximately 13.2 %. This research offers valuable insights into the toxic effects of swCNTs on dsDNA, contributing to a rationalization of the cancerous potential of swCNTs.

Cleavage of DNA Substrate Containing Nucleotide Mismatch in the Complementary Region to sgRNA by Cas9 Endonuclease: Thermodynamic and Structural Features

The non-ideal accuracy and insufficient selectivity of CRISPR/Cas9 systems is a serious problem for their use as a genome editing tool. It is important to select the target sequence correctly so that the CRISPR/Cas9 system does not cut similar sequences. This requires an understanding of how and why mismatches in the target sequence can affect the efficiency of the Cas9/sgRNA complex. In this work, we studied the catalytic activity of the Cas9 enzyme to cleave DNA substrates containing nucleotide mismatch at different positions relative to the PAM in the "seed" sequence. We show that mismatches in the complementarity of the sgRNA/DNA duplex at different positions relative to the protospacer adjacent motif (PAM) sequence tend to decrease the cleavage efficiency and increase the half-maximal reaction time. However, for two mismatches at positions 11 and 20 relative to the PAM, an increase in cleavage efficiency was observed, both with and without an increase in half-reaction time. Thermodynamic parameters were obtained from molecular dynamics results, which showed that mismatches at positions 8, 11, and 20 relative to the PAM thermodynamically stabilize the formed complex, and a mismatch at position 2 of the PAM fragment exerts the greatest stabilization compared to the original DNA sequence. The weak correlation of the thermodynamic binding parameters of the components of the Cas9/sgRNA:dsDNA complex with the cleavage data of DNA substrates containing mismatches indicates that the efficiency of Cas9 operation is mainly affected by the conformational changes in Cas9 and the mutual arrangement of sgRNA and substrates.

Dynamics of terminal fraying-peeling and hydrogen bonds dictates the sequential vs. cooperative melting pathways of nanoscale DNA and PNA triplexes

Peptide nucleic acids (PNAs) are charge-neutral synthetic DNA/RNA analogues. In many aspects of biology and biotechnology, the details of DNA and PNA melting reaction coordinates are crucial, and their associative/dissociative details remain inadequately understood. In the current study, we have attempted to gain insights into comparative melting pathways and binding affinity of iso-sequences of an 18-mer PNA-DNA-PNA triplex and the analogous DNA-DNA-DNA triplex, and DNA-DNA and PNA-DNA duplexes. It is intriguing that while the DNA-DNA-DNA triplex melts in two sequential steps, the PNA-DNA-PNA triplex melts in a single step and the mechanistic aspects for this difference are still not clear. We report an all-atom molecular dynamics simulation of both complexes in the temperature range of 300 to 500 K with 20 K intervals. Based on the trajectory analysis, we provide evidence that the association and dissociation are dictated by the differences in fraying-peeling effects from either terminus to the center in a zipper pattern among the PNA-DNA-PNA triplex and DNA-DNA-DNA triplexes. These are shown to be governed by the different characteristics of H-bonding, RMSD, and Free Energy Landscape (FEL) as analyzed by PCA, leading to the DNA-DNA-DNA triplex exhibiting sequential melting, while the PNA-DNA-PNA triplex shows cooperative melting of the whole fragment in a single-step. The PNA-DNA-PNA triplex base pairs are thermodynamically more stable than the DNA-DNA-DNA triplex, with the binding affinity of PNA-TFO to the PNA : DNA duplex being higher than that of DNA-TFO to the DNA : DNA duplex. The investigation of the association/dissociation of PNA-TFO to the PNA-DNA duplex has relevance and importance in the emerging effective applications of oligonucleotide therapy.

Modeling allosteric mechanisms of eukaryotic type II topoisomerases

Type II topoisomerases (TopoIIs) are ubiquitous enzymes that are involved in crucial nuclear processes such as genome organization, chromosome segregation, and other DNA metabolic processes. These enzymes function as large, homodimeric complexes that undergo a complex cycle of binding and hydrolysis of two ATP molecules in their ATPase domains, which regulates the capture and passage of one DNA double-helix through a second, cleaved DNA molecule. This process requires the transmission of information about the state of the bound nucleotide over vast ranges in the TopoII complex. How this information is transmitted at the molecular level to regulate TopoII functions and how protein substitutions disrupt these mechanisms remains largely unknown. Here, we employed extensive microsecond-scale molecular dynamics simulations of the yeast TopoII enzyme in multiple nucleotide-bound states and with amino acid substitutions near both the N and C termini of the complex. Simulation results indicate that the ATPase domains are remarkably flexible on the sub-microsecond timescale and that these dynamics are modulated by the identity of the bound nucleotides and both local and distant amino acid substitutions. Network analyses point toward specific allosteric networks that transmit information about the hydrolysis cycle throughout the complex, which include residues in both the protein and the bound DNA molecule. Amino acid substitutions weaken many of these pathways. Together, our results provide molecular level details on how the TopoII catalytic cycle is controlled through nucleotide binding and hydrolysis and how mutations may disrupt this process.

Structural and Dynamic Features of the Recognition of 8-oxoguanosine Paired with an 8-oxoG-clamp by Human 8-oxoguanine-DNA Glycosylase

8-oxoguanine (oxoG) is formed in DNA by the action of reactive oxygen species. As a highly mutagenic and the most common oxidative DNA lesion, it is an important marker of oxidative stress. Human 8-oxoguanine-DNA glycosylase (OGG1) is responsible for its prompt removal in human cells. OGG1 is a bifunctional DNA glycosylase with N-glycosylase and AP lyase activities. Aspects of the detailed mechanism underlying the recognition of 8-oxoguanine among numerous intact bases and its subsequent interaction with the enzyme's active site amino acid residues are still debated. The main objective of our work was to determine the effect (structural and thermodynamic) of introducing an oxoG-clamp in model DNA substrates on the process of 8-oxoG excision by OGG1. Towards that end, we used DNA duplexes modeling OGG1-specific lesions: 8-oxoguanine or an apurinic/apyrimidinic site with either cytidine or the oxoG-clamp in the complementary strand opposite to the lesion. It was revealed that there was neither hydrolysis of the N-glycosidic bond at oxoG nor cleavage of the sugar-phosphate backbone during the reaction between OGG1 and oxoG-clamp-containing duplexes. Possible structural reasons for the absence of OGG1 enzymatic activity were studied via the stopped-flow kinetic approach and molecular dynamics simulations. The base opposite the damage was found to have a critical effect on the formation of the enzyme-substrate complex and the initiation of DNA cleavage. The oxoG-clamp residue prevented the eversion of the oxoG base into the OGG1 active site pocket and impeded the correct convergence of the apurinic/apyrimidinic site of DNA and the attacking nucleophilic group of the enzyme. An obtained three-dimensional model of the OGG1 complex with DNA containing the oxoG-clamp, together with kinetic data, allowed us to clarify the role of the contact of amino acid residues with DNA in the formation of (and rearrangements in) the enzyme-substrate complex.

Comprehensive Comparison and Critical Assessment of RNA-Specific Force Fields

Molecular dynamics simulations play a pivotal role in elucidating the dynamic behaviors of RNA structures, offering a valuable complement to traditional methods such as nuclear magnetic resonance or X-ray. Despite this, the current precision of RNA force fields lags behind that of protein force fields. In this work, we systematically compared the performance of four RNA force fields (ff99bsc0χOL3, AMBERDES, ff99OL3_CMAP1, AMBERMaxEnt) across diverse RNA structures. Our findings highlight significant challenges in maintaining stability, particularly with regard to cross-strand and cross-loop hydrogen bonds. Furthermore, we observed the limitations in accurately describing the conformations of nonhelical structural motif, terminal nucleotides, and also base pairing and base stacking interactions by the tested RNA force fields. The identified deficiencies in existing RNA force fields provide valuable insights for subsequent force field development. Concurrently, these findings offer recommendations for selecting appropriate force fields in RNA simulations.

Modeling Allosteric Mechanisms of Eukaryotic Type II Topoisomerases

Type II topoisomerases (TopoIIs) are ubiquitous enzymes that are involved in crucial nuclear processes such as genome organization, chromosome segregation, and other DNA metabolic processes. These enzymes function as large, homodimeric complexes that undergo a complex cycle of binding and hydrolysis of two ATP molecules in their ATPase domains, which regulates the capture and passage of one DNA double-helix through a second, cleaved DNA molecule. This process requires the transmission of information about the state of the bound nucleotide over vast ranges in the TopoII complex. How this information is transmitted at the molecular level to regulate TopoII functions and how protein substitutions disrupt these mechanisms remains largely unknown. Here, we employed extensive microsecond scale molecular dynamics simulations of the yeast TopoII enzyme in multiple nucleotide-bound states and with amino acid substitutions near both the N- and C-terminals of the complex. Simulation results indicate that the ATPase domains are remarkably flexible on the sub-microsecond timescale and that these dynamics are modulated by the identity of the bound nucleotides and both local and distant amino acid substitutions. Network analyses point towards specific allosteric networks that transmit information about the hydrolysis cycle throughout the complex, which include residues in both the protein and the bound DNA molecule. Amino acid substitutions weaken many of these pathways. Together, our results provide molecular-level details on how the TopoII catalytic cycle is controlled through nucleotide binding and hydrolysis and how mutations may disrupt this process.

Rational Approach to Optimizing Conformation-Switching Aptamers for Biosensing Applications

The utilization of structure-switching aptamers (SSAs) has enabled the development of novel sensing platforms for the sensitive and continuous detection of molecules. De novo development of SSAs, however, is complex and laborious. Here we describe a rational approach to SSA optimization that simultaneously improves aptamer binding affinity and introduces target-dependent conformation-switching for compatibility with real-world biosensor applications. Key structural features identified from NMR and computational modeling were used to optimize conformational switching in the presence of target, while large-scale, microarray-based mutation analysis was used to map regions of the aptamer permissive to mutation and identify combinations of mutations with stronger binding affinity. Optimizations were carried out in a relevant biofluid to ensure a seamless transition of the aptamer to a biosensing platform. Initial proof-of-concept for this approach is demonstrated with a cortisol binding aptamer but can easily be translated to other relevant aptamers. Cortisol is a hormone correlated with the stress response that has been associated with various medical conditions and is present at quantifiable levels in accessible biofluids. The ability to continuously track levels of stress in real-time via cortisol monitoring, which can be enabled by the aptamers reported here, is crucial for assessing human health and performance.

Simple Adjustment of Intranucleotide Base-Phosphate Interaction in the OL3 AMBER Force Field Improves RNA Simulations

Molecular dynamics (MD) simulations represent an established tool to study RNA molecules. The outcome of MD studies depends, however, on the quality of the force field (ff). Here we suggest a correction for the widely used AMBER OL3 ff by adding a simple adjustment of the nonbonded parameters. The reparameterization of the Lennard-Jones potential for the -H8···O5'- and -H6···O5'- atom pairs addresses an intranucleotide steric clash occurring in the type 0 base-phosphate interaction (0BPh). The nonbonded fix (NBfix) modification of 0BPh interactions (NBfix0BPh modification) was tuned via a reweighting approach and subsequently tested using an extensive set of standard and enhanced sampling simulations of both unstructured and folded RNA motifs. The modification corrects minor but visible intranucleotide clash for the anti nucleobase conformation. We observed that structural ensembles of small RNA benchmark motifs simulated with the NBfix0BPh modification provide better agreement with experiments. No side effects of the modification were observed in standard simulations of larger structured RNA motifs. We suggest that the combination of OL3 RNA ff and NBfix0BPh modification is a viable option to improve RNA MD simulations.

AmberTools

AmberTools is a free and open-source collection of programs used to set up, run, and analyze molecular simulations. The newer features contained within AmberTools23 are briefly described in this Application note.

All-Atom Molecular Dynamics Simulations of Membrane-Spanning DNA Origami Nanopores

Building on the recent technological advances, all-atom molecular dynamics (MD) simulations have become an indispensable tool to study the molecular behavior at nanoscale. Molecular simulations have been used to characterize the structure, dynamics, and mechanical and electrical properties of DNA origami objects. In this chapter we describe a method to build all-atom model of lipid-spanning DNA origami nanopores and perform molecular dynamics simulations in explicit electrolyte solutions.

Polarizable Molecular Block Model: Toward the Development of an Induced Dipole Force Field for DNA

For flexible and highly ionized macromolecules such as DNA, it is important to correctly evaluate the intramolecular polarization in an induced dipole force field. In a proposed polarizable molecular block (PMB) model, a large molecule is divided into several molecular blocks. The atomic charges of the blocks are optimized by using the respective electrostatic potentials (ESPs) on the molecular surface. By using the capped hydrogen removal operation, the total charge of the blocks is controlled exactly to have an integer charge. The atomic polarizabilities of the blocks are optimized by using the respective polarized one-electron potentials that are the differences between ESPs with and without an external test charge. Induced dipole-charge interactions between the blocks are all included, but those interactions within the blocks are strictly excluded. All dipole-dipole interactions are included, but the damping functions are applied to the close dipole-dipole pairs. Several types of damping (simple scaling, exponential, linear, and Gaussian) are evaluated. The validity of the PMB model was verified by using trinucleotide duplexes which have A-, B-, and Z-DNA forms. The reference energies of trinucleotide duplexes including counterions (GGT3Na-ACC3Na, GAC3Na-GTC3Na, and GCG3Na-CGC3Na) are calculated using ωB97XD/aug-cc-pVDZ. All damping types reproduced well the reference interaction energies, dipole moments, and ESPs. Among them, the simple scaling with strong attenuation to 1-2 atomic pairs showed the highest stability against the polarization catastrophe. This study shows that it is possible to develop a high-quality polarizable force field by treating the intramolecular polarization on a block-by-block basis.

Structure- and Content-Dependent Efficiency of Cas9-Assisted DNA Cleavage in Genome-Editing Systems

Genome-editing systems, being some of the key tools of molecular biologists, represent a reasonable hope for progress in the field of personalized medicine. A major problem with such systems is their nonideal accuracy and insufficient selectivity. The selectivity of CRISPR-Cas9 systems can be improved in several ways. One efficient way is the proper selection of the consensus sequence of the DNA to be cleaved. In the present work, we attempted to evaluate the effect of formed non-Watson-Crick pairs in a DNA duplex on the efficiency of DNA cleavage in terms of the influence of the structure of the formed partially complementary pairs. We also studied the effect of the location of such pairs in DNA relative to the PAM (protospacer-adjacent motif) on the cleavage efficiency. We believe that the stabilization of the Cas9-sgRNA complex with a DNA substrate containing noncomplementary pairs is due to loop reorganization in the RuvC domain of the enzyme. In addition, PAM-proximal mismatches in the DNA substrate lower enzyme efficiency because the "seed" region is involved in binding and cleavage, whereas PAM-distal mismatches have no significant impact on target DNA cleavage. Our data suggest that in the case of short duplexes with mismatches, the stages of recognition and binding of dsDNA substrates by the enzyme determine the reaction rate and time rather than the thermodynamic parameters affected by the "unwinding" of DNA. The results will provide a theoretical basis for predicting the efficiency and accuracy of CRISPR-Cas9 systems at cleaving target DNA.

In-silico probing of AML related RUNX1 cancer-associated missense mutations: Predicted relationships to DNA binding and drug interactions

The molecular consequences of cancer associated mutations in Acute myeloid leukemia (AML) linked factors are not very well understood. Here, we interrogated the COSMIC database for missense mutations associated with the RUNX1 protein, that is frequently mis-regulated in AML, where we sought to identify recurrently mutated positions at the DNA-interacting interface. Indeed, six of the mutated residues, out of a total 417 residues examined within the DNA binding domain, evidenced reduced DNA association in in silico predictions. Further, given the prominence of RUNX1’s compromised function in AML, we asked the question if the mutations themselves might alter RUNX1’s interaction (off-target) with known FDA-approved drug molecules, including three currently used in treating AML. We identified several AML-associated mutations in RUNX1 that were calculated to enhance RUNX1’s interaction with specific drugs. Specifically, we retrieved data from the COSMIC database for cancer-associated mutations of RUNX1 by using R package “data.table” and “ggplot2” modules. In the presence of DNA and/or drug, we used docking scores and energetics of the complexes as tools to evaluate predicted interaction strengths with RUNX1. For example, we performed predictions of drug binding pockets involving Enasidenib, Giltertinib, and Midostaurin (AML associated), as well as ten different published cancer associated drug compounds. Docking of wild type RUNX1 with these 13 different cancer-associated drugs indicates that wild-type RUNX1 has a lower efficiency of binding while RUNX1 mutants R142K, D171N, R174Q, P176H, and R177Q suggested higher affinity of drug association. Literature evidence support our prediction and suggests the mutation R174Q affects RUNX1 DNA binding and could lead to compromised function. We conclude that specific RUNX1 mutations that lessen DNA binding facilitate the binding of a number of tested drug molecules. Further, we propose that molecular modeling and docking studies for RUNX1 in the presence of DNA and/or drugs enables evaluation of the potential impact of RUNX1 cancer associated mutations in AML.

Thermodynamic Swings: How Ideal Complex of Cas9-RNA/DNA Forms

Most processes of the recognition and formation of specific complexes in living systems begin with collisions in solutions or quasi-solutions. Then, the thermodynamic regulation of complex formation and fine tuning of complexes come into play. Precise regulation is very important in all cellular processes, including genome editing using the CRISPR-Cas9 tool. The Cas9 endonuclease is an essential component of the CRISPR-Cas-based genome editing systems. The attainment of high-specificity and -efficiency Cas9 during targeted DNA cleavage is the main problem that limits the practical application of the CRISPR-Cas9 system. In this study, we analyzed the thermodynamics of interaction of a complex's components of Cas9-RNA/DNA through experimental and computer simulation methods. We found that there is a small energetic preference during Cas9-RNA/DNA formation from the Cas9-RNA and DNA/DNA duplex. The small difference in binding energy is relevant for biological interactions and could be part of the sequence-specific recognition of double-stranded DNA by the CRISPR-Cas9 system.

Targeting the Essential Transcription Factor HP1043 of Helicobacter pylori: A Drug Repositioning Study

Antibiotic-resistant bacterial pathogens are a very challenging problem nowadays. Helicobacter pylori is one of the most widespread and successful human pathogens since it colonizes half of the world population causing chronic and atrophic gastritis, peptic ulcer, mucosa-associated lymphoid tissue-lymphoma, and even gastric adenocarcinoma. Moreover, it displays resistance to numerous antibiotics. One of the H. pylori pivotal transcription factors, HP1043, plays a fundamental role in regulating essential cellular processes. Like other bacterial transcription factors, HP1043 does not display a eukaryote homolog. These characteristics make HP1043 a promising candidate to develop novel antibacterial strategies. Drug repositioning is a relatively recent strategy employed in drug development; testing approved drugs on new targets considerably reduces the time and cost of this process. The combined computational and in vitro approach further reduces the number of compounds to be tested in vivo. Our aim was to identify a subset of known drugs able to prevent HP1043 binding to DNA promoters. This result was reached through evaluation by molecular docking the binding capacity of about 14,350 molecules on the HP1043 dimer in both conformations, bound and unbound to the DNA. Employing an ad hoc pipeline including MMGBSA molecular dynamics, a selection of seven drugs was obtained. These were tested in vitro by electrophoretic mobility shift assay to evaluate the HP1043–DNA interaction. Among these, three returned promising results showing an appreciable reduction of the DNA-binding activity of HP1043. Overall, we applied a computational methodology coupled with experimental validation of the results to screen a large number of known drugs on one of the H. pylori essential transcription factors. This methodology allowed a rapid reduction of the number of drugs to be tested, and the drug repositioning approach considerably reduced the drug design costs. Identified drugs do not belong to the same pharmaceutical category and, by computational analysis, bound different cavities, but all display a reduction of HP1043 binding activity on the DNA.

Toward Convergence in Folding Simulations of RNA Tetraloops: Comparison of Enhanced Sampling Techniques and Effects of Force Field Modifications

Atomistic molecular dynamics simulations represent an established technique for investigation of RNA structural dynamics. Despite continuous development, contemporary RNA simulations still suffer from suboptimal accuracy of empirical potentials (force fields, ffs) and sampling limitations. Development of efficient enhanced sampling techniques is important for two reasons. First, they allow us to overcome the sampling limitations, and second, they can be used to quantify ff imbalances provided they reach a sufficient convergence. Here, we study two RNA tetraloops (TLs), namely the GAGA and UUCG motifs. We perform extensive folding simulations and calculate folding free energies (ΔG_fold^°) with the aim to compare different enhanced sampling techniques and to test several modifications of the nonbonded terms extending the AMBER OL3 RNA ff. We demonstrate that replica-exchange solute tempering (REST2) simulations with 12-16 replicas do not show any sign of convergence even when extended to a timescale of 120 μs per replica. However, the combination of REST2 with well-tempered metadynamics (ST-MetaD) achieves good convergence on a timescale of 5-10 μs per replica, improving the sampling efficiency by at least 2 orders of magnitude. Effects of ff modifications on ΔG_fold^° energies were initially explored by the reweighting approach and then validated by new simulations. We tested several manually prepared variants of the gHBfix potential which improve stability of the native state of both TLs by ∼2 kcal/mol. This is sufficient to conveniently stabilize the folded GAGA TL while the UUCG TL still remains under-stabilized. Appropriate adjustment of van der Waals parameters for C-H···O5' base-phosphate interaction may further stabilize the native states of both TLs by ∼0.6 kcal/mol.

Loading of DOX into a tetrahedral DNA nanostructure: the corner does matter

With the rapid development of nanotechnology, various DNA nanostructures have been synthesized and widely used in drug delivery. However, the underlying mechanisms of drug molecule loading into the DNA nanostructure are still elusive. In this work, we systematically investigate the interactions of a tetrahedral DNA nanostructure (TDN) with the anti-cancer drug doxorubicin (DOX) by combining molecular docking and all-atom molecular dynamics simulations. It is found that there are five possible binding modes in the single TDN-DOX interactions, namely the outside-corner mode, the inside-corner mode, the major-groove mode, the minor-groove mode, and the intercalation mode, where the van der Waals (VDW) interaction and the electrostatic (ELE) interaction dominate in the case of unionized DOX and ionized DOX, respectively. Moreover, with the increase of the DOX number, some of the interaction modes may disappear and the inside-corner mode is the most energy-favorable mode. The present study enhances the molecular understanding of the role of TDN as the drug carrier, which may provide a useful guideline for the future design of DNA nanostructures.

Mechanism of RNA recognition by a Musashi RNA-binding protein

The Musashi RNA-binding proteins (RBPs) regulate translation of target mRNAs and maintenance of cell stemness and tumorigenesis. Musashi-1 (MSI1), long considered as an intestinal and neural stem cell marker, has been more recently found to be over expressed in many cancers. It has served as an important drug target for treating acute myeloid leukemia and solid tumors such as ovarian, colorectal and bladder cancer. One of the reported binding targets of MSI1 is Numb, a negative regulator of the Notch signaling. However, the dynamic mechanism of Numb RNA binding to MSI1 remains unknown, largely hindering effective drug design targeting this critical interaction. Here, we have performed extensive all-atom microsecond-timescale simulations using a robust Gaussian accelerated molecular dynamics (GaMD) method, which successfully captured multiple times of spontaneous and highly accurate binding of the Numb RNA from bulk solvent to the MSI1 protein target site. GaMD simulations revealed that Numb RNA binding to MSI1 involved largely induced fit in both the RNA and protein. The simulations also identified important low-energy intermediate conformational states during RNA binding, in which Numb interacted mainly with the β2-β3 loop and C terminus of MSI1. The mechanistic understanding of RNA binding obtained from our GaMD simulations is expected to facilitate rational structure-based drug design targeting MSI1 and other RBPs.

CAN: A new program to streamline preparation of molecular coordinate files for molecular dynamics simulations

Preparing molecular coordinate files for molecular dynamics (MD) simulations can be a very time-consuming process. Herein we present the development of a user-friendly program that drastically reduces the time required to prepare these molecular coordinate files for MD software packages such as AmberTools. Our program, known as charge atomtype naming (CAN), creates and uses a library of structures such as amino acid monomers to update the charge, atom type, and name of atoms in any molecular structure (mol2) file. We demonstrate the utility of this new program by rapidly preparing structural files for MD simulations for polypeptides ranging from small molecules to large protein structures. Both native and non-native amino acid residues are easily handled by this new program.

The Behavior of Triplet Thymine in a Model B-DNA Strand. Energetics and Spin Density Localization Revealed by ab initio Molecular Dynamics Simulations†

Among the naturally occurring nucleobases, thymine presents the lowest triplet state, hence it represents a hotspot for energy transfer and photosensitization. In turn, the population of the triplet state may lead to thymine dimerization and hence to the production of toxic DNA lesions and has been the subject of intensive theoretical and experimental investigations. Relying on QM/MM molecular dynamics simulations, we have sought to situate the energy of the lowest triplet state of thymine embedded in a B-DNA environment. The energy gap varies between 305 and 329 kJ mol^-1 when a single thymine is treated at the quantum chemistry level, depending on its position in the model double-stranded 16-bp oligonucleotide. The energy of triplet state decreases up to 300 kJ mol^-1 , due to polarization effects, when we consider coupled stacked nucleobases up to the inclusion of four nucleobases. Our value lies in good agreement with the energy inferred experimentally by Miranda and coworkers (270 kJ mol^-1 ), and our theoretical exploration opens the door to investigations toward other more complex and biologically relevant environments, such as thymines embedded in nucleosome core particles. Our investigations also provide a reference for further studies using semi-empirical approaches such as density functional-based tight-binding, allowing to further rationalize sequence effects.

DNA Aptamers Block the Receptor Binding Domain at the Spike Protein of SARS-CoV-2

DNA aptamers are versatile molecular species obtained by the folding of short single-stranded nucleotide sequences, with highly specific recognition capabilities against proteins. Here we test the ability of DNA aptamers to interact with the spike (S-)protein of the SARS-CoV-2 viral capsid. The S-protein, a trimer made up of several subdomains, develops the crucial function of recognizing the ACE2 receptors on the surface of human cells, and subsequent fusioning of the virus membrane with the host cell membrane. In order to achieve this, the S1 domain of one protomer switches between a closed conformation, in which the binding site is inaccessible to the cell receptors, and an open conformation, in which ACE2 can bind, thereby initiating the entry process of the viral genetic material in the host cell. Here we show, by means of state-of-the-art molecular simulations, that small DNA aptamers experimentally identified can recognize the S-protein of SARS-CoV-2, and characterize the details of the binding process. We find that their interaction with different subdomains of the S-protein can effectively block, or at least considerably slow down the opening process of the S1 domain, thereby significantly reducing the probability of virus-cell binding. We provide evidence that, as a consequence, binding of the human ACE2 receptor may be crucially affected under such conditions. Given the facility and low cost of fabrication of specific aptamers, the present findings could open the way to both an innovative viral screening technique with sub-nanomolar sensitivity, and to an effective and low impact curative strategy.

Molecular dynamics simulations of allosteric motions and competitive inhibition of the Zika virus helicase

The 2015 Zika outbreak sparked major global concern and emphasized the reality and dangers still posed by mosquito borne pathogens. While efforts have been made to develop a vaccine and other therapeutics, there is still a great demand for antiviral drugs targeting Zika and other flaviviruses. The non-structural protein 3 (NS3) helicase is a vital component of the viral replication complex, tasked with unwinding the viral dsRNA molecule into single strands. Given this critical function, the Zika virus helicase is a potential therapeutic target and the focus of many ongoing research efforts. Using a combination of drug docking and molecular dynamics simulations, we have identified a list of competitive helicase inhibitors targeting the ATP hydrolysis site and have discovered a potential allosteric site capable of distorting both of the protein's active sites.

Definition of the Binding Architecture to a Target Promoter of HP1043, the Essential Master Regulator of Helicobacter pylori

HP1043 is an essential orphan response regulator of Helicobacter pylori orchestrating multiple crucial cellular processes. Classified as a member of the OmpR/PhoB family of two-component systems, HP1043 exhibits a highly degenerate receiver domain and evolved to function independently of phosphorylation. Here, we investigated the HP1043 binding mode to a target sequence in the hp1227 promoter (P_hp1227). Scanning mutagenesis of HP1043 DNA-binding domain and consensus sequence led to the identification of residues relevant for the interaction of the protein with a target DNA. These determinants were used as restraints to guide a data-driven protein-DNA docking. Results suggested that, differently from most other response regulators of the same family, HP1043 binds in a head-to-head conformation to the P_hp1227 target promoter. HP1043 interacts with DNA largely through charged residues and contacts with both major and minor grooves of the DNA are required for a stable binding. Computational alanine scanning on molecular dynamics trajectory was performed to corroborate our findings. Additionally, in vitro transcription assays confirmed that HP1043 positively stimulates the activity of RNA polymerase.

Low-frequency vibrational modes in G-quadruplexes reveal the mechanical properties of nucleic acids

Low-frequency vibrations play an essential role in biomolecular processes involving DNA such as gene expression, charge transfer, drug intercalation, and DNA–protein recognition. However, understanding the vibrational basis of these mechanisms relies on theoretical models due to the lack of experimental evidence. Here we present the low-frequency vibrational spectra of G-quadruplexes (structures formed by four strands of DNA) and B-DNA characterized using femtosecond optical Kerr-effect spectroscopy. Contrary to expectation, we found that G-quadruplexes show several strongly underdamped delocalized phonon-like modes that have the potential to contribute to the biology of the DNA at the atomic level. In addition, G-quadruplexes present modes at a higher frequency than B-DNA demonstrating that changes in the stiffness of the molecule alter its gigahertz to terahertz vibrational profile.

Low-frequency vibrations play an essential role in biomolecular processes involving DNA such as gene expression, charge transfer, drug intercalation, and DNA–protein recognition.

Implicit Solvents for the Polarizable Atomic Multipole AMOEBA Force Field

Computational protein design, ab initio protein/RNA folding, and protein-ligand screening can be too computationally demanding for explicit treatment of solvent. For these applications, implicit solvent offers a compelling alternative, which we describe here for the polarizable atomic multipole AMOEBA force field based on three treatments of continuum electrostatics: numerical solutions to the nonlinear and linearized versions of the Poisson-Boltzmann equation (PBE), the domain-decomposition conductor-like screening model (ddCOSMO) approximation to the PBE, and the analytic generalized Kirkwood (GK) approximation. The continuum electrostatics models are combined with a nonpolar estimator based on novel cavitation and dispersion terms. Electrostatic model parameters are numerically optimized using a least-squares style target function based on a library of 103 small-molecule solvation free energy differences. Mean signed errors for the adaptive Poisson-Boltzmann solver (APBS), ddCOSMO, and GK models are 0.05, 0.00, and 0.00 kcal/mol, respectively, while the mean unsigned errors are 0.70, 0.63, and 0.58 kcal/mol, respectively. Validation of the electrostatic response of the resulting implicit solvents, which are available in the Tinker (or Tinker-HP), OpenMM, and Force Field X software packages, is based on comparisons to explicit solvent simulations for a series of proteins and nucleic acids. Overall, the emergence of performative implicit solvent models for polarizable force fields opens the door to their use for folding and design applications.

Effects of H2A.B incorporation on nucleosome structures and dynamics

The H2A.B histone variant is an epigenetic regulator involved in transcriptional upregulation, DNA synthesis, and splicing that functions by replacing the canonical H2A histone in the nucleosome core particle. Introduction of H2A.B results in less compact nucleosome states with increased DNA unwinding and accessibility at the nucleosomal entry and exit sites. Despite being well characterized experimentally, the molecular mechanisms by which H2A.B incorporation alters nucleosome stability and dynamics remain poorly understood. To study the molecular mechanisms of H2A.B, we have performed a series of conventional and enhanced sampling molecular dynamics simulation of H2A.B- and canonical H2A-containing nucleosomes. Results of conventional simulations show that H2A.B weakens protein-protein and protein-DNA interactions at specific locations throughout the nucleosome. These weakened interactions result in significantly more DNA opening from both the entry and exit sites in enhanced sampling simulations. Furthermore, free energy profiles show that H2A.B-containing nucleosomes have significantly broader free wells and that H2A.B allows for sampling of states with increased DNA breathing, which are shown to be stable on the hundreds of nanoseconds timescale with further conventional simulations. Together, our results show the molecular mechanisms by which H2A.B creates less compacted nucleosome states as a means of increasing genetic accessibility and gene transcription.

Combining Virtual Reality Visualization with Ensemble Molecular Dynamics to Study Complex Protein Conformational Changes

Molecular dynamics (MD) simulations are increasingly used to elucidate relationships between protein structure, dynamics, and their biological function. Currently, it is extremely challenging to perform MD simulations of large-scale structural rearrangements in proteins that occur on millisecond timescales or beyond, as this requires very significant computational resources, or the use of cumbersome "collective variable" enhanced sampling protocols. Here, we describe a framework that combines ensemble MD simulations and virtual reality visualization (eMD-VR) to enable users to interactively generate realistic descriptions of large amplitude, millisecond timescale protein conformational changes in proteins. Detailed tests demonstrate that eMD-VR substantially decreases the computational cost of folding simulations of a WW domain, without the need to define collective variables a priori. We further show that eMD-VR generated pathways can be combined with Markov state models to describe the thermodynamics and kinetics of large-scale loop motions in the enzyme cyclophilin A. Our results suggest eMD-VR is a powerful tool for exploring protein energy landscapes in bioengineering efforts.

UUCG RNA Tetraloop as a Formidable Force-Field Challenge for MD Simulations

Explicit solvent atomistic molecular dynamics (MD) simulations represent an established technique to study structural dynamics of RNA molecules and an important complement for diverse experimental methods. However, performance of molecular mechanical (MM) force fields (ff's) remains far from satisfactory even after decades of development, as apparent from a problematic structural description of some important RNA motifs. Actually, some of the smallest RNA molecules belong to the most challenging systems for MD simulations and, among them, the UUCG tetraloop is saliently difficult. We report a detailed analysis of UUCG MD simulations, depicting the sequence of events leading to the loss of the UUCG native state during MD simulations. The total amount of MD simulation data analyzed in this work is close to 1.3 ms. We identify molecular interactions, backbone conformations, and substates that are involved in the process. Then, we unravel specific ff deficiencies using diverse quantum mechanical/molecular mechanical (QM/MM) and QM calculations. Comparison between the MM and QM methods shows discrepancies in the description of the 5'-flanking phosphate moiety and both signature sugar-base interactions. Our work indicates that poor behavior of the UUCG tetraloop in simulations is a complex issue that cannot be attributed to one dominant and straightforwardly correctable factor. Instead, there is a concerted effect of multiple ff inaccuracies that are coupled and amplifying each other. We attempted to improve the simulation behavior by some carefully tailored interventions, but the results were still far from satisfactory, underlying the difficulties in development of accurate nucleic acid ff's.

Methods and Applications of In Silico Aptamer Design and Modeling

Aptamers are nucleic acid analogues of antibodies with high affinity to different targets, such as cells, viruses, proteins, inorganic materials, and coenzymes. Empirical approaches allow the design of in vitro aptamers that bind particularly to a target molecule with high affinity and selectivity. Theoretical methods allow significant expansion of the possibilities of aptamer design. In this study, we review theoretical and joint theoretical-experimental studies dedicated to aptamer design and modeling. We consider aptamers with different targets, such as proteins, antibiotics, organophosphates, nucleobases, amino acids, and drugs. During nucleic acid modeling and in silico design, a full set of in silico methods can be applied, such as docking, molecular dynamics (MD), and statistical analysis. The typical modeling workflow starts with structure prediction. Then, docking of target and aptamer is performed. Next, MD simulations are performed, which allows for an evaluation of the stability of aptamer/ligand complexes and determination of the binding energies with higher accuracy. Then, aptamer/ligand interactions are analyzed, and mutations of studied aptamers made. Subsequently, the whole procedure of molecular modeling can be reiterated. Thus, the interactions between aptamers and their ligands are complex and difficult to understand using only experimental approaches. Docking and MD are irreplaceable when aptamers are studied in silico.

Multipolar Description of Atom-Atom Electrostatic Interaction Energies in Single/Double-Stranded DNAs

Molecular force field simulation is an effective method to explore the properties of DNA molecules in depth. Almost all current popular force fields calculate atom-atom electrostatic interaction energies for DNAs based on the atomic charge and dipole or quadrupole moments, without considering high-rank atomic multipole moments for more accurate electrostatics. Actually, the distribution of electrons around atomic nuclei is not spherically symmetric but is geometry dependent. In this work, a multipole expansion method that allows us to combine polarizability and anisotropy was applied. One single-stranded DNA and one double-stranded DNA were selected as pilot systems. Deoxynucleotides were cut out from pilot systems and capped by mimicking the original DNA environment. Atomic multipole moments were integrated instead of fixed-point charges to calculate atom-atom electrostatic energies to improve the accuracy of force fields for DNA simulations. Also, the applicability of modeling the behavior of both single-stranded and double-stranded DNAs was investigated. The calculation results indicated that the models can be transferred from pilot systems to test systems, which is of great significance for the development of future DNA force fields.

In silico design and validation of high-affinity RNA aptamers targeting epithelial cellular adhesion molecule dimers

Nucleic acid aptamers hold great promise for therapeutic applications due to their favorable intrinsic properties, as well as high-throughput experimental selection techniques. Despite the utility of the systematic evolution of ligands by the exponential enrichment (SELEX) method for aptamer determination, complementary in silico aptamer design is highly sought after to facilitate virtual screening and increased understanding of important nucleic acid-protein interactions. Here, with a combined experimental and theoretical approach, we have developed two optimal epithelial cellular adhesion molecule (EpCAM) aptamers. Our structure-based in silico method first predicts their binding modes and then optimizes them for EpCAM with molecular dynamics simulations, docking, and free energy calculations. Our isothermal titration calorimetry experiments further confirm that the EpCAM aptamers indeed exhibit enhanced affinity over a previously patented nanomolar aptamer, EP23. Moreover, our study suggests that EP23 and the de novo designed aptamers primarily bind to EpCAM dimers (and not monomers, as hypothesized in previous published works), suggesting a paradigm for developing EpCAM-targeted therapies.

Strongly Bent Double-Stranded DNA: Reconciling Theory and Experiment

The strong bending of polymers is poorly understood. We propose a general quantitative framework of polymer bending that includes both the weak and strong bending regimes on the same footing, based on a single general physical principle. As the bending deformation increases beyond a certain (polymer-specific) point, the change in the convexity properties of the effective bending energy of the polymer makes the harmonic deformation energetically unfavorable: in this strong bending regime the energy of the polymer varies linearly with the average bending angle as the system follows the convex hull of the deformation energy function. For double-stranded DNA, the effective bending deformation energy becomes non-convex for bends greater than ~ 2° per base-pair, equivalent to the curvature of a closed circular loop of ~ 160 base pairs. A simple equation is derived for the polymer loop energy that covers both the weak and strong bending regimes. The theory shows quantitative agreement with recent DNA cyclization experiments on short DNA fragments, while maintaining the expected agreement with experiment in the weak bending regime. Counter-intuitively, cyclization probability (j-factor) of very short DNA loops is predicted to increase with decreasing loop length; the j-factor reaches its minimum for loops of ≃ 45 base pairs. Atomistic simulations reveal that the attractive component of the short-range Lennard-Jones interaction between the backbone atoms can explain the underlying non-convexity of the DNA effective bending energy, leading to the linear bending regime. Applicability of the theory to protein-DNA complexes, including the nucleosome, is discussed.

Hydrophobic unnatural base pairs show a Watson-Crick pairing in micro-second molecular dynamics simulations

Two unnatural hydrophobic nucleotides, d5SICS (2,6-dimethyl-2H-isoquiniline-1-thione) and dNaM (2-methoxy-3-methylnaphthalene), were previously replicated in vivo by a modified E. coli strand, however, a consistent structure for their pairing in terms of specific and selective directional interactions remains elusive, as data from spectroscopy experiments and simulations are inconsistent. The proposed d5SICS-dNaM pairing has been suggested to be a stacked configuration as suggested by NMR data; simulations have failed to reproduce this configuration and a Watson-Crick like pairing is observed. Previously, we focused on reproducing the d5SICS-dNaM Unnatural Base Pair (UBP) paring using an older (bsc0) AMBER force field, which was not able to correctly reproduce the experimental data. We present our efforts to reproduce the experimental pairing using the current version of the AMBER DNA force fields (OL15 and bsc1), two water models (OPC and TIP3P) and external electrostatic stabilization by Mg²⁺ ions. Opposite to previously reported simulations, a Watson-Crick-like pairing with no hydrogen bonds persists throughout all our results. Despite our efforts to replicate the reported stacked conformation, we cannot confirm its plausibility nor obtain a consistent structure that is independent of the neighboring nucleotides. Communicated by Ramaswamy H. Sarma.

Predicting Monovalent Ion Correlation Effects in Nucleic Acids

Ion correlation and fluctuation can play a potentially significant role in metal ion-nucleic acid interactions. Previous studies have focused on the effects for multivalent cations. However, the correlation and fluctuation effects can be important also for monovalent cations around the nucleic acid surface. Here, we report a model, gMCTBI, that can explicitly treat discrete distributions of both monovalent and multivalent cations and can account for the correlation and fluctuation effects for the cations in the solution. The gMCTBI model enables investigation of the global ion binding properties as well as the detailed discrete distributions of the bound ions. Accounting for the ion correlation effect for monovalent ions can lead to more accurate predictions, especially in a mixed monovalent and multivalent salt solution, for the number and location of the bound ions. Furthermore, although the monovalent ion-mediated correlation does not show a significant effect on the number of bound ions, the correlation may enhance the accumulation of monovalent ions near the nucleic acid surface and hence affect the ion distribution. The study further reveals novel ion correlation-induced effects in the competition between the different cations around nucleic acids.

Lessons learned in atomistic simulation of double-stranded DNA: Solvation and salt concerns [Article v1.0]

Nucleic acids are highly charged macromolecules sensitive to their surroundings of water, salt, and other biomolecules. Molecular dynamics simulations with accurate biomolecular force fields provide a detailed atomistic view into DNA and RNA that has been useful to study the structure and dynamics of these molecules and their biological relevance. In this work we study the Drew-Dickerson dodecamer duplex with the sequence d(GCGCAATTGCGC)₂ in three different salt concentrations and using different monvalent salt types to detect possible structural influence. Overall, the DNA shows no major structural changes regardless of amount or type of monovalent ions used. Our results show that only at very high salt conditions (5M) is a small structural effect observed in the DNA duplex, which mainly consist of narrowing of the grooves due to increased residence of ions. We also present the importance of sampling time to achieve a converged ensemble, which is of major relevance in any simulation to avoid biased or non-meaningful results.

5-Formylcytosine mediated DNA-protein cross-links block DNA replication and induce mutations in human cells

5-Formylcytosine (5fC) is an epigenetic DNA modification introduced via TET protein-mediated oxidation of 5-methyl-dC. We recently reported that 5fC form reversible DNA–protein conjugates (DPCs) with histone proteins in living cells (Ji et al. (2017) Angew. Chem. Int. Ed., 56:14130–14134). We now examined the effects of 5fC mediated DPCs on DNA replication. Synthetic DNA duplexes containing site-specific DPCs between 5fC and lysine-containing proteins and peptides were subjected to primer extension experiments in the presence of human translesion synthesis DNA polymerases η and κ. We found that DPCs containing histones H2A or H4 completely inhibited DNA replication, but the replication block was removed when the proteins were subjected to proteolytic digestion. Cross-links to 11-mer or 31-mer peptides were bypassed by both polymerases in an error-prone manner, inducing targeted C→T transitions and –1 deletions. Similar types of mutations were observed when plasmids containing 5fC-peptide cross-links were replicated in human embryonic kidney (HEK) 293T cells. Molecular simulations of the 11-mer peptide-dC cross-links bound to human polymerases η and κ revealed that the peptide fits well on the DNA major groove side, and the modified dC forms a stable mismatch with incoming dATP via wobble base pairing in the polymerase active site.

Improving the Performance of the Amber RNA Force Field by Tuning the Hydrogen-Bonding Interactions

Molecular dynamics (MD) simulations became a leading tool for investigation of structural dynamics of nucleic acids. Despite recent efforts to improve the empirical potentials (force fields, ffs), RNA ffs have persisting deficiencies, which hamper their utilization in quantitatively accurate simulations. Previous studies have shown that at least two salient problems contribute to difficulties in the description of free-energy landscapes of small RNA motifs: (i) excessive stabilization of the unfolded single-stranded RNA ensemble by intramolecular base-phosphate and sugar-phosphate interactions and (ii) destabilization of the native folded state by underestimation of stability of base pairing. Here, we introduce a general ff term (gHBfix) that can selectively fine-tune nonbonding interaction terms in RNA ffs, in particular, the H bonds. The gHBfix potential affects the pairwise interactions between all possible pairs of the specific atom types, while all other interactions remain intact; i.e., it is not a structure-based model. In order to probe the ability of the gHBfix potential to refine the ff nonbonded terms, we performed an extensive set of folding simulations of RNA tetranucleotides and tetraloops. On the basis of these data, we propose particular gHBfix parameters to modify the AMBER RNA ff. The suggested parametrization significantly improves the agreement between experimental data and the simulation conformational ensembles, although our current ff version still remains far from being flawless. While attempts to tune the RNA ffs by conventional reparametrizations of dihedral potentials or nonbonded terms can lead to major undesired side effects, as we demonstrate for some recently published ffs, gHBfix has a clear promising potential to improve the ff performance while avoiding introduction of major new imbalances.

Hexahydrated Mg2+ Binding and Outer-Shell Dehydration on RNA Surface

The interaction between metal ions, especially Mg²⁺ ions, and RNA plays a critical role in RNA folding. Upon binding to RNA, a metal ion that is fully hydrated in bulk solvent can become dehydrated. Here we use molecular dynamics simulation to investigate the dehydration of bound hexahydrated Mg²⁺ ions. We find that a hydrated Mg²⁺ ion in the RNA groove region can involve significant dehydration in the outer hydration shell. The first or innermost hydration shell of the Mg²⁺ ion, however, is retained during the simulation because of the strong ion-water electrostatic attraction. As a result, water-mediated hydrogen bonding remains an important form for Mg²⁺-RNA interaction. Analysis for ions at different binding sites shows that the most pronounced water deficiency relative to the fully hydrated state occurs at a radial distance of around 11 Å from the center of the ion. Based on the independent 200 ns molecular dynamics simulations for three different RNA structures (Protein Data Bank: 1TRA, 2TPK, and 437D), we find that Mg²⁺ ions overwhelmingly dominate over monovalent ions such as Na⁺ and K⁺ in ion-RNA binding. Furthermore, application of the free energy perturbation method leads to a quantitative relationship between the Mg²⁺ dehydration free energy and the local structural environment. We find that ΔΔG_hyd, the change of the Mg²⁺ hydration free energy upon binding to RNA, varies linearly with the inverse distance between the Mg²⁺ ion and the nearby nonbridging oxygen atoms of the phosphate groups, and ΔΔG_hyd can reach -2.0 kcal/mol and -3.0 kcal/mol for an Mg²⁺ ion bound to the surface and to the groove interior, respectively. In addition, the computation results in an analytical formula for the hydration ratio as a function of the average inverse Mg²⁺-O distance. The results here might be useful for further quantitative investigations of ion-RNA interactions in RNA folding.

Complex between a Multicrossover DNA Nanostructure, PX-DNA, and T7 Endonuclease I

Paranemic crossover DNA (PX-DNA) is a four-stranded multicrossover structure that has been implicated in recombination-independent recognition of homology. Although existing evidence has suggested that PX is the DNA motif in homologous pairing (HP), this conclusion remains ambiguous. Further investigation is needed but will require development of new tools. Here, we report characterization of the complex between PX-DNA and T7 endonuclease I (T7endoI), a junction-resolving protein that could serve as the prototype of an anti-PX ligand (a critical prerequisite for the future development of such tools). Specifically, nuclease-inactive T7endoI was produced and its ability to bind to PX-DNA was analyzed using a gel retardation assay. The molar ratio of PX to T7endoI was determined using gel electrophoresis and confirmed by the Hill equation. Hydroxyl radical footprinting of T7endoI on PX-DNA is used to verify the positive interaction between PX and T7endoI and to provide insight into the binding region. Cleavage of PX-DNA by wild-type T7endoI produces DNA fragments, which were used to identify the interacting sites on PX for T7endoI and led to a computational model of their interaction. Altogether, this study has identified a stable complex of PX-DNA and T7endoI and lays the foundation for engineering an anti-PX ligand, which can potentially assist in the study of molecular mechanisms for HP at an advanced level.

Challenges in the use of atomistic simulations to predict solubilities of drug-like molecules

Background: Solubility is a physical property of high importance to the pharmaceutical industry, the prediction of which for potential drugs has so far been a hard task. We attempted to predict the solubility of acetylsalicylic acid (ASA) by estimating the absolute chemical potentials of its most stable polymorph and of solutions with different concentrations of the drug molecule. Methods: Chemical potentials were estimated from all-atom molecular dynamics simulations.  We used the Einstein molecule method (EMM) to predict the absolute chemical potential of the solid and solvation free energy calculations to predict the excess chemical potentials of the liquid-phase systems. Results: Reliable estimations of the chemical potentials for the solid and for a single ASA molecule using the EMM required an extremely large number of intermediate states for the free energy calculations, meaning that the calculations were extremely demanding computationally. Despite the computational cost, however, the computed value did not agree well with the experimental value, potentially due to limitations with the underlying energy model. Perhaps better values could be obtained with a better energy model; however, it seems likely computational cost may remain a limiting factor for use of this particular approach to solubility estimation.    Conclusions: Solubility prediction of drug-like solids remains computationally challenging, and it appears that both the underlying energy model and the computational approach applied may need improvement before the approach is suitable for routine use.

Modeling and Simulation of Oligonucleotide Hybrids: Outlining a Strategy

Molecular dynamics simulations with a state-of-the-art force field provide an atomistic detailed description of the structural and thermodynamic features of biomolecules. Effects of chemical modifications and of the environment such as sequence, solvent, and ionic strength can explicitly be taken into account. Molecular simulation techniques can also provide insight in change in binding affinity, in protonation (pK_a shift) and tautomeric propensity due to changes in the environment or in the molecular system. The quality and reliability of a simulation depend on the quality of the force field and on the reproducibility of the data, and validation depends on the availability of suitable experimental data. Here, we describe the workflow to investigate oligonucleotide hybrids using molecular simulation including hardware and software information.

Twist and turn: a revised structural view on the unpaired bubble of class II CPD photolyase in complex with damaged DNA

Cyclobutane pyrimidine dimer (CPD) photolyases harness the energy of blue light to repair UV-induced DNA CPDs. Upon binding, CPD photolyases cause the photodamage to flip out of the duplex DNA and into the catalytic site of the enzyme. This process, called base-flipping, induces a kink in the DNA, as well as an unpaired bubble, which are stabilized by a network of protein-nucleic acid interactions. Previously, several co-crystal structures have been reported in which the binding mode of CPD photolyases has been studied in detail. However, in all cases the internucleoside linkage of the photodamage site was a chemically synthesized formacetal analogue and not the natural phosphodiester. Here, the first crystal structure and conformational analysis via molecular-dynamics simulations of a class II CPD photolyase in complex with photodamaged DNA that contains a natural cyclobutane pyrimidine dimer with an intra-lesion phosphodiester linkage are presented. It is concluded that a highly conserved bubble-intruding region (BIR) mediates stabilization of the open form of CPD DNA when complexed with class II CPD photolyases.

Molecular Modeling Applied to Nucleic Acid-Based Molecule Development

Molecular modeling by means of docking and molecular dynamics (MD) has become an integral part of early drug discovery projects, enabling the screening and enrichment of large libraries of small molecules. In the past decades, special emphasis was drawn to nucleic acid (NA)-based molecules in the fields of therapy, diagnosis, and drug delivery. Research has increased dramatically with the advent of the SELEX (systematic evolution of ligands by exponential enrichment) technique, which results in single-stranded DNA or RNA sequences that bind with high affinity and specificity to their targets. Herein, we discuss the role and contribution of docking and MD to the development and optimization of new nucleic acid-based molecules. This review focuses on the different approaches currently available for molecular modeling applied to NA interaction with proteins. We discuss topics ranging from structure prediction to docking and MD, highlighting their main advantages and limitations and the influence of flexibility on their calculations.

Characterization of Nine Cancer-Associated Variants in Human DNA Polymerase κ

Specialized DNA damage-bypass Y-family DNA polymerases contribute to cancer prevention by providing cellular tolerance to DNA damage that can lead to mutations and contribute to cancer progression by increasing genomic instability. Y-family polymerases can also bypass DNA adducts caused by chemotherapy agents. One of the four human Y-family DNA polymerases, DNA polymerase (pol) κ, has been shown to be specific for bypass of minor groove adducts and inhibited by major groove adducts. In addition, mutations in the gene encoding pol κ are associated with different types of cancers as well as with chemotherapy responses. We characterized nine variants of pol κ whose identity was inferred from cancer-associated single nucleotide polymorphisms for polymerization activity on undamaged and damaged DNA, their abilities to extend from mismatched or damaged base pairs at primer termini, and overall stability and dynamics. We find that these pol κ variants generally fall into three categories: similar activity to wild-type (WT) pol κ (L21F, I39T, P169T, F192C, and E292K), more active than WT pol κ (S423R), and less active than pol κ (R219I, R298H, and Y432S). Of these, only pol κ variants R298H and Y432S had markedly reduced thermal stability. Molecular dynamics (MD) simulations with undamaged DNA revealed that the active variant F192C and more active variant S423R with either correct or incorrect incoming nucleotide mimic WT pol κ with the correct incoming nucleotide, whereas the less active variants R219I, R298H, and Y432S with the correct incoming nucleotide mimic WT pol κ with the incorrect incoming nucleotide. Thus, the observations from MD simulations suggest a possible explanation for the observed experimental results that pol κ adopts specific active and inactive conformations that depend on both the protein variant and the identity of the DNA adduct.

Ionic liquids make DNA rigid

Persistence length of double-stranded DNA (dsDNA) is known to decrease with an increase in ionic concentration of the solution. In contrast to this, here we show that the persistence length of dsDNA increases dramatically as a function of ionic liquid (IL) concentration. Using all atom explicit solvent molecular dynamics simulations and theoretical models, we present, for the first time, a systematic study to determine the mechanical properties of dsDNA in various hydrated ILs at different concentrations. We find that dsDNA in 50 wt % ILs have lower persistence length and stretch modulus in comparison to 80 wt % ILs. We further observe that both the persistence length and stretch modulus of dsDNA increase as we increase the concentration of ILs. The present trend of the stretch modulus and persistence length of dsDNA with IL concentration supports the predictions of the macroscopic elastic theory, in contrast to the behavior exhibited by dsDNA in monovalent salt. Our study further suggests the preferable ILs that can be used for maintaining DNA stability during long-term storage.