Hypoglycaemia combined with mild hypokalaemia reduces the heart rate and causes abnormal pacemaker activity in a computational model of a human sinoatrial cell

Low blood glucose, hypoglycaemia, has been implicated as a possible contributing factor to sudden cardiac death (SCD) in people with diabetes but it is challenging to investigate in clinical studies. We hypothesized the effects of hypoglycaemia on the sinoatrial node (SAN) in the heart to be a candidate mechanism and adapted a computational model of the human SAN action potential developed by Fabbri et al., to investigate the effects of hypoglycaemia on the pacemaker rate. Using Latin hypercube sampling, we combined the effects of low glucose (LG) on the human ether-a-go-go-related gene channel with reduced blood potassium, hypokalaemia, and added sympathetic and parasympathetic stimulus. We showed that hypoglycaemia on its own causes a small decrease in heart rate but there was also a marked decrease in heart rate when combined with hypokalaemia. The effect of the sympathetic stimulus was diminished, causing a smaller increase in heart rate, with LG and hypokalaemia compared to normoglycaemia. By contrast, the effect of the parasympathetic stimulus was enhanced, causing a greater decrease in heart rate. We therefore demonstrate a potential mechanistic explanation for hypoglycaemia-induced bradycardia and show that sinus arrest is a plausible mechanism for SCD in people with diabetes.


Introduction
Blood glucose levels are elevated in both type 1 and type 2 diabetes. Insulin is used to therapeutically lower high glucose in individuals with type 1 and long-standing type 2 diabetes; however, this commonly results in iatrogenic hypoglycaemia, one major barrier in achieving and maintaining normal glucose control in diabetes [1]. Insulin causes glucose to fall, but in individuals with diabetes, insulin levels are reduced due to auto-immune damage to the pancreatic β cells (type 1 diabetes) or a combination of impaired insulin secretion and action (type 2 diabetes). In both types of diabetes, exogenous insulin is used to lower glucose therapeutically but due to the limitations of current insulin delivery, insulin replacement is unphysiological. It is challenging to achieve stable basal concentrations and those with insulin-treated diabetes are subject to both hyper-and hypoglycaemia. Attempts to mimic normal glucose levels through intensive insulin therapy increase the risk of severe hypoglycaemic episodes. The effects of hypoglycaemia on the cardiovascular system have only recently emerged and severe hypoglycaemia has been associated with increased mortality in randomized trials [2,3]. In one of the trials, intensive glucose therapy resulted in a sixfold increase in severe hypoglycaemia [2]. Over one-third of deaths were attributed to sudden cardiac death (SCD), raising the possible link between hypoglycaemia and cardiac arrhythmias.
Indeed, there is growing evidence in recent clinical and animal studies that hypoglycaemia might contribute to malignant cardiac arrhythmias leading to SCD. Severe hypoglycaemia in response to glucose therapy was associated with an increased risk of arrhythmic death in individuals with type 2 diabetes [3]. Spontaneous nocturnal hypoglycaemia was also implicated in the deaths of young people with type 1 diabetes. The term 'dead in bed syndrome' was introduced to describe these unexpected deaths [4], which were attributed to SCD. Subsequent studies reported a relatively large proportion of unexplained deaths in this population [5,6]. Although rare in absolute numbers, SCD is 10-fold higher in young people with type 1 diabetes compared to those without diabetes [7]. In observational studies, we and others reported an increased risk of arrhythmias, including slow heart rate, bradycardia, during spontaneous hypoglycaemia versus normoglycaemia (NG) [8][9][10][11]. These arrhythmias could be related to abnormal cardiac repolarization observed in spontaneous [12] and experimental hypoglycaemia studies, where hypoglycaemia was induced by an infusion of insulin [13][14][15]. Experimental hypoglycaemia in rats can cause lethal cardiac arrhythmias, although at glucose levels significantly lower compared to those in human studies [16].
Possible hypoglycaemia-induced proarrhythmic mechanisms include the direct effect of low extracellular glucose on the myocardial cell, which can block potassium channels that carry an outward current during repolarization [17]. This causes a prolongation in the cardiac action potential (AP) as well as lengthening of the QT interval duration in the ECG through a mechanism identical to that in proarrhythmic medications [18]. In addition, the counterregulatory sympathoadrenal activation in response to hypoglycaemia causes a release of adrenaline [19], which in turn lowers serum potassium concentration [20] and increases intracellular calcium concentration [21]. Combined, these can cause further QT interval prolongation and calcium overload and lead to proarrhythmogenic electrical instability [21]. Abnormal autonomic responses which are common in diabetes could cause further risk in these individuals.
The autonomic response to falling glucose levels includes sympathetic and parasympathetic activation, with the latter limiting the increase in heart rate during hypoglycaemia. However, bradycardia was also observed during spontaneous [22] and experimental hypoglycaemia [23]. In our observational studies, we detected increased rates of bradycardia during hypoglycaemia compared to NG in some individuals with diabetes [8,10,11]. Bradycardia was also common during experimental hypoglycaemia in rats [24]. These observations suggest that bradycardia might be an additional risk factor during hypoglycaemia, possibly by the development of early afterdepolarizations leading to ventricular tachycardia, particularly in the presence of hypokalaemia [25]. Our focus in this study was therefore on the effect of hypoglycaemia on the sinoatrial node (SAN), the heart's natural pacemaker.
Hypoglycaemia studies present challenges to investigators. While hypoglycaemia and cardiac arrhythmias are relatively common, SCD is rare, suggesting the involvement of additional factors, including variable and dysfunctional autonomic responses and impaired awareness of hypoglycaemia. Observational clinical hypoglycaemia studies are limited by an inability to measure electrolytes and catecholamines during spontaneous episodes of hypoglycaemia. Experimental studies on the other hand while more controlled require supra-physiological concentrations of insulin to induce hypoglycaemia, rarely encountered in clinical practice. In addition, a greater fall in potassium in these studies is not representative of what happens in people with diabetes in the 'free living' condition. Animal models offer an alternative approach but species differences in the hormonal responses to hypoglycaemia and ensuing cardiac electrophysiological responses mean there are limitations in drawing conclusions from these data.
In this study, we have therefore used a computational model to overcome some of these limitations. There are currently no computational models that describe the electrophysiological responses to hypoglycaemia, and in this study, we adapted a biophysically detailed computational model of a human SAN cell electrophysiology to examine the effect of hypoglycaemia on pacemaker rate.

