Updated stimulus current in several models

c/aguilar-2017.mmt

@@ -1,33 +1,36 @@
 [[model]]
 name: aguilar-2017
-version: 20240816
+version: 20240904
 mmt_authors: Michael Clerx
 desc: """
     Model of the human atrial action potential by Aguilar et al. [1], based on
     the model by Courtemanche et al. [2] and including the IKACh current
     described in [3].
-    
+
     This Myokit implementation is based on the CellML code for the Courtemanche
     model [4], but with IKur and IKAch equations from [1] and [2] added in. No
     source model was available to very against numerically.
-
+
+    The stimulus was set to 0.5 [ms] and approximately twice the threshold
+    value for depolarisation.
+
     [1] Aguilar, M., Feng, J., Vigmond, E. J., Comtois, P., & Nattel, S. (2017)
         Rate-Dependent Role of IKur in Human Atrial Repolarization and Atrial
         Fibrillation Maintenance. Biophysical Journal 112, 1997-2010.
         https://doi.org/10.1016/j.bpj.2017.03.022
-    
+
     [2] Courtemanche, M., Ramirez, R. J., & Nattel, S. (1998) Ionic mechanisms
         underlying human atrial action potential properties: insights from a
         mathematical model. American Journal of Physiology - Heart and
         Circulatory Physiology, 275, H301-H321.
         https://doi.org/10.1152/ajpheart.1998.275.1.H301
 
     [3] Kneller, J., Zou, R., Vigmond, E. J., Wang, Z., Leon, L. J., & Nattel,
-        S. (2002) Cholinergic Atrial Fibrillation in a Computer Model of a 
-        Two-Dimensional Sheet of Canine Atrial Cells With Realistic Ionic 
+        S. (2002) Cholinergic Atrial Fibrillation in a Computer Model of a
+        Two-Dimensional Sheet of Canine Atrial Cells With Realistic Ionic
         Properties. Circulation Research, 90:e73-e87.
-        https://doi.org/10.1161/01.RES.0000019783.88094.BA    
-    
+        https://doi.org/10.1161/01.RES.0000019783.88094.BA
+
     [4] https://models.cellml.org/exposure/0e03bbe01606be5811691f9d5de10b65
 
 """
