Comparison of Linear vs. Cyclic RGD Pentapeptide Interactions with Integrin αvβ3 by Molecular Dynamics Simulations

Abstract

Simple Summary

The integrin α_vβ₃-RGD motif interaction plays a key role in the progression of malignant tumor. Although two typical cyclic and linear RGD short peptides have been widely used in tumor diagnosis and therapy, little is known about the internal dynamic mechanism for different configurations of RGD peptides with different affinities interacting with the integrin α_vβ₃. Our results showed that the cyclic RGD peptide had a more stable configuration in binding to integrins α_vβ₃, which depended on the higher binding energy and higher static electrical energy, especially in the interaction between Asp^RGD-MIDAS. The steered molecular dynamics simulation showed a stronger interaction for the cyclic RGD-integrin α_vβ₃ system than the linear one, with a larger dissociation force (average peak force) and more time to dissociate. Our findings provide insights into the dynamics of integrin α_vβ₃ interactions with linear and cyclic RGD ligands and offer some new therapeutic approaches for the design and development of novel antitumor drugs.

Abstract

Integrin α_vβ₃ interacting with the short Arg-Gly-Asp (RGD) motif plays a critical role in the progression of several types of tumors. However, the effects of the RGD structure (cyclic or linear) with integrin α_vβ₃ at the atomic level remain poorly understood. Here, we performed association and dissociation dynamic simulations for integrin α_vβ₃ in complex with a linear or cyclic pentapeptide by steered molecular dynamics simulations. Compared with cyclic RGD, the linear RGD peptide triggers instability of the configurational changes, mainly resting with the RGD domain due to its flexibility. The main interaction energy between Mg²⁺ and cyclic RGD is much stronger than that of the linear RGD system by the well shield to lessen attacks by free water molecules. The force-dependent dissociation results show that it is easier for linear RGD peptides to leave the active site and much quicker than the cyclic RGD ligand, whereas it is harder to enter the appropriate active binding site in linear RGD. The Ser¹²³-Asp^RGD bond may play a critical role in the allosteric pathway. Our findings provide insights into the dynamics of α_vβ₃ interactions with linear and cyclic RGD ligands and contribute to the application of RGD-based strategies in preclinical therapy.

1. Introduction

Tumor angiogenesis, the formation of new blood vessels, is a critical process for tumor growth and metastasis [1,2,3]. Many types of adhesion molecules are involved in tumor angiogenesis [4]. Among cell adhesion molecules (CAMs), integrins are important CAMs that link the extracellular matrix (ECM) and the cytoskeleton and participate in adhesive events during various cancer stages, such as tumor growth, invasion, and metastasis [5]. Among all integrins, integrin α_vβ₃ plays a crucial role in angiogenesis and tumor metastasis, which is widely expressed in many kinds of human tumor biopsy samples but not in vessels in normal tissues [6,7]. As a receptor, integrin α_vβ₃ can specifically recognize one or more Arg-Gly-Asp tripeptide motifs, called RGD peptides. In the past decade, RGD peptides have been designed to target integrin for cancer therapy [8], and many radiolabeled RGD peptides have been used as α_vβ₃ -targeting tumor imaging agents [7,9]. In addition, RGD-based anticancer strategies have good prospects in the field of cancer diagnosis and therapy, even in the field of tissue regeneration, including cornea repair, artificial vascularization, and bone tissue regeneration [10,11]. Thus, the use of RGD-modified peptides to inhibit the activity of integrin α_vβ₃ is a promising strategy for tumor targeting.

Interaction studies of integrin α_vβ₃ with RGD peptides have greatly advanced our understanding of the dynamic binding process and molecular mechanism [12,13,14]. The activation mechanism of integrin α_vβ₃ binding to fibronectin has also been revealed to some extent recently [15]. Under different physiological conditions, integrin α_vβ₃ can form different dynamic conformations with distinct affinity [16]. When the RGD peptide binds to the active region of integrin α_vβ₃, a shallow crevice is located between the α_v and β₃ subunits, which changes the conformation of integrin and transduces signals from the ECM to the cytoplasm. There are two important interaction sites. One is the salt bridges formed between Arg^RGD residue of RGD with Asp²¹⁸ and Asp¹⁵⁰ residues of the β-propeller subunit, and the other is the Asp-MIDAS interaction, in which Asp^RGD carboxylate oxygen atoms of RGD coordinate with the metal ion at the metal ion-dependent adhesion site (MIDAS) of the βA domain. The breaking of the Asp-MIDAS ligand interaction corresponds to the major force peak as the largest barrier to unbinding [17]. The discovery of the structural mechanism of interaction between integrins and RGD-containing peptide ligands has contributed to the rational design of drugs that effectively inhibit integrin activation.

