Caffeine-induced Ca2+ oscillations in type I horizontal cell of carp retina: A mathematical model

Oscillations in intracellular free Ca2+ concentration ([Ca2+]i) have been observed in a variety of cell types. In the present study, we constructed a mathematical model to simulate the caffeine-induced [Ca2+]i oscillations based on experimental data obtained from isolated type I horizontal cell of carp retina. The results of model analysis confirm the notion that the caffeine-induced [Ca2+]i oscillations involve a number of cytoplasmic and endoplasmic Ca2+ processes that interact with each other. Using this model, we evaluated the importance of store-operated channel (SOC) in caffeine-induced [Ca2+]i oscillations. The model suggests that store-operated Ca2+ entry (SOCE) is elicited upon depletion of the endoplasmic reticulum (ER). When the SOC conductance is set to 0, caffeine-induced [Ca2+]i oscillations are abolished, which agrees with the experimental observation that [Ca2+]i oscillations were abolished when SOC was blocked pharmacologically, verifying that SOC is necessary for sustained [Ca2+]i oscillations.


Introduction
[Ca 2C ] i oscillations have been observed in a variety of cell types. 1 The present study is focused on caffeine-induced [Ca 2C ] i oscillations in type I horizontal cell (H1 HC) of carp retina. Horizontal cells are interneurons in vertebrate retinas, which receive glutamate input from photoreceptors. In carp H1 HC, [Ca 2C ] i increase is caused primarily by Ca 2C release from intracellular Ca 2C store (basically the endoplasmic reticulum, ER) via ryanodine receptor (RyR) activation 2,3 and Ca 2C influx from the extracellular environment, while Ca 2C removal depends on the activities of the sarco/endoplasmic reticulum Ca 2C -ATPase (SERCA), [2][3][4] Na C /Ca 2C exchanger (NCX) 5 and plasma membrane Ca 2C -ATPase (PMCA). 5 Ca 2C imaging study with carp H1 HC revealed that [Ca 2C ] i oscillations could be induced when caffeine was applied at various concentrations. 6 The application of caffeine activates RyR, leading to ER Ca 2C release and ER depletion. Subsequent ER refilling, which should be required for the maintenance of [Ca 2C ] i oscillations, is dependent on Ca 2C entry via storeoperated channel (SOC). Hence RyR and SOC, together with the Ca 2C removal processes are the major components for caffeine-induced [Ca 2C ] i oscillations. However, how these components interact with each other to give rise to [Ca 2C ] i oscillations observed in H1 HC is unclear.
In the present study, a computational model based on biophysical properties of H1 HC was constructed to explore the underlying mechanism(s) of caffeine-induced [

Results
Simulating caffeine-induced [Ca 2C ] i oscillations in isolated H1 HC Figure 1A-C shows the experimental observations of caffeineinduced [Ca 2C ] i oscillations of 3 representative H1 HCs in response to 3, 6, 10 mM caffeine (with identical caffeine application time of 90 s), respectively. 6 Statistical results of the experimental data show that, when the concentration of caffeine was increased from 3 to 10 mM, the time interval between the first and the second Ca 2C transients (DT 1-2 ; Fig. 1D) was decreased from 46.08 § 3.41 s to 28.37 § 2.30 s (P < 0.05, ANOVA, post hoc SNK test; Fig. 1E). The amplitude decrement during the Ca 2C oscillations was also concentration dependent. To quantify such amplitude decrement, the amplitude of the second transient was normalized against that of the first one (A 2 /A 1 ; Fig. 1D). When caffeine concentration was increased from 3 to 10 mM, the normalized amplitude of the second transient (A 2 /A 1 ) was decreased from 0.82 § 0.05 to 0.58 § 0.03 (P < 0.05, ANOVA, post hoc SNK test; Fig. 1F). 6 Caffeine-induced [Ca 2C ] i oscillations can be attributed to a number of factors (Fig. 6). Using the model cell, [Ca 2C ] i oscillations induced by 3, 6, 10 mM caffeine were simulated ( Fig. 2A-C) with parameter values listed in Table 3 and 4. The positive currents reflect Ca 2C efflux from cytoplasm to external environment or Ca 2C entering the ER.
Caffeine induces ER Ca 2C release by increasing the open probability of RyR (P RyR ), and P RyR is increased when the concentration of caffeine is increased. 7,8 To simulate the caffeine effect on [Ca 2C ] i , the RyR activation constant (K d,Ca in Eq. 10), which reflects the RyR's sensitivity to Ca 2C , was adjusted.
Model output simulating [Ca 2C ] i oscillations induced by 3 mM caffeine is shown in Figure 2A.  Figure 2B and C show model output induced by 6 and 10 mM caffeine, respectively. During 6 mM and 10 mM   (Fig. 2). By comparing [Ca 2C ] ER curves in Figure 2A-C, it can be noticed that u th is decreased from 85.43 mM to 66.74 mM when caffeine concentration is increased from 3 mM to 10 mM ( Table 1).
Caffeine effect at different concentrations is represented by using different K d,Ca values, which results in RyR current modulation by influencing P RyR (Eq. 10).  returns to the basal level, and a high caffeine concentration corresponds to a large P trough value and a quick P RyR increase (Fig. 3D). With a large P RyR , the second [Ca 2C ] i transient can be induce quickly. Therefore, when the concentration of caffeine is increased, the DT 1-2 interval is decreased (

Effect of SOC inhibition on caffeineinduced [Ca 2C ] i oscillations
Experimental data from isolated H1 HCs showed that prolonged application of 3 mM caffeine (7 min) induced sustained [Ca 2C ] i oscillations in H1 HCs ( Fig. 4A), 6 and such [Ca 2C ] i oscillations were abolished when 2-APB (SOC inhibitor) was co-applied ( Fig. 4B), 6 suggesting that SOCE is necessary for ER refilling.
Using our model, the above experimental results (Fig. 4) are simulated. The effect of 3 mM caffeine is simulated with K d,Ca D 0.16 mM, the effect of 2-APB was simulated by setting V SOC (the maximum SOC flux, Eq. 13) to 0. As shown in Figure 5A, prolonged application of caffeine ( ER level is recovered and reaches a certain threshold u th . When the concentration of caffeine is increased, the value of u th is decreased (Fig. 2, Table 1). The concentration-dependent caffeine effect on u th is consistent with the experimental results reported by Kong et al. 9 In their study, by recording [Ca 2C ] ER of HEK-293 cells expressing RyR, it was found that the threshold for Ca 2C release (u th ) was decreased when the concentration of caffeine was increased.
In carp retinal horizontal cell, Ca 2C entering the cytoplasmic is extruded into the extracellular space by both PMCA and NCX. 5 From the simulation output (Fig. 2, insets), it can be observed that during caffeine-induced [Ca 2C ] i oscillations, Ca 2C is extruded mainly by PMCA when [Ca 2C ] i is at lower levels; while at higher [Ca 2C ] i levels, NCX-mediated Ca 2C current is larger than the PMCA Ca 2C current. These results confirm the previously reported basic characteristics of PMCA and NCX, i.e., PMCA has higher affinity for Ca 2C than NCX, while NCX has higher turnover rate, 10,11 suggesting that PMCA and NCX cooperate to regulate Ca 2C signals over a wide [Ca 2C ] i range in H1 HC.
In addition to regulating the dynamic changes in [Ca 2C ] i , PMCA has also been shown to participate in the regulation of resting Ca 2C homeostasis. It was reported that in rod-driven bipolar cell of mouse retina, PMCA was involved in the maintenance of the resting [Ca 2C ] i levels. 12,13 For our model cell, after   Table 2) until a subsequent stimulus is applied. The model output shows that NCX is inactive at rest (I REST NCX 0.0 mA/cm 2 , Table 2) and Ca 2C extrusion occurs through PMCA (I REST PMCA 0.4 mA/cm 2 , Table 2). The simulation results suggest that in the resting state, extrusion of cytoplasmic Ca 2C depends mainly on PMCA in H1 HC.
It can also be noticed that the SOCE complex is activated at rest (I REST SOC 0.4 mA/cm 2 , Table 2). On the plasma membrane, Ca 2C extruded from the cell via PMCA is counterbalanced by Ca 2C entering the cell through SOC. On the ER membrane, the passive leakage of Ca 2C from the ER (I REST leak ¡0.5 mA/cm 2 , Table 2) is compensated by SERCA (I REST SERCA 0.1 mA/cm 2 s, Table 2) and SOCE complex-mediated Ca 2C uptake. Ca 2C homeostasis of the model cell is maintained by the activities of PMCA, SOCE complex, ER passive leak, and SERCA, keeping resting [Ca 2C ] i and [Ca 2C ] ER levels at 28 nM and 94 mM, respectively ( Table 2). The involvement of the SOCE complex in resting Ca 2C homeostasis has been reported in cells such as rat primary cortical neurons, mouse mammary epithelial cells, and HEK293 cells. [14][15][16][17] In these cells, 2 types of STIM proteins, i.e., STIM1 and STIM2, are expressed on the ER membrane. STIM1 triggers SOCE in response to store-depletion, while STIM2 activates Orai1 upon small fluctuations in [Ca 2C ] ER , keeping resting [Ca 2C ] i and [Ca 2C ] ER stable. Therefore, our simulation results suggest that in addition to ER refilling upon stimuli-induced store depletion (Fig. 5)

Methods
Ca 2C imaging data Experimental data were obtained from isolated type I horizontal cell of carp retina using fluo-3-based Ca 2C imaging technique. All the materials and methods have been previously reported. 6

Model description
All simulations were carried out using NEURON software. 18 A cylindrical neuron was constructed, with the [Ca 2C ] i dynamics being modeled at the whole-cell level.
The model system is schematically presented in Figure 6, which includes Ca 2C processes on the plasma membrane (PM), the ER, and the store-operated Ca 2C entry (SOCE) complex. Biophysical properties of the neuron membrane are spatially uniform.
The PM processes include the L-VGCC, NCX, and PMCA. L-VGCC mediates Ca 2C entry during membrane depolarization, while NCX and PMCA extrude Ca 2C from the cytosol into the extracellular space.
The ER processes involve 3 parallel components: Ca 2C release from the ER into the cytosol mediated by RyR, passive Ca 2C leakage from the ER into the cytosol, and Ca 2C uptake from the cytosol into the ER through SERCA.
SOCE is a process in which plasma membrane Ca 2C channel (SOC) is activated upon Ca 2C store depletion. 19 The SOCE complex is composed of proteins including STIM, SOC, and SERCA. STIM which resides on the ER membrane senses the ER Ca 2C content, 20,21 while Orai channels [22][23][24] and/or TRP channels 25,26 form the SOC on the PM. Upon store depletion, STIM and SOC translocate and cluster at the ER-PM junctions, leading to the formation of the STIM-SOC complexes and SOC activation. 23,27,28 Activated SOC mediates Ca 2C influx, and SERCA in the ER-PM junctions uptakes these Ca 2C to refill the depleted ER. [29][30][31][32][33] Besides, voltage-gated K C and Na C channel conductance also contribute to the regulation of [Ca 2C ] i by affecting the neuronal membrane potential, therefore, these processes are also included in our model.
Hence, in our model, Ca 2C dynamics can be described as: where [Ca 2C ] i is the free intracellular Ca 2C concentration, [Ca 2C ] ER is the free Ca 2C concentration inside the ER lumen; f cyt and f ER are Ca 2C buffering coefficients for cytoplasm and ER, respectively; vol cyt and vol ER are the volumes of the cell (excluding ER) and the ER, respectively; F is the Faraday constant; I All Ca represents the total Ca 2C currents across the cell membrane; J RyR represents ER Ca 2C release via RyR; J leak represents the passive Ca 2C leakage from the ER into the cytosol; J SERCA represents the Ca 2C uptake from the cytosol into the ER through SERCA; J SOC represents Ca 2C entering the ER through the SOCE complex.
The total Ca 2C currents across the plasma membrane are given by: where I VGCC , I NCX , and I PMCA denote the Ca 2C currents mediated by L-VGCC, NCX, and PMCA respectively. Ca 2C current through L-VGCC (I VGCC ) can be described as: 4 where g Ca is the maximum conductance of L-VGCC; m Ca and h Ca are the activation and inactivation variables, respectively; V m denotes the membrane potential; E Ca stands for the reversal potential of Ca 2C current; K Ca is the half-inactivation parameter; n is the Hill coefficient; a and b are forward and backward rate coefficients, respectively; t Ca represents the inactivation time constant.
Equation 3 shows that the Ca 2C dynamics are regulated by V m , while Ca 2C currents together with Na C and K C currents across the plasma membrane all contribute to the V m changes. Hence voltage-gated Na C conductance and 3 types of K C conductance (i.e., delay time rectifying K C current, outward rectifying K C current, and anomalous rectifying K C current), which are expressed in carp retinal HC, are also included in our model. 34 Current carried by voltage-gated Na C channel (I Nav ) can be described as: 34 where g Na is the maximum conductance of the voltage-gated Na C channel; m Na and h Na are the activation and inactivation variables, respectively; E Na is the reversal potential of Na C current.
Current carried by delay time rectifying K C channel (I kd ) can be described as: 34 where g Kd is the maximum conductance of the delay time rectifying K C channel; m Kd represents the activation variable; E K is the reversal potential of K C current.
Current carried by outward rectifying K C channel (I Kt ) can be described as: 34 where g Kt is the maximum conductance of the outward rectifying K C channel; m Kt and h Kt are the activation and inactivation variables, respectively. Current carried by anomalous rectifying K C channel (I Ka ) can be described as: 34 where g Ka is the maximum conductance of the anomalous rectifying K C channel; m Ka represents the activation variable.
On the plasma membrane, Ca 2C extrusion is mediated by NCX and PMCA. Current carried by NCX (I All NCX ) can be described as: 4 where where J PMCA represents the amount of Ca 2C efflux; A PMCA is the maximal pumping rate; K PMCA is the dissociation constant. RyR, passive ER Ca 2C leakage, and SERCA are involved in Ca 2C exchange between the cytoplasm and the ER. RyR-mediated Ca 2C flux from the ER into the cytosol is described as: [35][36][37] where K leak is the rate constant. ER Ca 2C uptake mediated by SERCA can be described as: 2,37 where V SERCA is the maximum SERCA flux; K SERCA1 and K SERCA2 are the activation and inactivation constants, respectively.
Upon store depletion, SOCE complex can be formed in the ER-PM junctions. Due to the tight coupling between the SOC and SERCA, the majority of Ca 2C entering the cytoplasm via SOC is pumped into the ER by SERCA. [29][30][31] Ca 2C entering the ER through the SOCE complex can thus be described as: 38 where V SOC is the maximum value of the SOC flux, and P SOC represents the fraction of STIM proteins bound to SOC, i.e., the fraction of activated SOC. P 1 SOC is the fraction of STIM proteins dissociated from the ER Ca 2C . K SOC is the affinity of STIM for [Ca 2C ] ER . t SOC represents the time constant for the adaptation of STIM to [Ca 2C ] ER changes. Parameter values for the above equations used in our model are listed in Table 3. The initial conditions (e.g., V m at t D 0) and geometric parameters of the model cell are listed in Table 4.

Disclosure of Potential Conflicts of Interest
No potential conflicts of interest were disclosed.

Funding
This work was supported by National Nature Science Foundation of China (grant no. 61375114). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. Funder's website: http://www. nsfc.gov.cn.