The interactions and recognition of cyclic peptide mimetics of Tat with HIV-1 TAR RNA: a molecular dynamics simulation study

The interaction of HIV-1 trans-activator protein Tat with its cognate trans-activation response element (TAR) RNA is critical for viral transcription and replication. Therefore, it has long been considered as an attractive target for the development of antiviral compounds. Recently, the conformationally constrained cyclic peptide mimetics of Tat have been tested to be a promising family of lead peptides. Here, we focused on two representative cyclic peptides termed as L-22 and KP-Z-41, both of which exhibit excellent inhibitory potency against Tat and TAR interaction. By means of molecular dynamics simulations, we obtained a detailed picture of the interactions between them and HIV-1 TAR RNA. In results, it is found that the binding modes of the two cyclic peptides to TAR RNA are almost identical at or near the bulge regions, whereas the binding interfaces at the apical loop exhibit large conformational heterogeneity. In addition, it is revealed that electrostatic interaction energy contributes much more to KP-Z-41 complex formation than to L-22 complex, which is the main source of energy that results in a higher binding affinity of KP-Z-41 over-22 for TAR RNA. Furthermore, we identified a conserved motif RRK (Arg-Arg-Lys) that is shown to be essential for specific binding of this class of cyclic peptides to TAR RNA. This work can provide a useful insight into the design and modification of cyclic peptide inhibitors targeting the association of HIV-1 Tat and TAR RNA.


Introduction
The interaction of HIV-1 trans-activator Tat protein with its cognate trans-activation response element (TAR) RNA plays a central role for viral transcription and replication (Biswas et al., 1995;Gatignol et al., 1991). HIV-1 TAR RNA, located in the 5′-untranslated region (UTR) of nascent viral transcripts, is characterized by two helical stems separated by a tri-nucleotide bulge and a hexa-nucleotide apical loop capping the terminal stem (Aboulela, Karn, & Varani, 1995;Ratnasabapathy et al., 1990) (Figure 1(A)). The association of Tat with TAR RNA triggers the recruitment of the essential cellular cofactor transcription elongation factor b (P-TEFb) to RNA polymerase II (Karn & Stoltzfus, 2012;Peterlin & Price, 2006). Consequently, the progressive transcription elongation is promoted, substantially facilitating viral gene expression and virion production.
A number of studies have focused on TAR RNA and Tat molecule systems through biochemical techniques. At the same time, the theoretical approaches provided significant information as well. Molecular dynamics (MD) simulations have been successfully applied to the interpretation and prediction of structurefunction relationships of biological macromolecules (Chen and Xiao, 2006;Cao and Wang, 2010;Su et al., 2011;Xu et al., 2011;Zhang and Liu, 2011). To our knowledge, Nifosi, Reyes, and Kollman (2000) first employed MD simulations to study the dynamic behavior of HIV-1 TAR RNA and its complex with argininamide. Subsequently, Reyes et al. (2001) investigated the differences of the structural dynamics and binding specificity between the two BIV TAR RNA complexes with Tat fragments of different lengths. Kulinski, Kulinska, Sarzynska, and Adamiak (2005) studied HIV-1 TAR RNA with the extensive MD simulation and confirmed the coexistence of alternative functional conformations of the apical loop. Mu and Stock (2006) presented a unique "minor groove" binding mode of a tripeptide to HIV-1 TAR RNA by combining docking and MD simulations. Musselman et al. (2007) elaborated the dynami-cal correlations within HIV-1 TAR RNA utilizing an isotropic reorientational eigenmode dynamics (iRED) analysis.
Apparently, these researches based on MD simulations concentrated on either the free TAR RNA or the TAR RNA associated with linear peptides or small molecules. With respect to the HIV-1 TAR RNA complexed with the aforementioned cyclic beta-hairpin peptidomimetics of Tat (cyclic peptide for short in the following parts), the relationships between the properties of structure, dynamics, and energetics and the high affinity and specificity recognition have not yet been clearly understood. In this work, we chose two representative cyclic peptides ( Figure 1B and Table 1), termed L-22 and KP-Z-41, respectively, as case studies (Davidson et al., 2009;Davidson, Patora-Komisarska, Robinson, & Varani, 2011). These two peptides were selected for their potent inhibitory activities and available complex structures with TAR RNA. Here, we delineated a detailed picture of the interactions between them and HIV-1 TAR RNA via all-atom MD simulations with explicit water solvent. Based on the two sets of simulation runs, we revealed the differences in binding mode and affinity of these two cyclic peptides to TAR RNA. Moreover, a conserved motif RRK (Arg-Arg-Lys) was identified which plays a central role in the high affinity and specificity recognition of TAR RNA by this class of cyclic peptide inhibitors.  -41 (top). In this article, the moiety stabilized by D Pro-L Pro is referred to as P-end ("posterior" end) and the opposite is termed A-end ("anterior" end) for simplicity. (C) L22/TAR RNA complex structure. (D) KP-Z-41/TAR RNA complex structure. The backbone of the bound peptide is shown with stick model and the bulge (U23 to U25) and apical loop (C30 to A35) of TAR RNA are highlighted in dark gray. Table 1. Sequence alignment between L-22 and KP-Z-41.

