Delay-induced uncertainty in the glucose-insulin system: Pathogenicity for obesity and type-2 diabetes mellitus

We have recently shown that physiological delay can induce a novel form of sustained temporal chaos we call delay-induced uncertainty (DIU) (Karamched et al. (Chaos, 2021, 31, 023142)). This paper assesses the impact of DIU on the ability of the glucose-insulin system to maintain homeostasis when responding to the ingestion of meals. We address two questions. First, what is the nature of the DIU phenotype? That is, what physiological macrostates (as encoded by physiological parameters) allow for DIU onset? Second, how does DIU impact health? We find that the DIU phenotype is abundant in the space of intrinsic parameters for the Ultradian glucose-insulin model—a model that has been successfully used to predict glucose-insulin dynamics in humans. Configurations of intrinsic parameters that correspond to high characteristic glucose levels facilitate DIU onset. We argue that DIU is pathogenic for obesity and type-2 diabetes mellitus by linking the statistical profile of DIU to the glucostatic theory of hunger.

1 Introduction

Clinical and laboratory practice throughout biomedicine and biochemistry proceeds from the assumption that the dynamics of measured quantities are predictable. For instance, a clinician administers medication to a patient based on the supposition that the medical intervention will not induce an unexpectedly erratic response. The presence of sustained temporal chaos would fundamentally undermine the assumption of predictability. Such chaos has been observed in certain classical physiological models (Abarbanel et al. (1993); Li and Yorke (2004); Mackey and Glass (1977); Glass et al. (1988); Glass and Malta (1990)).

We recently proposed a novel route through which physiological delay can induce sustained temporal chaos for concrete dynamical systems of interest in biomedicine (Karamched et al. (2021)). We termed the resulting chaos delay-induced uncertainty (DIU). We argued that DIU is relevant for glycemic management in the intensive care unit by exhibiting it for the Ultradian model, an archetypal model of glucose-insulin dynamics (Sturis et al. (1991); Drozdov and Khanina (1995)). Tools from the general theory of nonuniformly hyperbolic dynamical systems and the theory of rank-one maps yielded a precise characterization of the dynamical and statistical profiles of DIU. Clinicians may find DIU difficult to interpret because these profiles can be subtle.

DIU is potentially relevant for any physiological system wherein delayed regulatory feedback controls try to maintain healthy homeostasis. Examples include pulmonary and respiratory dynamics (Mackey and Glass (1977); Sottile et al. (2018)), cardiac dynamics (Christini and Glass (2002)), female endocrine dynamics (Graham et al. (2020); Urteaga et al. (2019)), and neurological dynamics (Stroh et al. (2020); Claassen et al. (2013); Hodgkin and Huxley (1952)). Indeed, the use of mathematical physiology within medicine has broad potential (Albers et al. (2018a); Zenker et al. (2007)).

This paper is a first attempt to assess the impact of DIU on the ability of the glucose-insulin system to maintain homeostasis when responding to the ingestion of meals. We work with the Ultradian model as we did before (Karamched et al. (2021)), but in a different regime. Our previous work focused on glycemic management in the intensive care unit (ICU) and therefore considered the regime wherein the intrinsic (unforced) system admits a glycemic oscillation (limit cycle). Here, we work in the regime wherein the intrinsic system admits a stable stationary state. In this regime, meals (glucose kicks) move trajectories away from the stationary point. After each kick, the glucose-insulin control system tries to efficiently return to the fixed point. We are therefore interested in how DIU impacts return to equilibrium.

In the context of return to equilibrium, the recipe for DIU has three ingredients. First, delay renders the unforced system excitable by weakening the stability of the stationary point. Second, shear is present near this stationary point. One can think of shear as velocity gradients. Third, external forcing (glucose kicking) interacts with shear during the relaxation phase between kicks. This interaction stretches and folds the phase space, creating hyperbolicity in the dynamics and producing sustained temporal chaos.

Here, we show that the physiological architecture of the glucose-insulin system possesses all three ingredients in the DIU recipe. We offer substantial evidence for the following two conjectures.

1) The DIU phenotype is abundant in the space of intrinsic parameters. In other words, a variety of physiological macrostates (as encoded by intrinsic parameters) lead to DIU emergence.

2) DIU is pathogenic for obesity and type-2 diabetes mellitus (T2DM).

This paper is a call to action—a first step toward verifying these conjectures.

