Validity of Linear Response Theory for Time ... - ACS Publications

Although linear response approach enables a more convenient way to calculate Stokes shift, its validity is unpredictable. Usually, the photon excitati...
0 downloads 0 Views 1MB Size
Download PDF
Article pubs.acs.org/JPCB

Validity of Linear Response Theory for Time-Dependent Fluorescence in Staphylococcus Nuclease Tanping Li* 700 Choppin Hall, Chemistry Department of Louisiana State University, Baton Rouge, Louisiana 70803, United States ABSTRACT: We report a theoretical study on the relaxation dynamics of Staphylococcus nuclease following photon excitation. Molecular dynamics simulations were implemented for both ground state and excited state surface of Trp140. Over the course of 40 ns ground state simulation, a structural transition was observed from one isomeric protein configuration to another. No obvious isomerization process was exhibited in the excited state simulations. Using linear response theory and direct nonequilibrium approach, time-dependent Stokes shift, as well as its timeinfinity value, was evaluated for both ground state isomers. Comparison between these methods exhibits severe disagreement. The nonequilibrium simulations show the similar relaxation dynamics as the excited state linear response approach, whereas they severely differ from the ground state linear response evaluation. Further examination reveals that the isomerization process breaks the characteristics of the energy gap histogram from Gaussian statistics. We also made the comparison between the experimental and simulation results. The significant inertial decay, in theory, is absent in the experiment. The ground state linear response evaluation provides a better agreement with the experimental results among all three approaches.

■

INTRODUCTION Time-resolved fluorescence experiments1−5 have been applied to explore solvation dynamics by employing optical probe on solute molecules. In these experiments, an electrostatic perturbation is introduced by exciting the chromophore of solute from its electronic ground state to the excited state, leading to a geometry change and charge redistribution. Environmental solvent molecules must respond to this perturbation and the system will relax to another new equilibrium state. This process upon photon excitation, clearly a nonequilibrium process, can be monitored by measuring the time-dependent Stokes shift in femtosecond-resolved spectroscopy,6−13 and is applied to investigate solvation dynamics on the time scale of femtosecond, picosecond, or longer. Molecular dynamics (MD) simulations have been performed to study Stokes shift extensively in order to reveal the intrinsic molecular mechanism of solvation dynamics.14−24 MD simulations with use of direct nonequilibrium ensemble average have achieved qualitative agreement with the experimental observations.14,16,19 A second theoretical approach, linear response theory, hypothesizes that the nonequilibrium relaxation dynamics of a system upon photon perturbation is determined by its equilibrium fluctuations as driven by either ground or excited state potential surface.25−29 Linear response theory is important not only because it evaluates solvation correlation function based on the system’s equilibrium fluctuations, which is fundamental in many biological activities, but also for the great advantage over nonequilibrium ensemble average. To calculate the Stokes shift relaxing on the nanosecond (ns) time scale, nonequilibrium simulation requires hundreds, or even thousands of ns trajectories for convergence, which is computationally © 2014 American Chemical Society

expensive even for regularly sized system. Linear response approach only involves a single equilibrium trajectory with the same length in principle, although a practically longer one might be required. Although linear response approach enables a more convenient way to calculate Stokes shift, its validity is unpredictable. Usually, the photon excitation causes a noticeable charge redistribution and size change of solute, giving rise to the switch of solvation energy by several times of kBT, where kB is the Boltzmann constant and T is temperature. The response of the system falls into the region of nonlinear behavior. However, it has been reported in numerous systems that linear response approximation, although not strictly applicable, still provides a reasonable description of the actual nonequilibrium dynamics.26−28,30 By employing the intrinsic residue tryptophan (Trp) as local optical probe, recent MD simulations on protein hydration dynamics especially discovered its myriad successes.14,15,17 It has been revealed that Gaussian statistics plays a critical role to validate the linear response approach, even though the system undergoes significant perturbation and the configurations far away from equilibrations are sampled.31 On the other hand, it is not uncommon to find counter examples where linear response theory breaks.32−38 Revealing the failure of linear response theory is important to understand solvation dynamics and its role in biological activities. Moskun et al. reported that an abrupt liquid-structure change triggered by the rapid rotation leads to the switch of the molecular relaxation Received: July 2, 2014 Revised: October 8, 2014 Published: October 9, 2014 12952

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959

The Journal of Physical Chemistry B

■

from linear to nonlinear response.33 By studying a simple Lennard-Jones solute, Schwartz and co-workers found the breakdown of linear response as the change of shape or size on solute is involved.34−36 The investigations of Ladanyi and coworkers and Turi et al. pointed the nonlinear response to the specific differences in hydrogen bonding patterns between the initial and final state of solute.39,40 One of the common features in these cases is the different equilibrium geometry of peripheral solvent molecules upon the solute’s state. It was also pointed out that the failure of linear response theory ties to the non-Gaussian statistics of energy gap, a quantity closely correlated with Stokes shift.31,41 Chandler and co-workers observed a nonstationary Gaussian statistics of energy gap during the nonequilibrium process, whereas the instability was assigned to explain the failure of linear response theory.41 In this article, we report theoretical studies on solvation dynamics of the enzyme Staphylococcus nuclease (SNase) following photoexcitation of the only tryptophan residue, Trp140, which has been well studied by femtosecond spectroscopy experiments.13 Figure 1 shows the X-ray structure of the

Article

METHODS

