Assessment of Functionals for TDDFT Calculations of One- and Two

Subscriber access provided by University of Rhode Island | University Libraries

Spectroscopy and Excited States

Assessment of functionals for TDDFT calculations of one- and two-photon absorption properties of neutral and anionic fluorescent proteins chromophores Dawid Grabarek, and Tadeusz Andruniow J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.8b00769 • Publication Date (Web): 28 Nov 2018 Downloaded from http://pubs.acs.org on November 30, 2018

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 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 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.

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 56 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

Journal of Chemical Theory and Computation

Assessment of functionals for TDDFT calculations of one- and two-photon absorption properties of neutral and anionic fluorescent proteins chromophores ´ ∗ Dawid Grabarek and Tadeusz Andruniow Advanced Materials Engineering and Modelling Group, Wroclaw University of Science and Technology, Wyb. Wyspianskiego 27, 50-370 Wroclaw, Poland E-mail: [email protected]

Abstract Performance of DFT functionals with different percentage of exact Hartree-Fock exchange energy (EX) is assessed for recovery of the CC2 reference one- (OPA) and two-photon absorption (TPA) spectra of fluorescent proteins chromophores in vacuo. The investigated DFT functionals, together with their EX contributions are: BLYP (0%), B3LYP (20%), B1LYP (25%), BHandHLYP (50%) and CAM-B3LYP (19% at short range and 65% at long range). Our test set consists of anionic and neutral chromophores as naturally occuring in the fluorescent proteins. For the first time, we compare TDDFT and CC2 methods for higher excited states, than the S1 state, exhibiting relatively large TPA intensity. Our TDDFT results for neutral chromophores reveal an increase of excitation energies as well as TPA and OPA intensities errors, compared to CC2-derived results, as the DFT functional contains less exact exchange. The long-range-corrected CAM-B3LYP functional performs the best, closely fol-

1

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 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

lowed by BHandHLYP, while BLYP usually significantly underestimates all investigated spectral properties, hence being the worst in reproducing the reference CC2 results. The hybrid B3LYP and B1LYP functionals can be roughly placed in between. We propose that TDDFT may underestimate the TPA intensities for neutral chromophores of fluorescent proteins due to underestimated oscillator strengths between some excited states. In case of anionic chromophores, we find that B3LYP and B1LYP functionals overcome others in terms of reproducing CC2 excitation energies. On the other hand, however, TPA intensity is usually significantly underestimated and in this respect CAM-B3LYP functional seems to be again superior. In contrast to the case of neutral chromophores, it seems that large magnitude of excited-state dipole moment or change in dipole moment upon excitation may be the driving force behind high TPA transition moments.

1 Introduction Since discovery of the green fluorescent protein (GFP) from Aequorea victoria jellyfish in 1970’s by Osamu Shimomura and co-workers, 1 fluorescent proteins (FPs) became an important tool in modern studies of biosystems in vivo at the molecular level. 2–7 In mid-90s Tsien and co-workers obtained novel, genetically modified fluorescent proteins based on GFP. 8 It was probably the first proof that the chromophore’s protein environment influences virtually all its spectral properties. 9–25 Since that, a pallete of FPs has been developed with absorption and emission spectra covering whole visual spectrum as well as near ultraviolet and near infrared. 3,5,6,26 Fluorescent proteins own their success due to several reasons: (i) the chromophore (Figure 1) – “the heart” of every fluorescent protein – is created in a series of autocatalyzed reactions 10,20,27–36 without a need of external cofactors except for molecular oxygen, (ii) non-toxicity and high-stability in host system, (iii) ability to create fusion proteins with host’s native proteins preserving their natural activity. One of the most recent and rapidly emerging directions of FPs’ development is their application for visualization techniques based on two-photon absorption (TPA) process 37–41 instead 2


Page 2 of 56

Page 3 of 56 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


