1. Introduction

Bioelectric phenomena predate all artificial forms of electricity. In bioelectricity endogenous electrically mediated signalling regulate cell, tissue, and organ-level modelling and behaviour. Cells and tissues of all types, in the main, use ion fluxes to communicate electrically. Bioelectric phenomena are usually described by volume conduction models. Such models describe the topology and conductivity of tissues wherein electric current flows and also the current sources in the tissue. The heart is typical of an autorhythmic volume conduction system, and it is therefore being considered here. It is here named autorhythmic volume conduction system (AVCS). Autorhytmicity is a unique attribute of cardiac muscle cells. The cells can generate an autonomous action potential (AP) to an enormous degree.

The passage of ionic current across the membrane of active excitable cells induces sources and sinks. In the AVCS the nodal structures, which include the sinoatrial node, atria, atrioventricular nodes (AVN), His bundles, Purkinje’s fibres act as electric sources and sinks. The heart is a site for the transmission of autonomous and automatic AP. The propagation of such APs is contingent on the efficient conduction of electrical current from an excited cell to an abutting quiescent cell. Sources and sinks are created when cells give charge to abutting quiescent cells. Almost every cell in the AVCS acts as a source and a sink. The source-sink balance is an essential issue in the cardiac system. In this regard, Unudurthi et al. [1] depicted the pivotal role of the sinoatrial node (SAN) in the integrity of cardiac source-sink balance.

Nikolaidou et al. [2] showed the elegant source-sink relationships that reside in the SAN region. The consummation of the conduction system presupposes the generating of the quantity of current that would be sufficient for the local sink. This is one possible way of avoiding source-sink mismatch which is largely implicated in the deleterious cardiac events. The choice of the term feasible in the caption of this work implies a mismatch-free dipole. In this regard, Boyle and Vigmond [3], Kleber and Rudy [4] stressed safety needs and therefore the quantification of the source-sink balance.

The question of how specialized cells prosecute, nay, control the conduction system given their susceptibilities requires great attention. For instance, how does the relatively small SAN cope with the electrical exigencies of the much larger atrial tissue? Such a source-sink mismatch is also seen in the Purkinje fiber–ventricular junction where a small Purkinje fiber (source) is coupled to a large mass of ventricular tissue (sink). How do the Purkinje cells satisfy the ventricular mass to cope with a source-sink mismatch? To be quite succinct, one may read what Joyner and Capelle [5] say as regards the earlier question, and see Morley et al. [6] for the latter question. We are yet to define the conduction network or, if you like, pathway: This is sequentially made up of the SAN, the atrioventricular node (AVN), and the His-Purkinje system. The anatomy of each of the above, which is not within the scope of this work, may be found in numerous Anatomy and Physiology works of literature including Monfredi et al. [7] and Kurian et al. [8]. Cardiac myocytes are connected end to end by intercalated disks. Adjacent to the discs is the gap junctions. Of physiological and electrical importance are the low resistance gap junctions. These junctions allow APs to spread from one myocyte to the next. They decisively determine how much depolarizing current passes from excited to quiescent regions of the network. Thus, they determine the speed and safety of the conduction process [9]. Gap junctions are known sites for electrical coupling between cells. Woodbury and Crill [10] supplied a treatise on electrical potential at the gap junction between two contiguous cells way back in 1970. Writing in contemporaneous times with the former, Heppner and Plonsey [11] produced a seminal work on the electrical interaction of cardiac cells. In real cardiac tissue, however electrical conduction takes place in both the intracellular and extracellular space in cardiac tissue, wherein there are variations in conductance and anisotropic ratios. Heppner and Plonsey [11] considered the intracellular and extracellular media as passive volume conductors of longitudinal resistances.

At present, mathematical treatment of the cardiac source-sink process is sparse in literature. A detailed understanding and application of graph theory [12,13,14] seem invaluable in the study of the AVCS. This work took this reasoning into account; it presents a mathematical approach to the understanding of the AVCS mechanism which is anchored on cell-to-cell source-sink electro-physiological flow. The issue of cell-to-cell flow here was treated as a network problem (see Rubido et al. [15,16], Wang and Lai [17], and Ahuja et al. [18]).

