Evolutionary Analysis As a Powerful Complement ... - ACS Publications

Aug 31, 2018 - International Clinical Research Center, St. Anne's University Hospital Brno, ... Institute of Thermodynamics and Thermal Process Engine...
3 downloads 0 Views 837KB Size
Subscriber access provided by University of South Dakota

Article

Evolutionary Analysis is a Powerful Complement to Energy Calculations for Protein Stabilization Koen Beerens, Stanislav Mazurenko, Antonin Kunka, Sérgio M. Marques, Niels Hansen, Milos Musil, Radka Chaloupkova, Jitka Waterman, Jan Brezovsky, David Bednar, Zbynek Prokop, and Jiri Damborsky ACS Catal., Just Accepted Manuscript • DOI: 10.1021/acscatal.8b01677 • Publication Date (Web): 31 Aug 2018 Downloaded from http://pubs.acs.org on September 4, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Catalysis

Evolutionary Analysis is a Powerful Complement to Energy Calculations for Protein Stabilization Koen Beerens⌘,†, Stanislav Mazurenko⌘, Antonin Kunka⌘, Sergio M. Marques⌘,§, Niels Hansen‡, Milos Musil⌘,$, Radka Chaloupkova⌘,§, Jitka Waterman¶, Jan Brezovsky⌘,§, David Bednar⌘,§, Zbynek Prokop⌘,§*, Jiri Damborsky⌘,§* ⌘Loschmidt

Laboratories, Department of Experimental Biology and Research Centre for Toxic Compounds in the Environment RECETOX, Masaryk University, Kamenice 5/A13, 625 00 Brno, Czech Republic ‡ Institute of Thermodynamics and Thermal Process Engineering, University of Stuttgart, D-70569 Stuttgart, Germany $ Department of Information Systems, Faculty of Information Technology, Brno University of Technology, 612 66 Brno, Czech Republic § International Clinical Research Center, St. Anne's University Hospital Brno, Pekarska 53, 656 91 Brno, Czech Republic ¶ Diamond Light Source, Harwell Science and Innovation Campus, Didcot, OX11 0DE, United Kingdom ABSTRACT: Stability is one of the most important characteristics of proteins employed as biocatalysts, biotherapeutics and biomaterials, and the role of computational approaches in modifying protein stability is rapidly expanding. We have recently identified stabilizing mutations in haloalkane dehalogenase DhaA using phylogenetic analysis but were not able to reproduce the effects of these mutations using force-field calculations. Here we tested four different hypotheses to explain the molecular basis of stabilization using structural, biochemical, biophysical and computational analyses. We demonstrate that stabilization of DhaA by the mutations identified using the phylogenetic analysis is driven by both entropy and enthalpy-contributions, in contrast to primarily enthalpy-driven stabilization by mutations designed by the force-field calculations. Comprehensive bioinformatics analysis revealed that more than half (53%) of 1,099 evolution-based stabilizing mutations would be evaluated as destabilizing by force-field calculations. Thermodynamic integration considers both folded and unfolded states and can describe the entropic component of stabilization, yet it is not suitable for predictive purposes due to computational demands. Altogether, our results strongly suggest that energetic calculations should be complemented by a phylogenetic analysis in protein stabilization endeavors. Keywords: Protein Stabilization, Thermostability, Evolutionary Analysis, Force-Field Calculations, Computational Tools, Entropy, Enthalpy, Thermodynamic integration

were not computationally predicted to be beneficial. Thus, despite the increasing success rates of computational methods5–10, a deeper understanding of the mechanistic basis for protein stability is needed to enable new and improved algorithms for computational tools. We recently reported a new computational method called FireProt6 that was used to direct the development of an improved variant of the haloalkane dehalogenase DhaA (Supporting Information Methods). This approach searches for mutations likely to stabilize a protein of interest via energy calculations of potential point mutants – called ‘energy-based mutations’ – and via phylogenetic analysis to identify residues that have drifted from more stable consensus sequences – called ‘evolution-based mutations’. FireProt then combines both of these calculations using smart filtering to allow the design of highly stabilized multiple-point mutants6. In our previous work, we experimentally characterized a triple-point mutant, DhaA101, containing three evolution-based mutations selected following the back-to-consensus approach. Even