Simulation Methods. Simulations for the protein were conducted with a double precision version of the GROMACS package43 and GROMACS96 force field,44 combined with the SPC/E water model. The nonbonded pair list was generated by using a 9 Å cutoff. Long-range electrostatic interactions were handled with the smoothed particle mesh Ewald (SPME) algorithm,45,46 with a real space cutoff length of 9 Å. The twin range cutoff scheme of 9 and 14 Å was applied for the LennardJones potential. All bond lengths were constrained by using the LINCS algorithm,47 allowing a 2 fs time step in the simulation. Periodic boundary conditions were implemented by using a truncated triclinic box of side length 61 Å and solvated with 4605 water molecules. The Nose-Hoover thermostat48−50 was used to maintain the system at 295 K. The initial configuration of SNase was taken from the crystal structure (PDB: 1SNO). Site charges of the indole chromophore in the ground state S0 came from the GROMOS96 force field; site charges of the La excited state of the indole ring, which is the fluorescing state of Trp, were obtained by applying the ab initio partial charge density differences calculated by Sobolewski and Domcke to the ground-state partial charges of the GROMOS96 force field as described in ref 51. Stokes Shift. Time-dependent Stokes shift is denoted as ΔEStokes(t) = ⟨ΔE(t) − ΔE(0)⟩, which is the shift of energy gap ΔE(t) relative to its value at time zero, the instant of photon excitation. For each configuration, ΔE(t) is denoted as the energy difference between the excited and ground state surface, where only the electrostatic interaction between indole and the rest of the system was involved in our simplified model. The angular brackets represent the ensemble average over the process upon photon probe. Note this nonequilibrium process requires the initial conditions sampled from the equilibrium fluctuation on the ground state surface followed by the relaxation dynamics driven by the excited state potential. Stokes shift can be divided into the contribution from a certain subsystem x by separating ⟨ΔE(t)⟩ into the corresponding component ⟨ΔEx(t)⟩, where x could be protein, water, or any portion of the system. In linear response approach, Stokes shift is estimated by evaluating the time correlation function ΔEStokes(t) = 1/kBT(⟨ΔE(t)ΔE(0)⟩ − ⟨ΔE(0)2⟩). In contrast to the nonequilibrium approach, the bracket indicates the average taken over equilibrium trajectories by either ground or excited state simulation. The decomposed contribution from a specific subsystem x is procured by applying the following equation: ΔEStokes(t) = 1/kBT(⟨ΔEx(t)ΔE(0)⟩ − ⟨ΔEx(0)ΔE(0)⟩). It is noted that a gas-phase transition between S0 and La of indole is involved for the energy gap ΔE(t), which is irrelevant for the calculation of Stokes shift. Energy Gap Histogram. Energy gap histograms Pg(ΔE) and Pe(ΔE) based on ground and excited state equilibrium simulations were applied to calculate the time-infinity Stokes shift. In the nonequilibrium approach, ΔEStokes(∞) is the difference of ⟨ΔE⟩ based on excited and ground state equilibrium fluctuations, where the ⟨ΔE⟩ can be achieved by applying the weighted average over the associated EGHs. In linear response theory, ΔEStokes(∞) = 1/kBT(⟨ΔE⟩2 − ⟨ΔE2⟩) = −(1/kBT)σ2, where σ2 is the variance of ground or excited state EGH. Furthermore, the energy gap histograms of ground and excited states are constrained with each other by the well-developed relationship Pe(ΔE) = e−β(ΔE−ΔA)Pg(ΔE),52,53 where ΔA supplies the associated free energy difference between the two surfaces.

Figure 1. Left: X-ray structure of Staphylococcus nuclease (PDB: 1SNO). Indole of Trp140 (yellow) is in close proximity to two positively charged residues K110, K133 (blue) and one negatively charged residue E129 (red) within the distance of 5 Å. Right: Local configuration at the periphery of Trp140. Residues E129 and K133 form a salt bridge with a distance of 3.5 Å.

protein, consisting of three α-helices and a five-strand β-barrel with a total of 149 amino acids.42 Trp140 has one edge exposed to the surface and inserts inside the protein to form a hydrophobic cluster at the C terminus, which was found to be essential for protein foldability, stability, and activity. Three charged residues K110, K133, and E129 exhibit within 5 Å of Trp140, whereas two neighboring residues E129 and K133 are closely coordinated with each other with a distance of 3.5 Å. By applying MD simulations, we examined the relaxation dynamics around Trp140 upon initial photoexcitation. Both linear response theory and direct nonequilibrium approach were performed to compute time-dependent Stokes shift and its time-infinity values. Comparison was made between these approaches. Energy gap histograms (EGH) sampled over ground and excited state simulations, as well as the coupling relationship between two histograms, were applied to analyze the origin of the discrepancy observed, which is the main focus of this paper. These studies further verify the role of non-Gaussian statistics, originating from protein isomerization, on the validity of linear response theory. We also made the comparison between the experimental measurements and the MD calculations. Finally, we present a unified picture that appears applicable to many other cases recently studied. 12953

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959

The Journal of Physical Chemistry B

Article

Figure 2. MD simulation of 40 ns performed on the ground state surface to generate (a) the time evolution of energy difference between the excited and ground state surface. A structural transition occurs at about 19 ns, and two configurations, defined as isomer 1 and isomer 2, were used to generate (b) the energy gap histograms and (c) the histogram sampled over the contribution from the pair E129−K133.

■

