DISCRETE AND CONTINUOUS doi:10.3934dcdsb.2012.17.1969 [610098]

DISCRETE AND CONTINUOUS doi:10.3934/dcdsb.2012.17.1969
DYNAMICAL SYSTEMS SERIES B
Volume 17, Number 6, September 2012 pp.1969–1990
GAP JUNCTIONS AND EXCITATION PATTERNS IN
CONTINUUM MODELS OF ISLETS
Pranay Goel
Indian Institute of Science Education and Research
Pune, Maharashtra 411021, India
James Sneyd
Department of Mathematics, University of Auckland
Private Bag 92019, Auckland, New Zealand
Abstract. We extend the development of homogenized models for excitab le
tissues coupled through “doughball” gap junctions. The ana lysis admits non-
linear Fickian fluxes in rather general ways and includes, in particular, calcium-
gated conductance. The theory is motivated by an attempt to u nderstand wave
propagation and failure observed in the pancreatic islets o f Langerhans. We re-
examine, numerically, the role that gap junctional strengt h is generally thought
to play in pattern formation in continuum models of islets.
1.Introduction. Cells in tissues such as in the heart and the pancreatic insulin-
secreting islets of Langerhans communicate almost exclusively via dir ect contact
coupling. Islets are typically composed of about a thousand β-cells that are elec-
trochemically coupled through gap junctions which form between po rtions of mem-
branes of adjacent cells. Gap junctions are important to cellular sig naling between
pancreatic beta-cells. The dominant mechanism of action is via intrac ellular elec-
trical fields that develop across the tissue. Gap junctions do allow s ome ions such
as calcium and small molecules to pass across them as well: in islet physio logy
the molecules of interest are likely to be some glycolytic intermediates and the nu-
cleotidesADPandATP.However,althoughcommunicationacrossth egapjunctions
is electro-chemical, these (slow) diffusion processes are likely to be le ss important
than the electrical conduction effects in regulating islet activity.
The roles that gap junctions play in controlling insulin secretion from is lets are
not yet fully understood. Common techniques of investigating gap j unction connec-
tivity between cells are dye-coupling experiments which estimate the permeability
of gap junctions to dye injected in one cell based on fluorescence m easurements in
its neighbors. However, electrical coupling through gap junctions is not well repre-
sented by such dye diffusion rates. To circumvent this difficulty elect rical conduc-
tance is directly measured between pairs of neighbors using a double two-electrode
voltage/current clamp or by a double whole-cell patch clamp, and es timating gap
junction conductance from an equivalent circuit representation o f the twin-cell sys-
tem (for details see for example [ 13]). This analysis is invariably limiting because
2000Mathematics Subject Classification. Primary: 74Q10, 35B27; Secondary: 92B25.
Key words and phrases. Homogenization, pancreatic islets, gap junctions, excite d activity,
partial wave propagation.
1969

