Ancillary Ligand in Ternary CuII Complexes Guides Binding Selectivity toward Minor-Groove DNA.

Copper-containing compounds known as Casiopeı́nas are biologically active molecules which show promising antineoplastic effects against several cancer types. Two possible hypotheses regarding the mode of action of the Casiopeı́nas have emerged from the experimental evidence: the generation of reactive oxygen species or the ability of the compounds to bind and interact with nucleic acids. Using robust molecular dynamics simulations, we investigate the interaction of four different Casiopeı́nas with the DNA duplex d(GCACGAACGAACGAACGC). The studied copper complexes contain either 4-7- or 5-6-substituted dimethyl phenanthroline as the primary ligand and either glycinate or acetylacetonate as the secondary ligand. For statistical significance and to reduce bias in the simulations, four molecules of each copper compound were manually placed at a distance of 10 Å away from the DNA and 20 independent molecular dynamics simulations were performed, each reaching at least 30 μs. This time scale allows us to reproduce expected DNA terminal base-pair fraying and also to observe intercalation/base-pair eversion events generated by the compounds interacting with DNA. The results reveal that the secondary ligand is the guide toward the mode of binding between the copper complex and DNA in which glycinate prefers minor-groove binding and acetylacetonate produces base-pair eversion and intercalation. The CuII complexes containing glycinate interact within the DNA minor groove which are stabilized principally by the hydrogen bonds formed between the amino group of the aminoacidate moiety, whereas the compounds with the acetylacetonate do not present a stable network of hydrogen bonds and the ligand interactions enhance DNA breathing dynamics that result in base-pair eversion.

Introduction

The family of copper based transition-metal compounds known as Casiopeínas, noting that copper is essential to our diet, has shown potential as chemotherapeutic drugs[1,2] with one of these compounds, Cas III-ia, currently in clinical phase I studies. With the general formula [Cu(N–N)(N–O)]NO3 and [Cu(N–N)(O–O)]NO3 where N–N denotes nonsubstituted and substituted 2,2′-bipyridine or 1,10-phenanthroline (identified as the primary ligand), N–O denotes α-aminoacidate, and O–O denotes acetylacetonate or salicylaldehydate (identified as the secondary ligand in the Cu-complex); this extensive family of compounds has shown significant biological activity both in vivo and in vitro.[3−6] Primary and secondary tags were assigned to the ligands because of their capacity to modulate the physicochemical properties of the Casiopeínas family of Cu compounds, mainly their redox potential.[7,8] The structural isomers shown in Figure , which belong to the Casiopeínas set, present antiproliferative activities against several human cancer cells.[9,10] Overall, experimental studies show that the biological mechanism of action for these copper complexes is the generation of reactive oxygen species and DNA binding interaction, independently or concomitantly. Experimental evidence of the interaction of copper compounds with nucleic acids was pioneered by Chikira and collaborators using spectroscopy techniques.[11,12] They also showed that copper-containing compounds have binding properties and nuclease activity in DNA.[13,14] In 2015, Bravo-Gómez and collaborators reported a study using six Casiopeínas and their interaction with calf thymus DNA. Using fluorescent displacement assays with two different dyes (ethidium bromide and SYBR Green), they report multiple modes of binding that depend on the aromatic (primary ligand) and secondary (glycinate or acetylacetonate) ligands that are coordinated to the metal center.[15]

Figure 1

Molecular structures of the isomers Cu[(4,7-dimethyl-1,10-phenanthroline)(glycinate)]2+ (Cas II-gly), Cu[(5,6-dimethyl-1,10-phenanthroline) (glycinate)]2+ (Cas VI-gly), Cu[(4,7-dimethyl-1,10-phenanthroline) (acetylacetonate)]2+ (Cas III-Ea), and Cu[(5,6-dimethyl-1,10-phenanthroline) (acetylacetonate)]2+ (Cas III-La).

