Isotope Substitution of Promiscuous Alcohol Dehydrogenase Reveals the Origin of Substrate Preference in the Transition State

Abstract The origin of substrate preference in promiscuous enzymes was investigated by enzyme isotope labelling of the alcohol dehydrogenase from Geobacillus stearothermophilus (BsADH). At physiological temperature, protein dynamic coupling to the reaction coordinate was insignificant. However, the extent of dynamic coupling was highly substrate‐dependent at lower temperatures. For benzyl alcohol, an enzyme isotope effect larger than unity was observed, whereas the enzyme isotope effect was close to unity for isopropanol. Frequency motion analysis on the transition states revealed that residues surrounding the active site undergo substantial displacement during catalysis for sterically bulky alcohols. BsADH prefers smaller substrates, which cause less protein friction along the reaction coordinate and reduced frequencies of dynamic recrossing. This hypothesis allows a prediction of the trend of enzyme isotope effects for a wide variety of substrates.

chloride (for the heavy enzyme). The enzyme was purified in 50 mM sodium phosphate buffer containing 5 mM 2mercabtoethanol, pH 7, using Q-sepharose anion exchange column then loaded into cibacron blue affinity column.
The enzyme was eluted from the final step using NAD + then stored at 4 °C. A Superdex G25 column was used to remove the cofactor NAD + before running kinetics. Purity and mass of the enzyme were confirmed by SDS-PAGE and LC-MS.

Staedy state Kinetics.
Kinetics has been performed on Jasco 2660 UV spectrometer in 25 mM sodium phosphate buffer, pH 7. The cofactor NAD + concentration has been determined at 260 using extension coefficient 17800 M -1 cm Then, $%& were calculated using non-linear fitting, where is the velocity $%& is the steady state turn over number, is the enzyme concentration, and are the cofactor NAD + and the substrate concentration, respectively, ,% is the dissociation constant of the E.NAD + binary complex and -0 and -. are the Michaelis constants for the cofactor NAD + and the substrate respectively.
For the enzyme kinetic isotope effect, the kinetics of the 'light' and 'heavy' enzymes have been performed at the same time and all the experiments has been repeated at least three times, then errors have been calculated form the data standered deviations.

Computational Details
The simulation model. The starting structure for dynamic simulations of ADH was obtained from the Protein Data Bank entry 1RJW [2] which codes for the crystal structure of the enzyme bounded to the inhibitor trifluoroethanol. The substrate isopropanol was introduced instead of the inhibitor while the cofactor NAD + was introduced manually based on the coordinates of the ligand-bound of other ADH with and entry 1ADC. [3] The PROPKA3 program [4] was employed to estimate the pKa values of the titratable protein residues to verify their protonation states at pH 7. To neutralize the system, 2 sodium counterions were placed in optimal electrostatic positions around the enzyme. Finally, the system was solvated using a cubic box of TIP3P water molecules with side lengths of 130.0 Å; water molecules with an oxygen atom within 2.8 Å of any heavy atom were removed.
The full systems were formed by 176206 atoms when isopropanol was the substrate, 176210 for benzyl alcohol and 176203 for ethanol. The model contained the protein (1356 residues for the whole tetramer; thus, 339 residues per subunit plus 8 zinc atoms) with 20664 atoms, the substrate and cofactor (70 atoms), 2 sodium ions and 51817 water molecules, 332 observed in the crystal structure and 51485 added during solvation (a total of 155451 atoms). After 500 minimization steps using a conjugate gradient method, 10 ns of molecular dynamics simulations using NAMD [5] software were carried out in order to equilibrate the system at 298 K. In the following simulations the whole system was divided into a QM part and a MM part to perform combined QM/MM calculations ( Figure 1). The quantum subsystem contained 61 atoms, including the full substrate, the zinc ion and parts of the cofactor (nicotinamide ring and the ribose), Cys38, Cys148 and His61. Four hydrogen 'link' atoms [6] were used to saturate the valence at the QM-MM boundary ( Figure 1). The quantum atoms were treated by the AM1 Hamiltonian, [7] modified using specific reaction parameters (denoted as AM1-SRP) developed previously for DHFR. [8] The protein atoms and the ions were described by OPLS-AA [9] force field while the water molecules were described by the TIP3P potential. [10] Cutoffs for the nonbonding interactions were applied using a switching function within a radius range of 13.0 to 9.0 Å. All the molecules further than 25 Å from the substrate were frozen, while the rest of the system was allowed to move. From the final structure, the substrate was substituted for benzyl alcohol and ethanol. After, 5000 steps of  Potentials of Mean Force (PMF). [11] From the final structures of the QM/MM simulations at each temperature, onedimensional PMFs, W CM , were computed using the antisymmetric combination of distances describing the hydride transfer, ξ = d C1Ht -d HtC4 , as the reaction coordinate. The umbrella sampling approach was used, [12] with the system restrained to remain close to the desired value of the reaction coordinate by means of the addition of a harmonic potential with a force constant of 2500 kJ mol -1 A -2 , which allows good overlap between windows. The reaction coordinate was then explored in a range from -1.93 to 1.36 Å and the total number of windows was 48. The probability distributions obtained from MD simulations within each individual window were combined by means of the weighted histogram analysis method (WHAM). [13] 100 ps of relaxation and 100 ps of production MD, with a time step of 0.5 fs, in the canonical ensemble (NVT, with reference temperatures at 293, 298, 303, 308, 313, 318 and 323 K) and the Langevin-Verlet integrator, [14] were used in the simulations.

