Quantum Calculations Indicate Effective Electron Transfer between

The ping-pong kinetic mechanisms of E. coli WrbA(8) and the Nqos,(9, 10) ..... of 10 that is predicted from the pKa of 6.2 for the N1 proton of free F...
2 downloads 0 Views 2MB Size
Subscriber access provided by ORTA DOGU TEKNIK UNIVERSITESI KUTUPHANESI

Article

Quantum calculations indicate effective electron transfer between FMN and benzoquinone in a new crystal structure of E. coli WrbA Oksana Degtjarik, Jiri Brynda, Olga Ettrichova, Michal Kuty, Dhiraj Sinha, Ivana Kuta Smatanova, Jannette Carey, Rudiger Ettrich, and David Reha J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.5b11958 • Publication Date (Web): 16 May 2016 Downloaded from http://pubs.acs.org on May 17, 2016

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 free 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 accessible to all readers and 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.

The Journal of Physical Chemistry B 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 27

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

The Journal of Physical Chemistry

Quantum Calculations Indicate Effective Electron Transfer between FMN and Benzoquinone in a New Crystal Structure of E. coli WrbA Oksana Degtjarik1,2,#, Jiři Brynda3, Olga Ettrichova1, Michal Kuty1,2, Dhiraj Sinha1,2, Ivana Kuta Smatanova1,2, Jannette Carey1,4, Rüdiger Ettrich1,2, David Řeha1,2,* 1Center

2Faculty

for Nanobiology and Structural Biology, Institute of Microbiology, Academy of Sciences of the Czech Republic, Zamek 136, CZ-373 33 Nove Hrady, Czech Republic

of Sciences, University of South Bohemia in Ceske Budejovice, Zamek 136, CZ-373 33 Nove Hrady, Czech Republic

3Institute

of Molecular Genetics and Institute of Organic Chemistry and Biochemistry, Academy of Sciences of the Czech Republic, Flemingovo nam. 2, 16610 Prague 6, Czech Republic

4Chemistry

#

Department, Princeton University, Princeton, New Jersey 08544-1009, USA

Current address: Department of Structural Biology, Weizmann Institute of Science, Rehovot, Israel

*

corresponding author: D. Řeha, [email protected]. Phone: +420389033801

Abstract Quantum mechanical calculations using the Marcus equation are applied to compare the electrontransfer probability for two distinct crystal structures of the Escherichia coli protein WrbA, an FMNdependent NAD(P)H:quinone oxidoreductase, with the bound substrate benzoquinone. The calculations indicate that the position of benzoquinone in a new structure reported here and solved at 1.33 Å resolution is more likely to be relevant for the physiological reaction of WrbA than a previously reported crystal structure in which benzoquinone is shifted by ~5 Å. Because the true electron-acceptor substrate for WrbA is not yet known, the present results can serve to constrain computational docking attempts with potential substrates that may aid in identifying the natural substrate(s) and physiological role(s) of this enzyme. The approach used here highlights a role for quantum mechanical calculations in the interpretation of protein crystal structures.

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 27

Introduction The FMN-dependent NAD(P)H:quinone oxidoreductase WrbA isolated from E. coli (Figure 1) is the founding member of an unusual family of flavodoxin-like proteins.1 WrbA-like proteins combine features found both in monomeric bacterial flavodoxins that use FMN as cofactor and in dimeric eukaryotic FAD-dependent oxidoreductases (Nqos).2 Like the flavodoxins, WrbAs adopt the canonical α/β twisted open-sheet fold and use FMN as physiological cofactor1, but like the Nqos they have a characteristic short sequence insertion after β2 that forms an additional beta-alpha secondary structural element lying outside the canonical fold.2-6 Like the Nqos, WrbAs populate a dimeric state in solution, but they also undergo a facile dynamic equilibrium to form active tetramers.2 Unexpectedly, the assembly of WrbA multimers relies not on their characteristic sequence insertion relative to the flavodoxins, but on secondary structural elements in common with the flavodoxin fold that in WrbA appear to be adapted for subunit interaction.5 The insertion instead is located at the two polar regions of the tetrahedral assembly, where it contributes to forming a hydrophobic channel leading from the protein surface to the active sites.2

Figure 1. Overall structure of WrbA-FMN-BQ. (a) Tetrameric assembly. Protein subunits are represented as cartoons named and colored by monomer (A, cyan; B, green; C, yellow; and D, pink). FMN and BQ are shown in ball-and-stick models with atomic colors (O, red; N, blue; P, orange) and black or orange carbons, respectively. (b) WrbA monomer. Monomer A from panel a zoomed in and rotated -95° about the horizontal axis and -35° about the vertical axis. Helices and strands are numbered as in Wolfova et al., 2009. Each cavernous active site is formed by residues from three WrbA subunits,2 indicating that ACS Paragon Plus Environment

Page 3 of 27

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

The Journal of Physical Chemistry

tetramers are the active species. The FMN isoalloxazine ring is stacked between aromatic residues similarly as in flavodoxins, although the topological origin of the active-site residues differs.5 Biochemical assays demonstrate reducing activity of WrbA toward short- and long-chain quinones7,8 including benzoquinone (BQ) that is used in a standard assay with NADH as electron donor. The pingpong kinetic mechanisms of E. coli WrbA8 and the Nqos9,10, together with computational docking experiments8, indicate that the electron donor substrate NADH and the electron acceptor substrate BQ do not occupy the active site simultaneously. The hydrophobic channel connecting the active site to the protein surface has been proposed as a possible conduit for long-chain, membrane-bound quinones.2 However, the natural substrate(s) and physiological role(s) of WrbA are unknown, and its ubiquity among bacteria and plants remains unexplained. One possible route to understanding the natural role of WrbA is to evaluate the suitability of potential physiological electron-acceptor substrates by docking them in the active site. A previous crystal structure of WrbA-FMN (i.e., holoWrbA) with BQ soaked in (PDB:3b6k)4 showed BQ in a location that appears to be too distant for efficient electron transfer from the isoalloxazine ring of FMN. If true, then that structure might not present a suitable active-site environment as the starting point for docking studies. This is because in docking applications proteins are generally treated as static, but real proteins are not static and structures can change between the free and ligand-bound states. For this reason a ligand-bound structure is the preferred starting point for docking studies (after ligand removal), provided the ligand is bound properly in the active site. New crystals of holoWrbA that diffract to 1.33 Å obtained under novel conditions are described in the present report, and structures are solved in the presence and absence of soaked-in BQ. The results reveal a position of BQ ~ 5 Å closer to the central ring of FMN isoalloxazine in the new crystal structure than in the earlier structure. To evaluate the suitability of the present and previous crystal structures as potential docking targets, the probability of electron transfer from FMN to a BQ molecule in the location found in each crystal structure was calculated in the present work using quantum mechanical (QM) and quantum mechanical/molecular mechanical (QM/MM) calculations to estimate the rate of electron transfer predicted for each structure by the Marcus equation.11 The approach used for the present calculations is borrowed from earlier work12 in which the effect of thermal fluctuations on the conductivity of DNA was studied. In that work pi-stacked DNA bases were treated as conductors, with charge transferred from one end of the DNA to the other along the helix. DNA bases were treated as isolated molecules for the purpose of the QM calculations, and a hopping mechanism was considered with electrons or electron holes transferred from each base to a consecutive base. The individual rates of charge-transfer steps between consecutive bases were ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 27

