Introduction

Eukaryotic protein kinases (EPKs) are plastic enzymes of paramount importance to signaling processes, catalyzing phosphoryl transfer reactions, or acting as scaffolds for other enzymes and/or binding partners ¹. Of all kinases, the catalytic subunit of PKA (PKA-C) was the first to be structurally characterized by X-ray crystallography.^2,3 In its inhibited state, PKA-C forms a heterotetrametric complex comprising two catalytic (C) and two regulatory (R) subunits.⁴ The canonical activation mechanism of PKA involves binding two cAMP molecules and disassembling the holoenzyme, which unleashes active PKA-C monomers that target signaling partners.⁵ In 1997, however, Scott and coworkers suggested that the holoenzyme does not disassemble under physiological conditions; rather, the intact holo-enzyme forms signaling islands localized by A-kinase anchoring proteins (AKAPs) in the proximity of substrates.⁶

X-ray crystallography studies revealed that PKA-C is a bilobal enzyme, with a dynamic N-terminal lobe comprising β-sheets and the αC helix and a more rigid C-terminal lobe, with mostly α-helices (Figure 1).^2,3,7 The N-lobe harbors the nucleotide-binding site, whereas the substrate binding cleft lays at the interface between the N-lobe and C-lobe. The three-dimensional structure of the enzyme features a highly conserved hydrophobic core decorated by catalytically important motifs, i.e., the Gly-rich loop, DFG-loop, activation loop, positioning loop, and magnesium and peptide positioning loops.⁸ In the catalytically active state, these motifs are all poised for phosphoryl transfer. However, the conformation of these motifs is necessary but not sufficient to define an active kinase. More recent studies revealed a critical role of the hydrophobic core, which is crossed by a catalytic (C) spine, a regulatory (R) spine, and shell residues (Figure 1).^9,10 The C spine comprises an array of hydrophobic residues and is assembled upon binding ATP, whereas the R spine is engaged when the activation loop is phosphorylated, which contributes to the positioning of the αC-helix.¹¹

Figure 1.

A distinct property of PKA-C is the binding cooperativity between ATP and substrate.¹² During the catalytic cycle, the kinase recognizes and binds substrates with positive binding cooperativity between ATP and unphosphorylated substrates, whereas a negative binding cooperativity between ADP and phosphorylated substrate characterizes the exit complex.¹³ The biological importance of binding cooperativity has been emphasized by our recent studies on disease-driven mutations of PKA-C,^14-16 which all feature disrupted cooperativity between nucleotides and protein kinase inhibitor (PKI) or typical substrates.^14-16 Since the recognition sequence of the substrate is highly homologous to the regulatory subunits,⁴ a loss of cooperativity of binding may affect not only substrate binding fidelity but also the regulation by the R subunit. However, the molecular determinants for the binding cooperativity between ATP and substrate and its role in PKA signalosome remain elusive to date.

Here, we combined NMR-restrained replica-averaged metadynamics (RAM)^17,18 and Markov State Model (MSM)¹⁹ to define the conformational landscape and the corresponding dynamics of PKA-C. We found that both the apo kinase and the nucleotide-bound kinase occupy three distinct basins: (i) a most populated ground state with constitutively active conformations competent for catalysis, (ii) a first high free energy basin representative of typical inactive states with a dislodged configuration of the αC helix, and (iii) a second high free energy basin with a disrupted hydrophobic array of residues at the core of the enzyme. Notably, the equilibrium between the most populated ground state and the other low-populated states agrees with previous NMR relaxation dispersion and CEST measurements.²⁰ To compound the existence of the second inactive basin, we mutated F100 at the tip of the αC-β4 loop into Ala (PKA-C^F100A) and characterized its binding thermodynamics and ligand binding response by isothermal titration calorimetry and NMR spectroscopy. We found that PKA-C^F100A phosphorylates canonical peptides and preserves the binding affinity for both nucleotide and substrate. However, PKA-C^F100A lacks binding cooperativity between nucleotide and substrate. NMR revealed that perturbing the hydrophobic packing around the αC-β4 loop interrupts the allosteric communication between the two lobes of the enzyme. These results further support the pivotal role of the αC-β4 loop in kinases function and may explain why single-site mutations or insertion mutations in the homologous kinase that stabilize this motif result in oncogenes and confer different drug sensitivity.²¹

Results

The free energy landscape of PKA-C charted by NMR-restrained replica-averaged metadynamics (RAM) simulations

To characterize the experimentally accessible conformational landscape of PKA-C in the µs-to-ms range, we performed NMR-restrained metadynamics simulations within the RAM framework.^17,18 We simulated the apo, binary (PKA-C/ATP), and ternary (PKA-C/ATP/PKI5-24) forms of PKA-C, using four replicas restrained using backbone chemical shifts (CS) to ensure a close agreement with conformational space explored by the kinase under the conditions used in this work. The enhanced sampling was achieved through bias-exchange metadynamics along different collective variables (CVs) to boost the conformational plasticity of the enzyme (Figure 2 - figure supplement 1).²² Back-calculation of the CS using the Sparta+ software²³ shows that the CS restraints improved the agreement between back-calculated and experimental CS. Specifically, we obtained an overall improvement of ∼0.2 ppm on the amide N atom and ∼0.1 ppm on the remainder backbone atoms (Figure 2 - figure supplement 2). The bias-exchange metadynamics allowed for each replica to span a broader conformational space relative to classical simulations (Figure 2 - figure supplement 3). The deposition of Gaussian biases required by the metadynamics approach converged after 300 ns, where the fluctuations along the first 3 CVs were less than 1 kcal/mol (Figure 2A, Figure 2 - figure supplement 4). The converged bias revealed that the apo kinase accesses multiple minima along each CV, and the conformational heterogeneity is significantly reduced upon ligand binding, especially in the ternary form (Figure 2A). The full free energy landscape was then reconstructed by sampling an extra 100 ns production phase with reduced biases along each CV. The free energy landscape shows how the population of the conformers is modulated by ligands within the NMR detection limit of sparsely populated states (∼ 0.5% or ΔG < 3.2 kcal/mol).

