Modelling endogenous insulin concentration in type 2 diabetes during closed-loop insulin delivery

Background Closed-loop insulin delivery is an emerging treatment for type 1 diabetes (T1D) evaluated clinically and using computer simulations during pre-clinical testing. Efforts to make closed-loop systems available to people with type 2 diabetes (T2D) calls for the development of a new type of simulators to accommodate differences between T1D and T2D. Presented here is the development of a model of posthepatic endogenous insulin concentration, a component omitted in T1D simulators but key for simulating T2D physiology. Methods We evaluated six competing models to describe the time course of endogenous insulin concentration as a function of the plasma glucose concentration and time. The models were fitted to data collected in insulin-naive subjects with T2D who underwent two 24-h visits and were treated, in a random order, by either closed-loop insulin delivery or glucose-lowering oral agents. The model parameters were estimated using a Bayesian approach, as implemented in the WinBUGS software. Model selection criteria were used to identify the best model describing our clinical data. Results The selected model successfully described endogenous insulin concentration over 24 h in both study periods and provided plausible parameter estimates. Model-derived results were in concordance with a clinical finding which revealed increased posthepatic endogenous insulin concentration during the control study period (P < 0.05). The modelling results indicated that the excess amount of insulin can be attributed to the glucose-independent effect as the glucose-dependent effect was similar between visits (P > 0.05). Conclusions A model to describe endogenous insulin concentration in T2D including components of posthepatic glucose-dependent and glucose-independent insulin secretion was identified and validated. The model is suitable to be incorporated in a simulation environment for evaluating closed-loop insulin delivery in T2D. Electronic supplementary material The online version of this article (doi:10.1186/s12938-015-0009-5) contains supplementary material, which is available to authorized users.

and efficacy of the closed-loop system applied in adolescents [2], adults [3] and pregnant women [4] with type 1 diabetes (T1D). Home trials have also shown favourable results of improved glucose control in participants with T1D utilizing closed-loop during their normal daily life [5,6]. Further clinical investigations are in the pipeline aiming to extend the application of closed-loop to benefit a larger population with type 2 diabetes (T2D). Recently reported clinical study tested the feasibility and safety of closed-loop in insulin-naïve T2D subjects [7] and showed promising results.
Clinical studies in humans are costly and time consuming. A virtual testing environment where a population of synthetic subjects with diabetes is tested in a virtual computer space provides a plausible alternative to assess the closed-loop performance in a more time-saving and cost-effective way. It also allows to test scenarios which might not be ethical to perform in real human subjects, so that it gives insight into the limit and a better design of the underlying control algorithm [8,9]. The simulator usually incorporates, as a key component, a glucoregulatory mathematical model and a number of virtual subjects represented by individualized model parameters. However, most of the present simulation models have been developed for T1D [9]. With the growing interest of bringing closed-loop to people with T2D [10], a simulation environment with a glucoregulatory model reflecting glucose-insulin interactions in T2D is desired. In order to create a virtual population of subjects with T2D, a simulation model developed for T1D needs to be supplemented with a model of endogenous insulin secretion and ensuing endogenous insulin concentration.
With this objective in mind, we aimed to develop a posthepatic insulin concentration model and evaluate it with insulin/glucose data collected during a clinical trial in T2D subjects [7]. Six competing models were proposed and the best model was selected on the basis of our model selection criteria. The model of choice assumes a linear relationship between glucose-dependent posthepatic insulin secretion and plasma glucose, and incorporates two two-segment linear functions representing glucose-independent posthepatic insulin secretion after breakfast and lunch, respectively. The model successfully described endogenous insulin concentration data in 11 subjects with T2D, who underwent two 24-hour visits and were treated by either closed-loop insulin delivery or glucose-lowering oral agents.

Subjects and experimental design
We utilised data obtained from a clinical study involving 11 subjects [6 male, age 59.7 (12.1) years, BMI 30.1 (3.9) kg/m 2 , diabetes duration 8.0 (6.2) years and HbA1c 8.3 (0.8)%, mean (SD)] with non-insulin treated T2D, who underwent two 24-hour visits (from 9:00 a.m. on Day 1 until 9:00 a.m. on Day 2), four to six weeks apart, in a random order, treated by either closed-loop insulin delivery or glucose-lowering oral agents [7]. The study was approved by the Cambridge Research Ethics Committee and subjects signed informed consent prior to the commencement of the study.
During closed-loop visits, subjects' routine diabetes therapy was discontinued the night before and replaced with model predictive control algorithm-driven subcutaneous insulin pump delivery of rapid-acting insulin analogue lispro (Eli Lilly, Indianapolis, IN), based on real-time continuous glucose monitoring. Meals were unannounced to the control algorithm and no prandial insulin bolus was given. During control visits, usual diabetes regimen was continued (metformin 92%, sulfonylureas 58%, DPP-4 inhibitors 33%). On both visits, subjects consumed matched 50-80 g carbohydrate meals at 0, 240 and 540 min relative to 9:00 a.m. on Day 1, and optional 15 g carbohydrate snacks. Blood samples were collected, every 15 minutes for measuring glucose concentration, and at the following time points: 0, 15