Model of the human sinoatrial node action potential
We used a computational model of the human SAN AP, developed by Fabbri et al. [26]. This model was in turn adapted from a computational rabbit SAN model, developed by Severi et al. [27] by incorporating experimental data from human SAN cells. These included updated cell capacitance and dimension data, a new formulation of the I f current, new steady-state activation and updated conductances of delayed rectifier K + currents (I Kr and I Ks ) and added I Kur current formulation. Parameters, for which human experimental data were not available, were fitted using automatic parameter optimization. The modulation of the pacemaker rate in response to autonomic stimulation was similarly adapted from the rabbit model, where acetylcholine-and isoprenaline-induced variations were introduced to mimic the activation of the sympathetic (increased rate) and parasympathetic (decreased rate) nervous system, respectively. Further details about the model are provided in the electronic supplementary material. We ran simulations in MATLAB (The MathWorks, Inc., Natick, MA, USA) vR2017a, using automatically generated code for the Fabbri model obtained from the CellML platform (www. cellml.org). For solving numerical differentiation, we used the 'ode15s' function with relative error tolerance 10 −7 , absolute error tolerance 10 −6 and a maximum step size of 10 −3 s. We ran the simulations for 60 s and used the last beat to calculate the features of AP waveforms and for figures. We confirmed the stability of the main features of the AP waveform (electronic supplementary material, figure S1).

Features of sinoatrial action potential waveform
We characterized the sinoatrial AP waveform using a set of parameters presented by Fabbri et al., and we calculated them based on the guidelines for standardization of AP parameters [28]. The following parameters were calculated: CL-cycle length of the pacemaker activity, defined as the duration between the peaks of two consecutive APs; PP-peak action potential, peak (maximum value) of the AP; APA-AP amplitude, defined as the difference between the maximum and minimum AP potentials; MDP-maximum diastolic potential preceding the PP, defined as the most negative repolarization potential; DDR-diastolic depolarization rate; DDR 100 -diastolic depolarization rate over the first 100 ms of diastolic depolarization; APD 90 -AP duration at 90% repolarization (figure 1).
The main parameter of our interest was CL. We explored the changes induced by hypoglycaemia and the main features of AP responsible for these changes.

