Three dimensional coupled reaction–diffusion modeling of calcium and inositol 1,4,5-trisphosphate dynamics in cardiomyocytes

Abstract

Nanoparticles have shown great promise in improving cancer treatment efficacy by changing the intracellular calcium level through activation of intracellular mechanisms. One of the mechanisms of the killing of the cancerous cell by a nanoparticle is through elevation of the intracellular calcium level. Evidence accumulated over the past decade indicates a pivotal role for the IP₃ receptor mediated Ca²⁺ release in the regulation of the cytosolic and the nuclear Ca²⁺ signals. There have been various studies done suggesting the role of IP₃ receptors (IP₃R) and IP₃ production and degradation in cardiomyocytes. In the present work, we have proposed a three-dimensional unsteady-state mathematical model to describe the mechanism of cardiomyocytes which focuses on evaluation of various parameters that affect these coupled dynamics and elevate the cytosolic calcium concentration which can be helpful to search for novel therapies to cure these malignancies by targeting the complex calcium signaling process in cardiomyocytes. Our study suggests that there are other factors involved in this signaling which can increase the calcium level, which can help in finding treatment for cancer. The cytosolic calcium level may be controlled by IP₃ signaling, leak, source influx of calcium (σ) and maximum production of IP₃ (V_P). We believe that the proposed model suggests new insight into finding treatment for cancer in cardiomyocytes through elevation of the cytosolic Ca²⁺ concentration by various parameters like leak, σ, V_P and especially by other complex cell signaling dynamics, namely IP₃ dynamics.

---

1 Introduction

The emerging field of nanomedicine aims to diagnose and treat various diseases with nanostructures. The design of such nanosystems for imaging and therapeutic applications requires a thorough understanding of the interactions between nanoparticles and biological systems. Cancer represents the 2nd largest cause of death in the United States, and one out of eight women is estimated to develop invasive breast cancer in her lifetime.¹ Nanoparticles have shown great promise in improving cancer treatment efficacy while reducing toxicity and treatment side effects.² Different nanoparticles can change the intracellular calcium level through activation of intracellular machinery. In cancer therapy, nanoparticles are used to kill the malignant cells.³ One of the mechanisms of the killing of the cancerous cell by nanoparticles is through elevation of the intracellular calcium level. One study reported that after 48 hours of nanoparticle treatment, the cytosolic calcium level increases 1.8 times in adenocarcinoma cells.⁴

Also, some attempts have been reported in the literature, for other different cells like hepatocytes,⁵ neurons,⁶ astrocytes,^7,8 fibroblasts,⁹ pancreatic acinar cells¹⁰ and oocytes,¹¹ to model the calcium diffusion and regulation in the presence of time-dependent influx for the smooth muscle cell to study the magnitude, and spatial and temporal characteristics of calcium transients in restricted diffusion spaces between the plasma membrane and intracellular organelles.¹² Gregory D. Smith et al.¹³ proposed a simple numerical model of calcium spark formation and detection in cardiomyocytes. According to their model, the elementary events of excitation–contraction coupling in the heart muscle are Ca²⁺ sparks, which arise from one or more ryanodine receptors in the endoplasmic reticulum (ER). Backx et al.¹⁴ constructed a model of propagating calcium-induced calcium release mediated by calcium diffusion in cardiomyocytes. The model was used to evaluate whether propagation of calcium transients and the range of propagation velocities observed experimentally (0.05–15 mm s⁻¹) could be predicted. Y. Hamam et al.¹⁵ proposed a model for calcium dynamics in cardiomyocytes. After the identification of the endoplasmic reticulum parameters, the SERCA pump and ryanodine channels activities, they made a comparison between experimental and calculated responses.

Consequently, it is of great interest to determine the mechanisms involved in these dynamics. In many cell types, these signaling processes depend on the inositol 1,4,5-trisphosphate (IP₃) signaling.^16–20 For both atrial and ventricular cells, a role for the cytosolic IP₃ has been reported for the modulation of the cytosolic Ca²⁺ transients in a variety of animal models.^21–25 IP₃R channel activity, with type-2 IP₃Rs as the most prevalent isoform in cardiomyocytes, depends on IP₃ concentration and Ca²⁺ concentration.^23,26,27 By the hydrolysis of phosphatidylinositol 4,5-bisphosphate (PIP₂) by phospholipase C (PLC), two distinct second messengers, diacylglycerol and IP₃, are produced that diffuse through the cytoplasm and bind to and open IP₃ receptors (IP₃Rs) on the endoplasmic reticulum (ER) membrane, leading to the release of Ca²⁺ from the ER. The released Ca²⁺ is then sent back into the ER by pumps like the SERCA pump in cardiomyocytes.²⁸ Recently, Li et al.²⁹ also uncaged IP₃ in cell cytoplasm from 2.1 to 5.3 nm AuNP-coated liposomes to activate cell signaling. They encapsulated IP₃ in the liposome complex which then binds with the IP₃R. The binding event then causes Ca²⁺ release from the endoplasmic reticulum (ER) into the cytosol. Also, El-Sayed et al.³⁰ demonstrated that AuNPs with integrin-targeting properties could inhibit cancer cell migration through affecting cytoskeletal proteins. Furthermore, the highest increase in the calcium signal was observed, and the calcium response was likely due to the inflow of Ca²⁺ as well as intracellular signaling via the IP₃ pathway.

