Coupled-Cluster Studies of Extensive Green Fluorescent Protein

Jan 22, 2015 - Department of Chemistry, University of Helsinki, P.O. Box 55 (A. I. Virtanens ... a major challenge for modern methods of theoretical c...
0 downloads 0 Views 2MB Size
Article pubs.acs.org/JPCB

Coupled-Cluster Studies of Extensive Green Fluorescent Protein Models Using the Reduced Virtual Space Approach Robert Send,*,†,∥ Carl-Mikael Suomivuori,‡,∥ Ville R. I. Kaila,*,§ and Dage Sundholm*,‡ †

BASF SE, Quantum Chemistry Group, GVM/M - B009, D-67056 Ludwigshafen, Germany Department of Chemistry, University of Helsinki, P.O. Box 55 (A. I. Virtanens plats 1), FI-00014 Helsinki, Finland § Department Chemie, Technische Universität München, Lichtenbergstraβe 4, D-85747 Garching, Munich, Germany ‡

S Supporting Information *

ABSTRACT: Accurate predictions of photoexcitation properties are a major challenge for modern methods of theoretical chemistry. We show here how approximate coupled-cluster singles and doubles (CC2) calculations in combination with the reduced virtual space (RVS) approach can be employed in studies of excited states of large biomolecular systems. The RVS-CC2 approach is used for accurately predicting optical properties of the p-hydroxybenzylidene-dihydroimidazolinone (p-HBDI) chromophore embedded in green fluorescent protein (GFP) models using quantum mechanical calculations in combination with large basis sets. We study the lowest excited states for the isolated and protein-embedded chromophore in two different protonation states, and show how omitting high-lying virtual orbitals in the RVS calculation of excitation energies renders large-scale CC2 studies computationally feasible. We also discuss how the error introduced by the RVS approach can be systematically estimated and controlled. The obtained CC2 excitation energies of 3.13−3.27 and 2.69−2.77 eV for the two protonation states of different protein models are in excellent agreement with the maxima of the experimental absorption spectra of 3.12−3.14 and 2.61−2.64 eV, respectively. Thus, the calculated energy splitting between the excited states of the two protonation states is 0.44−0.52 eV, which agrees very well with the experimental value of 0.48−0.51 eV. The calculations at the RVS-CC2 level on the protein models show the importance of using large QM regions in studies of biochromophores embedded in proteins.

1. INTRODUCTION

In both QM/MM and FDET methods, the separation of the molecular system into active and surrounding regions is, however, neither obvious nor trivial, when chemical bonds connecting the two regions have to be cut. Both the choice of the boundaries between the two regions and the size of the QM system have been shown to significantly affect the obtained results.16,17 Empirical force fields, which are most commonly employed in QM/MM calculations, describe the interaction of the classical surroundings with the active QM region through point charges and Lennard-Jones interactions, which are coupled to the former, e.g., by link atoms or frozen localized orbitals at the boundary regions.1 In photobiological systems, the charge distribution of the environment for the ground and excited states can, however, differ,18 which might lead to large effects in predicted optical properties. Polarization effects in the charge distribution of the surrounding protein can be taken into account in QM/MM calculations by using, e.g., the polarizable embedding (PE) scheme.6,7,19,20 In FDET calculations, the electron density of the surrounding protein is accurately considered, because the FDET scheme can be formulated such

The separation of molecular systems into a region treated using quantum mechanics (QM) for the active part of the system and a region considered using classical force fields (MM) for the surroundings is a general methodology for studying large molecular systems.1 Such a combination of QM and MM methods is an invaluable tool for studying complex photobiological systems, such as biochromophores embedded in their native protein surroundings.2−9 The QM/MM methods provide computationally inexpensive ways to treat the chromophore−protein interactions, which is crucial when aiming at a deeper understanding of photoreactions and spectral tuning mechanisms. The frozen-density embedding theory (FDET) scheme is an alternative approach to model the interactions between the active region and its molecular surroundings without significantly increasing the computational costs.10−13 In the FDET approach, the size of the quantum mechanical system can be increased relative to conventional quantum chemical calculations by dividing the system into an active fully optimized region embedded in a frozen electron density of the surroundings. The general FDET scheme has been successfully employed in modeling solvatochromic and protein induced spectral shifts.11,12,14,15 © XXXX American Chemical Society

Received: December 4, 2014 Revised: January 20, 2015

A

DOI: 10.1021/jp5120898 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

2. THE REDUCED VIRTUAL SPACE APPROACH The main ideas and motivations of the reduced virtual space (RVS) approach are • The CC2 excitation energies are calculated using linear response, suggesting that the most important contributions originate from low-lying virtual orbitals. The contributions from high-lying virtual orbitals are expected to be small due to the large energy denominator. • The CC2 excitation energies can be obtained at much lower computational costs, when the number of virtual orbitals is significantly reduced. • Large basis sets can be employed in CC2 calculations because the number of superfluous virtual orbitals, which increases with the size of the basis set, can be omitted in linear-response CC2 calculations of the excitation energies. • The truncation error introduced by the RVS approach in the calculation of the embedded biochromophore can be accurately estimated by performing a similar calculation on the isolated chromophore. The minimum requirement for obtaining an excitation energy of a protein model is thus one RVS-CC2 calculation on the entire system. The obtained results can be corrected and extrapolated to the full virtual space value by performing CC2 calculations on the chromophore in the gas phase. In a typical RVS-CC2 study, two CC2 calculations on the isolated chromophore are necessary, one with the same RVS threshold and one without the RVS approach. However, in general, the errors introduced by the RVS approach are much smaller than the uncertainty in the CC2 excitation energies. A RVS-CC2 study should thus begin with gas-phase calculations that determine the optimal RVS energy threshold for the studied chromophore. Alternatively, the computationally faster linearresponse time-dependent density functional theory (TDDFT) method can be used for estimating the RVS energy threshold. However, the RVS approach does not significantly speed up the TDDFT calculations, because atomic-orbital based algorithms are employed.22 On the other hand, the RVS idea could be used for rendering accurate calculations of excitation energies on larger molecules at higher-order coupled-cluster levels feasible. The RVS error can then be estimated at a lower correlation level, e.g., at the CC2 level, because high-order correlation contributions to excitation energies seem to be less sensitive to basis sets.23 Due to the reduced computational cost of the RVS-CC2 calculations, first- and second-sphere protein interaction effects can be included in the quantum chemical system, thus avoiding modeling nearby chromophore−protein interaction energies classically, which may be problematic when, e.g., polarization effects are significant. The RVS approach may thus provide a way to assess errors and biases in excitation energies calculated at different QM/MM levels by increasing the boundaries of the active quantum region. RVS calculations have so far been employed only in calculations of excitation energies, whereas RVS-CC2 calculations of gradients and other response properties might also be possible, but such calculations have not yet been published. The idea of introducing a cutoff threshold for canonical Hartree−Fock orbitals can certainly be refined, as there exists a long tradition of more elaborate approaches to virtual-space reduction in post-Hartree−Fock calculations.24−32

that the mutual ground-state polarization of the molecule and the surrounding molecule fragments is included. However, changes in the electron density of the surroundings upon excitation are generally not considered in FDET calculations.11 Although QM/MM and FDET calculations have been employed in many photobiological studies, the accuracy and reliability of the methods have not been thoroughly assessed using pure QM approaches, where the central chromophore− protein interactions are described quantum chemically, because benchmark calculations at QM levels are computationally expensive on extended systems. In addition, computational studies of excited-state properties often require high-level computational methods that consider electron-correlation effects. Filippi et al.4 recently presented careful benchmark studies on green fluorescent protein (GFP) models employing highlevel QM methods on the p-hydroxybenzylidene-dihydroimidazolinone (p-HBDI) chromophore, which was electrostatically embedded in the protein surroundings. Similar QM/MM calculations were also recently reported by Steindal et al.6 and Beerepoot et al.7 In the recent density functional theory (DFT) calculations using the QM/MM schemes, the excitation energies were calculated at the CAM-B3LYP level,21 employing ordinary electrostatic and polarized embedding (PE CAMB3LYP) schemes, respectively, to model the surrounding protein effects.4,6,7,19,20 The excitation energies obtained in the CAM-B3LYP/MM studies differed by ∼0.5 eV for the Aand B-forms of the protein, implying that the polarization effects might be significant, and thus important when aiming at a balanced description of the photophysical properties of the two GFP states. However, the recent CAM-B3LYP/MM studies employing the PE scheme showed that polarization effects red-shift the excitation energies for the A- and B-forms by only 0.1 eV, suggesting that the embedding scheme is not the main reason for the large deviation of the CAM-B3LYP excitation energies from experimental values.6,7 In this work, we employ approximate singles and doubles coupled-cluster (CC2) calculations in combination with the reduced virtual space (RVS) approach in extensive studies of the excitation energies of the two protonation states of the isolated p-HBDI chromophore as well as on p-HBDI embedded in a GFP model treated quantum mechanically using large cluster models. In the RVS-CC2 calculations, all virtual orbitals above a given orbital-energy threshold are neglected in the CC2 calculation, significantly lowering the computational costs. In order to rigorously benchmark the error introduced by the RVS-CC2 calculations, we systematically study the effect of the virtual-orbital cutoff energy as well as basis set effects on the isolated chromophores. The article is outlined as follows: After a description of the computational methodology and presentation of the p-HBDI chromophore and GFP protein models in sections 2 and 3, we discuss the applicability of the RVS approach on the isolated pHBDI chromophores in section 4. The RVS-CC2 results obtained for the protein models are presented in section 5. In section 6, we compare the calculated excitation energies to values deduced from experimental UV−vis spectra and to results obtained in QM/MM and FDET calculations. The performance of TDDFT excitation energy calculations of the two protein models is discussed in section 8. Our conclusions are summarized in section 9. B

DOI: 10.1021/jp5120898 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B The first RVS-CC2 studies on a rhodopsin and related cone pigment models yielded excitation energies in close agreement with experiment,18,33 which was not completely unexpected, as earlier CC2 calculations of the excitation energies of isolated retinal models also yielded values in good agreement with experiment and with results obtained at other ab initio levels of theory.33−37 These RVS-CC2 studies thus confirmed the applicability of the RVS approach. They showed that the computational costs are significantly reduced and that the RVS errors behave very similarly in the protein environment and in gas-phase calculations.18,33 Moreover, our recent study of the GFP suggested that RVS-CC2 performs well also for this system, yielding excitation energies in close agreement with experiment, as well as novel molecular insights about steric and electrostatic spectral tuning mechanisms.38 Filippi et al.4,8 recently illustrated problems which might arise when using small quantum chemical regions in QM/MM treatments of photobiological systems. In their GFP study,4 the p-HBDI chromophore was treated quantum chemically employing several DFT and ab initio methods, whereas the protein surroundings were modeled classically. Despite careful benchmarking of several QM methods, the experimental excitation energies were not reproduced. It was therefore suggested that polarizable force-field calculations would be necessary to reproduce the protein shift.4,39 However, the very recent GFP studies by Steindal et al.6 and by Beerepoot et al.7 indicate that the problem is more complicated. Due to the reduced computational cost, the RVS approach allows modeling of nearby environmental effects at the same level as the chromophore, thus avoiding possible problems originating from nonpolarizable protein environments. We also show how the RVS error can be controlled and benchmarked. Using the RVS methodology, we obtain at the CC2 level an excitation-energy separation of the two protonation states of GFP that reproduces the experimental value. The calculations suggest that CC2 calculations yield excitation energies for GFP that are in good agreement with the experimental values.40−42

PW91 functional for the nonadditive exchange correlation energy.50,65 The quantum system studied in the FDET calculations consisted of the 158 atom protein models constructed from the crystal structure of the wild-type GFP (PDB ID: 1GFL)66 with only the chromophore included in the active region. The FDET calculations were performed using the Amsterdam density functional (ADF) program, version 2013.1.67 3.2. RVS Excitation Energies. The RVS excitation energies were calculated at the CC2 level in combination with the TZVP basis sets. The 1s orbitals of C, N, and O were uncorrelated in the CC2 calculations, which is known as the frozen core approximation (FCA). In addition, all orbitals with orbital energies above the RVS-energy threshold are omitted in the CC2 calculation. The RVS error is defined as the difference between excitation energies obtained with and without the RVS approximation. Since the FCA introduces errors of less than 0.01 eV in the CC2 excitation energies, it is of minor importance for the discussion of the RVS error.33 In the present work, the reference values for the RVS error are obtained without the FCA. Extrapolated excitation energies for the GFP models are obtained by correcting the RVS excitation energies obtained for the protein model using the RVS error calculated for the isolated chromophore. 3.3. Isolated Chromophores. The molecular structures of the isolated p-hydroxybenzylidene-dihydroimidazolinone (pHBDI) chromophore and its anionic form (p-HBDI−) are shown in Figure 1A. The structures were optimized at the

3. COMPUTATIONAL DETAILS 3.1. Methods and Basis Sets. The excitation energies were calculated at the CC2 level using the resolution-of-theidentity (RI) approximation43−47 and at the TDDFT level using the B3LYP, PWLDA, BHLYP, and CAM-B3LYP functionals.21,22,48−54 The Karlsruhe split-valence polarization basis sets (def2-SVP) and triple-ζ polarization basis sets (def2TZVP) were employed.55,56 The aug-QZVP basis set used for the isolated chromophores is a combination of quadruple-ζ polarization basis sets (def2-QZVP) and diffuse basis functions taken from Dunning’s aug-cc-pVQZ basis sets.57 The def2 prefix is omitted in the following. We also performed DFT and CC2 calculations on the chromophore embedded in pointcharge surroundings. The Q-Chem program58 was used in the CAM-B3LYP calculations, whereas the other TDDFT and CC2 calculations were performed using Turbomole versions 6.2− 6.5.59,60 The excitation energies were also calculated at the B3LYP level employing the frozen-density embedding theory (FDET) approach.10,14,61 The frozen density of the environment was generated by performing DFT calculations at the BP86 level on the isolated environment using double-ζ Slater-type-orbital (STO) basis sets augmented with polarization functions (DZP).48,62−64 In the FDET calculations, the GGA97 functional was used for the nonadditive kinetic energy and the

Figure 1. (A) The molecular structures of p-hydroxybenzylideneimidazolinone in the neutral (p-HBDI) and anionic (p-HBDI−) forms. (B) The molecular structure of the 1EMB-I-H model for the two protonation states of the GFP photocycle: GFP-A, where the p-HBDI chromophore is protonated and Glu-222 deprotonated, and GFP-B, where the p-HBDI chromophore is deprotonated and Glu-222 is protonated. In the 1EMB-I-H model, the hydroxyl group of Thr-203 is directed toward the chromophore. The structures of the isolated chromophores and protein models were optimized at the MP2/TZVP and B3LYP/SVP levels of theory, respectively.

second-order Møller−Plesset perturbation theory (MP2) level using the TZVP basis sets and the RI approximation.44 Vertical CC2 excitation energies were calculated using the same basis sets. All orbitals were correlated in the calculations of the excitation energies. 3.4. Protein Models. Two molecular models of the green fluorescent protein (GFP) from Aequorea victoria were constructed on the basis of the crystal structure obtained from the Brookhaven protein data bank (PDB) (PDB ID: C

DOI: 10.1021/jp5120898 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B 1EMB),68 as described in ref 38. The models comprise residues Thr-62, Gln-69, Gln-94, Arg-96, His-148, Val-150, Thr-203, Ser-205, and Glu-222, in addition to the p-HBDI chromophore, four crystallographically observed water molecules, and part of the backbone of Phe-64 (atoms C, O, and Cα) and Val-68 (atoms N and Cα). The amino acid residues were cut at the Cβ atoms, which were fixed during structure optimization at the B3LYP/SVP level. The molecular structures of the GFP mutant model (1EMB-I-H), comprising 161 atoms, were optimized in two protonation states corresponding to the A- and B-forms of the protein. For GFP-A, the phenol group of the p-HBDI chromophore is protonated and Glu-222 is deprotonated, whereas, for GFP-B, the phenol group of the chromophore is deprotonated (p-HBDI−) and Glu-222 is protonated. In the structure optimization of the 1EMB protein models, we also find an alternative structure of the system for the GFPA and GFP-B states (1EMB-II). The most significant differences between the 1EMB-I and 1EMB-II structures are the location of one of the water molecules and the orientation of the molecular plane of His-148. In Figure 2, the molecular

however, possible that one of the two alternative structures might be better stabilized by the protein environment. It has been suggested that the hydroxyl group of the residue Thr-203 is directed away from the chromophore in the GFP-A form, whereas it forms a hydrogen bond to the chromophore in GFP-B.68 Both of these conformations are present in the 1EMB crystal structure,68 which is why we use it as a model in this work. To study the importance of the hydrogen-bond network around the chromophore, Thr-203 was modeled in both conformations for the GFP-A and GFP-B forms of the 1EMB-II model. In 1EMB-II-H, Thr-203 forms a hydrogen bond to the chromophore, whereas, in 1EMB-II, the hydroxyl group of Thr203 is directed away from the chromophore. The molecular structures of the two protonation states of the 1EMB-I and 1EMB-II-H protein models are shown in Figures 1b and 3, respectively. It should be noted that the protein structure deposited with the PDB ID: 1EMB is not the correct model of S65T mutated GFP, since the methyl group that replaces one of the hydrogens of Ser-65 is substituted in the incorrect position, leading to a structure that deviates from L-threonine. The 1EMB structure is though very similar to the wild-type GFP structure, the difference being one of the hydrogens of Ser-65 replaced by a methyl group. According to the X-ray structure, the hydrogen bonding network is not affected by the methyl substitution in 1EMB,68 whereas the S65T mutant with the correct chirality at Thr-65 (e.g., PDB ID: 4EUL, 2Y0G) has a significantly different hydrogen bonding network in the vicinity of Glu-222, leading to large changes in the optical properties.69,70 The 1EMB models were used mainly to study how the conformation of Thr-203 affects the excitation energies of the GFP chromophore. To assess how different protein models affect the calculated excitation energies, we also constructed the corresponding protein models based on the crystal structure of wild-type GFP (PDB ID: 1GFL),66 comprising 158 atoms, where the difference of three atoms relative to 1EMB originates from the methyl group at Thr-65. The hydroxyl group of the Thr203 residue is here assumed to be directed away from the chromophore in the GFP-A and GFP-B states of the 1GFL model. The B3LYP/SVP structures of the two protonation states of the 1GFL models, which were optimized similarly as for the 1EMB models, are shown in Figure 4. To probe the performance of earlier structural models of the GFP systems using the RVS/CC2 methodology, we also reoptimized the 158 atom models for the GFP-A and GFP-B states based on QM/MM-optimized structures provided by

Figure 2. Superposition of two alternative optimized structures (1EMB-I-H in cyan and 1EMB-II-H in purple) of the GFP-A state of the 1EMB protein model obtained using different starting geometries. The most significant structural differences are the location of one of the water molecules and the orientation of the His-148 residue.

structure of the A-state of 1EMB-I-H is superimposed on the structure of the corresponding state of the 1EMB-II-H model. Although both structures are feasible in the cluster models, it is,

Figure 3. Optimized molecular structures of the A (left) and B (right) states of the mutated 1EMB-II protein model; the hydroxyl group of Thr-203 points away from the chromophore. D

DOI: 10.1021/jp5120898 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B

Figure 4. Optimized molecular structures of the A (left) and B (right) states of the 1GFL model. In the wild-type protein model, Thr-203 points away from the chromophore.

Figure 5. Optimized molecular structures of the A (left) and B (right) states of the QM/MM-F model. In the wild-type QM/MM-F model, the Thr203 is pointing away from the chromophore in A and hydrogen bonded to the chromophore in B.

Table 1. Abbreviations Used for the Protein Modelsa

Filippi et al. (QM/MM-F).4 The QM/MM-F cluster models were optimized at the B3LYP/SVP level as described for the 1GFL and 1EMB and models above. The molecular structures of the two protonation states of the QM/MM-F protein model are shown in Figure 5. In the QM/MM-F-PC models, we employed the QM/MM optimized protein model as reported by Filippi et al.,4 by computing the excitation energies of the chromophore at the CC2 level, without reoptimizing the structure. The interactions with the complete protein and surrounding water molecules were described using point charges obtained from the CHARMM22 force field.71 We also performed QM(CC2)/ MM calculations with the same system size as the QM cluster models, i.e., where only the first solvation sphere of amino acids and water molecules is treated at the MM level. The abbreviations of the employed protein models are shown in Table 1. The Cartesian coordinates of the studied protein models are given as Supporting Information.

abbreviation 1EMB-I-H 1EMB-II-H 1EMB-II 1GFL QM/MM-F QM/MM-F-PC

explanation H-bond between Thr-203 and p-HBDI H-bond between Thr-203 and p-HBDI The OH group of Thr-203 points away from p-HBDI. The OH group of Thr-203 points away from p-HBDI. constructed from QM/MM-optimized structures QM/MM-optimized GFP protein model with point charges

ref 38, 68 PW, 68 PW, 68 PW, 66 PW, 4, 66 4, 66

a

1EMB-I-H refers to the structures of our previously optimized protein model,38 whereas the 1EMB-II structures have been optimized in this work (PW). QM/MM-F refers to the structure of a QM optimized model, whose initial QM/MM coordinates were taken from Filippi et al.4 QM/MM-F-PC refers to structures that were adopted from Filippi et al.4 In the QM/MM-F-PC calculations, the chromophore is described at the CC2 level and the remaining protein by classical point charges.

4. EXCITATION ENERGIES OF THE ISOLATED CHROMOPHORES The performance of CC2 and other quantum chemical methods for describing excited states of the isolated p-HBDI/ p-HBDI− chromophores has been benchmarked in previous studies.4,33,72 As shown in Table 2, the results suggest that the CC2 method accurately describes the lowest excited states of the isolated neutral and anionic p-HBDI chromophores.33 The experimental excitation energy for p-HBDI− is obtained from gas-phase photodissociation studies using ion-storage-ring experiments73 and from gas-phase measurements of the

electronic action spectra of the anion.74 Comparisons of calculated and experimental excitation energies are not straightforward for the isolated neutral chromophore. The experimental excitation energy for the neutral p-HBDI chromophore is obtained from spectroscopic measurements in solutions using different solvents and extrapolated to the vacuum dielectric limit.75 In the ion-storage-ring experiments, the molecules must have a net charge, which is why a protonated ethylamine group was chemically attached to the neutral p-HBDI chromophore.73 It was suggested that the E

DOI: 10.1021/jp5120898 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B Table 2. Basis Set Convergence of the Three Lowest Excitation Energies (EE in eV) of the Anionic and Neutral Forms of p-HBDIa basis set

EE1

TZVP aug-TZVP diff-aug-TZVP d-aug-TZVP QZVP aug-QZVP ext.

2.92 2.77 2.43 2.56 2.88 2.73

TZVP aug-TZVP QZVP aug-QZVP ext.

3.71 3.66 3.68 3.63

EE2 p-HBDI− 3.30 2.85 2.77 2.82 3.30 2.81 p-HBDI 3.90 3.87 3.89 3.86

EE3

Exp.

3.30 3.24 2.85 2.84 3.99 2.93

2.7674 2.5973

4.53 4.49 4.50 4.46

3.5175

Figure 6. Error of the RVS approach (RVS error) for the first singlet excitation energy given as a function of the RVS energy threshold. The results for the isolated chromophores are denoted as p-HBDI and pHBDI−. All values are given in eV.

Bright states are marked in boldface. For p-HBDI−, dark states appear below the bright state when diffuse functions are included in the basis set. diff-aug-TZVP refers to aug-TZVP with one function on a dummy atom in the center of mass with an exponent of 10−50, and d-augTZVP refers to d-aug instead of aug on all atoms except hydrogens. The data were obtained from ref 33 except for the diffuse basis set calculations for p-HBDI− performed in this study. a

We obtain a vertical detachment energy (VDE) of 2.43 eV for p-HBDI− at the CC2 level, using the approach proposed by Stanton and Gauss in the computation of ionization potentials.76 The VDE value calculated at the CC2 level agrees well with the value of 2.48 eV obtained at the EOM-IP-CCSD/ 6-31+G(d) level,77 whereas recent experiments suggest a larger VDE of 2.8(1) eV.78,79 The experimental value for the VDE indicates that the first excited state of p-HBDI− is bound,78,79 whereas the CC2 and EOM-IP-CCSD calculations suggest that the bright state of p-HBDI− lies above the ionization limit as the VDE is ∼0.4 eV below the calculated excitation energy of the bright state at 2.81 eV. Our calculations using very diffuse basis sets indicate the existence of a dark state below the lowest bright state, which most likely corresponds to the ionization limit in the employed basis set.

ethylamine substituent does not significantly affect the excitation energy of the chromophore, but CC2 calculations suggest that the ethylamine substitution red-shifts the excitation energy by 0.53 eV.33,72 The extrapolated CC2/aug-QZVP excitation energies for the lowest bright singlet states of p-HBDI and p-HBDI− are 3.63 and 2.81 eV, which agree well with the corresponding experimental excitation energies of 3.51 and 2.76 eV.74,75 An experimental excitation energy of 2.59 eV has also been obtained for p-HBDI−.73 In our previous benchmarking study, we incorrectly reported a CC2 excitation energy of 2.73 eV for the first bright state of p-HBDI−,33 which was obtained by adding contributions of diffuse functions at the TZVP level to the QZVP values. We find, however, that the order of the lowest excited states changes when diffuse basis functions are added. With additional Dunning-type diffuse functions included in the basis sets, we obtain a dark lowest excited state for pHBDI−, which lies less than 0.1 eV below the bright state at 2.81 eV. This result suggests that the basis set requirements for the first bright singlet state are fulfilled at the TZVP level used in this work. The difference to the extrapolated aug-QZVP level is 0.08 and 0.11 eV for p-HBDI and p-HBDI−, respectively. However, the basis set requirements for the p-HBDI − chromophore embedded in the protein can be expected to be more reminiscent of the neutral chromophore due to the positive charge of the surrounding protein residues. The dependence of the RVS error on the virtual orbital energy cutoff for the isolated chromophores is shown in Figure 6. The graph suggests that the RVS error is