A Dynamic Equilibrium of Three Hydrogen-Bond Conformers Explains

Sep 14, 2016 - We used the PROPKA server to determine the positions of the remaining protons.(31) The deprotonated pCA chromophore and the covalently ...
2 downloads 12 Views 2MB Size
Article pubs.acs.org/JCTC

A Dynamic Equilibrium of Three Hydrogen-Bond Conformers Explains the NMR Spectrum of the Active Site of Photoactive Yellow Protein Phillip Johannes Taenzler, Keyarash Sadeghian, and Christian Ochsenfeld* Chair of Theoretical Chemistry, Department of Chemistry, University of Munich (LMU), Butenandtstrasse 7, D-81377 Munich, Germany Center for Integrated Protein Science (CIPSM) at the Department of Chemistry, University of Munich (LMU), Butenandtstr. 5−13, D-81377 Muenchen, Germany S Supporting Information *

ABSTRACT: A theoretical study on the NMR shifts of the hydrogen bond network around the chromophore, paracoumaric acid (pCA), of photoactive yellow protein (PYP) is presented. Previous discrepancies between theoretical and experimental studies are resolved by our findings of a previously unknown rapid conformational exchange near the active site of PYP. This exchange caused by the rotation of Thr50 takes place in the ground state of PYP’s active site and results in three effectively energetically equal conformations characterized by the formation of new hydrogen bonds, all of which contribute to the overall NMR signals of the investigated protons. In light of these findings, we are able to successfully explain the experimental results and provide valuable insight into the behavior of PYP in solution. We further investigated related PYP mutants (T50V, E46Q, and Y42F), and found the same conformational exchange in E46Q and Y42F to be responsible for the experimentally observed NMR and UV/vis spectra.



INTRODUCTION Photoactive Yellow Protein (PYP) serves as a putative bacterial photoreceptor to Halorhodospira halophila, which enables the organism to avoid short-wavelength light via negative phototaxis.1 The active/signaling state of PYP is triggered by the photoinduced trans−cis isomerization of the negatively charged and covalently bound chromophore, p-coumaric acid (pCA). While extensive experimental and theoretical investigations2−13 to date were concentrated on the signaling state formation, the hydrogen bonding network involving the pCA chromophore in the inactive/dark-adapted state of PYP has been a matter of much debate and remains poorly understood. Numerous studies have dealt with the geometric arrangement of the hydrogen bond (HB) donors in PYP’s active site (shown schematically in Figure 1),14,15 namely the Glu46 and Tyr42 residues, as well as the strength and type of the HBs found between these residues and pCA.16−18 The Glu46 proton in particular has been the subject of considerable debate, as the formation of a so-called low barrier hydrogen bond (LBHB) with pCA was proposed.18 This idea is motivated by the crystal structure obtained via neutron diffraction experiments where (i) a short O(Glu46)-O(pCA) distance of 2.56 Å was observed and (ii) the Arg52 residue was assigned as a neutral species.18 Yamaguchi et al.18 proposed that in the absence of the positive © XXXX American Chemical Society

charge on the Arg52 residue, the formation of a very strong hydrogen bond stabilizes the embedded negative charge on the chromophore, which then resulted in the aforementioned shortening of the O−O distance. However, theoretical studies have questioned the deprotonation of Arg52.19 The formation of an LBHB in PYP is challenged by both experimental and computational NMR data. While the experimentally observed O(Glu46)−O(pCA) distance of 2.56 Å is within the 2.5−2.6 Å range expected of LBHBs,20−22 the experimentally assigned NMR peak of 15.2 ppm for the chemical shift of the Glu46 proton23,24 lies outside the characteristic 16−20 ppm range20,21,25 expected of LBHB protons. Although the previously computed 1H NMR shifts of 18.716 and 19.714 ppm for the Glu46-H in the neutron diffraction structure are consistent with the idea of LBHB formation in the crystalline phase, structural relaxation effects, however, led to conformations with more upfield NMR shifts in the range of 13.7−16.4 ppm.14,16 Overall, the experimental and theoretical data do not support the existence of an LBHB in PYP under solution conditions, hence questioning its relevance to the biological function of PYP. Received: April 9, 2016

A

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

range from 13.7 to 16.4 ppm, which could encompass the experimentally determined NMR shifts of both protons. Kanematsu and Tachikawa30 employed a combination of ONIOM and multicomponent QM approaches to compute the chemical shifts. Their study clearly shows that 1H NMR shifts are sensitive to how the system is treated at the QM level with regards to the choice of method. While their assigned peak difference of 2.3 ppm appears reasonably close to the experimental peak difference of 1.5 ppm, the assigned values of 15.5 (Tyr42−H) and 17.8 (Glu46−H) ppm30 are shifted downfield by over 2.5 ppm relative to the experiment. In this work, we demonstrate that the measured 1H NMR peaks assigned to the Glu46 and Tyr42 protons of PYP’s active site do not arise from a single conformation, as was assumed in the previously cited studies. Instead, a fast conformational exchange triggered by the rotation of the side chain of the neighboring Thr50 residue allows the active site to cycle through three conformers. Whereas the possibility of Thr50 inducing multiple conformations via this mechanism has been fully neglected in the previous computational investigations, our QM/MM calculations show that the geometric perturbations caused by its rotation couple with and consequently modulate the 1H NMR shifts of the Tyr42 and Glu46 protons. Furthermore, we demonstrate that over 400 atoms of the active site must be included in the quantum chemical calculations to ensure QM-size convergence. Finally, the QM/MM UV/vis spectra, computed at the coupled-cluster linear-response level of theory (LRC−CC2/MM), show that the vertical excitations are not affected by the Thr50 rotation, implying that all three conformations are relevant for the signaling state formation in PYP.

Figure 1. Schematic view of the active site of PYP.