calculated using the Marcus equation, and the effect of thermal fluctuations was incorporated by calculations on a set of geometries generated by MD simulation. The approximations introduced in that work are here adapted to treat BQ reduction, where stacking between the flavin ring and BQ may be considered analogous to the pi-stacked DNA bases. Electron transfer from the isoalloxazine ring to BQ is considered as one charge-transfer step analogous to one charge-transfer step between consecutive DNA bases. The WrbA reaction, FMNH2 oxidation to FMN and BQ reduction to hydroquinone (BQH2), involves proton transfer as well as electron transfer. Both transfers must be taken into account because it is uncertain at present which one is the rate-limiting step. Furthermore, the protonation state of reduced FMN is unknown under conditions of the WrbA enzyme assay, and may not be constant over the wide range of pH where WrbA is active.7 Therefore both protonated states of reduced FMN (FMNH- and FMNH2) were taken into account in this study, with an emphasis on FMNH2 as the generally assumed form of reduced FMN. To treat proton transfer from FMNH2 or FMNH- to BQ a method based on approximation of the most probable reaction coordinate is adopted here. To our knowledge the application of these computational methods to this kind of problem is novel. Although approximations are required in the computations to treat proton transfer, these limitations are the same for the two crystal structures. Thus a comparison of the predicted electron-transfer probability for the two crystallographic locations of BQ appears justified. The results predict a higher rate of electron transfer when BQ is in the position found in the new crystal structure. The electron-transfer probabilities for the two crystallographic locations differ more than could reasonably be attributed to the approximations required in the calculations. Thus the active site in the new structure should be a useful starting point for computational docking of potential substrates that can aid the search for the natural electron acceptor(s) and physiological role(s) of this enzyme. The computational methods introduced here may find application in other charge-transfer processes.

Materials and Methods

Protein preparation and crystallization WrbA protein was overexpressed and purified as described previously.6 Protein was dialyzed against 20 mM Tris-HCl, pH 7.5 at 4°C and concentrated using Amicon Ultra centrifugal filter units (3000

ACS Paragon Plus Environment

Page 5 of 27

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

The Journal of Physical Chemistry

nominal molecular weight limit; Merck Millipore, USA) to 8 mg ml-1 (~ 0.35 mM monomer). Prior to crystallization protein was mixed with 5mM FMN in 20 mM Tris-HCl, pH 7.5 to final FMN concentration 0.35 mM. Screening for new crystallization conditions was performed in MRC 96-well crystallization plates (Molecular Dimensions Limited, UK) by the sitting-drop vapor-diffusion method using Gryphon modular dispenser (Art Robbins Instruments, USA) with various commercial screens. Well-formed crystals of deep yellow color appeared within 3 days at room temperature in condition D4 of Morpheus screen (Molecular Dimensions Limited, UK) consisting of 0.1 M MES/imidazole buffer, pH 6.5, 12.5% (v/v) MPD, 12.5% (w/v) PEG 1000, 12.5% (w/v) PEG 3350, 0.02 M 1,6-hexanediol, 0.02 M 1-butanol, 0.02 M 1,2-propanediol, 0.02 M 2-propanol, 0.02 M 1,4-butanediol, 0.02 M 1,3-propanediol. The crystals reached their final dimensions in approximately 10 days.

Soaking and data collection

Due to very low solubility in water, dark brown crystals of 1,4-benzoquinone (BQ) were dissolved in 100% (v/v) ethanol, then diluted with water to reach 15% ethanol (v/v) prior to final dilution to 25 mM with 20 mM Tris-HCl, pH 7.5. Prior to data collection 1 µl of 25 mM BQ was added to a 4 µl drop with WrbA-FMN crystals and incubated for 10 min. The final concentration of BQ in the drop was 5 mM. Diffraction data were collected on beamline BL14.1 operated by the Joint Berlin MX Laboratory at the BESSY II electron-storage ring (Berlin-Adlershof, Germany)13 equipped with a PILATUS 6M detector (Dectris Ltd., Switzerland). Crystals were mounted in a LithoLoop (Molecular Dimensions Ltd., UK) and flash-cooled in liquid nitrogen without additional cryoprotection. Images (1200 and 1500 from WrbA-FMN and WrbA-FMN-BQ crystals, respectively) were collected at a wavelength of 0.91841 Å with an oscillation range of 0.1° and 0.5 s exposure time. Diffraction data were processed using the graphical user interface XDSAPP version 1.014,15. Solvent content was analyzed using MATTHEWS_COEF from the CCP4 package16.

Crystal structure determination and refinement The structure of WrbA in complex with FMN was solved by molecular replacement using the Balbes automated molecular replacement pipeline17 incorporated into the CCP4 suite16. The structure of WrbA (PDB:3b6i) was used as search model. The structure of WrbA with bound FMN and BQ was solved using the MOLREP program18 from the CCP4 suite with the coordinates of WrbA apoprotein ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 27

(PDB:2rg1)5 as search model. The positions of all atoms in both WrbA-FMN and WrbA-FMN-BQ models were refined in REFMAC19 applying standard geometrical restraints and NCS with anisotropic refinement of atomic displacement parameters for all protein, ligand, and water atoms. Manual model building and correction were performed in Coot20. The quality of the model was analyzed with the MolProbity server.21 Figures were prepared in PyMOL22 if not otherwise specified. Coordinates and structure factors for the crystal structures of WrbA-FMN and WrbA-FMN-BQ have been deposited in the RCSB Protein Data Bank with access codes 5F12 and 4YQE, respectively. Geometry optimization The geometry of the binding site containing both FMN and BQ was optimized using QM/MM calculations implemented in Qsite from the Schrodinger software package23,24. Binding site A from each crystal structure (the structure reported herein and the previously reported structure PDB:3b6k) was chosen for these calculations. For the QM region BQ, the isoalloxazine ring of FMN in its reduced form (FMNH2), and the side chain of Trp97 were selected. Selection of this QM region is based on the fact that these groups form a pi-stacked complex that dominates the positioning of BQ. All other parts of the system were treated as the MM region. Boundaries between the QM and MM regions were established by capping with hydrogen atoms, and a Gaussian grid was used to treat electrostatics. Density functional theory (DFT) with the M06-2X functional25 and D3 corrections of dispersion energy26 using the 6-31G* basis set were applied for the QM region. This method is well suited for calculations on pi-stacked structures27,28 and represents a compromise between computational accuracy and computational demand. To describe the MM region the OPLS2005 force field29 (including nonstandard residues like FMN tail, which is part of MM region) was used with implicit solvent and a 15 Å cutoff for non-bonding interactions. The whole QM region as well as all protein residues within 7Å of BQ was fully flexible, while the rest of the protein plus FMN and BQ in the other binding sites were restrained during geometry optimization to maintain the overall conformation of the protein.

Charge-transfer rate calculation

The probability of charge transfer (CT) between FMN and BQ was modeled by calculating the chargetransfer rate based on the Marcus equation: DA = ℏ

  TE

exp−ΔG+E  ⁄4  TE 

(1)

where kDA is CT rate, kB is Boltzman constant, h is Planck's constant, ∆G is the driving force for charge migration, Eλ is the reorganization energy, J is the electronic coupling constant, and T is temperature. ACS Paragon Plus Environment

Page 7 of 27

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

The Journal of Physical Chemistry

The energy Eλ represents the change of free energy necessary for the donor/acceptor pair to acquire a geometry conducive to charge transfer, as required by the Franck–Condon principle. Eλ also includes a contribution from reorganization of the surrounding environment including protein and solvent. The driving force for charge migration depends on the difference between the free energies of electron donor and acceptor. The electronic coupling constant J represents the orbital overlap between the initial and final state of the electron transfer reaction defined as: J=

(2)

