New Insights into the Structures of LigandâQuadruplex Complexes

J. Phys. Chem. B 2010, 114, 15301–15310

15301

New Insights into the Structures of Ligand-Quadruplex Complexes from Molecular Dynamics Simulations Jin-Qiang Hou, Shuo-Bin Chen, Jia-Heng Tan, Tian-Miao Ou, Hai-Bin Luo, Ding Li, Jun Xu, Lian-Quan Gu,* and Zhi-Shu Huang* School of Pharmaceutical Sciences, Sun Yat-sen UniVersity, Guangzhou 510006, China ReceiVed: July 19, 2010; ReVised Manuscript ReceiVed: October 16, 2010

G-quadruplexes are higher-order DNA and RNA structures formed from guanine-rich sequences, and they are attractive anticancer drug targets. Understanding the three-dimensional interactions between a G-quadruplex and its ligand in solution is the key to discovering a drug lead. Hence, from crystallographic or NMR structures, molecular dynamics studies have been performed on six ligand-quadruplex complexes. BRACO-19, BSU6039, daunomycin, RHPS4, MMQ1, and TMPyP4 are the six ligands that bind to the G-quadruplex structures in the studies. Based on molecular dynamics simulations and a series of computational analyses, the results suggest that ions move away from the external G-quartet to let the ligand bind to the quadruplex in aqueous solution. The ligand binding can increase the stability of the Hoogseen hydrogen bonds within the G-quartet. However, the G-quartet binding site can only fit one ligand molecule. The ligand can form hydrogen bonds at the loop or flank of the quadruplex. However, not all the interactions will stabilize the ligand-quadruplex complex in aqueous solution. These findings can assist in the design of selective and potent G-quadruplex ligands. TABLE 1: X-ray or NMR Structures Used in the MD Studiesa

Introduction Under appropriate cation conditions, guanine-rich singlestranded DNA or RNA can form four-stranded structures called G-quadruplexes.1-3 The G-quadruplex consists of π-stacked G-quartets, each of which is connected by a network of eight Hoogsteen hydrogen bonds. Possible G-quadruplex sequences are located in a number of human chromosomal regions, such as telomeres and oncogene promoters.4-8 G-quadruplexes play important roles in many biological events. The quadruplexes in telomeric DNA result in the inhibition of telomere elongation in cancer cells.9-11 The quadruplexes in many oncogene promoters including c-myc,12-15 c-kit,16,17 and bcl-218,19 can suppress gene expression and alter subsequent biological processes. Therefore, the induction of G-quadruplex formation, or stabilization of the quadruplex by a small molecule, opens an opportunity for cancer therapeutics.20-23 To discover a novel and effective G-quadruplex ligand, it is important to understand the interaction between a small molecule and its G-quadruplex receptor. Data from X-ray and NMR show that the ligands bind to the terminal G-quartets at the end of the G-quadruplex.24-26 The structural data from X-ray or NMR provide detailed structural information. However, dynamic conformational rearrangements and transition states may be more important.27,28 Molecular dynamics (MD) simulation may provide information regarding the rearrangements and transition states. A number of MD studies on G-quadruplexes,29-34 and ligand-quadruplex interactions,35-39 have been reported. There are still questions left unanswered. For examples, (1) why can three daunomycin molecules bind to the 5′- end of a G-quartet in a crystal structure, (2) is this complex stable in solution, (3) does the G-quartet surface accommodate all three ligands in * To whom correspondence should be addressed. Phone: 8620-39943056 (Z.-S.H.); 8620-39943055 (L.-Q.G.). Fax: 8620-39943056 (Z.-S.H. and L.-Q.G.). E-mail: [email protected] (Z.-S.H.); cesglq@ mail.sysu.edu.cn (L.-Q.G.).

PDB ID

sequence

source

ligand

experiment

ref

1l1h 1o0k 3ce5 1nzm 2jwq 2a5r

d(G4T4G4) d(TG4T) d(TAG3T2AG3T) d(T2AG3T) d (T2AG3T) Pu24Ib

Oxy-tel Oxy-tel h-tel h-tel h-tel c-myc

BSU6039 daunomycin BRACO-19 RHPS4 MMQ1 TMPyP4

X-ray(1.75 Å) X-ray(1.17 Å) X-ray(2.5 Å) NMR NMR NMR

40 41 42 43 44 45

a Three structures are from the human telomere (h-tel), two are from the Oxytricha nova telomere (Oxy-tel), and one is from the c-myc promoter region. b Pu24I: TGAGGGTGGIGAGGGTGGGGAAGG. I: 15N,13C-labeled guanine.

solution, (4) how do the ligands stabilize the G-quadruplex, and (5) what are the main driving forces for the ligand-quadruplex interactions? Previous studies have emphasized hydrogen binding related interactions; however, other interactions that stabilize the ligand-quadruplex complex in solution should also be investigated. So far, MD studies in explicit solvent have been carried out on three X-ray or three NMR ligand-quadruplex complexes, where ligands bind to the terminal G-quartets at the end of the quadruplex (Table 1). The DNA structures with different ions and ligands provide various structural data on G-quadruplexes. The data are fundamental for ligand design. The work aims to answer, through MD simulations, the above-mentioned questions. Six ligands, BRACO-19, BSU6039, daunomycin, RHPS4, MMQ1, and TMPyP4, are studied as shown in Figure 1. Method Data. More than 100 G-quadruplexes have been reported in the RCSB Protein Data Bank. However, many of them have no ligands, and some of them are almost the same, differing only in some metal ions. The data for this work are listed in Table 1.

10.1021/jp106683n  2010 American Chemical Society Published on Web 11/04/2010

15302

J. Phys. Chem. B, Vol. 114, No. 46, 2010

Hou et al.

Figure 1. Structures of ligands in a complex with G-quadruplex for MD studies.

TABLE 2: Models for MD Simulations environment model +

1l1h(K ) 3ce5(K+) 1o0k(Na+) 1o0k(K+) 1nzm(K+) 2jwq(K+) 2a5r(K+) a

structure bimolecular antiparallel bimolecular parallel tetramolecular parallel tetramolecular parallel tetramolecular parallel tetramolecular parallel intramolecular parallel

ligand

ligand net charge

BSU6039 BRACO-19 daunomycin daunomycin RHPS4 MMQ1 TMPyP4