Previous molecular dynamics (MD) simulations of 21 different compounds from the Casiopeínas family interacting with the Drew-Dickerson dodecamer sequence d(CGCGAATTCGCG)2 suggested five different modes of interaction between the compounds and the DNA.[16] These included (A) stacking on the terminal base pairs of the DNA chain, (B) minor-groove binding, (C) intercalation with base-pair eversion, (D) minor-groove binding with stacking of the terminal frayed bases, and (E) intercalation near the termini of the DNA duplex. Of special interest was the base-pair eversion binding mechanism in which the Cu complex initially binds to the minor groove of the DNA through the coordination of the phosphate’s oxygen atom to the metal center and stabilized by CH···π interactions with the backbone. After this, either the primary or the secondary ligands of the Cu complex push an AT base pair, disrupting the Watson–Crick pairing until eventually it is broken. The Cu complex continues forward into the double helix, pushing both nucleobases toward the major groove and then the compound finally inserts into the resulting cavity. The Cu complex remains in the intercalated position because of stacking interactions and electron depletion of the planar ligand because of charge transfer.[17] The described interaction was observed only in four compounds, specifically those that contained acetylacetonate as the secondary ligand, that is, [Cu(4,4-dimethyl-2,2′-bipyridine)(acetylacetonate)(H2O)]+, [Cu(1,10-phenanthroline)(acetylacetonate)(H2O)]+, [Cu(5-methyl-1,10-phenanthroline)(acetylacetonate)(H2O)]+, and [Cu(4,7-diphenyl-1,10-phenanthroline)(acetylacetonate)(H2O)]+ (noting that the names as they appear in the mentioned article are cas02, cas03, cas05α, and cas10, respectively).[16] This mode of interaction has been observed experimentally with rhodium and ruthenium compounds by the Barton group.[18] Additionally, Monari and collaborators also observed a single-base eversion binding mode using MD simulations of a palmatine photosensitizer interacting with a 16-mer[19] and a base-pair eversion mode between benzophenone interacting with a DNA 10-mer.[20]

The present work focuses on the DNA-binding properties of four distinct members of the Casiopeínas family (Figure ) to further elucidate the nature of (a) the formation of a DNA–copper-compound complex, (b) the detection of single or multiple modes of binding, (c) the observation of possible DNA-sequence selectivity, and (d) the influence of the aromatic and nonaromatic ancillary ligand and its possible role in the formation of the DNA–Cu complex. In the current work, a longer 18-mer double-stranded sequence, d(CGACGAACGAACGAACGC), was utilized to allow greater exploration of the nonterminal elements of the duplex (noting that interactions with the base pairs at each end of the duplex were the most commonly observed type of interaction in our previous work[16]). The structure and dynamics of this particular GAAC sequence have been extensively studied in converged and reproducible control sets of MD simulations by our group[21,22] and provides a proper balance of purines and pyrimidines that will allow us to discover possible sequence selectivity. The initial configuration for all the simulations consists of four ligand molecules at ∼10 Å from the geometric DNA center (Figure ).

Figure 2

Starting configuration structure for all systems. Four of the respective Cu compounds are manually placed around base pairs G9–C28 at a distance of 10 Å of the 18-mer DNA chain. In this image, the ligand Cas II-gly is shown around the GAAC sequence. Water molecules and ions, present in the simulations, were omitted for clarity.

Methods

MD Simulations

Structures for Cas II-gly, Cas VI-gly, Cas III-Ea, and Cas III-La were based on the X-ray crystal structures available[8,23] and were described with the general Amber force field[24] using additional parameters from Babu[25] and as previously described.[16,26] The double-stranded DNA sequence d(GCACGAACGAACGAACGC), referred to as GAAC, was generated using the nucleic acid builder.[27] The previously tested OL15[28] force field was used to describe the DNA which was included in a truncated octahedral box using the TIP3P[29] water model. Sodium ions were added to neutralize the charge and an excess of NaCl ions was added to reach ∼150 mM concentration (estimated based on the initial volume) using the Joung–Cheatham ion model.[30] The minimization and equilibration protocol was performed as per our previous work.[16] After initial equilibration of the system, 20 independent copies were run for 30 μs each at 300 K using the Berendsen thermostat for temperature control with a coupling value of 5.0.[31]

