Article pubs.acs.org/JPCB

Monte Carlo Modeling of Spectral Diffusion Employing Multiwell Protein Energy Landscapes: Application to Pigment−Protein Complexes Involved in Photosynthesis Mehdi Najafi and Valter Zazubovich* Department of Physics, Concordia University, Montreal H4B 1R6, Quebec, Canada S Supporting Information *

ABSTRACT: We are reporting development and initial applications of the light-induced and thermally induced spectral diffusion modeling software, covering nonphotochemical spectral hole burning (NPHB), hole recovery, and single-molecule spectroscopy and involving random generation of the multiwell protein energy landscapes. The model includes tunneling and activated barrier-hopping in both ground and excited states of a protein− chromophore system. Evolution of such a system is predicted by solving the system of rate equations. Using the barrier parameters from the range typical for the energy landscapes of the pigment−protein complexes involved in photosynthesis, we (a) show that realistic cooling of the sample must result in proteins quite far from thermodynamic equilibrium, (b) demonstrate hole evolution in the cases of burning, fixed-temperature recovery and thermocycling that mostly agrees with the experiment and modeling based on the NPHB master equation, and (c) explore the effects of different protein energy landscapes on the antihole shape. Introducing the multiwell energy landscapes and starting the hole burning experiments in realistic nonequilibrium conditions are not sufficient to explain all experimental observations even qualitatively. Therefore, for instance, one is required to invoke the modified NPHB mechanism where a complex interplay of several small conformational changes is poising the energy landscape of the pigment−protein system for downhill tunneling.

1. INTRODUCTION Spectral diffusion in pigment−protein complexes involved in photosynthesis1−12 and in other protein systems13−20 at low temperatures, and the underlying features of the protein energy landscapes have been a subject of numerous studies, including both theoretical and experimental ones. One of the questions that can be answered by these studies is to what extent proteins fit the frameworks developed for amorphous solids, such as glasses.4,5,14,15,20 The detailed information about the protein energy landscapes is also required to understand folding of the proteins as well as conformational changes involved in their functioning. The overall shape of the protein energy landscape is represented by the so-called “folding funnel”.21 A typical free energy difference between the unfolded and native states is on the order of 104 cm−1. However, the native state is not perfectly unique, and there are usually multiple conformational substates at the bottom of the funnel that are separated by much lower barriers. It also has been demonstrated that protein energy landscapes possess hierarchal character.3,14,15,18,22,23 Probing energy landscapes using optical methods involves either single molecule/single protein complex spectroscopy (SPCS),3−11 or subensemble techniques such as spectral hole burning (SHB), and particularly its nonphotochemical variety (NPHB).11,12,14,15,17,18,20 Both the highest spectral resolution and convenient slowing down of protein dynamics are achieved at cryogenic temperatures. Although the whole folding funnel in all its complexity cannot be probed at such temperatures, the © 2015 American Chemical Society

information that can be obtained pertains to the roughness of the protein energy landscape, likely on more than one hierarchal tiers.3 The landscape roughness affects the folding rates, for example. Protein energy landscapes can be calculated from the structural data in the course of molecular dynamics simulations, although the landscapes per se are infrequently reported explicitly. The parameters of the barriers on the energy landscapes obtained by optical methods are comparable to those featured in molecular dynamics simulations on a number of simple proteins.24,25 The reports on molecular dynamics simulations in photosynthesis-related systems are starting to emerge, but so far they did not focus on spectral diffusion and did not explicitly present the energy landscapes of the chromophore−protein system.26,27 On the other hand, these complexes are among the proteins most thoroughly studied by optical methods.1−12,28−31 In these complexes the pigments, serving as sensitive probes to their local environment, are built into protein by nature, in a broad variety of local environments, and one does not need to employ chemical manipulations or genetic engineering to artificially introduce the fluorescent probes into the protein. Received: March 22, 2015 Revised: May 20, 2015 Published: May 28, 2015 7911

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B

Figure 1. Representative energy landscapes of a pigment−protein system for energy landscapes dominated by randomness (A) and dominated by a parabola (B). The energies of the bottoms of the wells in frame A are also subject to distribution. The excited state landscapes are shifted up by an amount picked from a Gaussian distribution with the mean of 14560 cm−1 (686.8 nm) and fwhm of 180 cm−1. Excited state landscapes represented by the red curves are copies of the ground state landscapes (black) shifted and divided by a factor of 5.25, see text. Excited state landscapes represented by the green curves are copies of the ground state landscapes (blue) shifted and divided by a factor of 5.25. Also defined are the barrier height V, asymmetry Δ and the shift d along the generalized coordinate, associated with the conformational change. See text for details and distribution parameters.

change. As long as one is dealing with just one TLS per pigment, keeping track of the spectral shifts during both burning and recovery is possible even with relatively limited computational resources.35,36 However, this simple scenario is not always realistic even in low-temperature glasses where various departures from the classical TLS model were reported.34−39 SPCS results in pigment−protein complexes clearly demonstrate that the number of available conformational substates is larger than two per pigment molecule.3,5 On the other hand, this number must be limited. The evidence comes not only from SPCS, but also from the analysis of the spectral hole recovery.12,31 Dependence of the hole recovery (HR) rate on the fractional depth of the originally burnt hole is a clear signature of a significant degree of spectral memory. The functions describing the distributions of spectral shifts (“antihole” or “photoproduct” functions) in both glasses35,36 and proteins11,12 are phenomenological and not straightforwardly related to the features of the energy landscape. Finally, any HB model based on the master equation11,12,35,36 is applicable only to the ensembles of the pigment−protein systems and cannot be directly adapted to the interpretation of the SPCS results. Inclusion of the multiwell landscapes into the NPHB master equation results in a rapid increase of computational complexity.37 Here we will utilize an alternative approach to the NPHB modeling, using Monte Carlo simulation methods and the analysis of the evolution of the probability distributions that is based on solving the rate equations. Averaging the results for a large number of systems allows for modeling the evolution of the area or depth of the spectral hole as well as the antihole shape (but, at this stage, not yet the hole width) upon burning, fixed-temperature recovery and in thermocycling experiments. We will also address a related important issue, likely of interest to a broader audience than NPHB per se, namely how close the systems explored in the typical low-temperature spectroscopy experiments are to thermodynamic equilibrium, and the implications of performing SPCS and NPHB measurements on the far from equilibrium samples.