Given the importance of elucidating obesity pathogenesis (Schwartz et al. (2017)), the DIU pathogenicity conjecture is the primary contribution of this work. The two-part argument supporting it links the statistical distribution of glucose that DIU induces to the glucostatic theory (Chaput and Tremblay (2009); Mayer (1955)). First, when DIU is present, glucose level dips below its mean more frequently. Second, glucostatic theory asserts that such dips induce hunger. See Figure 1 for an illustration of this two-part argument. This conjectured form of obesity pathogenesis acts on long timescales (months and years). As we will show, DIU becomes more probable as intrinsic parameters move into regions of parameter space that correspond to elevated characteristic glucose levels. Development of early-stage obesity and T2DM would therefore act as a feedback mechanism by promoting DIU, leading to disease progression.

FIGURE 1

FIGURE 1. How DIU may be pathogenic for obesity and T2DM. When DIU is present, glucose level dips below its mean more frequently (bottom row). The glucostatic hypothesis asserts that such dips induce hunger. DIU would therefore impute elevated hunger frequency, leading to obesity on long timescales.

We work with the Ultradian model for two primary reasons, validity and flexibility. The model includes two major negative feedback loops describing effects of insulin on glucose use and glucose production. Both loops include glucose-based stimulation of insulin secretion. External forcing can include both meal ingestion and glucose infusion. The Ultradian model can be tuned so that the unforced system admits a limit cycle, as in (Karamched et al. (2021)), or a stationary state. Importantly, it has been used to accurately predict glucose dynamics in humans (Albers et al. (2017)).

2 The Ultradian model

In this section we describe the Ultradian glucose-insulin model (Sturis et al. (1991); Drozdov and Khanina (1995); Keener and Sneyd (1998)), the external forcing drive that we use for simulations, and intrinsic system parameters that we hypothesize can facilitate DIU onset.

The Ultradian model is a compartment model with three state variables: plasma glucose (G), plasma insulin (I_p), and interstitial insulin (I_i). See Figure 2 for the model schematic. These three state variables are coupled to a three-stage linear delay filter, producing a six-dimensional phase space. The model includes two major negative feedback loops describing effects of insulin on glucose use and glucose production. Both loops include glucose-based stimulation of insulin secretion. The Ultradian model includes physiologic delay, but the system is finite-dimensional because the delay assumes the form of a three-stage linear filter.

FIGURE 2

FIGURE 2. Schematic for the Ultradian model of glucose-insulin dynamics. Note the important delayed regulatory feedback between I_p and G.

The full model is given by

$\frac{\mathrm{d}{I}_p}{\mathrm{d}t}=f_1\left(G\right)-E\left(\frac{I_p}{V_p}-\frac{I_i}{V_i}\right)-\frac{I_p}{t_p}(\mathrm{1a})$

$\frac{\mathrm{d}{I}_i}{\mathrm{d}t}=E\left(\frac{I_p}{V_p}-\frac{I_i}{V_i}\right)-\frac{I_i}{t_i}(\mathrm{1b})$

$\frac{\mathrm{d}G}{\mathrm{d}t}=f_4\left(h_3\right)+I_G\left(t\right)-f_2\left(G\right)-f_3\left(I_i\right)G(\mathrm{1c})$

$\frac{\mathrm{d}{h}_1}{\mathrm{d}t}=\frac{1}{t_d}\left(I_p-h_1\right)(\mathrm{1d})$

$\frac{\mathrm{d}{h}_2}{\mathrm{d}t}=\frac{1}{t_d}\left(h_1-h_2\right)(1e)$

$\frac{\mathrm{d}{h}_3}{\mathrm{d}t}=\frac{1}{t_d}\left(h_2-h_3\right),(\mathrm{1f})$

where f₁(G) represents the rate of insulin production, f₂(G) represents insulin-independent glucose use, f₃(I_i)G represents insulin-dependent glucose use, and f₄ (h₃) represents delayed insulin-dependent glucose use. The functional forms of f₁, f₂, f₃, and f₄ are given by

$f_1\left(G\right)=\frac{R_m}{1+\exp \left(\frac{-G}{V_g{C}_1}+a_1\right)}(\mathrm{2a})$

$f_2\left(G\right)=U_b\left(1-\exp \left(\frac{-G}{C_2{V}_g}\right)\right)(\mathrm{2b})$

$f_3\left(I_i\right)=\frac{1}{C_3{V}_g}\left(U_0+\frac{U_m-U_0}{1+{\left(\kappa I_i\right)}^{-\beta }}\right)(\mathrm{2c})$

$f_4\left(h_3\right)=\frac{R_g}{1+\exp \left(\alpha \left(\frac{h_3}{C_5{V}_p}-1\right)\right)},(\mathrm{2d})$

with

$\kappa =\frac{1}{C_4}\left(\frac{1}{V_i}-\frac{1}{E{t}_i}\right).(3)$

Table 1 summarizes the meaning of each model parameter and provides the set of nominal parameter values.

TABLE 1

TABLE 1. Full list of intrinsic parameters for the Ultradian glucose-insulin model (Albers et al. (2017)). Note that IIGU and IDGU denote insulin-independent glucose utilization and insulin-dependent glucose utilization, respectively.

