Address correspondence to R.D. Hamer, Smith-Kettlewell Eye Research Institute, 2318 Fillmore Street, San Francisco, CA 94115. Fax: (415) 345-8455; email: russ{at}ski.org
![]() |
ABSTRACT |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Key Words: single-photon response reproducibility stochastic model multiple phosphorylation phototransduction rod responses
![]() |
INTRODUCTION |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Evidence for Reproducibility
Empirical measures have shown that the variability of vertebrate rod SPRs (assessed by any of several measures of amplitude or kinetics) is considerably lower than expected for a one-step R*-inactivation process. In this paper, we focus on four findings reported in the literature. The first finding is that the coefficient of variation of the amplitude of the SPRs (CVampl) is quite small (0.20.25: Baylor et al., 1979; Schneeweis and Schnapf, 1995
; Rieke and Baylor, 1998a
; Whitlock and Lamb, 1999
; Field and Rieke, 2002
). The second finding is that the ensemble variance of a set of dim-flash responses is very nearly proportional to square of the mean of the response over almost the entire response time course. This observation has been used to suggest that variability in the elementary response waveform must be small (Schnapf, 1983
; Schneeweis and Schnapf, 1995
; Rieke and Baylor, 1998a
). That argument has been rejected by Whitlock and Lamb (1999)
, and we will elaborate on their insight that this finding alone does not imply that the SPR has low variability over its entire time course. The third finding is the low coefficient of variation of the area (CVarea
0.3) under the SPR (Field and Rieke, 2002
). This measure has a natural pragmatic appeal, in that it captures variability in both the amplitude and time course of the SPR (Field and Rieke, 2002
), and we will show furthermore that there is a theoretical basis for choosing variability in SPR area to estimate the number of inactivation steps of R* over measures of variability of either amplitude or kinetics alone. The fourth measure is based on the finding that the variance of the SPR is at least an order of magnitude smaller than the square of the mean response, until some time after the mean response reaches its peak (Rieke and Baylor, 1998a
; Field and Rieke, 2002
). Field and Rieke (2002)
found that the relative times-to-peak of the SPR variance and SPR squared mean, as well as the width of the SPR variance waveform aided in discriminating between candidate models.
A New Stochastic Model of Phototransduction Accounts for the Variability/Reproducibility of Single-photon Responses, as well as Other Key Electrophysiological Data
Based, in part, on the ideas and biochemical data of Gibson et al. (2000), we develop in this paper a kinetic model of phototransduction that includes a detailed stochastic simulation of the activation and inactivation of rhodopsin, G-protein (G), and phosphodiesterase (PDE). Simulations of the stochastic activation of PDE have been done previously by Felber et al. (1996)
and Lamb (1994)
. Our modeling is distinguished from theirs in four main respects: (a) We model reaction kinetics in greater detail than the prior studies, including competitive, stochastic binding of R* with G-protein, arrestin (Arr), and rhodopsin kinase. (b) Rather than simulating the two-dimensional diffusional contact between molecules, we concentrate on the reactions of those molecules by assuming that the disc surface is well mixed. (c) Neither Lamb (1994)
nor Felber et al. (1996)
coupled the activated PDE to the downstream reactions generating photocurrent (including cGMP-hydrolysis, channel closure, Ca2+-feedback, and Ca2+-buffering). Our model does so, thus allowing direct comparisons between empirical electrophysiological responses and model responses under a wide range of conditions. (d) Lamb (1994)
did not simulate inactivation reactions.
Our model embodies the following four hypotheses: (a) that R* undergoes a series of sequential phosphorylation steps, mediated by rhodopsin kinase (Kühn and Wilden, 1982; Wilden and Kühn, 1982
; Aton et al., 1984
; Thompson and Findlay, 1984
; Palczewski et al., 1991
; Wilden, 1995
); (b) that competition for binding to R* occurs between three molecular species: inactive G-protein (G·GDP), rhodopsin kinase (RK), and Arr (Pfister et al., 1983
; Miller and Dratz, 1984
; Buczylko et al., 1991
; Pulvermuller et al., 1993
; Krupnick et al., 1997
; Gibson et al., 2000
); (c) that each sequential step of phosphorylation leads to a progressive decrease in affinity between R* and G·GDP, and concomitantly to a progressive increase in affinity between R* and Arr (Gibson et al., 2000
); and (d) that the affinity between R* and RK also ratchets down with each step of phosphorylation (Buczylko et al., 1991
; Gibson et al., 2000
). We will refer to our full model based on these premises as the sequential phosphorylation model.
Using this model, we simulate ensembles of SPRs and dim-flash responses using Monte-Carlo methods. The model accounts for the four quantitative measures of SPR variability/reproducibility delineated above. In addition, despite the fact that in the model phosphorylation is the dominant process inactivating R*, the simulations reproduce key experimental results that Rieke and Baylor (1998a) interpreted previously as evidence that phosphorylation could not be the dominant deactivation mechanism.
Moreover, without any additional parameter adjustments, the model reproduces the salient qualitative features of SPRs obtained from four genetic knockout and transgenic studies: three studies in which mouse rods had been genetically modified with specific pathways of R* inactivation disabled (Arr-/-, Xu et al., 1997; RK-/-, Chen et al., 1999
; CSM, Mendez et al., 2000
), and one study in which the mechanism of feedback synthesis of guanylate cyclase had been disabled (GCAPs-/-; Burns et al., 2002
). In general, we find that attempting to simulate the results of knockout and transgenic experiments can provide valuable insights into the kinds of mechanisms that must be present, or about candidate molecular schemes that can be ruled out even if, in principle, they could achieve the empirical SPR reproducibility in one or more of the four quantitative measures of reproducibility.
Theory for Mechanisms Underlying SPR Reproducibility
Finally, we present a theoretical analysis of mechanisms underlying SPR reproducibility. It reveals a heretofore unrecognized tradeoff between variability of SPR amplitude and variability of SPR kinetics that depends critically on R* lifetime in relation to the net kinetics of the downstream reactions in the cascade. We explain what the statistics of SPR area tell us about underlying molecular mechanisms. We show how the coefficient of variation of SPR area can be used, in the context of a certain class of models, to estimate a lower bound for the number of R* inactivation steps required to yield observed SPR reproducibility. Neither the CV of SPR amplitude nor the statistics of SPR duration variability alone can be used in this manner. An important conclusion of this theoretical analysis is that, in normal rods, across species, the kinetics of R* inactivation cannot be very different from the kinetics of the downstream reactions, including the inactivation of the activated transducinphosphodiesterase complex, in the dim-flash regime.
![]() |
MATERIALS AND METHODS |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
|
|
with the final, complete quench of R* activity occurring upon Arr binding to phosphorylated R*:
The notations "pre" and "post" in Eqs. distinguish RK-bound states of R* before and after phosphorylation.
Transducin is assumed to be activated by a conventional series of reactions (Eqs. ; see Lamb and Pugh, 1992
).
where G·GTP is the activated form of transducin.
Based on the biochemical results of Gibson et al. (2000), the affinity of rhodopsin for G-protein is assumed to decrease exponentially with increasing phosphorylations (Eq. 4), while its affinity for Arr is assumed to increase linearly with n (Eq. 5). Arrestin is able to bind R*, with increasing probability, at any time following the first phosphorylation (n
1), to terminate whatever R* activity remains in that state.
![]() | (4) |
![]() | (5) |
The exponential rate, , was set to 0.6 from an exponential fit to the data in Fig. 2 A in Gibson et al. (2000)
.
Choice of Phosphorylation Dependence for RRK Affinity
Although the explicit dependence of RRK binding affinity and rate on R* phosphorylation state has not been documented with the same detail as the RG and RArr affinities (Gibson et al., 2000), there is biochemical support for the notion that RRK affinity decreases systematically with the phosphorylation state of R*. For example, Buczylko et al. (1991)
found that phosphorylated RK has significantly lower affinity for phosphorylated R* than for unphosphorylated R*. Moreover, the data and analysis of Gibson et al. (2000)
shows that, given the opportunity for RK to interact with phosphorylated versus unphosphorylated rhodopsin, it preferentially interacts with unphosphorylated rhodopsin, consistent with the idea that RRK affinity decreases systematically as R* becomes phosphorylated.
We set the RRK affinity to have the same dependence on n as the affinity between R* and G-protein by varying kRK1.
![]() | (6) |
We treat the rate of incremental phosphorylation of the RRK complex by ATP as n dependent, in that it reduces to zero when all available sites on the R* carboxy terminus are occupied.
![]() | (7) |
In a model where RK and G bind R* competitively and Arr is only permitted to bind once R* is fully phosphorylated, this arrangement (Eqs. ) would cause, on average, an equal number of PDE* to be produced in each phosphorylation state, and it can be demonstrated analytically that this would minimize variability in the number of activated G-proteins produced per photoisomerization. In our implementation of the model, Arr is permitted to bind with increasing probability after the first phosphorylation (Fig. 4 D, inset). Consequently, derivation of the optimal phosphorylation-dependence for RRK affinity is quite complicated. However, we have estimated the ideal affinity profile by numerical optimization techniques, and have verified that, with model parameters that provide a good account of the data, the affinity profile embodied in Eq. 6 yields nearly ideal minimization of variability.1
We also examined the effect of having a different exponential phosphorylation dependence for the RG and RRK affinities. We found that (assuming the for RG affinity in Eq. 4 was fixed at the empirical estimate of 0.6), if
for RRK affinity in Eq. 6 varied between
0.3 and
0.7, the predicted variability of PDE* production would only vary by
5 percent. Finally, we examined other profiles for RRK affinity, such as linear decrease in affinity with n, or a flat profile (no dependence on phosphorylation state). These profiles led to a significant increase in the predicted variability in the number of PDE*s produced compared to the exponential profile used in the sequential phosphorylation model.
Because the phosphorylation dependence profile for RRK affinity we adopted yields good statistical performance of the model (e.g., low SPR variability), and because some of the details of such phosphorylation dependence have yet to be fleshed out experimentally, our choice of phosphorylation dependence for RRK affinity may be viewed as a testable biochemical prediction of our model.
Activation and inactivation of PDE (Fig. 1 B).
For the activation of PDE*, we assumed that activated transducin binds to the subunit of PDE (Eq.
), and that the inhibition of cGMP hydrolysis by the PDE
subunit is then relieved (Eq.
). For the inactivation of PDE*, we did not include explicit equations for acceleration of transducin (and, hence, PDE) inactivation due to interactions with the proteins RGS9-1 and Gß5 (Skiba et al., 2001
; Lishko et al., 2002
; Martemyanov and Arshavsky, 2002
; for review see He and Wensel, 2002
). Instead, we assumed that these reactions were fast, and we incorporated their effect into a single step of PDEinactivation (Eq.
).2 We have examined the effects of including such reactions, and found the effects on the simulations to be negligible.
where PDE*·G · GTP represents the activated form of the transducinPDE complex, referred to as PDE* throughout this paper, and PDE·G
· GDP represents the inactivated form (PDE). In all simulations, unless otherwise stated,
PDE was 3 s. The mean stochastic lifetime of R* (with all effects of multiple phosphorylation and Arr-binding included) was
2.8 s. However, the effective first-order time constant (defined as the first moment of the mean R* activity function) approximating the mean R* inactivation waveform was
1.3 s. Thus, PDEinactivation was (on average) the rate-limiting front-end reaction.
Deterministic model of "back-end" reactions (Fig. 2)
.
The front-end steps contribute considerable amplification, and Leskov et al. (2000) report that a single amphibian rod R* molecule activates
150 PDE*/s.3 The front-end parameters we used reproduce this PDE* activation rate up to the time of the first phosphorylation of R*. However, as phosphorylation of R* proceeds, the rate of PDE* production decreases (due to the phosphorylation-dependent decrease in RG affinity, Eq. 4). The net result is that each R* leads to the production of, on average,
220 PDE*.
|
![]() | (9) |
![]() | (10) |
![]() | (11) |
![]() | (12) |
We assumed that Ca2+ feedback occurs at guanylate cyclase only. Having cyclase feedback as the dominant feedback mechanism is consistent with recent evidence from GCAPs-/- rod responses (Burns et al., 2002), suggesting that Ca2+ feedback via recoverin and RK onto R* phosphorylation rate does not significantly affect the dim-flash response, though it does affect responses at high intensities.
Model Parameters
The values of the parameters for the sequential phosphorylation model (Eqs. ) are given in Table I. Parameters were set at modern literature estimates when these were known. In general, all parameter values for which estimates were available in the literature were within a factor of two of the literature values. A dark internal calcium concentration of 500 nM was assumed (Lagnado et al., 1992
; Gray-Keller and Detwiler, 1994
; Sampath et al., 1998
).
In order to capture the stochastic nature of the "front-end" of the cascade, in addition to R*'s interaction with Arr and RK, we implemented seven activation steps from photon absorption to PDE activation (Eqs. and
and
), based largely on the sequence delineated by Lamb and Pugh (1992)
. While empirical estimates were not available for each of these parameters, limits on some of these were estimated in Lamb and Pugh (1992)
(see legend to Table I for details). The front-end rate constants and the rate constant of cGMP hydrolysis per PDE hydrolytic subunit, ßsub, were adjusted as needed to achieve a close fit (by eye) to the ensemble mean of the Whitlock and Lamb (1999)
data, to generate a target number of PDE* activations per photoisomerization (200300), as well as to achieve the appropriate qualitative behavior in the ATP/GTP manipulation experiments of Rieke and Baylor (1998a)
.
The back-end parameters Kc, ßdark, fCa, k1, k2, and eT were optimized so that the model would generate flash responses that matched a representative set of empirical flash responses obtained from toad rods (from Rieke and Baylor, 1998b) over a wide range of flash intensities. More details about the choices of parameter values can be found in the legend to Table I.
Simulations and Monte-Carlo Methods
All aspects of the model and analyses were implemented using Matlab/Simulink (The Mathworks). The code is available upon request (contact first author).
Monte-Carlo simulations of 1,000 SPRs were run. Simulations were carried out using the Gillespie method (Gillespie, 1976, 1977
; Felber et al., 1996
). The probability that a molecular species takes a particular reaction pathway is given by the reaction rate for entering that pathway divided by the sum of rates of all available pathways. The dwell time or waiting time before taking one of the pathways is an exponentially distributed random variable with a mean given by the inverse sum of reaction rates for the possible pathways. For example, a free R* molecule's next event could be a binding with either G·GDP, RK, or Arr. To determine which of these happens, we generate a random number, r, uniformly distributed between 0 and 1. G·GDP binds if r < kG1
kG1 + kRK1 + kA, RK binds if kG1
kG1 + kRK1 + kA
r < kG1 + kRK1/kG1 + kRK1 + kA, otherwise Arr binds. The lag-time in the free R* state is determined by generating a second random number from an exponential distribution with mean equal to
kG1 + kRK1 + kA.
To simulate experimental data we generated a series of PDE* responses to Poisson numbers of photoisomerizations. Failures are all-zero vectors, and multiples are generated by summing the PDE* responses from randomly selected single-photon responses. The PDE* response is transformed to photocurrent through the deterministic back-end of the model (Eqs. ). Noise (both the continuous component of photoreceptor noise and instrument noise) was then added to the photocurrent response, as was done in Fig. 9 of Whitlock and Lamb (1999)
. Finally, a 2.5 Hz digital low-pass filter was applied to the responses (Whitlock and Lamb, 1999
).
Simulation of Other Conditions
In order to simulate alternate conditions (including effects of genetic modification of R* inactivation mechanisms), we made the following modifications to the sequential phosphorylation model and/or parameter values. Apart from the changes explicitly listed, all other parameters were unaltered. In each case, the mean number of photoisomerizations was 0.65, and we added noise as above.
Simulations of Transgenic and Genetic Knockout Data
The simulated genetically modified responses shown in Figs. 4 H, 7 F, and 8 F were generated using the means of 10 model SPRs in order to approximate the number of responses averaged in the corresponding experimental studies. These simulations represent the predictions for toad rods with the genetic manipulations used in four studies on mouse rods (Xu et al., 1997; Chen et al., 1999
; Mendez et al., 2000
; Burns et al., 2002
).
In order to simulate rhodopsin kinase knockout (RK-/-) responses (Chen et al., 1999), the parameter kRK1[RK] was set to zero across all n in Eq. 1a. Arr-/- (Xu et al., 1997
) was simulated by setting kA[Arr] = 0 across all n. In order to simulate the transgenic disabling of all phosphorylation sites on the carboxy terminus of rhodopsin (analogous to Mendez et al., 2000
, "completely substituted mutant", or CSM, in which serine and threonine residues were replaced by alanine), we set kRK3 = 0 across all n.
Altered Levels of ATP and GTP
In order to simulate the ATP and GTP manipulations used to generate the data in Fig. 14 of Rieke and Baylor (1998a), we scaled the parameters kRK3[ATP] and kG5[GTP] shown in Table I by 0.04 and 0.4, respectively, the same factors used in Rieke and Baylor (1998a)
(see Eqs.
and
, respectively).
Early Saturation Model
An early saturation model was implemented by restricting the amount of PDE locally available subsequent to a photoisomerization. R* shutoff was achieved in a single step by setting kA(1) to infinity, so that after the first phosphorylation, Arr-capping was automatic and instantaneous. This ensured that the only mechanism available to reduce SPR variability was the local depletion of PDE.
The local depletion of PDE was achieved by scaling kP1[PDE] (Eq. ) by the factor (1 - PDE*/PDEmax), where PDE* is the running count of PDE-subunit activations and PDEmax is the total number of PDE subunits locally available. In simulations, we found that in order to achieve empirical levels of SPR variability (CVarea) using only this local saturation mechanism, the number of PDEs locally available to R* (PDEmax) had to be restricted to
300. In contrast, for the full sequential phosphorylation model, PDEmax is set to infinity, so that there is no local saturation (i.e., the response to a single R* will not cause any local depletion of PDE).
Calcium Clamping
Ca2+ clamp was simulated by replacing the time-varying calcium variable, c, with the constant, cdark (so that ), which reduces the back-end deterministic equations (Eqs.
) to the single equation
![]() | (13) |
Data Analysis
In order to establish an empirical baseline for SPR variability, four analyses of SPR variability of Whitlock and Lamb's (1999) original data were performed. Two of these analyses were also presented in Whitlock and Lamb (1999)
, and two were new analyses not presented in their paper. The results analyzed here were those presented in their Figs. 13
, comprising responses of a toad rod to a set of 350 identical flashes that delivered on average
0.6 photoisomerizations per flash. The methods and procedures used in collecting Whitlock and Lamb's original electrophysiological data are described in detail in their paper (Lamb et al., 1986
; Whitlock and Lamb, 1999
).
|
The full data set consisted of 350 recording epochs of 10-s duration each, with the 20-ms (500 nm) flash presented 1 s into the recording epoch. Since the prestimulus portion of each recording comprises one tenth of the total epoch, we assume that at least 10% of each response should be at baseline levels, and extract the 10th percentile from the sorted current values for each record to form an estimate of the baseline drift over the experiment (1 h). This time-varying baseline drift waveform was low-pass filtered (-3 dB at 0.013 Hz) and subtracted from the raw responses.
The means of each of the 350, 1-s prestimulus intervals from the drift-corrected data were computed, and their distribution was fit with a Gaussian probability density function. Using the parameters of the Gaussian fit, responses whose pre-stimulus means were more than three standard deviations from the overall mean were excluded from further analysis ("bad-zero" records). 36 responses were excluded on this basis. The mean of the remaining 314 prestimulus means (a scalar DC offset for the entire experiment) was calculated and subtracted from each record, ensuring that all responses used for analysis start from zero, on average.
One additional response containing a large negative artifact was manually excluded, leaving 313 dim-flash response records for analysis.
Classification of Responses as Failures, SPRs, or Multiple-photon Responses (MPRs)
SPRs were identified from analysis of the histogram of response amplitudes. Amplitude was defined as the scaling factor providing the best (least squares) fit of the normalized mean response to each individual response. Fits were carried out over the time interval between t = 0 (time of the stimulus) and the time to peak (1.91 s) of the ensemble mean response (Field and Rieke, 2002). A sum of Gaussians model (Del Castillo and Katz, 1954
; Baylor et al., 1979
) was then fit to the histogram, yielding estimates of the mean number of photon absorptions, mean SPR amplitude, SPR amplitude variance, and additive background noise variance, for the underlying distribution. These four parameters were used to determine analytically a range of amplitudes where the probability that the response resulted from a single photon absorption was
50%. Responses were classified as either failures, SPRs, or MPRs on the basis of whether their amplitude fell below, within, or above this range. In simulations, this approach turned out to provide high sensitivity and positive predictive value (>95% for each).
Isolation of SPR Variance
To separate variability due to background noise from stochastic variability in the actual underlying response to photon absorption, we subtract the variance of responses classified as failures from the variance of those classified as SPRs. In any subsequent discussion of SPR variability, this correction is implied. This technique is used in the CV calculations for SPR amplitude and SPR area, as well as in the time-varying SPR variance plot described below.
SPR Area
In the first new analysis, we were interested in CVarea for the SPR (Field and Rieke, 2002). For each failure and SPR (identified by the method just described) we integrated the response from time zero until the end of the recording epoch at t = 9 s. CVarea was calculated as the square root of the variance difference between SPRs and failure areas, divided by the mean SPR area.
Squared Mean SPR Versus SPR Variance
Another new analysis was a comparison of the time-varying ensemble variance (SPR2) and squared mean (µSPR2), for identified SPRs only. This analysis was presented in Rieke and Baylor (1998a)
, and featured prominently in Field and Rieke (2002)
.
The other two variability measures we used had also been used in Whitlock and Lamb (1999): a histogram of dim-flash response amplitudes and calculation of CVampl, and a comparison of the variance of the full dim-flash ensemble (
dim2) with the squared mean (µdim2).
![]() |
RESULTS |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
In our analyses, we use four empirical measures of SPR variability/reproducibility to evaluate each model.
Variability of SPR Amplitude
A classical method of quantifying the variability of the SPR amplitude is to plot a histogram of the dim-flash amplitudes (Fig. 3 B) and fit a statistical (sum-of-Gaussians) model (black curve) to the histogram of dim-flash response amplitudes. As described in MATERIALS AND METHODS, we used the model fit to define amplitude boundaries with which to classify responses statistically as most likely resulting from zero, one, or multiple photon absorptions. The subset of the histogram resulting from SPRs identified by this method is shown as a red overlay. The CV of SPR amplitudes for the Whitlock and Lamb data estimated in this manner is 0.15, which is somewhat lower than some previously published values (0.20; Baylor et al., 1979
; Schneeweis and Schnapf, 1995
; Rieke and Baylor, 1998a
; Whitlock and Lamb, 1999
; Field and Rieke, 2002
). The relatively low CVampl is due to the definition of amplitude we used and our method of classification of responses. Larger values would be obtained using our classification scheme if amplitude were defined as in Whitlock and Lamb (1999)
, or if we had estimated CV from the parameters of the sum-of-Gaussians fit. The SPRs identified from the amplitude histogram were used in the analysis of SPR area (Fig. 3 D) and SPR variance over time (Fig. 3 E).
Dim-flash Variance Versus Dim-flash Squared Mean
A second measure of variability is based on comparing the light-evoked ensemble variance increase and the square of the ensemble mean (Schnapf, 1983; Rieke and Baylor, 1998a
). Fig. 3 C shows this classical comparison between the
dim2(red) and the µdim2 (blue) for the Whitlock and Lamb data. The variance increase was calculated by subtraction of the variance of the failures from the variance of all responses. As in Whitlock and Lamb (1999)
, our reanalysis of their data shows virtually the same relationship between
dim2 and µdim2.
This measure of SPR variability has some limitations. Whitlock and Lamb (1999) pointed out that, although a stereotyped elementary response will necessarily lead to a close match between the variance and squared mean responses, a close match between these does not necessarily imply a high degree of SPR reproducibility. This is because, for dim flashes, the variance is dominated by Poisson noise stemming from the quantal nature of light. Thus, the result in Fig. 3 C only shows that the SPR variability over time was not of sufficient magnitude to dominate the Poisson variability of the light stimulus.
In APPENDIX A we extend Whitlock and Lamb's insight about the limitations of this analysis, and provide a quantitative interpretation of the relationship between the squared mean response and variance of the response. Despite these caveats, this analysis can help to evaluate models since failure of proportionality between dim2 and µdim2 does imply nonstereotypic SPRs.
Variability of SPR Area
A third gauge of SPR variability is the variability of the response area for SPRs. The motivation for the area analysis comes from the fact that, as Field and Rieke (2002) pointed out, the area under the response waveform is expected to be a good gauge of overall response variability since it includes the effects of response variability throughout the entire response waveform. Apart from this pragmatic consideration, there is a theoretical basis for using the variability in response area (over other measures) in the analysis of mechanisms of SPR reproducibility. In DISCUSSION, we show why it is the CV of SPR area, and not the CV of either SPR amplitude or duration, that tracks the variability in integrated R* activity, and hence, can provide an estimate of the number of underlying R* inactivation steps.
A histogram of dim-flash response areas for the Whitlock and Lamb data is shown in Fig. 3 D. This analysis was not presented in the original Whitlock and Lamb (1999) paper. The red overlay in the histogram shows the subset of areas from the responses that had been classified as SPRs in Fig. 3 B. The CV for area for the SPRs derived in this fashion was 0.36, similar to values recently reported for mammalian rods (
0.3, Field and Rieke, 2002
). Using our sequential phosphorylation model, we found that after 1000 random additions of simulated noise to a single Monte-Carlo run of 350 trials (dim flashes), 95% of the estimates for CV of SPR area fell in the interval 0.30 to 0.44.
SPR Variance Versus SPR Squared Mean
The fourth measure of SPR reproducibility comes from an analysis of the time-dependent residual variability of the SPRs (Fig. 3 E). This analysis was introduced by Rieke and Baylor (1998a)(see their Fig. 5) and was featured prominently in a recent paper by Field and Rieke (2002)
. The panel depicts the time course of the noise-corrected SPR variance (red curve; SPR variance minus failure variance) and the square of the mean SPR (blue curve). The SPRs were categorized in the same manner as for the amplitude histogram in Fig. 3 B. As noted in Rieke and Baylor (1998a)
and Field and Rieke (2002)
,
SPR2 is approximately an order of magnitude smaller than µSPR2 until after the peak of the squared-mean response. In addition,
SPR2 peaks much later than µSPR2 (1.5 times), and is broader. Field and Rieke (2002)
emphasized that these features provide constraints on models of SPR reproducibility, aiding in the discrimination between models.
Aside from the variability/reproducibility measures described above, other data in the literature provide additional constraints on any candidate model.
Transduction Gain Manipulation by Alteration of Nucleotide Levels
Fig. 3, F and G, reproduces data from Fig. 14 of Rieke and Baylor (1998a) showing the effects of lowering transduction gain in the presence of normal (500 µM) or low (20 µM) ATP levels in dialyzed toad rod outer segments. With normal ATP (Fig. 3 F), lowering GTP by a factor of 2.5 caused the dim-flash response amplitude to decrease with no effect on the response kinetics. The kinetics were slowed by the GTP manipulation, however, if phosphorylation was substantially slowed by reducing ATP (Fig. 3 G). Rieke and Baylor (1998a)
interpreted these results to imply that neither phosphorylation nor arrestin-binding controlled the majority of rhodopsin's cumulative activity. We will show that this pattern of results under GTP and ATP manipulation does not rule out phosphorylation dominating R* inactivation, contrary to Rieke and Baylor's reasoning.
Genetic Knockout (KO) and Transgenic Manipulations
Implementation of a detailed stochastic model offers an opportunity to simulate genetic KO or transgenic substitution experiments targeting mechanisms expected to affect activation, inactivation or feedback in the phototransduction cascade. In the present study, we evaluate the ability of the model to predict the results from three recent experiments that genetically manipulated R* shutoff mechanisms, plus one experiment in which feedback synthesis of cGMP was disrupted. The results of these studies are reproduced in Fig. 3 H. The panel shows data obtained from mouse rods that had six major phosphorylation sites on the rhodopsin C-terminus disabled by substitution of alanine for the WT serine and threonine residues normally occurring at these sites (green: complete substitution, or CSM, rods; Mendez et al., 2000), that had Arr knocked out (red: Arr-/-; Xu et al., 1997
), or RK-/- (blue: Chen et al., 1999
), or GCAPs-/- (orange: Burns et al., 2002
). The WT responses from each of these studies are shown in Fig. 3 H as thin curves. The WT responses were scaled to the same relative peak amplitude (1.0), but the relationship to the corresponding genetically manipulated responses in each case was not altered.
Both the RK-/- and CSM responses rose at the same rate as the WT until 100 ms, and continued to rise until they reached a peak
2 times the WT peak amplitude. The SPRs in both CSM and RK-/- rods were step-like, shutting off abruptly at highly variable times. Histograms of the duration of SPRs from CSM and RK-/- rods were approximately exponential, with time constants of 5.1 and 3.3 s, respectively. The Arr-/- responses reached approximately the same peak amplitude as the WT responses, then exhibited a partial recovery with nearly normal kinetics. When viewed on a long time scale, the mean Arr-/- responses manifested the initial, relatively rapid recovery phase, and then settled into a long, slow recovery phase (with a mean recovery time constant = 51 s; results not shown here).4 The GCAPs-/- responses rose with
WT kinetics to peak at
4 times the WT amplitude at
300 ms.
The Sequential Phosphorylation Model Accounts for All the Data Examined
The sequential phosphorylation model dramatically reduces SPR variability relative to the expected behavior of a single-step inactivation model, bringing it to empirical values. The model yields the correct qualitative behavior in all tests, including reproduction of the Rieke and Baylor (1998a) transduction gain manipulation experiments (Fig. 4, F and G)
, as well as the response features from four genetic knockout and transgenic studies: three in rods in which R* inactivation mechanisms had been genetically disrupted (Xu et al., 1997
; Chen et al., 1999
; Mendez et al., 2000
), and one in which feedback synthesis of guanylate cyclase has been disabled (Burns et al., 2002
).
|
Low SPR Variability Matches Empirical Values
The relatively low response variability can be seen in raw model responses (Fig. 4 A), where the simulated response failures, SPRs and MPRs, can be distinguished. The distribution of amplitudes (Fig. 4 B) reproduces the behavior obtained empirically (Fig. 3 B). The subset of amplitudes classified as SPRs by the statistical method used in analysis of the Whitlock and Lamb data is shown as a red overlay. The solid blue curve in panel B depicts the histogram of the amplitudes of the actual SPRs (which, unlike in the electrophysiological data, can be identified perfectly in the simulated responses). This illustrates that our method of amplitude measurement and response classification identifies SPRs in the presence of (simulated) recording and photoreceptor noise with high statistical accuracy (sensitivity and positive predictive value both exceed 95%).
The CV for SPR amplitudes identified statistically is 0.16, nearly identical to the value we derived from the Whitlock and Lamb data (0.15; Fig. 3 B). The close match between the (blue) and the scaled
dim2(red) in Fig. 4 C indicates that SPR and MPR variability is low enough to allow the ensemble variance to be dominated by variability in the number of photon absorptions.
The distribution of SPR areas is Gaussian like (Fig. 4 D, red overlay) with a low CV of 0.38 that nearly matches the CVarea in the Whitlock and Lamb toad rod data (0.36; Fig. 3 D). This value is somewhat larger than the only other value reported in the literature (0.30 for mammalian CVarea; Field and Rieke, 2002
).
The CV for SPR area produced by the model (0.38) is close to the theoretical limit for an eight-step model (0.35; see Eq. , DISCUSSION), and reflects the contribution of approximately seven of the eight possible inactivation steps. The reason that all eight possible steps did not contribute is that, in our model, the RArr affinity increases monotonically with the number of phosphorylations and there is a finite probability of Arr capping as early as the first phosphorylation. Thus, Arr capping occurs before the maximum number of phosphorylations. On average, R* was capped when 6.1 of the 7 possible phosphorylations had occurred (Fig. 4 D, inset), corresponding to a total of
7 shutoff steps (including Arr-capping), which is in agreement with the observed CVarea (
=0.38). The fact that CVarea achieves the theoretical limit for the mean number of phosphorylations at shutoff indicates that the relative rates of G* activation and R* phosphorylation across phosphorylation states were nearly optimal.
The sequential phosphorylation model also reproduces the empirical relationship between µSPR2and SPR2 (compare panel E in Figs. 4 and 3). The
SPR2 waveform peaks at 1.6 times the time-to-peak of µSPR2, close to the peak shift when this analysis is applied to the Whitlock and Lamb data (1.5 times; Fig. 3 E).
Rieke and Baylor's Transduction Gain Manipulation Is Reproduced
The sequential phosphorylation model captures the transduction gain manipulation data from Fig. 14 of Rieke and Baylor (1998a). Fig. 4 F shows the results of a simulation under the control ATP condition. The decrease in GTP by a factor of 2.5 affects transduction gain (shown in inset) without significantly altering the kinetics of the response (shown by the larger normalized curves). However, when ATP is lowered by a factor of 2.5 (Fig. 4 G), the same GTP manipulation slowed the kinetics of the response in addition to decreasing the gain (inset). This finding is significant and will be discussed further in the DISCUSSION.
Qualitative Features of Transgenic and Genetic Knockout Data Are Reproduced
The sequential phosphorylation model reproduces the salient qualitative and quantitative features of the Arr-/-, CSM, RK-/- and GCAPs-/- genetic manipulations (Fig. 4 H). Allowing for the difference in timing between mammalian and amphibian rods, the model does well at capturing these features. These simulated responses represent predicted responses if the same genetic manipulations were performed in toad rods as were done in the mouse rods.
Analytical Expression for R* Activity and the Time Course of R* Inactivation
Because G-protein competes with RK and Arr, and because of reversibility in some of these reactions (Eqs. and
and
), G* activation rate (i.e., R* activity) is a complicated function of 10 front-end parameters. Thus, it is not feasible to identify a single parameter that controls the maximum G* activation rate or its dependence on n, the phosphorylation state of R*. However, we have solved the stochastic front-end equations to obtain an analytical expression for the mean steady-state activity of R* conditioned on the continued circulation of R* around the catalytic loop. This can be thought of as the rate of G* production at a late time (after the initial appearance of R* in a given phosphorylation state), averaged over all cases where neither capping nor further phosphorylation has yet occurred. The expression thus derived gives a theoretical rate of G* activation per R* (i.e.,
RG, in s-1) that is not physically observable, and therefore is not exactly the same as the chemical reaction rate of G* production. However, because R* typically circulates around the catalytic loop many times before being phosphorylated or capped, we can expect that the true maximum rate of G* production achieved will be close to this theoretical value. Although we do not show the full expression here, it can be approximated by the expression shown in Eq. 14.
Here, R* activity depends on 6 rate constants, plus the phosphorylation state of R*, n. Substituting for the rate constants in Table 1, and combining Eq. 4 with Eq. 14 we obtain
![]() | (15) |
Using our model, and simulations of genetic manipulations, we can illustrate the relative contributions of phosphorylation and Arr-binding to the inactivation process. To do this, we recorded the mean G-activation rate from 1,000 simulated SPRs under three conditions, the results of which are all illustrated in Fig. 5
: (a) WT (i.e., using our sequential phosphorylation model; blue curve), where both phosphorylation and Arr-capping contribute to R* inactivation. (b) CSM (green curve), depicting the case where both phosphorylation and Arr-quench are disabled. In the absence of phosphorylation, R* activity goes to a fixed, steady-state level close to the theoretical initial maximal R* activity of 146 G*/s (derived from Eq. 15, for n = 0). (c) Arr-/- (red curve) depicting the case where final Arr-quench is disabled, but phosphorylation is not. Thus, this curve depicts the mean reduction in R* activity due to phosphorylation per se.
|
These results are expected from general theoretical considerations if optimal variability reduction is to be achieved: when Arr quench of R* activity is only one final step out of eight potential inactivation steps, the multiple steps of phosphorylation preceding Arr-binding must dominate R* inactivation. Otherwise, the stochastic nature of the final Arr-binding would contribute a disproportionate amount of variability to the overall inactivation process.
![]() |
DISCUSSION |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Conclusions about the Sequential Phosphorylation Model
Previous authors have suggested that the processes of phosphorylation and arrestin binding may constitute multiple steps of R* inactivation that contribute to reducing SPR variability (e.g., Whitlock and Lamb, 1999; Gibson et al., 2000
; Mendez et al., 2000
; Field and Rieke, 2002
). Using a full stochastic, biochemical model, our simulations and analyses demonstrate quantitatively for the first time that multiple phosphorylation of R* (plus Arr-binding) can account for the SPR reproducibility observed in vertebrate rods. Sequential, phosphorylation-dependent ratcheting down of R* activity, and ratcheting up of inactivation rate, can change the distribution of R* lifetimes from approximately exponential (with CV
1) to a much less variable distribution, substantially reducing the variability of SPR amplitude, kinetics, and area below the levels that would otherwise occur if R* were shut off in a single, memoryless step (see APPENDIX B; compare Figs. 4 and 8).
We found that the sequential phosphorylation model accounts for the four measures of reproducibility and generates responses that exhibit the correct features in almost all details. The model also accounts for the signature qualitative features of four genetic experiments in mouse rods (see Fig. 4 H). The match to these data was achieved without any additional parameter adjustments, other than the simulation of the genetic manipulation. The model also reproduces the results of one of Rieke and Baylor's (1998a) transduction gain experiments that was thought to rule out phosphorylation (and/or Arr-binding) as the dominant mechanism in the deactivation of R*.
The Transduction Gain Manipulation Experiments of Rieke and Baylor (1998a) Do Not Rule Out R* Phosphorylation as the Dominant Mechanism Controlling R* Inactivation (and Hence SPR Recovery)
Rieke and Baylor (1998a) reasoned that the low GTP concentration (40% of control) in their transduction gain experiments would increase the amount of time a G-protein spent bound to R*. Because inactive G-protein, RK, and Arr are thought to compete in a mutually exclusive manner for R*, this should reduce the availability of R* for phosphorylation and arrestin capping. Consequently, if phosphorylation and Arr binding were responsible for a large portion of R* inactivation, one would expect both a reduction in amplitude and a slower response recovery in low GTP. They interpreted the result, that transduction gain was reduced and that response kinetics were unchanged, to mean that phosphorylation and arrestin binding do not make a major contribution to the inactivation kinetics of R*. They further tested this hypothesis by lowering ATP to slow phosphorylation so that it would be expected to control a significant fraction of R*'s cumulative activity. Under this condition, lowering GTP by the same amount did significantly slow the response (Fig. 3 G).
Our simulations reveal that Rieke and Baylor's nucleotide manipulation experiments do not, in fact, rule out R*'s inactivation being dominated by phosphorylation. In the sequential phosphorylation model, R* is inactivated entirely by phosphorylation followed by arrestin capping, with phosphorylation accounting for 66% of the total R* activity reduction, and the final Arr-quench accounting for the remaining 34% (Fig. 5). Yet, as discussed in the following section, when reversibility of some of the early reactions in the cascade is taken into account, we find that the model can reproduce the transduction gain results of Rieke and Baylor.
Reversibility of Early Reactions in the Cascade Allows the Sequential Phosphorylation Model to Reproduce the Transduction Gain Experiments of Rieke and Baylor (1998a)
Rieke and Baylor's reasoning would be correct if the reactions governing the interaction between R* and G-protein were unidirectional. The reasons that reversibility in these early reactions allows the sequential phosphorylation model to reproduce the Rieke and Baylor (1998a) nucleotide manipulation results may be understood in the following way.
Lowered GTP in the presence of normal ATP (Figs. 3 F and 4 F).
Reversibility of the interactions between R* and G-protein creates two possible outcomes for a G-protein molecule that has formed a complex with R*. It can dissociate from R* either in its active form, after replacement of the GDP by GTP (the catalytic route), or while still in the inactive GDP-bound form (see Eq. ).
Lowering the intracellular GTP concentration decreases the likelihood of G-proteins taking the catalytic route, so that a G·GDP molecule requires a greater number of R* encounters, on average, before it is activated. This, of course, reduces the rate of G·GTP formation, i.e. the gain of the response. In addition to decreasing the probability of catalysis, decreasing the GTP concentration will also increase the average time between the formation and dissociation of an RG·GDP complex (whether by the catalytic or non-catalytic pathway). The increase in lifetime of RG complexes will perforce slow phosphorylation and Arr-binding, due to the competitive nature of the interaction with R*.
Thus, Rieke and Baylor (1998a) were correct to conclude that, in addition to lowering the gain, the increase in lifetime of the RG complexes would lead to a slowing of inactivation kinetics, if phosphorylation and arrestin controlled R* shutoff. However, we found that if the rate constants in Eq.
are such that the net probability of proceeding along the noncatalytic path (i.e., back to free R* plus G·GDP) is made sufficiently high, the lifetime of the RG complex can be made to be much less sensitive to changes in GTP concentration than is the rate of catalysis. In particular, our parameters yield a probability of
0.3 that free R* will emerge from complex with G via the catalytic pathway, and a probability of 0.7 that it will proceed along the noncatalytic route. Under such conditions, our analysis shows that, with only phosphorylation and Arr binding to inactivate R*, a decrease in GTP concentration can cause a substantial decrease in gain but almost no slowing of kinetics (Fig. 4 F).
Lowered GTP in the presence of low ATP (Figs. 3 G and 4 G).
When the concentration of ATP is low, the slowing of kinetics elicited by the lowered GTP concentration will no longer be negligible. The lowered ATP concentration reduces the rate of phosphorylation of R* by RK, so that a given R* molecule will bind a substantially greater number of G-proteins in each of its phosphorylation states. Since each and every G-bound state of R* represents a timeout from the inactivation reactions, summing across a greater number of such timeouts magnifies the absolute difference in inactivation rate between the control and low GTP conditions.
Without reversibility.
When the RG reactions are assumed to be unidirectional, the RG complexes are unable to take the noncatalytic route, so that the GTP manipulation leading to a gain reduction will necessarily cause a comparable increase in the lifetime of the complex, and this in turn will slow the rate of phosphorylation and Arr capping, and hence the overall response kinetics. Thus, when reversibility is not included in these reactions, the only way of accounting for Rieke and Baylor's GTP/ATP manipulation data would be if R* inactivation occurred primarily by mechanisms that were not competitive with G-protein (e.g., Rieke and Baylor's putative multiple R* transitions).
For the sake of completeness we implemented such a modification to our stochastic model, and we found (as expected) that it failed to account for the results presented by Rieke and Baylor (1998a) in their Fig. 14 (unpublished data).
Stochastic Variability in Activation and Inactivation in the Cascade Is a Significant Source of Variability in Transduction Gain
Whitlock and Lamb (1999) analyzed variability in the rising phase of SPRs by fitting their discrete model to each response and extracting from each fit Lamb and Pugh's (1992)
amplification constant, A. They found that the CV of A was 0.13, and hypothesized that this variability in apparent transduction gain could be due to variability in the packing density of proteins on the 1,000 or so disc surfaces upon which the photoisomerizations occur across stimulus trials (see Calvert et al., 2001
).
We have applied this same analysis to the SPRs generated by the sequential phosphorylation model (by fitting Whitlock and Lamb's discrete model to the simulated SPRs in the absence of recording or photoreceptor noise), and find that, as in the physiological data, the simulated SPRs do not have stereotypical rising phasesthe CV for A is 0.19. Since the model does not have explicit implementation of packing density and diffusion kinetics, nor a mechanism simulating longitudinal variations in transduction gain (Schnapf, 1983), and since recording and photoreceptor noise were excluded from this analysis, this result demonstrates that a substantial proportion of the variability in photocurrent activation in the dim-flash regime may be due to the underlying stochastic variability in G-protein and PDE activation per se. The balance of contribution of these different mechanisms remains to be evaluated.
On the Number of Functional Phosphorylation Sites In Vivo
Biochemical studies in vitro indicate that seven or more phosphates per rhodopsin are incorporated (Kühn and Wilden, 1982; Wilden and Kühn, 1982
; Aton et al., 1984
; Thompson and Findlay, 1984
; Palczewski et al., 1991
; Wilden, 1995
). However, there are well-known technical difficulties in determining sufficiently quickly the number and identity of the rhodopsin sites phosphorylated under dim-flash conditions in vivo. For one thing, dephosphorylation events subsequent to physiologic quench can readily cause the number of phosphorylations to be underestimated (Ohguro et al., 1995
, 1996
; Hurley et al., 1998
; Kennedy et al., 2001
; Maeda et al., 2003
). Thus, in vivo studies have not been able to demonstrate the incorporation of more than one (Ohguro et al., 1995
, 1996
) or, at most, three or four phosphates per rhodopsin (Kennedy et al., 2001
). Moreover, these studies were conducted at intensities orders of magnitude above the single-photon level, where rhodopsin kinase may well be saturated, and therefore they may not accurately reflect the incorporation of phosphates under dim-flash conditions.
A recent study by Mendez et al. (2000) has provided evidence that, under physiological conditions, all the available phosphorylation sites are required for the normal kinetics of deactivation of rhodopsin. They used transgenic techniques to substitute various serine and threonine residues in the carboxy-terminal region with alanines, and they found that the rate of recovery of the dim-flash response increased systematically as the number of phosphorylation sites available was increased, implying that all the native phosphorylation sites were needed to support normal inactivation kinetics. While this in vivo study did not actually show that each of these sites underwent phosphorylation, it did show that the presence of all the native sites is important for normal response kinetics. Accordingly, the Mendez et al. (2000)
results are consistent with our assumption that all seven sites are available for phosphorylation, and furthermore we are not aware of any convincing evidence against the notion that, under single-photon conditions, as many as seven phosphorylations do indeed occur.
Three Phosphorylation Sites Are Not Sufficient to Account for SPR Reproducibility
Mendez et al. (2000) concluded that three phosphorylation sites are necessary for normal SPR reproducibility, and we agree with this, but they further concluded that more than three sites do not further improve SPR reproducibility, even though the rate of deactivation may be decreased. However, from theoretical considerations we can show that, in the absence of contributions from other mechanisms, the number of R* inactivation steps needed to achieve the observed SPR reproducibility must be greater than three. For example, our analyses of the sequential phosphorylation scheme show that if multiple, sequential phosphorylation (with Arr capping) is the mechanism that reduces SPR variability, then three phosphorylation sites cannot be sufficient to account for the observed variability in SPR area. The lowest CV for SPR area that three phosphorylations plus Arr-capping could support is 0.50, considerably higher than the empirical values (0.30, mammalian, Field and Rieke, 2002
; 0.36, amphibian, Fig. 3 D) or the value from simulations with our full model (0.38).
On the Relationship Between Variability of SPR Amplitude, Duration, and Area, and Its Dependence on the Rate-limiting Reactions in the Phototransduction Cascade
We will now show that there is not a straightforward relationship between the statistics of the random lifetime of a single activated rhodopsin molecule and the statistics of either dim-flash response amplitude or response kinetics alone. In particular, there is a theoretical tradeoff between variability in SPR amplitude and variability of SPR duration that depends on the kinetics of R* inactivation relative to the kinetics of reactions downstream to R*. Moreover, this tradeoff could be informative in determining what inactivation reactions might be rate-limiting in the recovery of rod responses in the dim-flash regime, since the relative CVs of SPR amplitude, duration, and area may indicate the extent to which R* is or is not rate-limiting in SPR recovery.
The variability of the photocurrent response in any biochemical kinetic scheme (in which non-linearities subsequent to PDE* production do not play a significant role) is determined by variation in the number of PDE* molecules produced, the time at which they are produced, and the variability in the lifetimes of individual PDE* molecules. In the sections below, we present a theoretical basis for the superiority of CV of SPR area as a measure of variability, and show why variability in SPR amplitude or kinetics alone are, in principle, less informative gauges of the variability.
The Theoretical Primacy of the Variability of SPR Area
Let N be a stochastic variable representing the cumulative R* activity during the SPR, i.e., the number of transducin molecules activated during R*'s lifetime, and hence the number of activated PDE* subunits. If SPR reproducibility were to derive purely from the regularization of R* activity by a sequence of n + 1 inactivation steps (e.g. n phosphorylations plus arrestin), then under optimal conditions, the lowest attainable CVN (i.e., N/µN) would be given by
![]() | (16) |
Below we show that CVarea of the SPR provides a good approximation of the CVN, and hence can be used to estimate a lower bound on the number of inactivation steps of R*.
Consider the net PDE* response to a flash to be the superposition of N discrete PDE* responses, each with unit amplitude, random onset time, and random (exponentially distributed) lifetime. The CV of the area under the PDE* waveform (APDE) will be determined entirely by the statistics of N. Assuming a linear dim-flash response subsequent to PDE* production, CVarea for the SPR will equal CV(APDE), which we show in APPENDIX C (Eq. C5) can be written as
![]() | (17) |
The first term under the square root can be thought of as the component of CVarea due to variability in the number of PDE* activations per R*. The second term under the square root represents variability, for a given mean number of PDE*, due to stochastic variation in individual PDE* lifetimes. In practice, we have empirical estimates for CVarea and µN, whereas the variability in N has not been directly measured, but can be estimated by a rearrangement of Eq. 17.
![]() | (18) |
Field and Rieke (2002), and our present analyses, find CVarea in the range of 0.260.35. Since the number of PDE* produced by a single activated R* during its lifetime is expected to be on the order of hundreds (e.g., Yee and Liebman, 1978
; Heck and Hofmann, 2001
), we expect the term µN-1 to be at least an order of magnitude smaller than CVarea2. Consequently, CVarea itself provides a close approximation to the variability in integrated R* activity, and therefore the minimum required number of inactivation steps (or, in this case, n phosphorylations) can be estimated from
![]() | (19) |
An Inherent Tradeoff Between Variability in SPR Amplitude and Kinetics Depends on What Rate-limits Recovery
To illustrate how variability in R* activity could manifest as variability in either the amplitude or duration of the response, we now consider two opposing, limiting-case scenarios using, for the sake of clarity, the following simple model: R* activity is assumed to be a rectangular pulse of fixed height and variable duration, and the SPR depends linearly and deterministically on this pulse. In the first scenario, the kinetics of the downstream cascade are assumed to be very fast compared to the typical R* inactivation rates. Here, the response is able to track the R* activity function with little lag, and the photocurrent response itself will also approximate a rectangular pulse. In the limiting-case, the SPR is just a scaled version of the R* activity, so there is no variation in SPR amplitude, and SPR duration will be equal to, and thus have the same variability as, the lifetime of R*.
In the second scenario, R* inactivation is very fast compared to downstream kinetics, so that with respect to the timescale of the response, R* activity is effectively an impulse. Here, the amplitude of the SPR will scale linearly with the random lifetime of R*. In particular, the coefficient of variation of the SPR amplitude will be the same as that for the random R* lifetime.
For example, if R* were shut off abruptly, following a series of x nonactivity-changing transitions with equally distributed waiting times (and there were no other mechanisms to reduce SPR variability), then the CV of R*'s lifetime would equal . In scenario 1, the CVs of SPR amplitude and duration would be zero and
, respectively, while in scenario 2, these values would be reversed. In contrast, the CV of SPR area (see Eq. 16) would equal
in both conditions; i.e., independent of the locus of the rate-limiting process in recovery.
In order to illustrate this tradeoff within the context of a full stochastic biochemical model, we ran Monte-Carlo simulations using seven different values of PDE, ranging from 0.75 s (R* recovery highly rate-limiting) to 12 s (PDE* recovery highly rate-limiting). The resulting CVs for SPR amplitude, duration, and area are shown in Fig. 6
A. The behavior of the full biochemical model clearly manifests the tradeoff behavior predicted from the limiting case scenarios: CVampl (open squares) is lowest when R* inactivation is rate-limiting, and monotonically increases as the downstream reactions become rate-limiting; conversely, CVdur (open triangles) is high and decreases monotonically as the rate-limitation switches from R* to downstream recovery kinetics. However, the corresponding CVarea (open circles) does not depend on what rate-limits recovery (CVarea is constant over all
PDE values).
|
The analysis above (see Eqs. 17 and 18) shows that for a given mean and variance of integrated R* activity, there will be a unique value of CVarea. However, as illustrated by the two tradeoff scenarios, neither the CV of SPR amplitude nor SPR duration are predetermined by the statistics of G* production per se; depending on the (average) speed of R* inactivation in relation to the kinetics of the downstream reactions, they may take on any value from near-zero to CVarea (which may be as high as 1 under some conditions, e.g., RK-/-).
Thus, without a priori knowledge of the relevant reaction rate constants, it is inappropriate to use either SPR amplitude or duration variability to estimate the minimum number of R* inactivation steps. This is illustrated in Fig. 6 B, where we have plotted the inferred number of R* inactivation steps corresponding to the CVs in Fig. 6 A. Only the number of steps inferred from the CVarea measures (open circles) match the actual number of R* inactivation steps used (X, dashed line).
In general, CVampl and CVdur measures will tend to overestimate the number of R* inactivation steps, and may do so severely. For example, a typical empirical CVampl of the single-photon response is 0.2 (Baylor et al., 1979
; Schnapf, 1983
; Baylor et al., 1984
; Schneeweis and Schnapf, 1995
; Rieke and Baylor, 1998a
; Whitlock and Lamb, 1999
; Field and Rieke, 2002
), which would correspond to an estimated minimum of 25 inactivation steps. The present analyses provide another striking example of the fallacy of this approach. The value of CVampl obtained from the sequential phosphorylation model was 0.16, which would correspond to an estimate of 39 inactivation steps, far greater than the number of steps actually used in the model. The value of CVarea, however, was 0.38, commensurate with the
7 inactivation steps that were used (at Arr-binding, an average of 6.1 phosphorylations occurred, plus Arr-binding = 7.1 steps, corresponding to a CVarea of 0.38; Fig. 4, B and D).
The Data and Analyses Suggest that the Kinetics of R* Inactivation and the Net Downstream Kinetics Cannot Differ Greatly
In practice, measures of CV of SPR amplitude and SPR duration have been reported to be of similar magnitudes, in the ranges 0.20.25 and 0.20.4, respectively (Baylor et al., 1979; Schnapf, 1983
; Baylor et al., 1984
; Schneeweis and Schnapf, 1995
; Rieke and Baylor, 1998a
; Whitlock and Lamb, 1999
; Field and Rieke, 2002
). These findings suggest that the rates of R* inactivation and downstream kinetics are not highly disparate in normal vertebrate rods.
What Kinds of Models Are Still Viable?
Four basic classes of models have been introduced to account for SPR reproducibility, and they are not mutually exclusive: late saturation, early (local) saturation, feedback, and multistep inactivation of R*. As we now discuss, none of the first three classes alone appear to be viable as the primary mechanism underlying SPR reproducibility.
Late Saturation (Depletion of a Species Subsequent to PDE* Production) Does Not Account for SPR Reproducibility
Single-photon responses might be rendered insensitive to variability in R* lifetime due to localized depletion of either cGMP or the number of available open cGMP-gated membrane channels. Some data from Field and Rieke (2002) argue against both of these possibilities in mammalian rods. They found that responses from local stimulation of the outer segment (12 µm, corresponding to
40 discs) were indistinguishable from responses obtained using diffuse illumination of the entire outer segment, and that, in both conditions, responses summed approximately linearly up to 3 R*. Photocurrent waveforms did not indicate any interaction between responses, even for two or three photo-isomerizations occurring in the same region of the outer segment. Since the PDE*s elicited by multiple photon absorptions within the local stimulation area were arguably competing for the same overall pool of cGMP (and cGMP-gated channels), these experiments indicated that activation of a single R* in normal (diffuse-light) conditions does not tax the limits of available channels or cGMP.
It is possible that the situation is different in rods of some lower vertebrates since Lamb et al. (1981) found that, even at very low intensities (
3 R*/flash), restricted stimuli generated smaller responses than diffuse stimuli eliciting an equal number of photoisomerizations in toad rods. However, Rieke and Baylor (1998a)
provided evidence against late saturation in amphibian rods by showing that manipulations that increased the amplitude of the SPR (e.g., lowering ATP to slow phosphorylation, or clamping Ca2+) did not generate responses with any evidence of saturation. Responses showed no clipping or other nonlinear behavior. In addition, when they lowered transduction gain by decreasing GTP (and, in some cases, also increasing GDP), the responses had the same shape as control responses when scaled to match in amplitude. If late saturation had played a significant role in shaping the responses in normal gain conditions, then the reduction in gain should have resulted in differently shaped responses (i.e., less distorted by saturation).
Early Saturation (Local Depletion of Transducin or PDE) Cannot be the Dominant Mechanism for SPR Reproducibility
In principle, responses could also be rendered insensitive to variability in integrated R* activity due to depletion of an early transduction intermediate (PDE or G-protein). This appears to be ruled out as a main mechanism by empirical data from at least four studies (Xu et al., 1997; Chen et al., 1999
; Mendez et al., 2000
; Rieke and Baylor, 1998a
), by the theoretical analyses in Field and Rieke (2002)
, and by our own analyses and simulations, as explained below.
First, a model in which SPR reproducibility is achieved solely by R* running out of G-protein (or G-proteins running out of PDE) in a local neighborhood on the disc surface would fail to account for the data from genetic experiments in which R* phosphorylation is disrupted. Assuming that in Arr-/-, RK-/-, and CSM rods rod structure is normal, and other proteins are unaffected by the genetic intervention, one would expect the same local saturation to occur in these rods as in WT rods. Thus, if saturation were responsible for the bulk of time-dependent effective R* activity reduction, one would expect the SPRs from Arr-/- (Xu et al., 1997), RK-/- (Chen et al., 1999
), and CSM (Mendez et al., 2000
) rods to exhibit recovery with kinetics similar to WT, contrary to what is observed (Fig. 3 H).
Second, the transduction gain manipulation experiments of Rieke and Baylor (1998a) argue against depletion of either G or PDE as the main variability-reducing mechanism for the same reasons as cited in the section above.
Third, in their theoretical analyses, Field and Rieke (2002) showed that local depletion of PDE predicts an incorrect relationship between the time to peak of the SPR variance and SPR mean (
SPR2 peaks too early), and predicts a waveform for SPR variance that is much less broad than the observed waveform.
Finally, we have implemented a stochastic biochemical model in which early saturation (local depletion of PDE) was the sole mechanism that reduced SPR variability (Fig. 7 ; see MATERIALS AND METHODS for model details). We systematically reduced the total pool of locally available PDE until the CVarea (0.32; Fig. 7 D) and CVampl (0.18; Fig. 7 B) matched empirical levels, as well as the CVs obtained from the full sequential phosphorylation model (Fig. 4). In addition, the variance of the ensemble dim-flash responses was close to the square of the ensemble mean (Fig. 7 C), and the individual simulated responses (Fig. 7 A) looked similar to the real data (Fig. 3 A).
|
The early saturation model also produces an inordinate number of small-amplitude, small-area SPRs that overlap substantially the distribution of response failures. This is not evident if the raw model responses are analyzed using the same statistical methods used to analyze real data (red overlays in Fig. 7, B and D). However, when we plot the distribution of veridical SPR amplitudes and areas (solid blue curves, Fig. 7, B and D), the substantial overlap with the corresponding distributions of response failures (gray distributions centered at zero on the abscissas of Fig. 7, B and D) can be seen. The variability of these veridical distributions is much higher than the CVs derived when the responses are analyzed as empirical data (CVampl = 0.45 vs. 0.18; CVarea = 0.57 vs. 0.32). More detailed discussion of this issue is presented in APPENDIX B, where we show the results of implementing a model with a 1-step inactivation of R* and no other variability reducing mechanisms.
Feedback Is Not the Agent that Regularizes Rhodopsin Lifetime or Activity
Feedback regulation of R* lifetime or activity, whether mediated by intracellular Ca2+ or some other messenger, seems to be ruled out as the main mechanism of SPR reproducibility by several lines of evidence. The evidence provided by Rieke and Baylor's (1998a) Ca2+-clamp and high-gain/low-gain experiments argue against any kind of feedback messenger (including Ca2+) being the sole mechanism for SPR reproducibility. Field and Rieke (2002)
examined the predictions of both Ca2+-mediated and non-Ca2+mediated feedback models. The non-Ca2+mediated model predicted SPR variance waveforms that were significantly more narrow than their data or than the predictions from a multi-step model. Moreover, the SPR variance waveforms peaked too soon relative to the SPR mean waveforms.
In addition to the experiments and analyses by Rieke and colleagues just cited, several other lines of evidence argue against Ca2+-mediated feedback as the main mechanism reducing R* lifetime variability. Briefly, some of the major arguments against a Ca2+-feedback regulation of SPR variability are:
(a) Burns et al. (2002) found that SPRs from GCAPs-/- mouse rods were nearly unaffected by application of the strong Ca2+ buffer, BAPTA. The logic here is that, in the absence of GCAPs, the main sites left for Ca2+ feedback are recoverin-RK and Ca2+ effects at the cyclic nucleotide gated channel (Hsu and Molday, 1993
, 1994
; Gordon et al., 1995
; Sagoo and Lagnado, 1996
; Haynes and Stotz, 1997
; Rebrik and Korenbrot, 1998
; Rebrik et al., 2000
). Since they found that BAPTA had no effect on the SPR amplitude or kinetics when cyclase feedback was taken out of the picture, the implication was that any Ca2+ feedback onto recoverin-RK in the GCAPs-/- rods was not strong enough to affect the responses in the single-photon regime. That such feedback does exist and is functional in GCAPs-/- rods was suggested by their finding that GCAPs-/- responses showed evidence of Ca2+-mediated effects at high intensities.
(b) At least two theoretical analyses argue against a Ca2+-feedback model. First, Field and Rieke (2002) simulated the effects of early Ca2+-feedback and found that the predicted
SPR2 waveform peaked too late relative to the µSPR2 in comparison to the empirical data or the predictions of their multi-step model. Second, we have simulated Ca2+-feedback onto R* lifetime via recoverin-RK with a deterministic version of the sequential phosphorylation model. We find that SPR amplitudes and kinetics are nearly unaffected by the Ca2+ changes induced in the single-photon regime, even when local Ca2+ changes are simulated by assuming that all the change in internal Ca2+ concentration is restricted to only 10% of the OS volume (Lamb et al., 1981
; Whitlock and Lamb, 1999
; unpublished data).
(c) The sequential phosphorylation model reproduces some key results from Whitlock and Lamb (1999) that were interpreted as evidence for Ca2+-mediated regulation of R* lifetime variability (via recoverin-RK). They fit a single-photon, discrete version of the Nikonov et al. (1998)
model to their SPR waveforms (Eqs. 2 and 3 in Whitlock and Lamb, 1999
) collected from control rods as well as from rods infused with BAPTA. The variability in their model was governed by a parameter, tlife, the stochastic lifetime of (all-or-none) R* activity. Whitlock and Lamb (1999)
found that both CVampl and the mean of the distribution of tlife values increased under BAPTA, suggesting (within the context of their model) that the slowing of Ca2+ dynamics had increased R* lifetime variability.
There is an alternative interpretation of the results of the Whitlock and Lamb's (1999) BAPTA experiments that derives from the tradeoff between amplitude and duration variability discussed in the previous section. The effect of BAPTA (or full Ca2+-clamp) would be to slow the net kinetics of the downstream reactions relative to R* inactivation by slowing (or stopping) the feedback synthesis of cGMP. Our tradeoff analysis would predict that, in addition to increasing mean tlife, BAPTA (or Ca2+-clamp) would increase CVampl, but would decrease the CV of tlife (a measure SPR duration). This is, in fact, what Whitlock and Lamb (1999)
observed. We confirmed this prediction in additional simulations with the sequential phosphorylation model (which has no Ca2+-feedback onto R* lifetime or activity). When full Ca2+-clamp was applied to the model, both CVampl and mean tlife increased, but the CV of tlife decreased (unpublished data).
The above interpretation provides a mechanistic explanation for Field and Rieke's (2002) caution that "...the inferred increase in rhodopsin lifetime [by Whitlock and Lamb] in BAPTA-loaded rods could have...been caused by a slowed Ca2+ feedback to cGMP synthesis (p. 743)."
(d) Some empirical evidence that Ca2+ is not affecting R* variability comes from Field and Rieke (2002), who found that the CV of SPR area in BAPTA-loaded rods was indistinguishable from that measured under control conditions. We have shown in the present paper that the CV of SPR area provides a gauge to underlying variability of integrated R* activity (Eqs. 17 and 18). If calcium feedback were contributing significantly to the reduction of variability of R* shutoff, then slowing down this feedback signal with BAPTA should have made single-photon responses more variable, and increased the CVarea, but it did not. The sequential phosphorylation model also reproduces this result (the CV of SPR area is unchanged under Ca2+-clamp; unpublished data).
Multistep R* Inactivation
Models based on multiple deactivation steps seem to the only viable remaining class of model.
The Rieke and Baylor (1998a) multistep model.
Rieke and Baylor (1998a) applied an array of experimental tests which, along with some modeling, they argued provided strong evidence against any form of saturation or feedback as the sole mechanism conferring SPR reproducibility. This left multistep rhodopsin shutoff as the only remaining viable mechanism. Although Rieke and Baylor (1998a)
acknowledged that the known processes of R* inactivationphosphorylation and Arr-cappingwould contribute to recovery and hence to the overall response variability, their experiments in which nucleotide levels were manipulated led them to conclude that these two processes could not dominate R* inactivation, a conclusion that our results directly dispute (see Fig. 4, F and G, and conclusions about the sequential phosphorylation model above). Their rejection of phosphorylation as the dominant inactivation mechanism (in conjunction with their arguments against saturation or feedback) is what led Rieke and Baylor (1998a)
to propose that 1020 phosphorylation-independent intrinsic conformational changes of R* occur to alter its catalytic activity, even though evidence for these steps in the known biochemistry and biophysics of rhodopsin is lacking. Because they did not have an explicit stochastic biochemical model, they could not test their hypotheses regarding the role of phosphorylation and Arr-binding, or the role of putative phosphorylation-independent R* states.
For the sake of completeness, we implemented a model that incorporated Rieke and Baylor's (1998a) putative conformational R* transitions within the context of a stochastic, biochemical model structure. In order for this model to capture dim-flash response behavior, as well as the Rieke and Baylor (1998a)
gain manipulation data as well as the knockout and transgenic data, we had to make a number of assumptions. First, we assumed that the putative rhodopsin transitions required a phosphorylation step, and would only occur after phosphorylation was complete. If phosphorylation and the transitions proceeded concurrently, the model could not reproduce the qualitative features of the RK-/- or CSM data; simulation of these responses would exhibit recovery not observed in the data. We implemented the model using only one phosphorylation step prior to the R* transitions, in order to capture the scenario suggested by Rieke and Baylor (1998b)
.
In order to be able to capture the Arr-/- data, we assumed that, in WT rods, after phosphorylation and the R* transitions were complete, there still remained 25% of R*'s activity. In this case, when Arr-/- was simulated, the model could exhibit the partial recovery to
50% baseline seen in the data. If, on the other hand, we assumed that in WT rods, each of the 1020 total inactivation steps (phosphorylation, R* transitions, and Arr included) ratcheted down R* activity by an equal amount, the simulated Arr-/- responses would recover nearly to baseline, in contrast to the data.
Under these constraints, we verified that a Rieke/Baylor-like multistep model can reproduce the observed properties of empirical SPRs, including the four variability measures, the transduction gain experiments of Rieke and Baylor (1998a), and the transgenic and KO data examined here (unpublished data). Despite this success, the model is not grounded in well-established biochemistry; there is no experimental evidence for a long sequence of postphosphorylation molecular transitions of rhodopsin. In contrast, the sequential phosphorylation model has substantial experimental (and now, theoretical) support.
Sequential phosphorylation model.
The sequential phosphorylation model is consistent with a broad array of electrophysiological and biochemical data (compare Figs. 3 and 4), and reproduces all the empirical results examined, including the nucleotide-manipulation experiment of Rieke and Baylor (1998a) that led them to reject such a scheme. There is thus no compelling need to invoke a long series of non-phosphorylation dependent state changes of R* in order to account for SPR reproducibility.
Concluding Remarks
One of the most difficult and enduring problems in the history of research on photoreceptor physiology has been to understand the single-photon response, and in particular its relatively high degree of reproducibility (Pugh, 1999; Pugh and Lamb, 2000
) despite the high variability inherent in all molecular reactions. Such reproducibility in the single-photon regime imposes strong constraints on any model of phototransduction.
We introduce a stochastic biochemical model, the sequential phosphorylation scheme, that is able to account for all aspects of observed SPR behavior examined here (Fig. 4), including four quantitative measures of reproducibility, and the transduction gain manipulations of Rieke and Baylor (1998a). The latter results and associated analyses call into question Rieke and Baylor's rejection of phosphorylation and arrestin-binding as the dominant mechanisms of R* inactivation, and suggest the testable biochemical hypothesis that one or more of the early reactions between R* and G-protein has a significant reverse component.
The paper also utilizes the analysis and simulation of knockout and transgenic experiments as an invaluable tool for testing and rejecting candidate models. The sequential phosphorylation model was able to reproduce salient response features from rods either with various R* inactivation mechanisms genetically disabled (RK-/-, CSM, Arr-/-), or with the mechanism of feedback synthesis of cGMP disrupted (GCAPs-/-). The simulations of the knockout and transgenic data illustrate vividly how some models may be able to reproduce the empirical levels of SPR reproducibility (CVarea), but may fail completely to reproduce fundamental qualitative features of a suite of extant (or simulated) knockout data, and can thus be rejected (e.g., see Fig. 7).
The sequential phosphorylation model has the virtue that it achieves all this while being based on a broad range of biochemical evidence. We believe that the model provides a solid foundation for further development and testing of models. It is clear that the phototransduction cascade is complicated enough so that stochastic modeling techniques of the sort used in this paper will play a central role in the evaluation of future hypotheses and data.
![]() | (14) |
![]() |
![]() |
APPENDIX A: |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
This idea can be stated in a quantitative manner as follows. The variance of the dim-flash response is related to the mean of the dim-flash response, the variance of the single-photon response, the mean of the single-photon response, and the mean number of photoisomerizations (), by
![]() | (A1) |
Note that the right side of the equation represents the scaling between µdim2 and dim2, and the term added to 1 in parentheses (
) is the square of the coefficient of variation of the single-photon response (i.e., CVSPR2
). Now with this equation in mind, the observation that the variance of the dim-flash response is approximately proportional to the square of the mean response over time means only the following: over the time course of the dim-flash response, the square of the coefficient of variation of the single photon response plus 1 does not vary much over time. It does not necessarily imply that the CVSPR(t) is small in magnitude.
Baylor et al. (1979), Schnapf (1983)
, and Rieke and Baylor (1998a)
assumed that the proportionality between
dim2 and µdim2 implied that the SPR had a stereotypic waveform, corresponding to an assumption that
SPR2=0 in the equation above (and, hence CVSPR2
will equal zero).
Whitlock and Lamb's (1999) insight may also be recast in terms of our equation: although the CVSPR2
must be substantially <1 to give the observed proportionality, the CVSPR(t) could still be significant and hence be associated with visibly non-stereotypic SPR waveforms (as their data analysis showed).
Our analysis extends the theoretical understanding of the problem by showing that, in principle, the CVSPR(t) could be of any magnitude as long as it was not time-varying, and the proportionality between dim2 and µdim2 would still hold.
In practice, data in the literature show that CVSPR(t) is time varying, but is small until late in the response (Field and Rieke, 2002; also, see Figs. 3 E and 4 E). Thus, the ensemble variance and ensemble squared mean tend to be closely matched, at least up to the peak of the squared ensemble mean response (Schnapf, 1983
; Schneeweis and Schnapf, 1995
; Xu et al., 1997
; Rieke and Baylor, 1998b
; Whitlock and Lamb, 1999
; Mendez et al., 2000
).
![]() |
APPENDIX B: |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
Fig. 8 shows the behavior of a stochastic model in which R* was forced to inactivate in a single step. This was achieved by setting kA(1) to infinity, so that after the first phosphorylation, Arr-capping was automatic and instantaneous. The result was that R* had fixed activity with approximately exponentially5 distributed lifetimes. In all other respects the model was identical to our implementation of the sequential phosphorylation model, including the simulation of recording and photoreceptor noise.
|
We know that the R* lifetimes associated with these responses must have an approximately exponential distribution, with a CV of 1. In this case, where R* activity is an all-or-none pulse, we have shown that R* lifetime is proportional to the number of PDE*s activated, and hence to SPR area as well (APPENDIX C). Consequently, the variability of the distribution of R* lifetimes must be reflected fully in the analogous distribution of SPR areas. We can thus predict in advance that, in the absence of noise, the CV of SPR area will be
1. This is confirmed by simulations without added noise, and where the SPRs were identified perfectly based on a priori knowledge from the model. The model generates an approximately exponential distribution of R* lifetimes with CV
1, which leads to the expected approximately exponential distribution of SPR area with CV of
1 (not depicted).
However, when recording and photoreceptor noise and failures are simulated, and the data are analyzed as was done for the Whitlock and Lamb (1999) data, the resulting distribution of SPR areas is not exponential, and CVarea is less than one (0.77; red overlay, Fig. 8 D). This result is not due to the added noise affecting the intrinsic variability of the SPRs. Rather, it is due to the empirical method of classification of responses (see MATERIALS AND METHODS) when we treat the model responses as real data in which one does not know in advance which responses are SPRs, MPRs, or failures. It turns out that this method leads to a biased classification of responses when one or more of the true underlying distributions of response amplitude categories is markedly asymmetric. The effect of this bias may be seen by comparing the amplitude histograms for the two cases, one in which the empirical classification method was used (red overlay in panel B), and one in which responses were classified perfectly (solid blue curve in B). The empirical classification identified SPRs incorrectly, so that the histogram of SPR amplitudes appeared Gaussian-like, with a low CV of 0.24. The CV of response areas calculated from this group of SPRs was also biased low (0.77; compare blue curve with red histogram in panel D). The group of actual SPR amplitudes (blue curve in B) is highly asymmetric and non-Gaussian, and when CVarea is calculated from these responses, it is close to the expected value of
1 (0.97; blue curve, panel D). The rounded onset of the area histogram is a result of the addition of noise to the model; without noise, the histogram is approximately exponential (not depicted).
We checked to see how much this factor might have biased the calculation of CVarea for the responses from the sequential phosphorylation model. In this case, the underlying distribution of response amplitudes was relatively symmetrical. Hence, the derived CVampl and CVarea did not change very much when perfect response identification was used instead of our empirical classification method (the CVarea increased modestly from 0.38 to 0.42). To the extent that the sequential phosphorylation model is a good representation of the underlying biology, this result reassures us that the empirical method did not significantly bias the estimate of CVarea from the Whitlock and Lamb data (Fig. 3 D).
From our understanding of the role of phosphorylation and Arr-binding in R* inactivation, we can predict the results of simulating the three genetic experiments (CSM, RK-/-, Arr-/-). In particular, since the single shutoff step is achieved by a single phosphorylation/Arr-binding event, we expect the simulated Arr-/- responses to parallel the simulated RK-/- and CSM responses. This, again, is confirmed by our modeling. Fig. 8 F shows that the simulated Arr-/- responses (red) peak at about twice the WT amplitude, and do not exhibit the partial recovery seen in the data (Fig. 3 H). In addition, although the simulated RK-/- and CSM responses are similar to the behavior seen in the data, a model in which R* terminates in a single step would not be able to reproduce any of the intermediate cases examined in the Mendez et al. (2000) experiments (where only one or two or three phosphorylation sites had been transgenically disabled).
Other behaviors cannot be predicted in advance. For example, due to the tradeoff that we now know exists between amplitude and duration variability, we cannot know a priori how the variability in R* lifetime will parse out between these two domains: we can only make a clear prediction about the variability of SPR area in the absence of noise, since this measure is not dependent on the relative "front-end" and "back-end" inactivation kinetics. In addition, the quantitative details of how response variability would manifest over time in the ensemble and SPR variances was not known in advance.
![]() |
APPENDIX C: |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|
![]() | (C1) |
![]() | (C2) |
When variability in N is taken into account, the overall expected value and variance of APDE are
![]() | (C3) |
![]() | (C4) |
Combining Eqs. C3 and C4, the CV of APDE can be derived
![]() | (C5) |
![]() |
FOOTNOTES |
---|
1 The choice to place all the phosphorylation dependence in the first of the RK reactions (kRK1, Eqs. 2a and 6) was made for computational simplicity. It is possible that other phosphorylation rate constants in Eqs.
![]() |
ACKNOWLEDGMENTS |
---|
Supported by NEI Grant EY115-13 and Smith-Kettlewell Eye Research Institute Fund 5809-0100 (R.D. Hamer), and an ARC Federation Fellowship to T.D. Lamb.
Olaf S. Andersen served as editor.
Submitted: 18 March 2003
Accepted: 19 August 2003
![]() |
REFERENCES |
---|
![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() ![]() |
---|