Despite numerous studies dealing with the assignment of the H NMR peaks in PYP, the experimental and computational data are still not compatible. Experimental 1H NMR shifts of 15.2 and 13.7 ppm, assigned to the Glu46 and Tyr42 protons,23 respectively, are separated by 1.5 ppm. However, none of the computational reports to date were successful in explaining the observed chemical shifts. For example, Saito and Ishikita14 took the crystal structures with the PDB codes 2ZOI18 (neutron diffraction structure) and 1OTB26 (X-ray structure), optimized them at the quantum mechanical/molecular mechanical (QM/ MM) level, and computed the proton NMR shifts. For the Glu46 proton, the authors assigned 14.5 (2ZOH) and 14.6 ppm (1OTB). For the Tyr42 hydrogen, the values obtained, 14.6 (2ZOH) and 14.0 ppm (1OTB), are (i) far too downfield and (ii) overlapping with those of Glu46. Nadal-Ferret et al.16 performed detailed investigations on the neutron diffraction structure of PYP18 by largely focusing on delocalization of the Glu46-H proton due to nuclear quantum effects. Such effects were also investigated in other biological systems featuring hydrogen bond networks, for example in ketosteroid isomerase,27 or stacked polyglutamine strands,28 where they were found to play crucial roles. Nuclear quantum effects were even found to impact the NMR shifts of protons in water by as much as 0.5 ppm,29 which the authors likened to the errors arising from the approximate treatment of the electronic structure.29 Although the calculations performed by Nadal-Ferret on the crystal structure of PYP were extensive, they were ultimately unable to replicate the experimental results of Sigala et al.16,23 The authors concluded that the crystal structure and the solvated protein differed sufficiently in structure to explain the differing results.16 The authors also computed the NMR spectrum of the PYP system, and sampled structures from a QM/MM (AM1/CHARMM) molecular dynamics simulation. However, they only reported NMR shifts for the Glu46 proton at snapshots featuring the shortest (2.61 Å), longest (2.71 Å), and mean (2.66 Å) observed O(Glu46)− O(pCA) distances.16 The reported NMR shifts of Glu46−H 1



COMPUTATIONAL DETAILS The crystal structure of PYP, obtained via neutron diffraction (PDB entry 2ZOI18), was used as the starting structure of our computational investigation. Only 87% of the protons were assigned in this structure. We used the PROPKA server to determine the positions of the remaining protons.31 The deprotonated pCA chromophore and the covalently attached Cys69 residues were parametrized as a single fragment using the ANTECHAMBER module of AMBER TOOLS.32 Arg52 was taken to be fully protonated (+1 charge). The system was solvated in a TIP3P water sphere of 35 Å radius centered at V120, and Na+ ions were added to neutralize the overall negative charge. The whole system was briefly equilibrated for 25 ps and was subsequently minimized at the force-field level and used as the starting point for subsequent QM/MM investigations. This was performed in the NVT ensemble, at a temperature of 300 K using Langevin dynamics (excluding hydrogens). Longer 5 ns MD simulations to investigate the impact of the surrounding environment were performed in the NVE ensemble with a temperature of 298 K. These simulations also included a spherical harmonic potential centered at V120 with a radius of 35 Å to maintain the shape of the system. All MD simulations were performed with a 1 fs time step. We employed two different QM regions for the QM/MM structural optimizations of the wild type (WT) and T50V and E46Q mutants, labeled here as QMS and QML, and a larger QM region (QMXL) for the Y42F mutant. To showcase sufficient convergence of the QM/MM NMR data, we chose seven QM regions of increasing size ranging from 70 to 518 atoms, which are denoted as R1−R7 (detailed structural information on the individual QM regions is available in the B

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Figure 2. Representation of the T (left), F (middle), and G (right) conformers. The residues pCA, Tyr42, Glu46, and Thr50, as well as existing hydrogen bonds are shown explicitly.

relaxed scan QM/MM calculation along the Glu46- and Tyr42pCA proton transfer coordinates. Along both paths, we observed the rotation of the Thr50 side chain, which was accompanied by a significant lowering of the relative energy. Free QM/MM minimization of the resulting structures led to the final T and G conformers. The stabilities of all three conformers were confirmed via numerical frequency analysis based on the QMS region and the def2-SVP basis set using B3LYP-D3/AMBER. In all three conformations, the Glu46 and Tyr42 residues remain protonated and retain their hydrogen bonding interactions with the chromophore. In all cases, the O(pCA)−O(Glu42/Tyr46) distances are in the 2.56−2.66 Å range, and the corresponding hydrogen bonds can be classified as short ionic.21,22 The impact of the mobile environment on the NMR shifts of the central protons in each conformer was investigated by sampling force-field MD simulations. The details are included in the Supporting Information. The relative QM/MM energies, presented in Table 1, indicate that the T, F, and G conformers are all energetically

Supporting Information). The remainder of the protein and solvent molecules outside these QM regions, containing up to 16 916 atoms, were treated at the MM level (AMBER forcefield). A persistent active region was chosen for all QM/MM structural optimizations, such that all residues (including water) having at least one atom within a 7.5 Å radius around the pCA chromphore were relaxed (787 atoms in total). All QM structural optimizations and NMR calculations were performed with the FermiONs++ program package developed in our group.33,34 The QM/MM structural optimizations were carried out at the B3LYP-D3/AMBER level,35,36 including the empirical dispersion correction proposed by Grimme.37 The def2-SVP and def2-TZVP basis sets were used to optimize the structures at the B3LYP-D3/AMBER level. The QM/MM NMR calculations were performed at the B97-2/AMBER38 level. A recent benchmark study showed this DFT functional to be highly reliable.39 Tetramethylsilane (TMS) was used as a reference to obtain the NMR shifts, and the pcS family of basis sets40 was used for the NMR calculations. The vertical excitation energies were calculated using the RI-CC2 implementation of Turbomole,41 with the aug-cc-pVDZ basis set. The mutants were generated by directly mutating the corresponding residues from the fully optimized (B3LYP-D, QML region, def2-TZVP basis) wild type conformers, followed by reoptimization of the resulting structures to ensure a consistent description.

Table 1. Relative B3LYP-D3/AMBER Energies (in kcal mol−1) of the T, F, and G Conformers Obtained with the QML and QMS Regions and the def2-SVP and def2-TZVP Basis Sets (All Values in kcal/mol and the T Conformer Taken As Reference)