2.1 Pulsatile glucose forcing drives

The term I_G(t) in Eq.1c represents the external nutritional drive. We call system (1) without this term the intrinsic system or unforced system. In this paper, we consider an idealized nutritional drive I_G(t) that consists of pulsatile kicks. This drive models meals that are eaten and digested instantaneously. That is, we assume that the nutritional content of each meal immediately affects the glucose state variable in the Ultradian system. The idealized nutritional drive is given by

$I_G\left(t\right)=\sum _{n=1}^{\infty }{A}_n\delta \left(t-T_n\right),(4)$

where δ(t) is the Dirac delta distribution (unit impulse), T_n is the time of meal n, and A_n is the amount of carbohydrate in meal n. Importantly, this pulsatile drive does not overwhelm the intrinsic dynamics. On the contrary, it can interact subtly with intrinsic shear to produce DIU, as we will see.

The form of I_G(t) in Eq. 4 induces the following dynamics. Between two consecutive kicks (T_n−1 < t < T_n), Ultradian dynamics evolve according to system (1) with I_G(t) = 0. At time T_n of meal n, the glucose state variable, G, undergoes the instantaneous change G↦G + A_n. That is, at time T_n we pause the flow generated by the intrinsic system and apply the diffeomorphism

$\left(I_p,I_i,G,h_1,h_2,h_3\right)\mapsto \left(I_p,I_i,G+A_n,h_1,h_2,h_3\right)(5)$

to the phase space. We call this diffeomorphism followed by flow of the intrinsic system cycle the kick-relaxation cycle.

In reality, meals produce glucose perturbations that are temporally localized but not instantaneous. Nevertheless, we have strong evidence that the emergence of DIU (or the absence of such emergence) is sensitive to neither the exact timing of the pulses nor to their shape. In previous work (Karamched et al. (2021)), we examined the emergence of DIU for the Ultradian model when the delay parameter t_d is tuned so that the intrinsic system admits a limit cycle (sustained oscillatory dynamics). There, we showed that DIU can emerge when each inter-meal time is drawn from an exponential distribution (Poissonian inter-meal timing) and when the drive 4) is replaced with square pulses of duration 30 min that arrive at 8 a.m., noon, and 6 PM. Here, we elect to work with drive 4) and consider only periodic pulsing (T_n = nT, where $T\in {\mathbb{R}}_{>0}$ is the inter-kick time) with constant kick amplitude (A_n = A for all $n\in {\mathbb{Z}}_{⩾0}$) in order to focus on how intrinsic parameters impact DIU emergence. Our previous work indicates that our new results for periodic pulsatile forcing will continue to hold for more complex forcing drives.

2.2 Key intrinsic parameters for DIU emergence

We hypothesize that intrinsic (unrelated to the forcing drive) parameters directly linked to G, the glucose state variable, play a key role in DIU onset. This hypothesis is partially inspired by recent work that established a positive correlation between mean glucose levels and glucose variance (Albers et al. (2018b)). Our numerical experiments examine the impact of the following parameters on DIU emergence.

• R_g - the uninhibited hepatic glucose production rate

• U_b - the maximal insulin-independent glucose usage rate

• U₀ - the basal insulin-dependent glucose usage rate

• α - the inhibition of hepatic glucose production

• a₁ - the basal glucose-based insulin inhibition

• C₁ - the sensitivity of insulin production to glucose

Importantly, each of these intrinsic parameters has a concrete physiological interpretation.

3 Methods

The maximal Lyapunov exponent as a diagnostic tool. We use the maximal Lyapunov exponent, Λ_max, as a DIU diagnostic: Λ_max > 0 indicates DIU whereas Λ_max < 0 indicates its absence. Computing Λ_max requires solving system (1). We do this in the following way. During the relaxation intervals (T_n−1, T_n) between kicks, we integrate the unforced differential equations using the MATLAB ode23s solver. At kick times T_n, we pause the differential equation solver and apply the diffeomorphism of phase space induced by the kick (see Eq. 5).

We compute the maximal Lyapunov exponent in the following way. We track two solutions to system (1), initially separated by d₀ = 10^–8. One of these solutions can be thought of as a base solution and the other as a perturbation. After the first kick-relaxation cycle, we compute the separation d₁ between the solutions and store the quantity log (d₁/d₀) in a vector. We then renormalize by rescaling the secondary orbit so that the distance between the solutions resets to d₀. We proceed in this manner for 10⁵ kick-relaxation cycles. This produces a vector containing 10⁵ values of log (d₁/d₀). Averaging over the vector produces Λ_max. The maximal Lyapunov exponent consequently quantifies the amount of expansion per kick-relaxation cycle.