where ΨA and ΨB are electronic wave functions of the system before (A) and after (B) charge transfer and Ĥ is the Hamiltonian of the system. Electron transport corresponds to electron transfer from the lowest unoccupied molecular orbital (LUMO) of system A to the LUMO of system B, where A is electron donor and B is electron acceptor. In case of hole transport12 the electron is transferred from the highest occupied molecular orbital (HOMO) of system B to the electron hole in the HOMO of system A, where A is the hole donor and B is the hole acceptor. The multi-electron wave function is described by a Slater determinant regardless of the position of the electron, and therefore the Slater determinants before and after charge transfer differ only in the occupancy of the respective LUMO and HOMO orbitals on A and B. Therefore using Slater rules: JAB=

(3)

where ψA , ψB are LUMO/HOMO orbitals of molecule A and B, respectively, and F is the approximate single-electron Fock matrix of the system12. The driving force (∆G) is estimated as the difference between site energies of molecules A and B (acceptor B minus donor A) associated with the HOMO or LUMO orbitals involved. The electronic coupling J was calculated in Gaussian 09 software package30 using the density functional theory (DFT) method with a wB97XD functional31 and cc-pVTZ basis set. The chosen functional was reported to give the best results for charge-transfer calculation32. First, the Fock matrix of the system (donor/acceptor pair) was calculated. Then the molecular orbitals of separate donor and acceptor molecules were calculated and the Fock matrix of the system was expressed by the molecular orbitals of separate donor and acceptor molecules. The coupling constant J then directly corresponds to the element of the Fock matrix (written in this form) defined by the involved LUMO or HOMO of the donor and acceptor.

Estimation of protein reorganization energy To estimate the reorganization energy of the protein environment surrounding the cofactor FMN and the ligand BQ that accompanies the redox reaction, the QM/MM-optimized geometry of each structure ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 8 of 27

in the initial redox state (reduced FMN, i.e., FMNH2, and oxidized BQ) was re-optimized for the final redox state (oxidized FMN and reduced BQ, i.e., BQH2) using the same method. This procedure generates the protein conformation for the situation after the redox reaction. In the next step, RMSD least-squares fit was used to overlay the cofactor and ligand in the initial redox state (FMNH2 and BQ) with their positions in the optimized final redox state (FMN and BQH2). Finally FMN and BQH2 were exchanged with the overlaid FMNH2 and BQ to create the situation representing the initial redox state surrounded by the protein conformation optimized in the final redox state (FMN and BQH2). Then single-point calculations were performed using QM/MM (Qsite, Schrodinger)23,24 for the protein conformation optimized with FMNH2-BQ ligands and the protein conformation optimized with FMNBQH2 ligands but containing FMNH2-BQ instead. In these calculations FMNH2-BQ was treated as the QM region and all other parts as the MM region. DFT with M06-2X functional25 and the 6-31G* basis set was used as the QM method and the OPLS2005 force field29 with a nonbonding cut-off of 15 Å and implicit solvent was applied for the MM region. The total QM/MM energy is the sum of the QM term, the MM term, and the QM/MM energy coupling. Because the geometries and redox states of FMNH2 and BQ are identical for the protein conformation optimized with FMNH2-BQ ligands and the protein conformation optimized with FMNBQH2 ligands but containing FMNH2-BQ instead, their respective QM energy terms must be the same as well. Thus the difference of QM/MM energies between the geometries originates from the MM energy term and/or the QM/MM energy coupling, and the difference of QM/MM energies is an estimate of the reorganization energy of the surrounding protein due to the redox reaction. In this way the reorganization energy was calculated for each crystal structure to estimate the difference in reorganization energies between them.

Estimation of proton-transfer probability Various protonation states of free FMN have been investigated using QM calculations.33 The results indicate that deprotonation of reduced FMN occurs preferentially at N1. The pKa of the proton at N1 of FMNH2 is 6.2.34 The activity of WrbA is highest in the pH range 5-8, and with a peak at pH ~6.5).7 Assuming a similar pKa value for FMNH2 bound to WrbA as for free FMN, both protonated forms of reduced FMN (FMNH2 and the anionic form FMNH- with protonated N5 atom) are likely to be physiologically relevant and must be taken into account when describing the enzymatic redox reaction. Spectroelectrochemical measurements with reduced FMN using a flavin-modified gold electrode35 indicate that reduced FMN in WrbA at pH 7.2 is mainly in the form of the monoprotonated anion FMNH-. This result at pH 7.2 is consistent with the ratio of [FMNH-] to [FMNH2] of 10 that is ACS Paragon Plus Environment

Page 9 of 27

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

The Journal of Physical Chemistry

predicted from the pKa of 6.2 for the N1 proton of free FMNH2. However, most kinetic studies with WrbA have been performed at pH 6.56,8, where the predicted ratio is only ~2. Furthermore, fully protonated FMNH2 is generally assumed in structural reports. For example, Driggers et al36 modeled the crystal structure of the WrbA homolog E.coli SsuE at pH 7.0 with FMNH2 based on absorption spectra. Although it is not possible to estimate the local pH in the protein, the cited results and the wide pH-activity range indicate that WrbA both FMNH2 must be taken into account in the calculations performed in the present work. The energy barrier of two-proton transfer from FMNH2 (or FMNH-) to benzoquinone for each crystal structure was estimated by semi-empirical calculations (PM7) implemented in MOPAC37,38 using the TRITON interface39,40 that enables calculations of energies along the reaction coordinate. This semi-empirical method is used to provide a qualitative estimate of the difference in probability of proton transfer for the two crystal structures. The PM7 method correctly describes pi-stacking interactions38. The system selected for the semi-empirical calculations consisted of BQ, FMNH2 or FMNH- with the phosphate group omitted to remove contributions from its charge, and protein residues within 5 Å of the bound ligand or cofactor. The system was capped using hydrogen atoms where covalent bonds were interrupted during the selection. To maintain the residue positions of the intact protein with most of the protein missing, position restraints were applied at each new N- or C-terminus of the protein created by the selection procedure, as well as at the C atom of the FMN ribosyl tail at the site of phosphate truncation. Two scenarios of proton transfer were considered, one in which both hydrogen atoms at N1 and N5 of FMNH2 are transferred simultaneously to two different oxygen atoms of BQ, and another in which the hydrogen atom at N5 of FMNH- is transferred to the closer oxygen of BQ. In the case of two-proton transfer the distance vector from the hydrogen atom at N1 of FMNH2 to one oxygen atom of BQ and the distance vector from the other oxygen of BQ to the hydrogen atom at N5 of FMNH2 were defined as the reaction coordinate. In structure 3b6k the distance from atom N5 of FMNH2 to the closer oxygen of BQ is 4.6 Å and the distance from atom N1 to the other oxygen of BQ is 8.8 Å, whereas in the present crystal structure the distances are similar, 3.9 Å to N5 and 3.5 Å to N1. As it is not known which proton is transferred first, both distance vectors are decreased simultaneously and proportionately to reach the target value of 0.85 Å. The vector with the shorter initial distance was decreased by 0.05 Å per step and is taken as the reaction coordinate. In the case of one-proton transfer the distance vector from the hydrogen atom at N5 of FMNH- to the closer oxygen of BQ was considered as the reaction coordinate, and was also decreased to reach the target value of 0.85 Å. At each distance the program optimizes the geometry by energy minimization using the PM7 method, producing one point of the energy profile. All iteratively calculated points taken together constitute the ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 10 of 27

energy profile along the reaction coordinate. The energies calculated in this way are strongly dependent on the environment and do not represent absolute values, but can be used only for comparison between the two structures. The product state geometries are almost identical in both cases. Thus energy of the product was assigned a value of zero. The final output therefore represents a relative energy profile for each structure.

Results and Discussion

