Computational Cardiac Modeling

1
Computational Cardiac Modeling Angelina
Altshuler and Iwen Wu August 25, 2006 PRIME 2006
Final Report Dr. David Abramson Monash
University Melbourne, Australia Dr. Andrew
McCulloch, Dr. Anushka Michailova, and Dr. Roy
Kerckhoffs University of California, San Diego

2
Background of Cardiac Modeling

• Developing and studying heart models is an
  important tool to understand health issues
• Ischemia
• arrhythmias
• implantation of pacemakers
• heart attacks
• By developing models that accurately model
  cellular level mechanisms, eventually entire
  organ models can be developed

3
Development of Cardiac Models

• Mathematical models have been developed from
  experimental data
• ODEs
• Initial conditions
• Ion movement in single cells (endo, mid, and epi
  myocytes)
• Increasing the complexity
• Add more channel to the models
• Add more metabolically important molecules
• Increase the amount of computations

Luo and Rudy model 1994
Michailova and McCulloch model, 2005
4
Connecting Cardiac Modeling
pacemaker
Model Implementation
Physiological Experiments and Data
Clinical and Engineering Applications
Computational Models
Physiological Validation
Cardiac Muscle
5
Ten Tusscher et al. ionic model of human
ventricular myocytes

• Mathematical model of action potentials in human
  ventricular cells that accurately reflects
  experimentally observed electrophysiological data
  
• Can reproduce electrophysiological behavior in
  Epicardial, Endocardial, and Midmyocardial cells

K. H. W. J. ten Tusscher et. Al., A Model for
Human Ventricular Tissue, AJP-Heart Circulation
Physiology. April 2004. 286H1573-H1589.
6
Coupled transmural heterogeneity of calcium
activity and its mechanical correlates

• Cordeiro et al. demonstrated the differences in
  cell shortening and Ca2 transients along with
  mechanical function across the three ventricular
  cell types in the canine left ventricle

John Cordeiro et. Al., Transmural heterogeneity
of calcium activity and mechanical function in
the canine left ventricle, AJ-Heart Circulation
Physiology, April 2004. 2861471-1479.
7
Coupled transmural heterogeneity of calcium
activity and its mechanical correlates

• Dr. Roy Kerckhoffs adapted the ten Tusscher et
  al. model to include the transmural heterogeneity
  of the Calcium transient in human ventricular
  myocytes
• The heterogeneous model was stabilized for
  calcium, but other ionic concentrations (Na and
  K) became unstable

Epicardial cell with Vmaxup 0.000425, Brel0.24,
Crel0.0086
8
Endocardial cell with Vmaxup0.0003, Brel0.24,
Crel0.0086
Midmyocardial cell with Vmaxup0.0003, Brel0.24,
Crel0.0086
9
Project Goals

• Develop more stable computational model for three
  different ventricular cell types
• Modify parameters Vmaxup, Brel, and Crel to
  stabilize Cai, CaSR, Nai, Ki concentrations for
  epicardial, endocardial, and midmyocardial cells
• Run parameter sweeps in Nimrod G to generate
  graphs at different parameter combinations
• Use the Simplex optimization algorithm in Nimrod
  O to find ideal set of parameters to stabilize
  the model

10
Nimrod G for parameter sweeps

• Perform parameter sweep over a specified range
  for parameters Vmaxup, Brel, Crel
• Vmaxup is the maximal pump current taking up
  calcium in the SR. Vmaxup from .0003-.000475,
  increment by .000025
• Brel is the CaSR half saturation constant of the
  calcium-induced calcium release (CICR) current.
  Brel from .2-.28, increment by .02
• Crel is the maximal CaSR-independent CICR
  current. Crel from .008-.009, increment by .0002
• Simulation from 0-50s at frequency of 1Hz
• Ran separate experiments for epicardial,
  endocardial, and midmyocardial cells
• Evaluate stability of ionic concentrations based
  on the graphs generated

11
Nimrod O for optimization

• Utilizes optimization algorithms to obtain best
  set of parameters to stabilize graphs
• First need a numeric value for stability
• Took regression of last ten periods in 50 sec
  experiment to get relative stability
• Sum squared regression values for all four
  concentrations and have Nimrod O try to minimize
  to get a value as close to zero as possible

