Fibrillatory-Like Arrhythmias

14 Pages·2005·4.12 MB·English

Biophysical Journal Volume 70 February 1996 581-594 581
Wavelet Formation in Excitable Cardiac Tissue: The Role of
Wavefront-Obstacle Interactions in Initiating High-Frequency
Fibrillatory-Like Arrhythmias
Josef M. Starobin, Yuri 1. Zilberter, Elizabeth M. Rusnak, and C. Frank Starmer
Departments of Medicine (Cardiology) and ComputerScience, Duke University Medical Center, Durham, North Carolina27710 USA
ABSTRACT High-frequency arrhythmias leading to fibrillation are often associated with the presence of inhomogeneities
(obstacles) in cardiac tissue and reduced excitability of cardiac cells. Studies of antiarrhythmic drugs in patients surviving
myocardial infarction revealed an increased rate of sudden cardiac death compared with untreated patients. These drugs
blockthecardiacsodium channel, therebyreducing excitability, which mayalterwavefront-obstacle interactions. Indiseased
atrial tissue, excitability is reduced by diminished sodium channel availability secondary to depolarized rest potentials and
cellular decoupling secondary to intercellular fibrosis. Excitability can also be reduced by incomplete recovery between
successive excitations. In all of these cases, wavefront-obstacle interactions in a poorly excitable medium may reflect an
arrhythmogenic processthatpermitsformation ofreentrantwavelets leadingtoflutter, fibrillation, andsudden cardiacdeath.
To probetherelationship between excitabilityand arrhythmogenesis, weexplored conditionsfornewwaveletformation after
collision of a plane wave with an obstacle in an otherwise homogeneous excitable medium. Formulating our approach in
terms ofthe balance between chargeavailable inthewavefrontand theexcitation chargerequirements ofadjacent medium,
wefound analyticallythecritical medium parametersthatdefined conditionsforwavefront-obstacle separation. Underthese
conditions, when a parent wavefront collided with a primitive obstacle, the resultant fragments separated from the obstacle
boundaries, subsequently curled, and spawned new "daughter" wavelets. We identified spatial arrangements of obstacles
such that wavefront-obstacle collisions leading to spawning of new wavelets could produce high-frequency wavelet trains
similarto fibrillation-like arrhythmias.
INTRODUCTION
High-frequency arrhythmias leading to fibrillation and sud- Konings et al., 1994). The complex interaction ofinhomo-
dencardiac deathareoften associated withinhomogeneities geneities and reduced excitability highlighted in these ex-
(obstacles) in cardiac tissue and reduced excitability of amples with the nonlinear dynamics of cardiac cells may
individual cells. However, the link between obstacles, re- reflectarrhythmogenic conditionsthatpromoteformationof
duced excitability, and arrhythmogenesis remains unclear. wavelets leading to reentrant arrhythmias, fibrillation, and
Obstacles and reduced excitability can arise from several sudden cardiac death.
different mechanisms. Studies of antiarrhythmic drugs in Reentrant activation (spiral wave) is observed in many
patients surviving myocardial infarction revealed an in- excitable preparations (Winfree, 1987, 1989; Zykov, 1988;
creased rate of sudden cardiac death compared with un- Davidenko, 1995). Variations of electrocardiographic po-
treated patients (CAST Investigators, 1989). Prior myocar- tentials (ECG) associated with reentrant cardiac tachycar-
dialinfarction mostprobably leftaresidueof"obstacles" or dias are synchronous with the rotation of a single spiral
poorly excitable myocardial regions, whereas the drugs (Krinsky, 1984; Markin et al., 1981; Pertsov et al., 1993).
utilized inthis study reducedexcitability byblockadeofthe Fibrillation, on the otherhand, often reflects multiple reen-
cardiac sodium channel. Studies ofatrial cells frompatients trant wavelets (Allessie et al., 1994; Konings et al., 1994;
with atrial disease revealed reduced excitability secondary Cox et al., 1994), leading to high-frequency oscillation of
to depolarized rest potentials (Ten Eick and Singer, 1979), the ECG potential that is not synchronized with rotation of
whereas studies of human diseased atrial tissue demon- a single spiral.
strated anincreased microfibrosis thatwouldhypothetically The initial work of Balakhovsky (1965) and Krinsky
lead to a loss of side-to-side cellular connections, thereby (1966) demonstrated that a temporary obstacle, caused by
increasing anisotropy, reducing excitability, and increasing prolonged refractoriness in a localized region of excitable
the size of inexcitable regions (Spach and Dolber, 1986).
medium, was adequatetorupture anincidentwavefrontin a
Excitability can also be attenuated by incomplete recovery
manner that led to spiral wave development. Pertsov and
associated withhighrates ofexcitation (Schalij etal., 1992;
co-workers (1983) extendedtheseresults innumerical stud-
ies of the collisions between a planar wavefront and a
permanent obstacle in an otherwise homogeneous medium.
Receivedforpublication 12April1995andinfinalform20October1995.
As the medium excitability was reduced (as determined by
Addressreprintrequests toC. FrankStarmer, Box 3181,DukeUniversity
Medical Center, Durham, NC 27710. Tel.: [phone]; Fax: 919-684- a model parameter similar to the Na conductance), they
8666; E-mail: [email protected]. found that above a critical excitability, the wave fragments
C 1996by the Biophysical Society remained "attached" to the obstacle borders and immedi-
0006-3495/96/02/581/14 $2.00 ately reattached proximal to the obstacle. Wave fragments
582 Biophysical Journal Volume 70 February 1996
derived from a wavefront colliding with an obstacle when (1993), and Agladze et al. (1994) all point to a critical role
the medium excitability was less than this critical value played by wavefront-obstacle separation in new wavelet
separated from the obstacle. In these studies, wavebreaks formation. To date, the relationship between medium pa-
followed by wavefront-obstacle separation could act as a rameters and conditions for wavefront-obstacle separation
sourceofnew spiral wavelets, as describedbyBalakhovsky hasbeenlimitedtomodelingdifferenttransienttime-depen-
(1965) and Krinsky (1966). dent regimes ofmedium excitability. In addition, there are
Kogan andcolleagues (1990) furtherexploredtherole of no studies that indicate how high-frequency arrhythmias
wavefront-obstacle separation in a series of numerical ex- associated with multiple wavelets might arise from a single
periments. As in the Krinsky model (1966), they showed wavefront propagating in a low-excitability medium.
thatiftheobstacle was "transient" (e.g.,physicallyremoved To explore these issues, we have defined three subprob-
afterwavefrontcollision), thentheresidual fragments could lems that lead from fractionation ofa wavefront to a high-
evolve into counter-rotating spiral wavefronts. frequency arrhythmia:
Using high-frequency stimulation to incrementally re- 1) How to spawn anew waveletafterwavefront-obstacle
duce medium excitability, Panfilov and Keener (1993) ex- separation;
tended the Pertsov et al. (1983, 1990) results and showed 2) How to insert the newly spawned wavelet behind the
that obstacle-wavefront separation and spiral wave devel- parent wavefront;
opment displayed "use-dependent" properties. More specif- 3) How to multiply wavelet formation.
ically, if the stimulation rate was fixed such that the inter- We explored multiple wavelet formation by first investi-
stimulus interval was less than the time necessary for full gating the behavior ofthe boundary layer between a single
recovery of excitability, they found that no separation oc- obstacle and a colliding wavefront as determined by the
cuffed with the first stimulus (wavefront propagated into balance of charge available in the wavefront and the exci-
rested, fully excitable medium); separation but no spiral tation charge requirements ofadjacent medium. Using this
development occurred with the next few successor stimuli approach we found analytical expressions relating critical
(due to the accumulated loss of excitability); and if the values of medium parameters to conditions for wavefront-
stimulation frequency was high enough, further successor obstacle separation. We approachedthequestionofwavelet
wavefronts fractionated, detached from the object, and formation and transformation ofa single wavelet into mul-
evolved into counter-rotating spiral waves (due to the ac- tiple wavelets leading to fibrillatory-like arrhythmias by
cumulated loss ofexcitability greater than a critical value). exploring spatial arrangements ofprimitive obstacles in an
In all ofthese studies, only single orpaired spiral processes otherwise homogeneous medium.
were studied.
The observation ofmultiple propagating wavelets during
fibrillation (Allessie et al., 1994, Konings et al., 1994; MODEL
Schuessler et al., 1995) confirmed an early hypothesis of
The normal cardiac action potential duration is approxi-
Moe (1962). Moeandco-workers (1959, 1962, 1964, 1965)
mately0.3 s, andtheexcitationwavepropagationvelocityis
proposed the multiple wavelet hypothesis, which explained approximately 0.5 m s-1. Such rapid propagation makes
these high-frequency arrhythmias as a result of repeated
multiple reentrant wavelets unlikely because the action po-
wavefront fractionation secondary to inhomogeneities in
tential wavelength (approximately 15 cm) can be compara-
refractory properties. This hypothesis postulated that high-
ble to the spatial extent ofthe heart. However, when excit-
frequency arrhythmias (e.g., fibrillation) were maintained
ability is diminished as a result ofreduced sodium channel
by multiple interfering wavelets thatmigrated around inho-
availability secondary to Na channel blockade, incomplete
mogeneous regions such as islets or strands of refractory
recovery from Na channel inactivation secondary to rapid
tissue. Although Moe's studies (Moe et al., 1964) demon-
reexcitation, orreducedcellularcoupling, thewavelength is
strated that multiple wavelets could be produced by certain
reduced so that the spatial requirements forreentry become
combinations of stimulus timing and distributed refractory
more probable. In some cases the spatial wavelength may
properties, the mechanism ofwavelet formation and main-
decrease to several millimeters (Spach et al., 1988).
tenancewasuncertain. Morespecifically, afterruptureofan
For our numerical and analytical investigations, we used
incident wavefront after collision with either an obstacle or
an excitable membrane formulation of the 2D FitzHugh-
a refractory region, the conditions that promoted or pre-
Nagumo (FitzHugh, 1961) model, describing wave motion,
vented the formation of independent wavelets were not
which is given by
explored.
Observations ofhigh-frequency arrhythmias inamedium au a82u a2u\
of reduced excitability (Leier et al., 1978; Buxton et al., Cmat =3x2+ ay2)
1984; Yamashita et al., 1994), combined with medium
inhomogeneities, may provide a substrate for formation of aV (OKU
multiple wavelets leading tohigh-frequency arrhythmias. In at TK
our opinion, the studies of Krinsky (1966), Pertsov et al. Here (Na('rK mimic the maximum sodium and potassium
(1983, 1990), Kogan et al. (1990), Panfilov and Keener conductances, respectively, andTKis a time constant of a
Starobin et al. Wavelet Formation in Cardiac Tissue 583
slow recovery process, INa = oNa(u) represents an inward like nature of the piecewise linear obstacle also approxi-
excitation current, andIK = Vrepresents an outward recov- mates micro-fibrotic regions observed in diseased or aged
ery current. To facilitate analyses, we used this system in atrial tissue(Spach andDolber, 1986). Wealsoinvestigated
dimensionless form (Tyson and Keener, 1988; Winfree, the interaction of a plane wave with a lattice of L-shaped
1991): obstacles. Here, we sought to approximate diffuse disease
and investigate the fate of newly spawned wavelets.
au a2u a2u
The lattice semiperiod, L, and the dimensions of the
(1)
at ax2+ay2+f( ) V obstacles were adjusted to investigate the fate of newly
av spawned wavelets. The boundary conditions for the inex-
at = E(yU- , (2) citable obstacles were "no-flux" conditions, au/an = 0,
where n is normal to the obstacle boundary. We used
where u(x,yt) is a dimensionless membrane potential, and periodic boundary conditions for the domain boundaries
V(x,y,t) is adimensionless slow recovery current. The scale along the lines of symmetry, Ai (Fig. 1 B).
ofuisthemaximumactionpotentialamplitude UO,thescale
of Vis given by crNaUO, the scale oftime is T0 = Cm/0rNa,
and the scale ofthe length Lo is given by (D/O'Na)1/2. Here
METHODS
Cm is amembrane capacity andDis adiffusion coefficient.
The scale factor, y, is the ratio of slow potassium to fast We solved the system (Eqs. 1 and 2) numerically, with
sodium conductances. The small parameter, E, is equal to boundary conditions described above, by using the implicit
Cm/(TkO(Na). difference scheme with second-order approximation on the
The functionJfu) is the piecewise linear function similar spacegridinterval, Ax, andfirst-orderapproximation onthe
to that used by Starobin et al. (1994) and represents the i/v time grid interval, At. The two-dimensional version ofthis
properties similar to those of a membrane sodium current method isbased onreplacing the multidimensional Laplace
(Fig. 1 A). The slope ofthe right- and left-hand branches of operator with a sequence of two one-dimensional second-
this function, A, determines the maximum current that is order Laplace differential operators, one for the X dimen-
available to depolarize adjoining regions of the excitable sion and one forthe Ydimension. This method is known as
medium. Pertsov et al. (1983) used a three-segment piece- the fractional-step difference method (Richtmayer, 1957).
wise linearjAu) and variedthe slope ofthe middle segment, Eachone-dimensionalproblemcanbesolvedindependently
which is an equivalent way ofaltering the current available in each direction by the implicit second-order approxima-
for propagation. tion difference scheme resulting in a computationally effi-
In anotherwise homogeneous medium, we exploredcon- cient method. This approach has two important advantages.
ditions forwavefront-obstacle separation resulting from the The first is that it allows a larger time step while maintain-
collision ofaplane wave with a simple obstacle. As a first ing unconditional numerical stability, andthe second is that
approximation, we have chosen an L-shaped obstacle (Fig. being locally one-dimensional, this method significantly
1 B) as a piecewise linear approximation of the proximal simplifies the approximation of the multidimensional
aspectofasolid (2D) obstacle that mightbeencountered in boundaryconditions whilemaintaining ahigh-ordernumer-
diseased ventricular tissue after infarction. The filament- ical approximation of differential equations.
A
FIGURE 1 (A) The null-cines of B
the systemofEqs. 1 and 2. Thefunc-
tion,fu),isapiecewiselinearfunction, 1t
wheretheslopeofeachlinearelement,
A,referstotherateofthefastexcitation incident
process.Withrespecttomembraneion wavefront
A+1+ L
channels, thefunctionfu)corresponds
V
to the Na current ilv relationship and L
the slope, A, determines themaximum
Na conductance and hence influences v=7 u
mediaexcitability.Theslopeyrefersto
V=f(u)
the rate ofthe slow recovery process. C Ao
(B) Therelationship between the inci- Vcrit
dentwavefrontandaperiodiclatticeof L
inexcitable filaments. Theyaxis coin- A_1
cideswiththedirection oftheincident
wavefront. Thexaxis isperpendicular F
to the direction ofthe incident wave- A_-L
front.LinesA+1andA 1correspondto
aperiodofthelatticethatisequalto2L.
584 Biophysical Journal Volume70 February 1996
The spatial numerical domain was arectangular grid that (e-f). After pinch-off, the parent wavefront fragments col-
did not exceed 1400 X 400 nodes in Yand X, respectively. lide and reattach. The residual excited regions attached to
The grid parameters were At = 0.2, Ax = Ay = 0.25. Such the obstacle propagate toward each other, eventually extin-
valuesofthegridparameters allowedustokeepatleastfour guishing propagation (f-h).
or five grid nodes per wavefront and, thus, to find a stable Fig. 3 illustrates a similar sequence of events where the
solution without significant influence of grid effects. Re- vertical component of the obstacle has been extended. In
sults ofthe two-dimensional studies were displayed with a this case, a large "pocket" ofexcited medium remains after
colorformatsothatpatternsofmembranepotentialcouldbe the"pinch-off' ofthecurled, daughtersegmentoftheparent
readily identified. wavefront (f). Note the polygon-like object formed by the
We studied the variations in response ofthe 2D lattice to daughter segments and the obstacle boundary. The velocity
single-plane wave stimuli for different values of the exci- vectors atallpoints along the daughterwave fragments that
tationparameter A andrelaxationparameter Eforconstant y define the boundary of this polygon are directed "inward"
= 7 and a = (M3 -M2)/(M2-M1) = 2.75 (see Fig. 1 A). toward the excitable medium. Because each daughter frag-
Varying A permitted us to adjust the amount of source ment propagates toward the interior of this region (f-j),
current available forpropagation andto simulate conditions escape to excitable medium outside this enclosed region is
associated with slow conduction secondary to Na channel not possible. The "attachment" of the wavefront to the
blockade, high-frequency stimulation, or depolarized rest obstacleboundary maintains "virtual" continuity, thus elim-
potentials. inating the possibility ofwavelet "birth." In other words, a
For analytical studies, we used the singular perturbation "classical" wavebreak in the Balakhovsky-Krinsky sense is
technique with anonasymptotic matching ofu(x,y,t)written notcreatedaftercollisionofthewavefront andtheobstacle.
interms ofafunctional series forthe smallparameter E(see Note that this argument is true only when the medium is
Eqs. 1 and 2) (Starobin et al., 1994). We have defined a homogeneous (except forthe presence ofinexcitable obsta-
safety factor forpropagation, SF, similar to that introduced cles). Ifthere are spatial differences inrefractory properties
by Rushton (1937) and used by Spach et al. (1981), and as used in the original Krinsky model (1966), then there
applied it to the boundary layer at the wavefront-obstacle may be conditions that permit "escape."
interface. Written in terms of the available wavefront In contrast, when the excitability of the medium is re-
"source" charge, QS, and the charge requirements (load) of duced, the "parent" wavefront separates from the obstacle
adjacent medium, QL' (Fig. 2 B, b)-and the right fragment end begins to curl
because ofthe increased load "seen" by the tip region. The
SF= Qs-QL. (3)
vertical element ofthe L obstacle acts to prevent curling of
From the resulting equation we estimated the critical value the left attached fragment, creating an asymmetry in spiral
development that reduces the spatial requirements for frag-
of A, Ac,zt, for a range of E where wavefront-obstacle sepa-
ration occurs, i.e., where SF = 0 with respect to the con- ment formation (c-e). In concert with curling, the tip re-
figuration ofthe wavefront and obstacle. gions ofboth parent fragments also propagate toward each
otherandeventually collide As thecurledregions ofthe
(f).
two parent fragments collide, a small "daughter" fragment
QUALITATIVE CONDITIONS FOR resulting from the asymmetric development of the two
WAVELET FORMATION parent fragments is "pinched off' (g) and propagates ap-
proximately perpendicular to the parent propagation vector
We hypothesize that wavefront-obstacle separation is nec- (h). This daughter fragment becomes the source of a new
essary for wavelet formation. To demonstrate this, we first wavelet. We will now concentrate our analysis on the con-
consider the likelihood of wavelet formation when the ditions for wavefront-obstacle separation before wavelet
wavefront remains attached to the obstacle after collision "birth."
(Figs. 2 and 3), i.e., the medium excitability is high.
Initially, the incident wavefront collides with the inexcit-
able obstacle (Fig. 2 A, a and b). Because the boundary
ANALYTICAL CONDITIONS FOR WAVEFRONT
conditions for the obstacle are "no flux" boundary condi-
SEPARATION FROM THE OBSTACLE
tions (i.e., aulan = O, where n is normal to the obstacle
boundary), the flow of charge at the wavefront-obstacle Let us consider the collision between a plane excitation
interface is directed parallel to the obstacle boundary, re- wave and an inexcitable obstacle. For analytical simplicity,
sulting in propagation parallel to the boundary. Ifthe prop- we will consider the element parallel to the x axis of the
agation vectorwere directed away fromtheobstaclebound- L-shaped obstacle, which is formedby astraight long (1>>
ary, separation would occur. Consequently, the wavefront Lf) and narrow (d < Lf) filament and aligned parallel to the
tips associated with the two parent wave fragments propa- wavefront (Fig. 4 A). Here Lf is the wavefront thickness.
gate toward each other following a path along the obstacle Therefore, one can assume that wavefront-obstacle separa-
boundary (Fig. 2 A, c-d). As the two wave fragments tion depends only on the behavior ofthe medium within a
approach each other, a small region will be "pinched off' small region that is on the order ofLf and located near the
Starobin et al. Wavelet Formation in CardiacTissue 585
A
e
a b d
..: :: :.
B
*a b
J*_
RIl
.e f .g h
FIGURE 2 Wavefront-obstacle interactions forsafetyfactorsnearzero. Thetemporalevolutionofthewavefrontis showninframes fromlefttoright,
top tobottom. The membrane potential is color coded: the orange background in all frames refers to the restpotential u(x,y,t) = ml; the red shadows
representtheareaofmembranerecovery;thegreenregionrepresentstheexcitationzones.(A)WhenthesafetyfactorispositivesothatA> Akct,thenthe
chargeinthewavefrontsegments adjacenttotheobstacleisadequatetoexciteallofthemediumbetweenthewavefronttipsandtheobstacle.Themodel
parameters were A = 0.9; E = 0.018; y = 7; a = 2.75. Undertheseconditions, thereisnowavefront-obstacle separation. Afterthetwowavefragments
rejoin(bottomsequenceofframes),the"pinched-offtail"diesbecausethevelocityvectorsoftheresidualwavefragments aredirectedtowardrefractory
media. (B)Thereversephenomenon: wavefront-obstacle separation andformation ofadaughterwavelet. Inthiscase, wereducedA = 0.738 = Acrt, so
thattheSF= 0.Withthischoiceofparameter,thereisinsufficientchargeinthewavefronttiptoexciteallofthemediumbetweenthetipandtheobstacle
bordersothatthewavefrontseparatesfromtheobstacleandbeginstocurl. Afterreattachmentofthetwoparentfragments (bottomsequenceofframes),
a daughter wavelet is "pinched off' from the parent wavefront and propagates in a direction that is approximately perpendicular to that ofthe parent
wavefront. Notethatthe reduction in Aresults in a shorterwavelength, indicating slowerpropagation.
corner ofthe filament. Under these conditions, the analysis As the excitability is further reduced, the tip velocity will
can be considered to be quasi-one-dimensional. become negative and a portion ofthe wavefront will with-
Earlier we showed that if the wavefront thickness Lf is draw from the obstacle boundary (Mikhailov and Zykov,
less than Lcit (the maximum wavefront thickness necessary 1991). The transition value of the excitability where a
forpropagationatzerovelocity), thereisinsufficientcharge wavefront separates from the obstacle is determined by the
available in the wavefront to support propagation in a ho- critical wavefront thickness, LfC > Lcrt. LfC is determined
mogeneous medium. Consequently, our analysis will as- when SF = 0, which reflects the balance of "charge"
sumemediumproperties suchthatLf> Lcnt(Starobinetal., available within the wavefront adjacent to the obstacle
1994). boundary and the amount of charge necessary to form an
Inahighlyexcitablemediumtheendsofawavefragment impulse that accelerates from C = 0 to C,rit. This approach
are able to excite adjacent medium such that they can turn hasbeendevelopedbyZeldovichetal. (1985) andmodified
aroundthecornersoftheobstacle, slipalongits surface, and by Starobin et al. (1994).
recombine to reestablish the incident wavefront. As excit- We assume that wavefront separation evolves within a
ability is reduced (a smaller safety factor), the fragment small region of order Lf,. Because we can neglect the
ends are progressively less able to excite adjacent medium, changes of wavefront thickness that are less than Lfc, we
thus reducing the coupling with the surface ofthe obstacle. assume that separation occurs in a square transition region
586 Biophysical Journal Volume 70 February 1996
a b c
d
..
UNW .. r.
No
m '01"
..
m:....V;..:...
i. ...
..I y
..l .
f g i i
h
FIGURE 3 Similar to Fig. 3, this sequence offrames illustrates the fate ofdaughter fragments when the media excitability is adequate to maintain
wavefront-obstacle attachment,i.e., A = 0.9;SF> 0.Inthisillustration, theverticalarmofthe"L"obstaclehasbeenextendedtodemonstratethenature
ofdaughterwaveletannihilationafter"pinch-off'fromtheparentwavefront.Theno-fluxboundaryconditionsattheobstacleboundary(au/dn =0)ensure
thatpropagation ofthedaughterfragmentsattheobstacleboundaryisparallelwiththeboundary.Thisconditionensuresthatthevelocityvectorsofeach
residualwavefragment(bottomsequenceofframes)pointtowardrefractorymedia. Asshown,thepolygonenclosedbythedaughterwaveletseventually
collapses.
withasideofLfC (Fig. 4B). Separationoccurs iftheportion Here Veq is the equilibrium level of a slow recovery
ofthetransitionregionadjacenttothecorneroftheobstacle current, and Vcit and Lcit are the critical values of a slow
has azero (inbothXand Ydirections) propagation velocity current and awavefront length associatedwith azero prop-
associated with a critical value of a wavefront thickness, agationvelocity,respectively. Accordingtotheassumptions
mentioned above, these values can be found from the ID
Lcrit.
Integrating Eqs. 1 and 2 over the transition region and theory described in (Starobin et al., 1994)
time interval At, one can find that the source amount of
charge, QS = (M3 - ml)At (Lfc)2/2, produced in the tran- b = -O.5Co + (0.25CO + A)'`2,
sition region (source charge) during the time interval At.
The factor of 2 arises from our assumption that the wave- b2= 0.5CO + (0.25CO + A)1'2
frontvelocity inthe Ydirectiondecreaseslinearly from Ccrit
to 0 (so we use Ccri/2), whereas in the X direction the /b bll(b,+2+ (b)-b2l(bl+b2)
velocity is 0. The load charge is equal to the sum of the G= +b K= (m3-ml)IG
tbl
charge flowing into thetwoleading wavefronts (this charge (5)
is approximated by the gradient over the distance, Lcit),
QLW = 2Lcrit (M3 - m1)At/Lcit, andthechargerequiredfor Lf= ln(M-' K) b b Lcrt = 2ln(M-' K)l X
wavefront acceleration from C = 0 to C = Ccnt ofthe two b,b2 C
developing wavefronts (approximated by square impulses
ofsize,Lun1)propagatinginbothXand Ydirections, QLR = a-i
2(Lcrt)2 (Veq - Vcrt)At. Under normal conditions ofprop- Veq = Vct-A(m2-iml) 2
agation, the time constant of the slow process is much
greaterthanthetimeconstantofthefastprocess, sothatthis Here C0 is a wavefront propagation velocity in a 1D excit-
termiszero(i.e.,themediumintheregionofthedeveloping able cable which (forfu) shown in Fig. 1 and E = 0) is
wavefronts is at Veq). However, during early impulsedevel- given by
opment when the wavefront velocity is near zero, the time
constants arecomparable (OstrovskiiandYakhno, 1975), so Co= (a -1)(A) (6)
thatthe "slow" current cannotbe ignored. We have conser-
vatively estimated this charge requirement in terms of the where
gradient, Veq - Vct. With these components, we can write
the equation for SF as a = (M3 -M2)/(M2 -MI)
SF = QS- QL = QS - (QLW + QLR) ThefactorMis asmallconstantthatdefinesthethreshold
4) ofu, relative to ml (ml + M) and M3 (M3-M), associated
SF = /2(M3-mi)Lf - 2(m3 -Ml)- 2(Vq -Vcn)L, with the wavefront boundaries.
Starobin et al. Wavelet Formation in CardiacTissue 587
A Taking into consideration that the factor G in Eq. 5 is
insensitive to the propagation velocity CO, and substituting
Incident Obstacle Eq. 5 into Eq. 4, we get
Wave
SF lC2 10 1 a - I
Lcrit a8 J t2r c + 1J
(M`(m3
= In - m,)/2)
In azero-order approximation ofE (singularlimit, E = 0)
I
the wave propagation velocity Cis equal to CO, but adding
the first-order term reduces the velocity: C = CO(1-40E)
(Mikhailov, 1990). It was shown by Starobin et al. (1994)
that the wavefront propagation velocity of the peak of the
wavefront, C(V*), is given by
Y > C(V*) = CO(V,)(I - C0E) V* = V, - Ey(m3 -M)T
x (8)
where CO(Veq) is the propagation velocity ofatriggerwave
where recovery is absent, V* represents the value of the
B recovery current at the peak of the wavefront, and ro =
Lf(CO)/CO is the front formation time. Using Eqs. 5 and 6,
one can find Co, keeping the linear term in Taylor's expan-
sion of C(V*) near CO(Veq),
Ccrit o= y(aA)-3/2(a + 1) (9)
a - 1
SubstitutingEqs. 8 and9intoEq. 7,wefindtheequationfor
SF in terms of medium parameters and the roots ofJ(u)
SF I I- (]A3 1)3 1)'o
a-- a
y(a
Lcnit L 2or' a-1J 8at(a-1) 8a2(a- 1)
(10)
Here we see that the transition between wavefront-obstacle
separation andwavefront-obstacle attachment (i.e., SF = 0)
=Ccrit/2 can be altered by varying any of the medium parameters.
Because of the link between A and the sodium current i/v
relationship in excitable cells, we have chosen to focus our
C(Xcrit)=Ccrit analyses on conditions where A = Acrit, where Acrit is aroot
ofthe equation SF = 0 (Fig. 5 A).
Because or > 1 and E > 0, Eq. 10 has a single positive
real root, Acit when SF =0. When A > Acrt then SF >0,
and propagation at the wavefront-obstacle boundary suc-
ceeds, whereas when A < kcntq then SF < 0 and local
FIGURE 4 Theboundarylayerbetweenaplanewavefrontandaninex- propagation fails at the obstacle-wavefront interface and
citable obstacle. (A) Orientation oftheplane waverelative totheinexcit- wavefront-obstacle separation occurs (Fig. 5 A). As we
ableobstacleoflengthIandthicknessd.(B)Areasofthewavefrontandthe indicated in the introduction, depolarization of the mem-
obstacle that are considered in the analysis of separation. The incident brane potential, Na channel blockade, and high-frequency
wavefront of thickness, LfC, supplies charge to the transition region of
similar thickness. The transition from wavefront-obstacle attachment to stimulation allreduce Nachannel availability, equivalentto
wavefront-obstacle separation occurs when SF = 0, i.e., when thecharge reducingA. ThisresultsinareductionoftheSFbecausethe
available fromtheincident wavefront isequal tothechargerequirements availablewavefront"source"charge, Qs,decreases. Fig. 5A
necessarytoexciteallofthemediawithinthetransitionregion.Theparent shows that the smallerthe value ofA, the smaller the SF. It
wavefront velocity Ccj, and the critical wavefront length LfC associated
is important to pointouthere that the safety factorwehave
withtransition fromwavefront-obstacle attachment to wavefront-obstacle
separation reflect properties ofthe wavefront when SF = 0 atthe wave- definedaboveisa"local"propertyandis determinedbythe
front-obstacle boundary. interactionbetween alocal sourceandloadpresentedbythe
obstacle-wavefront interface as indicatedinourapproxima-
tions.
588 Biophysical Journal Volume70 February 1996
A M equal to 0.293 (M = 4.688Ax2). The values of grid
intervals used for all numerical experiments were At = 0.2
and Ax = Ay = 0.25.
SF
With the value ofMdetermined, we evaluated the accu-
racy of Acrt derived from the roots of SF = 0 (Eq. 10) by
estimating numerically the critical value of A associated
with the transition from wavefront-obstacle attachment to
wavefront-obstacleseparationusingnumericalexperiments.
We set y = 7, selected a value of E and then numerically
solvedtheFHNequations (Eqs. 1 and2)fordifferentvalues
ofA using a binary search. We terminated the search when
the difference between two successive estimates of A asso-
ciated with wavefront-attachment and separation was
<0.01. Fig. 5 B illustrates the theoretical predictions (Eq.
10) and the numerically determined values (circles). These
x
experiments revealed good agreement with Eq. 10 for val-
ues of E < 0.04, comparable to physiologic conditions
reflecting a 25-fold difference between the excitation and
2.0
recovery time constants. For smaller E when CCrit < 1, our
Xcrit numerical results illustrate a departure from the theoretical
Wavefront-obstacleattachment curve because the linear approximation ofthe front forma-
1.5
tion time o- = Lf(CO)/CO in Eq. 7 is valid only for Cci =
0(1) (Starobin et al., 1994).
1.0
1.0 o~
WAVELET EVOLUTION AND INSERTION BEHIND
0.5 THE PARENT WAVEFRONT
Q0 Wavefront-obstacleseparation
Wavefront formation and insertion is a function of both
wavefront-obstacle interactions and wavefront-wavefront
0. 00 0.010 0.020 0.030 0.040 collisions. Here we will illustrate the role of wavefront-
£ wavefrontcollisions in theinsertionprocess. Becauseofthe
mirror symmetry of the L obstacles in the lattice, it is
FIGURE 5 (A) The relationship between A and SF/Lci,. Shown are helpful to follow the evolution ofpaired wavelets as shown
theoreticalplots ofthesafetyfactor(SF,Eq. 10))asafunctionofAand E in Fig. 6 for several lattice elements.
when -y = 7. Note that as E is increased (decreasing excitability), , After collision of the parent wavefront and the obstacle
defining the transition between wavefront-obstacle separation and wave-
front-obstacle attachment, also increases, indicating the increased charge
(Fig. 6, bc), attached (a) and unattached (u) wavefronts
requirementsnecessaryforpropagationintolessexcitablemedia. (B)The develop. Next the daughter wavelets (labeled d in frame d)
dependenceofthecritical slopeoffu) forwavefront-obstacleseparation, are "pinched off' as a result ofreattachment ofthe parent
Acit,as a function ofthe model parameter E while y = 7. The solid line wave fragments (as illustrated in detail in e-h ofFig. 2 B).
represents the analytical approximation of Ac*t(root ofSF = 0),ACa*tnd the These daughter wavelets propagate toward each other in a
circles represent numerically determined values. The curve of as a v
function of E separates the plane into a region where wavefront-obstacle diresction that is approximately perpendicular to the parent
separationoccurs (belowthecurve) andwherewavefront-obstacle attach- propagation vector (Fig. 6, d-e). After collision of the
ment is maintained (above the curve). The numerical experiments reveal daughterfragments at the lines ofsymmetry,A+1 andA
goodagreementwith theanalytical approximation. residual wavelets develop that propagate antegrade to th1e,
parentwavefront, andtheirend points begin to curl (Fig. 6,
f-g). The final eventin theinsertionprocess is thecollision
To compare the numerical solution of Eqs. 1 and 2 and shown in frame h resulting in reattachment of adjacent
analytical results given by Eq. 10, it is necessary to define fronts. With our choice of obstacle and lattice dimensions,
factorM (Eq. 7) designating the boundaries ofawavefront. new secondary daughter wavelets arise along theA_ -L,
The precision ofour numerical results is bounded by M = AO, andA+1 + L lines ofsymmetry (Fig. 6 h) that lead to a
rq(At + Ax2) = max(ij,n)IuiA,x jAy, nAt - Uexact(XyYt)I. We repetition ofthe insertion process as illustrated in Fig. 7.
used this quantity as a threshold for defining the beginning The growth ofthe two wave fragments shown in Fig. 6,
and the end of a wavefront. Here uexact(x,y,t) is the exact f-h, reflect what we have labeled as "insertion," that is, a
solution of Eqs. 1 and 2, and io, jo, and no are the grid new wave front has been inserted behind the parent wave,
parameters, with 1 si sio, 1 sj sjo, 1 sn sno. Our increasing the apparent frequency ofexcitation ofa region
numericalexperimentsfoundthat qisequalto 1.12,making in front of the parent. Moreover, it is now possible to
Starobin etal. Wavelet Formation in CardiacTissue 589
A-1-L A- A0 Art A-.1+L
I:
I I b d
1 1 a c
.rJ ? 'I
.i LI Li
.
f h
e g
FIGURE 6 Thedevelopment ofsecondary waves andinsertion behindtheincident wave (SF = 0). The orangebackgroundinallframes refers to the
restpotential u(x,y,t) = ml; thered shadows representtheareaofmembrane recovery; thegreenregionrepresents theexcitation zones. The inexcitable
obstacleis showninwhite. Timeincreases, lefttoright,toptobottom (a-h). Themedium is initially excitedbehindtheobstacle ina.Thebreakofthe
parent wavefront shown in a is followed by propagation ofan attached (a) wave fragment and unattached (u) wave fragment (b-c). The ends ofthe
unattached fragment begin to curl, daughter (d) fragments are "pinched off" and propagate toward the A+1 and A-, lines of symmetry (c-e).
Simultaneously,theparentwavefrontisrestored(d).Aftercollisionofthedaughterwavelets(e),residualsecondaryexcitationzones(newwavelets)begin
togrowaroundtheA+1 andA 1 linesofsymmetry (fandg).ThesewaveletspropagatebothantegradetotheparentwavefrontandtowardtheAOlineof
symmetry (g), wheretheycollide(h)andtheinsertionbehindtherestoredparentwavefront(f-h)iscompleted.ThemodelparameterswereA = 0.738 =
A,,ri; E = 0.018 and y = 7, At = 0.2.
visualize a mechanism for increasing the frequency of the behind the parent wavefront. Here we illustrate the role of
post-obstacle wave train. Ifthe post-obstacle train were the the lattice and obstacle dimensions in creating interference
resultofasinglerotating spiral (intheabsenceofobstacles), patterns that permit rapid additional wavelet formation and
then the time between successive wavefronts would be insertions.
equal to the time required for one spiral revolution. In Although we consider the events after a single incident
contrast, the time required forinsertion ofanew wavefront wavefront, our study approximates the clinical setting
using interfering wavelets is approximately one-fourth the where the excitation rate may be low (SA node rate of 1/s).
revolution time ofa single spiral. Aswewill see inthenext It is important to note that the combination of low-excit-
section, the lattice spacing can lead to two daughter wave- ability and low-frequency periodic stimulation results in a
letsarisingbetweeneachpairofobstacles (Fig.7c)thatwill large excitable gap (i.e., the wavelength is much less than
further double the insertion rate. the distance between wavefronts initiated by stimulation),
sothatinsertion ofadditional wavefronts isfeasible. Recent
studies of regional excitability during atrial fibrillation re-
MULTIPLE WAVEFRONT INSERTIONS BEHIND
vealed that the reentry circuits do contain excitable gaps,
THE INCIDENT WAVEFRONT
allowing regional control of atrial fibrillation by local ex-
The final stage ofour study was to demonstrate how high- ternal stimulation (Kirchhofet al., 1993).
rate excitation of a region of excitable medium in front of In our model, the collision of wavelets near the lines of
the lattice is continued in the absence of additional stimu- symmetry results in fragments that can potentially be in-
lation. This is accomplished by rapidly creating and insert- serted behind the parent wavefront. The wavelet multipli-
ing new wavelets behind the original incident wavefront. cation results from repeated wavelet spawning along the
Above, we illustrated insertion of the first new wavefront A+1 + L,AO, andA 1 + L lines ofsymmetry as illustrated
590 Biophysical Journal Volume 70 February 1996
AI-L .A-, AIJ A,, A*a+1.
__1.7 ".P
* .I........l. ft.
!9:. :;i.!
.:e f:..
............ .:..
................
.:...:.:...'...............................'.:,-.+3
..:.. ... ,: { ;
..W:... .'
......
FIGURE 7 This sequenceofframes (lefttoright, toptobottom) is acontinuation oftheeventsshown inFig. 6andhighlights the sequenceofevents
leading tomultiple wavefront insertion. Note thatin Fig. 6, thewavelets thatare eventually insertedbehind theparentwavefront are located attheA+1
andA 1 linesofsymmetry (Fig. 6f). Afterinsertion andcollision ofthecurledregions (Fig. 6g-h), residual fragments remainattheA+, + L,AOand
A 1 + Llinesofsymmetry, whichbecomethenewwaveletgenerators. Fig. 7illustrates thedevelopmentofthesefragments. Notethatfirst,theresidual
fragments propagate retrograde to the parent and inserted wavefronts (b). These wavelets collide with the obstacles, resulting in paired new zones of
excitation(c)thatpropagatetowardtheA_l andA+1 linesofsymmetry(d).Thewaveletsabovetheobstaclecollidewhilethewaveletsbelowtheobstacle
grow andcollide (e). Now theprocess shown inFig. 6isrepeated. Theresidual daughterwavelets (developedboth infrontofandbehind the obstacle)
propagateantegradetotheparentwavefront(fandg)andeventuallyjoin,completingtheinsertionprocess. Notethatinthisillustration,theobstacleand
lattice semiperiod dimensions are suchthattwo new wave fragments areformedthatareeventually insertedbehind the parentwave.
in Fig. 6, c-f, and Fig. 7, b-e. Fig. 7 a illustrates the first wavelet to develop subsequent to "pinch-off' (Fig. 6, c-d).
repetitionofresidual waveletdevelopmentneartheobstacle Similarly, lattice semiperiods < Lmin provides insufficient
attheA+1 + L,AO, andA1- Llines ofsymmetry (shown space for wavelet development such that subsequent to
in Fig. 6 h). These wavelets propagate retrograde to the collision, an inserted wavelet is born (Fig. 6J).
parentwavefrontandcollide withtheobstacle (Fig. 7, b-c). The rate of insertion is also dependent on Acrit Fig. 9
Inthis case, twodaughterfragments areborn (Fig. 7 c), one shows therateofinsertion perdimensionless timeunit. The
thatpropagates behind the obstacle and one thatpropagates larger Acrit (indicating large e; Fig. 5, A and B), the shorter
within the obstacle and escapes. As these two pairs ofnew thetimeintervalrequiredforpairedfragments initiatednear
fragments grow, they collide with theirsymmetrical partner symmetry linesA+ andA to reach theAo symmetry line
1 1
aroundtheA+1 andA 1 linesofsymmetry (Fig. 7,d-e) and (Fig. 6, e-h; Fig. 7, e-h). Similarly, the fragment initiated
the residual wavelets propagate antegrade to the parent near theAo line ofsymmetry (Fig. 6, b-e; Fig. 7, a-e) can
wavefront, developing into two newly inserted wavefronts evolve rapidly toward theA+ andA -lines of symmetry,
1
(Fig. 7,f-h). and hence, the largerthe number ofadditional plane waves
In this example, wavelet multiplication and insertion is available to be inserted behind the restored incident wave-
dependent ontheobstacle andlattice dimensions. Shown in front.
Fig. 8 aretherelationships betweenthecritical slopeofftu),
Acr,t, the minimal lattice semiperiod (Lmn,), and the mini-
mum horizontal length of the obstacle ( ,m) for wavelet DISCUSSION
formation and subsequent interference leading to insertion.
ACrit was computed as afunction ofE when y = 7. Obstacle Atrial and ventricular fibrillation are significant clinical
lengths < lmi. provide insufficient area for the pinched-off problems, yet our understanding ofthe mechanisms ofini-