RESULTS AND DISCUSSION Conformational Exchange Induced by the Rotation of Thr50. While the neutron diffraction structure of PYP (PDB code 2ZOI18) suggests a structural motive in which the Thr50 side chain appears to be hydrogen bonded with the phenol group of Tyr42 (denoted in this work as the T conformation), two additional conformations with differing hydrogen bonding networks were found as part of this work. These conformations are not only thermally accessible but also significantly influence the 1H NMR shifts of the Glu42 and Tyr46 hydrogens. These additional two conformers, denoted here as F and G, can be structurally characterized by the orientation of the Thr50 side chain. In the F conformation, the Thr50 side chain is hydrogen bonded with the phenolate oxygen of the pCA chromophore. In the G conformation, the side chain of Thr50 is pointing away from the chromophore and acts as a hydrogen bond donor to the carbonyl backbone of the Glu46 residue (see Figure 2 for a schematic presentation of all three conformers). We identified the F conformer via free QM/MM relaxation of the force-field minimized structure. From F, we performed a

structure

QMS/def2SVP

QMS/def2TZVP

QML/def2SVP

QML/def2TZVP

T F G

0.0 3.7 2.2

0.0 2.0 1.7

0.0 3.3 2.0

0.0 1.2 1.6

competitive. If the small QMS region (86 atoms) and the def2SVP basis set are employed for the QM/MM structural optimization, the G and F conformers are 2.2 and 3.7 kcal mol−1 above T, respectively. However, theses differences are reduced to 1.6 (G) and 1.2 kcal mol−1 (F) above T once the larger QM region QML (160 atoms) and the larger basis set def2-TZVP are used (Table 1). The QM/MM calculations reveal that the T, F, and G conformers are connected via small energetic barriers of approximately 4 kcal mol−1 (B3LYP-D3/ AMBER and def2-SVP basis set) and are therefore kinetically accessible at room temperature. This was investigated by computing the minimum energy paths connecting these conformers by employing the nudged elastic band (NEB) method. The corresponding energy profile is shown in Figure 3 (bottom). These results suggest F and G to be only intermediary states during the conformational exchange. However, as mentioned earlier, by employing a larger QM C

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

Figure 3. Top: B79-2/AMBER NMR peaks (pcS-2 basis sets) for the Glu46-H (red) and Tyr42-H (blue) calculated along the NEB paths. Bottom: B3LYP-D3/AMBER NEB paths connecting the T, F, and G minima (def2-SVP basis set, all values in kcal mol−1). The red and blue dotted lines show the experimental NMR shifts of Glu46-H and Tyr42-H, respectively.

Simultaneously, the chemical shift of Glu46−H moves slightly downfield to 15.3 ppm. As the G conformer is approached, the weak hydrogen bond between Thr50 and pCA is broken, causing both remaining hydrogen bonds to strengthen despite retaining effectively constant distances of 2.60 Å and 2.63 Å for Glu46−H and Tyr42−H, respectively. This causes both chemical shifts to move downfield to 16.0 ppm (Glu46−H) and 12.7 ppm (Tyr42−H). In returning to the T conformer, both chemical shifts approach and pass their experimental values (13.7 and 15.223) and, as mentioned previously, overlap. The key observation from our QM/MM NMR data is that, in contrast to previous studies, no single conformation is able to reproduce the experimental 15.2 and 13.7 ppm chemical shifts23 for Glu46 and Tyr42 protons, respectively; the Thr50 rotation induces three energetically competitive conformations, all of which contribute to the 1H NMR spectrum of PYP. In order to ensure that the obtained chemical shifts are converged with respect to the QM region size chosen for the QM/MM calculation, we carried out a benchmark study, which is described further below. Role of the Environment and the Coupling of the Chemical Shifts. The QM/MM 1H NMR data presented here highlight the key role played by the rotation of the Thr50 side chain. As the side chain of Thr50 cycles through the T, F, and G conformations, the electronic structure of the active site, particularly the pCA-Tyr42 and pCA-Glu46 hydrogen bonds, is altered such that the chemical shifts of the Tyr42 and Glu46 protons fluctuate in the range of 2.6 and 1.0 ppm, respectively. In order to explore the influence of the PYP environment on the discussed chemical shifts, we performed a set of NMR calculations in which the surroundings were systematically removed (see Table. 2). Comparing the QM-only and QM/ MM NMR data (QM region = R1, basis set = pcS-1) for all three conformations, the Tyr42 proton appears more sensitive toward the environment (maximum difference between QM and QM/MM: 0.57 ppm for Tyr42 vs 0.10 ppm for Glu46). Upon removing the Thr50 residue from the R1 (QM only) set of atoms in the F conformation, a downfield change of 0.16

region (QML) and basis set (def2-TZVP) for the QM/MM structure optimization, a lowering of the relative QM/MM energies for both the F and G conformers relative to the T state was observed (Table 1). Nevertheless, even by taking the maximum barrier of 4 kcal mol−1 obtained from the NEB calculations as an approximate value of the upper bound, the rate of the T-F-G-T conformational exchange can be estimated to be on the order of 109 s−1 via Eyring’s equation. With such a fast rate of interconversion, we do not expect the experimental NMR spectrum to be able to differentiate between the contribution of each individual conformation. Instead, we propose that all three conformations contribute proportionally to the 1H NMR spectrum of PYP, as discussed further below. Conformational Exchange Modulates the 1H NMR Spectrum of PYP. The QM/MM NMR data clearly demonstrate that the orientation of the side chain of the Thr50 residue has a substantial impact on the chemical shifts of the Glu46 and Tyr42 protons. The computed chemical shifts for both protons, along the paths connecting the T, F, and G conformers, are plotted in Figure 3 (top). For the T conformer, the hydrogen bond distance, defined by the donor and acceptor oxygen separation, is 2.62 Å for the Glu46−H, and 2.54 Å for the Tyr42−H hydrogen bond. The additional hydrogen bond between Thr50 and Tyr42 weakens Tyr42−H, resulting in the two discussed protons having nearly identical chemical shifts of approximately 15.0 ppm. As mentioned earlier, Saito et al.14 also computed very similar chemical shifts for both protons (14.5 and 14.6 ppm, respectively). We compared the geometric arrangement of our T conformer with the minimized structure presented by Saito et al.14 and found them to have identical structural motifs, despite different QM/MM setups employed in both studies. As the Thr50-rotation path progresses toward the F conformer, the Glu46−H distance remains constant at 2.61 Å, while Tyr42−H lengthens to 2.62 Å. The formation of the additional weak hydrogen bond between Thr50 and pCA (OThr50−OpCA 2.92 Å) and the lengthening of Tyr42−H weakens this hydrogen bond considerably, which causes a significant upfield move of its chemical shift to 12.4 ppm. D

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation

NMR chemical shifts of the remaining protons were computed. As expected, and shown in Table 2, the obtained chemical shifts show substantial differences of up to 0.52 ppm across all conformations, which is a clear sign for the coupling of the Tyr42 and Glu46 protons, as experimentally demonstrated by Sigala et al.23 These findings show that the coupling between the chemical shift of the two protons is also influenced by the three distinct conformations, which are accessible at room temperature and under solution conditions. Sigala et al. concluded23 after studying the NMR shifts in mutants, and after isotopic substitution, that Glu46−H and Tyr42−H are coupled hydrogen bonds, such that any perturbation in one hydrogen bond affects the other. Their experimental findings and our computed data can now be combined to understand the behavior of PYP’s active site in the electronic ground state. QM-Size Convergence of the QM/MM NMR Data. The benchmark QM/MM NMR calculations of the PYP system confirm our earlier experience with various large biological systems,42−45 as well as those of other groups,46,47 where it was shown that the computed QM/MM properties are highly sensitive toward the size of the QM regions. Whereas previous computational studies of PYP14,16,30 used relatively small QM regions (fewer than 200 QM atoms14,16,30) to calculate the chemical shifts of key protons in PYP, our calculations indicate that more than 400 atoms are needed to reach converged (differences below 0.1 ppm) 1H NMR shifts with respect to the QM system size. In our benchmark study, we utilized all three QM/MM optimized structures (presented in Table 1) and

Table 2. Investigation of the Coupling of E46-H (Glu46-H) and Y42-H (Tyr42-H) by the Successive Elimination of Hydrogen Bond Forming Residues around pCAa T QM selection pCA+E46+Y42+T50 (QM/MM) pCA+E46+Y42+T50 pCA+E46+Y42 pCA+E46 pCA + Y42 pCA+E46 + T50 pCA + Y42+T50 a

F

E46-H Y42-H

G

E46-H Y42-H

E46-H Y42-H

13.9

13.6

14.3

11.0

15.2

11.5

13.8 13.8 14.4

14.1 14.0

14.2 14.4 14.9

11.6 12.0

15.1 14.9 15.3

12.1 12.0

14.3 14.3

12.5 14.7

14.3

12.4 15.5

12.0

12.4

The NMR data (in ppm) was obtained at the B97-2/pcS-1 level.

ppm in the chemical shift of the Glu46 proton is observed, whereas the downfield change in the Tyr42 proton is much larger with 0.42 ppm. In the G conformation, where the Thr50 is pointing away from the chromophore, very small changes in both chemical shifts upon removing the Thr50 are observed. This is expected, as the interaction with the third hydrogen bond (with Thr50 as the hydrogen donor) is no longer present in this conformation. In the T conformation, the removal of Thr50 from the R1 set of QM atoms has a small downfield influence of 0.17 ppm on the Tyr42 proton, while the effect on the chemical shift of the Glu46 proton, no longer interacting with Thr50 in the T conformation, is negligible with 0.03 ppm. After reducing the system down to pCA, Tyr42, and Glu46, the

Figure 4. Convergence of the B97-2/AMBER NMR chemical shifts of the Glu46 (H46: red) and Tyr42 (H42: blue) protons (relative to TMS in ppm; pcS-1 basis set) with respect to the QM regions R1−R7 in the T, F, and G conformations optimized at the B3LYP-D3/AMBER level of theory (QML region, def2-TZVP basis set). Schematic representation of the hydrogen bonding networks and the corresponding B3LYP-D3/AMBER relative energies are also shown for comparison. E

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation computed the chemical shifts of the discussed protons by performing single point B97-2/AMBER NMR calculations with the QM regions R1−R7 (70−517 atoms). The corresponding QM/MM NMR data for the conformations optimized using the QML and def2-TZVP basis set are shown in Figure 4. The NMR shifts of the Try42 and Glu46 protons move continuously downfield in all conformations and converge to within less than 0.1 ppm after including 425 atoms in the QM region. Interestingly, the computed proton shifts across the R1−R7 QM regions show downfield variations in the 0.5−0.7 ppm range, which highlights the sensitivity of the obtained values with respect to the size of system treated quantum chemically. Conversely, the peak difference for each conformation remains almost constant for all QM regions and varies by at most 0.1 ppm. Consequently, the QM-size converged peak difference reflects a similar trend already observed with the smaller QM region, QMS (see Figure 3); for the F and G conformers, converged peak differences are 3.6 and 3.2 ppm, respectively, implying a large separation between the two proton shifts, while for the T conformer, the corresponding value is 0.4 ppm, an indication of a close overlap of the two chemical shifts. The results of the convergence study can be combined with the energetic calculations, thereby further improving the quality of the calculated NMR peaks. The Boltzmann weighted 1H NMR shifts across all three conformations were calculated using the NMR shifts obtained from the largest QM region (R7), and the relative energies of each conformation. The obtained average peak differences are sensitive with respect to the choice of the QM region and basis set employed for the structural optimization. This is readily apparent in the case where the small QM region, QMS, and basis set def2-SVP are chosen. Here, the Boltzmann-averaged peak difference is 0.2 ppm. Furthermore, depending on the QM region and basis set, the Tyr42 proton shifts can vary by up to 1.7 ppm while the corresponding change in the Glu46 proton chemical shift is at most 0.5 ppm. This difference in the sensitivity of the chemical shifts can be traced back to the disproportionally larger influence of the conformational exchange on the Tyr42 proton. However, if the largest setup is employed (def2-TZVP basis and QML region), a Boltzmann-averaged peak difference of 0.9 ppm was found, which is in excellent agreement with the experimentally observed value of 1.5 ppm (see Table 3).

