Research Article Cite This: ACS Catal. 2018, 8, 2375−2384
pubs.acs.org/acscatalysis
Unraveling Allosteric Mechanisms of Enzymatic Catalysis with an Evolutionary Analysis of Residue−Residue Contact Dynamical Changes Phuoc Jake Vu,† Xin-Qiu Yao,† Mohamed Momin, and Donald Hamelberg* Department of Chemistry, Georgia State University, Atlanta, Georgia 30303-2515, United States S Supporting Information *
ABSTRACT: The evolution of protein conformational dynamics contains important information about protein function and regulation. Here, we describe an approach to dynamical-evolution analysis based on multiple microsecond molecular dynamics simulations and residue− residue contact analysis. We illustrate our approach by comparing three human cyclophilin isoforms, cyclophilin A, D, and E, which belong to a family of enzymes catalyzing peptidyl-prolyl cis−trans isomerization. Our results reveal that despite distinct overall equilibrium conformations between cyclophilins under substrate-f ree conditions, functional dynamical changes resembling substrate-binding and catalytic processes tend to be conserved. Key residues displaying either concerted or specific dynamical changes among isoforms during the reactions are identified, which delineate two distinct allosteric pathways for cyclophilin function consistent with recent nuclear magnetic resonance experiments. A sequence-based coevolution analysis is also employed for further understanding dynamical consequences. Our results collectively provide a framework where both common and specific functional mechanisms of a protein family can be elucidated. KEYWORDS: evolution, allosteric regulation, enzyme dynamics, residue−residue contact, cyclophilin, molecular dynamics
■
INTRODUCTION Protein internal motions or dynamics have been increasingly recognized to play a crucial role in protein function and regulation.1−3 Typical protein dynamics span a broad range of spatiotemporal scales, from subtle backbone and side chain fluctuations (∼10−9 s) and domain motions (∼10−6 s) to largeamplitude conformational changes normally observed in biological molecular machines (10−3−102 s). Notable examples in which dynamics determine function include the dynamical rearrangements in enzymes to facilitate catalytic turnover,4 the conformational changes in transporters to pump small molecules in and out of cell membrane,5 the force-producing structural changes in molecular motors,6 and the prevailing roles of dynamics in the allosteric regulation during signal transduction.7−9 Despite numerous efforts in inspecting protein dynamics with both experimental and computational methods, such as nuclear magnetic resonance (NMR)10 and molecular dynamics (MD) simulation,11 how protein molecules harness thermal fluctuations for function remains elusive. Protein dynamics reflect the underlying energy landscape, which is predominantly determined by protein sequence. In this perspective, evolution modified protein dynamics and function by altering the energy landscape. During evolution, certain patterns of protein dynamics must be conserved to retain the core protein function that arose early. Similar to how multiple sequence alignment screens out variable residues to © XXXX American Chemical Society
detect structurally or functionally critical sites, comparing protein dynamics at aligned residues across a protein family identifies dynamically conserved sites that are intimately related to the core function. Meanwhile, variable dynamics in synergy with variable sequences across family members underlie the subtle functional diversity, and hence, comparative analysis of dynamics also helps to foster our understanding of protein functional specificity (e.g., selective substrate binding or distinct kinetics). Indeed, this new paradigm of an evolution-based approach to delineate the complex sequence−structure− dynamics−function relationship has gained increasing interests. By comparing the crystallographic structures of several homologous proteins under distinct functional states, Babu and colleagues recently derived the universal mechanisms governing the activation of heterotrimeric guanine nucleotidebinding proteins (G proteins) and G protein-coupled receptors.12,13 With X-ray crystallography and NMR spectroscopy, Wright and colleagues compared the dynamics of human and Escherichia coli dihydrofolate reductase (DHFR), an enzyme that catalyzes the NADPH-dependent reduction of dihydrofolate to tetrahydrofolate, and identified the residue determinants underlying the distinct catalytic efficiency and Received: December 12, 2017 Revised: January 29, 2018 Published: February 1, 2018 2375
DOI: 10.1021/acscatal.7b04263 ACS Catal. 2018, 8, 2375−2384
Research Article
ACS Catalysis robustness between the human and the bacterial enzymes.14 In addition to the comparisons between extant proteins, the structures and dynamics of ancestral proteins have been modeled with a computational ancestor reconstruction method, which empowers a direct tracking of the evolution of protein dynamics.15−17 These recent advances demonstrate how a study of protein dynamics in the context of evolution provides unprecedented insights into functional mechanisms. The knowledge obtained complements that derived from a largescale evolutionary analysis aimed at inferring sequence and structural traits underlying conserved and divergent functions of enzyme superfamilies18,19 and can be further leveraged to develop new algorithms for rational protein design and drug discovery. However, a general approach to the evolutionary analysis of conformational dynamics for most biological systems is still lacking. Cyclophilin A (CypA) is a ubiquitous protein where function is strongly coupled with dynamics.4,20 As a peptidyl-prolyl cistrans isomerase (PPIase), cyclophilin A catalyzes the interconversion between the cis (with the peptidyl-prolyl torsion angle ω = 0°) and the trans (ω= ± 180°) states of the prolyl peptide bond in Xaa-Pro motifs, where Xaa represents any amino acid (Figure 1A,B). Humans have 17 cyclophilin
1A). The three isoforms of cyclophilin have the same PPIase activity; however, they function in different subcellular locations: CypA is generally found in the cytosol, CypD in the mitochondria matrix, and CypE in the nucleus.22 We recently examined the dynamical properties, represented by the breaking and formation of residue−residue contacts, of CypA under various substrate-binding and mutational conditions using a combined approach of MD simulations and NMR experiments.23 In particular, we found an interesting “dynamic cluster” containing residues showing substantial dynamical changes upon substrate binding located ∼15 Å away from the active site. This finding reveals a new allosteric mechanism in CypA, but its generality across the cyclophilin family is not known. In this work, we develop a new approach for general evolutionary analysis of protein conformational dynamics based on multiple microsecond-long MD simulations and the similar contact analysis method previously described.23 We illustrate our approach by comparing the dynamics between human CypA, CypD, and CypE derived from MD simulations (totaling 26 μs). For each cyclophilin isoform, three functional states collectively representing substrate-binding and catalytic isomerization processes are examined. This allows us to evaluate the conservation of dynamics across isoforms with respect to distinct enzymatic processes. Dissecting the contact dynamics further enables us to identify the key residues determining the common and isoform-specific dynamical changes. We also perform a sequence coevolution analysis to identify correlated amino acid substitutions between residues in the cyclophilin family and find that dynamically conserved contacts do not always represent the most highly coevolving residue pairs, suggesting that dynamics is the missing link and sequence and structure alone cannot fully describe function. The consensus groups of residues derived from the dynamical conservation analysis are consistent with recent NMR experiments.24
■
RESULTS AND DISCUSSION Protein Backbone and Residue Side Chain Dynamics under the Substrate-Free State Are Distinct among CypA, CypD, and CypE. Fluctuation analysis and principal component analysis (PCA) of MD simulation trajectories reveal that backbone dynamics are different between CypA, CypD, and CypE. Multiple long-time (2.2−2.7 μs) MD simulations were performed under substrate-f ree conditions for each cyclophilin isoform. Snapshots of the latter 2-μs of each simulation trajectory were analyzed. Residue-wise averaged root-mean-square fluctuation (RMSF) of backbone atoms derived from the simulations shows that, while the atomic fluctuations at most residues look similar across isoforms, the fluctuations located at the β4-β5 loop (the loop between the fourth and the fifth β-strands), the β5-β6 loop, the β6-β7 loop, and the α2-β8 loop (the loop between the second α-helix and the eighth β-strand) are apparently different (Figure 2A). In general, CypA and CypE are more flexible than CypD. PCA was then performed on the Cartesian coordinates of backbone atoms from the simulations to examine backbone conformational distributions at equilibrium (Figure 2B). It reveals that CypE has the broadest distribution and samples multiple conformational states in the subspace spanned by the top two principal components (i.e., PC 1 and PC 2, which collectively capture nearly 50% of total atomic mean displacements or variance), indicating the overall highest flexibility of CypE among the isoforms. CypA has a distribution with
Figure 1. Backbone and active site structures are highly similar across CypA, CypD, and CypE. (A) The substrate-free crystallographic structures of CypA (white; PDB: 3K0M), CypD (red; PDB: 4O8H), and CypE (green; PDB: 3UCH) are superimposed and are represented as a cartoon. The modeled substrate (see Methods) is displayed as cartoon (yellow) with “Gly-Pro” motif shown as licorice and colored by atom types. The enlarged view shows the substratebinding pocket represented as white transparent surface. Side chains of active site residues are displayed as sticks and are color-coded the same as backbone. (B) Schematic cis−trans isomerization of the “Gly-Pro” motif. Ts, transition state.
isoforms (Figure S1), among which CypA is the best characterized, both experimentally and computationally. Cyclophilins are known targets for immunosuppressant, and their function is critical to many important cellular processes, including protein folding and signal transduction.21 Cyclophilin D (CypD) and cyclophilin E (CypE) are the family members closest to CypA, with high sequence identity (68−75%) and almost identical structures to CypA (the backbone root-meansquare deviation, or RMSD, is within 0.50−0.65 Å; see Figure 2376
DOI: 10.1021/acscatal.7b04263 ACS Catal. 2018, 8, 2375−2384
Research Article
ACS Catalysis
error estimate of f23) distribute all over the cyclophilin molecule for all the comparisons (Figure 2C−E), indicating that side chain dynamics are distinct among the isoforms. In summary, although in both sequence and structure cyclophilin isoforms are very similar, they possess different conformational ensembles at equilibrium at least under the substrate-f ree state. Dynamical Changes of Residue−Residue Contacts during Substrate Binding Are Conserved between CypA and CypE but Are Not Conserved between CypA/E and CypD. Although the conformational dynamics under an individual state (i.e., the substrate-free state) are very different between the cyclophilin isoforms, dynamical changes from the substrate-free to the substrate-bound state display a certain extent of similarity. Multiple 2.4-μs additional simulations were performed under substrate-bound conditions where the peptidyl-prolyl torsion angle (ω) was in the cisconformation (termed “cis-bound state”). The probability difference of contact formation during simulations between the substrate-f ree and the cis-bound states was used to characterize the dynamical changes during substrate binding. Analysis of the CypA simulations reveals a site ∼15 Å away from the active site showing substantial dynamical changes (|df = f Y − f X| ≥ 0.1, where f X and f Y are contact probabilities calculated for the two compared states) of residue contacts (Figure 3A). This observation resembles the dynamic cluster identified in the previous simulation study of CypA,23 although the substrate employed in the present work is five residues longer than that used in the previous study.23 This consistency indicates that the revealed dynamical changes are an intrinsic dynamical characteristic of CypA independent from the identity of the bound substrate. Intriguingly, overall similar patterns of dynamical changes are observed between CypA and CypE (Figure 3A,C). A large portion of contacts that are either more often formed (with an increase of contact probability) or more often broken (decrease of contact probability) upon substrate binding in CypA are shown to have the same direction of changes in CypE (Figure 3E). These contacts with the same trends of changes from one state to the other between isoforms are defined as “dynamically conserved contacts”. To quantitatively measure the overall similarity of dynamical changes between isoforms, we developed a dynamical conservation index (DCI), which is defined by the percentage of dynamically conserved contacts. In the calculation of DCI, contacts showing small dynamical changes (absolute contact probability difference |df|