For these reasons, although experimental methods^17,31,32 are being developed to measure [IP₃] directly, the question of the underlying oscillatory mechanism often remains unresolved. Many investigations on calcium dynamics in cardiomyocytes have not included the role of ryanodine receptors (RyR) and IP₃R.^33–36 whereas there are various investigations to understand the role of ryanodine receptors (RyR) in cardiomyocytes^37–39 but the importance of IP₃ dynamics in Ca²⁺ in the cardiomyocyte is less explored. But there are pieces of literature^40–42 that suggest the impact of IP₃ signaling on Ca²⁺ signaling in other cells. Also, there are few works showing the important role of IP₃ signaling in Ca²⁺ signaling in the cardiomyocyte.^26,28,43 Motivated by the above and the intriguing possibility of cancer treatment with nanoparticles by elevation of the cytosolic calcium level, we have explored another biological factor inside the cardiomyocyte, IP₃, which strongly influences the cytosolic calcium ions concentrations. Earlier we have developed one^44–46 and two dimensional models⁴⁷ for the coupled dynamics but to make the study more realistic and have better insights into various processes involved in calcium regulation it is necessary to develop a model in a three-dimensional cylindrical-shaped cardiomyocyte. In this paper, an attempt has been made in this direction to propose a three-dimensional model to study spatiotemporal calcium and IP₃ distribution in cardiomyocytes involving various parameters and their role in controlling calcium concentration in the cytosol for the better understanding and exploring of various treatments for cancer.

2 Methods

2.1 The mathematical model

Experimental investigations have been reported in the literature for interdependent Ca²⁺ and IP₃ dynamics in oocytes.⁴⁸ Also, few works have been reported for interdependent Ca²⁺ and IP₃ dynamics in cardiomyocytes. However, elementary Ca²⁺ release events and functional properties of IP₃ are available in the literature.^17,49,50

Apart from this, the reaction–diffusion equations for Ca²⁺ and IP₃ dynamics are available in the literature.⁴⁸ Thus, the functional properties of IP₃ release and the dynamics with Ca²⁺ in cardiomyocytes can be used as a basis to develop a coupled model of Ca²⁺ and IP₃ dynamics in a cell. This coupling will be based on the relationships between the processes of Ca²⁺ and IP₃ dynamics in a cell. In the present study the relationship between reaction–diffusion of Ca²⁺, as well as IP₃, has been explored to propose a coupled model for Ca²⁺ and IP₃ in cardiomyocytes (Fig. 1).

Fig. 1 Schematic representation of the Ca²⁺ and IP₃ dynamics model. In this model, various fluxes are involved such as Ca²⁺ is released from the ER via the IP₃ receptor (J_IPR), and reuptake occurs via SERCA ATPase pumps (J_SP). A leak is used to maintain basal Ca²⁺ levels (J_Le). IP₃ is produced at the ER membrane via Ca²⁺-dependent activation of PLC (J_P). The terms responsible for IP₃ degradation are these two fluxes J_Pp (Phosphatase) and J_K (Kinase).

The model incorporates three mechanisms for Ca²⁺ release and reuptake: the IP₃ receptor/Ca²⁺ channel, a small leak to maintain resting Ca²⁺ levels and SERCA pumps. The reaction–diffusion equations for Ca²⁺ and IP₃ are given by,⁴⁸

where [Ca²⁺]_C and [Ca²⁺]_E are the concentrations of Ca²⁺ in the cytosol and ER, respectively. D_C and D_E are the diffusion coefficients of Ca²⁺ in the cytosol and ER. Also, D_I is the diffusion coefficient of IP₃. J_P, J_K and J_Pp are fluxes of IP₃ production, Kinase and Phosphatase respectively in terms of total cell volume. All the flux terms are given by,⁴⁸

J_IPR = V_IPRd³e³([Ca²⁺]_E − [Ca²⁺]_C), (4)

J_Le = V_Le([Ca²⁺]_E − [Ca²⁺]_C), (5)

In this formulation, F_C and F_E are the volume fractions, relative to total cell volume, of the cytosol and ER, respectively (with F_C + F_E = 1), so the fluxes are also in terms of total cell volume. IP₃ receptor dynamics are treated using the Li–Rinzel simplification of the eight-state De Young–Keizer model. In both formulations, the IP₃ receptor is comprised of four identical subunits, each of which can bind an IP₃, an activating Ca²⁺, and an inactivating Ca²⁺, with the conducting state corresponding to any three of the subunits having IP₃ and activating Ca²⁺ bound.⁴⁸

and

The Ca²⁺ dependence of the 3-kinase degradation pathway is described by a Hill function,⁴⁸

The Li–Rinzel model, however, utilizes time scale arguments to eliminate two of the three binding reactions, replacing them with the equilibrium equation,⁴⁸