4 Results

We have designed our numerical experiments to support two primary conjectures. These conjectures are the animating force behind this paper.

1) The DIU phenotype is abundant in the space of intrinsic parameters. In other words, a variety of physiological macrostates (as encoded by intrinsic parameters) lead to DIU emergence.

2) DIU is pathogenic for obesity and T2DM.

4.1 Numerical experiments: Design, rationale, and expectations

Tuning of intrinsic parameters. To support the conjecture that the DIU phenotype is abundant in the space of intrinsic parameters, we begin by setting the intrinsic parameters in the unforced Ultradian model to the nominal values listed in Table 1. Crucially, the delay timescale t_d acts as a bifurcation parameter for the intrinsic system. There exists a value $t_d^{*}$ at which the intrinsic system undergoes a supercritical Hopf bifurcation. The intrinsic system admits a stable stationary point (I_p,eq, I_i,eq, G_eq, h_1,eq, h_2,eq, h_3,eq) for $t_d<t_d^{*}$ (homeostasis) that gives birth to a stable limit cycle (glycemic oscillation) for $t_d>t_d^{*}$. For our numerical experiments, we set t_d to the nominal value 10.5 min, a value strictly less than $t_d^{*}$, thereby placing the intrinsic system in the stable stationary point regime. This is the appropriate regime for our current study because we are interested in how the dynamical variables relax to homeostatic levels between glucose kicks.

Using the nominal values of the intrinsic parameters as a starting point, we look for DIU along six one-dimensional slices of parameter space. We select a parameter from the list given in Section 2.2 and then vary this parameter while holding all other intrinsic parameters fixed.

Testing for DIU onset. Having set the intrinsic parameters, we test for DIU onset by tuning the external pulsatile forcing drive (4). For the sake of simplicity, we select a kick amplitude A and set A_n = A for all $n\in {\mathbb{Z}}_{⩾0}$. We work with periodic kicks, so we set T_n = nT, where T is the time between consecutive kicks. The forcing drive 4) for the experiments is therefore given by

$I_G\left(t\right)=A\sum _{n=1}^{\infty }\delta \left(t-nT\right).(6)$

To test for DIU onset, we compute the maximal Lyapunov exponent Λ_max as a function of T.

Expectations. DIU may or may not emerge as T increases, depending on the dynamics of the intrinsic flow near the stationary point. If contraction to the stationary state is strong and shear near the stationary state is weak, DIU will not emerge. The maximal Lyapunov exponent Λ_max will indicate this by remaining negative as T increases. In fact, Λ_max will decrease as T increases because the phase space has more time to contract between kicks as T increases.

On the other hand, if contraction to the stationary state is weak and shear near the stationary state is strong, then DIU can emerge as T increases. This can happen because when T is large, shear has a long time to act between kicks. Shear causes the phase space to stretch and fold, thereby producing DIU. In our experiments, a transition from Λ_max < 0 to Λ_max > 0 as T increases indicates that DIU has emerged.

4.2 The DIU phenotype is abundant in the space of intrinsic parameters

Figure 3 illustrates how Λ_max varies with T as we individually tune each of the six parameters identified in Section 2.2. Each panel corresponds to tuning a single parameter while holding all other intrinsic parameters fixed at the nominal values. Importantly, DIU emerges in every one of the six panels when we tune the selected parameter so as to increase characteristic glucose levels in the intrinsic dynamics.

FIGURE 3

FIGURE 3. The DIU phenotype is abundant in the space of intrinsic parameters. Plots show the maximal Lyapunov exponent Λ_max as a function of inter-kick time T for the time-T map induced by the Ultradian system (1) with T-periodic pulsatile forcing (6). DIU is present when Λ_max > 0 and absent when Λ_max < 0. As T increases, DIU emerges when intrinsic parameters are tuned so as to increase characteristic glucose levels. Intrinsic parameters are set to the nominal values in Table 1 except for the single intrinsic parameter that is tuned in each panel: (A) R_g; (B) U_b; (C) U₀; (D) α; (E) a₁; (F) C₁. Kick amplitude: A = 10 mg/dl.

Figure 4 confirms the expected link between strength of contraction to the stationary point, shear near the stationary point, and DIU emergence. For Figure 4, we replace the periodic pulsatile forcing used to generate Figure 3 with a forcing signal that consists of three kicks (meals). After the final kick, the glucose variable converges to the equilibrium level G_eq as t → ∞. The panels in Figure 4 indicate that our experiments have captured two behaviors. Either we see rapid convergence to G_eq (as in Figure 4D (top)), or we see slow convergence to G_eq by way of a damped oscillation (as in Figure 4D (bottom)). Notice that in each panel of Figure 4, we tune the same parameter that we tune in the corresponding panel of Figure 3, while holding all other intrinsic parameters fixed at the nominal values.