MD was performed using the pmemd.cuda GPU code from the AMBER 16 and AMBER 18[32,33] suite of programs. Originally, simulations for Cas II-gly, Cas VI-gly, Cas III-Ea, and Cas III-La were previously performed in an orthorhombic box with only 5 Å of solvent padding (using the same simulation protocols with simulations on the 20 μs time scale). This was done in an attempt to reduce the computational cost required for the simulations. However, upon analysis and visualization of the trajectory, it was observed that the DNA rotated to interact with its periodic image with either one, two, or three of the ligands binding and bridging the two DNA ends effectively creating an infinite stack or aligned DNA duplex (Figure S1). When the periodic images align up, this effectively stabilizes the terminal end binding, reducing the ability of compounds to bind elsewhere. This leads to reduced atomic density at some of the minor-groove-binding sites, which otherwise should all be equivalent because of the repeating sequence (Figure S2). This is due to the effective lowered concentration of available ligands and requisite increased sampling time to sample all the bound modes. This artifact of end-to-end stacking has been seen previously in minimally solvated systems of DNA in orthorhombic boxes where rotation on the nanosecond scale leads to infinite stacking and also in systems containing multiple DNA duplexes (Cheatham, personal communication). In order to eliminate this artifact, MD simulations were performed using a truncated octahedral box of larger diameter (with more water). No subsequent direct periodic interactions of the duplex or double-end trapping of ligands were observed.

Trajectory analysis was performed using CPPTRAJ v18.00.[34,35] Quantum chemical calculations were performed using the D.09 version of Gaussian.[36] To manage the significant volume of MD simulation data, the analysis was performed concatenating each of the 20 individual copies into an aggregated, water-stripped, 600 μs trajectory. The complex of the copper compound with DNA-binding energies were computed using the MM-PBSA methodology[37] using MMPBSA.py.[38] Scripts for performing the analyses are provided in the Supporting Information. Grid density analysis was performed first orienting the DNA and ligands into a common reference frame based on RMS fitting the DNA with a 0–4 μs average structure, and then, a 90 × 90 × 145 grid with 0.5 Å spacing was superimposed on the DNA and the atomic density for the ligands, independently, are binned into the grid elements and visualized at a density isosurface value of 30 molecules/Å3 using VMD.[39]

Results and Discussion

Overall, we have produced 2.4 ms of MD simulation-generated trajectory data, which required careful considerations to process and manage the data. Our first question is that where do the Cu complexes find a stable binding site within the DNA? To better visualize the general binding of the molecules to the GAAC sequence, a 3D molecular grid atomic density analysis is shown in Figure for each of the four molecules present in the simulation. What this does is put the DNA into a common reference frame where for each ligand its atom density is binned into a grid throughout the MD trajectory and the grid is visualized at various isocontours of density. Our previous work showed five types of binding interactions of these compounds with a 12-mer DNA: (a) stacking of the compound on the terminal base pairs, (b) minor-groove binding, (c) base-pair eversion, (d) minor-groove binding with stacking of one of the terminal frayed bases, and (e) intercalation near the end of the DNA chain (Figure of ref (10)). It is important to remember that each of the 20 simulations for each of the four tested molecules (Cas II-gly, Cas VI-gly, Cas III-Ea, and Cas III-La) may have four Cu complexes actively interacting with the DNA. For compounds Cas II-gly and Cas VI-gly, which have glycinate as the secondary ligand, we observe Cu-complex accumulation forming stacking interactions at the terminal base pairs, minor-groove binding, and a few very short-lived intercalation events. Stacking interactions at the DNA terminal base pair were also observed in the grid density visualizations for the compounds Cas III-Ea and Cas III-La that have acetylacetonate as the secondary ligand; however, an increase in intercalation events and minimal binding to the DNA minor groove were found with these compounds compared with Cas II-gly and Cas VI-gly. Terminal base-pair stacking was the most common binding mode in our previous studies and the main motivation to include four copper complexes in this study. Two of the Cu compounds interact as expected with DNA through terminal base-pair stacking leaving the other two Cu compounds the freedom to explore other binding modes. In the case of Cas II-gly and Cas VI-gly, the stacking at the terminal base pairs stabilizes the WC pairing and reduces fraying events, which is observed in the RMS fluctuation plot of Figure . There is a slight increase in the terminal base-pair fluctuations for Cas III-Ea and Cas III-La which suggests a weaker interaction in this particular binding mode. Some degree of sequence selectivity can be inferred by the grid density analysis, pointing toward a preference of the Cu compounds to interact within the A:T regions of the GAAC sequence. This could be explained by the fact that A:T base pairs located at the floor of the minor groove readily attract the cationic Cu complexes because of the greatest negative potential compared with other nucleobase pairs.[40,41] Further evidence of the A:T preference is observed when the grid density analysis for each of the ligands is overlaid on top of each other. From this overlay representation, we observe terminal base-pair stacking and minor-groove binding for both Cas II-gly and Cas VI-gly. Cas II-gly shows some evidence of intercalation at the A17-T22 region and a more diffuse binding through the minor groove. Cas VI-gly shows a more localized minor-groove binding and evidence of multiple eversion/intercalation within the A11-T26 base-pair region. These results agree with the proposed mechanism of action for Cu(phenanthroline)2 derivatives by hydrogen abstraction from C1′, C4′, and C5′ of DNA’s deoxyribose located within the minor groove.[42] The oxidative attack by Cu(phenanthroline)2 derivatives is favored by the interaction with the minor-groove floor by intercalation. The fact that multiple binding modes are observed despite only effectively two free compounds suggests that the binding is reversible and/or statistically sampled in this set of MD simulations. Cas II-gly and Cas VI-gly also show tight binding through the minor groove, although multiple events where the Cu compound is between base steps (not necessarily intercalated) in several regions of the GAAC sequence are evident and will be discussed in detail.

