Modeling the dynamics of glucose, insulin, and free fatty acids with time delay: The impact of bariatric surgery on type 2 diabetes mellitus

The role of free fatty acids (FFA) on Type 2 diabetes mellitus (T2DM) progression has been studied extensively with prior studies suggesting that individuals with shared familial genetic predisposition to metabolic-related diseases may be vulnerable to dysfunctional plasma FFA regulation. A harmful cycle arises when FFA are not properly regulated by insulin contributing to the development of insulin resistance, a key indicator for T2DM, since prolonged insulin resistance may lead to hyperglycemia. We introduce a hypothesis-driven dynamical model and use it to evaluate the role of FFA on insulin resistance progression that is mathematically constructed within the context of individuals that have genetic predisposition to dysfunctional plasma FFA. The dynamics of the nonlinear interactions that involve glucose, insulin, and FFA are modeled by incorporating a fixed-time delay with the corresponding delay-differential equations being studied numerically. The results of computational studies, that is, extensive simulations, are compared to the known minimal ordinary differential equations model. Parameter estimation and model validation are carried out using clinical data of patients who underwent bariatric surgery. These estimates provide a quantitative measure that is used to evaluate the regulation of lipolysis by insulin action measured by insulin sensitivity, within a metabolically heterogeneous population (non-diabetic to diabetic). Results show that key metabolic factors improve after surgery, such as the effect of insulin inhibition of FFA on insulin and glucose regulation, results that do match prior clinical studies. These findings indicate that the reduction in weight or body mass due to surgery improve insulin action for the regulation of glucose, FFA, and insulin levels. This reinforces what we know, namely, that insulin action is essential for regulating FFA and glucose levels and is a robust effect that can be observed not only in the long-term, but also in the short-term; thereby preventing the manifestation of T2DM.

bring glucose levels back to the normal range. IR can also develop in the liver, known as hepatic insulin resistance, which is defined as impaired insulin signalling affecting the liver. Clinical studies show that insulin resistance in skeletal muscle and adipose tissue are common in individuals with T2DM [6][7][8][9], and thus IR is considered a strong predictor for the development of T2DM.
Several studies support the hypothesis that prolonged presence of free fatty acids (FFA) circulating in the bloodstream or in organs (e.g., skeletal muscle, pancreas, or liver) contributes to the development of T2DM by impacting insulin signalling, glucose transport, insulin secretory function, or other mechanisms leading to disruptions in glucose-insulin regulation [10][11][12]. FFA is used as fuel for the liver, resting skeletal muscle, renal cortex, and myocardium [10]. During starvation or exercise, the demand for FFA as a fuel source increases, where after an overnight fast lipid oxidation can account for over 70% of total body energy expenditure [10]. When FFA is needed as a source of energy, then triglycerides (TG) are broken down into FFA and released into the bloodstream via lipolysis. Hence, FFA can be approximated by TG. Since TG is associated with cardiovascular diseases [13], then high TG levels coexisting with normal or impaired fasting glucose further increase risk of metabolic-related diseases, such as T2DM. FFA is stored in adipocytes (adipose cells) from TG which are made up of three FFA molecules and one glycerol molecule. As of today, several hypotheses have been proposed characterizing the mechanistic role of FFA on the progression of insulin resistance (summarized in Table 1). It is widely accepted that elevated FFA levels promote insulin resistance in skeletal muscle and liver [10,14,15] and yet the mechanisms underlying the progression of T2DM based on obesity levels are not well understood. Here we formulate a mathematical framework based on what is currently known about the role of FFA on IR progression as a starting point, or baseline model, to be expanded upon in future work to test different hypotheses through computational studies.