Figure 1: Investigated chromophores models. The neutral structures are gathered in the top and middle rows. The anionic ones in the bottom row. of one-photon absorption (OPA). The former was first theoretically predicted by Maria G¨oppertMayer in 1931. 42 TPA-based visualization allows for utilization of less invasive and phototoxic light source (red and infrared) for which the biomaterials (e.g. tissues) are more transparent, i.e. allowing for deeper visualization. What is more, quadratic dependence of the TPA processes on light intensity offers three-dimensional imaging with high resolution and less out-of-focus bleaching. 43,44 Fluorescent proteins’ one-photon spectral properties were subject of extensive experimental 13,14,16–18,20,22–26,45,46 and theoretical 9,11,15,19,21,47–59 investigations, while two-photon absorption properties of these markers seem to be significantly lagged and only recently have drawn significant attention from both theoreticians 60–69 and experimentalists. 70–75 It seems that molecular modelling will play particularly significant role in investigating TPA properties of FPs. This is because, currently applied experimental techniques for TPA intensity measurements struggle with other nonlinear processes that accompany the TPA measurements, 76 i.e. providing TPA cross-section (σ T P A ) values differing even by two orders of magnitude for the same fluorescent protein as may be concluded from literature data for the enhanced green fluorescent protein (EGFP). 37,74,77–79 Even

3



more, efficiency of these processes may vary between chromophores, as suggested by results of Patterson et al. 80 i.e. hindering comparison of σ T P A values obtained for FPs with different chromophores and/or different chromophore’s protein environment (amino acid residues). Finally, the non-linear FPs’ chromophores photolysis (photobleaching) 73 depends in a complex manner on laser parameters and currently cannot be predicted quantitatively. 81–83 The computational chemistry approaches offer a way to bypass the experimental issues as well as provide consistent results if one theoretical ansatz is utilized. Recent works by Kongsted and co-workers 64,65 seem to be promising and tempting in terms of rational engineering and obtaining novel fluorescent proteins for TPA microscopy (i.e. with high σ T P A value). Currently, however, literature data on computational methodology impact on two-photon absorption cross-section values for fluorescent proteins chromophores are incomplete and not systematic. This raises a question about a methodology that is on one hand accurate and on the other computationally cheap in a sense of required CPU and memory resources. Possibly, the time-dependent density functional theory (TDDFT) could provide a good compromise. However, the question remains which DFT functional should be chosen for TPA properties calculations? In 2014, Salem and Brown 61 attempted to assess the performance of TDDFT utilizing hybrid PBE0, B3LYP, hybrid long-range corrected CAM-B3LYP and long-range corrected LC-BLYP DFT functionals for two-photon and one-photon absorption intensities and excitation energies (∆E) calculations against experimental data and approximated second-order coupled cluster (CC2) calculations results. They have investigated chromophores differing in both electronic structure and electric charge. The authors conclude that long-range corrected DFT functionals (CAM-B3LYP, LC-BLYP) usually provide larger excitation energies and σ T P A values compared to two other hybrid functionals for the S0 → S1 transition. The discrepancy seems to strongly depend on the chromophore’s electronic structure and may be up to few tenths of electronovolts (eV) and few tens of G¨oppert-Mayers (GM, where 1 GM is 10−50 cm4 ·s/photon) for the excitation energy and two-photon absorption cross-section, respectively. In case of the S0 → S1 σ T P A , again all DFT functionals underestimate this value compared to the CC2 method about 3–10 times, although in all cases, the CAM-B3LYP and LC-

4


Page 4 of 56

Page 5 of 56 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