+3 +3 +1 +1 +1 0 +4

ionsa

watersb

time scale (ns)

+

13 062 14 787 13 668 13 575 14 643 12 981 15 620

10 10 32 26 10 10 10

K (19) K+(19) Na+(12) K+(12) K+(22) K+(24) K+(20)

Type and numbers of ions in the models. b Numbers of water atoms.

MD Simulations. The ligands’ partial charges were computed using the HF/6-31G* basis set from GAUSSIAN 03,46 and refined by RESP calculation in the antechamber module of the Amber 10 program.47 Other parameters were adopted from the Generalized Amber Force field (GAFF).48 Consecutive K+/Na+ ions are vertically aligned along the axis in the central core of the complex. The inner ions are kept in the positions as observed in the crystal structure. The coordinates of inner ions are not available in NMR structures; K+ is placed manually in the central core of the complexes. Every ligandreceptor system was then solvated in a truncated octahedron box of TIP3P water molecules with a 10.0 Å buffer along each dimension. Additional positively charged counterions were added in the system to neutralize the negative charge on the DNA backbone. The counterions were automatically placed, using the Coulombic potential, into the most negative locations by the LEAP program. These settings are consistent with the experimental conditions. The structure of the daunomycinquadruplex complex (PDB ID: 1O0K) was measured in a Na+ environment. Since the K+ ion is more biologically relevant, another model, 1O0K(K+), was built from the crystal structure by replacing the inner Na+ with K+. All models are listed in Table 2. The AMBER ff03 force field49 was applied for G-quadruplexes, ions, and water molecules. Simulations were performed with the sander program of Amber 10. Periodic boundary conditions were applied to avoid the edge effect. The particle mesh Ewald (PME) method was employed to calculate long-range electrostatic interactions.50 The hydrogen bonds were

constrained using SHAKE algorithm.51 For the nonbonded interactions, a residue-based cutoff of 10 Å was used. Temperature regulation was achieved by Langevin coupling with a collision frequency of 1.0. The final systems were subjected to initial minimization to equilibrate the solvent and countercations. The G-quadruplex complex and inner ions were initially fixed with force constants of 100 kcal · mol-1 and were then gradually reduced in subsequent runs. About 2000 steps were used for each stage. The first 1000 steps were carried out using the steepest descent algorithm, and the remaining steps were performed with a conjugate gradient algorithm. The system was then slowly heated from 0 to 300 K within 100 ps. After the heating process, a further 100 ps of equilibration at 300 K was carried out to obtain a stable density. A force constant of 50 kcal mol-1 was maintained for the complex and inner ions during the heating and equilibration. Afterward, an unconstrained production phase was initiated and continued for 10-30 ns in an NPT ensemble at 1 atm and 300 K. The time step used for the MD simulations was set to 2.0 fs and the trajectory files were collected every 1 ps for the subsequent analysis. All trajectory analysis was done with the Ptraj module in the Amber 10.0 software and examined visually using the VMD software.52 Binding Free Energy Analysis. The binding free energies (∆Gbind) were calculated using the MM-PBSA approach53,54 inside the AMBER program. The dielectric constants were set to 1 for the solute, and 80 for the surrounding solvent molecules. The K+ radius was set to 2.025 Å.55 A total of 100 snapshots

Structures of Ligand-Quadruplex Complexes

J. Phys. Chem. B, Vol. 114, No. 46, 2010 15303

Figure 2. Schematic representation of G-quadruplex conformations and ligands (in orange) binding sites.

were taken from the last 2 ns trajectory with an interval of 20 ps. Prior to the analysis, all counterions and water molecules (except for the K+/Na+ ions present in the central channel of the G-quadruplex) were stripped from the trajectories. For each snapshot, the free energy was calculated for each molecular species (complex, ligand, and quadruplex) using the following equation:36,56

∆Gbind ) Gcomplex - [Gligand + Gquadruplex]

(1)

where Gcomplex, Gligand, and Gquadruplex were the free energies for the complex, quadruplex, and ligand, respectively. Each of them was estimated with the following equation:

∆Gbind ) ∆Εmm + ∆Gsolv - T∆S

(2)

where ∆Emm was the change in the molecular mechanics free energy upon complex formation in the gas phase, ∆Gsolv was the solvation free energy, and -T∆S was the entropy term. The molecular mechanics free energy was calculated as the following:

∆Emm ) ∆Eelec + ∆Evdw + ∆Eini

(3)

where ∆Eelec was the Coulomb interaction, ∆Evdw was the van der Waals interaction, and ∆Eini was the sum of the bond, angle, and dihedral energies; in this case, ∆Eini ) 0. The solvation free energy was calculated in (4):

∆Gsolv ) ∆GPB + ∆Gnp

(4)

where ∆GPB was the polar contribution to solvation, which was obtained by solving the Poisson-Boltzmann equation57 for the MM-PBSA method. ∆Gnp was the nonpolar solvation term, which was computed in (5):

∆Gnp ) γ∆SASA + b

(5)

where γ was the surface tension that was set to 0.0072 kcal/ (mol Å2) and b was a constant that was set to 0.58 SASA is the solvent accessible surface area (Å2) that was estimated using the MOLSURF algorithm. The solvent probe radius was set to 1.4 Å to define the dielectric boundary around the molecular surface. The calculation for the entropic contribution was performed with the NMODE module, which computes vibrational, rotational, and translational entropies. A total of 100 snapshots were collected from the last 2 ns of the trajectory at 20 ps intervals. Each snapshot was minimized using the conjugate gradient method for 10 000 steps with a distance-dependent dielectric constant of 4 interatomic distances. The terminating criterion