FIGURE 4

FIGURE 4. Glucose trajectories generated by the Ultradian system (1) with forcing that consists of three glucose kicks spaced 100 min apart. After the final kick, the glucose level converges to the equilibrium value G_eq. Convergence is either rapid (top of each panel) or via a slow damped oscillation (bottom of each panel). Intrinsic parameters are set to the nominal values in Table 1 except for the single intrinsic parameter that is tuned in each panel: (A) R_g; (B) U_b; (C) U₀; (D) α; (E) a₁; (F) C₁. Kick amplitude: A = 10 mg/dl.

Comparing Figures 3, 4 shows that without exception, the geometry of the glucose trajectory predicts whether or not DIU will emerge. If we observe rapid convergence to G_eq, as in Figure 4D (top) for instance, then DIU does not emerge. If, however, we observe slow convergence to G_eq by way of a damped oscillation, as in Figure 4D (bottom) for instance, then DIU emerges.

Figure 5 illustrates the DIU dynamical profile and acts as a companion to Figure 4. Each glucose trajectory in Figure 5 results from forcing with periodic pulsatile kicks 6) and corresponds to a companion glucose trajectory in Figure 4 (produced by applying only three kicks). When contraction toward the equilibrium glucose level G_eq is strong (Figures 4A–F (top)), driving with periodic pulsatile kicks produces rhythmic behavior (Figures 5A–F (top)). When periodic pulsatile kicks produce DIU, glucose trajectories exhibit sustained temporal chaos (Figures 5A–F (bottom)).

FIGURE 5

FIGURE 5. Sustained temporal chaos associated with DIU. Plots show glucose trajectories produced by the Ultradian system (1) with T-periodic pulsatile forcing (6). Each trajectory in Figure 5 corresponds to a companion trajectory in Figure 4. When DIU is absent, T-periodic pulsatile forcing results in a rhythmic glucose signal (top of each panel). When DIU is present, we observe sustained temporal chaos (bottom of each panel). Intrinsic parameters are set to the nominal values in Table 1 except for the single intrinsic parameter that is tuned in each panel: (A) R_g; (B) U_b; (C) U₀; (D) α; (E) a₁; (F) C₁. Forcing parameters: A = 10 mg/dl, T = 100 min.

4.3 DIU is pathogenic for obesity and T2DM

We have established that the DIU phenotype is abundant in the space of intrinsic parameters for the Ultradian model. But why does this matter? Delayed regulatory feedback pathways are common in mathematical physiology. Since DIU emerges in a natural way for the Ultradian model, it may appear in a variety of physiological models. When present, DIU can profoundly impact medical practice because medicine proceeds from the assumption that the outcome of an intervention can be predicted when the state of the patient at the time of intervention is known. Sustained temporal chaos undercuts this assumption. See (Karamched et al. (2021)) for an assessment of the impact of DIU on glycemic management in the intensive care unit.

Here, we conjecture that DIU is pathogenic for obesity and T2DM. This conjecture is based on how the statistical signature of DIU links to the glucostatic theory. The glucostatic theory asserts that drops in blood glucose levels induce hunger and therefore energy intake (Chaput and Tremblay (2009); Mayer (1955)). If such drops are frequent in time and sizable in magnitude, excess energy intake could result.

Figure 6 shows that DIU induces frequent, sizable drops in blood glucose levels! We start with all intrinsic parameters set at the nominal values and we then tune R_g, the uninhibited hepatic glucose production rate. Figure 6 shows glucose distributions for the time-T map induced by Ultradian dynamics 1) with T-periodic pulsatile forcing (6). That is, each histogram gives the distribution of

$\left\{G\left(nT\right):n\in {\mathbb{Z}}_{⩾0}\right\}(7)$

for a different value of R_g. We set T = 100 min. When R_g = 120 mg/min (Figure 6A), a value for which DIU is absent, the glucose distribution is essentially a Dirac measure concentrated at the mean. (Blue indicates the mean of the glucose distribution and orange indicates the distribution itself.) However, when R_g = 180 mg/min (Figure 6B), a value for which DIU is present, the glucose distribution is approximately Gaussian. This is as it should be—the mathematical theory behind DIU predicts Gaussian statistics when DIU is present. Notice that the variance of the approximately Gaussian distribution is large. This means that the glucose level frequently drops well below its mean. In light of glucostatic theory, this observation directly supports the conjecture that DIU is pathogenic for obesity and T2DM.

FIGURE 6