INTRODUCTION Proteins are used in an ever-expanding list of biotechnological applications, ranging from household products and technical industries, to food and animal feed, to fine chemicals and biopharmaceuticals1. Not only do proteins provide high specificity and activity but they are also environmentally friendlier than typical chemical synthesis protocols. However, their application can be hampered by limited stability, since production-line conditions including high temperatures, extreme pH, or the presence of organic solvents or proteases are often far from the natural conditions for which proteins were evolved2. Protein engineering can be used to improve natural proteins via directed evolution3 or computational prediction of hotspots4, as well as prediction of single-point5 or multiplepoint mutations6. Recently, the computational approaches have been increasingly applied for protein stabilization as they allow fast and focused redesign of existing proteins, requiring only limited resources and experimental effort. However, some of the beneficial mutations observed in these studies 1

ACS Paragon Plus Environment

ACS Catalysis 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 10

Table 1: Stability and Functional Characteristics of the DhaA Variants Variant

Mutations

Tm (°C) b

DhaAwt

-

a

∆Tm (°C) b a

50.4 ± 0.5

-

58.4 ± 0.2

+8.0

Specific activity (nmol.s-1.mg-1) c 18.0 ± 0.1

Evolution-based triple-point variant d DhaA101

E20S+F80R+A155P

49.3 ± 0.6

Evolution-based single-point variants d DhaA123

E20S

57.6 ± 0.2

+7.2

26.7 ± 0.5

DhaA124

F80R

49.9 ± 0.3

−0.5

27.3 ± 1.3

DhaA125

A155P

51.5 ± 0.2

+1.1

25.1 ± 1.3

65.1 ± 0.1

+14.7

5.5 ± 0.1

Energy-based multiple-point variant e DhaA112

C128F+T148L+A172I+C176F+ D198W+V219W+C262L+D266F

Combined multiple-point variant d,e DhaA115

E20S+F80R+C128F+T148L+

73.5 ± 0.1

+23.1

5.6 ± 0.1

DhaA115 monomer

A155P+A172I+C176F+D198W+

73.4 ± 0.1

+23.0

ND

DhaA115 dimer

V219W+C262L+D266F

73.4 ± 0.1

+23.0

ND

a

b

c

not applicable; determined by differential scanning calorimetry; activity determined with 1-iodohexane at 37oC and pH 8.6; designed based on phylogenetic analysis6; e designed based on force-field calculations6; Tm - apparent thermal transition midpoint or melting temperature; ND, not determined. d

stability prediction due to incorrect modeling of the evolutionbased mutations by the energy-based approach, (H2) unanticipated stabilization due to changes in the oligomeric state, (H3) changes in formation of protein-ligand complexes and (H4) entropy-driven stabilization (Figure 1). In order to investigate these hypotheses, we have determined the crystal structure of DhaA101 and examined the biophysical properties of this protein variant and its deconstructed single-point mutants. Additionally, more accurate molecular modeling approaches were applied to perform free energy predictions to explore the role of non-proteinaceous components in the studied system. Overall, our data highlight the strong added value of computational tools, such as FireProt, which combine energy-based force-field calculations with evolutionary analysis6. This combination indirectly considers both enthalpic and entropic factors for protein stabilization, allowing bigger stability improvements in contrast to cases where only force-field calculations are applied, naturally biasing predictions towards the enthalpy changes. Moreover, evolution-based methods have a good potential to uncover stabilizing mutations with additional beneficial features, such as increased activity11, as we have also observed. RESULTS Deconstructing contributions to DhaA101 stability. We first constructed the single-point mutants corresponding to the three mutations in DhaA101 to determine whether the triplepoint mutant might have unique properties that would explain the observed stabilization. Each of the single mutations had previously been predicted to be mildly stabilizing at best by FoldX6; additional calculations using various force fields showed that these approaches consistently failed to reproduce the stabilizing effects of the combined mutations (Table S1). To determine whether individual mutations might display a different outcome from their combination, we followed protein unfolding upon thermal denaturation of the single-point mutants by Differential Scanning Calorimetry (DSC), previously shown to be an appropriate method for this purpose12,13. In

