A dynamic model of the type-2 inositol trisphosphate receptor

Computational Model of Complex Calcium Dynamics: Store Operated Ca2+ Channels and Mitochondrial Associated Membranes in Pancreatic Acinar Cells

This proposed model explores the intricate Ca2+ dynamics within the pancreatic acinar cells (PACs) by emphasizing the role of store-operated Ca2+ entry (SOCE) and the mitochondrial-associated membranes (MAMs) in the secretory region (apical) of the PACs. Traditionally, Ca2+ releases from the endoplasmic reticulum (ER) via calcium-induced calcium release (CICR). It has been shown to be important in regulating functions such as secretion of digestive enzymes in PACs. However, this model posits that upon the depletion of Ca2+ in the ER, the signaling protein stromal interaction molecule (STIM1) is activated. Activated STIM1, then facilitates the opening of Orai channels, allowing Ca2+ influx through the store-operated calcium channels (SOCCs). The model highlights the complexity of the Ca2+ dynamics, and the importance of SOCE and MAMs in the PACs Ca2+ homeostasis. The numerical and bifurcation analysis illustrate how changes in agonist concentrations can lead to the diverse Ca2+ oscillation patterns, such as thin to broader oscillations, sinusoidal patterns, and baseline fluctuations, driven by the feedback mechanisms involving Ca2+ and inositol 1,4,5 trisphosphate (IP3). This understanding could have broader implications for cellular physiology and the development of therapies targeting Ca2+ signaling pathways.

The bioelectric mechanisms of local calcium dynamics in cancer cell proliferation: an extension of the A549 in silico cell model

Introduction

Advances in molecular targeting of ion channels may open up new avenues for therapeutic approaches in cancer based on the cells' bioelectric properties. In addition to in-vitro or in-vivo models, in silico models can provide deeper insight into the complex role of electrophysiology in cancer and reveal the impact of altered ion channel expression and the membrane potential on malignant processes. The A549 in silico model is the first computational cancer whole-cell ion current model that simulates the bioelectric mechanisms of the human non-small cell lung cancer cell line A549 during the different phases of the cell cycle. This work extends the existing model with a detailed mathematical description of the store-operated Ca2+ entry (SOCE) and the complex local intracellular calcium dynamics, which significantly affect the entire electrophysiological properties of the cell and regulate cell cycle progression.

Methods

The initial model was extended by a multicompartmental approach, addressing the heterogenous calcium profile and dynamics in the ER-PM junction provoked by local calcium entry of store-operated calcium channels (SOCs) and uptake by SERCA pumps. Changes of cytosolic calcium levels due to diffusion from the ER-PM junction, release from the ER by RyR channels and IP3 receptors, as well as corresponding PM channels were simulated and the dynamics evaluated based on calcium imaging data. The model parameters were fitted to available data from two published experimental studies, showing the function of CRAC channels and indirectly of IP3R, RyR and PMCA via changes of the cytosolic calcium levels.

Results

The proposed calcium description accurately reproduces the dynamics of calcium imaging data and simulates the SOCE mechanisms. In addition, simulations of the combined A549-SOCE model in distinct phases of the cell cycle demonstrate how Ca2+ - dynamics influence responding channels such as KCa, and consequently modulate the membrane potential accordingly.

Discussion

Local calcium distribution and time evolution in microdomains of the cell significantly impact the overall electrophysiological properties and exert control over cell cycle progression. By providing a more profound description, the extended A549-SOCE model represents an important step on the route towards a valid model for oncological research and in silico supported development of novel therapeutic strategies.

Mathematical Modeling of PI3K/Akt Pathway in Microglia

The motility of microglia involves intracellular signaling pathways that are predominantly controlled by changes in cytosolic Ca2+ and activation of PI3K/Akt (phosphoinositide-3-kinase/protein kinase B). In this letter, we develop a novel biophysical model for cytosolic Ca2+ activation of the PI3K/Akt pathway in microglia where Ca2+ influx is mediated by both P2Y purinergic receptors (P2YR) and P2X purinergic receptors (P2XR). The model parameters are estimated by employing optimization techniques to fit the model to phosphorylated Akt (pAkt) experimental modeling/in vitro data. The integrated model supports the hypothesis that Ca2+ influx via P2YR and P2XR can explain the experimentally reported biphasic transient responses in measuring pAkt levels. Our predictions reveal new quantitative insights into P2Rs on how they regulate Ca2+ and Akt in terms of physiological interactions and transient responses. It is shown that the upregulation of P2X receptors through a repetitive application of agonist results in a continual increase in the baseline [Ca2+], which causes the biphasic response to become a monophasic response which prolongs elevated levels of pAkt.

Likelihood-ratio test statistic for the finite-sample case in nonlinear ordinary differential equation models

Likelihood ratios are frequently utilized as basis for statistical tests, for model selection criteria and for assessing parameter and prediction uncertainties, e.g. using the profile likelihood. However, translating these likelihood ratios into p-values or confidence intervals requires the exact form of the test statistic's distribution. The lack of knowledge about this distribution for nonlinear ordinary differential equation (ODE) models requires an approximation which assumes the so-called asymptotic setting, i.e. a sufficiently large amount of data. Since the amount of data from quantitative molecular biology is typically limited in applications, this finite-sample case regularly occurs for mechanistic models of dynamical systems, e.g. biochemical reaction networks or infectious disease models. Thus, it is unclear whether the standard approach of using statistical thresholds derived for the asymptotic large-sample setting in realistic applications results in valid conclusions. In this study, empirical likelihood ratios for parameters from 19 published nonlinear ODE benchmark models are investigated using a resampling approach for the original data designs. Their distributions are compared to the asymptotic approximation and statistical thresholds are checked for conservativeness. It turns out, that corrections of the likelihood ratios in such finite-sample applications are required in order to avoid anti-conservative results.

Calcium oscillation on homogeneous and heterogeneous networks of ryanodine receptor

Calcium oscillation is an important calcium homeostasis, imbalance of which is the key mechanism of initiation and progression of many major diseases. The formation and maintenance of calcium homeostasis are closely related to the spatial distribution of calcium channels on endoplasmic reticulum, whose complex structure was unveiled by recent observations with superresolution imaging techniques. In the current paper, a theoretical framework is established by abstracting the spatial distribution of the calcium channels as a nonlinear biological complex network with calcium channels as nodes and Ca^{2+} as edges. A dynamical model for a ryanodine receptor (RyR) is adopted to investigate the effect of spatial distribution on calcium oscillation. The mean-field model can be well reproduced from the complete graph and dense Erdös-Rényi network. The synchronization of RyRs is found important to generate a global calcium oscillation. Below a critical density of the Erdös-Rényi or BaraBási-Albert network, the amplitude and interspike interval decrease rapidly with the end of disappearance of oscillation due to the desynchronization. The clique graph with a cluster structure cannot produce a global oscillation due to the failure of synchronization between clusters. A more realistic geometric network is constructed in a two-dimensional plane based on the experimental information about the RyR arrangement of clusters and the frequency distribution of cluster sizes. Different from the clique graph, the global oscillation can be generated with reasonable parameters on the geometric network. The simulation also suggests that existence of small clusters and rogue RyRs plays an important role in the maintenance of global calcium oscillation through keeping synchronization between large clusters. Such results support the heterogeneous distribution of RyRs with different-size clusters, which is helpful to understand recent observations with superresolution nanoscale imaging techniques. The current theoretical framework can also be extent to investigate other phenomena in calcium signal transduction.