Figure 3

Grid atomic density visualization for each of the copper compounds present at the start of the simulations represented with the same isodensity value (using a solid representation in each case). Each column represents the grid analysis based on one of the four Cu complexes present in the simulation. The overlay structure is with the volume data for each Cu complex superimposed (using a mesh representation). The CPPTRAJ analysis script used is available in the Supporting Information, and the image was generated with the isosurface plugin available in VMD.[39]

Figure 4

Normalized RMSD population using a copper-compound-free GAAC sequence (black). The solid line is considering all residues, the dashed line represents inner residues (3–16 and 21–34). All calculations used a 3.8 μs average structure as a reference. Right, root-mean-square fluctuations.

Structural Influence of the Copper Complexes to the DNA

From Figure , we can infer the overall structural deviation that the presence of the Cu compounds has on the GAAC sequence (root mean square population). The control simulation with no compounds has a maximum in the root-mean-square deviation (RMSD) histogram of ∼2.3 Å. When the Cu compounds with glycinate are present, the deviation is ∼2.2 Å, suggesting less deviation from the reference structure which is due to the inhibition of terminal base-pair fraying because of the ligand binding at the termini and the extended interaction of the molecules in the minor groove, which in turn, slightly increases the values of the helical twist (Table ). This increased twist and reduction of the minor-groove width due to the presence of small molecules have been observed previously in similar MD simulations.[43] The slight increase in fluctuations around residues 11 and 25 for Cas VI-gly is due to one simulation in which Cas VI-gly has an intercalation event. Although unexpected, this observed event explains the increase in fluctuations that can be confirmed via C1–C1 distance analysis (Figure S3). In contrast, the presence of the Cu compounds with acetylacetonate, Cas III-Ea and Cas III-La, shows an increased number of distorted structures with higher RMSD values (∼1–1.5 Å from the reference structure) and greater internal base-pair fluctuations because of the base-pair eversion and intercalation events.

Table 1

Average Intra- and Interhelical Parameters for the Studied Systems (Values in Angstroms and Degrees)a

	reference	Cas III-Ea	Cas III-La	Cas II-gly	Cas VI-gly
shear	0.05	–0.05 ± 4	0.00 ± 1	0.06 ± 0.5	0.02 ± 0.8
stretch	0.03	–0.24 ± 4	0.02 ± 0.9	0.01 ± 0.4	0.02 ± 0.6
stagger	0.02	0.13 ± 4	0.08 ± 1	0.06 ± 0.5	0.09 ± 0.7
buckle	3.75	3.80 ± 26	4.10 ± 14	3.50 ± 12	3.60 ± 14
propeller	–11.0	–10.4 ± 28	–10.6 ± 13	–10.3 ± 10	–10.1 ± 11
opening	2.04	1.40 ± 23	1.90 ± 10	2.40 ± 6	2.50 ± 8
X-displacement	–0.53	–0.32 ± 1	–0.36 ± 0.9	–0.20 ± 0.8	–0.08 ± 0.9
Y-displacement	0.05	0.04 ± 1	0.02 ± 0.7	0.01 ± 0.6	0.02 ± 0.6
inclination	5.6	4.6 ± 21	4.5 ± 8	2.9 ± 6	2.5 ± 8
tip	0.7	0.4 ± 29	0.4 ± 14	0.5 ± 9	0.4 ± 11
axial bend	1.6	3.1 ± 5	1.9 ± 1	1.9 ± 1	2.0 ± 1
shift	–0.05	–0.06 ± 1	–0.02 ± 1	–0.02 ± 0.9	–0.02 ± 1
slide	–0.01	0.11 ± 1	0.03 ± 0.8	0.02 ± 0.7	0.07 ± 0.7
rise	3.29	3.26 ± 1	3.37 ± 0.6	3.33 ± 0.4	3.38 ± 0.6
tilt	–0.42	–0.20 ± 22	–0.20 ± 9	–0.30 ± 6	–0.20 ± 9
roll	3.4	2.5 ± 24	2.8 ± 12	1.7 ± 8	1.5 ± 10
twist	34.6	33.1 ± 22	34.5 ± 8	35.3 ± 7	35.0 ± 7
minor groove	6.7	6.2 ± 2	6.1 ± 2	5.5 ± 2	5.6 ± 1
major groove	10.9	11.4 ± 2	11.6 ± 2	11.6 ± 1	11.7 ± 2