In preclinical studies, various RGD peptides or peptidomimetics are designed as selective integrin inhibitors [9,18,19]. The steric conformation of the peptide and the structural features of the RGD ligand have become important because these factors influence the affinity between RGD and integrin α_vβ₃. Recent studies have shown that cyclic RGD peptides commonly help improve the binding properties of RGD peptides, whereas linear RGD peptides are easily susceptible to chemical degradation [9,20,21]. Therefore, the cyclic RGD is more stable and more active. However, differences in the dynamic properties of integrin α_vβ₃ interacting with linear or cyclic RGD peptides remain unknown at the molecular, and even at the atomic, level.

Integrins have been studied as therapeutic targets in many diseases due to their involvement in modulating various vital physiological and pathological processes, including proliferation, survival, differentiation, and migration. Recently, several crystal structures of integrins with or without ligands have been analyzed [22]. Although the details of the molecular interaction mechanism of these structures are not clear, molecular dynamics (MD) simulation may have the ability to derive these insights at the molecular level, which will help better design drugs targeting integrin [23,24]. In this study, we used computational MD and SMD to investigate the binding divergence of different RGD structures (cyclic or linear) interacting with integrin α_vβ₃. The corresponding structural stability and energy fluctuation profiles were obtained. These results provide insights into the differences in the dynamic properties of RGD-containing peptides with the integrin objective, which will help better understand the field of RGD-mediated drug delivery and imaging constructs.

2. Materials and Methods

2.1. System Setup

Two simulated systems of integrin α_vβ₃ binding with cyclic or linear RGD were set up for the MD simulations (Figure 1). The simulation for the cyclic system was started from the crystal structure 1L5G (resolution: 3.20 Å) containing the complex of integrin α_vβ₃ headpiece with the cyclic RGD peptide [25]. The linear MD system was then obtained from the same crystal structure but only the RGD sequence was linearized (Figure 1D,E). To reduce the system size, we only used the headpiece of integrin α_vβ₃, including the βA domain of the β₃ subunit (residues 110 to 353) and the β-propeller domain of the α_v subunit (residues 1 to 438) (Figure 1A,B). Before MD simulations, the complex was solvated in a 126 × 107 × 81 ų water box together with Na⁺ and Cl⁻ ions to neutralize the system at a 150 mM ionic concentration, resulting in 104,246 atoms (Figure 1C). Due to the lack of Mn²⁺ force field parameters, a universal method in which Mn²⁺ is replaced by Mg²⁺, metal ion-binding sites were occupied by Mg²⁺ ions in the integrin α_vβ₃ headpiece in our simulations.

2.2. Simulation Procedure and Parameters

All MD simulations were performed with the program NAMD 2.6, using periodic boundary conditions and the CHARMM22 force field [26,27]. The particle-mesh Ewald method was used to calculate the full electrostatic calculations and van der Waals interactions were evaluated using a smooth cutoff (12 Å). Visualization, molecular graphics, and analyses of simulations, including root mean square deviation (RMSD), the distance between two atoms or the centers of mass of the helix, and solvent accessible surface area (SASA) with a 1.4 Å probe radius, were measured in Tcl within the program VMD [28]. A hydrogen bond was defined when the donor–acceptor distance was less than 0.35 nm, and the donor–hydrogen–acceptor angle was less than 30°. Occupancy was evaluated by the proportion of the bond survival time in the simulation time. In all figures, the red and blue curves indicate the linear and cyclic RGD-containing liganded systems, respectively.