Comparison of the Simulated Response of Three in Silico Human Stem Cell-Derived Cardiomyocytes Models and in Vitro Data Under 15 Drug Actions

Objectives: Improvements in human stem cell-derived cardiomyocyte (hSC-CM) technology have promoted their use for drug testing and disease investigations. Several in silico hSC-CM models have been proposed to augment interpretation of experimental findings through simulations. This work aims to assess the response of three hSC-CM in silico models (Koivumäki2018, Kernik2019, and Paci2020) to simulated drug action, and compare simulation results against in vitro data for 15 drugs.

Methods: First, simulations were conducted considering 15 drugs, using a simple pore-block model and experimental data for seven ion channels. Similarities and differences were analyzed in the in silico responses of the three models to drugs, in terms of Ca²⁺ transient duration (CTD₉₀) and occurrence of arrhythmic events. Then, the sensitivity of each model to different degrees of blockage of Na⁺ (I_Na), L-type Ca²⁺ (I_CaL), and rapid delayed rectifying K⁺ (I_Kr) currents was quantified. Finally, we compared the drug-induced effects on CTD₉₀ against the corresponding in vitro experiments.

Results: The observed CTD₉₀ changes were overall consistent among the in silico models, all three showing changes of smaller magnitudes compared to the ones measured in vitro. For example, sparfloxacin 10 µM induced +42% CTD₉₀ prolongation in vitro, and +17% (Koivumäki2018), +6% (Kernik2019), and +9% (Paci2020) in silico. Different arrhythmic events were observed following drug application, mainly for drugs affecting I_Kr. Paci2020 and Kernik2019 showed only repolarization failure, while Koivumäki2018 also displayed early and delayed afterdepolarizations. The spontaneous activity was suppressed by Na⁺ blockers and by drugs with similar effects on I_CaL and I_Kr in Koivumäki2018 and Paci2020, while only by strong I_CaL blockers, e.g. nisoldipine, in Kernik2019. These results were confirmed by the sensitivity analysis.

Conclusion: To conclude, The CTD₉₀ changes observed in silico are qualitatively consistent with our in vitro data, although our simulations show differences in drug responses across the hSC-CM models, which could stem from variability in the experimental data used in their construction.

Quantifying the dose-dependent impact of intracellular amyloid beta in a mathematical model of calcium regulation in xenopus oocyte

Alzheimer’s disease (AD) is a devastating illness affecting over 40 million people worldwide. Intraneuronal rise of amyloid beta in its oligomeric forms (iAβOs), has been linked to the pathogenesis of AD by disrupting cytosolic Ca²⁺ homeostasis. However, the specific mechanisms of action are still under debate and intense effort is ongoing to improve our understanding of the crucial steps involved in the mechanisms of AβOs toxicity. We report the development of a mathematical model describing a proposed mechanism by which stimulation of Phospholipase C (PLC) by iAβO, triggers production of IP₃ with consequent abnormal release of Ca²⁺ from the endoplasmic reticulum (ER) through activation of IP₃ receptor (IP₃R) Ca²⁺ channels. After validating the model using experimental data, we quantify the effects of intracellular rise in iAβOs on model solutions. Our model validates a dose-dependent influence of iAβOs on IP₃-mediated Ca²⁺ signaling. We investigate Ca²⁺ signaling patterns for small and large iAβOs doses and study the role of various parameters on Ca²⁺ signals. Uncertainty quantification and partial rank correlation coefficients are used to better understand how the model behaves under various parameter regimes. Our model predicts that iAβO alter IP₃R sensitivity to IP₃ for large doses. Our analysis also shows that the upstream production of IP₃ can influence Aβ-driven solution patterns in a dose-dependent manner. Model results illustrate and confirm the detrimental impact of iAβOs on IP₃ signaling.

Mechanisms Underlying Spontaneous Action Potential Generation Induced by Catecholamine in Pulmonary Vein Cardiomyocytes: A Simulation Study

Cardiomyocytes and myocardial sleeves dissociated from pulmonary veins (PVs) potentially generate ectopic automaticity in response to noradrenaline (NA), and thereby trigger atrial fibrillation. We developed a mathematical model of rat PV cardiomyocytes (PVC) based on experimental data that incorporates the microscopic framework of the local control theory of Ca²⁺ release from the sarcoplasmic reticulum (SR), which can generate rhythmic Ca²⁺ release (limit cycle revealed by the bifurcation analysis) when total Ca²⁺ within the cell increased. Ca²⁺ overload in SR increased resting Ca²⁺ efflux through the type II inositol 1,4,5-trisphosphate (IP₃) receptors (InsP₃R) as well as ryanodine receptors (RyRs), which finally triggered massive Ca²⁺ release through activation of RyRs via local Ca²⁺ accumulation in the vicinity of RyRs. The new PVC model exhibited a resting potential of −68 mV. Under NA effects, repetitive Ca²⁺ release from SR triggered spontaneous action potentials (APs) by evoking transient depolarizations (TDs) through Na⁺/Ca²⁺ exchanger (AP_TDs). Marked and variable latencies initiating AP_TDs could be explained by the time courses of the α1- and β1-adrenergic influence on the regulation of intracellular Ca²⁺ content and random occurrences of spontaneous TD activating the first AP_TD. Positive and negative feedback relations were clarified under AP_TD generation.

SemGen: a tool for semantics-based annotation and composition of biosimulation models

SUMMARY

As the number and complexity of biosimulation models grows, so do demands for tools that can help users understand models and compose more comprehensive and accurate systems from existing models. SemGen is a tool for semantics-based annotation and composition of biosimulation models designed to address this demand. A key SemGen capability is to decompose and then integrate models across existing model exchange formats including SBML and CellML. To support this capability, we use semantic annotations to explicitly capture the underlying biological and physical meanings of the entities and processes that are modeled. SemGen leverages annotations to expose a model's biological and computational architecture and to help automate model composition.

AVAILABILITY AND IMPLEMENTATION

SemGen is freely available at https://github.com/SemBioProcess/SemGen.

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

A mathematical model for the effects of amyloid beta on intracellular calcium

The accumulation of Alzheimer's disease (AD) associated Amyloid beta (Aβ) oligomers can trigger aberrant intracellular calcium (Ca2+) levels by disrupting the intrinsic Ca2+ regulatory mechanism within cells. These disruptions can cause changes in homeostasis levels that can have detrimental effects on cell function and survival. Although studies have shown that Aβ can interfere with various Ca2+ fluxes, the complexity of these interactions remains elusive. We have constructed a mathematical model that simulates Ca2+ patterns under the influence of Aβ. Our simulations shows that Aβ can increase regions of mixed-mode oscillations leading to aberrant signals under various conditions. We investigate how Aβ affects individual flux contributions through inositol triphosphate (IP3) receptors, ryanodine receptors, and membrane pores. We demonstrate that controlling for the ryanodine receptor's maximal kinetic reaction rate may provide a biophysical way of managing aberrant Ca2+ signals. The influence of a dynamic model for IP3 production is also investigated under various conditions as well as the impact of changes in membrane potential. Our model is one of the first to investigate the effects of Aβ on a variety of cellular mechanisms providing a base modeling scheme from which further studies can draw on to better understand Ca2+ regulation in an AD environment.

Mathematical Techniques for Understanding Platelet Regulation and the Development of New Pharmacological Approaches

