**6. Interaction of molecules with DNA**

196 Bioinformatics

2004].

**5. Selective formation of complexes** 

ligands with this purpose [Kaye, 2011].

31G\*\* calculation of atomic charges [Sakakura et al., 2011].

There are several examples of molecular modelling studies on complexes between cyclic receptors and ammonium ions, calixarenes [Choe & Chang, 2002] and crown ethers being the most used. As an example, it is noteworthy to mention the theoretical studies on calix[4]crown-5 and a series of alkyl ammonium ions [Park et al., 2007], having shown that the energy of complex formation depends on the number of amine groups in the alkyl chain as well as on the number of methylene groups between the primary and secondary amine groups, results that agree with experimental measurements. Although the calculations are performed under quite different conditions of vacuum compared with the experimental conditions of the phase system of chloroform-water, the binding properties of calixarenetype compounds towards alkyl ammonium ions have been successfully simulated,

providing general and useful explanations for the molecular recognition behaviour.

Complex formation of compounds containing benzene rings with ammonium cations has also been theoretically studied using many computational techniques, including *ab initio* calculations [Kim et al., 2000]. It has been shown that two types of NH-aromatic π and CHaromatic π interactions, which are important in biological systems, are responsible for binding, and that charged hydrogen bonds versus cation-π interaction is the origin of the high affinity and selectivity of novel receptors for NH4+ over K+ ions [Oh et al., 2000]. Organic molecules complexed with metal cations have also been studied by MM2 molecular modelling [Mishra, 2010]. The search for metal ion selectivity is of interest in the field of biomimetic models of metalloenzymes and molecular modelling helps in the design of new

