Nothing Special   »   [go: up one dir, main page]
Next Article in Journal
Antimalarial Activity of Aqueous Extracts of Nasturtium (Tropaeolum majus L.) and Benzyl Isothiocyanate
Previous Article in Journal
Screening Antioxidant Components in Different Parts of Dandelion Using Online Gradient Pressure Liquid Extraction Coupled with High-Performance Liquid Chromatography Antioxidant Analysis System and Molecular Simulations
You seem to have javascript disabled. Please note that many of the page functionalities won't work as expected without javascript enabled.
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Molecules
Volume 29
Issue 10
10.3390/molecules29102317
molecules-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Molecular Mechanism of Phosphorylation-Mediated Impacts on the Conformation Dynamics of GTP-Bound KRAS Probed by GaMD Trajectory-Based Deep Learning

by
Jianzhong Chen
1,2,*,
Jian Wang
1,
Wanchun Yang
1,
Lu Zhao
1,
Juan Zhao
1 and
Guodong Hu
2,*
1
School of Science, Shandong Jiaotong University, Jinan 250357, China
2
Shandong Key Laboratory of Biophysics, Institute of Biophysics, Dezhou University, Dezhou 253023, China
*
Authors to whom correspondence should be addressed.
Molecules 2024, 29(10), 2317; https://doi.org/10.3390/molecules29102317
Submission received: 23 April 2024 / Revised: 9 May 2024 / Accepted: 12 May 2024 / Published: 15 May 2024
Browse Figures
Figure 1
<p>Molecular structure: (<b>A</b>) the GTP-bound WT KRAS, in which secondary structure and phosphorylation group are labeled, (<b>B</b>) the geometry position of the phosphorylated residues relative to GTP, (<b>C</b>) phosphorylation with two net negative charges, and (<b>D</b>) GTP. In this Figure, KRAS is shown in cartoon modes while GTP and phosphorylation-related residues are shown in ball-stick modes.</p> "> Figure 2
<p>Classification and saliency map of residue contact gradients: (<b>A</b>) classification of the GTP-bound WT, pY32, pY64, and pY137 KRAS, (<b>B</b>–<b>E</b>) the saliency map of residue contact gradients for the GTP-bound WT, pY32, pY64, and pY137 KRAS, and (<b>F</b>) key structural domains revealed by the DL. The gradient of each residue contact is shown in a 0.3 (white) to 0.6 (black) color scale.</p> "> Figure 3
<p>Free energy profiles and representative structures of the GTP-pY32 KRAS: (<b>A</b>) FEL, (<b>B</b>) superimposition of initial optimized structure with the EV1, EV2, and EV3 structures, (<b>C</b>) structural superimposition of GTP and magnesium ions (Mg<sup>2+</sup>) in the initially optimized structure and the EV1, EV2, and EV3 structures, and (<b>D</b>–<b>F</b>) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, and EV3 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.</p> "> Figure 4
<p>Free energy profiles and representative structures of the GTP-pY64 KRAS: (<b>A</b>) FEL, (<b>B</b>) superimposition of initial optimized structure with the EV1, EV2, EV3, and EV4 structures, and (<b>C</b>–<b>F</b>) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, EV3, and EV4 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.</p> "> Figure 5
<p>Free energy profiles and representative structures of the GTP-pY137 KRAS: (<b>A</b>) FEL, (<b>B</b>) alignment of initial optimized structure with the EV1, EV2, EV3, and EV4 structures, and (<b>C</b>–<b>F</b>) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, EV3, and EV4 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.</p> "> Figure 6
<p>Dynamics indexes: (<b>A</b>) the probability distribution of GTP RMSDs, (<b>B</b>) RMSFs of KRAS, (<b>C</b>) the probability distribution of the distances between the mass center of all Cα atoms in helix α2 away from that of helix α3, and (<b>D</b>) the probability distribution of the distances between the mass center of all Cα atoms in helix α3 away from that of helix α4.</p> "> Figure 7
<p>Concerted motions of the GTP-KRAS revealed by the PC1 from the PCA: (<b>A</b>) the GTP-WT KRAS, (<b>B</b>) the GTP-pY32 KRAS, (<b>C</b>) the GTP-pY64 KRAS, and (<b>D</b>) the GTP-pY137 KRAS. KRAS was shown in cartoon mode and GTP is shown in ball-stick modes. The blue, green and red indicate the P-loop, SW1 and SW2, respectively.</p> "> Figure 8
<p>Collective movement of the GTP-KRAS revealed by the PC2 from the PCA: (<b>A</b>) the GTP-WT KRAS, (<b>B</b>) the GTP-pY32 KRAS, (<b>C</b>) the GTP-pY64 KRAS, and (<b>D</b>) the GTP-pY137 KRAS. KRAS was shown in cartoon mode and GTP is shown in ball-stick modes. The blue, green and red indicate the P-loop, SW1 and SW2, respectively.</p> "> Figure 9
<p>The probability of the χ angle for the sidechains of three phosphorylation-related residues: (<b>A</b>) residue 32, (<b>B</b>) residue 64, (<b>C</b>) residue 137, and (<b>D</b>) geometric positions of key residues.</p> "> Figure 10
<p>Key interactions of GTP with KRAS and their corresponding distance distributions. (<b>A</b>) Relative position of GTP in the binding pocket of KRAS, (<b>B</b>) geometry positions of key interactions, (<b>C</b>) the distances of the π-π interaction between GTP and F28, (<b>D</b>) the distances of the salt bridge interaction between the phosphorus atom PB of GTP and the nitrogen atom NZ of K16, (<b>E</b>) the distances of the salt bridge interaction between the phosphorus atom PG of GTP and the nitrogen atom NZ of K16, and (<b>F</b>) the distances of the salt bridge interaction between the guanine group of GTP and the carbonyl group of D119.</p> "> Figure 11
<p>Key interactions of magnesium ions (Mg<sup>2+</sup>) with GTP and KRAS together with their distance distributions: (<b>A</b>) the geometry position of interactions, (<b>B</b>) the distances between Mg<sup>2+</sup> and the oxygen atom OG of S17, (<b>C</b>) the distances of Mg<sup>2+</sup> away from the oxygen atom O of T35, (<b>D</b>) the distances between Mg<sup>2+</sup> and the oxygen atom O of T35, (<b>E</b>) the distances of Mg<sup>2+</sup> away from the mass center of oxygen atoms of OD1 and OD2 in D57, (<b>F</b>) the distances of Mg<sup>2+</sup> away from oxygen atom O2B of GTP, and (<b>F</b>) the distances between Mg<sup>2+</sup> and oxygen atom OG of GTP.</p> "> Figure 12
<p>Workflow of deep leaning from GaMD simulations: (<b>A</b>) the GTP-bound KRAS with three phosphorylation sites, (<b>B</b>) the initialized GTP-bound KRAS, (<b>C</b>) conformational ensembles recorded in three independent GaMD trajectories, (<b>D</b>) images extracted though the MDTraj program, (<b>E</b>) convolution neural networks, (<b>F</b>) the saliency maps calculated through backward propagation, (<b>G</b>) key residue contacts identified by deep learning, and (<b>H</b>) free energy landscapes used for revealing the phosphorylation-mediated effect on free energy profiles of the GTP-bound KRAS.</p> ">
Versions Notes

Abstract

:
The phosphorylation of different sites produces a significant effect on the conformational dynamics of KRAS. Gaussian accelerated molecular dynamics (GaMD) simulations were combined with deep learning (DL) to explore the molecular mechanism of the phosphorylation-mediated effect on conformational dynamics of the GTP-bound KRAS. The DL finds that the switch domains are involved in obvious differences in conformation contacts and suggests that the switch domains play a key role in the function of KRAS. The analyses of free energy landscapes (FELs) reveal that the phosphorylation of pY32, pY64, and pY137 leads to more disordered states of the switch domains than the wild-type (WT) KRAS and induces conformational transformations between the closed and open states. The results from principal component analysis (PCA) indicate that principal motions PC1 and PC2 are responsible for the closed and open states of the phosphorylated KRAS. Interaction networks were analyzed and the results verify that the phosphorylation alters interactions of GTP and magnesium ion Mg2+ with the switch domains. It is concluded that the phosphorylation pY32, pY64, and pY137 tune the activity of KRAS through changing conformational dynamics and interactions of the switch domains. We anticipated that this work could provide theoretical aids for deeply understanding the function of KRAS.

1. Introduction

RAS proteins, a group of small GTPases, play vital roles in intracellular signal transduction pathways and are responsible for cell growth, proliferation, differentiation, and survival [1,2,3,4]. In fact, the RAS proteins primarily include four isoforms in human cells, namely HRAS, NRAS, KRAS4A, and KRAS4B. They function as a molecular switch by cycling between the GTP-bound active state and the GDP-bound inactive one [5,6,7,8]. Recently, increasing attention has been paid to these RAS proteins due to their close relations to a variety of human cancers. Previous studies have indicated that mutant RAS proteins with a constitutively active state trigger tumorigenesis by hyperactivating downstream signaling pathways. Among HRAS, NRAS, and KRAS, KRAS is involved in the most frequently mutated oncogene in human cancers, accounting for approximately 25% of lung, 40% of colorectal, and 95% of pancreatic cancers [9]. Overall, KRAS oncogenes make up about 82% of the Ras nonsynonymous mutations found in human oncogenes [10]. Therefore, KRAS has been given special attention and is of interest in understanding treatment strategies for human cancers.
KRAS has similar structures to HRAS and NRAS, primarily composed of the effector lobe and allosteric lobe. The effector lobe, including residues 1–86, is located at the N-terminal region of the catalytic domain, which forms three common secondary structures, namely the P-loop (9–18), switch domain 1 (SW1: residue 28–40), and switch domain 2 (SW2: residue 59–75), as depicted in Figure 1A. The allosteric lobe, including residues 87–166, consists of an α-helical dimerization interface and nucleotide base binding motifs (Figure 1A). In terms of function, the effector lobe participates in the binding of KRAS to GTPase activating proteins (GAPs) and guanosine exchange factors (GEFs), while the allosteric lobe is responsible for relaying information through protein interactions [11,12,13,14,15,16]. Point mutations and residue modifications can directly affect the activity of KRASB by disturbing the conformational dynamics of the switch domains SW1 and SW2 [17,18,19,20,21,22,23,24,25,26]. Furthermore, point mutations in codons 12, 13, or 61 convert RAS into an active oncoprotein by impeding the GTPase activity of RAS and changing the conformational states of the switch domains [9,27,28,29,30,31,32,33,34,35,36]. Based on the roles of conformational dynamics for SW1 and SW2, clarifying the molecular mechanism of the alterations in conformational behavior of KRAS caused by mutations and residue modifications remains a significant and outstanding aim because of its involvement in a variety of cancers.
Phosphorylation, one of the residue modifications, is a common occurrence in proteins within specific environments, which can impact the functionality of proteins. Recent studies have focused on exploring the effects of phosphorylation on the activity of RAS proteins [37,38]. Kano et al. demonstrated that KRAS phosphorylation through Src on residues Y32 and Y64 changes the conformation of the SW1 and SW2, which stalls multiple steps of the GTPase cycle and disturbs binding of KRAS to effectors [39]. Wang et al. proposed that phosphorylation of Y32 and Y64 influences the catalytic transformations between GTP- and GDP-bound states, as mediated by the intrinsic or GAP/SOS (Son of sevenless) and affects its signal transduction to Raf [2]. Khaled et al. employed all-atom molecular dynamics simulations (MD) and Markov state models to conduct a comparative analysis of the oncogenic G12D KRAS mutant and its phosphorylated variant at Y32. Their results revealed that the G12D KRAS with Y32 phosphorylation exhibits greater flexibility in the SW1 compared to the G12D KRAS [40]. The work of Qiu et al. highlighted the significance of phosphorylation at multiple residues in regulating both wild-type (WT) and mutated RAS proteins, shedding light on the mechanistic insights into associated signaling pathways [41]. Although a deeper understanding of RAS phosphorylation can provide valuable insights for future therapeutic strategies, the current knowledge regarding the phosphorylation-mediated impacts on KRAS activity remains limited. Therefore, further investigation into the influence of phosphorylation on the conformational dynamics of KRAS is essential for the future design of anticancer therapies.
Among the current tools for probing the conformational dynamics of targets, conventional molecular dynamics (cMD) [42,43,44,45,46,47,48], Gaussian accelerated molecular dynamics (GaMD) simulations [49,50,51], and analysis of binding free energy [52,53,54,55,56,57] have been extensively used to reveal the molecular mechanisms and free energy profiles underlying conformational changes. Compared to cMD simulations, GaMD simulations more easily overcome energy barriers in protein systems. Furthermore, GaMD simulations have obtained successes in probing alterations in the conformational dynamics of targets [58,59,60,61,62,63,64]. Despite the great advantages of cMD and GaMD, there is a significant challenge in efficiently extracting dynamic information and functional differences of targets from a large MD trajectory. To address this issue, cMD and machine learning (ML) were integrated to deeply decode the function of G protein-coupled receptors (GPCRs) from the MD trajectory [65,66]. Miao’s group proposed a trajectory-based deep learning (DL) approach, called GaMD, DL, and free energy profiling workflow (GLOW), to successfully decipher the molecular mechanisms underlying the activation and allosteric modulation of G protein-coupled receptors [67,68]. The aforementioned methods provide reliable tools for gaining insights into the phosphorylation-mediated effects on the conformational dynamics of KRAS.
To achieve our goal, we selected the GTP-bound WT KRAS and the GTP-bound phosphorylated KRAS at residues Y32, Y64, and Y137 to investigate the influence of phosphorylation on KRAS conformations. The phosphorylation sites and sidechain changes are shown in Figure 1B,C, along with the structures of GTP displayed in Figure 1D. It is evident that a phosphate group on the sidechain of Tyr can have a greater effect due to its size and negative charge [37]. Phosphorylation of Y32 and Y64 significantly affects the binding of KRAS to its effectors. Phosphorylation of Y137 induces allosteric changes in protein conformation and effector binding [2,39], providing a mechanism for effector-initiated modulation of RAS signaling. However, the molecular mechanism underlying the effects of phosphorylation currently remains unclear, which is why we selected these four systems. In this work, multiple independent GaMD (MI-GaMD) simulations were conducted to enhance the conformational sampling of KRAS and DL was employed to identify the phosphorylation-sensitive sites. Additionally, principal component analysis (PCA) [69,70,71] and calculations of free energy landscapes (FELs) were performed to reveal changes in conformational dynamics and free energy profiles induced by phosphorylation at these three positions. We anticipate that this study will provide valuable insights into the regulation of KRAS activity.