Peptide
Residue The corresponding residues between L-22 and KP-Z-41 are aligned to represent similarity in position upon binding HIV-1 TAR RNA.

Molecular dynamics simulations
The two initial complex structures (Figure 1(C) and (D)) were derived from Protein Data Bank (PDB entries: 2KDQ for L-22/TAR RNA complex (Davidson et al., 2009) and 2KX5 for KP-Z-41/TAR RNA complex (Davidson, Patora-Komisarska, Robinson, & Varani, 2011). Each complex system contains an ensemble of 10 structure models resolved by nuclear magnetic resonance spectroscopy (NMR), in which the structure with the minimum cumulative root-mean-square displacement (RMSD) from the rest of the ensemble was used as the starting coordinate. In addition, the respective monomer coordinates were obtained by splitting the selected complex models. For each type of peptide-RNA complex, three separate MD simulations involving the complex, receptor TAR RNA, and ligand cyclic peptide were conducted independently. The two sets of simulations were setup based on the same procedures as follows. The AMBER10 program package (Case et al., 2005) and the force field parm99SB (Perez et al., 2007) were chosen for simulations. Every structure was immersed in a truncated octahedral box filled with the TIP3P water molecules, keeping a minimum distance of 10 Å between the solute and each face of the box. The proper type and number of counterions (sodium or chloride ions) were added to neutralize the charged solute molecule(s). The solvated systems were subjected to a thorough energy minimization and then were gradually heated to 300 K with weak harmonic constraints over a period of 1 ns. The production run without constraints was extended up to 40 ns at constant temperature (300 K) and constant pressure (1 atm) with an integration time step of 2 fs. The trajectory of production phase was saved every 1 ps, thus yielding 40,000 snapshots for each system. The long-range electrostatic interactions were treated with particle mesh Ewald method (Darden, York, & Pedersen, 1993), while covalent bonds involving hydrogen atoms were constrained by the SHAKE algorithm (Ryckaert, Ciccotti, & Berendsen, 1977).

Binding free energy calculation and free energy decomposition
The binding free energy calculations were carried out based on molecular mechanics Poisson-Boltzmann surface area (MM-PBSA) and molecular mechanics generalized Born (MM-GBSA) approaches, using the pbsa program as implemented in AMBER10. In principle, the binding free energy can be evaluated on the basis of either the "single trajectory protocol" or the "separate trajectory protocol" . Compared with the former protocol, the latter is computationally more expensive, but it is more accurate due to the conformational energy of adaption taken into account explicitly . In this study, the separate trajectory protocol was adopted in view of the high flexibility of RNA. The binding free energies of cyclic peptides to TAR RNA were estimated according to the following equation: where the absolute free energies of the complex ( " G complex ), peptide ( " G peptide ), and RNA ( " G RNA ) were averaged over 3000 snapshots taken at 10 ps intervals from the last 30 ns equilibration MD trajectories. The contributions of conformational entropy were also considered through the normal mode analysis. More details about the binding free energy calculation are included in Supplementary material. In addition, the free energy decomposition on a per-residue basis was performed with the MM-GBSA approach.

Assessment of structural stability and flexibility during simulations
Conformational changes of L-22 complex, KP-Z-41 complex, and their respective isolated monomers over the course of simulations were monitored by the massweighted heavy atom RMSDs in reference to the starting structures as a function of time. From Figure 2, the RMSD for each system has an updrift during the early MD simulation. It can be partly attributed to the sudden removal of constraints on the solute molecule when the simulation progressed into the production phase from the gradual heating stage (see Materials and methods). Noticeably, those free monomers show a larger fluctuation in RMSD as compared to the two complexes within the first $5 ns of the simulations, which is mainly because the monomer simulations were started from the conformations separated from the complexes. On the whole, all these simulation runs tend to reach a plateau since $10 ns. Hence, all the following analyses were performed on the last 30 ns equilibrated MD trajectories.
To track the extent of variation in flexibility of individual residues/nucleotides upon binding, the per-residue/nucleotide root-mean-square fluctuations (RMSFs) of heavy atoms were computed over respective simulation runs (see Figure 3). As shown in Figure 3(A) and (B), the overall flexibility of each peptide is reduced upon TAR RNA binding. Rather intriguingly, the fluctuations occurring at the P-end moiety of each bound peptide are much larger than those observed for the corresponding part in the free one. Upon examination, we found that the A-end freedom degrees of both peptides are largely constrained due to extensive contacts with the major groove binding pockets; by contrast, their P-ends deviating from the major groove are relatively more free to oscillate owe to the lack of structural constrains imposed by TAR (see Supplementary material Table SI). This can account for the exceptional fluctuations in the P-end of each bound peptide.
According to the per-nucleotide RMSF plots for apo-TAR RNA (free state) and holo-TAR RNA (bound state) in Figure 3(C) and (D), it can be seen that TAR RNA overall becomes much less flexible upon binding of peptide. For the two regions apical loop and bulge, they exhibit large fluctuations in both free and bound states, which can be expected due to lack of hydrogen bonding constrains within them. Comparing the fluctuations of the two holo-TAR RNAs, three nucleotide sites should be specially pointed out: G33 and A35 in the apical loop and U23 in the bulge. Both G33 and A35 show the opposite fluctuating tendencies, while U23 shares the lowest fluctuations. In the following two sections concerning the analyses of hydrogen bonding and binding interfaces of complexes, we will give an account for the origins of these unusual fluctuations at length.

Comparative analyses of hydrogen bonding between complexes
To reveal similarities and differences in hydrogen bonding pattern, we carried out comparative analyses of the hydrogen bonds (H-bonds) with an occupancy > 45% over both complex trajectories. Table 2 presents all these intermolecular (upper part) and partial intra-RNA Hbonds (lower part). From Table 2, it is found that the guanidinium groups of Arg3 and Arg5 in L-22 make stable H-bonds with the bases of G26 and G28, respectively. Similarly, the stable H-bond pairs also form between the side chains of Arg5 and Arg7 in KP-Z-41 and the bases of G26 and G28, respectively. The Hbonding of G28 to arginine observed here is in agreement with the NMR data (Davidson et al., 2009(Davidson et al., , 2011. The H-bonds between G26 and arginine were not imme-diately established in the NMR experiments, but this nucleobase is in proximity to the residue side chain in both NMR ensembles, and analogous H-bonding interactions (G11 and Arg3) are present in the BIV TAR RNA complex with a peptidomimetic BIV-2 (just one residue replaced compared to L-22) (Leeper, Athanassiou, Dias, Robinson, & Varani, 2005). Thus, it is likely that this type of bis-arginine (bis-Arg) among these cyclic peptides is associated with specific nucleotide recognition. Coincidently, both side chains of Thr4 in L-22 and the corresponding Gln6 in KP-Z-41 form H-bonds with the same phosphate oxygen of bulge nucleotide C24. Arg1 in L-22 and Arg3 in KP-Z-41 make contacts with the phosphate of A22 in the lower helix, and the latter is additionally H-bonded to the phosphate of G21. These H-bondings with the nucleotides in the bulge and lower stem of TAR RNA accord well with the intermolecular NOE measurements (Davidson et al., 2009(Davidson et al., , 2011. Arg15 in KP-Z-41, located opposite Arg3, interacts with the phosphate of C37 in the upper helix. It should be noted that the equivalent residue to Arg15 in KP-Z-41 does not exist in L-22 (see Table 1). At the apical loop, Arg9 of L-22 makes an H-bond with the ribose hydroxyl group of G33, while Arg11 of KP-Z-41 forms contacts with phosphate oxygens of G34, which are consistent with the NMR data (Davidson et al., 2009(Davidson et al., , 2011. In contrast to the lack of any H-bonds mediated by Arg10 in KP-Z-41 complex, the corresponding Arg8 in L-22 complex is able to form five H-bonds with nucleotides C30, G33, and G34 in the apical loop (see below Figure 4(B)), which confirms the observed low flexibility of Arg8 upon binding (Figure 3(A)). As to RNA, the phosphate group of loop nucleotide A35 makes the stable Hbond with the backbone N atom of Cys13 in KP-Z-41 (see below Figure 4(C)), leading into the diminished mobility of A35 observed in KP-Z-41 complex run rather than in L-22 complex run (Figure 3(C) and (D)). Finally, summing up the H-bonds at interface, we can see for KP-Z-41 complex, five residues Arg5, Gln6, Gly9, Arg11, and Cys13 are involved in H-bonding with RNA apical loop or bulge and for L-22 complex, the residues are Arg3, Thr4, Arg8, and Arg9. It is these extensive H-bonding interactions that collectively drive TAR RNA to wrap around the cyclic peptides. Additionally, we investigated the H-bonds of base pairs at or around the bulge (see lower part of Table 2). As indicated by the high occupancy, the Hoogsteen H-bond alignment between U23 and A27 and the Watson-Crick base pair between A27 and U38 are well maintained in both complex runs, in good accordance with NMR NOE data (Davidson et al., 2009(Davidson et al., , 2011. This suggests that the base triple U23:A27-U38 (see Figure S1) is defined unequivocally. The base triple is considered as an essential structural requirement for the bound conformation of TAR RNA characterized by coaxial stacking between two adjacent helices (Davidson et al., 2009;Ye, Kumar & Patel, 1995). Evidently, the participation of U23 in the formation of this base triple greatly reduces the flexibility of this bulge nucleotide, thereby causing an unusually low fluctuation (Figure 3(C) and (D)). In addition, we have noticed the two base pairs (A22-U40 and G26-C39) around the bulge that are much more relaxed in the free TAR RNA (Aboulela, Karn, & Varani, 1995;Musselman, Al-Hashimi, & Andricioaei, 2007) are also well defined. From Figure S1, the establishment of base triple necessarily entails the approaching of U23 to A27, during which the 5′-end backbone of U23 through U25 bulges upward and thus, the A22 gets close to U40. Meanwhile, the G26 near the 3′-end of the bulge is moderately collapsed to fine-tune its orientation against the C39. As a result, both base pair planes of A22-U40 and G26-C39 become more coplaner. Intriguingly, unlike KP-Z-41 complex, the apical loop nucleotide G33 in L-22 complex can make the steady Watson-Crick base pair with C29 in the upper helix at the expense of disrupting C29-G36 base paring, whilst making several H-bonds with Arg8 and Arg9 (see below Figure 4(B)). Thus, the dynamic level of this RNA loop nucleotide G33 is reduced remarkably upon binding of L-22 (Figure 3(C)), in marked contrast to that observed for KP-Z-41 complex (Figure 3(D)).

Comparison of peptide-RNA binding interfaces
To further dissect the differences between the two peptide-RNA binding interfaces, we performed a structural comparison between the two lowest energy  Table 2. conformations obtained from the respective complex trajectories. Figure 4(A) illustrates the structure superimposition between TAR RNAs bound to L-22 and KP-Z-41, respectively, with a pairwise backbone RMSD of 3.5 Å (as compared to 4.6 Å computed for the two initial NMR structures). From Figure 4(A), both the RNA structures globally converge well toward each other, including the region at the interhelical junction (marked with rectangle). The most distinct difference between them lies in the apical loop moiety (marked with ellipse). Hence, we focused on the structural organization around the loop (see Figure 4(B) and (C)). As shown in Figure 4 (B), Arg8 in L-22 is sandwiched between G32 and G33 to form cation-pi interactions while the Arg9 is stacked below G34. Interestingly, the corresponding Arg10 and Arg11 in KP-Z-41 are also involved in analogous stacking interactions with G32 and G34, respectively (Figure 4 (C)); however, the orientations of both bases against the stacked residues are distinct from those observed in L-22 complex (Figure 4(B)). In addition, it can be seen in Figure 4(C) that the G33 at the apex of loop is away from the upper helix major groove and most solvent-exposed. This observation shows structural evidence for the large fluctuation of G33 in KP-Z-41/TAR complex (Figure 3 (D)). Similarly, the nucleotide A35 in L-22/TAR complex is also exposed to the solution (Figure 4(B)) and thereby exhibiting a high flexibility (Figure 3(C)). With regard to the holo-TAR RNA with L-22, it is the very juxtaposition of Arg8 guanidinium with bases of G34 and C30 via H-bonds as well as the base pairing between G33 and C29 that leads into the apparent distortion of apical loop in a "zigzag" course ( Figure 4(B)), which is in contrast with the "bow" shape observed in KP-Z-41 complex (Figure 4(C)). Together with the Hbonding analyses above, it is suggested that the structural arrangements between L-22 complex and KP-Z-41 complex are nearly identical at the interhelical junction, whereas the local conformations at the apical loop exhibit large heterogeneity.

Binding free energies and residual energy decomposition
To obtain insights into the contributions of different energy terms to the affinity of cyclic peptides for HIV-1 TAR RNA, the binding free energies for L-22/TAR RNA and KP-Z-41/TAR RNA were computed with both MM-PBSA and MM-GBSA methods (for comparison). Table 3 contains various contributions to the binding free energy as calculated using the "separate trajectory protocol" mentioned in Materials and methods. Overall, both GB and PB models yielded qualitatively very similar results. The calculated total binding free energies are remarkably favorable for the formation of both L-22 complex (À33.59 and À46.23 kcal/mol for PB and GB models, respectively) and KP-Z-41 complex (À60.09 and À66.55 kcal/mol, respectively). As a result, KP-Z-41 is predicted to bind to target TAR RNA with a higher affinity than L-22, which is in qualitative agreement with the experimental data (Davidson et al., 2009(Davidson et al., , 2011. Average over 3000 (300 in the case of entropy calculations) snapshots. c Standard error of mean value in parenthesis.
As reported in Table 3, the molecular mechanics electrostatic energy (ÁELE) strongly favors the association of cyclic peptides with TAR RNA, while the electrostatic solvation terms (ÁPB CAL and ÁGB CAL ) significantly oppose complex formation. Both are rather substantial in magnitude but compensate each other, overall resulting in favorable total electrostatic contributions. From PB model, the absolute values of total electrostatic energies (ÁPB ELE ) reach to $17 kcal/mol for L-22/TAR RNA and $52 kcal/mol for KP-Z-41/TAR RNA; from GB model, the corresponding values (ÁGB ELE ) are $30 and $58 kcal/mol for these two complexes. The van der Waals interaction (ÁVDW) energy, the major contributor, is favorable to the binding free energy by 65-70 kcal/mol. Also, the nonpolar solvation term (ÁSUR NP ) contributes favorably to the binding affinity by 11-12 kcal/mol. The contributions from internal energy (ÁINT), however, are unfavorable by 8-9 kcal/mol, indicating that the conformational changes of the binding partners lead to internal strain upon complex formation. In addition, the changes in the conformation entropy (ÀT ÁS) disfavor both peptides binding to TAR RNA by up to 57-60 kcal/mol. Apparently, taking the entropic contribution into account (more reasonable) shifts the binding free energy to much more positive values. This is as expected given that the complex formation is accompanied with a considerable loss of conformational flexibility especially for TAR RNA.
Adding together the above terms (marked in italic in Table 3), ÁVDW +ÁINT +ÁPB ELE (ÁPB ELE ) +ÁSUR NP -TÁS TOT , yield the total gains in free energy of cyclic peptides binding to TAR RNA. Comparing the magnitudes of all these energy terms between the two complexes, we found that ÁVDW, ÁINT, ÁSUR NP , and TÁS TOT make comparable contributions to the binding free energies of the two complexes. By contrast, the total electrostatic terms (ÁPB ELE and ÁGB ELE ) differ markedly between complexes, according to both PB and GB models. It is evident that the total electrostatic contribution favors the association of HIV-1 TAR RNA with KP-Z-41 much more than that with L-22. This is the main source of energy that results in a higher affinity of KP-Z-41 over L-22 for target TAR RNA. Figure S2 illustrates the electrostatic potentials mapped on the surfaces of the peptides from their respective lowest energy complex conformations. In comparison with L-22 ( Figure S2A), KP-Z-41 ( Figure S2B) apparently forms a wider area of positive potential (blue) at the RNA binding interface and this area was examined to match well with the highly electronegative RNA via electrostatic complementarity. The extra positive potential regions at KP-Z-41 are mostly owing to the two distal residues Arg1 and Arg15 at P-end ( Figure S2B). These can account for the more favorable electrostatic contribution to KP-Z-41 complex formation as compared to L-22 complex.
In order to identify the key residues important for association, the free energy decomposition on residual basis was performed through MM-GBSA method. As shown in Figure 5, the basic resides of each peptide substantially contribute to the binding to TAR RNA, with an energy value <À7 kcal/mol. Overall, the charged resides of KP-Z-41 make almost equivalent contributions (around À9 kcal/mol) to the binding free energy, whereas the binding energy per residue differs markedly between charged residues in L-22. In particular, Arg3 and Arg8 in L-22 make the most energetic contributions, almost twice as much as those of other basic residues. It is unsurprising since these two resides participate in the Hbonding network (Table 2) and stacking interactions with adjacent nucleotides (Figure 4(B)). With respect to KP-Z-41, it is a point worth noting that its two P-end residues Arg1 and Arg15 are considerable contributors although relatively away from the RNA phosphate backbone, which is consistent with the observation that the electrostatic force can extend up to 11 Å from the RNA phosphate group (Bahadur, Kannan & Zacharias, 2009;Garcia-Garcia & Draper, 2003;Law et al., 2006). Given the lack of the corresponding residues in L-22, it is indicated that the enhanced electrostatic contribution within KP-Z-41 complex (as mentioned above) can be partly ascribed to the two residues Arg1 and Arg15 in the extended P-end.

Discussion
Essential role of RRK motif in cyclic peptide/TAR RNA recognition With respect to this class of lead peptides, it is interesting to explore whether there exists a certain consensus motif that is important for potent and specific binding to TAR RNA. This is also useful for the design of cyclic peptide inhibitors. On the ground of our MD studies and available experimental data, here, we carefully identified Figure 5. Per-residue contributions of L22 (light gray) and KP-Z-41 (dark gray) to the binding free energies. Bars represent mean and error bars represent standard error of the mean. The residues of the two peptides are aligned to show similarity in position upon binding to TAR RNA according to Table 1. a conserved motif Arg-Arg-Lys (RRK, see Figure 6) that is shown to be critical for the high affinity and specificity association with TAR RNA. The residue Lys (Lys6 in L-22 and Lys8 in KP-Z-41, see Figure 6(C) and (F)) at A-end serves to guide the correct location of peptide relative to TAR RNA, while the bis-Arg sandwich (Arg3-Arg5 in L-22 and Arg5-Arg7 in KP-Z-41, see Figure 6(A) and (D)) adjacent to that mentioned Lys chiefly contributes to the specific recognition. This motif is proposed based on the following analyses.
On the one hand, the bis-Arg sandwich is involved in specific base recognition, which can be explained in terms of hydrogen bonding and electrostatic features at the binding interface. As revealed by the H-bond analysis (Table 2), for each complex the bis-Arg side chains of A-end exclusively make stable H-bonds with the nucleobases of G26 and G28, respectively. Other basic residues, however, either make contacts with RNA phosphate backbones or form inconsistent H-bonds with the nucleobases. Moreover, we carried out electrostatic potential calculations on around the bis-Arg sandwich and observed two negative potential pockets above and below the base plane of U23 (Figure 6(B) and (E)), respectively. It is clear that the bis-Arg sandwich helps to establish the base triple by buttressing the conserved nucleotide U23 in the flexible bulge. Additionally, it is noteworthy that the two electrostatic potential pockets are dominantly formed by electronegative groups not from the backbone phosphates, but from the bases in the major groove near the bulge (Figure 6(B) and (E)). In previous studies, the analogous bis-Arg signature is also found important for the binding of Tat-derived nonclosed peptides (Calabro, Daugherty & Frankel, 2005;Puglisi et al., 1995;Ye, Kumar, & Patel, 1995) and other cyclic peptides Leeper et al., 2005) to BIV TAR RNA. Taken together, these results demonstrate that the bis-Arg sandwich is of essence in conferring recognition specificity.
On the other hand, the residue Lys at A-end is related to the proper positioning of the cyclic peptide relative to TAR RNA. To our knowledge, all known cyclic peptides Davidson et al., 2009Davidson et al., , 2011Leeper et al., 2005) adopt the opposite orientation toward TAR RNA compared to nonclosed peptides (Aboulela et al., 1995;Puglisi, Chen, Blanchard, & Frankel, 1995;Ye et al., 1995). On this point, we explain that the residue Lys at A-end plays a primary role in orienting the correct face of cyclic peptides toward TAR RNA. As illustrated in Figure 6(C) and (F), for each peptide the positively charged terminus of Lys at A-end penetrates into an electrostatic potential pocket composed of the phosphate groups of nucleotides C24 through G28 at the interhelical junction. Thus, the Lys is regarded as a key residue for the initial electrostatic attraction (before binding) and neutralization (upon binding) with TAR RNA, and this kind of role regarding Lys apparently  Table 2. (B) Surface potential map around the two electrostatic pockets as outlined in (A). One is formed between nucleotides A22 and U23 and the electronegative groups of U23 and A22 bases and the phosphate backbone between them constitute the charge cluster. The other is above U23, composed of electronegative groups on U23 (O2, N3, O4), A27 (N6, N7), and G28 (O6, N7). (C) The positively charged side chain of Lys6 at L-22 A-end is oriented toward a negative potential pocket formed by the phosphate groups of U24 through C28. (D), (E), and (F) The Arg5-Arg7-Lys8 combination found in KP-Z-41/TAR RNA complex with the same drawing and coloring schemes as in (A), (B), and (C), respectively. Electrostatic potentials were color-coded on the molecular surfaces from red (À20 kT) to blue (+20 kT). All the figures were prepared using the two lowest energy structures from the respective complex trajectories. cannot be conferred by those hydrophobic residues at Pend. Moreover, the experiments found that the substitutions of A-end amino acid lysine with some other residues affect the binding affinity to TAR RNA to a different extent Davidson et al., 2009Davidson et al., , 2011.

Conclusion and perspective
In this work, we focused on the two representative cyclic peptide inhibitors namely L-22 and KP-Z-41 and delineated detailed characteristics of the interactions and recognition between them and HIV-1 TAR RNA via MD simulations. It is found that the TAR RNA base triple (U23: A27-U38) near the bulge is well defined over the course of both complex simulation runs. The binding modes of these two peptides to TAR RNA are almost identical at the interhelical junction, whereas the binding interfaces at the apical loop exhibit large conformational heterogeneity. In addition, we revealed that the electrostatic interaction energy contributes much more to KP-Z-41 complex formation than to L-22 complex, which is partly attributed to the two P-end residues Arg1 and Arg15 in KP-Z-41. This is the main source of energy that results in a higher affinity of KP-Z-41 than L-22 for TAR RNA. Furthermore, we identified a conserved motif RRK that is shown to be essential for the high affinity and selectivity recognition of TAR RNA by this class of lead peptides.
The results of this work could provide some implications for the improvement in activity of Tat peptidomimetic inhibitors. The basic RRK combination identified here is crucial for the maintenance of conserved contacts of peptides with TAR RNA; hence, this motif should be exploited for future development of cyclic peptides. Meanwhile, its location on the beta-hairpin scaffold needs to be optimized to achieve desired potency. Another point worth noting is that, compared with the intimate contacts of RRK motif with the area near the bulge, the interactions around the apical loop are relatively weak, as high fluctuations still observed in this region. This implies that the additional binding affinity for TAR could be obtained by introducing the residues with larger side chains (e.g. Trp, Phe, and His) in peptides to form more extensive interactions with RNA loop.