Figure 1: Overview of studied haloalkane dehalogenase DhaA variants and schematic representation of different working hypotheses (H1-H4). DhaA101 contains three evolution-based mutations predicted by phylogenetic analysis, DhaA112 contains eight energy-based mutations predicted by force-field calculations and DhaA115 combines all eleven mutations.

though the impact of these mutations had been predicted by various computational tools as destabilizing, neutral or slightly stabilizing at best, the triple-point mutant displayed apparent melting temperature eight degrees higher than the wild type. Here we explore potential explanations for the discrepancy between our computational and experimental data through consideration of multiple hypotheses, including: (H1) failed 2

ACS Paragon Plus Environment

Page 3 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Catalysis three evolution-based mutations, we determined the crystal structure of DhaA101 to look for unpredicted changes in the backbone and side chain orientations (Figure 2a). The structure of DhaA101 was solved to the resolution of 0.99 Å by molecular replacement (Table S2). The search model used for the structure determination was DhaA14 (PDB 3G9X)14. The final structural model contains residues 4–296 and shows that the enzyme exists as a monomer in the crystal. The overall structure of DhaA101 resembles that of DhaAwt, namely the typical α/β hydrolase core domain and a helical cap domain. Circular dichroism (CD) spectra also indicated no secondary structure changes of the mutants6. Based on the superimposition of the crystal structures and the Rosetta predicted structure of DhaA101, it was clear that the mutations do not induce significant changes in the protein’s backbone and that the mutant’s backbone was thus correctly modeled (Figure 2c). The RMSD values of the backbone atoms for the overlap of both crystal structures and the overlay of both DhaA101 structures are 0.104 Å and 0.152 Å, respectively. Comparison of the modeled and crystal structures of DhaA101 further reveals that two out of three mutations (E20S and A155P) were correctly predicted (Figure 2d,e). Only the side chain of the newly introduced arginine (F80R) was partially misplaced in the model (Figure 2f). The Cζ-Cζ distance d between the new Arg and the naturally present Arg204 in the model was 6.6 Å, whilst in the determined structure these are stacked closely together, with a Cζ-Cζ distance of only 3.8 Å. In addition, the angle between the planes formed by the guanidinium groups of these arginines (θ) is 35.1° in the crystal structure compared to 52.3° in the modeled structure. Hence, the geometry of the arginine pair observed in the crystal structure is close to parallel (d ≤ 4.0 Å but not θ ≤ 30°)15, whereas that in the model is neither parallel nor perpendicular (neither 60°C ≤ θ ≤ 90° nor 5.0 Å ≤ d ≤ 6.0 Å)15. The repulsions between the two positively charged side chains of the arginines are most likely the reason for the partially incorrect positioning of the newly introduced arginine side chain in the model (Figure 2f). We concluded that the modeled and crystal structures showed an overall very good agreement, suggesting that the source of the discrepancy must lay elsewhere. Determination of the crystal structures of the individual mutants would likely not reveal any additional valuable information for our study. Therefore, we decided to focus on the other hypotheses. Ligands and ions do not impact protein stability. In examining the crystal structure of DhaA101, we further noted strong electron density for potential ligands in the active site and at two positions on the protein surface. Three ligands were subsequently modeled in the structure of DhaA101 (Figure 2a). The ligand bound in the enzyme active site was interpreted as 2-(N-morpholino)-ethanesulfonic acid (MES) and occupies two alternative conformations A and B with the occupancy of 0.6 and 0.4, respectively (Figure 2b). The MES position in the active site is stabilized by interactions with halidestabilizing amino acids (Trp107 and Asn41), nucleophile (Asp106) and Phe149. Phe149 adopts two different conformations depending on conformations of MES in the enzyme active site. Both ligands located on the protein surface were identified as di(hydroxyethyl)ether (PEG), designated as PEG1 and PEG2 (Figure S1). The co-crystallization of these ligands with DhaA101 could play a role in stabilizing the protein if these molecules were bound by DhaA101 more generally. However, all three ligands were components of the