S9
Calculation of the recrossing transmission coefficient. Grote-Hynes (GH) theory can be applied to describe the evolution of the system along the reaction coordinate at the TS. In particular, the transmission coefficient can be obtained as the ratio between the reactive frequency and the equilibrium barrier frequency: [15] (S1) with the equilibrium frequency derived from a parabolic fit of the PMF around the free energy maximum and the reactive frequency w r is obtained via the GH equation: [16] (S2) z TS (t) is the friction kernel obtained at the TS, assuming that recrossings take place in the proximity of this dynamic bottleneck: were saved at each simulation step. We previously tested that the GH approach gives transmission coefficients in very good agreement with those obtained from activated trajectories initiated at the TS ensemble. [19], [20] The recrossing transmission coefficients γ(T,ξ) were calculated using eq. S1 for the light and heavy versions of ADH.
Because the heavy enzyme has larger mass some of its internal motions are slower. Then, this version of the enzyme presents a higher friction and thus a smaller value of the transmission coefficient. The transmission coefficients of the two versions were found to be statistically different.

Temperature dependence of the recrossing transmission coefficient
As γ < 1, the first term R·ln(γ) in Eq. S4 is negative because of the contribution of substrate and protein motions, leading to an increase of free energy in barrier crossing. These degrees of freedom must reach a particular value at the transition state, such that the entropy of the system is reduced with respect to the equilibrium description. is negative and thus thermal activation contributes to an increase of the mass-dependent entropic barrier (SI , Table S6).
For isopropanol, the transmission coefficient in the heavy enzyme is slightly lower than that in the light counterpart at all the temperatures analyzed (SI Figure S8). This is likely due to the decrease in the frequencies of the protein motions and an increase in the friction on the reaction coordinate. The transmission coefficients for both light and heavy enzymes decrease smoothly with temperature; the resulting slopes A@ A?
are therefore similar in magnitude (SI, Figure S8 and Table S6). The measured enzyme KIE is only slightly larger than unity and also largely temperature independent ( Figure 2C). In contrast, when benzyl alcohol is used, a significantly different observation was made.
The recrossing coefficients γ in both the light and heavy enzymes decrease with temperature, but the light enzyme shows a steeper decrease (i.e. a greater magnitude of A@ A? , SI, Figure S8 and Table S6). Particularly, at low temperatures the recrossing coefficients computed in the light enzyme are significantly larger than those of the heavy enzyme (SI, Figure S8). Accordingly, the computational enzyme KIE ( Figure 2C) is temperature dependent and large at low temperatures.

Ensemble Averaged Variational Transition State Theory. Acording to the Ensemble-Averaged Variational Transition
State Theory (EA-VTST), [21] the theoretical estimation of the rate constant can be obtained as shown in Eq 1. of the manuscript. is the quasiclassical activation free energy at the transition state, obtained from the classical mechanical (CM) PMF (W CM (T, ξ)) as: where ΔW vib (T, ξ*) corrects W CM (T, ξ*) to account for quantized vibrations orthogonal to the reaction coordinate along which the PMF is defined, ξ at the maximum of the PMF, ξ*; ΔW vib,R (T) corrects W CM (T, ξ R ) for quantized vibrations at the reactant side minimum of the PMF, ξ R , and G CM R,T,F is a correction for the vibrational free energy of the reactant mode that correlates with motion along the reaction coordinate. [21a] To correct the classical mechanical PMF, W CM , normal mode analyses were performed for the quantum region atoms.
To perform these calculations we localized 15 TS structures starting from different configurations of the corresponding simulation windows. After tracing the minimum energy path, we optimized 15 reactant structures and obtained the Hessian matrix for all the stationary structures.
[22] The final quantum mechanical vibrations correction to the quasi-classical activation free energy was obtained as an average over these structures.
The tunneling transmission coefficients were calculated with the small-curvature tunneling (SCT) approximation, which includes reaction-path curvature appropriate for enzymatic hydride transfers. The final tunneling contribution (see Tables S4 and S5) is obtained as the average over the reaction paths of 15 TS structures, as made in our previous works. [20,23] QC act