The Organization of the Sinoatrial Node Microvasculature Varies Regionally to Match Local Myocyte Excitability

Abstract The cardiac cycle starts when an action potential is produced by pacemaking cells in the sinoatrial node. This cycle is repeated approximately 100 000 times in humans and 1 million times in mice per day, imposing a monumental metabolic demand on the heart, requiring efficient blood supply via the coronary vasculature to maintain cardiac function. Although the ventricular coronary circulation has been extensively studied, the relationship between vascularization and cellular pacemaking modalities in the sinoatrial node is poorly understood. Here, we tested the hypothesis that the organization of the sinoatrial node microvasculature varies regionally, reflecting local myocyte firing properties. We show that vessel densities are higher in the superior versus inferior sinoatrial node. Accordingly, sinoatrial node myocytes are closer to vessels in the superior versus inferior regions. Superior and inferior sinoatrial node myocytes produce stochastic subthreshold voltage fluctuations and action potentials. However, the intrinsic action potential firing rate of sinoatrial node myocytes is higher in the superior versus inferior node. Our data support a model in which the microvascular densities vary regionally within the sinoatrial node to match the electrical and Ca2+ dynamics of nearby myocytes, effectively determining the dominant pacemaking site within the node. In this model, the high vascular density in the superior sinoatrial node places myocytes with metabolically demanding, high-frequency action potentials near vessels. The lower vascularization and electrical activity of inferior sinoatrial node myocytes could limit these cells to function to support sinoatrial node periodicity with sporadic voltage fluctuations via a stochastic resonance mechanism.

The cardiac cycle starts when an action potential is produced by pacemaking cells in the sinoatrial node. This cycle is repeated approximately 100,000 times in humans and 1 million times in mice per day, imposing a monumental metabolic demand on the heart, requiring efficient blood supply via the coronary vasculature to maintain cardiac function. Although the ventricular coronary circulation has been extensively studied, the relationship between vascularization and cellular pacemaking modalities in the sino-atrial node is poorly understood. Here, we tested the hypothesis that the organization of the sino-atrial node micro-vasculature varies regionally, reflecting local myocyte firing properties. We show that vessel densities are higher in the superior versus inferior sino-atrial node. Accordingly, sino-atrial node myocytes are closer to vessels in the superior versus inferior regions. Superior and inferior sino-atrial node myocytes produce stochastic subthreshold voltage fluctuations and action potentials. However, the intrinsic action potential firing rate of sino-atrial node myocytes is higher in the superior versus inferior node. Our data support a model in which the micro-vascular densities vary regionally within the sino-atrial node to match the electrical and Ca 2+ dynamics of nearby myocytes, effectively determining the dominant pacemaking site within the node. In this model, the high vascular density in the superior sino-atrial node places myocytes with metabolically demanding, high frequency action potentials near vessels. The lower vascularization and electrical activity of inferior sino-atrial node myocytes could limit these cells to function to support sino-atrial node periodicity with sporadic voltage fluctuations via a stochastic resonance mechanism.
The function of the sino-atrial node (SAN)  is to produce the action potentials (AP) that initiate each heartbeat. These APs are generated by clusters of pacemaker cells firing in unison 1 and propagate via gap junctions to neighboring SAN myocytes and, eventually, surrounding atrial myocytes 2 . SAN myocytes that are firing APs at the highest frequency become the dominant leading pacemaker site.
The dominant pacemaking site in the SAN is not static 20,22,23 . Instead, it dynamically shifts within the SAN in response to physiological stimuli, including activation of the autonomic nervous system 2 . A recent study on ex-vivo rat and human SAN tissue supports the concept that two distinct pacemaker sites co-inhabit the SAN, one near the superior section of the node and another further below in the inferior section 24 . Depending on the degree of sympathetic and parasympathetic nerve activity, the leading pacemaker site in the SAN shifts between superior and inferior regions. With high sympathetic drive, the superior SAN functions as the dominant pacemaking center. The pacemaker center shifts towards the inferior SAN during strong parasympathetic activation. Accordingly, a shift in pacemaker origin has been proposed to trigger changes in heart rate 24 .
Irrespective of pacemaker location, SAN myocytes require energy in the form of adenosine triphosphate (ATP) for normal function. For example, significant amounts of ATP are consumed for the maintenance of Na + and K + gradients across the sarcolemma by the Na + /K + pump, critical for Na + /Ca 2+ exchanger function. The SR Ca 2+ pump, responsible for refilling the internal Ca 2+ store during each beat, also requires a significant amount of energy to operate in muscle 25 . Although SAN myocytes do not have as many sarcomeres as atrial and ventricular myocytes, they still consume ATP during crossbridge cycling. Furthermore, the constant generation of cytosolic adenosine monophosphate (cAMP), fulfilled by the plasma membrane-associated enzyme adenylate cyclase, is an energetically demanding process that consumes significant amounts of ATP 26  To sustain the energetic requirements of SAN myocytes, the SAN artery, a branch of the right coronary artery, delivers oxygenated blood to the node. During diastole, when the compressive forces of the myocardium are at their lowest, blood moves through SAN artery and capillary network to deliver oxygen to pacemaker cells 31 . Deoxygenated blood exits the SAN tissue via venules and specialized Thebesian veins where blood directly enters the right atrial lumen 32 . In humans, the SAN artery has multiple anatomical variations 33 and damage to this artery can be pro-arrhythmogenic 34 . Despite the high metabolic demand of the SAN, a detailed understanding of the vascular network that supports the node for aerobic respiration is limited.
In this study, we used a combination of imaging approaches to generate high-resolution 3dimensional (3D) maps of SAN vascular anatomy and its relationship to SAN myocytes across the entire node. Analysis of these images show that vascular density is high and myocyte-tovessel distances are relatively short in the superior region of the node. However, vascular density decreases, sharply increasing myocyte-to-vessel distances towards the inferior region of the node. Using whole-cell patch-clamp techniques and Ca 2+ imaging of isolated SAN myocytes from densely (i.e., superior) and sparsely (i.e., inferior) vascularized sections of the SAN, we discovered that superior-derived SAN myocytes are intrinsically capable of firing higher frequency APs versus inferior-derived SAN myocytes. Based on these data, we posit that blood supply could be a determinant in sculpting the firing properties of SAN myocytes, effectively determining the dominant pacemaking site within the SAN.