Insulin and glucose measurement
Samples were centrifuged immediately with plasma kept on ice and stored at -80°C until assayed. Endogenous plasma insulin was measured by AutoDELFIA immunoassay (Perkin Elmer Life Sciences, Wallac Oy, Turku, Finland; inter-assay CV 3.1% at 29 pmol/l, 2.1% at 79.4 pmol/l, 1.9% at 277 pmol/l, 2.0% at 705 pmol/l) which has zero cross reactivity with insulin lispro administered during the closed-loop visits. Plasma glucose was measured by YSI2300STAT Plus Analyser (YSI, Fleet, Hampshire UK).

Model description
We developed six models of increasing complexity to describe endogenous plasma insulin concentration [I ENDO (t)] as a function of plasma glucose concentration [G(t)] and time. A schematic representation of the six competing models is shown in Figure 1. We incorporate the glucose-dependent model components (in Model 1 to 6) and additionally the glucose-independent components (in Model 4 to 6). Models' mathematical formulations are provided in Additional file 1.
Model 1, the base model, assumes a linear relationship between the rate of posthepatic insulin secretion I s (t) (mU/min) and measured plasma glucose concentration G(t) (mM). This assumption originates from previous studies which demonstrated linearity between the prehepatic insulin/C-peptide secretion and the plasma glucose concentration [11,12]. I ENDO (t) (mU/l), the measured endogenous plasma insulin concentration and model output, is assumed proportional to I s (t), through the inverse of the product of insulin metabolic clearance rate MCR I (l/kg/min) and subject's body weight W (kg). This assumption is based on the relatively short plasma insulin half-life (~5 min) compared to the sampling frequency in our study (15 to 60 min), allowing us to assume an instantaneous equilibration between posthepatic insulin appearance in plasma and plasma insulin [13]. The individual parameter values of MCR I are taken from our previous study [13]. In Eq. (1), Additional file 1, G b (mM) is the fasting plasma glucose concentration; M I (mU/min/mM) is the posthepatic glucose sensitivity, representing the effect of unit change in blood glucose concentration on posthepatic insulin secretion and M 0 (mU/min/mM) is the basal glucose sensitivity, representing effect of fasting plasma glucose on posthepatic insulin secretion. values (see Eq. (2), Additional file 1). In Eq. (2), Additional file 1, M I turns into four parameters with subscripts b, l, d and f, representing the four time intervals. In the same fashion, M 0 adopts four subscripts, with M 0,b being the estimated parameter from which values of M 0,l , M 0,d and M 0,f can be derived. This derivation is to satisfy the assumption that the time course of G(t), Is(t) in a two-dimensional space is continuous throughout the study period. t l , t d , t f and t end are at 240, 540, 960 and 1440 min, respectively; these values are assigned according to meal times.
Model 3 expands on Model 2 by further dividing the post-breakfast period into two sub-periods with an additional breakpoint at t b = 90 min, highlighting the difference of the first meal of the day after an overnight fast.
The weighted residuals of Model 3 ( Figure 2) show postprandial underestimation of plasma insulin level indicating that underlying factors other than glucose stimulation influence the posthepatic insulin secretion. Thus, in Model 4, a two-segment piecewise linear function representing additional posthepatic insulin secretion after breakfast I add, b (t) is added. This additional flux assumes value 0 at t = 0, peaks at t=t peak,b = 30 min and falls gradually to 0 at t=t end . The peak value I peak,b (mU/min) is estimated.
The weighted residuals of Model 4 ( Figure 2) shows improved post-breakfast model fit but underestimation of insulin concentration still exist after lunch and dinner time. Therefore, in Model 5, another two-segment piece-wise linear function I add,l (t) is added after lunch time. Similar to the post-breakfast flux in Model 3, the post-lunch flux assumes value 0 at t = t l = 240 min, peaks at t=t peak,l =270 min and drops to 0 at t=t end . The peak value I peak,l (mU/min) is estimated.
Finally, in addition to Model 5, Model 6 assumes a derivative action of G(t)on the posthepatic insulin secretion as a number of studies claimed derivative action of glucose concentration on prehepatic insulin secretion [11]. Id(t) (mU/min) is the posthepatic insulin secretion triggered by the positive rate of change of glucose concentration; M d (mU/mM) is the dynamic glucose responsivity representing the ability of glucose rate of change to stimulate posthepatic insulin secretion (Eq. (6), Additional file 1).