1970 PRANAY GOEL AND JAMES SNEYD
it involves isolated pairs of cells and it is not straightforward to exten d the mea-
surement in the context of an entire islet. It has therefore been in teresting to study
gap junctional conductivity using mathematical modeling of networ ks of coupled
excitable cells that mimic islet activity [ 29,1].
The“strength”ofintercellularcommunicationingeneralisakeypar ameterwhen
studying organized or emergent behavior in tissue. To investigate h ow insulin se-
cretion might synchronize under stimulation Maki and Keizer [ 21] modeled glucose
uptake, metabolism and glucose-dependent insulin secretion in a per ifused bed of
HIT-15 cells in a continuous stirred tank reactor, and demonstrat ed organized in-
sulin oscillations. Later Keener [ 18] revisited this analysis to compare the rates of
flow in the reactor against insulin diffusion between cells and concluded that for re-
alisticestimates ofexperimentalparametersandthe insulin diffusion coefficient (the
communication variable) the model does not show oscillations. Gener ally speaking
though, in terms of communication strengths within an islet, given a s ufficiently
strong gap junctional connectivity modeling studies indicate that b eta-cell oscilla-
tions will all synchronize. This can be taken to mean that gap junctio ns ensure that
even if the delivery of blood induces a gradient of stimulation in the islet (blood
vessels are typically profuse in the islet; although that may not be ca se in diabetes)
the entire islet responds uniformly to stimulation. Although numerou s experiments
have reported islets cells burst in unison, it is not clear if junctional c oupling is
indeed this strong in islets; typical experiments on islet synchrony a re done in vitro
and synchrony may be a function of the perfusion system and unifo rm stimula-
tion. Rocheleau et al. [ 28] have carefully tried to reproduce the heterogeneity in
the microenvironment of an islet: they captured a single islet within a m icrofluidic
chamber and embedded two opposite ends of the islet in two different glucose condi-
tions; one of the glucose concentrations was substimulatory and t he other provided
stimulation for oscillations. They found cases for which, interesting ly, waves prop-
agated from the stimulated side but penetrated only partially into th e islet, that is,
the unstimulated side did not oscillate. The authors claim that since th ey see waves
penetrate the islet partially, and sometimes even observed indepen dent oscillations,
this implies gap junctional communication is weak; in other words, gap junctions
do not necessarily function to synchronize beta-cells in an islet. Ben ninger et al. [ 3]
have investigated this behavior in a model consisting of beta-cells on a lattice with
resistive coupling along edges, they conclude that such wave failure can occur below
a certain percolation threshold arising from a sparse expression of gap junctions in
the matrix.
However, the matter is far from settled. From a theoretical point of view the
choice of modeling framework has been known to influence the propa gation proper-
ties of waves in networks of coupled excitable cells. For a discrete lat tice of Nagumo
cells coupled resistively to their nearest neghbours Keener [ 17] showed that wave
propagationcan fail for sufficiently small coupling strength. Such la ttice differential
equations (LDEs) have now been investigated extensively; it is know n that wave
failure issues depend subtly on the specific nonlinearities in the model cells. For ex-
ample, the authors in [ 14] obtained conditions on the nonlinearities such that wave
propagation failure is guaranteed, while Hupkes et al. [ 15] recently demonstrated a
class of nonlinearities for which propagation failure cannot occur. O n the one hand,
most studies of networks representing islets have typically been do ne using LDE
models. On the other hand, however, theoretical work on LDEs te nds to employ
the Nagumo model as the canonical example, and while these argume nts do focus

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1971
on essentially the nonlinearity associated with bistability in the model c ells it is not
clear if these results extend directly to the more elaborate islet mod els such as those
used in [3].
We remarked that most of the modeling done in islet networks has bee n of the
LDE kind. LDEs themselves can be obtained via a discretization, of a P DE system
for example. However, the LDEs in islet modeling have tended not to b e obtained
in this fashion; for example, the coupling constants between cells in [ 29] were bio-
physically justified against the paired-cell recordings of gap junct ional conductance
measured in [ 25]. An alternate modeling approach towards gap junctions-coupled
tissue is to develop an entirely continuum description. In this framew ork, indi-
vidual conductances between pairs of cells are discarded in favor o f a conductance
tensor varying across the entire tissue. This, in principle, has the t heoretical ad-
vantage of directly engaging experiments on measuring conductan ce across entire
islets, without having to obtain those estimates extrapolated from twin-cell record-
ings. Such an approach was pioneered for cardiac tissues in the 70s : the so-called
“bidomain equations” for syncytial myocardia, which are PDEs that describe evo-
lution on macroscopic scale, i.e. at the scale of the entire tissue, rat her than of
cells. In [ 22,19] the bidomain equations were formally derived using homogeniza-
tion techniques. However, although the general framework has e xisted for some
time now, islet equations have usually been written in LDE format rath er than in
the styleofthe bidomain equations(although see[ 30] foran exampleofa continuum
model). In a series of recent papers Goel and colleagues have inves ted in a contin-
uum framework of equations appropriate for modeling the marcros cale dynamics of
gap junction-coupled excitable cells, islets in particular: in [ 10] we developed the
islet tridomain equations. The islet tridomain equations were themselv es derived
based on our earlier work that obtained homogenized calcium models in [11].
Homogenization techniques have largely been limited to the “syncytia l mem-
brane viewpoint”; in these models gap junctions are simply holes in a co ntinuous
membrane stretching across cells. Another model of gap junction s that may bet-
ter represent the physiology of gap junctional plaques is what we h ave termed the
“doughball gap junction model” in [ 23]: in this view gap junctions are thought of as
thin conductorsbetween cells by meansof which cells are electricallya nd chemically
coupled. Doughball junctions were used in [ 12] to derive macroscopic conductivities
in myocardia for example, while we have computed conductivities more suited to
islet geometry [ 23].
Another aspect of the problem that we are interested in is the role c alcium plays
in gating the conductivity of gap junctions. Several authors have commented that
calcium coupling across gap junctions is of little importance in determin ing syn-
chronicity of oscillations in islets (see, for example, [ 31] for a recent discussion).
Both theoretical as well as experimental studies have focussed o n (electrical) con-
ductivity as the dominant coupling mechanism. Indeed, it is hard to ma ke a case for
a diffusive calcium coupling in the presence of electrical coupling; elect rical coupling
is by far the stronger mechanism of the two in determining the net isle t activity,
to which calcium is then slaved. However, an intriguing possibility emerg es when
calcium is allowedto gate conductance. Calcium is known to be able to ex ert gating
control over gap junction permeability in general [ 24,20], and this may be the case
in islets as well. In general, high calcium – whether in the physiological ra nge or
beyond depends on the cell type – closes gap junctions. (Junction al sensitivity to
calcium, however, varies tremendously in different cells, and is influen ced strongly

1972 PRANAY GOEL AND JAMES SNEYD
by pH as well; little seems to be known about those details in islets.) It is t here-
fore important to obtain how calcium control of the gap junctions e nters into the
homogenized equations. For syncytial models this is a straightforw ard extension of
the usual techniques, as in [ 10] for example. To the best of our knowledge this has
not so far been done in the context of the doughball-junction mode l.
The plan of the paper is as follows. In section 2we extend the derivation of ex-
citable cells coupled through doughball1gap junctions using homogenization tech-
niques [23] to incorporatea calcium-dependent conductive flux acrossjunc tions. We
showthat the macroscopicmodel obtained with the doughballmode l anda calcium-
dependent conductivity follows a reaction-diffusionequation similar t o that used for
cardiac and islet bidomain (syncytial) models:
∇x·(˜σ∇xvm) =β(cm∂vm
∂t+iion), (1)
∂c
∂t=f(c,vm), (2)
wherevmis the membrane potential difference and βis the ratio of the total mem-
brane surface area to the fractional volume occupied by the cells in the periodic
cube, however the homogenized conductivity tensor, ˜ σ, now emerges as a function
of calcium ˜ σ≡˜σ(c). (We shall drop the tilde in subsequent references to σ.) The
main difference between the homogenized syncytial and doughball m odels lies in
the computation of an intermediate “cell problem” which is used to co mpute the
homogenized diffusion and conductivities. At the end of the derivatio n we discuss
an extension to a more general class of nonlinear junctional fluxes using a numerical
example. We then apply the newly equations developed equations to n umerically
study aspects of wave propagation in simple islet models in section 4. In particular
we ask if, in principle, wave failure can be shown to occur in these – con tinuum –
islet equations.
2.Homogenization.
2.1.The doughball gap junction plaque. Most modeling in gap junctions-
coupled tissues takes the syncytial view of cell coupling in which gap j unctions
are modeled as holes in a syncytial membrane; macroscopic equation s of this type,
suitable for modeling the dynamics of islets, were obtained in [ 10] . In reality,
however,a gap junction is made up of thousands of closelyclustere d proteins (called
“connexons”) that rivet together portions of the membranes of adjacent cells; each
connexon has a pore that allows a direct exchange of ions across it. Together, this
forms a plaque that mediates electrical and diffusive coupling throug h the tissue.
In previous work [ 23] we have modeled the plaque as a thin plate with conductivity
and permeability properties and derived the corresponding homoge nized equations.
Figure1shows a schematic depicting the essential geometry of two neighbo ring
beta-cells in the islet separated by gap junctional plaques.
In [23] we considered a constant conductivity, g, of the plaques. In this paper
we focus on the homogenized equations that correspond to a calciu m modulation
of the conductive flux across a gap junction (see below, in particula r Eq. (12)). In
the doughball view: a calcium-dependent closure of the connexon p ores results in a
decreased permeability to calcium, and concomitantly, a decreased the conductivity
of the plaque as well. In other words in the doughball model conduct ivity emerges
1So-called because in this model the islet resembles a stack o f dough balls stuck to each other!

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1973
naturally as a function of calcium concentrations in the cells. We point out that the
corresponding microscopic interpretation in the syncytial model is to assume that
the (effective) size of the junctional pore shrinks with increasing c alcium. It would
be useful to further compare the two models in this sense; we leave this to future
work.
The homogenized model in [ 23] can be viewed as essentially an adaptation of the
modelsderivedin[ 2](fordetailssee, forexample,[ 10]) usingformaltechniques; such
models have since been extensively supported by rigorous proof. T he calculations
below are formal as well, but to the best of our knowledge the homog enization
models described below have yet to be supported with rigorous proo f.
Figure 1. Twoneighboringbeta-cellsintheislet. Γ gisthebound-
ary that represents the gap junctional plaques between cells, an d
Γmis the membrane boundary of the cell(s).
2.2.The microscopic equations. We model the gap junction plate as a thin uni-
form disk of resistivity ρ(conductivity σ= 1/ρ), of thickness l(see also [ 23]). Using
Ohm’s law the current flux jis proportional to electric field across the conductor:
j=σE (3)
=σ∆V
l, (4)
where ∆Vis the voltage jump across the gap junctional plate. Expressing th e gap
junction thickness as some fraction, γ, of the edge length, ε, of the periodic cube
surrounding each cell, the flux therefore scales as
j=g
ε∆V, (5)
withg≡σ/γ. In the sequel we ignore the geometry (thickness) of the junctio nal
plate and simply account for the gap junction using a flux jump condit ion (5) on
the boundary between neighboring cells, Γ gin figure 1.