Figure 2.

According to these simulations, apo PKA-C populates preferentially a ground state and five readily accessible low-populated excited states (Figure 2B, Figure 2 – supplementary table 1). The nucleotide-bound PKA-C (binary form) features a similar ground state and a broad higher energy basin (Figure 2C, Figure 2 – supplementary table 1). Finally, the ternary complex occupies a narrow dominant ground state (Figure 2D, Figure 2 – supplementary table 1). This free energy landscape obtained from the RAM simulations is consistent with the qualitative picture previously inferred from our NMR spin relaxation experiments,²⁴ while providing a detailed structural characterization of the excited states.

MSM reveals the conformational transitions of PKA-C from ground to high free energy (excited) states

To explore the conformational transitions of the kinase upon ligand binding, we performed additional unbiased sampling to build a Markov State Model (MSM). MSMs are commonly used to describe the dynamic transitions of macromolecules in terms of (1) the probabilities of occupation of a specific set of states and (2) the transition probabilities of moving between these states. In practice, a MSM is typically created by combining thousands of short unbiased simulations. ^25,26 Following this strategy, we performed several short simulations (10 – 20 ns) using thousands of the low free energy conformations (ΔG < 3.2 kcal/mol) chosen randomly from the three forms of PKA-C as starting structures. The conformational ensembles were clustered into microstates and seeded to start a second round of adaptive sampling (see Methods). The iterative process was repeated three times to assure convergence and yielded a total of 100 μs trajectories for both the apo and binary forms, whereas for the less dynamic ternary complex, we collected trajectories of 60 μs. Once we reached a sufficient sampling, we built a MSM including L95, V104, L106, M118, M120, Y164, and F185 to investigate the dynamic transitions of the hydrophobic R spine and shell residues (Figure 3 – figure supplement 1). These residues are ideal reporters of the dynamic processes governing the activation/deactivation of the kinase.²⁷ To compare the free energy landscape of different complexes, we first projected the conformational ensembles of three forms and existing crystal structures along the first two time-lagged independent components (tICs) of the apo form, which were obtained by a time-lagged independent component analysis (tICA) analysis (see Methods). These tICs represent the directions of the slowest motion of the kinase and visualize the conformational transitions of V104, L95, and F185 (Figure 3 – figure supplement 1). Apo PKA-C accesses three major basins. The broadest basin represents the ground state (GS) and represents the conformations of the kinase captured by essentially all crystal structures (Figure 3A). Additionally, there are two distinct excited states: the first excited state (ES1) features a disrupted hydrophobic packing of L95, V104, and F185. This conformation features an inactive state with an orientation of the αC-helix typical of the inhibited states as found for the PKA-C bound to regulatory subunits RIα and RIIβ. The second excited state (ES2), to our knowledge, was never captured in crystal structures. The ES2 state displays a flipped configuration of the V104 side chain and a rearrangement of the αC-β4 loop, with the dihedral angles of F100 and F102 adopting a gauche⁺ configuration. This orientation of the F100 and F102 aromatic rings breaks the hydrophobic packing of the αC-β4 loop with the C-lobe, causing steric contacts between F100 and V104 (Figure 3A). In contrast, the active GS ensemble features a trans configuration of the F100 and F102 side chains that stabilizes the hydrophobic interactions with αE- and αJ-helices in the C-lobe (Figure 3A). Upon binding ATP, the conformational space span by the kinase becomes narrower, and the conformers populate mostly the GS, with a small fraction in the ES1 state (Figure 3B). This is consistent with the role of the nucleotide as an allosteric effector, enhancing the affinity of the enzyme for the substrate. In the ternary form (ATP and PKI-bound), PKA-C populates only the GS consistent with the competent conformation observed in the first ternary complex structure. In this case, the αC-β4 loop of the enzyme is locked in a well-defined configuration, as shown in Figure 3A.

Figure 3.

In the apo form, the αC-β4 loop is quite dynamic due to transient hydrophobic interactions between F100 and V104 as well as W222 and the APE motif (A206 and P207) (Figure 4A). The binding of both nucleotide and PKI increases the rigidity of residues near F100 and V104 such as V103 and F185 (Figure 4A). In addition, several electrostatic interactions essential for catalysis (D166-N171, K168-T201, and Y204-E230), which are transient in the apo PKA-C, become more persistent (Figure 4B).

Figure 4.

The MSM makes it possible to use kinetic Monte Carlo sampling to characterize the slow transition between different states.²⁸ In the apo form, the GS features F100 and F102 in trans, a configuration that stabilizes the interactions with the αE and αJ helices and, together with the nucleotide, locks the αC-β4 loop to elicit an active kinase conformation (Figure 3A). The GS to ES1 transition features the disruption of the K72-E91 salt bridge, whereas the GS to ES2 transition involves a 120° flip of the F100 aromatic group that interacts with V104, a conformation found only in the uncommitted apo enzyme (Figure 3A). The GS to ES2 transition involves a concerted disruption of the D166-N171 and K168-T201 electrostatic interactions, essential for catalysis. Also, this event destabilizes the packing between W222 and the APE motif (A206 and P207) required for substrate recognition (Figure 3C). All these conformational transitions suggest that ES1 and ES2 represent destabilized states of the kinase.