The variable e represents the fraction of subunits not yet inactivated by Ca²⁺ given by,⁴⁸

whereis the equilibrium value, and τ is the inactivation time scale. The analysis of the complete model is not possible using basic phase plane techniques. However, in a whole-cell model (i.e. where the diffusion terms are eliminated), the [Ca²⁺]_E equation can be eliminated using the conservation relation for the total cellular Ca²⁺ concentration,⁴⁸

[Ca²⁺]_T = F_C[Ca²⁺]_C + F_E[Ca²⁺]_E. (14)

The e nullcline is the surface e = e_∞, which after substitution is,⁴⁸

a relatively simple, monotonically decreasing function of Ca²⁺. This allows the model to be reduced to the three variables [Ca²⁺]_C, [IP₃], and e.⁴⁸ This model has been employed to perform a numerical study of Ca²⁺ and IP₃ dynamics in cardiomyocytes. The calcium diffusion and IP₃ diffusion in cardiomyocytes for a three-dimensional unsteady-state case can be rewritten as,

The height of the cardiomyocyte is 100 μm and width of the cytosol of the cardiomyocyte is 8 μm. Therefore, it is assumed that 0 ≤ r ≤ 4 μm, 0 ≤ θ ≤ 2π, 0 ≤ z ≤ 4, where z = 0 means 0 μm, z = 1 means 25 μm, z = 2 means 50 μm, z = 3 means 75 μm and z = 4 means 100 μm. The three-dimensional finite difference discretization of the cardiomyocyte is shown in Fig. 2.

Fig. 2 Discretization of the cytosol of the cell, where the yellow point at the top (at r = 4, θ = π, z = 2) and the green point at the bottom (at r = 4, θ = 0, z = 2) denote the point sources of Ca²⁺ and IP₃, respectively.

The initial and boundary conditions governing the Ca²⁺ and IP₃ diffusion processes are given by,^34,51,52

(i) Initial condition

[Ca²⁺]_Cr,θ,z,t=0 = 0.1 μM, (18)

[IP₃]_r,θ,z,t=0 = 0.16 μM. (19)

(ii) Boundary condition

Point sources of calcium and IP₃ are assumed at r = 4, θ = π, z = 2 and r = 4, θ = 0, z = 2 respectively. Thus, the boundary conditions at t > 0 are given below. Here, the source influx can be expressed as

here, an influx of free Ca²⁺ is taken at the rate σ by Faraday’s law,¹³where I is the amplitude of the element of Ca²⁺ release, F is Faraday’s constant and z is the valence of the Ca²⁺ ion.

At the rest of the boundary points there is no flux, thus the boundary conditions are as follows,

It is assumed that at the boundary far away from the source of calcium, the calcium concentration is maintained at background calcium concentration which can be written as follows,

It is assumed that at the boundary far away from the source of IP₃, the IP₃ concentration is maintained at background concentration which can be written as follows,

2.2 Solution

The model eqn (16)–(30) are solved numerically using the Forward Time Centered Space Approach. Using the Forward Time Centered Space Approach, the above eqn (16) and (17) take the following form,where ‘c’ denotes the calcium concentration in the cytosol, ‘c_E’ denotes calcium concentration in the ER and ‘d’ denotes the IP₃ concentration, which are a function of (r, θ, z, t). Here, ‘h’ represents the radial step, ‘l’ represents the angular step, ‘q’ represents the height step and ‘k’ represents the time step. The indexes of the radial, angular and height are ‘i’, ‘j’ and ‘m’ respectively and n represent the index of time. V_a is the maximum rate constant at low Ca²⁺ (3-kinase), V_b is the maximum rate constant at high Ca²⁺ (3-kinase) and V_Pp is the maximum rate constant (phosphatase). The initial and boundary conditions governing the Ca²⁺ and IP₃ diffusion processes are given by eqn (18)–(29), which are rewritten as given below.

(i) Initial condition

At time t = 0 s, we have following initial conditions,

c⁰_i,j,m = 0.1 μM, (33)

d⁰_i,j,m = 0.16 μM. (34)

(ii) Boundary condition

The above equation is not valid at (r = 0, 0 ≤ θ ≤ 2π, 0 ≤ z ≤ 4), near the sources (r = 4, θ = 0, π, z = 2) and far away from the sources (r = 4, θ = 0, π, 1 ≤ z ≤ 3) and (0 ≤ r ≤ 4, 0 ≤ θ ≤ 2π, z = 0, 4); therefore the approximation at these nodes at t > 0 is given by,

d_4,0,2ⁿ = 0.1882(kn)⁶ + 1.3121(kn)⁵ + 3.5391(kn)⁴ + 4.5312(kn)³ + 2.5893(kn)² + 0.3648(kn) + 0.1691 ≤ 3 (36)

c_4,0,mⁿ = 0.1 μM; at (1 ≤ z ≤ 3), (41)

c_i,j,mⁿ = 0.1 μM; at (0 ≤ r ≤ 4, 0 ≤ θ ≤ 2π, z = {0,4}), (42)

d_4,π,mⁿ = 0.16 μM at (1 ≤ z ≤ 3), (43)