Table 4. RI-CC2/AMBER Excitation Energies in eV of PYP in the T, F, and G Conformers Calculated with the R1 and R2 (in parentheses) QM Regions, and the aug-cc-pVDZ Basis Set

Glu46-H [ppm]

Tyr42-H [ppm]

Δδ [ppm]

QMS/def2-SVP QMS/def2-TZVP QML/def2-SVP QML/def2-TZVP exptl.23

14.9 14.8 15.3 14.6 15.2

14.7 14.2 14.7 13.7 13.7

0.2 0.6 0.6 0.9 1.5

S0 → S1 [eV]a

oscillator strengthb

T F G

2.98 (2.98) 2.98 (2.99) 2.92 (2.92)

1.39 1.38 1.42

a The experimental UV/vis absorption maximum is 2.78 eV.9,11 bThe oscillator strength in length representation corresponds to R1.

(see Computational Details for further information). The π → π* local transition within the pCA chromophore was found to correspond to the lowest excited state at the RI-CC2 level, along with a high oscillator strength. The nature of the excited state was verified by inspection of the participating molecular orbitals. The vertical QM/MM excitations computed for both R1 and R2 QM regions (the latter also containing the Thr50) show very small variations across the conformers. The excitation energies of the T, F, and G conformers are within 0.2 eV of the experimental value.9,11 Considering the estimated 0.3 eV intrinsic error of the CC2 method itself,5 the computed values can be considered to be in agreement with the experiment. The occupied and virtual orbitals involved in the vertical excitation are, as expected, the phenolate oxygen and the carbon−carbon double bond just before the thioester linkage between pCA and Cys69. The geometric perturbations caused by the conformational exchange do not appear to differently affect the electronic transition, which in turn leads to almost identical vertical excitation energies for the T, F, and G conformers. Future investigations may elucidate the direct role the conformational exchange plays in the photocycle of PYP, such as the possible stabilization of intermediary states.



PYP MUTANTS To investigate whether the presented conformational exchange participates in other variants of PYP, we investigated three mutants, T50V, E46Q, and Y42F. Sigala et al. used these mutants as part of their investigation into the coupling of Glu46−H and Tyr42−H,23 which in addition to spectroscopic studies48,49 provides ample reference values to validate the presented mechanism. As mentioned in the computational details, the mutants were generated by direct mutation of the listed residues from the fully optimized T, F, and G conformers to ensure a consistent description of all structures. Geometric descriptions are based on the largest optimizations (B3LYP-D with QML and def2-TZVP for T50V and E46Q, QMXL and def2-TZVP for Y42F). The Boltzmann weighted NMR shifts based on the final QM size converged NMR shifts (R7) are reported for all mutants. T50V Mutant. As predicted, the active site of the T50V mutant exhibits only a single stable conformation due to the absence of Thr50. The fully optimized structure shows that both the pCA-Glu46 and pCA-Tyr42 hydrogen bonds are somewhat shortened compared to the WT to 2.59 and 2.63 Å, respectively. The Glu46−H and Tyr42−H chemical shifts are 15.6 and 12.2 ppm, with a peak difference of 3.4 ppm, which are in excellent agreement with the experimental values of 15.8 and 12.4 ppm reported by Sigala et al.23 The calculated vertical excitation energy based on the R2 region is 2.94 eV compared to the experimental UV/vis absorption maximum of 2.71 eV.49

Table 3. Weighted Average NMR Shifts of the R7 QM Region of the T, F, and G Conformers Using the Coefficients of a Boltzmann Distribution at T = 300 K optimization setup

structure

Influence of the Conformational Exchange on the UV/ vis Spectrum. Considering PYP’s role as a bacterial photosensor, we investigated the extent of the influence the conformational exchange has on the modulation of the UV/vis absorption spectrum. For this purpose, we performed single point vertical excitation energy calculations at the RI-CC2/ AMBER level of theory (aug-cc-pVDZ basis set) for all three conformations. The computed data are presented in Table 4 F

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation Considering the 0.3 eV intrinsic error of RI-CC2,5 this result can be taken to be accurate. E46Q Mutant. We predicted E46Q to possess three stable conformers, corresponding to the rotation of Thr50 to face different acceptor groups. The conformer labels of the WT are retained. The calculations confirm the predictions. As shown in Table 5, the E46Q mutant is able to cycle through thermodynamically accessible conformers, which contribute to the observed NMR spectrum of these mutants.

Table 7. Relative B3LYP-D3/AMBER Energies in kcal mol −1 of the F, Fp, Gp, and G Conformers of the Y42F Mutant Obtained with the QMXL Region and the def2-TZVP Basis Set Y42F

Table 5. Relative B3LYP-D3/AMBER Energies in kcal mol of the T, F, and G Conformers of the E46Q Mutant Obtained with the QML Regions and the def2-TZVP Basis Set T

F

G

0

2.3

2.2

The structural effects of the mutation on the hydrogen bond network are significant. The glutamine group is a much weaker hydrogen donor than the original Glu46, which causes the pCA-Tyr42 hydrogen bond to contract to 2.54, 2.62, and 2.60 Å in the T, F, and G conformers, respectively, to compensate for the absence of the stabilizing Glu46−H. Likewise, the pCAThr50 hydrogen bond is slightly shortened in the F conformation, resulting in an OpCA−OThr50 distance of 2.91 Å. The Gln46−H shows a Boltzmann-averaged chemical shift of 10.5 ppm, while the corresponding Tyr42−H chemical shift is 14.8 ppm, with a peak difference of 4.3 ppm. Sigala et al. reported an experimental value of 15.2 ppm for the Tyr42−H proton in this mutant.23 The calculated vertical excitation energies are shown in Table 6. As with the WT, the calculations show marked deviations Table 6. RI-CC2/AMBER Excitation Energies in eV of the E46Q Mutant in the T, F, and G Conformers Calculated with the R1 QM Regions, and the aug-cc-pVDZ Basis Set mutant