FIGURE 6. Support for the conjecture that DIU is pathogenic for obesity and T2DM. Distributions of the glucose variable (7) for the time-T map induced by the Ultradian system (1) with T-periodic pulsatile forcing (6). Blue bar indicates mean. (A) When R_g = 120 mg/min, DIU is absent and the glucose distribution concentrates at the mean. (B) When R_g = 180 mg/min, DIU is present. Consistent with the mathematical structure of the DIU profile, the glucose distribution is approximately Gaussian. All of the other intrinsic parameters are set to the nominal values in Table 1. Forcing parameters: A = 10 mg/dl, T = 100 min.

4.4 DIU emerges for generic pulsatile meal drives

Our results do not depend on the precise form of the pulsatile forcing that appears in (6). The forcing need not be periodic, and it need not consist of δ-pulses. DIU should emerge for a generic pulsatile forcing drive as long as the forcing interacts with intrinsic shear¹.

To support this claim, we have varied the intrinsic parameter R_g to test for DIU emergence after replacing 6) with

$I_G\left(t\right)=A\sum _{n=1}^{\infty }\mathrm{\Theta }\left(t-m_n\right)e^{-\upsilon \left(t-m_n\right)},(8)$

where A > 0 denotes meal amplitude, Θ(t) is the Heaviside function, υ > 0 is a constant, and m_n denotes the time of meal n. For this set of experiments, meals are consumed daily at 8 a.m., noon, and 6 p.m. Figure 7A shows Λ_max as a function of meal amplitude A for three values of R_g. For two of the three values of R_g, the top Lyapunov exponent becomes positive as A increases, indicating DIU onset. For Figure 7bc, we replace (8) with a single day of meals (three meals) in order to show that the nature of return to equilibrium correlates with DIU onset. Rapid return to equilibrium correlates with the absence of DIU (Figure 7B), while slow, oscillatory return to equilibrium correlates with the presence of DIU (Figure 7C).

FIGURE 7

FIGURE 7. DIU phenotype for a realistic nutritional driver. (A) We have replaced (6) with the exponential-type drive in (8). Meals are consumed daily at 8 a.m., noon, and 6 p.m. Plot shows Λ_max as a function of meal amplitude A for three values of R_g (B,C) We replace (8) with a single day of meals (three meals). Rapid return to equilibrium correlates with the absence of DIU, while slow, oscillatory return to equilibrium correlates with the presence of DIU. Here, A = 50 mg/(dL ⋅ min) and υ = 1/120 min⁻¹.

5 Discussion

We have found that DIU is abundant in the space of parameters for the Ultradian glucose-insulin model. Such DIU could result in obesity and T2DM if induced low-glucose excursions produce excess hunger frequently enough, but much work remains to verify the conjecture that DIU is pathogenic for obesity and T2DM. Crucially, DIU and the theory behind it must be anchored to data. Methods for DIU detection directly from data should be developed for the clinical and self-care settings. The impact of DIU on the techniques by which models are fit to data should be assessed.

We have assumed in this paper that the intrinsic parameters in the Ultradian model do not vary over time. On long timescales, however, DIU may affect physiological state. At the modeling level, this would correspond to DIU causing intrinsic model parameters to drift (perhaps slowly) over time. Such drift might enhance the pathogenicity of DIU through a feedback mechanism: When DIU is present, intrinsic parameters may slowly drift into a region of parameter space that is even more favorable for DIU. A mathematical investigation of this phenomenon would involve developing a theory of DIU for nonstationary dynamical systems.

We have shown here that the DIU phenotype is abundant in the space of intrinsic parameters for the Ultradian model. An important next step will be to precisely characterize the DIU phenotype in terms of physiological architecture. Such a characterization may reveal the most essential physiological mechanisms that lead to DIU onset. Mathematically speaking, we must quantify shear near stationary states of flows. Shear near limit cycles has received considerable attention (Ott and Stenlund (2010); Wang and Young (2003)). Shear near stationary states, though, has only been quantified in dimension two (Ott (2008)).

The rigorous mathematical theory behind DIU is known as the theory of rank-one maps. This theory has been developed for finite-dimensional dynamical systems (Wang and Young (2001, 2008; 2013)). The Utradian model is finite-dimensional as a dynamical system because the delay in the Ultradian model takes the form of a three-stage linear filter. The theory of rank-one maps therefore characterizes the sustained temporal chaos that we see in the Ultradian model. However, models that include explicit delays—systems of nonlinear delay differential equations—permeate mathematical physiology. Models that include explicit delays are infinite-dimensional when viewed as dynamical systems. Important infinite-dimensional analogs of the Ultradian model have been studied (Li et al. (2006); Li and Kuang (2007)). The theory of rank-one maps must be extended to infinite-dimensional dynamical systems in order to analyze delay differential equations in the DIU context. See (Lu et al. (2013)) for an approach that combines the existing theory of rank-one maps with invariant manifold techniques.