Figure 2: Crystal structure of DhaA101. (a) Overall structure of DhaA101 (PDB ID 5FLK) shown in cartoon representation and the three ligands identified in the crystal structure shown as sticks. (b) Expanded view of the DhaA101 active site. 2Fo-Fc electron density map contoured at 1σ is shown for the active site residues with bound MES ligand. The residues of catalytic pentad and MES are represented as sticks. Two alternative conformations of MES and Phe149 illustrate how the side-chain affects binding modes of MES inside the active site. (c-f) Comparison of the crystal structure of DhaA101 with the crystal structure of DhaAwt and the theoretical model of DhaA101. (c) Superimposition of the crystal structures of DhaAwt (PDB ID 4E46) and DhaA101 (PDB ID 5FLK). The root-mean-square deviation (RMSD) of the alignment of the crystal structures is 0.104 Å. The Cα of mutations are represented as magenta spheres. (d-f) Expanded view of the overlays of the introduced mutations E20S (d), A155P (e) and F80R (f) and their close surroundings. (c-f) The crystal structure of DhaAwt is given in grey with the respective mutated residues in magenta sticks. The modeled and crystal structures of DhaA101 are shown in dark blue and cyan, respectively. The mutated and interacting residues are shown in sticks. contrast to force-field calculations, DSC data showed stabilization for two out of three evolution-based mutants (Table 1). These data thus confirmed that the prediction failure is already occurring at the level of the individual mutations, encouraging us to undertake more detailed investigations of these sites and their contributions in the structure-function relationships of DhaA101. The tertiary structure of DhaA101 is correctly predicted. To test the first hypothesis that incorrect modeling of the mutant structure led to the miss-assignment of the roles of the 3

ACS Paragon Plus Environment

ACS Catalysis 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 10

Figure 3: Occupancy density of sodium ions over molecular dynamics simulations with DhaAwt and two DhaA variants. Ion densities displayed above the structures of DhaAwt (a), the single-point mutant DhaA123 (E20S) mutant (b), and the triple mutant DhaA101 (E20S+F80R+A155P) (c). The blue occupancy surfaces correspond to the 0.025 isovalue and the cyan ones to the 0.005 isovalue.

crystallization solution and not present in the buffers used for the activity and thermostability assays, suggesting they are non-native and highly unlikely to influence the observed stability. Broadening our consideration of small molecules that might influence protein stability further, we noted that solvated ions had not been taken into account during our initial force-field calculations and could potentially be relevant factors that influence protein stability16. Indeed, the distribution of the charges on the protein surface is a crucial determinant of protein-ion interactions. E20S and F80R mutations removed a negative surface charge and introduced a positive surface charge, respectively. The occupancy of cations calculated for 200 ns-long molecular dynamics (MD) simulations revealed significant differences of the regions with the highest density of sodium ions compared to the wild type (Figure 3a) and the evolution-based triple-point mutant DhaA101 (Figure 3c), over the entire simulations. Deconvoluted single-point mutants were also analyzed and MD trajectories revealed significant differences due to E20S mutation (Figure 3b), whereas no changes where observed for the other two mutations (F80R and A155P). The largest change in sodium ion density was observed for the E20S mutation, which had high sodium ion density around the wild type glutamate, but which was strongly reduced in the presence of the serine mutation, even showing a long-distance influence. Similarly, the introduced positive charge of F80R also significantly reduced the high density of sodium ions in its region (Figure 3c). On the other hand, the occupancy of the chloride ions is much more dispersed than that of sodium ions. Despite a slight chloride ion density increase around the R80 residue in the mutant, this effect is negligible compared to the effect of sodium ion density. These observations were corroborated by the radial distribution functions of the ions around the proteins (Figure S2), and suggest that the different interactions with the sodium ions in solution might be related to the observed stabilization. DSC thermograms were, therefore, also determined at increasing concentrations of sodium chloride. The results obtained did not reveal any clear dependence of the stability of those proteins on the ion concentrations (Table S3), thus rejecting this hypothesis.