2. Graph Theoretics of the Source-Sink Network

Let G = (V, E) be a directed graph (digraph), where V and E are the node (vertex) set and segment (edge) set respectively. Let N = N(s, t) be a network, with two distinct vertices- a source s and a sink t, together with a non-negative real-valued function c defined on its arc set A. The vertex s corresponds to a production centre, or in the present case the impulse generating node, and the vertex t to a consumer, or the adjoining impulse receiving node. The function c is the capacity function of N and its value on a segment l is the capacity of l. The capacity of a segment may be seen as the maximum rate at which a substance can be transported through it.

An (s, t)-flow in N is a real-valued function f defined on a set S such that

(1)

where v represents a generalized vertex, I is the set of all intermediate vertices. Condition (1) is known as the conservation condition. A feasible flow must satisfy, in addition, the capacity constraint:

(2)

The upper bound in condition (2) imposes the natural restriction that the rate of flow along a segment cannot exceed the capacity of the segment.

Definition 1. [13] Let f be a feasible flow in a network N. A source-to-sink segment in a graph G is an f-augmenting path P such that for each l ϵ E(P):

(3)

Suppose τ(l) = c (l) – f (l) when l is forward on P, and let τ(l) = f(l) when l is backward on P. The tolerance of P is

(4)

Definition 2: A directed graph G = (V, E) is strongly connected if and only if, for each pair of nodes u, v ∈ V, there exists a path from u to v.

2.1. Graph Metrics

Graph metrics reveal the diverse characteristic of a defined network. Some salient ones are as follows:

2.1.1. Integration Metric

This estimates the relative ease with which cell regions communicate. If v_i and v_j are two adjoining nodes then the Shortest Path Length (SPL) between them is the shortest path that reaches v_j beginning from v_i. It is the shortest number of edges between v_i and v_j; for a weighed graph SPL is the sum of all the edge weights in the shortest path between v_i and v_j, calculated as

(5)

where G_{i⇔ j} is the shortest path between v_i and v_j.

In measuring functional integration global efficiency is the choice metric.

Definition 3 (Bryan et al. [19]) Let P_n denote the path on vertices v₁,v₂,…,v_n with edges v₁v₂, v₂v₂ …, v_n−1v_n. The distance d(v_i, v_j) between distinct vertices v_i and v_j is |i − j|. Hence the efficiency between different vertices v_i and v_j is

The global efficiency of the graph is calculated as

(6)

2.1.2. Segregation Metric

This measures the presence of clusters or modules within the network. Modularity is the degree to which a network may be partitioned into prominently delineated and non-overlapping groups. We may not bother much about the intricacies of this metric since the network in this study is not as complex as to require its use.

2.1.3. Measure of Centrality

Degree centrality is calculated by counting the neighbours of each vertex. It is given by

(7)

with v_k ∈ V, V the set of nodes, n = |V| and deg(v_k) the degree of node v_k.

Closeness centrality is calculated by summing distances from a vertex to each other. It is encoded as

(8)

with v_k ∈ V, n = |V| and dis(v_i, v_j) the distance from node v_i to node v_j,.

For a vertex v, betweenness centrality is defined as [20]

(9)

where σ_st is the total number of shortest paths from node and σ_st(v) is the number of those paths that pass through v (not where v is an end point).

Describing betweenness centrality of a node v in proportionality to the number of occurrences of itself on all geodesics of a graph, the evaluation is [21]:

(10)

where g_jk is the number of geodesics from node x_j to node v_k, g_jk (i) the number of geodesics from x_j to v_k containing v_i; the double sum is evaluated on all pairs (j,k) such that j ≠ i ≠ k and j < k.

2.1.4. Measure of Resilience

This is a metric of the vulnerability of the network as regards functional integration.