When assessing the impact of DIU on a given physiological system, one should ask the following questions. Are we interested in precisely predicting the temporal evolution of individual orbits, or do we care more about the statistics of the system? What are the relevant timescales? For the glucose-insulin system, we have now studied two contrasting settings. In the ICU context, we showed that DIU can disrupt single-orbit prediction on short timescales (Karamched et al. (2021)). In the present paper, we have argued that over long timescales, DIU-induced glucose statistics may be pathogenic for obesity and T2DM.

Data availability statement

The original contributions presented in the study are included in the article/supplementary materials, further inquiries can be directed to the corresponding author/s. The code used to generate all figures is available here: https://github.com/Bargo727/DIU.

Author contributions

BK performed the simulations. BK and WO wrote the first draft of the manuscript. All authors contributed to subsequent revisions.

Funding

This work has been supported by the National Science Foundation under grant DMS 1816315 (WO) and by the National Institutes of Health under grants LM006910 (GH), R01DK052431 (RL), P30DK26687 (RL), and LM012734 (DA).

Acknowledgments

BK would like to thank his wife, Hajra Habib, for her unending support.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Footnotes

¹Here we compute Λ_max by averaging over time intervals of length 12 h.

References

Abarbanel H. D., Brown R., Sidorowich J. J., Tsimring L. S. (1993). The analysis of observed chaotic data in physical systems. Rev. Mod. Phys. 65, 1331–1392. doi:10.1103/RevModPhys.65.1331

CrossRef Full Text | Google Scholar

Albers D. J., Elhadad N., Claassen J., Perotte R., Goldstein A., Hripcsak G. (2018b). Estimating summary statistics for electronic health record laboratory data for use in high-throughput phenotyping algorithms. J. Biomed. Inf. 78, 87–101. doi:10.1016/j.jbi.2018.01.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Albers D. J., Levine M., Gluckman B., Ginsberg H., Hripcsak G., Mamykina L. (2017). Personalized glucose forecasting for type 2 diabetes using data assimilation. PLoS Comput. Biol. 13, e1005232. doi:10.1371/journal.pcbi.1005232

PubMed Abstract | CrossRef Full Text | Google Scholar

Albers D. J., Levine M. E., Stuart A., Mamykina L., Gluckman B., Hripcsak G. (2018a). Mechanistic machine learning: How data assimilation leverages physiologic knowledge using bayesian inference to forecast the future, infer the present, and phenotype. J. Am. Med. Inf. Assoc. 25, 1392–1401. doi:10.1093/jamia/ocy106

PubMed Abstract | CrossRef Full Text | Google Scholar

Chaput J., Tremblay A. (2009). The glucostatic theory of appetite control and the risk of obesity and diabetes. Int. J. Obes. 33, 46–53. doi:10.1038/ijo.2008.221

PubMed Abstract | CrossRef Full Text | Google Scholar

Christini D. J., Glass L. (2002). Introduction: Mapping and control of complex cardiac arrhythmias. Chaos 12, 732–739. doi:10.1063/1.1504061

PubMed Abstract | CrossRef Full Text | Google Scholar

Claassen J., Perotte A., Albers D., Kleinberg S., Schmidt J. M., Tu B., et al. (2013). Nonconvulsive seizures after subarachnoid hemorrhage: Multimodal detection and outcomes. Ann. Neurol. 74, 53–64. doi:10.1002/ana.23859

PubMed Abstract | CrossRef Full Text | Google Scholar

Drozdov A., Khanina H. (1995). A model for ultradian oscillations of insulin and glucose. Math. Comput. Model. 22, 23–38. doi:10.1016/0895-7177(95)00108-E

CrossRef Full Text | Google Scholar

Glass L., Beuter A., Larocque D. (1988). Time delays, oscillations, and chaos in physiological control systems. Math. Biosci. 90, 111–125. doi:10.1016/0025-5564(88)90060-0

CrossRef Full Text | Google Scholar

Glass L., Malta C. P. (1990). Chaos in multi-looped negative feedback systems. J. Theor. Biol. 145, 217–223. doi:10.1016/s0022-5193(05)80127-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Graham E., Elhadad N., Albers D. (2020). Reduced model for female endocrine dynamics: Validation and functional variations. arXiv Prepr. arXiv:2006.05034.

Google Scholar

Hodgkin A. L., Huxley A. F. (1952). A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 117, 500–544. doi:10.1113/jphysiol.1952.sp004764

PubMed Abstract | CrossRef Full Text | Google Scholar

Karamched B., Hripcsak G., Albers D., Ott W. (2021). Delay-induced uncertainty for a paradigmatic glucose–insulin model. Chaos 31, 023142. doi:10.1063/5.0027682

PubMed Abstract | CrossRef Full Text | Google Scholar