Direct correspondence between the conformationally excited states identified by MD simulations and NMR

NMR relaxation dispersion and CEST experiments performed on the apo PKA-C revealed the presence of conformationally excited states for several residues embedded into the hydrophobic core of the enzyme.²⁰ Based on the free energy landscape and the conformational transitions in the core of the kinase, we reasoned that the inactive states identified by the simulations would correspond to a disruption of the hydrophobic packing near the αC-β4 loop. To test this hypothesis, we calculated the values of the change in CS (Δω) for the methyl groups of L103, V104, I150, L172, and I180 from the MD simulations and compared them with those obtained from the fitting of the CPMG dispersion curves and CEST profiles.²⁰ We first sampled more than 500 snapshots as representative structures of the ES2 and GS states (Figure 3 – figure supplement 2), and computed the distribution of the methyl ¹³C chemical shifts on these sites using ShiftX2.²⁹ These calculations yield Δω values of 0.87 ± 0.02 and 0.85 ± 0.02 ppm between ES2 and GS for Val104 and Ile150, respectively. These values are in good agreement with the experimental Δω values (1.10 ± 0.06 and 0.94 ± 0.09 ppm) obtained from the fitting of the CPMG dispersion curves (Figure 5A-B). The CS difference obtained for the other three sites also yielded a good agreement, with the same direction of the chemical shift changes found experimentally (Figure 5 – figure supplement 1). Remarkably, the fitting of the calculated vs. experimental Δω values gives a slope of 0.86 with R² of 0.82 (Figure 4C) supporting that the excited state observed experimentally may correspond to the calculated ensemble. In addition, MSM estimates a population of ES2 of 6 ± 2%, consistent with the 5 ± 1% reported in NMR experiments.²⁷ These results support that the alternate conformations identified by NMR represent inactive states and may correspond to the disruption of the hydrophobic anchor of the αC-β4 loop via F100 and F102.

Figure 5.

The F100A mutation disrupts the allosteric network of the kinase

The above analysis suggests that F100 located at the αC-β4 loop is critical for the GS to ES1 transition path of the kinase. Therefore, we modeled the F100A mutant and analyzed its dynamic trajectories. First, we performed a short equilibration using classical MD simulations starting from the coordinates of the X-ray structure of the ternary complex (PDB ID: 4WB5). During 1 μs of MD simulations, the αC-β4 loop of the F100A mutant undergoes a significant motion as manifested by the increased values of the backbone rmsd and the conformational transition (flipping) of the F102 side chain (Figure 6A). This region in the WT PKA-C adopts a stable β-turn in the WT, with a persistent H-bond between the backbone oxygen of F100 and the amide hydrogen of L103. In the simulation of F100A, H-bond is formed more frequently between A100 and F102, and this region adopts a γ- turn conformation. Such local rearrangement disrupts the hydrogen bond between N99 and Y156 of αE, altering the anchoring of the αC-β4 to the αE helix (Figure 6B).

Figure 6.

In the F100A mutant, the local changes caused by the αC-β4 loop structural transitions propagate and alter the response of the nucleotide binding as shown on the rmsd changes for key hydrophobic motifs such as the R- and C- spines and shell residues (Figure 7A). While the nucleotide binding decreases the average RMSD for these hydrophobic motifs for WT PKA-C, it stabilizes only the C-spine and fails to drive the R-spine and shell residues to an intermediate state competent for substrate binding (Figure 7B). The latter can be attributed to the perturbation of hydrophobic packing of L95, L106 of R-spine, and V104 of the shell residues close to the mutation site. Using the lowest principal components, we also analyzed the global dynamic response to ATP binding. Not only does F100A change the breathing mode of the two lobes (PC1), but it also alters the shearing motion (PC2) of the kinase, emphasizing the importance of this allosteric mutation on the internal communication across the hydrophobic core (Figure 7C-E).

Figure 7.

In fact, these altered motions have global repercussions on the allosteric networks as detected by mutual information for rotameric states (Figure 8). For PKA-C^WT, there are numerous correlations observable within both lobes of the enzyme, involving the Gly-rich loop, αC-β4 loop, activation loop, and the C-terminal tail These dynamic correlations between the two lobes are the hallmark for the dynamically-committed state. In contrast, the mutual information for F100A mutants shows significantly less correlations in the N-lobe and the correlations within the C-lobe are redistributed.

Figure 8.

Taken together, these simulations on F100A suggest that the local perturbation of the αC-β4 loop diminishes the structural connection between the two lobes, disrupting the correlated shearing motion highlighted in other simulation studies. Although the single F100A mutation does not drive the kinase into a completely dysfunctional state (ES2), it is sufficient to abolish the binding cooperativity of the enzyme.

The F100A mutation disrupts nucleotide-substrate binding cooperativity in PKA-C