@@ -66,7 +69,7 @@ dot(V) = -(I_ion + I_diff + stimulus.I_stim)
     in [mV]
     label membrane_potential
 I_ion = (
-        + ina.I_Na
+        + ina.INa
         + ik1.IK1
         + ito.Ito
         + ikur.IKur
@@ -113,7 +116,7 @@ pace = 0
 [phys]
 R = 8.3143 [J/mol/K]
     in [J/mol/K]
-    desc: 
+    desc:
 T = 310 [K]
     in [K]
     desc: Temperature
@@ -179,7 +182,7 @@ KQ10 = 3
 #
 [sodium]
 use phys.F, geom.V_i, geom.Cm
-dot(Nai) = (-3 * inak.INaK - (3 * inaca.INaCa + ib.IbNa + ina.I_Na)) * Cm / (V_i * F)
+dot(Nai) = (-3 * inak.INaK - (3 * inaca.INaCa + ib.IbNa + ina.INa)) * Cm / (V_i * F)
     in [mM]
     label Na_i
 
@@ -274,10 +277,10 @@ dot(j) = (inf - j) / tau
     tau = 1 / (alpha + beta)
         in [ms]
     inf = alpha / (alpha + beta)
-g_Na = 7.8 [nS/pF]
+gNa = 7.8 [nS/pF]
     in [nS/pF]
     desc: Maximal INa conductance
-I_Na = g_Na * m^3 * h * j * (V - nernst.ENa)
+INa = gNa * m^3 * h * j * (V - nernst.ENa)
     in [A/F]
     label I_Na
 
@@ -457,7 +460,7 @@ ICaL = gCaL * d * f * fCa * (V - ECaL)
 
 #
 # Acetylcholine-activated potassium current
-# 
+#
 # From [3], without changes.
 #
 [ikach]
@@ -585,7 +588,7 @@ dot(v) = (inf - v) / tau
         in [ms]
     inf = 1 - (1 + exp(-(Fn - 0.2 * c1) / c2)) ^ (-1)
 dot(w) = (inf - w) / tau
-    tau = 6 [ms*mV] * if(abs(V - 7.9 [mV]) < 1e-6 [mV], 
+    tau = 6 [ms*mV] * if(abs(V - 7.9 [mV]) < 1e-6 [mV],
             2 / 13 [mV],
             (1 - exp(-(V - 7.9 [mV]) / 5 [mV])) / ((1 + 0.3 * exp(-(V - 7.9 [mV]) / 5 [mV])) * (V - 7.9 [mV]))
         )

c/annotations.md

@@ -9,8 +9,10 @@
 - `stimulus_current`, a current injected with an electrode to pace the cell.
 - `membrane_capacitance`, the modelled cell's capacitance.
 
-- `ICaL`, the L-type calcium current (all species and compartments combined).
-- `gCaL`, a variable that can be scaled up and down to scale `ICaL`.
+- `I_CaL`, the L-type calcium current (all species and compartments combined).
+- `g_CaL`, a variable that can be scaled up and down to scale `I_CaL`.
+- `I_Kr`, the rapid delayed rectifier potassium current.
+- `g_Kr`, a variable that can be scaled up or down to scale `I_Kr`.
 
 ## Model meta properties
 

c/bai-2018.mmt

@@ -1,21 +1,24 @@
 [[model]]
 name: bai-2018
-version: 20240813
+version: 20240904
 mmt_authors: Michael Clerx
 display_name: Bai et al., 2018
 desc: """
-    The human left-atrial AP model by Bai et al. [1], based on the 
+    The human left-atrial AP model by Bai et al. [1], based on the
     ventricular model by Ten Tusscher and Panfilov [2].
 
     This implementation is based on the author-provided CellML [3], but
     includes the unit fixes described in [4]. The implementation was verified
     numerically by comparing the calculates derivatives to the values from the
     author-provided CellML. Differences from [2] are indicated throughout.
 
+    The stimulus was set to 0.5ms and approximately twice the threshold for
+    depolarisation.
+
     References:
 
-    [1] Bai, J. Gladding, P. A., Stiles, M. K., Fedorov, V. V. & Zhao, J. 
-        (2018) Ionic and cellular mechanisms underlying TBX5/PITX2 
+    [1] Bai, J. Gladding, P. A., Stiles, M. K., Fedorov, V. V. & Zhao, J.
+        (2018) Ionic and cellular mechanisms underlying TBX5/PITX2
         insufficiency-induced atrial fibrillation: Insights from mathematical
         models of human atrial cells. Scientific Reports 8:15642
         https://doi.org/10.1038/s41598-018-33958-y
@@ -29,7 +32,7 @@ desc: """
 
     [4] Clerx, M., Barral, Y., Agrawal, A. & Mirams, G. R. (2020)
         TNNP unit fixes
-        https://models.cellml.org/workspace/604/file/de5a4e600b57c8b9b0bd477695648c82351bc6dd/TNNP_unit_fixes.pdf   
+        https://models.cellml.org/workspace/604/file/de5a4e600b57c8b9b0bd477695648c82351bc6dd/TNNP_unit_fixes.pdf
 
     """
 # Initial values
@@ -95,7 +98,7 @@ i_ion = (+ ina.INa
 [stimulus]
 i_stim = engine.pace * amplitude
     in [A/F]
-amplitude = -80 [A/F]
+amplitude = -29 [A/F] * 2
     in [A/F]
 
 #
@@ -114,6 +117,7 @@ Vsr = 1094 [um^3]
 Cm = 185 [pF]
     in [pF]
     desc: Cell capacitance
+    label cell_capacitance
 
 #
 # Physical constants
@@ -521,7 +525,7 @@ IK_total = ik1.IK1 + ito.Ito + ikr.IKr + iks.IKs + ipk.IpK + ikur.IKur + stimulu
 
 [[protocol]]
 # Level  Start    Length   Period   Multiplier
-1        10       0.5      1000.0   0
+1        50       0.5      1000.0   0
 
 [[script]]
 import matplotlib.pyplot as plt

c/bartolucci-2020.mmt

@@ -18,6 +18,9 @@ desc: """
     Note that the CellML version uses an elevated external calcium of 2.7,
     while the paper and the matlab code use the default value of 1.8.
 
+    The stimulus was set to 0.5ms at approximately twice the threshold for
+    depolarisation.
+
     References:
 
     [1] Bartolucci C, Passini E, Hyttinen J, Paci M & Severi S (2020).
@@ -127,13 +130,10 @@ i_ion = (+ ina.INa
          + inal.INaL
          + ito.Ito
          + ical.ICaL
-         + ical.ICaNa
-         + ical.ICaK
          + ikr.IKr
          + iks.IKs
          + ik1.IK1
-         + inaca.INaCa
-         + inacass.INaCa_ss
+         + inacass.INaCa_total
          + inak.INaK
          + inab.INab
          + ikb.IKb
@@ -482,17 +482,19 @@ PhiNa =  V * FFRT * (0.75 * Na_ss * exp(V * FRT) - 0.75 * Na_o) / (exp(V * FRT)
 PhiK = V * FFRT * (0.75 * K_ss * exp(V * FRT) - 0.75 * K_o) / (exp(V * FRT) - 1)
     in [mC/L]
 # Current
-ICaLnp = PCa * PhiCa * (O + OCa)
+ICaLCa_np = PCa * PhiCa * (O + OCa)
+    in [A/F]
+ICaLNa_np = PNa * PhiNa * (O + OCa)
     in [A/F]
-ICaNanp = PNa * PhiNa * (O + OCa)
+ICaLK_np = PK * PhiK * (O + OCa)
     in [A/F]
-ICaKnp = PK * PhiK * (O + OCa)
+ICaLCa = (icalp.ICaLCa_p * camk.f + ICaLCa_np * (1 - camk.f))
     in [A/F]
-ICaK = (icalp.ICaKp * camk.f + ICaKnp * (1 - camk.f))
+ICaLNa = (icalp.ICaLNa_p * camk.f + ICaLNa_np * (1 - camk.f))
     in [A/F]
-ICaL = (icalp.ICaLp * camk.f + ICaLnp * (1 - camk.f))
+ICaLK = (icalp.ICaLK_p * camk.f + ICaLK_np * (1 - camk.f))
     in [A/F]
-ICaNa = (icalp.ICaNap * camk.f + ICaNanp * (1 - camk.f))
+ICaL = ICaLCa + ICaLNa + ICaLK
     in [A/F]
 
 #
@@ -536,11 +538,11 @@ dot(CCa) = beta * OCa + psi * I2Ca - (omega + alpha) * CCa + r_up * C - r_down *
 OCa = 1 - CCa - I1Ca - I2Ca - C - I1 - I2 - O
 # Current
 fp = 1.1
-ICaLp = fp * ical.PCa * ical.PhiCa * (O + OCa)
+ICaLCa_p = fp * ical.PCa * ical.PhiCa * (O + OCa)
     in [A/F]
-ICaNap = fp * ical.PNa * ical.PhiNa * (O + OCa)
+ICaLNa_p = fp * ical.PNa * ical.PhiNa * (O + OCa)
     in [A/F]
-ICaKp = fp * ical.PK * ical.PhiK * (O + OCa)
+ICaLK_p = fp * ical.PK * ical.PhiK * (O + OCa)
     in [A/F]
 
 #
@@ -822,6 +824,8 @@ JncxCa = E2 * k2 - E1 * k1
 INaCa_ss = 0.2 * inaca.fNaCa * inaca.gNaCa * allo * (JncxNa + 2 * JncxCa)
     in [A/F]
     desc: Sodium/Calcium exchange current into the T-Tubule subspace
+INaCa_total = inaca.INaCa + INaCa_ss
+    in [A/F]
 
 #
 # INaK: Sodium/potassium ATPase current
@@ -1074,7 +1078,7 @@ INa_tot = ina.INa + inal.INaL + inab.INab + 3 * inaca.INaCa + 3 * inak.INaK
 dot(Na_i) = -INa_tot * AFC / vmyo + diff.JdiffNa * vss / vmyo
     desc: Intracellular potassium concentration
     in [mM]
-INa_ss_tot = ical.ICaNa + 3 * inacass.INaCa_ss
+INa_ss_tot = ical.ICaLNa + 3 * inacass.INaCa_ss
     in [A/F]
 dot(Na_ss) = -INa_ss_tot * AFC / vss - diff.JdiffNa
     in [mM]
@@ -1094,7 +1098,7 @@ IK_tot = (
     - 2 * inak.INaK
 )
     in [A/F]
-IK_ss_tot = ical.ICaK
+IK_ss_tot = ical.ICaLK
     in [A/F]
 dot(K_i) = -(IK_tot + stimulus.i_stim) * AFC / vmyo + diff.JdiffK * vss / vmyo
     desc: Intracellular potassium concentration
@@ -1137,7 +1141,7 @@ dot(Ca_i) = buff * (-ICa_tot * AFC / (2 * vmyo) - serca.Jup * vsr / vmyo + serca
     in [mM]
     desc: Intracellular calcium concentration
     buff = 1 / (1 + cmdnmax * kmcmdn / (kmcmdn + Ca_i)^2 + trpnmax * kmtrpn / (kmtrpn + Ca_i)^2)
-ICa_ss_tot = ical.ICaL - 2 * inacass.INaCa_ss
+ICa_ss_tot = ical.ICaLCa - 2 * inacass.INaCa_ss
     in [A/F]
 dot(Ca_ss) = buff * (-ICa_ss_tot * AFC / (2 * vss) + ryr.Jrel * vsr / vss - diff.Jdiff)
     in [mM]

c/courtemanche-1998.mmt

@@ -1,35 +1,38 @@
 [[model]]
 name: courtemanche-1998
-version: 20240903
+version: 20240904
 mmt_authors: Michael Clerx
 desc: """
     Model of the human atrial action potential by Courtemanche et al. [1].
-    
+
     Adapted from the CellML implementation by Penny Noble [2], and checked
     against the published equations.
-    
+
     For some constants, the CellML implementation has more digits than the
     paper. In addition, the CellML model corrects a sign error in the equation
     for ikur.ui.tau.beta (as can be seen by comparing the graph for ikr.ui.tau
     with that in the paper). It seems like the CellML model may have been based
     on the original code, so this Myokit implementation follows the CellML
     file. Differences with the equations in the paper are highlighted with
     comments in the code below.
-    
+
     This implementation makes some cosmetic changes, and three larger ones:
-    
+
     1. Currents are all normalised, as in the paper.
     2. The stimulus was set to 0.5ms and approximately twice the threshold
        value for depolarisation at 1Hz pacing.
     3. To ensure that [K]i stabilised, the stimulus current was added to the
        equation for [K]i.
-
+
+    In this Myokit implementation, the stimulus was set to 0.5 [ms] and
+    approximately twice the threshold value for depolarisation.
+
     [1] Courtemanche, M., Ramirez, R. J., & Nattel, S. (1998) Ionic mechanisms
         underlying human atrial action potential properties: insights from a
         mathematical model. American Journal of Physiology - Heart and
         Circulatory Physiology, 275, H301-H321.
         https://doi.org/10.1152/ajpheart.1998.275.1.H301
-    
+
     [2] https://models.cellml.org/exposure/0e03bbe01606be5811691f9d5de10b65
 
 """
@@ -99,7 +102,7 @@ amplitude = 2 * -4618 [pA]
     desc: """
         Under baseline conditions, at 1Hz with a 0.5ms stimulus, -4618 pA is
         approximately the smallest value that triggers three successive APs.
-        
+
         The paper mentions a 1Hz stimulus of 2ms and 2000pA, used in figure 14:
         "Figure 14 shows a model AP generated during stimulation at 1,000 ms in
         response to 2-ms pulses of 2-nA amplitude (twice diastolic threshold)."
@@ -123,7 +126,7 @@ pace = 0
 [phys]
 R = 8.3143 [J/mol/K]
     in [J/mol/K]
-    desc: 
+    desc:
 T = 310 [K]
     in [K]
     desc: Temperature
@@ -577,7 +580,7 @@ dot(v) = (inf - v) / tau
         in [ms]
     inf = 1 - (1 + exp(-(Fn - 0.2 * c1) / c2)) ^ (-1)
 dot(w) = (inf - w) / tau
-    tau = 6 [ms*mV] * if(abs(V - 7.9 [mV]) < 1e-6 [mV], 
+    tau = 6 [ms*mV] * if(abs(V - 7.9 [mV]) < 1e-6 [mV],
             2 / 13 [mV],
             (1 - exp(-(V - 7.9 [mV]) / 5 [mV])) / ((1 + 0.3 * exp(-(V - 7.9 [mV]) / 5 [mV])) * (V - 7.9 [mV]))
         )
@@ -673,7 +676,7 @@ Ca_CSQN = CSQN_max * CaRel / (CaRel + Km_CSQN)
 
 [[protocol]]
 # Level  Start    Length   Period   Multiplier
-1.0      50       0.5      1000     0
+1        50       0.5      1000     0
 
 [[script]]
 import matplotlib.pyplot as plt