Effect of low extracellular glucose-experimental data
In an experimental study, Zhang et al. demonstrated a direct effect of low extracellular glucose concentration (Gluo) on the rapid delayed rectifier potassium current (I Kr ) [17]. Using a whole-cell patch clamp, they showed an impairment in the function of the corresponding human ether-a-go-go-related gene (hERG) channel, proportional to the reduction in Gluo. Zhang et al. used human embryonic kidney cells (HEK293 cells), and in this study, we assume that electrophysiological properties and responses of the hERG channels are similar in the heart, where hERG is highly expressed. Zhang et al. demonstrated two separate mechanisms by which reduced Gluo alters the function of the hERG channel: a reduction in I Kr current density and a shift in its current-voltage (I-V) relationship. We digitized their graphs to numerically quantify the extent of changes. We further describe the effect size for both mechanisms and how we incorporated them in our model.

Reduced I Kr current density
Zhang et al. presented experimental data for I Kr I-V relationships at various Gluo: Gluo = 5 (NG), 2.5, 1 and 0 mmol l −1 . They demonstrated a 10-20% reduction in current density at Gluo = 2.5 mmol l −1 compared to 5 mmol l −1 , a 20-30% reduction at 1 mmol l −1 and a 50% reduction at 0 mmol l −1 . The current density decreased proportionally with decreased Gluo, but the reduction was not constant across the voltage range. In our model, we incorporated the drop in current density by decreasing the conductance (g Kr ) of I Kr . For easier interpretation, g Kr in our model represents normalized conductance, a multiplier of maximum conductance (4.2 nS in the Fabbri model). We decreased g Kr from 1 (maximum conductance, no effect of glucose) to 0 to investigate the overall response of the model. We then varied g Kr within the [1, 0.7] range (100 to 70% of maximum conductance) to consider only physiologically relevant glucose concentrations, based on the above data.

Shift in the I Kr current-voltage relationship (IVshift)
Zhang et al. demonstrated that the I Kr I-V relationship shifted towards more negative potentials with reduced Gluo. Specifically, they reported a −10 mV shift at Gluo = 0 versus 5 mmol l −1 . There was a similar −10 mV shift in normalized conductance. I-V relationships were not shown for other Gluo, but they performed additional experiments under conditions with corrected osmolarity at Gluo = 1 mmol l −1 . In these experiments, they measured a −3 mV shift in normalized conductance at Gluo = 1 versus 5 mmol l −1 and for the purpose of our study, we assumed that the shift in I-V curves is similar. This is since the shifts in the I-V relationship and conductance were similar in the experiment without correction for osmolarity. To model physiologically relevant effects of hypoglycaemia, we varied IVshift within [−3, 0] mV range but we explored shifts below −3 mV to test the overall response of the model.
We extracted all relevant data by digitizing the graphs, and we measured the shifts at 50% I-V curve or conductance amplitude. In the Fabbri model, I Kr is formulated as a function of membrane potential by combining the activation and inactivation kinetics. To model the shift in the I-V curve, we incorporated IVshift into the membrane potential variables of all equations describing the activation and inactivation relationships and their corresponding time constants.

Latin hypercube sampling
We used Latin hypercube sampling for simulations where two or more input parameters were randomly selected from a predefined range. The sampling was performed by using the 'lhsdesign' algorithm in MATLAB. The variable input parameters included g Kr , IVshift and extracellular potassium concentration (Ko).

Outline of the study
We investigated the individual effect of each parameter g Kr , IVshift and Ko on CL and features of the AP waveform. For each parameter, the modelling was performed for physiologically relevant values as well as outside this range to test the behaviour of the model and to identify input values that produce meaningful APs. We also looked at the combined effects of these parameters when their values were within physiologically relevant ranges. Using Latin hypercube sampling, we varied g Kr (g Kr is a multiplier of maximum conductance) between 0.7 and 1, IVshift between −3 and 0 mV and Ko between 3.0 and 4.0 mmol l −1 . Ko at 4.0 mmol l −1 represents normal Ko while Ko = 3.0 mmol l −1 indicates moderate hypokalaemia. Finally, we checked the combined effect of these parameters with additional activation of either the sympathetic or parasympathetic autonomic system. The current model allows only one of them to be activated at any time.