d_i,j,mⁿ = 0.16 μM; at (0 ≤ r ≤ 4, 0 ≤ θ ≤ 2π, z = {0,4}), (44)

The resulting system (31)–(44) provides simultaneous algebraic equations for unknown node concentrations c_i,j,mⁿ and d_i,j,mⁿ. The Gaussian elimination method has been employed to solve the resulting equations.

3 Results and discussion

To start with, we analyzed the cardiomyocyte parameters for biological systems to confirm the accuracy of the computational analysis used in this study. Table 1 presents these properties in a cardiomyocyte. The parameters used for the computation of numerical results are given in Table 1.⁴⁸ To facilitate our discussion, all calculations utilize the standard model parameters in Table 1, except where noted otherwise.

Table 1 The experimental values of different parameters in the cardiac cardiomyocyte

Parameter	Value^48	Parameter	Value^48
VIPR	8.5 s^−1	Va	0.001 μM s^−1
VLe	0.01 s^−1	Vb	0.005 μM s^−1
VSP	0.65 μM s^−1	VPp	0.02 μM s^−1
KIP3	0.15 μM	VP	0.075 μM s^−1
KP	0.4 μM	FC	0.83
KAt	0.8 μM	λ	30
KSP	0.4 μM	DC	16 μm^2 s^−1
KIn	1.9 μM	DI	283 μm^2 s^−1

---

Firstly, we discuss the variations of IP₃ and Ca²⁺ in the cardiomyocyte as a function of angle, radius, height and time (spatio-temporal variation).

Oscillatory [Ca²⁺] dynamics

The variation of the calcium profile in the cardiomyocyte with respect to time is shown in Fig. 3. Here, the Ca²⁺ source is assumed to be situated at r = 4, θ = π and z = 2 (i.e. 50 μm). With the increase in time, the calcium concentration near the source also increases and reaches the maximum concentration at r = 4 μm. Also, the maximum calcium concentration was found at r = 4, θ = π, z = 2 because the source of Ca²⁺ is located there. Moreover, a small variation occurs but visible away from r = 4 and symmetric variations in both sides around θ = π. Also, the calcium concentration increases with increasing time and at time t = 0 it satisfied the boundary condition because there is no variation of calcium concentration at the source of calcium.

Fig. 3 Ca²⁺ concentration profile in the cardiomyocyte at different points of time (t = 0 s to 0.03 s with intervals of 0.01 s) in the cytosol at the z = 2 (i.e. 50 μm) position.

Initially, the calcium concentration is high at the source because the source channel opens and starts releasing Ca²⁺, the calcium concentration in the cytosol increases very fast and then it decreases due to the removal of calcium out of the cell by the pump as we go away from the source. Away from the source, the calcium concentration reaches its steady-state of 0.1 μM to maintain the Ca²⁺ concentration in the cytosol.

Oscillatory [IP₃] dynamics

Now, we discuss the variation of calcium concentration which depends on the spatiotemporal function. Fig. 4 shows the radial and angular distributions of Ca²⁺ concentration in the cytosol of the cardiomyocyte at different positions z and z is varying from 0 to 4 with intervals of 1 in the cell at fixed time t = 0.03 s. It can be observed from the graph that the Ca²⁺ concentration near the position of the source of Ca²⁺, i.e. at z = 2, is lower than that at z = 1 and z = 3 as the released Ca²⁺ from the source moves through the cytosol to bind with many proteins like calmodulin to activate various cellular responses like contraction, metabolism etc. Also, we can see that the calcium concentrations are very small at the end of the cytosol of the cardiomyocyte (z = 0 and z = 4) because it is away from the calcium source. Moreover, z = 2 means at the calcium source, and the Ca²⁺ concentration is high near the source, and away from the source, its concentration decreases very slowly as shown in Fig. 4 at z = 2. According to this, calcium signaling shows a wave-like nature at the source of calcium. Also, away from the source, i.e. at z = 0 and z = 4, IP₃ concentration again decreases gradually to attain its background concentration 0.1 μM.

Fig. 4 Ca²⁺ concentration profile in the cardiomyocyte at different positions z in the cytosol at t = 0.03 s.

Furthermore, we have systematically investigated the effect of the IP₃ source on each parameter of the spatiotemporal function and its variation is shown in Fig. 5. Fig. 5 shows the changes in IP₃ concentration at different time intervals. Here, the IP₃ source is assumed to be situated at r = 4, θ = 0 and z = 2 (i.e. 50 μm). With the increase in time, the IP₃ concentration near the source also increases. As the time increases, the concentration of IP₃ increases and shows an oscillatory type of behaviour in the radial as well as the angular directions. Also, the maximum variation is seen at time t = 0.03 s at the source of IP₃ and oscillations start from θ = π where maximum values are reached either side of the source (θ = 0, 2π). Initially, the IP₃ concentration decreases as we go away from the source then it increases to match the boundary concentration i.e. 3 μM.

Fig. 5 IP₃ concentration profile in the cardiomyocyte at different points of time in the cytosol at the z = 2 (i.e. 50 μm) position.