RESULTS Equilibrium Simulations. Ground state MD simulation of SNase, extending over 40 ns after initial equilibration, reveals distinct structural fluctuation occurring on the nanosecond time scale. Figure 2a, a plot of energy gap, shows a sudden change after 19 ns, indicating a structure transition that probably originates from a dynamical hopping between local substates. We refer to the structure of the first 19 ns as isomer 1, and the reminder of the trajectory as isomer 2. Energy gap histograms sampled over two isomeric structures exhibit an obvious shift, as displayed in Figure 2b. The histogram based on the joint contribution from E129 and K133 to energy gap shows a similar shift between two isomers (Figure 2c). Further examination on local structure reveals that neighboring residues E129 and K133 form a stable pair by salt bridge interaction with a distance of 3.5 Å. The pair demonstrates a significant displacement during the structural transition. The interaction with the indole ring is stronger in isomer 1 configurations than in isomer2, and is responsible for most of the energy shift in Figure 2a. On the other hand, two excited state MD trajectories were extended over 40 ns after the equilibration relaxation, with the initial configurations taken from isomer 1 or isomer 2 configurations. The resulting energy gap fluctuations are shown in Figure 3a. There is no obvious isomerization process observed until the end of the maximal length of our simulations. The corresponding EGH curves (Figure 3b) almost overlap with each other, implying that the two sets of equilibrium fluctuations sampled over excited state surface are equivalent, in spite of the distinct initial isomeric

configurations the simulations take. The photon excitation produces a noticeable charge redistribution of the indole, where the dipole moment increases from 2.4 D in the S0 to 4.96 D in the La state.15 The strong interaction of the excited state indole with the charged pair E129−K133 enhances the local rigidity around W140. Thus, the isomeric fluctuations exhibited on the ground state surface do not show on the excited state surface in SNase. Time-Dependent Stokes Shift. The relaxation dynamics upon photon excitation of Trp140 was examined by measuring time-dependent Stokes shift, as well as the specific contribution from protein, surrounding water molecules, and the pair E129− K133. For isomer 1, both linear response approach and direct nonequilibrium MD simulations were employed to calculate the Stokes shift up to 400 ps. For nonequilibrium simulations, 150 initial configurations were uniformly sampled over a 14 ns time interval during which the system is clearly in isomer 1 configurations, as indicated in Figure 2a. The ground state linear response was accumulated over the same time period. The excited state linear response was accumulated over one of the excited trajectories in Figure 3a whose initial configuration was taken from isomer 1 configurations. As shown in Figure 4, all of the calculations exhibit the relaxation dynamics on three time scales. An ultrafast decay on the femtosecond time scale was observed to lead to a severe inertial drop, and is followed by a picosecond decay originating from the local reorientational and translational motions of water and protein. At the slow time scale of tens of picoseconds, protein and water contribute in compensate manner as revealed in our published papers.14,15 Severe discrepancies are exhibited between different approaches. The nonequilibrium simulations demonstrate qualitatively similar relaxation dynamics to the excited state linear response, however seriously different from ground state linear response evaluations. Up to 400 ps, the total Stokes shifts by nonequilibrium approach, excited state and ground state linear response are −29.8, −32.2, and −22.6 kJ/mol, respectively. Closer examination reveals that the pair E129−K133 results in most of the discrepancy, which is −20.4 kJ/mol and −22.8 kJ/ mol for nonequilibrium simulations and excited state linear response, and only −3.4 kJ/mol for ground state linear response. Figure 5 shows the Stokes shift relaxation dynamics for the isomer 2 configurations, as well as the decomposed contribution from protein, water, and the pair E129−K133. For nonequilibrium simulations, 200 initial configurations were sampled over 15 ns interval designated as isomer 2 substate in Figure 2a. The ground state linear response was accumulated over the same time interval, and the excited state linear response was accumulated over the equilibrated excited trajectory in Figure 3a whose initial configuration was taken from isomer 2

Figure 3. Two MD simulations of 40 ns, initiated from isomer 1 or isomer 2 configurations, were performed on excited stte surface to generate (a) the time evolution of energy gap between the excited and ground state surface and (b) the resulting histogram. The results are presented by a straight line for isomer 1, and a straight line with circles for isomer 2. The fluctuation of isomer 2 is shifted to the time interval of 40−80 ns in panel a to make a visual clarification. No obvious isomeric process is exhibited. 12954

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959

The Journal of Physical Chemistry B

Article

Figure 4. Time-resolved Stokes shift (total) for isomer 1, as well as the specific contributions from protein, water, and the pair E129−K133 by (a) nonequilibrium (NE) simulations, (b) excited state linear response (LRE), and (c) ground state linear response (LRG) approaches. Nonequilibrium average and the excited state linear response approach show similar relaxation dynamics, which are different from the ground state linear response. Most of the discrepancy is attributed to the contribution of the pair E129−K133.

Figure 5. Time-resolved Stokes shift for isomer 2, as well as the specific contribution from protein, water, and the pair E129−K133 by (a) nonequilibrium simulations (NE), (b) excited state linear response (LRE), and (c) ground state linear response (LRG) approaches. Nonequilibrium average and excited state linear response show significant decay due to the contribution of E129−K133. The excited state linear responses for isomer 1 (black dotted line) and isomer 2 (black straight line) are in close proximity with each other in panel b.

Table 1. Time Infinity Stokes Shift Dataa