Although significant progress has been achieved, full understanding of the spectral diffusion phenomena is elusive. Multiple pieces of the puzzle are available, but the whole picture is yet to emerge. For instance, it is a matter of debate if tunneling or barrier-hopping or both are behind the LH2 lowtemperature SPCS observations.6 The onset of the barrierhopping in the excited state of the pigment−protein systems was observed at temperatures as low as 10−13 K in CP4312 and Cytochrome b6f31 using NPHB. It appears that it was implicitly assumed in many SPCS papers that observed spectral diffusion is thermally induced and is occurring anyway, whether one observes it by means of SPCS (i.e., frequently excites the chromophore electronically) or not. However, analysis of the SHB yields and photon budgets of the SPCS experiments indicates that the majority of spectral shifts must be lightinduced, i.e., represent nonphotochemical hole burning (NPHB) on a single-molecule level. SPCS and NPHB sometimes give contradictory information on the distributions of the spectral line shifts for the same system.11 Additional hole recovery channels, whose existence may contradict the presence of spectral memory (see below), are required to explain the results of thermocycling experiments in CP43 12 and Cytochrome b6f31 and so on. Here we are going to explore, by means of computer modeling, how certain realistic energy landscapes, more complex than a simple two-level system32,33(TLS), are expected to affect the SHB and some SPCS observations. In other words, we are extending the TLS-based NPHB model to multiwell systems. In general, NPHB involves the rearrangement of the immediate environment of the pigment molecule triggered by optical excitation. The original NPHB model assumed that the pigment is in interaction with just one structural element of the amorphous solid (or a protein) capable of conformational changes between two nearly identical substates (TLS) and that its spectral line is shifting between two frequencies.34−36 This scenario also implies the presence of perfect “spectral memory”the hole recovers as burnt pigment−protein systems return to the preburn conformation. The transition energy of a pigment is changing by several wavenumbers to several tens of wavenumbers as a result of the conformational 7912

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B

or barrier-hopping. (This value is more realistic than 7.6 × 1012 s−1 from,11,12,35,36,40 that is continued to be used for the sake of easier comparison of the parameters between different amorphous solids and proteins. Fitting the same experimental HGK curves employing Ω0 = 1012 s−1 instead of Ω0 = 7.6 ×1012 s−1 simply requires decreasing the mean of the λ-distribution by 1.04.) Hopping rates are Ω0 exp(−ΔE/kBT), where ΔE is the energy difference between the top of the barrier and the bottom of the well from which the hop is initiated. For the sake of simplicity, we assume that both m and d = 0.1 nm are constant, and md2 = 1.4 × 10−46 kg·m2. In other words, the spacing between neighboring wells along the generalized coordinate is constant, and wells of the excited state landscape are located at the same values of the generalized coordinate as the wells of the ground state landscape. The equal spacing assumption is based on an argument that similar entities (e.g., protons7,31 or, alternatively, small protein side groups12) are likely involved in conformational changes on a particular tier of the protein energy landscape. Invoking distributions of m and d values will further affect the λ-distribution and likely make it closer to a Gaussian even if the barrier height V distribution per se is not Gaussian. The simplification involving no shift between the ground and excited states along the generalized coordinate is fairly standard in NPHB modeling;35,36 additionally, the absence of shifts between excited and ground states is consistent with us not including the linear electron−phonon coupling in the current version of the model (quadratic electron−phonon coupling is involved in phonon scattering responsible for phonon-assisted tunneling). Two types of the energy landscapes were considered. In the first model, the barrier heights and depths of the wells are both simply drawn from the (separate) Gaussian distributions (in the following, “flat energy landscape”), Figure 1A. In the second model, both ground and excited state landscapes are composed of a parabola (representing the bottom of a well on a higherorder tier of the protein energy landscape3,14,15,18,22,23) and a collection of randomly generated barriers (in the following, “parabolic energy landscape”), Figure 1B. Each frame of Figure 1 depicts two examples of the energy landscapes in the ground (black, blue) and excited (red, green) states of the chromophore−protein system. The parameters of the barrier distributions were chosen to yield λ-distribution with the parameters resembling those reported for photosynthetic pigment−protein complexes, i.e., λ0 = 10.7 and σλ = 1.25 in the excited state,11,12,40 where λ0 and σλ are the mean and the standard deviation of the λ-distribution. The magnitudes of the barriers in the excited and ground states are correlated. Barriers are 5.25 times higher in the ground state than in the excited state, which results in a 2.3-fold difference in the mean of the λ−distribution. The numbers here and below originate from refs 12 and 31; barriers in the excited state must be significantly lower than in the ground state, as noticeable NPHB occurs despite nanosecond excited state lifetimes, while the recovery of a spectral hole takes hours or even days at cryogenic temperatures. In the “flat” model, the Gaussian distribution of ground state barriers is peaked at 1187 cm−1 and has the width of 283 cm−1. The distribution for the bottoms of the wells is a Gaussian peaked at 20 cm−1 and with the width of 10 cm−1. Therefore, the average difference of the depths of the neighboring wells (the equivalent of the TLS asymmetry) is about 10 cm−1. The distribution of a sum of two independent normally distributed parameters is also normally distributed. This barrier distribution corresponds to the ground state λ-