Each system was first subjected to energy minimization for two consecutive 50,000 conjugate gradient steps: first with backbone atoms fixed, and second with all atoms free. The energy-minimized structures were then gradually heated from 0 to 310 K and subsequently equilibrated for 10 ns under constant pressure and temperature conditions. The temperature was maintained at 310 K and controlled with Langevin dynamics, and the pressure was maintained at 1 atm using the Langevin piston method.

Constant velocity-steered molecular dynamics simulations (SMD) were performed to accelerate integrin unbinding and started from snapshots sampled after 10 ns of equilibration. We used a time step of 2 fs, a uniform dielectric constant of 1.0, a 12 Å cut-off for non-bonded interactions, and a scaling factor of 1–4 interactions of 1.0. During the SMD simulations, the RGD ligand was pulled away from the binding shallow crevice located between the α_v and β₃ subunits of the integrin α_vβ₃ headpiece. The C or N termini of the βA and β-propeller subunits were fixed as original points. The pulling potential moved with a constant velocity v, and the Cα atom of RGD⁵⁰⁰⁵ was steered as the SMD atom. In constant-velocity simulations, the time dependence of the external force is F = k(vt − Δx), where Δx is the displacement along the pulling direction at time t = 0, v is the pulling velocity, and k is the spring constant. The rupture force presents the maximum of the force spectrum. Here, we set the pulling speed as v = 0.02 Å⋅ps⁻¹, 0.04 Å⋅ps⁻¹, and 0.06 Å⋅ps⁻¹ and the spring constant as k = 7 kcal⋅mol⁻¹⋅Å⁻². The SMD simulation was performed in an NVE ensemble at 310 K and lasted for at least 400,000 steps.

In the binding process, the ligand was placed at a distance of 25 Å from the target binding gorge first, and then solvated with a TIP3P water box, which was neutralized at a 150 mM ionic concentration, resulting in a system of 87,098 atoms in a 102 × 94 × 92 ų water box. After the system minimized and equilibrated, we chose a pushing speed of 0.1 Å⋅ps⁻¹ and a spring constant of 2 kcal⋅mol⁻¹⋅Å⁻², and next pushed the RGD ligand into integrin α_vβ₃ from the direction of the Asp sidechain oxygen OD2 of the ligand to the MIDAS site.

3. Results

3.1. Configurational Changes for RGD-Integrin α_vβ₃ Interaction during Equilibration

To reveal the structural stability of different structures of RGD-containing liganded integrins, we first analyzed the root mean square deviations of different domains of integrin α_vβ₃ systems during equilibration (Figure 2). RMSD represents the root mean square deviation of Cα atoms. The time–RMSD profiles of the global structure were smooth and then leveled off at ~2 Å, indicating that equilibrium had been reached (Figure 2A). The time–RMSD profiles of the global were smoother and remained at a similar level in both cyclic and linear RGD systems, whereas small fluctuations of the β-propeller and βA domain were observed in the linear system (Figure 2B,C). In particular, the trends of RMSD at the RGD-containing segment were obviously different (Figure 2D), suggesting that the structural difference is commonly located in the RGD region. Large-amplitude fluctuations may play a vital role in inhibiting the ligand from stable binding to the active-site site. To make clear the reason for the large-amplitude fluctuations, time series of the Cβ atom distance between Arg^RGD and Val^RGD were calculated (Figure 2E). Obviously, the fluctuations of the distance in the linear system were large, ranging within 6~11 Å, and the fluctuation trend was consistent with the trends of RMSD at the RGD-containing segment along the simulation time. Meanwhile, from the trajectory of the linear RGD ligand as time changed from dark to light grey, we observed that the distance at the ends of the linear RGD residues increased gradually and a little bit of rotation happened in the side chain of the benzene ring (Figure 2F). Thus, the linear RGD form induced local destabilization of the system due to its flexibility, whereas according to the rigidity of the ring structure, the cyclic peptides were more stable.

Integrin can form strong noncovalent bonds with RGD-containing peptides that bind to a shallow crevice rather than into a deep binding pocket, which is not well shielded from attacks by water molecules. To explore the influence of such differences on the binding surface, we tested the SASA value of the buried surface. It was found that the buried SASA value in the linear RGD complex was approximately 420 Ų, whereas the buried SASA value in the cyclic RGD system was markedly decreased (~340 Ų) (Figure 3A). It was speculated that due to the flexibility of the linear RGD, the steric conformation became more relaxed, and the contact area was larger, thereby blocking the access of free water molecules to the most critical RGD-α_vβ₃ interaction surface.