Model selection and parameter estimation
Parameters of the six competing models were identified on clinical data obtained from the closed loop periods and the results were used in model selection criteria. The parameters of the model of choice were subsequently identified on data obtained from the control period.
In the first instance, all six models were numerically identified using glucose/insulin data obtained from the closed-loop period, using a Bayesian approach. WinBUGS (Bayesian inference Using Gibbs Sampling for Windows) software version 1.4 (MRC Biostatistics Unit, Cambridge, UK) [14] was used for model identification. The measurement error in I ENDO (t) was assumed to be uncorrelated, normally distributed with zero mean, and with a coefficient of variation at 6%. The plasma glucose concentration represented by G(t), being the model input, was assumed error free.
The initial model identification results were used to select the best model on the basis of the following criteria: (i) the deviance information criterion (DIC) which is a measurement of model parsimony, taking into account the goodness of fit and the model complexity [15], (ii) the weighted residual profiles documenting overall model ability to describe the data, and (iii) the presence of zero estimated parameters which indicate increased model complexity leading to over-fitting. The model of choice was then identified on data obtained from the eleven control periods.
The individual parameter sets of the best model were identified using the Bayes Theorem:  M 0,f , I peak,b and I peak,l . The prior information on the model parameters is represented by their a priori distributions p(θ)which are updated by the likelihood p(y|θ), yielding parameters' a posteriori distributions p(θ|y). The integral ∫p(θ)p(y|θ)dθ acts as a scale factor.
In the present study, vague prior information was assigned to each parameter. For each study period, it took WinBUGS 550 seconds to run 100,000 iterations and produce posterior parameter distributions. A standard 64-bit desktop PC (OPTIPLEXTM 7010, DELLTM Computers Ltd) was used. The first 10,000 iterations were discarded as burn-ins to allow the initial stabilization of the Markov chain [16].

Statistical analysis
Results are shown as median (interquartile range), if not differently indicated. The Wilcoxon signed-rank test was used to assess the between-sample difference due to the non-normal data. SPSS software version 21 [17] was used to carry out the statistical analysis where P values less than 0.05 were considered statistically significant. Table 1 compares the area-under-curve (AUC ins ) of measured endogenous plasma insulin concentration between visits. A significant increase in the control period AUC ins can be observed during the three postprandial time intervals as well as during the fasting conditions (P < 0.05) suggesting an increase in posthepatic insulin secretion throughout the control period. Figure 2 illustrates weighted residual profiles obtained with the six competing models. Model 5 and Model 6 result in comparable profiles and provide the best fit to the data. Table 2 shows results of the DIC analysis. As expected, the deviance D decreased and the 'effective number of parameters' p D increased with the increasing model complexity. Model 6 outperforms the other five models in terms of the lowest DIC value. However, following model identification, we noted that the parameter M d , exclusive to Model 6, was estimated to be zero in 9 out of the 11 tested data sets. Model 6 was therefore rejected and Model 5 became the model of choice as it adequately described the data and held the second lowest, after Model 6, DIC value.     satisfactory model behaviour during both study periods and a positive correlation between plasma glucose and endogenous plasma insulin levels. Table 3     However, I peak,l holds a higher median estimate during the control period compared to the closed-loop visit (P < 0.05). Table 4 compares model-derived AUC ins and its glucose-dependent and glucoseindependent components between visits. The AUC ins is significantly larger during the control period (P < 0.05), with the excess amount only attributable to the glucoseindependent component (P < 0.05). The glucose-dependent AUC ins show no discordance between visits (P > 0.05).