1974 PRANAY GOEL AND JAMES SNEYD
2.3.The microscopic problem. Letuandvbe two potentials, and c1andc2
be calcium concentrations defined alternately in neighboring cells. Fo r simplicity of
presentation we will ignore the potential in the extracellular space, ve, entirely; for
full detailsofthe correspondingderivationswereferthe readert o similarderivations
in [11,10]. The concentrations satisfy the diffusion equation in the two domain s,
−∇·(Aε∇uε) = 0, (6)
−∇·(Aε∇vε) = 0, (7)
−∇·(Dε
c∇cε
1) = 0, (8)
−∇·(Dε
c∇cε
2) = 0, (9)
in Ω1,2respectively, where AεandDε
cthe conductivity and diffusion tensors respec-
tively.
The boundary conditions on the membrane boundary Γ gare:
−Aε∇uε·η=Aε∇vε·η, (10)
−Aε∇uε·η=1
εg(cε
1,cε
2) (uε−vε), (11)
−Dε
c∇cε
1·η=Dε
c∇cε
2·η= 0. (12)
An example of g(cε
1,cε
2) is
g(cε
1,cε
2)≡g/parenleftbiggcε
1+cε
2
2/parenrightbigg
. (13)
This model of the calcium-dependence of the gap junction is the simp lest possible,
and is only a caricature. In general, one supposes gwill be a complicated function
ofc1andc2;gis presumably symmetrical, but otherwise arbitrary. The analysis
below will reveal that, quite generally, different g(c1,c2) can be studied in much
the same way (see Section 2.3.2). However, for concreteness we will focus on the
model (13).
On the membrane boundary Γ mthe boundary conditions are
−Aε∇uε·η=ε(cm∂(uε−vε
e)
∂t+iion), (14)
−Aε∇vε·η=ε(cm∂(uε−vε
e)
∂t+iion), (15)
−Dε
c∇cε
1·η=ε jCa(cε
1,uε−vε
e), (16)
−Dε
c∇cε
2·η=ε jCa(cε
2,vε−vε
e), (17)
on the membranes of cells 1 and 2, respectively. ηis the outward normal to Ω 1,
jCais the calcium influx current across the cell membrane, and we abuse notation
slightly, using Γ mjointly for the boundaries of the twodomains, Ω 1a well as Ω 2.
We assume the ansatz
uε(x) =u0/parenleftbigg
x,x
ε/parenrightbigg
+ε u1/parenleftbigg
x,x
ε/parenrightbigg
+ε2u2/parenleftbigg
x,x
ε/parenrightbigg
+…, (18)
vε(x) =v0/parenleftbigg
x,x
ε/parenrightbigg
+ε v1/parenleftbigg
x,x
ε/parenrightbigg
+ε2v2/parenleftbigg
x,x
ε/parenrightbigg
+…, (19)
cε
1(x) =c0
1/parenleftbigg
x,x
ε/parenrightbigg
+ε c1
1/parenleftbigg
x,x
ε/parenrightbigg
+ε2c2
1/parenleftbigg
x,x
ε/parenrightbigg
+…, (20)
cε
2(x) =c0
2/parenleftbigg
x,x
ε/parenrightbigg
+ε c1
2/parenleftbigg
x,x
ε/parenrightbigg
+ε2c2
2/parenleftbigg
x,x
ε/parenrightbigg
+…, (21)

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1975
and
∇ → ∇ x+1
ε∇y (22)
wherey=x/ε, we can write the operator Aε=∇·(A(y)∇) in the form
Aε=1
ε2A0+1
εA1+A2, (23)
where
A0=∇y·(A(y)∇y), (24)
A1=∇y·(A(y)∇x)−∇x·(A(y)∇y), (25)
A2=∇x·(A(y)∇x). (26)
We next proceed to match the equations and the boundary conditio ns at the
various orders of ε. For clarity of presentation, in what follows we suppress the
equations for calcium and only show the voltage equations; the deriv ation of the
calcium equations follow similarly.
2.3.1.Equations at lowest order. The equations at order O(1
ε2) in Ω1,2are
−∇y·(A∇yu0) = 0, (27)
−∇y·(A∇yv0) = 0. (28)
together with the boundary conditions at O(1
ε):
−A∇yu0·η=A∇yv0·η, (29)
−A∇yu0·η=g(c0) (u0−v0),on Γg, (30)
and
−A∇yu0·η= 0 (31)
−A∇yv0·η= 0,on Γm. (32)
The only periodic solutions are u0≡u0(x,t) andv0≡v0(x,t). Therefore by
(30)
u0(x,t) =v0(x,t). (33)
2.3.2.Equations at the next order. The equations at order O(1
ε) in Ω1,2are
−∇y·(A∇xu0+A∇yu1) = 0, (34)
−∇y·(A∇xv0+A∇yv1) = 0. (35)
The boundary conditions at O(1) are obtained from O(ε) terms in the Taylor ex-
pansion:
g(cε
1,cε
2) (uε−vε) =g(c0
1,c0
2) (u0−v0)
+ (εu1+ε2u2+…)g(c0
1,c0
2)
−(εv1+ε2v2+…)g(c0
1,c0
2)
+ (εc1
1+ε2c2
1+…)∂g
∂c1(c0
1,c0
2) (u0−v0)
+ (εc1
2+ε2c2
2+…)∂g
∂c2(c0
1,c0
2) (u0−v0)
+….