Partial oligomerization occurs but does not influence stability. Our initial energetic calculations on the evolutionbased mutations were carried out using monomeric structures. However, if the introduced mutations at the surface of the protein are involved in the formation of multimeric structures, then our energetic predictions using monomers could be misleading. Therefore, we examined the oligomeric state of DhaA101 as well as two mutants that contained eight energybased mutations: (i) DhaA112 designed by the energy-based protocol of the FireProt tool6 and (ii) DhaA115 combining three evolution-based and eight energy-based mutations (Table 1). Analysis by native polyacrylamide gel electrophoresis demonstrated that the mutants were mainly present as monomers, though a limited amount of higher order oligomer formation was also observed (Figure S3a). The dimeric form of DhaA101 was not present in sufficient quantities to allow isolation, but the dimeric and monomeric forms of DhaA115 were successfully separated by gel permeation chromatography, as was confirmed by native polyacrylamide gel electrophoresis (Figure S3b,c). The lack of a rapid equilibrium between the monomeric and oligomeric states allowed the DSC analysis of the separated monomer and dimer fractions to determine whether the minor dimer population (~8%) was skewing the observed melting temperature (Figure S3b). Stability analysis of the separated monomer and dimer fractions via DSC revealed the same melting temperatures for the different oligomerization states (Table 1). For the DhaA115 dimer, small additional heat capacity was observed before the main peak, most likely resulting from dimer dissociation (Figure S3d). Therefore, we ruled out oligomerization as the source of the computational mismatch. Calorimetry identifies changes in entropy of evolutionbased mutants. Finally, we have carried out the analysis of calorimetric data to investigate the potential role of entropy in the stabilization. The deconvolution of the DSC data allowed the determination of the apparent melting temperature (Tm) and the calorimetric enthalpy change (∆Hcal) of each variant (Figure 4). The introduction of the energy-based mutations (DhaA112) to the wild type substantially increased both the ∆Hcal and Tm, indicated by an increased area under the DSC thermogram peak and its shift to higher temperatures, respec4

ACS Paragon Plus Environment

Page 5 of 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ACS Catalysis

tively (Figure 4a), and thus serves as a positive control for enthalpic stabilization. The introduction of the evolution-based mutations (DhaA101), however, had only a minor effect on the area under the peak while still shifting it to higher temperatures. Combination of the two sets of mutations (DhaA115) resulted in the combination of the observed effects, i.e., further increase in Tm with the change in ∆Hcal similar to that of the energy-based mutations alone (DhaA112). A similar picture appears when analyzing the thermograms of single and double-point precursors of DhaA101: a mixed direction of changes in Tm and enthalpy differences was observed upon the introduction of the corresponding mutations (Figures 4b and S4). These observations illustrate that changes in Tm induced by evolutionary mutations are poorly correlated with the corresponding changes in experimentally obtained enthalpy, which is in strong contrast to energy-based mutations. The above observations imply a significant role of entropy for the evolution-based mutations. Calculation of the entropy change from DSC thermograms is difficult. Privalov and Dragan suggested using the Kirchhoff's relation ∂∆H/∂Τ=∆Cp and ∂∆S/∂T=∆Cp/T, where ∆Cp stands for the difference between the heat capacities of the folded and unfolded protein states17. Following their procedure, we calculated the apparent entropy changes (∆Scal) by integration of the melting curves (Figure S4). The evolutionary-based mutations decrease ∆Scal (DhaAwt versus DhaA101; DhaA112 versus DhaA115), whereas an increase in ∆Scal is observed upon introduction of the energy-based mutations as expected due to the common enthalpy-entropy compensation18 (DhaAwt versus DhaA112; DhaA101 versus DhaA115). Thus, the result of the direct integration of the DSC curves supports our conclusions regarding the role of entropy in the stabilization of DhaA by the evolution-based mutations. Nevertheless, the calculated values of the entropy changes are only approximates since we noted a significant sensitivity of the method17 to the baseline and peak border determination, and also different degrees of reversibility. A small decrease in reversibility of unfolding was observed after introduction of the evolutionarybased mutations (86% for DhaAwt and 65% for DhaA101), whereas a more significant decrease in reversibility was observed upon introduction of the energy-based mutations (90 % were filtered out. Maximum of 200 sequences were randomly selected for the construction of the alignment. Back-to-consensus mutations were identified in the alignment by both simple consensus and frequency ratio approaches. For the validation set, only the records with Gibbs free energies changes of folding 0.5 kcal.mol-1 were included, omitting neutral mutations and in this way increasing the reliability of the analysis. In total, we have extracted 1,328 mutations, from which 256 (19.3 %) mutations were stabilizing and 1,072 (80.7 %) mutations were destabilizing. A total number of 1,099 potentially stabilizing mutations was identified for these 103 proteins by back-to-consensus analysis, from which 515 (46.9 %) were predicted as stabilizing by the energy-based approach, while in more than a half of the cases (53.1 %), the back-to-consensus and force-field methods disagreed in the identification of the stabilizing mutations (Table S7, Supporting Data). The experimental and back-to-consensus datasets have overlapped in nine of the cases. Whilst all nine mutations were correctly assigned as stabilizing by the evolutionary analysis, only six of the mutations were predicted as stabilizing by the force-field calculations. This bioinformatics analysis highlights the complementary information derived from the energy-based and the evolution-based analyses. DISCUSSION Our combined structural, biophysical and computational analysis allowed resolution of the factors that led to the failure of our energy-based calculations to predict the outcome of evolution-based mutations. The first hypothesis – that the mutant structure was incorrectly modeled – was ruled out by good agreement between the model and crystal structure of both the protein backbone and the positions of the E20S and A155P side chains. Only the side chain from the newly introduced arginine (F80R) was partially incorrectly positioned in the model, most likely due to prediction of a repulsive interaction between the new arginine R80 and the neighboring R204. We can conclude that the modeled structure resembled the crystal structure well and that the minor differences observed cannot explain the incorrect force-field predictions. 5