Set limits for parameters and Nimrod O chooses
combinations within the parameter space to use,
going in the direction that gives the smallest
summed regression value
12
Nimrod O Results
iwenw_at_valdore epicardial50o2 cat
nimrodo.runtime.error.log grep best sgt
(0.000425, 0.28, 0.008232) cost 1.37101e-11 best
sgt Reflected point better than best sgt
(0.000405, 0.27, 0.007932) cost 1.22961e-11 best
sgt (0.000405, 0.27, 0.007932) cost 1.22961e-11
best sgt (0.000405, 0.27, 0.007932) cost
1.22961e-11 best sgt (0.000405, 0.27, 0.007932)
cost 1.22961e-11 best sgt (0.000405, 0.27,
0.007932) cost 1.22961e-11 best iwenw_at_valdore
endocardial50o cat nimrodo.runtime.error.log
grep best sgt (0.000425, 0.28, 0.008232) cost
1.92313e-11 best sgt Reflected point between
second worst and best sgt (0.000425, 0.28,
0.008232) cost 1.92313e-11 best sgt (0.00044,
0.275, 0.008332) cost 1.53008e-11 best sgt
Reflected point between second worst and best sgt
(0.00044, 0.275, 0.008332) cost 1.53008e-11 best
sgt (0.00044, 0.275, 0.008332) cost 1.53008e-11
best sgt (0.00044, 0.275, 0.008332) cost
1.53008e-11 best
iwenw_at_valdore midmyocardial50o cat
nimrodo.runtime.error.log grep best sgt
(0.000425, 0.25, 0.008232) cost 3.19246e-11 best
sgt (0.00044, 0.265, 0.008382) cost 2.68118e-11
best sgt Reflected point between second worst and
best sgt (0.00044, 0.265, 0.008382) cost
2.68118e-11 best sgt Reflected point between
second worst and best sgt (0.00044, 0.265,
0.008382) cost 2.68118e-11 best sgt Reflected
point between second worst and best sgt (0.00044,
0.265, 0.008382) cost 2.68118e-11 best sgt
(0.00044, 0.265, 0.008382) cost 2.68118e-11 best
sgt Reflected point between second worst and best
sgt (0.00044, 0.265, 0.008382) cost 2.68118e-11
best sgt (0.000439028, 0.268194, 0.0082945) cost
2.59326e-11 best sgt (0.000439028, 0.268194,
0.0082945) cost 2.59326e-11 best sgt Reflected
point better than best sgt (0.000442654,
0.266636, 0.00833756) cost 2.58487e-11 best sgt
(0.000442654, 0.266636, 0.00833756) cost
2.58487e-11 best
13
Visualization of results
CaSR- slope 2
Cai- slope 1
Nai- slope 3
Ki- slope 4
Special thanks to Donny Kurniawan
14
Visualization of results

• With two inputs x y
• Level curves- curves joining all the points with
  the same output value (contour map)
• Generalize to 3D- x,y,z
• Level surfaces- surface joining all the points
  with the same output value
• Make curves different colors and transparent to
  show different surfaces
• Example- Ki (slope 4)
• Themed h- Vmaxup, w- Crel, l- Brel
• Stable low bottom left corner

15
Project Plan

• The UCSD MATLAB version of the LabHEART model
  (Puglisi and Bers, Am J Physiol, 2001), which
  calculates ionic currents, Ca-transients, and
  action potentials in epicardial, mid-myocardial,
  and endocardial rabbit ventricular myocytes has
  been found unstable at high frequencies over long
  periods of time
• Use Nimrod/G to sweep over a large parameter
  space to identify which model parameters to
  adjust in order to improve the model stability
  and relate these results to the experimental data
  among the three cell sub-types at high frequencies

16
Puglisi and Bers Model
17
Cardiac Modeling with LabHEART

• Parameters Used to Stabilize Model
• G_Cab (maximum I_Cab conductance)
• G_Nab (maximum I_Nab conductance)
• Run the code from 1-10 hertz increasing the
  frequency 1 hertz every 30 seconds for a total
  simulated time of 5 minutes
• Run experiments in endocardial and epicardial
  myocyte cells because there is comparable data
• Focus on the calcium transient and action
  potential calculations for each cell to determine
  stability of the model