Simultaneously, the major interactions between integrin and RGD peptides in the equilibrium processes were investigated. It was found that the salt bridge between Arg^RGD and Asp²¹⁸ remained intact throughout the cyclic RGD simulation, whereas the occupancies were significantly decreased in the linear system (Figure 3B,G). Hydrogen bond counting showed that compared with the cyclic RGD peptide, several hydrogen bonds (Asp²¹⁸-Arg^RGD, Ser¹²³-Asp^RGD, Arg²¹⁶-Asp^RGD) had significantly higher occupancies in the linear RGD system (Table 1). Quantitatively, the number of hydrogen bonds between integrin and RGD peptides was approximately 6 and 2–4 in the linear/cyclic RGD liganded system, respectively (Figure 3C). Moreover, the linear RGD peptides could easily form intramolecular hydrogen bonds with themselves, which could not be formed in the cyclic RGD sequences (Figure 3D). Thus, it was speculated that the form of intramolecular and intermolecular hydrogen bonds may block the binding of linear RGD peptides to integrin.

The occupancies of Ser¹²³-Asp^RGD were 0.21% and 66.81% in the cyclic RGD and linear RGD systems, respectively (Table 1). We found that one carboxylate oxygen of Ser¹²³ came into contact with the MIDAS Mg²⁺ ion, and another side chain was directly coordinated with ADMIDAS in the cyclic RGD, whereas the Ser¹²³ side-chain group formed two salt bridges with MIDAS ion and Asp^RGD in the linear RGD (Figure 3E,F). In the linear RGD, Ser¹²³ bound to the Asp of RGD, except in the cyclic RGD, indicating that the cyclic RGD was far away from Ser¹²³, thus presenting a smaller binding surface.

We further analyzed the dynamic trajectories during equilibration to estimate whether the different conformational features of the ligand caused a local structural change in integrin α_vβ₃. Vogel et al. reported that ligand binding could induce the activating integrin α_vβ₃ conformational change via the formation of the T-junction between the middle of the α1 helix and the top of the α7 helix [13]. Similarly, we measured the distance between Leu¹³⁴ of the α1 helix and Leu³³³ on the β6 strand to identify the formation of T-junctions. This distance is shown in Figure 4A. The comparison between the linear and cyclic RGD revealed that there was a considerable decrease in distance in the linear RGD, whereas the cyclic RGD remained nearly identical throughout the simulation time. Corresponding to a decrease in the C_β-atom distance, the bidirectional allosteric signal process was accompanied by an increase in distance between the β1 and α1 and β6–α7 loops. A similar increase in the separation between the β1–α1 and β6–α7 loops was observed in both systems. However, the increase in distance in the linear region was obvious (Figure 4B).

These findings were further illustrated (Figure 4C) by the alignment of both the linear and cyclic RGD and the corresponding domains from the unliganded α_vβ₃ integrin after 10 ns of equilibration. Restoration of the β1–α1 loop was found in both the linear and cyclic liganded structures. In addition, another characteristic helical structural alteration was the inward movement of the α1 helix. Moreover, the modifications were more apparent in the linear RGD-occupied integrin α_vβ₃ due to the firm bond of Ser¹²³-Asp^RGD (occupancy 66.81% vs. 0.21%) than in the cyclic RGD. The inward shift of the α1 helix has been reported to promote βA/hybrid domain hinge opening along the allosteric pathway [12,29,30]. Whether it is easier for the linear RGD-containing peptides to trigger the allosteric pathway of integrin activation requires further verification combined with other domains, such as the hybrid domain and transmembrane domains.

3.2. Interaction Energy of RGD Peptide Binding to Integrin α_vβ₃