Mathematical and computational modeling is currently in the process of becoming an accepted tool in the arsenal of methods utilized for the investigation of complex biological systems. For some problems in the field, like cellular metabolic regulation, neural impulse propagation, or cell cycle, progress is already unthinkable without use of such methods. Mathematical models of platelet signaling, function, and metabolism during the last years have not only been steadily increasing in their number, but have also been providing more in-depth insights, generating hypotheses, and allowing predictions to be made leading to new experimental designs and data. Here we describe the basic approaches to platelet mathematical model development and validation, highlighting the challenges involved. We then review the current theoretical models in the literature and how these are being utilized to increase our understanding of these complex cells.

A computational model of intracellular calcium oscillations in urinary bladder smooth muscle cells

This paper presents a computational model for calcium (Ca²⁺) oscillations in detrusor smooth muscle (DSM) cells. The proposed model simulates the temporal profile of oscillations by incorporating various cellular and subcellular components. The cellular components include calcium influx via membrane and the plasma membrane calcium ATPase (PMCA) pump. The subcellular components include ryanodine receptors (RyRs), inositol 1, 4, 5 trisphosphate receptors (IP₃Rs) and the sarco/endoplasmic reticulum calcium ATPase (SERCA) pump. The use of both cellular and subcellular components provides a better estimation of the origin and factors affecting these oscillations. Moreover, our work correlates these computational findings with associated physiology of the smooth muscle cell that aids our understanding of intracellular calcium oscillations and its inception in DSM. A deeper insight into calcium signalling in DSM cells is expected to provide a firmer basis for understanding the mechanical contractile activity.

Comparison of gating dynamics of different IP3R channels with immune algorithm searching for channel parameter distributions

The gating properties of the inositol 1, 4, 5-trisphosphate (IP₃) receptor (IP₃R) are determined by the binding and unbinding capability of Ca²⁺ ions and IP₃ messengers. With the patch clamp experiments, the stationary properties have been discussed for Xenopus oocyte type-1 IP₃R (Oo-IP₃R1), type-3 IP₃R (Oo-IP₃R3) and Spodoptera frugiperda IP₃R (Sf-IP₃R). In this paper, in order to provide insights about the relation between the observed gating characteristics and the gating parameters in different IP₃Rs, we apply the immune algorithm to fit the parameters of a modified DeYoung-Keizer model. By comparing the fitting parameter distributions of three IP₃Rs, we suggest that the three types of IP₃Rs have the similar open sensitivity in responding to IP₃. The Oo-IP₃R3 channel is easy to open in responding to low Ca²⁺ concentration, while Sf-IP₃R channel is easily inhibited in responding to high Ca²⁺ concentration. We also show that the IP₃ binding rate is not a sensitive parameter for stationary gating dynamics for three IP₃Rs, but the inhibitory Ca²⁺ binding/unbinding rates are sensitive parameters for gating dynamics for both Oo-IP₃R1 and Oo-IP₃R3 channels. Such differences may be important in generating the spatially and temporally complex Ca²⁺ oscillations in cells. Our study also demonstrates that the immune algorithm can be applied for model parameter searching in biological systems.

Calcium Dynamics of Ex Vivo Long-Term Cultured CD8+ T Cells Are Regulated by Changes in Redox Metabolism

T cells reach a state of replicative senescence characterized by a decreased ability to proliferate and respond to foreign antigens. Calcium release associated with TCR engagement is widely used as a surrogate measure of T cell response. Using an ex vivo culture model that partially replicates features of organismal aging, we observe that while the amplitude of Ca2+ signaling does not change with time in culture, older T cells exhibit faster Ca2+ rise and a faster decay. Gene expression analysis of Ca2+ channels and pumps expressed in T cells by RT-qPCR identified overexpression of the plasma membrane CRAC channel subunit ORAI1 and PMCA in older T cells. To test whether overexpression of the plasma membrane Ca2+ channel is sufficient to explain the kinetic information, we adapted a previously published computational model by Maurya and Subramaniam to include additional details on the store-operated calcium entry (SOCE) process to recapitulate Ca2+ dynamics after T cell receptor stimulation. Simulations demonstrated that upregulation of ORAI1 and PMCA channels is not sufficient to explain the observed alterations in Ca2+ signaling. Instead, modeling analysis identified kinetic parameters associated with the IP3R and STIM1 channels as potential causes for alterations in Ca2+ dynamics associated with the long term ex vivo culturing protocol. Due to these proteins having known cysteine residues susceptible to oxidation, we subsequently investigated and observed transcriptional remodeling of metabolic enzymes, a shift to more oxidized redox couples, and post-translational thiol oxidation of STIM1. The model-directed findings from this study highlight changes in the cellular redox environment that may ultimately lead to altered T cell calcium dynamics during immunosenescence or organismal aging.

Modulation and pre-amplification of PAR1 signaling by ADP acting via the P2Y12 receptor during platelet subpopulation formation

BACKGROUND

Two major soluble blood platelet activators are thrombin and ADP. Of these two, only thrombin can induce mitochondrial collapse and programmed cell death leading to phosphatidylserine (PS) exposure required for blood clotting reactions acceleration. Thrombin can also greatly potentiate collagen-induced PS exposure. However, ADP acting through the P2Y12 receptor was shown to increase the PS-exposing (PS+) platelets fraction produced by thrombin or thrombin-plus-collagen via an unknown mechanism.

METHODS

We developed a comprehensive multicompartmental computational model of platelet PAR1-and-P2Y12 calcium signal transduction that included cytoplasmic signaling, dense tubular system and mitochondria. To test model predictions, flow cytometry experiments with washed, annexin V-labeled platelets were performed.

RESULTS

Stimulation of thrombin receptor PAR1 in the model induced cytoplasmic calcium oscillations, calcium uptake by mitochondria, opening of the permeability transition pore and collapse of the mitochondrial membrane potential. ADP stimulation of P2Y12 led to cAMP decrease that, in turn, caused changes in phospholipase C phosphorylation by protein kinase A, increase in cytoplasmic calcium level and, consequently, PS+ platelet formation. ADP addition before stimulation of PAR1 produced much greater increase of the PS+ fraction because cAMP concentration had time to go down prior to calcium oscillations; this prediction was also tested and confirmed experimentally.

CONCLUSION

These results suggest a mechanism of ADP-dependent PS exposure regulation and show a likely mode of action that could be important for the PS exposure regulation in thrombi, where ADP is released before thrombin formation.

The role of IP3 receptor channel clustering in Ca2+ wave propagation during oocyte maturation

During oocyte maturation, the calcium-signaling machinery undergoes a dramatic remodeling resulting in distinctly different calcium-release patterns on all organizational scales from puffs to waves. The dynamics of the Ca(2+) release wave in mature as compared to immature oocytes are defined by a slower propagation speed and longer duration of the high Ca(2+) plateau. In this chapter, we use computational modeling to identify the changes in the signaling machinery, which contribute most significantly to the alterations observed in Ca(2+) wave propagation during Xenopus oocyte maturation. In addition to loss of store-operated calcium entry and internalization of plasma membrane pumps, we propose that spatial reorganization of the IP3 receptors in the plane of the ER membrane is a key factor for the observed signaling changes in Ca(2+) wave propagation.

Compartmentalized calcium signaling triggers subpopulation formation upon platelet activation through PAR1

A computational model of PAR1-stimulated platelet calcium signaling is developed to analyze the formation of platelet subpopulations. This occurs via a mitochondria-dependent decision-making mechanism. This is a stochastic phenomenon caused by a small number of PARs.

Systems modeling of Ca(2+) homeostasis and mobilization in platelets mediated by IP3 and store-operated Ca(2+) entry

