Open Access Article
This Open Access Article is licensed under a
Creative Commons Attribution 3.0 Unported Licence

Q61 mutant-mediated dynamics changes of the GTP-KRAS complex probed by Gaussian accelerated molecular dynamics and free energy landscapes

Jianzhong Chen*a, Qingkai Zenga, Wei Wanga, Qingquan Hua and Huayin Bao*b
aSchool of Science, Shandong Jiaotong University, Jinan, 250357, China. E-mail: chenjianzhong1970@163.com; jzchen@sdjtu.edu.cn
bSchool of Pharmacy, Shandong University of Traditional Chinese Medicine, Jinan, 250355, China. E-mail: baohuayin@163.com

Received 28th October 2021 , Accepted 4th January 2022

First published on 11th January 2022


Abstract

Understanding the molecular mechanism of the GTP-KRAS binding is significant for improving the target roles of KRAS in cancer treatment. In this work, multiple replica Gaussian accelerated molecular dynamics (MR-GaMD) simulations were applied to decode the effect of Q61A, Q61H and Q61L on the activity of KRAS. Dynamics analyses based on MR-GaMD trajectory reveal that motion modes and dynamics behavior of the switch domain in KRAS are heavily affected by the three Q61 mutants. Information of free energy landscapes (FELs) shows that Q61A, Q61H and Q61L induce structural disorder of the switch domain and disturb the activity of KRAS. Analysis of the interaction network uncovers that the decrease in the stability of hydrogen bonding interactions (HBIs) of GTP with residues V29 and D30 induced by Q61A, Q61H and Q61L is responsible for the structural disorder of the switch-I and that in the occupancy of the hydrogen bond between GTP and residue G60 leads to the structural disorder of the switch-II. Thus, the high disorder of the switch domain caused by three current Q61 mutants produces a significant effect on binding of KRAS to its effectors. This work is expected to provide useful information for further understanding function and target roles of KRAS in anti-cancer drug development.


1 Introduction

RAS proteins HRAS, KRAS and NRAS are considered as significant family members of small guanosine triphosphatases (GTPases) and they can be functionally utilized as a molecular switch to regulate cell growth and differentiation in intracellular signaling pathways.1,2 The GTPase activating proteins (GAPs) accelerate the hydrolysis reaction of guanosine triohosphte (GTP) into guanosine diphosphate (GDP), which leads to an inactive state of RAS proteins.3 Guanosine exchange factors (GEFs) catalyze the transformation of GDP toward GTP and induce an active state of the GTP-bound RAS proteins.3 Multiple signaling pathways playing a key role in cell proliferation and survival can be activated through interactions of GTP-bound RAS proteins with effectors, namely RAF, phosphoinositide 3-kinase (PI3K) and Ral guanine nucleotide dissociation stimulator (RalGDS).4 The previous reports suggested that the active RAS proteins mediated by residue mutations promote tumorigenesis through hyperactivating the downstream signal pathways.4–8 Among mutations of three RAS isoforms, KRAS oncogenes account for ∼82% of RAS mutations relating with human cancer based on information from COSMIC.9 Therefore, KRAS has been an interesting focus on design of clinically available drugs toward human cancer treatment.

By now, a number of molecular structures of GTP- and GDP-bound KRAS have been determined by X-ray and NMR experiments from different work groups,10–15 which provide significant structural basis for understanding function of KRAS. According to function of structural domains, the conserved catalytic domain of KRAS is divided into two functional lobes, namely the effector lobe (residues 1–86) and the allosteric lobe (residues 87–166). The effector lobe is located at the N-terminal region of the catalytic domain, including P-loop (residue 10–17), switch-I (residues 30–40) and switch-II (residues 60–76), and these secondary structures form the main binding sites (Fig. 1A).16 The allosteric lobe is involved in a α-helical dimerization interface and binding motifs of GTP or GDP's nucleotide base and functions as a relaying of information through the protein.16–18 In the above mentioned structures, the switch-I and switch-II generate functional interactions with GEFs and GAPs or effectors of KRAS, indicating that the conformational transformation between switch-I and switch-II highly affects the activity of KRAS.19–24 Hence it is of high significance to deeply explore molecular mechanism underlying conformational changes of KRAS for understanding how to tune the activity of KRAS.


image file: d1ra07936k-f1.tif
Fig. 1 Molecular structures: (A) the GTP-bound KRAS shown by using the crystal structure (PDB ID: 6MNX), in which KRAS, GTP, mutated site and magnesium ion (Mg2+) are displayed in cartoon, stick, dot and ball modes, respectively, moreover the P-loop, switch-I and switch-II exhibited in green, cyan and orange, (B) Q61 and mutants A61, H61 and L61, (C) binding pocket of GTP depicted in surface styles and (D) GTP shown in stick patterns.

The previous works suggested that mutations significantly affect the conformational alterations of the switch-I and switch-II in KRAS and change binding activity of KRAS to its regulators and effectors.25–33 From residue substitutions at different sites, G12, G13 and Q61 mutations respectively account for 89%, 9% and 1% of the observed mutants of KRAS.34 In addition, the effect of mutations at the other sites on the GTP- and GDP-KRAS binding were also investigated by different work groups.29,35–37 Based on significant roles of the P-loop in function of KRAS, many experimental and theoretical works have focused on insights into impacts of G12 and G13 mutations on the conformational changes of KRAS and the results provide molecular mechanisms underlying the regulation of mutations from these two sites on the activity of KRAS.38–50 Although the number of works correlating with Q61 mutations in KRAS is less than that of G12 and G13 mutations, Q61 mutations also play a vital role in regulations on binding activity of KRAS to its effectors and regulators. For instance, the X-ray crystal structure of the GTP-associated Q61H determined by Zhou et al. unveiled that a hyperdynamic switch-II induces a more stable interaction with switch-I and verified that the enhanced RAF activity stems from not only the absent intrinsic GTP hydrolysis activity but also the improved affinity for RAF.51 The investigation from Yan et al. suggested that the GMPPNP-bound Q61L KRAS generates a stronger binding to NF1 than the GMPPNP-bound Q61R one and provided a new mechanistic insight about how SPRED1 interacts with neurofibromin.52 The work of Buhrman et al. on transformation efficiency of Q61 mutations clarified a molecular mechanism of intrinsic hydrolysis in the RAF/RAS complex and elucidated structural features in different Q61 mutants correlating with their potency of cell transformation.53 Despite these successful studies, dynamics information involved in effect of Q61 mutations on the activity of KRAS is insufficient. Therefore, it is highly requisite to further explore Q61 mutant-mediated conformational changes of KRAS for deeply understanding function of KRAS.

Up to now, multiple simulation methods, including conventional molecular dynamics (cMD),54–62 replica-exchange molecular dynamics (REMD),63,64 accelerated molecular dynamics (aMD)65–69 and Gaussian accelerated molecular dynamics (GaMD)70–78 etc., have been proposed to probe conformational changes of receptors because of ligand bindings and residue mutations. Compared to cMD, MR-GaMD simulations can obtain more rational conformational sampling on receptors.71,76 Recently, MD simulations were also applied to successfully investigate impacts of Q61 mutations on the activity of KRAS. The NMR S2 calculated by Kumar et al. with the DynaMine suggested that Q61H interferes with water molecule coordination and affects hydrolysis of GTP.24 The results from MD simulations of Vatansever et al. indicated that Q61H on the switch-II does not generate obvious influences on the fluctuations of the switch-II in active KRAS but evidently affects the flexibility of inactive protein.40 Lu et al. performed MD simulations on KRAS of multiple mutations and their studies verified that Q61H is more prone to alter the KRAS4B–GTP conformation to the active state than the G12C mutation.34 Mutations at Q61 play vital roles in regulations on bindings of KRAS to the effectors and regulators. Thus, it is still high significance to further explore dynamics changes of KRAS by using rational conformational sampling technology for efficiently understanding the role of KRAS in drug targets toward cancer treatments.

In this work, three mutations Q61A, Q61H and Q61L at the GTP-bound KRAS, observed at the experimental works,10,53,79 were selected to decipher influences of Q61 mutants on binding activity of KRAS. Q61 is a residue without charge and located at the switch-II. The sidechain of Q61 is replaced by a small hydrophobic alkyl, a big hydrophobic alkyl and a hydrophobic ring in Q61A, Q61L and Q61H, respectively (Fig. 1B). Binding pocket and molecular structure of GTP were separately displayed in Fig. 1C and D. Insights into impacts of the changes in sizes of hydrophobic sidechains of Q61 mutants on conformational alterations of KRAS are helpful for understanding molecular mechanism for tuning binding activity of KRAS. To reach our current aims, multiple replica GaMD (MR-GaMD) simulations were employed to enhance conformational sampling of KRAS. Principal component analysis (PCA), construction of free energy landscapes (FELs) and interaction network analysis were performed to explore Q61 mutant-mediated regulation on conformational changes of KRAS. This work is also expected to provide useful information for deeply understanding activity and function of KRAS.