2. Results and Discussion

2.1. Phosphorylation-Mediated Difference in Domain Contacts Revealed by Deep Learning

MD trajectory-based DL has been used to explore molecular mechanisms underlying the activation and allosteric modulation of G protein-coupled receptors [67,68]. Residue contacts were calculated to produce images for the DL. Classification of the GTP-bound WT, pY32, pY64, and pY137 KRAS were achieved by the DL and the results were shown in Figure 2. The overall accuracy and loss reached on the validation set after 25 epochs are 0.9982 and 0.0046, respectively (Supporting Information Figure S1). This high accuracy should be due to the high efficiency of conformational sampling performed by GaMD simulations. From the 8000 snapshots for the validation set of each system, the DL accurately classified most of four current systems (Figure 2A), including 7984 snapshots of the GTP-WT KRAS, 7988 snapshots of the GTP-pY32 KRAS, 7972 snapshots of the GTP-pY64 KRAS, and 7991 snapshots of the GTP-pY137 KRAS individually. As shown in Figure 2A, only a marginal portion of the snapshots were inaccurately classified. In detail, ~0.01% (1 snapshots), ~0.04% (3 snapshots), and ~0.15% (12 snapshots) of the GTP-WT KRAS were, respectively, predicted as the GTP-pY32, pY64, and pY137 KRAS. The ~0.03% (3 snapshots), 0.04% (3 snapshots), and 0.08% (6 snapshots) of the GTP-pY32 KRAS were inaccurately classified as the GTP-WT, pY64, and pY137 KRAS, respectively. The ~0.16% (13 frames), ~0.05% (4 frames), and ~0.13% (11 frames) of the GTP-pY64 KRAS were inaccurately predicted as the GTP-WT, pY32, and pY137 KRAS, respectively, while the ~0.05% (4 frames), 0.03% (2 frames), and ~0.04% (3 frames) of the GTP-pY137 KRAS were erroneously recognized as the GTP-WT, pY32, and pY64 KRAS, separately. These results indicate that the information revealed by our current DL is reliable.
The pixel-attributed residue contact gradient maps of the most populated KRAS structures were depicted at Figure 2B,E. On the whole, the characteristic residue contacts of the GTP-WT KRAS are located between the SW1 and SW2, the P-loop and SW2, the helixes α3 and α4, as well as the helix α3 and SW2 (Figure 2B,F). The characteristic residue contact of the GTP-pY32 KRAS is primarily situated between the SW1 and the P-loop (Figure 2C,F). Compared to the GTP-WT KRAS, the phosphorylation pY32 leads to the disappearance of the characteristic residue contact between the SW1 and SW2 as well as α3 and α4 together with the α3 and SW2. The characteristic residue contacts of the GTP-pY64 KRAS were detected between the P-loop and SW1 as well as the P-loop and SW2 together with the α3 and SW2 (Figure 2D,F). By comparison with the GTP-WT KRAS, the phosphorylation pY64 results in the disappearance of the characteristic residue contact between the SW1 and SW2 as well as the α3 and α4. The characteristic residue contacts of the GTP-pY137 were located between the helixes α3 and α4 (Figure 2E,F). By referencing GTP-WT KRAS, the phosphorylation pY137 leads to the disappearance of the characteristic residue contacts between the SW2 and P-loop as well as the SW1 and SW2 together with the α3 and SW2.
Through the characteristic residue contacts learned by our DL, the phosphorylation-mediated impacts are involved in the structure domains SW1, SW2, α3, and α4 of KRAS. It is worth noting that the SW1 and SW2 participate in binding of KRAS to its effectors or regulators while the α3 and α4 are involved in allosteric regulations on the activity of KRAS [72]. The tyrosyl phosphorylation of KRAS disturbs the GTPase cycle through conformational alteration of the SW1 and SW2 [39], which agrees with our findings learned by the DL. The phosphorylation pY32 highly affects the flexibility of the SW1 and SW2 conformation and tunes the activity of KRAS [40]. The RAS phosphorylation at Y137 allosterically changes protein conformation and effector binding, revealing a mechanism for effector-initiated modulation of the RAS signaling [37]. These previous works support our current findings well.

2.2. Free Energy Profiles Affected by Phosphorylation

The previous DL suggests that the residue contacts between the switch domains are heavily affected by the phosphorylation of three residues: Y32, Y64, and Y137. The distance of the Cα atom of residue 32 in the SW1 away from that of Q61 in the SW2 was selected as an RC. The root-mean-square deviations (RMSDs) of backbone atoms were calculated and the results were provided in Figure S2. It was noted that RMSDs can reflect the total structural fluctuation of KRAS, which is selected as another RC. Thus, these two RCs were used to build FELs to reveal the energy basis of the phosphorylation-mediated impacts on the conformations of KRAS. The FELs and the corresponding representative structures are depicted in Figure 3, Figure 4, Figure 5 and Figure S3–S6.
For the GTP-bound WT KRAS, only an energy valley (EV1) was captured by GaMD simulation and it is situated at the (RMSD, distance) of (2.8, 19.2) Å (Figure S3A), implying that GaMD simulations do not lead to conformational rearrangement of KRAS. To compare the structural difference between the crystal structure and the low energy structure EV1, the initial optimized structure was superimposed with the structure EV1 (Figure S3B). It is observed that only the SW1 produces a small deviation while GTP and magnesium ion (Mg2+) are aligned well (Figure S3C). The distance between Y32 and Q61 in the structure EV1 is only greater than 0.8 Å compared to that in the initial optimized structure (Figure S3D,E). These results suggest that the GTP-bound WT KRAS maintains structural stability through the entire GaMD simulation.
With respect to the GTP-bound pY32 KRAS, three energy valleys (EV1–EV3) were identified by GaMD simulations (Figure 3A). The EV1, EV2, and EV3 are, respectively, situated at the (RMSD, distance) of (2.8, 18.1), (2.6, 22.7), and (2.6, 26.1) Å (Figure 3A and Figure S4). The superimposition of the initially optimized structure with the representative structures EV1–EV3 shows that the SW1 evidently deviates from each other among different structures and has highly disordered states while the SW2 produces a small deviation (Figure 3B). On the contrary, GTP and magnesium ions (Mg2+) are aligned well among different structures (Figure 3C), indicating that they are stable during GaMD simulations. Compared to the GTP-bound WT KRAS, pY32 leads to big structural deviations of SW1 relative to the initially optimized structure and induces more disordered conformational states of the SW1 (Figure 3B). The switch domains of the EV1 structure form a closed state (Figure 3D), those of the EV2 structure show a semi-open one (Figure 3E), and those of the EV3 structure have an open state (Figure 3F). Based on these results, the pY32 results in the conformational transition of the switch domains between the closed and open states. It is well known that SW1 is involved in the binding of KRAS to its effectors. Therefore, pY32 affects the KRAS-effector and GTP-KRAS binding, which regulates the activity of KRAS. In addition, two sodium ions (Na+) appear near the phosphate group of GTP (Figure 3C).
With regard to the GTP-bound pY64 KRAS, four energy valleys EV1–EV4 were captured by GaMD simulations and they are located at the (RMSD, distance) of (3.4, 21.5), (2.7, 19.1), (2.7, 13.2), and (3.1, 11.8) Å, respectively (Figure 4A). The alignment of the initially optimized structure with the structures EV1-EV4 verifies that the switch domains SW1 and SW2 not only produce evident deviations but also possess highly disordered states (Figure 4B). By comparison, the disorder extent of SW2 is higher than that of SW1 (Figure 4B). However, GTP and magnesium ions (Mg2+) agree well with each other among the initially optimized structures and four representative structures (Figure S5A), reflecting that these two ligands are stably maintained in the binding pocket. The switch domains of the structure EV1 shows an open conformation (Figure 4C) while those of the structures EV2–EV4 have a closed conformation (Figure 4D,F); in particular, the structure EV4 verifies that pY64 leads to the conformational transition between the closed and open states of the switch domains. Among four representative structures, the switch domains of the EV4 are the most compact (Figure S5). Compared to the GTP-bound WT KRAS, pY64 induces more energy states and more disordered switch domains, which exerts significant influences on the KRAS-effector binding and the activity of KRAS.
In the case of the GTP-bound pY137 KRAS, four energy valleys EV1–EV4 were detected by GaMD simulations and they are individually trapped at the (RMSD, distance) of (2.8, 18.9), (2.6, 23.7), (3.5, 29.6), and (4.0, 33.1) Å (Figure 5A and Figure S6B–E). The alignment of four representative structures EV1-EV4 with the initially optimized structure demonstrates that the switch domains SW1 and SW2 obviously deviate from each other and are located in a highly disordered state; moreover, the disorder extent of the SW1 is much higher than that of SW2 (Figure 5B). Different from the switch domains, GTP and magnesium ions (Mg2+) are aligned well among the structures EV1–EV4 and the initial structure (Figure S6A), implying that pY137 hardly affects the stability of GTP and Mg2+ through the entire GaMD simulations. The switch domains of the structure EV1 show a slightly closed state (Figure 5C) while those of the structures EV2-EV4 are located in an open state (Figure 5D,F). In particular, the switch domains of the structure EV4 are completely separated and do not form any overlap (Figure 5F). Thus, pY137 induces the conformational transition between the closed and open states of the switch domains. By comparison with the GTP-bound WT KRAS, pY137 not only induces more energy states but also leads to more disordered states of the switch domains. Thus, pY137 certainly impacts the binding of KRAS to its effector and highly disturbs the activity of KRAS.
Based on our current results, the phosphorylation of pY32, pY64, and pY137 yield significant influences on conformations and functions of KRAS: (1) the phosphorylation induces more energy states and changes free energy profile of KRAS, (2) the phosphorylation leads to more disordered states of the switch domains and induces conformational transition between the closed and open states of the switch domains, and (3) the phosphorylation hardly affects the structural stability of GTP and magnesium ions during GaMD simulation. GaMD simulations capture the conformational transition between the open and closed states of the switch domains, which shows the high efficiency of GaMD conformation sampling. The switch domains participate binding of KRAS to its effectors or regulators. Therefore, the highly disordered states of the switch domains caused by the phosphorylation necessarily generate significant impacts on the activity of KRAS. Differently, pY32 and pY137 regulate the activity of KRAS by primarily changing the conformation of the SW1 while pY64 tunes the activity by mainly altering the conformation of the SW2. The work of Kano et al. verified that KRAS phosphorylation via Src on Y32 and Y64 alters the conformation of the SW1 and SW2 and impairs binding to effectors, supporting our current findings [39]. The study by Bunda et al. indicated that the phosphorylation of pY32 alters the positioning or conformation of the SW1 so that it no longer accommodates Raf binding [38], agreeing with our current results.