Molecular modelling has been used to suggest possible contributions of carrier effectivity and selectivity to complex formation in accordance with experimental results [Chipot et al., 1996; Ilioudis et al., 2005]. Our research group has evaluated the possible cation-receptor interactions involved in the complexes with ammonium and metal cations of selective carriers using the Amber force field with appropriate parameters developed by us. The complexation energies obtained are in reasonable agreement with experimental values, taking into account that complexation/decomplexation processes have a great influence on transport rates and are not equally favoured in cyclic and acyclic carriers [Campayo et al.,

Both binding and selectivity in binding can be understood through the combined efforts of several non-covalent interactions, such as hydrogen bonding, electrostatic interactions, hydrophobic interactions, cation-π interactions, π-π stacking interactions and steric complementarity [Späth & König, 2010]. Formation of complexes is also possible in the case of neutral ligands. For example, the interactions between cholesterol and cyclodextrins have been theoretically studied to investigate their 1:1 and 1:2 complexes [Castagne et al., 2010], while the formation of stable complexes between trehalose and benzene compounds have been investigated by the general Amber force field (GAFF) and Gaussian 03 for MP2/6G-

Anthracycline antibiotics such as doxorubicin and its analogues have been in common use as anticancer drugs for around half a century. There has been intense interest in the DNAbinding sequence specificity of these compounds in recent years, with the hope of identifying a compound that can modulate gene expression or exhibit reduced toxicity. Cashman and Kellog have studied models of binding for doxorubicin and derivatives [Cashman & Kellog, 2004], looking for sequence specificity and the effects of adding aromatic or aliphatic ring substituents or additional amino or hydroxyl groups. They performed a hydropathic interaction analysis using the HINT program (a Sybyl program module, Tripos Inc.) and four double base pair combinations. Interaction of some intercalators with two double DNA base pairs have also been studied with the density functional based tight binding (DFTB) method [Riahi et al., 2010], despite DFT methods being known to be inherently deficient in calculating stacking interactions, and the Amber force field and then AM1 to dock the intercalator between DNA base pairs [Miri et al., 2004].

Studies on sequence-selectivity of DNA minor groove binding ligands have shown that the most reliable results for AT-rich DNA sequences are obtained when MD simulations are performed in explicit solvent, when the data are processed using the MM-PB/SA approach, and when normal mode analysis is used to estimate configurational entropy changes [Shaikh et al., 2004; Wang & Laughton, 2009]. Use of the GB/SE model with a suitable choice of parameters adequately reproduces the structural and dynamic characteristics in explicitly solvated simulations in approximately a quarter of the computational time, although limitations become apparent when the thermodynamic properties are evaluated [Sands & Laughton, 2004]. Water molecules taking part in the complexation have been studied using

the MMX force field and then PMB for gas phase optimisation, followed by re-optimisation in aqueous phase with the PM3 method using the AMSOL package [Silva & Jayasundera, 2002]. Optimisation geometries with AM1 and the use of implicit solvent have been taken into account when considering intercalation versus insertion into the minor or major groove [Bendic & Volanschi, 2006].

Using Molecular Modelling to Study Interactions Between Molecules with Biological Activity 199

ΔU = EinterES + ΔES + ΔEE + EdesolvS + EdesolvE (3)

estimate ligand binding affinities and selectivities, with various levels of accuracy and computational expense: free energy perturbation (FEP), thermodynamic integration (TI), lineal response (LR), and molecular mechanics Poisson-Boltzman surface area (MM/PBSA). Identification of conformational preferences and binding site residues, as well as structural and energetic characterisation, is possible using MD simulations [Anzini et al., 2011; Dastidor et al., 2008; Stoika et al., 2008]. It is also possible to estimate conformational energy penalties for adopting the bioactive conformation identified by using a pharmacophore

A model based on van der Waals intermolecular contribution from Amber and electrostatic interactions derived from the Poisson-Boltzman equation has been used to predict the change in the apparent dissociation constant for a series of six enzyme-substrate complexes during COMBINE analysis [Kmunicek et al., 2001]. In COMBINE analysis, binding energies are calculated for the set of enzyme-substrate complexes using the molecular mechanics force field. The total binding energy, ΔU, may be assumed to be the sum of five terms: the intermolecular interaction energies between the substrate and each enzyme residue, EinterES, the change in the intramolecular energy of the substrate upon binding to the enzyme, ΔES, the change in the intramolecular energy of the enzyme upon binding, ΔEE, the desolvation

When the substrate is a rather small molecule, there is no evidence for large differences in the structure of the enzyme when different substrates are bound and so the second and third are neglected. This method identifies the amino acid residues responsible for modulating

Molecular modelling of proteins is sometimes directed towards homology modelling, enabling progress in understanding the mechanisms of action despite the lack of detailed information on the 3D structure of a protein. Molecular dynamic simulations are usually used to test the stability of the complete structure derived from homology modelling

Molecular docking examples can be used to compare relative stabilities of the complexes, but not calculate binding affinities, since changes in entropy and solvation effects are not taken into account [Pastorin et al., 2006; Tschammer et al., 2011]. In any case, docking calculations are common studies on novel drugs, Autodock being one of the most used docking programs [see for example Venskutonyte et al., 2011]. Docking programs treat enzymes and substrates as rigid entities, but flexible docking is also possible, if several different protein conformations extracted from molecular dynamic simulations are used

In our laboratory, molecular modelling has been tentatively used to study the trypanosomicidal activity of some phthalazine derivatives. Results obtained with Amber force field implemented in HyperChem 8.0 plus our own necessary parameters, and with AutoDock 4.2 using the PDB structure for *T. cruzi* Fe-SOD enzyme, were in accordance with

energy of the substrate, EdesolvS, and the desolvation energy of the enzyme, EdesolvE.

model [Frølund et al., 2005].

enzyme activity [Kmunicek et al., 2005].

[Srinivas et al., 2006].

[Roumen et al., 2010].

Calculating the curvature radius of molecular DNA structures has been reported [Slickers et al., 1998] as a new method for understanding the dependence of binding affinity on ligand structure, assuming that strong binders should have a shape complementary to the DNA minor groove. A method for predicting sequence selectivity and minor groove binding, based on MD simulations on DNA sequences with and without the bound ligand, to obtain an approximate free energy of binding has been proposed [Wang & Laughton, 2010].

Amber force field, developing the necessary parameters, has also been used together with electrostatic potential-derived (ESP) charges and explicit solvent molecules to study bisintercalation into DNA. The targeted molecular dynamics (tMD) approach has been considered for comparing the relative energetic cost involved in creating the intercalation sites and also studying the mechanisms of action [Braña et al., 2004]. It has been found that the electrostatic contribution is a critical characteristic of binding selectivity [Marco et al., 2005]. Reports on duplex and triplex formation of oligonucleotides by stacking aromatic moieties in the major groove, using Amber force field and the GB/SA solvation model in molecular dynamic simulations, can be found in the literature [Andersen et al., 2011]. Studies on docking using GOLD [Kiselev et al., 2010] to optimise the starting structures with the MMFF94 force field have also been performed.

Most of the published molecular modelling studies use two double base pairs or more than eight double base pairs to represent DNA. In our opinion, molecular modelling of DNA intercalation complexes should be done using at least the two base pairs of the intercalation site and an additional base pair at the two strand ends to maintain DNA shape and avoid distortion leading to inaccurate results. That means four base pairs for monointercalation studies and five or six base pairs for bisintercalation ones should be used. Using these DNA models, our studies on the mono and bisintercalation of benzo[*g*]phthalazine derivatives strongly suggest the possibility of bisintercalation and the important role played by an Nmethyl group in stabilising the DNA complex of one of the compounds, throwing some light over the experimental results obtained [Rodríguez-Ciria et al., 2003]. The possibility of bisintercalation for a 1,4-disubstituted piperazine has been studied on duplexes of five and six base pairs, obtaining much better results in the case of five base pairs, in accordance with the theoretical calculations of binding mode not conforming to the neighbouring exclusion principle proposed by different authors [Veal et al., 1990].