Based on the results from the MD simulations, we hypothesized that the disruption of the coupling between the N-lobe and C lobe would affect the binding cooperativity between nucleotide (ATPγN) and pseudosubstrate inhibitor (PKI5-24) binding. We first evaluated the catalytic efficiency for the wild-type and F100A mutant by carrying out steady-state coupled enzyme assays³⁰ using the standard substrate, Kemptide.³¹ F100A showed a slight increase in KM and Vmax compared to PKA- C^WT, resulting in a reduction of the catalytic efficiency (kcat/KM = 0.50 ± 0.04 and 0.41 ± 0.08 for PKA-C^WT and PKA-C^F100A, respectively; Supplementary table 2). We then performed isothermal titration calorimetry (ITC)³² to obtain ΔG, ΔH, -TΔS, Kd, and determine the cooperativity coefficient (σ) for ATPγN and PKI5-24 binding. We first analyzed the binding of ATPγN to the apo PKA-C^F100A and, subsequently, the binding of PKI5-24 to the ATPγN-saturated PKA-C^F100A (Supplementary Table 3). We found that PKA-C^WT and PKAC^F100A have similar binding affinities for ATPγN (Kd = 83 ± 8 μM and 73 ± 2 μM, respectively). However, in the apo form F100A showed a 3-fold higher binding affinity for the pseudosubstrate relative to PKA-C^WT (Kd = 5 ± 1 μM and 17 ± 2 μM, respectively - Supplementary table 3). Upon saturation with ATPγN, PKA-C^F100A displayed a 12- fold reduction in binding affinity for PKI5-24, resulting in a σ of ∼3, a value significantly lower than the wild-type (σ greater than 100).

NMR mapping of nucleotide/pseudosubstrate binding response

To understand the atomic details of the disrupted binding cooperativity for PKA-C^F100A, we used solution NMR spectroscopy and mapped the amide fingerprints of PKA-C^F100A using [¹H, ¹⁵N]-WADE-TROSY-HSQC experiments ²⁰ upon adding nucleotide (ATPγN) and pseudosubstrate (PKI5-24). Specifically, we monitored the ¹H and ¹⁵N chemical shift perturbations (CSPs, Δδ) of the amide fingerprints for the PKA-C^F100A saturated with ATPγN and in complex with ATPγN/PKI5-24 and we compared them with PKA-C^WT (Figure 9 – figure supplement 1). ATPγN binding causes the dramatic broadening of several amide resonances throughout the structure of PKA-C^F100A, suggesting the presence of an intermediate conformational exchange ³³. The fingerprints of both PKAC^WT and PKA-C^F100A show similar CSP profiles upon binding ATPγN (Figure 9 – figure supplement 1B), although the extent of the changes is significantly attenuated for PKA-C^F100A. The reduction of CSP is apparent for residues in the N-lobe (β2- β3 region), around the mutation (β4), and at the C-terminal tail. A similar pattern is observed upon binding PKI5-24 to the ATPγN-saturated PKA-C^F100A (Figure 9). A substantial decrease in CSP is observed for the C-lobe near the domains critical for substrate-binding (i.e., αE, αF, and αG).

Figure 9.

To determine the global response to ligand binding for WT and F100A, we analyzed the chemical shifts of the amide fingerprints of the two proteins using the COordiNated ChemIcal Shifts bEhavior (CONCISE).³⁴ CONCISE describes the global response of the protein fingerprint resonances by providing the probability density (population) of each state along the conformational equilibrium for binding phenomena that follow linear chemical shift trajectories. For PKA-C^WT, both nucleotide and pseudosubstrate binding shift the overall population of the residues from an open to an intermediate and a fully closed state (Figure 10 – figure supplement 1). A similar trend is observed for PKA-C^F100A, though the probability densities are broader, indicating that the amide resonances follow a less coordinated response.³⁴ Also, the maximum probability density for the closed state shows that the ternary complex of the mutant is slightly more open than the corresponding wild-type (Figure 10 – figure supplement 1). Overall, the shape of the probability distributions of the amide chemical shifts suggests that several residues do not respond in a coordinated manner and that PKI5-24 binding shifts the conformation of the kinase toward a partially closed state, which may explain the loss in binding cooperativity as previously observed.^14-16

To analyze these coordinated chemical shift changes in detail and define the allosteric network of the kinase upon binding ligands and substrate, we utilized the CHEmical Shift Covariance Analysis (CHESCA),^15,16,20 a statistical method that identifies correlated responses of residues to a specific perturbation. CHESCA works under the assumption that pairwise correlated chemical shift changes between residues identify possible allosteric networks.^35,36 We found that coordinated structural rearrangements, as identified by CHESCA, are a strong indication of binding cooperativity in PKA-C.^15,16,20 Therefore, we compared the CHESCA maps for PKA-C^WT and PKAC^F100A for four different states: apo, ATPγN-, ADP-, and ATPγN/PKI5–24-bound. The CHESCA matrix of PKA-C^F100A exhibits sparser and more attenuated (i.e., lower value of correlation coefficient) correlations relative to PKA-C^WT (Figure 10). Although many inter-lobe correlations are still present for the mutant, several correlations in specific structural domains such as the αG-, αH-, and αI-helices are absent or attenuated. For instance, the F100A mutation does not display correlations between the αA-helix and the C-terminal tail. We also utilized the CHESCA maps to assess the allosteric communication among the communities as defined by McClendon et al.³⁷ Using the community definition, the CHESCA map for PKA-C^WT shows a strong correlation across the enzyme, especially for structurally adjacent communities and at the interface between the two lobes (see for instance the correlations among ComA, ComB, ComC, ComE, and ComH). The community CHESCA for F100A shows that the cross-talk between the nucleotide-binding (ComA) and positioning of αC-helix (ComB) is preserved, as well as the coupling between the R-spine assembly (ComC), and the activation loop (ComF). Except for ComD, the correlations among these communities in the N-lobe and ComF1, ComG, and Com H are completely abolished. Note that ComF1, ComG, and ComH are involved in the substrate binding and are instrumental in docking the R-subunits. Overall, the CHESCA analysis shows that the reduced degree of cooperativity we determined thermodynamically for PKA-C^F100A corresponds to a decrease in coordinated structural changes upon ligand binding, and these changes affect the C-lobe communities involved in the substrate binding and protein/protein interactions.