BLYP functionals gave results in the best agreement with CC2 data. For the transitions to higherenergy excited states, the long-range corrected functionals provide larger excitation energies by up to ca. 1 eV and larger σ T P A , even by two orders of magnitude, compared to PBE0 and B3LYP functionals. The higher excited states were not investigated at CC2 level of theory, even though according to TDDFT results they may show distinctively high TPA cross-section as originally shown theoretically by Nifos`ı and Luo. 60 Beerepoot et al. 68 investigated performance of TDDFT utilizing CAM-B3LYP functional against the CC2 methodology in recovering two-photon absorption spectrum of the yellow fluorescent protein (YFP) chromophore in the gas phase within 250–450 nm (2.76–4.97 eV) spectral region, i.e. for transitions to excited states higher in energy than the S1 state. The authors conclude that TDDFT/CAM-B3LYP and CC2 methods are in a qualitative agreement in terms of description of intermolecular charge-transfer transitions and that the former approach underestimates the σ T P A value approximately by an order of magnitude. In their other work on chosen neutral chromophores of FPs, Beerepoot et al. 66 benchmarked the TDDFT/CAM-B3LYP and CC2 results against equation-of-motion second-order coupled cluster methodology (EOM-CCSD) 84 for the excitation energy, σ T P A and OPA oscillator strength (f ) values for transition to the S1 state. The authors found that both approximated and full second-order CC methods are in a good agreement in terms of all investigated absorption properties. The former approach overestimates the σ T P A values by up to 40% for the lowest A’-symmetry excited state. On the other hand, the CAM-B3LYP functional underestimated the two-photon absorption cross-section by a factor of 1.9–3.7, i.e. in a quantitative agreement with results obtained by Salem and Brown. 61 The good performance of the CAM-B3LYP functional against approximated third-order coupled cluster (CC3) method for evaluating the σ T P A was also shown by Paterson et al. 85 and Knippenberg et al. 86 for small organic and inorganic molecules. In the former work, it was also shown that BLYP and B3LYP functionals perform visibly worse in this respect, showing importance of both exact exchange energy and long-range correlation effects for quantitative recovery of TPA intensities. Analogously, Friese et al. 67 concluded that the σ T P A values for both charge-transfer and non-charge-transfer transitions

5



for benzene and azobenzene derivatives, AF240 dye as well as large molecular tweezer complex exhibiting states with intermolecular charge transfer obtained with CC2 and TDDFT/CAM-B3LYP methods are in a good agreement. To summarize, based on analysis of literature data for FPs’ chromophores, TDDFT/CAMB3LYP method seems to provide reliable TPA cross-section values for the lowest excited state, usually exhibiting quite low σ T P A values, as compared to second-order coupled cluster methods. We note, however, that most of these reports 61–63 do not assess performance of TDDFT with respect to coupled cluster methods in calculating σ T P A for transitions to higher-lying excited states. Even though they reveal high TPA cross-section values as firstly shown theoretically (TDDFT/B3LYP) by Nifos`ı and Luo 60 in vacuo and later on confirmed by experiments on a variety of fluorescent proteins. 71 It is also not clear which DFT functional will be the best choice for simultaneous simulation of TPA and OPA spectra of FPs’ chromophores. Aim of the present work is to investigate impact of several DFT functionals, among them gradient-corrected (GGA) BLYP, 87,88 hybrid B3LYP, 89 B1LYP, 90 and BHandHLYP, 91 as well as long-range corrected hybrid CAM-B3LYP 92 on two-photon absorption cross-section, one-photon absorption oscillator strength and vertical excitation energies to the S1 and higher excited states. Our choice of DFT functionals is dictated by increasing percentage of the exact HF exchange in energy functional formula, i.e. we will be able to assess importance of this parameter for accuracy of obtained values (see Methodology section for details). The reference method will be the CC2 93 approach which, as discussed in previous paragraph, seems to be a reliable methodology for TPA and OPA properties calculations for FPs’ chromophores 61,66–68 as well as for other small- and medium-sized organic molecules in terms of excitation energies and OPA oscillator strengths 94–99 compared to higher-order coupled cluster or CASPT2 (complete active space second-order perturbation theory) approaches. We will also investigate how the size of chosen Dunning-style correlation-consistent basis sets 100,101 impacts ∆E, σ T P A and f values within TDDFT/B3LYP, TDDFT/CAM-B3LYP and CC2 methods. In particular, we were interested in the effect of a number of diffuse and polarization functions on obtained values. According to Paterson et al., 85 a

6


Page 6 of 56

Page 7 of 56 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