Animals
Male wildtype C57BL/6J mice (The Jackson Laboratory) aged between 6-14 weeks were used in this study. Mice were euthanized with a single, intraperitoneally administered lethal dose of sodium pentobarbital (250 mg/kg). All mice were maintained, and experiments conducted, in accordance with the University of California, Davis Institutional Animal Care and Use Committee guidelines.

Immunohistochemistry image acquisition and analysis
Imaging of fixed SAN tissues was performed using an Olympus FluoView FV3000 confocal laser scanning microscope with an Olympus APON (40X, numerical aperture = 1.3) oilimmersion lens and acquired using FluoView software. Images were acquired at 0.5 µm/pixel with a field of view of 320 µm. Image volumes consisted of between 80-120 z-planes. In order to image the entirety of each SAN preparation, a mosaic of 3-dimensional image volumes were acquired. All image analysis was performed using a combination of Imaris (9.7, Oxford Instruments), ImageJ/Fiji, and Python (3.7). Both vascular (Cluster of differentiation 31-positive cells, CD31 + ) and myocyte (HCN4-positive, HCN4 + ) image volume color channels were thresholded to generate 3-dimensional binary masks. When necessary, CD31 + positive endocardium was manually removed from the image volumes using Imaris. Vascular and myocyte fractional volumes were calculated from 3-dimensional binary masks with custom written Python scripts. With fractional volume = (number of pixels in the 3-dimensional mask) / (total number of pixels in the volume). The 3-dimensional binary masks were used to calculate the 3-dimensional distance from each HCN4 + pixel to the nearest CD31 + pixel (a Euclidean distance transform, using the SciPy Python package). 3-dimensional masking of HCN4 + and CD31 + volumes was repeated using either Imaris or Python and obtained similar results. For 3dimensional masking in Python, The Allen Institute for Cell Science -Segmentation package (https://github.com/AllenInstitute/aics-segmentation) was used. For 3-dimensional vascular segmentation, the total vessel segments and length per field of view were calculated using Imaris. Electrophysiological signals were recorded using a Multiclamp 700B amplifier controlled with pCLAMP 10 software (Molecular Devices). Cells were initially patched in voltage-clamp and upon the formation of a giga-ohm seal, a short period of time was required (1-5 minutes) until amphotericin B permeabilized the membrane, allowing electrical access to the cell. Accessed cells were transitioned to current-clamp mode to record membrane voltage. All recordings were acquired at 10 kHz.
All current-clamp recordings were analyzed using custom Python scripts for the detection and analysis of SAN myocytes AP parameters as previously described 38 . The take-off potential was defined as the membrane potential when the first derivative of voltage with respect to time (dV/dt) reached 10% of its maximum value. The take-off potential was used as the time of an AP and to calculate both AP frequency (Hz) and inter-AP-interval (ms). The maximum depolarizing potential was defined as the most negative membrane potential between APs. The diastolic duration was defined as the interval between the maximum depolarizing potential and take off potential. The early diastolic depolarization rate was estimated as the slope of a linear fit between 10% and 50% of the diastolic duration and the early diastolic duration was the corresponding time interval.
To detect subthreshold membrane fluctuations, we used the "threshold search" function of Clampfit. A subthreshold voltage fluctuation was defined as an event that had an amplitude larger than the mean + 3 standard deviations of baseline voltage and a duration of at least 1 ms, but failed to exceed the experimentally determined AP take-off potential (i.e., AP threshold) of -29 mV.
After plating single SAN myocytes on the recording chamber, cells were loaded with the membrane-permeable acetoxymethyl-ester form of Fluo-4 (Fluo-4 AM, Invitrogen) for measurement of [Ca 2+ ] i . Fluo-4 (5 µM) was mixed with Tyrode solution and introduced to plated cells during Ca 2+ reintroduction (described above). Cells were loaded for 10 minutes, followed by a de-esterification period of 10 minutes. Temporal fluorescence fluctuations caused by [Ca 2+ ] transients were imaged using the line-scan mode on an Olympus FluoView 3000 confocal microscope with an Olympus APON (60X, numerical aperture = 1.49) oil-immersion lens. Fluo-4 was excited with a 473 nm solid-state laser. Line-scan images were analyzed using ImageJ. Background subtracted fluorescence signals were normalized by dividing fluorescence at each point (F) with the baseline fluorescence (F 0 ). Ca 2+ sparks were identified automatically using the SparkMaster plugin of ImageJ 39 .

Statistics
Parametric, normally distributed data are presented as mean  standard error of the mean. The coefficient of variation was calculated by dividing the standard deviation by the mean of the data set. In the case of nonparametric data, median and range values are reported. Group comparisons were made using a student's T-test or ANOVA followed by a Tukey multicomparison test for parametric data. For non-parametric data, Mann-Whitney tests were implemented. A p value of less than 0.05 was considered significant. The asterisk symbol is used in the figures to illustrate the significant difference between groups where NS (i.e., not significant) represents p >0.05; * represents p ≤ 0.05; ** represents p ≤ 0.01; and *** represents p ≤ 0.001.

The distribution of vessels and SAN myocyte varies along the SAN
We investigated the distribution of vasculature and pacemaking myocytes along the SAN. To do this, we fixed and double labeled SAN tissue with an HCN4 antibody to reveal SAN myocytes and a CD31 (i.e., cluster of differentiation 31 protein) antibody as a marker of endothelial cells (n = 6 SANs). We then used laser scanning confocal microscopy to image sub-micron resolution 3dimensional stacks of HCN4-positive (HCN4 + ) SAN myocytes and CD31-positive (CD31 + ) endothelial cells within the entire SAN ( Figure 1a).
Our initial inspection of these image volumes suggested regional variations in the density of both HCN4 + SAN myocytes ( Figure 1b) as well as CD31 + vessels ( Figure 1c). These regional heterogeneities are more evident in orthogonal images of SAN myocytes and endothelial cells from the superior and inferior sections of the SAN (Figure 1d). Note that in the superior SAN, vessels are embedded within SAN myocytes and become progressively segregated in the inferior SAN.
We zoomed into regions of the superior and inferior SAN to get a closer look at myocytes and vessels in these regions of the node. These images further suggest that myocyte and vessel density is highest in the superior SAN ( Figure 1 e/f), progressively decreasing in the inferior SAN ( Figure 1 h/i). These regional variations in myocyte distribution seems to be a unique feature of the SAN. Images of vessels in the superior vena cava ( Figure 1g) and right auricle ( Figure 1j) surrounding the SAN show that vessel density is homogeneous.

The distance between SAN myocytes and vessels varies along the SAN
We performed a detailed analysis of HCN4 + myocytes and CD31 + vessels along the entire length of the SAN ( Figure 2). To determine vascular and myocyte densities, we generated 3dimensional masks of CD31 + endothelial cells and HCN4 + myocytes (see Materials and Methods). Masks of SAN myocytes and vessel masks were used to determine the distance of SAN myocytes to their nearest vessel. An example of this workflow is shown in Figure 2a.
To quantify the density of myocytes and vessels from the superior to inferior regions of the SAN, we calculated the fractional volume occupied by the CD31 + and HCN4 + masks and generated spatial heatmaps encompassing the entire SAN ( Figure 2b). We found the fractional volumes of both HCN4 + and CD31 + labelling decreases from superior to inferior regions.
What would the functional implications of this difference in vessel and myocyte density between superior and inferior SAN regions have for resident SAN myocytes? To address this, we calculated the 3-dimensional distance from each HCN4 + pixel to its nearest CD31 + positive pixel (a 3-dimensional Euclidean distance transform). This is an important parameter determining myocyte-vessel exchange and communication 31,40,41 . Our analysis revealed increasing distances between myocytes and their neighboring vessels traversing from the superior to the inferior portions of the SAN (Figure 2c). In some regions of the inferior node, myocyte-vessel distance could exceed 10-15 μm.
We set out to develop a measure that would allow us to quantitatively segregate superior and inferior regions of the SAN. Fitting the mean myocyte-vessel distance with a sigmoidal function ( Figure 2c, blue line) yields a separation of these two regions. We used this metric to separate superior and inferior regions for the remaining analysis.
To further reveal differences between the superior and inferior regions of the SAN, we performed additional analyses across a number of fixed-tissue SAN preparations and segregating each SAN imaging dataset into corresponding superior and inferior regions. A cumulative distribution histogram plotted for a single SAN preparation illustrates the difference in proximity of myocytes to vessels (top, Figure 2d). Indeed, we found that superior SAN myocytes were significantly closer to neighboring vessels (bottom, Figure 2d; superior 2.54 ± 0.4 µm; inferior 4.69 ± 1.8 µm, p = 0.016, n = 6 SANs).
Next, we quantified the number of vessel segments and total vessel length and compared the superior and inferior regions within each fixed SAN tissue. We averaged the number of vessel segments and total vessel length between different fields of view between each region and performed a pairwise analysis. This analysis reveals a significantly higher number of vessel segments in the superior versus inferior regions (Figure 2e; superior 1101 ± 126 segments; inferior 269 ± 87 segments, p = 0.0052, n = 4 SANs, fields of view = 9 superior and 12 inferior regions). Finally, we find that the superior region has a significantly higher total vessel length versus the inferior (Figure 2f; superior 12.6 ± 1.2 mm; inferior 5.7 ± 1.3 mm, p = 0.0001, n = 4 SANs, fields of view = 9 superior and 12 inferior regions).
Taken together, these results indicate a progressive decrease in both vessel and SAN myocyte densities along the primary axis of the SAN from the superior to the inferior venae cava. The functional effect is that SAN myocytes in the superior versus inferior region are in closer proximity to their neighboring vessels. This could impact the range of physiological signatures of SAN myocytes with superior cells receiving substantially more metabolic support via greater access to blood supply.

Regional variations in vascular density along the SAN are not due to gross vascular topological features
To determine if the observed differences in vascular density between the superior and inferior regions of the SAN were due to gross vascular topological features, we examined the branching patterns of the main SAN artery as it bifurcated into arterioles and then capillaries.
The main SAN artery arborizes in multiple locations as it descends the node giving rise to intermediate arterioles and finally capillaries (Figure 1c, f, i). Interestingly, the diameter of the primary SAN artery (first order; 1 o ) was larger in the superior versus inferior regions (Figure 2g; superior 30.62 ± 4.3 µm, inferior 23.37 ± 3.9 µm, p = 0.02, n = 3 SANs). Given the decrease in caliber of the primary SAN artery from superior to inferior regions, we investigated if the secondary, tertiary, and quaternary vascular branches also varied between the superior and inferior regions. To do this, we segregated 2 nd order (2 o  Together, this analysis indicates that the gross vascular topological features does not differ between the superior and inferior SAN and in turn does not contribute to regional variations in vascular density.

Regional variations in AP frequency SAN
Our data describing the SAN vasculature and its relationship to SAN myocytes have important functional implications. For example, the intrinsic electrical signaling modalities of SAN myocytes in the superior could be fundamentally similar but vary in vivo in response to differential blood supply. Alternatively, the firing properties of SAN myocytes themselves may vary along the SAN. One plausible hypothesis is that the superior SAN is populated by myocytes that have a higher intrinsic capacity to fire APs than the inferior SAN. In this model, regional variations in AP firing patterns in vivo would likely be due to differences in intrinsic excitability as well as blood supply.
As a first step to test these hypotheses, we developed protocols to identify, dissect, and enzymatically dissociate myocytes from the superior and inferior SAN. We used our vessel and myocyte density data (see Figures 1 and 2) to anatomically define the superior and inferior SAN. We chose a region of approximately 2-3 mm 2 for the superior SAN and similar area for the inferior SAN and dissociated cells from each region individually.
We performed whole-cell perforated-patch current-clamp recordings in individually isolated superior and inferior SAN myocytes ( Figure 3). Figure 3a shows single AP records from representative superior and inferior SAN myocytes. These records show key differences in the waveforms of the APs of superior and inferior SAN myocytes, including early diastolic duration and early diastolic duration rates. A summary of AP parameters is included in Table 1.
We found that the frequency of spontaneous APs was significantly higher in the superior versus inferior SAN myocytes (superior: 4.04 ± 0.36 Hz, n = 1374 APs from from 10 cells; inferior: 1.44 ± 0.33 Hz, n = 319 APs from 7 cells; p = 0.0028). Although the take off potential of APs (mV) was similar in these cells, the early diastolic depolarization rate was nearly 2.3-fold faster in superior than inferior SAN myocytes. This is reflected in the fact that the early diastolic duration in superior myocytes was only about 0.16 times that in inferior SAN myocytes. These data indicates that, in addition to having higher spontaneous firing frequencies, superior SAN myocytes exhibit a number of intrinsic membrane properties making them more excitable than their inferior counterparts.
In addition to APs, we observed subthreshold voltage fluctuations in superior and inferior SAN myocytes (Figure 3b). The differences between superior and inferior SAN myocytes, could be due to a combination of factors, which include regional variations in input resistance and the magnitude of spontaneous currents during diastole in these cells.
We next asked if the occurrence of subthreshold voltage fluctuations followed any notable periodicity or were random by examining the distribution of the timing between these signals. This analysis revealed that the inter-subthreshold voltage fluctuations intervals in both superior and inferior SAN follow a normal distribution and could each be fit with a Gaussian function (superior: center at 333 ms and width of 234.7 ms; inferior: center at 414 ms and width of 217 ms) (Figure 3d). Plotting each inter-subthreshold voltage fluctuations (Inter-subthreshold voltage fluctuations i ) interval against the following one (Inter-subthreshold voltage fluctuations i+1 ) (i.e., "joint interval plot" 42 ) showed a large degree of dispersion of the data, suggesting independence in the timing between successive subthreshold voltage fluctuations ( Figure 3e). These data are consistent with the hypothesis that subthreshold membrane fluctuations occurrence is stochastic.
A similar analysis was performed for APs to quantify the degree of periodicity of AP firing in superior and inferior SAN myocytes. As part of this analysis, we also determined mean and standard deviation of the inter-AP-interval and used them to calculate the coefficient of variation (Coefficient of variation = standard deviation/mean) (Figure 4). The rationale for calculating the coefficient of variation was to provide a quantitative measure of inter-AP-interval dispersion (i.e., periodicity). The higher the coefficient of variation, the lower the periodicity of APs. This analysis revealed that APs in the superior region were more periodic with a significantly smaller coefficient of variation than inferior cells (superior inter-AP-interval coefficient of variation = 0.46 ± 0.14, n = 1374 APs from 10 cells; inferior inter-AP-interval coefficient of variation 0.96 ± 0.18, n = 319 APs from 7 cells; p = 0.02). It is important to note that a coefficient of variation of 1 (i.e., mean = standard deviation) is typically associated with stochastic events. Thus, inferior SAN myocytes seem to fire APs randomly.
To further quantify the diversity in periodic APs between superior and inferior SAN myocytes, we generated joint inter-AP-interval plots. SAN myocytes with a regular tonic firing mode have a joint-ISI plot with a tight cluster of intervals, suggesting a high degree of periodicity and dependency in the timing of successive APs (Figure 4a). The mean inter-AP-interval was 243.5 ms, standard deviation was 103.7 ms, and coefficient of variation in these cells was 0.37. SAN myocytes with irregular AP firing patterns had joint inter-AP-interval plots with a higher degree of data dispersion than in regularly firing tonic cells (Figure 4b). Accordingly, for these cells, the mean, standard deviation, and coefficient of variation of the inter-AP-interval were 1003.2 ms, 1222.3 ms, and 1.23, respectively. The combination of high highly dispersed joint inter-APinterval plot and high coefficient of variation value suggest that myocytes with irregular firing patterns are stochastic.
In sharp contrast to regular and irregular AP firing cells, SAN myocytes that generated APs in bursts had an 'L-shaped' distribution. This is due to alternating long periods of AP silence with brief periods of higher frequency APs. The mean, standard deviation, and coefficient of variation of the inter-AP-interval of 'burst' firing cells was 703.6 ms, 1251.4 ms, and 1.62, respectively.
Our analysis indicates that the proportion of cells with specific intrinsic-excitability phenotypes varied between the superior and inferior SAN ( Figure 5). The pie charts in Figure 5a illustrate that the superior SAN is populated with fewer cells that are electrically silent or produce subthreshold voltage fluctuations and more cells capable of producing APs than inferior SAN. Importantly, most cells in the superior SAN are regularly firing tonic APs while inferior SAN myocytes are mostly irregular, random AP firing modality (Figure 5b).

SR Ca 2+ release is necessary for AP firing but has no direct influence on subthreshold voltage fluctuations in superior and inferior SAN myocytes
We tested the hypothesis that SR Ca 2+ release contributes to differential electrical activity in superior and inferior SAN myocytes. To do this, we recorded membrane potential before and after the application of the SR Ca 2+ pump inhibitor thapsigargin (5 μM) 43 . We found that thapsigargin eliminated APs in all superior (n = 5 cells) and inferior (n = 4 cells) SAN myocytes. Interestingly, inhibiting the SR Ca 2+ pump differentially altered the resting potential of these cells. Superior SAN myocytes tended to depolarize (observed in 4 out of 5 cells) by +20 ± 9 mV (Figure 6a), while inferior cells kept their resting potential unaltered or hyperpolarized by a few millivolts (change in membrane potential = -1.5 ± 1 mV) (Figure 6b). SAN myocytes was unchanged (i.e., control = 6.7 ± 1.6 versus thapsigargin = 6.0 ± 2.2 mV, p = 0.46) by exposure to thapsigargin.
Removing thapsigargin from the bath restored control AP firing frequency in superior and inferior SAN myocytes, on average, within 107 ± 11 s and 116 ± 14 s, respectively. This is interesting, as thapsigargin is virtually a non-reversible inhibitor of SR Ca 2+ pump 43 , suggesting a relatively high turnover rate of SR Ca 2+ pump protein in these cells. These data suggest that SR Ca 2+ release is necessary for AP firing, but not subthreshold voltage fluctuations in superior and inferior SAN myocytes.

Diverse global and local Ca 2+ signaling modalities in superior and inferior SAN myocytes
Next, we visualized superior and inferior SAN myocytes loaded with the fluorescent Ca 2+ indicator fluo-4-AM (Figure 7). We show the spatially-averaged, whole-cell [Ca 2+ ] i records from representative superior ( Figure 7a) and inferior (Figure 7b) (Figure 7d).
To determine the rate of clearance of intracellular Ca 2+ , we determined the time to 50% of the [Ca 2+ ] i transient (Figure 7e). We found that the full duration at half maximum amplitude of the global [Ca 2+ ] i transient of superior SAN cells (median = 291 ms; mean = 337 ± 81 ms) was longer than that of inferior SAN myocytes (median = 193 ms; mean = 196 ± 32 ms; p = 0.04).
Next, we examined diastolic Ca 2+ sparks before diastolic depolarization in both superior and inferior SAN myocytes. Figure 8a  Ca 2+ spark rate ranged from 0.6 to 10.9 events/100 µm*s and 0.3 to 40.7 events/100 µm*s in superior and inferior SAN myocytes, respectively (Figure 8d). The median Ca 2+ spark rate was 2.8 (mean = 4.7 ± 1.5; n = 7 cells) in superior and 2.1 events/100 µm*s (mean = 8.8 ± 2.9; n = 16 cells) in inferior SAN myocytes (p = 0.19). These Ca 2+ spark amplitude data indicate that SR Ca 2+ release sites vary in activity and amplitude between regions of the node, likely reflecting differences in ryanodine receptor cluster sizes and/or SR Ca 2+ load (i.e., driving force).
The spatial spread (i.e., full width at half maximum amplitude; Figure 8e) and duration (i.e., full duration at half maximum amplitude; Figure 8f) was similar in superior and inferior SAN myocytes (p = 0.67). Together, this analysis suggests that the geometries as well as clearance mechanisms of Ca 2+ spark are similar in superior and inferior SAN myocytes.
Coupling of global and local Ca 2+ signals to changes in membrane potential in superior and inferior SAN myocytes A comparison of our electrophysiological and Ca 2+ data suggests that there is a partial alignment of these signals. For example, inferior SAN myocytes have a single population of subthreshold voltage fluctuations (Figure 3c), but a bimodal Ca 2+ spark amplitude distribution ( Figure 8c). Furthermore, although superior and inferior SAN myocytes have a subgroup of Ca 2+ spark with similar amplitudes, their subthreshold voltage fluctuations are significantly different. Similarly, while there were more superior cells with higher frequencies of whole-cell [Ca 2+ ] i transients than inferior cells, we did not detect the same level of segregation or similarity in frequencies in APs (Figure 5b). This could be due to some single [Ca 2+ ] i transients caused by a train of APs. Together, these data suggest a lack of 1:1 correspondence of APs and wholecell [Ca 2+ ] i transients as well as potential variations in coupling strength between Ca 2+ sparks and changes in membrane potential between superior and inferior SAN myocytes.
To address this issue, we simultaneously recorded membrane potential and [Ca 2+ ] i from superior and inferior SAN myocytes to examine the relationship between the occurrence of Ca 2+ sparks and the resultant impact on membrane potential depolarization. We performed dualrecordings from SAN myocytes by simultaneously monitoring [Ca 2+ ] i with confocal line-scans of cells loaded with Fluo-4 paired with whole-cell current clamp recording (Figure 9a). Visually inspecting these dual-recordings of Ca 2+ sparks and whole-cell membrane potential we did not see pronounced effects on membrane potential before, during or after individual Ca 2+ sparks ( Figure 9b). To test if the occurrence of Ca 2+ sparks had an impact on membrane potential depolarization, we quantified the relationship between Ca 2+ spark amplitude and the resultant membrane potential depolarization and found no significant correlation (Figure 9c; superior r 2 = -0.068, n = 4 cells, p = 0.93; inferior r 2 = 0.453, n = 6 cells, p = 0.37).
If individual Ca 2+ spark produce little to no membrane potential changes between APs, is the pattern of Ca 2+ sparks leading to an AP different between the superior and inferior regions? To examine this, we asked if the number of Ca 2+ sparks or the slope of Ca 2+ sparks approaching an AP was different between superior and inferior SAN cells (Figure 9d). Interestingly, we found significantly more Ca 2+ sparks leading to an AP in the superior versus the inferior SAN ( Figure  9e; superior 4.5 ± 2.8 Ca 2+ sparks, inferior 3.4 ± 1.7 Ca 2+ sparks, p = 0.006). Likewise, we found a significantly steeper slope in Ca 2+ spark accumulation before each AP in superior versus inferior SAN myocytes (Figure 9f; superior 2.0 ± 2.1 Ca 2+ spark slope; inferior 1.2 ± 1.1 Ca 2+ spark slope, p = 0.007).
This analysis indicates that while the coupling of individual Ca 2+ sparks to membrane potential depolarization is weak throughout the SAN, there is a regional difference in the number and rate of accumulation of Ca 2+ sparks before an AP. With the superior region receiving a stronger driving force of Ca 2+ sparks leading to AP generation. This could be a mechanism contributing to the higher spontaneous AP frequency of superior versus inferior myocytes.

Discussion
We generated high-resolution 3-dimensional maps of the vasculature and SAN myocytes across the entirety of the mouse SAN. Our data show that vascular and myocyte densities vary regionally along the SAN. The superior SAN is highly vascularized, creating a redundant network whereby all pacemaking myocytes are near at least one vessel. This affords a high surface area for exchange and minimizes the diffusional distances in a region of the node where myocytes have highest intrinsic firing rates. In contrast, vascular and myocyte densities are lowest in the inferior SAN. The inferior SAN is also where myocytes with the lowest inherent AP firing rate are located. Thus, the SAN seems to be designed so that vascular density varies regionally to match the electrical signaling modalities of the myocytes that populate this heterogenous tissue. An important implication of these findings is that regional variations in blood supply could be a key factor determining the site of AP origin in the SAN as well as its firing frequency.
Based on our detailed vascular and myocyte anatomical and functional data, we propose a new, objective identification of the superior and inferior mouse SAN. Our data indicate that SAN myocyte-to-vessel distance is relatively short and remains relatively similar from the initiation point of the SAN artery until about 3.25 mm towards the inferior section of the node. We call this region the superior SAN. We propose that the starting point of the inferior SAN is set at the lowest point of the superior node, where there is a relatively sharp, nearly 1.5-fold increase in myocyte-to-vessel distance. Using this framework, the superior and inferior SAN account for about 60% and 40% of the mouse node, respectively.
An interesting observation in our study was that while the diameter of the SAN artery decreased nearly 20% from the superior to inferior SAN, the diameter of 2 o -4 o vessels did not. Thus, differences in vascular density are not due to fundamental differences in the structure of microvessels, rather their number.
Our work is consistent with a recent study from the Efimov lab 24 involving rat hearts. They found that, although the superior and inferior sections of the node have the capacity to fire action potentials, the superior node has the fastest firing rate and hence functions as the dominant pacemaking center in of the SAN. The rat inferior SAN was capable of becoming the dominant pacemaking center, but only when the superior node was electrically silenced. Although the concept of a functional heterogeneous SAN has been proposed before, our study, in combination with Brennan et al 24 , supports a model where heterogeneity encompasses all elements of nodal function, i.e., vascular, neuronal, and myocyte, but in two distinct regions of the node, at least in rodent hearts.
What is the physiological significance of the reciprocal association of high excitability with high vascular density in the SAN? Our study, in combination with others, provide a potential answer to this difficult question. A reasonable hypothesis is that vessel density and hence blood supply scales up with tissue metabolic needs. It is estimated that a neuron consumes about 7-9 million ATP molecules during a single AP 44,45 . At low firing rates (i.e., 0.1 Hz), ATP levels in neurons remain relatively stable, suggesting a balance between ATP synthesis and hydrolysis. However, increasing AP frequency (>10 Hz) leads to rapid ATP depletion within seconds of the onset of high frequency stimulation as ATP consumption exceeds production 46 .
Our electrophysiological data suggest that the intrinsic AP frequency is about 3.9-fold higher in superior versus inferior SAN myocytes. Although similar ATP consumption studies like those in neurons have not been performed in SAN myocytes, assuming that energetic cost of firing an AP is likely similar in these cells, our data suggest that in terms of APs alone, ATP consumption could be 3.9-fold higher in superior compared to inferior SAN myocytes. Note, however, that a relatively large fraction the inferior SAN myocytes do not produce APs. Thus, the energetic differential between superior and inferior SAN myocytes could be even larger than 3.9-fold in these cells.
Although neurons and SAN myocytes have the capacity to spontaneously fire APs, their energetic needs are likely quite different. As in neurons, these include ion transport across the sarcolemma and intracellular organelles such as the endoplasmic reticulum 25 . Additional ATP is consumed by myosin-actin crossbridge cycling and the constant generation of cytosolic cAMP by adenylate cyclase 26 . The latter is significant as our data show that even the fastest firing superior SAN myocytes produce APs at a lower rate than the heart rate of a mouse, which is about [12][13]48 . This disparity between in vivo and isolated single cell firing rates reflects strong input from the sympathetic arm of the autonomic nervous system, which increases firing rate through the activation of cAMP/protein kinase A signaling. Indeed, recent studies by Brennan et al 24 and Hanna et al 49 suggested that sympathetic innervation and drive is larger in the superior compared to the inferior SAN. Thus, a combination of higher AP firing, Ca 2+ cycling, and contraction rates as well as cAMP production may conspire to produce a much larger energetic load in superior than inferior SAN myocytes. In combination with our vascular mapping, this implies that  as originally proposed Krogh 50 for skeletal muscle  micro-vessel density scales up with the energetic needs of SAN myocytes.
Although performing a detailed analysis of the electrophysiological properties of superior and inferior SAN myocytes was beyond the scope of this study, our data provide insights into this interesting issue. We found that the coupling strength between Ca 2+ sparks and membrane potential is relatively weak, but similar in superior and inferior SAN myocytes. However, the number of Ca 2+ sparks immediately before an AP was larger in superior than inferior SAN myocytes. Thus, increased SR Ca 2+ release could contribute to faster diastolic depolarization rates leading to AP threshold in superior myocytes and hence higher firing rates. In principle, these regional variations in Ca 2+ clock function could be the result of differences in the activity of Ca 2+ channels, SR Ca 2+ pump, ryanodine receptors, and/or Na + /Ca 2+ exchanger. Future experiments should determine whether regional differences in the function and expression of these proteins as well as HCN 51 and Ca 2+ currents contribute to differential electrical firing patterns across the SAN.
An interesting finding in this study is that SR Ca 2+ pump inhibition with thapsigargin eliminated AP firing, but not subthreshold voltage fluctuations in superior and inferior SAN myocytes. In combination with our finding that Ca 2+ sparks are weakly coupled to changes in membrane potential, this suggests that SR Ca 2+ release may be necessary, but not sufficient to trigger an AP in superior and inferior SAN myocytes. Furthermore, our data suggest that subthreshold voltage fluctuations are unlikely produced by Ca 2+ spark-activation of Na + /Ca 2+ exchanger currents. In a follow up study, we will investigate the molecular identity of the proteins underlying subthreshold voltage fluctuations. In this model, all SAN myocytes function like a bi-stable system, alternating from the maximum diastolic potential to a full AP. SAN myocyte in the sparsely vascularized inferior section of the node fire stochastic subthreshold voltage fluctuations and sporadic APs. These events do not lead to periodic pacemaking on their own. However, when coupled to a more periodic voltage oscillators such as superior SAN myocytes, subthreshold voltage fluctuations and rare APs integrate (i.e., resonance) and increase the probability that superior SAN myocytes reach AP threshold. Thus, inferior SAN myocytes likely increase the strength and periodicity of superior SAN activation and hence pacemaking activity. Stochastic resonance models have also been implicated in the periodic activation of APs in neurons 53,54 . A key new insight provided by our data is that the capacity of specific regions of the node to operate as periodic oscillators or stochastic signal generators could be critically dependent on vessel density and hence blood flow. Thus, blood supply could effectively determining the dominant pacemaking site within the SAN.

Acknowledgements
This work was supported by grants from the US National Institutes of Health to LFS (1R01HL144071 and 1OT2OD026580) and the American Heart association to RHC (16SDG27130006).
LG was supported by training grant 5T32HL086350.

Disclosures
The authors declare no financial interests in the publication of this work or any other conflict of interests.

Data Availability Statement
The data underlying this article will be shared on reasonable request to the corresponding author.