Keener J., Sneyd J. (1998). Mathematical physiology, interdisciplinary applied mathematics, vol. 8. New York: Springer-Verlag, xx+766.

Google Scholar

Li J., Kuang Y. (2007). Analysis of a model of the glucose-insulin regulatory system with two delays. SIAM J. Appl. Math. 67, 757–776. doi:10.1137/050634001

CrossRef Full Text | Google Scholar

Li J., Kuang Y., Mason C. (2006). Modeling the glucose-insulin regulatory system and ultradian insulin secretory oscillations with two explicit time delays. J. Theor. Biol. 242, 722–735. doi:10.1016/j.jtbi.2006.04.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Li T. Y., Yorke J. A. (2004). “Period three implies chaos,” in The theory of chaotic attractors (Berlin, Germany: Springer), 77–84.

CrossRef Full Text | Google Scholar

Lu K., Wang Q., Young L. S. (2013). Strange attractors for periodically forced parabolic equations. Mem. Am. Math. Soc. 224, 1. doi:10.1090/S0065-9266-2012-00669-1

CrossRef Full Text | Google Scholar

Mackey M. C., Glass L. (1977). Oscillation and chaos in physiological control systems. Science 197, 287–289. doi:10.1126/science.267326

PubMed Abstract | CrossRef Full Text | Google Scholar

Mayer J. (1955). Regulation of energy intake and the body weight: The glucostatic theory and the lipostatic hypothesis. Ann. N. Y. Acad. Sci. 63, 15–43. doi:10.1111/j.1749-6632.1955.tb36543.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Ott W., Stenlund M. (2010). From limit cycles to strange attractors. Commun. Math. Phys. 296, 215–249. doi:10.1007/s00220-010-0994-y

CrossRef Full Text | Google Scholar

Ott W. (2008). Strange attractors in periodically-kicked degenerate Hopf bifurcations. Commun. Math. Phys. 281, 775–791. doi:10.1007/s00220-008-0499-0

CrossRef Full Text | Google Scholar

Schwartz M. W., Seeley R. J., Zeltser L. M., Drewnowski A., Ravussin E., Redman L. M., et al. (2017). Obesity pathogenesis: An Endocrine Society scientific statement. Endocr. Rev. 38, 267–296. doi:10.1210/er.2017-00111

PubMed Abstract | CrossRef Full Text | Google Scholar

Sottile P. D., Albers D., Higgins C., Mckeehan J., Moss M. M. (2018). The association between ventilator dyssynchrony, delivered tidal volume, and sedation using a novel automated ventilator dyssynchrony detection algorithm. Crit. Care Med. 46, e151–e157. doi:10.1097/CCM.0000000000002849

PubMed Abstract | CrossRef Full Text | Google Scholar

Stroh J., Bennett T., Kheyfets V., Albers D. (2020). Estimating intracranial pressure via low-dimensional models: Toward a practical tool for clinical decision support at multi-hour timescales. bioRxiv.

Google Scholar

Sturis J., Polonsky K., Mosekilde E., Van Cauter E. (1991). Computer model for mechanisms underlying ultradian oscillations of insulin and glucose. Am. J. Physiol. 260, E801–E809. doi:10.1152/ajpendo.1991.260.5.E801

PubMed Abstract | CrossRef Full Text | Google Scholar

Urteaga I., Bertin T., Hardy T. M., Albers D. J., Elhadad N. (2019). Multi-task Gaussian processes and dilated convolutional networks for reconstruction of reproductive hormonal dynamics. arXiv Prepr. arXiv:1908.10226.

Google Scholar

Wang Q., Young L. S. (2013). Dynamical profile of a class of rank-one attractors. Ergod. Th. Dynam. Sys. 33, 1221–1264. doi:10.1017/S014338571200020X

CrossRef Full Text | Google Scholar

Wang Q., Young L. S. (2003). Strange attractors in periodically-kicked limit cycles and Hopf bifurcations. Commun. Math. Phys. 240, 509–529. doi:10.1007/s00220-003-0902-9

CrossRef Full Text | Google Scholar

Wang Q., Young L. S. (2001). Strange attractors with one direction of instability. Commun. Math. Phys. 218, 1–97. doi:10.1007/s002200100379

CrossRef Full Text | Google Scholar

Wang Q., Young L. S. (2008). Toward a theory of rank one attractors. Ann. Math. (2), 167349–167480. doi:10.4007/annals.2008.167.349

CrossRef Full Text | Google Scholar

Zenker S., Rubin J., Clermont G. (2007). From inverse problems in mathematical physiology to quantitative differential diagnoses. PLoS Comput. Biol. 3, e204. doi:10.1371/journal.pcbi.0030204

PubMed Abstract | CrossRef Full Text | Google Scholar