Structure solution The present work was motivated by the discovery of a new crystal of holoWrbA that grows under conditions very different from all prior WrbA crystals. Deep yellow crystals of WrbA-FMN indicating the presence of oxidized FMN were obtained in the following crystallization condition: 0.1 M MES/imidazole, pH 6.5, 12.5% (v/v) 2-methyl-2,4-pentanediol, 12.5% (w/v) PEG 1000, 12.5% (w/v) PEG 3350, and 0.02 M each of six alcohols (1,6-hexanediol, 1-butanol, 1,2-propanediol, 2-propanol, 1,4-butanediol, 1,3-propanediol). HoloWrbA with BQ (WrbA-FMN-BQ) was obtained by adding BQ to a final concentration of 5 mM to a droplet containing HoloWrbA crystals and soaking for 10 min prior to commencing data collection. Upon solving the new structures the position of BQ was observed to differ substantially from an earlier reported structure with soaked-in BQ (PDB:3b6k)4. Although the origin of this difference is unknown, it could reflect differences in the crystallization conditions, which differ from all previous crystals of WrbA, and/or in the crystal packing, which differs from all previous crystals of WrbA with bound substrates. Complete data collection and processing statistics are reported in Table 1. The structure of WrbA-FMN was solved to 1.5 Å resolution in space group P41212. The asymmetric unit contains two protein monomers (A and B) and has a Matthews coefficient of 1.77 Å3 Da-1 with solvent content of 30.45%.41 The WrbA primary structure consists of 197 residues starting from Ala1.1 The crystal structure shows clear electron density for residues 1-152 and 156-197 of chain A and 1-144 and 155197 of chain B. Residues 153-155 of chain A and 145-154 of chain B, both located in surface loop β5aβ5b, had ambiguous electron density and were omitted from the final structure. High resolution allowed building the side chains of the following residues in two conformations: Ser26, Gln92, Ser111, Ser112, and Thr162 of chain A and Met17, Ser99, Ser111, Ser112, and Ser179 of chain B. The final model has Rwork = 0.178 and Rfree = 0.205. 97.3% of residues are in Ramachandran favored regions

ACS Paragon Plus Environment

Page 11 of 27

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

The Journal of Physical Chemistry

with one outlier Arg171 of chain B located in a surface loop. The crystal structure of WrbA-FMN-BQ was determined in space group P41212 at 1.3 Å resolution. Similarly to the WrbA-FMN structure, residues 144-154 of chain A and 153-155 of chain B had poorly defined electron density and were omitted from the final structure. High resolution allowed building the side chains of the following residues in two conformations: Arg19, Glu32, Thr54, Gln92, Ser99, Ser111, Ser112, Ser114, Thr162, Ser179, and Leu191 of chain A and Lys27, Asp29, Glu32, Met42, Arg78, Gln92, Ser111, Ser112, and Thr162 of chain B. The final model has Rwork = 0.126 and Rfree = 0.167 with 99.2% of the residues in Ramachandran favored regions and 0.8% in allowed regions.

Table 1. Data collection and refinement statistics. Data collection and processing Diffraction source Wavelength (Å) Detector Rotation range per image (°) Total rotation range (°) Space group Unit-cell parameters (Å, °) Mosaicity (°) Resolution (Å) Total reflections Unique reflections Completeness (%) Multiplicity Rmeas (%) Overall B factor from Wilson plot (Å2) Refinement Resolution (Å) Completeness (%) No. of reflections Final Rwork (%) Final Rfree(%) No. of non-H atoms Protein Ligand Water Total R.m.s. deviations Bonds (Å) Angles (°) Mean B factors (Å2) Protein Ligand Water

WrbA-FMN BESSY, BL14.1 0.91841 Pilatus 6M 0.1 120 P41212 a = b = 60.75, c = 169.02, α = β = γ = 90 0.15 50.0-1.5 (1.59-1.5) 432796 51696 (8187) 99.9 (99.4) 8.37 11.7 (2.21) 10.9 (84.4) 24.3

WrbA-FMN-BQ BESSY, BL14.1 0.91841 Pilatus 6M 0.1 150 P41212 a = b = 61.18, c = 169.99, α = β = γ = 90 0.08 43.26-1.33 (1.411.33) 799147 (127436) 141939 (22871) 99.9 (99.8) 5.63 18.02 (2.36) 5.4 (72.6) 20.8

42.96-1.50 99.85 49068 17.8 20.5

42.50-1.33 99.92 71491 12.64 16.73

2847 62 206 3115

2912 78 224 3214

0.02 2.0 21.5 21.5 16.9 30.3

0.013 1.66 23.3 22.4 23.4 35.3

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

Model statistics Ramachandran plot Favored regions (%) Outliers (%) MolProbity score Clashscore PDB code

97.3 1.02 3.64 [96th percentile] 1.26 [96th percentile] 5F12

Page 12 of 27

99.24 0 3.19 [97th percentile] 1.11 [98th percentile] 4YQE

Overall structure Crystals of WrbA-FMN and WrbA-FMN-BQ belong to the same space group and have very similar unit cell parameters (Table 1) and their structures are very similar. The two monomers in the asymmetric unit are related to each other by non-crystallographic two-fold symmetry, and two dimers build up a homotetrameric complex with 222 symmetry as in previously reported structures of WrbA. The overall structure of each monomer exhibits architecture typical for the WrbA family5: five flavodoxin β–α repeats with two characteristic insertions: β6-α6 following β2, and α’ following β5a. The monomers are highly similar to each other, superimposing with RMSD of ~ 0.45 Å over 181 Cα atoms. The most significant difference between the two monomers of each asymmetric unit, and also between the two structures with and without BQ, is observed in the β4-α4 loop situated close to the active site (Figure 2a). In molecule A part of the loop corresponding to residues Thr115-Gly117 is bent toward the active site, and in molecule B this segment is bent in the opposite direction. In both monomers this loop participates in oligomerization and coordination of the isoalloxazine ring. In both WrbA-FMN and WrbA-FMN-BQ crystals the two Thr115-Gly117 segments have similar alternative conformations in each monomer (not shown), indicating that the observed differences between monomer conformations do not reflect BQ binding. A positive peak in the Fourier difference map suggested the presence of Met10 in oxidized form as observed previously in a high-resolution structure of holoWrbA (PDB:3zho)6. Thus, Met10 was modeled as methionine sulfoxide with good fit in both chains.

ACS Paragon Plus Environment

Page 13 of 27

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

The Journal of Physical Chemistry

Figure 2. FMN and BQ. Secondary structure elements of the protein are shown as a cartoon; FMN and BQ are shown as skeletal models in atomic colors with carbons colored to correspond with the monomer. (a) Superimposition of the active sites of WrbA monomers A (cyan) and B (yellow) showing the alternative positions of the β4-α4 loop. Residues 165-171 of the β5b-α5 loops are hidden for clarity. (b) Electron density (gray mesh) of BQ in the 2Fo-Fc map contoured at 1.0 σ. Active site Three of the four subunits contribute residues to each of the four active sites of the tetramer. The active site is formed mainly by loops at the C-termini of the β-strands. FMN has well-defined electron density. Positive and negative peaks in the Fo-Fc difference map of the FMN isoalloxazine ring imply propeller twist along the length of the ring, as reported previously (PDB:3zho)6. After removing restraints on ring planarity these difference peaks were eliminated. In the WrbA-FMN-BQ structure planar electron density located above the isoalloxazine ring was interpreted as a BQ molecule and gave acceptable fit (Figure 2b), although the quality of the electron density map differs between the two monomers. Negative electron density after refinement and B-factors higher than the average for the structure suggest partial occupancy of BQ in this position. Therefore, the occupancy of BQ was set to 0.6 in both monomers. Modeling of water, PEG, ethylene glycol, or each of the alcohols present in the crystallization cocktail resulted in less good fit. Additional electron density in the Fo-Fc difference map connected with the density of BQ is present in monomer A and could indicate partial occupancy of some unknown molecule at this location. Similar observations have been reported previously for the structure PDB:3b6k.4 This additional density was impossible to interpret, and thus it remained unmodeled. Three of the four monomers contribute to coordination of BQ in each active site. BQ is within hydrogen-bonding distance of His132 of the neighboring monomer in the dimer and Tyr142 of the ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 14 of 27