From the crystal structure of the cyclic pentapeptide in integrin α_vβ₃, we observed that Asp^RGD coordinated with a metal ion located at the MIDAS; Arg^RGD formed two salt bridges with Asp²¹⁸ and Asp¹⁵⁰ in the α_v subunit. This phenomenon indicates that cooperative interaction between RGD and integrin α_vβ₃ occurs not only through the neighboring polar or charged amino acids, but also through the Mg²⁺ metal ion. To elucidate the effect of different integrin domains interacting with RGD, we analyzed the electrostatic energies and van der Waals interactions between RGD and integrin based on a series of molecular dynamics simulations, and the interaction energies and the key interactions are shown in Figure 5. Our results show that the kinetic energy barrier between the metal ion Mg²⁺ and cyclic RGD was approximately −420 kcal/mol, which is higher than that observed for the linear RGD system (−350 kcal/mol) (Figure 5A). Among the three Mg²⁺ ions at MIDAS, ADMIDAS, and LIMBS, compared to the linear RGD, the interaction energies of the cyclic RGD residues interacting with the Mg²⁺ at MIDAS appeared to be much higher (Figure 5C). The fluctuation of interaction energies was consistent with the distance of the mass center of two oxygen atoms of Arg^RGD and MIDAS ions (Figure 5D). The most significant change was that only one of the two carboxylic oxygen atoms from Asp^RGD formed contact with the MIDAS ion during the equilibration of the linear system (Figure 5F), whereas both carboxylic oxygen atoms of Asp^RGD remained in contact with the MIDAS ion in the cyclic RGDfV (Figure 5E). Therefore, even though the total interaction energy of RGD-α_vβ₃ was at a similar level in both systems, the larger energy fluctuation observed in linear RGD may have been the cause of structural instability in the linear system (Figure 5B).

3.3. Force Induced Unbinding of Liganded Integrin α_vβ₃

To map the unbinding dynamics of different RGD structures on the interaction mechanism with integrin α_vβ₃ at the atomic level, external forces were applied to the Cα atom of the RGD⁵⁰⁰⁵ of the ligand to facilitate its binding with integrin α_vβ₃. The pulling speed v = 0.02 Å⋅ps⁻¹ was set based on others’ publication [31]. The rupture force profiles in the dissociative trajectories are shown in Figure 6A. At the same time, we also set up two other speeds, v = 0.04 Å⋅ps⁻¹ (Figure 6B) and v = 0.06 Å⋅ps⁻¹ (Figure 6C). The rupture force of pulling cyclic RGDfV out of the binding site was much larger than that of the linear RGD pentapeptide and took more time to dissociate (Figure 6D,E). The difference in the rupture force profiles implied that it was easier for the linear RGDfV ligand to unbind to the active site and much quicker to dissociate than the cyclic RGDfV ligand. The results show good agreement with the experimental evidence that linear RGD-containing peptides have a short circulation half-life and are more susceptible to chemical degradation [20,32,33,34].

3.4. Characterization of Binding Pathway

After the unbinding simulation of the integrin α_vβ₃-RGD complex, we simulated the binding process between integrin α_vβ₃ with linear RGDfV and cyclic RGDfV. Steered MD simulations were performed to enforce the binding of the RGD-containing ligand to the binding site, as they were unable to effectively and suitably ligate in free MD simulations. For this purpose, the pushing direction by loading force was chosen along the vector pointing from the side-chain oxygen atom OD2 of Asp^RGD to the MIDAS ion, the key RGD-integrin-binding contact. The force-induced binding test was applied at a constant pulling velocity of 0.1 Å⋅ps⁻¹, and the spring constant was set to 2 kcal/(mol·Å²). Ten separate SMD simulations were conducted for each system. However, only two of these trajectories entered the active site within a short simulation time, whereas eight of these simulations successfully recovered the bound pose in the binding site in the cyclic RGDfV.

The interaction features were examined to assess the behavior of the binding process (Figure 7). Compared with the cyclic RGDfV, the distance between the side-chain oxygen atom OD2 of Asp^RGD and the headpiece of integrin α_vβ₃, which was the most vital bond in the complex, was larger in the linear system (Figure 7A), implying that the interaction was much weaker for the linear system. The formation of hydrogen bonds occurred continuously while moving along the binding site (Figure 7B). The final number of hydrogen bonds was greater than that of cyclic RGD. The more hydrogen bonds, the stronger the binding forces. The linear RGD peptide may have arrived at its binding site slowly. The interaction energy profile during the binding process is shown in Figure 7C. The energies increased along both RGD peptides, pushing into the active binding site. The lowest valley was produced in the cyclic RGD system. The corresponding time was also much shorter. The results indicate that it was easier for the cyclic RGD-containing peptide to bind to the active site than the linear RGD peptides.

4. Discussion

This study attempted to explore the effect of the RGD structure (cyclic or linear) on integrin α_vβ₃. Our simulations provided insight into the mechanism to explain why cyclic RGD improves tumor-targeting efficacy at the atomic level. Since the amino acids flanking the RGD sequence could alter the mechanical stability and other properties, cyclic RGDfV and linear RGDfV pentapeptides with the same residues interacting with integrin α_vβ₃ were set up for comparison. Our data indicate that the linear RGD peptides have a higher flexibility (Figure 2D), and larger amplitude fluctuations may play a negative role in the binding or unbinding process of the α_vβ₃-RGD interaction. Because of the flexible structure of the linear RGD, except for the linear RGD that can form a combination that was reported in [32], it was easier to form hydrogen bonds with integrin α_vβ₃ than with cyclic RGD (Figure 3C). However, the linear RGD could interact with many residues at the ligand-binding site. Therefore, its ability to accurately locate the critical site was significantly reduced. Cyclic RGD had strong localization ability, and both the binding energy and binding efficiency with Mg²⁺ were outstanding (Figure 5A and Figure 7A). Thus, when considering the synthesis of RGD-targeting drugs, cyclic RGD peptides should be considered for their stability and effectiveness, and it is true that many drugs have been designed in recent years based on cyclic RGD, including Cilengitide [35,36].

Most interestingly, we found that the Ser¹²³-Asp^RGD interaction could improve the formation of α1 the α7 T-junction. The formation of T-junctions may promote the βA/hybrid domain hinge opening in the integrin activation process. The residue Ser¹²³, located at the top of α1 helix and close to the ligand-binding site, was conducive to the inward and upward movement of the α1 helix interacting with Asp^RGD. The Ser¹²³-Asp^RGD bond may play a critical role in the allosteric pathway. The findings of our study reveal that in linear RGD, the occupancy of the Ser¹²³-Asp^RGD interaction is larger than that of cyclic RGD (Table 1). Considering the high sensitivity to minor structural perturbations along the allosteric pathway, allosteric activation requires further verification.

Previously, Yu et al. reported that the targeted recognition between RGD peptide and integrin α_vβ₃ was mainly driven by electrostatic interactions [32]. Indeed, we were able to confirm that the interaction energy of the complex mainly involved electrostatic interactions between the residues in RGD and the MIDAS ion in integrin α_vβ₃. However, van der Waals forces are also non-bonding forces; therefore, they are included in the energy values. Yu et al. mentioned that the interaction between RGD and ions at the ion-binding site is the mechanism of the interaction between RGD and the ion-binding site, which is consistent with the interaction between RGD and Mg²⁺ mentioned in our study. We also showed that the interaction between ions and RGD was stronger in cyclic RGD than that in linear RGD. The evidence comes from the energy level (Figure 5A) and the accuracy and velocity of interaction (Figure 7A,B). The explanation of the distinct differences is whether the carboxylic oxygen of Asp^RGD remains in contact with the MIDAS ion. It was shown that the oxygen atoms of Asp^RGD contributed significantly to the electrostatic effect on the MIDAS cationic ligands. Furthermore, the cyclic RGD resulted in a much greater decrease in the buried SASA than in the linear RGD ligand (Figure 3A). It was possible to stabilize the interaction by shielding the critical Asp^RGD-MIDAS interaction from attacks by free water molecules in the cyclic RGD. This would enhance the mechanical stability of protein–protein interaction [17]. Our simulations provide a valid explanation for the experimental result that cyclic peptides are more stable at the atomic level.

In this study, the binding and unbinding properties of cyclic RGDfV and linear RGDfV interacting with integrin α_vβ₃ were obtained. The rupture force of the cyclic RGD RGDfV-α_vβ₃ complex to dissociate was much larger than that of the linear RGD pentapeptide, and the corresponding time at the peak force was also greater than the time required to pull the linear RGDfV out of the active binding site (Figure 6). The difference in forced dissociation tests suggests that the barrier along the unbinding pathway in the cyclic RGD was much larger than that in the linear RGD system. The results indicate that the binding affinity between cyclic RGDfV and integrin α_vβ₃ was higher than that of linear RGDfV peptides. This is consistent with the experimental evidence that cyclic RGD is commonly employed to improve the binding properties of RGD peptides because of the rigidity of the ring structure [10,37]. For the RGD pentapeptide entry process, the linear RGD-containing peptides easily formed hydrogen bonds because of their flexibility, and hydrogen bonds negatively impacted the ability of linear peptides to bind integrin α_vβ₃ accurately. Thus, they took much more time to bind to the active site than cyclic RGD (Figure 7). Additionally, the interaction energy in binding integrin α_vβ₃ exhibited greater potential energy in cyclic peptides in a short simulation time. This distinct difference implies that it is easier for cyclic peptides to bind to the active site than to linear peptides. This could be the reason that most cyclic RGD-containing peptides have a longer circulation survival life, which results in ideal effects on treatment [38].

Recently, researchers have not only found that cyclic RGDs have important pharmaceutical value, but that polycyclic RGDs are also gaining increasing attention. For example, tetrameric cyclo (DKP-RGD) ligands have the potential to improve tumor targeting for diagnosis and therapy [39]. Therefore, there is still much room for research on the medicinal value of RGD.

5. Conclusions

Owing to the application value of integrin–RGD interaction in the design of oncology drugs, an increasing number of researchers are focusing on the mechanism of interaction between integrins and RGD. In our study, we found that the linear RGD peptides triggered the instability of the configurational changes due to their flexibility, and the linear RGD peptides formed their own combination more easily than cyclic RGD. The main interaction energy between Mg²⁺ and cyclic RGD was much stronger than that of the linear RGD system due to the well shield lessening attacks by free water molecules. In addition, our results show that it was easier for linear RGDfV to leave the active site and more quickly than the cyclic RGDfV ligand, whereas it was more difficult to enter the appropriate active binding site in linear RGD. The Ser¹²³-Asp^RGD bond may play a critical role in the allosteric pathway. In conclusion, our study may provide new ideas for new drug creation and development in the treatment of tumors.

Author Contributions

Q.L. and N.L. were responsible for the overall design and investigation. N.L. and S.Q. performed data analysis. N.L. S.Q., Y.F., J.W. and Q.L. were responsible for manuscript writing and participated in the discussion of the results. All authors discussed the basic structure of the manuscript. All authors have read and agreed to the published version of the manuscript.

Funding

This work was supported by the National Natural Science Foundation of China (No. 31870928 (Q.L.) and No. 12072117(J.W.)) and the Natural Science Foundation of Guangdong Province, China (No. 2021A1515010040).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

No new data were created or analyzed in this study.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Table

Figure 1. System setup. (A) A schematic of the extracellular domain of integrin αvβ3. (B) Secondary structure of the integrin headpiece fragments in a NewCartoon representation. β–propeller and βA domains are indicated in cyan and orange, respectively. RGD peptide ligand is shown in red VDWs; the T-junction structure is represented in a green NewCartoon, with β6, α7, β1, and α1 from left to right; the blue sphere is an alpha carbon atom of residues. (C) The RGD ligated with integrin αvβ3 headpiece in a water box used for equilibration. (D) The cyclic RGD peptide is shown in ResID. (E) The linear RGD peptide is shown in ResID. The cyclo-RGDfV and linear RGDfV peptides have the same meaning as the cyclic and linear RGD peptides in this study, respectively.

Figure 2. The simulation time curves of RMSD for RGD liganded integrin αvβ3 during equilibration for the global systems (A), the β-propeller regions (residues 1–438) (B), the βA regions (residues 110–353) (C), and the RGD peptide segments (D). Each time–RMSD profiles was averaged from three independent runs. The unliganded αvβ3 integrin regions are shown in black, the cyclic RGD and its complex with integrin αvβ3 are shown in blue, and the linear RGD and its complex with integrin αvβ3 are shown in red. (E) Time series of the Cβ-atom distance between ArgRGD and ValRGD. (F) The trajectory of the linear RGD ligand as time changed from dark to light grey. Image of all 50 frames shown at once, smoothed with a 100-frame window.