Reference values are from 60 μs GAAC simulation with no Cu compound. For the DNA–Cu compound complex, the entire aggregated trajectory of 400 μs was used. Only the internal base pairs are considered in all the cases (residues 3–16 and 21–34).

Cas III-Ea and Cas II-La intercalation events generate a slight decrease in the helical twist, unwinding the DNA (Table ). The decrease in the average helicoidal value is consistent with experimental observations of known intercalator agents such as ethidium bromide.[44] Visual inspection of the simulations with these Cu compounds presents multiple long-lived base-pair eversion events as the main source of structural deviation. A selected example of the eversion binding mode is depicted in Figure . The process starts (Figure A) with the Cu compound interacting with the DNA backbone, mainly through interactions between the copper atom and the deoxyribose O4′. We have previously observed this interaction between Cu–phosphate as the first step of Casiopeína–DNA binding formation;[2,26,45] as the simulation progresses, the Cu compound is “locked” within the minor groove using the weak interaction between the Cu center and the O4′ of the deoxyribose present in the n-1 position as a pivot, allowing the Cu compound to dynamically move within the walls of the minor groove (Figure B,C). During the course of the simulation, the DNA breathing causes spontaneous, short-lived events where the opening of the WC pairing between A:T is increased. This small A:T base-pair opening increase, combined with the presence of the Cu compound within the minor groove, pushes both dA and dT toward the major groove (Figure C,D). The opening allows the Cu compound to move inside the hydrophobic pocket to form a strong interaction with the top and bottom base pairs from the resulting cavity, strongly stabilized by π···π stacking interactions formed with the primary ligand (Figure D,E). Alternative configurations of the base-pair eversion binding mode are presented in Figure S4.

Figure 5

Selected frames showing the base-pair eversion mechanism between Cas III-Ea and a G:C pair (top) and an A:T pair (bottom). The Cu compound is depicted in green; hydrogen atoms have been hidden for clarity.

Evidence from multiple sets of experiments with aromatic planar molecules reveals that intercalation causes DNA unwinding because of the separation of the base pairs to make room for the intercalating molecule.[46] Considering this, we can observe that for the Cu compounds that have acetylacetonate, there is a small decrease in the twist values (Table ) suggesting unwinding of the DNA because of the presence of the Cu compound in the cavity formed by the A:T base pair. The previously mentioned result contrasts with the results obtained with Cu compounds that have glycinate, whose binding within the minor groove reduces the minor-groove width and slightly increases the twist value.

H-Bonding Analysis

For each of the four Cu compounds present in the simulation, a hydrogen-bond analysis was performed (see the Supporting Information. See CPPTRAJ script 3 in the Supporting Information). This information is presented in Table which shows the atoms involved in producing a hydrogen bond, as well as the fraction of time that the interaction persists while the molecule stays in that particular binding site. It is clear that both Cas II-gly and Cas VI-gly present similar interaction profiles with the GAAC sequence, as already observed visually in the previous grid atomic density analysis. These two compounds show consistent hydrogen bonds with the pyrimidine O2 atoms, which are accessible at the floor of the minor groove of DNA, that provide an anchor point on which to form hydrogen bonds that are being recognized by the amino group. For Cas II-gly, the interaction is through the H7 hydrogen which belongs to the glycinate amino group; this happens also with the Cas VI-gly H11 hydrogen from the amino part. Overall, Cas II-gly and Cas VI-gly interact with cytosine O2 through the glycinate moiety and, to a lesser extent, with thymine O2 and adenine N3. The same happens with the N3 nitrogen atom of the adenine base, although in lower population. Contrary to the robust hydrogen-bond matrix formed by the glycinate derivatives, compounds with the acetylacetonate present low populations of H-bond-type interactions within the minor groove. It is important to mention that in order to classify two atoms as an H bond, we have considered a simple geometric criterion: any contact within 3.0 Å length and inside 135° will be marked as a hydrogen bond. With this criterion, it readily detects bonds between the GAAC sequence and Cas II-gly and Cas VI-gly, whereas for Cas III-Ea and Cas III-La, only detects “possible” candidates that are very short-lived and are between C–H···H atoms. Nevertheless, H···H interactions should not be discarded as a stabilizing contributor because it has proven to be very important in different crystal systems, mainly of nonpolar molecules.[47]