ACS Paragon Plus Environment

ACS Catalysis 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 10

Figure 4: Thermostability analysis of DhaA variants. (a) The calorimetric curves of DhaA variants show increasing apparent melting temperatures (Tm,app, marked with ♦) upon introduction of evolutionary (red), energy-based mutations (green) and their combination (purple) compared to the wild type (blue). The area under the calorimetric curve, corresponding to calorimetric enthalpy (∆Hcal), is large for variants containing energy-based mutations, and comparable to the wild type in case of DhaA101 which carries evolutionary-based mutations. (b) ∆Hcal calculated via integration of the peak area under the DSC thermograms using linear baseline subtraction versus Tm,app for single and double-point evolutionary variants. Although the single-point mutants demonstrate mixed effects in terms of enthalpy changes, the double-point variants follow the pattern of DhaA101: an increase in melting temperatures with minor increase or even decrease of enthalpy change as compared to the wild type. The wild type Tm,app level is indicated by the dotted red line.

late29 and are generally approximated or neglected30, often resulting in incorrect predictions of entropy-driven mutations6. The entropic effects require extensive conformational sampling, while tested computational methods work well for rigid systems that involve changes in an enthalpic term. Hence, we compared the calorimetric enthalpy levels of the different DhaA variants. It is clear from the data that the evolutionarybased mutations have only limited effect on the enthalpy difference, whereas the energy-based mutations clearly show enthalpy-driven stabilization. This suggests an entropic nature of the stabilization due to the mutations inferred by phylogenetic analysis. Although approximate, entropy calculations revealed that the entropy change curves vary in both the respective position and the slope for the evolutionary triple mutant. Moreover, DhaA variants with single and doublepoint evolutionary mutations also showed a more pronounced entropy-related effect as compared to the energy-based mutants. The entropy change curves of the energy-based and combined mutants were shifted to higher values compared to the wild type, which is in agreement with the enthalpy-entropy compensation31. Nonetheless, accurate estimation of the entropy of large systems currently remains very challenging and the results can be biased by available sampling29. TI was the computational method providing the most consistent result with the experimental observations. However, recent large-scale screening studies show that especially for charge-changing mutations the predictive capabilities of alchemical free energy simulations are less reliable32,33. For complex systems such as the one studied here the results are most likely influenced by both sampling issues and force field bias. Implementations of alchemical transformations making use of graphics processing units34 will greatly enhance the accessible time scales of MD simulations. While this reduces the sampling problem in freeenergy simulations, the reduction of the force field bias remains a topic of active research in the coming years, making these calculations demanding if several force fields are to be

Our second hypothesis – that oligomerization was influencing stability – can also be ruled out. We do observe a small fraction of oligomeric protein, presumably because the change in surface charge of the more stable variants could remove potential repulsive interaction between monomers that prevent the wild type from oligomerization. However, the thermostability analysis of the isolated dimer and monomer fractions of DhaA115 revealed no improved stability due to dimerization. The small additional heat capacity peak observed for the DhaA115 dimer most likely results from dimer dissociation. Since only a limited fraction (