configurations. Similarly as in isomer 1, all three calculations show the inertial drop on the femtosecond time scale, a fast decay on the picosecond time scale, and a slow decay on the tens of picoseconds time scale. The nonequilibrium simulations demonstrate qualitatively similar relaxation dynamics as the excited state linear response, and deviate severely from the ground state linear response results. Up to 400 ps, the total Stokes shifts are −41.5 and −34.5 kJ/mol for nonequilibrium evaluation and excited state linear response, and −19.5 kJ/mol for ground state linear response. The pair E129−K133 contributes significantly in nonequilibrium simulations (−41.8 kJ/mol) and excited state linear response approach (−27.4 kJ/ mol), which is mild for ground state linear response (−9.1 kJ/ mol). For comparison purposes, two excited state linear response calculations of isomer 1 and isomer 2 are plotted together in Figure 5b. The curves are in close proximity to each other, which is attributed to the equivalent equilibrium fluctuations of two excited state simulations. Time-Infinity Stokes Shift. The time-infinity Stokes shifts of both nonequilibrium and linear response approaches are calculated by analyzing the EGH. The results are listed in Table 1. Figure 6 shows that the isomerization process during ground state simulation leads to an obvious shift between two EGH curves. Since two excited state simulations demonstrate almost equivalent equilibrium statistics and EGH curves, the EGH of isomer 1 is used in our calculation. It is noticed that the EGH curves of two ground state isomers are much narrower than the excited state one. Correspondingly, the linear response evaluations are less for ground state isomers than the excited state result, which are −22.6 kJ/mol (isomer 1), −17.4 kJ/mol

isomer 1 (0.45) isomer 2 (0.55) CGS

NE

LRG

LRE

−27.5 −38.1 −33.3

−22.6 −17.4 −32.5

−36.5

a Time infinity Stokes shift determined by ground state linear response (LRG), excited state linear response (LRE), and nonequilibriun (NE) average for isomer 1, isomer 2, and the combined ground state (CGS). All energies are in units of kJ/mol. By applying Pe(ΔE) = e−β(ΔE−ΔA)[c1Pg1(ΔE) + c2Pg2(ΔE)], the energy gap histograms of isomer 1 and isomer 2 were fitted with Pg1(ΔE) = (0.093 + 0.0045ΔE 2

2

+ 0.000062ΔE2)e−[(ΔE+12) /(1.2 × 10 )] and Pg2(ΔE) = (0.053 + 2

0.0055ΔE + 0.00032ΔE2)e−[(ΔE+1.0) /78] with the percentage c1 and c2 in parentheses.

(isomer 2), and −36.5 kJ/mol (excited state). Further examination on local structure reveals that the system is trapped in either isomer 1 or isomer 2 configuration by rigid alignment between the pair E129−K133 and the indole ring during ground state simulation, and is restricted from inspecting as many configurations as excited state simulation samples. Nonequilibrium Stokes shift is −27.5 kJ/mol for isomer 1 and −38.1 kJ/ mol for isomer 2. The value is less for isomer 1 than isomer 2, as demonstrated obviously in Figure 6. The difference is mostly attributed to the extensive structural rearrangement of the pair upon electric probe in isomer 2 configurations, as displayed in Figure 5a. It is interesting to revisit the ground state linear response and nonequilibrium approaches considering the fact that the ground 12955

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959

The Journal of Physical Chemistry B

Article

relaxation dynamics along similar time scales, although deviations exist between the various approaches. However, the magnitude of the Stokes shift is severely different between the experimental and MD observations. The experimental measurement is 10.1 kJ/mol, much smaller than all of the MD simulations. Further examination reveals that the significant inertial drop (less than 1 ps) in the MD simulations is absent in the experimental result. For all of the simulations, except the nonequilibrium approach in isomer 2, the inertial decays yield more than 50% of the total energy. For example, the total Stokes shift in the nonequilibrium simulation is 29.8 kJ/mol for isomer 1 by fitting Figure 4a, whereas the contribution from the inertial decay is 16.7 kJ/mol. Similarly for isomer 2, the inertial decay contributes 15.1 kJ/mol out of the total decay (41.5 kJ/mol) by fitting Figure 5a. These discrepancies would be resolved if the large inertial response seen in the simulations were removed, either by the improvement of the force field, such as by an incorporation of electronic polarization,54,55 or by the inclusion of electronically nonadiabatic effects that might suppress the inertial regime. Further work is required to achieve a quantitative agreement between theory and experiment. As all three approaches show significant inertial decay, the ground state linear response evaluation provides a better agreement with the experimental observations. By subtracting the ultrafast inertial contribution, the decay by the linear response approach is 6.2 kJ/mol for isomer 1, and 8.3 kJ/mol for isomer 2, which is comparable to the experimental result (10.1 kJ/mol). Furthermore, the associated time scales are also close to each other, with 5.8 and 115.2 ps for isomer 1, 4.2, and 78 ps for isomer 2, and 5.1 and 153 ps for experiment. On the other hand, the results from the excited state linear response and the nonequilibrium approaches deviate from the experiment by either the magnitude of the decay, even excluding the inertial drop, or the time constants.

Figure 6. Energy gap histograms of ground state isomers and excited state (ES) from the MD simulations (noisy) and the fitting (smooth). The time-infinity Stokes shift by ground state linear response, excited state linear response, and direct nonequilibrim are displayed by arrows for isomer 1, which is not shown for isomer 2 for visual clarification. The nonequilibrium Stokes shift is more in isomer 2 than isomer 1. The combined ground state energy gap histogram is shown by dashed line curve with squares.