monomer from the symmetry-related dimer. Distances between BQ atom O5 and the closest atoms of His132 and Tyr142 are 2.80/2.99 and 3.11/3.15 Å for A/B monomers, respectively (Figure 3). The position of BQ is nearly identical in both monomers (Figure 2a); only the position of BQ in monomer A is depicted in Figure 3 and further described. BQ is situated perfectly in plane with the FMN pyrimidine ring and the Trp97 side chain of the neighboring monomer of the dimer, forming a stacked structure with distances from the BQ ring center of 3.8 Å to FMN and 3.3 Å to Trp97. The position of the BQ ring appears suitable for electron transfer from FMN isoalloxazine with distances N1(FMN)C2(BQ), 3.37 Å; N3(FMN)-C4(BQ), 3.14 Å; and N5(FMN)-C6(BQ), 3.34 Å. In structure PDB:3b6k4 BQ is rotated ~45° along its O2-O5 axis and shifted ~5 Å toward the dimethylbenzene end of FMN isoalloxazine (Figure 3).

Figure 3. BQ positions. Zoomed-in view of one active site in the present WrbA-FMN-BQ crystal structure model showing monomer A with FMN and BQ (cyan) and nearby regions of monomers B (green) and D (pink). Distances (in Å) discussed in the text between atoms of BQ and WrbA are indicated by black dashed lines. Residues 165-171 of monomer A are hidden for clarity. The two positions of BQ in the previously reported structure (PDB:3b6k)4 are shown in grey and labeled BQ'. QM calculations Geometry optimization As the crystal structures do not necessarily represent minima on the potential energy surface, the

ACS Paragon Plus Environment

Page 15 of 27

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

The Journal of Physical Chemistry