for minimization was set to 0.01 kJ/mol Å gradient rms. Graphic models were prepared using the Pymol version 0.9,59 Discovery Studio (Accelrys, San Diego, CA, U.S.), or the Xmgrace program (http://plasma-gate.weizmann.ac.il/Grace/). Results The structures of ligand-quadruplexes were measured under different conditions as shown in Figure 2. The structure of the BRACO-19-quadruplex complex (PDB ID: 3ce5) indicates that the ligand stacked at the end of the G-quartet and formed hydrogen bonding interactions with quadruplex and water molecules.42 The structure of the BSU6039-quadruplex complex (PDB ID: 1l1h) indicates that the acridine moiety of the BSU6039 molecule binds to the end of G-quartet and forms hydrogen bonds with thymine bases in the loop.40 The structure of the daunomycin-quadruplex complex (PDB ID: 1o0k) demonstrates that three daunomycin ligands jointly bind to the end of the G-quartet. The three daunosamine substituents occupy three of the four quadruplex grooves.41 The NMR structures of the MMQ1-quadruplex (PDB ID: 2jwq) and RHPS4-quadruplex (PDB ID: 1nzm) show that two molecules of the ligand separately bind at the 5′ and 3′ ends of the tetramolecular parallel G-quadruplex.43,44 Another NMR structure of TMPyP4quadruplex complex (PDB ID: 2a5r) shows that a TMPyP4 molecule binds to the 5′ end of the parallel G-quadruplex.45 These binding models consist of G-quadruplexes, ions, and ligands. To understand the dynamic properties of ligandquadruplexes, their stabilities are analyzed as follows. 1. Structural Stability of the G-Quadruplex. 1.1. Stability of G-Quartets. The root mean-square deviation (RMSD) values were used to measure the conformational stability of a structure. To measure the stabilities of G-quartets during the MD simulations, RMSD values for heavy atoms in G-quartets were monitored for all MD trajectories and compared against those for the corresponding crystal or NMR structures as shown in Figure 3. By inspecting RMSD curves in Figure 3, one can conclude that the G-quartets are stable over the course of the simulation. The average RMSD values of G-quartet heavy atoms are between 0.47 and 0.66 Å (Table 3). In a G-quadruplex structure, each G-quartet is held together by eight Hoogsteen hydrogen bonds (N2-H · · · N7 and N1H · · · O6), which are the main forces for keeping quartets stable. To study the impact of the ligand binding to the stability of the G-quartet, the hydrogen bond occupancies of the G-quartet were calculated along all MD trajectories (Figure 4), and Table 4 shows the average occupancy of eight hydrogen bonds for each G-quartet. The results indicate that the outer G-quartets (Gquartet 1, G-quartet 3, or G-quartet 4) are also stable. The average hydrogen bond occupancies are more than 97.55%. The inner G-quartets (G-quartet 2 or G-quartet 3) are less stable (the average hydrogen bond occupancies range from 85.55 to 99.10%). The N2-H · · · N7 bonds are more stable than the N1-H · · · O6 bonds in all G-quartet models (Figure 4). This may

15304


Hou et al.

Figure 3. Time dependence of the RMSD of G-quartet heavy atoms (black) and all ligand atoms (lig1, lig2, and lig3 are shown in red, green, and blue, respectively).

TABLE 3: Average RMSD Values (Å) of the Heavy Atoms in the G-Quartet and the Average Mass-Weighted RMSD Values of the Ligands of All Atoms Calculated from the MD Simulationsa model +

3ce5(K ) 1l1h(K+) 1o0k(Na+) 1nzm(K+) 2jwq(K+) 2a5r(K+)

G-quartet

Lig1

0.48 0.56 0.64 0.47 0.56 0.66

1.16 (0.14) 1.20 (0.11) 1.01 (0.28) 0.33 1.23 (0.16) 0.70

Lig2

Lig3

0.79 (0.23) 0.34 1.47 (0.16)

1.11 (0.32)

a Numbers in parentheses are the average RMSD values for the cyclic carbons in ligands.

be due to the strong attractive force between the inner ions and the O6 atoms of the guanine base. N1-H · · · O6 hydrogen bonds are weakened because O6 atoms are moved closer to the inner ions. An ab initio study60 on a G-quadruplex without cations indicated that the G-quartet was stabilized by bifurcated N1-H · · · O6 and N2-H · · · O6 hydrogen bonds. Chowdhury and colleagues reported their MD simulations29 on a 7-mer Gquadruplex and suggested that normal Hoogsteen type hydrogen bonds in the G-quartet were restored in the presence of Na+ ions at the ligand binding sites. This study reveals that both N1-H · · · N7 and N2-H · · · O6 bonds (Figure 5) are formed with lower occupancies in the presence of K+ ions. The average occupancy of eight hydrogen bonds in the ligandbinding G-quartet (G-quartet 1) is greater in comparison with that of a ligand-free G-quartet (G-quartet 3 or G-quartet 4). This suggests that the ligand-bounded G-quartet is more stable than the ligand-free G-quartet. Therefore, ligand binding has a positive impact on G-quartet stabilization, and helps the overall G-quadruplex stabilization. 1.2. OWerall Stability of the G-Quadruplexes. To investigate the overall stability of a G-quadruplex, RMSD values of backbone atoms (including P, O3′, O5′, C3′, C4′, and C5′ atoms) were monitored in all MD trajectories by comparing against those for the crystal or NMR structures (Supporting Information, Figure S7). A jump in RMSD was observed at the 2 ns equilibration for all complexes. The jump was due to the

backbone relaxation of the starting structure. The complex structure was stabilized after 5 ns of equilibration, where the RMSD value converged to a steady value. The backbones of G-quadruplexes did not significantly move away from starting structures in model 3ce5(K+) (containing two propeller loops and two flanking nucleotides), 1l1h(K+) (containing two diagonal loops), and 2a5r(K+) (containing three propeller loops and one diagonal loop). More backbone changes took place in the model 1nzm(K+) and in 2jwq(K+) (both contain four flanking sequences). This indicates that flanking nucleotides are flexible in a G-quadruplex. The averaged structures were obtained by averaging 50 MD simulation snapshots, which were taken from the last 2 ns on the MD trajectory with an interval of 20 ps. Each averaged structure was superimposed with the starting structure, as shown in Figure 6. Figure 6 demonstrates that the averaged structures and experimental structures fit well. The conformations at loop domains have greater differences due to the limitations of force fields.32,55 Overall, the trajectory structures have little fluctuation in comparison with those of the experimental structures. The nucleotides at the loop and flank of a G-quadruplex become rigid when ligands bind. In other words, ligand binding induces a conformational change of the loop and flanking nucleotides. 2. Ions Behavior. Metal ions play an important role in stabilizing the G-quadruplex structure and induce the formation of topologies of G-quadruplex. Visual inspection of the trajectories reveals that metal ions between the two G-quartets were stable during the entire simulation. The structure of the daunomycin-quadruplex complex (PDB ID: 1o0K) was determined in a Na+ ionic environment, others were in a K+ ionic environment. Since K+ is more biological related, a model 1o0k(K+) was built by replacing Na+ with K+. The trajectories show that the G-quartets are stable within the first 500 ps of simulation (Figure 7). At 500 ps after starting the simulation, the K+ ion between the G-quartet3 and the G-quartet4 moved from the ion channel to the entrance, and was close to the negatively charged phosphate. Consequently, the G-quartets became unstable (Supporting Information, Figure S8). Previous MD studies29,61



Figure 4. Hydrogen bond occupancies within every G-quartet during MD simulations. Occupancy is provided as a percentage of the investigated time period.

TABLE 4: Average Occupancy (%) of Eight Hydrogen Bonds in Each G-Quartet Calculated along Entire MD Trajectories G-quartet1 G-quartet2 G-quartet3 G-quartet4

3ce5(K+)

1l1h(K+)

1o0k(Na+)

1nzm(K+)

2jwq(K+)

2a5r(K+)

99.59 98.39 99.41

99.89 85.55 89.75 99.67

99.50 91.06 88.41 97.98

99.74 98.16 99.24

99.82 99.10 99.88

99.75 96.83 97.55

suggested that Na+ ions could enter the empty central channel of the G-quartets. To find out whether K+ could return to the central channel in our model, this simulation time was increased to 26 ns, and K+ did not return. This result suggested that metal ions are important for the stabilization of ligand-quadruplex complexes. The result is consistent with previous experimental studies. For model 1l1h(K+), an ion below the G-quartet layers moved rapidly from the initial position, and entered the solvent in one ns. A similar phenomenon was also observed by Spackova and

Figure 5. Hydrogen bonding interactions in a representative G-quartet, indicating normal Hoogseen hydrogen bonds (N2-H · · · N7 and N1-H · · · O6) and other hydrogen bonds (N2-H · · · O6) and N1-H · · · N7).