T

F

G

ref.

E46Q

2.93

2.94

2.89

2.6948

Fp

Gp

G

0

2.6

2.3

1.6

including the Gly29, Ala30, Ala45, Gly47, Ile49, Phe96, Val120, and Val122 residues in the QM region, an artificial instability causing a deprotonated Glu46 to move toward pCA, and recapture its proton is eliminated. This QM region is denoted QMXL (270 atoms). The Glu46−H hydrogen bond is strengthened in this mutant over the WT, shrinking to 2.56 Å in the F, and to 2.52 Å in the G conformers. The weaker hydrogen bond between Thr50 and pCA is likewise strengthened, resulting in a final OpCA−OThr50 distance of 2.87 Å. The Glu46−H bond shortens to 2.45 Å in the Fp conformer, due to a partial deprotonation of Glu46 (OGlu46−H distance of 1.32 Å) and a partial protonation of pCA (H−OpCA distance of 1.13 Å). The hydrogen bond between Thr50 and pCA is retained, although it lengthens to 2.97 Å. The Glu46−H bond lengthens to 2.51 Å in the Gp conformation, with pCA remaining protonated, and Glu46 deprotonated. Thr50 faces the backbone of Glu46. As the Glu46−H bond is significantly shorter in the Y42F mutant over the WT, the corresponding NMR shifts tend to be further downfield. In the F conformer, the final QM size converged NMR shifts of Glu46−H is 17.0 ppm. As the proton migrates to pCA in the Fp conformer, the NMR shifts move far downfield to 20.9 ppm and return upfield to 18.2 ppm in the Gp conformer. The NMR shifts of the G conformer are very similar to the Gp conformer with 18.3 ppm, indicating that the proton experiences a similar chemical environment in both, and that pCA and Glu46 are equivalent donor/acceptor groups. The Boltzmann weighted NMR shift of the Glu46−H proton in this mutant is 17.1 ppm, which agrees well with the experimental reference value of 16.7 ppm.23 In addition to a slightly red-shifted absorption maximum (458 nm) relative to the WT, the Y42F mutant shows a spectral shoulder at 391 nm in its UV/vis absorption spectrum.49,50 The experimental NMR spectrum likewise displays an anomaly, in that the assigned Glu46−H peak centered at 16.7 ppm is very broad.23 These findings suggest that the spectral shoulder is due to the presence of an additional active site conformation, commonly denoted as the intermediate spectral form,51 in which pCA appears protonated.49,52,53 The calculated vertical excitation energies of the chromophore in the F, Fg, Gp, and G conformers are shown in Table 8 The data show two distinct absorption maxima, with the conformers corresponding to the deprotonated chromophore being responsible for the excitation maximum, and the other

−1

E46Q

F

from the experimental references. However, they fall well within the experimentally observed bands and display the expected minimal redshift with respect to the WT of 0.03−0.05 eV relative to the experimentally observed redshift of 0.09 eV.9,11,48 The calculated values can therefore be understood to be reliable and describe the same electronic transition as the WT. The Y42F Mutant. Lacking Tyr42, and hence one possible hydrogen bond acceptor for Thr50, only two active site conformations due to the rotation of Thr50 are predicted (F and G) to exist. In light of spectroscopic studies,49,52−54 differing protonation states of pCA were taken into consideration. Two additional structures were found. In the first, Thr50 forms a weak hydrogen bond with pCA, while the Glu46 proton is shifted considerably toward pCA. This structure is denoted as “Fp.” In the second, Thr50 forms a hydrogen bond with the backbone of Glu46, while the Glu46 proton is shifted even closer to pCA. This structure is denoted as “Gp.” The relative energies of these four thermodynamically accessible conformers are listed in Table 7. The QM region had to be enlarged to include the residues surrounding Glu46 in order to properly describe the system. By

Table 8. RI-CC2/AMBER Excitation Energies in eV of the Y42F Mutant in the F, Fp, Gp, and G Conformers Calculated with the R2 QM Regions (Excluding Phe42 and Ile31), and the aug-cc-pVDZ Basis Set F, Fp

G, Gp

ref.a

2.92, 3.15

2.86, 3.17

2.71, 3.1749

a

The experimental reference lists the main absorption peak (2.71 eV), and the spectral shoulder (3.17 eV).

G

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation conformers with the protonated chromophore being responsible for the spectral shoulder, with the small energetic differences between all four suggesting an equilibrium. The stability of the Fp and Gp conformers not only confirms the experimental study by Joshi et al.,52 which found the intermediate spectral form to possess a protonated chromophore, and a deprotonated Glu46, but also expands upon the interpretation, as the rotation of Thr50 results in two additional contributing conformers.