state energy gap histogram should be sampled over both isomer 1 and isomer 2 configuration ensembles, which virtually provides the counterpart of experimental measurements. Hence we propose a combined ground state energy gap distribution by Pg(ΔE) = c1Pg1(ΔE) + c2Pg2(ΔE). Pgi(ΔE) represents the EGH of isomer i. ci is the related percentage with the restrain c1 + c2 = 1, since only two isomers exhibit up to the current length of our simulation. Pgi(ΔE) and ci were obtained by simultaneous fitting on restrain Pe(ΔE) = e−β(ΔE−ΔA)[c1Pg1(ΔE) + c2Pg2(ΔE)]. The extracted parameters are shown in Table 1. Pgi(ΔE) histograms were fitted with a polynomial multiplied Gaussian function. The percentage is 0.45 for isomer 1 and 0.55 for isomer 2. It should be noted that the parameters significantly depend on the sampling of the overlapped tails between curves Pg1, Pg2, and Pe, where statistics is notoriously difficult. By analyzing the combined ground state EGH shown as a dashed line curve with squares in Figure 6, time-infinity Stokes shift is recalculated as −32.5 kJ/ mol for ground state linear response, and −33.3 kJ/mol for nonequilibrium approach (Table 1). Comparison between Experiments and Simulations. By applying the femtosecond-resolved fluorescence technique, the experiments examined the relaxation dynamics of the Stokes shift.13 The comparisons with the MD simulation, including the direct nonequilibrium average and the ground and excited state linear response approaches of isomer 1 and isomer 2, are displayed in Table 2. Experiments found the biphasic relaxation dynamics with 5.1 and 153 ps. The MD simulations procured the

■

DISCUSSION AND CONCLUSIONS The linear response theory is widely applied to explore solvation dynamics. Its performance has been evaluated by using the comparison with the nonequilibrium simulations, as well as with the experimental results. The validity, in terms of the former, was reported to correlate with Gaussain statistics of energy gap in the literature. By examining time-dependent Stokes shift and its time-infinity value in SNase, our results point the observed discrepancy to an isomerization process on the ground state surface. Figure 7 elucidates a unified mechanism to further verify the role of Gaussian statistics, which appears applicable not only to the case in this context, but also to many systems recently studied.

Table 2. Stokes Shift Obtained from the Experiment and the MD Simulations, Including Nonequilibrium (NE) Average, Excited State Linear Response (LRE), and Ground State Linear Response (LRG) Evaluations for Isomer 1 and Isomer 2a isomer 1 cg; τg c1; τ1 c2; τ2 s∝

isomer 2

NE

LRE

LRG

NE

LRE

LRG

exptlb

16.7; 0.12 9.6; 9.9 3.5; 143.8 29.8

16.9; 0.13 4.2; 2.0 11.1; 22.1 32.2

16.4; 0.12 4.9; 5.8 1.3; 115.2 22.6

15.1; 0.14 16.5; 23.8 9.9; 236 41.5

20.5; 0.11 10.7; 6.4 3.3; 151 34.5

11.2; 0.12 3.8; 4.2 4.5; 78 19.5

4.7; 5.1 5.4; 153 10.1

2

Stokes shift was fitted with cge−(t/τg) + c1e−(t/τ1) + c2e−(t/τ2) − s∝ with cg + c1 + c2 − s∝ = 0. cg and τg represent the magnitude and time scale of a Gaussian-type inertial decay. The amplitudes and time constants are in units of kJ/mol and ps, respectively. bExperimental data were obtained from ref 13.

a

12956

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959

The Journal of Physical Chemistry B

Article

may also lead to the failure of linear response theory. This isomerization process in proteins, originating from the structural transition between local substates, has been observed extensively in the literature.57−60 The system is dynamically heterogeneous on the time scale of nanoseconds or even much longer. Thus, the employment of linear response theory on fluorescence spectroscopy might be restricted by the inherent local multiple states on energy landscape, by which the system could easily be trapped in one of the isomeric configurations during the MD simulations. It should be noted that isomerization process in proteins may not always result in the failure of linear response theory. Actually, linear response approach was found to agree well with nonequilibrium results in some systems, although the isomerization process apparently exhibits. In our previous publication on apomyoglobin, the linear response evaluations for the isomeric substates on the ground state surface provide a reasonable approximation on their nonequilibrium process.14 The mechanism can be illustrated by Figure 7c. Isomerization occurs on both the ground and excited state potential surface. The coupling relationship requires that the combined EGH of ground and excited state are constrained by ce1Pe1(ΔE) + ce2Pe2(ΔE) = e−β(ΔE−ΔA)[cg1Pg1(ΔE) + cg2Pg2(ΔE)]. Here, x = g, e represents ground or excited state surface, i = 1, 2 is the index of isomer, and cxi is the percentage of isomer i on surface x. For each ground state isomer i, Gaussian statistics on EGH is assumed to 2 2 be Pgi(ΔE) = 1/Nie−(ΔE−ΔEgi) /2σi . The excited state homology is required to maintain a shifted Gaussian form with the same 2 2 variance Pei(ΔE) = 1/Nie−(ΔE−ΔEei) /2σi , with the constrain ΔEei − ΔEgi = −σi2β. Therefore, the time-infinity Stokes shifts via both ground and excited state linear response approach agree with the nonequilibriun value for each isomer. It is interesting to notice that the ground and excited state EGH curves for isomer i preserve the same Gaussian variance. However, the percentage cgi and cei does not maintain the equality, and can be determined by fitting the coupling relationship. Although the ground state linear response approach may provide different results from the nonequilibrium simulations due to the conformational transition near the tryptophan, its performance should not be denied. It has been implemented extensively to examine the fluorescence Stokes shift, and thus helps to interpret the inherent mechanism of solvation dynamics from the molecular point of view. The direct comparisons with the experimental results show a reasonable agreement not only in the SNase, but also in the other systems in the literature.14,17,19 For example, the ground state linear response approach was reported to consistently reproduce the experimental findings of the fluorescent probe Hoechst 33258 for both cases: free in solution and bound to DNA.19 Giving that deviations commonly appear between the experiments and theories, the observed agreement might provide a valuable insight into explaining and then diminishing the discrepancies.

