Deciphering the structural and dynamic effects of SHP2-E76 mutations: mechanistic insights into oncogenic activation

Abstract

The tyrosine phosphatase known as SHP2 is a cytoplasmic protein and encodes by proto-oncogene PTPN11. This protein is essential for the regulation of cell growth, differentiation, programed cell death, and survival. This regulation is achieved through the release of intramolecular autoinhibition and the modulation of several signaling pathways, including the signaling cascade of Ras-MAPK. Mutations in SHP2 are frequently associated with human malignancies and neurodevelopmental disorders (NDDs). Specifically, a germline mutation (E76D) in SHP2 is linked to neurodevelopmental disorders, such as Noonan syndrome, while somatic mutations (E76G and E76A) and altered SHP2 expression are implicated in several forms of leukemia. These mutations disrupt the closed conformation, which normally keeps SHP2 in an inactive, auto-inhibited state, thereby enhancing phosphatase activity and activating SHP2, leading to a gain-of-function effect. However, the structural and functional implications of these disease-related mutants are not well elucidated. Therefore, in this study, we investigate the structural mechanisms underlying three distinct gain-of-function SHP2 mutations (E76D, E76G, and E76A) through the application of molecular dynamics (MD) simulations, focusing on how a single amino acid mutation at the same position result in different disease phenotypes, either cause cancer or NDDs. Notably, Patients with Noonan Syndrome have an increased risk of developing cancer, suggesting a potential link between these diseases and their mutations. MD simulation was employed to elucidate this mechanism, examining four distinct states: Apo-state (E76), M1-state (E76D), M2-state (E76G), and M3-state (E76A). The dynamics and conformational changes of SHP2 in both its Apo-state and mutant states (M1, M2, and M3) were compared. Our findings indicate that both cancer-related and NDD-related mutations destabilize the N-SH2 and PTP interface, facilitating SHP2 activation. However, the cancer-associated mutations induce more severe disruption at the N-SH2 and PTP interface than the NDD mutations. Additionally, dynamic analyses revealed that mutations at the interface (M1, M2, and M3) not only alter the native folded conformation of SHP2 but also significantly enhance the C-distance between the N-SH2 and PTP domains. Overall, this study provides a comprehensive understanding of the structural dynamics of SHP2 at the atomic level, revealing how mutations disrupt its auto-inhibition and increase PTP activity, providing valuable insights into the molecular mechanisms driving both cancer and neurodevelopmental disorders.

Introduction

The implications of mutations on both structure and function commonly found in oncogenes and tumor suppressor genes are still not fully understood. The initial proto-oncogene, PTPN11 (MIM 176876), located at 12q24.13, consists of 16 exons that cover a genomic area of approximately 92 kb [1]. This gene encodes a non-receptor cytoplasmic protein known as SRC homology phosphatase 2 (SHP2), which serves as a signal-enhancing element in the signaling pathways of the regulatory components of tissue homeostasis. In the regulation of various cellular processes, SHP2 plays a crucial role, including proliferation, differentiation, survival, metabolism, motility, and migration [2]. Biochemical and genetic studies have demonstrated that SHP2 is expressed in a wide range of tissues and functions upstream of the RAS signaling pathway. Additionally, it transcriptionally influences several downstream signaling pathways, including NF-κB, RAS-ERK1/2, RAS-RAF-MAP kinase, PI3K-AKT, JAK-STAT, apoptosis, and immune checkpoint (BTLA) pathways; however, its involvement in these pathways is dependent on the receptor and/or the specific cellular context [3,4,5,6]. Numerous germline mutations in the SHP2 gene are responsible for various neurodevelopmental disorders (NDDs), including 50% of cases of Noonan syndrome and 90% of cases of LEOPARD syndrome [7, 8]. In contrast, gain-of-function somatic mutations are linked to several forms of leukemia, being identified in 35% of juvenile myelomonocytic leukemia cases, 10% of myelodysplastic syndromes, 7% of B-precursor acute lymphoblastic leukemia cases, and 5% of acute myeloid leukemia cases, with lower frequencies observed in other hematological malignancies and solid tumors [9,10,11,12,13]. Additionally, it has been found the mutations in SHP2 occur in multiple myeloma as well as in solid tumors affecting different organs [14, 15]. Furthermore, cancer-related mutations in SHP2 have been shown to induce leukemia in murine models [16].

The structural organization of the classical NRPTP, SHP2, which has a molecular weight of 92 kDa and consists of 593 amino acids, features two N-terminal regulatory SH2 domains: the N-SH2 domain (residues 2–104) and the C-SH2 domain (residues 112–215), as illustrated in Fig. 1. Additionally, it contains a catalytic PTP domain (residues 220–525) located at the amino and carboxyl terminals of the SHP2 protein. The structure also includes a C-terminal tail that possesses four tyrosine phosphorylation sites (Tyr542, Tyr546, Tyr580, and Tyr584) along with a [17, 18] proline-rich region.

Fig. 1

a Schematic representation of SHP2 activation. In the Apo-state, SHP2 remains in a closed autoinhibited conformation. Peptide binding can either stabilize this autoinhibited state or induce a conformational transition of the N-SH2 domain to the α-state, weakening N-SH2–PTP interactions. This transition exposes the catalytic site, leading to SHP2 activation. b Location of PTPN11 pathogenic variants mapped to exons and SHP-2 functional domains. The top bar represents PTPN11 exons, with the coding region in pink. Blue boxes indicate the number of patients with each pathogenic variant. Below, the SHP-2 protein structure is shown, highlighting the N-SH2, C-SH2, and PTP domains [19, 20]