1976 PRANAY GOEL AND JAMES SNEYD
Inserting this into Equation ( 11) we find that since u0=v0(andc0
1=c0
2), the only
O(ε) term that survives in the expansion is g(c0)(u1−v1), i.e. with conductance
gthat is a function only of xbut not of y. The resulting boundary conditions are
therefore:
−A(∇xu0+∇yu1)·η=A(∇xv0+∇yv1)·η, (36)
−A(∇xu0+∇yu1)·η=g(c0) (u1−v1),on Γg, (37)
and
−A(∇xu0+∇yu1)·η= 0, (38)
−A(∇xu0+∇yu1)·η= 0, on Γm. (39)
(34) and (37) are linear in y. We introduce the solutions χ(y) andξ(y) of the
following system – the cell problem –
−∇·(A(∇yχ+I)) = 0 , (40)
−∇·(A(∇yξ+I)) = 0 , (41)
on Ω1,2, with boundary conditions on Γ g:
−A(∇yχ+I)·η=A(∇yξ+I)·η, (42)
−A(∇yχ+I)·η=g(u0,v0)(χ−ξ), (43)
and
−A(∇yχ+I)·η= 0, (44)
−A(∇yξ+I)·η= 0, (45)
on Γm, whereχ(y) andξ(y) are 1-periodic in y. Notice that twocells are contained
in the unit cube, and are coupled; we shall therefore refer to this a s the “twin-cell
problem”.
Then we can write
u1=χ∇xu0+u1, (46)
v1=ξ∇xv0+v1, (47)
whereu1≡u1(x,t) andv1≡v1(x,t), andu1=v1. The cell problem is coupled in
the variables χandξ;χandξare unique up to a constant, and from Eq. ( 43) this
constant is the same for both χandξ.
Cell problems can similarly be set up for the calcium equations, and it tu rns out
that due to ( 12),˜Dc= 0 .
2.3.3.The macroscopic equations. To get the macroscopic equations we integrate
the equations at order O(1)
−∇x·(A∇xu0+A∇yu1)−∇y·(A∇xu1+A∇yu2) = 0,
−∇x·(A∇xv0+A∇yv1)−∇y·(A∇xv1+A∇yv2) = 0,(48)
on Ω1∪Ω2together with the boundary conditions at O(ε) on Γ gand Γm. On Γ m
the boundary conditions are:
−A(∇xu1+∇yu2)·η=cm∂v0
m
∂t+i0
ion,
−A(∇xv1+∇yv2)·η=cm∂v0
m
∂t+i0
ion, (49)

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1977
The membrane potential v0
m= (u0−v0
e) = (v0−v0
e) because the fields in the two
domains agree to lowest order, and further, i0
ionis dependent on v0
mand involves
gating variables and calcium etc. at lowest order.
WedonotexplicitlycomputetheboundaryconditionsonΓ g: fromthedivergence
theorem it is seen that the contribution to the integral from the se cond terms on
the left hand sides of ( 48) is zero on Γ g.
The contribution to the second terms on the left hand sides of ( 48) from integra-
tion over Γ mis obtained from ( 49). The macroscopic equations (for the model ( 13))
are therefore:
γi∇x·(˜σ∇xv) =λm(cm∂vm
∂t+iion), (50)
∂c
∂t=jCa(c,vm) (51)
wherev≡u0=v0andvm=v−v0
e,γiis the fractional volume occupied by the
cells in the periodic cube and λmis the total membrane surface area/integraltext
Γm.
The homogenized conductivity tensor entering the macroscopic eq uation (51) is
computed from the cell problem solutions as
˜σ=/integraldisplay
Ω1A(∇yχ+I)+/integraldisplay
Ω2A(∇yξ+I). (52)
We note that for clarity of the essential argument we have ignored the dynamics
of the extracellular potential ve. If we include it in the derivation, and use isotropy
of the conductivities of the intra- and extracellular media, we can ob tain a single
equation for the membrane potential difference vm(the so-called “monodomain
reduction”, see [ 10] for example), that is, we recoveran equation ofthe form Eq.( 1).
3.Nonlinear gap junctions. In the previous section we have derived the homog-
enizedequationsforthe caseofacalcium-controlledjunctional flu x: Fickian(linear)
in voltage with a calcium dependence, g(c), modulating the conductivity amplitude
that is quite generally nonlinear (in c). In the macroscopic equations voltage ap-
pears as a “reaction-diffusion” field and calcium enters passively thr ough ODEs.
On the other hand, diffusion of calcium in the islet and its permeability ac ross gap
junctions is an important aspect of the physiology in itself; it is not diffi cult to
develop the equations if the calcium were modeled as a reaction-diffus ion field as
well, since the homogenized diffusion coefficient for calcium can be calcu lated by
the solution of a “twin-cell problem” similar to what we describe for th e voltage
equations.
Note that in either of these scenarios, the nonlinearity of the junc tional flux is
contained in the calcium dependence of conductivity. The equations above were
not derived for a situation in which the amplitude prefixing the Fickian t erm is a
functionof voltageitself(andsimilarquestionsholdforcalciumaswell). Indeed, itis
interestingto askwhat the homogenizedequationswould be when pla que junctional
fluxes are nonlinear in voltage, for example. To the best of our know ledge such
models have not yet been derived from a homogenization procedure . It appears
that the derivation we have just described may continue to hold for such nonlinear
fluxes as well, with the alteration that the cell problem is computed wit h the new
nonlinear flux conditions instead.
Instead of the presenting the corresponding derivation, we next show a numerical
example in 1D that homogenization might hold for reaction-diffusion fie lds with

1978 PRANAY GOEL AND JAMES SNEYD
nonlinear gap junctions. For simplicity of presentation, we concent rate on a single
diffusion field (which we arbitrarilyreferto as calcium) with reactionsa nd nonlinear
gap junctional fluxes on the cell boundaries. We compute the full m icroscopic
solution, which we show agrees well with the homogenized version.
The microscopic problem is
c1t=c1xx−0.1c12(53)
c2t=c2xx−0.1c22(54)
withc1andc2definedalternatelyonlinesegmentsoflength0 .05inthe unitinterval.
Boundary fluxes are equal to
−c1x=c2x=1
εg(c1,c2) (c1−c2) (55)
=1
ε1
1+c1+c2(c1−c2) (56)
on the boundaries on x= 0.05n, n= 1,2,…∈(0,1) where the “cell length”
ε= 0.1, together with c1x=−1 atx= 0 and c2= 0 atx= 1. The calcium-
dependent gin the fluxes above expresses that permeability decreases inverse ly
with the “average” calcium across the gap junction, ( c1+c2)/2. Initial values are
taken as c1= 1 and c2= 1 uniformly.
The corresponding homogenized equation for c≡c0
1=c0
2is
ct=D(c)cxx−c2(57)
withcx=−1 atx= 0 and c= 0 atx= 1 and initial condition c= 1 uniformly.
D(c) is computed from the following twin-cell problem:
−∇·(∇χ+1) = 0 on [0 ,1/2] (58)
−∇·(∇ξ+1) = 0 on [1 /2,1] (59)
with boundary conditions
−(∇χ+1)·η=g(c0
1,c0
2)(χ−ξ) (60)
−(∇ξ+1)·η=g(c0
1,c0
2)(ξ−χ) (61)
solved for the 1-periodic functions χandξ, withg(c0
1,c0
2) = 1/(1+c0
1+c0
2).
In [23] we computed the variation of homogenized conductivity tensors w ith
constant conductance, g, for cell geometries in one and higher dimensions. In one
dimension we showed analytically that σ(normalized to the conductivity of the free
medium) varied with gsigmoidally, σ=g/(2+g); in higher dimensions numerical
simulations suggest the relationship σ=g/(¯g+g) in general, with the constant
¯gvarying with the geometry and dimension. A calcium-dependent homo genized
conductivity can be obtained by substituting these relations into th e corresponding
constant- gproblem: For example, if the conductance of the junctional plate is
g(c) = 1/(1+c) (i.e. the flux is g(c) [v],vbeing the intracellular potential), in the
1Dcasethehomogenizedconductivitywouldbe σ(c) =1/(1+c)
2+1/(1+c). Thehomogenized
diffusion coefficient for this problem is therefore1/(1+2c0
1)
2+1/(1+2c0
1). In other words, the
homogenized diffusion coefficient in Eq. ( 57) is
D(c) =1
4c+3(62)
Figure2shows good agreement between the two solutions with as few as 20
partitions of the domain.

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1979
Figure 2. Concentration profiles for the microscopic model ( 54)
(broken lines) compared to the corresponding solution of the ho-
mogenized equation, Eq. ( 57) (solid lines). Solutions are seen to
increase from t=0.5 to 1 to t=10 towards the steady state solution .
4.Application to wave propagation in islets. We have derived a rather gen-
eralformat ofthe homogenized equationssuitable for studying no t only the effective
conductivity of islet tissue, but also its calcium dependence. We next wish to use
these equations to understand how conductivity affects spatiote mporal patterns
such as synchrony and wave propagation in the islet. Recall from th e Introduction
that an important consideration is to examine the issue of the stren gth of gap junc-
tions: discrete models indicate that if the junctional conductivity is large then we
expect synchrony, but if the conductivity is low then wave failure ca n possibly take
place; we would like to ask if such behavior can be seen in the continuum equations.
We first comment on the current state of modeling of beta-cells or is lets. Modern
models of glucose-stimulated islet bursting are now far improved bey ond the simple
Fitzhugh-Nagumo or Morris-Lecar style of modeling that were repr esented in the
original Chay-Keizer models. We believe that the Dual Oscillator Mode l (DOM) [ 4]
is the current state-of-the-art in islet modelling today. We have re cently presented
a detailed analysis dissecting the properties of the DOM using geomet ric slow-
fast arguments [ 9]. In that paper we show that the DOM naturally “splits” into
components, an electrical one that is composed of voltage and calc ium dynamics
and another, a metabolic oscillator that drives the dynamics along th e bifurcation
diagram of the electrical subsystem. However, the DOM is rather c omplex and not
yet well understood in terms of coupled networks of beta-cells (ev en with LDEs).
It is also very expensive to compute with. Note that the electrical p ortion of the
DOM oscillator is directly influenced by electrical coupling but not the g lycolytic
component (for details see [ 4]). For a first pass at the problem, therefore, we choose
touseMorrisLecardynamicstorepresentisletbursting. Weleavet hemorecomplex
considerations of the DOM for future work.

