Modeling elec odi usion in ca diac in e cala ed disc nanos uc u es
Ena I ano ic1, Ma ina Vannuci1,2, Nicolae Moise3, Se h H. Weinbe g3 and Jan P. Kuce a1
1 Depa men o Physiology, Uni e si y o Be n, Be n, Swi ze land
2 G adua e School o Cellula and Biomedical Sciences, Uni e si y o Be n, Be n, Swi ze land
3 Depa men o Biomedical Enginee ing, The Ohio S a e Uni e si y, Columbus, OH, USA
This is he MATLAB sou ce code pe mi ing o eplica e he esul s p esen ed in ou a icle en i led
”Elec odi usion in ca diac in e cala ed disc nanos uc u es al e s cell-cell ac ion po en ial
ansmission ia ephap ic coupling: a model s udy” ha was accep ed o publica ion in he Jou nal o
Physiology“.
Download he ZIP ile “CodeAndMeshes.zip” and decomp ess all iles in o a dedica ed olde . Make
su e ha all iles a e loca ed in he same olde . De ails ega ding he ope a ion o he sc ip s and
unc ions a e p esen ed b ie ly below. Mo e de ails a e p o ided as commen s a he beginning o he
sc ip s/ unc ions as well as h oughou he code.
The code was de eloped wi h he 2019a e sion o MATLAB and should also un wi h newe e sions.
Use in e ms o he C ea i e Commons CC-BY license (A ibu ion 4.0 In e na ional).
Main sc ip s wi h hei auxilia y unc ions
Main_MATLAB.m
Main_COMSOL.m
Main_OSU.m
The esul s shown in Figu es 1, 2, 4, 5 and 6 can be eplica ed using Main_MATLAB.m. To eplica e
Figu e 3, use Main_COMSOL.m. To ob ain he esul s wi h he o uous mesh shown in Figu e 7, use
Main_OSU.m. O no e, he di e ence be ween hese iles esides mainly in he gene a ion o he ini e
elemen mesh. Main_MATLAB.m gene a es he mesh wi hin he code, while Main_COMSOL.m and
Main_OSU.m load an ex e nal ile (MeshCOMSOL.ma and MeshOSU.ma , espec i ely). The mesh
da a in MeshCOMSOL.ma we e gene a ed using COMSOL Mul iphysics ( e sion 6.3) while he mesh
da a in MeshOSU.ma we e gene a ed using dedica ed code a Ohio S a e Uni e si y (OSU).
A he beginning o each sc ip , a numbe o pa ame e s and lags a e de ined (as MATLAB a iables,
see Table 1) and he use mus adjus hem acco dingly as desc ibed below.
Table 1: Va iables (no e ha no all a e a ailable in he h ee di e en main sc ip s)
Va iable
Possible Values
De ini ion
De aul
alue
NaC
0, 1
Na+ channel dis ibu ion: uni o m (0) o
clus e ed (1)
1
PATTERN
1, 2, 3, 4, 5, 6
Na+ channel clus e dis ibu ion pa e n (in he
la in e cala ed disc)
1: Small cen al clus e
2: La ge cen al clus e
3: 5 cen al clus e s (3 d ow o Figu e 6)
4: 7 clus e s (4 d ow o Figu e 6)
5: 13 clus e s (5 h ow o Figu e 6)
6: 3 g oups o 5 clus e s (6
h
ow o Figu e 6)
1
Cle Wid h_nm
Be ween 10 and
100
Cle wid h in nm ou side he pe inexus o gap
junc ion plaque
30
gap ac o
Be ween 0 and 1
ac ion o gap junc ional coupling ela i e o
no mal (e.g., 0.01 o 1%)
0
Dscaling
1
Scaling ac o o simula e educed di usion due
o ex acellula p o eins in he im a ound Na+
channel clus e s
1
DYNCONC
0, 1
Cons an (0) o dynamically changing (1) ion
concen a ions
1
BLKOUTINA
0, 1
Ou wa d Na+ cu en : ac i e (0) o blocked (1)
0
GJp
0, 1
Gap junc ion dis ibu ion: uni o m (0) o
clus e ed (1)
0
Pe iWid h_nm*
Be ween 10 and
Cle Wid h_nm
Pe inexus wid h in nm
10
To ep oduce he simula ion esul s o each Figu e, se hese a iables as ollows:
Figu e 1: Open Main_MATLAB.m. Fo a uni o m Na+ channel dis ibu ion (le column o Figu e 1), se
NaC=0. To eg oup all Na+ channels in he cen e o he in e cala ed disc (middle and igh columns o
Figu e 1), se NaC=1. Choose he adius o he clus e by se ing PATTERN=2 (middle column) o
PATTERN=1 ( igh column, de aul ). Lea e he emaining a iables a hei de aul alues
(BLKOUTINA=0, Dscaling=1). To eplica e he esul s shown in panel B, se Cle Wid h_nm=40 and
gap ac o =0. Ob ain he black aces wi h DYNCONC=0 and he ed aces wi h DYNCONC=1. A he
end o he simula ion, he delay be ween he ups okes is shown in he MATLAB command window.
To ob ain he indi idual da a poin s used o gene a e panel C, di e en simula ions can be un wi h
di e en alues o Cle Wid h_nm and gap ac o .
Figu e 2: Open Main_MATLAB.m. Se NaC=1 and PATTERN=1 o ob ain a mesh wi h a small cen al
Na+ channel clus e . Mo eo e , se Dscaling=1, GJp=0, gap ac o =0 and Cle Wid h_nm=50. Simila o
Figu e 1, se DYNCONC=0 and DYNCONC=1 o esul s wi h cons an and dynamic ion concen a ions.
To ob ain simula ion esul s whe e ou wa d INa is blocked se BLKOUTINA=1. A he end o he
simula ion, he delay be ween he ups okes is shown in he MATLAB command window. To ob ain
he indi idual da a poin s used o gene a e panel C, di e en simula ions can be un wi h di e en
alues o Cle Wid h_nm and gap ac o .
Figu e 3: Open Main_COMSOL.m. To eplica e he esul s in panel B se gap ac o =0 and
Cle Wid h_nm=50. Le column: DYNCONC=0 and BLKOUTINA=0. Middle column: DYNCONC=1 and
BLKOUTINA=0. Righ column: DYNCONC=1 and BLKOUTINA=1. Va y Pe iWid h_nm be ween 10 and 50
nm in s eps o 10 nm o ep oduce he di e en cu es. The simula ions analyzed in panel C can be
ep oduced by se ing Cle _wid h_nm o alues be ween 10 and 100 nm (in s eps o 10 nm),
Pe iWid h_nm o alues be ween 10 and Cle _wid h_nm (in s eps o 10 nm) and se ing DYNCONC
and BLKOUTINA acco dingly.
Figu e 4: Open Main_MATLAB.m. Se NaC=1 and PATTERN = 1 o ob ain a mesh wi h a small cen al
Na+ channel clus e . To eplica e he esul s in panel B se gap ac o =0, Cle Wid h_nm=50 and a y
Dscaling om 1 (de aul ) o 0.2 in s eps o 0.1. Se DYNCONC and BLKOUTINA as in Figu e 3. The
simula ions analyzed in panel C can be ep oduced by se ing Cle _wid h_nm o alues be ween 10
and 100 nm (in s eps o 10 nm), Dscaling o alues be ween 0.2 and 1 (in s eps o 0.1), and se ing
DYNCONC and BLKOUTINA acco dingly.
Figu e 5: Open Main_MATLAB.m. Se NaC=1 and PATTERN=4 o ob ain a mesh wi h mul iple Na+
channel clus e s. To eplica e he esul s, se Cle Wid h_nm=30, gap ac o =0.01, DYNCONC=1
BLKOUTINA=0 and a y Dscaling om 1 (de aul ) o 0.2 in s eps o 0.1.
Figu e 6: Open Main_MATLAB.m. Choose a Na+ channel dis ibu ion pa e n (co esponding o he
di e en ows o Figu e 6) and se PATTERN acco dingly (see Table 1). Se DYNCONC=1 and
BLKOUTINA=0. To ob ain he schema ics shown in panel B, se gap ac o = 0 and o ob ain he esul s
in panel C, se gap ac o =0.01. Repea he simula ions o di e en alues o Cle Wid h_nm ( om 10
o 100 nm in s eps o 10 nm) and Dscaling ( om 0.2 o 1 in s eps o 0.1) o econs uc he da a. A he
end o each simula ion, he delay be ween he ups okes is shown in he MATLAB command window.
Figu e 7: Open Main_MATLAB.m o eplica e he esul s o a la in e cala ed disc mesh wi h uni o m
Na+ channel dis ibu ion (NaC=0, Mesh 1). Open Main_OSU.m o eplica e he esul s o he o uous
in e cala ed disc mesh, as ollows. Mesh 2: se NaC=0, GJp=0 and Dscaling=1. Mesh 3: se NaC=1,
GJp=0 and Dscaling=1. Mesh 4: se NaC=1, GJp=1 and Dscaling=1. Mesh 5: se NaC=1, GJp=1 and
Dscaling=0.2, 0.3 o 0.4. To ob ain he esul s depic ed in panel B, se gap ac o =0.01, Cle Wid h_nm
= 20, DYNCONC=1 and BLKOUTINA=0.
Then, execu e he sc ip . The p og am displays he mesh (and he bounda y nodes) in a dedica ed
igu e. The p og am also displays he memb ane po en ial o he in e cala ed disc o cell 1 and cell 2,
he sodium cu en densi ies o bo h cell memb anes, he ex acellula po en ial (in mV) in he cle
as well as ion concen a ions in he cle (ion species: sodium, po assium, calcium) as a colo map. A
he end o he simula ion, he p og am p oduces a igu e simila o he ones p esen ed in he a icle.
Fu he mo e, he p og am ou pu s in he MATLAB command p omp he delay be ween he wo
ups okes o he in acellula po en ials. The sc ip Main_OSU.m also gene a es mo ies, as de ailed
below.
Mo ie S1: Open Main_OSU.m. Se he pa ame e s as ollows: NaC=0, GJp=0, Dscaling=1,
gap ac o =0.01, Cle Wid h_nm = 20, DYNCONC=1 and BLKOUTINA=0 ( hese a e he pa ame e s o
Mesh 2). Execu e he sc ip , which by de aul will gene a e a mo ie ile “Mo ie.mpg” in he wo king
di ec o y.
Mo ie S2: Open Main_OSU.m. Se he pa ame e s as ollows: NaC=1, GJp=1 and Dscaling=0.2,
gap ac o =0.01, Cle Wid h_nm = 20, DYNCONC=1 and BLKOUTINA=0 ( hese a e he pa ame e s o
Mesh 5). Execu e he sc ip , which by de aul will gene a e a mo ie ile “Mo ie.mpg” in he wo king
di ec o y.
The unc ions used by he h ee main sc ip s a e b ie ly desc ibed below. Fo u he de ails, see he
commen s in he co esponding MATLAB code.
Func ions ela ed o he ini e elemen me hod used in his wo k
The pa ial di e en ial equa ions we e disc e ized in space using he ini e elemen me hod wi h
iangula elemen s and piecewise linea basis unc ions, esul ing in sys ems o o dina y di e en ial
equa ions.
_p_lin i.m
Gi en a mesh o iangles, his unc ion compu es se e al auxilia y a iables and a ays (e.g., iangle
a eas, con ibu ion o indi idual elemen s ( iangles) o he s i ness/mass ma ices and load ec o s)
o be used subsequen ly o assemble he a ious ele an s i ness and mass ma ices as well as he
load ec o s used in he sys ems o o dina y di e en ial equa ions. These a ays and a iables a e
packed as ields in o a s uc u e “S” e u ned by he unc ion. These a ays a e o ma ed in a manne
o op imize he gene a ion (in e ms o memo y and speed) o he ele an ma ices in spa se o ma
(using he MATLAB unc ion “spa se”). These ma ices and load ec o s a e hen compu ed om he
s uc u e “S” by he ollowing unc ions.
_k_lin i_0.m
Gene a es s i ness (conduc ance) ma ix (K, in spa se o ma ) om he s uc u e “S” ou pu by
“ _p_lin i.m”, assuming ha he s i ness pa ame e (e.g., ex acellula conduc ance) is cons an
wi hin e e y iangle.
_m_lin i_0.m
Gene a es mass (capaci ance) ma ix (M, in spa se o ma ) om he s uc u e “S” ou pu by
“ _p_lin i.m”, assuming ha he s i ness pa ame e (e.g., memb ane capaci ance) is cons an wi hin
e e y iangle.
_ _lin i_0.m
Gene a es load ec o (F) om he s uc u e “S” ou pu by “ _p_lin i.m”, assuming ha he load
pa ame e (e.g., sodium cu en conduc ance pe uni a ea) is cons an wi hin e e y iangle.
_k_lin i_1.m
Gene a es s i ness (conduc ance) ma ix (K, in spa se o ma ) om he s uc u e “S” ou pu by
“ _p_lin i.m”, assuming ha he s i ness pa ame e (e.g., ex acellula conduc ance) is a linea
unc ion wi hin e e y iangle.
_k_lin i_1p .m
Va ia ion o _k_lin i_1.m. Gene a es s i ness (conduc ance) ma ix (K, in spa se o ma ) om he
s uc u e “S” ou pu by “ _p_lin i.m”, assuming ha he s i ness pa ame e (e.g., ex acellula
conduc ance) is a linea unc ion wi hin e e y iangle and he e is ano he s i ness pa ame e ha is
uni o m o e each iangle.
_m_lin i_1.m
Gene a es mass (capaci ance) ma ix (M, in spa se o ma ) om he s uc u e “S” ou pu by
“ _p_lin i.m”, assuming ha he s i ness pa ame e (e.g., memb ane capaci ance) is a linea unc ion
wi hin e e y iangle.
_ _lin i_1.m
Gene a es load ec o (F) om he s uc u e “S” ou pu by “ _p_lin i.m”, assuming ha he load
pa ame e (e.g., sodium cu en conduc ance pe uni a ea) is a linea unc ion wi hin e e y iangle.
Func ion used o dila e a meshed objec
_Dila e.m
Func ion used o iden i y he su ounding o objec s in a wo-dimensional ini e elemen mesh o
iangles, e.g., he pe iphe al im o "dense ex acellula p o eins" wi h a "wid h" o N iangles (de aul
alue: 2 iangles) a ound Na+ channel clus e s.
Func ion used o modeling he Na+ cu en
mhjLR1_li shi z.m
This unc ion compu es and e u ns ga ing pa ame e s o he h ee ga es m, h and j o he sodium
cu en acco ding o Luo and Rudy (Ci c Res 1991, 68:1501-1526) wi h modi ica ions by Li shi z and
Rudy (Biophys J 2009, 97:1265-1276). This unc ion was al eady used in p e ious p ojec s and uploaded
p e iously on Zenodo (I ano ic and Kuce a, J Physiol 2021, 599: 4779-4811, DOI:
10.5281/zenodo.5226268; Cells 2022, 11, 3477, DOI: 10.5281/zenodo.7271839).