Moreover, the variation of IP₃ concentration which is the function of the radial as well as the angular part also depends on the height of the cylindrical cardiomyocyte. Fig. 6 shows the radial and angular distributions of the IP₃ concentration in the cytosol of the cardiomyocyte at different positions z in the cell at t = 0.03 s. It can be observed from the graph that the IP₃ concentration near the position of the source of IP₃, i.e. at z = 2, is lower than that at z = 1 and z = 3 as the released IP₃ from the source moves through the cytosol to bind with the calcium release channel IP₃R in order to allow the channel to open so that free Ca²⁺ ions can be released from the endoplasmic reticulum to balance the Ca²⁺ concentration in the cytosol and also to allow various cellular responses. Also, away from the source, i.e. at z = 0 and z = 4, the IP₃ concentration again decreases gradually to attain its background concentration of 0.16 μM.

Fig. 6 IP₃ concentration profile in the cardiomyocyte at different positions z in the cytosol at t = 0.03 s.

Influence of [IP₃] in the [Ca²⁺] dynamics

Now, we discuss the ratio of Ca²⁺ and IP₃ concentrations as a function of time at fixed height z = 2 and different radial and angular positions as shown in Fig. 7. As time increases at a different position in the cytosol, the ratio of c = [Ca²⁺] and d = [IP₃] decreases very fast, after that it starts increasing and after a small decrease, it increases again and attains a steady-state. The maximum variation is seen at (r = 4, θ = 0, z = 2) which is the location of the source of IP₃. Physiologically this behaviour of interdependent Ca²⁺ and IP₃ is justified. At the center of the cytosol, i.e. at z = 2, maximum interdependent dynamics and non-linearity in these signaling patterns are observed which is shown in Fig. 7. This is due to the positioning of the two source channels of Ca²⁺ and IP₃. The source of Ca²⁺ is located at (r = 4, θ = π, z = 2) and the source of IP₃ is located at (r = 4, θ = 0, z = 2) which are just opposite to each other. Initially, d = [IP₃] is higher as compared to c = [Ca²⁺]_C due to the initial conditions of both the signals. At time t = 0 s, initial conditions are c⁰_i,j,m = 0.1 μM and d⁰_i,j,m = 0.16 μM. Also, the equation for the source of IP₃ increases rapidly with time as compared to the source of Ca²⁺. Therefore, a rapid decrease in the ratio of c = [Ca²⁺]_C and d = [IP₃] is seen in the first few milliseconds. But after that, it is observed that there is a sharp peak in the ratio c/d, after this point of time the ratio c/d starts increasing. At this point of time, the IP₃ binds with the IP₃ receptor (IP₃R), opens the IP₃R channel and facilitates the outflow of Ca²⁺ ions from the ER to the cytosol. As a result, the amount of Ca²⁺ ions increases in the cytosol. Therefore, the ratio c/d increases. This increase in the cytosolic Ca²⁺ concentration will help in finding treatment for cancer.⁴ These released free Ca²⁺ ions from the ER are not all breaking PIP2 to form IP₃ ions, there are many proteins available in the cytosol which bind with the released free Ca²⁺ ion to stimulate various cell responses. Also, the SERCA pump is present which pumps out the extra free Ca²⁺ ions from the cytosol to the ER to maintain a balanced amount of free Ca²⁺ ions inside the cytosol. This binding of free Ca²⁺ ions with proteins and pumping out process lead to a decrease in the amount of free cytosolic Ca²⁺ ions as compared to the cytosolic IP₃. This is the reason behind the decrease in ratio c/d. Near 0.6 s again the ratio c/d increases and attains an equilibrium state. Physiologically, when the amount of free Ca²⁺ ions in the cytosol decreases then the leak (Le) located at the ER membrane starts releasing a small amount of Ca²⁺ ions into the cytosol to maintain the balance and initiate the whole process of intracellular cell signaling cell signaling. It is observed that after 0.6 s, at all the positions equilibrium between c = [Ca²⁺]_C and d = [IP₃] is attained.

Fig. 7 Ratio of c = [Ca²⁺]_C and d = [IP₃] with respect to time in the cytosol of the cardiomyocyte at different positions (r = 4 μm, θ = π, z = 2), (r = 2 μm, θ = π, z = 2), (r = 0 μm, θ = 0, z = 2), (r = 2 μm, θ = 0, z = 2) and (r = 4 μm, θ = 0, z = 2).

Effect of leak on Ca²⁺ dynamics

Fig. 8 shows the behaviour of the ratio of concentrations as a function of time, and we have as well included the effects of leak. The main function of leak helps to facilitate the flow of Ca²⁺ ions in the cytosol from the ER which increases the Ca²⁺ concentration in the cytosol and initiates the complete Ca²⁺ dynamic processes. The ratio concentration of Ca²⁺ and IP₃ initially decreases with time and then immediately the ratio concentration starts increasing at a certain time, showing an oscillatory pattern and then it goes to a steady-state with the presence of leak. In the case of absent leak, the behaviour is the same but there is a decrease in ratio concentration which is shown in Fig. 8(a)–(e).