2. MODEL The model described here is predicting the evolution of the distributions of probability to find a spectral line (zero-phonon line, ZPL) at various wavelengths (or frequencies) upon lightand thermally induced spectral diffusion. These distributions are ideologically similar to the first cumulant distributions deduced from the SPCS data.5 Each {pigment inside protein}in the following, a “system”is described with two correlated energy landscapes, one for the excited state of the pigment in interaction with the environment and one for the ground state (see Figure 1). The landscapes for each of M systems are randomly generated, drawing barrier heights and well depths from the (Gaussian-shaped12) distributions. At each step of the evolution of a system, its state is described by a set of probabilities to be in a given well, Pi,a where a = e or g for ground and excited states, respectively, and i = 0···N − 1, where N is the total number of wells in the excited or the ground states of the system. ΣPi,a = 1. The energy difference ℏωB between the baselines for the ground and excited states (red arrow) is subject to a separate Gaussian distribution, with the width equal to the typical inhomogeneous bandwidth in photosynthetic complexes, i.e., about 100−200 cm−1, see Figure 1. Hole burning and subsequent hole recovery involve changes to the set of probabilities Pi,a. The evolution of this set of probabilities is determined by solving the rate matrix, containing the rates of transitions between different wells of the protein energy landscape. Technical issues related to solving rate matrices numerically are addressed in the Supporting Information. Note that in this half-continuous half-discrete approach each system in fact represents the whole subclass of similar systems, which, despite featuring identical energy landscapes, may be in different conformational substates. Both phonon-assisted tunneling and barrier-hopping are included. Tunneling rates are Ω0 exp(−2λ)(n(Δ) + 1) for tunneling involving transition to a deeper well and an increase of the energy of the scattered phonon (“downhill tunneling”) and Ω0 exp(−2λ)n(Δ) for “uphill” tunneling; λ = d(2mV)1/2/ℏ, where λ, d, and V are tunneling parameter, displacement along the generalized coordinate, and barrier height, respectively, and m is the effective mass of the entity rearranging during conformational change. n(Δ) = 1/(exp(Δ/kBT) − 1) is the population number for phonons being scattered, Δ is the asymmetry, the energy difference between the neighboring wells; see Figure 1. We define it as always positive; the sign of it is taken care of when we determine if given act of tunneling is uphill or downhill. The above expression for λ, widely used in papers addressing spectral diffusion and other low-temperature properties of amorphous solids,32−36 implies that the barriers are rectangular. In the case of parabolic barriers, λ is proportional to the barrier height. Exact shape of the barriers and the scaling between V and λ are undoubtedly important for precise quantitative modeling of the experimental data. However, for Gaussian barrier distributions and rectangular barriers employed here the λ-distribution is still fairly close to Gaussian, i.e. it is a product of a Gaussian term and a term linear in λ.12 From the viewpoint of the current manuscript that does not include linear electron−phonon coupling and does not attempt exact quantitative fitting of the experimental spectra, the issue of the barrier shape is not critically important. Ω0 = 1012 s−1 is the attempt frequency, determining how often the system is approaching the barrier it may cross by tunneling 7913

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B distribution with λ0 = 24.6 and σλ = 2.9. In the case of the parabolic model, there was no additional distribution for the bottoms of the wells, and the mean and standard deviation values of the Gaussian barrier distribution were 1207 and 284 cm−1. The coefficient of the parabola is also assumed to be smaller in the excited state, based on the argument that the barriers have to be lower in the excited state also on the next, higher hierarchal tier of the energy landscape. The ground state parabola coefficient kg ≈ 5000 cm−1/nm2 can be estimated based on the average antihole shift with respect to the hole, which in turn is close to the most likely ground state asymmetry ∼10 cm−1, i.e., the energy difference between the deepest, most populated wells. The probability distributions at the start of the lowtemperature optical experiment are a result of the cooling process. In most cases the sample is relatively slowly (within an hour or two) cooled down from room temperature to 5 K in the dark. To simulate this process, we start with the Boltzmann distribution (Pi,g ∼ exp(−Ei,g/kBT) for 300 K, where Ei,g is the bottom of the ith well of the ground state landscape). The preburn distribution Pi,g at 5 K is found via multistep procedure imitating gradual cooling. Once the preburn distribution Pi,g is established, we start simulating the hole burning. For every system a well is chosen so that the energy difference between ground and excited states for this well matches the fixed burn frequency (same for all systems) within 30 GHz (1 cm−1) or 10 GHz. This well may be different for different individual systems, even in the parabolic landscape model, as the energy difference between the bottoms of the ground and excited states is subject to a distribution. Let us assume that in a particular system the resonant well is the well number j, and the respective probability is Pj,g,before. This is the probability that given system will be in resonance with the narrow-band laser excitation to begin with. As a result of absorption of a photon, the system finds itself in the excited electronic state. The evolution in the excited state is started with Pj,e = 1, and all other Pi≠j,e = 0. System is given 3.5 ns (excited state lifetime) to evolve in the excited state. Rate equations are used to determine the respective time constants and the final excited state probabilities Pi,e. The system then returns to the ground state. All evolved Pi.e are multiplied by Pj,g,before (see above) and results are added to Pi,g for all i, with one exception: the Pj.e is multiplied by preburn Pj.g,before and placed instead of the preburn Pj.g,before. This procedure results in a ground-state distribution of Pi,g still normalized to one, but modified compared to the preburn distribution. The probability to find the system in resonance with the laser (in well j) has been reduced, the probabilities to find the system in other wells have increased somewhat. This modified ground-state probability distribution can be used as (i) the initial state for modeling of subsequent acts of burning, (ii) a starting point for the modeling of recovery, at burn temperature or at different temperatures, or (iii) for calculating the hole spectra. Additional hole burning involves repeating the single burning step described above, but with the initial ground-state distribution being the result of the previous burning step. The Pj,g after the previous step becomes the new Pj,g,before. Note that the probability of a photon being absorbed is determined by absorption cross section and the actual photon flux that are usually too small to supply a photon to every system every several nanoseconds. Thus, in reality there is a delay between two consecutive acts of burning. During this delay, the system can experience some recovery due to relaxation in the ground state. Rigorous modeling of hole-

growth kinetics (HGK) in this context involves a sequence of alternating acts of burning and fixed-temperature recovery. However, as shown in the Supporting Information, this correction may be ignored at 5 K and for realistic illumination intensities. In general, hole recovery at fixed temperature (either in-between the acts of burning, or after burning is finished) starts with postburn population probability distribution Pi,g,after, i = 0···N−1 and the system is allowed to evolve for some time, in ground state (in the dark) with respective rates. In the case of thermocycling,12,14,31 one has to simulate multiple relatively short recovery steps sequentially, according to the temperature change profile resembling those actually employed in the experiment.12,31 In both cases, Pi,g i = 0···N − 1 after one step serves as the initial condition for the next step. The hole spectrum is calculated in the following way: For every system separately, the difference between the preburn and postburn ground state distributions is converted into the difference in absorption spectra. For resonant well j (j obviously may differ from system to system), the Pj,g, af ter − Pj,g,before is negative, which corresponds to the resonant hole. For other wells Pi,g, af ter − Pi,g,before is positive, which corresponds to the antihole. When the results for all systems are being summarized and binned, they are added with the weight equal to preburn Pj,g,before (before any burning) for this particular system. The results of our hole burning and hole recovery simulations are shown in the following sections for an ensemble consisting of 5000 systems. Such simulations could be performed in reasonable time using a PC with an i7 processor and 8 GB of RAM. 50000 systems were employed in the simulations of the absorption spectra. It should be noted that the real number of pigments in a typical bulk sample is on the order of 1015 per approximately 40 mm3 illuminated sample volume. However, since for each system we are considering evolution of the probability distributions rather than discrete jumps, one could argue that each of the 5000 systems in turn represents a whole subclass of the real pigment−protein systems with similar energy landscapes but occupying different substates on that landscape. On the other hand, vast majority of the available 1015 systems do not possess a well in resonance with the burn frequency (see Discussion on the ensemble absorption spectra below), and would not participate in hole burning anyway.

3. RESULTS 3.1. Cooling the Sample Down. The first question which we want to address is if the samples are in thermal equilibrium prior to a typical SHB or SPCS experiment, given the landscape parameters matching the results from12 and reviewed above. The preburn condition of a system (equilibrium vs lack thereof) must affect the burning and recovery processes, as well as SPCS observations. For example, the out-of equilibrium system may exhibit, on average, higher SHB rate than an equilibrated system, since the probability to be in the deepest well will be smaller and the probability of burning via downhill tunneling will be higher. But if the system is not in equilibrium after cooling downhow far is it from the equilibrium? To answer these questions the cooling process was simulated by first assigning the Boltzmann distribution at room temperature (RT) to each system in the ensemble. It appears that the initial distribution is not very important since above 100 K the system with barriers described above equilibrates very quickly. However, for realistic cooling speeds (≈300 K/2 h), both the 7914

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B

Figure 2. Representative probability distributions Pwell number for flat (A) and parabolic (B) energy landscape for a single system. Red/blue curves: practically indistinguishable distributions for 5 and 25 K resulting from gradual cooling. Light green curve: result of gradual cooling to 50 K. Magenta: result of the gradual cooling to 100 K. Light blue: 200 K. Dark green: indistinguishable distribution after 5 min at 300 K and Boltzmann distribution at 300 K. Black: Boltzmann distribution at 5 K. The histogram in frame B represents the population probability distribution in the ensemble of 5000 systems resulting from realistic cooling.

