Analysis of the electronic structure of the primary electron donor of photosystem I of Spirodela oligorrhiza by photochemically induced dynamic nuclear polarization (photo-CIDNP) solid-state nuclear magnetic resonance (NMR)

Abstract The electron donor in photosystem I (PSI), the chlorophyll dimer P700, is studied by photochemically induced dynamic nuclear polarization (photo-CIDNP) magic angle spinning (MAS) nuclear magnetic resonance (NMR) on selectively 13C and uniformly 15N labeled PSI core preparations (PSI-100) obtained from the aquatic plant duckweed (Spirodela oligorrhiza). Light-induced signals originate from the isotope-labeled nuclei of the cofactors involved in the spin-correlated radical pair forming upon light excitation. Signals are assigned to the two donor cofactors (Chl a and Chl a ') and the two acceptor cofactors (both Chl a ). Light-induced signals originating from both donor and acceptor cofactors demonstrate that electron transfer occurs through both branches of cofactors in the pseudo- C2 symmetric reaction center (RC). The experimental results supported by quantum chemical calculations indicate that this functional symmetry occurs in PSI despite similarly sized chemical shift differences between the cofactors of PSI and the functionally asymmetric special pair donor of the bacterial RC of Rhodobacter sphaeroides. This contributes to converging evidence that local differences in time-averaged electronic ground-state properties, over the donor are of little importance for the functional symmetry breaking across photosynthetic RC species.


Determination of the level of isotope labeling in Synechocystis sp. PCC 6803 and duckweed by LC-MS.
.1: 13 C isotope incorporation determined by LCMS for Chl a isolated from duckweed leaves grown on unlabeled substrate (A) and with the 13 C 4-ALA precursor in the medium (B).

Model Setup
The desired models were extracted from the crystal structure of Photosystem I in plants (PDB entry 2WSC (Amunts et al., 2010), provided by the Protein Data Bank (Berman et al., 2000). The protein and non-protein part of the models was processed separately. The Reduce program (Word et al., 1999) was used to add missing hydrogen atoms of co-factors. This procedure was not accurate enough to add hydrogen atoms of water molecules, why they were added manually with Avogadro (Hanwell et al., 2012). Missing hydrogens of the protein were added using CHARMM22 topology files (MacKerell et al., 1998;MacKerell et al., 2004) with the Automatic PSF Builder within the VMD (Humphrey et al., 1996) program. Additionally, cut protein bonds were saturated with neutral groups -C(O)CH3 and -NH2 as N-and C-terminus, respectively. The values of bond distances for the C-and N-termini were manually adjusted. The distances of corresponding groups in the crystal structure were added to a modified topology file (see Fig. S2.1). Binding pocket models of PSI were created by specifying radii of 3.2 and 3.4 Å around each atom of the co-factor of interest. All surrounding co-factors, water molecules, and amino acid residues with at least one atom within these radii were included explicitly into the models. Addition of missing hydrogens and cut bond saturation were conducted in a similar fashion as for smaller models considered in this work (see above). The binding pocket models are indicated by abbreviations corresponding to the chosen radii "r32" and "r34", while the isolated co-factor models are abbreviated "iso".

Geometry setup
To assess the quality of the calculated chemical shifts, the structure of PA was extracted from the crystal structure, truncated and processed according to Sec. 2.1. The molecular structure was then fully optimized using the Turbomole v7.4.1 (Ahlrichs et al., 1989;TURBOMOLE, 2019) program package. This model will be called "isoOpt" in the following. For the optimization the def2-SVP basis set (Schäfer et al., 1992) and PBE0 functional (Perdew et al., 1996a;Adamo and Barone, 1999) were used. Additionally, the resolution-of-theidentity approximation in conjunction with the auxiliary Coulomb fitting basis (Schäfer et al., 1992;Schäfer et al., 1994) was enabled. Dispersion corrections were taken into account with the D3 correction with Becke-Johnson damping (Grimme et al., 2010;Grimme et al., 2011). Secondly, the structure of PA was extracted from the crystal structure and optimized using the DFTB3 (Gaus et al., 2012) method, which is described in detail in the following. This model will be called "Pa extract".
For all other geometry optimizations, the DFTB3 (Gaus et al., 2012) method within the AMS-DFTB module from the ADF 2019 package (Amsterdam; Velde et al., 2001) was used. The "Third-Order Parametrization for Organic and Biological Systems" (3ob) (Gaus et al., 2013;Kubillus et al., 2015) parameters from the corresponding Slater-Koster file were used. The optimization was carried out within two steps, where in the first step the coordinates corresponding to hydrogen atoms were optimized, while all other nuclear coordinates were kept fix. In the second step for each Chl co-factor the atoms C-3 1 , C-3 2 and their attached hydrogen atoms were optimized while all other coordinates were kept fix. This step was done, to encounter the poor quality of this type of double bond in the crystal structure, where atoms C-3, C-3 1 and C-3 2 are oriented linearly. The general two-step procedure ensured to keep the relative arrangement of co-factors and environment residues, while the co-factor's geometry is partly relaxed in the presence of the protein pocket.

NMR Calculations
For the "iso" models 15 N and 13 C nuclear magnetic shieldings were calculated using the KT2 (Keal and Tozer, 2003) exchange-correlation functional and a TZP (Perdew et al., 1992) basis set from the ADF 2019 package (Amsterdam; Velde et al., 2001) library. The numerical quality of the density fit and grid construction procedures were set to "good". For binding pocket models, labeled as "r32" or "r34", calculations were carried out within a subsystem DFT approach (Jacob and Visscher, 2006;Jacob and Neugebauer, 2014;Wesolowski et al., 2015) using the TZP (van Lenthe and Baerends, 2003) basis set and the PW91 (Perdew and Wang, 1991;Perdew et al., 1992) exchange-correlation functional with the conjoint (Lee et al., 1991) kinetic-energy functional PW91k (Lembarki and Chermette, 1994). Mutual relaxations of subsystem densities were accounted for using 3 freeze-and-thaw (FaT) cycles (Wesolowski and Weber, 1996). 15 N chemical shifts were calculated with respect to the ammonia shieldings, while 13 C chemical shifts were calculated with respect to tetramethylsilane (TMS). Therefore, both molecules were optimized using the KT2 (Keal and Tozer, 2003) exchange-correlation functional and a TZP (Perdew et al., 1992) basis set from the ADF 2019 package (Amsterdam; Velde et al., 2001) library. 15 N as well as 13 C nuclear magnetic shieldings were calculated, where the chemical shift was calculated by δ " = σ %&' -σ " , where *+, andare the chemical shielding of the reference and atom of interest, respectively. Ring current effects of other subsystems were considered by calculating nuclear independent chemical shifts (NICS) as it was done in Jacob and Visscher (2006).
For comparison of the NMR calculations of the "isoOpt" and "Pa extract" models the ADF program with the PBE functional (Perdew et al., 1996b) and triple-ζ-basis set was used.

Calculation of the Environment effect
To calculate the effect of the protein environment on the chemical shifts of the co-factors the following equation was applied, Δ = 0CoF A/B 78 9 − 0CoF </= iso 9. (1) In the equation above 0CoF A/B 78 9 is the chemical shift of a particular co-factor CoF A/B within a protein environment of 3.4 Å, where 0CoF </= iso 9 is the chemical shift of a particular isolated cofactor. The difference Δ is, thus, the effect of the protein environment on the chemical shifts of CoF.

Validation of the NMR Calculations
In this section results are shown that assess the quality of the calculated chemical shifts. Therefore, 15 N Chemical shifts obtained experimentally by Boxer et al. (1974) and by the authors of this work, which agree very well, are compared to shifts calculated with ADF for the "isoOpt" and "Pa extract" model. The data are shown in Fig. S2.3. The calculated chemical shifts are overestimated by about 10 to 20 ppm compared to the experimental value, but correctly predict the trend of the chemical shifts (see also Fig. S2.4). Additionally, the calculated shifts of NII and NIII for "Pa extract" are in very good agreement with the experimental values, while NI and NIV deviate from the experimental results. The reader should note here, that the results for NIV strongly deviate in all four shown lines (see Figs. S2.3 and S2.4). The deviations in the calculated NMR shifts strongly depend on the geometry of the calculated models. The values calculated for "Pa extract" show good agreement with the experiment especially for nitrogen atoms 2 and 3, while the "isoOpt" values better represent the general trend of the shifts. This shows that a possible error source for the calculated NMR shifts arises from the use of the static crystal structure rather than averaging over conformations accessible during the protein dynamics. This could be assessed by performing a short molecular dynamics simulation and calculating NMR shifts for an ensemble of structures, which is, however, beyond the scope of this work.