Discussion
We present a model describing the time course of 24-h endogenous plasma insulin concentration in people with T2D. The model was identified on insulin/glucose data obtained from 11 insulin-naïve T2D subjects who underwent two 24-hour visits treated by either closed-loop insulin delivery or glucose-lowering oral agents. Both clinical and model-derived results suggested an increase in the amount of posthepatic endogenous insulin secretion during the control period, namely the oral agents-treated period. The modelling results further suggested that the surplus insulin secretion was glucoseindependent. We consider the model suitable to be incorporated in a simulation environment for testing closed-loop inulin delivery performance in T2D.
In the model development process, six competing models of increasing complexity have been evaluated. All models adopted the following assumptions: a linear relationship exists between plasma glucose concentration and posthepatic insulin secretion rate; MCR I is identical between study visits; the equilibration is instantaneous between posthepatic insulin appearance in plasma and plasma insulin, reflecting a short plasma insulin half-life relative to the frequency of plasma insulin measurements. The latter assumption was adopted in our previous study during which individual MCR I values were estimated [13].
In order to select the model best representing our data sets, three criteria were adopted. Based on the DIC criterion (Table 1), Model 6 would have been chosen as the most parsimonious. However, the parameter estimates of M d converged to zero in 82% of the cases implying that Model 6 collapsed into Model 5 and that the contribution of the derivative action of glucose concentration on posthepatic insulin secretion was overall negligible. Therefore, Model 6 was rejected and Model 5 became our model of choice. In terms of model parameter estimation, we found that, unlike M I,b , M I,l and M I,d , which represent posthepatic glucose sensitivities after breakfast, lunch and dinner and which had non-zero estimates, the median estimates of M I,f were zero for both visits. This may suggest that the linear relationship between plasma glucose and plasma insulin may be less pronounced during the fasting period, probably due to the less variable glucose and insulin levels, compared to postprandial periods. A more sophisticated model structure may be required to describe the 'plasma glucose-endogenous plasma insulin' relationship overnight. Additional file 1 shows plots of measured plasma glucose concentration versus endogenous plasma insulin concentration. Given the assumption of linearity, we should have observed a unique linear dependency between measurements for each of the time intervals. However, in a few subjects, namely in subjects 2, 3, 4, 6, 8 during the closedloop period and subjects 2, 4, 6, 8, 10 during the control period, the endogenous plasma insulin concentration remains high or continues to increase despite a pronounced drop in the plasma glucose level, especially after the first meal of the day. This observation encouraged us to incorporate the glucose-independent insulin secretion component.
The bioactive incretins, namely glucagon-like peptide-1 (GLP-1) and glucose-dependent insulinotropic polypeptide (GIP), are hormones released from the small intestine in response to oral meal intake, acting as enhancers to endogenous insulin secretion. Woerle et al. demonstrated that incretins play a relatively more important role in amplifying insulin secretion in people with T2D compared to healthy subjects [18]. We therefore hypothesize that the glucose-independent insulin secretion functions in our model may account for, at least in part, the insulinotropic effect of incretin hormones, and may correspond to the 'potentiation' of insulin secretion after meal proposed by Mari et al. [11].
The excess amount of glucose-independent insulin secretion during the control period revealed by the modelling process may be attributed to the effect of glucoselowering agents. As suggested in literature, metformin (taken by 92% of the studied subjects during the control period) enhances GLP-1 action [19], sulfonylureas (taken by 58% of the studied subjects) stimulates endogenous insulin secretion [20] and DPP-4 inhibitors (taken by 33% of the studied subjects) improves insulin secretion and prolongs GLP-1 action [21].
Insulin and C-peptide are co-secreted equimolarly by the pancreatic beta-cells. After travelling through the liver, about 50% of insulin is extracted before entering the plasma while C-peptide is not affected. The prehepatic insulin secretion profile can, therefore, be reconstructed using a model of C-peptide kinetics [22] where the beta-cell function can be assessed quantitatively [11,12,23]. However, the objective of this study was to describe posthepatic insulin concentration with a relatively simple model. Thus, plasma C-peptide concentration was not measured in the present study.
The model can be implemented in a simulation environment designed specifically for accelerating the development of closed-loop insulin delivery systems for T2D. The endogenous insulin secretion characteristics of the virtual subjects with T2D will be represented by individual model parameters identified in the present study. The model may also be used in pharmaceutical development. In a possible application, model parameters identified using data collected before and after the ingestion of a tested glucose-lowering agent could be compared. The comparison may help to evaluate the efficacy of a new drug to influence, for instance, the posthepatic glucose sensitivity M I . The model's ability to quantitatively evaluate effects of altered glucose turnover on endogenous insulin secretion, when patients' usual glucose-lowering oral agent is replaced with subcutaneous insulin therapy, may also be of clinical relevance [24].