CONCLUSION AND OUTLOOK Using linear-scaling QM/MM techniques,55,56 we computed the 1H NMR shifts for the Tyr42 and Glu46 side chains within the hydrogen bond network around the pCA chromophore of PYP. Our investigation revealed that the hydrogen bond network undergoes a previously unknown fast conformational exchange, which is triggered by the rotation of a third nearby amino acid, Thr50. This exchange can induce thermally competitive conformations, denoted in this work as T, F, and G, which influence the electronic structure, and consequently the 1H NMR shifts of the Tyr42 and Glu46 side chain protons. Taking this conformational exchange into account allows for a proper description of PYP’s active site, which in turn allows for a successful explanation of the experimentally observed NMR spectrum. Via systematic combination of the QM/MM conformational energetics and the QM/MM NMR convergence data, a peak difference of 0.9 ppm between the Glu46 and Tyr42 protons can be assigned, which is in good agreement with the experimental value of 1.5 ppm. Our data support, and fully explain, the mechanism of the coupling between the Glu46−H and Tyr42−H hydrogen bonds proposed by the experimental work of Sigala et al.23 by quantifying to what extent perturbations in the chemical environment around one hydrogen bond affect the other. We further applied our proposed mechanism to mutants of PYP. We found the conformational exchange involving the rotation of Thr50 to play an essential role regarding the NMR properties of the E46Q and Y42F mutants, while the T50V mutant, which lacks Thr50, serves as an illustrative example of the accuracy of the employed quantum chemical methods with regard to large systems. Furthermore, the calculated vertical excitation energies of the wild type and these mutants allow for reasonable comparison with experiments, and also explain the abnormal spectral shoulder in the UV/vis absorption spectrum of Y42F.49,50 We also performed a QM/MM NMR benchmark study to ensure convergence of the NMR data with respect to the size of the QM region. Our data show that more than 400 atoms of PYP’s active site need to be included in the QM region of the calculations in order to obtain 1H shifts which are consistent to within 0.1 ppm.



mutants, FF MD simulations on the conformers, investigation of possible BSSE effects (PDF) largest QM/MM optimizations (QML/def2-TZVP) for the T conformer (PDB) largest QM/MM optimizations (QML/def2-TZVP) for the F conformer (PDB) largest QM/MM optimizations (QML/def2-TZVP) for the G conformer (PDB)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We acknowledge financial support by the DFG funding initiatives SFB749 (C7), the Excellence Cluster EXC114 (CIPSM), and by the Volkswagen Stiftung within the funding initiative “New Conceptual Approaches to Modeling and Simulation of Complex Systems.”



REFERENCES

(1) Sprenger, W. W.; Hoff, W. D.; Armitage, J. P.; Hellingwerf, K. J. J. Bacteriol. 1993, 175, 3096−3104. (2) Kamiya, M.; Saito, S.; Ohmine, I. J. Phys. Chem. B 2007, 111, 2948−2956. (3) Hellingwerf, K. J.; Hendriks, J.; Gensch, T. J. Phys. Chem. A 2003, 107, 1082−1094. (4) Borgstahl, G. E. O.; Williams, D. R.; Getzoff, E. D. Biochemistry 1995, 34, 6278−6287. (5) Gromov, E. V.; Burghardt, I.; Koeppel, H.; Cederbaum, L. S. J. Am. Chem. Soc. 2007, 129, 6798−6806. (6) Groenhof, G.; Schaefer, L. V.; Boggio-Pasqua, M.; Grubmueller, H.; Robb, M. A. J. Am. Chem. Soc. 2008, 130, 3250−3251. (7) Isborn, C. M.; Götz, A. W.; Clark, M. A.; Walker, R. C.; Martínez, T. J. J. Chem. Theory Comput. 2012, 8, 5092−5106. (8) Groenhof, G.; Bouxin-Cademartory, M.; Hess, B.; de Visser, S. P.; Berendsen, H. J. C.; Olivucci, M.; Mark, A. E.; Robb, M. A. J. Am. Chem. Soc. 2004, 126, 4228−4233. (9) Stahl, A. D.; Hospes, M.; Singhal, K.; van Stokkum, I.; van Grondelle, R.; Groot, M. L.; Hellingwerf, K. J. Biophys. J. 2011, 101, 1184−1192. (10) Changenet-Barret, P.; Plaza, P.; Martin, M. M.; Chosrowjan, H.; Taniguchi, S.; Mataga, N.; Imamoto, Y.; Kataoka, M. Chem. Phys. Lett. 2007, 434, 320−325. (11) Imamoto, Y.; Kataoka, M.; Tokunaga, F. Biochemistry 1996, 35, 14047−14053. (12) Groot, M. L.; van Wilderen, L. J. G. W.; Larsen, D. S.; van der Horst, M. A.; van Stokkum, I. H. M.; Hellingwerf, K. J.; van Grondelle, R. Biochemistry 2003, 42, 10054−10059. (13) Schotte, F.; Cho, H.; Kaila, V.; Kamikubo, H.; Dashdorj, N.; Henry, E.; Graber, T.; Henning, R.; Wulff, M.; Hummer, G.; Kataoka, M.; Anfinrud, P. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 19256− 19261. (14) Saito, K.; Ishikita, H. Biochemistry 2012, 51, 1171−1177. (15) Thellamurege, N. M.; Si, D.; Cui, F.; Li, H. J. Chem. Phys. 2014, 140, 174115. (16) Nadal-Ferret, M.; Gelabert, R.; Moreno, M.; Lluch, J. M. J. Am. Chem. Soc. 2014, 136, 3542−3552. (17) Hirano, K.; Sato, H. Chem. Phys. 2013, 419, 163−166. (18) Yamaguchi, S.; Kamikubo, H.; Kurihara, K.; Kuroki, R.; Niimura, N.; Shimizu, N.; Yamazaki, Y.; Kataoka, M. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 440−444. (19) Saito, K.; Ishikita, H. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 167−172.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.6b00359. Structural information regarding the chosen QM regions, ZPE corrections for the relative energies of the conformers, NMR and UV/vis vertical excitation energy calculations on the crystal structure of PYP, additional NMR calculations on the presented conformers and PYP H

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Article