Resting platelets maintain a stable level of low cytoplasmic calcium ([Ca(2+)]cyt) and high dense tubular system calcium ([Ca(2+)]dts). During thrombosis, activators cause a transient rise in inositol trisphosphate (IP3) to trigger calcium mobilization from stores and elevation of [Ca(2+)]cyt. Another major source of [Ca(2+)]cyt elevation is store-operated calcium entry (SOCE) through plasmalemmal calcium channels that open in response to store depletion as [Ca(2+)]dts drops. A 34-species systems model employed kinetics describing IP3-receptor, DTS-plasmalemma puncta formation, SOCE via assembly of STIM1 and Orai1, and the plasmalemma and sarco/endoplasmic reticulum Ca(2+)-ATPases. Four constraints were imposed: calcium homeostasis before activation; stable in zero extracellular calcium; IP3-activatable; and functional SOCE. Using a Monte Carlo method to sample three unknown parameters and nine initial concentrations in a 12-dimensional space near measured or expected values, we found that model configurations that were responsive to stimuli and demonstrated significant SOCE required high inner membrane electric potential (>-70 mV) and low resting IP3 concentrations. The absence of puncta in resting cells was required to prevent spontaneous store depletion in calcium-free media. Ten-fold increases in IP3 caused saturated calcium mobilization. This systems model represents a critical step in being able to predict platelets' phenotypes during hemostasis or thrombosis.

An allosteric model of the inositol trisphosphate receptor with nonequilibrium binding

The inositol trisphosphate receptor (IPR) is a crucial ion channel that regulates the Ca(2+) influx from the endoplasmic reticulum (ER) to the cytoplasm. A thorough study of the IPR channel contributes to a better understanding of calcium oscillations and waves. It has long been observed that the IPR channel is a typical biological system which performs adaptation. However, recent advances on the physical essence of adaptation show that adaptation systems with a negative feedback mechanism, such as the IPR channel, must break detailed balance and always operate out of equilibrium with energy dissipation. Almost all previous IPR models are equilibrium models assuming detailed balance and thus violate the dissipative nature of adaptation. In this article, we constructed a nonequilibrium allosteric model of single IPR channels based on the patch-clamp experimental data obtained from the IPR in the outer membranes of isolated nuclei of the Xenopus oocyte. It turns out that our model reproduces the patch-clamp experimental data reasonably well and produces both the correct steady-state and dynamic properties of the channel. Particularly, our model successfully describes the complicated bimodal [Ca(2+)] dependence of the mean open duration at high [IP3], a steady-state behavior which fails to be correctly described in previous IPR models. Finally, we used the patch-clamp experimental data to validate that the IPR channel indeed breaks detailed balance and thus is a nonequilibrium system which consumes energy.

Modelling the transition from simple to complex Ca²⁺ oscillations in pancreatic acinar cells

A mathematical model is proposed which systematically investigates complex calcium oscillations in pancreatic acinar cells. This model is based on calcium-induced calcium release via inositol trisphosphate receptors (IPR) and ryanodine receptors (RyR) and includes calcium modulation of inositol (1,4,5) trisphosphate (IP3) levels through feedback regulation of degradation and production. In our model, the apical and the basal regions are separated by a region containing mitochondria, which is capable of restricting Ca2+ responses to the apical region. We were able to reproduce the observed oscillatory patterns, from baseline spikes to sinusoidal oscillations. The model predicts that calcium-dependent production and degradation of IP3 is a key mechanism for complex calcium oscillations in pancreatic acinar cells. A partial bifurcation analysis is performed which explores the dynamic behaviour of the model in both apical and basal regions.

Translating intracellular calcium signaling into models

The rich experimental data on intracellular calcium has put theoreticians in an ideal position to derive models of intracellular calcium signaling. Over the last 25 years, a large number of modeling frameworks have been suggested. Here, I will review some of the milestones of intracellular calcium modeling with a special emphasis on calcium-induced calcium release (CICR) through inositol-1,4,5-trisphosphate and ryanodine receptors. I will highlight key features of CICR and how they are represented in models as well as the challenges that theoreticians face when translating our current understanding of calcium signals into equations. The selected examples demonstrate that a successful model provides mechanistic insights into the molecular machinery of the Ca²⁺ signaling toolbox and determines the contribution of local Ca²⁺ release to global Ca²⁺ patterns, which at the moment cannot be resolved experimentally.

Calcium and IP3 dynamics in cardiac myocytes: experimental and computational perspectives and approaches

Calcium plays a crucial role in excitation-contraction coupling (ECC), but it is also a pivotal second messenger activating Ca²⁺-dependent transcription factors in a process termed excitation-transcription coupling (ETC). Evidence accumulated over the past decade indicates a pivotal role of inositol 1,4,5-trisphosphate receptor (IP₃R)-mediated Ca²⁺ release in the regulation of cytosolic and nuclear Ca²⁺ signals. IP₃ is generated by stimulation of plasma membrane receptors that couple to phospholipase C (PLC), liberating IP₃ from phosphatidylinositol 4,5-bisphosphate (PIP₂). An intriguing aspect of IP₃ signaling is the presence of the entire PIP₂-PLC-IP₃ signaling cascade as well as the presence of IP₃Rs at the inner and outer membranes of the nuclear envelope (NE) which functions as a Ca²⁺ store. The observation that the nucleus is surrounded by its own putative Ca²⁺ store raises the possibility that nuclear IP₃-dependent Ca²⁺ release plays a critical role in ETC. This provides a potential mechanism of regulation that acts locally and autonomously from the global cytosolic Ca²⁺ signal underlying ECC. Moreover, there is evidence that: (i) the sarcoplasmic reticulum (SR) and NE are a single contiguous Ca²⁺ store; (ii) the nuclear pore complex is the major gateway for Ca²⁺ and macromolecules to pass between the cytosol and the nucleoplasm; (iii) the inner membrane of the NE hosts key Ca²⁺ handling proteins including the Na⁺/Ca²⁺ exchanger (NCX)/GM1 complex, ryanodine receptors (RyRs), nicotinic acid adenine dinucleotide phosphate receptors (NAADPRs), Na⁺/K⁺ ATPase, and Na⁺/H⁺ exchanger. Thus, it appears that the nucleus represents a Ca²⁺ signaling domain equipped with its own ion channels and transporters that allow for complex local Ca²⁺ signals. Many experimental and modeling approaches have been used for the study of intracellular Ca²⁺ signaling but the key to the understanding of the dual role of Ca²⁺ mediating ECC and ECT lays in quantitative differences of local [Ca²⁺] in the nuclear and cytosolic compartment. In this review, we discuss the state of knowledge regarding the origin and the physiological implications of nuclear Ca²⁺ transients in different cardiac cell types (adult atrial and ventricular myocytes) as well as experimental and mathematical approaches to study Ca²⁺ and IP₃ signaling in the cytosol and nucleus. In particular, we focus on the concept that highly localized Ca²⁺ signals are required to translocate and activate Ca²⁺-dependent transcription factors (e.g., nuclear factor of activated T-cells, NFAT; histone deacetylase, HDAC) through phosphorylation/dephosphorylation processes.

Systems biology of platelet-vessel wall interactions

Blood systems biology seeks to quantify outside-in signaling as platelets respond to numerous external stimuli, typically under flow conditions. Platelets can activate via GPVI collagen receptor and numerous G-protein coupled receptors (GPCRs) responsive to ADP, thromboxane, thrombin, and prostacyclin. A bottom-up ODE approach allowed prediction of platelet calcium and phosphoinositides following P2Y₁ activation with ADP, either for a population average or single cell stochastic behavior. The homeostasis assumption (i.e., a resting platelet stays resting until activated) was particularly useful in finding global steady states for these large metabolic networks. Alternatively, a top-down approach involving Pairwise Agonist Scanning (PAS) allowed large data sets of measured calcium mobilization to predict an individual's platelet responses. The data was used to train neural network (NN) models of signaling to predict patient-specific responses to combinatorial stimulation. A kinetic description of platelet signaling then allows prediction of inside-out activation of platelets as they experience the complex biochemical milieu at the site of thrombosis. Multiscale lattice kinetic Monte Carlo (LKMC) utilizes these detailed descriptions of platelet signaling under flow conditions where released soluble species are solved by finite element method and the flow field around the growing thrombus is updated using computational fluid dynamics or lattice Boltzmann method. Since hemodynamic effects are included in a multiscale approach, thrombosis can then be predicted under arterial and venous thrombotic conditions for various anatomical geometries. Such systems biology approaches accommodate the effect of anti-platelet pharmacological intervention where COX1 pathways or ADP signaling are modulated in a patient-specific manner.

The most precise computations using Euler's method in standard floating-point arithmetic applied to modelling of biological systems

The explicit Euler's method is known to be very easy and effective in implementation for many applications. This article extends results previously obtained for the systems of linear differential equations with constant coefficients to arbitrary systems of ordinary differential equations. Optimal (providing minimum total error) step size is calculated at each step of Euler's method. Several examples of solving stiff systems are included.

Comparison of models for IP3 receptor kinetics using stochastic simulations

Inositol 1,4,5-trisphosphate receptor (IP₃R) is a ubiquitous intracellular calcium (Ca²⁺) channel which has a major role in controlling Ca²⁺ levels in neurons. A variety of computational models have been developed to describe the kinetic function of IP₃R under different conditions. In the field of computational neuroscience, it is of great interest to apply the existing models of IP₃R when modeling local Ca²⁺ transients in dendrites or overall Ca²⁺ dynamics in large neuronal models. The goal of this study was to evaluate existing IP₃R models, based on electrophysiological data. This was done in order to be able to suggest suitable models for neuronal modeling. Altogether four models (Othmer and Tang, 1993; Dawson etal., 2003; Fraiman and Dawson, 2004; Doi etal., 2005) were selected for a more detailed comparison. The selection was based on the computational efficiency of the models and the type of experimental data that was used in developing the model. The kinetics of all four models were simulated by stochastic means, using the simulation software STEPS, which implements the Gillespie stochastic simulation algorithm. The results show major differences in the statistical properties of model functionality. Of the four compared models, the one by Fraiman and Dawson (2004) proved most satisfactory in producing the specific features of experimental findings reported in literature. To our knowledge, the present study is the first detailed evaluation of IP₃R models using stochastic simulation methods, thus providing an important setting for constructing a new, realistic model of IP₃R channel kinetics for compartmental modeling of neuronal functions. We conclude that the kinetics of IP₃R with different concentrations of Ca²⁺ and IP₃ should be more carefully addressed when new models for IP₃R are developed.

A data-driven model of a modal gated ion channel: the inositol 1,4,5-trisphosphate receptor in insect Sf9 cells

The inositol 1,4,5-trisphosphate (IP(3)) receptor (IP(3)R) channel is crucial for the generation and modulation of intracellular Ca(2+) signals in animal cells. To gain insight into the complicated ligand regulation of this ubiquitous channel, we constructed a simple quantitative continuous-time Markov-chain model from the data. Our model accounts for most experimentally observed gating behaviors of single native IP(3)R channels from insect Sf9 cells. Ligand (Ca(2+) and IP(3)) dependencies of channel activity established six main ligand-bound channel complexes, where a complex consists of one or more states with the same ligand stoichiometry and open or closed conformation. Channel gating in three distinct modes added one complex and indicated that three complexes gate in multiple modes. This also restricted the connectivity between channel complexes. Finally, latencies of channel responses to abrupt ligand concentration changes defined a model with specific network topology between 9 closed and 3 open states. The model with 28 parameters can closely reproduce the equilibrium gating statistics for all three gating modes over a broad range of ligand concentrations. It also captures the major features of channel response latency distributions. The model can generate falsifiable predictions of IP(3)R channel gating behaviors and provide insights to both guide future experiment development and improve IP(3)R channel gating analysis. Maximum likelihood estimates of the model parameters and of the parameters in the De Young-Keizer model yield strong statistical evidence in favor of our model. Our method is simple and easily applicable to the dynamics of other ion channels and molecules.

Stochastic simulations of pattern formation in excitable media

We present a method for mesoscopic, dynamic Monte Carlo simulations of pattern formation in excitable reaction-diffusion systems. Using a two-level parallelization approach, our simulations cover the whole range of the parameter space, from the noise-dominated low-particle number regime to the quasi-deterministic high-particle number limit. Three qualitatively different case studies are performed that stand exemplary for the wide variety of excitable systems. We present mesoscopic stochastic simulations of the Gray-Scott model, of a simplified model for intracellular Ca2+ oscillations and, for the first time, of the Oregonator model. We achieve simulations with up to 10(10) particles. The software and the model files are freely available and researchers can use the models to reproduce our results or adapt and refine them for further exploration.

Multi-scale data-driven modeling and observation of calcium puffs

The spatiotemporal dynamics of elementary Ca(2+) release events, such as "blips" and "puffs" shapes the hierarchal Ca(2+) signaling in many cell types. Despite being the building blocks of Ca(2+) patterning, the mechanism responsible for the observed properties of puffs, especially their termination is incompletely understood. In this paper, we employ a data-driven approach to gain insights into the complex dynamics of blips and puffs. We use a model of inositol 1,4,5-trisphosphate (IP(3)) receptor (IP(3)R) derived directly from single channel patch clamp data taken at 10 μM concentration of IP(3) to simulate calcium puffs. We first reproduce recent observations regarding puffs and blips and then investigate the mechanism of puff termination. Our model suggests that during a puff, IP(3)R s proceed around a loop through kinetic states from "rest" to "open" to "inhibited" and back to "rest". A puff terminates because of self-inhibition. Based on our simulations, we rule out the endoplasmic reticulum (ER) Ca(2+) depletion as a possible cause for puff termination. The data-driven approach also enables us to estimate the current through a single IP(3)R and the peak Ca(2+) concentration near the channel pore.

Termination of Ca²+ release for clustered IP₃R channels

In many cell types, release of calcium ions is controlled by inositol 1,4,5-trisphosphate () receptor channels. Elevations in concentration after intracellular release through receptors (R) can either propagate in the form of waves spreading through the entire cell or produce spatially localized puffs. The appearance of waves and puffs is thought to implicate random initial openings of one or a few channels and subsequent activation of neighboring channels because of an “autocatalytic” feedback. It is much less clear, however, what determines the further time course of release, particularly since the lifetime is very different for waves (several seconds) and puffs (around 100 ms). Here we study the lifetime of signals and their dependence on residual microdomains. Our general idea is that microdomains are dynamical and mediate the effect of other physiological processes. Specifically, we focus on the mechanism by which binding proteins (buffers) alter the lifetime of signals. We use stochastic simulations of channel gating coupled to a coarse-grained description for the concentration. To describe the concentration in a phenomenological way, we here introduce a differential equation, which reflects the buffer characteristics by a few effective parameters. This non-stationary model for microdomains gives deep insight into the dynamical differences between puffs and waves. It provides a novel explanation for the different lifetimes of puffs and waves and suggests that puffs are terminated by inhibition while unbinding is responsible for termination of waves. Thus our analysis hints at an additional role of and shows how cells can make use of the full complexity in R gating behavior to achieve different signals.