number of diffuse functions as well as cardinality of the basis set may have significant impact on the excitation energy and TPA intensity values. To achieve our goals, we selected neutral and anionic chromophores naturally occuring in FPs such as, e.g. green, 1 cyan 8 and blue 8,102 fluorescent proteins, mBlueberry, 25 mFruits series 18 (Figure 1; see the Methodology section for details). We believe that our critical and systematic assessment of TDDFT functionals for one- and two-photon absorption properties for anionic and neutral FPs’ chromophores in vacuo will allow for a more reliable choice of methodology for more computationally demanding calculations of spectral properties employing electrostatic or polarizable embedding of the chromophore inside the protein environment. The paper is organized as follows. In section 2, we present methodology details including conversion of TPA transition moments [atomic units] into TPA cross-section [GM]. Next, in section 3, we present and discuss our results. Concluding remarks are gathered in section 4.

2 Methodology 2.1

Models and geometry

Figure 1 presents chromophores models investigated herein. As all chromophores are covalently bound to the FP amino acid chain in their natural environment, we saturated the cut bonds with hydrogen atoms. The green fluorescent protein (GFP) chromphore model is represented by both anionic and neutral forms, denoted as GFP-a and GFP-n, respectively. Such a choice is justified as GFP chromophore is believed to occur naturally in both forms in wild-type GFP 103–105 or solely in the anionic one in mutants, e.g. enhanced GFP (EGFP) 24 or yellow fluorescent protein (YFP). 24 As the chromophores of the cyan fluorescent protein (CFP) 8 and mBlueberry (BLB) 25 do not contain moieties that could be easily protonated/deprotonated, we describe them as neutral structures which to the best of our knowledge is not discussed in the literature. The orange fluorescent protein (OFP) and red fluorescent protein (RFP) chromophores are widely acknowledged to have deprotonated 7



Page 8 of 56

hydroxyl moiety of hydroxybenzyl group, i.e. they are anionic. The former is found in mOrange 18 and the latter in DsRed 106–109 and many of its derivatives, notably mFruits series of FPs. 18 The protonation state of chromophore found in blue fluorescent protein (BFP) 8,102 is possibly the most challenging to establish among investigated FPs’ chromophores as there are four possibilities of the imidazole ring protonation. According to theoretical work by Rubio and co-workers, 21 the anionic and cationic imidazole ring can be rejected to naturally occur in the protein environment, at least in the ground electronic state, leaving room for the neutral forms with one of the nitrogen atoms protonated. As it is not clear which structure is more realistic to be found in the BFP, we perform calculations for both of them, denoting them as BFP1 and BFP2 models. The ground-state geometry for each chromophore was obtained with MP2 110 method in def2TZVPP 111,112 basis set (Cartesian coordinates are provided in the Supporting Information) with Gaussian09d01 113 computational package.

2.2

Basis sets

To assess impact of the basis set size on the excitation energies and both OPA and TPA intensities, we have utilized correlation-consistent polarizable valence X-tuple-ζ (cc-pVXZ, where X=D, T, Q) 100 as well as their diffuse functions-augmented derivatives (aug-cc-pVXZ, where X=D, T, Q and X=5 for GFP-a) 101 and doubly-augmented one (d-aug-cc-pVDZ) basis sets. All of the listed basis sets were used for CC2, including auxiliary basis set for resolution of identity (RI) 114 approximation, as well as the TDDFT calculations utilizing B3LYP and CAM-B3LYP DFT functionals. For other DFT functionals (BLYP, B1LYP, BHandHLYP), we have used only aug-cc-pVDZ and aug-cc-pVTZ basis sets. For all basis sets, spherical functions were used. Our analysis of the basis set size impact was carried out for GFP-n and GFP-a only as representatives of the neutral and anionic FP chromophore, respectively. For other chromophores, calculations were performed only in aug-cc-pVDZ and aug-cc-pVTZ – see the Results and Discussion section.

8


Page 9 of 56 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


2.3

TDDFT calculations