Fig. 8 The effect of leak (Le) on the ratio of c = [Ca²⁺]_C and d = [IP₃] with respect to time at various positions i.e., radius (r) and angle (θ) in the cardiomyocyte. (a) (r = 4 μm, θ = π, z = 2), (b) (r = 2 μm, θ = π, z = 2), (c) (r = 0 μm, θ = 0, z = 2), (d) (r = 2 μm, θ = 0, z = 2) and (e) (r = 4 μm, θ = 0, z = 2) in the cytosol of the cardiomyocyte.

Also, we have simulated the other compositions of the radial and angular parts of a cardiomyocyte inside the cytosol which is the function of ratio concentration and time. We found that the behaviour is similar to Fig. 8(a)–(e) but at z = 2 the graph shows a better oscillatory type nature due to the positioning of both the sources is at z = 2. This shows that the cytosolic Ca²⁺ concentration can be increased by increasing the effect of leak in the cardiomyocyte.

Effect of σ on Ca²⁺ dynamics

Fig. 9 shows the behaviour of the ratio concentration as a function of time at different amounts of the Ca²⁺ source influx. The same behaviour is found as in Fig. 8 but if we increase the amount of Ca²⁺ source influx, the ratio concentrations of Ca²⁺ and IP₃ continuously increase. But in the case of r = 4 μm and θ = 0, there is no change in the ratio as at this position which is far away from the source of calcium it is assumed in the boundary condition that c = [Ca²⁺]_C is fixed i.e. 0.1 μM. Therefore, the effect of an increase in Ca²⁺ source influx (σ) is negligible here. This analysis also suggests that the cytosolic Ca²⁺ concentration can be increased by increasing the source influx of calcium (σ) in the cardiomyocyte which will help in finding various treatments for cancer.

Fig. 9 The effect of the source influx of Ca²⁺ (σ) on the ratio of c = [Ca²⁺]_C and d = [IP₃] with respect to time at various positions i.e., radius (r) and angle (θ) in the cardiomyocyte. (a) θ = π and r = 4 μm, (b) θ = π and r = 2 μm, (c) θ = 0 and r = 0 μm, (d) θ = 0 and r = 2 μm, (e) θ = 0 and r = 4 μm in the cytosol of the cardiomyocyte.

Effect of V_P on Ca²⁺ dynamics

Lastly, we have simulated the ratio concentration with respect to time at various values of the maximum rate of production of IP₃ as depicted in Fig. 10. Here, we can see that the variation in ratio concentration initially decreases and then it increases very rapidly showing an oscillatory behaviour after that it goes to a steady-state. But in the case of r = 4 μm and θ = 0, the ratio concentration is the same at different V_P because it is a boundary point which is away from the source so it is not affected by the maximum rate production of IP₃. From Fig. 9 and 10 we can easily see that the Ca²⁺ and IP₃ concentrations affect each other non-linearly. But Fig. 9 and 10 suggest that the maximum rate of production of IP₃ (V_P) can affect more as compared to the Ca²⁺ source influx (σ) because if we are increasing the maximum rate of production of IP₃, the ratio peaks are shifted towards higher values. This means that the Ca²⁺ concentration in the cytosol is much more affected by IP₃ parameters which indicates that experimental researchers can explore various ways to kill the cancerous cell by increasing the cytosolic calcium concentration with the help of controlling IP₃ dynamics in cardiomyocytes as IP₃ signaling plays a very important role in regulation of the calcium concentration in cardiomyocytes.

Fig. 10 The effect of the maximum IP₃ production rate (V_P) on the ratio of c = [Ca²⁺]_C and d = [IP₃] with respect to time at various positions i.e., radius (r) and angle (θ) in the cardiomyocyte. (a) θ = π and r = 4 μm, (b) θ = π and r = 2 μm, (c) θ = 0 and r = 0 μm, (d) θ = 0 and r = 2 μm, (e) θ = 0 and r = 4 μm in the cytosol of the cardiomyocyte.

[Ca²⁺] dynamics in cancer therapy

Dynamics of free intracellular calcium concentration ([Ca²⁺]) are a crucial control mechanism in many cell types.⁵³ The temporal and spatial information encoded by these dynamic signals controls many processes, including contraction of the cardiac tissues, information processing in the brain, gene expression, release of digestive enzymes by the liver, differentiation and cell movement.⁵⁴ Imbalanced calcium signaling has been associated with major pathologies, such as cardiovascular disease, stroke and cancer.^55,56 The results found here are in good agreement with the experimental as well as theoretical work^26,28,43,57 and suggest new insights into finding treatment for cancer in cardiomyocytes through elevation of the cytosolic Ca²⁺ concentration by various parameters like leak, σ, V_P and especially by other complex cell signaling dynamics, namely the IP₃ dynamics.

4 Conclusion