Journal of Chemical Theory and Computation (20) Cleland, W. W.; Kreevoy, M. M. Science (Washington, DC, U. S.) 1994, 264, 1887−1890. (21) Frey, P. A.; Whitt, S. A.; Tobin, J. B. Science (Washington, DC, U. S.) 1994, 264, 1927−1930. (22) Frey, P. In Isotope Effects in the Characterization of Low Barrier Hydrogen Bonds; Kohen, A., Limbach, H.-H., Eds.; CRC Press: Boca Raton, FL, 2006; p 975. (23) Sigala, P. A.; Tsuchida, M. A.; Herschlag, D. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 9232−9237. (24) Pool, T. J.; Oktaviani, N. A.; Kamikubo, H.; Kataoka, M.; Mulder, F. A. A. Biomol. NMR Assignments 2013, 7, 97−100. (25) Mildvan, A.; Massiah, M.; Harris, T.; Marks, G.; Harrison, D.; Viragh, C.; Reddy, P.; Kovach, I. J. Mol. Struct. 2002, 615, 163−175. (26) Anderson, S.; Crosson, S.; Moffat, K. Acta Crystallogr., Sect. D: Biol. Crystallogr. 2004, 60, 1008−1016. (27) Wang, L.; Fried, S. D.; Boxer, S. G.; Markland, T. E. Proc. Natl. Acad. Sci. U. S. A. 2014, 111, 18454−18459. (28) Rossi, M.; Fang, W.; Michaelides, A. J. Phys. Chem. Lett. 2015, 6, 4233−4238. (29) Ceriotti, M.; Cuny, J.; Parrinello, M.; Manolopoulos, D. E. Proc. Natl. Acad. Sci. U. S. A. 2013, 110, 15591−15596. (30) Kanematsu, Y.; Tachikawa, M. J. Chem. Phys. 2014, 141, 185101. (31) Søndergaard, C. R.; Olsson, M. H. M.; Rostkowski, M.; Jensen, J. H. J. Chem. Theory Comput. 2011, 7, 2284−2295. (32) Case, D. A.; Darden, T. A.; Cheatham, T. E.; Simmerling, C. L.; Wang, J.; Duke, R. E.; Luo, R.; Crowley, M.; Walker, R. C.; Zhang, W.; Merz, K. M.; Wang, B.; Hayik, S.; Roitberg, A.; Seabra, G.; Kolossváry, I.; Wong, K. F.; Paesani, F.; Vanicek, J.; Wu, X.; Brozell, S. R.; Steinbrecher, T.; Gohlke, H.; Yang, L.; Tan, C.; Mongan, J.; Hornak, V.; Cui, G.; Mathews, D. H.; Seetin, M. G.; Sagui, C.; Babin, V.; Kollman, P. A. Amber 11; University of California: San Francisco. (33) Kussmann, J.; Ochsenfeld, C. J. Chem. Phys. 2013, 138, 134114. (34) Kussmann, J.; Ochsenfeld, C. J. Chem. Theory Comput. 2015, 11, 918−922. (35) Stephens, P.; Devlin, F.; Chabalowski, C.; Frisch, M. J. Phys. Chem. 1994, 98, 11623−11627. (36) Becke, A. J. Chem. Phys. 1993, 98, 1372. (37) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. J. Chem. Phys. 2010, 132, 154104. (38) Wilson, P. J.; Bradley, T. J.; Tozer, D. J. J. Chem. Phys. 2001, 115, 9233. (39) Flaig, D.; Maurer, M.; Hanni, M.; Braunger, K.; Kick, L.; Thubauville, M.; Ochsenfeld, C. J. Chem. Theory Comput. 2014, 10, 572−578. (40) Jensen, F. J. Chem. Theory Comput. 2008, 4, 719−727. (41) TURBOMOLE, V6.3; University of Karlsruhe and Forschungszentrum Karlsruhe GmbH: Karlsruhe, Germany, 1989−2007; TURBOMOLE GmbH: Karlsruhe, Germany, 2007. Available from http://www.turbomole.com. (42) Blank, I. D.; Sadeghian, K.; Ochsenfeld, C. Sci. Rep. 2015, 5, 10369. (43) Sadeghian, K.; Flaig, D.; Blank, I. D.; Schneider, S.; Strasser, R.; Stathis, D.; Winnacker, M.; Carell, T.; Ochsenfeld, C. Angew. Chem., Int. Ed. 2014, 53, 10044−8. (44) Flaig, D.; Beer, M.; Ochsenfeld, C. J. Chem. Theory Comput. 2012, 8, 2260−2271. (45) Sumowski, C. V.; Schmitt, B. B. T.; Schweizer, S.; Ochsenfeld, C. Angew. Chem., Int. Ed. 2010, 49, 9951−9955. (46) Liao, R.-Z.; Thiel, W. J. Chem. Theory Comput. 2012, 8, 3793− 3803. (47) Hu, L.; Söderhjelm, P.; Ryde, U. J. Chem. Theory Comput. 2011, 7, 761−777. (48) Philip, A. f.; Eisenman, K. T.; Papadantonakis, G. A.; Hoff, W. D. Biochemistry 2008, 47, 13800−13810. (49) Brudler, R.; Meyer, T. E.; Genick, U. K.; Devanathan, S.; Woo, T. T.; Millar, D. P.; Gerwert, K.; Cusanovich, M. A.; Tollin, G.; Getzoff, E. D. Biochemistry 2000, 39, 13478−13486. (50) Mihara, K.; Hisatomi, O.; Imamoto, Y.; Kataoka, M.; Tokunaga, F. J. Biochem. 1997, 121, 876−880.

(51) Mizuno, M.; Kamikubo, H.; Kataoka, M.; Mizutani, Y. J. Phys. Chem. B 2011, 115, 9306−9310. (52) Joshi, C. P.; Otto, H.; Hoersch, D.; Meyer, T. E.; Cusanovich, M. A.; Heyn, M. P. Biochemistry 2009, 48, 9980−9993. (53) El-Mashtoly, S. F.; Unno, M.; Kumauchi, M.; Hamada, N.; Fujiwara, K.; Sasaki, J.; Imamoto, Y.; Kataoka, M.; Tokunaga, F.; Yamauchi, S. Biochemistry 2004, 43, 2279−2287. (54) Meyer, T. E.; Devanathan, S.; Woo, T.; Getzoff, E. D.; Tollin, G.; Cusanovich, M. A. Biochemistry 2003, 42, 3319−3325. (55) Kussmann, J.; Beer, M.; Ochsenfeld, C. WIREs Comput. Mol. Sci. 2013, 3, 614−636. (56) Kussmann, J.; Ochsenfeld, C. J. Chem. Phys. 2007, 127, 054103− 054119.

I

DOI: 10.1021/acs.jctc.6b00359 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX