-
Notifications
You must be signed in to change notification settings - Fork 2
/
RatAPtest01.ode
212 lines (201 loc) · 10.8 KB
/
RatAPtest01.ode
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
# RatAPtest01.ode
# Rat action potential is simulated (bondarenko_model_2004)
# (cf. Bondarenko et al. Am J Physiol Heart Circ Physiol 2004;287: H1378-H1403)
# (cf. Wang et al. Toxicol Sci 2008;106:454-463)
# implemented by Dr. Sheng-Nan Wu
# Constant values
Number Temp=298.0, R=8.314, Fara=96.5
Number E_CaL=63.0, E_Cl=-40.0
Number Cm=1.0, Vmyo=2.584E-5, VNSR=2.098E-6, VJSR=1.2E-7
# Initial values
Initial V=-82.4202, Cai=0.115001
Initial Nai=14237.1, Ki=143720.0
Initial C_Na1=2.79132E-4, C_Na2=0.020752, O_Na=7.13483E-7
Initial IF_Na=1.53176E-4, I1_Na=6.73345E-7, I2_Na=1.55787E-9
Initial IC_Na2=0.0113879, IC_Na3=0.34278
Initial ato_f=0.00265563, ito_f=0.999977
Initial ato_s=4.17069E-4, ito_s=0.998543
Initial nKs=2.62753E-4, aur=4.17069E-4, iur=0.998543
Initial aKss=4.17069E-4, iKss=1.0
Initial C_K1=9.92513E-4, C_K2=6.41229E-4
Initial O_K=1.75298E-4, I_K=3.19129E-5
Initial Cass=0.115001, CaJSR=1299.5, CaNSR=1299.5
Initial P_RyR=0.0
Initial LTRPN_Ca=11.2684, HTRPN_Ca=125.29
Initial O=9.30308E-19, O1=1.49102E-5, O2=9.51726E-11
Initial P_C2=1.6774E-4
Initial C2=1.24216E-4, C3=5.78679E-9, C4=1.19816E-13
Initial I1=4.97023E-19, I2=3.45847E-14, I3=1.85106E-14
# Stimulus protocol
Par r0=1, period=300, pulse=10
Par tf=0, tp=5, tstart=10
ts=t-tstart
rstar = r0 + pulse*(heav(mod(ts,period)-tf)-heav(mod(ts,period)-(tf+tp)))
Par Nao=140000.0
Par LTRPN_tot=70.0
Par kf=0.023761
Par kb=0.036778
Par Km_CSQN=800.0
Par i_pCa_max=1.0
Par CMDN_tot=50.0
Par Kpc_max=0.23324
Par Km_CMDN=0.238
Par K_mCa=1380.0
Par Km_Nai=21.0
Par v1=4.5, v2=1.7E-5, v3=0.45
Par k_plus_c=0.0090
Par k_minus_ltrpn=0.196
Par g_Cab=3.67E-4
Par k_plus_b=0.00405
Par k_plus_a=0.006075
Par Km_up=0.5
Par K_mNa=87500.0
Par Km_Cl=10.0
Par CSQN_tot=15000.0
Par k_minus_htrpn =3.2E-5
Par Km_pCa=0.5
Par k_plus_ltrpn=0.0327
Par Acap=1.534E-4
Par g_Na=13.0, g_CaL=1.51729, g_Ks=0.000575, g_Kr=0.0078, g_Kur=0.0016
Par g_ClCa=10.0
Par Vss=1.485E-9
Par k_NaCa=992.8
Par k_plus_htrpn=0.00237
Par Ko=5400.0, Cao=1800.0
Par tau_xfer=8.0
Par k_sat =0.1
Par n=4.0
Par g_Kto_f=0.4067, g_Kto_s=0.01, g_Nab=0.0026
Par m=3.0
Par Km_Ko=1.5
Par g_Kss=0.05
Par HTRPN_tot=140.0
Par i_NaK_max=0.88
Par i_CaL_max=7.0
Par tau_tr=20.0
Par k_minus_c=8.0E-4
Par k_minus_b=0.965
Par k_minus_a=0.07125
Par Kpc_half=20.0
Par eta=0.35
Par Kpcb=5.0E-4
ass = (1.0 / (1.0 + exp( - (0.12987012987012986 * (22.5 + V)))))
E_K = (ln((Ko / Ki)) * R * Temp / Fara)
i_Ks = (g_Ks * (nKs^2.0) * (V - E_K))
E_Na = (ln((((0.9 * Nao) + (0.1 * Ko)) / ((0.9 * Nai) + (0.1 * Ki)))) * R * Temp / Fara)
i_Kr = (g_Kr * O_K * (V - (ln((((0.98 * Ko) + (0.02 * Nao)) / ((0.98 * Ki) + (0.02 * Nai)))) * R * Temp / Fara)))
i_Kur = (g_Kur * aur * iur * (V - E_K))
sigma = (0.14285714285714285 * (-1.0 + exp((1.4858841010401188E-5 * Nao))))
f_NaK = (1.0 / (1.0 + (0.1245 * exp( - (0.1 * V * Fara / (R * Temp)))) + (0.0365 * sigma * exp( - (Fara * V / (R * Temp))))))
J_leak = (v2 * (CaNSR - Cai))
i_NaCa = (k_NaCa * ((exp((eta * V * Fara / (R * Temp))) * (Nai^3.0) * Cao) - (exp(((-1.0 + eta) * V * Fara / (R * Temp))) * (Nao^3.0) * Cai)) / ((K_mNa^3.0) + (Nao^3.0)) / (K_mCa + Cao) / (1.0 + (k_sat * exp(((-1.0 + eta) * V * Fara / (R * Temp))))))
O_ClCa = (0.2 / (1.0 + exp( - (0.12820512820512822 * (-46.7 + V)))))
i_ClCa = (g_ClCa * O_ClCa * (V - E_Cl) * Cai / (Cai + Km_Cl))
i_Nab = (g_Nab * (V - E_Na))
tau_Kss = (13.17 + (39.3 * exp( - (0.0862 * V))))
beta_i1 = (9.5E-4 * exp((0.14285714285714285 * (33.5 + V))) / (1.0 + (0.051335 * exp((0.14285714285714285 * (33.5 + V))))))
beta = (0.05 * exp( - (0.07692307692307693 * (12.0 + V))))
i_Kss = (g_Kss * aKss * iKss * (V - E_K))
i_K1 = (0.2938 * Ko * (V - E_K) / (210.0 + Ko) / (1.0 + exp((0.0896 * (V - E_K)))))
i_NaK = (i_NaK_max * f_NaK * Ko / (1.0 + ((Km_Nai / Nai)^0.66667)) / (Ko + Km_Ko))
J_rel = (v1 * (O1 + O2) * (CaJSR - Cass) * P_RyR)
alpha_n = (4.81333E-6 * (26.5 + V) * (1.0 - exp( - (0.128 * (26.5 + V)))))
alpha_i= (0.090821 * exp((0.023391 * (5.0 + V))))
alpha_a = (0.18064 * exp((0.03577 * (30.0 + V))))
i_pCa = (i_pCa_max * (Cai^2.0) / ((Km_pCa^2.0) + (Cai^2.0)))
alpha_i1 = (1.52E-4 * exp( - (0.14285714285714285 * (13.5 + V))) / (1.0 + (0.0067083 * exp( - (0.14285714285714285 * (33.5 + V))))))
# i_stim = - (10.0 * (t > 10.0) * (t < 15.0))
beta_Na13 = (0.22 * exp( - (0.04926108374384236 * (-7.5 + V))))
beta_Na12 = (0.2 * exp( - (0.04926108374384236 * (-2.5 + V))))
beta_Na11 = (0.1917 * exp( - (0.04926108374384236 * (2.5 + V))))
alpha = (0.4 * exp((0.1 * (12.0 + V))) * (1.0 - (0.75 * exp( - (0.0025 * ((20.0 + V)^2.0)))) + (0.7 * exp( - (0.1 * ((40.0 + V)^2.0))))) / (1.0 + (0.12 * exp((0.1 * (12.0 + V))))))
Bi = ((1.0 + (CMDN_tot * Km_CMDN / ((Km_CMDN + Cai)^2.0)))^(-1.0))
gamma = (Kpc_max * Cass / (Kpc_half + Cass))
P_C1= (1.0 - (P_C2 + O1 + O2))
C1= (1.0 - (O + C2 + C2 + C3 + C4 + I1 + I2 + I3))
beta_n= (9.53333E-5 * exp( - (0.038 * (26.5 + V))))
beta_i= (0.006497 * exp( - (0.03268 * (5.0 + V))))
beta_a= (0.3956 * exp( - (0.06237 * (30.0 + V))))
i_Kto_s = (g_Kto_s * ato_s * ito_s * (V - E_K))
i_Kto_f = (g_Kto_f * (ato_f^3.0) * ito_f * (V - E_K))
alpha_Na3 = (7.0E-7 * exp( - (0.12987012987012986 * (7.0 + V))))
beta_Na5 = (0.02 * alpha_Na3)
beta_Na4 = alpha_Na3
beta_Na3 = (0.0084 + (2.0E-5 * (7.0 + V)))
alpha_Na13 = (3.802 / ((0.1027 * exp( - (0.08333333333333333 * (2.5 + V)))) + (0.25 * exp( - (0.006666666666666667 * (2.5 + V))))))
alpha_Na2 = (1.0 / (0.393956 + (0.188495 * exp( - (0.06024096385542168 * (7.0 + V))))))
beta_Na2 = (alpha_Na13 * alpha_Na2 * alpha_Na3 / (beta_Na13 * beta_Na3))
BJSR = ((1.0 + (CSQN_tot * Km_CSQN / ((Km_CSQN + CaJSR)^2.0)))^(-1.0))
tau_iur= (1200.0 - (170.0 / (1.0 + exp((0.17543859649122806 * (45.2 + V))))))
E_CaN= (ln((Cao / Cai)) * R * Temp / (2.0 * Fara))
i_Cab= (g_Cab * (V - E_CaN))
tau_ti_s= (270.0 + (1050.0 / (1.0 + exp((0.17543859649122806 * (45.2 + V))))))
Bss= ((1.0 + (CMDN_tot * Km_CMDN / ((Km_CMDN + Cass)^2.0)))^(-1.0))
alpha_Na12= (3.802 / ((0.1027 * exp( - (0.06666666666666667 * (2.5 + V)))) + (0.23 * exp( - (0.006666666666666667 * (2.5 + V))))))
alpha_Na11= (3.802 / ((0.1027 * exp( - (0.058823529411764705 * (2.5 + V)))) + (0.2 * exp( - (0.006666666666666667 * (2.5 + V))))))
tau_aur= (2.058 + (0.493 * exp( - (0.0629 * V))))
beta_a1= (6.89E-5 * exp( - (0.04178 * V)))
beta_a0= (0.047002 * exp( - (0.0631 * V)))
J_xfer= ((Cass - Cai) / tau_xfer)
tau_ta_s= (2.058 + (0.493 * exp( - (0.0629 * V))))
J_up= (v3 * (Cai^2.0) / ((Km_up^2.0) + (Cai^2.0)))
alpha_Na5 = (1.0526315789473684E-5 * alpha_Na2)
alpha_Na4 = (0.0010 * alpha_Na2)
iss = (1.0 / (1.0 + exp((0.17543859649122806 * (45.2 + V)))))
Kpcf = (13.0 * (1.0 - exp( - (0.01 * ((14.5 + V)^2.0)))))
C_K0 = (1.0 - (C_K1 + C_K2 + O_K + I_K))
J_tr = ((CaNSR - CaJSR) / tau_tr)
alpha_a1= (0.013733 * exp((0.038198 * V)))
alpha_a0 = (0.022348 * exp((0.01176 * V)))
J_trpn = ( - ((k_minus_htrpn * HTRPN_Ca) + (k_minus_ltrpn * LTRPN_Ca)) + (k_plus_htrpn * Cai * (HTRPN_tot - HTRPN_Ca)) + (k_plus_ltrpn * Cai * (LTRPN_tot - LTRPN_Ca)))
C_Na3 = (1.0 - (O_Na + C_Na1 + C_Na2 + IF_Na + I1_Na + I2_Na + IC_Na2 + IC_Na3))
i_Na= (g_Na * O_Na * (V - E_Na))
i_CaL= (g_CaL * O * (V - E_CaL))
V'= - ((i_CaL + i_pCa + i_NaCa + i_Cab + i_Na + i_Nab + i_NaK + i_Kto_f + i_Kto_s + i_K1 + i_Ks + i_Kur + i_Kss + i_Kr + i_ClCa - rstar) / Cm)
Cai'=(Bi * ( - (J_up + J_trpn + (( - (2.0 * i_NaCa) + i_Cab + i_pCa) * Acap * Cm / (2.0 * Vmyo * Fara))) + J_leak + J_xfer))
Cass'=(Bss * ((J_rel * VJSR / Vss) - ((J_xfer * Vmyo / Vss) + (i_CaL * Acap * Cm / (2.0 * Vss * Fara)))))
CaJSR'=(BJSR * (J_tr - J_rel))
CaNSR'=(((J_up - J_leak) * Vmyo / VNSR) - (J_tr * VJSR / VNSR))
P_RyR'=( - (0.04 * P_RyR) - (0.1 * exp( - (0.0015432098765432098 * ((-5.0 + V)^2.0))) * i_CaL / i_CaL_max))
LTRPN_Ca'=((k_plus_ltrpn * Cai * (LTRPN_tot - LTRPN_Ca)) - (k_minus_ltrpn * LTRPN_Ca))
HTRPN_Ca'=((k_plus_htrpn * Cai * (HTRPN_tot - HTRPN_Ca)) - (k_minus_htrpn * HTRPN_Ca))
O1'=( - ((k_minus_a * O1) + (k_plus_b * (Cass^m) * O1) + (k_plus_c * O1)) + (k_plus_a * (Cass^n) * P_C1) + (k_minus_b * O2) + (k_minus_c * P_C2))
O2'=((k_plus_b * (Cass^m) * O1) - (k_minus_b * O2))
P_C2'=((k_plus_c * O1) - (k_minus_c * P_C2))
O'=( - ((4.0 * beta * O) + (gamma * O)) + (alpha * C4) + (Kpcb * I1) + (0.0010 * ((alpha * I2) - (Kpcf * O))))
C2'=( - ((beta * C2) + (3.0 * alpha * C2)) + (4.0 * alpha * C1) + (2.0 * beta * C3))
C3'=( - ((2.0 * beta * C3) + (2.0 * alpha * C3)) + (3.0 * alpha * C2) + (3.0 * beta * C4))
C4'=( - ((3.0 * beta * C4) + (alpha * C4) + (gamma * Kpcf * C4)) + (2.0 * alpha * C3) + (4.0 * beta * O) + (0.01 * ((4.0 * Kpcb * beta * I1) - (alpha * gamma * C4))) + (0.0020 * ((4.0 * beta * I2) - (Kpcf * C4))) + (4.0 * beta * Kpcb * I3))
I1'=( - (Kpcb * I1) + (gamma * O) + (0.0010 * ((alpha * I3) - (Kpcf * I1))) + (0.01 * ((alpha * gamma * C4) - (4.0 * beta * Kpcf * I1))))
I2'=( - (gamma * I2) + (0.0010 * ((Kpcf * O) - (alpha * I2))) + (Kpcb * I3) + (0.0020 * ((Kpcf * C4) - (4.0 * beta * I2))))
I3'=( - ((4.0 * beta * Kpcb * I3) + (Kpcb * I3)) + (0.0010 * ((Kpcf * I1) - (alpha * I3))) + (gamma * I2) + (gamma * Kpcf * C4))
Nai'= - (Acap * Cm * (i_Na + i_Nab + (3.0 * i_NaK) + (3.0 * i_NaCa)) / (Vmyo * Fara))
C_Na2'=( - ((beta_Na11 * C_Na2) + (alpha_Na12 * C_Na2) + (beta_Na3 * C_Na2)) + (alpha_Na11 * C_Na3) + (beta_Na12 * C_Na1) + (alpha_Na3 * IC_Na2))
C_Na1'=( - ((beta_Na12 * C_Na1) + (alpha_Na13 * C_Na1) + (beta_Na3 * C_Na1)) + (alpha_Na12 * C_Na2) + (beta_Na13 * O_Na) + (alpha_Na3 * IF_Na))
O_Na'=( - ((beta_Na13 * O_Na) + (alpha_Na2 * O_Na)) + (alpha_Na13 * C_Na1) + (beta_Na2 * IF_Na))
IF_Na'=( - ((beta_Na2 * IF_Na) + (alpha_Na3 * IF_Na) + (alpha_Na4 * IF_Na) + (beta_Na12 * IF_Na)) + (alpha_Na2 * O_Na) + (beta_Na3 * C_Na1) + (beta_Na4 * I1_Na) + (alpha_Na12 * IC_Na2))
I1_Na'=( - ((beta_Na4 * I1_Na) + (alpha_Na5 * I1_Na)) + (alpha_Na4 * IF_Na) + (beta_Na5 * I2_Na))
I2_Na'=((alpha_Na5 * I1_Na) - (beta_Na5 * I2_Na))
IC_Na2'=( - ((beta_Na11 * IC_Na2) + (alpha_Na12 * IC_Na2) + (alpha_Na3 * IC_Na2)) + (alpha_Na11 * IC_Na3) + (beta_Na12 * IF_Na) + (beta_Na3 * IC_Na2))
IC_Na3'=( - ((alpha_Na11 * IC_Na3) + (alpha_Na3 * IC_Na3)) + (beta_Na11 * IC_Na2) + (beta_Na3 * C_Na3))
Ki'= - (Acap * Cm * (i_Kto_f + i_Kto_s + i_K1 + i_Ks + i_Kss + i_Kur + i_Kr + (2.0 * i_NaK)) / (Vmyo * Fara))
ato_f'=((alpha_a * (1.0 - ato_f)) - (beta_a * ato_f))
ito_f'=((alpha_i1 * (1.0 - ito_f)) - (beta_i1 * ito_f))
ato_s' = ((ass - ato_s) / tau_ta_s)
ito_s' = ((iss - ito_s) / tau_ti_s)
nKs'=((alpha_n * (1.0 - nKs)) - (beta_n * nKs))
aur'=((ass - aur) / tau_aur)
iur'=((iss - iur) / tau_iur)
aKss'=((ass - aKss) / tau_Kss)
iKss'=0.0
C_K2'=( - ((kb * C_K2) + (alpha_a1 * C_K2)) + (kf * C_K1) + (beta_a1 * O_K))
C_K1'=( - ((beta_a0 * C_K1) + (kf * C_K1)) + (alpha_a0 * C_K0) + (kb * C_K2))
O_K'=( - ((beta_a1 * O_K) + (alpha_i * O_K)) + (alpha_a1 * C_K2) + (beta_i * I_K))
I_K'=((alpha_i * O_K) - (beta_i * I_K))
aux ina=i_Na
aux ical=i_CaL
aux ikur=i_Kur
# Numerical and plotting parameters for xpp
@ maxstor=5000000, bounds=10000000, total=1100, xp=t, yp=V, trans=860
@ meth=Gear, dt=0.05, toler=0.01, xlo=860, xhi=1100, ylo=-100, yhi=60
done