Figure 7. Mechanism on the validity of linear response theory to evaluate the time-infinity Stokes shift. (a) Linear response theory equal to the nonequilibrium average with Gaussian statistics. Ground state (GS) and excited state (ES) EGHs have the same variance. (b) Linear response theory fails for the ground state isomeric substate GS-isomer 1 and GS-isomer 2. (c) Linear response theory equals the nonequilibrium average, where the isomerization process exhibits in both the ground state (GS-isomer 1 and GS-isomer 2) and excited state surface (ESisomer 1 and ES-isomer 2). For each isomer, the ground and excited state counterpart retain the same variance, yet not the same percentage.

Figure 7a presents an ideal example where linear response approach equals that of the nonequilibrium process even as the system undergoes a significant perturbation like electric probe. Gaussian statistics is supposed in our picture, which is actually a reasonable assumption as illustrated in the work by Skinner and co-workers.56 With this assumption, ground state EGH has the 2 2 form Pg(ΔE) = 1/Ne−(ΔE−ΔEg) /2σ , where N is the normalization, σ2 is the variance and ΔEg defines the center of the curve. By applying the coupling relationship Pe(ΔE) = e−β(ΔE−ΔA)Pg(ΔE), excited EGH is found to follow a similar Gaussian formulism 2 2 Pe(ΔE) = 1/Ne−(ΔE−ΔEe) /2σ with the same variance as the ground state. The distance between the centers of two Gaussian curves, virtually representing the nonequilibrium Stokes shift of time infinity, is ΔEe − ΔEg = −σ2β. Therefore, both ground and excited state linear response theory correctly generate the timeinfinity Stokes shift, as long as energy gap histograms satisfy the Gaussian distribution. As far as the time-dependent Stokes shift is concerned, Laird et al. has proofed rigorously that excited state linear response evaluation exactly equals the nonequilibrium dynamics once the corresponding ΔE fluctuation obeys the Gaussian process, although it is not clear yet for ground state linear response evaluation.31 Figure 7b displays the breakdown of the equal variance of EGH characterized by Gaussian statistics. In contrast to Figure 7a, the isomerization process on the ground state surface turns a restriction on the sampling of energy gap for the individual substate, leading to the marked difference on the EGH variances between ground and excited state. Thus, the Stokes shifts by linear response approach, measuring the system’s response upon the perturbation based on its equilibrium fluctuation, are much less for ground state isomers than the excited state evaluation, and inconsistent with nonequilibrium values. The isomerization process on the excited state surface, although not observed yet,

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The calculations reported here were supported by HPC@LSU computing resources and the Ohio Supercomputer Center. The 12957

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959

The Journal of Physical Chemistry B

Article

(21) Furse, K. E.; Lindquist, B. A.; Corcelli, S. A. Solvation Dynamics of Hoechst 33258 in Water: An Equilibrium and Nonequilibrium Molecular Dynamics Study. J. Phys. Chem. B 2008, 112, 3231−3239. (22) Furse, K. E.; Corcelli, S. A. Effects of Long-Range Electrostatics on Time-Dependent Stokes Shift Calculations. J. Chem. Theory Comput. 2009, 5, 1959−1967. (23) Pal, S.; Maiti, P. K.; Bagchi, B.; Hynes, J. T. Multiple Time Scales in Solvation Dynamics of DNA in Aqueous Solution: The Role of Water, Counterions, and Cross-Correlations. J. Phys. Chem. B 2006, 110, 26396−26402. (24) Golosov, A. A.; Karplus, M. Probing Polar Solvation Dynamics in Proteins: A Molecular Dynamics Simulation Analysis. J. Phys. Chem. B 2007, 111, 1482−1490. (25) Schwartz, B. J.; Rossky, P. J. Aqueous Solvation Dynamics with a Quantum-Mechanical Solute - Computer-Simulation Studies of the Photoexcited Hydrated Electron. J. Chem. Phys. 1994, 101, 6902−6916. (26) Maroncelli, M.; Fleming, G. R. Computer-Simulation of the Dynamics of Aqueous Solvation. J. Chem. Phys. 1988, 89, 5044−5069. (27) Stratt, R. M.; Maroncelli, M. Nonreactive Dynamics in Solution: The Emerging Molecular View of Solvation Dynamics and Vibrational Relaxation. J. Phys. Chem. 1996, 100, 12981−12996. (28) Carter, E. A.; Hynes, J. T. Solvation Dynamics for an Ion-Pair in a Polar-Solvent - Time-Dependent Fluorescence and Photochemical Charge-Transfer. J. Chem. Phys. 1991, 94, 5961−5979. (29) Ladanyi, B. M.; Perng, B. C. Solvation Dynamics in DipolarQuadrupolar Mixtures: A Computer Simulation Study of Dipole Creation in Mixtures of Acetonitrile and Benzene. J. Phys. Chem. A 2002, 106, 6922−6934. (30) Bader, J. S.; Chandler, D. Computer-Simulation of Photochemically Induced Electron-Transfer. Chem. Phys. Lett. 1989, 157, 501−504. (31) Laird, B. B.; Thompson, W. H. On the Connection between Gaussian Statistics and Excited-State Linear Response for TimeDependent Fluorescence. J. Chem. Phys. 2007, 126, 211104. (32) Bragg, A. E.; Cavanagh, M. C.; Schwartz, B. J. Linear Response Breakdown in Solvation Dynamics Induced by Atomic ElectronTransfer Reactions. Science 2008, 321, 1817−1822. (33) Moskun, A. C.; Jailaubekov, A. E.; Bradforth, S. E.; Tao, G. H.; Stratt, R. M. Rotational Coherence and a Sudden Breakdown in Linear Response Seen in Room-Temperature Liquids. Science 2006, 311, 1907−1911. (34) Aherne, D.; Tran, V.; Schwartz, B. J. Nonlinear, Nonpolar Solvation Dynamics in Water: The Roles of Electrostriction and Solvent Translation in the Breakdown of Linear Response. J. Phys. Chem. B 2000, 104, 5382−5394. (35) Bedard-Hearn, M. J.; Larsen, R. E.; Schwartz, B. J. Understanding Nonequilibrium Solute and Solvent Motions through Molecular Projections: Computer Simulations of Solvation Dynamics in Liquid Tetrahydrofuran (Thf). J. Phys. Chem. B 2003, 107, 14464−14475. (36) Bedard-Hearn, M. J.; Larsen, R. E.; Schwartz, B. J. Hidden Breakdown of Linear Response: Projections of Molecular Motions in Nonequilibrium Simulations of Solvation Dynamics. J. Phys. Chem. A 2003, 107, 4773−4777. (37) Tao, G. H.; Stratt, R. M. The Molecular Origins of Nonlinear Response in Solute Energy Relaxation: The Example of High-Energy Rotational Relaxation. J. Chem. Phys. 2006, 125, 114501. (38) Vuilleumier, R.; Tay, K. A.; Jeanmairet, G.; Borgis, D.; Boutin, A. Extension of Marcus Picture for Electron Transfer Reactions with Large Solvation Changes. J. Am. Chem. Soc. 2012, 134, 2067−74. (39) Skaf, M. S.; Ladanyi, B. M. Molecular Dynamics Simulation of Solvation Dynamics in Methanol-Water Mixtures. J. Phys. Chem. 1996, 100, 18258−18268. (40) Turi, L.; Minary, P.; Rossky, P. J. Non-Linear Response and Hydrogen Bond Dynamics for Electron Solvation in Methanol. Chem. Phys. Lett. 2000, 316, 465−470. (41) Geissler, P. L.; Chandler, D. Importance Sampling and Theory of Nonequilibrium Solvation Dynamics in Water. J. Chem. Phys. 2000, 113, 9759−9765. (42) Truckses, D. M.; Somoza, J. R.; Prehoda, K. E.; Miller, S. C.; Markley, J. L. Coupling between Trans/Cis Proline Isomerization and