The enzymatic activity of SHP2 as phosphatase is crucial for its various biological functions, with the phosphorylation site at the C-terminus playing a regulatory role in certain signaling pathways. The X-ray crystallography of SHP2 indicates that in its resting state, SHP2 is inactive and adopts a closed, auto-inhibited structure. In this state, the N-SH2 domain interacts with the PTP domain at catalytic site, effectively obstructing the active site [21]. This contact involves the insertion of a βD’-D’E-βE loop, comprising residues N58, T59, G60, D61, Y62, and A72, into the catalytic cleft formed by residues C459, D464, R465, Q506, and G503. Consequently, this arrangement diminishes solvent accessibility and prevents substrate entry into the catalytic site [22]. The reverse side of the N-SH2 domain interacts with proteins that contain phosphorylated tyrosine residues, such as immune inhibitory receptors and scaffolding adaptors like Grb2-associated binder proteins (Gabs) and fibroblast growth factor receptor substrate [23]. This interaction prompts SHP2 to adopt an open conformation, resulting in increased phosphatase activity and the exposure of the catalytic pocket to substrates [23, 24]. Conversely, oncogenic mutations in SHP2 disrupt the auto-inhibited conformation, leading to heightened basal activity even in the absence of stimulation by tyrosine phosphorylated ligands [25]. Such binding events or induced point mutations interfere with the interaction between the N-SH2 domain and the PTP, thereby activating the phosphatase. Consequently, the N-SH2 domain functions as a conformational switch, capable of either binding to and inhibiting the phosphatase or binding to phosphotyrosyl peptides (pY) to activate the enzymes.