co-workers.30 The results suggest that ions move away from the external G-quartet to let the ligand bind to the quadruplex in aqueous solution. 3. Interactions of Ligands with G-Quadruplexes. 3.1. Dynamics BehaWior of Ligand Binding to G-Quadruplex. To investigate the ligand conformational changes, mass-weighted RMSDs (mwRMSD) between simulated ligand conformations and the initial conformation (crystal or NMR structure) were calculated and depicted in Figure 3. Little mwRMSD fluctuations were observed in models 1nzm(K+) and 2a5r(K+) (Figure 3D,F). Greater mwRMSD fluctuations were observed in models 3ce5(K+), 1l1h(K+), and 2jwq(K+) (Figure 3A,B,E). The side chains at the ligand cyclic scaffold were responsible for the fluctuations. In Table 3, the average mwRMSD values for the ligand with flexible side chains were between 1.16 and 1.47 Å (BSU6039, 1.20 Å; BRACO-19, 1.16 Å; MMQ1, 1.23, 1.47 Å). The ligand without side chains had smaller mwRMSD fluctuations (TMPyP4, 0.70 Å; RHPS4, 0.33, 0.34 Å). When the side chains were discounted, the average RMSD values were between 0.11 and 0.16 Å (BSU6039, 0.11 Å; BRACO-19, 0.14 Å; MMQ1, 0.16, 0.16 Å), which were comparable with values for one of the non-side-chain groups. These results indicated that the side chains were more flexible than the cyclic scaffolds. These results are thus consistent with previous reports of 2 ns MD studies on BSU6039-quadruplex crystal structure.40 Daunomycin (lig2) with the daunosamine sugar moiety was stable in the Model 1o0k(Na+) simulation (Figure 3C). The average mwRMSD was 0.79 Å (Table 3), which was compa-

15306


Hou et al.

Figure 6. Averaged structure (yellow) from the last 2 ns trajectories superimposed with the X-ray or NMR structures (green). Purple and red spheres are metal ions in the starting structures and in the averaged structures. Red arrows indicate ligand binding sites.

Figure 7. Snapshots from the MD simulation for Model 1o0k(K+). The G-quadruplex is in green, and the ligand is in yellow. Red arrows indicate the direction of cation movement.

rable to the RMSD value of cyclic scaffold atoms of 0.23 Å. This is because daunomycin’s chiral carbon C13 is in the S configration, which makes it impossible for the daunosamine sugar moiety to move from one side of the aromatic plane to another. This structural feature decreases the mobility of the ligands. The mwRMSD values of the other two daunomycin molecules, lig3 and lig1, demonstrated peak fluctuations at the eighth and the thirteenth nanosecond (Figure 3C blue and red). The translocation of these two daunomycin molecules from the G-quartet surface to the nearby groove region was observed. Daunomycin lig2 migrated to occupy more G-quartet surface. This suggested that not all of the three daunomycin molecules stably docked to the G-quartet in the aqueous solution. To better understand fluctuations in the conformation and orientation of ligands, ten conformation snapshots were taken from the last 2 ns trajectories, with an interval of 200 ps. Five structures were randomly chosen and superposed with the X-ray or NMR structures (Figure 8). For model 1o0k(Na+) (Figure 8C), the pose of lig1 was significantly changed and was not binding to the G-quartet surface. Due to the steric clashes, Lig3 was pushed away by lig2, and lost binding to the G-quartet surface. Therefore, lig2 gained more room to effectively bind to the G-quartet. This confirms that not all three daunomycin

molecules bind to the G-quartet in the same way, and that the G-quartet surface allows only one ligand to effectively bind. For the other models, the conformation and orientation of the ligands were similar to those from crystal or NMR structures. 3.2. Energetic Analyses of the LigandsQuaduplex Complexes. The binding free energies were computed by means of the MM-PBSA. The results are listed in Table 5. The binding free energies for all the models ranged from -6.67 to -49.49 kcal · mol-1. The calculated energies could be overestimated, i.e., the model 2a5r(K+) (∆Gbind ) -49.49 kcal · mol-1). This was probably due to the limitation of the MM-PBSA for highly charged species (TMPyP4 possesses four positive charges). Nevertheless, the individual energy components allowed a better understanding of ligand-quadruplex complex formation and their driving forces. For all complexes, the electrostatic energy (∆Eelec) and the van der Waals energy (∆Evdw) favorably contributed to the binding free energies. The nonpolar solvation energy (∆Enp) also made a favorable contribution. This implied that electrostatic energy, van der Waals forces, and nonpolar solvation energy were the main driving forces for the formation of ligand-quadruplex complexes. However, the solvent electrostatic energy (∆EPB) and entropic effect (T∆S) made unfavor-