1980 PRANAY GOEL AND JAMES SNEYD
ML1ML2
VK-84-84
VL-60-60
VCa120120
¯gK88
¯gL22
¯gCa44.4
V1-1.2-1.2
V21818
V3122
V417.4 30
φ0.0667 0.04
C2020
Table 1. Parameter values for the Morris Lecar models, ML1 and
ML2, based on [ 8].
For simplicity we work on a finite interval [0 ,100µm] with no-flux boundary
conditions: this is in keeping with the experimental design in [ 28] in which the
essential dynamics is more or less along a single direction between the two glucose
concentrations at opposite ends. We stimulate the islet with a glucos e differential
between two polar ends of the islet by we applying a (linear) gradient o fIappacross
the interval, with Iapp=Imaxatx= 0 and Iapp= 0 atx= 100.
We investigate synchrony and waves, in particular as parameterize d the homog-
enized conductivity of the tissue.
4.1.Model parameters. Morris Lecar dynamics have frequently been used to
modelburstenvelopes,see[ 31]forarecentexample. TheMorrisLecarcellequations
we will use below are:
CdV
dt=σ∂2V
∂x2−¯gCam∞(V)(V−VCa)−¯gKc(V−VK)−¯gL(V−VL)+I(63)
dc
dt=φ[c∞(V)−c]
τc(V)(64)
with
m∞(V) = 0.5/bracketleftbigg
1+tanh/parenleftbiggV−V1
V2/parenrightbigg/bracketrightbigg
(65)
c∞(V) = 0.5/bracketleftbigg
1+tanh/parenleftbiggV−V3
V4/parenrightbigg/bracketrightbigg
(66)
τc(V) = 1/cosh/parenleftbiggV−V3
2V4/parenrightbigg
(67)
In absence ofthe laplacianterm (i.e. in the ODE models) the twoparam etersets,
models ML1 and ML2, correspond to Type 1 and Type 2 excitability. We recall
(see, for example, [ 27,8]) that Type 1 and Type 2 models both have the following
general features: at very low and very high stimuli the cell is at equ ilibrium, and
oscillatory in between. The major difference in excitability in the two mo dels is: in
Type 1 models the loss from equilibria to oscillations (near the low stimulu s side) is
via a saddle-node on an invariant circle (i.e. a SNIC bifurcation), while f or a Type

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1981
2 model it is through a (subcritical) Hopf bifurcation. Near the high s timulus side,
the loss from oscillations to equilibria is through a (subcrticial) Hopf bif urcation in
either model.
4.2.Synchrony and waves in islets. In the Rocheleau experiment an islet was
captured in a microfluidic chamber between two different glucose con centrations at
opposite ends. This setup is probably somewhat artificial as blood ve ssels typi-
cally innervate islets profusely. Nonetheless, this is a useful exper imental strategy
to characterize coupling strength in the islet. Further, the exact nature of micro-
circulation in the islet is not yet fully understood: on the one hand, ex periments
imaging the cytoarchitecture of human islets have revealed beta ce lls are situated
along vascular innervations of the islet [ 6] which indicates that stimulation is likely
to be uniform in the islet; on the other hand, calcium recordings in [ 6] did not
find the typical synchronized beta-cell oscillatory activity that is u sually thought
to be important to pulsatile insulin secretion. Additionally, opposing vie ws seem to
exist regarding whether a direction in blood flow establishes an differe ntial cellular
perfusion order within an islet [ 5,16], or not [ 6]. It is therefore useful to study
theoretical models in which one side of the islet senses a higher conce ntration of
glucose than the other. In our stylized view of the islet, and in keepin g with the
Rocheleau experiment, we impose a linear glucose gradient across th e interval.
Fig.3shows the effect of varying σin the islet. These are computed with ML2
dynamics with glucose drive Imax=200. At high σ, Fig.3(a), the entire islet is
synchronized. As σis lowered to 1 ×103and 1×102a noticeable phase difference
can be seen between the two ends of the islet. At σ= 1×102note that while
the frequencies of oscillation are the same across the islet the duty cycle near the
low glucose end is significantly smaller than at the high glucose end. At σ= 10
wave propagation is slow enough to result in anti-synchronous oscilla tions between
the two ends of the islet, and a sharp drop in amplitude is apparent ne ar the low
glucose end as well.
In generalone assumes that “high”conductivity values are those that correspond
tosynchrony; forwavestopassacrossthe isletconductivitymus tbe lowerthanthat.
4.3.Low conductivity can explain wave failure. At even lower conductivity,
σ= 1, the wave fails to propagate across the islet completely, Fig. 4.
It is interesting to note from the figure that the wave velocity is not uniform
throughout: The wave is seen to slow down as it propagates into the islet. The
amplitude decreases steadily as the wave penetrates the islet, and about four-fifths
of the way in it decays completely to rest.
4.4.Model dependence of partial wave propagation. We comment that sim-
ilar results, especially wave propagation failure, can be obtained with ML1 as well.
Fromgeneralphysicalconsiderationsit doesnot seem toosurpris ingthat synchrony,
waves, and propagation failure can be obtained in either model ML1 o r ML2; how-
ever, a detailed analysis of the instability that leads, in particular, to wave failure,
merits further investigation which we leave for future work. We rem ark here that
propagationfailure studiesin LDEs aredone with homogeneousdrive , in contrastto
thepresentsetup, inwhichthedriveisnonuniform. Itwouldtheref orebeinteresting
to compare wave failure in LDEs with heterogeneous drive. Another consideration
that remains to be checked is whether similar behavior can also be est ablished with
the more recent comprehensive models such as the dual-oscillator m odel in [4] that
include metabolic components as well is as yet an open question.