author thanks Prof. Sherwin Singer (OSU) for stimulating discussions and constructive comments.

■

REFERENCES

(1) Ware, W. R.; Lee, S. K.; Brant, G. J.; Chow, P. P. Nanosecond Time-Resolved Emission Spectroscopy - Spectral Shifts Due to SolventExcited Solute Relaxation. J. Chem. Phys. 1971, 54, 4729−4737. (2) Chakraba, S.; Ware, W. R. Nanosecond Time-Resolved Emission Spectroscopy of 1-Anilino-8-Naphthalene Sulfonate. J. Chem. Phys. 1971, 55, 5494−5498. (3) Castner, E. W.; Maroncelli, M.; Fleming, G. R. Subpicosecond Resolution Studies of Solvation Dynamics in Polar Aprotic and Alcohol Solvents. J. Chem. Phys. 1987, 86, 1090−1097. (4) Maroncelli, M.; Fleming, G. R. Picosecond Solvation Dynamics of Coumarin-153 - the Importance of Molecular Aspects of Solvation. J. Chem. Phys. 1987, 86, 6221−6239. (5) Nagarajan, V.; Brearley, A. M.; Kang, T. J.; Barbara, P. F. TimeResolved Spectroscopic Measurements on Microscopic Solvation Dynamics. J. Chem. Phys. 1987, 86, 3183−3196. (6) Zhong, D. P.; Pal, S. K.; Zhang, D. Q.; Chan, S. I.; Zewail, A. H. Femtosecond Dynamics of Rubredoxin: Tryptophan Solvation and Resonance Energy Transfer in the Protein. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 13−18. (7) Pal, S. K.; Peon, J.; Zewail, A. H. Biological Water at the Protein Surface: Dynamical Solvation Probed Directly with Femtosecond Resolution. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 1763−1768. (8) Peon, J.; Pal, S. K.; Zewail, A. H. Hydration at the Surface of the Protein Monellin: Dynamics with Femtosecond Resolution. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 10964−10969. (9) Zhao, L.; Pal, S. K.; Xia, T. B.; Zewail, A. H. Dynamics of Ordered Water in Interfacial Enzyme Recognition: Bovine Pancreatic Phospholipase a(2). Angew. Chem., Int. Ed. 2004, 43, 60−63. (10) Pal, S. K.; Peon, J.; Bagchi, B.; Zewail, A. H. Biological Water: Femtosecond Dynamics of Macromolecular Hydration. J. Phys. Chem. B 2002, 106, 12376−12395. (11) Qiu, W. H.; Zhang, L. Y.; Kao, Y. T.; Lu, W. Y.; Li, T. P.; Kim, J.; Sollenberger, G. M.; Wang, L. J.; Zhong, D. P. Ultrafast Hydration Dynamics in Melittin Folding and Aggregation: Helix Formation and Tetramer Self-Assembly. J. Phys. Chem. B 2005, 109, 16901−16910. (12) Qiu, W. H.; Zhang, L. Y.; Okobiah, O.; Yang, Y.; Wang, L. J.; Zhong, D. P.; Zewail, A. H. Ultrafast Solvation Dynamics of Human Serum Albumin: Correlations with Conformational Transitions and Site-Selected Recognition. J. Phys. Chem. B 2006, 110, 10540−10549. (13) Qiu, W. H.; Kao, Y. T.; Zhang, L. Y.; Yang, Y.; Wang, L. J.; Stites, W. E.; Zhong, D. P.; Zewail, A. H. Protein Surface Hydration Mapped by Site-Specific Mutations. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 13979− 13984. (14) Li, T. P.; Hassanali, A. A. P.; Kao, Y. T.; Zhong, D. P.; Singer, S. J. Hydration Dynamics and Time Scales of Coupled Water-Protein Fluctuations. J. Am. Chem. Soc. 2007, 129, 3376−3382. (15) Li, T. P.; Hassanali, A. A.; Singer, S. J. Origin of Slow Relaxation Following Photoexcitation of W7 in Myoglobin and the Dynamics of Its Hydration Layer. J. Phys. Chem. B 2008, 112, 16121−16134. (16) Hassanali, A. A.; Li, T. P.; Zhong, D. P.; Singer, S. J. A Molecular Dynamics Study of Lys-Trp-Lys: Structure and Dynamics in Solution Following Photoexcitation. J. Phys. Chem. B 2006, 110, 10497−10508. (17) Nilsson, L.; Halle, B. Molecular Origin of Time-Dependent Fluorescence Shifts in Proteins. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 13867−13872. (18) Sen, S.; Andreatta, D.; Ponomarev, S. Y.; Beveridge, D. L.; Berg, M. A. Dynamics of Water and Ions near DNA: Comparison of Simulation to Time-Resolved Stokes-Shift Experiments. J. Am. Chem. Soc. 2009, 131, 1724−1735. (19) Furse, K. E.; Corcelli, S. A. The Dynamics of Water at DNA Interfaces: Computational Studies of Hoechst 33258 Bound to DNA. J. Am. Chem. Soc. 2008, 130, 13103−13109. (20) Furse, K. E.; Corcelli, S. A. Molecular Dynamics Simulations of DNA Solvation Dynamics. J. Phys. Chem. Lett. 2010, 1, 1813−1820. 12958

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959