2 Methods

2.1 Modeling simulated systems

The crystal structure 6MNX was extracted from protein data bank (PDB) to obtain the initial structure of the GTP-bound Q61H KRAS.51 Two missing residues Y64 and S65 in the A section of 6MNX were repaired with the program MODELLER.80 To keep atomic coordinate consistency, residue H61 in the repaired structure was separately changed into Q61, A61 and L61 to obtain the initial structures of the GTP-bound WT, Q61A and Q61L KRAS. A Web server H++ 3.0 (http://biophysics.cs.vt.edu/H++) was used to check and assign rational protonated states to residues of KRAS.81 All crystal water molecules and magnesium ion (Mg2+) were left at the initial structures of four simulated systems. The Leap module in Amber 20 (ref. 82 and 83) was utilized to obtain force field parameters of four simulated systems based on the following procedure: (1) forming the chemical bonds between heavy atoms and hydrogen atoms at the crystal structures, (2) producing force field parameters of GTP by aid of the work of Meagher et al.,84 (3) generating force field parameters of the WT and mutated KRAS with ff14SB force field,85 (4) building an octahedral water box of 12 Å buffer with the TIP3P model86 to solve the GTP-KRAS complexes and assigning force field parameters to water molecules and (5) neutralizing the simulated systems with the appropriate sodium ions (Na+) in a salt environment of 0.15 M NaCl and yielding force field parameters of Na+ and Mg2+ by using the Aqvist force field.87

2.2 MR-GaMD simulation

The program pmemd.cuda inlayed at Amber 20 was employed to carry out all the following simulations.88 3000-step steepest descent minimization, 5000-step conjugate gradient minimization, 2 ns heating process of 0 to 310 K at the NVT condition and 2 ns equilibrium process of 310 K at the NPT condition were sequentially run to deeply relax every simulated system. Then, the relaxed structure was utilized to execute 200 ns cMD simulation on each simulated system. Two structures randomly extracted from the previous cMD simulation were used as initial conformation to restart two new 200 ns cMD simulations. Three ending structures from the above cMD simulations were adopted to perform three independent MR-GaMD simulations and each simulation was run for 1.6 μs.

GaMD simulations can efficiently leap over energy barrier of simulated systems and enhance conformational samplings with the following potential function described by eqn (1) and (2)

 
image file: d1ra07936k-t1.tif(1)
 
image file: d1ra07936k-t2.tif(2)
where image file: d1ra07936k-t3.tif represents potential function of simulated systems and image file: d1ra07936k-t4.tif is revised as image file: d1ra07936k-t5.tif when image file: d1ra07936k-t6.tif is lower than a threshold energy E. Two parameters E and k can be updated based on the following principals
 
image file: d1ra07936k-t7.tif(3)
 
image file: d1ra07936k-t8.tif(4)
in which when E is lower than an energy E = Emax, k0 is obtained from the eqn (5) as below
 
image file: d1ra07936k-t9.tif(5)
but when E is set as an upper bound image file: d1ra07936k-t10.tif, k0 is got from the eqn (6)
 
image file: d1ra07936k-t11.tif(6)
where Vmax, Vmin and Vavg correspond to the minimum, maximum, and averaged potential energy of each system, respectively. The parameter σ represents standard deviation of potential energies and the parameter σ0 is artificially set as an upper limit to accurately reweight the original energy profile of simulated systems. In our work, 4.8 μs MR-GaMD simulations, composed of three-replica GaMD simulations of 1.6 μs, were run on the GTP-bound WT, Q61A, Q61H and Q61L KRAS to rationally sample conformations of KRAS. Three-replica GaMD trajectories were connected into a single integrated MR-GaMD trajectory to facilitate the post-processing analysis. A program PyReweighting89 provided by Miao was employed to decode energy profile of KRAS and probe conformational changes of KRAS. Through all simulations involved in this work, chemical bonds linking hydrogen with heavy atoms were restrained by utilizing the SHAKE algorithm.90 The Langevin dynamics,91 assigned as a collision frequency of 2.0 ps−1, was used for regulation on temperatures of four simulated systems. The particle mesh Ewald (PME) was wielded to estimate electrostatic interactions between atoms and the cutoff value for treating non-bonded interactions was set as 12 Å.

2.3 Principal component analysis

To recognize correlated motions relating with functional significance, PCA was used in our current work.92–95 The first step of PCA is to build a covariance matrix C by using eqn (7)
 
C = 〈(qi−〈qi〉)(qj−〈qj〉)T (7)
where qi and qj are the Cartesian coordinates of the ith and jth Cα atoms in KRAS while 〈qi〉 and 〈qj〉 indicate their averaged positions. The average is estimated by superimposing the single integrated MR-GaMD trajectory with a defined reference structure to remove overall translations and rotations with a least-square fit procedure.96 The second step of PCA is to diagonalize the symmetric matrix C into a diagonal one A with an orthogonal coordinate transformation matrix T based on the following equation
 
A = TTCijT (8)
from which the diagonal elements represent the eigenvalues λi and the columns signify the eigenvectors characterizing the motion direction of movement relative to 〈qi〉. The third step of PCA is to describe a single correlated displacement of a group of atoms in a multidimensional space using the eigenvector and explain the amplitude of the motion along an eigenvector, which can rationally reflect conformational changes of the domain in KRAS. In this work, PCA was carried out by utilizing the CPPTRAJ module in Amber 20.97 The software VMD98 was employed to visualize the data from PCA, plot pictures and unveil impacts of Q61 mutants on conformational alterations of KRAS.

2.4 Cross-correlation analysis

To check the extent of correlated movements of structural domain in KRAS, cross-correlation matrix, Cij, was calculated by using the x, y and z coordinates of the Cα atoms in KRAS based on the following equation99
 
Cij=(Δri × Δrj)/(〈Δri2〉〈Δrj2〉)2 (9)
where Δri and Δrj are the displacement vectors of atoms i and j relative to their averaged positions while the angle brackets indicate ensemble averages over the structures recorded at the single integrated MR-GaMD trajectory. The element values of cross-correlation matrix are located at a range from −1 to 1. The Cij of positive values describe the positively correlated motions between atoms i and j while the Cij of negative values reflect the anti-correlated movements of atoms i relative to j. The color-coded patterns were adopted to characterize the extent of correlated motions. In our current study, the module CPPTRAJ in Amber 20 was wielded to execute cross-correlation analysis.

2.5 Analysis of interaction network

Identification of interaction network plays an important role in probing ligand-receptor binding mechanisms. Hydrogen bonding interactions (HBIs) and hydrophobic interactions are two favorable factors for bindings of ligands to receptors. In this work, the module CPPTRAJ inlayed in Amber 20 was applied to analyze HBIs between GTP and KRAS. A protein−ligand interaction profiler (PLIP) server100 was employed to recognize interaction networks and types of GTP with KRAS. Meanwhile the software PyMOL (https://www.pymol.org.) was utilized to depict hot interaction spots of GTP with KRAS.

3 Results and discussion

3.1 Structural fluctuations and flexibilities of KRAS

Root-mean-square deviations (RMSDs) were computed by using the coordinates of backbone atoms from KRAS recorded at the single integrated MR-GaMD trajectory to understand changes in structural fluctuations of KRAS due to Q61 mutations. Frequency distributions of RMSDs were plotted in Fig. 2A. RMSDs of the GTP-bound Q61A, Q61H and Q61L KRAS are increased by 0.22, 0.39 and 0.10 Å compared to that of the GTP-bound WT KRAS, respectively, moreover the distribution shape of RMSDs of three mutated KRAS moves toward the right. Therefore, three mutations Q61A, Q61H and Q61L strengthen structural fluctuations of KRAS relative to the WT one and the changes in structural stability certainly affect interactions of KRAS with GEFs and GAPs or effectors.
image file: d1ra07936k-f2.tif
Fig. 2 Structural fluctuations and flexibilities: (A) frequency distributions of backbone RMSDs, (B) RMSF difference of Cα atoms between the WT and mutated KRAS, (C) structural domains corresponding to obvious changes of RMSF values and (D) frequency distributions of molecular surface areas.

Difference in root-mean-square fluctuations (RMSFs) of the Cα atoms between the mutated and WT KRAS was calculated by using the equation ΔRMSF = RMSFmutatedRMSFWT to explore impacts of Q61 mutants on structural flexibility of KRAS (Fig. 2B). As a result, Q61A, Q61H and Q61L yield obvious effect on the switch-I and switch-II of KRAS. According to Fig. 2C, the switch-I consists of two parts (S1-L and S1-β) while the-switch II is composed of two parts (S2-L and S2-α2). Q61A evidently enhances the structural flexibility of the switch-I relative to the WT KRAS. By comparison with the WT KRAS, Q61H increases the structural flexibility of S1-L in the switch-I but decreases that of S1-β, on the contrary, Q61L weakens the flexibility of S1-L but strengthens that of S1-β (Fig. 2B and C). Similar to the switch-I, Q61A highly strengthens the flexibility of the switch-II relative to the WT KRAS. Compared to the WT KRAS, Q61H slightly reduces the flexibility of S2-L in the switch-II but evidently enhances that of S2-α2, differently, Q61L strengthens the flexibility of S2-L but softly weakens that of S2-α2 (Fig. 2B and C). Meanwhile, Q61A, Q61H and Q61L increase the structural flexibility of a loop L2 linking α3 with a β-strand relative to the WT KRAS. In addition, Q61A evidently enhances the structural flexibility of P-loop and α3-L1 compared to the WT KRAS. According to structural information, the switch-I and switch-II are involved in functional interactions of KRAS with GEFs and GAPs or effectors,21 thus the alterations in the flexibility of these two switches necessarily affect the activity of KRAS, which is also supported by the previous works.34,40

Molecular surface areas (MSAs) of KRAS were also computed to probe influences of Q61 mutants on total flexibility of KRAS and frequency distribution of MSAs was displayed in Fig. 2D. The MSAs of the GTP-bound Q61A, Q61H and Q61L KRAS are increased by ∼72, 150 and 221 Å2 by comparison with that of the WT KRAS, respectively, moreover the distribution shape of MSAs in three mutated KRAS move toward the right. Hence it is concluded that Q61A, Q61H and Q61L enhance the total flexibility of KRAS relative to that of the WT KRAS.

Based on the aforementioned information, Q61A, Q61H and Q61L obviously increase the structural fluctuation of KRAS by referencing the WT state. Meanwhile, three mutants Q61A, Q61H and Q61L not only evidently alter the local flexibility of the switch-I, switch-II and loop L2 in KRAS but also strengthen the total flexibility of KRAS. These structural changes inevitably generate significant influences on binding of KRAS to GEFs and GAPs or effectors.

3.2 Changes in dynamics behavior of KRAS affected by Q61 mutations

Cross-correlation analysis was executed to uncover effect of Q61 mutations on correlated motions between residues of KRAS and the results were visually display in ESI Fig. S1. Color-coded styles were adopted to reflect the extent of correlated motions. It is observed that Q61A, Q61H and Q61L produce different influences on correlated movements of KRAS. In the WT state, the switch-II yields obvious anticorrelated motions (blue) relative to the P-loop but the switch-I hardly produces correlated movements relative to the P-loop (Fig. S1A). In the meantime, a slight anticorrelated motion (light blue) also appear between the switch-II and switch-I. In addition, the structure α3-L1 not only generates a strong positive correlated motion (yellow) relative to the P-loop but also produces a weak anticorrelated movement relative to the switch-II (Fig. S1A and 2C).

By comparison with the WT KRAS, Q61A enhances the anticorrelated motion between the switch-II and the P-loop but Q61H and Q61L highly weaken this anticorrelated motion (Fig. S1B–D). Compared to the WT state of KRAS, Q61A obviously strengthens the anticorrelated motions (dark blue) of the switch-I relative to the P-loop and the α3-L1 relative to the switch-II, however Q61H and Q61L hardly affect these two anticorrelated motions. Besides, Q61A extremely enhances the positive correlated movement (yellow and red) between the α3-L1 and the P-loop, but Q61H and Q61L do not bring obvious effect on this positive correlated motion (Fig. S1B–D and 2C).

In general, eigenvalues obtained from PCA are utilized to reflect total motion intensity of proteins. For this study, function of eigenvalues vs. eigenvector indexes were depicted in Fig. S2. In fact, the first several bigger eigenvalues describe mainly concerted motions of structural domain in KRAS. As shown in Fig. S2, the first five eigenvalues account for 60.1, 86.5, 82.8 and 83.3% of the total movements of the GTP-associated WT, Q61A, Q61H and Q61L KRAS, respectively. By referencing the WT KRAS, the first eigenvalue of three mutated KRAS is extremely improved compared to that of the WT KRAS. Thus Q61A, Q61H and Q61L highly enhance the total motion intensity of KRAS relative to the WT one.

Through the above mentioned eigenvalues, the first eigenvector embodies main motion behavior of KRAS. Because of so, the first eigenvector from PCA is visually exhibited in Fig. 3 by aid of the software VMD and the optimized structure. It is found that the switch-I and switch-II show stronger movement by comparison with the other structural domain of KRAS. In the GTP-associated WT KRAS, two switches move toward a completely opposite direction and are close to each other, which results in a more compact switch domain, furthermore the switch-II displays the strongest movement among the structural domain (Fig. 3A). Compared to the WT state of KRAS, Q61A not only obviously enhances the motions of the P-loop and switch-II but also induces a tendency of the P-loop, switch-I and switch-II to move outwards, which leads to a more incompact binding pocket of GTP (Fig. 3B). By referencing the WT KRAS, Q61H changes the motion direction of the switch-II and strengthens its motion intensity, moreover makes the switch-I move toward the left, which also causes an untighter binding pocket of GTP (Fig. 3C). In the case of the GTP-bound Q61L KRAS, two switches move toward an opposite direction and are far away from each other, which forms an completely opening tendency of the switch domain (Fig. 3D).


image file: d1ra07936k-f3.tif
Fig. 3 Concerted motions of structural domains in KRAS revealed by the first eigenvector from principal component analysis: (A) GTP-associated WT KRAS, (B) GTP-associated Q61A KRAS, (C) GTP-associated Q61H KRAS and (D) GTP-associated Q61L KRAS.

Based on the above analyses, Q61A, Q61H and Q61L not only alter motion modes of the switch-II relative to the P-loop and switch-I but also obviously change motion intensity and dynamics behavior of the switch domain of KRAS, which yields significant influences on binding activity of KRAS. The Pantsar's work verified that conformational changes of the switch domain heavily affect the activity of KRAS.16 The μs-scaling MD simulation from Sayyed-Ahmad et al. indicated that mutation-mediated differences in dynamics and interaction networks can disturb binding of KRAS to GEFs and GAPs or effectors,101 which basically agrees with our current findings.

3.3 Difference in free energy landscapes induced by mutations

To understand changes in free energies of KRAS induced by Q61 mutants, the RMSDs of backbone atoms and the distance of Y32 away from E63 were used as reaction coordinates to construct FELs of KRAS. RMSDs can efficiently reflect structural fluctuation of KRAS through the entire MR-GaMD simulation. Residues Y32 and E63 are situated at the switch-I and switch-II, respectively, and the alterations in the distance between them can rationally embody conformational transformation of the switch domain of KRAS. Therefore, we selected RMSDs and the distance between Y32 and E63 as reaction coordinates to build FELs. The constructed FELs and the extracted structural information from MR-GaMD simulations were displayed at Fig. 4–6 and S3–S6.
image file: d1ra07936k-f4.tif
Fig. 4 Free energy landscapes and structures of the GTP-bound Q61A KRAS: (A) showing free energy landscape constructed by using RMSD of backbone atoms and the distance of Y32 away from E63 as reaction coordinates, in which E1, E2, E3 and E4 indicate energy valleys detected by MR-GaMD simulations, (B), (C), (D) and (E) corresponding to representative structures located at energy valleys E1, E2, E3 and E4, respectively, and (F) displaying superimpositions of GTP, residue Y32 and ions in four representative structures with that in the crystal structure (6MNX). The unit of the distances is Å and that of free energy is kcal mol−1.

image file: d1ra07936k-f5.tif
Fig. 5 Free energy landscape and structures: (A) depicting free energy landscape built by utilizing RMSDs of backbone atoms and the distance of Y32 away from E63, (B) and (C) indicating the representative structures situated at energy valleys E1 and E2, and (D) suggesting superimposition of residue Y32, GTP and ions in two representative structures with that in the crystal structure (6MNX). The unit of the distances is Å and that of free energy is kcal mol−1.

image file: d1ra07936k-f6.tif
Fig. 6 Free energy landscape and structures of the GTP-associated Q61L KRAS: (A) indicating free energy landscape built by employing RMSDs of backbone atoms and the distance of Y32 away from E63 as reaction coordinates, in which E1, E2 and E3 signify energy valleys identified by MR-GaMD simulations, (B), (C) and (D) displaying representative structures falling into energy valleys E1, E2 and E3, separately. The unit of the distances is Å and that of free energy is kcal mol−1.

For the GTP-bound WT KRAS, only an energy valley (E1) is detected by MR-GaMD simulations (Fig. S3A). This result indicates that the conformations of the GTP-bound WT KRAS are mostly populated at an energetic space. The distance of Y32 away from E63 in the crystal structure (6MNX) is 19.3 Å and that in the structure E1 is ∼18.8 Å, moreover the distances of G13 away from E63 and T35 shows a small change between the crystal structure and the structure E1 (Fig. S3B and C). The alignment between these two structures shows that the most regions of KRAS, GTP, magnesium ion Mg2+ and an important residue Y32 are stable in the GTP-bound WT state, only the switch-II has an obvious deviation.

In the case of the GTP-bound Q61A KRAS, MR-GaMD simulations capture four different energy valleys E1, E2, E3 and E4 (Fig. 4A), indicating that the conformations of the GTP-bound Q61A KRAS are mainly distributed at four energetic spaces. The distances of Y32 away from E63 in four energetic states E1, E2, E3 and E4 are 12.9, 16.5, 21.1 and 14.1 Å (Fig. 4B–E), respectively, which suggests that the switch domain of KRAS in energetic state E1 is the tightest while that of KRAS in energetic state E3 is the most incompact. Meanwhile the distances of G13 away from E63 and Y32 in four different energetic states also yield different changes, which also reflects alterations in sizes of binding pockets. Superimposition of four representative structures located at energy valleys E1–E4 with the crystal structure (6MNX) exhibits that the switch-I and switch-II are highly disordered (Fig. S4), implying that Q61A produces significant effect on structural flexibility of the switch domain. More interestingly, a sodium ion Na+ appears at the structures E1 and E2 and two sodium ions Na+ are detected in the structure E4, which possibly changes electrostatic environment around the binding pocket. Superimposition of GTP, residue Y32 and ions in four representative structures with that in the crystal structure (6MNX) is depicted in Fig. 4F. The result suggests that Q61A hardly affects stability of GTP and Mg2+ through the entire MR-GaMD simulations but evidently alters conformations of Y32 in four energetic states by comparison with the crystal structure. Besides, Fig. 4F also shows that sodium ions (Na+) appear at different sites near the phosphate group of GTP and change electrostatic environment around GTP, but these sodium ions are unstable. The work of Xu et al. also indicated that conformational arrangement of the switch domain induced by Q61A also disturbs binding of KRAS to effectors, which is in good agreement with our current findings.10

With regard to the GTP-bound Q61H KRAS, MR-GaMD simulations identify two energetic valleys E1 and E2 and this result signifies that the conformations of the GTP-bound Q61H KRAS are primarily populated at two energetic spaces (Fig. 5A). In the energetic state E1, although the distance between Y32 and E63 of the Q61H KRAS is 19.1 Å and similar to that of the crystal structure (6MNX), but the distances of G13 away from E63 and Y32 in the Q61H KRAS have obvious difference compared to that in the crystal structure (Fig. 5B and S4B). The distance of Y32 away from E63 in the energetic state E2 is ∼16.6 Å, which results in a more compact switch domain by referencing the GTP-bound WT KRAS (Fig. 5C). Alignment of the structures E1 and E2 with the crystal structure exhibits that Q61H generates evident influences on conformational alterations of the switch domain, especially for the switch-II (Fig. S5). According to superimposition of GTP, ions and Y32 in the structures E1 and E2 of the GTP-bound Q61H KRAS with that in the WT one (Fig. 5D), Q61H hardly disturbs stability of GTP and Mg2+ but obviously alters conformations of Y32. Besides although sodium ions (Na+) appear near the phosphate group of GTP in the Q61H KRAS, these sodium ions have different sites in two structures E1 and E2 (Fig. 5D), hence these sodium ions certainly change electrostatic environment around GTP. Compared to the GTP-bound Q61A KRAS, the GTP-associated Q61H one has more stable conformations. The study from Zhou et al. unveiled that a hyperdynamic switch-II allows for a more stable interaction with the switch-I in the GTP-associated Q61H KRAS and their results also verified that an enhanced RAF activity stems from not only absent intrinsic GTP hydrolysis activity but also increased affinity for RAF.51 MD simulations performed by Lu et al. suggested that conformational transformation of the switch domain because of Q61H affects intrinsic GTPase activity.34 The Vatansever's work found that Q61H leads to distinctly strong deviations in the switch-II and decreases the number of conformations of the switch region in the inactive KRAS.40 These previous findings basically support our current work.

As for the GTP-bound Q61L KRAS, three energetic valleys E1, E2 and E3 are recognized by MR-GaMD simulations (Fig. 6A), demonstrating that conformations of the GTP-bound Q61L KRAS are primarily distributed at three energetic spaces. As shown in Fig. 6B–D, the distances between Y32 and E63 in three structures E1, E2 and E3 are 17.7, 20.1 and 16.4 Å, separately, thus the structures E1 and E3 form a tighter switch domain while the structure E2 has a more incompact one compared to the WT KRAS. Three structures E1, E2 and E3 are aligned together with crystal structure (Fig. S6A) and the result suggests that Q61L exerts an apparent effect on the conformation of the switch domain of KRAS but barely affects stability of GTP and Mg2+. Meanwhile a sodium Na+ appears at the structure E3 (Fig. 6D), which possibly brings an impact on electrostatic environment near GTP. Superimposition of GTP, ions and Y32 in the structures E1–E3 with that in crystal structure verifies that conformations of Y32 are heavily affected by Q61L but GTP and Mg2+ are highly stable through the entire MR-GaMD simulations (Fig. S6B). The crystal structure of the GTP-bound Q61L KRAS determined by Hunter et al. suggested that the switch-I and switch-II are disordered and KRAS exhibits a low level of intrinsic activity,79 which agrees well with the information revealed by our current study.

The aforementioned details reveal that Q61A, Q61H and Q61L highly affect structural stability of the switch domain and make the switch domain locate at a highly disordered state. In fact, the switch-I and switch-II take part in functional interactions of KRAS with GEFs and GAPs or effectors,21 hence high disorder of the switch domain of course affects the activity of KRAS. Multiple previous studies uncovered that conformational changes of the switch domain induced by mutations produce key influences on the intrinsic GTPase activity,53,102–104 which is in basic consistence with our results. Our work also finds that Q61A, Q61H and Q61L yield significant effect on conformations of Y32, which tunes the intrinsic activity of KRAS. It has also been reported that residue Y32 situated at the switch-I plays a significant role in intrinsic hydrolysis of KRAS,26,79,105,106 thus changes in the orientation of the sidechain in Y32 directly disturb the activity of KRAS. More interestingly, our work verifies that Q61A, Q61H and Q61L induces a possibility of sodium ions (Na+) appearing in different energetic structures, which possibly yields an important impact on electrostatic environment around binding sites. This interesting phenomenon is also observed at the previous studies,30,46 supporting our current findings.

3.4 Alterations in interaction network of GTP with KRAS caused by Q61 mutants

Identification of ligand–receptor interaction network is highly helpful for understanding molecular mechanism on binding of ligands to receptors. The PLIP sever was applied to analyze the interaction network of GTP with KRAS and the results were displayed at Fig. 7, 8, S7 and S8. HBIs of GTP with KRAS were also explored by using the CPPTRAJ module in Amber 20 and the information of HBIs were given in Table 1. HBIs are determined through a standard of an acceptor⋯donor distance of <3.5 Å and acceptor⋯H-donor angle of >120°. The occupancy of hydrogen bonds listed in Table 1 reflects stability of HBIs though the entire MR-GaMD simulations.
Table 1 Hydrogen bonding interactions of GTP with the WT and mutated KRAS
Hydrogen bonds aOccupancy (%)
a Occupancy (%) is defined as the percentage of simulation time that a specific hydrogen bond exists.
Residue GTP WT Q61A Q61H Q61L
G13–N–H O3B 93.1 44.2 51.2 91.3
V14–N–H O1B 25.9 28.2 25.1 26.3
G15–N–H O1B 98.1 98.3 97.8 97.9
K16–N–H O1B 99.9 99.9 99.9 99.9
S17–N–H O2B 99.6 98.7 99.6 99.6
A18–N–H O1A 99.7 99.6 99.6 99.7
V29–N–H O2′-HO′2 65.8 45.3 40.1 36.8
D30–O O3′-H3T 54.2 30.2 34.1 33.6
T35–N–H O1G 3.5 4.2 29.2 15.3
G60–N–H O3G 92.1 56.7 78.7 72.7
N116–ND2–HD2 N7 82.6 83.1 80.8 80.6
K117–N–H O4′ 40.1 45.4 44.8 44.7
S145–OG–HG N1 66.9 66.2 37.8 43.4
A146–N–H O6 75.6 72.3 73.0 76.6
K147–N–H O6 73.1 74.1 75.3 71.0


The structure E1 extracted from the single integrated MR-GaMD trajectory of the GTP-bound WT KRAS was used to analyze interaction network with the PLLP sever (Fig. 7). It is found from Fig. 7A that GTP yields HBIs with six residues G13, V14, G15, K16, S17 and A18 situated at the P-loop and except for V14, the occupancy of HBIs are higher than 93.1% (Table 1), indicating that HBIs of GTP with G13, G15, K16, S17 and A18 are extremely stable. According to Table 1, GTP forms three hydrogen bonds with residues V29, D30 and T35 located at the switch-I and their occupancy is lower than 65.8% and the hydrogen bond between GTP and V29 does not appear at the structure E1 (Fig. 7A). It is also observed that GTP produces a HBI with residue G60 in the switch-II and its occupancy is 92.1%, implying that this hydrogen bond is highly stable through the entire MR-GaMD simulations in the GTP-bound WT state of KRAS (Table 1 and Fig. 7A). In addition, GTP generate HBIs with residues N116, S145, A146 and K147 with the occupancy higher than 66.9% (Table 1 and Fig. 7A), showing that these hydrogen bonds have a higher stability. Although GTP also forms a hydrogen bond with K117, its occupancy is only 21.1%, hence this hydrogen bond shows a weak stability during MR-GaMD simulations. Based on Fig. 7A, the phosphate group of GTP produces two salt bridge interactions with K16 and the adenine group of GTP forms a salt bridge interaction with D119. These three salt bridges play a significant role in stabilization of binding of GTP to KRAS. Meanwhile, the PLIP server also detects electrostatic interactions of magnesium ion Mg2+ with GTP, residues S17 and T35 (Fig. 7B), thus the appearance of Mg2+ between GTP and S17 is helpful for maintaining stability of the P-loop conformation.


image file: d1ra07936k-f7.tif
Fig. 7 Interaction network in the GTP-bound WT KRAS: (A) interactions of GTP with key residues, including hydrogen bonds and salt bridge and (B) electrostatic interactions of magnesium ion Mg2+ with residues and GTP. The unit of the distances is Å.

For the GTP-bound Q61A, Q61H and Q61L KRAS, the most compact and incompact structures taken from the single integrated MR-GaMD trajectory were adopted to analyze the interaction network of GTP with three mutated KRAS using the PLIP server (Fig. 8, S7 and S8). It is noted that Q61A, Q61H and Q61L KRAS share almost same binding spots of GTP as the WT KRAS. As a common feature, three salt bridge interactions of GTP with K16 and D119 can be identified in all energetic states of the GTP-bound WT and mutated KRAS (Fig. 7A, 8A, B, S7A, B, S8A and B). To understand effect of Q61A, Q61H and Q61L on these salt bridge interactions, the distances of the positive charge in K16 away from the phosphorus atoms PB and PG of GTP and that of the negative charge in D119 away from the adenine group of GTP were computed by means of the single integrated MR-GaMD trajectory (Fig. 9A–C). The results exhibit that Q61A, Q61H and Q61L hardly produce obvious influences on the salt bridge between D119 and the adenine group of GTP (Fig. 9C), but Q61A slightly strengthens the salt bridge interaction between K16 and the phosphorus atom PB of GTP compared to three other states and Q61L enhances the salt bridge interaction between K16 and the phosphorus atom PG of GTP by referencing three other states (Fig. 9A and B). A π–π interaction between the phenyl group of F28 and the adenine group of GTP was recognized in the compact and incompact states of the GTP-associated Q61A KRAS (Fig. 8A and B) ant the incompact states of the GTP-bound Q61H and Q61L KRAS (Fig. S7B and S8B). To understand stability of this π–π interaction at four simulated systems, the distance between the phenyl group of F28 and the adenine group of GTP was estimated and the results were depicted in Fig. 9D. As a result, the distance of this π–π interaction is distributed at 5.1 Å in four simulated systems and the currently studied Q61 mutants hardly generate evident impacts on stability of this π–π interaction. Although the π–π interaction between F28 and GTP does not appear at the structure E1 extracted from MR-GaMD simulations of the GTP-bound WT KRAS (Fig. 7A), this π–π interaction is still stably maintained through MR-GaMD simulation (Fig. 9D). In addition, a cation–π interaction of K117 with the adenine of GTP is recognized at the incompact state of the GTP-bound Q61H and Q61L KRAS (Fig. S7B and S8B).


image file: d1ra07936k-f8.tif
Fig. 8 Interaction network in the GTP-bound Q61A KRAS: (A) interactions of GTP with residues in compact state of the switch domain, (B) interactions between GTP and residues in uncompact state of the switch domain, (C) interactions of magnesium ion Mg2+ and sodium ion Na+ with GTP and residues and (D) interactions of magnesium ion Mg2+ with GTP and residues.

image file: d1ra07936k-f9.tif
Fig. 9 Key interactions between GTP and several residues in the WT and mutated KRAS: (A) salt bridge interaction between K16 and phosphorus atom PB of GTP, (B) salt bridge interactions of K16 with phosphorus atom PG of GTP, (C) salt bridge interaction between D119 and GTP and (D) π–π stacking interaction of F28 with GTP.

As shown in Fig. 7B, 8C, D, S7C, D, S8C and D, the PLIP server also detects electrostatic interactions of magnesium ion Mg2+ with the oxygen atoms OG of S17, OG1 of T35, O2B and O2G of GTP, but the interaction of Mg2+ with T35 does not appear in the compact state of the GTP-bound Q61A KRAS (Fig. 8C). The distances of these oxygen atoms away from Mg2+ were calculated and their frequency distribution were plotted at Fig. 10. The information demonstrates that Q61A, Q61H and Q61L do not affect stability of Mg2+ through the entire MR-GaMD simulations. Thus electrostatic interactions of Mg2+ with GTP and S17 is helpful for stabilizing the interaction of GTP with the P-loop. Besides, sodium ions Na+ appear near the phosphate group of GTP at different states of the GTP-associated Q61A, Q61H and Q61L KRAS (Fig. 8C, S7C, D and S8C), and these sodium ions (Na+) interact with different oxygen atoms in the phosphate group of GTP, which changes electrostatic environment around GTP and affects conformational changes of KRAS.


image file: d1ra07936k-f10.tif
Fig. 10 Interactions of magnesium ion Mg2+ with GTP and residue: (A) interaction between Mg2+ and oxygen atom O2B of GTP, (B) interaction of Mg2+ with oxygen atom O2G of GTP, (C) interaction between Mg2+ and oxygen atom OG of S17 and (D) interaction of Mg2+ with oxygen atom OG1 of GTP.

The distinct difference in interaction network of GTP with the WT, Q61A, Q61H and Q61L KRAS primarily stems from HBIs of GTP with residues V29 and D30 in the switch-I. Compared to the GTP-associated WT KRAS, the occupancy of the hydrogen bond between GTP and V29 is decreased by 20.5, 25.8 and 29.0% because of Q61A, Q61H and Q61L, respectively, while that of the hydrogen bond between GTP and D30 is reduced by 24.0, 20.1 and 20.6%, separately (Table 1). Therefore the decrease in stability of hydrogen bonds between GTP and residues V29 and D30 is responsible for the high disorder of the switch-I observed from the previous FEL analyses. The occupancy of the hydrogen bond between GTP and G60 is reduced by 35.4, 13.4 and 19.4% due to Q61A, Q61H and Q61L by comparison with the WT KRAS, respectively, which plays a key role in the high disorder of the switch-II seen from the previous FEL analyses (Table 1). Apart from effect of Q61A and Q61H on the HBI of GTP with G13, three Q61 mutants hardly disturb stability of HBIs of GTP with V14, G15, K16, S17 and A18 in the P-loop (Table 1), which rationally explains the high order of the P-loop observed from previous structural superimposition. In addition, Q61A, Q61H and Q61L do not yield evident impacts on stability of HBIs of GTP with N116, K117, S145, A146 and K147 located at the non-switch regions compared to the WT KRAS (Table 1) and this result verifies a key factor maintaining structural order of the non-switch domain.

According to the aforementioned results, the following conclusions can be drawn: (1) salt bridge interactions of K116 and D119 with GTP and electrostatic interactions of Mg2+ with GTP, S17 and T35 play a key role in maintaining stability of the GTP-KRAS binding, (2) the decrease in stability of HBIs of GTP with V29 and D30 induced by Q61A, Q61H and Q61L is responsible for structural disorder of the switch-I and that in the occupancy of the hydrogen bond between GTP and G60 leads to structural disorder of the switch-II, and (3) the appearance of sodium ions near the phosphate group of GTP in the mutated KRAS alters electrostatic environment around GTP and affects the activity of KRAS.

4 Conclusions

KRAS has been regarded as an important drug target for cancer treatment and its conformational changes caused by ligand bindings and mutations are requisite for drug design. In this work, 4.8-μs MR-GaMD simulations were performed on the GTP-bound WT, Q61A, Q61H and Q61L KRAS to explore effect of Q61 mutants on the activity of KRAS. The calculated RMSDs of backbone atoms show that Q61A, Q61H and Q61L enhance total structural fluctuations of KRAS and the estimated RMSF difference of the Cα atoms between the WT and mutated KRAS suggest that Q61A, Q61H and Q61L obviously alter structural flexibility of the switch domain in KRAS, which affects the activity of KRAS. The information stemming from cross-correlation analysis and PCA verifies that three Q61 mutants affect motion modes of the switch domain and evidently change dynamics behavior of the switch domain in KRAS, which in turn regulates binding of KRAS to effectors.

FELs were constructed by utilizing RMSDs of backbone atoms and the distance of Y32 away from E63 as reaction coordinates. The results display that Q61A, Q61H and Q61L lead to high disorder of the switch domain in KRAS and reveal energetic basis with regard to conformational changes of KRAS. Moreover, MR-GaMD simulations detect the appearance of sodium ions (Na+) around the phosphate group of GTP in different structures falling into energy valleys, which changes electrostatic environment around GTP and certainly affects conformational changes of KRAS.

The analysis of interaction network reveals that the decrease in stability of HBIs of GTP with the switch-I and switch-II is responsible for the structural disorder of the switch domain. Dynamic analysis suggests that the flexibility change of the loop L2 is involved in the regulation on allosteric sites of KRAS and affects the alteration in allosteric binding pocket, which provides a hint for design and development of inhibitors tuning the activity of KRAS in the future.

Conflicts of interest

The authors declare that they have no known competing financial interests.

Acknowledgements

This work is supported by the National Natural Science Foundation of China (11504206), Shandong Provincial Natural Science Foundation (ZR2021MA069, ZR2017MA040 and ZR2020ME231), and Key Research and Development Project of Shandong province (No. 2019GGX102050).

References

  1. J. Spiegel, P. M. Cromm, G. Zimmermann, T. N. Grossmann and H. Waldmann, Nat. Chem. Biol., 2014, 10, 613–622 CrossRef CAS PubMed.
  2. J. Ma and M. Karplus, Proc. Natl. Acad. Sci. U. S. A., 1997, 94, 11905–11910 CrossRef CAS PubMed.
  3. J. L. Bos, H. Rehmann and A. Wittinghofer, Cell, 2007, 129, 865–877 CrossRef CAS PubMed.
  4. S. P. Mo, J. M. Coulson and I. A. Prior, Biochem. Soc. Trans., 2018, 46, 1325–1332 CrossRef CAS PubMed.
  5. M. Eren, N. Tuncbag, H. Jang, R. Nussinov, A. Gursoy and O. Keskin, J. Phys. Chem. B, 2021, 125, 5210–5221 CrossRef CAS PubMed.
  6. G. A. Hobbs, C. J. Der and K. L. Rossman, J. Cell Sci., 2016, 129, 1287–1292 CAS.
  7. A. E. Karnoub and R. A. Weinberg, Nat. Rev. Mol. Cell Biol., 2008, 9, 517–531 CrossRef CAS PubMed.
  8. I. A. Prior, P. D. Lewis and C. Mattos, Cancer Res., 2012, 72, 2457–2467 CrossRef CAS PubMed.
  9. J. G. Tate, S. Bamford, H. C. Jubb, Z. Sondka, D. M. Beare, N. Bindal, H. Boutselakis, C. G. Cole, C. Creatore, E. Dawson, P. Fish, B. Harsha, C. Hathaway, S. C. Jupe, C. Y. Kok, K. Noble, L. Ponting, C. C. Ramshaw, C. E. Rye, H. E. Speedy, R. Stefancsik, S. L. Thompson, S. Wang, S. Ward, P. J. Campbell and S. A. Forbes, Nucleic Acids Res., 2018, 47, D941–D947 CrossRef PubMed.
  10. S. Xu, B. N. Long, G. H. Boris, A. Chen, S. Ni and M. A. Kennedy, Acta Crystallogr., Sect. D: Struct. Biol., 2017, 73, 970–984 CrossRef CAS PubMed.
  11. J. M. Ostrem, U. Peters, M. L. Sos, J. A. Wells and K. M. Shokat, Nature, 2013, 503, 548–551 CrossRef CAS PubMed.
  12. K. Bum-Erdene, D. Liu, G. Gonzalez-Gutierrez, M. K. Ghozayel, D. Xu and S. O. Meroueh, Proc. Natl. Acad. Sci. U. S. A., 2020, 117, 7131–7139 CrossRef CAS PubMed.
  13. A. C. Nelson, T. J. Turbyville, S. Dharmaiah, M. Rigby, R. Yang, T.-Y. Wang, J. Columbus, R. Stephens, T. Taylor, D. Sciacca, G. Onsongo, A. Sarver, S. Subramanian, D. V. Nissley, D. K. Simanshu and E. Lou, J. Biol. Chem., 2020, 295, 9335–9348 CrossRef CAS PubMed.
  14. E. J. Poulin, A. K. Bera, J. Lu, Y.-J. Lin, S. D. Strasser, J. A. Paulo, T. Q. Huang, C. Morales, W. Yan, J. Cook, J. A. Nowak, D. K. Brubaker, B. A. Joughin, C. W. Johnson, R. A. DeStefanis, P. C. Ghazi, S. Gondi, T. E. Wales, R. E. Iacob, L. Bogdanova, J. J. Gierut, Y. Li, J. R. Engen, P. A. Perez-Mancera, B. S. Braun, S. P. Gygi, D. A. Lauffenburger, K. D. Westover and K. M. Haigis, Cancer Discovery, 2019, 9, 738–755 CrossRef CAS PubMed.
  15. G. Pálfy, I. Vida and A. Perczel, Biomol. NMR Assignments, 2020, 14, 1–7 CrossRef PubMed.
  16. T. Pantsar, Comput. Struct. Biotechnol. J., 2020, 18, 189–198 CrossRef CAS PubMed.
  17. S. Lu, H. Jang, S. Muratcioglu, A. Gursoy, O. Keskin, R. Nussinov and J. Zhang, Chem. Rev., 2016, 116, 6607–6665 CrossRef CAS PubMed.
  18. S. Lu, H. Jang, S. Gu, J. Zhang and R. Nussinov, Chem. Soc. Rev., 2016, 45, 4929–4952 RSC.
  19. G. W. K. Moore, S. E. L. Howell, M. Brady, X. Xu and K. McNeil, Nat. Commun., 2021, 12, 1 CrossRef CAS PubMed.
  20. J. Chen, W. Wang, L. Pang and W. Zhu, Phys. Chem. Chem. Phys., 2020, 22, 21238–21250 RSC.
  21. P. A. Boriack-Sjodin, S. M. Margarit, D. Bar-Sagi and J. Kuriyan, Nature, 1998, 394, 337–343 CrossRef CAS PubMed.
  22. D. Rabara, T. H. Tran, S. Dharmaiah, R. M. Stephens, F. McCormick, D. K. Simanshu and M. Holderfield, Proc. Natl. Acad. Sci. U. S. A., 2019, 116, 22122–22131 CrossRef CAS PubMed.
  23. Y. Li, Y. Zhang, F. Großerüschkamp, S. Stephan, Q. Cui, C. Kötting, F. Xia and K. Gerwert, J. Phys. Chem. Lett., 2018, 9, 1312–1317 CrossRef CAS PubMed.
  24. S. U. Kumar and C. G. Priya Doss, Comput. Biol. Med., 2021, 135, 104597 CrossRef CAS PubMed.
  25. J. C. Hunter, D. Gurbani, S. B. Ficarro, M. A. Carrasco, S. M. Lim, H. G. Choi, T. Xie, J. A. Marto, Z. Chen, N. S. Gray and K. D. Westover, Proc. Natl. Acad. Sci. U. S. A., 2014, 111, 8895–8900 CrossRef CAS PubMed.
  26. S. Lu, D. Ni, C. Wang, X. He, H. Lin, Z. Wang and J. Zhang, ACS Catal., 2019, 9, 7188–7196 CrossRef CAS.
  27. D. K. Simanshu, D. V. Nissley and F. McCormick, Cell, 2017, 170, 17–33 CrossRef CAS PubMed.
  28. A. C. Kaushik, Y.-J. Wang, X. Wang and D.-Q. Wei, Briefings Bioinf., 2020, 22, bbaa149 CrossRef PubMed.
  29. J. Lu, A. K. Bera, S. Gondi and K. D. Westover, Biochemistry, 2018, 57, 324–333 CrossRef CAS PubMed.
  30. M. Khaled, A. Gorfe and A. Sayyed-Ahmad, J. Phys. Chem. B, 2019, 123, 7667–7675 CrossRef CAS PubMed.
  31. Y. Lee, C. Phelps, T. Huang, B. Mostofian, L. Wu, Y. Zhang, K. Tao, Y. H. Chang, P. J. S. Stork, J. W. Gray, D. M. Zuckerman and X. Nan, eLife, 2019, 8, e46393 CrossRef PubMed.
  32. S. E. Acuner, F. Sumbul, H. Torun and T. Haliloglu, Biophys. J., 2021, 120, 866–876 CrossRef CAS PubMed.
  33. J. Zeng, J. Weng, Y. Zhang, F. Xia, Q. Cui and X. Xu, J. Phys. Chem. B, 2021, 125, 8805–8813 CrossRef CAS PubMed.
  34. S. Lu, H. Jang, R. Nussinov and J. Zhang, Sci. Rep., 2016, 6, 21949 CrossRef CAS PubMed.
  35. P. Y. Ting, C. W. Johnson, C. Fang, X. Cao, T. G. Graeber, C. Mattos and J. Colicelli, FASEB J., 2015, 29, 3750–3761 CrossRef CAS PubMed.
  36. E. Denayer, A. Parret, M. Chmara, S. Schubbert, A. Vogels, K. Devriendt, J.-P. Frijns, V. Rybin, T. J. de Ravel, K. Shannon, J. Cools, K. Scheffzek and E. Legius, Hum. Mutat., 2008, 29, 232–239 CrossRef PubMed.
  37. S. Lukman, B. J. Grant, A. Gorfe, G. H. Grant and J. A. McCammon, PLoS Comput. Biol., 2010, 6, e1000922 CrossRef PubMed.
  38. Y. Matsuo, P. M. Campbell, R. A. Brekken, B. Sung, M. M. Ouellette, J. B. Fleming, B. B. Aggarwal, C. J. Der and S. Guha, Mol. Cancer Res., 2009, 7, 799–808 CrossRef CAS PubMed.
  39. Y. Hu, W. Lu, G. Chen, P. Wang, Z. Chen, Y. Zhou, M. Ogasawara, D. Trachootham, L. Feng, H. Pelicano, P. J. Chiao, M. J. Keating, G. Garcia-Manero and P. Huang, Cell Res., 2012, 22, 399–412 CrossRef CAS PubMed.
  40. S. Vatansever, B. Erman and Z. H. Gümüş, Comput. Struct. Biotechnol. J., 2020, 18, 1000–1011 CrossRef CAS PubMed.
  41. M. R. Mehaffey, C. L. Schardon, E. T. Novelli, M. B. Cammarata, L. J. Webb, W. Fast and J. S. Brodbelt, Chem. Sci., 2019, 10, 8025–8034 RSC.
  42. A. A. Gorfe, B. J. Grant and J. A. McCammon, Structure, 2008, 16, 885–896 CrossRef CAS PubMed.
  43. M. E. Welsch, A. Kaplan, J. M. Chambers, M. E. Stokes, P. H. Bos, A. Zask, Y. Zhang, M. Sanchez-Martin, M. A. Badgley, C. S. Huang, T. H. Tran, H. Akkiraju, L. M. Brown, R. Nandakumar, S. Cremers, W. S. Yang, L. Tong, K. P. Olive, A. Ferrando and B. R. Stockwell, Cell, 2017, 168, 878–889.e829 CrossRef CAS PubMed.
  44. J. Chen, L. Wang, W. Wang, H. Sun, L. Pang and H. Bao, Comput. Biol. Med., 2021, 135, 104639 CrossRef CAS PubMed.
  45. J. Chen, S. Zhang, W. Wang, L. Pang, Q. Zhang and X. Liu, J. Chem. Inf. Model., 2021, 61, 1954–1969 CrossRef CAS PubMed.
  46. P. Prakash, J. F. Hancock and A. A. Gorfe, Proteins, 2015, 83, 898–909 CrossRef CAS PubMed.
  47. S. Vatansever, B. Erman and Z. H. Gümüş, Sci. Rep., 2019, 9, 11730 CrossRef PubMed.
  48. C. I. Nnadi, M. L. Jenkins, D. R. Gentile, L. A. Bateman, D. Zaidman, T. E. Balius, D. K. Nomura, J. E. Burke, K. M. Shokat and N. London, J. Chem. Inf. Model., 2018, 58, 464–471 CrossRef CAS PubMed.
  49. C. W. Johnson, Y.-J. Lin, D. Reid, J. Parker, S. Pavlopoulos, P. Dischinger, C. Graveel, A. J. Aguirre, M. Steensma, K. M. Haigis and C. Mattos, Cell Rep., 2019, 28, 1538–1550.e1537 CrossRef CAS PubMed.
  50. A. M. Levinson, J. H. McGee, A. G. Roberts, G. S. Creech, T. Wang, M. T. Peterson, R. C. Hendrickson, G. L. Verdine and S. J. Danishefsky, J. Am. Chem. Soc., 2017, 139, 7632–7639 CrossRef CAS PubMed.
  51. Z.-W. Zhou, C. Ambrogio, A. K. Bera, Q. Li, X.-X. Li, L. Li, J. Son, S. Gondi, J. Li, E. Campbell, H. Jin, J. J. Okoro, C.-X. Xu, P. A. Janne and K. D. Westover, Cancer Res., 2020, 80, 3719–3731 CAS.
  52. W. Yan, E. Markegard, S. Dharmaiah, A. Urisman, M. Drew, D. Esposito, K. Scheffzek, D. V. Nissley, F. McCormick and D. K. Simanshu, Cell Rep., 2020, 32 Search PubMed.
  53. G. Buhrman, G. Wink and C. Mattos, Structure, 2007, 15, 1618–1629 CrossRef CAS PubMed.
  54. F. Yan, X. Liu, S. Zhang, J. Su, Q. Zhang and J. Chen, Int. J. Mol. Sci., 2018, 19, 2496 CrossRef PubMed.
  55. A. H. Mazurek, Ł. Szeleszczuk and T. Gubica, Int. J. Mol. Sci., 2021, 22, 9422 CrossRef CAS PubMed.
  56. W. Xue, F. Yang, P. Wang, G. Zheng, Y. Chen, X. Yao and F. Zhu, ACS Chem. Neurosci., 2018, 9, 1128–1140 CrossRef CAS PubMed.
  57. T. Hou and R. Yu, J. Med. Chem., 2007, 50, 1177–1188 CrossRef CAS PubMed.
  58. J. Chen, J. Wang, B. Yin, L. Pang, W. Wang and W. Zhu, ACS Chem. Neurosci., 2019, 10, 4303–4318 CrossRef CAS PubMed.
  59. M.-J. Yang, X.-Q. Pang, X. Zhang and K.-L. Han, J. Struct. Biol., 2011, 173, 57–66 CrossRef CAS PubMed.
  60. R. Wang and Q. Zheng, Langmuir, 2020, 36, 13817–13832 CrossRef CAS PubMed.
  61. Z. Sun, Z. Gong, F. Xia and X. He, Chem. Phys., 2021, 548, 111245 CrossRef CAS.
  62. Z. Sun, Z. Huai, Q. He and Z. Liu, J. Chem. Inf. Model., 2021, 61(12), 6107–6134 CrossRef CAS PubMed.
  63. J. Chen, C. Peng, J. Wang and W. Zhu, Proteins, 2018, 86, 1294–1305 CrossRef CAS PubMed.
  64. Y. Sugita and Y. Okamoto, Chem. Phys. Lett., 1999, 314, 141–151 CrossRef CAS.
  65. B. J. Grant, A. A. Gorfe and J. A. McCammon, PLoS Comput. Biol., 2009, 5, e1000325 CrossRef PubMed.
  66. L. C. T. Pierce, R. Salomon-Ferrer, C. A. F. de Oliveira, J. A. McCammon and R. C. Walker, J. Chem. Theory Comput., 2012, 8, 2997–3002 CrossRef CAS PubMed.
  67. X. Peng, Y. Zhang, Y. Li, Q. Liu, H. Chu, D. Zhang and G. Li, J. Chem. Theory Comput., 2018, 14, 1216–1227 CrossRef CAS PubMed.
  68. L. Duan, X. Guo, Y. Cong, G. Feng, Y. Li and J. Z. H. Zhang, Front. Chem., 2019, 7, 540 CrossRef CAS PubMed.
  69. J. Chen, B. Yin, W. Wang and H. Sun, ACS Chem. Neurosci., 2020, 11, 1811–1826 CrossRef CAS PubMed.
  70. Y. Miao, V. A. Feher and J. A. McCammon, J. Chem. Theory Comput., 2015, 11, 3584–3595 CrossRef CAS PubMed.
  71. J. Wang and Y. Miao, J. Phys. Chem. B, 2019, 123, 6462–6473 CrossRef CAS PubMed.
  72. J. Wang, P. R. Arantes, A. Bhattarai, R. V. Hsu, S. Pawnikar, Y.-m. M. Huang, G. Palermo and Y. Miao, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2021, 11, e1521 Search PubMed.
  73. J. Wang and Y. Miao, J. Chem. Phys., 2020, 153, 154109 CrossRef CAS PubMed.
  74. X. An, Q. Bai, Z. Bing, H. Liu, Q. Zhang, H. Liu and X. Yao, ACS Chem. Neurosci., 2020, 11, 628–637 CrossRef CAS PubMed.
  75. J. Chen, S. Zhang, W. Wang, H. Sun, Q. Zhang and X. Liu, ACS Chem. Neurosci., 2021, 12, 2591–2607 CrossRef CAS PubMed.
  76. Y. Miao and J. A. McCammon, Proc. Natl. Acad. Sci. U. S. A., 2016, 113, 12162–12167 CrossRef CAS PubMed.
  77. J. Chen, W. Wang, H. Sun, L. Pang and H. Bao, Comput. Biol. Med., 2021, 134, 104485 CrossRef CAS PubMed.
  78. Q. Shao and W. Zhu, J. Chem. Theory Comput., 2017, 13, 4240–4252 CrossRef CAS PubMed.
  79. J. C. Hunter, A. Manandhar, M. A. Carrasco, D. Gurbani, S. Gondi and K. D. Westover, Mol. Cancer Res., 2015, 13, 1325–1335 CrossRef CAS PubMed.
  80. B. Webb and A. Sali, Comparative Protein Structure Modeling Using MODELLER, Current Protocols in Bioinformatics, John Wiley & Sons, Inc., 2014, p. 2015, 2016.2011-2015.2016.2032 Search PubMed.
  81. R. Anandakrishnan, B. Aguilar and A. V. Onufriev, Nucleic Acids Res., 2012, 40, W537–W541 CrossRef CAS PubMed.
  82. R. Salomon-Ferrer, D. A. Case and R. C. Walker, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2013, 3, 198–210 CAS.
  83. D. A. Case, T. E. Cheatham III, T. Darden, H. Gohlke, R. Luo, K. M. Merz Jr, A. Onufriev, C. Simmerling, B. Wang and R. J. Woods, J. Comput. Chem., 2005, 26, 1668–1688 CrossRef CAS PubMed.
  84. K. L. Meagher, L. T. Redman and H. A. Carlson, J. Comput. Chem., 2003, 24, 1016–1025 CrossRef CAS PubMed.
  85. J. A. Maier, C. Martinez, K. Kasavajhala, L. Wickstrom, K. E. Hauser and C. Simmerling, J. Chem. Theory Comput., 2015, 11, 3696–3713 CrossRef CAS PubMed.
  86. W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey and M. L. Klein, J. Chem. Phys., 1983, 79, 926–935 CrossRef CAS.
  87. J. Aqvist, J. Phys. Chem., 1990, 94, 8021–8024 CrossRef.
  88. R. Salomon-Ferrer, A. W. Götz, D. Poole, S. Le Grand and R. C. Walker, J. Chem. Theory Comput., 2013, 9, 3878–3888 CrossRef CAS PubMed.
  89. Y. Miao, W. Sinko, L. Pierce, D. Bucher, R. C. Walker and J. A. McCammon, J. Chem. Theory Comput., 2014, 10, 2677–2689 CrossRef CAS PubMed.
  90. J.-P. Ryckaert, G. Ciccotti and H. J. C. Berendsen, J. Comput. Phys., 1977, 23, 327–341 CrossRef CAS.
  91. J. A. Izaguirre, D. P. Catarello, J. M. Wozniak and R. D. Skeel, J. Chem. Phys., 2001, 114, 2090–2098 CrossRef CAS.
  92. R. M. Levy, A. R. Srinivasan, W. K. Olson and J. A. McCammon, Biopolymers, 1984, 23, 1099–1112 CrossRef CAS PubMed.
  93. A. Amadei, A. B. M. Linssen and H. J. C. Berendsen, Proteins, 1993, 17, 412–425 CrossRef CAS PubMed.
  94. M. Laberge and T. Yonetani, Biophys. J., 2008, 94, 2737–2751 CrossRef CAS PubMed.
  95. J. Chen, W. Wang, H. Sun, L. Pang and B. Yin, J. Comput.-Aided Mol. Des., 2020, 34, 1289–1305 CrossRef CAS PubMed.
  96. A. D. McLachlan, J. Mol. Biol., 1979, 128, 49–79 CrossRef CAS PubMed.
  97. D. R. Roe and T. E. Cheatham, J. Chem. Theory Comput., 2013, 9, 3084–3095 CrossRef CAS PubMed.
  98. W. Humphrey, A. Dalke and K. Schulten, J. Mol. Graphics, 1996, 14, 33–38 CrossRef CAS PubMed.
  99. T. Ichiye and M. Karplus, Proteins, 1991, 11, 205–217 CrossRef CAS PubMed.
  100. S. Salentin, S. Schreiber, V. J. Haupt, M. F. Adasme and M. Schroeder, Nucleic Acids Res., 2015, 43, W443–W447 CrossRef CAS PubMed.
  101. A. Sayyed-Ahmad, P. Prakash and A. A. Gorfe, Proteins, 2017, 85, 1618–1632 CrossRef CAS PubMed.
  102. J. Zhao, L. Wang, H. Bao and J. Chen, Chem. Phys. Lett., 2021, 781, 138995 CrossRef CAS.
  103. P. Prakash, A. Sayyed-Ahmad and A. A. Gorfe, PLoS Comput. Biol., 2012, 8, e1002394 CrossRef CAS PubMed.
  104. X. Wang, J. Chem. Inf. Model., 2021, 61(10), 5212–5222 CrossRef CAS PubMed.
  105. M. V. Milburn, L. Tong, A. M. deVos, A. Brünger, Z. Yamaizumi, S. Nishimura and S.-H. Kim, Science, 1990, 247, 939–945 CrossRef CAS PubMed.
  106. J. A. Parker, A. Y. Volmar, S. Pavlopoulos and C. Mattos, Structure, 2018, 26, 810–820.e814 CrossRef CAS PubMed.

Footnote

Electronic supplementary information (ESI) available. See DOI: 10.1039/d1ra07936k

This journal is © The Royal Society of Chemistry 2022