Phototransduction in retinal cones: Analysis of parameter importance

Colin Klaus; Giovanni Caruso; Vsevolod V. Gurevich; Heidi E. Hamm; Clint L. Makino; Emmanuele DiBenedetto

doi:10.1371/journal.pone.0258721

Abstract

In daylight, cone photoreceptors in the retina are responsible for the bulk of visual perception, yet compared to rods, far less is known quantitatively about their biochemistry. This is partly because it is hard to isolate and purify cone proteins. The issue is also complicated by the synergistic interaction of these parameters in producing systems biology outputs, such as photoresponse. Using a 3-D resolved, finite element model of cone outer segments, here we conducted a study of parameter significance using global sensitivity analysis, by Sobol indices, which was contextualized within the uncertainty surrounding these parameters in the available literature. The analysis showed that a subset of the parameters influencing the circulating dark current, such as the turnover rate of cGMP in the dark, may be most influential for variance with experimental flash response, while the shut-off rates of photoexcited rhodopsin and phosphodiesterase also exerted sizable effect. The activation rate of transducin by rhodopsin and the light-induced hydrolysis rate of cGMP exerted measurable effects as well but were estimated as relatively less significant. The results of this study depend on experimental ranges currently described in the literature and should be revised as these become better established. To that end, these findings may be used to prioritize parameters for measurement in future investigations.

Citation: Klaus C, Caruso G, Gurevich VV, Hamm HE, Makino CL, DiBenedetto E (2021) Phototransduction in retinal cones: Analysis of parameter importance. PLoS ONE 16(10): e0258721. https://doi.org/10.1371/journal.pone.0258721

Editor: Stephan C.F. Neuhauss, University Zürich, SWITZERLAND

Received: May 24, 2021; Accepted: October 5, 2021; Published: October 28, 2021

Copyright: © 2021 Klaus et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Data Availability: The software which underlies this study has been made freely available on GitHub (doi:10.5281/zenodo.4731156). The data for this study will be available on Dryad with citation Klaus, Colin (2021), Data for: Phototransduction in retinal cones: Analysis of parameter importance, Dryad, Dataset, https://doi.org/10.5061/dryad.6djh9w11c.