• Check if model calculates alternans at high
  frequencies (9 and 10 hertz) by extracting
  maximum peak values of action potential and
  calcium transient calculations

18
What are Alternans?

• Calcium Transient
• Instability in intracellular calcium cycling
  which causes action potential duration alternans
  by affecting the calcium sensitive membrane
  currents
• Action Potential
• Gating kinetics that create a steep action
  potential duration curve

Spatially Discordant Alternans in Cardiac Tissue
Role of Calcium Cycling. Circ. Res.
200699520-527.
19
Nimrod/G Distributed Computing for Parametric
Modeling

• Distributes jobs to many nodes
• Parallel processing of jobs
• Performs large parameter sweeps
• Advantages to using Nimrod With Cardiac Modeling
• Faster since jobs run in parallel
• Starts experiments automatically
• Automatically runs all parameter combinations

20
The Advantage of Using Nimrod/G for Testing
Cardiac Models

• Parameters G_Cab.001-.006 with step .0005 and
  G_Nab.00005-.00015 with step .00001 creates 110
  possible parameter combinations
• Without Nimrod/G
• Must reset each parameter combination after the
  previous experiment finishes
• With Nimrod/G
• Set parameters once and run multiple jobs at a
  time
• When one job finishes, queued jobs begin to run

21
Stabilizing the Endocardial Calcium Transient and
Action Potential Calculations Using Nimrod/G

• Endocardial rabbit myocyte simulation (G_Cab.001
  and G_Nab.00005)
• The model is stable and physiologically correct
  from 1-7 hertz
• At 8 hertz (210 seconds), the calcium transient
  values do not increase enough to be
  physiologically stable
• Although the model is not physiologically stable
  at 9-10 hertz, the model is able to predict
  calcium and action potential alternans

22
Conclusions for Endocardial Rabbit Myocyte
Simulations

• Values of G_Cab0.001 and G_Nab0.00005 produces
  the most stable (1-7 hertz), physiologically
  accurate graphs
• Alternans values at high frequencies (maximum
  peaks extracted from the added MATLAB code
• Maximum values shown are taken from the middle
  calculations at 9 and 10 hertz since values tend
  to slightly increase with time

23
Stabilizing the Endocardial Calcium Transient and
Action Potential Calculations from 1-10 Hertz

• epicardial myocyte simulation (G_Cab.0005 and
  G_Nab.00007)
• The model is stable and physiologically correct
  from 1-8 hertz
• At 9 hertz (240 seconds), the calcium
  transients do not increase enough to be
  physiologically stable, yet the action potential
  calculations are still follow the decreasing
  trend
• The model predicts calcium and action potential
  alternans at 9 and 10 hertz

24
Conclusions for Epicardial Rabbit Myocyte
Simulations

• Values of G_Cab.0005 and G_Nab.00007
• Alternans are calculated at high frequencies
  (maximum peaks extracted from the added MATLAB
  code
• Maximum values shown are taken from the middle
  calculations at 9 and 10 hertz since values tend
  to slightly increase with time

25
Alternans at 10 Hertz in the Action Potential
Calculations for Epicardial Cells
G_Cab.0005 and G_Nab.00007
26
Further Work

• Investigate the alternans observed at high
  frequencies in more detail
• Run more experiments
• Change different parameters
• Test the mechanisms regulating EC-coupling in
  rabbit myocytes at high rates
• Splice the code to make the jobs smaller but
  dependent upon one another
• Slow nodes will only run fractions of the project
  
• Data is already stored so an experiment could be
  picked up in the middle if there were a crash
• Jobs can move around the nodes to decrease the
  time spent running a long job on a slow node

27
Acknowledgements

• Monash University
• Dr. David Abramson
• Colin Enticott, Slavisa Garic, Tom Peachy
• University of California, San Diego
• Dr. Andrew McCulloch
• Dr. Anushka Michailova
• Dr. Roy Kerckhoffs
• PRIME Coordinators
• Peter Arzberger
• Linda Feldman
• Teri Simas
• Gabrielle Wienhausen
• The PRIME program was made possible by
• University of California San Diego
• National Science Foundation
• The California Institute for Telecommunications
  and Information Technology

28
Australia
29
Sydney and New Zealand