Figure 10

Discussion

A distinct feature of PKA-C is the binding cooperativity between nucleotide and substrate that originate from the allosteric coupling between the nucleotide-binding pocket and the interfacial region between the two lobes that harbors the substrate binding cleft. Structural and dynamic NMR data suggested that a well-tuned coupling between the two lobes of the kinase is required for efficient substrate binding. Also, additional NMR and functional studies showed that mutations in the activation loop linked to Cushing’s syndrome reduce drastically substrate binding affinity, and more importantly, reduce the communication between the ligand binding pockets.^14-16 Interestingly, a mutation (E31V) distal from the active site and linked to the progression of Cushing’s syndrome has a similar effect, suggesting a possible allosteric modulation of the kinase substrate recognition.¹⁶ NMR chemical shift perturbation data suggested that these mutations are connected to allosteric nodes that, once perturbed, radiate their effects in the periphery of the enzyme and prevent an efficient dynamic coupling between the two lobes of the enzyme. Interestingly, these mutations do not prevent Kemptide phosphorylation. Rather, they cause a loss of substrate fidelity with consequent aberrant phosphorylation of downstream substrates.^38-43 Additionally, thermodynamics and recent NMR studies using different nucleotides and inhibitors demonstrated that it is possible to control substrate binding affinity by changing the chemistry of the ligand at the ATP binding pocket. Altogether, these studies suggest that phosphorylation reaction and binding synergy between ATP and substrates may be controlled independently.

Our NMR data combined with RAM simulations and MSM enabled us to comprehensively map the free energy landscape of PKA-C in various forms. We found that the active kinase unleashed from the regulatory subunits occupies a broad energy basin (GS) that corresponds to the conformation of the ternary structure of PKA-C with ATP and pseudosubstrate (PKI5-24) that exemplifies a catalytic competent state, poised for phosphoryl transfer.⁴⁴ We also found two orthogonal conformationally excited states ES1 and ES2. While the ES1 state corresponds to the inactive kinase conformations, the ES2 state was never observed in the crystallized structures. Our previous CEST NMR measurements suggested the presence of a sparsely populated state that, at that time, we were unable to characterize structurally. These new simulations and MSM show that the transition from GS to ES2 is due to a disruption of hydrophobic packing, featuring a conformational rearrangement for the αC-β4 loop, which causes a partial disruption of the hydrophobic R-spine. These structural changes interrupt the allosteric coupling between the two lobes, as shown by mutual information analysis. A single mutation (F100A), suggested by our simulations, promotes the flip of the αC-β4 loop and reproduces the hypothesized structural uncoupling between the two lobes of PKA-C. We experimentally tested the effects of the F100A mutations and found that it prevents the enzyme from adopting a conformation competent for substrate binding, resulting in a drastic reduction of the cooperativity between ATP and nucleotide. These NMR data further support our working model, showing that the inter-lobe communication is interrupted and the binding response of the F100A kinase is attenuated based on chemical shift perturbation. Furthermore, our study further emphasizes the active role of the hydrophobic interior of the kinase and shows that small alterations in the spines and shell sequences may lead to a dysfunctional kinase by either preventing phosphorylation (see V104G and I150A mutations)²⁷ or by disrupting the binding cooperativity as for F100A.

The αC-β4 loop is a regulatory element present in all EPKs and its importance has been stressed in bioinformatics studies⁴⁵ and supported by computational work, showing that F100 and F102 are at intersections of various communities and constitute a critical hydrophobic node, anchoring the N- to the C-lobe.³⁷ More importantly, studies on EGFR and ErbB2 kinases lead to the hypothesis that the αC-β4 loop may act as a molecular brake ^46,47 or an autoinhibitory switch.⁴⁸ Therefore, it is not surprising that activating mutations and in-frame insertions in the αC-β4 loop are frequently found in kinase-related cancers, such as somatic oncogenic mutations such as P101S, P101L, and L103F.⁴⁹ Finally, an elegant study by Kannan and coworkers emphasized the role of the αC-β4 loop in dimerization of EGFR.²¹ These researchers found that in-frame insertions rigidify and activate the kinase in a length dependent manner. More importantly, these human mutations display a gradual response to drugs, a factor that may exploited for designing mutant-selective inhibitors of EGFR.⁵⁰ The data presented here show that it is possible to abolish the binding cooperativity of a kinase by turning the dial in the opposite direction, i.e., increasing the dynamics of the αC-β4 loop and disconnecting the allosteric network between the N- and C-lobes. The identification of this new, partially inactivating pathway provides further understanding of how to control the dynamics and function of kinases.

Material and methods

Replica-averaged metadynamics (RAM) simulations

I System setup

We used the crystal structure of the wild type PKA-C (PDB ID: 1ATP) as the template and added the missing residues 1-14 at the N terminus. The protonation state of histidine residues followed our previous settings.⁵¹ The protein was solvated in a rhombic dodecahedron solvent box with the tree-point charge TIP3P model⁵² for water extended approximately 10 Å away from the surface of the proteins. Counter ions (K⁺ and Cl^-) were added to ensure electrostatic neutrality corresponding to an ionic concentration of ∼ 150 mM. All protein covalent bonds were constrained with the LINCS algorithm.⁵³ and long-range electrostatic interactions are treated with the particle-mesh Ewald method with a real-space cut-off of 10 Å.⁵⁴ Parallel simulations on the apo form, the binary form with one Mg²⁺ ion and one ATP, and the ternary form with two Mg²⁺ ions, one ATP and one PKI5-24 are performed simultaneously using GROMACS 4.6⁵⁵ with CHARMM36a1 force field.⁵⁶ For the two mutants, F100A and V104G, the corresponding residues were mutated through the mutagenesis wizard of PYMOL.