whole ensemble and the individual systems find themselves trapped fairly far from equilibrium at low temperatures (roughly, below 50 K) due to very fast decrease of the barrier-hopping rate. In Figure 2, the results of simulating cooling of a single system in 5 K steps within 2 h (2 min per step) for both types of energy landscapes are compared with Boltzmann distributions at RT and 5 K (the burn temperature). The higher-temperature probability distribution in the case of a flat energy landscape (Figure 2A), dominated by randomness, is almost uniform. Since the relaxation rates become much smaller as the system cools down, the probability distribution stays rather close to uniform (populating deeper wells becomes only slightly more probable). Note that the deepest well could be located anywhere in the flat landscape model. In the example shown in Figure 2A the deepest wells are #9, #1, and #20, featuring the highest probability spikes in the case of the 5 K equilibrium. This equilibrium is not reached via realistic cooling. The average population of the deepest well in the case of realistic cooling is less than 8%. However, the population of the deepest well is not a useful characteristic of a system in the flat model. Instead, we can consider the proportion of the wells that have both nearest neighbor wells shallower than them. In the case of a parabolic energy landscape (Figure 2B), most of the probability is always distributed within just several wells near the bottom of the parabola. As the sample is cooling down, the probability to find the system in the deepest well will increase. However, the lowtemperature population distributions resulting from the realistic cooling process still deviate significantly from the 5 K Boltzmann distribution. Importantly, large fraction of the systems will be trapped in the wells other than the deepest one, and would exhibit downhill tunneling and effective SHB down to the lowest temperatures. Summarizing, both the individual systems and the whole ensemble end up far away from equilibrium after a reasonable cooling procedure in either energy landscape model. Also, even allowing the ensemble a long time (couple of years) to relax at low (burn) temperature still does not result in a perfect equilibration (not shown). Freezing of the nonequilibrium situation will likely be even more pronounced when the sample is quickly plunged into liquid helium. In refs 15, 17, 18, and 20, the effects of waiting time prior to burning on the small-step spectral diffusion

(spectral hole broadening) were explored. Here we will demonstrate that initial conditions must affect the hole growth and hole recovery as well. Note that the degree of equilibration does have an effect also on the whole absorption spectrum of the ensemble. Figure 3 compares the absorption spectra which

Figure 3. 5 K absorption spectra of 50 000-system ensembles in case of equilibrium probability distributions (black) and nonequilibrium distribution resulting from realistic cooling from 300 to 5 K within 2 h (red); flat, randomness-dominated energy landscape. Arrow indicates the burn frequency in subsequent simulations.

would follow from our model in the case of equilibrium situation and in the case of a nonequilibrium situation resulting from realistic cooling. Both spectra appear nearly Gaussian (again we point out that linear electron−phonon coupling is not included in this model), although the spectrum for the equilibrium situation is somewhat blue-shifted. This is a manifestation of the higher likelihood, in the equilibrium case, to find the systems in the deepest wells of their respective ground states, that is positively correlated with the higher transition energies (features of the excited state landscape are smaller than those of the ground state landscape). The magnitude of the absorption spectrum shift is a reflection of the average landscape asymmetry (or, more precisely, the difference of the asymmetries of the ground and excited states), and, in principle, low-resolution measurements of the absorption spectrum evolution over a very long period of time may be used to determine this average asymmetry 7915

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B

Figure 4. Frame A: Comparison of the HGK curves calculated for 5 and 20 K starting from equilibrium situation at indicated burn temperature (black, red) and from nonequilibrium situation resulting from realistic cooling (blue, green). The correction for recovery between the acts of burning was on. Frame B: HGK curves for 5 K (black, same as blue curve in frame A), 7 (blue), 9 (red), 11 (light green), 13 (magenta), 15 (light blue), 17 (brown), and 20 K (dark green), calculated for realistic cooling conditions.

in the HGK measurements. Thus, in the following, we will not utilize this correction, allowing for a significant increase in the calculation speed. Figure 4A presents the dependence of HGK (reflecting the burning rate) on the initial conditions, namely on how far from equilibrium the ensemble was prior to burning, as well as on burn temperature. First, as expected, the hole burning is slower when performed on the sample in equilibrium conditions. In equilibrium the deepest wells are more populated (see Figure 2), and burning, if any, more likely involves uphill tunneling rather than downhill tunneling. As the temperature increases, the ratio of uphill and downhill tunneling rates is decreased, and the probability to find the system in wells other than the deepest one is increased. Thus, the difference between equilibrium and nonequilibrium situation is diminished (see green and red curves in Figure 4A). Second, the hole burning rate is clearly increasing with temperature. HGK curves calculated for realistic cooling conditions and additional temperatures between 5 and 20 K are shown in Frame B. In ref 12, the slowdown of burning with temperature was reported and attributed to the increase of the homogeneous line width and electron−phonon coupling. These effects are not included in this simple model, and result in a decrease of the peak absorption cross-section. In Frame B of Figure 4, the hole burning is accelerated by a factor of about five between 5 and 15 K. This is comparable to the increase of the homogeneous line width observed in CP4341 and Cytochrome b6f31 over the same temperature range. However, explaining the slowdown of burning with temperature by a factor of about ten12,31 is impossible with the current model. Thus, in experiments the burning appears to be dominated by the downhill tunneling to a degree that is impossible in either equilibrium situation or nonequilibrium situation resulting from a realistic cooling process. See the Discussion on NPHB mechanism34 below for more details. Overall, the HGK are multiexponential and are in qualitative agreement with experiments and modeling based on the NPHB master equation,11,12,31,35,36. One has to remember that at the present stage our software is incapable of properly modeling the evolution of the width of a narrow hole. Thus, here the hole width does not change, the hole depth is proportional to the hole area and care should be taken when comparing the experimental HGK data with those simulated here. Since in reality the holes broaden with the increase of the depth,42 the