Bariatric surgery as a treatment strategy
Obesity, characterized by excessive body fat, is considered a risk factor for diabetes since among those individuals diagnosed, 84.7% turned out to be overweight, with 56.9% obese [1]. In many cases, patients with severe obesity and diagnosed with T2DM can have improved metabolic health following bariatric surgery, which is a treatment strategy for patients with severe obesity. These procedures may include: gastric banding (such as adjustable and nonadjustable bands), gastric bypass (such as Rouxen-Y variations or any other procedure combined with gastric bypass), gastroplasty (such as vertical banded gastriplasty), biliopancreatic diversion or duodenal switch (such as various modifications), or a combination of these methods [24]. In a review of 134 studies [24], it was found that 76.8% patients who underwent bariatric surgery improved diabetes-related outcomes, and a total of 85.4% patients had either resolved or improved diabetes-related outcomes. Bariatric surgery patients with diabetes have demonstrated recovery of acute insulin response, also referred to as the "first phase" of insulin secretion [25], decreases of inflammatory indicators (C-reactive protein and interleukin 6) which are associated with cardiovascular and metabolic health [26], improvement in insulin sensitivity [27,28], significant changes in the response of enteroglucagon, defined as a peptide hormone that helps to facilitate the absorption of nutrients in the presence of glucose and fats, to glucose levels [29], significant reduction in ghrelin, also referred to as the "hunger hormone" which triggers appetite [30], and significant improvement in beta-cell function [31]. Hence, a reduction of body fat following bariatric surgery leads to improvements in the patient's metabolic health.
In this paper, we consider the case of a prototypic individual with genetic susceptibility to T2DM, that is, an individual who self-identified as having a family member with T2DM or other metabolic-associated diseases, and thus, could be susceptible defined by its inability to regulate FFA or by a disruption on the regulation of FFA. We refer to this process as the "harmful cycle hypothesis" and define the mechanism as follows. A dysfunctional regulation of FFA that leads to an increase demand for insulin, a demand that is stimulated by both glucose and FFA, which in turn, puts the individuals at higher risk of hyperglycemia when beta-cell compensation fails. In the long-term, insulin resistance is observed when insulin does not effectively regulate plasma glucose and FFA. We introduce a model, a hypothesisdriven model, built off of models previously used to study the joint dynamics of insulin, glucose, and FFA. These dynamics are studied quantitatively and qualitatively within a framework that includes model validation. Parameters are fitted to data obtained from a heterogeneous sample of patients ranging from non-diabetic and non-obese to diabetic and obese so that the model results to may provide insights into the physiological factors that govern glucose homeostasis. By fitting the model to clinical data of patients who underwent bariatric surgery, we explore to what extent does a reduction of body fat following surgery alter a patient's metabolic health by recovering insulin sensitivity, beta-cell function, and other related parameters. This paper is organized as follows: In section 2 the mathematical models are described including the classic minimal model (in section 2.1) of insulin, glucose, and FFA as described in the work of [32] and an explicit time-delay model of insulin, glucose, and FFA dynamics (in section 2.2) incorporating the harmful cycle hypothesis with corresponding analytical results, in section 4 the results of the model validation and parameter estimates are shown to compare both the minimal model and explicit time delay model, and in section 5 the conclusions and future work is described.

Mathematical modeling framework
The insulin sensitivity is defined through the clamp protocol. Although the clamp protocol is the gold standard for assessing insulin sensitivity, not only do the subjects suffer pain from the procedure, but also the test is labor-intensive and financially expensive. The intravenous glucose tolerance test (IVGTT) has been clinically considered the most accurate protocol next to the clamp protocol to determine insulin sensitivity and glucose effectiveness through the approach of mathematical modeling with curve-fitting [34][35][36][37][38]. The data provided through an IVGTT offer rich information and offer a more realistic picture of a subjects metabolic portrait of insulin sensitivity and glucose effectiveness by differentiating glucose production from the liver and the exogenous glucose administered during the study protocol [39]. Furthermore, the IVGTT protocol has been extended to gain a better understanding of the dynamical regulations underlying FFA, insulin and glucose [32,40], in which, as the original IVGTT, subjects fast overnight and then are given a bolus of glucose infusion intravenously (e.g., 0.33 g/kg of body weight or 0.3 g/kg body weight of a 50% solution), which is administered into the antecubital vein in approximately 2 minutes [39,41]. FFA, plasma glucose and insulin levels are sampled over the duration of the test. The short dynamics captured by the mathematical model begins with a rise in plasma glucose due to bolus infusion, which triggers pancreatic beta-cells to quickly secrete insulin into the bloodstream. Insulin mediates glucose removal, also referred to as insulin-dependent removal, which in turn, lowers plasma glucose to basal level and then the demand for insulin is inhibited, i.e. negative feedback. Meanwhile, FFA production is inhibited by insulin when glucose supply is high (see Figure 1). On the other hand, some prior studies suggest that insulin inhibition of FFA is weak in individuals genetically predisposed to metabolic-related diseases. Furthermore, it is known that FFA may enhance basal and glucose-stimulated insulin secretion among individuals including those with diabetes [10]. We hypothesize and validate that insulin ineffectively regulates FFA and higher FFA may reduce glucose transport, which leads to a harmful cycle promoting hyperglycemia and developing insulin resistance in the long-term.