geometry of the binding site in chain A of each structure was optimized by QM/MM (M06-2X-D3/631G*//OPLS2005; for details see Methods) prior to QM calculation of charge transfer. The resulting optimized structures and their comparison with the crystal structures are depicted in Figure S1 (Supporting Information). For the new structure optimization results in only a small change of the position of BQ, whereas in structure 3b6k a larger rearrangement is observed. The difference can be quantified by the RMSD of BQ heavy atoms before and after optimization. The isoalloxazine ring of FMN was used as the reference when overlaying the structures before and after minimization. The resulting RMSD is 1.0 Å for the new crystal structure and 2.7 Å for structure 3b6k (Figures S1A, S1B). In the new structure optimization results in a small rotation of BQ along its long axis (approximately 15 degrees clockwise), a shift of less than 1 Å toward the central isoalloxazine ring of FMN in the plane parallel to FMN, and a tilt of ~5 degrees of the BQ plane (Figure S1A). Changes in FMN position are within the expected dynamic behavior and can be considered negligible (RMSD 0.4 Å for FMN isoalloxazine heavy atoms when the entire protein is used as the reference for overlaying the structures before and after optimization). The changes in structure 3b6k upon minimization are much larger in comparison. The BQ plane is tilted about 40 degrees, bringing it to a position almost parallel to the FMN isoalloxazine ring, and BQ is shifted about 2.3 Å from its initial position to just above the dimethylbenzene ring of isoalloxazine (Figure S1B). FMN also changes position significantly. The isoalloxazine ring plane is tilted by approximately 15 degrees, making it more parallel to the BQ plane, and FMN is shifted about 0.7 Å closer to BQ (Figures S1C, S1D). The internal geometry of FMN is also changed in both structures upon minimization. The isoalloxazine ring becomes deformed and the hydrogen atoms (present in the reduced form only) on nitrogens N1 and N5 are shifted out of the ring plane, indicating that these nitrogen atoms have sp3 hybridization. These changes were quantified using the backbone of the entire protein as the reference, resulting in RMSD values for FMN isoalloxazine heavy atoms before and after optimization of 2.7Å in structure 3b6k and 0.69Å in the new structure (Figures S1C, S1D). These results suggest that the new crystal structure is closer to the local energy minimum than the previously reported structure 3b6k, where the geometry changes extensively during minimization.

Charge-transfer calculations Quantum mechanical calculations based on the Marcus equation were employed to describe the probability of charge transfer (CT) between BQ and FMN in the present and previous crystal structures. A simplified model was used containing only BQ and the isoalloxazine ring of FMN (Figure 4) because the sugar-phosphate tail of FMN can allow partial redistribution of negative charge from FMN ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 16 of 27

onto BQ. This simplification appears justified as electron transfer occurs exclusively between the molecular orbitals of the flavin and BQ rings. The simplified model was prepared for two different geometries, one based on the new crystal structure and the other based on chain A of the previous crystal structure PDB:3b6k, which is the chain with BQ located closer to its position in the new crystal structure. These two positions of BQ are depicted in Figure 4. This simplified model was also prepared for the geometry-optimized structures produced in the previous step by QM/MM calculations as described above. To include the effect of the protein on the QM wave function all protein atoms were described by point charges (atomic charges of OPLS2005 force field) for the case of the simplified model based on optimized geometries. Further simplifications were made in describing electron transfer between FMN and BQ. Because the exact mechanism of the redox reaction is not known, several possible electron transfers routes were considered. First, electron transfer between reduced FMN isoalloxazine (FMNH2) and oxidized BQ was considered, and interactions between the HOMO of FMNH2 and LUMO of BQ were calculated, representing the situation when a single electron transfer is the first step of the redox reaction. Alternatively, electron transfer between FMNH- (with the hydrogen at N9 of FMN) and BQH+ (with the hydrogen on the BQ oxygen located closer to N9 of FMN) was considered, and interactions between the HOMO of FMNH- and LUMO of BQH+ were calculated, representing the situation when electron transfer occurs following one-proton transfer. Additionally, this alternative also describes the situation when reduced FMN is in the deprotonated form FMNH- as discussed in Methods, representing electron transfer from FMNH- to BQ that has received a proton from solvent water. Finally, electron transfer between FMN2- and BQH22+ was considered, and interactions between the HOMO of FMN2and LUMO of BQH22+ were calculated, representing the situation when electron transfer occurs following two-proton transfer or one-proton transfer from FMNH-. The intramolecular geometries of BQ and FMN do not differ significantly between the present and previous crystal structures; thus the driving force (∆G) that depends on the energy of the HOMO/LUMO orbitals of the donor (isoalloxazine) and acceptor (BQ) can be assumed to be approximately the same for both crystal structures. This approximation was evaluated by estimating ∆G for geometrically-optimized structures originating from each crystal structure in the presence of MM point charges representing the protein. The estimation was based on LUMO-LUMO interactions between oxidized forms of both FMN and BQ. ∆G estimated for the optimized structure originating from the new crystal structure is -7.8 kcal/mol, and for the previous crystal structure is -8.3kcal/mol. The difference in ∆G between the two optimized structures (0.5 kcal/mol or ~ 6% of 7.8) can be considered negligible. QM/MM calculations carried out as described in Methods to estimate the ACS Paragon Plus Environment

Page 17 of 27

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

The Journal of Physical Chemistry

reorganization energy of the surrounding protein residues upon reduction indicate that the reorganization energy is approximately the same for both crystal structures as well. In these calculations the reorganization energy of the protein in the present crystal structure is estimated to be 20.1 kcal/mol, and the reorganization energy calculated for structure 3b6k is 18.7 kcal/mol. The internal reorganization energy associated with changes to the geometries of FMN and BQ upon reduction will be the same for both crystal structures. Thus the difference in total reorganization energy of the two crystal structures (20.1 - 18.7 = 1.4 kcal/mol or ~ 7% of 18.7) can be considered negligible. With these assumptions, the electronic coupling constant J that depends on the positions of donor and acceptor can be compared for the two crystal structures as the parameter of the Marcus equation that approximates the probability of CT between FMN and BQ.

Figure 4. Model systems used for charge-transfer calculations. The present and previous crystal structures were aligned on their FMN isoalloxazine ring systems, which superimpose nearly identically; thus only the ring system of the present structure is shown. The positions of BQ from the present (cyan) and previous (grey) crystal structures are shown. Top, view orthogonal to the plane of FMN; lower left, view from the long side of FMN; lower right, view from the short side of FMN.

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 18 of 27

Table 2. Calculated electronic coupling constants, J (absolute values, kcal/mol) redox partners

FMNH2 - BQ FMNH- - BQH+ FMN2- - BQH22+

structure

PDB:4yqe PDB:3b6k PDB:4yqe PDB:3b6k PDB:4yqe PDB:3b6k

crystal structure geometry 6.9 0.7 29.6 13.1 6.2 2.1

QM/MMoptimized geometry 24.7 17.1 53.3 39.8 17.0 9.9

QM/MM-optimized geometry including protein charges 17.9 12.8 48.5 36.4 10.6 5.8

Table 2 presents the absolute values of coupling constants calculated for three possible chargetransfer processes for three geometries of each crystal structure: between reduced FMN and oxidized BQ; between FMNH- and BQH+; and between FMN2- and BQH22+; and in each case using the two crystal structures directly; using QM/MM-optimized geometries for the local active-site environment; and using QM/MM-optimized geometries for the local active-site environment including protein charges. Accounting for the fact that the coupling constant J is squared in the Marcus equation, for the non-optimized structures (column three) the charge-transfer rate between FMNH2 and BQ is ~ 100 times higher in structure 4yqe compared to structure 3b6k [(6.9/0.7)2 = 97], ~ 5 times higher for transfer between FMNH- and BQH+ [(29.6/13.1)2 = 5.1], and ~ 10 times higher for transfer between FMN2- and BQH22+ [(6.2/2.1)2 = 9]. All coupling constants are larger for the geometry-optimized structures; this is the expected outcome of optimization, which brings donor and acceptor to a minimum on the potential energy surface. As with the crystal structure geometries, charge-transfer rates after QM/MM geometry optimization are higher for the new structure than for structure 3b6k for all three redox processes considered, although the magnitudes are smaller, ranging from ~ 2-fold higher for transfer between FMNH2 and BQ (or FMNH- and BQH+) to ~ 3-fold higher for transfer between FMN2- and BQH22+. Finally, inclusion of protein charges produces intermediate values of the coupling constants, with charge-transfer rates ranging from ~ 2-fold higher for transfer between FMNH2 and BQ (or FMNH- and BQH+) to ~ 3.5-fold higher for transfer between FMN2- and BQH22+. Thus in all cases the charge-transfer rate is predicted to be substantially larger for the new crystal structure than for the previous one. The calculations indicate that the protein environment affects the absolute values of the chargetransfer coupling constants, as expected, but has only a minor effect on the ratio of charge-transfer rates between the present structure and structure 3b6k. QM/MM optimization reduces the difference between ACS Paragon Plus Environment

Page 19 of 27

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

The Journal of Physical Chemistry

the present and previous crystal structures and leads to an increase of the absolute values of J for both structures. This can be understood from the structural results of optimization, which led to a far larger change of BQ geometry for structure 3b6k than for the new structure (RMSD 2.7 Å vs. 1.0 Å), shifting BQ into a more favorable position for electron transfer. The charge-transfer calculations thus support the conclusion that the new crystal structure is closer to the optimal geometry (i.e., to a local minimum on the potential energy surface) than structure 3b6k. Even though the geometry of structure 3b6k changes substantially upon QM/MM optimization, and becomes more similar to the new structure, the charge-transfer rate calculated for this new local minimum is still ~two to three times lower than for the optimized new structure. In summary, the charge-transfer rates calculated here indicate that in all cases the new structure is predicted to be considerably more favorable for electron transfer than the previous structure. The more favorable rates calculated for the present structure suggest that the position of BQ found in the present structure is more likely to reflect its position during electron transfer. These results suggest that the new position found for BQ in the present structure may aid computational efforts to dock potential substrates and thus point to the physiological role of WrbA.

Estimation of proton-transfer probability The above analysis and its application to estimate the most favorable geometry of BQ binding for enzymatic reduction rests on the assumption that the first electron transfer from FMNH2 to BQ (or FMNH- to BQH+) is the rate-limiting step of the WrbA redox reaction. In fact the rate-limiting step cannot be predicted with confidence because the reaction mechanism is not known with certainty. Twoelectron transfer is not absolutely certain for WrbA; rather, what can be stated is that transfer of the second electron is apparently very fast relative to the first.35 The typical spectroscopic signatures of one-electron transfer were not detected for WrbA, nor for human QR.10 Nevertheless, one electron must be transferred first. Therefore proton transfer cannot be excluded as the rate-limiting step. Thus, the energy barriers for two-proton transfer from FMNH2 to BQ and for one-proton transfer from FMNH- to BQ were estimated for both crystal structures using semi-empirical (PM7) calculations. The hydrogen atoms to be transferred in the case of two-proton transfer are those on N1 and N5 of FMNH2, and in the case of one-proton transfer it is on N5 of FMNH-. The initial position of BQ in the two crystal structures shows that BQ is closer to both these hydrogen atoms in the new crystal structure than in structure 3b6k, as seen in Figure 3. To reach a conformation favorable for reaction the BQ molecule in the new structure must rotate clockwise along its long axis by ~ 40 degrees. Already the initial unrestrained minimization with the PM7 semi-empirical method results in this conformation, ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 20 of 27

which is shown in Figure S2A (Supporting Information). On the other hand, the initial position of BQ in structure 3b6k is much farther from the hydrogen atoms on both N1 and N5 of FMNH2. In particular, the position of the oxygen located on the side of BQ farthest from FMNH2 is extremely unfavorable for proton transfer. In order to reach a conformation favorable for the two-proton transfer reaction BQ must rotate ~ 90 degrees clockwise and shift ~ 2.3 Å toward the position of BQ found in the present crystal structure. Unrestrained optimization does not result in this conformation for crystal structure 3b6k. Only upon applying distance restraints is the BQ molecule of structure 3b6k brought into a position (Figure S2B) similar to the optimized position found for the new crystal structure after optimization. In the case of one-proton transfer from N5 of FMNH- to the closest oxygen of BQ, the required shift of BQ by ~ 1.2 Å toward the N5 atom of FMNH- in structure 3b6k is smaller compared with two-proton transfer. The energy profiles along the reaction coordinate for two-proton transfer starting from the optimized BQ position of each crystal structure are plotted in Figure 5. The reaction coordinate represents transfer of the hydrogen atom on N5 of FMN isoalloxazine to the closest oxygen atom of BQ, which is the shortest distance for hydrogen transfer between the two reaction partners (see Methods for details). The energy profile for two-proton transfer starting from the optimized new structure is short (~ 2.3 Å to ~ 1.0 Å) and smooth, with a reaction barrier of ~ 65 kcal/mol as the transition state is approached at ~ 1.0 Å. The reaction barrier height is defined as the difference between the energy of the reactive conformation at ~ 2.3 Å on the reaction coordinate and the highest energy as the transition state is approached at ~ 1.0 Å. In contrast, the energy profile for two-proton transfer starting from optimized structure 3b6k has a much longer path (~ 4.5 Å to ~ 1.1 Å) that presents several steep energy barriers. The early barriers of ≥ 120 kcal/mol at ≥ 4.0 Å are those described above that are not overcome without applying distance restraints. These early barriers separate the initial conformation from the reactive conformation that is reached between 2.5 and 2.0 Å. Thereafter the transition state is again approached smoothly as it is with the new crystal structure, with a final reaction barrier of ~90 kcal/mol that is considerably smaller than the early barriers.

ACS Paragon Plus Environment

Page 21 of 27

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

The Journal of Physical Chemistry

Figure 5. Reaction energy profiles and transition-state structures for two-proton transfer. The relative energy profile along the reaction coordinate was calculated using the semi-empirical PM7 method. The horizontal axis represents the reaction coordinate from reactant on the left to product on the right using the distance between the hydrogen atom on N5 of FMN isoalloxazine and the closest oxygen atom of BQ. TS, transition state. The inset structures show the transition-state geometries of FMN and BQ for each corresponding energy profile, with FMN isoalloxazine rings superimposed. Gray, PDB:3b6k; cyan, PDB:4yqe.

The reaction energy barrier of ~90 kcal/mol is almost 50% higher than the barrier starting from the optimized new crystal structure of ~ 65 kcal/mol. This difference reflects a number of factors, including the slight difference evident in the transition-state structures in the overlaid view shown in the Figure 5 inset. However the transition state discussed here is only a relatively crude approximation due to the semi-empirical methods used in the present work, and it does not necessarily represent a true saddle point on the potential energy surface. These factors probably account for the fact that the reaction energy is reached at slightly, and probably negligibly, different distances as the transition state is approached (1.0 Å for structure 3b6k and 1.1 Å for the new structure). Nevertheless, the two reaction energy values contain all the same approximations and sources of uncertainty and thus the large

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 22 of 27

difference in their values is likely to be meaningful. The energy profiles along the reaction coordinate for one-proton transfer starting from the optimized BQ position of each crystal structure are plotted in Figure S3 (Supporting Information). The reaction coordinate represents transfer of the hydrogen atom on N5 of FMN isoalloxazine to the closest oxygen atom of BQ. The products have considerably higher energy (~ 30 kcal/mol for 4yqe and ~ 40 kcal/mol for 3b6k) than the reactants. This results reflects the fact that monoprotonated BQ is less stable than oxidized BQ or fully reduced BQH2, and the expectation that in a following step monoprotonated BQ forms fully reduced BQH2 by accepting a proton from solvent water that is not present in the calculations. Similarly as for two-proton transfer, the energy profile for one-proton transfer starting from the optimized new structure is short (~ 1.8 Å to ~ 1.0 Å) and smooth, with a reaction barrier of ~ 50 kcal/mol as the transition state is approached at ~ 1.0 Å. Similarly as for twoproton transfer, the energy profile of one-proton transfer starting from optimized structure 3b6k has a longer path (~ 3.0 Å to ~ 0.9 Å) that presents several energy barriers. The largest barrier of ~75 kcal/mol at 1.3 Å is associated with the deformation of BQ prior to the transition state, which is then approached smoothly with a reaction barrier that is ~ 10kcal/mol smaller than the preceding one. The transition-state geometry of structure 3b6k and the new structure differ significantly (Figure S3 inset), because in the case of one proton transfer, structure 3b6k does not have to undergo large conformational changes thus its transition-state geometry is closer to the original crystal structure. The much lower reaction energy predicted for the new crystal structure reported here in both two-proton (from FMNH2) and one-proton (from FMNH-) transfer suggests that this structure is more likely to be functionally relevant than the earlier 3b6k structure. This interpretation applies even if proton transfer is the rate-limiting step for WrbA catalysis. Furthermore, the results imply that along the reaction path starting from structure 3b6k the BQ molecule is likely to adopt a position similar to the one observed for the reaction path starting from the new crystal structure.

Conclusions QM and QM/MM methods can fill many roles in crystal structure analysis, e.g., when the interpretation of crystallographic data is limited by resolution, by the static nature of the data, or by the absence of hydrogen atoms. Indeed, another recent analysis of WrbA used QM/MM methods to determine the oxidation state of the FMN isoalloxazine ring6. Despite very high resolution (1.2 Å) the oxidation state could not be determined directly because hydrogen atoms are not detected in x-ray crystal structures. Instead, the oxidation state was investigated using QM/MM calculations, which are uniquely suited to the task, as classical MM force fields do not adequately describe the subtle ACS Paragon Plus Environment

Page 23 of 27

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

The Journal of Physical Chemistry

differences among oxidation states. The present work illustrates that charge-transfer calculations can be used to evaluate the CT probabilities of alternative binding poses for BQ that are observed in crystal structures. To the best of the authors' knowledge this application of QM methods is novel. It may find applications in other charge-transfer processes.

Supporting Information Figures S1-S3 depicting optimized structures, reactive conformations and energy profile for one-proton transfer

Acknowledgements The authors thank I. Iermak and T. Prudnikova for help with collecting the data from WrbA-FMN crystals and the MX user-support team at BESSYII, Helmholtz-Zentrum, Berlin, for their assistance during data collection. Access to instruments and other facilities was supported by the Czech research infrastructure for systems biology C4SYS (project no LM2015055). Support from the Czech Science Foundation (P207/10/1934) and joint Czech - US National Science Foundation International Research Cooperation (OISE08-53423) is acknowledged.

References

1.

Grandori, R.; Khalifah, P.; Boice, J. A.; Fairman, R.; Giovanielli, K.; Carey, J. Biochemical

characterization of WrbA, founding member of a new family of multimeric flavodoxin-like proteins. J Biol. Chem. 1998, 273, 20960-20966. 2.

Carey, J.; Brynda, J.; Wolfova, J.; Grandori, R.; Gustavsson, T.; Ettrich, R.; Smatanova, I. K.

WrbA bridges bacterial flavodoxins and eukaryotic NAD(P)H:quinone oxidoreductases. Protein Sci. 2007, 16, 2301-2305. 3.

Gorman, J.; Shapiro, L. Crystal structures of the tryptophan repressor binding protein WrbA

and complexes with flavin mononucleotide. Protein Sci. 2005, 14, 3004-3012. 4.

Andrade, S. L.; Patridge, E. V.; Ferry, J. G.; Einsle, O. Crystal structure of the NADH:quinone

oxidoreductase WrbA from Escherichia coli. J. Bacteriol. 2007, 189, 9101-9107. 5.

Wolfova, J.; Smatanova, I. K.; Brynda, J.; Mesters, J. R.; Lapkouski, M.; Kuty, M.; Natalello,

A.; Chatterjee, N.; Chern, S. Y.; Ebbel, E. et al. Structural organization of WrbA in apo- and holoprotein crystals. BBA-Proteins Proteom. 2009, 1794, 1288-1298.

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

6.

Page 24 of 27

Kishko, I.; Carey, J.; Reha, D.; Brynda, J.; Winkler, R.; Harish, B; Guerra, R.; Ettrichova, O.;

Kukacka, Z.; Sheryemyetyeva, O. et al. 1.2 A resolution crystal structure of Escherichia coli WrbA holoprotein. Acta Crystallogr. D 2013, 69, 1748-1757. 7.

Patridge, E. V.; Ferry, J. G. WrbA from Escherichia coli and Archaeoglobus fulgidus is an

NAD(P)H:quinone oxidoreductase. J. Bacteriol. 2006, 188, 3498-3506. 8.

Kishko, I.; Harish, B.; Zayats, V.; Reha, D.; Tenner, B.; Beri, D.; Gustavsson, T.; Ettrich, R.;

Carey, J. Biphasic kinetic behavior of E. coli WrbA, an FMN-dependent NAD(P)H:quinone oxidoreductase. PloS one 2012, 7, e43902. 9.

Ernster, L. DT-Diaphorase-a historical review. Chem. Scripta 1987, 27, 1-13.

10.

Li, R. B.; Bianchet, M. A.; Talalay, P.; Amzel, L. M. The 3-dimensional structure of NAD(P)H-

quinone reductase, a flavoprotein involved in cancer chemoprotection and chemotherapy - mechanism of the 2-electron reduction. Proc. Natl. Acad. Sci. U.S.A. 1995, 92, 8846-8850. 11.

Marcus, R. A.; Sutin, N. Electron transfers in chemistry and biology. Biochim. Biophys. Acta

1985, 811, 265-322. 12.

Řeha, D.; Barford, W.; Harris, S. A multi-scale method for the calculation of charge transfer

rates through the Π-stack of DNA: application to DNA dynamics. Phys. Chem. Chem. Phys. 2008, 10, 5436-5444. 13.

Mueller, U.; Darowski, N.; Fuchs, M. R.; Forster, R.; Hellmig, M.; Paithankar, K. S.; Pühringer,

S.; Steffien, M.; Zocher, G.; Weiss, M. S. Facilities for macromolecular crystallography at the Helmholtz-Zentrum Berlin. J. Synchrotron Radiat. 2012, 19, 442-449. 14.

Krug, M.; Weiss, M. S.; Heinemann, U.; Mueller, U. XDSAPP: a graphical user interface for

the convenient processing of diffraction data using XDS. J. Appl. Crystallogr. 2012, 45, 568-572. 15.

Kabsch ,W. Xds. Acta Crystallogr. D 2010, 66, 125-132.

16.

Winn, M. D.; Ballard, C. C.; Cowtan, K. D.; Dodson, E. J.; Emsley, P.; Evans, P. R.; Keegan,

R. M.; Krissinel, E. B.; Leslie, A. G. W.; McCoy, A. et al. Overview of the CCP4 suite and current developments. Acta Crystallogr . D 2011, 67, 235-242. 17.

Long, F.; Vagin, A. A.; Young, P.; Murshudov, G. N. BALBES: a molecular-replacement

pipeline. Acta Crystallogr . D 2008, 64, 125-132. 18.

Vagin, A.; Teplyakov, A. Molecular replacement with MOLREP. Acta Crystallogr. D 2010, 66,

22-25. 19.

Murshudov, G. N.; Skubak, P.; Lebedev, A. A.; Pannu, N. S.; Steiner, R. A.; Nicholls, R. A.;

Winn, M. D.; Long, F; Vagin, A. A. REFMAC5 for the refinement of macromolecular crystal structures. Acta Crystallogr. D 2011, 67, 355-367. ACS Paragon Plus Environment

Page 25 of 27

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

20.

The Journal of Physical Chemistry

Emsley, P.; Lohkamp, B.; Scott, W. G.; Cowtan, K. Features and development of Coot. Acta

Crystallogr. D 2010, 66, 486-501. 21.

Chen, V. B.; Arendall, W. B.; Headd, J. J.; Keedy, D. A.; Immormino, R. M.; Kapral, G. J.;

Murray, L. W.; Richardson, J. S.; Richardson, D. C. MolProbity: all-atom structure validation for macromolecular crystallography. Acta Crystallogr. D 2010, 66, 12-21. 22.

Delano, W. L. The PyMOL Molecular Graphics System 2002.

23.

Small-Molecule Drug Discovery Suite 2015-1: QSite, version 6.6, Schrödinger, LLC, New

York, NY, 2015. 24.

Murphy, R. B.; Philipp, D. M.; Friesner, R. A. A mixed quantum mechanics/molecular

mechanics (QM/MM) method for large-scale modeling of chemistry in protein environments. J. Comp. Chem. 2000, 21, 1442-1457. 25.

Zhao, Y.; Truhlar, D. G. The M06 suite of density functionals for main group thermochemistry,

thermochemical kinetics, noncovalent interactions, excited states, and transition elements: two new functionals and systematic testing of four M06-class functionals and 12 other functionals. Theor. Chem. Acc. 2008, 120, 215-41. 26.

Grimme, S. ;Antony, J.; Ehrlich, S.; Krieg, S. A consistent and accurate ab initio

parametrization of density functional dispersion correction (dft-d) for the 94 elements H-Pu', J. Chem. Phys. 2010, 132, 154104. 27.

Hohenstein, E. G.; Chill, S. T.; Sherrill, C. D. Assessment of the performance of the M05-2X

and M06-2X exchange-correlation functionals for noncovalent interactions in biomolecules. J. Chem. Theory Comput. 2008, 4, 1996–2000. 28.

Josaa, D.; Rodriguez-Oteroa, J.; Cabaleiro-Lagob, E. M.; Rellan-Piñeiro, M. Analysis of the

performance of DFT-D, M05-2X and M06-2X functionals for studying π..π interactions. Chem. Phys. Lett. 2013, 557, 170-175. 29.

Banks, J. L.; Beard, H .S.; Cao, Y.; Cho, A. E.; Damm, W.; Farid, R.; Felts, A. K.; Halgren, T.

A.; Mainz, D. T.; Maple, J. R. et al. Integrated modeling program, applied chemical theory (IMPACT). J. Comp. Chem. 2005, 26, 1752. 30.

Frisch, M. J. et al. Gaussian 09, revision A.1; Gaussian, Inc., Wallingford, CT, 2009.

31.

Chai, J. D.; Head-Gordon, M. Long-range corrected hybrid density functionals with damped

atom-atom dispersion corrections. Phys. Chem. Chem. Phys. 2008, 10, 6615-6620. 32.

Sutton, C.; Sears, J. S.; Coropceanu V.; Brédas J. L. Understanding the density functional

dependence of DFT-calculated electronic couplings in organic semiconductors. J. Phys. Chem. Lett. 2013, 4, 919-924. ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

33.

Page 26 of 27

Zheng, Y. J.; Ornstein, R. L. A theoretical study of the structures of flavin in different oxidation

and protonation states. J. Am. Chem. Soc. 1996, 118, 9402-9408. 34.

Hemmerich, P.; Veeger, C.; Wood, H. C. S. Progress in chemistry and molecular biology of

flavins and flavocoenzymes. Angew. Chem. 1965, 77, 699-716.

35.

Nöll, G.; Kozma, E.; Grandori, R.; Carey, J.; Schödl, T.; Hauska, G.; Daub J.

Spectroelectrochemical investigation of a flavoprotein with a flavin-modified gold electrode. Langmuir 2006, 22, 2378-2383. 36.

Driggers, C. M.; Dayal, P. V.; Ellis, H. R.; Karplus, P. A. Crystal structure of Escherichia coli

SsuE: defining a general catalytic cycle for FMN reductases of the flavodoxin-like superfamily Biochemistry 2014, 53, 3509-3519. 37.

MOPAC2012, Stewart J.J.P., Stewart computational chemistry; Colorado Springs, CO, USA.

38.

Stewart, J. J. P. Optimization of parameters for semiempirical methods VI: more modifications

to the NDDO approximations and re-optimization of parameters, J. Mol. Model. 2013, 19, 1-32. 39.

Černohorský, M.; Kutý M.; Koča J. A Multidimensional driver for quantum chemistry program

MOPAC. Computers & Chemistry 1996, 21, 35-44. 40.

Damborský, J.; Prokop, M.; Koča, J. TRITON: Graphic software for rational engineering of

enzymes. Trends Biochem. Sci. 2001, 26, 71-73. 41.

Matthews, B.W. Solvent content of protein crystals. J. Mol. Biol. 1968, 33, 491-497.

ACS Paragon Plus Environment

Page 27 of 27

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

The Journal of Physical Chemistry

Graphical table of contents:

ACS Paragon Plus Environment