(although this is obviously not practical). Arrow indicates the burn frequency of 14560 cm−1 employed in subsequent simulations. This is also the mean of the distribution of the gaps between ground and excited states baselines, which has the width of 180 cm−1. Two remarks are in order here, one philosophical and another technical. A narrow spectral hole represents a much more drastic departure from the equilibrium than the nonequilibrium result of the realistic cooling. Thus, spectral hole recovery is expected to appear dominated by the return to the preburn state. However, there is no fundamental difference between different kinds of evolution toward equilibrium. Thus, any recovery observed in experiment can include contributions that are independent of the return of the systems to exact preburn configuration, and are the manifestations of the general drift toward equilibrium. In other words, this model does not imply the presence of the perfect “spectral memory”. Second, the absorption spectra reported in Figure 3 are the weighted sums of the spectra of all possible systems, whether they possess a pair of the ground and excited state wells with particular transition energy (corresponding to the burn frequency) or not. On the other hand, for subsequent simulations of hole-burning related phenomena we selected only a subensemble of 5000 systems which do possess a state with the particular transition energy. The absorption spectrum of such a subensemble may deviate from the Gaussian shape (not shown). However, as long as we are concerned with the difference spectra, the results are not affected by this choice of the subensemble: the systems with no wells resonant with the laser do not burn and therefore cannot contribute to our observations in any wayin the form of the hole spectra per se, or the hole growth or hole recovery kinetics data. 3.2. Hole Burning. The hole burning is performed after the sample has been cooled down. To compare the results of spectral hole burning simulations with the experiment, one can normalize to one the sum of the probabilities to find each system in the resonant well prior to burning (this sum is proportional to the preburn absorption). Subsequent evolution of the sum of the resonant well population probabilities will result in a hole growth kinetics (HGK) curve. As justified in the Supporting Information, correction for recovery in between individual acts of burning is negligible for phonon fluxes corresponding to tens of microwatts per square centimeter used 7916

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B HGK predicted here will be somewhat faster than in the experiment or in the simulations based on the NPHB master equation11,12,31,35,36for the same parameters. Additionally, the model reported here does not yet include the linear electron− phonon coupling. In the first approximation, introducing linear electron−phonon coupling with Huang−Rhys factor S may be taken into account just by reducing the peak absorption crosssection of the zero-phonon line by a factor of exp(−S). However, this approach would be sufficient only for modeling the HGK far from saturation, and it still would not yield any phonon sidebands in the hole spectra. Note that this correction will have to be introduced for all temperatures, and therefore this will not resolve the problem with the slowdown of burning observed in the experiment. Calculated hole spectra for both flat and parabolic energy landscape models are depicted in Figure 5. No significant

Figure 6. Fixed temperature (5 K, burn temperature) recovery of spectral holes of two different initial fractional depths. Flat randomness-dominated landscape; the burning is performed in nonequilibrium conditions resulting from a realistic cooling procedure. Black: 58%-deep hole. Blue: 37%-deep hole. Red curve represents the recovery of a 58% hole that was produced after the sample was allowed to equilibrate for a month at 5 K.

imental, theoretical or simulation data on the dependence of the hole recovery on the initial burn conditions. Red curve in Figure 6 represents the recovery of a 58%-deep hole that was burned after the ensemble was allowed to equilibrate for a month. This hole recovered slightly faster than the hole produced right after realistic 2 h cooling. This is in agreement with larger contribution of the uphill tunneling to burning and larger contribution of the downhill tunneling closer to equilibrium. Note that full recovery of the holes from Figure 6 would require ∼108 s or about 3 years! 3.4. Thermocycling. The result of simulating a thermocycling experiment similar to those reported in refs 12 and 31 is shown in Figure 7. Horizontal axis presents time rather than

Figure 5. 30%-deep spectral holes produced with flat (red) and parabolic (black) energy landscape models at 5 K. Resolution 10 GHz. Hole burning starts in nonequilibrium conditions resulting from realistic cooling procedure.

qualitative differences can be seen. In the case of the “flat” landscape model, the NPHB antihole is concentrated closer to the burn frequency. In the case of the parabolic model, larger spectral shifts become possible. Only a slight preference for red NPHB-induced spectral shifts is observed. Phonon sidebands are not included in the model yet. 3.3. Hole Recovery at Fixed (Burn) Temperature. To determine the parameters of the energy landscape in the ground state, one can study the recovery of the spectral hole at fixed burn temperature (e.g., 5 K) and/or upon thermocycling. First, in agreement with the experiments and with the simpler TLS model with the perfect spectral memory,12,31 the rate of the recovery process depends on the initial fractional hole depth. Deeper holes recover slower than the shallower ones, see Figure 6. This effect is present for both the flat randomnessdominated landscape and the parabolic one. As stated in the introduction, it can be interpreted as a manifestation of “spectral memory” or hole being recovered as a result of burnt pigment−protein systems returning to the preburn state. In the light of the above arguments, it is not surprising that the recovery of a hole also depends on the conditions under which the hole was burned, namely on how far the system was from equilibrium at the time of burning. While dependence of the hole broadening (small-shift spectral diffusion) on the equilibration state of the system (“waiting time”) has been reported previously,15,17,18,20 there appears to be no exper-

Figure 7. Thermocycling data for a 37%-deep hole simulated using our MLS model (blue curve). “Flat” randomness-dominated landscape; realistic cooling conditions. Black diamonds: experimental data for a 37%-deep CP43 hole.12

temperature, and the temperature followed the T(t) profile from the actual experiment. Some temperatures are explicitly indicated in the figure. The first 3 h was spent at burn temperature (5 K). Thermocycling to somewhat higher temperatures initially did not result in an acceleration of the calculated hole recovery. However, calculated recovery accelerates drastically between 30 and 40 K. This behavior 7917

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B

Figure 8. Probability distributions for a single 21-well system in (A) equilibrium; (B) preburn state resulting from realistic cooling; (C) after the first laser scan and (D) after 2000 laser scans. Note different scales of the vertical axis in different frames.

about 3 × 107 photons per second are absorbed under these conditions if the laser is in resonance with the ZPL, which corresponds to 30 ns intervals between individual acts of excitation while the system is in resonance with the laser. The population probability distribution evolution upon repeated wavelength scanning can be modeled as follows: While the wavelength of the laser is changing during the scan over the absorption band, different wells are chosen accordingly as the resonant ones. These wells obviously are accessed not in the order of the well numbers. The scanning speed is slow enough for multiple individual acts of excitation (∼104) to occur while the laser is in resonance with a particular well. Then the average value of the postburn probability after all acts of burning in this well is considered as the post burn probability for the purpose of calculating the observable (optical signal at given wavelength). The average value is chosen because in practice, the signal at each wavelength is measured during the scan process and thus is a result of averaging. For instance, some single molecule line will be detected even if a conformational change has occurred during the time the ZPL of that molecule was in resonance and by the end of this time interval the molecule was no longer in resonance with the laser. (Here it is important to distinguish between probability distributions and observed spectra. The probabilities P should add up to one at any given moment, but the spectra do not have to be normalized, since one may chase one and the same protein−chromophore system through several wells and excite it at several different wavelengths within one and the same laser scan.) As the wavelength is changing, another well will come into resonance with the laser, and the SHB-style part of the modeling will be repeated. Otherwise, if the system does not possess a resonant well at the next wavelength, the system will spend respective time recovering. Below, we report the simulation of a hypothetical singlecomplex experiment on a chlorophyll-containing sample. The scanned range is located between 14506 and 14622 cm−1 (corresponding to 689.4 and 683.9 nm, respectively). The