2.2. Impulse Transmission on the Conduction Network

The graph G = (V, E) here is treated as analogous to an electric network. Kirchhoff’s flow network is ideal in the treatment of minimal cost in transportation, including biological steady-state systems. We assume that the flow network obeys Ohm’s law. Thus,

(11)

In the above equation, the left-hand side represents voltage drop across the segment (i→j) connecting nodes at i and j, with a source situated at node s and a sink situated at node t. The segment (i→j) has a resistance R_ij against the flux f_ij which is determined by the conductance of the segment l. By conservation condition (11), the net flow at any node in the network must vanish. Therefore,

(12)

except at the source s and the sink t. The net flow at node i is governed by the Kirchhoff law:

(13)

where the bandwidth W_ij of the segment ij is inversely proportional to the resistance R_ij. W_ij ≈1/R_ij is the matrix representation of the network arrangement, (see [22] in which such arrangement was treated).

2.2.1. Constructing the Adjacency Matrix of the AVCS

Figure 1(a). Electrical conduction system of the heart [23]

Figure 1(b). Node-edge schematic of the CCS [24]

Definition 4. The adjacency matrix of a directed graph G is the number of edges that issue from vertex v_i and go into vertex v_j. Thus, if D = (d_ij) is the adjacency matrix of G, then

(14)

The adjacency matrix of the AVCS is constructed using Fig. 1(b) or equation (14) as (see Nzerem and Ugorji [24])

(15)

In (15) above v_i on the row indicates the initial node of an arc e_j while each v_j on the column indicates the terminal vertex of the arc. From (14) we have

(16)

Each of the nodes represents a part of the ACVS, to wit- Sinoatrial node (SAN), v_1, atrioventricular node (AVN), v₂, the bundle of HIS bifurcation point, v₃, the left bundle of HIS branch v₅, and the right bundle of HIS branch v₆. It is instructive to note that

(17)

For any isolated node v_k,

(18)

For a networked system, the total flow at any node obeys Kirchhoff conservation law:

(19)

where A is the adjacency matrix of the network. Another defining component of the network is the Laplacian matrix, M, of the fundamental network, which is obtainable from MV = F, whose elements are given by

(20)

Such a matrix was detailed in [24].

3. Ionic Flow and Action Potential

Like biological fluid-structure interaction problems [25], the cell-to-cell flow consists of the cardiac cellular structures and the ionic substances in current. However, unlike the Eulerian formulation that characterizes fully developed fluid flow, the present flow is amenable to Kirchhoff’s laws. The treatment of general cardiac APs is somewhat herculean due to variations in the constitutive ions and APs of the various cardiac myocytes. It would seem rather plausible to treat conduction in cell-to-cell segmentations, thereby accommodating their idiosyncratic APs (see Qu et al. [26]). In the AVCS, the electrical impulse is generated in the SAN, which spreads quickly through the atria to the atrioventricular node (AVN). (Note that the SAN is electrically coupled to the atria). The impulse travels from the AVN to the His-Purkinje system. Eventually, the Purkinje fibres transmit the impulse to the ventricular muscles at numerous discrete locations. The problem of source-to-sink distribution involves triggering an AP in a quiescent cell by an abutting active cardiac cell. The architecture of a cardiac strand, with its many individual cells of irregular shapes, is a mathematical task. Here each cell is assumed cylindrical and non-tapering.

3.1. Bidomain Electrical Flow

Bidomain structure defines a model of the heart tissue comprising two interpenetrating domains representing cardiac cells and the surrounding space, the intracellular and the extracellular domain [27,28]. The continuous topology of the intracellular space enhances ion travel from the inside of one cell to another via gap junctions without getting into the extracellular space, nor would ion travelling in the extracellular region require entry to the cell. Both spaces have electrical autonomous properties, therefore they form an anisotropic medium.

Let the associated internal and external potentials of the domains be and . Then the transmembrane potential is . The current associated with the internal and external potentials are and respectively. By Kirchhoff law

(21)

(22)