2.3. Dynamics Behavior of KRAS Influenced by Phosphorylation

To understand the influences of the phosphorylation on structural stability of GTP, the RMSDs of heavy atoms for GTP were calculated relative to the initial structure (Figure S7). It is observed that GTP shows a stable fluctuation in four current systems. The RMSDs of GTP in the WT, pY32, pY64, and pY137 KRAS are populated at 1.5, 0.8, 1.0, and 0.9 Å (Figure 6A), respectively, indicating that pY32, pY64, and pY137 increase the structural stability of GTP in the binding pocket of KRAS. Root-mean-square fluctuations (RMSFs) of KRAS were estimated by using the coordinates of the Cα atoms (Figure 6B). The phosphorylation of pY32 and pY137 weakens the structural flexibility of the SW1 relative to the WT KRAS while pY64 slightly increases that of the SW1. The phosphorylation of pY64 and pY137 strengthens the structural flexibility of the SW2 compared to the WT KRAS but pY32 slightly reduces that of the SW2 (Figure 6B). By comparison with the WT KRAS, the phosphorylation pY64 increases the structural flexibility of the loop L3 linking the α3 and β5 while pY137 weakens that of this loop (Figure 1A and Figure 6B). Different from pY32 and pY64, pY137 obviously enhances the structural flexibility of the α4 (Figure 1A and Figure 6B). In fact, a common mechanistic basis inherent in the high flexibility of the switch regions relates to the conformational transition between two states, which has been revealed by the previous study [2,73].
To probe the effect of the phosphorylation on allosteric regulation, the distances between the mass centers of all Cα atoms in the helix α2 and α3 were computed (Figure S8A) and the results show a stable fluctuation. The distances between α2 and α3 are primarily distributed at 14.9, 15.9, 15.7, and 15.7 Å in the GTP-bound WT, pY32, pY64, and pY137 KRAS, respectively, (Figure 6C), verifying that pY32, pY64, and pY137 increase the distance between α2 and α3. In addition, the distances between the mass centers of all Cα atoms in the helix α3 and α4 were also calculated (Figure S8B) and the results indicate that this distance fluctuates at a range of 11.8–19.8 Å during GaMD simulations of four systems. The distances between α3 and α4 are distributed at 13.9 Å in the GTP-bound WT, pY32, and pY64 KRAS while this distance is populated at 14.5 Å with a wider distribution range in the GTP-bound pY137 KRAS (Figure 6D), suggesting that pY137 increases the distance between α3 and α4. Two helixes α3 and α4 form an interface of allosteric regulation; thus, the changes in distances of the α3 away from α2 and α4 induced by the phosphorylation affect the allosteric activity of KRAS.
To explore phosphorylation-mediated impacts on the dynamics behavior of KRAS, PCA was performed on the coordinates of the Cα atoms recorded at the SGT. The function of eigenvalues over eigenvectors is depicted in Figure S9, which is used to characterize structural fluctuation along the eigenvectors. The first eigenvalues of the pY32, pY64, and pY137 KRAS are obviously increased relative to the WT KRAS, hence the phosphorylation of pY32, pY64, and pY137 enhances the structural fluctuation of KRAS along the first eigenvector. It is also observed that the first two eigenvalues of the four systems are much higher than the other eigenvalues (Figure S8). To check the roles of the first two PCs (PC1 and PC2) in conformational dynamics, they were visualized by means of the VMD1.9.3 software [74] and the results corresponding to PC1 and PC2 were separately depicted in Figure 7 and Figure 8. The switch domains display well-concerted motions in the first two PCs; furthermore, the phosphorylation greatly affects the dynamics behavior of the switch domains. In the PC1 of the GTP-bound WT KRAS, the SW1 and SW2 have an opposite motion tendency and tend to be close to each other (Figure 7A). In the PC2 of the GTP-bound WT KRAS, the SW1 also has a motion tendency to be close to the SW2 (Figure 8A). This information agrees with the previous results that no conformation transition occurs between the closed and open states of the switch domains in the GTP-WT KRAS. In the PC1 of the GTP-pY32 KRAS, the switch domains have a tendency to open (Figure 7B) but it yields a closed tendency in the PC2 (Figure 8B). Therefore, in the GTP-pY32 KRAS, the PC1 is responsible for the open state of the switch domains while the PC2 mediates the closed state. Different from the GTP-pY32 KRAS, the switch domains of the GTP-pY64 and pY137 KRAS tend to be closed in the PC1 (Figure 7C,D) but they have an opening tendency in the PC2 (Figure 8C,D). It is concluded that the PC1 drives the closed conformation of the switch domains in the GTP-pY64 and pY137 KRAS while the PC2 is responsible for the open conformation of the switch domains. These results are in good agreement with the phosphorylation-mediated conformation transition between the open and closed states of the switch domains revealed by the previous analysis of FELs.
In summary, pY32, pY64, and pY137 highly influence the dynamics behavior and activity of KRAS: (1) the phosphorylation stabilizes the structure of GTP in the binding pocket of KRAS and affects the GTP-KRAS interactions, (2) the phosphorylation changes the relative geometric position of the helix α3 to α2 and α4, which affects the allosteric regulation of KRAS, and (3) the phosphorylation alters concerted motions of the switch domains and two principal components PC1 and PC2 are responsible for conformational transition between the open and closed states of the switch domains, which disturbs binding of KRAS to effectors or regulators. Our current findings are in basic agreement with the previous studies [2,40,75].

2.4. Dihedral Angle of Phosphorylated Residues

The phosphorylation not only enlarges the sidechain of the phosphorylated residues but also brings two net negative charges for the sidechain. To clarify the molecular mechanism underlying the phosphorylation-mediated effect, the atoms N, C, CB, and CG of three residues 32, 64, and 137 were utilized to calculate the dihedral angles through the CPPTRAJ program in Amber. The probability distributions of these dihedral angles were plotted in Figure 9A–C and key residues near the phosphorylation sites were depicted in Figure 9D. Indeed, the phosphorylation changes the rotation and distribution of the sidechains of the phosphorylated residues.
The sidechain of residue 32 partly rotates from −166.9° to 96.1° to −157.1° due to pY32 and pY64 compared to the WT KRAS while that of residue 32 partly rotates from −166.9° to 94.7° because of pY137 (Figure 9A); thus, the phosphorylation changes the conformation distribution of the sidechain for residue 32. Structurally, Y32 is located near D33 in the WT KRAS (Figure 9D), hence the negative charges of pY32 yielding strong electrostatic repulsive interactions with the carbonyl group of D33. The negative charges of pY64 and pY137 also form long-range electrostatic repulsive interactions with the carbonyl group of D33. These electrostatic repulsive interactions drive the conformational arrangement of the sidechain for residue 32. The phosphorylation of Y32 leads to an electrostatic repulsion against the negatively charged D38 and D57 from the nucleotide-binding groove, which changes the SW1 positioning or conformation so that it no longer accommodates Raf binding [38].
The dihedral angles of the sidechain for residue 64 are mainly distributed at ~−157.4°, −39.5°, and 85.7° in four systems (Figure 9B). Compared to the WT KRAS, pY32 leads to a decrease in the distribution of the sidechain for residue 64 in −39.5° and 85.7° and an increase in the distribution in −157.4°. The phosphorylation of pY64 and pY137 leads to partial rotations from −157.4° to 85.7° by comparison with the WT KRAS (Figure 9B). Structurally, the sidechain of Y64 is next to the guanidino group of R68 with a net positive charge in the WT KRAS (Figure 9D). The negative charges of pY64 produce strong electrostatic attraction interactions with the guanidino group of R68 in the SW2, which drives the rotation of the sidechain of pY64 (Figure 9B). The negative charges of pY32 and pY137 also generate long-range electrostatic attraction interactions with the guanidino group of R68, which induces partial rotation of the sidechain for Y64. These electrostatic attraction interactions are responsible for the conformational arrangement of the sidechain for residue 64. The previous work indicated that the conformational arrangement of residue 64 caused by the phosphorylation possibly affects the association of RAS and RAN mediated by R68 and R76 with GEFs [76], which agrees with our current result.
The dihedral angles of the sidechain for Y137 are populated at −54.5° in the GTP-bound WT, pY32, and pY137 KRAS while the dihedral angle of pY137 is situated at two peaks, −145.7° and −64.5°, in the GTP-bound pY137 KRAS (Figure 9C). This result indicates that pY137 results in a rotation of the residue 137 sidechain from −64.5° to −147.5° compared to the WT KRAS. In the geometry position, Y137 in the helix α4 is located near R97 in the α3 in the WT KRAS. Thus, the negative charges of the sidechain in pY137 yield strong electrostatic attraction interactions with the guanidino group of R97 (Figure 9D), which makes the sidechain of pY137 rotate toward the gunnidio group of R97. This rotation of pY137 forms a steric hindrance effect with the α3, which increases the distance between the α3 and α4. The phosphorylation of Y137 triggers allosteric changes in GTP-bound RAS that enhance its affinity for RAF1-RBD and severely impair or even abolish the ability of RAS to hydrolyze GTP [37].
In summary, the phosphorylation of tyrosyl not only enlarges the sidechain of residue but also brings two net negative charges. This residue modification easily forms a steric hindrance effect with its neighboring secondary structure and also produces electrostatic attraction or repulsive interactions with the charged residues, which drives the conformational rearrangement of the phosphorylated residues. Thus, electrostatic interactions and steric hindrance effects caused by the phosphorylation of tyrosyl affect the binding of KRAS to effectors and allosteric changes in GTP-bound KRAS [37].

2.5. Interaction Networks Affected by Phosphorylation