Reduced conductance of the hERG channel
We reduced g Kr (g Kr is a multiplier of maximum conductance) from 1 (maximum conductance) to 0 in steps of 0.05. The model failed to generate spontaneous AP activity at g Kr < 0.11. The relationship between g Kr and CL is presented in figure 2a. CL shortens with decreased g Kr and the relationship is approximately linear down to g Kr = 0.4, followed by a plateau. The minimum CL at g Kr = 0.11 is 0.47 s (−49%). The vertical blue lines indicate the [0.7, 1] range, which represents our estimated reduction in g Kr at clinically relevant reduced glucose concentrations. CL at g Kr = 0.7 is 0.75 s (−18%), APD is prolonged (+21%), the MDP potential is elevated (+8%) and the diastolic depolarization rates DDR and DDR 100 are increased (+32% and +27%) (table 1). Examples of the corresponding AP waveforms are presented in figure 2b.

Negative shift in I Kr current-voltage relationship
We introduced a negative shift (towards more negative potentials) in the I-V relationship (IVshift), starting from 0 mV in steps of −0.1 mV. CL exponentially increases with negative IVshift (figure 3a) and the longest valid CL is obtained at −5.4 mV (CL = 5.94 s). The blue vertical lines indicate the [−3, 0] mV range which represents our estimated shift at clinically relevant glucose concentrations. CL at IVshift = −3 mV equals 1.2 s (+30%) and the prolongation is mostly caused by decreased MDP (−8%) and depolarization rates DDR (−43%) and DDR 100 (−15%) (table 1). The corresponding AP waveforms are presented in figure 3b.

Combined effects of reduced g Kr and negative IVshift
We combined the effects of reduced g Kr and negative IVshift by simultaneously sampling g Kr from [0.7, 1] and IVshift from [−3, 0] mV ranges. We used Latin hypercube sampling to generate 300 pairs of random and independent input variables for 300 evaluations of the model. The resulting CLs are presented in figure 4. Red data indicate CL prolongation and blue data CL shortening versus baseline (CL = 0.92 s) (figure 4a). CLs are longest at baseline g Kr (g Kr = 1) combined with maximum IVshift (IVshift = −3 mV) and are shortest at g Kr = 0.7 combined with baseline IVshift (0 mV). The grey dashed line indicates a simultaneous linear change in both parameters. CLs along this line are approximately constant and are slightly longer than baseline CL (approx. 0.93 s, +1%). In figure 4b, a projection of CL values to the g Kr axis shows a gradual shortening of CL with decreased g Kr , and for each g Kr , CL is prolonged with the negative IVshift. The maximum and minimum CLs obtained with these parameters are 1.2 s (+30%) and 0.75 s (−18%), respectively (table 1).

Features of the action potential waveform
The individual effects of reduced g Kr and negative IVshift produce CL changes in opposite directions (figures 2 and 3) and their combined effects cancel each other out to some degree, producing a slight CL prolongation. A similar opposite effect is obtained with other AP features, including APA, MDP and DDR. For example, baseline MDP equals −62.0 mV (table 1) and the combined effects of g Kr = 0.7 and IVshift = −3 mV slightly raise MDP to −60.3 mV (+3%). Individual changes in g Kr and IVshift, however, produce a bigger difference compared to baseline with minimum/maximum MDP values −66.8 mV (−8%) and −56.9 mV (+8%), respectively. On the other hand, PP, APD 90 and DDR 100 are predominantly affected by one of the input parameters but not by the other. For these AP features, the combined changes in g Kr and IVshift produce values outside the range of their individual effects. For example, APD 90 is prolonged at g Kr = 0.7 (APD 90 = 0.17 s, +21%) but with added IVshift = −3 mV, there is a further prolongation (APD 90 = 0.18 s, +29%). Maximum values of PP, Table 1. Features of the AP waveform (rows 4 and onwards) calculated using a combination of input parameters (rows 1 to 3) representing NG, LG and a combination of LG and low potassium. Baseline values of input parameters (italics), their extreme values and their combinations are considered. We changed the baseline Ko = 5.4 mmol l −1 in the Fabbri model to Ko = 4.0 mmol l −1 .

NG
LG LG + low Ko

Effect of low extracellular potassium
With the other input parameters at baseline (NG), we lowered Ko from 5 mmol l −1 in steps of −0.1 mmol l −1 . The relationship between Ko and CL is presented in figure 6a in blue colour. The lowest Ko that produces periodic pacemaker activity is 0.6 mmol l −1 (CL = 1.54 s, +67%). CL is slightly prolonged with the decrease in Ko from 4 mmol l −1 to 3 mmol l −1 (blue vertical lines) and reaches a local maximum of 0.94 s (+2%) at 3.3 mmol l −1 . The prolongation is mainly caused by more negative MDP (−6%) and decreased DDR 100 (−12%) (figure 5c and table 1).