Figure 8. Dynamics snapshots of the ligands taken from the last 2 ns trajectories (blue), which were superposed with X-ray or NMR structures (orange).

TABLE 5: Computed Binding Free Energies for Ligand-Quadruplex Complexes (kcal/mol) 1o0k(Na+) +

∆Eelec ∆Evdw ∆Eini ∆Gmm ∆Gnp ∆GPB ∆Gsolv ∆Gelec,tot ∆Gtot -T∆S ∆Gbind

+

3ce5(K )

1l1h(K )

Lig1

Lig2

Lig3

-1379.62 (25.47) -60.87 (2.88) 0.00 (0.00) -1440.49 (26.63) -7.43 (0.29) 1399.89 (27.63) 1392.46 (27.45) 20.28 (4.55) -48.02 (4.26) 33.52 (6.33) -14.5

-1327.41 (22.79) -61.63 (2.76) 0.00 (0.00) -1389.04 (23.69) -6.38 (0.24) 1331.37 (21.43) 1324.99 (21.29) 3.96 (5.27) -64.05 (5.37) 26.23 (6.55) -37.83

-340.23 (11.76) -30.87 (2.19) 0.00 (0.00) -371.10 (12.05) -3.81 (0.16) 347.25 (12.07) 343.43 (12.01) 7.02 (2.91) -27.66 (2.88) 18.63 (3.41) -9.03

-402.34 (20.55) -44.84 (2.86) 0.00 (0.00) -447.18 (19.83) -5.12 (0.17) 419.49 (18.17) 414.37 (18.12) 17.15 (5.48) -32.81 (4.66) 18.13 (4.59) -14.68

-337.67 (7.88) -35.90 (2.70) 0.00 (0.00) -373.57 (8.28) -4.18 (0.14) 348.43 (8.21) 344.24 (8.16) 10.76 (2.72) -29.32 (2.85) 19.53 (3.28) -9.79

1nzm(K+) ∆Eelec ∆Evdw ∆Eini ∆Emm ∆Gnp ∆GPB ∆Gsolv ∆Gelec,tot ∆Gtot -T∆S ∆Gbind

2jwq(K+)

Lig1

Lig2

Lig1

Lig2

2a5r(K+)

-416.20 (7.97) -43.98 (4.12) 0.00 (0.00) -460.18 (9.57) -5.05 (0.37) 432.27 (8.04) 427.22 (7.93) 16.07 (2.62) -32.96 (4.13) 17.73 (5.21) -15.23

-509.17 (5.85) -59.12 (1.71) 0.00 (0.00) -568.29 (6.30) -5.75 (0.13) 526.25 (6.18) 520.49 (6.17) 17.07 (2.51) -47.80 (2.86) 14.68 (5.47) -33.12

-11.22 (3.26) -63.32 (4.89) 0.00 (0.00) -74.54 (5.06) -7.00 (0.30) 38.65 (3.21) 31.65 (3.13) 27.43 (3.45) -42.89 (4.28) 21.29 (4.69) -21.6

-6.02 (3.14) -66.05 (2.89) 0.00 (0.00) -72.07 (4.61) -6.04 (0.22) 28.51 (4.14) 22.47 (4.22) 22.50 (4.98) -49.60 (4.30) 15.58 (5.09) -34.02

-2078.22 (40.77) -92.02 (3.29) 0.00 (0.00) -2170.24 (42.11) -8.79 (0.25) 2101.38 (41.34) 2092.59 (41.19) 23.16 (3.09) -77.65 (3.51) 28.16 (6.48) -49.49

∆Eelec is the electrostatic interaction calculated by the MM force field. ∆Evdw is the van der Waals contribution from MM. ∆Eini is the internal energy. ∆Emm is the total molecular-mechanical energy (∆Eele + ∆Evdw + ∆Eini). ∆Gnp is the nonpolar contribution to the solvation energy. ∆GPB is the electrostatic contribution to the solvation energy calculated by the PB approach. ∆Gsolv is the total solvation energy (∆Gnp+∆GPB). ∆Gelec,tot is the total electrostatic energy (∆Gelec + ∆GPB). ∆Gtot is the total energy without solute entropic contribution (∆Emm + ∆Gsolv). -T∆S is solute entropic contribution, where T ) temperature and S ) sum of translational, rotational, and vibrational entropies. ∆Gbind is the total energy with solute entropic contribution (∆Gtot - T∆S).

able contributions to all complexes. More positively charged ligands result in more negative electrostatic energy (∆Eelec). For models 2jwq(K+) and 1nzm(K+), ligands MMQ1 and RHPS4 were both in complexes with the intermolecular parallelstranded DNA quadruplex d(T2AG3T)4, which contains the human telomeric repeat. A direct comparison of the energy components of these two complexes enabled a better understanding of the structure-activity relationship. Two van der Waals energies in the MMQ1-quadruplex complex (-63.32, -66.05 kcal · mol-1) are more negative than those of the RHPS4-quadruplex complex (-43.98, -59.12 kcal · mol-1), although both ligands possessed five aromatic rings. This was

because the crescent-like scaffold in MMQ1 molecule stacking more effectively with G-quartet. For the model 2jwq(K+), the free energy of MMQ1 binding to the lig2 site was -34.02 kcal · mol-1. This means that the ligand more favorably binds to the lig2 site rather than to the lig1 site (∆Gbind ) -21.6 kcal · mol-1). For model 1nzm(K+), the free energy of RHPS4 binding to the lig2 site and the lig1 site are -33.12 and -15.23 kcal · mol-1, respectively. This indicates that the ligands prefer to bind to the lig2 site. This is because the lig2 pocket is deeper inside the G-quadruplexes; the lig1 site is more exposed to the solvent. The result suggests that the two ends of the quadruplex have different environments,

15308


Hou et al.

TABLE 6: Hydrogen Binding Data for the Four Modelsa model +

3ce5(K )

1l1h(K+)

1o0k(Na+) 2jwq(K+)