Minimal model of glucose, insulin, and FFA
The minimal model was introduced in 1979 and 1980 by Bergman, Cobelli and their colleagues. It was the first model to define two significant indices, the glucose effectiveness index and the insulin sensitivity index, which quantify two clinically and physiologically relevant features [34,35], and to date, continues to be used in clinical settings or improved in modeling studies. Insulin kinetics, including both first phase and second phase insulin secretion, reflects pancreatic responsiveness. Insight into an individual's glucose tolerance or intolerance are obtained via the estimation of pancreatic responsiveness, glucose disappearance, and insulin sensitivity. The minimal model has been extended and is now widely used in various experimental settings. The software MINMOD based on the minimal model is used by clinicians and researchers who are interested in quantifying insulin sensitivity and beta-cell responsiveness [42]. However, fewer mathematical models have been proposed that link insulin, glucose, and FFA [32,40,43,44]. The model developed by Chow and his colleagues [40] captures the interactions between remote insulin X(t), glucose G(t), and FFA F(t). Glucose enters the body intravenously and then it is removed from immediate use by other tissues at a constant glucose effectiveness rate S G or via insulinmediated removal modeled by the interaction term S I XG, where S I represents insulin sensitivity. A proportion of insulin is available for use at rate c X , while the remainder is removed either by natural degradation or by the kidney and liver. The maximal lipolysis rate is given by l 0 + l 2 . Insulin inhibition of lipolysis is denoted by X 2 with the exponent A, and the clearance rate of FFA is denoted by c f . The system of equations describing these dynamics is given by where I(t) represents the insulin concentration in the body over the time t. In minimal model and its siblings in ordinary differential equation systems, the known physiological delay of insulin secretion into the body in response to the rise of glucose level is incorporated using compartment-split techniques, where the insulin compartment is split into two equations I(t) and X′(t) leading to higher dimensional systems of ordinary differential equations [34,35,41,[45][46][47]. The Model (2.1) is next modified via the incorporation of an explicit fixed timedelay that is explicitly linked to insulin secretion.