To investigate the phosphorylation-mediated effect on the interaction networks of GTP with KRAS, the protein–ligand interaction profiler (PLIP) web server [77,78] was employed to identify GTP-KRAS interaction modes (Figure 10A,B and Figure S10). Moreover, the CPPTRAJ program was adopted to recognize the hydrogen bonding interactions (HBIs) of GTP with KRAS (Table 1). The distances involved in the salt bridge, π-π, and electrostatic interactions were computed and the corresponding results are illustrated in Figure 10B,E. Moreover, the distances of Mg2+-mediated electrostatic interactions with GTP and residues were also estimated and the relevant information was presented in Figure 11 and Figure S11.
According to Figure 10A, Figure 11B, and Figure S10, the guanine group of GTP forms a π-π interaction with the phenyl of F28. The distances between the mass center of the guanine group and that of phenyl are populated at 5.1 Å in four systems (Figure 11C), indicating that pY32, pY64, and pY137 hardly affect this π-π interaction. The phosphorus atom (PG) of GTP yields a salt bridge interaction with the nitrogen atom (NZ) of K16 (Figure 11B,C). The distances of PG away from NZ are located at 3.64 Å in the GTP-bound WT and pY32 KRAS (Figure 11D). This distance is distributed at 3.65 and 4.76 Å in the GTP-bound pY64 KRAS while it is populated at 3.54 and 4.68 Å in the GTP-bound pY137 (Figure 11D), showing that pY64 and pY137 weaken the stability of this salt bridge interaction. The phosphorus atom (PB) of GTP forms a salt bridge with the nitrogen atom (NZ) of K16 (Figure 11B). The distances between PB and NZ are distributed at 3.58, 3.64, 3.60, and 3.32 Å in the GTP-bound WT, pY32, pY64, and pY137 KRAS (Figure 11E), respectively, indicating that pY32 and pY64 hardly influence this salt bridge while pY137 slightly strengthens this salt bridge. The guanine group of GTP generates a salt bridge interaction with the carbonyl of D119 (Figure 11B) and the distances between the mass center of three nitrogen atoms (N1-N3) of GTP and that of oxygen atoms (OD1 and OD2) are populated at 3.49 Å (Figure 11F), verifying that the phosphorylation does not affect this salt bridge interaction.
Based on Table 1 and Figure S10, the phosphate group of GTP forms HBIs with six residues from the P-loop, including G13, V14, G15, K16, S17, and A18. Compared to the WT KRAS, pY64 obviously decreases the occupancy of hydrogen bonds between GTP and S17 while pY137 evidently reduces the occupancy of HBIs of GTP with V14 and S17 (Table 1). GTP produces HBIs with residues V29 and D30 from the SW1 and their occupancy is lower than 30.6% (Table 1), suggesting that these hydrogen bonds are highly unstable through the GaMD simulation. By comparison with the WT KRAS, the phosphorylation decreases the occupancy of HBIs of GTP with V29 and D30 in the SW1, which explains the reason why the phosphorylation leads to more disordered states of the switch domains than the WT KRAS. GTP produces HBIs with residues N116, D119, S145, A146, and K147 (Table 1 and Figure S10) but the phosphorylation hardly affects the stability of these hydrogen bonds.
Magnesium ion plays an important role in the GTP-KRAS binding and the function of KRAS. To reveal the effect of the phosphorylation on the stability of Mg2+, the distances involved in the Mg2+-mediated electrostatic interactions with residues and GTP were calculated (Figure 11 and Figure S11). The distances of Mg2+ away from the oxygen atom O of S17 fluctuate from 1.83 to 3.16 Å in the GTP-bound WT and pY137 KRAS (Figure 11A and Figure S11A) and they are located at 2.12 Å (Figure 11B); thus, pY137 scarcely influences the interaction of Mg2+ with the OG of S17. The distances of Mg2+ away from the OG of S17 transform between two different states in the pY32 and pY64 KRAS (Figure S11A) and their probabilities have two position sites (Figure 11B), suggesting that pY32 and pY64 slightly weaken the interaction between Mg2+ and the OG of S17. The distances of Mg2+ away from the oxygen atom O of T35 fall into a fluctuation range of 3.33–21.39 Å (Figure S11B), showing that the position of the SW1 relative to Mg2+ is greatly changeable. On the whole, pY137 enhances the electrostatic interaction of Mg2+ with the O of T35 (Figure 11C). The distances between Mg2+ and the mass center of OD1 and OD2 in D57 fluctuate from 1.48 to 3.49 Å in the GTP-bound WT, pY32, and pY64 KRAS while this distance falls into a range from 3.49 to 8.29 Å in the GTP-bound pY137 KRAS (Figure S11C). The distances of Mg2+ away from the mass center of OD1 and OD2 of D57 are distributed near 2.88 Å in the GTP-bound WT, pY32, and pY64 KRAS but it is populated at 4.01 and 4.91 Å in the GTP-bound pY137 KRAS (Figure 11D), suggesting that pY137 weakens the electrostatic interaction of Mg2+ with the carbonyl of D57. Based on Figure S11D and Figure 12E, the phosphorylation hardly produces the effect on the electrostatic interaction of Mg2+ with the oxygen atom O2B of GTP. The distances of Mg2+ away from the oxygen atom O2G of GTP transform between two different states in the GTP-bound pY32 and pY64 KRAS (Figure S11E). This distance is distributed at 1.86 Å in the GTP-bound WT and pY137 KRAS (Figure 11F). Although the first distance between Mg2+ and the O2G of GTP in the GTP-bound pY32 and pY64 KRAS is also located at 1.86 Å, its second value is distributed at 3.97 and 4.19 Å, separately (Figure 11F). As a result, pY32 and pY64 weaken the electrostatic interaction of Mg2+ with the O2G of GTP.
Through the above analyses, the phosphorylation-mediated influences on interaction networks are summarized: (1) the phosphorylation changes the salt bridge interaction between the nitrogen atom NZ and the phosphorus atom PB of GTP, (2) the phosphorylation reduces the occupancy of HBIs between GTP and the SW1, and (3) the phosphorylation produces a significant effect on the electrostatic interactions of Mg2+ with T35 and D57 located at the SW1 and SW2. The Mg2+-free state not only plays an important role in the nucleotide exchange of oncogenic KRAS mutants but also results in characteristic differences in the switch regions of the WT and mutated KRAS [79,80]. In the work of Hu et al., the 1D- and 2D-potential mean forces (PMFs) constructed by using the Mg2+-T35 distance as the RCs indicate that Mg2+ participates in the conformation transformation of RAS proteins [81]. Thus, the phosphorylation-mediated impacts on the HBIs of GTP with the SW1 and Mg2+ with T35 and D57 in the switch regions can regulate the activity of KRAS.

3. Materials and Methods

3.1. Scheme of Operating Calculations

MI-GaMD simulations, DL, and PCA were integrated to recognize key residue contacts and reveal phosphorylation-mediated effects on conformations of the GTP-bound KRAS and the overall scheme is shown in Figure 12. The operating procedure is as follows: (1) the initial atomics coordinates of systems were extracted from the protein data bank (PDB) and simulation systems were constructed using the Amber program to produce force field parameters, (2) three independent GaMD simulations were performed to relax conformations and collect conformational ensembles, (3) the MDTraj program was adopted to convert conformational ensembles into images for DL, (4) the images were randomly divided in the training set and validation set to perform image classification based on the two-dimensional (2D) convolutional neural network (CNN), (5) key residue contacts were identified through classic gradient-based pixel attribution [82], and (6) the reaction coordinates (RCs) were obtained using key residue contacts to construct FELs and reveal the energy basis for the phosphorylation-mediated effects.

3.2. Constructions of Simulated Systems

The complete structure of the GTP-bound WT KRAS is unavailable in the protein data bank (PDB). The GTP-bound WT KRAS was obtained by deleting GDP and incomplete KRAS from superimposed structures of the GDP-KRAS (PDB entry: 5W22) [83] and the GTP-KRAS (PDB entry: 5VQ2) [83]. To keep atomic coordinate consistence, residues Y32, Y64, and Y137 in the above obtained GTP-WT bound KRAS were, respectively, phosphorylated as pY32, pY64, and pY137 by using the Leap module in Amber22 to produce the GTP-bound pY32, pY64, and pY137 KRAS [84,85]. The magnesium ion (Mg2+) in the crystal structure was kept at the initial model. The protonated states of residues in KRAS were examined by using the program H++ 3.0 [86] and the rational protonation was assigned to each residue. All force field parameters used in our current study were assigned through the Leap module in Amber22 according to the procedure as follows: (1) the missing hydrogen atoms in the crystal were lined to their corresponding heavy atoms, (2) the parameters of KRAS were produced with the ff19SB force field [87] and the parameters of three phosphorylated resides were obtained from the work of Homeyer et al. [88], (3) the parameters of GTP was taken from the work of Meagher et al. [89], (4) an octahedral periodic box of water with a buffer of 12.0 Å was adopted to solve each of the GTP-bound WT, pY32, pY64, and pY137 KRAS, in which the force field parameters of water molecules were obtained from the TIP3P model [90] and (5) each simulated system was neutralized by adding sodium ion (Na+) in 0.15 M NaCl salt environment, in which in which the parameters of Mg2+, Na+, and Cl− ions were derived from the study of Joung and Cheatham [91,92]. The corresponding information of the constructed systems was provided in Table S1.

3.3. Multiple Independent Gaussian Accelerated Molecular Dynamics

With the expectation for relieving bad contacts between atoms caused by initialization of four current GTP-KRAS systems, each system endured two-step minimizations consisting of 5000-cycle steepest descent minimization and 10,000-cycle conjugate gradient minimization. The optimized systems were slowly heated from 0 to 310 K within 2 ns in the canonical ensemble (NVT) with a weak harmonic restriction of 2 kcal mol−1‧Å2 on heavy atoms from the GTP-bound WT and phosphorylated KRAS. Then, four systems were further equilibrated at 310 K under the isothermal–isobaric ensemble (NPT). Subsequently, the 10-ns NPT conduction was executed to keep the density of the system in 1.01 g/cm3. Finally, independent 300-ns cMD production was implemented at the NVT with periodic boundary conditions and the particle mesh Ewald method (PME) [93], during which the initial atomic velocities of each structure were assigned by means of the Maxwell distribution. Four well-equilibrated systems were used as starting points for four independent GaMD simulations.
GaMD simulations adopt a harmonic boost potential to reduce free energy barriers in the biomolecules and enhance conformational samplings of systems. In GaMD simulations, if the potential energy V r ⃑ of the system is lower than a threshold energy E, the V r ⃑ is revised as V * ( r ⃑ ) according to the following Equations (1) and (2)
V * r ⃑ = V r ⃑ + ∆ V ( r ⃑ )
∆ V r ⃑ = 0 ,                                                         V ( r ⃑ ) ≥ E 1 2 k E − V r ⃑ 2 ,       V r ⃑ < E
in the above equations, the parameter k represents the harmonic force constant, while the parameters E and k can be regulated by following the enhanced sampling principles defined in Equations (3) and (4) as below
V m a x ≤ E ≤ V m i n + 1 k
k = k 0 1 V m a x − V m i n
where if E is set as the lower bound E = V m a x , then k 0 is determined with the Equation (5)
k 0 = min ⁡ ( 1.0 ,     σ 0 σ V · V m a x − V m i n V m a x − V a v g )
on the contrary, if E is set as the upper bound E = V m i n + 1 k , then k 0 is derived from the Equation (6)
k 0 = ( 1.0 − σ 0 σ V ) · ( V m a x − V m i n V a v g − V m i n )
of which three energy parameters V m a x , V m i n , and V a v g indicate the maximum, minimum, and averaged potential energies of the simulated systems extracted from the previous cMD simulations, separately. The parameter σV corresponds to the standard deviation of the system potential energies and the σ0 is a user-determined upper limit for accurately reweighting. In our current study, 4-µs GaMD simulations, composed of four independent GaMD simulations of 1 µs, were run on the GTP-bound WT, pY32, pY64, and pY137 KRAS. To facilitate the DL and post-processing analysis, four independent GaMD trajectories were connected into a single GaMD trajectory (SGT) and the CPPTRAJ module inlayed in Amber 22 was applied to extract data for insight into the function of KRAS [94]. A program, PyReweighting, contributed by Miao et al. was employed to accurately reweigh and recognize the original free energy profiles of our current KRAS systems [95]. All cMD and GaMD simulations use the SHAKE algorithm to constrain chemical bonds connecting hydrogen atoms to heavy ones [96]. The Langevin thermostat with a collision frequency of 2.0 ps−1 was utilized to regulate the temperatures of four KRAS-related systems [97]. The PME method with a 12-Å cutoff was utilized to estimate the non-bonded interactions. The program pmemd.cuda implemented in Amber 20 [98,99] was employed to run all simulations. PCA was performed by using the CPPTRAJ module [94] to uncover the phosphorylation-mediated effect on dynamics behavior of KRAS and the details have been provided in our previous works [27,100].

3.4. Deep Learning