Calcium signals are important for a host of cellular processes such as neurotransmitter release, cell contraction and gene expression. While the principles of activation and spreading of calcium signals have been largely understood, it is much less clear how their spatio-temporal appearance is shaped. This issue is of high relevance since the spatio-temporal signature is thought to carry the information content. In our paper we study the dynamical mechanisms that determine the time course of calcium release from receptor channels. We use a stochastic channel description combined with a recently developed model for the distribution of released calcium in a microdomain. The simulations uncover a complex control mechanism, which allows for the tuning of release from short frequent puffs to extended and less frequent wave-like release. Unexpectedly, the model predicts that for wave-like release the dissociation of from the receptors leads to termination of the calcium signal. This effect relies on a well-known gating property of R channels, which earlier has been regarded as superfluous in studies for groups of channels. Our results also provide a missing link to understand cellular response to calcium-binding proteins and present a novel mechanism for information processing by R channels.

Modelling the effect of specific inositol 1,4,5-trisphosphate receptor isoforms on cellular Ca2+ signals

BACKGROUND INFORMATION

Oscillations of cytosolic Ca2+ are well-known to rely on the regulatory properties of the InsP3R (inositol 1,4,5-trisphosphate receptor). Three isoforms of this channel have been identified. They differ in their regulatory properties by Ca2+ and InsP3. Experiments in different cell types clearly indicate that the relative amounts of each isoform affect the time course of Ca2+ changes after agonist stimulation. In the present study, we investigate whether different steady-state curves for the open probability of the InsP3Rs as a function of Ca2+ imply different dynamical behaviours when these receptors are present in a cellular environment. We therefore describe by a specific phenomenological model the three main types of curves that have been reported: (i) the classical bell-shaped curve, (ii) the bell-shaped curve that is shifted towards higher Ca2+ concentrations when InsP3 is increased, and (iii) a monotonous increasing function of cytosolic Ca2+.

RESULTS

We show that, although these types of curves can be ascribed to slight differences in the channel regulation by Ca2+ and InsP3, they can indicate important variations as to the receptor role in cellular Ca2+ control. Thus the receptor associated with the classical bell-shaped curve appears to be the most robust Ca2+ oscillator. If the steady-state curve is supposed to be a monotonous increasing function of cytosolic Ca2+, the modelled receptor cannot sustain Ca2+ oscillations in the absence of Ca2+ exchanges with the extracellular medium. When the bell-shaped curve is shifted towards higher Ca2+ concentrations with increasing InsP3 levels, the model predicts that the receptor is less robust to changes in density; this receptor, however, provides a finer control of the steady-state level of Ca2+ when varying the InsP3 concentration.

CONCLUSIONS

Our model allows us to propose an explanation for the experimental observations about the effect of selectively expressing or down-regulating InsP3R isoforms, as well as to make theoretical predictions.

Endoplasmic reticulum remodeling tunes IP₃-dependent Ca²+ release sensitivity

The activation of vertebrate development at fertilization relies on IP₃-dependent Ca²⁺ release, a pathway that is sensitized during oocyte maturation. This sensitization has been shown to correlate with the remodeling of the endoplasmic reticulum into large ER patches, however the mechanisms involved are not clear. Here we show that IP₃ receptors within ER patches have a higher sensitivity to IP₃ than those in the neighboring reticular ER. The lateral diffusion rate of IP₃ receptors in both ER domains is similar, and ER patches dynamically fuse with reticular ER, arguing that IP₃ receptors exchange freely between the two ER compartments. These results suggest that increasing the density of IP₃ receptors through ER remodeling is sufficient to sensitize IP₃-dependent Ca²⁺ release. Mathematical modeling supports this concept of ‘geometric sensitization’ of IP₃ receptors as a population, and argues that it depends on enhanced Ca²⁺-dependent cooperativity at sub-threshold IP₃ concentrations. This represents a novel mechanism of tuning the sensitivity of IP₃ receptors through ER remodeling during meiosis.

Intercellular calcium waves in glial cells with bistable dynamics

A two-dimensional model is proposed for intercellular calcium (Ca(2 +)) waves with Ca(2 +)-induced IP(3) regeneration and the diffusion of IP(3) through gap junctions. Many experimental observations in glial cells, i.e. responding to local mechanical stimulation, glutamate application, mechanical stimulation followed by ACh application, and glutamate followed by mechanical stimulation, are reproduced and classified by the model. We show that a glial cell model with bistable dynamics, i.e. a Ca(2 +) oscillation state coexisting with a fixed point, can cause a prolonged plateau of Ca(2 +) signals in the cells nearby the stimulated cell when the cell network responds to the local mechanical stimulation.

Designing dynamical output feedback controllers for store-operated Ca²+ entry

Store-operated calcium entry (SOCE) has been proposed as the main process controlling Ca²+ entry in non-excitable cells. Although recent breakthroughs in experimental studies of SOCE have been made, its mathematical modeling has not been developed. In the present work, SOCE is viewed as a feedback control system subject to an extracellular agonist disturbance and an extracellular calcium input. We then design a dynamic output feedback controller to reject the disturbance and track Ca²+ resting levels in the cytosol and the endoplasmic reticulum (ER). The constructed feedback control system is validated by published experimental data and its global asymptotic stability is proved by using the LaSalle's invariance principle. We then simulate the dynamic responses of STIM1 and Orai1, two major components in the operation of the store-operated channels, to the depletion of Ca²+ in the ER with thapsigargin, which show that: (1) Upon the depletion of Ca²+ in the ER, the concentrations of activated STIM1 and STIM1-Orai1 cluster are elevated gradually, indicating that STIM1 is accumulating in the ER-PM junctions and that the cytosolic portion of the active STIM1 is binding to Orai1 and driving the opening of CRAC channels for Ca²+ entry; (2) after the extracellular Ca²+ addition, the concentrations of both STIM1 and STIM1-Orai1 cluster decrease but still much higher than the original levels. We also simulate the system responses to the agonist disturbance, which show that, when a sequence of periodic agonist pulses is applied, the system returns to its equilibrium after each pulse. This indicates that the designed feedback controller can reject the disturbance and track the equilibrium.

Two slow calcium-activated afterhyperpolarization currents control burst firing dynamics in gonadotropin-releasing hormone neurons

Gonadotropin-releasing hormone (GnRH) neurons release GnRH in a pulsatile manner to control fertility in all mammals. The mechanisms underlying burst firing in GnRH neurons, thought to contribute to pulsatile GnRH release, are not yet understood. Using minimally invasive, dual electrical-calcium recordings in acute brain slices from GnRH-Pericam transgenic mice, we find that the soma/proximal dendrites of GnRH neurons exhibit long-duration (approximately 10 s) calcium transients that are perfectly synchronized with their burst firing. These transients were found to be generated by calcium entry through voltage-dependent L-type calcium channels that was amplified by inositol-1,4,5-trisphosphate receptor-dependent store mechanisms. Perforated-patch current- and voltage-clamp electrophysiology coupled with mathematical modeling approaches revealed that these broad calcium transients act to control two slow afterhyperpolarization currents (sI(AHP)) in GnRH neurons: a quick-activating apamin-sensitive sI(AHP) that regulates both intraburst and interburst dynamics, and a slow-onset UCL2077-sensitive sI(AHP) that regulates only interburst dynamics. These observations highlight a unique interplay between electrical activity, calcium dynamics, and multiple calcium-regulated sI(AHP)s critical for shaping GnRH neuron burst firing.