donor

acceptor

occupancy (%)

distance (å)b

angle (deg)

T24@O4 WAT WAT T12@O2 T6@O2 lig1@O14 T7@O1P T8@O1P G17@O1P lig2@O10 A3@N3

lig1@H31 lig1@N21 lig1@H90 lig1@N7 T24@H3 T24@N3 lig1@H17 lig1@N19 lig1@H70 lig1@N T6@H3 T6@N3 lig1@H69 29@N27 lig1@H13 29@N1 lig2@H23 lig2@N3′ G2@H22 G2@N2 lig1@H16 lig1@N3

90.40 99.70 99.15 32.60 99.94 99.74 87.72 31.16 65.64 57.57 52.40

2.92 (0.13) 2.94 (0.12) 3.01(0.14) 3.08 (0.20) 2.84 (0.11) 2.94 (0.14) 2.84 (0.13) 2.95 (0.19) 2.95 (0.18) 3.08 (0.18) 3.18 (0.16)

24 (14) 16 (9) 18 (9) 40 (13) 18 (8) 19 (9) 23 (13) 26 (13) 28 (14) 39 (12) 20 (11)

a The listed donor and acceptor pairs satisfy the criteria for the hydrogen bond over 30.0% of the time during the simulations. b The average distance between the hydrogen acceptor atom and the hydrogen donor atom in the investigated time period.

and the same ligand has different binding affinities when it is docked at different ends. In model 1o0k(Na+), daunomycin is stably bound to the G-quartet at the lig2 position for all simulations. However, daunomycin unsuccessfully bound to the G-quartet at the lig1 or lig3 position. The data from Table 5 support the observation (daunomycin binding free energy at lig2 position is -14.68 kcal · mol-1, and at lig1 and lig3 position are -9.03 and -9.79 kcal · mol-1, respectively). This concludes that the G-quartet surface could only accommodate one ligand binding in aqueous solution. 3.3. Hydrogen Bonding. An interaction was considered to be a hydrogen bond if the distance between the hydrogen donor and acceptor was less than 3.5 Å, and the binding angle was greater than 120.0°. The hydrogen binding data from the four models have been listed in Table 6. For model 3ce5(K+), BRACO-19 bound at the 3′ end of G-quartet and formed hydrogen bonds with the quadruplex and water molecules. Six hydrogen bonds were formed between the ligand and the receptor in the crystal structure; five were watermediated hydrogen bonds, as shown in Figure 9A. As shown in Table 6, the hydrogen bond between N21 of BRACO-19 and O4 of T24 is strong and stable in the aqueous solution with a mean distance of 2.9 Å (3.0 Å in the crystal structure), and its occupancy is 90.4%. N7 of BRACO-19 and N3 of T24 are mediated by a water molecule, and the binding is stable during the simulation. The mean distance of N7 · · · water and N3 · · · water are 2.9 and 3.0 Å, respectively (2.7 and 2.9 Å in the crystal structure, respectively). Other water molecules, however, could not maintain their starting positions during the simulations (Figure 9A). A new hydrogen bond was formed between N19 of BRACO-19 and O2 of T12 with a distance of 3.1 Å in the simulation, although the binding score was low. Figure 9B depicts the hydrogen bond distance versus the simulation time in the MD trajectory, and indicates that most hydrogen bonds are still stable. For model 1l1h(K+), BSU6039 bound at one end of the G-quartet kept the pose in the loop. The average hydrogen bond distance between the N of BSU6039 and the O2 of T6 was 2.8 Å (3.0 Å in crystal structure) (Table 6), and the hydrogen bond was stable during MD simulation with an occupancy of 99.9%. The mean hydrogen bond distance between O14 of BSU6039 and N3 of T6 was 2.9 Å (2.7 Å in crystal structure) with an occupancy of 99.7%. Two new hydrogen bonds were also observed in the simulation with lower occupancies. For model 1o0k(Na+), the average hydrogen bond distance between N3′ of daunomycin and phosphate oxygen of G17 was 3.0 Å (2.7 Å in crystal structure) with an occupancy of 65.6%.

A newly formed hydrogen bond between O10 of daunomycin and N2 of G2 was observed with an occupancy of 57.6%. For model 2jwq(K+), one hydrogen bond was formed between the MMQ1 molecule and G3, with an occupancy of 52.4%. There was no hydrogen bond observed in models 1nzm(K+) and 2a5r(K+). The hydrogen bonds increase the binding affinity of ligands and G-quadruplexes, and enhance the selectivity. The thymine base seems easier to form hydrogen bonds in comparison with others. Discussion and Conclusions Six ligand-quadruplex complex models have been simulated by molecular dynamics method for a total of 108 ns. The dynamic interactions between ligands and G-quadruplexes are understood. Ionic changes will destabilize the G-quadruplexes. It is hard for ions to enter the ion channel of the G-quadruplexes from outside. This observation is consistent with the previous MD study. Metal ion binding can change ligand binding status with the G-quartet. Therefore, most of NMR structures of complexes (in solution) are bound with more than one ligand, and most of the X-ray structures of complexes are bound with one ligand at the end of the G-quartet. The ligand prefers to bind at the surface of the G-quartet. Usually, only one ligand is allowed in solution. Once a ligand is stably bound, it increases the G-quartet stability of Hoogseen hydrogen bonds; thus, the overall G-quadruplex is stabilized. With the increasing number of quadruplexes being identified in the genome, the success of gene targeting therapy will depend not only on the small molecule selectivity for the quadruplex against duplex DNAs but also crucially on the selectivity for the targeted quadruplex over other quadruplexes. So far, the binding energy analysis indicates that electrostatic energy dominates ligand-quadruplex interactions. However, this cannot improve ligand selectivity for the targeted G-quadruplex. This is because nucleic acids possess negative net charges. The positively charged ligands will bind to them without selectivity. The π-π stacking produces stronger van der Waals energy, which cannot improve the selectivity either. Hydrogen bond analysis indicates that hydrogen bonds between ligands, loop, and flanking nucleotides (especially for adenine base) of G-quadruplexes are stable in solution. Binding site rearrangement may introduce a different G-quartet surface environment, which may enhance the selective targeting of a particular G-quadruplex. To design a selective ligand, the electrostatics, van der Waals interactions, and specific arrangement of hydrogen bonds should be taken into account together. Thymine bases in a quadruplex seem to be easier for forming hydrogen bond interactions with the ligand, in comparison with the other