DL was utilized to decode the phosphorylation-mediated impacts on conformational dynamics of KRAS and our entire scheme, including system preparation, GaMD simulations, image extraction, DL, and construction of FELs, was depicted in Figure 2. The residue contact map in each GaMD trajectory snapshot was calculated by utilizing the Python packages MDTraj1.9.8.dev0 [101] and contact map explorer [101], in which a contact definition of ≤4.5 Å between any Cα atoms was adopted. Then, the obtained 170 × 170 residue contacts were converted into the 170 × 170-pixel grayscale images that can be identified by a two-dimensional (2D) convolutional neural network (CNN). Then, 80% and 20% of images randomly selected from 160,000 images extracted from the SGT by the MDTraj were, respectively, used for training and validation. The current 2D-CNN was constructed by using the Keras3.0 module [102] and the Python PyTorch package (https://pytorch.org/). The current CNN model was mainly composed of four convolutional layers of 3 × 3 kernel size, with 32, 32, 64, and 64 filters, individually, followed by three fully connected layers, the first two of which included 512 and 128 filters with a dropout rate of 0.5 each. The fully connected layer was the classification layer for discriminating four different KRAS. Apart from the classification layers, the“ReLu” activation function was utilized in all layers of the 2D-CNN, from which the“softmax” activation function was employed. A maximum pooling layer with 2 × 2 kernel size was followed after each convolutional layer. Lastly, the backpropagation by vanilla gradient-based pixel attribution [82] was applied to compute the saliency map of residue contact gradients so as to identify the function difference of KRAS caused by phosphorylation, in which residue contact was mapped onto the most populated structural cluster of each KRAS-related systems. We declare that our program was rewritten with the PyTorch package based on the work of Miao’s group [67] and a previous work has been performed [103].

3.5. Construction of Free Energy Landscapes

Key residue contacts were detected by the DL according to three following criterions: (1) the contact gradients were higher than 0.75 in the saliency maps, (2) the residue contacts were calculated by using all Cα atoms, and (3) obvious changes in contacts between structure domains can be captured. The distances between the Cα atoms of key residues recognized by the DL and root-mean-square deviations (RMSDs) of backbone atoms from KRAS were used as reaction coordinates (RCs) to reweigh GaMD simulations and construct FELs.
In reweighting of GaMD simulations, the reweighted free energy F A = − k B T l n ( ρ A ) is computed as
F A = F * A − ∑ k = 1 2 β k k ! C k + F C
where F * A = − k B T l n p * A represents the modified free energy arising from GaMD simulations, F C indicates a constant, and β = k B T . The probability distribution p * A of selected RCs from GaMD simulations can be reweighted to recover the canonical ensemble distribution ρ A . All calculations in the free energy reweighting were realized by using the program PyReweighting developed by Miao et al. and the detail for the reweighting procedure has been clarified in the work of Miao et al. [95].

4. Conclusions

The phosphorylation of different sites produces significant impacts on the activity of KRAS. Four independent GaMD simulations, each for running 1 μs, were conducted to probe the influences of the phosphorylation of pY32, pY64, and pY137 on conformational dynamics of the GTP-bound KRAS. The CNN-based DL was performed to identify significant functional domains and the results reveal that the switch domains are involved in obvious difference in structural contacts. The RCs learned by the DL were used to construct FELs and the results indicate that the phosphorylation not only induces more energy states than the WT KRAS but also leads to more disordered states of the switch domains. The information from the PCA suggests that the phosphorylation alters structural fluctuation in the switch domains. The interaction networks uncover that the phosphorylation changes not only the HBIs of the GTP with the SW1 but also the electrostatic interactions of Mg2+ with SW1 and SW2. The SW1 and SW2 participate in binding of KRAS to its effectors or regulators; hence, it is known that the alterations in conformational dynamics of the switch domains certainly affect the function of KRAS. This work is also anticipated to provide useful theoretical aids for deeply understanding the function of KRAS.

Supplementary Materials

The following supporting information can be downloaded at https://www.mdpi.com/article/10.3390/molecules29102317/s1. Table S1: The corresponding information for construction of simulation systems. Figure S1: Learning curves of the training and validation datasets for four GTP-bound KRAS systems; Figure S2: RMSDs of backbone atoms from the GTP-bound WT, pY32, pY64 and pY137 KRAS; Figure S3: Free energy profiles and representative structures of the GTP-WT KRAS; Figure S4: Geometric positions of key residues in representative structures of the GTP-pY32 KRAS; Figure S5 Geometric positions of key residues in representative structures of the GTP-pY64 KRAS; Figure S6: Geometric positions of key residues in representative structures of the GTP-pY137 KRAS; Figure S7: The time course of RMSDs for heavy atoms of GTP; Figure S8: The distances between the mass centers of all Cα atoms in the helix α2 and α3 and the distances between the mass centers of all Cα atoms in the helix α3 and α4; Figure S9: The function of eigenvalues as eigenvector indexes; Figure S10 Geometric position of key residues relative to GTP; Figure S11 The time course of distances for Mg2+-mediated electrostatic interactions: (A) the distances between Mg2+ and the oxygen atom OG of S17, (B) the distances of Mg2+ away from the oxygen atom O of T35, (C) the distances of Mg2+ away from the mass center of oxygen atoms of OD1 and OD2 in D57, (D) the distances of Mg2+ away from oxygen atom O2B of GTP and (E) the distances between Mg2+ and oxygen atom OG of GTP.

Author Contributions

Conceptualization, J.C. and G.H.; methodology, J.C. and G.H.; software, J.W., W.Y. and L.Z.; validation, J.W., W.Y. and L.Z.; formal analysis, J.C. and J.Z.; investigation, J.C.; resources, J.W.; data curation, W.Y.; writing—original draft preparation, J.C.; writing—review and editing, J.C.; visualization, J.C.; supervision, J.Z.; project administration, J.C.; funding acquisition, J.C. All authors have read and agreed to the published version of the manuscript.

Funding

This work was supported by the high-caliber talent of a Tuojiang scholar from Shandong Jiaotong University (No. TJXZ202204) and Natural Science Foundation of Shandong Province Grant (ZR2021MA069) and Key Research and Development Project of Shandong province (No. 2019GGX102050).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data are contained within the article and Supplementary Materials.

Conflicts of Interest

The authors declare no conflicts of interest.

References

  1. Crespo, P.; León, J. Ras proteins in the control of the cell cycle and cell differentiation. Cell. Mol. Life Sci. 2000, 57, 1613–1636. [Google Scholar] [CrossRef]
  2. Wang, Y.; Ji, D.; Lei, C.; Chen, Y.; Qiu, Y.; Li, X.; Li, M.; Ni, D.; Pu, J.; Zhang, J.; et al. Mechanistic insights into the effect of phosphorylation on Ras conformational dynamics and its interactions with cell signaling proteins. Comput. Struct. Biotec. 2021, 19, 1184–1199. [Google Scholar] [CrossRef]
  3. Drosten, M.; Dhawahir, A.; Sum, E.Y.M.; Urosevic, J.; Lechuga, C.G.; Esteban, L.M.; Castellano, E.; Guerra, C.; Santos, E.; Barbacid, M. Genetic analysis of Ras signalling pathways in cell proliferation, migration and survival. EMBO J. 2010, 29, 1091–1104. [Google Scholar] [CrossRef]
  4. Brown, W.; Wesalo, J.; Tsang, M.; Deiters, A. Engineering Small Molecule Switches of Protein Function in Zebrafish Embryos. J. Am. Chem. Soc. 2023, 145, 2395–2403. [Google Scholar] [CrossRef]
  5. Lu, S.; Jang, H.; Muratcioglu, S.; Gursoy, A.; Keskin, O.; Nussinov, R.; Zhang, J. Ras Conformational Ensembles, Allostery, and Signaling. Chem. Rev. 2016, 116, 6607–6665. [Google Scholar] [CrossRef]
  6. Lu, S.; Jang, H.; Gu, S.; Zhang, J.; Nussinov, R. Drugging Ras GTPase: A comprehensive mechanistic and signaling structural view. Chem. Soc. Rev. 2016, 45, 4929–4952. [Google Scholar] [CrossRef]
  7. Vetter, I.R.; Wittinghofer, A. The Guanine Nucleotide-Binding Switch in Three Dimensions. Science 2001, 294, 1299–1304. [Google Scholar] [CrossRef]
  8. Zeng, J.; Chen, J.; Xia, F.; Cui, Q.; Deng, X.; Xu, X. Identification of functional substates of KRas during GTP hydrolysis with enhanced sampling simulations. Phys. Chem. Chem. Phys. 2022, 24, 7653–7665. [Google Scholar] [CrossRef]
  9. Kano, Y.; Cook, J.D.; Lee, J.E.; Ohh, M. New structural and functional insight into the regulation of Ras. Semin. Cell Dev. Biol. 2016, 58, 70–78. [Google Scholar] [CrossRef]
  10. Tate, J.G.; Bamford, S.; Jubb, H.C.; Sondka, Z.; Beare, D.M.; Bindal, N.; Boutselakis, H.; Cole, C.G.; Creatore, C.; Dawson, E.; et al. COSMIC: The Catalogue of Somatic Mutations in Cancer. Nucleic Acids Res. 2018, 47, D941–D947. [Google Scholar] [CrossRef] [PubMed]
  11. Bos, J.L.; Rehmann, H.; Wittinghofer, A. GEFs and GAPs: Critical Elements in the Control of Small G Proteins. Cell 2007, 129, 865–877. [Google Scholar] [CrossRef] [PubMed]
  12. Moore, G.W.K.; Howell, S.E.L.; Brady, M.; Xu, X.; McNeil, K. Anomalous collapses of Nares Strait ice arches leads to enhanced export of Arctic sea ice. Nat. Commun. 2021, 12, 1. [Google Scholar] [CrossRef] [PubMed]
  13. Boriack-Sjodin, P.A.; Margarit, S.M.; Bar-Sagi, D.; Kuriyan, J. The structural basis of the activation of Ras by Sos. Nature 1998, 394, 337–343. [Google Scholar] [CrossRef] [PubMed]
  14. Rabara, D.; Tran, T.H.; Dharmaiah, S.; Stephens, R.M.; McCormick, F.; Simanshu, D.K.; Holderfield, M. KRAS G13D sensitivity to neurofibromin-mediated GTP hydrolysis. Proc. Natl. Acad. Sci. USA 2019, 116, 22122–22131. [Google Scholar] [CrossRef] [PubMed]
  15. Xiong, Y.; Zeng, J.; Xia, F.; Cui, Q.; Deng, X.; Xu, X. Conformations and binding pockets of HRas and its guanine nucleotide exchange factors complexes in the guanosine triphosphate exchange process. J. Comput. Chem. 2022, 43, 906–916. [Google Scholar] [CrossRef]
  16. Berta, D.; Gehrke, S.; Nyíri, K.; Vértessy, B.G.; Rosta, E. Mechanism-Based Redesign of GAP to Activate Oncogenic Ras. J. Am. Chem. Soc. 2023, 145, 20302–20310. [Google Scholar] [CrossRef] [PubMed]
  17. Narayan, B.; Kiel, C.; Buchete, N.-V. Classification of GTP-dependent K-Ras4B active and inactive conformational states. J. Chem. Phys. 2023, 158, 091104. [Google Scholar] [CrossRef] [PubMed]
  18. Chao, F.-A.; Dharmaiah, S.; Taylor, T.; Messing, S.; Gillette, W.; Esposito, D.; Nissley, D.V.; McCormick, F.; Byrd, R.A.; Simanshu, D.K.; et al. Insights into the Cross Talk between Effector and Allosteric Lobes of KRAS from Methyl Conformational Dynamics. J. Am. Chem. Soc. 2022, 144, 4196–4205. [Google Scholar] [CrossRef] [PubMed]
  19. Chen, J.; Zhang, S.; Wang, W.; Pang, L.; Zhang, Q.; Liu, X. Mutation-Induced Impacts on the Switch Transformations of the GDP- and GTP-Bound K-Ras: Insights from Multiple Replica Gaussian Accelerated Molecular Dynamics and Free Energy Analysis. J. Chem. Inf. Model. 2021, 61, 1954–1969. [Google Scholar] [CrossRef]
  20. Bao, H.; Wang, W.; Sun, H.; Chen, J. Probing mutation-induced conformational transformation of the GTP/M-RAS complex through Gaussian accelerated molecular dynamics simulations. J. Enzym. Inhib. Med. Chem. 2023, 38, 2195995. [Google Scholar] [CrossRef]
  21. Bao, H.; Wang, W.; Sun, H.; Chen, J. The switch states of the GDP-bound HRAS affected by point mutations: A study from Gaussian accelerated molecular dynamics simulations and free energy landscapes. J. Biomol. Struct. Dyn. 2024, 42, 3363–3381. [Google Scholar] [CrossRef] [PubMed]
  22. Spoerner, M.; Herrmann, C.; Vetter, I.R.; Kalbitzer, H.R.; Wittinghofer, A. Dynamic properties of the Ras switch I region and its importance for binding to effectors. Proc. Natl. Acad. Sci. USA 2001, 98, 4944–4949. [Google Scholar] [CrossRef]
  23. Hocker, H.J.; Cho, K.-J.; Chen, C.-Y.K.; Rambahal, N.; Sagineedu, S.R.; Shaari, K.; Stanslas, J.; Hancock, J.F.; Gorfe, A.A. Andrographolide derivatives inhibit guanine nucleotide exchange and abrogate oncogenic Ras function. Proc. Natl. Acad. Sci. USA 2013, 110, 10201–10206. [Google Scholar] [CrossRef] [PubMed]
  24. Sakamoto, K.; Kamada, Y.; Sameshima, T.; Yaguchi, M.; Niida, A.; Sasaki, S.; Miwa, M.; Ohkubo, S.; Sakamoto, J.-i.; Kamaura, M.; et al. K-Ras(G12D)-selective inhibitory peptides generated by random peptide T7 phage display technology. Biochem. Biophys. Res. Commun. 2017, 484, 605–611. [Google Scholar] [CrossRef] [PubMed]
  25. Poulin, E.J.; Bera, A.K.; Lu, J.; Lin, Y.-J.; Strasser, S.D.; Paulo, J.A.; Huang, T.Q.; Morales, C.; Yan, W.; Cook, J.; et al. Tissue-Specific Oncogenic Activity of KRASA146T. Cancer Discov. 2019, 9, 738–755. [Google Scholar] [CrossRef] [PubMed]
  26. Hu, Z.; Marti, J. In silico drug derivatives for KRAS-G12D: Free-energy surfaces in aqueous solution by well-tempered metadynamics simulations. Mol. Phys. 2024, e2316883. [Google Scholar] [CrossRef]
  27. Mao, C.; Huang, Y.-F.; Yang, Z.-Y.; Zheng, D.-Y.; Chen, J.-Z.; Tang, J.-L. KRAS p.G13D mutation and codon 12 mutations are not created equal in predicting clinical outcomes of cetuximab in metastatic colorectal cancer. Cancer 2013, 119, 714–721. [Google Scholar] [CrossRef]
  28. Hu, Y.; Lu, W.; Chen, G.; Wang, P.; Chen, Z.; Zhou, Y.; Ogasawara, M.; Trachootham, D.; Feng, L.; Pelicano, H.; et al. K-rasG12V transformation leads to mitochondrial dysfunction and a metabolic switch from oxidative phosphorylation to glycolysis. Cell Res. 2012, 22, 399–412. [Google Scholar] [CrossRef]
  29. Mehaffey, M.R.; Schardon, C.L.; Novelli, E.T.; Cammarata, M.B.; Webb, L.J.; Fast, W.; Brodbelt, J.S. Investigation of GTP-dependent dimerization of G12X K-Ras variants using ultraviolet photodissociation mass spectrometry. Chem. Sci. 2019, 10, 8025–8034. [Google Scholar] [CrossRef]
  30. Hunter, J.C.; Manandhar, A.; Carrasco, M.A.; Gurbani, D.; Gondi, S.; Westover, K.D. Biochemical and Structural Analysis of Common Cancer-Associated KRAS Mutations. Mol. Cancer Res. 2015, 13, 1325–1335. [Google Scholar] [CrossRef]
  31. Chen, S.; Shu, L.; Zhao, R.; Zhao, Y. Molecular dynamics simulations reveal the activation mechanism of mutations G12V and Q61L of Cdc42. Proteins 2022, 90, 1376–1389. [Google Scholar] [CrossRef] [PubMed]
  32. Buhrman, G.; Wink, G.; Mattos, C. Transformation Efficiency of RasQ61 Mutants Linked to Structural Features of the Switch Regions in the Presence of Raf. Structure 2007, 15, 1618–1629. [Google Scholar] [CrossRef] [PubMed]
  33. Cruz-Migoni, A.; Canning, P.; Quevedo, C.E.; Bataille, C.J.R.; Bery, N.; Miller, A.; Russell, A.J.; Phillips, S.E.V.; Carr, S.B.; Rabbitts, T.H. Structure-based development of new RAS-effector inhibitors from a combination of active and inactive RAS-binding compounds. Proc. Natl. Acad. Sci. USA 2019, 116, 2545–2550. [Google Scholar] [CrossRef] [PubMed]
  34. Parker, J.A.; Volmar, A.Y.; Pavlopoulos, S.; Mattos, C. K-Ras Populates Conformational States Differently from Its Isoform H-Ras and Oncogenic Mutant K-RasG12D. Structure 2018, 26, 810–820.e814. [Google Scholar] [CrossRef] [PubMed]
  35. Matsumoto, S.; Taniguchi-Tamura, H.; Araki, M.; Kawamura, T.; Miyamoto, R.; Tsuda, C.; Shima, F.; Kumasaka, T.; Okuno, Y.; Kataoka, T. Oncogenic mutations Q61L and Q61H confer active form-like structural features to the inactive state (state 1) conformation of H-Ras protein. Biochem. Biophys. Res. Commun. 2021, 565, 85–90. [Google Scholar] [CrossRef] [PubMed]
  36. Johnson, C.W.; Lin, Y.-J.; Reid, D.; Parker, J.; Pavlopoulos, S.; Dischinger, P.; Graveel, C.; Aguirre, A.J.; Steensma, M.; Haigis, K.M.; et al. Isoform-Specific Destabilization of the Active Site Reveals a Molecular Mechanism of Intrinsic Activation of KRas G13D. Cell Rep. 2019, 28, 1538–1550.e1537. [Google Scholar] [CrossRef]
  37. Ting, P.Y.; Johnson, C.W.; Fang, C.; Cao, X.; Graeber, T.G.; Mattos, C.; Colicelli, J. Tyrosine phosphorylation of RAS by ABL allosterically enhances effector binding. FASEB J. 2015, 29, 3750–3761. [Google Scholar] [CrossRef] [PubMed]
  38. Bunda, S.; Heir, P.; Srikumar, T.; Cook, J.D.; Burrell, K.; Kano, Y.; Lee, J.E.; Zadeh, G.; Raught, B.; Ohh, M. Src promotes GTPase activity of Ras via tyrosine 32 phosphorylation. Proc. Natl. Acad. Sci. USA 2014, 111, E3785–E3794. [Google Scholar] [CrossRef]
  39. Kano, Y.; Gebregiworgis, T.; Marshall, C.B.; Radulovich, N.; Poon, B.P.K.; St-Germain, J.; Cook, J.D.; Valencia-Sama, I.; Grant, B.M.M.; Herrera, S.G.; et al. Tyrosyl phosphorylation of KRAS stalls GTPase cycle via alteration of switch I and II conformation. Nat. Commun. 2019, 10, 224. [Google Scholar] [CrossRef]
  40. Khaled, M.; Gorfe, A.; Sayyed-Ahmad, A. Conformational and Dynamical Effects of Tyr32 Phosphorylation in K-Ras: Molecular Dynamics Simulation and Markov State Models Analysis. J. Phys. Chem. B 2019, 123, 7667–7675. [Google Scholar] [CrossRef]
  41. Qiu, Y.; Wang, Y.; Chai, Z.; Ni, D.; Li, X.; Pu, J.; Chen, J.; Zhang, J.; Lu, S.; Lv, C.; et al. Targeting RAS phosphorylation in cancer therapy: Mechanisms and modulators. Acta Pharm. Sin. B 2021, 11, 3433–3446. [Google Scholar] [CrossRef] [PubMed]
  42. Sun, Z.; Gong, Z.; Xia, F.; He, X. Ion dynamics and selectivity of Nav channels from molecular dynamics simulation. Chem. Phys. 2021, 548, 111245. [Google Scholar] [CrossRef]
  43. Wang, Y.; Yang, F.; Yan, D.; Zeng, Y.; Wei, B.; Chen, J.; He, W. Identification Mechanism of BACE1 on Inhibitors Probed by Using Multiple Separate Molecular Dynamics Simulations and Comparative Calculations of Binding Free Energies. Molecules 2023, 28, 4773. [Google Scholar] [CrossRef] [PubMed]
  44. Xue, W.; Yang, F.; Wang, P.; Zheng, G.; Chen, Y.; Yao, X.; Zhu, F. What Contributes to Serotonin–Norepinephrine Reuptake Inhibitors’ Dual-Targeting Mechanism? The Key Role of Transmembrane Domain 6 in Human Serotonin and Norepinephrine Transporters Revealed by Molecular Dynamics Simulation. ACS Chem. Neurosci. 2018, 9, 1128–1140. [Google Scholar] [CrossRef] [PubMed]
  45. Sun, J.; Liu, X.; Zhang, S.; Li, M.; Zhang, Q.; Chen, J. Molecular insights and optimization strategies for the competitive binding of engineered ACE2 proteins: A multiple replica molecular dynamics study. Phys. Chem. Chem. Phys. 2023, 25, 28479–28496. [Google Scholar] [CrossRef] [PubMed]
  46. Sun, H.; Li, Y.; Shen, M.; Tian, S.; Xu, L.; Pan, P.; Guan, Y.; Hou, T. Assessing the performance of MM/PBSA and MM/GBSA methods. 5. Improved docking performance using high solute dielectric constant MM/GBSA and MM/PBSA rescoring. Phys. Chem. Chem. Phys. 2014, 16, 22035–22045. [Google Scholar] [CrossRef] [PubMed]
  47. Gupta, A.; Chaudhary, N.; Aparoy, P. MM-PBSA and per-residue decomposition energy studies on 7-Phenyl-imidazoquinolin-4(5H)-one derivatives: Identification of crucial site points at microsomal prostaglandin E synthase-1 (mPGES-1) active site. Int. J. Biol. Macromol. 2018, 119, 352–359. [Google Scholar] [CrossRef] [PubMed]
  48. Wang, X.-F.; Sun, J.; Wang, X.-L.; Tian, J.-K.; Tian, Z.-W.; Zhang, J.-L.; Jia, R. MD investigation on the binding of microphthalmia-associated transcription factor with DNA. J. Saudi Chem. Soc. 2022, 26, 101420. [Google Scholar] [CrossRef]
  49. Miao, Y.; Feher, V.A.; McCammon, J.A. Gaussian Accelerated Molecular Dynamics: Unconstrained Enhanced Sampling and Free Energy Calculation. J. Chem. Theory Comput. 2015, 11, 3584–3595. [Google Scholar] [CrossRef]
  50. Wang, J.; Miao, Y. Peptide Gaussian accelerated molecular dynamics (Pep-GaMD): Enhanced sampling and free energy and kinetics calculations of peptide binding. J. Chem. Phys. 2020, 153, 154109. [Google Scholar] [CrossRef]
  51. Wang, J.; Arantes, P.R.; Bhattarai, A.; Hsu, R.V.; Pawnikar, S.; Huang, Y.-M.M.; Palermo, G.; Miao, Y. Gaussian accelerated molecular dynamics: Principles and applications. WIREs Comput. Mol. Sci. 2021, 11, e1521. [Google Scholar] [CrossRef] [PubMed]
  52. Sun, Z.; Huai, Z.; He, Q.; Liu, Z. A General Picture of Cucurbit[8]uril Host–Guest Binding. J. Chem. Inf. Model. 2021, 61, 6107–6134. [Google Scholar] [CrossRef]
  53. Yang, F.; Wang, Y.; Yan, D.; Liu, Z.; Wei, B.; Chen, J.; He, W. Binding Mechanism of Inhibitors to Heat Shock Protein 90 Investigated by Multiple Independent Molecular Dynamics Simulations and Prediction of Binding Free Energy. Molecules 2023, 28, 4792. [Google Scholar] [CrossRef] [PubMed]
  54. Sun, H.; Li, Y.; Tian, S.; Xu, L.; Hou, T. Assessing the performance of MM/PBSA and MM/GBSA methods. 4. Accuracies of MM/PBSA and MM/GBSA methodologies evaluated by various simulation protocols using PDBbind data set. Phys. Chem. Chem. Phys. 2014, 16, 16719–16729. [Google Scholar] [CrossRef] [PubMed]
  55. Xue, W.; Wang, P.; Tu, G.; Yang, F.; Zheng, G.; Li, X.; Li, X.; Chen, Y.; Yao, X.; Zhu, F. Computational identification of the binding mechanism of a triple reuptake inhibitor amitifadine for the treatment of major depressive disorder. Phys. Chem. Chem. Phys. 2018, 20, 6606–6616. [Google Scholar] [CrossRef] [PubMed]
  56. Hou, T.; Yu, R. Molecular Dynamics and Free Energy Studies on the Wild-type and Double Mutant HIV-1 Protease Complexed with Amprenavir and Two Amprenavir-Related Inhibitors: Mechanism for Binding and Drug Resistance. J. Med. Chem. 2007, 50, 1177–1188. [Google Scholar] [CrossRef] [PubMed]
  57. Wang, Q.; Luo, S.; Xiong, D.; Xu, X.; Zhao, X.; Duan, L. Quantitative investigation of the effects of DNA modifications and protein mutations on MeCP2-MBD-DNA interactions. Int. J. Biol. Macromol. 2023, 247, 125690. [Google Scholar] [CrossRef] [PubMed]
  58. Wang, J.; Miao, Y. Mechanistic Insights into Specific G Protein Interactions with Adenosine Receptors. J. Phys. Chem. B 2019, 123, 6462–6473. [Google Scholar] [CrossRef] [PubMed]
  59. Miao, Y.; McCammon, J.A. Graded activation and free energy landscapes of a muscarinic G-protein–coupled receptor. Proc. Natl. Acad. Sci. USA 2016, 113, 12162–12167. [Google Scholar] [CrossRef]
  60. Chen, J.; Zeng, Q.; Wang, W.; Sun, H.; Hu, G. Decoding the Identification Mechanism of an SAM-III Riboswitch on Ligands through Multiple Independent Gaussian-Accelerated Molecular Dynamics Simulations. J. Chem. Inf. Model. 2022, 62, 6118–6132. [Google Scholar] [CrossRef]
  61. Wang, J.; Lan, L.; Wu, X.; Xu, L.; Miao, Y. Mechanism of RNA recognition by a Musashi RNA-binding protein. Curr. Res. Struct. Biol. 2022, 4, 10–20. [Google Scholar] [CrossRef]
  62. Célerse, F.; Inizan, T.J.; Lagardère, L.; Adjoua, O.; Monmarché, P.; Miao, Y.; Derat, E.; Piquemal, J.-P. An Efficient Gaussian-Accelerated Molecular Dynamics (GaMD) Multilevel Enhanced Sampling Strategy: Application to Polarizable Force Fields Simulations of Large Biological Systems. J. Chem. Theory Comput. 2022, 18, 968–977. [Google Scholar] [CrossRef]
  63. An, X.; Bai, Q.; Bing, Z.; Liu, H.; Yao, X. Insights into the molecular mechanism of positive cooperativity between partial agonist MK-8666 and full allosteric agonist AP8 of hGPR40 by Gaussian accelerated molecular dynamics (GaMD) simulations. Comput. Struct. Biotec. 2021, 19, 3978–3989. [Google Scholar] [CrossRef]
  64. Chen, J.; Wang, W.; Sun, H.; He, W. Roles of Accelerated Molecular Dynamics Simulations in Predictions of Binding Kinetic Parameters. Mini-Rev. Med. Chem. 2024; Online ahead of print. [Google Scholar] [CrossRef]
  65. Plante, A.; Shore, D.M.; Morra, G.; Khelashvili, G.; Weinstein, H. A Machine Learning Approach for the Discovery of Ligand-Specific Functional Mechanisms of GPCRs. Molecules 2019, 24, 2097. [Google Scholar] [CrossRef]
  66. Plante, A.; Weinstein, H. Ligand-Dependent Conformational Transitions in Molecular Dynamics Trajectories of GPCRs Revealed by a New Machine Learning Rare Event Detection Protocol. Molecules 2021, 26, 3059. [Google Scholar] [CrossRef]
  67. Do, H.N.; Wang, J.; Bhattarai, A.; Miao, Y. GLOW: A Workflow Integrating Gaussian-Accelerated Molecular Dynamics and Deep Learning for Free Energy Profiling. J. Chem. Theory Comput. 2022, 18, 1423–1436. [Google Scholar] [CrossRef] [PubMed]
  68. Do, H.N.; Wang, J.; Miao, Y. Deep Learning Dynamic Allostery of G-Protein-Coupled Receptors. JACS Au 2023, 3, 3165–3180. [Google Scholar] [CrossRef]
  69. Amadei, A.; Linssen, A.B.M.; Berendsen, H.J.C. Essential dynamics of proteins. Proteins 1993, 17, 412–425. [Google Scholar] [CrossRef] [PubMed]
  70. Ichiye, T.; Karplus, M. Collective motions in proteins: A covariance analysis of atomic fluctuations in molecular dynamics and normal mode simulations. Proteins 1991, 11, 205–217. [Google Scholar] [CrossRef] [PubMed]
  71. Yan, F.; Liu, X.; Zhang, S.; Su, J.; Zhang, Q.; Chen, J. Molecular Dynamics Exploration of Selectivity of Dual Inhibitors 5M7, 65X, and 65Z toward Fatty Acid Binding Proteins 4 and 5. Int. J. Mol. Sci. 2018, 19, 2496. [Google Scholar] [CrossRef]
  72. Buhrman, G.; Kumar, V.S.S.; Cirit, M.; Haugh, J.M.; Mattos, C. Allosteric Modulation of Ras-GTP Is Linked to Signal Transduction through RAF Kinase. J. Biol. Chem. 2011, 286, 3323–3331. [Google Scholar] [CrossRef] [PubMed]
  73. Shima, F.; Ijiri, Y.; Muraoka, S.; Liao, J.; Ye, M.; Araki, M.; Matsumoto, K.; Yamamoto, N.; Sugimoto, T.; Yoshikawa, Y.; et al. Structural Basis for Conformational Dynamics of GTP-bound Ras Protein. J. Biol. Chem. 2010, 285, 22696–22705. [Google Scholar] [CrossRef] [PubMed]
  74. Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual molecular dynamics. J. Mol. Graph. 1996, 14, 33–38. [Google Scholar] [CrossRef] [PubMed]
  75. Eren, M.; Tuncbag, N.; Jang, H.; Nussinov, R.; Gursoy, A.; Keskin, O. Normal Mode Analysis of KRas4B Reveals Partner Specific Dynamics. J. Phys. Chem. B 2021, 125, 5210–5221. [Google Scholar] [CrossRef] [PubMed]
  76. Rojas, A.M.; Fuentes, G.; Rausell, A.; Valencia, A. The Ras protein superfamily: Evolutionary tree and role of conserved amino acids. J. Cell Biol. 2012, 196, 189–201. [Google Scholar] [CrossRef] [PubMed]
  77. Adasme, M.F.; Linnemann, K.L.; Bolz, S.N.; Kaiser, F.; Salentin, S.; Haupt, V.J.; Schroeder, M. PLIP 2021: Expanding the scope of the protein–ligand interaction profiler to DNA and RNA. Nucleic Acids Res. 2021, 49, W530–W534. [Google Scholar] [CrossRef] [PubMed]
  78. Salentin, S.; Schreiber, S.; Haupt, V.J.; Adasme, M.F.; Schroeder, M. PLIP: Fully automated protein–ligand interaction profiler. Nucleic Acids Res. 2015, 43, W443–W447. [Google Scholar] [CrossRef] [PubMed]
  79. Pálfy, G.; Menyhárd, D.K.; Ákontz-Kiss, H.; Vida, I.; Batta, G.; Tőke, O.; Perczel, A. The Importance of Mg2+-Free State in Nucleotide Exchange of Oncogenic K-Ras Mutants. Chem. Eur. J. 2022, 28, e202201449. [Google Scholar] [CrossRef] [PubMed]
  80. Hu, Z.; Marti, J. Discovering and Targeting Dynamic Drugging Pockets of Oncogenic Proteins: The Role of Magnesium in Conformational Changes of the G12D Mutated Kirsten Rat Sarcoma-Guanosine Diphosphate Complex. Int. J. Mol. Sci. 2022, 23, 13865. [Google Scholar] [CrossRef]
  81. Hu, F.; Wang, Y.; Zeng, J.; Deng, X.; Xia, F.; Xu, X. Unveiling the State Transition Mechanisms of Ras Proteins through Enhanced Sampling and QM/MM Simulations. J. Phys. Chem. B 2024, 128, 1418–1427. [Google Scholar] [CrossRef]
  82. Kotikalapudi, R. Keras-Vis; GitHub: San Francisco, CA, USA, 2017. [Google Scholar]
  83. Xu, S.; Long, B.N.; Boris, G.H.; Chen, A.; Ni, S.; Kennedy, M.A. Structural insight into the rearrangement of the switch I region in GTP-bound G12A K-Ras. Acta Crystallogr. D Struct. Biol. 2017, 73, 970–984. [Google Scholar] [CrossRef] [PubMed]
  84. Case, D.A.; Cheatham, T.E., III; Darden, T.; Gohlke, H.; Luo, R.; Merz, K.M., Jr.; Onufriev, A.; Simmerling, C.; Wang, B.; Woods, R.J. The Amber biomolecular simulation programs. J. Comput. Chem. 2005, 26, 1668–1688. [Google Scholar] [CrossRef] [PubMed]
  85. Salomon-Ferrer, R.; Case, D.A.; Walker, R.C. An overview of the Amber biomolecular simulation package. WIREs Comput. Mol. Sci. 2013, 3, 198–210. [Google Scholar] [CrossRef]
  86. Anandakrishnan, R.; Aguilar, B.; Onufriev, A.V. H++ 3.0: Automating pK prediction and the preparation of biomolecular structures for atomistic molecular modeling and simulations. Nucleic Acids Res. 2012, 40, W537–W541. [Google Scholar] [CrossRef] [PubMed]
  87. Tian, C.; Kasavajhala, K.; Belfon, K.A.A.; Raguette, L.; Huang, H.; Migues, A.N.; Bickel, J.; Wang, Y.; Pincay, J.; Wu, Q.; et al. ff19SB: Amino-Acid-Specific Protein Backbone Parameters Trained against Quantum Mechanics Energy Surfaces in Solution. J. Chem. Theory Comput. 2020, 16, 528–552. [Google Scholar] [CrossRef]
  88. Homeyer, N.; Horn, A.H.C.; Lanig, H.; Sticht, H. AMBER force-field parameters for phosphorylated amino acids in different protonation states: Phosphoserine, phosphothreonine, phosphotyrosine, and phosphohistidine. J. Mol. Model. 2006, 12, 281–289. [Google Scholar] [CrossRef]
  89. Meagher, K.L.; Redman, L.T.; Carlson, H.A. Development of polyphosphate parameters for use with the AMBER force field. J. Comput. Chem. 2003, 24, 1016–1025. [Google Scholar] [CrossRef]
  90. Jorgensen, W.L.; Chandrasekhar, J.; Madura, J.D.; Impey, R.W.; Klein, M.L. Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 1983, 79, 926–935. [Google Scholar] [CrossRef]
  91. Joung, I.S.; Cheatham, T.E., III. Determination of Alkali and Halide Monovalent Ion Parameters for Use in Explicitly Solvated Biomolecular Simulations. J. Phys. Chem. B 2008, 112, 9020–9041. [Google Scholar] [CrossRef]
  92. Joung, I.S.; Cheatham, T.E., III. Molecular Dynamics Simulations of the Dynamic and Energetic Properties of Alkali and Halide Ions Using Water-Model-Specific Ion Parameters. J. Phys. Chem. B 2009, 113, 13279–13290. [Google Scholar] [CrossRef]
  93. Essmann, U.; Perera, L.; Berkowitz, M.L.; Darden, T.; Lee, H.; Pedersen, L.G. A smooth particle mesh Ewald method. J. Chem. Phys. 1995, 103, 8577–8593. [Google Scholar] [CrossRef]
  94. Roe, D.R.; Cheatham, T.E., III. PTRAJ and CPPTRAJ: Software for Processing and Analysis of Molecular Dynamics Trajectory Data. J. Chem. Theory Comput. 2013, 9, 3084–3095. [Google Scholar] [CrossRef]
  95. Miao, Y.; Sinko, W.; Pierce, L.; Bucher, D.; Walker, R.C.; McCammon, J.A. Improved Reweighting of Accelerated Molecular Dynamics Simulations for Free Energy Calculation. J. Chem. Theory Comput. 2014, 10, 2677–2689. [Google Scholar] [CrossRef] [PubMed]
  96. Ryckaert, J.-P.; Ciccotti, G.; Berendsen, H.J.C. Numerical integration of the cartesian equations of motion of a system with constraints: Molecular dynamics of n-alkanes. J. Comput. Phys. 1977, 23, 327–341. [Google Scholar] [CrossRef]
  97. Izaguirre, J.A.; Catarello, D.P.; Wozniak, J.M.; Skeel, R.D. Langevin stabilization of molecular dynamics. J. Chem. Phys. 2001, 114, 2090–2098. [Google Scholar] [CrossRef]
  98. Salomon-Ferrer, R.; Götz, A.W.; Poole, D.; Le Grand, S.; Walker, R.C. Routine Microsecond Molecular Dynamics Simulations with AMBER on GPUs. 2. Explicit Solvent Particle Mesh Ewald. J. Chem. Theory Comput. 2013, 9, 3878–3888. [Google Scholar] [CrossRef]
  99. Götz, A.W.; Williamson, M.J.; Xu, D.; Poole, D.; Le Grand, S.; Walker, R.C. Routine Microsecond Molecular Dynamics Simulations with AMBER on GPUs. 1. Generalized Born. J. Chem. Theory Comput. 2012, 8, 1542–1555. [Google Scholar] [CrossRef]
  100. Kulkarni, A.M.; Kumar, V.; Parate, S.; Lee, G.; Yoon, S.; Lee, K.W. Identification of New KRAS G12D Inhibitors through Computer-Aided Drug Discovery Methods. Int. J. Mol. Sci. 2022, 23, 1309. [Google Scholar] [CrossRef]
  101. McGibbon, R.T.; Beauchamp, K.A.; Harrigan, M.P.; Klein, C.; Swails, J.M.; Hernández, C.X.; Schwantes, C.R.; Wang, L.-P.; Lane, T.J.; Pande, V.S. MDTraj: A Modern Open Library for the Analysis of Molecular Dynamics Trajectories. Biophys. J. 2015, 109, 1528–1532. [Google Scholar] [CrossRef]
  102. Gulli, A.; Pal, S. Deep Learning with Keras; Packt Publishing Ltd.: Birmingham, UK, 2017. [Google Scholar]
  103. Wang, J.; Yang, W.; Zhao, L.; Wei, B.; Chen, J. Binding Mechanism of Inhibitors to BRD4 and BRD9 Decoded by Multiple Independent Molecular Dynamics Simulations and Deep Learning. Molecules 2024, 29, 1857. [Google Scholar] [CrossRef]
Molecules 29 02317 g001
Figure 1. Molecular structure: (A) the GTP-bound WT KRAS, in which secondary structure and phosphorylation group are labeled, (B) the geometry position of the phosphorylated residues relative to GTP, (C) phosphorylation with two net negative charges, and (D) GTP. In this Figure, KRAS is shown in cartoon modes while GTP and phosphorylation-related residues are shown in ball-stick modes.
Figure 1. Molecular structure: (A) the GTP-bound WT KRAS, in which secondary structure and phosphorylation group are labeled, (B) the geometry position of the phosphorylated residues relative to GTP, (C) phosphorylation with two net negative charges, and (D) GTP. In this Figure, KRAS is shown in cartoon modes while GTP and phosphorylation-related residues are shown in ball-stick modes.
Molecules 29 02317 g001
Molecules 29 02317 g002
Figure 2. Classification and saliency map of residue contact gradients: (A) classification of the GTP-bound WT, pY32, pY64, and pY137 KRAS, (B–E) the saliency map of residue contact gradients for the GTP-bound WT, pY32, pY64, and pY137 KRAS, and (F) key structural domains revealed by the DL. The gradient of each residue contact is shown in a 0.3 (white) to 0.6 (black) color scale.
Figure 2. Classification and saliency map of residue contact gradients: (A) classification of the GTP-bound WT, pY32, pY64, and pY137 KRAS, (B–E) the saliency map of residue contact gradients for the GTP-bound WT, pY32, pY64, and pY137 KRAS, and (F) key structural domains revealed by the DL. The gradient of each residue contact is shown in a 0.3 (white) to 0.6 (black) color scale.
Molecules 29 02317 g002
Molecules 29 02317 g003
Figure 3. Free energy profiles and representative structures of the GTP-pY32 KRAS: (A) FEL, (B) superimposition of initial optimized structure with the EV1, EV2, and EV3 structures, (C) structural superimposition of GTP and magnesium ions (Mg2+) in the initially optimized structure and the EV1, EV2, and EV3 structures, and (D–F) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, and EV3 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.
Figure 3. Free energy profiles and representative structures of the GTP-pY32 KRAS: (A) FEL, (B) superimposition of initial optimized structure with the EV1, EV2, and EV3 structures, (C) structural superimposition of GTP and magnesium ions (Mg2+) in the initially optimized structure and the EV1, EV2, and EV3 structures, and (D–F) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, and EV3 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.
Molecules 29 02317 g003
Molecules 29 02317 g004
Figure 4. Free energy profiles and representative structures of the GTP-pY64 KRAS: (A) FEL, (B) superimposition of initial optimized structure with the EV1, EV2, EV3, and EV4 structures, and (C–F) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, EV3, and EV4 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.
Figure 4. Free energy profiles and representative structures of the GTP-pY64 KRAS: (A) FEL, (B) superimposition of initial optimized structure with the EV1, EV2, EV3, and EV4 structures, and (C–F) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, EV3, and EV4 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.
Molecules 29 02317 g004
Molecules 29 02317 g005
Figure 5. Free energy profiles and representative structures of the GTP-pY137 KRAS: (A) FEL, (B) alignment of initial optimized structure with the EV1, EV2, EV3, and EV4 structures, and (C–F) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, EV3, and EV4 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.
Figure 5. Free energy profiles and representative structures of the GTP-pY137 KRAS: (A) FEL, (B) alignment of initial optimized structure with the EV1, EV2, EV3, and EV4 structures, and (C–F) geometric positions of the P-loop, SW1, and SW2 in the EV1, EV2, EV3, and EV4 structures, in which KRAS was shown in surface modes. The PMF is scaled in kcal/mol.
Molecules 29 02317 g005
Molecules 29 02317 g006
Figure 6. Dynamics indexes: (A) the probability distribution of GTP RMSDs, (B) RMSFs of KRAS, (C) the probability distribution of the distances between the mass center of all Cα atoms in helix α2 away from that of helix α3, and (D) the probability distribution of the distances between the mass center of all Cα atoms in helix α3 away from that of helix α4.
Figure 6. Dynamics indexes: (A) the probability distribution of GTP RMSDs, (B) RMSFs of KRAS, (C) the probability distribution of the distances between the mass center of all Cα atoms in helix α2 away from that of helix α3, and (D) the probability distribution of the distances between the mass center of all Cα atoms in helix α3 away from that of helix α4.
Molecules 29 02317 g006
Molecules 29 02317 g007
Figure 7. Concerted motions of the GTP-KRAS revealed by the PC1 from the PCA: (A) the GTP-WT KRAS, (B) the GTP-pY32 KRAS, (C) the GTP-pY64 KRAS, and (D) the GTP-pY137 KRAS. KRAS was shown in cartoon mode and GTP is shown in ball-stick modes. The blue, green and red indicate the P-loop, SW1 and SW2, respectively.
Figure 7. Concerted motions of the GTP-KRAS revealed by the PC1 from the PCA: (A) the GTP-WT KRAS, (B) the GTP-pY32 KRAS, (C) the GTP-pY64 KRAS, and (D) the GTP-pY137 KRAS. KRAS was shown in cartoon mode and GTP is shown in ball-stick modes. The blue, green and red indicate the P-loop, SW1 and SW2, respectively.
Molecules 29 02317 g007
Molecules 29 02317 g008
Figure 8. Collective movement of the GTP-KRAS revealed by the PC2 from the PCA: (A) the GTP-WT KRAS, (B) the GTP-pY32 KRAS, (C) the GTP-pY64 KRAS, and (D) the GTP-pY137 KRAS. KRAS was shown in cartoon mode and GTP is shown in ball-stick modes. The blue, green and red indicate the P-loop, SW1 and SW2, respectively.
Figure 8. Collective movement of the GTP-KRAS revealed by the PC2 from the PCA: (A) the GTP-WT KRAS, (B) the GTP-pY32 KRAS, (C) the GTP-pY64 KRAS, and (D) the GTP-pY137 KRAS. KRAS was shown in cartoon mode and GTP is shown in ball-stick modes. The blue, green and red indicate the P-loop, SW1 and SW2, respectively.
Molecules 29 02317 g008
Molecules 29 02317 g009
Figure 9. The probability of the χ angle for the sidechains of three phosphorylation-related residues: (A) residue 32, (B) residue 64, (C) residue 137, and (D) geometric positions of key residues.
Figure 9. The probability of the χ angle for the sidechains of three phosphorylation-related residues: (A) residue 32, (B) residue 64, (C) residue 137, and (D) geometric positions of key residues.
Molecules 29 02317 g009
Molecules 29 02317 g010
Figure 10. Key interactions of GTP with KRAS and their corresponding distance distributions. (A) Relative position of GTP in the binding pocket of KRAS, (B) geometry positions of key interactions, (C) the distances of the π-π interaction between GTP and F28, (D) the distances of the salt bridge interaction between the phosphorus atom PB of GTP and the nitrogen atom NZ of K16, (E) the distances of the salt bridge interaction between the phosphorus atom PG of GTP and the nitrogen atom NZ of K16, and (F) the distances of the salt bridge interaction between the guanine group of GTP and the carbonyl group of D119.
Figure 10. Key interactions of GTP with KRAS and their corresponding distance distributions. (A) Relative position of GTP in the binding pocket of KRAS, (B) geometry positions of key interactions, (C) the distances of the π-π interaction between GTP and F28, (D) the distances of the salt bridge interaction between the phosphorus atom PB of GTP and the nitrogen atom NZ of K16, (E) the distances of the salt bridge interaction between the phosphorus atom PG of GTP and the nitrogen atom NZ of K16, and (F) the distances of the salt bridge interaction between the guanine group of GTP and the carbonyl group of D119.
Molecules 29 02317 g010
Molecules 29 02317 g011
Figure 11. Key interactions of magnesium ions (Mg2+) with GTP and KRAS together with their distance distributions: (A) the geometry position of interactions, (B) the distances between Mg2+ and the oxygen atom OG of S17, (C) the distances of Mg2+ away from the oxygen atom O of T35, (D) the distances between Mg2+ and the oxygen atom O of T35, (E) the distances of Mg2+ away from the mass center of oxygen atoms of OD1 and OD2 in D57, (F) the distances of Mg2+ away from oxygen atom O2B of GTP, and (F) the distances between Mg2+ and oxygen atom OG of GTP.
Figure 11. Key interactions of magnesium ions (Mg2+) with GTP and KRAS together with their distance distributions: (A) the geometry position of interactions, (B) the distances between Mg2+ and the oxygen atom OG of S17, (C) the distances of Mg2+ away from the oxygen atom O of T35, (D) the distances between Mg2+ and the oxygen atom O of T35, (E) the distances of Mg2+ away from the mass center of oxygen atoms of OD1 and OD2 in D57, (F) the distances of Mg2+ away from oxygen atom O2B of GTP, and (F) the distances between Mg2+ and oxygen atom OG of GTP.
Molecules 29 02317 g011
Molecules 29 02317 g012
Figure 12. Workflow of deep leaning from GaMD simulations: (A) the GTP-bound KRAS with three phosphorylation sites, (B) the initialized GTP-bound KRAS, (C) conformational ensembles recorded in three independent GaMD trajectories, (D) images extracted though the MDTraj program, (E) convolution neural networks, (F) the saliency maps calculated through backward propagation, (G) key residue contacts identified by deep learning, and (H) free energy landscapes used for revealing the phosphorylation-mediated effect on free energy profiles of the GTP-bound KRAS.
Figure 12. Workflow of deep leaning from GaMD simulations: (A) the GTP-bound KRAS with three phosphorylation sites, (B) the initialized GTP-bound KRAS, (C) conformational ensembles recorded in three independent GaMD trajectories, (D) images extracted though the MDTraj program, (E) convolution neural networks, (F) the saliency maps calculated through backward propagation, (G) key residue contacts identified by deep learning, and (H) free energy landscapes used for revealing the phosphorylation-mediated effect on free energy profiles of the GTP-bound KRAS.
Molecules 29 02317 g012
Table 1. Hydrogen bonding interactions of GTP with KRAS identified by the CPPTRAJ program.
Table 1. Hydrogen bonding interactions of GTP with KRAS identified by the CPPTRAJ program.
a Hydrogen Bondsb Occupancy(%)
ResiduesGTPWTpY32pY64pY137
G13-N-HO3B89.288.189.387.2
V14-N-HO1B20.120.319.613.1
G15-N-HO1B99.098.397.699.5
K16-N-HO1B99.999.998.899.9
S17-N-HO2B99.596.287.388.8
A18-N-HO1A98.696.499.299.8
V29-OO2′-HO’230.616.415.717.4
D30-OO2′-HO’227.218.218.411.1
N116-ND2-HD21N790.087.288.791.1
D119-OD1N1-H191.193.393.190.4
D119-OD2N1-H176.378.375.777.1
S145-OG-HGO659.461.359.657.6
A146-N-HO661.563.167.563.4
K147-N-HO6 84.486.584.982.3
a Hydrogen bonds are analyzed by an acceptor···donor distance of <3.5 Å and acceptor···H-donor angle of >120°. b Occupancy (%) is defined as the percentage of simulation time that a specific hydrogen bond exists.
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Chen, J.; Wang, J.; Yang, W.; Zhao, L.; Zhao, J.; Hu, G. Molecular Mechanism of Phosphorylation-Mediated Impacts on the Conformation Dynamics of GTP-Bound KRAS Probed by GaMD Trajectory-Based Deep Learning. Molecules 2024, 29, 2317. https://doi.org/10.3390/molecules29102317

AMA Style

Chen J, Wang J, Yang W, Zhao L, Zhao J, Hu G. Molecular Mechanism of Phosphorylation-Mediated Impacts on the Conformation Dynamics of GTP-Bound KRAS Probed by GaMD Trajectory-Based Deep Learning. Molecules. 2024; 29(10):2317. https://doi.org/10.3390/molecules29102317

Chicago/Turabian Style

Chen, Jianzhong, Jian Wang, Wanchun Yang, Lu Zhao, Juan Zhao, and Guodong Hu. 2024. "Molecular Mechanism of Phosphorylation-Mediated Impacts on the Conformation Dynamics of GTP-Bound KRAS Probed by GaMD Trajectory-Based Deep Learning" Molecules 29, no. 10: 2317. https://doi.org/10.3390/molecules29102317

APA Style

Chen, J., Wang, J., Yang, W., Zhao, L., Zhao, J., & Hu, G. (2024). Molecular Mechanism of Phosphorylation-Mediated Impacts on the Conformation Dynamics of GTP-Bound KRAS Probed by GaMD Trajectory-Based Deep Learning. Molecules, 29(10), 2317. https://doi.org/10.3390/molecules29102317

Article Metrics

Back to TopTop