1982 PRANAY GOEL AND JAMES SNEYD
050100150200250300350400−60−40−2002040
(a)σ= 1×105.050100150200250300350400−60−40−2002040
(b)σ= 1×103.
050 100 150 200250300350400−60 −40 −20 020 40
x=100 x=0
(c)σ= 1×102.
(d)σ= 10.
Figure 3. Synchrony and traveling wave solutions in the islet for
varying values of σ. The graphs show voltage against time plotted
at three points: the two endpoints, and at x= 50.
4.5.Oscillator death: The role of driving stimulus and conducti vity.We
next examine the role of conductivity relative to the driving stimulus in synchrony
and wave propagation in the islet.
Figure5is computed with ML1 at lowerglucose driverelative to the synchrono us
case in Sec. 4.2,Imax=50. At the lower value of conductivity σ= 1 we find partially
propagatingwavesasbefore. Thisiscouldleadtotheassumptionth atalargervalue
ofσmight therefore recover synchrony. That is, if larger gap junctio nal strengths
lead to synchrony then a (low) glucose concentration that is sufficie nt to produce
excitation at lower σought to result in a fully synchronous state at high junctional
strength. This is a naive argument, as shown in Panel (b): if conduc tivity is raised
toσ= 1×102and higher, wave propagation ceases completely; and a standing
wave of subthreshold voltage is seen instead.
It seems this might be due to the phenomenon of oscillator death: in g eneral
it is well known that for oscillating systems with diffusive coupling, a suffi ciently
large coupling strength can sometimes lead to oscillations being quenc hed in the
system [26,7]. Thus the increasein junctional strength does not immediately imply ,
at least in the context of the reaction-diffusion models described he re, a guarantee
of synchronized oscillations.

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1983
Figure 4. Partially penetrating waves at σ= 1.
050 100 150 200250300350400−60 −40 −20 020 40
x=100 x=0; x=5
(a)σ= 1.050100150200250300350400−60−40−2002040
(b)σ= 1×102.
Figure 5. Oscillator death at stimulus drive, Imax=50. The plots
are displayed at three points: the two endpoints, and at x= 5.
Oscillator death is seen for a conductivity σ= 1×102and higher,
while at the lower σ= 1 partially propagating waves are recovered.
We point out in theoretical studies with coupled (modified) Morris Lec ar islet
(ODE) models Tsaneva-Atanasova et al. [ 31] noted that a large diffusivity of cal-
cium through gap junctions can result in oscillator death. In [ 31] oscillator death

1984 PRANAY GOEL AND JAMES SNEYD
occurs via a pitchfork (or saddle-node, in the case of heterogene ous oscillators) bi-
furcation, which the authors speculate is a generic feature of rela xation oscillators
coupled through the slow variable. In contrast, in our simulations we find that
even in the absence of calcium diffusion, electrical coupling alone can le ad to par-
tially propagating waves at low gap junction conductivity and to oscilla tor death at
high conductivity. A similar behavior can also be seen with ML2 (simulatio ns not
shown).
4.6.Calcium-gated gap junctions and recovery from quenched osc illa-
tions.So far we have not explicitly invoked a calcium dependence of conduct ivity
in our simulations. We next consider interesting examples of role that calcium can
play if it modulates islet conductivity.
A calcium modulation of conductivity likely to be most effective in islets in h igh
glucose conditions. We speculate that since islet oscillations are thou ght to be
important to insulin secretion, a calcium modulation of gap junctions m ight be
required precisely to ensure islet pulsatility. We show with numerical e xamples
an islet that is saturated in high glucose conditions is able to restore p ulsatility
by modulating conductivity in response to elevated calcium concentr ations. Thus,
calcium gating of gap junctions in islets may effectively extend the phy siological
range over which the islet is able to respond to high glucose with oscillat ions and
secretion.
The previous section showed that when the high-end glucose stimulu s is weak,
oscillatordeath can result from a high conductivity: the unstimulate d regionsof the
islet couple strongly to the stimulated regions and overwhelm them. A n alternate
scheme of quenched oscillations is the following scenario: the islet is un able to
oscillate if the high-end glucose is so strong, i.e. at suprathreshold le vels, that beta-
cells close to this end do not oscillate at this concentration but are sa turated (i.e.
at steady-state high calcium and voltage levels). Figure 6, panel (a), shows voltage
saturating across an ML2 islet for which Imax=450; both, the high-end glucose as
well as the conductivity, are high here, σ= 1×104. Notice that the entire islet
has seized up, due to a suprathresold stimulus inducing oscillator dea th from high
conductivity as described above; this can be seen by noting that if t he conductivity
is lowered to σ= 4×102oscillations can be recovered, as seen in the top right
panel of the figure. In this figure the conductivity is treated as a p arameter, and
not influenced by calcium; we next demonstrate two examples of calc ium-dependent
conductivity rules that result in this recovery behavior.
In Figure 7we allow conductivity to vary with calcium according to
σ(c) = 100+2×103
1+exp(( c−0.55)/0.1). (68)
Ifσis held constant at the maximal value, 2 .1×103, this islet has stationary
behavior, i.e. oscillations are eliminated (simulations not shown). For t he islet
with a calcium-dependent σ, it is seen that the conductivity adjusts to calcium,
becoming low when calcium is high and high when calcium is low; the net resu lt is
that oscillations are possible in the islet. Notice that the calcium is slight ly higher
near the high-glucose end than the opposite end; this is in keeping wit h a lower
conductivity at the high glucose end.

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1985
0 50 100 150 200 250 300 350 400−60−40−20020406080
(a) v at σ= 1×104.0 50 100 150 200 250 300 350 40000.10.20.30.40.50.60.70.8
(b) c at σ= 1×104.
050100150200250300350400−60−40−20020406080
(c) v at σ= 4×102.05010015020025030035040000.10.20.30.40.50.60.70.8
(d) c at σ= 4×102.
Figure 6. Steady-state and travelling wave solutions in the islet
for varying values of σ. The top panels show transmembrane volt-
age and the lower panels the corresponding calcium. The plots are
displayed at three points: the two endpoints, and at half way in
between.
It is also possible that gap junctions do not respond instantaneous ly to calcium;
in that case, we simulate a time-dependent change in σas:
σ∞(c) = 100+3×103
1+exp(( c−0.55)/0.02),
dσ
dt= (σ∞(c)−σ)/τσ, (69)
whereτσ= 40 and σ(x,t= 0) = 3 .1×103. As before, if σis held constant at
the maximal value, 3 .1×103, oscillations are eliminated (simulations not shown).
The oscillations resulting from a σthat varies as Eq. ( 69) are shown in Figure 8.
The recovered oscillations are very similar between Figs. 7and8, although it is
interesting to note the sharper sawtooth waveform of σin Fig.8.
Thus, an islet in a hyperglycemic environment without a calcium modulat ion of
conductivity, is unable to respond with adequate oscillations and sec retion. How-
ever, once we assume islet gap junctions are sensitive to a [Ca] ithreshold and