The Journal of Physical Chemistry B

Article

Protein Stability in Staphylococcal Nuclease. Protein Sci. 1996, 5, 1907− 1916. (43) Lindahl, E.; Hess, B.; Van Der Spoel, D. Gromacs 3.0: A Package for Molecular Simulation and Trajectory Analysis. J. Mol. Model. 2001, 7, 306−317. (44) Scott, W. R. P.; Hunenberger, P. H.; Tironi, I. G.; Mark, A. E.; Billeter, S. R.; Fennen, J.; Torda, A. E.; Huber, T.; Kruger, P.; van Gunsteren, W. F. The Gromos Biomolecular Simulation Program Package. J. Phys. Chem. A 1999, 103, 3596−3607. (45) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald - an N.Log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (46) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577−8593. (47) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. Lincs: A Linear Constraint Solver for Molecular Simulations. J. Comput. Chem. 1997, 18, 1463−1472. (48) Nose, S. A Molecular-Dynamics Method for Simulations in the Canonical Ensemble. Mol. Phys. 1984, 52, 255−268. (49) Nose, S. A Unified Formulation of the Constant Temperature Molecular-Dynamics Methods. J. Chem. Phys. 1984, 81, 511−519. (50) Hoover, W. G. Canonical Dynamics - Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695−1697. (51) Sobolewski, A. L.; Domcke, W. Ab Initio Investigations on the Photophysics of Indole. Chem. Phys. Lett. 1999, 315, 293−298. (52) Bennett, C. H. Efficient Estimation of Free-Energy Differences from Monte-Carlo Data. J. Comput. Phys. 1976, 22, 245−268. (53) Tachiya, M. Relation between the Electron-Transfer Rate and the Free-Energy Change of Reaction. J. Phys. Chem. 1989, 93, 7050−7052. (54) Bader, J. S.; Berne, B. J. Solvation Energies and Electronic Spectra in Polar, Polarizable Media: Simulation Tests of Dielectric Continuum Theory. J. Chem. Phys. 1996, 104, 1293−1308. (55) Kumar, P. V.; Maroncelli, M. Polar Solvation Dynamics of Polyatomic Solutes - Simulation Studies in Acetonitrile and Methanol. J. Chem. Phys. 1995, 103, 3038−3060. (56) Stephens, M. D.; Saven, J. G.; Skinner, J. L. Molecular Theory of Electronic Spectroscopy in Nonpolar Fluids: Ultrafast Solvation Dynamics and Absorption and Emission Line Shapes. J. Chem. Phys. 1997, 106, 2129−2144. (57) Frauenfelder, H.; Parak, F.; Young, R. D. Conformational Substates in Proteins. Annu. Rev. Biophys. Biophys. Chem. 1988, 17, 451− 479. (58) Frauenfelder, H.; Alberding, N. A.; Ansari, A.; Braunstein, D.; Cowen, B. R.; Hong, M. K.; Iben, I. E. T.; Johnson, J. B.; Luck, S.; Marden, M. C.; et al. Proteins and Pressure. J. Phys. Chem. 1990, 94, 1024−1037. (59) Frauenfelder, H.; Sligar, S. G.; Wolynes, P. G. The Energy Landscapes and Motions of Proteins. Science 1991, 254, 1598−1603. (60) Young, R. D.; Frauenfelder, H.; Johnson, J. B.; Lamb, D. C.; Nienhaus, G. U.; Philipp, R.; Scholl, R. Time-Dependence and Temperature-Dependence of Large-Scale Conformational Transitions in Myoglobin. Chem. Phys. 1991, 158, 315−327.

12959

dx.doi.org/10.1021/jp506599d | J. Phys. Chem. B 2014, 118, 12952−12959