Explicit time delay model of glucose, insulin, and FFA
A time delay for the insulin secretion in the body is a key physiological factor in the endocrine regulation of insulin, glucose and FFA, since without insulin is essential for lowering glucose levels and suppressing the production of FFA during the hyperglycemic state [7]. This delay can be incorporated explicitly in a system of delay differential equations in order to investigate more realistic intrinsic phenomena in the biological process [39,41]. More recent models incorporating an explicit time delay provide more accurate quantification of insulin sensitivity and glucose effectiveness since these models are more robust [39,41]. Mathematical analysis and numerical simulations [39,48,49] show that these models generate results that can match the observations generated by clinical studies [50,51,55], both in the short-and long-term dynamics (see work by [39,48,[52][53][54]). The above model is an extension of the insulin and glucose model studied in [39] and the equation representing FFA was adapted from the work in [40]. We hypothesize that when insulin ineffectively regulates FFA, then higher FFA may reduce glucose transport, leading to a harmful cycle promoting hyperglycemia and thus contributing to the development of insulin resistance in the long-term.
Next, we consider the interplay of glucose denoted by G(t), insulin denoted by I(t), and FFA denoted by F(t). The parameter (S i I b + S g )G b represents the rate of constant average hepatic glucose input in the short dynamics, that is, in the beginning of the IVGTT protocol. The elevated glucose is either immediately used from other cells at a rate S g or by insulinmediated removal by the interaction term S i GI at the rate S i representing insulin sensitivity.
Both glucose-stimulated and FFA-stimulated insulin production follow the dynamics of Hill's function, β , respectively. Here τ represents the explicit time delay for glucose-stimulated insulin secretion. Insulin has a natural degradation rate d i . For FFA, g 0 + g 1 represents the maximal lipolysis rate. The concentration of insulin inhibition of lipolysis is I 2 with the exponent k , and the FFA clearance rate is denoted d f . The model is illustrated in Figure 2.
Hence, this system is given by with positive initial conditions, where the parameters β and γ represent the Hills function coefficient and σ and α represent the values of half-saturation, respectively. The descriptions of the parameters for Model (2.2) are shown in Table 3.

Mathematical treatment of explicit time delay model
Basic properties of the model. In this paper we only show some biologically relevant results of Model (2.2). First we show that the basal levels of glucose, insulin and FFA are the only equilibrium point of the model. Proof. In fact, it is obvious that (G b , I b , F b ) is a steady state of Model (2.2). Suppose Model (2.2) has another steady state (G * , I * , F * ), a solution of the following system We will show that G* = G b , I* = I b and F* = F b . If G* > G b , then (3.1) implies that I* < I b , and thus F* > F b by (3.3). Therefore G* > G b , F* > F b and (3.2) lead to a contradiction I* > I b .
If G* < G b , then (3.1) implies that I* > I b , and (3.3) implies that F* < F b . Again G* < G b , F* < F b and (3.2) lead to a contradiction I* < I b . Assume that G(t) is non-positive for some t, then there must exist a t 0 > 0 such that G(t 0 ) = 0 and G(t) > 0 for 0 ≤ t < t 0 . Moreover, then dG t 0 dt ≤ 0, which contradicts the following Therefore G(t) > 0 for all t > 0. Similarly, for I(t), assume that ∃ t 1 > 0 such that I(t 1 ) = 0 and I(t) > 0 for 0 ≤ t < t 1 . Then, Therefore I(t) is positive for t > 0. Finally, the similar and standard approach as above ensures F(t) > 0 for t > 0.
For the boundedness, clearly G′(t) ≤ (S i I b + S g )G b -S g G, and thus 0 ≤ G ≤ ( , and I(t) < (σ 1 + σ 2 )/d i for t > 0 by standard treatments.
Remark 1. Theorems 1 and 2 assure that the system of equations in Model (2.2) is well posed, that is, it supports positive bounded solutions under any positive initial condition.

Equilibria and local stability
The steady state is obtained by setting equations in Model 2.2 equal to 0 shown in Eqs (3.1)-(3.3). Rearranging terms and substituting G* and F* into Eq (3.2), yields the equilibrium point implicitly in terms of I*: In order to determine the number of roots in the system, let us consider: It is easy to see that y 1 (I*) and y 2 (I*) must intersect once, which indicates that (3.6) has at least one positive root. Further, we can show that y 1 (I*) < 0, and therefore the system (2.2) has a unique steady state. Let us define z 1 I* = . Substituting z 1 (I*) and z 2 (I*) into (3.6) yields, Next, obtaining d dI* y 1 I* gives, It is easy to see that d dI * y 1 I * < 0, and since σ 2 << σ 1 , then this statement is true and there is one unique steady state. This conclusion is consistent with findings of other glucoseinsulin regulation models in the literature and demonstrate that our model is well-posed [32,34,35,[39][40][41]43,44].
We now turn to study the stability of the unique steady state. It is straightforward that the characteristic equation of Model (2.2) is given by Case with no delay. We analyze the local stability of a positive equilibrium point E b for our system with no time delay by evaluating Eq (3.7) with τ = 0, which gives Applying Routh-Hurwitz Stability Criterion [59] for a cubic polynomial, that is, for the cubic polynomial: where all a i are positive. The Routh array is so the condition that all roots have negative real parts is a 1 a 2 > a 0 a 3 . Therefore, in this case, the equilibrium point is asymptotically stable if Case with delay. Next, we investigate the stability for the case τ > 0. Notice that the equilibrium point E b is stable for τ = 0. If there is some τ > 0 such that E b is unstable, then the characteristic Eq (3.7) must have a pair of pure imaginary roots ±wi with w > 0 [48,49,56,57]. Thus the characteristic equation becomes After algebraic rearrangement, Then we have, which leads to Let u = w 2 > 0, then According to Descartes' Rules, Eq (3.10) has no positive root if c 2 , c 1 , c 0 > 0. Numerically, we test these conditions and observe that the conditions are met for all groups (see Table 4). Therefore we summarize the above results in the following.

Results of model validation and parameter estimation
To determine whether Model (2.2) could capture the metabolic profile of an individual, we fit the model to the clinical data in order to estimate the physiologically relevant parameters. Data were obtained from a study [58], where an IVGTT was given to N = 38 patients (14 men and 24 women) with severe obesity who underwent bariatric surgery and varied metabolically. The control group was neither obese nor diabetic and did not undergo bariatric surgery. It consisted of 6 men and 6 women (N = 12). The treatment group consisted of 14 men and 24 women (N = 38) and included severely obese individuals divided into 3 groups: First, the normal fasting glucose (NFG) tolerance group (N = 9); second, the impaired fasting glucose (IFG) group (N = 17); and third, the T2DM group (N =12). Participants in the control and NFG groups were in similar age groups: 36.7 ± 1.9 and 35.9 ± 3.4 years, respectively. Similarly, the IFG and T2DM groups were closer in age: 45 Plasma glucose, insulin, and FFA levels were measured by the IVGTT before and seven months after undergoing bariatric surgery (N = 38). The control group did not undergo bariatric surgery and were sampled at baseline and at follow-up in seven months. The IVGTT study protocol was performed after a 10 to 12 hours fast. Baseline blood samples were obtained before glucose administration at 15, 10, and 5 min time marks. At time 0 glucose was administered (50% dextrose; 11.4 g/m2 body surface area) in less than one min. Blood samples were obtained after 2, 3,4,5,6,8,10,12,14,16,19,22,25,30,40,50,60,70, 80, 100, 120, 140, 160, and 180 min to measure glucose, insulin, and FFA concentration levels [58]. Here the minimal model (Model (2.1)) and the explicit time delay model (Model (2.2)) are fit to a dataset [58] and model parameters were estimated using the function fmincon for nonlinear programming in MATLAB R2019a. The data are the average values of multiple individual data. The data was extracted from the paper using "Plot Digitizer." Initial conditions from the data are considered for time t = 0.

Validation of the explicit time delay model
The performance of each model on predicting the qualitative trends of the clinical data were determined through many computational studies and numerical simulations. The estimated model parameters for the explicit time delay model is shown in Table 5 and simulations are shown in Figure 3. Qualitatively, comparing the plasma glucose levels pre-and post-bariatric surgery, show a significant improvement, where the NFG, IFG, and T2DM data nearly resembles the glucose trends of the control group. The overall glucose trends for each group pre-and post-surgery were captured by Model (2.2) (shown in Figure 3). The constant rate of glucose effectiveness (S g ) and constant insulin sensitivity rate (S i ) had insignificant changes in the control group but significantly improved (e.g. increased) for the NFG, IFG, and T2DM groups.
The plasma insulin levels significantly improved comparing pre-and post-bariatric surgery results, where the NFG, IFG, and T2DM data closely match the overall trends of the control group. Model (2.2) captures the overall trends for insulin (shown in Figure 3). The constant insulin degradation rate (d i ) decreased from post-compared to pre-surgery for all groups and remained constant for the control group. The maximum secretion rate (σ 1 ) decreased postsurgery for IFG and T2DM but not for NFG nor control groups. The constant rate of FFAstimulated insulin secretion (σ 2 ) increased for all groups except for T2DM.
We also evaluated changes in the estimated explicit time delay parameter τ. As expected, the time delay in the control group pre-and post-surgery were similar since they did not have surgery. The estimated time delay decreased after surgery in the NFG and IFG groups, which supports prior findings that have observed improvements in acute insulin response and beta-cell function after bariatric surgery [25,31]. In contrast, an overall increased time delay was observed in patients with T2DM. One possible explanation is that the treatment strategy implemented for changes and/or the presence of other clinical health factors may explain variations in these metabolic parameters.
The variance for plasma FFA levels significantly reduced comparing pre-and post-bariatric surgery results, where the NFG, IFG, and T2DM data nearly overlap. The model results qualitatively match the overall trends, except for the NFG group (shown in Figure 3). The maximal lipolysis rate (g 0 + g 1 ) reduced significantly from pre-to post-surgery for the IFG and T2DM groups. The insulin inhibition (I 2 ) rate also decreased post-surgery compared to pre-surgery in the NFG, IFG, and T2DM groups, whereas no changes in FFA clearance rate (d f ) were observed.

Comparison of minimal model and explicit time delay model
Model validation and parameter estimation were completed for Model (2.1) and Model (2.2) with the estimated model parameters for each model shown in Tables 5 and 6, respectively. The numerical results are shown in Figure 3. Both models qualitatively predicted the overall trends for the control, NFG, IFG, and T2DM groups of plasma glucose. In this case, both models fit the data well. To determine the overall goodness of fit, the Akaike information criterion (AIC) was calculated (see Table 7). The AIC allows us to determine the overall quality of the statistical model for a given set of data and estimated parameter values. The delay model performed better for the NFG (Total AIC = 1014.47), IFG (Total AIC = 811.09), and T2D (Total AIC = 799.56) group settings in comparison to the minimal model (see Table 7). In summary, the AIC for all groups combined was lower for the delay model (AIC = 2625.12) than in the minimal model (AIC = 2910.85).

Discussion
We investigated the efficacy of insulin suppression on lipolysis and assessed the hypothesis that FFA-stimulated insulin secretion might play a vital role on the progression of insulin resistance. While the role of FFA-stimulated insulin secretion on the progression of insulin resistance is known, we formulated a mathematical model that incorporates FFA-stimulated insulin secretion explicitly and were able to quantify differences in this effect pre-and postsurgery. In Table 5, we found that ˙1 and estimated pre-surgery is much greater than the estimated values post-surgery in the impaired fasting glucose (IFG) and T2DM settings. Hence, our findings indicate that the effect of FFA-stimulated insulin secretion on insulin and glucose levels in patients with IFG and T2DM resembled the estimates obtained for the control group indicating that the bariatric surgery patients could recover moderate insulin secretion activity stimulated by FFA. Similarly, the efficacy of insulin inhibition of FFA. Similarly, I 2 , which represents the efficacy of insulin inhibition of FFA production, on the regulation of FFA, insulin and glucose in patients with NFG, IFG, or T2DM resembled the control group after surgery. These findings indicate that the reduction in weight or body mass due to surgery improve insulin action. This reinforces what we know, that insulin action is essential for regulating plasma FFA and glucose levels. Our results highlight the robust effect of insulin, that is, it can be observed not only in the short-term dynamics during an IVGTT protocol but also in the long-term in longitudinal epidemiological studies.
Overall, it seems that the inclusion of the explicit time delay in Model (2.2) was able to capture the overall qualitative trends well except for the NFG group when compared to Model (2.1). Particularly, Models (2.1) and (2.2) captured the peak and steady states for glucose concentration levels pre-and post-surgery in all groups. Both models captured the initial peak in insulin levels as well as the steady states for the control group. Though the steady states and rapid decline of insulin levels predicted by the models were not exactly equivalent to the averaged data for the NFG group, both models were able to predict the rise and decline in insulin with (Model (2.1)) and without (Model (2.2)) time delay. As shown for the IFG group, the explicit time delay in Model (2.2) captured the two peaks in insulin, which has been shown in previous studies [39,41]. Further, both models captured overall rise, decline, and steady state for insulin levels in the IFG and T2D groups. In contrast, the models fit FFA levels over time less precisely. Among the control group, the initial high basal levels and decrease in FFA were predicted by both models. Among the IFG, NFG, and T2D groups the fit to the FFA data were less precise. However, both models had an initial decline and rise of FFA levels until reaching steady states. Overall, obtaining the final steady states for FFA were difficult to compute across all groups due to the complexity of the model as well as large variance in study participants that could be driven by several biological factors, such as age, gender, stage of diabetes or diabetes-related conditions, and other factors.
The results show that the incorporation of an explicit time delay in the model led to better approximations of the qualitative dynamics when compared to those generated by the minimal model. In previous studies [22,38], a key feature of τ is not only the biological interpretation given by quantitatively estimating the fixed time delays across patients and comparing their differences, but also its ability to capture different qualitative behaviors (or peaks) that may give insight into beta-cell function among patients with and without impaired fasting glucose. In this study, it was shown that the explicit time delay model could capture the two peaks for the averaged data of the impaired fasting glucose group, but the interpretation of τ for a sample average of this group and between groups is limited due to the data being averaged since the proposed model is intended to model the dynamics of a single person, and not a group of individuals. In our future work, we plan to analyze individual patient data and explore the heterogeneity in estimates of τ across different groups based on metabolic health, gender, and age.
The novelty of our work is proposing a dynamic model of glucose, insulin, and FFA with the inclusion of an explicit fixed time delay. Additionally, our model was constructed to directly study the harmful cycle that is activated when FFA is not regulated properly by insulin and, in turn, FFA continues to circulate in the bloodstream, eventually promoting insulin resistance. Though we used averaged data for our model validation, our computational studies demonstrated that the model qualitatively captured the glucose, insulin, and FFA levels trends over time after an IVGTT protocol. However, the use of average data instead of individual patient data limits the conclusions of this model, and hindered our ability to measure insulin sensitivity for each individual, which is needed to understand and assess insulin suppression by FFA. Individual patient data is needed to better assess models' results since the use of averages eliminates the variance that needs to be quantified in each group. Further work is needed to fit more data to evaluate this observation in depth at the patient level. That is, the variabilities in glucose effectiveness, insulin suppression, and insulin sensitivity cannot be assessed at the individual level using average data between groups, and therefore, patient level data is needed. More specifically, the dynamics of insulin, glucose, and FFA need to be better understood for individuals within and across metabolic groups (e.g., nondiabetic or prediabetic, type 2 diabetic) in order to understand the conditions in which FFA is suppressed effectively and in which cases it is not, and how bariatric surgery impacts the FFA dynamics. Lastly, the physiological parameters adjusted to fit the data matched findings from the literature on bariatric surgery, an indicator that the model captures some of the observed phenomena. Future work would require individual patient data and possibly a new model that better captures the impact of bariatric surgery on insulin regulation of glucose and FFA. Furthermore, while the results of our analyses were similar to the work shown in [40], our work provided new insight by the inclusion of τ in our model as well as other terms to explore the effect of FFA on insulin action. Additionally, models modified to include specific mechanisms corresponding to different bariatric surgery types would also provide greater insight into the roles of FFA and how it is impacted by bariatric surgery.

Figure 1.
A schematic diagram illustrating insulin, glucose, and FFA regulation after a meal. Adapted from [33].

Figure 2.
A schematic diagram illustrating the mathematical model of glucose, insulin, and FFA adapted from previous work shown in [40]. Numerical simulations for Model (2.1) and (2.2) are fit to the data for parameters summarized in Tables 2 and 3, respectively. A description of estimated values for both models can be found in Evidence supporting the mechanistic role of FFA on diabetes progression.

Organ Observation [Reference]
Skeletal Muscle High FFA disrupt the insulin signaling process [11,[16][17][18] Higher FFA interfere with the action of insulin to skeletal muscle or hinder insulin signaling, and reduce glucose transport [14,19,20] Liver High FFA levels increase hepatic glucose production in diabetes [18,21] Excessive endogenous glucose production increase when FFA levels influx rise in the liver from lipolysis of visceral adipose depots [10,14,18,22] Pancreas Prolonged high FFA levels impair insulin secretory function and have toxic effects (e.g., "lipotoxicity hypothesis") on pancreatic beta-cells [12,14] Adipose Suppressed inhibitory effect of insulin on lipolysis increases FFA levels [18,21] Increased release of FFA from adipocytes can induce IR [18] Imbalanced production of adipokines (or cell signalling proteins secreted by adipose tissue) promote IR [18,23] Murillo et al.