1986 PRANAY GOEL AND JAMES SNEYD
(a) v
x=100 x=0
(b) c
x=0 x=100
(c)σ
Figure 7. Calcium-dependent conductivity computed from
Eqs. (68) shows synchronized voltage and calcium oscillations in
the islet, top panels; conductance responds instantaneously to lo –
cal calcium and hence oscillates as well, bottom panel. The plots
are displayed at three points: the two endpoints, and at half way
in between.
respond by shutting down gap junctions, we find that with a plastic c onductivity
the islet can potentially oscillate. We note that there are probably ad aptive mech-
anisms in the cell that can enable tuning of the junctional paramete rs to a suitable
calcium threshold that allows for oscillations. However, although ada ptation rules
that enable such parameter control are hotly investigated theor etically, they are not
well understood in general.
An islet that is able to monitora high glucose condition viaelevated calciu m, and
in turn closes gap junctions and lowers the effective conductivity of the tissue, can
thus rescue itself from oscillator death and enhance the effective r ange over which
it is able to secrete.
5.Summary. In the present study we describe fresh progress made in address ing
the long-standing, and invariably difficult, problem of islet synchrony and wave
propagation (failure), and the attendant questions of the stren gth of gap junctions
ingeneral. Wehaveintroducedacompletelynewframeworkofcontin uumequations
suitable for modeling this problem.

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1987
(a) v
x=0
x=100
(b) c
x=0 x=100
(c)σ
Figure 8. Calcium-dependent conductivity computed from
Eqs. (69) evolves to generate synchronized voltage and calcium
oscillations in the islet, top panels. Notice that conductance oscil-
lates aswell, and oscillationamplitude varieswith position, bottom
panel. The plots are displayed at three points: the two endpoints,
and at half way in between.
In general we have been interested in extending homogenization de scriptions of
tissues electrochemically coupled through gap junctions; the pres ent work has fea-
tures that make this framework best representative of islet mode ling. Firstly, we
work with the doughball gap junction. Historically, homogenization is such tis-
sues has largely focussed upon the syncytial model; in syncytial mo dels “coupling
strength” is a function of geometry, i.e. pore size determines the e ffective electri-
cal conductivity (and permeability to ions such as calcium) of the tiss ue. In the
doughball model of cellular coupling each cell is partitioned from its ne ighbor by
gap junctional plaques, and intracellular fluxes (both electrical an d chemical) are
determined by the conductivity and permeability across this bounda ry. In other
words, junctional “strength” measures these properties of th e plaques boundaries
(and in principle these properties can be nonlinear). Thus this type o f junction is
better suited to accommodate gating laws of greater sophisticatio n than the simple
pore of syncytial models.

1988 PRANAY GOEL AND JAMES SNEYD
Secondly, the homogenized equations we derive represent a large c lass of models
in which an associated concentration field can modulate a linear condu ction flux
across the gap junction in quite general (i.e. nonlinear) ways. In pr evious work [ 23]
we studied how the macroscopic electrical and chemical equations d epended on the
strength, g, of Fickian fluxes. The effective conductivity and permeability tenso rs
varied with dimension and geometry, and are parameterized by g. Here we have
extended the analysis beyond the constant gcase. We are specifically interested in
two situations:
i. Electrical flux (across the gap junction) is of the type g(c)[v] , i.e. calcium-
dependent. In general, g(c) can be an arbitrary nonlinear function of calcium.
The calcium permeability is taken to be zero for simplicity, but poses no real
restriction. We obtain the complete derivation in this case.
ii. This derivation is also attractive in that it appears it may extend to n onlin-
earities in the constitutive laws for junctional conductance as well (Section 3).
Moreover, if electrical coupling is as described above and in addition, gap junc-
tions are calcium-permeable: Linear fluxes (i.e. p[c]) are easily handled via
twin-cell problems as earlier; and we speculate that nonlinear fluxes , where
p≡p(c), can be handled in much the same way, i.e. with the boundary condi-
tions used to solve the twin-cell problems then appearing nonlinearly . We give
one numerical example that confirms this by comparing an exact micr oscopic
simulation to its homogenized version; the solution to the nonlinear (t win-cell)
problem is facilitated by first solving the linear problem as a function of p(as
is done numerically in [ 23]), then rescaling through the nonlinearity p(c).
We point out that while these features are in general also true for c ardiac models,
we do not discuss differences in how calcium engages voltage respect ively in the two
cases.
We have next numerically studied wave propagation in pancreatic islet s as a
key motivation of the analysis. The homogenized equations themselv es have been
derived in a general setting to enable them to accommodate a large c lass of voltage,
calcium and even metabolite (islet) dynamics. The Morris-Lecar burs ting models
used for the simulations are admittedly simplistic; nonetheless they a re a useful
guide to anchor ideas regarding the fuller problem represented in, s ay, the DOM.
Further, while we have derived the equations anticipating different s cenarios of
increasing complexity, we restrict ourselves in the numerics to the s implest case, of
studying how σparameterizes the behavior of the islet: This is, in part, because
very little seems to be known of the nonlinearities associated with gap junction
dynamics in islets, but also because from the physiological point of vie w we are
largely interested in the question of the “strength” of islet coupling . We leave
a fuller discussion of the role of junctional nonlinearities, and assoc iated calcium
dynamics for example, for later work.
Wave propagation is a function of the excitable properties of individu al cells and
the driving stimulus, but also of the effective conductivity of the tiss ue. Indeed, the
case of wave failure as it penetrates the islet indicates that the con ductivity tensor
in islets may be several orders of magnitude smaller than previously r eported for
the entire pancreas, and needs to be investigated further, both theoretically and
experimentally.
Acknowledgments. We areindebted to AvnerFriedmanwhofirst helped usdown
this road; the scenery has had us occupied for numerous years no w.