As suggested by Chaochu Cui et al.,⁵⁸ alterations in intracellular Ca²⁺ signaling for cancer therapy has become an emerging research area. Thus, keeping in mind this ideology, the physiological processes like reaction–diffusion of calcium and IP₃ are taken into consideration. One of the mechanisms of the killing of the cancerous cell by nanoparticles is through elevation of the intracellular calcium level. Our study suggests that there are other factors involved in this signaling which can increase the calcium level, which can help in finding a cure for cancer. Thus, a possible treatment for cancer may be studied and the preventive step of increasing the calcium concentration can be studied mathematically. In this piece of research, the effect of IP₃ dynamics in the presence of leak, source influx of calcium (σ) and maximum production of IP₃ (V_P) on the cytosolic calcium dynamics is presented. Several graphs depicting the physiological phenomenon of calcium and IP₃ dynamics are plotted here to give a better insight into the interdependence of these signaling processes. In all the ratio graphs (see Fig. 7–10), the nonlinear relationship between these two signaling processes is visible. Here, leak, source influx of calcium (σ) and maximum production of IP₃ (V_P) are taken into consideration and are found to play important roles in increasing the amount of cytosolic calcium concentration with the help of IP₃ signaling that would be helpful in finding treatment for cancer. Hence, the impact of IP₃ signaling in calcium diffusion is found to be significant. The cytosolic calcium level may be controlled by IP₃ signaling, leak, source influx of calcium (σ) and maximum production of IP₃ (V_P). The present work focuses on the evaluation of the various parameters that affect these coupled dynamics and elevate the cytosolic calcium concentration which can be helpful to search for novel therapies to treat these malignancies by targeting complex calcium signaling processes in cardiomyocytes.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

The authors are greatly thankful to the Department of Biotechnology, New Delhi, India for providing the bioinformatics infrastructure facility for carrying out this work.

References