Modeling and analysis of the molecular basis of pain in sensory neurons

Intracellular calcium dynamics are critical to cellular functions like pain transmission. Extracellular ATP plays an important role in modulating intracellular calcium levels by interacting with the P2 family of surface receptors. In this study, we developed a mechanistic mathematical model of ATP-induced P2 mediated calcium signaling in archetype sensory neurons. The model architecture, which described 90 species connected by 162 interactions, was formulated by aggregating disparate molecular modules from literature. Unlike previous models, only mass action kinetics were used to describe the rate of molecular interactions. Thus, the majority of the 252 unknown model parameters were either association, dissociation or catalytic rate constants. Model parameters were estimated from nine independent data sets taken from multiple laboratories. The training data consisted of both dynamic and steady-state measurements. However, because of the complexity of the calcium network, we were unable to estimate unique model parameters. Instead, we estimated a family or ensemble of probable parameter sets using a multi-objective thermal ensemble method. Each member of the ensemble met an error criterion and was located along or near the optimal trade-off surface between the individual training data sets. The model quantitatively reproduced experimental measurements from dorsal root ganglion neurons as a function of extracellular ATP forcing. Hypothesized architecture linking phosphoinositide regulation with P2X receptor activity explained the inhibition of P2X-mediated current flow by activated metabotropic P2Y receptors. Sensitivity analysis using individual and the whole system outputs suggested which molecular subsystems were most important following P2 activation. Taken together, modeling and analysis of ATP-induced P2 mediated calcium signaling generated qualitative insight into the critical interactions controlling ATP induced calcium dynamics. Understanding these critical interactions may prove useful for the design of the next generation of molecular pain management strategies.

Inositol trisphosphate receptor and ion channel models based on single-channel data

The inositol trisphosphate receptor (IPR) plays an important role in controlling the dynamics of intracellular Ca(2+). Single-channel patch-clamp recordings are a typical way to study these receptors as well as other ion channels. Methods for analyzing and using this type of data have been developed to fit Markov models of the receptor. The usual method of parameter fitting is based on maximum-likelihood techniques. However, Bayesian inference and Markov chain Monte Carlo techniques are becoming more popular. We describe the application of the Bayesian methods to real experimental single-channel data in three ion channels: the ryanodine receptor, the K(+) channel, and the IPR. One of the main aims of all three studies was that of model selection with different approaches taken. We also discuss the modeling implications for single-channel data that display different levels of channel activity within one recording.

Toward a predictive model of Ca2+ puffs

We investigate the key characteristics of Ca(2+) puffs in deterministic and stochastic frameworks that all incorporate the cellular morphology of IP(3) receptor channel clusters. In the first step, we numerically study the Ca(2+) liberation in a three-dimensional representation of a cluster environment with reaction-diffusion dynamics in both the cytosol and the lumen. These simulations reveal that Ca(2+) concentrations at a releasing cluster range from 80 to 170 microM and equilibrate almost instantaneously on the time scale of the release duration. These highly elevated Ca(2+) concentrations eliminate Ca(2+) oscillations in a deterministic model of an IP(3)R channel cluster at physiological parameter values as revealed by a linear stability analysis. The reason lies in the saturation of all feedback processes in the IP(3)R gating dynamics, so that only fluctuations can restore experimentally observed Ca(2+) oscillations. In this spirit, we derive master equations that allow us to analytically quantify the onset of Ca(2+) puffs and hence the stochastic time scale of intracellular Ca(2+) dynamics. Moving up the spatial scale, we suggest to formulate cellular dynamics in terms of waiting time distribution functions. This approach prevents the state space explosion that is typical for the description of cellular dynamics based on channel states and still contains information on molecular fluctuations. We illustrate this method by studying global Ca(2+) oscillations.

An investigation of models of the IP3R channel in Xenopus oocyte

We consider different models of inositol 1,4,5-trisphosphate (IP(3)) receptor (IP(3)R) channels in order to fit nuclear membrane patch clamp data of the stationary open probability, mean open time, and mean close time of channels in the Xenopus oocyte. Our results indicate that rather than to treat the tetrameric IP(3)R as four independent and identical subunits, one should assume sequential binding-unbinding processes of Ca(2+) ions and IP(3) messengers. Our simulations also favor the assumption that a channel opens through a conformational transition from a close state to an active state.

A kinetic model of the inositol trisphosphate receptor based on single-channel data

In many cell types, the inositol trisphosphate receptor is one of the important components controlling intracellular calcium dynamics, and an understanding of this receptor is necessary for an understanding of calcium oscillations and waves. Based on single-channel data from the type-I inositol trisphosphate receptor, and using a Markov chain Monte Carlo approach, we show that the most complex time-dependent model that can be unambiguously determined from steady-state data is one with three closed states and one open state, and we determine how the rate constants depend on calcium. Because the transitions between these states are complex functions of calcium concentration, each model state must correspond to a group of physical states. We fit two different topologies and find that both models predict that the main effect of [Ca(2+)] is to modulate the probability that the receptor is in a state that is able to open, rather than to modulate the transition rate to the open state.

A computational analysis of localized Ca2+-dynamics generated by heterogeneous release sites

We investigate the role of heterogeneous expression of IP3R and RyR in generating diverse elementary Ca2+ signals. It has been shown empirically (Wojcikiewicz and Luo in Mol. Pharmacol. 53(4):656-662, 1998; Newton et al. in J. Biol. Chem. 269(46):28613-28619, 1994; Smedt et al. in Biochem. J. 322(Pt. 2):575-583, 1997) that tissues express various proportions of IP3 and RyR isoforms and this expression is dynamically regulated (Parrington et al. in Dev. Biol. 203(2):451-461, 1998; Fissore et al. in Biol. Reprod. 60(1):49-57, 1999; Tovey et al. in J. Cell Sci. 114(Pt. 22):3979-3989, 2001). Although many previous theoretical studies have investigated the dynamics of localized calcium release sites (Swillens et al. in Proc. Natl. Acad. Sci. U.S.A. 96(24):13750-13755, 1999; Shuai and Jung in Proc. Natl. Acad. Sci. U.S.A. 100(2):506-510, 2003a; Shuai and Jung in Phys. Rev. E, Stat. Nonlinear Soft Matter Phys. 67(3 Pt. 1):031905, 2003b; Thul and Falcke in Biophys. J. 86(5):2660-2673, 2004; DeRemigio and Smith in Cell Calcium 38(2):73-86, 2005; Nguyen et al. in Bull. Math. Biol. 67(3):393-432, 2005), so far all such studies focused on release sites consisting of identical channel types. We have extended an existing mathematical model (Nguyen et al. in Bull. Math. Biol. 67(3):393-432, 2005) to release sites with two (or more) receptor types, each with its distinct channel kinetics. Mathematically, the release site is represented by a transition probability matrix for a collection of nonidentical stochastically gating channels coupled through a shared Ca2+ domain. We demonstrate that under certain conditions a previously defined mean-field approximation of the coupling strength does not accurately reproduce the release site dynamics. We develop a novel approximation and establish that its performance in these instances is superior. We use this mathematical framework to study the effect of heterogeneity in the Ca2+-regulation of two colocalized channel types on the release site dynamics. We consider release sites consisting of channels with both Ca2+-activation and inactivation ("four-state channels") and channels with Ca2+-activation only ("two-state channels") and show that for the appropriate parameter values, synchronous channel openings within a release site with any proportion of two-state to four-state channels are possible, however, the larger the proportion of two-state channels, the more sensitive the dynamics are to the exact spatial positioning of the channels and the distance between channels. Specifically, the clustering of even a small number of two-state channels interferes with puff/spark termination and increases puff durations or leads to a tonic response.