where I_m represents the current through the cell membrane dividing the two regions, B_m is the ratio of the surface area and the cell volume. The transmembrane current has the form [29]

(23)

Current flows through the boundary of extracellular space, but this is not so at boundary of intracellular space [30]. The boundary conditions for the bidomain model is of the form

(24)

where is the cardiac domain, n is the outer unit normal vector to the boundary . Since the flow of net current must be zero, we have

(25)

Thus, charge accumulation does not exist in a feasible flow since any entering current must be exiting.

3.2. Potential Field and Current in the Gap Region

Here we consider the potential field in the inter-gap region between two abutting cells. Two abutted identical cylindrical fibres of radius r_o, one assumed active and the other assumed inactive, are separated by a distance δ.

Define the potential field ϵ(r, z) in the region, 0 < z < δ; 0 ≤ r ≤ r_o. It is assumed that: the fibres are draped with non-excitable membranes, with specific resistance R_d (Ω-cm²); the end-cap resistances R_d are entirely passive and do not depend on voltage; the flow of current depends only on the voltage drop across the gap. In the specified region, the potential ϵ will satisfy Laplace’s equation . With the cylinder assumed axisymmetric we have

(26)

With the membranes assumed thin and situated at z = 0 and z = δ the associated boundary conditions are [10]:

With the rest of the present section emphasizing [10], conditions (iii) and (iv) are imposed to have a solution of equation (24) is in the form

(27)

where k_n are the roots of Applying boundary conditions (i) and (ii) gives

and

The coefficients A_n and B_n by multiplying both sides of the above two equations by rJ₀(k_mr), integrating with respect to r over [0,r₀] and applying the orthogonality property of the Bessel functions. Substitute the resulting expressions into equation (25) and get

(28)

where

The gap region accommodates current density with given by equation (28) above. Three notable currents come to play. They are the currents crossing the active and quiescent disc surfaces and the current across the cylindrical wall at r = r₀ (0 < z < δ). Denote these currents by IA, I_Q, and I_W respectively. The total current flowing into the gap from IA is

The current flowing into the quiescent cell is

The residual current, I_W, across the cylindrical gap boundary is

Substitute equation (26) into the above three equations and interchange the order of integration to get

(29)

(30)

(31)

where

3.3. Source-Sink Analysis

Source-sink interaction, which is achieved through excitation, is a sine qua non to the proper functioning of the SAN and other cardiac cells. The excitability of a cell can be expressed loosely as the ease with which a response could be prompted. It is often described as the least possible current necessary to depolarize the membrane to the threshold potential [31]. Electric sources of bioelectric origin are generated by the passage of current across the membrane of active (excitable) cells. How depolarizing ‘source’ current generated by the SAN drives depolarization and activation of the contiguous atrial tissue (current ‘sink’) remains uncertain. Bioelectric sources can refer to surface/volume distributions of two types of source elements, that is, the monopole and/or dipole. The dipole is a combination of a source and sink of equal strength. In the conduction system, the self-excitatory SAN is the leading source of action potential. (Note: From Fig. 1(b) or using (17) we have)

The zero indegree indicates the inherent self-excitatory property Of the SAN. The SAN and the AVN constitute a dipole. In a similar mode, the AVN and the HIS bundle constitute a dipole. In effect, every activated sink becomes a source. Suppose f is a feasible flow and [S, T] is a source/sink cut, then the net flow out of S and net flow into T equal val(f). Assume flow from SAN to HIS bundle is feasible (as in a physiological flow), we have

(32)

with U as a set of nodes in the section of the network being considered [13].

Assume a cylindrical conducting medium wherein a point source (the SAN) is located. In this regard, we may require the equation of the field of the point source inside the medium, with co-ordinates . Thus, the potential of a unit source located at some point inside the tube, which is bounded above and below by the planes z = −L and z = L and on the sides by the cylinder = a may be sought. As L approaches infinity the field of the point source inside the conducting cylinder is given by [32]

(33)