closely resembles predictions of the model based on the NPHB master equation12,31 (TLS, perfect spectral memory) and differs significantly from the experimental observations (black circles). Thus, introducing multilevel energy landscapes and realistic imperfections of the spectral memory does not explain more gradual recovery actually observed in thermocycling experiments on CP4312 and Cytochrome b6f.31 Crossing between experimental and calculated curves at about 39 K indicates that md2 = 1.4 ×10−46 kg·m2 is somewhat too large, in agreement with ref 12. 3.5. Quasi-SPCS. The advantage of our current approach over NPHB master equation is in keeping track of all changes of the well population probabilities, which makes the SPCS simulations possible as well. The procedure of converting the evolution of the probability distributions into the discrete line jumps is well-known: the probability to find the system in one of the wells is set to one and the rest are set to zero; the probability of transition between the wells is repeatedly compared to a uniformly distributed random number between zero and one and if the transition probability is larger than the random number it is assumed that the transition has occurred; the probability to find the system in the new well is set to one, and so on. However, one can obtain useful information comparable to published SPCS results5 directly from the probability distributions. The procedure described here yields the distributions of probability to find the spectral line of a pigment−protein system at various wavelengths in the course of the SPCS experiment, which are similar to the first cumulant distributions of Köhler and co-workers.5 As the current model does not include small spectral shifts affecting the apparent line width or the spectral hole width, our results should be compared to the probabilities of finding the system in the different higher-tier states from refs 3 and 5, separated by about 10 cm−1. In SPCS experiments, the absorption band of the protein complex is scanned with relatively high-intensity light (∼5−10 W/cm2) multiple (∼103) times with a speed of ∼40 cm−1/s.5 With the peak absorption cross-section of ∼10−12 cm2, 7918

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B laser intensity and scan rate are 5 W/cm2 and 44 cm−1/s, respectively.5 Flat randomness-dominated energy landscape is employed. Figure 8 compares probability distributions corresponding to thermodynamic equilibrium (A) and before the beginning of the optical experiment after realistic cooling (B). Note significantly different vertical scales in different frames. Comparison between the distribution after the first scan (C) and after 2000 scans (D) and equilibrium distribution (A) reveals that repeated scanning of the laser gradually drives the distribution toward the thermal equilibrium.

4. DISCUSSION: EXTENSION OF THE SHU−SMALL NPHB MECHANISM TO PROTEINS As discussed above, the slowdown of the hole burning with the increase of temperature observed in experiments12,31 can only be explained assuming that the NPHB yield is very weakly dependent on temperature (at least up to a certain point, 10− 13 K12,31). This is possible only if the NPHB process is dominated by the downhill tunneling, i.e. the tunneling associated with the increase of energy of the scattered phonon. According to Figure 5, however, increase of the apparent NPHB rate is expected for both equilibrium and realistic nonequilibrium starting conditions. Thus, one cannot explain the experimental observations by the realistically nonequilibrium starting conditions of the multiwell system or by the imperfect spectral memory. One should still require the system to get poised for the downhill tunneling in the excited state prior to tunneling, and that it has to be poised for the downhill tunneling also in the ground electronic state for the efficient hole recovery. (If a particular system ends up significantly closer to thermal equilibrium as a result of burning, it is not expected to recover.) This in turn requires permanent opposite tilts of the ground and excited states of the landscape, and strictly the red shift of the antihole absorption, in disagreement with the experiment. Similar problems arise with interpreting spectral hole burning results in glasses.34−36 To address these problems, Shu and Small proposed a NPHB mechanism that involves interactions between “extrinsic” (containing the chromophore) and “intrinsic” TLS of the glassy solid that result in the outsidein cascade of the tunneling events and migration of the free volume toward the chromophore.34 The sequence of the events in the NPHB process modified in view of Shu and Small’s arguments is depicted in Figure 9, with only a double-well system shown for clarity (The TOC graphic shows the threewell version, but not all the stages of the process). Before burning, the system is somewhat more likely in the lowerenergy minimum of the ground state of the TLS/pigment system due to the relaxation during cooling process. At stage 1, absorption of a resonant photon brings the pigment to the excited state. At the second step, excitation of the pigment molecule triggers the transitions in the “intrinsic” TLS and diffusion of the free volume toward the chromophore, changing the tilt in the “extrinsic” TLS, and poising it for downhill tunneling in the excited state. Once the pigment returns to the ground electronic state (stage 3), the tilt of the excited state of the TLS/pigment system is reversed. The tilt of the ground state may or may not be reversibly switched together with the tilt of the excited state, it is only required that the tilt after burning is the same as prior to burning. Stage 4 represents hole recovery. The particular TLS arrangement shown in Figure 9 results in a red shift of the postburn absorption (Eeduct > Eproduct) but for a different combination of the ground and excited state TLS asymmetries blue shift is possible; also note that the actual

Figure 9. Shu−Small NPHB mechanism.34 The small TLS between stages 1 and 2 symbolize interactions with the surroundings (described in amorphous solids as “intrinsic TLS”, as opposed to the “extrinsic TLS” involving the pigment).