Figure 3. The buried solvent accessible surface area and the major interactions during equilibration. Each curve presented an average of data from three independent runs. (A) The buried solvent accessible surface area (SASA) during equilibration. (B) The distance between Asp218 and ArgRGD, (C) Time series of the number of hydrogen bonds for receptor–ligand interaction. (D) Time series of the number of hydrogen bonds for the RGD interactions. (E,F) Coordination spheres of the three important metal ions in the (E) cyclic RGD liganded system and (F) linear RGD liganded system. Atoms are shown with a CPK representation (carbon atoms are cyan, oxygen atoms are red, and hydrogen atoms are white), and a portion of integrin is shown as a yellow ribbon. The ADMIDAS (purple), MIDAS (green), and LIMBS (yellow) ions are shown in a VDW representation.

Figure 4. (A) Time series of the Cβ-atom distance between Leu134 of the α1 helix and Leu333 on the β6 strand. (B) Time series of the distance between the centers of mass of the β1-α1 loop and the β6-α7 loop. (C, top and bottom) The superposition of the equilibrated complex and the unliganded structure (silver). Each conformation was an average of the last 10 ps simulation times with three runs during equilibrating. (C, middle) Both of the equilibrated complexes were closely aligned with the conformation of the unliganded structure (silver). The selected regions are shown in a NewCartoon representation. Red arrows identify the region where the β1-α1 loop met the β6-α7 loop, and the region where Leu134 of the α1 helix contacted with Leu333 of the β6 strand and α7 helix during T-junction formation.

Figure 5. The interaction energies of RGD peptides targeted with integrin αvβ3. Each curve presented an average of data from three independent runs. (A) Interaction energies of three Mg2+ ions contacting RGD peptides. (B) Interaction energies of the RGD peptide domain and integrin αvβ3. (C) Interaction energies of MIDAS contacting RGD peptides. The term of interaction energy included electrostatic and van der Waals. (D) Time series of the distance between the mass center of two oxygen atoms of ArgRGD and MIDAS ions. (E) The binding sites of the cyclic RGD pentapeptide and the MIDAS ion. (F) Binding site of the linear RGD pentapeptide and the MIDAS ion.

Figure 6. Force-induced dissociation of linear vs. cyclic RGD from integrin αvβ3 by steered molecular dynamics. The force profiles were obtained from five independent simulations of pulling the ligand away from integrin αvβ3 by pulling speed v = 0.02 Å⋅ps−1 (A), v = 0.04 Å⋅ps−1 (B), and v = 0.06 Å⋅ps−1 (C). Average peak force (D) and average time required to reach peak force (E) are indicated. Error bars: standard deviation from the mean; asterisks: statistical significance following the t-test (* p < 0.05; *** p < 0.01).

Figure 7. Binding pathway simulation. The distance from the oxygen atom OD2 of AspRGD to the MIDAS ion (A), numbers of hydrogen bonds (B), and interaction energies between the RGD peptide and integrin αvβ3 (C) vs. time in the process of pushing RGD-containing peptides into the active binding site. Each curve presents an average of data from two independent runs. The term “interaction energies” refers to the total number of electrostatic and VDW interactions.

Table 1

The occupancies of interactions at the RGD peptide/integrin interfaces of both simulation systems.

	Residues	Cyclic	Linear	Residues	Cyclic	Linear
Hydrogen Bond	αv-TYR178&Arg5001	6.05	0.07	β3-TYR122&ASP5003	6.75	2.40
	αv-ALA213&Arg5001	7.40	3.87	β3-SER123&ASP5003	0.21	66.81
	αv-GLN214&Arg5001	0.20	17.57	β3-ASN215&ASP5003	11.20	1.99
	αv-ASP218&Arg5001	97.86	83.25	β3-ARG216&ASP5003	19.55	42.13
	αv-ASP218&GLY5002	0.00	16.73	β3-LYS125&VAL5005	0.00	8.95
	αv-LYS253&GLY5002	0.00	9.00	β3-ARG214&VAL5005	0.00	10.56
Salt Bridge	αv-ASP218&ARG5001	93.10	49.07	αv-ASP150&ARG5001	0.20	0.00