Simplified notation in ICaL components of ord models

c/bartolucci-2020.mmt

@@ -231,9 +231,9 @@ ENa = phys.RTF * log(extra.Na_o / sodium.Na_i)
 EK = phys.RTF * log(extra.K_o / potassium.K_i)
     in [mV]
     desc: Reversal potential for Potassium currents
-PNaK = 0.01833
+PKNa = 0.01833
     desc: Permeability ratio K+ to Na+
-EKs = phys.RTF * log((extra.K_o + PNaK * extra.Na_o) / (potassium.K_i + PNaK * sodium.Na_i))
+EKs = phys.RTF * log((extra.K_o + PKNa * extra.Na_o) / (potassium.K_i + PKNa * sodium.Na_i))
     in [mV]
     desc: Reversal potential for IKs
 EKshift = 8 [mV]

c/ohara-2011.mmt

@@ -1,6 +1,6 @@
 [[model]]
 name: ohara-2011
-version: 20240814
+version: 20240815
 mmt_authors: Michael Clerx
 display_name: O'Hara et al., 2011
 desc: """
@@ -77,13 +77,13 @@ ito.is            = 1
 ito.ap            = 0
 ito.ifp           = 1
 ito.isp           = 1
-ical.d            = 0
+ical.d            = 0.2
 ical.ff           = 1
 ical.fs           = 1
 ical.fcaf         = 1
 ical.fcas         = 1
 ical.jca          = 1
-ical.nca          = 0
+ical.nca          = 0.2
 ical.ffp          = 1
 ical.fcafp        = 1
 ikr.xf            = 0
@@ -139,7 +139,7 @@ amplitude = -80 [A/F]
 #
 [cell]
 mode = 1
-    desc: The type of cell. Endo = 0, Epi = 1, M = 2
+    desc: The type of cell. Endo = 0, Epi = 1, Mid = 2
 L = 0.01 [cm] : Cell length
     in [cm]
 r = 0.0011 [cm] : Cell radius
@@ -205,13 +205,13 @@ K_o = 5.4 [mM] : Extracellular K+ concentration
 [rev]
 ENa = phys.RTF * log(extra.Na_o / sodium.Na_i)
     in [mV]
-    desc: Reversal potential for Sodium currents
+    desc: Reversal potential for sodium currents
 EK = phys.RTF * log(extra.K_o / potassium.K_i)
     in [mV]
-    desc: Reversal potential for Potassium currents
-PNaK = 0.01833
+    desc: Reversal potential for potassium currents
+PKNa = 0.01833
     desc: Permeability ratio K+ to Na+
-EKs = phys.RTF * log((extra.K_o + PNaK * extra.Na_o) / (potassium.K_i + PNaK * sodium.Na_i))
+EKs = phys.RTF * log((extra.K_o + PKNa * extra.Na_o) / (potassium.K_i + PKNa * sodium.Na_i))
     in [mV]
     desc: Reversal potential for IKs
 
@@ -392,128 +392,97 @@ Ito = fto * gto * (V - rev.EK) * ((1 - camk.f) * a * i + camk.f * ap * ip)
 use membrane.V
 use extra.Ca_o, extra.K_o, extra.Na_o
 use calcium.Ca_ss, potassium.K_ss, sodium.Na_ss
-vf = V * phys.FRT
-vff = V * phys.FFRT
-    in [C/mol]
-evf = exp(vf)
-evf2 = exp(2 * vf)
 # Activation
-sd = 1 / (1 + exp((V + 3.94 [mV]) / -4.23 [mV]))
-    desc: Steady-state value for activation gate of ICaL
-td = 0.6 [ms] + 1 [ms] / (exp(-0.05 [1/mV] * (V + 6 [mV])) + exp(0.09 [1/mV] * (V + 14 [mV])))
-    desc: Time constant for activation gate of ICaL
-    in [ms]
-dot(d) = (sd - d) / td
-    desc: Activation gate of ICaL channel
-# Inactivation
-sf = 1 / (1 + exp((V + 19.58 [mV]) / 3.696 [mV]))
-    desc: Steady-state value for inactivation gate of ICaL
-tff = 7 [ms] + 1 [ms] / (0.0045 * exp((V + 20 [mV]) / -10 [mV]) + 0.0045 * exp((V + 20 [mV]) / 10 [mV]))
-    desc: Time constant for fast inactivation of ICaL
+dot(d) = (inf - d) / tau
+    inf = 1 / (1 + exp((V + 3.94 [mV]) / -4.23 [mV]))
+    tau = 0.6 [ms] + 1 [ms] / (exp(-0.05 [1/mV] * (V + 6 [mV])) + exp(0.09 [1/mV] * (V + 14 [mV])))
+        in [ms]
+# Voltage-dependent inactivation (fast and slow)
+f_inf = 1 / (1 + exp((V + 19.58 [mV]) / 3.696 [mV]))
+ff_tau = 7 [ms] + 1 [ms] / (0.0045 * exp((V + 20 [mV]) / -10 [mV]) + 0.0045 * exp((V + 20 [mV]) / 10 [mV]))
     in [ms]
-tfs = 1000 [ms] + 1 [ms] / (3.5e-5 * exp((V + 5 [mV]) / -4 [mV]) + 3.5e-5 * exp((V + 5 [mV]) / 6 [mV]))
-    desc: Time constant for fast inactivation of ICaL
+fs_tau = 1000 [ms] + 1 [ms] / (3.5e-5 * exp((V + 5 [mV]) / -4 [mV]) + 3.5e-5 * exp((V + 5 [mV]) / 6 [mV]))
     in [ms]
-dot(ff) = (sf - ff) / tff
-    desc: Fast inactivation of ICaL
-dot(fs) = (sf - fs) / tfs
-    desc: Slow inactivation of ICaL
-Aff = 0.6 : Fraction of ICaL channels with fast inactivation
-Afs = 1 - Aff : Fraction of ICaL channels with slow inactivation
-f = Aff * ff + Afs * fs
-    desc: Inactivation of ICaL
-# Ca-dependent inactivation
-sfca = sf
-    desc: Steady-state value for Ca-dependent inactivation of ICaL
-tfcaf = 7 [ms] + 1 [ms] / (0.04 * exp((V - 4 [mV]) / -7 [mV]) + 0.04 * exp((V - 4 [mV]) / 7 [mV]))
-    desc: Time constant for fast Ca-dependent inactivation of ICaL
+dot(ff) = (f_inf - ff) / ff_tau
+dot(fs) = (f_inf - fs) / fs_tau
+Aff = 0.6 : Fraction of channels with fast inactivation
+Afs = 1 - Aff : Fraction of channels with slow inactivation
+f = Aff * ff + Afs * fs : Total voltage-dependent inactivation
+# Calcium-dependent inactivation of non-phosphorylated channels
+fcaf_tau = 7 [ms] + 1 [ms] / (0.04 * exp((V - 4 [mV]) / -7 [mV]) + 0.04 * exp((V - 4 [mV]) / 7 [mV]))
     in [ms]
-tfcas = 100 [ms] + 1 [ms] / (0.00012 * exp(V / -3 [mV]) + 0.00012 * exp(V / 7 [mV]))
-    desc: Time constant for slow Ca-dependent inactivation of ICaL
+fcas_tau = 100 [ms] + 1 [ms] / (0.00012 * exp(V / -3 [mV]) + 0.00012 * exp(V / 7 [mV]))
     in [ms]
-dot(fcaf) = (sfca - fcaf) / tfcaf
-    desc: Fast Ca-dependent inactivation of ICaL
-dot(fcas) = (sfca - fcas) / tfcas
-    desc: Slow Ca-dependent inactivation of ICaL
+dot(fcaf) = (f_inf - fcaf) / fcaf_tau
+dot(fcas) = (f_inf - fcas) / fcas_tau
 Afcaf = 0.3 + 0.6 / (1 + exp((V - 10 [mV]) / 10 [mV]))
     desc: Fraction of ICaL channels with fast Ca-dependent inactivation
 Afcas = 1 - Afcaf
     desc: Fraction of ICaL channels with slow Ca-dependent inactivation
 fca = Afcaf * fcaf + Afcas * fcas
-    desc: Ca-dependent inactivation of ICaL
+    desc: Ca-dependent inactivation of non-phosphorylated ICaL
 # Recovery from Ca-dependent inactivation
-tjca = 75 [ms] : Time constant of recovery from Ca-dependent inactivation
+jca_tau = 75 [ms]
     in [ms]
-dot(jca) = (sfca - jca) / tjca
+dot(jca) = (f_inf - jca) / jca_tau
     desc: Recovery from Ca-dependent inactivation
-# Inactivation of phosphorylated ICaL
-tffp = 2.5 * tff
-    in [ms]
-    desc: Time constant for fast inactivation of phosphorylated ICaL
-dot(ffp) = (sf - ffp) / tffp
-    desc: Fast inactivation of phosphorylated ICaL
+# Voltage-dependent inactivation of phosphorylated channels
+dot(ffp) = (f_inf - ffp) / (2.5 * ff_tau)
 fp = Aff * ffp + Afs * fs
-    desc: Inactivation of phosphorylated ICaL
-# Ca-dependent inactivation of phosphorylated ICaL
-tfcafp = 2.5 * tfcaf
-    in [ms]
-    desc: Time constant for fast Ca-dependent inactivation of phosphorylated ICaL
-dot(fcafp) = (sfca - fcafp) / tfcafp
-    desc: Fast Ca-dependent inactivation of phosphorylated ICaL
+    desc: Voltage-dependent inactivation of phosphorylated ICaL
+# Calcium-dependent inactivation of phosphorylated channels
+dot(fcafp) = (f_inf - fcafp) / (2.5 * fcaf_tau)
 fcap = Afcaf * fcafp + Afcas * fcas
-    desc: Ca-dependent inactivation of phosphorylated ICaL
+    desc: Ca-dependent inactivation of phosphorylated channels
 # Fraction of channels in Ca-dependent inactivation mode
+Kmn = 0.002 [mM] in [mM]
+k2n = 1000 [1/ms] in [1/ms]
+km2n = jca * 1 [1/ms] in [1/ms]
+anca = 1 / (k2n / km2n + (1 + Kmn / Ca_ss)^4)
 dot(nca) = anca * k2n - nca * km2n
-    desc: Fraction of channels in Ca-depdent inactivation mode
-    Kmn = 0.002 [mM]
-        in [mM]
-    k2n = 1000 [1/ms]
-        in [1/ms]
-    km2n = jca * 1 [1/ms]
-        in [1/ms]
-    anca = 1 / (k2n / km2n + (1 + Kmn / Ca_ss)^4)
-# Total currents through L-type calcium channels
-PhiCaL = if(abs(vf) < 1e-6,
-            2 * phys.F * (Ca_ss - 0.341 * Ca_o),
-            4 * vff * (Ca_ss * evf2 - 0.341 * Ca_o) / (evf2 - 1))
+    desc: Fraction of channels in Ca-dependent inactivation mode
+# Permeabilities
+vf = V * phys.FRT
+vff = V * phys.FFRT
+    in [C/mol]
+evf = exp(vf)
+evf2 = exp(2 * vf)
+PhiCa = if(abs(vf) < 1e-6,
+           2 * phys.F * (Ca_ss - 0.341 * Ca_o),
+           4 * vff * (Ca_ss * evf2 - 0.341 * Ca_o) / (evf2 - 1))
     in [mC/L]
-PhiCaNa = if(abs(vf) < 1e-6,
-             0.75 * phys.F * (Na_ss - Na_o),
-             0.75 * vff * (Na_ss * evf - Na_o) / (evf - 1))
+PhiNa = if(abs(vf) < 1e-6,
+           0.75 * phys.F * (Na_ss - Na_o),
+           0.75 * vff * (Na_ss * evf - Na_o) / (evf - 1))
     in [mC/L]
-PhiCaK = if(abs(vf) < 1e-6,
-            0.75 * phys.F * (K_ss - K_o),
-            0.75 * vff * (K_ss * evf - K_o) / (evf - 1))
+PhiK = if(abs(vf) < 1e-6,
+          0.75 * phys.F * (K_ss - K_o),
+          0.75 * vff * (K_ss * evf - K_o) / (evf - 1))
     in [mC/L]
-PCa_base = 0.0001 [L/ms/F]
-    in [L/ms/F]
-    label gCaL
+PCa_base = 0.0001 [L/F/ms]
+    in [L/F/ms]
 PCa = PCa_base * piecewise(cell.mode == 0, 1, cell.mode == 1, 1.2, 2.5)
-    in [L/ms/F]
-PCap = 1.1 * PCa
-    in [L/ms/F]
-PCaNa = 0.00125 * PCa
-    in [L/ms/F]
-PCaK = 3.574e-4 * PCa
-    in [L/ms/F]
-PCaNap = 0.00125 * PCap
-    in [L/ms/F]
-PCaKp = 3.574e-4 * PCap
-    in [L/ms/F]
-g = d * (f  * (1 - nca) + jca * fca  * nca)
-    desc: Conductivity of non-phosphorylated ICaL
-gp = d * (fp * (1 - nca) + jca * fcap * nca)
-    desc: Conductivity of phosphorylated ICaL
-ICaL = (1 - camk.f) * PCa * PhiCaL * g + camk.f * PCap * PhiCaL * gp
+    in [L/F/ms]
+PNa = 0.00125 * PCa
+    in [L/F/ms]
+PK = 3.574e-4 * PCa
+    in [L/F/ms]
+Pp = 1.1
+    desc: Permeation multiplication factor for phosphorylated channel    
+# Conductivity (weighted sum of gates)
+g = d * ((f  * (1 - nca) + jca * fca * nca) * (1 - camk.f) +
+         (fp * (1 - nca) + jca * fcap * nca) * camk.f * Pp)
+# Current
+ICaLCa = g * PCa * PhiCa
     in [A/F]
     desc: Ca-component of L-type calcium current
-ICaNa = (1 - camk.f) * PCaNa * PhiCaNa * g + camk.f * PCaNap * PhiCaNa * gp
+ICaLNa = g * PNa * PhiNa
     in [A/F]
     desc: Na-component of L-type calcium current
-ICaK = (1 - camk.f) * PCaK * PhiCaK * g + camk.f * PCaKp * PhiCaK * gp
+ICaLK = g * PK * PhiK
     in [A/F]
     desc: K-component of L-type calcium current
-ICaL_total = ICaL + ICaK + ICaNa
+ICaL_tot = ICaLCa + ICaLK + ICaLNa
     in [A/F]
     label ICaL
 
@@ -805,7 +774,7 @@ Kxkur = 292 [mM]    # Original
 Knap = 224 [mM]     # Original
     in [mM]
     desc: Dissociation constant relating [NaPi] and [Na]i
-# Table 1 parameters in Smith & Cramin
+# Table 1 parameters in Smith & Crampin
 k1p = 949.5 [1/s]       # Refit from 1050
     in [1/s]
 k1m = 182.4 [1/s/mM]    # Refit from 172.1
@@ -848,38 +817,38 @@ Knao = Knao0 * exp((1 - delta) * V * phys.FRT / 3)
     in [mM]
 # Equations 30: Forward (clockwise) rates
 a1 = k1p * (Na_i / Knai)^3 / ((1 + Na_i / Knai)^3 + (1 + K_i / Kki)^2 - 1)
-    in [Hz]
+    in [1/s]
 a2 = k2p
-    in [Hz]
+    in [1/s]
 a3 = k3p * (K_o / Kko)^2 / ((1 + Na_o / Knao)^3 + (1 + K_o / Kko)^2 - 1)
-    in [Hz]
+    in [1/s]
 a4 = k4p * MgATP / Kmgatp / (1 + MgATP / Kmgatp)
-    in [Hz]
+    in [1/s]
 # Equations 31: Backward (anticlockwise) rates
 b1 = k1m * MgADP
-    in [Hz]
+    in [1/s]
 b2 = k2m * (Na_o / Knao)^3 / ((1 + Na_o / Knao)^3 + (1 + K_o / Kko)^2 - 1)
-    in [Hz]
+    in [1/s]
 b3 = k3m * P * H / (1 + MgATP / Kmgatp)
-    in [Hz]
+    in [1/s]
 b4 = k4m * (K_i / Kki)^2 / ((1 + Na_i / Knai)^3 + (1 + K_i / Kki)^2 - 1)
-    in [Hz]
+    in [1/s]
 # Equation 33
 P = eP / (1 + H / Khp + Na_i / Knap + K_i / Kxkur)
     in [mM]
 # Terms used to calculate the steady-state occupancies of the four states
 # Based on the "diagram method" described in [2], these terms are the sums of
 # the reaction rates of all (non-cyclical) paths leading to each state
 x1 = a4 * a1 * a2 + b1 * b4 * b3 + a2 * b4 * b3 + b3 * a1 * a2
-    in [Hz^3]
+    in [s^-3]
     note: Corrected from the original code (b1 in second term)
 x2 = b2 * b1 * b4 + a1 * a2 * a3 + a3 * b1 * b4 + a2 * a3 * b4
-    in [Hz^3]
+    in [s^-3]
 x3 = a2 * a3 * a4 + b3 * b2 * b1 + b2 * b1 * a4 + a3 * a4 * b1
-    in [Hz^3]
+    in [s^-3]
 x4 = b4 * b3 * b2 + a3 * a4 * a1 + b2 * a4 * a1 + b3 * b2 * a1
-    in [Hz^3]
-# Cycle rate (obtaining by writing any one of the four flux equations: they all
+    in [s^-3]
+# Cycle rate (obtained by writing any one of the four flux equations: they all
 # give the same result so that the pump count is conserved).
 r = (a1 * a2 * a3 * a4 - b1 * b2 * b3 * b4) / (x1 + x2 + x3 + x4)
     in [1/s]
@@ -964,7 +933,7 @@ dot(Jrelnp) = (inf - Jrelnp) / tau
             in [ms]
     inf = base * if(cell.mode == 2, 1.7, 1)
         in [mM/ms]
-        base = -1 [mM/ms/mV] * a_rel * ical.ICaL / (1 + (1.5 [mM] / Ca_jsr)^8)
+        base = -1 [mM/ms/mV] * a_rel * ical.ICaLCa / (1 + (1.5 [mM] / Ca_jsr)^8)
             in [mM/ms]
 btp = 1.25 * bt
     in [ms]
@@ -978,7 +947,7 @@ dot(Jrelp) = (inf - Jrelp) / tau
             in [ms]
     inf = base * if(cell.mode == 2, 1.7, 1)
         in [mM/ms]
-        base = -1 [mM/ms/mV] * a_relp * ical.ICaL / (1 + (1.5 [mM] / Ca_jsr)^8)
+        base = -1 [mM/ms/mV] * a_relp * ical.ICaLCa / (1 + (1.5 [mM] / Ca_jsr)^8)
             in [mM/ms]
 Jrel = (1 - camk.f) * Jrelnp + camk.f * Jrelp
     desc: SR Calcium release flux via Ryanodine receptor
@@ -1027,7 +996,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]
@@ -1047,7 +1016,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
@@ -1090,7 +1059,7 @@ dot(Ca_i) = buff * (-ICa_tot * AFC / (2 * vmyo) - serca.Jup * vnsr / vmyo + diff
     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 * vjsr / vss - diff.Jdiff)
     in [mM]