Table 2

Hydrogen-Bond Analysis Using the Entire Aggregated Trajectoriesa

Cas III-Ea				Cas II-gly
acceptor	donorH	donor	%	acceptor	donorH	donor	%
Compound 1				Compound 1
DG_33@H8	CAS_37@H3	CAS_37@C4	0.02	DA_7@N3	CAS_37@H7	CAS_37@N8	3.7
DC_32@H41	CAS_37@H2	CAS_37@C4	0.02	DT_23@O2	CAS_37@H7	CAS_37@N8	2.6
DT_26@H1′	CAS_37@H10	CAS_37@C12	0.02	DC_32@O2	CAS_37@H7	CAS_37@N8	1.9
DG_33@H8	CAS_37@H1	CAS_37@C4	0.02	DC_24@O2	CAS_37@H7	CAS_37@N8	1.8
DC_32@H41	CAS_37@H1	CAS_37@C4	0.02	DT_31@O2	CAS_37@H6	CAS_37@N8	1.7
Compound 2				Compound 2
DG5_1@H5′	CAS_38@H32	CAS_38@C33	0.19	DC_24@O2	CAS_38@H7	CAS_38@N8	6.6
DG_35@H2″	CAS_38@H36	CAS_38@C37	0.03	DT_23@O2	CAS_38@H7	CAS_38@N8	5.5
DC3_18@H2″	CAS_38@H32	CAS_38@C33	0.01	DC_32@O2	CAS_38@H7	CAS_38@N8	3.0
DC3_18@H2″	CAS_38@H5	CAS_38@C8	0.01	DA_15@N3	CAS_38@H7	CAS_38@N8	2.6
DC3_18@H2″	CAS_38@H25	CAS_38@C26	0.01	DT_31@O2	CAS_38@H7	CAS_38@N8	2.3
Compound 3				Compound 3
DT_22@H2″	CAS_39@H36	CAS_39@C37	0.04	DC_28@O2	CAS_39@H7	CAS_39@N8	2.6
DA_15@H2	CAS_39@H25	CAS_39@C26	0.03	DC_24@O2	CAS_39@H7	CAS_39@N8	2.6
DT_22@H5′	CAS_39@H30	CAS_39@C31	0.02	DC_4@O2	CAS_39@H7	CAS_39@N8	2.1
DC3_18@H2″	CAS_39@H36	CAS_39@C37	0.01	DT_23@O2	CAS_39@H7	CAS_39@N8	1.9
DC3_18@H2″	CAS_39@H9	CAS_39@C12	0.01	DT_27@O2	CAS_39@H7	CAS_39@N8	1.9
Compound 4				Compound 4
DT_22@H1′	CAS_40@H9	CAS_40@C12	0.05	DA_7@N3	CAS_40@H7	CAS_40@N8	4.8
DT_22@H1′	CAS_40@H11	CAS_40@C12	0.05	DC_32@O2	CAS_40@H7	CAS_40@N8	3.0
DT_22@H1′	CAS_40@H10	CAS_40@C12	0.05	DT_31@O2	CAS_40@H6	CAS_40@N8	2.3
DT_22@H1′	CAS_40@H15	CAS_40@C16	0.03	DC_4@O2	CAS_40@H7	CAS_40@N8	2.2
DC3_36@H2″	CAS_40@H36	CAS_40@C37	0.02	DT_23@O2	CAS_40@H7	CAS_40@N8	2.2

Each copper compound corresponds to one of the four molecules included in the simulations. The acceptor column refers to the residue number and involved atom, and the donorH column refers to the proton involved in the hydrogen bonding. The Frac column refers to the percentage of total frames in which a particular bond is present. Only the top five populated bonds for each compound are shown. Refer to Figure S7 for atom labels.