The excitation energies, one-photon oscillator strengths and two-photon transition moments were calculated with the TDDFT methodology employing BLYP (0 %), 87,88 B3LYP 89 (20 %), B1LYP 90 (25 %), BHandHLYP 91 (50 %) and CAM-B3LYP (19 % at short range + 65 % at long range) 92 DFT functionals in Dalton 2015.0 115,116 suite of programmes. In parenthesis, we provide contribution of the Hartree-Fock exchange energy to the DFT functional formula. In the course of the TDDFT calculations, no core orbitals were kept frozen. In order to avoid deletion of molecular orbitals (MOs) due to linear dependencies, we have decreased the default threshold of eigenvalues from 10−6 to 10−10 . One-photon oscillator strengths are provided in length representation, while two-photon absorption cross-sections were calculated based on TPA of linearly polarized light with parallel polarization and two photons of equal energy using the quadratic response function formalism. Table S18 in Supporting Information provides number of excited states calculated with chosen DFT functionals and basis sets for GFP-n and GFP-a. For other chromophores, we calculated and analyzed 25 excited states. One needs to be aware that in the manifold of 25 excited states, described by single-reference TDDFT approach, some may display non-negligible multireference character leading to poor description of spectral properties. In course of all TDDFT calculations, we did not employ symmetry (C1 point group was used). Ground- and excited-state electron densities and dipole moments were calculated with Gaussian09-d01 113 programme for all basis sets and functionals. We calculate the dipole changes P 0.5 as ∆µ = ( a (µa,ES − µa,GS )2 ) . The a denotes Cartesian coordinates (a = x, y and z), while GS and ES stand for ground- and excited state, respectively.

2.4 CC2 calculations CC2 93 calculations of excitation energies, OPA oscillator strengths, ground- or excited-state electron densities as well as dipole moments were performed in Turbomole 6.6 117 computational package. To obtain the CC2 σ T P A values, we have utilized ricc2 module 67 which has recently been implemented in Turbomole 7.2 118 release, although we used pre-released version provided by the 9



Authors. The Authors 67 have introduced numerical Laplace transformation of the coupled-cluster ground-state amplitudes to avoid storage of large arrays, i.e. one needs to choose the threshold for Laplace decomposition. We have chosen it to be 10−6 and the same value for the solution of the response equations. Resolution-of-identity 114 approximation was employed for more efficient treatment of the two-electron integrals giving rise to the RI-CC2 method. Nanda and Krylov 84 stress that RI-approximated CCSD methodology gives the TPA cross-section values in an excellent agreement with the non-approximated approach 85 while computational cost is significantly reduced in the former case. As for TDDFT calculations, no core orbitals were frozen in course of the calculations and the one-photon oscillator strengths were calculated in length representation. To reduce computational cost by not calculating A“ states, we employ the Cs point group symmetry for all CC2 calculations and present results only for the A’ symmetry excited states. In fact, our preliminary CC2/aug-cc-pVDZ results reveal that 5 lowest-lying (∆E within a range of 3.79–5.54 eV) A“ states are dark for both OPA (f not exceeding 0.001) and TPA (σ T P A in range 0.00–0.65 GM) processes for GFP-n. Analogously, this is the case for GFP-a structure (∆E within 2.82–3.83 eV, σ T P A 0.05–2.43 GM and f 0.000–0.002). For GFP-n and GFP-a, number of calculated A‘ excited states is given in Table S18 in the Supporting Information. For other chromophores, always 8 lowest-lying excited states were calculated (except OFP – 10 states). The dipole moment changes were calculated according to the same protocol as for TDDFT results. Validity of single-reference CC2 results was checked by running D1 119 and D2 120 diagnostics (D1 < 0.12 and D2 < 0.24 for all investigated models and basis sets; see Tables S25 and S27 in the Supporting Information) which all fall within the trust regions (D1 < 0.15 and D2 < 0.25). 121 Moreover, the double-excitation character of the response state, as measured by %T2 diagnostic, 122 reveals a weight of double excitations within the recommended region (

Assessment of Functionals for TDDFT Calculations of One- and Two

Recommend Documents