Figure 9. Superposition of crystal structure (green) and average structure calculated from the last 2 ns trajectories (orange). Water molecules binding to ligands are shown in red spheres. Green and black dashed lines represent hydrogen bonds observed in the X-ray and average structure, respectively. The hydrogen bond lengths observed in the X-ray and average structure are shown in green and black, respectively. The numbers in the parentheses are hydrogen bond scores, which are given as percentages.

bases. Therefore, the quadruplex with structurally distinctive loops containing more thymine bases can be a druggable target. Introducing chiral centers into the ligands may enhance affinity and selectivity. Recently, a chiral supramolecular62 and a chiral cyclic helicene63 have been reported to possess enantioselective recognition of the telomere G-quadruplex and enantioselective inhibition of telomerase. Ligands can induce G-quadruplex changes in the loop or flanking nucleotides, which are different from that of the ligand-free structure. This suggests that the loops are flexible and can adopt conformational changes that result in the formation of a compacted secondary structural arrangement. This will effectively compete with qudruplex binding or unwinding proteins. Acknowledgment. We thank the Natural Science Foundation of China (Grants U0832005, 90813011, 20772159, 30801436),

and the Science Foundation of Guangzhou (2009A1-E011-6) for financial support of this study. Supporting Information Available: Tables of ligand charges and figures showing the time dependence of the RMSD. This material is available free of charge via the Internet at http:// pubs.acs.org. References and Notes (1) Burge, S.; Parkinson, G. N.; Hazel, P.; Todd, A. K.; Neidle, S. Nucleic Acids Res. 2006, 34, 5402–5415. (2) Dai, J.; Carver, M.; Yang, D. Biochimie 2008, 90, 1172–1183. (3) Davis, J. T. Angew. Chem., Int. Ed. 2004, 43, 668–698. (4) Blackburn, E. H. Nature 1991, 350, 569–573. (5) Balagurumoorthy, P.; Brahmachari, S. K. J. Biol. Chem. 1994, 269, 21858–21869.

15310


(6) Huppert, J. L.; Balasubramanian, S. Nucleic Acids Res. 2005, 33, 2908–2916. (7) Todd, A. K.; Johnston, M.; Neidle, S. Nucleic Acids Res. 2005, 33, 2901–2907. (8) Huppert, J. L.; Balasubramanian, S. Nucleic Acids Res. 2007, 35, 406–413. (9) De Cian, A.; Lacroix, L.; Douarre, C.; Temime-Smaali, N.; Trentesaux, C.; Riou, J. F.; Mergny, J. L. Biochimie 2008, 90, 131–155. (10) Harrison, R. J.; Cuesta, J.; Chessari, G.; Read, M. A.; Basra, S. K.; Reszka, A. P.; Morrell, J.; Gowan, S. M.; Incles, C. M.; Tanious, F. A.; Wilson, W. D.; Kelland, L. R.; Neidle, S. J. Med. Chem. 2003, 46, 4463– 4476. (11) Tan, J. H.; Ou, T. M.; Hou, J. Q.; Lu, Y. J.; Huang, S. L.; Luo, H. B.; Wu, J. Y.; Huang, Z. S.; Wong, K. Y.; Gu, L. Q. J. Med. Chem. 2009, 52, 2825–2835. (12) Siddiqui-Jain, A.; Grand, C. L.; Bearss, D. J.; Hurley, L. H. Proc. Natl. Acad. Sci. U. S. A. 2002, 99, 11593–11598. (13) Ou, T.-M.; Lu, Y.-J.; Zhang, C.; Huang, Z.-S.; Wang, X.-D.; Tan, J.-H.; Chen, Y.; Ma, D.-L.; Wong, K.-Y.; Tang, J. C.-O.; Chan, A. S.-C.; Gu, L.-Q. J. Med. Chem. 2007, 50, 1465–1474. (14) Ma, Y.; Ou, T. M.; Hou, J. Q.; Lu, Y. J.; Tan, J. H.; Gu, L. Q.; Huang, Z. S. Bioorg. Med. Chem. 2008, 16, 7582–7591. (15) Agarwal, T.; Roy, S.; Chakraborty, T. K.; Maiti, S. Biochemistry 2010, 49, 8388–8397. (16) Bejugam, M.; Sewitz, S.; Shirude, P. S.; Rodriguez, R.; Shahid, R.; Balasubramanian, S. J. Am. Chem. Soc. 2007, 129, 12926–12927. (17) Gunaratnam, M.; Swank, S.; Haider, S. M.; Galesa, K.; Reszka, A. P.; Beltran, M.; Cuenca, F.; Fletcher, J. A.; Neidle, S. J. Med. Chem. 2009, 52, 3774–3783. (18) Dexheimer, T. S.; Sun, D.; Hurley, L. H. J. Am. Chem. Soc. 2006, 128, 5404–5415. (19) Wang, X. D.; Ou, T. M.; Lu, Y. J.; Li, Z.; Xu, Z.; Xi, C.; Tan, J. H.; Huang, S. L.; An, L. K.; Li, D.; Gu, L. Q.; Huang, Z. S. J. Med. Chem. 2010, 53, 4390–4398. (20) Hurley, L. H. Nat. ReV. Cancer 2002, 2, 188–200. (21) Neidle, S.; Parkinson, G. Nat. ReV. Drug. DiscoV. 2002, 1, 383– 392. (22) Ou, T. M.; Lu, Y. J.; Tan, J. H.; Huang, Z. S.; Wong, K. Y.; Gu, L. Q. ChemMedChem. 2008, 3, 690–713. (23) Wong, H. M.; Payet, L.; Huppert, J. L. Curr. Opin. Mol. Ther. 2009, 11, 146–155. (24) Neidle, S.; Parkinson, G. N. Biochimie 2008, 90, 1184–96. (25) Campbell, N. H.; Patel, M.; Tofa, A. B.; Ghosh, R.; Parkinson, G. N.; Neidle, S. Biochemistry 2009, 48, 1675–1680. (26) Neidle, S. Curr. Opin. Struct. Biol. 2009, 19, 239–250. (27) Li, J.; Correia, J. J.; Wang, L.; Trent, J. O.; Chaires, J. B. Nucleic Acids Res. 2005, 33, 4649–4659. (28) Dailey, M. M.; Hait, C.; Holt, P. A.; Maguire, J. M.; Meier, J. B.; Miller, M. C.; Petraccone, L.; Trent, J. O. Exp. Mol. Pathol. 2009, 86, 141– 150. (29) Chowdhury, S.; Bansal, M. J. Phys. Chem. B 2001, 105, 7572– 7578. (30) Spackova, N.; Berger, I.; Sponer, J. J. Am. Chem. Soc. 2001, 123, 3295–3307. (31) Stefl, R.; Cheatham, T. E., III; Spackova, N.; Fadrna, E.; Berger, I.; Koca, J.; Sponer, J. Biophys. J. 2003, 85, 1787–1804. (32) Fadrna, E.; Spackova, N.; Stefl, R.; Koca, J.; Cheatham, T. E., III; Sponer, J. Biophys. J. 2004, 87, 227–242. (33) Cavallari, M.; Calzolari, A.; Garbesi, A.; Di Felice, R. J. Phys. Chem. B 2006, 110, 26337–26348. (34) Sponer, J.; Spackova, N. Methods. 2007, 43, 278–290. (35) Haider, S.; Parkinson, G. N.; Neidle, S. Biophys. J. 2008, 95, 296– 311.