II Standard MD simulations

Each system was minimized using the steepest decent algorithm to remove the bad contacts, and then gradually heated to 300 K at a constant volume over 1 ns, using harmonic restraints with a force constant 1000 kJ/(mol *Ų) on heavy atoms of both proteins and nucleotides. Over the following 12 ns of simulations at constant pressure (1 atm) and temperature (300 K), the restraints were gradually released. The systems were equilibrated for an additional 20 ns without positional restraints. A Parrinello-Rahman barostat⁵⁷ was used to keep the pressure constant, while a V-rescale thermostat ⁵⁸ with a time step of 2 fs was used to keep the temperature constant. Each system was simulated for 1.05 µs, with snapshots recorded every 20 ps.

III Replica exchange (REX) simulations

Following the standard MD simulations, parallel replica exchange (REX) simulations on the apo, binary and ternary form of PKA-C were set up. For each REX simulations, four replicas were used, and the initial structures were randomly chosen from the µs-scale unbiased simulations. Chemical shifts of PKA-C for N, CA, CO, CB, HN from NMR experiments were imposed as restraints based on the penalty function

where Α is the force constant, k runs over all residues of the protein, l denotes the different back-bone atoms, and m runs over the four replicas. is computed using CamShift that is a module of ALMOST-2.1.⁵⁹ The force constant was gradually increased from 0 (unbiased) to 20 (maximum restraints for production) over 50 ns. All other settings are the same as the previous unbiased simulations. REX simulations were carried out with GROMACS 4.6,⁵⁵ with the replica exchange controlled by the module PLUMED 2.1.1.⁶⁰ These simulations were further conducted ∼100 ns for each replica of the three forms.

IV Replica-averaged metadynamics (RAM) simulations

The subsequent RAM simulations were started from the final structures of REX simulations. The CS restraints were imposed in the same way as the REX simulations. Four CVs are chosen to increase the conformational plasticity around the catalytic cores (detailed in Figure S1): (CVI) the ψ angles of backbone of all the loops that are not in contact with ATP (Back-far), (CVII) the ϕ angles of backbone of all the loops that are in contact with ATP (Back-close), (CVIII) the χ1 angles of side-chains of all the loops that are in contact with ATP (Side-close), (CVIV) the radius of gyration calculated over the rigid part (i.e., residues 50-300) of the protein (rgss). Gaussian deposition rate was performed with an initial rate of 0.125 kJ/mol/ps, where the σ values were set to 0.5, 0.2, 0.2, and 0.01 nm for the four CVs, respectively. The RAM simulations were also carried out with GROMACS 4.6 in conjunction with PLUMED 2.1 and ALMOST 2.1, and continued for ∼400 ns for each replica, with exchange trails every 1 ps.

V Reconstruction of free energy surface (FES) from the RAM simulations

After about 300 ns, the sampling along the first 3 CVs reached convergence, allowing a reliable reconstruction of the corresponding FES. The production run was continued for another 100 ns to sample enough conformations. These conformations were first clustered into microstates using the regular spatial method (cut-off radius of 0.13), and then the free energy of each state is reweighted according to the deposited potential along each CV, which can be obtained from the analysis module of METAGUI.⁶¹ To visualize the distribution of these microstates and their relative energies differences, we further performed principal component analysis to project the microstates represented in the 3-dimennsional CV space into 2-dimensional space spanned by PC1 and PC2. Then we can plot these microstates in a 3-dimensional space spanned by PC1, PC2 and free energy differences ΔG.

VI Independent validation of chemical shifts with SPARTA+

During the REX and RAM simulations, chemical shifts were computed via CamShift in ALMOST 2.1. As an independent validation for the efficacy of the bias, we further calibrated chemical shifts of 2000 snapshots with MDTraj⁶² and SPARTA+.²³

VII Adaptive sampling

The first round of adaptive sampling started from the snapshots of low-energy microstates obtained from the previous step, i.e., 1200 structures for the apo form, 400 structures for the binary form and 200 structures for the ternary form. The initial velocities were randomly generated to satisfy the Maxwell distribution at 300K. For the apo form, a 10 ns simulation was performed for each run, whereas for the binary, each simulation lasted 30 ns, resulting in a total of 12 μs trajectories for both the apo and binary forms. To obtain converged free energy landscape, a total of three rounds of adaptive sampling was started from the 400 microstates that was obtained by K-mean clustering of all snapshots of previous ensembles. Therefore, a total of 100 µs trajectories and 100,0000 snapshots (100 ps per frame) were collected for both the apo and binary form after three rounds of adaptive sampling, and a total of 60 µs trajectories were collected for the ternary form.

VIII Markov state model (MSM) and time-lagged independent component analysis (tICA)

The Cartesian coordinates of key hydrophobic residues, include R-spine residues, L95, L106, Y164 and F185, and the shell residues, V104, M118 and M120, were chosen as the metrics to characterize the conformational transition of the hydrophobic core of PKA-C. Specifically, each snapshot was first aligned to the same reference structure by superimposition of αE (residues 140-160) and αF helices (residues 217-233), and represented by the deviation of Cartesian coordinates of the key residues. The representation in this metric space was further reduced to 10-dimension vectors using time-lagged independent component analysis (tICA)⁶³ at a lag time of 1 ns. All the snapshots were clustered into 400 microstates with K-mean clustering. A MSM was built upon the transition counts between these microstates.

IX Kinetic Monte Carlo trajectory of PKA-C in different forms