GAP JUNCTIONS AND EXCITATION IN HOMOGENIZED ISLETS 1989
REFERENCES
[1]Gerda de Vries and Arthur Sherman, Beyond synchronization: Modulatory and emergent
effects of coupling in square-wave bursting , in “Bursting: The Genesis of Rhythm in the
Nervous System” (eds. Stephen Coombes and Paul C. Bressloff) , World Scientific, London,
(2005), 243–272.
[2] J.-L. Auriault and H. I. Ene, Macroscopic modelling of heat transfer in composites with
interfacial thermal barrier , Int. J. Heat Mass Transfer, 37(1994), 2885–2892.
[3] R. K. Benninger, M. Zhang, W. S. Head, L. S. Satin and D. W. P iston,Gap junction coupling
and calcium waves in the pancreatic islet , Biophys. J., 95(2008), 5048–5061.
[4] R. Bertram, L. Satin, M. Zhang, P. Smolen and A. Sherman, Calcium and glycolysis mediate
multiple bursting modes in pancreatic islets , Biophys. J., 87(2004), 3074–3087.
[5] F. C. Brunicardi, J. Stagner, S. Bonner-Weir, H. Wayland , R. Kleinman, E. Livingston,
P. Guth, M. Menger, R. McCuskey, M. Intaglietta, A. Charles, S. Ashley, A. Cheung, E. Ipp,
S. Gilman, T. Howard and E. Passaro, Microcirculation of the islets of Langerhans, Long
Beach Veterans Administration Regional Medical Education Center Symposium, Diabetes,
45(1996), 385–392.
[6] O. Cabrera, D. M. Berman, N. S. Kenyon, C. Ricordi, P. O. Be rggren and A. Caicedo, The
unique cytoarchitecture of human pancreatic islets has imp lications for islet cell function ,
Proc. Natl. Acad. Sci. U.S.A., 103(2006), 2334–2339.
[7] B. Ermentrout and K. Bar-Eli, Oscillation death , Scholarpedia, 3(2008), 5371.
[8]G. B. Ermentrout and D. H. Terman, “Mathematical Foundation s of Neuroscience,” Interdis-
ciplinary Applied Mathematics, 35, Springer, New York, 2010.
[9]P. Goel and A. Sherman, The geometry of bursting in the dual oscillator model of panc reatic
beta-cells , SIAM Journal on Applied Dynamical Systems, 8(2009), 1664–1693.
[10]P. Goel, A. Sherman and A. Friedman, Multiscale modeling of electrical and intracellular
activity in the pancreas: The islet tridomain equations , Multiscale Modeling & Simulation, 7
(2009), 1609–1642.
[11]P. Goel, J. Sneyd and A. Friedman, Homogenization of the cell cytoplasm: The calcium
bidomain equations , Multiscale Modeling & Simulation, 5(2006), 1045–1062.
[12]P. E. Hand, B. E. Griffith and C. S. Peskin, Deriving macroscopic myocardial conductivities
by homogenization of microscopic models , Bull. Math. Biol., 71(2009), 1707–1726.
[13] E. Hertzberg, ed., “Gap Junctions,” Advances in Molecu lar and Cell Biology, Elsevier Science,
2000.
[14]A. Hoffman and J. Mallet-Paret, Universality of crystallographic pinning , Journal of Dynam-
ics and Differential Equations, 22(2010), 79–119.
[15]H. J. Hupkes, D. Pelinovsky and B. Sandstede, Propagation failure in the discrete nagumo
equation , Proceedings of the American Mathematical Society, 139(2011), 3537–3551.
[16] T. Kanno, S. O. Gopel, P. Rorsman and M. Wakui, Cellular function in multicellular system
for hormone-secretion: Electrophysiological aspect of st udies on alpha-, beta- and delta-cells
of the pancreatic islet , Neurosci. Res., 42(2002), 79–90.
[17]J. P. Keener, Propagation and its failure in coupled systems of discrete e xcitable cells , SIAM
Journal on Applied Mathematics, 47(1987), 556–572.
[18] J. P. Keener, Diffusion induced oscillatory insulin secretion , Bull. Math. Biol., 63(2001),
625–641.
[19] W. Krassowska and J. C. Neu, Effective boundary conditions for syncytial tissues , IEEE
Trans Biomed Eng, 41(1994), 143–150.
[20] A. Lazrak and C. Peracchia, Gap junction gating sensitivity to physiological internal calcium
regardless of pH in Novikoff hepatoma cells , Biophys. J., 65(1993), 2002–2012.
[21] L. W. Maki and J. Keizer, Mathematical analysis of a proposed mechanism for oscillat ory
insulin secretion in perifused HIT-15 cells , Bull. Math. Biol., 57(1995), 569–591.
[22] J. C. Neu and W. Krassowska, Homogenization of syncytial tissues , Crit. Rev. Biomed. Eng.,
21(1993), 137–199.
[23] G. P. and S. J., A comparison of effective conductivity in two models of gap ju nction coupling
in tissues, , submitted.
[24] C. Peracchia, Chemical gating of gap junction channels; roles of calcium, pH and calmodulin ,
Biochim. Biophys. Acta, 1662(2004), 61–80.

1990 PRANAY GOEL AND JAMES SNEYD
[25] M. Perez-Armendariz, C. Roy, D. C. Spray and M. V. Bennet t,Biophysical properties of
gap junctions between freshly dispersed pairs of mouse panc reatic beta cells , Biophys. J., 59
(1991), 76–92.
[26]A. Pikovsky, M. Rosenblum and J. Kurths, “Synchronization: A Universal Concept in Non-
linear Sciences,” Cambridge Nonlinear Science Series, 12, Cambridge University Press, Cam-
bridge, 2001.
[27] J. Rinzel and G. B. Ermentrout, Analysis of neural excitability and oscillations , “Methods in
Neuronal Modeling,” MIT Press, Cambridge, MA, USA, (1989), 135–169.
[28] J. V. Rocheleau, G. M. Walker, W. S. Head, O. P. McGuinnes s and D. W. Piston, Microfluidic
glucose stimulation reveals limited coordination of intra cellular Ca2+ activity oscillations in
pancreatic islets , Proc. Natl. Acad. Sci. U.S.A., 101(2004), 12899–12903.
[29] A. Sherman and J. Rinzel, Model for synchronization of pancreatic beta-cells by gap j unction
coupling , Biophys. J., 59(1991), 547–559.
[30] C. L. Stokes and J. Rinzel, Diffusion of extracellular K+ can synchronize bursting osci llations
in a model islet of Langerhans , Biophys. J., 65(1993), 597–607.
[31] K. Tsaneva-Atanasova, C. L. Zimliki, R. Bertram and A. S herman, Diffusion of calcium and
metabolites in pancreatic islets: Killing oscillations wi th a pitchfork , Biophys. J., 90(2006),
3434–3446.
Received October 2011; revised March 2012.
E-mail address :pgoel@iiserpune.ac.in
E-mail address :sneyd@math.auckland.ac.nz