Hou et al. (36) Agrawal, S.; Ojha, R. P.; Maiti, S. J. Phys. Chem. B 2008, 112, 6828–6836. (37) Cavallari, M.; Garbesi, A.; Di Felice, R. J. Phys. Chem. B 2009, 113, 13152–13160. (38) Yang, D. Y.; Chang, T. C.; Sheu, S. Y. J. Phys. Chem. A 2007, 111, 9224–9232. (39) Aixiao, L.; Francois, M.; Florent, B.; Michel, D.; Baoshan, W.; Xiang, Z.; Ping, W. Eur. J. Med. Chem. 2009, 45, 983–991. (40) Haider, S. M.; Parkinson, G. N.; Neidle, S. J. Mol. Biol. 2003, 326, 117–125. (41) Clark, G. R.; Pytel, P. D.; Squire, C. J.; Neidle, S. J. Am. Chem. Soc. 2003, 125, 4066–4067. (42) Campbell, N. H.; Parkinson, G. N.; Reszka, A. P.; Neidle, S. J. Am. Chem. Soc. 2008, 130, 6722–6724. (43) Gavathiotis, E.; Heald, R. A.; Stevens, M. F. G.; Searle, M. S. J. Mol. Biol. 2003, 334, 25–36. (44) Hounsou, C.; Guittat, L.; Monchaud, D.; Jourdan, M.; Saettel, N.; Mergny, J. L.; Teulade-Fichou, M. P. ChemMedChem 2007, 2, 655–666. (45) Phan, A. T.; Kuryavyi, V.; Gaw, H. Y.; Patel, D. J. Nat. Chem. Biol. 2005, 1, 167–173. (46) Frisch, M. J. T., G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A. et al. Gaussian 03, Revision E.01; Gaussian, Inc.: Wallingford, CT, U.S., 2004. (47) Case, D. A.; Cheatham, T. E.; Darden, T.; Gohlke, H.; Luo, R.; Merz, K. M.; Onufriev, A.; Simmerling, C.; Wang, B.; Woods, R. J. J. Comput. Chem. 2005, 26, 1668–1688. (48) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. J. Comput. Chem. 2004, 25, 1157–1174. (49) Duan, Y.; Wu, C.; Chowdhury, S.; Lee, M. C.; Xiong, G.; Zhang, W.; Yang, R.; Cieplak, P.; Luo, R.; Lee, T.; Caldwell, J.; Wang, J.; Kollman, P. J. Comput. Chem. 2003, 24, 1999–2012. (50) Darden, T.; York, D.; Pedersen, L. J. Chem. Phys. 1993, 98, 10089– 10092. (51) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327–341. (52) Humphrey, W.; Dalke, A.; Schulten, K. J. Mol. Graph. 1996, 14, 33–38, 27-28. (53) Fogolari, F.; Brigo, A.; Molinari, H. Biophys. J. 2003, 85, 159– 166. (54) Kollman, P. A.; Massova, I.; Reyes, C.; Kuhn, B.; Huo, S.; Chong, L.; Lee, M.; Lee, T.; Duan, Y.; Wang, W.; Donini, O.; Cieplak, P.; Srinivasan, J.; Case, D. A.; Cheatham, T. E. Acc. Chem. Res. 2000, 33, 889–897. (55) Hazel, P.; Parkinson, G. N.; Neidle, S. Nucleic Acids Res. 2006, 34, 2117–2127. (56) Zeng, J.; Li, W.; Zhao, Y.; Liu, G.; Tang, Y.; Jiang, H. J. Phys. Chem. B 2008, 112, 2719–2726. (57) Gilson, M. K.; Sharp, K. A.; Honig, B. J. Comput. Chem. 1987, 9, 327–335. (58) Still, W. C.; Tempczyk, A.; Hawley, R. C.; Hendrickson, T. J. Am. Chem. Soc. 1990, 112, 6127–6129. (59) DeLano, W. L. The PyMOL Molecular Graphics System; DeLano Scientific: San Carlos, CA, U.S., 2002. (60) Gu, J.; Leszczynski, J.; Bansal, M. Chem. Phys. Lett. 1999, 311, 209–214. (61) Pagano, B.; Mattia, C. A.; Cavallo, L.; Uesugi, S.; Giancola, C.; Fraternali, F. J. Phys. Chem. B 2008, 112, 12115–12123. (62) Yu, H.a J.; Zhao, C. Q.; Chen, Y.; Fu, M. L.; Ren, J. S.; Qu, X. G. J. Med. Chem. 2010, 53, 492–498. (63) Shinohara, K.; Sannohe, Y.; Kaieda, S.; Tanaka, K.; Osuga, H.; Tahara, H.; Xu, Y.; Kawase, T.; Bando, T.; Sugiyama, H. J. Am. Chem. Soc. 2010, 132, 3778–3782.

JP106683N

New Insights into the Structures of LigandâQuadruplex Complexes

Recommend Documents

New Insights into the Structures of LigandâQuadruplex Complexes