Combined effects of low glucose concentration and potassium concentration
We first combined the effects of low Gluo and Ko by fixing g Kr and IVshift at the maximum of their ranges: g Kr = 0.7, IVshift = −3 mV and by reducing Ko in a stepwise manner. Second, we randomly sampled all three parameters (g Kr , IVshift and Ko) from their predefined ranges.   By combining the effects of randomly sampled input parameters, CL stays within the range obtained from parameters at the extreme ends of their predefined ranges (table 1).

Autonomic modulation
We investigated the effect of autonomic modulation by separately including in the model either sympathetic or parasympathetic activity. In both scenarios, we looked at the combined responses of autonomic activity, LG and low potassium. To model LG, we first fixed g Kr and IVshift at the extreme ends of their ranges. Second, we explored all possible outcomes by randomly sampling all input parameters (g Kr , IVshift and Ko).

Fixed g Kr and IVshift
At NG and Ko = 4 mmol l −1 , an addition of sympathetic activity shortens the CL from 0.92 s to 0.69 s (−25%). In LG, Ko = 4 + Symp.Act.
LG, Ko = 3 + Parasymp.Act.   The AP waveforms at Ko = 3, 3.5 and 4 mmol l −1 are presented in figure 7b. CL is prolonged with reduced Ko, but the prolongation is greater at LG versus NG (red data above blue), with the biggest difference at Ko = 3 mmol l −1 (CL = 1.73 s versus 1.50, +15%). Some combinations of input parameters with added parasympathetic activity do not produce periodic pacemaker activity (data not presented, see below).
To identify input parameters which in combination with parasympathetic activity produce long CL or no pacemaker activity, we divided the parameter space as shown in figure 9a. Red data indicate CL prolongation and blue data CL shortening compared to NG with added parasympathetic activity (CL = 1.27 s). Green data indicate CL > 2 s and purple data indicate invalid CL with no periodic pacemaker activity. At Ko = 4 mmol l −1 , all simulations result in valid periodic activity, with some CLs longer than 2 s (see also figure 8b). These are roughly confined to the area with g Kr > 0.9 and IVshift < −2 mV. At Ko = 3 mmol l −1 , this area extends to g Kr > 0.8 and IVshift < −1 mV (green circles). At Ko = 3.7 mmol l −1 , the model stops producing pacemaker activity around g Kr = 1 and IVshift = −3 mV. At Ko = 3 mmol l −1 , this area extends to g Kr > 0.88 and IVshift < −1.8 mV. Figure 9b shows examples of APs over 40 s simulations. Baseline AP (CL = 0.92 s) is shown in the top row, followed by CL = 2 s in the second row and sinoatrial pauses in rows 3 to 5. The penultimate row shows the longest CL (around 25 s), and the last row shows an example with no periodic pacemaker activity.