The abovementioned information can be used to hypothesize that the mentioned differences in the hydrogen-bonding network between the secondary ligands of the Cu compound promote different interactions with DNA. The “poor” hydrogen-bond network formed by the acetylacetonate derivatives with the DNA minor groove allows the first recognition process but subsequently, searching for a more stable energy situation, the Cu compounds push the base-pair floor of the minor groove to start the eversion process. On the other hand, the glycinate moiety of both Cas II-gly and Cas VI-gly presents a more energetically favored binding environment because of the created hydrogen bonds with the atoms present on the floor of the minor groove; therefore, the molecule is not shifting positions that could potentially produce intercalation or base-pair eversion. Furthermore, it has been observed experimentally that only purine N3, pyrimidine O2, guanine N2, and deoxyribose O4′ are involved in minor-groove interactions[48] which are the target atoms for Cas II-gly and Cas VI-gly to form hydrogen bonding, as quantified by the data presented in Table .

Cu Compounds with Glycinate Show Long-Lived Base-Pair Eversion and Intercalation Events

The increased hydrogen-bond contacts present for Cas II-gly and Cas VI-gly stabilize the compound within the minor groove for most of the sampled time, whereas Cas III-Ea and Cas III-La present multiple events of base-pair eversion and intercalation. To further confirm this behavior, we plotted the distance over time between the C1′ carbons for each base pair (excluding the terminal base pairs) in the GAAC sequence (Figure ) using all 20 copies. The average distance between the C1′–C1′ atoms in a Watson–Crick base pair is 10.1 Å. Values above it correspond to events where the distance between the base pairs is increased because of the base-pair eversion effect. The compound Cas III-La presents events where the C1′–C1′ distance is increased from the canonical value of 10.3 Å to values ranging from ∼14 to 16 Å. In each of those cases, the interaction of the Cas III-La with DNA is similar to the one described previously (Figure ). None of the trajectories showed an event where the Cu compound would exit the cavity toward either of the grooves or to the solvent. From the distance versus time information, we notice that the interaction is effectively irreversible in the sampled time scale. As observed in the grid analysis, the interaction is driven primarily to A:T sites, with two copies out of the 20 interacting with a G:C site, one which is depicted in Figure . Contrary to what is observed for Cas III-La, Cas II-gly base-pair eversion events are quite infrequent. In these short-lived events, the intercalated Cu compound moves back toward the minor groove and the shifted nucleotides move back to reform the Watson–Crick pairing. This behavior is consistent for the four Cu compounds included for each simulation (a similar plot for Cas III-Ea and Cas VI-gly is included, as shown in Figure S3).

Figure 6

Distance (Å) vs time plots between the C1′ atoms for all the base pairs. Each of the rows represents a single base pair in the GAAC sequence having the 5′ termini at the top. Each line in a single row represents an independent copy out of the 20 calculated. Y-Axis represents distance ranging from 9 to 17 Å. Only Cas III-La and Cas II-gly are presented, please refer to the Supporting Information for the Cas III-Ea and Cas VI-gly plots.

Base-Pair Eversion Binding Mode is the Most Energetically Favored

Binding energy analysis for each of the four compounds present in each simulation was analyzed using a subsection of frames extracted from all the sampled data. Four regions are observed representing four stable binding modes (Figure ). A highly populated and highly packed region is observed around −1.0 kcal/mol (region D) that represents the Cu compounds interacting with the GAAC sequence at both ends of the DNA chain. This binding mode probes to be the least energetically favored, although it is highly populated. Region C represents the minor-groove binding of Cu compounds that contain acetylacetonate between −15 and −20 kcal/mol and region B represents the minor-groove binding of Cu compounds that contain glycinate between −29 and −38 kcal/mol. The result suggests that Cas II-gly and Cas VI-gly are more tightly bound to the minor groove of the GAAC sequence because of the presence of the glycinate ligand, providing a large number of contacts that stabilize the interaction. Region A represents all the events on which Cas III-Ea and Cas III-La are interacting with the DNA via base-pair eversion and intercalation events. This mode of binding is then the most energetically stable showing values between −33 and −55 kcal/mol. We hypothesize that extended simulations will increase the population of this region, with the corresponding population decrease in regions B and C. If extended simulations were calculated, the population of C and D should decrease and more population of regions A and B should increase. The energy distribution also suggests that the position of the methyl groups in the primary ligand of the Cu compounds has little influence on the binding mode with the GAAC sequence.