The existing body of literature has indicated that various mutations, including E76D, E76K, T73I, E76V, D61Y, E76Q, E76G, S502P, and E76A, are linked to neurofibromatosis (germline origin) and leukemia (somatic origin), both of which are associated with cancer. Consequently, the investigations of structure and function of pro-oncogenic SHP2 variants have demonstrated a relationship between the location of amino acid mutations within the sequence of protein and the increased phosphatase activity. The majority of mutations related to cancer are found at the auto-inhibition interface between the N-SH2 and PTP domains, generally leading to the activation of the phosphatase [26]. The germline missense mutation from glutamic acid to aspartic acid at position 76 of SHP2 (E76D) is a relatively conserved mutation, identified in 50% of patients with (NDDs), such as Noonan syndrome (MIM#163950) [27]. The first human gain-of-function (GOF) mutation in SHP2 was reported early in this decade among individuals with Noonan syndrome [7]. Noonan syndrome (NS) is a clinically diverse condition associated with craniofacial anomalies, cardiac defects, reduced stature, and learning difficulties. The majority of mutations documented in NS families are located in exons 3 and 8, which correspond to the N-SH2 and PTP domains, respectively. A locus for Noonan syndrome on the long arm of chromosome 12 has been mapped in some, but not all, NS families [28], leading to the identification of mutations in the candidate gene PTPN11 within the chromosome 12q22-qter region [7]. Patients diagnosed with (NDDs), face an elevated risk of developing juvenile myelomonocytic leukemia (JMML), a severe myeloproliferative disorder characterized by the excessive production of myelomonocytic cells that infiltrate various organs, including the spleen, liver, and skin [29,30,31].

Additionally, the missense somatic gain-of-function mutation from glutamic acid to glycine at position 76 of SHP2 (E76G) is implicated in acute myeloid leukemia (AML) [32]. The potential causal relationship between activating mutations in PTPN11 phosphatase and the pathogenesis of acute leukemias remains poorly defined. Consequently, the significance of SHP2 mutations in AML has not been fully elucidated. However, Chen et al. discovered gain-of-function PTPN11 mutations in 6.6% of 91 pretreatment AML samples. The E76G mutations are situated in exon 3, which encodes the N-SH2 domain. The SHP2-E76G mutation results in diminished interactions between the N-SH2 and PTP domains, leading to the destabilization of the auto-inhibited structure. Consequently, patients with acute myeloid leukemia (AML) are more likely to experience over activation, which may also result in juvenile myelomonocytic leukemia (JMML). Additionally, another variant of the somatic SHP2 gain-of-function mutation, E76A, located in exon 3, is associated with myelodysplastic syndrome (MDS). Myelodysplastic syndromes represent a collection of heterogeneous bone marrow disorders characterized by the production of insufficient healthy blood cells. The E76A mutation occurs at the N-SH2 domain interface and enhances SHP2 activity. Patients diagnosed with MDS are at an increased risk of progressing to AML. It is noteworthy that mutations found in leukemia and (NDDs), occur at the same amino acid positions; however, the types of substitutions are rarely identical. Numerous mutations have been identified at position E76, all of which are situated within exon-3, responsible for encoding the N-SH2 domain of SHP2. Consequently, the E76 position is regarded as a mutational hotspot for leukemia, characterized by various residue substitutions, and serves as a primary target for drug development. Saturation analysis pertaining to cancer indicate that the PTPN11 gene is anticipated to be functionally critical and represents a promising candidate for the development of cancer therapeutics. Despite the existing knowledge regarding SHP2 variations linked to pathology, the structural and functional implications of many mutations remain inadequately understood, and the structural comparison among different pro-oncogenic mutations is still largely incomplete. Furthermore, the literature clearly indicates that the aforementioned mutations (E76D, E76G, and E76A) interfere with the residue-residue interactions at the interface between the N-SH2 and PTP domains. This disruption enhances the phosphatase activity of SHP2 by eliminating the autoinhibition mechanism, thereby resulting in gain-of-function (GOF) activity. Although these mutations share similarities by occurring at the same position, they all contribute to increased SHP2 activity (GOF), elevate the RAS/RAF/ERK/MAPK signaling pathways, and disrupt the βD’-D’E-βE loop residues. This raises questions about the conformational or structural changes that may lead to various disease conditions. To investigate these matters, we will conduct molecular dynamics simulations.

In this study we examine the structural and mechanistic features associated with three mutations (SHP2E76D, SHP2E76G, and SHP2E76A) of SHP2 utilizing MD simulations. This research aims to identify the conformational or dynamic alterations that significantly affect the physiological properties of these mutations. Molecular dynamics simulations have been extensively employed to analyze structural variations at the atomic level and have been used in elucidating the consequences of missense mutations at the interface between the N-SH2 and PTP domains. It has been established that the extent of mutation activity correlates with the severity of associated diseases. Notably, Patients with NDDs have an increased risk of developing cancer, suggesting a potential link between these diseases and their mutations. Furthermore, the phosphatase activity of SHP2 is further augmented due to the distance between the interface residues of the N-SH2 and PTP domains, which significantly increases. The resulting variations in CA distance between the PTP and N-SH2 domains indicate that certain paired residues either lose (F7/Q256, E76/S502, and E76/R265) or exhibit weakened (E76/N308) hydrogen interactions, which ultimately facilitates the opening of the catalytic site for substrate binding. The dynamic mechanical properties of polymorphic SHP2 mutations remain ambiguous. Consequently, MD simulations of both the mutant states (SHP2E76D, SHP2E76G, and SHP2E76A) and the Apo state (SHP2 + E76) were conducted to elucidate the comprehensive physiological mechanisms of the SHP2 protein. These findings imply that pathological mutations in SHP2 exert their influence through various mechanisms, either by increasing the distance between interface residues or by altering the βD’-D’E-βE loop, which is directly associated with SHP2 activation and disease severity. This also establishes a framework for understanding the contributions of these mutations to neurodevelopmental disorders (NDDs) and human malignance. Our results provide a direct comparison of the effects of each mutation on the auto-inhibition of the phosphatase and offer mechanistic insights that enhance the understanding of mutation-specific SHP2 therapies.

Material and methods

System preparation

The X-ray crystallographic structure of SHP2WT (PDB ID: 2SHP) along with two mutant variants, SHP2E76D and SHP2E76A (PDB IDs: 6CMR and 6IHZ, respectively), were downloaded from the Protein Data Bank and served as the starting coordinates for this study. However, the three dimensional coordinates of the SHP2E76G mutant is not yet available in the protein data bank. Consequently, the structure of the SHP2E76G mutant was created by introducing the mutation at a designated site within the crystallographic structure of SHP2WT using PYMOL[33]. The systematic exemplification was demonstrated by Apo (E76), M1 (E76D), M2 (E76G), and M3 (E76A). To ensure high-quality structures for further studies, all crystal structures were meticulously prepared using UCSF Chimera [34] and CHARMM-GUI [35]. Protein preparation is essential to eliminate unwanted molecules, resolve steric clashes, and optimize the structural integrity of the system. Using Chimera, non-essential water molecules, ligands, and redundant heteroatoms were removed, and hydrogen atoms were added to ensure proper valency. The protonation states of ionizable residues were assigned using the CHARMM-GUI PDB Reader & Manipulator, which optimizes the protonation states based on physiological pH conditions. Additionally, missing heavy atoms and side chains were reconstructed to ensure structural completeness. The prepared structures were then energy-minimized to remove steric hindrances and stabilize the system before further molecular dynamics simulations. The finalized structures were subjected to molecular dynamics simulations using AMBER to assess their conformational stability and dynamic properties [36].

All-atom MD simulation

MD simulations is a powerful technique for extensive study of dynamical, structural, and energetic insights into bio-macromolecules, target proteins, and their interactions with drugs in therapeutic contexts. All-atom MD simulations of both the Apo-state and Mutant-state were performed using the AMBER [36]. Four distinct systems were established for the MD simulations: Apo-state, M1-state, M2-state, and M3-state. The leaprc.ff14SB force field was utilized through Antechamber to parameterize the proteins accordingly [37]. Additionally, the histidine residues were protonated by pdb4amber script, thereby adjusting the protein system automatically before executing the AMBER-tLEaP module [38]. The hydrogen atoms were added to the crystal lattice by tLEaP. Subsequently, the system was neutralized by adding the counter, leading to the creation of topologies and parameter files. The solvation of all the systems within a truncated octahedral box containing TIP3P water molecules was performed with a cut-off buffer of 10.0 Å. To address long-range electrostatic interactions the Particle Mesh Ewald (PME) method was employed [39]. In all MD simulations the ff14SB force field was applied. The SHAKE algorithm, with a tolerance of 10^-5 Å, was implemented to constrain all covalent bonds involving hydrogen atoms [40]. The solvated system was minimized over 20, 000 steps using steepest descent method, followed by heating up 400 ps and an equilibration of 200 ps within the NVT ensemble. Temperature and pressure were coupled with a time constant of 1.0 ps, utilizing isotropic position scaling and a relaxation time of 2.0 ps, in accordance with Langevin’s algorithm [41]. Subsequent analyses, including root mean square deviation (RMSD), root mean square fluctuation (RMSF), dynamic cross-correlation mapping (DCCM), distance calculations, solvent-accessible surface area, and radius of gyration, were performed utilizing the CPPTRAJ modules integrated within Amber version 2018. Structural and visual assessments, along with data plotting, were executed using PyMol V1.3 and Origin version 9.1 data analysis software, respectively [42].

Dynamic cross-correlation matrix (DCCM) analysis

The dynamic properties of the SHP2 protein in both the Apo-state and the Mutant-state can be examined through MD simulations to provide insights into the correlated movements across all states. We constructed and evaluated the DCCM to explore the fluctuations and movements of the backbone of the Cα atoms in each state [43, 44]. This analysis offers comprehensive information regarding the conformational alterations resulting from mutations on the overall dynamics of the protein. The degree of correlated motion is represented by the magnitude of the respective correlation coefficient. The cross-correlation coefficient for the displacement of any two atoms, i and j, is determined using the equation presented below:

$$Cij\, = \,\langle \Delta ri. \, \Delta rj\rangle (\langle \Delta ri{2}\rangle \langle \Delta rj{2}\rangle ){ 1}/{2}$$

Cij denotes the deviation from the average position of the ith atom, with the angle brackets indicating the time average across the entire trajectory. The Cij values can be organized into a matrix format, which can then be represented as a three-dimensional DCCM [45]. The elements within the cross-correlation matrix range from -1 to 1. Strongly positive Cij values suggest a significant correlated movement between residues i and j, whereas negative Cij values indicate that the two residues are moving in opposite directions, reflecting anticorrelated motion. The DCCM matrix was generated utilizing the CPPTRAJ package in Amber 14, and the resulting matrices were plotted and analyzed using Origin software [42].

Results and discussion

SHP2 structural feature

The 3D coordinates of SHP2 (PDB: 2SHP, including residues 1-527) is characterized by a closed auto-inhibited conformation, which features two tandem SH2 domains (N-SH2 and C-SH2) alongside a catalytic (PTP) domain. In this arrangement, the N-SH2 domain engages in an intramolecular interaction with the PTP domain. The structure reveals that both SH2 domains are composed of multiple β sheets centrally located, flanked by a single α helix on each side, illustrated in (Fig. 2). Conversely, the PTP domain exhibits a mixed architecture, incorporating both α helices and β sheets, with ten β sheets encircling the αE helix.

Fig. 2

Crystallographic structure of the SH2 domain, featuring a central β-sheet (green) flanked by two α-helices, highlighting its structural organization and functional architecture [46]

The interaction between the N-SH2 and PTP domains forms a polar interface, which is crucial for modulating the phosphatase activity of SHP2, either enhancing (gain of function) or diminishing (loss of function) it. Specifically, the βD’-D’E-βE region within the N-SH2 domain interacts with phosphotyrosyl (pTyr) residues located in the catalytic pocket of the PTP domain, effectively inhibiting the catalytic site. Although the C-SH2 domain connects to either the N-SH2 or PTP domains without significant interaction, its orientation remains stable. The overall structure maintains a closed configuration, resulting in the binding sites for pTyr peptides being fully exposed on the side of the SH2 domains that is opposite to the PTP domain [17, 21].

The structural reorganization of SHP2 resulting from various disease-associated mutations is not unexpected. Subsequently, mutations in SHP2 have been linked to Noonan syndrome (NS) and neoplastic disorders; however, the biochemical mechanisms underlying these effects remain ambiguous. Certain mutations in SHP2 diminish the interaction at the interface between the N-SH2 and PTP domains, which disrupts the formation of the pTyr substrate recognition surface and the process of dephosphorylation, thereby leading to a loss of catalytic activity, exemplified as a loss-of-function (LOF) defect [47]. Conversely, many mutations associated with cancer, such as SHP2-E76D, E76G, and E76A, induce an open conformation. In this state, the C-SH2 domain rotates 120 degrees relative to the PTP domain, resulting in the N-SH2 domain shifting to the opposite side of the catalytic pocket and fully exposing the catalytic site. Thus, the rotation of the C-SH2 domain prompts a rearrangement of the N-SH2 domain, transforming its conformation from closed to open and enhancing the catalytic activity of SHP2, which is characterized as a gain-of-function (GOF) defect [48]. There is significant interest in not only comprehending how frequently observed mutations in oncogenes and tumor suppressors influence protein structure and function, but also in determining the implications of these mutations on the overall physiology of cancer. To gain insights into the effects of gain-of-function (GOF), we conducted an investigation into the structure and function of three cancer variants (E76D, E76G, and E76A) utilizing essential molecular dynamics (MD) simulations. Furthermore, we will examine the dynamic conformational mechanisms at the atomic level to identify the specific alterations occurring in the catalytic loops (i.e., βD’-D’E-βE loop) by comparing them to the Apo-state, as well as how the severity of these mutations correlates with the activation state of SHP2. Consequently, this study aims to elucidate the intrinsic effects of the missense SHP2 mutations located at position 76 (E76D, E76G, and E76A) in alleviating the autoinhibited state of SHP2 through various analyses conducted post-MD simulation.

To assess the structural alteration induced by the E76 mutations, we generated 3D visualizations of the backbone structures for each mutant system (SHP2^E76D, SHP2^E76G, and SHP2^E76A) and superimposed them onto SHP2^WT (Fig. 3). This comparison highlights any significant conformational deviations caused by the mutations. The resulting structure showed that the overall protein fold remained intact, but subtle structural variations were observed near the mutation site. These differences may affect local interactions, flexibility, and potentially influence the protein dynamic behavior. Although no major conformational shifts were detected, the observed deviations could play a role in altering protein stability and function. To further investigate these effects, detailed analyses including RMSD, RMSF, Rg, SASA, distance analysis, and DCCM were performed.

Fig. 3

Superimposed structures of SHP2^WT, and its mutants (SHP2^E76D, SHP2^E76G, and SHP2^E76A)

The stable and flexible dynamic behavior of SHP2^WT and SHP2^Mutants

The SHP2 is widely expressed and is essential for normal developmental processes. It plays a critical role in cellular signaling by regulating protein functions. Deregulation, such as mutations in SHP2 or variations in its expression levels, can lead to various forms of leukemia and related disorders, establishing SHP2 as the first identified proto-oncogene. Notably, most cancer-related mutations in SHP2 occur within the N-SH2/PTP domain, which disrupts the equilibrium between its active and inactive conformations due to allosteric activation and structural characteristics. We conducted a total of 400 ns of MD simulations in an explicit water environment. For the evaluation of dynamics stability and MD simulation quality, we monitored the RMSD of the Cα atoms, comparing simulated snapshots to the original crystal structure throughout the entire MD trajectory, as illustrated in Fig. 4. The curve representing smaller deviations signifies a high level of stability, whereas the opposite curve suggests reduced stability. The RMSD of all states in relation to the original structures demonstrates that a total of 400 ns of MD simulation time is sufficient to reach equilibration at 310 K. The RMSd results indicate that the SHP2 protein exhibits different conformations across all three mutation states (E76D, E76G, and E76A) when compared to the SHP2 Apo-state. The RMSD plots reveal that the conformations of the Apo-state attain equilibrium at approximately above 250 ns, subsequently decreasing with smaller deviations and maintaining stability throughout the simulation period, as illustrated in Fig. 4.

Fig. 4

Illustration of the RMSD graph for Cα atoms (Å) relative to their initial conformation vs simulation time (ns) for all four states during 400 ns MD simulation. The RMSD demonstrate a notable shift after 100-150 ns, i.e., both M2-state and M3-state gradually increase, which indicate the unstability. At the same time, in the case of M1-state, the RMSd curve show smaller increase

The RMSD results for all three mutant states exhibited significant differences. The conformational dynamics of SHP2 across these mutant states may be attributed to the distortion and residual movement of the βD’-D’E-βE loop in various directions, which encodes the N-SH2 domain of SHP2. In the initial phase of the simulation (0-150 ns), the M1-statet (blue) exhibits lower RMSD values compared to the Apo-state, indicating a stable conformation with minimal fluctuations ranging from 2 Å to 3 Å. This suggests that the M1-state may initially stabilize the protein structure. However, after 200 ns, the M1-state displays increased fluctuations, with RMSD values reaching up to 4 Å. This shift indicates a transition to a more flexible or less stable conformation over time. Such behavior could be attributed to the mutation-induced disruption of critical interactions that maintain SHP2's autoinhibited state, leading to increased conformational flexibility, as illustrate in Fig. 4. The gradual increase in the backbone RMSD may be attributed to the weakening interactions of the interface residues, as the mutated residues are situated on the surface between the N-SH2 and PTP domains. Therefore, these observations support the notion that mutations contribute to the destabilization of SHP2, with heightened fluctuations indicating increased catalytic activity, as seen in the initial gain-of-function SHP2 mutants leading to the normal state.

Similarly, the dynamic analysis of the M2-state and M3-state exhibit distinct fluctuation behaviors when compared to the M1-state. At the beginning of the simulation (0–50 ns), both systems show a rapid increase in RMSD, which is expected as the protein adjusts from its initial energy-minimized state to a more relaxed conformation. Around 100–150 ns, the Apo-state of SHP2 (black) stabilizes within the range of ~ 2.5 Å to 4.2 Å, maintaining a relatively stable conformation with moderate fluctuations as shown in Fig. 4. While compare to the Apo-state, the M2- state (Purple) stabilizes earlier, with lower RMSD values (~ 2.0 Å to 3.5 Å), suggesting a more compact structure. This suggests that the Apo-state remains more flexible and dynamic, whereas the M2-state stabilizes an active conformation, leading to reduced fluctuations. As the simulation progresses, the Apo-state exhibits higher fluctuations (~ 4.2 Å) between 200 and 300 ns, indicating dynamic domain rearrangements. In contrast, M2-statet maintains stability within a narrower RMSD range (~ 3.5 Å), suggesting that the mutation leads to a structurally more rigid conformation. By 350–400 ns, both systems appear to have reached equilibrium, as indicated by smaller fluctuations in RMSD. The Apo-state still shows slight variations (~ 3.8 Å to 4.2 Å), while the M2-state remains more stable (~ 3.2–3.5 Å). These observations suggest that the 400 ns simulation is sufficient for reaching equilibrium, as no major conformational changes occur beyond this point.

Similarly, M3-state E76A (red) initially follows a similar stabilization trend but begins to show greater fluctuations after 100 ns, reaching RMSD values above 4.5 Å and peaking around 5.0 Å towards the final stage of the simulation. These fluctuations suggest that M2-state introduces greater structural flexibility, possibly promoting a more dynamic, destabilized state compared to both Apo-state and mutant states. Such as between 200 and 300 ns, the Apo-state maintains its RMSD within the range of 3.2–4.2 Å, indicating a relatively stable structure. However, the M3-state exhibits a gradual increase in RMSD and surpasses the Apo-state in fluctuations beyond 300 ns, reaching values above 4.5 Å. This indicates that the M3-state increases structural instability at later stages of the simulation, potentially promoting a more dynamic or extended conformation. By 350–400 ns, both systems appear to have reached equilibrium, although the M3-state continues to show higher fluctuations (~ 3.5 Å to 5.0 Å), suggesting a more flexible conformation compared to Apo-state. Thus the final 400 ns simulation reveals M3-state as the most dynamic variant, suggesting that it may further disrupt the autoinhibited state of SHP2, potentially leading to hyperactivation. By the end of 400 ns, all three systems appear to have reached equilibrium, as no major RMSD fluctuations are observed. This suggests that the simulation duration is sufficient to capture equilibrium behavior. The key findings indicate that M2-state stabilizes a conformational state that favors activation, whereas M3-state introduces additional flexibility, making it potentially more prone to activation and oncogenic signaling. The results align with previous studies highlighting that the mutation disrupts SHP2 autoinhibited conformation, potentially leading to heightened phosphatase activity and oncogenic potential.

To gain deeper insights into the influence of individual residues, which may shed light on the flexibility of residues in both the Mutant and Apo states, we conducted an analysis of the Cα root-mean-square fluctuations for all side-chain atoms of the proteins in these two states. To enhance the precision of our analysis, we extracted the coordinates of the minimal energy structure from the equilibrium phase. This structure was subsequently aligned and utilized as a reference for calculating the RMSf. The observed fluctuations exhibited a negative correlation with the stability of the residues; specifically, greater fluctuations in residues indicated a lack of stability, and conversely, smaller fluctuations suggested greater stability.

The RMSf analysis provides insights into the local flexibility of individual residues within Apo-state (black) and its mutant state during the 400 ns simulation. The results of RMSf for the cancerous mutant states exhibited significant differences compared to the Apo state, demonstrating pronounced fluctuations, particularly in the βD’-D’E-βE loop region (Fig. 5). In contrast, the Apo-state displayed a stable and consistent fluctuation pattern within the βD’-D’E-βE loop, indicating the closed autoinhibited conformation of SHP2, which is critical for the molecular switching mechanism. The M1-state reveals increased flexibility in particular area compared to the Apo-state. Notably, higher fluctuation were observed in the residue regions 150–200 with peaks up to 6.5 Å, while a sharp fluctuation increase were demonstrated at the region of 228–300, with pronounced peaks reaching up to 7.2 Å, as shown in Fig. 5. These findings suggest that the E76D mutation disrupts stabilizing interactions within SHP2, leading to enhanced flexibility in these regions.

Fig. 5

The RMSF plots for each residue of the Apo-state, M1-state, M2-state, and M3-state obtained from 400 ns MD simulation. Residue fluctuations monitored for the Cα atoms of the protein over the entire trajectory. Each panel represent variants of SHP2 in comparison with SHP2^Apo

While the RMSf analysis of the M2-state and M3-state exhibit higher fluctuations compared to the Apo-state, particularly in specific regions such as βD’-D’E-βE loop critical for SHP2's regulatory function.​ The M2 mutant state generally shows higher fluctuations at several positions, particularly around residues 150–225 and 275–325, with peaks exceeding 6.0 Å, suggesting increased flexibility in these regions (Fig. 5). Notably, at the residue region of 300, some sharp peaks were observed reaching up to 7.0 Å. These regions correspond to key functional domains of SHP2, including the regulatory PTP domain and the C-terminal tail, essential for enzymatic activity and conformational regulation. The increased flexibility in these areas suggests that the M2-mutant state disrupts the autoinhibited conformation of SHP2, leading to an open, active state, consistent with its gain-of-function behavior.​ The similar pattern of fluctuation were also observed for the M3-state, but the M3-state display significantly a more higher peaks as compared to the Apo-state, starting form 4.0 Å and reaching peak above 7.0 Å as shown in Fig. 5. The most pronounced fluctuations are observed, particularly in 150–225 and 228–315 residue regions, with peaks exceeding 6.0 Å, indicating increased flexibility in these regions. This enhanced flexibility suggests that the M3-state also disrupts the autoinhibited conformation of SHP2, leading to increased phosphatase activity and oncogenic potential.

Overall, the RMSF analysis highlights that the most significant fluctuations were noted at the regulatory regions. This indicates that the cancerous mutations (M2 and M3) at position 76 result in the loss of essential interactions that tether the N-SH2 and PTP domains together, thereby disrupting the autoinhibited conformation of SHP2 and promoting an active state by increasing the flexibility of these domains. Thus our MD analysis results reveal that the corresponding regions exhibits dynamic fluctuation behavior across all three mutant states and the Apo state, with stable behavior in the Apo state, less stability in the M1 state, and high instability in the M2 and M3 states. This suggests that certain residues in the corresponding regions of the Apo and M1 states maintain strong interactions with the catalytic cleft of the PTP domain. Conversely, in the M2 and M3 mutations, the residue-residue interactions at the interface were weakened, further increasing the solvent accessibility of the catalytic site, thereby providing clear evidence of the active state of the SHP2 protein. This structural insight provides a deeper understanding of the molecular mechanisms underlying SHP2-related pathologies and may inform the development of targeted therapies.

Exploring Cα-distance analysis of the D′E loop to activate SHP2

To gain a deeper understanding of the interaction effects between the interface residues of the N-SH2 and PTP domains, as well as to clarify the role of the βD’-D’E-βE loop in the activation of SHP2, we conducted a comprehensive analysis of the Cα distances among specific residue pairs, particularly those within the N-SH2 (βD’-D’E-βE loop) and PTP domains. The interface residues play a crucial role in the catalytic function of SHP2. Consequently, any alteration in the Cα distance between the N-SH2 and PTP domains can significantly impact the physiological behavior of the SHP2 protein, in accordance with the "Molecular switching mechanism." This mechanism maintains a close interaction that inhibits the PTP activity of SHP2 by obstructing substrate access, or conversely, it can lead to an increased Cα distance that enhances PTP activity. Additionally, missense mutations can disrupt the formation of the pTyr substrate recognition surface and hinder dephosphorylation, thereby reducing catalytic activity. While many SHP2 mutants exhibit basal activation, they vary considerably in the degree of this activation, which can be influenced by whether the Cα distance between interface residues increases, decreases, or is completely lost, providing direct evidence of loss-of-function (LOF) or gain-of-function (GOF) activity in SHP2. Moreover, there is no definitive biochemical distinction among mutations associated with neurodevelopmental disorders (NDDs), leukemia, or both conditions.

The results of the dynamics discussed above indicate that the βD’-D’E-βE loop exhibited distinct conformations across different states. In the Apo-state, the βD’-D’E-βE loop maintained stability throughout the molecular dynamics simulation, supporting the notion of D′E loop occupancy (switch) over the catalytic site of PTP. This observation suggests a “closed” autoinhibited conformation of the SHP2 protein across all dynamic states. In contrast, the M1-state displayed a behavior somewhat akin to the Apo-state; however, minor fluctuations or instabilities in the residual movements pointed towards an activated state of SHP2. Furthermore, the dynamic analyses of the M2-state and M3-state demonstrated distinct behaviors, with both domains slightly separating from one another. This minor displacement in the D′E loop increased the solvent-accessible surface area (SASA) of the catalytic site of PTP, favoring conformations with an “open” switch and allowing substrate access to the catalytic site, thereby activating the SHP2 protein to further propagate the RAS-RAF signaling pathway. Overall, the dynamic results clearly illustrate the allosteric effects of inducing point mutations, which not only manipulate the D′E loop in relation to the catalytic site but also enhance the Cα-distance between the N-SH2 and PTP domains.

The Cα-distance between the interface residues was systematically calculated, with specific examples including D1 (F7/Q256), D2 (E76/S502), D3 (E76/R265), and D4 (E76/N308). The findings indicate that the residue − residue distance analysis reflects atomic interactions such as hydrogen bonds, Van der Waals forces, ionic bonds, and hydrophobic interactions, demonstrating a consistent pattern in the stabilization of the SHP2 protein against activation (refer to Fig. 6). Overall, the results suggest that while some interface residues exhibited stronger interactions, others diminished, and some maintained weak interactions throughout the molecular dynamics (MD) simulation. In the Apo-state, all interface residues displayed strong interactions and remained stable during the simulation period (see Fig. 6). Conversely, in the M1-state, the distance between F7 (N-SH2 domain) and Q256 (PTP domain) (D1) ranged from 11.0 Å to 12.5 Å, with a noticeable increase over time in the mutant, suggesting weakened domain interaction stability. Similarly, the E76-S502 (D2) interaction, critical for the N-SH2 and C-SH2 domain interface, exhibited a more higher fluctuations between 13.0 Å and 14.0 Å. While both Apo-state and mutant proteins maintained this interaction, the M1-state showed increased distance between certain pairs of residues. Such the E76-N308 (D4) interaction, another crucial regulatory interface, displayed larger fluctuations, although interactions persisted without complete loss, as illustrated in (Fig. 6).

Fig. 6

The dynamic Cα-Distance Analysis of the Interface pair-wise residues, and were represented by D1 (F7/Q256), D2 (E76/S502), D3 (E76/R265), and D4 (E76/N308) respectively. The less the distance, the more stable the binding should be, and vice versa

In contrast, the cancer-associated mutant states (M2 and M3) exhibited distinct distance patterns compared to the Apo-state and M1-state. The M2 mutant showed a decreased distance between certain residue pairs relative to M1 and M3, suggesting a tighter domain conformation. However, in both M2 and M3 mutants, the distance between F7/Q256 (D1) increased after 150 ns, reaching 8.0 Å to 10.0 Å, further indicating weakened domain interaction stability. The E76-S502 (D2) interaction showed differing trends between M2 and M3: while M2 exhibited a slightly reduced distance (8.5 Å to 10.2 Å), M3 demonstrated greater fluctuations (9.5 Å to 11.0 Å). Similarly, consistently higher distances were observed between key interface residues, particularly E76/R265 (D3) and E76/N308 (D4), which ranged from 13.8 Å to 19.0 Å (Fig. 6).

These findings highlight that the mutations significantly alter interdomain flexibility, particularly at the N-SH2/PTP and N-SH2/C-SH2 interfaces, potentially disrupting SHP2’s allosteric regulation. Notably, the E76/S502 interaction formed a robust hydrogen bond in the wild-type protein, but the Cα-distance analysis revealed that this bond was compromised in all mutants due to an increased bond distance. This disruption promotes N-SH2 domain inhibition, preventing substrate access to the catalytic domain's active site. Consequently, the mutations caused a significant displacement of the N-SH2 domain from the PTP domain, though the overall interdomain distance remained largely unchanged. This variation may be attributed to extensive water-mediated hydrogen bonding, which connects E76 to distant residues within the PTP catalytic pocket. The perturbation of E76 likely further destabilizes this interaction, underscoring the structural and functional consequences of these mutations.

The interface mutant adjust SASA of the catalytic site

The SASA refers to the surface area that is traced by the center of a probe sphere as it rolls over a molecule, where the atoms are represented as spheres of varying radii. SASA serves as a metric to quantify the fraction of a protein's surface that interacts with solvent molecules, and it can be utilized to predict the degree of conformational changes that occur during the binding process, whether between domains or between proteins and ligands. This surface is defined by an imaginary solvent sphere surrounding a protein, in conjunction with the van der Waals contact surface of the molecule. The method for calculating SASA was developed by Miller and colleagues in 1987, allowing for the assessment of polar, nonpolar, and total molecular surfaces. In the context of protein folding and stability research, SASA has consistently been regarded as a crucial factor. In this study, we examined the SASA measurements for the interface residues of the SHP2 protein, which include the N-SH2 and PTP domains, to elucidate the significant movements of the βD’-D’E-βE loop across various states (Apo-state, M1, M2, and M3 states). Additionally, we aimed to investigate whether mutations affect the stability or folding of the SHP2 protein by altering its SASA. The results indicated notable differences in SASA for the SHP2 protein across all systems, with a particularly significant increase observed in the βD’-D’E-βE loop region, providing clear evidence regarding the activation state of SHP2. In the Apo-state, the protein's native structure is compact, characterized by strong interactions among residues, resulting in low solvent accessibility with minimal fluctuations throughout the simulation as illustrated in Fig. 7. The findings demonstrate that the loop (βD’-D’E-βE) of the N-SH2 domain is tightly bound to the catalytic cleft of the PTP domain, indicating a buried surface with low SASA in relation to the solvent. In the M1 state, a significant increase was observed, suggesting the release or partial release of the loop from the catalytic cleft, allowing it to move freely and become more exposed to the solvent (see Fig. 7).

Fig. 7

Solvent Accessible Surface Area of the Interdomain interface residues between N-SH2 and PTP domains

In all cancerous mutations (M2, M3), the solvent-accessible surface area (SASA) consistently increased throughout the simulation period, indicating a more solvent exposed conformation. Initially, the solvent accessibility was low, up to 100 ns, but exhibited sudden oscillations between both domains after 150 ns, as illustrated in Fig. 7. The M2-state (purple) shows increased SASA values, suggesting partial destabilization or loosing of the protein structure, which may enhance solvent accessibility and lead to altered functional properties. In a similar way, the M3-state (red) display an even greater increase in SASA value, implying a more pronounced structural perturbation compared to Apo-state and other mutant state. The increase solvent exposure observed in both mutants suggest that E76G and E76A substations impact the overall protein folding and stability. The differences in SASA may occur when the mutated residues is solvent-accessible, and the nature of the amino acid substitution also influences the extent of structural alterations. The E76A mutation, which replaces glutamic acid with a non-polar alanine, may disrupt hydrogen bonding interactions more severely than the glycine substitution, leading to greater unfolding or increased flexibility. While in the absence of mutations, hydrophobic forces maintain the protein's structure by sequestering hydrophobic groups away from water into a solvent-protected hydrophobic core. This structural perturbation could potentially impact SHP2’s functional interactions with other biomolecules, altering its activity and regulatory mechanisms. Overall, the SASA analysis highlights that all the mutations increase solvent accessibility, suggesting potential destabilization or conformational changes in the protein structure.

The compactness of the SHP2 protein was further assessed using the Radius of Gyration (Rg) across all states. The Rg serves as a metric for evaluating the structural compactness of biomolecules. Given the observed fluctuations and deviations in the backbone, a more detailed analysis of the overall compactness in both the Apo-state and Mutant-state was warranted. The Rg analysis revealed distinct patterns of compactness for all states, as illustrated in Fig. 8. Notably, the Apo-state exhibits a peak in Rg distribution around 26.0 Å, with a relatively narrow spread. This suggests that the Apo protein maintains a stable and compact conformation throughout the simulation, which is consistent with its native fold. The limited fluctuations in Rg values indicate that the Apo-state remains structurally intact with minimal unfolding or expansion. In contrast, the M1-state displayed a slight reduction in compactness as compared to the Apo-state, with a notable peak around 25.9 Å over the course of the MD simulation, suggesting an expansion of the closed/native structure of SHP2, particularly in the interface region. The altered compactness due to the E76D mutation may affect SHP2 interaction with substrates and regulators, potentially leading to dysregulated signaling pathways. Such perturbations are associated with developmental disorders like Noonan syndrome.

Fig. 8

Superposed plot for radius of gyration (Rg), in all four state

The cancerous mutants (M2 and M3) exhibited a marked difference in Rg values compared to the Apo-state, as depicted in Fig. 8. Such as M2-state shows a shift in Rg distribution towards lower values, with a peak around 25.8 Å. This decrease in Rg suggests that the M2-state results in a slightly more compact structure compared to the Apo-state. The increased frequency of lower Rg values implies that the protein might be experiencing structural tightening or localized conformational changes that reduce its overall dimensions. This could be due to altered intramolecular interactions or increased rigidity within certain regions of the protein. Similarly, the M3-state exhibits a broader Rg distribution with a shift towards higher values, peaking slightly above 26.0 Å and extending beyond 26.4 Å. This indicates that the M3-state leads to increased structural flexibility and expansion. This mutation likely disrupts hydrogen bonding or electrostatic interactions, leading to a looser protein structure with increased solvent exposure.

These findings indicate that the induced cancerous mutations lead to a decrease in the overall compactness of the SHP2 protein, revealing instability among the residues, which are no longer in close contact and are free to move. The M1 and M2 states induces a more compact conformation, while the M3-state promotes a more expanded and flexible structure. These findings align with the SASA analysis, where M3-state demonstrated greater solvent exposure, further supporting the notion that this mutation destabilizes the protein overall fold. The Apo-state remains in a balanced structural state, whereas the mutants display distinct shifts in compactness, potentially affecting their functional dynamics and interactions with other molecules. Overall, the Rg analysis highlights how point mutations at E76 differentially affect the global structural properties of SHP2. These alterations may have significant implications for SHP2 regulatory mechanisms, potentially influencing its role in cellular signaling pathways.

Mutations weaken the correlated motion of the interface residues

The Dynamic Cross-Correlation Matrix (DCCM) analysis was employed to examine the correlation motion patterns both within and between monomers across all systems. The figures provided represent the DCCM plots for Apo-state and mutant states, allowing us to compare how these mutations affect the protein’s internal motions. The x-axis and y-axis of the plots both represent the residue index of the SHP2 protein. Each axis spans from residue 1 to approximately residue 530, covering the entire length of the protein sequence. Every point (i, j), in the matrix corresponds to the correlation between the fluctuations of residue i and residue j over the course of the simulation. If a pair of residues move in a coordinated manner, they will have a positive correlation (green regions), while residues moving in opposite directions will show a negative correlation (brown/black regions). The diagonal line running from the bottom-left to the top-right represents the correlation of each residue with itself, which is always 1.0 (perfect correlation). The color scale on the right side of the plots quantifies the degree of correlation as illustrated in (Fig. 9).

Fig. 9

The DCCM shows the collective atom fluctuations for all four states with the correlated motions of protein residues. The (red-yellow) color represents the positive correlation, white color represents the local displacement, and the (black-dodger blue) color represents the negative correlation

The DCCM analysis of all the mutant and Apo-state were compared, to find correlated motions between different residues in a protein. The M1-state demonstrate a more structured and continuous correlation pattern, indicating that M1-state preserves the natural dynamic behavior of SHP2 better than Apo-state. While there are still some anti-correlated regions, they are less extensive, and the overall correlation network is more stable. This suggests that replacing glutamate (E) with aspartate (D), which retains a negative charge, has a less disruptive effect on protein dynamics, illustrated in Fig. 9. The M1 mutation preserves much of the native correlation pattern, suggesting that the structural and dynamic effects of this mutation are milder.

In contrast the cancerous mutants (M2, M3) illustrated diverse correlation motion as compared to the M1 and Apo-state. The M2-state correlation pattern appears weaker overall, as indicated by the lower intensity of green regions along the diagonal. This suggests a loss of coordinated motion, meaning that structural stability and inter-domain communication may be disrupted. The fragmented pattern in the lower-left quadrant (residues 1–200) suggests weaker intra-domain interactions compared to Apo-state M1-state, shown in Fig. 9. The presence of fewer strongly correlated (green) and anti-correlated (brown) regions, indicating that structural changes caused by M2-state may lead to a more flexible but less functionally stable protein. Similarly, the M3-state shows a correlation map that is more fragmented, with disrupted regions of both positive and negative correlations compared to both Apo-state and mutant states. The anti-correlated regions (brown/black patches) are more pronounced in M3-state distributed across the structure, particularly in the 200–300 residue range, suggesting that this mutation significantly alters the natural movement of SHP2, as given in Fig. 9. The green regions are less continuous, meaning that coordinated motions are weakened, which could impact the protein's structural stability and function. There are noticeable patches of anti-correlated regions (brown/black) suggesting disrupted long-range interactions. Overall, the results from the DCCM analysis highlight the importance of residue E76 in SHP2 function. The E76A mutation causes a significant disruption in residue-residue correlations, likely due to the loss of electrostatic interactions. This could affect allosteric regulation, where conformational changes in one part of the protein influence its catalytic activity or interactions with binding partners. In contrast, the other mutant preserves much of the native correlation pattern, suggesting that the structural and dynamic effects of this mutation are milder. These findings are relevant in the context of SHP2-associated diseases such as Noonan Syndrome and certain cancers, where mutations can alter the protein’s regulatory function. Thus the following study provide a deeper understanding of how specific mutations influence protein dynamics, which could be valuable for biomedical research and drug development.

Conclusion

SHP2, a classical non-receptor protein tyrosine phosphatase (NRPTP), plays a crucial role in cell proliferation, differentiation, survival, and apoptosis by modulating major signaling pathways such as Ras/Raf/MAPK, PI3K/Akt, and JAK/STAT. Dysregulation of SHP2 through gain-of-function (GOF) or loss-of-function (LOF) mutations has been implicated in Noonan syndrome (NS) and various cancers, including leukemia. This study explores the structural and dynamic consequences of three GOF SHP2 mutations (E76D, E76G, and E76A) using molecular dynamics (MD) simulations. These mutations, located at the interface between the N-SH2 and PTP domains, disrupt the autoinhibited conformation of SHP2, leading to an increase in PTP activity and persistent activation. The E76D mutation, associated with Noonan syndrome, induces moderate structural changes, while the E76G and E76A mutations, linked to cancer, cause significant disruption of key interdomain interactions. Our dynamic analysis reveals that the cancer-associated mutants (E76G and E76A) exhibit higher fluctuations and increased Cα-distances between interface residues, weakening key hydrogen bonds (E76/N308, F7/Q256, E76/S502, and E76/R265). This structural destabilization releases the βD’-D’E-βE loop, promoting an open, active conformation that enhances SHP2 activity. In contrast, the E76D mutation (M1-state) shows relatively lower fluctuation, maintaining a degree of stability while still disrupting the interface interactions. These findings provide compelling evidence that oncogenic E76G and E76A mutations drive hyper activation of SHP2, contributing to pathological signaling in cancer. The insights gained from this study emphasize the structural basis of SHP2 activation and offer potential targets for therapeutic intervention in diseases driven by SHP2 dysregulation.