Discussion
In this first study investigating the effects of hypoglycaemia on the SAN electrophysiology, we demonstrate the following findings. First, in response to hypoglycaemia, the reduced conductance of the hERG channel and the shift in its I-V relationship produce opposite effects on the pacemaker rate. When combined, they cancel each other out, resulting in a small decrease in heart rate. Second, mild hypokalaemia causes a slight decrease in heart rate. When combined with hypoglycaemia, there is a marked further decrease in heart rate. Third, hypoglycaemia combined with sympathetic activity enhances the effect of sympathetic activity by further increasing the heart rate compared to NG. However, when hypoglycaemia is combined with hypokalaemia and sympathetic activity, the effect of sympathetic activity is diminished, causing a smaller increase in heart rate compared to that at NG with added sympathetic activity. Fourth, hypoglycaemia enhances the effect of parasympathetic activity by causing a greater decrease in heart rate compared to NG. Hypoglycaemia combined with hypokalaemia and parasympathetic activity causes a further marked reduction in heart rate, resulting in sinoatrial pauses, sinus arrest and block of pacemaker activity.
Combined, the above findings indicate that hypoglycaemia tends to decrease the pacemaker rate of the SAN. During experimental hypoglycaemia, where glucose was reduced in a controlled environment using insulin infusion, hypoglycaemia counterregulation caused a significant increase in circulating adrenaline and noradrenaline concentrations together with a significant decrease in potassium in healthy people [15] and in people with type 2 [15] and type 1 diabetes [14]. Despite the significant increase in adrenaline and noradrenaline, there was a relatively small (about 3 bpm) increase in heart rate during hypoglycaemia compared to NG [14,15]. It is therefore possible that decreased glucose and potassium levels counteract the effects of the autonomic nervous system at the level of the SAN. Robinson et al. have reported the results of experimental studies in healthy subjects where they restored normal concentrations of potassium during both NG and hypoglycaemia [14]. They have shown that heart rate did not change during NG after potassium has been added (Ko = 3.6 versus 4 mmol l −1 ). On the other hand, there was an approximate 7 bpm increase in heart rate (61 versus 68 bpm, +11%) when potassium was added (Ko = 3.2 versus 3.9 mmol l −1 ) during hypoglycaemia. These data are in line with our findings that decreased potassium does not significantly affect the heart rate during NG but causes a decrease in heart rate when combined with decreased glucose. For comparable values of Ko, our model shows a 2% higher heart rate during NG and a 9% higher heart rate during hypoglycaemia. In our study, we demonstrate a potential mechanistic explanation for hypoglycaemia-induced bradycardia, which has been observed during clinical episodes in people with type 1 [30] and type 2 diabetes [22]. We have undertaken several observational studies, where we compared rates of cardiac arrhythmia during spontaneous hypoglycaemia versus NG. We observed increased incidences of bradycardia in people with type 1 diabetes [10], with type 2 diabetes [8] and with type 2 diabetes after discharge from ICU [11]. However, increased rates of bradycardia were only detected in a few individuals. We speculated that these effects are idiosyncratic and that mechanisms other than those involved in hypoglycaemia counterregulation might contribute to the development of bradycardia [10,11,31]. These could include autonomic dysfunction and/or abnormal and variable autonomic responses such as those observed during experimental hypoglycaemia, where the early parasympathetic withdrawal was followed by parasympathetic reactivation during sustained hypoglycaemia [15]. Nocturnal hypoglycaemia could also increase the risk as sleep and prone position diminish the sympathoadrenal responses during the night [32]. Our findings that decreased blood glucose and potassium diminish the effect of the sympathetic activity and enhance the parasympathetic activity might further exacerbate the risk of bradycardia in the above cases.
To validate the responses of our model during hypoglycaemia, we used the heart rate differences between matched spontaneous hypoglycaemia and NG episodes from our observational study in people with type 1 diabetes [10]. We chose only nocturnal episodes due to diminished sympathoadrenal responses which in turn also cause a smaller reduction in Ko. Our computational model shows an approximate −1% decrease in heart rate during hypoglycaemia at Ko = 4 mmol l −1 . In our observational study, there was a nonsignificant +2% increase in heart rate during hypoglycaemia versus NG with a mean difference of 1.5 bpm (95%CI −1.3, 4.4 bpm). In a similar study, Koivikko et al. [33] reported a +2% increase. Our computational data are within the confidence intervals of those from the above observational studies and the difference in trends could be partly due to counterregulatory sympathetic activation.
There are several conditions which could combine with the effects of hypoglycaemia to further promote decreased pacemaker rate at the SAN. Channelopathies, such as mutations in the HCN4, SCN5A and KCNQ1 genes, which encode ion currents I f , I Na and I Ks , can lead to reduced pacemaker rate. In addition, remodelling of the I f depolarization current in trained athletes, which causes decreased resting heart rate and potential training-induced bradycardia [34] could lead to increased risk during hypoglycaemia. These cases could be further explored in the computational model of the human SAN as well as in observational or experimental hypoglycaemia studies.
In our model, we show a gradual increase in pacemaker rate with decreased g Kr (figure 2). This contradicts experimental data from animals, where selective blockage of I Kr resulted in decreased rates in rabbit [35] and guinea pig [36] SAN. Due to the lack of human experimental data to confirm our findings, we explored the reduction of g Kr in a computational model of the rabbit SAN. We chose the Severi model [27], which is the parent model of the one used in this study. The AP waveforms and the relationship between g Kr and CL are presented in figure 10. CL is proportionally prolonged with reduced g Kr (figure 10a) and the pacemaker activity is abolished with complete block of I Kr (g Kr = 0, figure 10b). CL prolongation is accompanied by APD prolongation, the elevation of MDP and decreased DDR. This is in line with experimental data from the rabbit SAN [35]. In our human model, a reduction in g Kr is accompanied by APD prolongation and elevation of MDP, as in the rabbit model. In contrast with the rabbit model, DDR increases in the human model (table 1). One of the reasons for differences between the models is the roughly 2.5 times higher resting pacemaker rate in rabbits (CL = 0.92 s in human versus 0.36 s in rabbit). While APD prolongation is similar in both models, this makes up a much bigger proportion of CL in the rabbit AP. In figure 10b, APD prolongations are similar to the prolongations of CL, while the opposite effects of elevated MDP and decreased DDR on CL seem to cancel each other out. By contrast, the duration royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20210612 of diastolic depolarization is much longer in the human model, making the CL more susceptible to changes in DDR.
To further investigate the differences between DDR responses in both models, we compared the major currents responsible for diastolic depolarization at g Kr = 1 versus g Kr = 0.5. In the rabbit model, the amplitude of I f current is relatively high versus the other currents and its amplitude decreases by about 50% with reduced g Kr (electronic supplementary material, figure S2). Similarly, the amplitude of I NaCa markedly decreases. By contrast, the amplitude of the I f current is much smaller in absolute terms as well as relative to the other currents in the human model (electronic supplementary material, figure S3). In addition, the amplitude of I NaCa slightly increases with reduced g Kr . These differential responses in currents, combined with the higher resting rate in rabbits could explain the differences between the g Kr /CL relationship in humans versus rabbits. It is also possible that the decreased pacemaker rate in response to reduced g Kr in rabbits contributes to the reduction of heart rate and occurrence of bradycardia during experimental hypoglycaemia in rats [24].
The key strength of our study is that we explored for the first time the electrophysiological effects of hypoglycaemia on the SAN. Furthermore, to mimic the physiological conditions in humans, we added the stimuli of hypokalaemia as well as sympathetic and parasympathetic stimulation. This approach overcomes some of the technical, biological and ethical constraints of in vivo human and animal studies. One of the limitations of the study is the lack of detailed experimental data on the effect of Gluo on cardiac cell electrophysiology. Zhang et al. demonstrated the main effects of LG on the hERG channel [17], but experimental data at finer glucose increments would greatly improve the understanding of the balance between increased heart rate, caused by reduced conductance of the channel versus slowing of the rate caused by the shift in its I-V properties and how these relate to clinically relevant glucose levels. To the best of our knowledge, there is no evidence of any other cardiac ion channel being affected by hypoglycaemia.
Fabbri et al. optimized their model to reproduce baseline electrophysiological characteristics from a limited number of human experimental data and validated it against baseline heart rates observed in common ion channel mutations. Given the lack of experimental data, the model is not optimized and validated for variations in heart rate and altered extracellular ion concentrations. Loewe et al. [37] extended the Fabbri model to work at variable Ko concentrations by including a new formulation of I Kr , where g Kr is dependent on Ko. Their model showed a positive relationship between heart rate and Ko and they found that this was in disagreement with data from the general population, which showed a negative linear relationship. Our model, as the extended Fabbri model, shows a positive relationship between heart rate and Ko. This could mean that our model overestimates the CL at Ko = 3 mmol l −1 during NG as well as with LG and in combination with parasympathetic activity. This behaviour of our model could be inherited from the parent model. Indeed, existing human atrial models show divergent and sometimes unphysiological responses to changing extracellular Ko [38] and the corresponding electrophysiological conditions remain to be explored. Finally, the Fabbri model offers selective activation of either sympathetic or parasympathetic activity. In reality, both branches of the autonomic system are continuously active, and they complement each other to modulate the pacemaker rate. The responses of their combined activities and how these relate to their individual effects, observed in this study, remain to be investigated.

Conclusion
Computational modelling has the potential to clarify electropathophysiological effects of hypoglycaemia. We have shown that sinus arrest, a case of extreme bradycardia, is a plausible mechanism for SCD, possibly more probable than tachycardia during the night. Investigation of confounding factors that lead to the development of abnormal cardiac responses, as well as modelling of different cardiac cells at different levels may help to identify individuals with diabetes who are at an increased risk of sudden death due to sinus arrest during hypoglycaemia. This has important clinical implications because equivalent studies are not feasible or are challenging to perform in a clinical setting.
Data accessibility. Numerical data, source code of the model and the code to reproduce the figures are available in a Github repository: https://github.com/AlanBernjak/Bernjak_etal_Human_SAN_model_ hypoglycaemia.