Figure 7

Normalized population of the binding energy distribution. Each line corresponds to one of the four Cu compounds present in each simulation. To speed up the calculation, a total of 2000 frames were extracted using all 20 copies for each system.

Proving an Alternate Binding Position for Cas II-gly and Cas VI-gly

As discussed, one of the main observations of the simulations is the fact that the coordination compound preferably binds to the minor groove of the GAAC sequence. The glycinate part of both Cas II-gly and Cas VI-gly orients the amino group toward the floor of the minor groove, and interactions consisting of hydrogen bonding and CH···π stabilize the DNA–copper compound complex. To test this interaction for possible artifacts, we performed MD simulations starting from one of the structures with the Cu compound inside the minor groove rotated 180° compared with the most stable configuration observed so that the carboxylate anion faces toward the floor of the minor groove. We then used the same MD protocol as described in the methodology section and run five independent copies for each Cas II-gly and Cas VI-gly for ∼1.5 μs each. In all the cases, the bound compound detaches from the minor groove and explores different conformations. Figure S5 shows a histogram of the distance between the center of mass of the GAAC sequence at the starting point of the ligand and the ligand itself, for each copy. The lower distance values would correspond to the Cu compound remaining in the starting position, which is not the case in either of the five independent copies; instead, the electron density accumulation on the minor-groove floor promotes the fact that the shift of the secondary ligand detaches from the minor groove and explores the DNA until a final conformation is reached. The conformation for both Cas IV-gly and Cas VI-gly corresponds to the same type of interaction as observed, with the amino part of the glycinate molecule forming hydrogen-bond interactions with the cytosine O2 oxygen atom. It was not necessary to perform this experiment with Cas III-Ea and Cas III-La because those compounds are symmetrical.

Copper Compounds Exhibit Sequence Selectivity

From the grid analysis presented, there is some hint that the Cu compounds present some degree of DNA sequence selectivity. The Cu compounds interacted more frequently with A:T dense regions of the sequence in comparison with the G:C base-pair regions. In order to verify these observations, MD of a separate control sequence, d(GCATAAACAGGTCTGCGC), was performed following the same methods. Modified from the mini-ABC sequence library, this sequence was chosen because it contains features commonly found in genomic DNA including the motifs: TG/CA, GG/CC, AT/TA, and AAA.[39] Once the trajectories converged, a grid analysis (Figure S6) similar to the one performed with the GAAC sequence was performed. These results were consistent with the previous observations. Compounds preferentially interacted with the A:T base pairs while G:C base pairs had much less-frequent interactions, indicating a certain level of sequence selectivity. These results are in agreement with the previously described sequence-dependent, but not nucleotide-specific, DNA cleavage activity for Cu-phen derivatives.[42]

Conclusions

The use of MD simulations allowed us to study the binding modes of four copper compounds with two 18-mer double-stranded DNAs. To avoid possible bias in the starting structure used to perform the simulations, four copies of the copper compound were positioned 10 Å away from the center of the duplexes without any restraints. The results show that these compounds bind to DNA through terminal base-pair stacking, minor-groove binding, or intercalation/base-pair eversion mechanism. It is clear that the secondary ligand of the Cu compounds drives the binding mode with DNA. The glycinate ligand for Cas II-gly and Cas VI-gly provides a unique directionality for the binding to DNA as observed in the simulations. This is caused by specific hydrogen-bond interactions between the amino group of the glycinate with the pyrimidine O2 atom and to a lower extent the N3 nitrogen atom of adenine. Hydrogen bonds that help in this directionality are generated between the oxygen atoms of the carboxylate group and hydrogen in the sugar backbone of the DNA. This direction preference was confirmed through simulations where the carboxylate part of the glycinate was switched inside the minor groove and resulted in the Cu compound unbinding from the DNA. The Cu compounds Cas III-Ea and Cas III-La which possess an acetylacetonate secondary ligand showed increased events of base-pair eversion that are both long-lived and consistent within the different copies we calculated. We suggest that the DNA minor-groove binding events initially observed could reach a more stable conformation through extended sampling time, which could lead to higher base-pair eversion events. We have not detected any influence of the positions of the methyl groups present in the aromatic ligand of the Cu compounds. Analysis focusing on studying any relevant interaction between these methyl groups and the DNA failed to show distances that could potentially form any type of binding.