1. R. L. Siegel, K. D. Miller and A. Jemal, Ca-Cancer J. Clin., 2015, 65, 5–29.
2. T. A. Brocato, E. N. Coker, P. N. Durfee, Y.-S. Lin, J. Townson, E. F. Wyckoff, V. Cristini, C. J. Brinker and Z. Wang, Sci. Rep., 2018, 8, 7538.
3. A. Onodera, K. Yayama, H. Morosawa, Y. Ishii, Y. Tsutsumi and Y. Kawai, Biochem. Biophys. Rep., 2017, 9, 330–334.
4. J. Wang, X. Fang and W. Liang, ACS Nano, 2012, 6, 5018–5030.
5. T. Höfer, Biophys. J., 1999, 77, 1244–1256.
6. J. W. Shuai, S. Nadkarni, P. Jung, A. Cornell-Bell and V. Trinkaus-Randall, Adv. Mol. Cell Biol., 2003, 31, 689–706.
7. B. K. Jha, N. Adlakha and M. N. Mehta, International Journal of Modeling, Simulation, and Scientific Computing, 2013, 4, 1250030.
8. S. Tewari and K. R. Pardasani, IAENG Int. J. Appl. Math., 2010, 40, 108–112.
9. M. Kotwani and N. Adlakha, Int. J. Comput. Mater. Sci. Eng., 2017, 6, 1750004.
10. N. Manhas and K. R. Pardasani, J. Bioenerg. Biomembr., 2014, 46, 403–420.
11. P. A. Naik and K. R. Pardasani, Alexandria J. Med., 2016, 52, 43–49.
12. N. Singh and N. Adlakha, Mathematical Modelling and Scientific Computing with Applications, Springer Nature Singapore Pte Ltd, 2019, ch. 6.
13. G. D. Smith, J. E. Keizer, M. D. Stern, W. J. Lederer and H. Cheng, Biophys. J., 1998, 75, 15–32.
14. P. H. Backx, P. P. De Tombe, J. H. Van Deen, B. J. Mulder and H. E. Ter Keurs, J. Gen. Physiol., 1989, 93, 963–977.
15. Y. Hamam, F. Pecker, F. Rocaries, H. Lorino, C. Pavoine and R. Natowicz, Simul. Pract. Theory, 2000, 8, 3–15.
16. M. J. Berridge, Biochim. Biophys. Acta, Mol. Cell Res., 2009, 1793, 933–940.
17. A. P. Dawson, Curr. Biol., 1997, 7, R544–R547.
18. T. J. Hund, A. P. Ziman, W. J. Lederer and P. J. Mohler, J. Mol. Cell. Cardiol., 2008, 45, 159–161.
19. G. Ullah, P. Jung and K. Machaca, Cell Calcium, 2007, 42, 556–564.
20. A. Tanimura, T. Morita, A. Nezu and Y. Tojyo, J. Med. Investig., 2009, 56, 357–361.
21. A. V. Zima and L. A. Blatter, J. Physiol., 2004, 555, 607–615.
22. A. Proven, J. Cell Sci., 2006, 119, 3363–3375.
23. T. L. Domeier, A. V. Zima, J. T. Maxwell, S. Huke, G. A. Mignery and L. A. Blatter, Am. J. Physiol. Heart Circ. Physiol., 2008, 294, H596–H604.
24. D. Harzheim, M. Movassagh, R. S.-Y. Foo, O. Ritter, A. Tashfeen, S. J. Conway, M. D. Bootman and H. L. Roderick, Proc. Natl. Acad. Sci. U. S. A., 2009, 106, 11406–11411.
25. J.-C. Kim, M.-J. Son, K. P. Subedi, Y. Li, J. R. Ahn and S.-H. Woo, Prog. Biophys. Mol. Biol., 2010, 103, 59–70.
26. J. Kockskämper, A. V. Zima, H. L. Roderick, B. Pieske, L. A. Blatter and M. D. Bootman, J. Mol. Cell. Cardiol., 2008, 45, 128–147.
27. R. H. Michell, C. K. Kirk, L. M. Jones, C. P. Downes and J. A. Creba, Philos. Trans. R. Soc., B, 1981, 296, 123–138.
28. F. Hohendanner, A. D. McCulloch, L. A. Blatter and A. P. Michailova, Front. Pharmacol., 2014, 5, 1–15.
29. X. Li, Z. Che, K. Mazhar, T. J. Price and Z. Qin, Adv. Funct. Mater., 2017, 27, 1605778.
30. M. R. K. Ali, Y. Wu, Y. Tang, H. Xiao, K. Chen, T. Han, N. Fang, R. Wu and M. A. El-Sayed, Proc. Natl. Acad. Sci. U. S. A., 2017, 114, E5655–E5663.
31. B. Ciapa, D. Pesando, M. Wilding and M. Whitaker, Nature, 1994, 368, 875–878.
32. G. Dupont and C. Erneux, Cell Calcium, 1997, 22, 321–331.
33. A. Hatano, J.-i. Okada, T. Washio, T. Hisada and S. Sugiura, Biophys. J., 2011, 101, 2601–2610.
34. K. B. Pathak and N. Adlakha, J. Med. Imaging Health Inform., 2015, 5, 683–688.
35. K. Pathak and N. Adlakha, Adv. Comput. Sci. Technol., 2017, 10, 11–23.
36. K. Pathak and N. Adlakha, Alexandria J. Med., 2016, 52, 261–268.
37. D. A. Eisner, J. L. Caldwell, K. Kistamás and A. W. Trafford, Circ. Res., 2017, 121, 181–195.
38. A. Michailova, F. DelPrincipe, M. Egger and E. Niggli, Biophys. J., 2002, 83, 3134–3151.
39. J. S. Sham, L. Cleemann and M. Morad, Proc. Natl. Acad. Sci. U. S. A., 1995, 92, 121–125.
40. C. E. Adkins and C. W. Taylor, Curr. Biol., 1999, 9, 1115–1118.
41. G. W. De Young and J. Keizer, Proc. Natl. Acad. Sci. U. S. A., 1992, 89, 9895–9899.
42. B. D. Stewart, C. E. Scott, T. P. McCoy, G. Yin, F. Despa, S. Despa and P. M. Kekenes-Huskey, Cell Calcium, 2018, 71, 65–74.
43. P. Lipp, M. Laine, S. C. Tovey, K. M. Burrell, M. J. Berridge, W. Li and M. D. Bootman, Curr. Biol., 2000, 10, 939–942.
44. N. Singh and N. Adlakha, Netw. Model. Anal. Health Inform. Bioinform., 2019, 8, 18.
45. N. Singh and N. Adlakha, J. Anal. Comput., 2019, 1–4.
46. N. Singh and N. Adlakha, Communications in Mathematical Biology and Neuroscience, 2019, 1–17.
47. N. Singh and N. Adlakha, Math. Biol. Bioinf., 2019, 14, 290–305.
48. J. Wagner, C. P. Fall, F. Hong, C. E. Sims, N. L. Allbritton, R. A. Fontanilla, I. I. Moraru, L. M. Loew and R. Nuccitelli, Cell Calcium, 2004, 35, 433–447.
49. G. Handy, M. Taheri, J. A. White and A. Borisyuk, J. Comput. Neurosci., 2017, 42, 257–273.
50. A. Harootunian, J. Kao, S. Paranjape and R. Tsien, Science, 1991, 251, 75–78.
51. S. A. Brown, F. Morgan, J. Watras and L. M. Loew, Biophys. J., 2008, 95, 1795–1812.
52. C. C. Fink, B. Slepchenko, I. I. Moraru, J. Watras, J. C. Schaff and L. M. Loew, Biophys. J., 2000, 79, 163–183.
53. N. L. Allbritton, T. Meyer and L. Stryer, Science, 1992, 258, 1812–1814.
54. M. J. Berridge, J. Physiol., 1997, 499, 291–306.
55. E. Przybytkowski, M. Behrendt, D. Dubois and D. Maysinger, FEBS J., 2009, 276, 6204–6217.
56. N. Rizaner, R. Onkal, S. P. Fraser, A. Pristerá, K. Okuse and M. B. Djamgoz, Eur. Biophys. J., 2016, 45, 735–748.
57. M. Cooling, P. Hunter and E. J. Crampin, Biophys. J., 2007, 93, 3421–3433.
58. C. Cui, R. Merritt, L. Fu and Z. Pan, Acta Pharm. Sin. B, 2017, 7, 3–17.

---