Long trajectories were generated using a kinetic Monte Carlo method based on the MSM transition probability matrix of the three forms of PKA-C. Specifically, the discrete jumps between the 100 microstates were sampled for 60 us. And then random conformations were chosen for each state from all the available snapshots. Time series of various order parameters were analyzed subsequently.

Protein expression and purification

The recombinant human Cα subunit of cAMP-dependent protein kinase with the Phe to Ala mutation in position 100 (PKA-C^F100A) was generated from the human PKA-Cα wild-type using Quik-Change Lightning mutagenesis kit (Agilent genomics). The key resource table lists the PCR primers used to modify the pET-28a expression vector encoding the wild-type human PKA-Cα gene (PRKACA – uniprot P17612) ^14-16. The unlabeled and uniformly ¹⁵N-labeled PKA-C^F100A mutant was expressed and purified following the same protocols used for the wild-type protein.⁶⁴ Briefly, transformed E. coli BL21 (DE3) pLyss cells (Agilent) were cultured overnight at 30 °C in LuriaBertani (LB) medium. The next morning, the cells were transferred to fresh LB medium for the overexpression of the unlabeled protein or to M9 minimal medium supplied with ¹⁵NH4Cl (Cambridge Isotope Laboratories Inc.) as the only nitrogen source for the labeled protein overexpression. In both cases, protein overexpression was induced with 0.4 mM of β-D-thiogalactopyra-noside (IPTG) and carried out for 16 hours at 20 °C. The cells were harvested by centrifugation and resuspended in 50 mM Tris-HCl, 30 mM KH2PO4, 100 mM NaCl, 5 mM 2-mercaptoethanol, 0.15 mg/mL lysozyme, 200 μM ATP, DNaseI, 1 tablet of cOmplete ULTRA EDTA-free protease inhibitors (Roche Applied Science) (pH 8.0) and lysed using French press at 1,000 psi. The cell lysate was cleared by centrifugation (60,000 × g, 4 °C, 45 min), and the supernatant was batch-bound with Ni²⁺-NTA agarose affinity resin (Thermo Scientific). The his-tagged PKA-C^F100A was eluted with 50 mM Tris-HCl, 30 mM KH2PO4, 100 mM NaCl, 5 mM 2-mercaptoethanol, 0.5 mM phenylmethylsulfonyl fluoride (PMSF) (pH 8.0) supplied with 200 mM of imidazole. The tail of poly-His was cleaved using a stoichiometric amount of recombinant tobacco etch virus (TEV) protease in 20 mM KH2PO4, 25 mM KCl, 5 mM 2-mercaptoethanol, 0.1 mM PMSF (pH 6.5), over-night at 4 °C. The different phosphorylation states of PKA-C^F100A were separated using a cation exchange column (HiTrap Q-SP, GE Healthcare Life Sciences) using a linear gradient of KCl in 20 mM KH2PO4 at pH6.5.⁶⁵ The purified protein isoforms were then stored in phosphate buffer containing 10 mM dithiothreitol (DTT), 10 mM MgCl2, and 1 mM NaN3 at 4 °C. The protein purity was assessed by sodium dodecyl sulfate-polyacrylamide gel electrophoresis (SDS–PAGE).

Peptide synthesis

The Kemptide (LRRASLG) and PKI5–24 (TTYADFIASGRTGRRNAIHD) peptides were synthesized using a CEM Liberty Blue microwave synthesizer using standard Fmoc chemistry. All peptides were cleaved with Reagent K (82.5% TFA, 5% phenol, 5% thioanisole, 2.5% ethanedithiol, and 5% water) for 3 h and purified using a semipreparative Supelco C18 reverse-phase HPLC column at 3 mL/min.

The purified peptides were concentrated, lyophilized, and stored at −20 °C. Molecular weight and quantity were verified by MALDI-TOF and/or amino-acid analysis (Texas A&M University).

Isotherm titration calorimetry (ITC) measurements

PKA-C^F100A was dialyzed into 20 mM MOPS, 90 mM KCl, 10 mM DTT, 10 mM MgCl2, and 1 mM NaN3 (pH 6.5) and concentrated using conical spin concentrator (10 KDa membrane cut-off, Millipore) to a solution at 80-100 μM, as confirmed by A280 = 55,475 M⁻¹ cm⁻¹. Approximately 300 μL of protein was used for each experiment, with 50 μL of 2 mM ATPγN and/or 1 mM PKI5–24 in the titrant syringe. All measurements were performed at 300 K in triplicates with a low-volume NanoITC (TA Instruments). The binding was assumed to be 1:1, and curves were analyzed with the NanoAnalyze software (TA Instruments) using the Wiseman isotherm ³²

where d[MX] is the change in total complex relative to the change in total protein concentration, d[Xtot] is dependent on r (the ratio of Kd relative to the total protein concentration), and Rm (the ratio between total ligand and total protein concentration). The heat of dilution of the ligand into the buffer was considered for all experiments and subtracted.

The free energy of binding was determined from:

where R is the universal gas constant and T is the temperature at measurement (300 K). The entropic contribution to binding was calculated using:

The degree of cooperativity (σ) was calculated as:

where Kd apo is the dissociation constant of PKI5–24 binding to the apo-enzyme, and Kd Nucleotide is the corresponding dissociation constant for PKI5–24 binding to the nucleotide-bound the kinase.

Enzyme assays

Steady-state coupled enzyme activity assays using Kemptide as substrate were performed under saturating ATP concentrations and spectrophotometrically at 298 K, as reported by Cook et al.³⁰ The values of Vmax and KM were obtained from a nonlinear fit of the initial velocities to the Mich-aelis-Menten equation.