energy of absorbable photon (before burning, stage 1 in Figure 9, or after burning, stage 4 may not be straightforwardly related to the energy differences at stage 2. In pigment−protein complexes involved in photosynthesis, the distinction between extrinsic and intrinsic TLS fades as the pigment is a natural part of the complex and not an artificially introduced dopant. The exact nature of the tunneling systems is not known. In general, either protons or small side groups of the protein (−CH3, −OH) were considered as the candidates for the tunneling entities in proteins.7,12,31 One can also consider the possible effects of surface TLS or some aspects of the protein dynamics following the dynamics of the surrounding amorphous host (buffer/glycerol mixture)i.e., “slaving”.43 Thus, at the moment one cannot answer the question about the exact mechanism of interaction between several tunneling entities, beyond requiring that such a mechanism must involve diffusion of the free volume toward the chromophore while the chromophore is in the excited electronic state. We may also refer the reader to an excellent review of various cooperative effects and extended definition of allostery, discussing the interplay of entropic and enthalpic contributions.44 Although44,45 are devoted mainly to the changes involving the protein backbone, similar logic could be extended to the movements of the side groups as well. Other manifestations of cooperative effects in photosynthetic proteins at low temperatures include switching between different energy transfer pathways in Photosystem I8 and the drastic reduction of the attempt frequency proposed for Cytochrome b6f,31 and, in retrospect, probably for CP4312 as well.

5. CONCLUSIONS We reported the development and initial applications of the spectral diffusion modeling software involving random generation of multiwell protein energy landscapes and including both temperature-dependent tunneling and barrierhopping in both the ground and excited states of a protein/ chromophore system. In other words, our model includes 7919

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

The Journal of Physical Chemistry B

■

realistic proportions of light- and thermally induced spectral diffusion. The model, also possessing the realistic degree of spectral memory, adequately captures the key features of the evolution of the pigment−protein system(s), including initial cooling down of the sample, spectral hole burning (NPHB) and hole recovery. The results resemble those obtained using models based on the NPHB master equation, meaning that there are no fundamental differences between the NPHB observations involving TLS and MLS. Moreover, modeling of the thermocycling process is also in agreement with the simulations of ref 12 and 31 and not with the experimental results, indicating that neither MLS nor realistic nonequilibrium population probability distributions resolve the problem of the mismatch between the barrier distributions manifesting in fixedtemperature recovery and thermocycling experiments.12,31 Typical cooling procedures result in proteins quite far from thermodynamic equilibrium (for the protein energy landscape tier involved in spectral hole burning), but neither nonequilibrium starting conditions nor increased number of wells can eliminate the need for the Shu−Small type NPHB mechanism34 involving migration of the free volume toward the chromophore. It is important to keep in mind the sample’s cooling history while performing and analyzing the data produced in SHB and SPCS experiments. Determining the effects of different protein energy landscapes on the antihole shape is particularly important for determining how lightinduced site energy shifts (hole burning) alter the picture of excitonic interactions in PS complexes, and for modeling their various resonant and nonresonant hole spectra. Initial simulations imitating typical SPCS experiments suggest that in such experiments each pigment−protein system experiences frequent light-induced spectral shifts, and has a finite probability to be eventually observed in all available wells. Thus, sufficiently long-duration SPCS experiments likely yield the approximate number of available conformational substates, which is relatively small (but larger than two). Our calculations also show that in SPCS experiments repeated scanning of the whole spectrum with the laser drives the system toward thermodynamic equilibrium. Further improvements of the model will include incorporating electron−phonon coupling and phonon sidebands. The modeling framework outlined here can easily be extended to energy landscapes with different topologies, including the multidimensional landscapes, as well as different relationships between the ground and excited state landscapes. One can also envision a hybrid approach when the energy landscapes derived in the course of molecular dynamics simulations are used for spectral diffusion modeling. Doing it properly would require understanding to what extent the dynamics of the chromophore pocket is affected by the dynamics of protein, surrounding solvent, lipid membrane, or detergent micelle.

■

Article

AUTHOR INFORMATION

Corresponding Author

*(V.Z.) E-mail: [email protected]. Telephone: 1514-848-2424 ext. 5050. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS Authors express gratitude to NSERC and Concordia University. We also thank Mr. Jan Stubben, DAAD Rise in North America undergraduate scholarship holder, who participated in the initial stages of this project.

■

REFERENCES

(1) Berlin, Y.; Burin, A.; Friedrich, J.; Köhler, J. Spectroscopy of Proteins at Low Temperature. Part I: Experiments with Molecular Ensembles. Phys. Life Rev. 2006, 3, 262−292. (2) Berlin, Y.; Burin, A.; Friedrich, J.; Köhler, J. Low Temperature Spectroscopy of Proteins. Part II: Experiments with Single Protein Complexes. Phys. Life Rev. 2007, 4, 64−89. (3) Hofmann, C.; Aartsma, T. J.; Michel, H.; Köhler, J. Direct Observation of Tiers in the Energy Landscape of a Chromoprotein: A Single-molecule Study. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 15534− 15538. (4) Baier, J.; Richter, M. F.; Cogdell, R. J.; Oellerich, S.; Köhler, J. Do Proteins at Low Temperature Behave as Glasses? A Single-Molecule Study. J. Phys. Chem. B 2007, 111, 1135−1138. (5) Baier, J.; Richter, M. F.; Cogdell, R. J.; Oellerich, S.; Köhler, J. Determination of the Spectral Diffusion Kernel of a Protein by SingleMolecule Spectroscopy. Phys. Rev. Lett. 2008, 100, 018108. (6) Oikawa, H.; Fujiyoshi, S.; Dewa, T.; Nango, M.; Matsushita, M. How Deep Is the Potential Well Confining a Protein in a Specific Conformation? A Single-Molecule Study on Temperature Dependence of Conformational Change between 5 and 18 K. J. Am. Chem. Soc. 2008, 130, 4580−4581. (7) Brecht, M.; Studier, H.; Radics, V.; Nieder, J. B.; Bittl, R. Spectral Diffusion Induced by Proton Dynamics in Pigment-Protein Complexes. J. Am. Chem. Soc. 2008, 130, 17487−17493. (8) Brecht, M.; Radics, V.; Nieder, J. B.; Bittl, R. Protein Dynamicsinduced Variation of Excitation Energy Transfer Pathways. Proc. Natl. Acad. Sci. U.S.A. 2009, 106, 11857−11861. (9) Tietz, C.; Jelezko, F.; Gerken, U.; Schuler, S.; Schubert, A.; Rogl, H.; Wrachtrup, J. Single Molecule Spectroscopy of the LightHarvesting Complex II of Higher Plants. Biophys. J. 2001, 81, 556− 562. (10) Krüger, T. P. J.; Novoderezhkin, V. I.; Ilioaia, C.; van Grondelle, R. Fluorescence Spectral Dynamics of Single LHCII Trimers. Biophys. J. 2010, 98, 3093−3101. (11) Grozdanov, D.; Herascu, N.; Reinot, T.; Jankowiak, R.; Zazubovich, V. Low-Temperature Protein Dynamics of the B800 Molecules in the LH2 Light-Harvesting Complex: Spectral Hole Burning Study and Comparison with Single Photosynthetic Complex Spectroscopy. J. Phys. Chem. B 2010, 114, 3426−3438. (12) Najafi, M.; Herascu, N.; Seibert, M.; Picorel, R.; Jankowiak, R.; Zazubovich, V. Spectral Hole Burning, Recovery, and Thermocycling in Chlorophyll-Protein Complexes: Distributions of Barriers on the Protein Energy Landscape. J. Phys. Chem. B 2012, 116, 11780−11790. (13) Burin, A. L.; Berlin, Yu. A.; Patashinski, A. Z.; Ratner, M. A.; Friedrich, J. Fluctuations of Conformational States in Biological Molecules: Theory for Anomalous Spectral Diffusion Dynamics. Physica B 2002, 316-317, 321−323. (14) Köhler, W.; Friedrich, W.; Scheer, H. Conformational Barriers in Low-temperature Proteins and Glasses. Phys. Rev. A 1988, 37, 660− 662. (15) Schlichter, J.; Friedrich, J. Glasses and Proteins: Similarities and Differences in their Spectral Diffusion Dynamics. J. Chem. Phys. 2001, 114, 8718−8721.

ASSOCIATED CONTENT

S Supporting Information *

Brief description of the numerical methods used to determine the evolution of a chromophore−protein system, estimation of times between the individual acts of excitation under realistic conditions of SHB experiments, as well as some results, and discussion pertaining to possible corrections for recovery in between individual acts of burning. The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.5b02764. 7920

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921

Article

The Journal of Physical Chemistry B

(38) Barkai, E.; Naumov, A. V.; Vainer, Yu. G.; Bauer, M.; Kador, L. Lévy Statistics for Random Single-Molecule Line Shapes in a Glass. Phys. Rev. Lett. 2003, 91, 075502. (39) Naumov, A. V.; Vainer, Yu. G.; Kador, L. Does the Standard Model of Low-Temperature Glass Dynamics Describe a Real Glass? Phys. Rev. Lett. 2007, 98, 145501. (40) Herascu, N.; Najafi, M.; Amunts, A.; Pieper, J.; Irrgang, K.-D.; Picorel, R.; Seibert, M.; Zazubovich, V. Parameters of the Protein Energy Landscapes of Several Light-Harvesting Complexes Probed via Spectral Hole Growth Kinetics Measurements. J. Phys. Chem. B 2011, 115, 2737−2747. (41) Jankowiak, R.; Zazubovich, V.; Rätsep, M.; Matsuzaki, S.; Alfonso, M.; Picorel, R.; Seibert, M.; Small, G. J. The CP43 Core Antenna Complex of Photosystem II Possesses Two Quasi-Degenerate and Weakly Coupled Qy-Trap States. J. Phys. Chem. B 2000, 104, 11805−11815. (42) Herascu, N.; Ahmouda, S.; Picorel, R.; Seibert, M.; Jankowiak, R.; Zazubovich, V. Effects of the Distributions of Energy or Charge Transfer Rates on Spectral Hole Burning in Pigment-Protein Complexes at Low Temperatures. J. Phys. Chem. B 2011, 115, 15098−15109. (43) Fenimore, P. W.; Frauenfelder, H.; McMahon, B. H.; Parak, F. G. Slaving: Solvent Fluctuations Dominate Protein Dynamics and Functions. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 16047−16051. (44) Tsai, C.-J.; del Sol, A.; Nussinov, R. Allostery: Absence of a Change in Shape Does not Imply that Allostery is not at Play. J. Mol. Biol. 2008, 378, 1−11. (45) Tsai, C.-J.; Ma, B.; Nussinov, R. Folding and Binding Cascades: Shifts in Energy Landscapes. Proc. Natl. Acad. Sci. U.S.A. 1999, 96, 9970−9972.

(16) Skinner, J. L.; Friedrich, J.; Schlichter, J. Spectral Diffusion in Proteins: A Simple Phenomenological Model. J. Phys. Chem. A 1999, 103, 2310−2311. (17) Zollfrank, J.; Friedrich, J.; Vanderkooi, J. M.; Fidy, J. Conformational Relaxation of a Low-temperature Protein as Probed by Photochemical Hole Burning: Horseradish Peroxidase. Biophys. J. 1991, 59, 305−312. (18) Fritsch, K.; Friedrich, J.; Parak, F.; Skinner, J. L. Spectral Diffusion and the Energy Landscape of a Protein. Proc. Natl. Acad. Sci. U.S.A. 1996, 93, 15141−15145. (19) Leeson, D. T.; Wiersma, D. A.; Fritsch, K.; Friedrich, J. The Energy Landscape of Myoglobin: An Optical Study. J. Phys. Chem. B 1997, 101, 6331−6340. (20) Gafert, J.; Pschierer, H.; Friedrich, J. Proteins and Glasses: A Relaxation Study in the Millikelvin Range. Phys. Rev. Lett. 1995, 74, 3704−3707. (21) Wolynes, P. G.; Luthey-Schulten, Z.; Onuchic, J. N. Fast-folding Experiments and the Topography of Protein Folding Energy Landscapes. Chem. Biol. 1996, 3, 425−432. (22) Frauenfelder, H.; Sligar, S. G.; Wolynes, P. G. The Energy Landscapes and Motions of Proteins. Science 1991, 254, 1598−1603. (23) Fenimore, P. W.; Frauenfelder, H.; McMahon, B. H.; Young, R. D. Proteins are Paradigms of Stochastic Complexity. Physica A 2005, 351, 1−13. (24) Tavernelli, I.; Cotesta, S.; Di Iorio, E. E. Protein Dynamics, Thermal Stability, and Free-Energy Landscapes: A Molecular Dynamics Investigation. Biophys. J. 2003, 85, 2641−2649. (25) Okazaki, K.; Koga, N.; Takada, S.; Onuchic, J. N.; Wolynes, P. G. Multiple-basin Energy Landscapes for Large-amplitude Conformational Motions of Proteins: Structure-based Molecular Dynamics Simulations. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 11844−11849. (26) Harris, B. J.; Cheng, X.; Frymier, P. All-atom Molecular Dynamics of a Photosystem I/Detergent Complex. J. Phys. Chem. B 2014, 118, 11633−11645. (27) Chandler, D.; Strümpfer, J.; Sener, M.; Scheuring, S.; Schulten, K. Light Harvesting by Lamellar Chromatophores in Rhodospirillum photometricum. Biophys. J. 2014, 106, 2503−2510. (28) van Amerongen, H.; Valkunas, L.; van Grondelle, R. Photosynthetic Excitons; World Scientific: Singapore, 2000. (29) Purchase, R.; Völker, S. Spectral Hole Burning: Examples from Photosynthesis. Photosynth. Res. 2009, 101, 245−267. (30) Jankowiak, R.; Reppert, M.; Zazubovich, V.; Pieper, J.; Reinot, T. Site Selective and Single Complex Laser-Based Spectroscopies: A Window on Excited State Electronic Structure, Excitation Energy Transfer, and Electron-Phonon Coupling of Selected Photosynthetic Complexes. Chem. Rev. 2011, 111, 4546−4598. (31) Najafi, M.; Herascu, N.; Shafiei, G.; Picorel, R.; Zazubovich, V. Conformational Changes in Pigment-Protein Complexes at Low Temperatures - Spectral Memory and Possibility of Cooperative Effects. J. Phys. Chem. B 2015, DOI: 10.1021/acs.jpcb.5b02845. (32) Anderson, P. W.; Halperin, B. I.; Varma, C. M. Anomalous Lowtemperature Thermal Properties of Glasses and Spin Glasses. Philos. Mag. 1972, 25, 1−9. (33) Phillips, W. A. Tunneling States in Amorphous Solids. J. Low. Temp. Phys. 1972, 7, 351−360. (34) Shu, L.; Small, G. J. On the Mechanism of Nonphotochemical Hole Burning of Optical Transitions in Amorphous Solids. Chem. Phys. 1990, 141, 447−455. (35) Reinot, T.; Dang, N. C.; Small, G. J. A Model for Persistent Hole Burned Spectra and Hole Growth Kinetics that Includes Photoproduct Absorption: Application to Free Base Phthalocyanine in Hyperquenched Glassy Ortho-dichlorobenzene at 5 K. J. Chem. Phys. 2003, 119, 10404−10414. (36) Dang, N. C.; Reinot, T.; Reppert, M.; Jankowiak, R. Temperature Dependence of Hole Growth Kinetics in AluminumPhthalocyanine- Tetrasulfonate in Hyperquenched Glassy Water. J. Phys. Chem. B 2007, 111, 1582−1589. (37) Reinot T. Private communication 7921

DOI: 10.1021/acs.jpcb.5b02764 J. Phys. Chem. B 2015, 119, 7911−7921