The modified equation of SAN membrane current is given by [23]

(34)

(35)

where:

A be the intra-cellular region consisting of the SAN

i_Ca,L, i_CaT are L- and T-type Ca²⁺ currents;

i_K,r, i_K,s are rapid and slow delayed rectifying K⁺ currents;

i_f (the funny current) is the hyperpolarization-activated current given by

Na⁺ is the sodium ion;

i_st is the sustained current;

i_to is the transient outward current.

C_m is the membrane capacitance.

The canonical form of the equation of transmembrane current is (23)

Considering equations (29) and (23) and taking the SAN as the active cell (source), the transiting current, which would drive depolarization and activation of the contiguous cell is

(36)

where I_ion is the total ionic currents exiting the SAN.

Equation (36) holds well for the SAN as a source, which recruits excitation through the atria to the adjoining AVN cells. The time-dependent membrane voltage may be furnished by the non-linear two-dimensional Fitzhugh-Nagumo equations, which describe the behaviour of an excitable cell. The equations read:

(37)

(38)

where V* encodes the membrane voltage, W is a quantity associated with the recovery of the cell after it fires, ϵ, γ, and β are non-negative dimensionless parameters. Equation (38) may be seen as the current equation; the capacitive charging term, proportional to dV*/dt, on the left and the nonlinear channel currents on the right. The form of equation (37) shows that the recovery of the cell is denoted as a gate opening equation with voltage-dependent rates. The amplitude of ϵ corresponding to the inverse of a time constant determines the rate with which W changes relative to V. The analysis of the system of equations (37) and (38 is by the method of phase plane, which may be found in [33].

3.4. Node-Edge Capacity and Maximal Flow

A source-sink mismatch is a deleterious liability to numerous impulsive systems. Some biological systems manage to overcome source-sink imbroglio through some anatomical topology. Bartos et al. [34] noted that the ability of SAN to overcome the source-sink mismatch to activate the surrounding atrial myocardium stem from the distinctive anatomy of gap junction proteins. The capacity of an edge c_ij may be expressed as the bandwidth of the edge, being the maximum flux that the edge can deliver devoid of cramming or damage. Thus, the capacity bandwidth of an edge is inversely proportional to the edge resistance, we write

(39)

where λ is a capacity parameter.

Assume for convenience first that the AVCS network has n nodes that consist of a single source SAN and a single sink AVN. We tag the source node 1 and the sink node n. Let c_ij denote the capacity of arc (i,j). For a flow in this network, let x_ij denote the quantity flowing from node i to node j along an arc (i,j). Then x_ij must satisfy

(40)

One finds that

(41)

A formulation of the maximal flow problem is

(42)

Subject to

(43)

The maximal flow problem furnishes a linear programming problem, which may be solved using some programming methods (see Kolman and Beck [35]).

The cost implication of the flow may be required. The work by Hung and Chien [36] may be employed to find such cost. In this regard, let E be the path from vertex u to vertex v through edges e_j, j=1, 2, ..., (r+1), and vertices u_j, j=1, 2, ..., h. Thus,

(44)

The cost of passing a unit of service of kind l (in this case, current), l = 1, 2,…,w, through the path P, encoded by Je_l(P), is given by

(45)

where Je_l(e) is the edge cost function, Jv_l the node switch cost function. An issue of great importance is how to maximize flow at the least cost. This however is not within the purview of the present work. If Je_l(e) = ∞, it is illegal to pass service on a path e. If Jv_i = ∞, it is an illegal edge to pass a service from e_i through v_i to edge e_j (i ≠ j). In the AVCS segment, where d⁺(v₁) ≅ 3 and d⁻(v₁)= 0, the edge cost of SAN-AVN transmission involves maintaining three (3) transmission lines to ensure that the individual (and by extension the total) capacity c_ij of the edges is not compromised.

3.4.1. AVCS Global Efficiency

The question of AVCS global efficiency may be difficult to furnish from the standard equation (6). This arises from the fact that for each of the expressed nodes v_k, other than the SAN, d⁻(v_k) = 1. It is of note that, by Definition 2, the AVCS graph is not strongly connected.

The centrality concept measures the importance of nodes in a graph network. Here the graph is the G_AVCS whose nodes are not interconnected. Centrality may be from the standpoint of the impact of a node on other nodes. Following [37], let v* be the node with the highest degree centrality in G. Let be the |Y| node connected graph that maximizes the quantity (with y* being the node with highest degree centrality in X):

(46)

where C_D(.) is the centrality degree of a given node. Therefore, the degree centralization of the graph G is

(47)

H is maximized when graph X contains one central node to which all other nodes are connected. For any graph G,

(48)

Degree centrality is the most widely used measure of centrality. Since it measures the number of direct neighbours of a node, it explains the importance of the node in a given network. In the brain network, it measures the influence of a brain region on adjacent brain regions [35]. In the AVCS network analysis, the degree centrality of a cardiac nodal cell shall measure the influence of a nodal cell on other nodal cells. In the AVCS we have,

SAN:

AVN:

(49)

HIS Bundle:

(50)

It may matter to consider the system based on weighted degree centrality on the out-degree of the neighbours. Thus, the Katz centrality measure

(51)

holds, A is the adjacency matrix of the graph G with eigenvalues lambda, β is a constant initial weight given to each vertex so that vertices with zero in-degree (or out-degree) are included in calculations (see Gera [38]).

4. Summary and Conclusions

Sources and sinks are typical of flow phenomena. A source is a place/point of inception of flow whereas a sink is a place/point of accumulation of flow. Thus, a sink is a negative source. In nature, there is no presupposition of all-time source and sink of equal strength (dipole). In the bioelectric flow considered, sources and sink mismatch exist without prejudice to the physiological vivacity of the immediate and remote cells involved in the mechanism of flow. The bioelectric conductor considered here is the AVCS. As a rule, the driver of cell-to-cell AP is ionic current. Although the heart is both functionally and anatomically complex, humans, and perhaps animated life with cardiac structures, are compensated by the non-strongly connected nature of its electric flow network. A strongly connected network system is more or less a culprit in the event of deleterious cascading failures. In the event of SAN failure, the AVN is a known surrogate pacemaker, which can initiate the firing of AP to the contiguous cell, the His bundle. In the absence of autonomic nervous stimulation, the rate of ventricular contraction is set at 40–60 bpm by the AVN. It is instructive to note that the serial surrogate responsibility of the AVCS cells is largely compensatory and as well averting.

Source and sink system is a network issue. This work employed the benefits of graph theory in analyzing the AVCS. An essential requirement of a network is the efficiency of the edge and nodal structures. Global efficiency is a measure of cumulative effect of the edge-nodal system, especially concerning strongly a connected weighted graph. Given the nature of connectedness of the AVCS, cardiac efficiency is often defined in terms of the mechanics of the heart rather than the electrical system (see for instance [39]). Some interesting points are worthy of note/ emphasis:

(i). The AVCS consists of pairs of sources and sinks of unequal strength.

(ii). Source and sink mismatch of the AVCS is anatomical, but it may be implicated in pathophysiological conditions. (However, in an anatomical state, it would not be an independent marker of a cardiac event.)

(iii). The non-strongly connected nature of the AVCS is an ethereal blessing, as it slows down the rate of cascading failures following any segmental/nodal disorder. This accounts for the relative enduring life of the critical organs in the event of mild to serious AVCS disorder.

(iv). The cost of achieving maximal flow is not of necessity an economic one; it requires maintaining the node-edge vivacity.

Since autochthonous electrical excitation resides in the SAN, which is the source of AP transmission, and since the other conduction nodes are surrogate pace-makers, it is reasonable to believe that the resolution of conduction anomaly can be found by in-depth analysis of source and sink system associated with the AVCS.

ACKNOWLEDGEMENTS

The authors are grateful to the reviewers for their meaningful input that helped in the revision of this work.