Open Access Article
Julian E.
Fuchs
a,
Inés G.
Muñoz
b,
David J.
Timson
c and
Angel L.
Pey
*d
aInstitute of General, Inorganic and Theoretical Chemistry, Faculty of Chemistry and Pharmacy, University of Innsbruck, Innsbruck, Austria
bCrystallography and Protein Engineering Unit, Structural Biology and Biocomputing Programme, Spanish National Cancer Research Centre (CNIO), Madrid, Spain
cSchool of Pharmacy and Biomolecular Sciences, The University of Brighton, Brighton, UK
dDepartment of Physical Chemistry, Faculty of Sciences, University of Granada, Granada, Spain. E-mail: angelpey@ugr.es
First published on 13th June 2016
Theoretical and experimental evidence has shown that protein function, regulation and degradation are intrinsically linked to the dynamic and fluctuating nature of protein ensembles. However, the effect of missense mutations on catalytic performance are often interpreted from conformational analyses derived from X-ray crystallography, molecular dynamics and modeling, while effects on conformational fluctuations at the active site as the source of catalytic defects are rarely investigated. Here, we explore the role of conformational fluctuations in the catalytic efficiency of WT and three missense mutations in the UDP-galactose 4′-epimerase (GALE) protein causing type III galactosemia. Using comprehensive molecular dynamics simulations and small angle X-ray scattering we correlate low NAD+ binding affinity in some mutants with an increased population of non-competent conformations for NAD+ binding. Proteolysis studies combined with thermodynamic calculations reveal that mutations affecting GALE catalytic performance favor larger conformational fluctuations at the N-terminal domain and NAD+ binding site, shifting the equilibrium towards non-binding competent states in the native ensemble. Therefore, we provide a novel ensemble-based thermodynamic mechanism to explain catalytic defects caused by missense mutations that links large and transient conformational fluctuations and loss of catalytic efficiency and substrate/coenzyme affinity. In the context of this mechanism, we propose that allosteric ligands aimed at modulating these transient conformational fluctuations might correct catalytic defects in inherited metabolic diseases, representing a different approach to current small ligand therapies which target the low stability, but not catalytic defects, of mutations.
Despite the emergent importance of native state dynamics on protein functionality, its role in mutation-associated catalytic defects in inherited metabolic diseases has been largely unexplored. Here, we investigate the relationship between dynamic alterations and catalytic defects in three catalytic mutations of UDP-galactose 4′-epimerase (GALE, EC 5.1.3.2) causing type III galactosemia. GALE catalyzes the interconversion between UDP-galactose and UDP-glucose using NAD+ as cofactor.25 Two common GALE disease-causing mutations (p.N34S and p.G90E) decrease NAD+ binding affinity, selectively decrease protein unfolding cooperativity and facilitate proteolytic cleavage at the loop formed by residues 34–45 near the NAD+ binding site.12 The impaired catalytic performance in these two mutants may arise from structural fluctuations of the NAD+ binding site that energetically penalize coenzyme binding. The p.V94M mutant greatly reduces GALE catalytic efficiency,25,26 possibly due to binding of the substrate in a non-competent fashion for catalysis,27 without largely affecting protein stability and flexibility.12 We have combined molecular dynamic (MD) simulations and small angle X-ray scattering (SAXS) to ascertain the effects of these mutations on the native conformation and dynamics. Furthermore, we have experimentally characterized the overall conformation and relative population of partially unfolded states under native conditions using partial proteolysis in combination with thermodynamic analyses. Our results allow to propose a thermodynamic mechanism that explains how mutational defects on NAD+ binding occur due to an increment of the population of non-competent binding states and changes in their structure, and as result, allosteric ligands that bind to and stabilize the N-terminal domain of GALE might be used to overcome functional defects in type III galactosemia. We propose that these results may generally apply to inherited metabolic diseases in which point mutations, distant from the enzyme active site, influence binding or catalysis by similar mechanisms.
28) that was prepared for simulations using MOE's protonate3D and mutate function.29
MD simulations were performed using the GPU code of pmemd in Amber12.30 Amino acids were parameterized using the Amber force field 99SB with ILDN corrections,31 whilst the NAD+ co-factors were described according to parameters from Walker and Pavelites.32,33 We employed an extensive equilibration protocol34 before performing unrestrained sampling simulations for 100 ns in NpT ensemble at 300 K and 1 bar. We stored 1000 equal-spaced snapshots to trajectory for later analysis. Shorter sampling periods (i.e. 10 ps) yielded essentially the same results (data not shown).
We checked the simulations for energetic and structural stability and found maximum Cα RMSD values to be below 3.5 Å over all 692 protein residues. Resulting conformational ensembles were analyzed using ptraj and cpptraj from AmberTools.35 We extracted residue-wise B-factors derived from root mean squared fluctuations of Cα atoms following a global alignment. Additionally, we calculated thermodynamic entropies from the distributions of protein backbone torsion angles and summarized the three individual degrees of freedom to residue-wise and total entropies.36 To study the behavior of the 34–45 loop, we performed an alignment to the whole respective protein subunits and recorded occurring loop movements. Subsequently, we clustered the loop structures according to their Cα RMSD into two states using a hierarchical agglomerative clustering approach. Respective cluster centroids were extracted as representative structures and compared in terms of structure and cluster occupancy. Values for dimeric systems are reported as average over both individual sub-units. Structures were visualized using Pymol.37
K. During the experiment the samples were exposed for 300 s in 10 s acquisition blocks using a sample to detector distance of 3.9 m and X-ray wavelength of 1 Å. The data were analyzed, buffer-subtracted, scaled, and merged using the Scatter software package (http://www.bioisis.net). This software was also used to check possible radiation damage of the samples by visual inspection of the Guinier region as a function of exposure time during data collection. RG and Dmax values were calculated with PRIMUS and GNOM, and shape estimation was carried out with DAMMIF/DAMMIN. Real-space scattering profiles of atomic models were calculated using CRYSOL and the final ab initio models were superimposed with the high-resolution structure using the program SUPCOMB. The program OLIGOMER was used to calculate the ratio of the species present in the solutions. All the SAXS programs are included in the ATSAS programs package.38
The apparent rate constants for proteolysis (kp) were determined using a single exponental function. For the estimation of equilibrium m values, kp was corrected for the inhibitory effect of urea on thermolysin degradation using the values previously reported.39 The corrected rate constants (k′p) were used to determine the equilibrium m values between the native and cleavable state (X) by determining their urea concentration dependence using the following equation:
![]() | (1) |
To estimate the effect of mutations on the thermodynamic stability of the native vs. the cleavable state, we have assumed that mutations do not affect the intrinsic proteolytic step (X → F), and thus, changes in the overall proteolysis rate after normalization using protease concentration (kWT and kmut) reflect changes in the equilibrium constant between N and X (N ↔ X;41), and therefore in their free energy changes (ΔΔG), following eqn (2):
ΔΔGmut–WT = −RT ln(kmut/kWT) | (2) |
MD simulations were used to characterize the conformational microstates present in native GALE ensembles at atomic resolution. GALE monomeric enzymes are rather unstable over 100 ns simulation time, and therefore, subsequent analyses were carried out using the functional GALE dimer only. We have focused our analyses on the 34–45 loop where major conformational transitions are found (reaching RMSD values up to 8 Å in apo GALE systems; see Fig. S1† for a 2D RMSD plot). By contrast, NAD+-bound GALE complexes show drastically reduced loop mobility (Fig. S1†), supporting its stabilizing effect towards proteolysis.12 The aggregate simulation time of 1.6 μs of the dimer (i.e. 3.2 μs sampling time for individual loops) make us confident on the observed trends for loop dynamics, despite individual simulation (100 ns) is at the lower limit for investigating changes in loop movements.45 In this regard, enhanced sampling methods46 could be helpful to provide further insight into the loop trends described here.
We performed a cluster analysis to investigate plausible loop conformations during molecular dynamics simulations (Fig. 1A–D). MD trajectories were clustered into two states: the major one representing a ground state (close to that of the crystal structure) and the minor excited (i.e. high-energy) state. The relative populations of the major and minor states depend on the cut-off values used in the clustering procedure, even though the qualitative trends described here for the effect of mutations and NAD+ binding are robust.
In the absence of NAD+, the minor conformation is more populated in the p.N34S and p.G90E mutants (Fig. 1E), defective in NAD+ binding, while generally the thermolysin cleavage site Ala38–Phe39 is more solvent exposed in these high-energy states (Fig. 1A and C). In the presence of NAD+, a much lower population of the high-energy states is found, except for p.V94M (Fig. 1E), an in all cases the conformation of the ground and excited states is more similar (Fig. 1B and D). These analyses support that two conformational states exist in equilibrium in the absence of NAD+, with the more populated displaying higher affinity for NAD+ due to higher structural complementarity that likely prevents from thermolysin cleavage (Fig. 1A–D). Therefore, the lower affinity for NAD+ in p.N34S and the inability of p.G90E to bind it, can be partly explained by a shift in the conformational equilibrium towards non-binding competent states, remodeling the protein energy landscape in such a way that cofactor binding is thermodynamically penalized.4,47
According to this proposed mechanism, p.N34S and p.G90E could be more susceptible to proteolytic attack at the 34–45 loop since they increase the population of cleavable (excited) non-binding competent states, while NAD+ binding reduces the population of cleavable states explaining its protective effect.12 Alternatively, accelerated proteolysis and low binding affinity could also be explained by enhancement of local dynamics at the NAD+ binding site that impedes its correct spatial complementarity towards the cofactor. To distinguish between these two mechanisms, we have analyzed local dynamics at the residue level by MD simulations (Fig. 1F). Interestingly, we find that the flexibility around the cleavage site (loop 34–45) does not correlate with the experimentally determined proteolysis susceptibility (Fig. 1F and ref. 12), implying that proteolysis susceptibility is not directly linked to local dynamics in this case. Additionally, p.N34S and p.G90E facilitate thermally induced partial unfolding12 which supports a population-shift mechanism. Addition of NAD+ decreases the local dynamics around the cleavage site, which likely reflects the shift towards binding competent states (Fig. 1E).
![]() | ||
| Fig. 2 Structural models of GALE WT and p.G90E in solution by SAXS. (A and B) Experimental scattering curves (dots) and calculated (lines) SAXS data for GALE WT (A) and p.G90E (B) at different protein concentrations. (C and D) P(r) functions for GALE WT (C) and p.G90E (D). The tail reveals the more elongated dimer shape of the p.G90E mutant. (E and F) Ab initio shape reconstructions of GALE WT and p.G90E were superimposed on the crystal structure (PDB entry 1EK6). | ||
In combination with our results using MD simulations (Fig. 1), these SAXS analyses support the existence of conformational fluctuations (i.e. high-energy states) in GALE in solution, which are enhanced by the p.G90E mutation.
![]() | ||
| Fig. 3 Urea dependence of proteolysis rate constants and thermodynamic analysis. (A–D) proteolysis rate constants (kp) at different urea concentrations in the absence or presence of NAD+; (E–H) corrected proteolysis rate constants (k′p) by urea inhibition of thermolysin as reported by ref. 39. The slopes of these plots are used to determined equilibrium m values. (I and J) Equilibrium m values (I) for the native-to-cleavable state (X) equilibrium determined by proteolysis and the corresponding fraction of accessible surface area (ASA) unfolded in the X state (J). (K) Effect of mutations on the change in free energy (ΔΔG(mut–WT)) between native and cleavable states. Error bars are fitting errors or those propagated from them. | ||
The sensitivity to thermolysin degradation correlates well with the effect of urea on proteolysis kinetics (Fig. 3I–K): the higher the sensitivity to proteolysis of a GALE variant, the higher the effect of urea. For WT and p.V94M, we found fairly similar equilibrium m values, on average a 0.55 kcal mol−1 M−1, which implies that ∼15% of the native structure (∼50 residues) is unfolded in the cleavable state X (Fig. 3J). In the case of p.N34S and p.G90E, the equilibrium m value is around 1.4 kcal mol−1 M−1, and thus, about ∼40% of the native structure (∼130 residues) is unfolded in the cleavable state for these two mutants. Therefore, p.N34S and p.G90E populate a more extensively unfolded cleavable state under native conditions, and since proteolysis occur at Ala38–Phe39,12 likely involves the N-terminal domain that contains the NAD+ binding site. Moreover, a direct comparison of proteolysis rate constants normalized by protease concentration provides a measure of the thermodynamic destabilization induced by p.N34S and p.G90E on the N ↔ X equilibrium (eqn (2)), which amounts to 1.5–2 kcal mol−1, while the effect of p.V94M is marginal (Fig. 3K). These results show that p.N34S and p.G90E shift the N ↔ X towards the cleavable state by a factor of 22.5 ± 4.5 and 14.1 ± 3.1. Overall, in the context of the population-shift mechanism, the lower affinity for NAD+ in p.N34S and p.G90E is explained by the stabilization of non-competent binding states and changes in their structure under native conditions.
Addition of NAD+ has a stabilizing effect towards proteolysis (Fig. 3 and S3†), although this effect is variant-dependent. In p.G90E, only a 20% decrease in the proteolysis rate constant is found, while this effect is of 40% in WT, 60% in p.N34S and 130% in p.V94M. The stabilization agrees well with their corresponding affinities for NAD+.12 Importantly, the stabilization exerted by NAD+ does not cause large changes in the equilibrium m values (Fig. 3I), indicating that NAD+ binding does not alter the structure of the cleavable state but instead simply shifts the conformational equilibrium towards more protease-resistant NAD+ binding competent states.
It must be noted that conformational fluctuations at the NAD+ binding site occur even in GALE variants with WT-like binding affinity (Fig. 1 and 3). Non-binding competent excited states can be used by nature to provide the appropriate binding affinity for enzyme function at physiologically relevant ligand concentrations, allowing conformational excursions to states with different degrees of unfolding.6,47 In this context, p.N34S and p.G90E are deleterious because they shift the conformational equilibrium towards more unstructured partially unfolded conformations in the native state ensemble thus energetically penalizing NAD+ binding (Fig. 4A and B). These excited states seem to be populated to significant levels, as seen by the ability of our MD simulations to detect them even at the low microsecond sampling time scale (Fig. 1), although variants with low binding affinity for NAD+ clearly promote the population of excited states. This population-shift is reflected in our SAXS analyses as an expansion in the average ensemble in solution (Fig. 2), thus overall supporting this conformational-selection mechanism for the binding of NAD+ to GALE. In addition, proteolysis experiment revealed an increase in population of these excited states in p.N34S and p.G90E as well as a change in their overall structure, showing a higher degree of unfolding (Fig. 3). Importantly, since the population of these excited states is inversely proportional to NAD+ binding affinity, it is plausible that allosteric effectors would be able to bind to WT-like excited states destabilized by catalytic mutations such as p.N34S and p.G90E, thus partially restoring the binding affinity for NAD+ (Fig. 4C).
![]() | ||
| Fig. 4 Pictorial representation of the effect of catalytic mutations on the structure of excited states and NAD+ binding affinity and its potential correction by allosteric effectors. Panels (A) and (B) show the equilibrium between native GALE WT (A) and p.G90E (B) and the excited, cleavable states as inferred from native proteolysis experiments. In panel (C), the potential effect of an allosteric ligand binding to an unstructured region of the excited state of a catalytic mutant leading to correction of its NAD+ binding affinity by a coupled binding-folding process. For sake of simplicity, only the GALE monomer is displayed (based on PDB:1EK5). | ||
So far, we have identified a thermodynamic signature for mutations affecting NAD+ binding to GALE. Consistently, this signature is not found in p.V94M, a mutant strongly reducing catalytic efficiency with mild effects on the apparent affinities for the substrate and cofactor.12,25 Previous MD simulations have suggested little effect of p.V94M on the interaction with substrate and cofactor.27 However, p.V94M seems to severely affect the ability of GALE active site to bring together the NAD+ cofactor, the substrate and the catalytic Tyr157, leading to non-optimal ternary complex for catalysis.25 In the case of p.N34S and p.G90E, their low binding affinity for NAD+ (Fig. 4 and ref. 12) is caused by a shift in the conformational equilibrium of the apo-state towards NAD+ non-binding competent states, and in addition, p.G90E could develop unfavorable electrostatic interactions that may also contribute to its apparent inability to bind the cofactor and its dramatically low activity.12,25
Footnote |
| † Electronic supplementary information (ESI) available. See DOI: 10.1039/c6ra05499d |
| This journal is © The Royal Society of Chemistry 2016 |