Funding: The authors CK, GC, VG, HH, CM, and ED are partially supported by National Science Foundation (https://nsf.gov/) grant DMS1812601. CK is supported at The Mathematical Biosciences Institute. MBI receives funding through the National Science Foundation Grant DMS1440386. VG receives support from the National Institutes of Health (https://www.nih.gov/) grant R01EY011500. HH receives support from the National Institutes of Health grants R01HL133923, R01NS111749. CM receives support from the National Institutes of Health grant R01EY031702. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Competing interests: The authors have declared that no competing interests exist.

Introduction

Diurnal vertebrates are mostly active in fairly high light levels, when visual perception is dominated by cone photoreceptors, which are significantly less light sensitive than rods [1]. This is particularly true for modern humans using artificial lights to enable cone vision. In fact, out of more than 90,000,000 photoreceptors in the human retina, approximately 100,000 cones concentrated in the virtually rod-free fovea are used for the tasks requiring high spatial acuity, such as reading [1]. The structure of rods, in which a narrow connecting cilium is located between the outer segment containing visual transduction machinery and the rest of the cell, made preparation of fairly pure rod outer segments feasible many decades ago. This was performed simply by breaking the cell at the connecting cilium and then using density gradients to separate outer segments from other retinal components (reviewed in [2]). A high abundance of rods, constituting more than 90% of all photoreceptors in the retinas of most model species, resulted in high yields of rod-specific proteins, which allowed their biochemical characterization. Extensive sets of rod phototransduction parameters are now available for several species, including mouse [3], which are very similar to human. Hence, most biochemically detailed models of visual transduction described rods [4–10]. Few model animals have a cone-dominated retina; ground squirrel and tree shrew are notable exceptions [11, 12]. As a consequence, the preparation of relatively pure cone outer segments suitable for biochemical characterization of transduction components is often not possible although progress in cone purification techniques has been made, for example, with carp [13, 14]. While the general structure of the signaling cascade and its shutoff mechanisms are similar in rods and cones, cones use many distinct phototransduction proteins including critical components of the cascade: photopigment, G-protein transducin, effector phosphodiesterase, cyclic nucleotide gated channel (reviewed in [15]). Thus, one cannot rely on known rod parameters to model cone responses.

At present the characterization of cone-specific proteins is woefully incomplete, and there is no model species for which all the values necessary for quantitative modeling have been experimentally determined. Here we analyzed the impact of variation of individual parameters on the predictions of our space-resolved model of cone signaling [16]. Parameter ranges were coarsely suggested by interspecies measurements and, in the face of such uncertainty, we used global sensitivity analysis (GSA) to identify parameters that are most influential as measured by variance in model output and therefore first priority for future investigations. It should be emphasized that parameter influence was measured with respect to literature uncertainty which may not yet coincide with that of biological significance. In contrast to local sensitivity analysis, which depends on a particular choice of fitting parameters, GSA measures parameter importance when these vary over prescribed ranges. To the authors’ knowledge this is the first time GSA, by Sobol indices [17, 18], has been applied to phototransduction. Our mathematical model has proved particularly useful for predicting the functional behavior of rod photoreceptors [6, 8, 19–21]. Now GSA with a revised version of this model, adapted to cones [16], has evidenced that a subset of parameters which determine the dark-adapted state of the photoreceptor are the most influential, over the presented uncertainty ranges, for reproducing experimental cone flash responses. The most significant parameter was found to be the turnover rate of cGMP in the dark, β_dark. The second most significant was a newly derived parameter that quantifies how nearly the photoreceptor is biochemically tuned towards the impossibility of a dark-adapted state, a_min. This latter quantity has strong physical meaning and originates from the need for balance between cGMP synthesis by guanylyl cyclase and its hydrolysis by PDE in the dark, concurrently with the required balance between Ca²⁺ influx-efflux. In particular, sigmoidal Hill and Michaelis-Menten expressions create the possibility for maximal and minimal synthesis or hydrolysis rates (similarly for Ca²⁺ influx-efflux rates) to overwhelm the other if biochemical parameters are not properly constrained (See Eqs 9 and 10). While this finding is retrospectively intuitive, this may be the first time this constraint has been presented as potentially significant for phototransduction so that the biological range of a_min should be better quantified. Additional parameters which also affect the dark current and were found to be significant were the saturated exchanger current, , and the maximum synthesis rate of cGMP by GCAP-activated guanylyl cyclase, α_max. The shut-off rates of light-activated rhodopsin, k_R, and phosphodiesterase, k_E, were also significant and influential. The influence of the activation rate of transducin by photoexcited rhodopsin, ν_RG, and the hydrolytic efficiency of the activated PDE dimer, k_cat/K_m, upon model variance with experiment was also appreciable but comparatively less.

Materials and methods

Initial ranges for model parameters were based on values reported for several different species. While working within a single animal model is highly preferred, there is no complete, experimentally determined parameter set in the literature for any one species. To mitigate this situation, ranges were chosen to contain values that reproduced trends in experimental flash response without violating known parameter constraints. This selection was performed by stochastic optimization and a Metropolis-Hastings Markov Chain Monte Carlo (MCMC) random walk [22–24] over parameter space for a stationary distribution whose modes minimized root-mean-square (rms) error. Next GSA was performed to weight parameters by their relative importance. This was technically performed by the method of Sobol indices [17, 18, 25, 26]. In particular, this method demonstrated which parameters contributed the most to variance, when they were varied over their prescribed uncertainty ranges, with the examined flash response. For completeness, local sensitivity analysis is also reported.

Proposing first ranges for parameter values of the cone outer segment

Geometry.

Table 1 reports the geometric parameters for the mouse cone outer segment with experimental ranges from [27]. Some features of the sliver, the cytoplasmic volume that surrounds the closed section of disks and is encased by plasma membrane, were not known for mouse, so values were taken from striped bass [28] and frog [29].

Download:

Table 1. Geometric parameters for the mouse cone outer segment (COS).

https://doi.org/10.1371/journal.pone.0258721.t001

G-protein/effector cascade.

Parameters for disc membrane proteins, with their experimental ranges from the literature, were collected in Table 2. The volumic concentrations of [R]_vol = 3 mM, [G]_vol = 0.21 mM, and [PDE]_vol = 15 μM were reported for frog and carp in [30]. In [31] it was measured that [G]_vol in carp cone was 0.6x that in rod. These values were converted into surface densities through multiplication by the volume-to-surface conversion factor [6, 19]. For the geometric parameters in Table 1, η = 5.5 nm. The surface density of PDE on a cone disc differs considerably from the rod range of [500 μM⁻², 1000 μM⁻²]. For k_R, estimates for mouse rod were given [3]. The parameter range for k_E was based on the mouse rod value [3] and the observation that k_E for cones is ∼ 2.3x the mouse rod value [32]. This range is similar to the value k_E = 18.5 s⁻¹ obtained by numerical fit and reported for striped bass [33], while in carp a GTP hydrolysis rate as much as ∼25x higher than rod was reported [31]. This higher estimate led to an alternative upper bound of k_E ∼ 150 s⁻¹.

Download:

Table 2. G-protein and effector-related parameters.

https://doi.org/10.1371/journal.pone.0258721.t002

Catalytic activity.

Table 3 collected the activation and hydrolysis parameters used in the model. These parameters inform the model through the equations below. As before, is the volume-surface conversion factor. (1) (2) (3) (4) (5)

Download:

Table 3. Activation and hydrolysis parameters.

https://doi.org/10.1371/journal.pone.0258721.t003

The ranges of B_cG and for mouse rod are given above [3]. The parameter value for k_cat/K_m used by [33] in the analysis of striped bass cone was derived from measurements from frog rod [34]. N_Av is Avogadro’s number. The PDE dark activity was reported as [0.3%, 4.7%] of the maximal PDE activity, 18 cGMP/R*/s, in carp cone [14]. The dark rate of cGMP hydrolysis by PDE, β_dark, was estimated by multiplying these factors by the concentration of R* in carp, [R*] = 3 mM [30], and then setting this resulting value to β_dark [cG]_dark. The value [cG]_dark = 3 μM was also taken (see Table 5 for a discussion of [cG]_dark).

The activation rate ν_RG was experimentally measured for carp in [31, 35] as ν_RG ≈ 33 s⁻¹. Values for the rate of PDE activation by R*, estimated from modeling of striped bass and goldfish cones [33, 36, 37], were considered in order to estimate the range for ν_RG. Finally, the effectiveness of transducin in carp cone to activate PDE was reported as one-tenth of its effectiveness in rod [14, 30]. This led to the estimate k_GE = 0.1 μm² s⁻¹ since in amphibian rod this value was reported as unity [38, 39].

Guanylyl cyclase activity.

Parameters for guanylyl cyclase activity with their reported experimental ranges were given in Table 4. These parameters inform the model through the equation (6)

Download:

Table 4. Guanylyl cyclase (GC) activity parameters.

α_max is the maximum cGMP synthesis rate in the absence of Ca²⁺ and α_min is the synthesis rate at saturating Ca²⁺ concentration. These activities were measured in the absence of bicarbonate.

https://doi.org/10.1371/journal.pone.0258721.t004

The ratio α_max/α_min = ∞ was effectively adopted in [37] by their choice of mathematical model, since cyclase activity vanishes as [Ca²⁺] → ∞ in that framework. The measurements for α_max were given for striped bass [37] and implicitly for carp in [30]. There the volumic concentrations of guanylyl cyclase and the Ca²⁺ sensing GCAPS were reported as [GC]_vol = 72 μM and [GCAP]_vol = 33 μM. The value α_min in carp cone was estimated as α_min = (72 μM GC)(1.7 cGMP formed/1GC/s) where the latter was the reported activity rate [30]. Similarly α_max was estimated assuming that all available GCAP was bound to GC and its reported activity rate: Together, these considerations estimated that the ratio α_max/α_min ∼ 2. For mouse cones, GCAP1 normally dominates GCAP2 for binding of GC [40]. GCAP1 binding affinities were reported in [41].

Ionic current parameters.

Table 5 reports the parameters for ionic current with experimental ranges given in the literature for striped bass primarily [33, 37, 42–44]. These inform the model through the equations (7)

Download:

Table 5. Parameters for ionic currents of cone outer segments.

https://doi.org/10.1371/journal.pone.0258721.t005

Here stands for the Faraday constant and Σ_S is the surface area over which ion channels are distributed. Unless otherwise stated, channels were uniformly distributed over the sliver’s lateral boundary. The value K_ex = 19 nM was estimated numerically by [42] and is more than an order of magnitude smaller than the range 0.9 − 1.6 μM given for mouse rod. It is also evident from Fig 1 that, for the mouse cone examined in [45], J_dark ∼ 25.75 − 27 pA while circulating dark current for striped bass was reported as J_dark = 27.3 ± 10.5 pA [37, 42]. Finally, K_cG = 20 μM was taken from mouse rod [3]. For striped bass, K_cG was reported in [37, 43] to depend sigmoidally on Ca²⁺ and across the range 105 μM − 316 μM. Since this range was used in tandem with a high, computed dark cGMP concentration, [cG]_dark = 27.9 μM, we used the mouse rod value of [3] as other authors have reported [cG]_dark to be similar between rods and cones [30, 46]. The Ca²⁺-buffer, B_Ca²⁺, is reported for mouse rod; however, this value is similar to that in [37] except there a functional relationship is used instead of a single value.

Download:

Fig 1. Modeling flash responses of a mouse cone.

Black traces are flash responses recorded from a Gnat1^−/− mouse cone for an estimated range of 40—6,000 photoisomerizations with a half-maximal intensity that produced 940 photoisomerizations [45]. The colored traces show model predictions for indicated flash intensities. All simulations use the single set of parameters in Table 10. This set was found by stochastically minimizing the rms error between the model and experimental response solely for the 940-photoisomerizations trace while constraining the parameter values to satisfy known experimental constraints.

https://doi.org/10.1371/journal.pone.0258721.g001

Diffusion coefficients.

The diffusion coefficients for the G-protein machinery and second messengers (Table 6) were taken as those reported for mouse rod [3]. However, membrane diffusion is likely to be slower in cones than in rods because of the lower unsaturated fatty acid content in the cone outer segments [47]. Indeed, the diffusion coefficient for visual pigment was found to be approximately 16% lower in catfish cones than in amphibian red rods [48].

Download:

Table 6. Diffusion coefficients for cascade components in the membrane and for second messengers in the cytoplasm.

https://doi.org/10.1371/journal.pone.0258721.t006

Choice of ranges for Markov Chain Monte Carlo.

The experimental ranges in Tables 3–6 are the basis for a Metropolis-Hastings Markov Chain Monte Carlo optimization scheme for estimating parameter values from experimental data. The MCMC stationary distribution penalized some parameter values if they did not belong to the ranges in Table 8 while other parameters were fully constrained to belong to the ranges in Table 7.

Download:

Table 7. Strict constraints imposed on many parameters that influence the dark current and the values derived from those parameters.

The range for dark current was centered about the experimental flash response of [45], also shown in Fig 1.

https://doi.org/10.1371/journal.pone.0258721.t007

Download:

Table 8. Parameter ranges within which there was no penalty imposed by the Metropolis-Hastings search.

The geometric parameters were held fixed at R_b = 0.6 μm, R_t = 0.4 μm, H = 13.4 μm, ϵ₀ = 16.8 nm, ν = 0.65, ω₀ = π, and the Ca²⁺ diffusion coefficient was held constant at D_Ca²⁺ = 15 μm² s⁻¹. [R]_σ was omitted since flash response depended only on the initial population of R*, and was otherwise independent of its surface density.

https://doi.org/10.1371/journal.pone.0258721.t008

Choice of ranges for global sensitivity analysis

The parameter ranges used for Markov Chain Monte Carlo in Tables 7 and 8 were slightly revised, once the best fit was obtained, to reflect the optimized parameter values (Table 10). The Sobol, global sensitivity analysis was conducted using ranges in Table 9, mostly enlarged from Table 8.

Download:

Table 9. Ranges over which parameters were varied when conducting Sobol sensitivity analysis.

https://doi.org/10.1371/journal.pone.0258721.t009

Download:

Table 10. Parameter values for a mouse cone found by minimizing the rms error between experiment and model predictions for a flash producing 940 photoisomerizations according to the Metropolis-Hastings random walk.

https://doi.org/10.1371/journal.pone.0258721.t010

Numerical implementation by finite elements

We used the homogenized, finite element model of cone phototransduction published in [16] with adjustments to handle the case of continuously distributed illumination such as for cone’s native bright light settings. The software library used in this work is available on GitHub [49], and the simulation data produced by it are available through Dryad repository [50]. See also S1 Appendix.

Identifying a consistent set of parameters

To choose a consistent set of parameters in murine cones, the root-mean-square (rms) error between the experimental flash response reported in [45] and model prediction was minimized for 940 photoisomerizations (Fig 1). This flash intensity was chosen because it was the value reported by [45] for attaining a half-maximal response. Some parameter values were constrained to a subset, K, of parameter space defined by the ranges in Table 7. Since intervals for many biochemical parameters were unknown, exploration of parameter space for other parameter values outside of anticipated ranges was only penalized rather than prohibited (Table 8). Let e(x) stand for the error between the experimental flash response and the model prediction for a choice of parameter values x. A Metropolis-Hastings random-walk was used to construct a Markov chain with stationary probability distribution over parameter space with density π(x) satisfying the proportionality relation (8) for user-selected values γ, β > 0. (See also S1 Appendix). The interval I_i was the range reported in Table 8 for the i^th parameter. Sampling the Markov chain then preferentially explored regions of parameter space with high probability [22, 23]. By construction this was where the error, e(x), was small.

Choice of functionals for sensitivity analysis

Sensitivity analysis was performed by measuring the change in certain quantities of the flash response, hereafter called functionals, that were elicited by changes in the parameter values. The following functionals were used to measure various aspects of the photoresponse: I_act was defined by first identifying the current drop for times t in the interval [0, 10 ms] with the quadratic function and then setting I_act = A. E_act was similar but for the total amount of E* across the outer segment at any instant of time, where E* = 2[PDE*] is the concentration of active catalytic subunits of PDE. Each PDE had two such subunits, which in the model could be activated independently. I_drop was defined as the maximum drop in flash response attained at any instant of time expressed in proportion to the dark current. E_peak was the total amount of E* aggregated across the outer segment at any instant of time. A functional for the recovery phase of E* was also included. E_rec was found by fitting the current drop over the time interval [0.135 s, 0.5 s] with an exponential ce^−αt and then taking E_rec = α. The analog for current drop was not included because an exponential does not correctly model cone overshoot. J_dark and T_peak were more simply the circulating dark current and time-to-peak of the current drop. J_over was the greatest magnitude of overshoot (nonpositive) current values exhibited by the flash response expressed in proportion to the dark current. In addition to these we also introduced an error functional, denoted L². This functional quantified the rms error between the experimental flash response reported by [45] and the model prediction for the flash producing 940 isomerizations.

Local sensitivity analysis

Once a choice of parameter values was made, , the local sensitivity of the model for a functional y at x* was computed with a gradient-based method [3] by the quantity Q_i measures the instantaneous relative change in the functional y with respect to the relative change in parameter x_i based at the point x*.

Global sensitivity analysis

In contrast to local sensitivity analysis, which is based on a specific choice of parameter values, global sensitivity analysis examines the statistical variance of a functional, when its input parameters vary independently, in the statistical sense, over given ranges. Using the GSA method of Sobol indices [17, 18, 51], a functional was decomposed into components whose own respective variances quantified interactions between subgroups of parameters. This analysis precisely determined the percentage of the functional’s total statistical variance which was due to interactions between the freely chosen, prefixed subset of parameters. In this work, we considered the indices S_y and where y was either a single parameter x_i or was a collection of two parameters x_i, x_j. The index S_y was the percentage of total variance that could be explained by the parameters in y alone while was the percentage of total variance that could be explained using only parameters not in y. Equivalently, was the fraction of variance that was due to the parameters in y interacting with all other parameters when they were randomly varied at once. It follows that . It further holds that as S_y approaches 1, the functional tends to only depend on the parameters in y. Similarly as approaches 0, the functional tends to not depend on the parameters in y. Owing to the inherent randomness in parameters that define biological systems, GSA is an apt tool for quantifying the significance of parameters upon model output. (See also S1 Appendix).

The Sobol indices were estimated by Monte Carlo integration through a scheme first presented in [25]. See also [26]. In total, 6.8 million Monte Carlo trials were performed at the Ohio Supercomputer Center [52], so that each statistic needed for Sobol analysis was estimated using 100, 000 trials. Confidence intervals were constructed assuming that sufficiently many trials had been conducted so that the corresponding sample means were normally distributed. The validity of this assumption was investigated using a bootstrap sampling procedure (e.g. [53]) of the empirical samples. This procedure along with convergence rates and confidence intervals are given in S1 Appendix. Since the Sobol indices were ratios of two Monte Carlo estimated quantities, after the numerators and denominators were given 95% confidence intervals, the indices were then estimated up to 90% confidence using also the order-preserving properties of division.

Renormalizing α_min to a_min.

The Sobol method required that the considered parameters be independently sampled. However, the existence of a dark-adapted steady state required a first-principles balance of fluxes between the CNG channel and exchanger currents as well as a balance between cGMP synthesis and hydrolysis in the dark: (9) (10) In particular, Eqs 9 and 10 could not be satisfied for arbitrary parameter choices. A criterion for the existence of the two dark concentrations [Ca²⁺]_dark and [cG]_dark was the constraint [3] (11)

This constraint constituted a dependence between the constituent parameters, and a change of variables was required to reestablish a parameter set compatible with independent sampling. Once ranges for all parameters except α_min were fixed, the monotonicity of Eq 11 ensured that a steady state would exist so long as α_min belonged to the interval [0, ξ]. Here where χ satisfies the equality (12) We defined the unitless parameter (13) No further information on the biological uncertainty of this parameter was assumed, and so it was taken to be uniformly distributed over the interval [0, 1]. The value of a_min represents the relative position of α_min within the interval [0, ξ] over which the dark steady state exists. In particular, as a_min approaches either 0 or 1, α_min approaches extremes where the existence of a steady state becomes impossible. Through this change of variables, the new collection of parameters with α_min replaced by a_min could be sampled independently for Sobol analysis and the model evaluated through the substitution .

Results

Consistent parameter set

The Markov chain was sampled to obtain the best fit with values shown in Table 10 for a Gnat1^−/− cone (Fig 1). Knockout of rod transducin in the Gnat1^−/− mouse rendered rods nonfunctional [54], thereby precluding any intrusion of rod responses in the mouse cone recordings (e.g., [55]).

This fit was performed for a response to a flash producing 940 photoisomerizations, which was reported by [45] to be at half-maximal flash response. Parameter values were subject to the constraints in Table 7 and penalized if they fell outside the ranges in Table 8. Because the distribution, Eq 8, was chosen to prioritize exploration of quality fitting regions of parameter space, and was not, for example, used in a Bayesian framework [56, 57], the chain was not necessarily sampled until a stationary distribution was attained. For convenience, the deposited software library may be used to continue the MCMC chain when desired. Modeling results for striped bass cone flash responses subject to the same biochemical ranges of Table 9 are included in S1 Appendix as supplementary materials.

Local sensitivity of parameters

Table 11 reports the local sensitivity of functionals y(x) at x = x* with respect to the coordinate x_i as described in Methods. The partial derivatives were numerically computed by increasing the x_i parameter by a relative 5% when forming the numerical difference quotient.

Download:

Table 11. Local sensitivity indices for a flash of 940 isomerizations uniformly distributed throughout the outer segment.

https://doi.org/10.1371/journal.pone.0258721.t011

Global sensitivity of parameters

Owing to the number of model parameters and the uncertainty in their experimental values, global sensitivity analysis was performed to assess which uncertainties had the biggest effect upon model output. Parameter values were uniformly sampled across the ranges given in Table 9. For each choice of parameters, the model simulated the response to the flash of 940 isomerizations. Associated functionals that quantified individual components of the cascade were computed along with their statistical variance. The Sobol method then assigned percentages of that variance that were due to individual and subcollections of parameters. This allowed the parameters to be ranked by their effect across all ranges of uncertainty (Table 9). Hereafter, these percentages are reported as numbers in the interval [0, 1].

Findings.

Tables 12 and 13 and Figs 2–4 report global sensitivity of functionals using the Sobol method described in Methods (100,000 samples per estimate) for ranges of parameters in Table 9. The eight most influential (pairs of) parameters are shown for the given functional. The Monte Carlo Sobol trials along with additional convergence rates and confidence intervals are available at Dryad [50].

Download:

Fig 2. Sobol indices for functionals quantifying E*.

The dot at the center of a circle is the Sobol index obtained by Monte Carlo evaluation (100,00 samples). The blue bars define a 90% confidence interval. Plots show the eight most influential parameters ordered from most significant to least significant. (a) Pairwise sensitivity indices for E* activation. (b) Pairwise sensitivity indices for E* recovery. (c) Pairwise sensitivity indices for peak E* production. (d) Total sensitivity indices for peak E* production.

https://doi.org/10.1371/journal.pone.0258721.g002

Download:

Fig 3. Sobol indices for functionals quantifying the drop in current due to flash response, the rms error between simulation and experiment, and the dark current.

Blue bars define a 90% confidence interval (100,000 samples). Plots show the eight most influential parameters. Confidence intervals could be off-center of the estimated Sobol index, because the indices were ratios of two Monte Carlo estimated quantities. (a) Single sensitivity indices for the current drop. (b) Total sensitivity indices for the current drop. (c) Total sensitivity indices for the rms error between model prediction and experiment. (d) Total sensitivity indices for the circulating dark current.

https://doi.org/10.1371/journal.pone.0258721.g003

Download:

Fig 4. Sobol indices for functionals quantifying the time-to-peak of the current drop and its overshoot.

Blue bars define a 90% confidence interval (100,000 samples). Plots show the eight most influential parameters. Confidence intervals could be off-center of the estimated Sobol index, because the indices were ratios of two Monte Carlo estimated quantities. (a) Single sensitivity indices for the time-to-peak. (b) Total sensitivity indices for the time-to-peak. (c) Single sensitivity indices for the overshoot. (d) Total sensitivity indices for the overshoot.

https://doi.org/10.1371/journal.pone.0258721.g004

Download:

Table 12. Single sensitivity indices.

As an index approached 1, the functional became dependent only on that parameter. Most values shown are close to 0, which indicated that nonlinear interactions between parameters dominated the cone flash response. While the theoretical value of the Sobol index must fall in the interval [0, 1], small negative values sometimes occurred above as an artifact of the Monte Carlo approximation. These should be regarded as approximately 0. Confidence intervals are given in S1 Appendix.

https://doi.org/10.1371/journal.pone.0258721.t012

Download:

Table 13. Total sensitivity indices.

As an index approached 0, the functional became essentially independent of that parameter. A large index value indicated that the considered parameter contributed to significant nonlinear interactions with other model parameters, so that ignoring it would amount to that index’s loss, as a proportion, of the total variance. Some parameters that were negligible, e.g. m_cyc, may have been so because their prescribed uncertainties were smaller than other parameters. Confidence intervals are given in S1 Appendix.

https://doi.org/10.1371/journal.pone.0258721.t013

Sensitivity indices for E* production are given in Fig 2. Monte Carlo estimated that k_R and k_E alone accounted for approximately 70% of the variance in peak E* production while ν_RG also played an influential but lesser role (Fig 2c). Similarly, total sensitivity indices found that these parameters were the most influential for peak E* production (Fig 2d). The insignificance of [PDE]_σ for peak E* production indicated that shut-off of R* and E* happened rapidly, before exhausting the population of available E. However, [PDE]_σ and k_GE were influential for the activation and recovery time courses (Fig 2a and 2b). Their simultaneous occurrence was expected as they appear similarly in the equations for G* and E* production.

Sensitivity indices for the current drop due to flash response, the rms error between simulation and experiment, and the dark current are given in Fig 3. Monte Carlo estimated that the five most influential parameters for goodness of fit were β_dark, a_min, , k_R, and k_E, with each implicated in at least 20% of the variance, and β_dark implicated in as much as 90% (Fig 3c). Moreover, β_dark was also the most influential parameter of the current drop and the dark current (Fig 3a, 3b and 3d). Generally, the total sensitivity indices for the dark current functional exhibited some similarity between its ranking of parameters and those for goodness of fit. This suggested that accurately reproducing the dark current was one of the more important criteria for reproducing experimental results.

Sensitivity indices for the time-to-peak of the current drop and its overshoot are plotted in Fig 4. Monte Carlo estimated that the four most influential parameters for time-to-peak were k_R, k_E, k_GE, and [PDE]_σ (Fig 4a and 4b). The first two were the rates of shut-off for R* and G*, and their importance was expected. The significance of the second two likely derived from these parameters postponing G* complexing with E*, which in the model was required before G* could decay into inactivated G. Monte Carlo also estimated that an overshoot, where during the flash response the total current temporarily exceeded the dark current value, could not be attributed to any one parameter (Fig 4c). Therefore, an overshoot occurred as a consequence of fundamentally nonlinear interactions between several model parameters. Total sensitivity indices for the overshoot indicated that parameters determining the dark current and time-to-peak were most implicated in its variance (Fig 4d). The uncertainty in these estimates was large, which was partly a consequence of an overshoot occurring infrequently.

Discussion

The large number of biochemical and geometric parameters together with the relative scarcity of experimental data for cone photoreceptors challenges visual transduction models of increasing biological sophistication. In response, many models in the literature, concerning rods as well as cones, have made simplifying spatial homogeneity assumptions for describing the signaling cascade or aggregated several distinct components of the cascade into simpler phenomenological terms [28, 36, 37, 39, 58–60]. Others have even taken a strictly phenomenological approach to reduce model uncertainty [61]. Faced with such complexity, we retained all the features of our spatio-temporal model [16], and employed a statistical approach to parameter sensitivity analysis that emphasized global behavior (statistical distributions), instead of restricting consideration to local behavior (pointwise derivatives). To better understand the relative importance of parameters across such uncertainty in the literature, we used Sobol indices to quantify which parameters were the most influential for the considered experimental cone flash response of [45]. The advantage of the Sobol method was that it allowed us to encode parameter uncertainty as probability distributions and then decompose the statistical variance of a functional into fractions attributed to any grouping of parameters. In the absence of a compelling reason to do otherwise, we described parameter uncertainty with uniform distributions over simple intervals. While the estimated Sobol indices did depend on the intervals assigned to the parameters, we considered a choice of interval as more robust than a choice of pointwise value. Moreover, this approach could track interactions between parameters while derivative-based sensitivity measures could not.

After reporting known parameter ranges from the literature, these ranges were validated by stochastic optimization, insofar as a parameter set was found that reproduced behaviors of experimental flash responses. The Sobol analysis (Tables 12 and 13 and Figs 2–4) showed the uncertainties in β_dark and a_min to be the most implicated in the error between model prediction and experiment. Among the eight parameters estimated by Monte Carlo to be the most influential for the L² error functional, four of these were also influential in determining the dark current as measured by their indices: β_dark at 94%, a_min at 55%, at 30%, and α_max at 15%. The other four could be attributed to shut-off of R*, k_R at 29%, and E*, k_E at 20%, the activation rate of G* by R*, ν_RG at 7%, and the light-induced hydrolysis of cGMP by E*, k_cat/K_m at 7%. No single parameter by itself could account for more than an estimated 6% of the model’s error with experiment (Table 12). Consequently, the flash response inextricably depended on nonlinear interactions between model parameters. This conclusion appeared true for the overshoot in the flash response as well. While the uncertainties associated to J_over’s indices were larger than other functionals, which may have been a consequence of an overshoot not occurring in all Monte Carlo trials, the eight parameters estimated as most significant could be loosely characterized as parameters influencing the dark current and parameters influencing the time-to-peak of the flash response. In particular, the value of [PDE]_σ was shown to be important for time-to-peak but not important for the peak amount of E* produced. In the latter case, this was expected to be a consequence of the rapid shut-off of R* and G* before the available population of E could be exhausted. In the former case, its significance was expected to be a consequence of the model requiring that G* complex with E before it could decay to its inactivated state. In practice, the implications would be that [PDE]_σ may be most important when it is used to estimate the dark turnover rate of cGMP, β_dark. Once that turnover rate is fixed, the significance of [PDE]_σ’s value may be much less. Similar conclusions may hold for the geometric parameters ν and ϵ₀, when these are used to estimate a unit conversion between the surface hydrolysis rate of cGMP at the bounding membrane disc and an equivalent volumic rate of cGMP turnover in the thin interdiscal space.

It was not surprising that the geometric parameters R_b, R_t, H, ν, and ϵ₀ exhibited relatively little influence on the photoresponse given that the Sobol analysis modeled β_dark as a parameter independent of ν and ϵ₀, that the flashes modeled in this work were of uniform intensity throughout the photoreceptor, and that the prescribed uncertainty associated with these parameters was much more narrow than their biochemical counterparts. In [16] it was observed in silico that mouse rod biochemistry simply expressed in a fish cone morphology could exhibit a single photon response with an overshoot. In the present work, these parameters were uninvolved in the overshoot for two reasons. First, the ranges of R_b and R_t did not encompass the size of the striped bass cone photoreceptor used in [16]. Second, the single photon response was the most spatially localized case while flashes of uniform intensity were the most homogeneous. We expected that spatial localization created greater opportunity for transduction machinery to become asynchronous as one component may have outpaced another. On the other hand, the parameter ω₀ exhibited a slightly stronger influence on overshoot as indicated by Table 13. This was expected to be a consequence of the sliver angle determining the width across which ion channels were localized on the lateral side.

Our parameter set of Table 10 was broadly consistent with that of other recent modeling efforts. For example [59] found that the dark hydrolysis of cGMP is ∼4x greater in cone than rods and an even higher ratio was given by [36]. Our ratio between β_dark selected in [3] for mouse rod and the value obtained here was ∼ 3x. [59] and [36] reported ∼8–17x faster PDE decay rates in cones compared to rods. In the present work, this rate was ∼13.7x faster. [59] also inferred that the rate of transducin activation by rhodopsin, ν_RG, is likely much smaller in cones than rods and an even higher ratio was given by [36]. Although, the value presented in Table 10 was comparable to values reported for mouse rod [3], the obtained rate for activation of PDE* by transducin (ν_GE ≈ 38 s⁻¹) was much smaller than that derived for mouse rod (ν_GE ≈ 500–1000 s⁻¹) [8, 62]. Since the current response depended on the E* population, a low value for ν_GE could compensate for ν_RG. A somewhat smaller reduction in the activation of PDE* by R* was estimated in red-sensitive cones, but interestingly, not in green-sensitive cones, by [36]. Here and with regard to other parameters, the presently considered model accounted for additional features of the cascade, so some differences compared to other models were to be expected.

The identification of parameters given in Table 10 cannot substitute for future experimental findings and must be refined as more becomes known. A natural future step would be a Bayesian inference of parameter values, and it is expected that such analysis could improve upon the model’s presented experimental predictions. At present, the GSA variance estimates may be strongly influenced by the large uncertainty ranges for some of the parameters, such as β_dark (Tables 8 and 9). This range for β_dark is between 1–2 orders of magnitude larger than that estimated by several other authors. For example, [36] found the dark cGMP turnover rate to be about 10 times faster in cones than rods. As future studies reduce these uncertainties, this approach would be better able to reveal how biological variation in parameters affects the cone responses, e.g., across different types of cones, during light adaptation and in diseased states.

Conclusion

We end by suggesting, on the basis of Fig 3, the following prioritization of parameters for future measurement and refinement. The turnover rate of cGMP in the dark, β_dark, has highest priority. The saturated exchanger current, , has second highest priority. After these, the rates of R* and E* shutoff, k_R and k_E, have priority. While it would certainly be very valuable to quantify a_min directly, this parameter may not be immediately accessible to experiment but instead would have to be estimated from other measurements. This, for example, would be possible in conjunction with Bayesian inference.

Supporting information

S1 Appendix. Supplementary materials.

https://doi.org/10.1371/journal.pone.0258721.s001

(PDF)

References

1. Lamb TD. Why rods and cones? Eye (Lond). 2016;30(2):179–85. pmid:26563661
- View Article
- PubMed/NCBI
- Google Scholar
2. Daemen FJ. Vertebrate rod outer segment membranes. Biochim Biophys Acta. 1973;300(3):255–88. pmid:4587475
- View Article
- PubMed/NCBI
- Google Scholar
3. Shen L, Caruso G, Bisegna P, Andreucci D, Gurevich VV, Hamm HE, et al. Dynamics of mouse rod phototransduction and its sensitivity to variation of key parameters. IET Sys Biol. 2010;4:12–32. pmid:20001089
- View Article
- PubMed/NCBI
- Google Scholar
4. Pugh EN Jr, Lamb TD. Amplification and kinetics of the activation steps in phototransduction. Biochim Biophys Acta. 1993;1141:111–149. pmid:8382952
- View Article
- PubMed/NCBI
- Google Scholar
5. Hamer RD, Nicholas SC, Tranchina D, Liebman PA, Lamb TD. Multiple steps of phosphorylation of activated rhodopsin can account for the reproducibility of vertebrate rod single-photon responses. J Gen Physiol. 2003;122:419–444. pmid:12975449
- View Article
- PubMed/NCBI
- Google Scholar
6. Andreucci D, Bisegna P, Caruso G, Hamm HE, DiBenedetto E. Mathematical model of the spatio-temporal dynamics of second messengers in visual transduction. Biophys J. 2003;85:1358–1376. pmid:12944255
- View Article
- PubMed/NCBI
- Google Scholar
7. Reingruber J, Holcman D. The dynamics of phosphodiesterase activation in rods and cones. Biophys J. 2008;94:1954–1970. pmid:18065454
- View Article
- PubMed/NCBI
- Google Scholar
8. Caruso G, Bisegna P, Andreucci D, Lenoci L, Gurevich VV, Hamm HE, et al. Identification of key factors that reduce the variability of the single photon response. Proc Natl Acad Sci USA. 2011;108(19):7804–7807. pmid:21518901
- View Article
- PubMed/NCBI
- Google Scholar
9. Gross OP, Pugh EN Jr, Burns ME. Spatiotemporal cGMP dynamics in living mouse rods. Biophys J. 2012;102(8):1775–1784. pmid:22768933
- View Article
- PubMed/NCBI
- Google Scholar
10. Invergo BM, Dell’Orco D, Montanucci L, Koch KW, Bertranpetit J. A comprehensive model of the phototransduction cascade in mouse rod cells. Mol Biosyst. 2014;10(6):1481–1489. pmid:24675755
- View Article
- PubMed/NCBI
- Google Scholar
11. Li W. Ground squirrel—A cool model for a bright vision. Semin Cell Dev Biol. 2020;106:127–134. pmid:32593518
- View Article
- PubMed/NCBI
- Google Scholar
12. Sajdak BS, Salmon AE, Cava JA, Allen KP, Freling S, Ramamirtham R, et al. Noninvasive imaging of the tree shrew eye: Wavefront analysis and retinal imaging with correlative histology. Exp Eye Res. 2019;185:107683. pmid:31158381
- View Article
- PubMed/NCBI
- Google Scholar
13. Fukagawa T, Takafuji K, Tachibanaki S, Kawamura S. Purification of cone outer segment for proteomic analysis on its membrane proteins in carp retina. PLoS ONE. 2017;12(3):e0173908. pmid:28291804
- View Article
- PubMed/NCBI
- Google Scholar
14. Tachibanaki S, Tsushima S, Kawamura S. Low amplification and fast visual pigment phosphorylation as mechanisms characterizing cone photoresponses. Proc Natl Acad Sci USA. 2001;98(24):14044–14049. pmid:11707584
- View Article
- PubMed/NCBI
- Google Scholar
15. Ingram NT, Sampath AP, Fain GL. Why are rods more sensitive than cones? J Physiol (Lond). 2016;594(19):5415–26. pmid:27218707
- View Article
- PubMed/NCBI
- Google Scholar
16. Klaus C, Caruso G, Gurevich VV, DiBenedetto E. Multi-scale, numerical modeling of spatio-temporal signaling in cone phototransduction. PLoS ONE. 2019;14(7):e0219848. pmid:31344066
- View Article
- PubMed/NCBI
- Google Scholar
17. Sobol IM. Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates. Math Comput Simulation. 2001;55(1-3):271–280.
- View Article
- Google Scholar
18. Saltelli A, Ratto M, Andres T, Campolongo F, Cariboni J, Gatelli D, et al. Global sensitivity analysis. The primer. John Wiley & Sons, Ltd., Chichester; 2008.
19. Caruso G, Khanal H, Alexiades V, Rieke F, Hamm HE, DiBenedetto E. Mathematical and computational modeling of spatio-temporal signaling in rod phototransduction. IEE Proc Syst Biol. 2005;152:119–137.
- View Article
- Google Scholar
20. Caruso G, Bisegna P, Shen L, Andreucci D, Hamm HE, DiBenedetto E. Modeling the role of incisures in vertebrate phototransduction. Biophys J. 2006;91:1192–1212. pmid:16714347
- View Article
- PubMed/NCBI
- Google Scholar
21. Caruso G, Bisegna P, Lenoci L, Andreucci D, Gurevich VV, Hamm HE, et al. Kinetics of rhodopsin deactivation and its role in regulating recovery and reproducibility in wild type and transgenic mouse photoresponse. PLoS Comput Biol. 2010;6(12):1–15. pmid:21200415
- View Article
- PubMed/NCBI
- Google Scholar
22. Robert CP, Casella G. Monte Carlo statistical methods. 2nd ed. Springer Texts in Statistics. Springer-Verlag, New York; 2004.
23. Brooks S, Gelman A, Jones GL, Meng XL, editors. Handbook of Markov chain Monte Carlo. Chapman & Hall/CRC Handbooks of Modern Statistical Methods. CRC Press, Boca Raton, FL; 2011.
24. Varadhan SRS. Probability theory. vol. 7 of Courant Lecture Notes in Mathematics. New York University, Courant Institute of Mathematical Sciences, New York; American Mathematical Society, Providence, RI; 2001.
25. Saltelli A. Making best use of model evaluations to compute sensitivity indices. Comput Phys Commun. 2002;145(2):280–297.
- View Article
- Google Scholar
26. Saltelli A, Annoni P, Azzini I, Campolongo F, Ratto M, Tarantola S. Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index. Comput Phys Commun. 2010;181(2):259–270.
- View Article
- Google Scholar
27. Carter-Dawson LD, Lavail MM. Rods and cones in the mouse retina. J Comp Neurol. 1979;188:245–262. pmid:500858
- View Article
- PubMed/NCBI
- Google Scholar
28. Holcman D, Korenbrot JI. Longitudinal diffusion in retinal rod and cone outer segment cytoplasm: the consequence of cell structure. Biophys J. 2004;86:2566–2582. pmid:15041693
- View Article
- PubMed/NCBI
- Google Scholar
29. Fetter RD, Corless JM. Morphological components associated with frog cone outer segment disc margins. Invest Ophthalmol Vis Sci. 1987;28(4):646–657. pmid:3493999
- View Article
- PubMed/NCBI
- Google Scholar
30. Kawamura S, Tachibanaki S. Explaining the functional differences of rods versus cones. Wiley Interdiscip Rev: Membr Transp Signal. 2012;1(5):675–683.
- View Article
- Google Scholar
31. Tachibanaki S, Yonetsu S, Fukaya S, Koshitani Y, Kawamura S. Low activation and fast inactivation of transducin in carp cones. J Biol Chem. 2012;287(49):41186–41194. pmid:23045532
- View Article
- PubMed/NCBI
- Google Scholar
32. Majumder A, Pahlberg J, Muradov H, Boyd KK, Sampath AP, Artemyev NO. Exchange of cone for rod phosphodiesterase 6 catalytic subunits in rod photoreceptors mimics in part features of light adaptation. J Neurosci. 2015;35(24):9225–9235. pmid:26085644
- View Article
- PubMed/NCBI
- Google Scholar
33. Holcman D, Korenbrot JI. The limit of photoreceptor sensitivity: molecular mechanisms of dark noise in retinal cones. J Gen Physiol. 2005;125:641–660. pmid:15928405
- View Article
- PubMed/NCBI
- Google Scholar
34. Leskov IB, Klenchin VA, Handy JW, Whitlock GG, Govardovskii VI, Bownds MD, et al. The gain of rod phototransduction: reconciliation of biochemical and electrophysiological measurements. Neuron. 2000;27:525–537. pmid:11055435
- View Article
- PubMed/NCBI
- Google Scholar
35. Koshitani Y, Tachibanaki S, Kawamura S. Quantitative aspects of cGMP phosphodiesterase activation in carp rods and cones. J Biol Chem. 2014;289(5):2651–2657. pmid:24344136
- View Article
- PubMed/NCBI
- Google Scholar
36. Astakhova L, Firsov M, Govardovskii V. Activation and quenching of the phototransduction cascade in retinal cones as inferred from electrophysiology and mathematical modeling. Mol Vis. 2015;21:244–263. pmid:25866462
- View Article
- PubMed/NCBI
- Google Scholar
37. Korenbrot JI. Speed, sensitivity, and stability of the light response in rod and cone photoreceptors: Facts and models. Prog Retin Eye Res. 2012;31:442–466. pmid:22658984
- View Article
- PubMed/NCBI
- Google Scholar
38. Detwiler PB, Ramanathan S, Sengupta A, Shraiman BI. Engineering aspects of enzymatic signal transduction: photoreceptors in the retina. Biophys J. 2000;79(6):2801–2817. pmid:11106590
- View Article
- PubMed/NCBI
- Google Scholar
39. Lamb TD, Pugh EN Jr. A quantitative account of the activation steps involved in phototransduction in amphibian photoreceptors. J Physiol (Lond). 1992;449:719–758. pmid:1326052
- View Article
- PubMed/NCBI
- Google Scholar
40. Vinberg F, Peshenko IV, Chen J, Dizhoor AM, Kefalov VJ. Guanylate cyclase-activating protein 2 contributes to phototransduction and light adaptation in mouse cone photoreceptors. J Biol Chem. 2018;293(19):7457–7465. pmid:29549122
- View Article
- PubMed/NCBI
- Google Scholar
41. Peshenko IV, Olshevskaya EV, Savchenko AB, Karan S, Palczewski K, Baehr W, et al. Enzymatic properties and regulation of the native isozymes of retinal membrane guanylyl cyclase (RetGC) from mouse photoreceptors. Biochemistry. 2011;50(25):5590–5600. pmid:21598940
- View Article
- PubMed/NCBI
- Google Scholar
42. Korenbrot JI. Speed, adaptation, and stability of the response to light in cone photoreceptors: The functional role of Ca-dependent modulation of ligand sensitivity in cGMP-gated ion channels. J Gen Physiol. 2012;139(1):31–56. pmid:22200947
- View Article
- PubMed/NCBI
- Google Scholar
43. Rebrik TI, Kotelnikova EA, Korenbrot JI. Time course and Ca2+ dependence of sensitivity modulation in cyclic GMP-gated currents of intact cone photoreceptors. J Gen Physiol. 2000;116:521–524. pmid:11004202
- View Article
- PubMed/NCBI
- Google Scholar
44. Ohyama T, Hackos DH, Frings S, Hagen V, Kaupp UB, Korenbrot JI. Fraction of the dark current carried by Ca(2+) through cGMP-gated ion channels of intact rod and cone photoreceptors. J Gen Physiol. 2000;116:735–754. pmid:11099344
- View Article
- PubMed/NCBI
- Google Scholar
45. Ingram NT, Sampath AP, Fain GL. Voltage-clamp recordings of light responses from wild-type and mutant mouse cone photoreceptors. J Gen Physiol. 2019;151(11):1287–1299. pmid:31562185
- View Article
- PubMed/NCBI
- Google Scholar
46. Yau KW. Phototransduction mechanism in retinal rods and cones. The Friedenwald Lecture. Invest Opthalmol Visual Sci. 1994;35:9–32. pmid:7507907
- View Article
- PubMed/NCBI
- Google Scholar
47. Agbaga MP, Merriman DK, Brush RS, Lydic TA, Conley SM, Naash MI, et al. Differential composition of DHA and very-long-chain PUFAs in rod and cone photoreceptors. J Lipid Res. 2018;59(9):1586–1596. pmid:29986998
- View Article
- PubMed/NCBI
- Google Scholar
48. Gupta BD, Williams TP. Lateral diffusion of visual pigments in toad (Bufo marinus) rods and in catfish (Ictalurus punctatus) cones. J Physiol (Lond). 1990;430:483–496. pmid:2128335
- View Article
- PubMed/NCBI
- Google Scholar
49. Klaus C, Caruso G. klauscj68/Homogenized-Cone-Sensitivity-Analysis: Homogenized Cone Sensitivity Analysis; 2021. https://doi.org/10.5281/zenodo.4731156
50. Klaus C. Data for: Phototransduction in retinal cones: Analysis of parameter importance. Dryad, Dataset; 2021. https://doi.org/10.5061/dryad.6djh9w11c
51. Saltelli A. What is sensitivity analysis? In: Saltelli A, Chan K, Scott EM, editors. Sensitivity analysis. Wiley Ser Probab. Stat. John Wiley & Sons, Ltd., Chichester; 2000. p. 3–13.
52. Ohio Supercomputer Center; 1987. http://osc.edu/ark:/19495/f5s1ph73.
53. Kulesa A, Krzywinski M, Blainey P, Altman N. Sampling distributions and the bootstrap. Nature Methods. 2015;12(6):477–478. pmid:26221652
- View Article
- PubMed/NCBI
- Google Scholar
54. Calvert PD, Krasnoperova NV, Lyubarsky AL, Isayama T, Nicoló M, Kosaras B, et al. Phototransduction in transgenic mice after targeted deletion of the rod transducin alpha-subunit. Proc Natl Acad Sci USA. 2000;97(25):13913–13918. pmid:11095744
- View Article
- PubMed/NCBI
- Google Scholar
55. Schneeweis DM, Schnapf JL. Photovoltage of rods and cones in the macaque retina. Science. 1995;268(5213):1053–1056. pmid:7754386
- View Article
- PubMed/NCBI
- Google Scholar
56. Box GEP, Tiao GC. Bayesian inference in statistical analysis. Addison-Wesley Publishing Co., Reading, Mass.-London-Don Mills, Ont.; 1973.
57. Gelman A, Carlin JB, Stern HS, Dunson DB, Vehtari A, Rubin DB. Bayesian data analysis. 3rd ed. Texts in Statistical Science Series. CRC Press, Boca Raton, FL; 2014.
58. Chen J, Woodruff ML, Wang T, Concepcion FA, Tranchina D, Fain GL. Channel modulation and the mechanism of light adaptation in mouse rods. J Neurosci. 2010;30(48):16232–16240. pmid:21123569
- View Article
- PubMed/NCBI
- Google Scholar
59. Reingruber J, Ingram NT, Griffis KG, Fain GL. A kinetic analysis of mouse rod and cone photoreceptor responses. J Physiol (Lond). 2020;598(17):3747–3763. pmid:32557629
- View Article
- PubMed/NCBI
- Google Scholar
60. Rotov AY, Astakhova LA, Firsov ML, Govardovskii VI. Origins of the phototransduction delay as inferred from stochastic and deterministic simulation of the amplification cascade. Mol Vis. 2017;23:416–430. pmid:28744093
- View Article
- PubMed/NCBI
- Google Scholar
61. Clark DA, Benichou R, Meister M, Azeredo da Silveira R. Dynamical adaptation in photoreceptors. PLoS Comput Biol. 2013;9(11):e1003289. pmid:24244119
- View Article
- PubMed/NCBI
- Google Scholar
62. Bisegna P, Caruso G, Andreucci D, Shen L, Gurevich VV, Hamm HE, et al. Diffusion of the second messengers in the cytoplasm acts as a variability suppressor of the single photon response in vertebrate photoresponse. Biophys J. 2008;94:3363–3383. pmid:18400950
- View Article
- PubMed/NCBI
- Google Scholar

[ref1] 1. Lamb TD. Why rods and cones? Eye (Lond). 2016;30(2):179–85. pmid:26563661
View Article
PubMed/NCBI
Google Scholar

[2] View Article

[3] PubMed/NCBI

[4] Google Scholar

[ref2] 2. Daemen FJ. Vertebrate rod outer segment membranes. Biochim Biophys Acta. 1973;300(3):255–88. pmid:4587475
View Article
PubMed/NCBI
Google Scholar

[6] View Article

[7] PubMed/NCBI

[8] Google Scholar

[ref3] 3. Shen L, Caruso G, Bisegna P, Andreucci D, Gurevich VV, Hamm HE, et al. Dynamics of mouse rod phototransduction and its sensitivity to variation of key parameters. IET Sys Biol. 2010;4:12–32. pmid:20001089
View Article
PubMed/NCBI
Google Scholar

[10] View Article

[11] PubMed/NCBI

[12] Google Scholar

[ref4] 4. Pugh EN Jr, Lamb TD. Amplification and kinetics of the activation steps in phototransduction. Biochim Biophys Acta. 1993;1141:111–149. pmid:8382952
View Article
PubMed/NCBI
Google Scholar

[14] View Article

[15] PubMed/NCBI

[16] Google Scholar

[ref5] 5. Hamer RD, Nicholas SC, Tranchina D, Liebman PA, Lamb TD. Multiple steps of phosphorylation of activated rhodopsin can account for the reproducibility of vertebrate rod single-photon responses. J Gen Physiol. 2003;122:419–444. pmid:12975449
View Article
PubMed/NCBI
Google Scholar

[18] View Article

[19] PubMed/NCBI

[20] Google Scholar

[ref6] 6. Andreucci D, Bisegna P, Caruso G, Hamm HE, DiBenedetto E. Mathematical model of the spatio-temporal dynamics of second messengers in visual transduction. Biophys J. 2003;85:1358–1376. pmid:12944255
View Article
PubMed/NCBI
Google Scholar

[22] View Article

[23] PubMed/NCBI

[24] Google Scholar

[ref7] 7. Reingruber J, Holcman D. The dynamics of phosphodiesterase activation in rods and cones. Biophys J. 2008;94:1954–1970. pmid:18065454
View Article
PubMed/NCBI
Google Scholar

[26] View Article

[27] PubMed/NCBI

[28] Google Scholar

[ref8] 8. Caruso G, Bisegna P, Andreucci D, Lenoci L, Gurevich VV, Hamm HE, et al. Identification of key factors that reduce the variability of the single photon response. Proc Natl Acad Sci USA. 2011;108(19):7804–7807. pmid:21518901
View Article
PubMed/NCBI
Google Scholar

[30] View Article

[31] PubMed/NCBI

[32] Google Scholar

[ref9] 9. Gross OP, Pugh EN Jr, Burns ME. Spatiotemporal cGMP dynamics in living mouse rods. Biophys J. 2012;102(8):1775–1784. pmid:22768933
View Article
PubMed/NCBI
Google Scholar

[34] View Article

[35] PubMed/NCBI

[36] Google Scholar

[ref10] 10. Invergo BM, Dell’Orco D, Montanucci L, Koch KW, Bertranpetit J. A comprehensive model of the phototransduction cascade in mouse rod cells. Mol Biosyst. 2014;10(6):1481–1489. pmid:24675755
View Article
PubMed/NCBI
Google Scholar

[38] View Article

[39] PubMed/NCBI

[40] Google Scholar

[ref11] 11. Li W. Ground squirrel—A cool model for a bright vision. Semin Cell Dev Biol. 2020;106:127–134. pmid:32593518
View Article
PubMed/NCBI
Google Scholar

[42] View Article

[43] PubMed/NCBI

[44] Google Scholar

[ref12] 12. Sajdak BS, Salmon AE, Cava JA, Allen KP, Freling S, Ramamirtham R, et al. Noninvasive imaging of the tree shrew eye: Wavefront analysis and retinal imaging with correlative histology. Exp Eye Res. 2019;185:107683. pmid:31158381
View Article
PubMed/NCBI
Google Scholar

[46] View Article

[47] PubMed/NCBI

[48] Google Scholar

[ref13] 13. Fukagawa T, Takafuji K, Tachibanaki S, Kawamura S. Purification of cone outer segment for proteomic analysis on its membrane proteins in carp retina. PLoS ONE. 2017;12(3):e0173908. pmid:28291804
View Article
PubMed/NCBI
Google Scholar

[50] View Article

[51] PubMed/NCBI

[52] Google Scholar

[ref14] 14. Tachibanaki S, Tsushima S, Kawamura S. Low amplification and fast visual pigment phosphorylation as mechanisms characterizing cone photoresponses. Proc Natl Acad Sci USA. 2001;98(24):14044–14049. pmid:11707584
View Article
PubMed/NCBI
Google Scholar

[54] View Article

[55] PubMed/NCBI

[56] Google Scholar

[ref15] 15. Ingram NT, Sampath AP, Fain GL. Why are rods more sensitive than cones? J Physiol (Lond). 2016;594(19):5415–26. pmid:27218707
View Article
PubMed/NCBI
Google Scholar

[58] View Article

[59] PubMed/NCBI

[60] Google Scholar

[ref16] 16. Klaus C, Caruso G, Gurevich VV, DiBenedetto E. Multi-scale, numerical modeling of spatio-temporal signaling in cone phototransduction. PLoS ONE. 2019;14(7):e0219848. pmid:31344066
View Article
PubMed/NCBI
Google Scholar

[62] View Article

[63] PubMed/NCBI

[64] Google Scholar

[ref17] 17. Sobol IM. Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates. Math Comput Simulation. 2001;55(1-3):271–280.
View Article
Google Scholar

[66] View Article

[67] Google Scholar

[ref18] 18. Saltelli A, Ratto M, Andres T, Campolongo F, Cariboni J, Gatelli D, et al. Global sensitivity analysis. The primer. John Wiley & Sons, Ltd., Chichester; 2008.

[ref19] 19. Caruso G, Khanal H, Alexiades V, Rieke F, Hamm HE, DiBenedetto E. Mathematical and computational modeling of spatio-temporal signaling in rod phototransduction. IEE Proc Syst Biol. 2005;152:119–137.
View Article
Google Scholar

[70] View Article

[71] Google Scholar

[ref20] 20. Caruso G, Bisegna P, Shen L, Andreucci D, Hamm HE, DiBenedetto E. Modeling the role of incisures in vertebrate phototransduction. Biophys J. 2006;91:1192–1212. pmid:16714347
View Article
PubMed/NCBI
Google Scholar

[73] View Article

[74] PubMed/NCBI

[75] Google Scholar

[ref21] 21. Caruso G, Bisegna P, Lenoci L, Andreucci D, Gurevich VV, Hamm HE, et al. Kinetics of rhodopsin deactivation and its role in regulating recovery and reproducibility in wild type and transgenic mouse photoresponse. PLoS Comput Biol. 2010;6(12):1–15. pmid:21200415
View Article
PubMed/NCBI
Google Scholar

[77] View Article

[78] PubMed/NCBI

[79] Google Scholar

[ref22] 22. Robert CP, Casella G. Monte Carlo statistical methods. 2nd ed. Springer Texts in Statistics. Springer-Verlag, New York; 2004.

[ref23] 23. Brooks S, Gelman A, Jones GL, Meng XL, editors. Handbook of Markov chain Monte Carlo. Chapman & Hall/CRC Handbooks of Modern Statistical Methods. CRC Press, Boca Raton, FL; 2011.

[ref24] 24. Varadhan SRS. Probability theory. vol. 7 of Courant Lecture Notes in Mathematics. New York University, Courant Institute of Mathematical Sciences, New York; American Mathematical Society, Providence, RI; 2001.

[ref25] 25. Saltelli A. Making best use of model evaluations to compute sensitivity indices. Comput Phys Commun. 2002;145(2):280–297.
View Article
Google Scholar

[84] View Article

[85] Google Scholar

[ref26] 26. Saltelli A, Annoni P, Azzini I, Campolongo F, Ratto M, Tarantola S. Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index. Comput Phys Commun. 2010;181(2):259–270.
View Article
Google Scholar

[87] View Article

[88] Google Scholar

[ref27] 27. Carter-Dawson LD, Lavail MM. Rods and cones in the mouse retina. J Comp Neurol. 1979;188:245–262. pmid:500858
View Article
PubMed/NCBI
Google Scholar

[90] View Article

[91] PubMed/NCBI

[92] Google Scholar

[ref28] 28. Holcman D, Korenbrot JI. Longitudinal diffusion in retinal rod and cone outer segment cytoplasm: the consequence of cell structure. Biophys J. 2004;86:2566–2582. pmid:15041693
View Article
PubMed/NCBI
Google Scholar

[94] View Article

[95] PubMed/NCBI

[96] Google Scholar

[ref29] 29. Fetter RD, Corless JM. Morphological components associated with frog cone outer segment disc margins. Invest Ophthalmol Vis Sci. 1987;28(4):646–657. pmid:3493999
View Article
PubMed/NCBI
Google Scholar

[98] View Article

[99] PubMed/NCBI

[100] Google Scholar

[ref30] 30. Kawamura S, Tachibanaki S. Explaining the functional differences of rods versus cones. Wiley Interdiscip Rev: Membr Transp Signal. 2012;1(5):675–683.
View Article
Google Scholar

[102] View Article

[103] Google Scholar

[ref31] 31. Tachibanaki S, Yonetsu S, Fukaya S, Koshitani Y, Kawamura S. Low activation and fast inactivation of transducin in carp cones. J Biol Chem. 2012;287(49):41186–41194. pmid:23045532
View Article
PubMed/NCBI
Google Scholar

[105] View Article

[106] PubMed/NCBI

[107] Google Scholar

[ref32] 32. Majumder A, Pahlberg J, Muradov H, Boyd KK, Sampath AP, Artemyev NO. Exchange of cone for rod phosphodiesterase 6 catalytic subunits in rod photoreceptors mimics in part features of light adaptation. J Neurosci. 2015;35(24):9225–9235. pmid:26085644
View Article
PubMed/NCBI
Google Scholar

[109] View Article

[110] PubMed/NCBI

[111] Google Scholar

[ref33] 33. Holcman D, Korenbrot JI. The limit of photoreceptor sensitivity: molecular mechanisms of dark noise in retinal cones. J Gen Physiol. 2005;125:641–660. pmid:15928405
View Article
PubMed/NCBI
Google Scholar

[113] View Article

[114] PubMed/NCBI

[115] Google Scholar

[ref34] 34. Leskov IB, Klenchin VA, Handy JW, Whitlock GG, Govardovskii VI, Bownds MD, et al. The gain of rod phototransduction: reconciliation of biochemical and electrophysiological measurements. Neuron. 2000;27:525–537. pmid:11055435
View Article
PubMed/NCBI
Google Scholar

[117] View Article

[118] PubMed/NCBI

[119] Google Scholar

[ref35] 35. Koshitani Y, Tachibanaki S, Kawamura S. Quantitative aspects of cGMP phosphodiesterase activation in carp rods and cones. J Biol Chem. 2014;289(5):2651–2657. pmid:24344136
View Article
PubMed/NCBI
Google Scholar

[121] View Article

[122] PubMed/NCBI

[123] Google Scholar

[ref36] 36. Astakhova L, Firsov M, Govardovskii V. Activation and quenching of the phototransduction cascade in retinal cones as inferred from electrophysiology and mathematical modeling. Mol Vis. 2015;21:244–263. pmid:25866462
View Article
PubMed/NCBI
Google Scholar

[125] View Article

[126] PubMed/NCBI

[127] Google Scholar

[ref37] 37. Korenbrot JI. Speed, sensitivity, and stability of the light response in rod and cone photoreceptors: Facts and models. Prog Retin Eye Res. 2012;31:442–466. pmid:22658984
View Article
PubMed/NCBI
Google Scholar

[129] View Article

[130] PubMed/NCBI

[131] Google Scholar

[ref38] 38. Detwiler PB, Ramanathan S, Sengupta A, Shraiman BI. Engineering aspects of enzymatic signal transduction: photoreceptors in the retina. Biophys J. 2000;79(6):2801–2817. pmid:11106590
View Article
PubMed/NCBI
Google Scholar

[133] View Article

[134] PubMed/NCBI

[135] Google Scholar

[ref39] 39. Lamb TD, Pugh EN Jr. A quantitative account of the activation steps involved in phototransduction in amphibian photoreceptors. J Physiol (Lond). 1992;449:719–758. pmid:1326052
View Article
PubMed/NCBI
Google Scholar

[137] View Article

[138] PubMed/NCBI

[139] Google Scholar

[ref40] 40. Vinberg F, Peshenko IV, Chen J, Dizhoor AM, Kefalov VJ. Guanylate cyclase-activating protein 2 contributes to phototransduction and light adaptation in mouse cone photoreceptors. J Biol Chem. 2018;293(19):7457–7465. pmid:29549122
View Article
PubMed/NCBI
Google Scholar

[141] View Article

[142] PubMed/NCBI

[143] Google Scholar

[ref41] 41. Peshenko IV, Olshevskaya EV, Savchenko AB, Karan S, Palczewski K, Baehr W, et al. Enzymatic properties and regulation of the native isozymes of retinal membrane guanylyl cyclase (RetGC) from mouse photoreceptors. Biochemistry. 2011;50(25):5590–5600. pmid:21598940
View Article
PubMed/NCBI
Google Scholar

[145] View Article

[146] PubMed/NCBI

[147] Google Scholar

[ref42] 42. Korenbrot JI. Speed, adaptation, and stability of the response to light in cone photoreceptors: The functional role of Ca-dependent modulation of ligand sensitivity in cGMP-gated ion channels. J Gen Physiol. 2012;139(1):31–56. pmid:22200947
View Article
PubMed/NCBI
Google Scholar

[149] View Article

[150] PubMed/NCBI

[151] Google Scholar

[ref43] 43. Rebrik TI, Kotelnikova EA, Korenbrot JI. Time course and Ca2+ dependence of sensitivity modulation in cyclic GMP-gated currents of intact cone photoreceptors. J Gen Physiol. 2000;116:521–524. pmid:11004202
View Article
PubMed/NCBI
Google Scholar

[153] View Article

[154] PubMed/NCBI

[155] Google Scholar

[ref44] 44. Ohyama T, Hackos DH, Frings S, Hagen V, Kaupp UB, Korenbrot JI. Fraction of the dark current carried by Ca(2+) through cGMP-gated ion channels of intact rod and cone photoreceptors. J Gen Physiol. 2000;116:735–754. pmid:11099344
View Article
PubMed/NCBI
Google Scholar

[157] View Article

[158] PubMed/NCBI

[159] Google Scholar

[ref45] 45. Ingram NT, Sampath AP, Fain GL. Voltage-clamp recordings of light responses from wild-type and mutant mouse cone photoreceptors. J Gen Physiol. 2019;151(11):1287–1299. pmid:31562185
View Article
PubMed/NCBI
Google Scholar

[161] View Article

[162] PubMed/NCBI

[163] Google Scholar

[ref46] 46. Yau KW. Phototransduction mechanism in retinal rods and cones. The Friedenwald Lecture. Invest Opthalmol Visual Sci. 1994;35:9–32. pmid:7507907
View Article
PubMed/NCBI
Google Scholar

[165] View Article

[166] PubMed/NCBI

[167] Google Scholar

[ref47] 47. Agbaga MP, Merriman DK, Brush RS, Lydic TA, Conley SM, Naash MI, et al. Differential composition of DHA and very-long-chain PUFAs in rod and cone photoreceptors. J Lipid Res. 2018;59(9):1586–1596. pmid:29986998
View Article
PubMed/NCBI
Google Scholar

[169] View Article

[170] PubMed/NCBI

[171] Google Scholar

[ref48] 48. Gupta BD, Williams TP. Lateral diffusion of visual pigments in toad (Bufo marinus) rods and in catfish (Ictalurus punctatus) cones. J Physiol (Lond). 1990;430:483–496. pmid:2128335
View Article
PubMed/NCBI
Google Scholar

[173] View Article

[174] PubMed/NCBI

[175] Google Scholar

[ref49] 49. Klaus C, Caruso G. klauscj68/Homogenized-Cone-Sensitivity-Analysis: Homogenized Cone Sensitivity Analysis; 2021. https://doi.org/10.5281/zenodo.4731156

[ref50] 50. Klaus C. Data for: Phototransduction in retinal cones: Analysis of parameter importance. Dryad, Dataset; 2021. https://doi.org/10.5061/dryad.6djh9w11c

[ref51] 51. Saltelli A. What is sensitivity analysis? In: Saltelli A, Chan K, Scott EM, editors. Sensitivity analysis. Wiley Ser Probab. Stat. John Wiley & Sons, Ltd., Chichester; 2000. p. 3–13.

[ref52] 52. Ohio Supercomputer Center; 1987. http://osc.edu/ark:/19495/f5s1ph73.

[ref53] 53. Kulesa A, Krzywinski M, Blainey P, Altman N. Sampling distributions and the bootstrap. Nature Methods. 2015;12(6):477–478. pmid:26221652
View Article
PubMed/NCBI
Google Scholar

[181] View Article

[182] PubMed/NCBI

[183] Google Scholar

[ref54] 54. Calvert PD, Krasnoperova NV, Lyubarsky AL, Isayama T, Nicoló M, Kosaras B, et al. Phototransduction in transgenic mice after targeted deletion of the rod transducin alpha-subunit. Proc Natl Acad Sci USA. 2000;97(25):13913–13918. pmid:11095744
View Article
PubMed/NCBI
Google Scholar

[185] View Article

[186] PubMed/NCBI

[187] Google Scholar

[ref55] 55. Schneeweis DM, Schnapf JL. Photovoltage of rods and cones in the macaque retina. Science. 1995;268(5213):1053–1056. pmid:7754386
View Article
PubMed/NCBI
Google Scholar

[189] View Article

[190] PubMed/NCBI

[191] Google Scholar

[ref56] 56. Box GEP, Tiao GC. Bayesian inference in statistical analysis. Addison-Wesley Publishing Co., Reading, Mass.-London-Don Mills, Ont.; 1973.

[ref57] 57. Gelman A, Carlin JB, Stern HS, Dunson DB, Vehtari A, Rubin DB. Bayesian data analysis. 3rd ed. Texts in Statistical Science Series. CRC Press, Boca Raton, FL; 2014.

[ref58] 58. Chen J, Woodruff ML, Wang T, Concepcion FA, Tranchina D, Fain GL. Channel modulation and the mechanism of light adaptation in mouse rods. J Neurosci. 2010;30(48):16232–16240. pmid:21123569
View Article
PubMed/NCBI
Google Scholar

[195] View Article

[196] PubMed/NCBI

[197] Google Scholar

[ref59] 59. Reingruber J, Ingram NT, Griffis KG, Fain GL. A kinetic analysis of mouse rod and cone photoreceptor responses. J Physiol (Lond). 2020;598(17):3747–3763. pmid:32557629
View Article
PubMed/NCBI
Google Scholar

[199] View Article

[200] PubMed/NCBI

[201] Google Scholar

[ref60] 60. Rotov AY, Astakhova LA, Firsov ML, Govardovskii VI. Origins of the phototransduction delay as inferred from stochastic and deterministic simulation of the amplification cascade. Mol Vis. 2017;23:416–430. pmid:28744093
View Article
PubMed/NCBI
Google Scholar

[203] View Article

[204] PubMed/NCBI

[205] Google Scholar

[ref61] 61. Clark DA, Benichou R, Meister M, Azeredo da Silveira R. Dynamical adaptation in photoreceptors. PLoS Comput Biol. 2013;9(11):e1003289. pmid:24244119
View Article
PubMed/NCBI
Google Scholar

[207] View Article

[208] PubMed/NCBI

[209] Google Scholar

[ref62] 62. Bisegna P, Caruso G, Andreucci D, Shen L, Gurevich VV, Hamm HE, et al. Diffusion of the second messengers in the cytoplasm acts as a variability suppressor of the single photon response in vertebrate photoresponse. Biophys J. 2008;94:3363–3383. pmid:18400950
View Article
PubMed/NCBI
Google Scholar

[211] View Article

[212] PubMed/NCBI

[213] Google Scholar

Figures

Abstract

Introduction

Materials and methods

Proposing first ranges for parameter values of the cone outer segment

Geometry.

G-protein/effector cascade.

Catalytic activity.

Guanylyl cyclase activity.

Ionic current parameters.

Diffusion coefficients.

Choice of ranges for Markov Chain Monte Carlo.

Choice of ranges for global sensitivity analysis

Numerical implementation by finite elements

Identifying a consistent set of parameters

Choice of functionals for sensitivity analysis

Local sensitivity analysis

Global sensitivity analysis

Renormalizing αmin to amin.

Results

Consistent parameter set

Local sensitivity of parameters

Global sensitivity of parameters

Findings.

Discussion

Conclusion

Supporting information

S1 Appendix. Supplementary materials.

References

Renormalizing α_min to a_min.