Steady-state kinetic modeling constrains cellular resting states and dynamic behavior

A defining characteristic of living cells is the ability to respond dynamically to external stimuli while maintaining homeostasis under resting conditions. Capturing both of these features in a single kinetic model is difficult because the model must be able to reproduce both behaviors using the same set of molecular components. Here, we show how combining small, well-defined steady-state networks provides an efficient means of constructing large-scale kinetic models that exhibit realistic resting and dynamic behaviors. By requiring each kinetic module to be homeostatic (at steady state under resting conditions), the method proceeds by (i) computing steady-state solutions to a system of ordinary differential equations for each module, (ii) applying principal component analysis to each set of solutions to capture the steady-state solution space of each module network, and (iii) combining optimal search directions from all modules to form a global steady-state space that is searched for accurate simulation of the time-dependent behavior of the whole system upon perturbation. Importantly, this stepwise approach retains the nonlinear rate expressions that govern each reaction in the system and enforces constraints on the range of allowable concentration states for the full-scale model. These constraints not only reduce the computational cost of fitting experimental time-series data but can also provide insight into limitations on system concentrations and architecture. To demonstrate application of the method, we show how small kinetic perturbations in a modular model of platelet P2Y₁ signaling can cause widespread compensatory effects on cellular resting states.

Cells respond to extracellular signals through a complex coordination of interacting molecular components. Computational models can serve as powerful tools for prediction and analysis of signaling systems, but constructing large models typically requires extensive experimental datasets and computation. To facilitate the construction of complex signaling models, we present a strategy in which the models are built in a stepwise fashion, beginning with small “resting” networks that are combined to form larger models with complex time-dependent behaviors. Interestingly, we found that only a minor fraction of potential model configurations were compatible with resting behavior in an example signaling system. These reduced sets of configurations were used to limit the search for more complicated solutions that also captured the dynamic behavior of the system. Using an example model constructed by this approach, we show how a cell's resting behavior adjusts to changes in the kinetic rate processes of the system. This strategy offers a general and biologically intuitive framework for building large-scale kinetic models of steady-state cellular systems and their dynamics.

Markov chain Monte Carlo fitting of single-channel data from inositol trisphosphate receptors

In many cell types, the inositol trisphosphate receptor (IPR) is one of the important components that control intracellular calcium dynamics, and an understanding of this receptor (which is also a calcium channel) is necessary for an understanding of calcium oscillations and waves. Recent advances in experimental techniques now allow for the measurement of single-channel activity of the IPR in conditions similar to its native environment, and these data can be used to determine the rate constants in Markov models of the IPR. We illustrate a parameter estimation method based on Markov chain Monte Carlo, which can be used to fit directly to single-channel data, and determining, as an intrinsic part of the fit, the times at which the IPR is opening and closing. We show, using simulated data, the most complex Markov model that can be unambiguously determined from steady-state data and show that non-steady-state data is required to determine more complex models.

Mathematical model of mouse embryonic cardiomyocyte excitation-contraction coupling

Excitation-contraction (E-C) coupling is the mechanism that connects the electrical excitation with cardiomyocyte contraction. Embryonic cardiomyocytes are not only capable of generating action potential (AP)-induced Ca(2+) signals and contractions (E-C coupling), but they also can induce spontaneous pacemaking activity. The spontaneous activity originates from spontaneous Ca(2+) releases from the sarcoplasmic reticulum (SR), which trigger APs via the Na(+)/Ca(2+) exchanger (NCX). In the AP-driven mode, an external stimulus triggers an AP and activates voltage-activated Ca(2+) intrusion to the cell. These complex and unique features of the embryonic cardiomyocyte pacemaking and E-C coupling have never been assessed with mathematical modeling. Here, we suggest a novel mathematical model explaining how both of these mechanisms can coexist in the same embryonic cardiomyocytes. In addition to experimentally characterized ion currents, the model includes novel heterogeneous cytosolic Ca(2+) dynamics and oscillatory SR Ca(2+) handling. The model reproduces faithfully the experimentally observed fundamental features of both E-C coupling and pacemaking. We further validate our model by simulating the effect of genetic modifications on the hyperpolarization-activated current, NCX, and the SR Ca(2+) buffer protein calreticulin. In these simulations, the model produces a similar functional alteration to that observed previously in the genetically engineered mice, and thus provides mechanistic explanations for the cardiac phenotypes of these animals. In general, this study presents the first model explaining the underlying cellular mechanism for the origin and the regulation of the heartbeat in early embryonic cardiomyocytes.

Modeling Ca2+ feedback on a single inositol 1,4,5-trisphosphate receptor and its modulation by Ca2+ buffers

The inositol 1,4,5-trisphosphate receptor/channel (IP(3)R) is a major regulator of intracellular Ca(2+) signaling, and liberates Ca(2+) ions from the endoplasmic reticulum in response to binding at cytosolic sites for both IP(3) and Ca(2+). Although the steady-state gating properties of the IP(3)R have been extensively studied and modeled under conditions of fixed [IP(3)] and [Ca(2+)], little is known about how Ca(2+) flux through a channel may modulate the gating of that same channel by feedback onto activating and inhibitory Ca(2+) binding sites. We thus simulated the dynamics of Ca(2+) self-feedback on monomeric and tetrameric IP(3)R models. A major conclusion is that self-activation depends crucially on stationary cytosolic Ca(2+) buffers that slow the collapse of the local [Ca(2+)] microdomain after closure. This promotes burst-like reopenings by the rebinding of Ca(2+) to the activating site; whereas inhibitory actions are substantially independent of stationary buffers but are strongly dependent on the location of the inhibitory Ca(2+) binding site on the IP(3)R in relation to the channel pore.

A molecular signaling model of platelet phosphoinositide and calcium regulation during homeostasis and P2Y1 activation

To quantify how various molecular mechanisms are integrated to maintain platelet homeostasis and allow responsiveness to adenosine diphosphate (ADP), we developed a computational model of the human platelet. Existing kinetic information for 77 reactions, 132 fixed kinetic rate constants, and 70 species was combined with electrochemical calculations, measurements of platelet ultrastructure, novel experimental results, and published single-cell data. The model accurately predicted: (1) steady-state resting concentrations for intracellular calcium, inositol 1,4,5-trisphosphate, diacylglycerol, phosphatidic acid, phosphatidylinositol, phosphatidylinositol phosphate, and phosphatidylinositol 4,5-bisphosphate; (2) transient increases in intracellular calcium, inositol 1,4,5-trisphosphate, and G(q)-GTP in response to ADP; and (3) the volume of the platelet dense tubular system. A more stringent test of the model involved stochastic simulation of individual platelets, which display an asynchronous calcium spiking behavior in response to ADP. Simulations accurately reproduced the broad frequency distribution of measured spiking events and demonstrated that asynchronous spiking was a consequence of stochastic fluctuations resulting from the small volume of the platelet. The model also provided insights into possible mechanisms of negative-feedback signaling, the relative potency of platelet agonists, and cell-to-cell variation across platelet populations. This integrative approach to platelet biology offers a novel and complementary strategy to traditional reductionist methods.