using ModelingToolkit
using OrdinaryDiffEq, SteadyStateDiffEq, DiffEqCallbacks
using Plots
using CSV
using DataFrames
using Dates
using CaMKIIModel
using CaMKIIModel: second
Plots.default(lw=1.5)Caffeine Effects
Caffeine increase RyR opening sensitivity to luminal and subspace calcium.
In this model, we decrease the mid saturation sub-SR calcium concentration for the opening rate.
sys = build_neonatal_ecc_sys(simplify=true, reduce_iso=true, reduce_camk=true)
tend = 500second
prob = ODEProblem(sys, [], tend)
stimstart = 100second
stimend = 300second
alg = KenCarp47()
function add_coffee_affect!(integrator)
integrator.ps[sys.RyRsensitivity] = 10
endadd_coffee_affect! (generic function with 1 method)
@unpack Istim = sys
callback = build_stim_callbacks(Istim, stimend; period=1second, starttime=stimstart)SciMLBase.CallbackSet{Tuple{}, Tuple{SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}, SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}}}((), (SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}(100000:1000:300000, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim)), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim), DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}(100000:1000:300000, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim)), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ()), SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}(100000.5:1000.0:300000.5, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim)), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim), DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}(100000.5:1000.0:300000.5, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim)), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ())))
Add caffeine at t = 200 econd
callback_caf = CallbackSet(build_stim_callbacks(Istim, stimend; period=1second, starttime=stimstart), PresetTimeCallback(200.0second, add_coffee_affect!))SciMLBase.CallbackSet{Tuple{}, Tuple{SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}, SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}, SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##277".add_coffee_affect!)}, typeof(Main.var"##277".add_coffee_affect!), DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##277".add_coffee_affect!)}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}}}((), (SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}(100000:1000:300000, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim)), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim), DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}(100000:1000:300000, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim)), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ()), SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}(100000.5:1000.0:300000.5, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim)), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim), DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}(100000.5:1000.0:300000.5, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim)), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ()), SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##277".add_coffee_affect!)}, typeof(Main.var"##277".add_coffee_affect!), DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##277".add_coffee_affect!)}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##277".add_coffee_affect!)}([200000.0], true, SciMLBase.INITIALIZE_DEFAULT, Main.var"##277".add_coffee_affect!), Main.var"##277".add_coffee_affect!, DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##277".add_coffee_affect!)}([200000.0], true, SciMLBase.INITIALIZE_DEFAULT, Main.var"##277".add_coffee_affect!), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ())))
Single-dose caffeine
@time sol = solve(prob, alg; callback)
@time sol_caf = solve(prob, alg; callback=callback_caf) 10.154164 seconds (23.35 M allocations: 1.454 GiB, 2.17% gc time, 83.38% compilation time)
4.511128 seconds (8.11 M allocations: 864.038 MiB, 6.60% gc time, 61.92% compilation time)
retcode: Success
Interpolation: 3rd order Hermite
t: 8773-element Vector{Float64}:
0.0
0.019333123284849377
0.09098129546758026
0.20735582943983927
0.4561191027686269
1.047256292840616
2.0342667657223767
3.4517766411484665
6.733792423912831
12.45011708102324
⋮
400587.1730516022
410672.0397456905
421716.44761840976
434008.49190990673
450151.06267123023
462974.23931692366
480604.5572769363
498234.87523694895
500000.0
u: 8773-element Vector{Vector{Float64}}:
[150952.75035000002, 13838.37602, -68.79268, 0.12113, 0.12113, 0.12113, 0.12113, 0.12113, 0.12113, 0.12113 … 0.26081, 0.00977, 0.00188, 0.09243, 0.22156, 0.966, 0.00702, 830.0, 830.0, 0.0026]
[150952.7504077194, 13838.375711519551, -68.79732946186226, 0.1211319185255175, 0.1211301421535992, 0.12113000811850358, 0.12113000037367234, 0.12113000001389383, 0.12113000000038378, 0.12113000000000405 … 0.2608055799112056, 0.009769882842069916, 0.001879826176081169, 0.09242766669618695, 0.2215639938807744, 0.9660018276402699, 0.007020208561323646, 829.9999654017583, 829.9994975880693, 0.0025985659590295374]
[150952.75062244057, 13838.374513299192, -68.81607005265552, 0.12115248970210647, 0.12113624253259858, 0.12113143375871602, 0.12113027926811751, 0.12113004704056851, 0.12113000693481112, 0.12113000089551408 … 0.26078941009382306, 0.009769459049732313, 0.0018789042446514422, 0.09241901973901424, 0.2215787945110049, 0.9660085993318799, 0.00702008986027517, 829.9998348635061, 829.9976422525165, 0.002593310524681686]
[150952.7509740238, 13838.372555328333, -68.84681913285226, 0.12118951551963696, 0.12115723910142134, 0.12114109135582982, 0.12113403405679991, 0.12113131833063846, 0.12113038966346883, 0.12113010483012127 … 0.2607638707001999, 0.009768805322556195, 0.0018766436937193197, 0.09240497563401104, 0.22160283236040595, 0.9660195931792828, 0.007017234889558555, 829.9996133960318, 829.9946506157202, 0.0025849478335103772]
[150952.7517373371, 13838.368384399068, -68.91207316777846, 0.12125418995384042, 0.12120713813483534, 0.12117519662647318, 0.12115494411090462, 0.12114295934422188, 0.12113633873395535, 0.12113292110288301 … 0.2607122383588733, 0.00976754294689657, 0.001869919626294057, 0.0923749578279845, 0.2216542075240971, 0.9660430731013478, 0.007003014305234393, 829.9990964774663, 829.9883445519497, 0.0025677078010499796]
[150952.75361533565, 13838.358559292059, -69.06429960994139, 0.12136450208851986, 0.12130632552358439, 0.12125925929399749, 0.12122229800810431, 0.12119414515426226, 0.12117335985880448, 0.1211584926966588 … 0.26060507831375035, 0.009765228548607825, 0.0018505337694469416, 0.09230364339951613, 0.22177624785138977, 0.9660987611411453, 0.006946470397388756, 829.9976227528975, 829.9738229154025, 0.0025298083713907]
[150952.75694915713, 13838.34240329337, -69.31024285744344, 0.1214919651325425, 0.12142927821043346, 0.1213742637991964, 0.12132672249013038, 0.12128628091507475, 0.12125242557826563, 0.12122454121844704 … 0.26047114911179464, 0.009763403076884485, 0.0018178174419014655, 0.09218462776986637, 0.2219798989175973, 0.9661914345949096, 0.006833754172875444, 829.9944193668179, 829.9509077855928, 0.002474807978215895]
[150952.76215620188, 13838.319723120678, -69.64637448767434, 0.12162015426876964, 0.12155743299711909, 0.12149988983034504, 0.12144759471003595, 0.12140052118634267, 0.12135855488207911, 0.12132150441596078 … 0.2603648729749554, 0.009765046552800477, 0.0017738156907532213, 0.09201383579531924, 0.22227216398019892, 0.9663239343003682, 0.0066755120626567415, 829.9883533251192, 829.9205247545176, 0.0024106893964811306]
[150952.77598869594, 13838.2694938995, -70.35063788657429, 0.12180447948880603, 0.12174768443532381, 0.12169340184085989, 0.12164180995533579, 0.12159304920072045, 0.12154722283384199, 0.12150439833798565 … 0.26041507667541997, 0.009786689200123853, 0.00168558990949223, 0.09161908329520764, 0.22294828828976576, 0.9666286957812559, 0.006355039713677853, 829.9689941817725, 829.8589858006806, 0.002310655727734336]
[150952.80529264064, 13838.18882346682, -71.35304509818073, 0.1219637836440816, 0.12191937753328323, 0.1218757903829034, 0.12183318573340736, 0.12179171693004404, 0.12175152667167365, 0.12171274673876364 … 0.26104417804603025, 0.009873892913308327, 0.0015688568862426027, 0.09093422465365385, 0.22412560303273427, 0.9671554512230358, 0.005926150085676066, 829.923524981244, 829.7704141396295, 0.002228995040224513]
⋮
[151617.6676592487, 13261.115089921832, -69.79226031208901, 0.1325222509745926, 0.13252225738575027, 0.13252225946706916, 0.1325222571327741, 0.13252225029362644, 0.1325222388567477, 0.13252222272543157 … 0.0012318407145399332, 0.0016437171132340357, 0.0025443412202716138, 0.0016562010175567839, 0.999712289565461, 0.9997103653553531, 0.00656466248745671, 51.350329001024605, 50.044323004799274, 0.115407907233062]
[151662.98368801433, 13215.809140524345, -69.70313014241552, 0.13198341629455648, 0.13198342106338237, 0.13198342172031022, 0.13198341818362508, 0.131983410368313, 0.1319833981858924, 0.1319833815442357 … 0.0012489507778964324, 0.0016546371159752405, 0.0025612404919185556, 0.0016775215870001295, 0.9997066419289161, 0.9997047470874393, 0.006604554777213319, 51.643412591601084, 50.347438135877425, 0.11383532954450551]
[151710.30673930817, 13168.495907932493, -69.60892118803808, 0.13142872799590954, 0.1314287313223141, 0.13142873076453182, 0.13142872624513433, 0.13142871768356684, 0.13142870499598855, 0.1314286880951032 … 0.001267289894171816, 0.0016662599533203907, 0.002579227556422485, 0.0017003472952818626, 0.9997005540036995, 0.9996986884554673, 0.006646996318156452, 51.94957514605457, 50.66389625800978, 0.1122294374515424]
[151760.2768320601, 13118.535431464285, -69.5081814029278, 0.13085191887946276, 0.130851920740174, 0.1308519189553221, 0.130851913451993, 0.13085190415432818, 0.13085189098337457, 0.13085187385692496 … 0.0012872048952438473, 0.0016787760463155051, 0.002598596191368288, 0.0017250944283408766, 0.9996939010960414, 0.9996920724410305, 0.006692668893964975, 52.27265495674743, 50.99764934817281, 0.11057349973713322]
[151821.82285620715, 13057.000208053652, -69.38225539818865, 0.13015436298139466, 0.13015436331034388, 0.1301543602854011, 0.1301543538392183, 0.13015434390172795, 0.13015433040000418, 0.13015431325811516 … 0.0013125371907681568, 0.0016945515225815888, 0.002623008775847405, 0.0017565205192994384, 0.9996853754371718, 0.9996835988371647, 0.006750196341662394, 52.669984323265304, 51.407840682047784, 0.10859009372114425]
[151867.62161195502, 13011.208765438236, -69.28718604073126, 0.129644636814347, 0.12964463619644956, 0.12964463243996532, 0.12964462548169214, 0.12964461525587573, 0.12964460169407915, 0.1296445847250447 … 0.0013319950194807383, 0.0017065584238387316, 0.002641589270677801, 0.0017806175373668359, 0.9996787803273957, 0.9996770470658016, 0.006793952882277252, 52.96497273704035, 51.712196201042566, 0.10715401804989967]
[151926.42556655762, 12952.413327034952, -69.16335640675119, 0.12900206903827027, 0.12900206742140347, 0.1290020629391372, 0.12900205553356603, 0.12900204514444594, 0.1290020317090752, 0.12900201516216758 … 0.0013577735754901903, 0.001722325323410439, 0.0026659882439735477, 0.001812489345798519, 0.9996699802073966, 0.9996683088789485, 0.006851380782062976, 53.342515523529606, 52.101509287684145, 0.10535972793783961]
[151980.74692857257, 12898.099009846683, -69.04716026244081, 0.12842057830989534, 0.12842057596093734, 0.1284205709955261, 0.12842056336061755, 0.1284205530010255, 0.12842053985931212, 0.1284205238756721 … 0.0013824159994341742, 0.0017372497768260962, 0.0026890832862942935, 0.0018429012121774467, 0.9996615026901544, 0.9996598950775745, 0.00690569860319312, 53.68968063181806, 52.459289332119724, 0.10375143088643606]
[151985.95219894202, 12892.894377572118, -69.03593427867983, 0.12836545533469, 0.12836545289178009, 0.1283654478538549, 0.12836544016835133, 0.12836542978058196, 0.12836541663362658, 0.1283654006682173 … 0.0013848250878191657, 0.0017386985654080562, 0.0026913251579741165, 0.001845871912028587, 0.9996606720946383, 0.9996590681190823, 0.006910970123733874, 53.7228354417192, 52.49344894514705, 0.10359972019584344]
i = (sys.t / 1000, sys.vm)
plot(sol, idxs=i, title="Action potential", lab="Ctl", tspan=(198second, 205second))
plot!(sol_caf, idxs=i, lab="Caf", tspan=(198second, 205second), ylabel="Voltage (mV)", xlabel="Time (s)")
i = (sys.t / 1000, sys.Cai_sub_SR * 1000)
plot(sol, idxs=i, title="Calcium transient (During caffeine addition)", lab="Ctl", tspan=(198second, 205second))
plot!(sol_caf, idxs=i, tspan=(198second, 205second), lab="Caf", ylabel="Subspace calcium (nM)", xlabel="Time (s)")
i = (sys.t / 1000, sys.PO1RyR)
plot(sol, idxs=i, title="RyR open (During caffeine addition)", lab="Ctl", tspan=(198second, 205second))
plot!(sol_caf, idxs=i, tspan=(198second, 205second), lab="Caf", ylabel="Open probability", ylims=(0, 1), xlabel="Time (s)")
i = (sys.t / 1000, sys.Cai_sub_SR * 1000)
plot(sol, idxs=i, title="Calcium transient (After caffeine addition)", lab="Ctl", ylabel="Subspace calcium (nM)", tspan=(198second, 205second))
plot!(sol_caf, idxs=i, lab="Caf", xlabel="Time (s)", tspan=(198second, 205second))
i = (sys.t / 1000, sys.CaJSR)
plot(sol, idxs=i, title="SR Calcium (During caffeine addition)", lab="Ctl", ylabel="SR calcium (μM)", tspan=(198second, 205second))
plot!(sol_caf, idxs=i, tspan=(198second, 205second), lab="Caf", ylims=(0, 850), xlabel="Time (s)")
i = (sys.t / 1000, sys.Jrel)
plot(sol, idxs=sys.Jrel, title="Ca flux", lab="Ctl (Jrel)", tspan=(198second, 205second))
plot!(sol_caf, idxs=sys.Jrel, lab="Caf (Jrel)", tspan=(198second, 205second), ylabel="μM/ms", xlabel="Time (s)")
i = (sys.t / 1000, sys.CaMKAct * 100)
plot(sol, idxs=i, title="Active CaMKII", lab="Ctl")
plot!(sol_caf, idxs=i, lab="Caf", ylabel="CaMKII activity (%)", xlabel="Time (s)")
Caffeine and electrophysiology
- Add caffeine in the beginning of the simulation.
- Add caffeine and nifedipine in the beginning of the simulation (nifedipine blocks 90% of L-type calcium channel).
sys = build_neonatal_ecc_sys(simplify=true, reduce_iso=true, reduce_camk=true)
tend = 205second
stimstart = 30second
stimend = 120second
alg = KenCarp47()
@unpack Istim = sys
callback = build_stim_callbacks(Istim, stimend; period=1second, starttime=stimstart)SciMLBase.CallbackSet{Tuple{}, Tuple{SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}, SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}}}((), (SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}(30000:1000:120000, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim)), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim), DiffEqCallbacks.PresetTimeFunction{StepRange{Int64, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}}(30000:1000:120000, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#20#21"{Float64, Float64, Symbolics.Num}(-80.0, 0.5, Istim)), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ()), SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}, DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}(30000.5:1000.0:120000.5, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim)), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim), DiffEqCallbacks.PresetTimeFunction{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, typeof(SciMLBase.INITIALIZE_DEFAULT), CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}}(30000.5:1000.0:120000.5, true, SciMLBase.INITIALIZE_DEFAULT, CaMKIIModel.var"#22#23"{Float64, Float64, Symbolics.Num}(0.0, 0.5, Istim)), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ())))
prob = ODEProblem(sys, [], tend)
prob_caf = ODEProblem(sys, [sys.RyRsensitivity => 10], tend)
gCaL = prob.ps[sys.GCaL]
prob_nif_caf = ODEProblem(sys, [sys.RyRsensitivity => 10, sys.GCaL => 0.1 * gCaL], tend)ODEProblem with uType Vector{Float64} and tType Int64. In-place: true Initialization status: FULLY_DETERMINED Non-trivial mass matrix: false timespan: (0, 205000) u0: 75-element Vector{Float64}: 150952.75035000002 13838.37602 -68.79268 0.12113 0.12113 0.12113 0.12113 0.12113 0.12113 0.12113 ⋮ 0.00977 0.00188 0.09243 0.22156 0.966 0.00702 830.0 830.0 0.0026
ssalg = DynamicSS(alg)
sprob_caf = SteadyStateProblem(prob_caf)
sssol = solve(sprob_caf, ssalg; abstol=1e-10, reltol=1e-10)
@time sol = solve(prob, alg; callback)
@time sol_caf = solve(remake(prob_caf, u0=sssol.u), alg; callback)
@time sol_nif_caf = solve(remake(prob_nif_caf, u0=sssol.u), alg; callback) 0.776985 seconds (63.30 k allocations: 208.387 MiB, 5.83% gc time)
0.860369 seconds (285.97 k allocations: 224.830 MiB, 3.54% gc time, 13.84% compilation time)
0.443770 seconds (46.90 k allocations: 122.731 MiB, 2.51% gc time)
retcode: Success
Interpolation: 3rd order Hermite
t: 2780-element Vector{Float64}:
0.0
0.1842530093836164
0.9882003294241174
4.0712734063454485
11.468877042666529
25.163064632955443
49.586038028951926
83.36346955211738
137.50303891933373
206.18214075191028
⋮
144181.8999375682
149370.41327000043
154558.92660243265
161070.99389492723
169013.16056859927
178032.10045064535
187357.12845898993
196682.1564673345
205000.0
u: 2780-element Vector{Vector{Float64}}:
[152640.04602394556, 12238.844347207129, -67.47370318939181, 0.12238742933122324, 0.12238742928703969, 0.12238742923269894, 0.1223874291679882, 0.12238742909268607, 0.12238742900656217, 0.12238742890937654 … 0.0017638899722423758, 0.001952362450087965, 0.003021926989795491, 0.002307204591585509, 0.9995225437942479, 0.999522529858612, 0.007684751290720895, 57.619305316284475, 56.49698307150776, 0.08793706442329717]
[152640.04602675923, 12238.844344472389, -67.47608385405596, 0.12238150246352823, 0.12238487899828919, 0.12238646228031345, 0.12238710396683712, 0.12238733097534742, 0.12238740218353872, 0.12238742223223555 … 0.0017638903432949167, 0.0019523622796301684, 0.003021791443131714, 0.002307204498140665, 0.9995225437995142, 0.9995225299226123, 0.007684556194472399, 57.61930533073416, 56.49698308631932, 0.08793706437212857]
[152640.0460412888, 12238.844335628484, -67.4863419588738, 0.1223615384582426, 0.12236820385901406, 0.12237352433837904, 0.12237764766092506, 0.12238074283456951, 0.12238299054015957, 0.12238456892777363 … 0.0017639002894421544, 0.001952344562701918, 0.003019913817156804, 0.002307201868610228, 0.9995225439620977, 0.9995225317016053, 0.007681109442438214, 57.619305393782035, 56.49698315094613, 0.0879370641488662]
[152640.04613117408, 12238.844369009661, -67.52375733772861, 0.12231878096368617, 0.12232726067935705, 0.12233499971402545, 0.12234201738668696, 0.12234833829539732, 0.12235399250991266, 0.12235901490648916 … 0.0017640459438024142, 0.001951923697210742, 0.003011514625036301, 0.0023071595625756254, 0.9995225466102455, 0.999522560183886, 0.007661875426907316, 57.619305635604306, 56.4969833987625, 0.0879370633001058]
[152640.04655034712, 12238.84489100134, -67.60266725864268, 0.1222491074415269, 0.1222590594164575, 0.12226851979829308, 0.12227749456196159, 0.12228599082846867, 0.12229401691843596, 0.1223015823424691 … 0.0017648365861398052, 0.0019490721282138946, 0.0029938018661324487, 0.0023068702445830512, 0.9995225647721525, 0.9995227533412281, 0.007620380270868266, 57.619304101796516, 56.49698422891565, 0.08793612653190575]
[152640.04795075575, 12238.847244104358, -67.7148791835132, 0.12214343863327022, 0.12215492461892645, 0.12216604256656051, 0.12217679480255758, 0.12218718403807635, 0.1221972133656129, 0.1222068862540012 … 0.0017667245465694675, 0.001939974898076939, 0.002968746195410057, 0.0023057740579724637, 0.9995226337805907, 0.999523474332407, 0.007561731153973638, 57.6192067202659, 56.497012425406844, 0.0879012449672452]
[152640.05182398023, 12238.85439858778, -67.83680356837378, 0.12197763253750511, 0.121990404262627, 0.12200289319934955, 0.1220150987719831, 0.12202702049826185, 0.12203865798116212, 0.1220500109007837 … 0.001768804431439548, 0.0019217398393596324, 0.002941688101330225, 0.002302645361648815, 0.9995228323121723, 0.9995254812830298, 0.007498489491083083, 57.618450458683725, 56.497428326669294, 0.08768483274019537]
[152640.05872927356, 12238.867184311712, -67.9145349932821, 0.12176796424585437, 0.12178112375242402, 0.12179404269592281, 0.12180671884972011, 0.12181914995295658, 0.12183133370425427, 0.12184326775510433 … 0.0017684333983344995, 0.0019036036640327674, 0.0029245250738852047, 0.002297124747454202, 0.9995231890396084, 0.9995288950919317, 0.007458460520247774, 57.616899372972775, 56.499171522987425, 0.08722104990388099]
[152640.07121698445, 12238.888862142012, -67.95943296168402, 0.12145262459741218, 0.12146536388339307, 0.12147788147770658, 0.12149017449090486, 0.12150223994927498, 0.12151407478954794, 0.12152567585327985 … 0.001764199196622062, 0.001890588364224957, 0.002914671403603113, 0.002287499913165783, 0.9995238338263256, 0.9995345163031448, 0.007435512847136509, 57.61582107315886, 56.50479722352017, 0.08638377123655408]
[152640.0878821591, 12238.91431161045, -67.98524694945176, 0.12107343199180916, 0.12108532476451642, 0.12109700275291573, 0.12110846294866884, 0.12111970224890595, 0.12113071745100744, 0.12114150524707776 … 0.001757303378289173, 0.0018848891864662522, 0.0029090700526981646, 0.0022753872135243816, 0.9995246924404243, 0.99954103152635, 0.00742243650902569, 57.61889479280782, 56.5165398633969, 0.08533600528765778]
⋮
[152693.26001180694, 12185.48959099249, -68.08121163577833, 0.1156954443400513, 0.11569543825540965, 0.11569543218291896, 0.11569542612208146, 0.11569542007234727, 0.11569541403311275, 0.11569540800371819 … 0.001606291931051056, 0.001866398102651513, 0.002888865096970763, 0.0021164059913932384, 0.9989428859588773, 0.999580689265147, 0.007374124200127033, 62.65301516833919, 61.646623035724645, 0.07224338392822158]
[152692.94148220762, 12185.807508555421, -68.08187526627945, 0.11570107513707062, 0.11570107605716645, 0.1157010770247509, 0.1157010780414748, 0.11570107910900512, 0.11570108022902607, 0.11570108140323979 … 0.0016062295439797291, 0.0018663052786664313, 0.002888712549611433, 0.002116221839583373, 0.9994207726430377, 0.9995812295762702, 0.007373656704163249, 62.652458675696856, 61.645939316201584, 0.07225709490663301]
[152692.63002585256, 12186.118979693392, -68.08262780538278, 0.11570351709935411, 0.11570351667281244, 0.11570351631835661, 0.11570351603724614, 0.11570351583080105, 0.11570351570040462, 0.11570351564750629 … 0.0016059679351072206, 0.001866200914687765, 0.002888561122081898, 0.0021159810717152985, 0.9995408991472481, 0.9995811845158319, 0.007373406325148475, 62.65099191716847, 61.64443345089232, 0.07226260561216306]
[152692.2488699209, 12186.500186041383, -68.08345919087942, 0.11570774868745277, 0.11570774925176999, 0.11570774987615821, 0.11570775056191906, 0.11570775131040643, 0.1157077521230291, 0.11570775300125352 … 0.0016056889143966266, 0.001866085919964952, 0.0028883840071190605, 0.002115725677701799, 0.9995739059972749, 0.9995813278251715, 0.007373028548984384, 62.64765994394472, 61.641028299376416, 0.07227200665284517]
[152691.79744968357, 12186.951730803436, -68.08448651690412, 0.11571206600727083, 0.11571206612738491, 0.11571206630850304, 0.1157120665518736, 0.11571206685879545, 0.1157120672306206, 0.11571206766875704 … 0.0016054969563990644, 0.0018659450217724633, 0.0028881631937818925, 0.0021154234133918108, 0.9995804418148716, 0.9995814667714464, 0.007372462098276309, 62.64397091329136, 61.637266756729396, 0.07228161028593971]
[152691.30458436458, 12187.444715802152, -68.0855859291069, 0.11571709813433517, 0.11571709855422545, 0.11571709903148317, 0.11571709956725161, 0.11571710016272127, 0.11571710081913254, 0.11571710153777814 … 0.0016052240705574288, 0.001865793093181328, 0.0028879343959607967, 0.002115099374294442, 0.9995813883988347, 0.9995822273840506, 0.00737214460101173, 62.63973269635586, 61.63294136957926, 0.07229330523873762]
[152690.81450019413, 12187.934862119164, -68.08666137749051, 0.11572237181242624, 0.11572237300519034, 0.11572237425005179, 0.1157223755481652, 0.115722376900734, 0.11572237830901266, 0.1157223797743094 … 0.001604996637039195, 0.0018656451079984047, 0.0028876976129570153, 0.0021147745454111275, 0.999581642306344, 0.9995818264594976, 0.0073714159608775365, 62.635758052625675, 61.628877623064476, 0.07230423327507342]
[152690.34446290907, 12188.40505801524, -68.08774071150177, 0.11572672704469153, 0.11572672750297519, 0.11572672801173878, 0.11572672857210893, 0.1157267291852564, 0.11572672985239814, 0.11572673057480008 … 0.0016046841235744075, 0.0018654943066299724, 0.002887466244914795, 0.002114440465045214, 0.9995816010171785, 0.9995816922358677, 0.007370842669770941, 62.63175302797355, 61.6248010418821, 0.07231367700297381]
[152689.941462889, 12188.808146934349, -68.088640103516, 0.11573074721780516, 0.11573074752004323, 0.11573074787247822, 0.11573074827619276, 0.11573074873231155, 0.11573074924200386, 0.11573074980648565 … 0.001604456481001354, 0.0018653697052767896, 0.0028872740098933813, 0.002114166076709945, 0.9995816977395495, 0.9995816885158855, 0.007370435174449336, 62.6284128429165, 61.62139324907382, 0.07232248450078153]
i = (sys.t / 1000, sys.vm)
tspan = (100second, 102second)
plot(sol, idxs=i, title="Action potential", lab="Ctl"; tspan)
plot!(sol_caf, idxs=i, lab="Caf"; tspan)
plot!(sol_nif_caf, idxs=i, lab="Caf + Nif", tspan=tspan, ylabel="Voltage (mV)", xlabel="Time (s)")
savefig("caf-ap.pdf")"/home/runner/work/camkii-cardiomyocyte-model/camkii-cardiomyocyte-model/docs/caf-ap.pdf"
i = (sys.t / 1000, sys.Cai_sub_SR * 1000)
tspan = (100second, 102second)
plot(sol, idxs=i, title="Calcium transient", lab="Ctl"; tspan)
plot!(sol_caf, idxs=i, lab="Caf"; tspan)
plot!(sol_nif_caf, idxs=i, lab="Caf + Nif", ylabel="Subspace calcium (nM)", xlabel="Time (s)"; tspan)
savefig("caf-cat.pdf")"/home/runner/work/camkii-cardiomyocyte-model/camkii-cardiomyocyte-model/docs/caf-cat.pdf"
i = (sys.t / 1000, sys.CaJSR)
tspan = (100second, 102second)
plot(sol, idxs=i, title="SR Calcium", lab="Ctl", ylabel="SR calcium (μM)"; tspan)
plot!(sol_caf, idxs=i, lab="Caf"; tspan)
plot!(sol_nif_caf, idxs=i, lab="Caf + Nif", ylims=(0, 850), xlabel="Time (s)"; tspan)
CaMKII activities
Simulations
i = (sys.t / 1000, sys.CaMKAct * 100)
plot(sol, idxs=i, title="Simulation", lab="Ctl")
plot!(sol_caf, idxs=i, lab="Caf")
plot!(sol_nif_caf, idxs=i, lab="Caf + Nif", ylabel="CaMKII active fraction (%)", xlabel="Time (s)")
savefig("caf-camkact.pdf")"/home/runner/work/camkii-cardiomyocyte-model/camkii-cardiomyocyte-model/docs/caf-camkact.pdf"
Experiments
chemicaldf = CSV.read(joinpath(@__DIR__, "data/CaMKAR-chemical.csv"), DataFrame)
ts = Dates.value.(chemicaldf[!, "Time"]) ./ 10^9
ctl = chemicaldf[!, "Ctrl Mean"]
ctl_error = chemicaldf[!, "Ctrl SD"] ./ sqrt.(chemicaldf[!, "Ctrl N"])
caf = chemicaldf[!, "caffeine 20mM Mean"]
caf_error = chemicaldf[!, "caffeine 20mM SD"] ./ sqrt.(chemicaldf[!, "caffeine 20mM N"])
plot(ts, ctl, yerr=ctl_error, lab="Control", color=:blue, markerstrokecolor=:blue)
plot!(ts, caf, yerr=caf_error, lab="Caffeine 20mM", color=:red, markerstrokecolor=:red)
plot!(xlabel="Time (s)", ylabel="CaMKII activity (A.U.)", title= "Experiment")
savefig("caf-exp.pdf")"/home/runner/work/camkii-cardiomyocyte-model/camkii-cardiomyocyte-model/docs/caf-exp.pdf"
Runtime information
using InteractiveUtils
InteractiveUtils.versioninfo()Julia Version 1.12.4
Commit 01a2eadb047 (2026-01-06 16:56 UTC)
Build Info:
Official https://julialang.org release
Platform Info:
OS: Linux (x86_64-linux-gnu)
CPU: 4 × AMD EPYC 7763 64-Core Processor
WORD_SIZE: 64
LLVM: libLLVM-18.1.7 (ORCJIT, znver3)
GC: Built with stock GC
Threads: 4 default, 1 interactive, 4 GC (on 4 virtual cores)
Environment:
JULIA_CPU_TARGET = generic;icelake-server,clone_all;znver3,clone_all
JULIA_CONDAPKG_OFFLINE = true
JULIA_CONDAPKG_BACKEND = Null
JULIA_CI = true
LD_LIBRARY_PATH = /opt/hostedtoolcache/Python/3.14.2/x64/lib
JULIA_NUM_THREADS = auto
using Pkg
Pkg.status()Project CaMKIIModel v0.7.0
Status `~/work/camkii-cardiomyocyte-model/camkii-cardiomyocyte-model/Project.toml`
[336ed68f] CSV v0.10.15
[a93c6f00] DataFrames v1.8.1
[459566f4] DiffEqCallbacks v4.11.0
[f6369f11] ForwardDiff v1.3.1
[682c06a0] JSON v1.3.0
[23fbe1c1] Latexify v0.16.10
[98b081ad] Literate v2.21.0
[2fda8390] LsqFit v0.15.1
⌅ [961ee093] ModelingToolkit v10.31.2
[77ba4419] NaNMath v1.1.3
[1dea7af3] OrdinaryDiffEq v6.105.0
[91a5bcdd] Plots v1.41.4
[2913bbd2] StatsBase v0.34.9
[9672c7b4] SteadyStateDiffEq v2.8.0
[ea8e919c] SHA v0.7.0
Info Packages marked with ⌅ have new versions available but compatibility constraints restrict them from upgrading. To see why use `status --outdated`
This notebook was generated using Literate.jl.