NMR spectroscopy

NMR measurements were performed on a Bruker Avance NEO spectrometer operating at a ¹H Larmor frequency of 600 MHz equipped with a cryogenic probe or on a Bruker Avance III 850 MHz spectrometer equipped with a TCI cryoprobe. The NMR experiments were recorded at 300K in 20 mM KH2PO4 (pH 6.5), 90 mM KCl, 10 mM MgCl2, 10 mM DTT, 1 mM NaN3, 5% D2O, and 0.1% 4-benzene sulfonyl fluoride hydrochloride (AEBSF, Pefabl–c - Sigma-Aldrich). Concentrations for samples were 0.15 mM of uniformly ¹⁵N-labeled PKA-C^F100A, as determined by A280 measurements, 12 mM ATPγN or ADP was added for the nucleotide-bound form, and 0.3 mM PKI5-24 for the ternary complex. [¹H, ¹⁵N]-WADE-TROSY-HSQC pulse sequence ⁶⁶ was used to record the amide fingerprint spectra of PKA-C^F100A in the apo, nucleotide-bound (ADP- or ATPγN-bound – binary form), and ternary complex (PKAC^F100A/ATPγN/PKI5-24). All [¹H, ¹⁵N]-WADE-TROSY-HSQC experiments were acquired with 2048 (proton) and128 (nitrogen) complex points, processed using NMRPipe⁶⁷ and visualized using NMRFAM-SPARKY ⁶⁸ and POKY.⁶⁹ Combined chemical shift perturbations (CSP) were calculated using ¹H and ¹⁵N chemical shifts according to:

in which Δδ is the CSP; ΔδH and ΔδN are the differences of ¹H and ¹⁵N chemical shifts, respectively, between the first and last point of the titration; and 0.154 is the scaling factor for nitrogen.⁷⁰

COordiNated ChemIcal Shift bEhavior (CONCISE)

The normal distributions reported in the CONCISE plot were calculated using principal component analysis for residues whose chemical shifts responded linearly to ligand binding.³⁴ In this work, we use the ¹H and ¹⁵N chemical shifts derived from the [¹H, ¹⁵N]-WADE-TROSY-HSQC experiments for the apo, ADP, ATPγN, and ATPγN/PKI5–24 bound forms of PKA-C.

CHEmical Shift Covariance Analysis (CHESCA)

This analysis was used to identify and functionally characterize allosteric networks of residues eliciting concerted responses to, in this case, nucleotide and pseudosubstrate. To identify inter-residue correlations, four states were used: apo, ATPγN-bound, ADP-bound, and ATPγN/PKI5–24. The identification of inter-residue correlations by CHESCA relies on agglomerative clustering (AC) and singular value decomposition (SVD). Pairwise correlations between chemical shift variations experienced by different residues were calculated to identify networks. When plotted on a correlation matrix, this allows for the identification of regions that are correlated to one another. For each residue, the maximum change in chemical shift was calculated in both the ¹H (x) and ¹⁵N (y) dimensions (Δδx,y). The residues included in the CHESCA analysis were the ones that satisfied the following: Δδx,y > ½ ΔνxA,yA+½ ΔνxB,yB, where A and B correspond to two different forms analyzed (note that there is no dependence on which two forms satisfied this statement), and Δν denotes the linewidth. Correlation scores were used to quantify the CHESCA correlation of a single residue or a group of residues with another group. Correlation scores were evaluated for both (a) a single residue and (b) the full protein. The generalized expression for evaluating either case is:

where RiJ is the CHESCA correlation matrix and i and j denote (a) a single residue and all other assigned residues in the protein, respectively, or (b) both represent all the assigned residues in the protein. For all the analysis a Rcutoff of 0.98 was used.

Community CHESCA analysis is a chemical shift-based correlation map between functional communities within the kinase. Each community is a group of residues ³⁷ associated with a function or regulatory mechanism. To represent community-based CHESCA analysis, we lowered the correlation cut-off such that Rcutoff > 0.8.

Suppose communities X and Y have nA and nB number of assigned residues respectively, the correlation score between A and B is defined as,

where RiJ is the CHESCA correlation coefficient between residue i (belongs to community A) and residue j (belongs to community B). Rcutoff is the correlation value cutoff. RA,B can take values from 0 (no correlation between residues in A and B) to 1 (all residues in A has correlation > cutoff with all residues in B).

Acknowledgements

This work was supported by the National Institute of Health GM 100310 and HL 144130 to GV. The authors would like to acknowledge the Minnesota Supercomputing Institute for MD calculations. YW would like to thank Guangdong Pearl River Talent Program (2021QN02Y618) and National Natural Science Foundation of China (22007069, 92269102) for part of the MD analysis carried out at the Shenzhen Bay Laboratory Supercomputing Centre.

Competing interest

The authors declare no financial and non-financial competing interests.

Author contributions

Y.W., C.O., and G.V. designed research. Y.W. performed the MD simulations and analysis. C.O. prepared the kinase samples, executed all NMR and ITC experiments, and analyzed the results. C.W. performed the kinase activity assay and prepared the peptides for the assay. C.O., Y.W. contributed to the manuscript, preparing all the text and figures. M.V.S., K.N.H., C.C., J.G., and M.V. contributed to the critical analysis of the data and writing of the manuscript. D.A.B. contributed to the critical analysis of the data and writing the manuscript. S.S.T. contributed to the writing of the manuscript. G.V. designed and directed the experiments and contributed to the writing of the manuscript. All authors have given approval to the final version of the manuscript.

To whom correspondence should be addressed. Email: vegli001@umn.edu.