Incorporating Unannounced Meals and Exercise in Adaptive Learning of Personalized Models for Multivariable Artificial Pancreas Systems

Abstract

Background:

Despite the recent advancements in the modeling of glycemic dynamics for type 1 diabetes mellitus, automatically considering unannounced meals and exercise without manual user inputs remains challenging.

Method:

An adaptive model identification technique that incorporates exercise information and estimates of the effects of unannounced meals obtained automatically without user input is proposed in this work. The effects of the unknown consumed carbohydrates are estimated using an individualized unscented Kalman filtering algorithm employing an augmented glucose-insulin dynamic model, and exercise information is acquired from noninvasive physiological measurements. The additional information on meals and exercise is incorporated with personalized estimates of plasma insulin concentration and glucose measurement data in an adaptive model identification algorithm.

Results:

The efficacy of the proposed personalized and adaptive modeling algorithm is demonstrated using clinical data involving closed-loop experiments of the artificial pancreas system, and the results demonstrate accurate glycemic modeling with the average root-mean-square error (mean absolute error) of 25.50 mg/dL (18.18 mg/dL) for six-step (30 minutes ahead) predictions.

Conclusions:

The approach presented is able to identify reliable time-varying individualized glucose-insulin models.

Improved glycemic control in individuals with type 1 diabetes mellitus (T1DM) via the use of artificial pancreas (AP) systems, which integrate a continuous glucose monitoring (CGM) sensor, continuous subcutaneous insulin infusion (CSII) pump, and insulin dosing control algorithm, reduces the risks of immediate life-threatening conditions, such as severe hypoglycemia and ketoacidosis, and long-term health complications, such as cardiovascular disease, nephropathy, neuropathy, and retinopathy.^1-13 A fully automated AP system eliminating the need for patients with T1DM to enter user inputs for meal and exercise announcements represents a substantial step toward achieving better insulin delivery systems. Although AP systems are shown to be effective compared to conventional multiple daily insulin injections or sensor augmented pump therapy, improvements in the insulin dosing algorithms are needed if the AP systems are to compensate for meals and exercise without requiring patients to manually interact with the systems. Model-based predictive control has emerged as an effective insulin dosing algorithm for AP systems as it inherently considers future projections of the glycemic dynamics for predictive hypo- and hyperglycemia alarms and manipulating insulin delivery when glucose concentrations are forecast to deviate from the desired outcomes. The ability of the predictive controllers is dependent on accurate models of the metabolic system that can effectively encode a comprehensive understanding of the physiology of T1DM, handle the long dynamical effects of the insulin action, and characterize the effects of various diurnal disturbances on glycemic dynamics. Furthermore, many mathematical models of the glucose and insulin dynamics employ real-time feedback solely from the CGM signal, which limits their ability to automatically accommodate meals and exercise without manual interventions by people with T1DM.

One of the challenges to achieve a fully automated and reliable AP is the lack of an accurate model to represent the dynamic changes in physiology under various conditions. The future glucose measurements are difficult to predict accurately as glycemic dynamics vary substantially due to the effects of numerous factors such as current and historical glucose trends, carbohydrate consumption, previously administered insulin, exercise or physical activity levels, and concentrations of certain hormones. Moreover, metabolic processes vary substantially among subjects and temporally within people due to the diversified lifestyles and erratic routines of individuals. The diversity of physiology and behaviors in people causes the glucose-insulin dynamics and insulin sensitivities of individuals to vary over time, while large perturbations such as meals and physical activity cause significant excursions in glucose concentrations that may not be described accurately by generalized models. Hence, the necessity in AP systems of personalized glucose-insulin models rather than models with generalized parameters that do not reflect the dynamic characteristics of subjects in different situations.^14-19

The glycemic models proposed in literature can be categorized as either physiological or data-driven models, often with carbohydrate intake and infused insulin as inputs. Physiological models, usually based on compartmental models, consist of simultaneous differential equations describing the insulin and glucose metabolism with a number of physiological parameters to be identified. Detailed physiological models are very useful in simulations, though model complexity and computational load impede their real-time implementation in AP systems. In contrast to physiological models, data-driven models with relatively simpler structures that effectively characterize the intricate relationships among the measured variables generally require less computational effort. ¹⁵ Among the various empirical modeling approaches proposed to predict future glycemic measurements, subspace-based system identification methods are capable of efficiently identifying linear state-space models from multi-input, multi-output sampled data of a dynamic system. Nevertheless, a fixed empirical model may not predict glycemic measurements well in all scenarios and across all subjects due to the significant variability in glycemic dynamics. Therefore, the empirical models need to be appropriately adapted on-line to characterize the current dynamics of the individuals and make accurate short-term predictions of glucose concentration measurements. For this purpose, adaptive system identification approaches are proposed to determine linear, time-varying models and effectively characterize the evolving glycemic dynamics, thus allowing the adaptive models to be valid over a diverse range of daily conditions.

Despite the predictive ability of the adaptive system identification approach, the discrete basal changes and acute bolus impulses present a challenging configuration for empirical modeling techniques. Further complicating the model identification procedure is the fact that administered insulin (basal or bolus) gradually accumulates in the bloodstream and is eventually utilized by the body. One of the factors that prolongs the utilization of administered insulin is the significant time delays involved in the diffusion and absorption of the subcutaneously injected insulin analogues. For instance, the plasma insulin response to subcutaneously infused fast-acting insulin has a time lag of approximately 15 minutes, peak effect at about 45-90 minutes, and an overall effective duration of 4-6 hours. The protracted effects of the administered insulin may cause glucose concentration measurements to continue rising in response to carbohydrate consumption even though sufficient insulin is already administered. Under such postprandial circumstances where insulin effects are not observed in CGM measurements for lengthy periods of time after insulin infusion, methods relying entirely on data-driven models may lead to inappropriate generalizations with regards to the effects of insulin. This well-recognized artifact of modeling glycemic measurements can be addressed by filtering the administered insulin dose into a newly constructed variable that readily accommodates the effects of previously administered insulin. It is well recognized that filtering the input variables has a significant effect on the predictive performance of the identified model. Although such schemes are shown to improve glycemic prediction accuracy, the filtering algorithms utilized, whether based on numerical signal processing techniques or derived from compartment models describing the underlying physiological phenomena, are often time-invariant. The application of fixed time-invariant filters thus may be suboptimal and may diminish the potential for improvement, especially when employed simultaneously with adaptive system identification techniques. Utilizing PIC estimates from an adaptive and personalized PIC estimator as the newly established input for the identification algorithm provides a filtered variable that improves prediction ability. In contrast to the conventional insulin on board approach that estimates the amount of active insulin in the body through static appropriations of action profiles and decay curves, the estimated PIC is a physically meaningful and quantifiable measure of insulin the bloodstream that may benefit control design as well. The estimated plasma insulin concentration (PIC) can be considered as a more appropriate input variable in the identification of a dynamic model for predicting future glucose concentrations because it directly quantifies the insulin in the bloodstream that affects the glycemic evolution.^20-24

Besides the effects of administered insulin, another limitation of current AP systems is the requirement of manual meal and exercise announcements. Several studies have incorporated unannounced meals through the estimation of time-varying parameters or analysis of glucose trends.^25-30 Concerning the automation of AP systems, additional physiological variables related to physical activity are also considered to automatically accommodate exercise.^31-34 Despite this progress, automatically handling of unannounced meals and exercise in adaptive and personalized glycemic models for AP systems is not sufficiently studied.

Motivated by the above considerations, in this work, an adaptive and personalized compartment model that translates the abrupt bolus and discrete basal changes into estimates of PIC is integrated with an adaptive system identification approach to characterize the transient dynamics of glycemic measurements. The adaptations by the recursive system identification assist in handling stochastic disturbances like the effects of meal consumption and physical activities. Furthermore, the system identification algorithm is extended to improve reliability for use in fully automated artificial pancreas systems. The proposed approach involves the modification of the recursive predictor-based subspace identification (PBSID) algorithm to incorporate constraints on the fidelity and accuracy of the identified models, correctness of the sign of the input-to-output gains, and the integration of heuristics to ensure the stability of the recursively identified models. To achieve this, the proposed adaptive models also include estimated meal effects and the reported exercise related biometric variables by BodyMedia SenseWear armband as additional inputs to automatically accommodate unannounced meals and exercise. The proposed adaptive and personalized modeling approach considering the effects of unannounced meals and exercise on the transient glycemic dynamics is applied to 15 clinical data sets involving closed-loop experiments of the AP systems.

Methods

The predictor-based subspace identification (PBSID) approach is used to identify a recursively updated state-space model characterizing the glucose-insulin dynamics. The PBSID algorithm estimates a vector autoregressive with exogenous inputs model and uses it to construct matrices related to the state-space model parameters. A low-rank approximation then identifies a state sequence, and the state-space model parameters are recursively updated. The PBSID technique is further extended to provide stable, time-varying, and individualized models for glycemic predictions by using estimates of PIC and meal effects, and measured physiological signals. The estimates of PIC and meal effects are obtained by-using an adaptive and personalized PIC estimator designed based on the Hovorka’s glucose-insulin dynamic model. In this section, a review of the PBSID algorithm for the identification of linear, time-varying state-space models is provided, followed by a brief overview of the adaptive and personalized PIC estimator. Figure 1 shows a flowchart of the proposed personalized and adaptive modeling algorithm.

Figure 1.

Flowchart of the proposed personalized and adaptive modeling algorithm.

Adaptive Modeling Approach

In this work, a recursive subspace-based system identification approach is used to build linear, time-varying state-space models of the form:

$${\mathrm{x}}_{\mathrm{k}+1}={\mathrm{A}}_{\mathrm{k}}{\mathrm{x}}_{\mathrm{k}}+{\mathrm{B}}_{\mathrm{k}}{\mathrm{u}}_{\mathrm{k}}+{\omega }_{\mathrm{k}}$$

$${\mathrm{y}}_{\mathrm{k}}={\mathrm{C}}_{\mathrm{k}}{\mathrm{x}}_{\mathrm{k}}+{\mathrm{D}}_{\mathrm{k}}{\mathrm{u}}_{\mathrm{k}}+{\nu }_{\mathrm{k}}$$

(1)

where $\mathrm{x}\in R^{n_x}$ , $\mathrm{u}\in R^{n_u}$ and $\mathrm{y}\in R^{n_y}$ denote vectors of the state, input, and output variables, respectively, with $\omega$ and $\nu$ denoting the vectors of process disturbance and measurement noise (respectively) that are assumed to be zero-mean white Gaussian sequences. The input variables considered are the insulin concentration in the bloodstream, energy expenditure, and the estimated meal effects, and the output variable is the CGM measurement. The system matrices ${\mathrm{A}}_{\mathrm{k}},$ ${\mathrm{B}}_{\mathrm{k}},$ ${\mathrm{C}}_{\mathrm{k}},$ and ${\mathrm{D}}_{\mathrm{k}}$ are the state transition, input, output, and direct feedthrough matrices of appropriate dimensions (respectively) that are adaptively estimated at each sampling instance using subspace identification techniques. A detailed description of the proposed technique is provided in Appendix A, as well as in Hajizadeh et al ¹⁴ and Houtzager et al. ³⁵ The identification procedure employs constraints on the estimated model coefficient matrices to ensure that the appropriate gains are determined relating the insulin to the glucose outputs. In another consideration of the fidelity of the identified models, constraints are also imposed to restrict the identified model to be stable. Consider the system matrices ${\Theta }_{\mathrm{k}}=\left[{\mathrm{A}}_{\mathrm{k}},{\mathrm{B}}_{\mathrm{k}},{\mathrm{C}}_{\mathrm{k}},{\mathrm{D}}_{\mathrm{k}},\right]$ to constrain for fidelity and stability of the identification. Then the constrained optimization problem becomes

$${\Theta }_{\mathrm{k}}^{*}=\underset{{\Theta }_{\mathrm{k}}\in {\Omega }_{\Theta }}{\mathrm{argmin}}{\displaystyle \sum }_{\mathrm{i}=1}^{\mathrm{n}}{\Vert {\mathrm{y}}_{\mathrm{k}}-{\widehat{\mathrm{y}}}_{\mathrm{k}}\Vert }_{\Upsilon }^2$$

(2)

where a weighted quadratic norm of the deviation of the predicted CGM outputs is penalized from the actual CGM measurements while ensuring that the estimated model coefficients satisfy the imposed constraints ${\Theta }_{\mathrm{k}}\in {\Omega }_{\Theta }$ for model fidelity and stability purposes. The constrained optimization problems are solved using the computationally efficient quadratic programming solver, which seeks to minimize a scalar function in multiple variables within a region specified by constraints and bounds. A detailed description of the proposed technique is provided in Appendices B and C.

One modification compared to traditional glycemic models is the use of PIC as an input to the system identification procedure. The insulin present in the bloodstream is derived from the insulin subsystem of Hovorka’s model, which translates the abrupt basal and discrete bolus inputs into PIC estimates with time-varying model parameters simultaneously estimated using a nonlinear observer.^36,37 The observer is of the form

where ${\widehat{\mathrm{x}}}_{\mathrm{k}}$ denotes the estimated state vector, with PIC an estimated state variable, ${\tilde{\mathrm{u}}}_{\mathrm{k}}$ is the infused (basal and bolus) insulin, ${\mathrm{y}}_{\mathrm{k}}$ is the CGM measurements providing feedback correction for the estimation of states and parameters, and ${\widehat{\theta }}_{\mathrm{k}}$ are the model parameters that are estimated simultaneously with the state variables to characterize the time-varying physiological processes of people. The estimated parameters of Hovorka’s model include the time to maximum insulin absorption ( ${\mathrm{t}}_{\max ,\mathrm{I}}$ ), the fractional elimination rate of insulin from the bloodstream ( ${\mathrm{k}}_{\mathrm{e}}$ ), and the gut absorption rate ( ${\mathrm{U}}_{\mathrm{G}}$ ) representing the rate of appearance of glucose in the blood from consumed meals. The on-line estimation of the insulin compartment parameters characterizes the transient pharmacodynamics of subcutaneously administered insulin. The estimation of the glucose appearance rate allows for quantification of the effects of unannounced meals on the glycemic dynamics without requiring user specified information on the amounts of carbohydrates consumed. The easily attainable demographic information of subjects, such as weight, height, BMI, age, total daily insulin and duration with T1DM, is exploited to personalize the initialization of the PIC estimator. To this end, the readily available demographic information is used to identify a relationship between the demographic variables and the model parameters to be individually initialized to each patient. Data-driven approaches, such as partial least squares (PLS) regression models, are widely used when the exact underlying mathematical relationship between two sets of data, demographic inputs and the output parameters, is not explicitly and mathematically formalized through fundamental physiological models. PLS is a multivariate regression method for modeling the relationship between two groups of data consisting of numerous noisy and correlated variables while appropriately handling potentially incomplete measurements with missing data. The PLS based regression relationship is derived from data collected from several experiments. The latent variables of the PLS model are identified to maximize the prediction performance of the model. This is achieved by finding components that maximize the covariance between the independent and the dependent variables. The independent variables matrix consists of the demographic information of the patients, such as bodyweight, height, BMI, total daily insulin dose, and the dependent variables matrix is defined to be the parameters of the PIC estimator that were determined for each patient using actual clinical experimental PIC measurements. Following the initialization of the parameter values using the PLS regression models, the proposed estimation technique can rapidly converge to provide good estimates of the PIC in a timely manner by using only the CGM measurements and infused insulin data. A detailed description of the proposed technique is provided in Appendix D, as well as in Hajizadeh et al.^22-24

$${\widehat{\mathrm{x}}}_{\mathrm{k}+1}={\mathrm{f}}^{\mathrm{obs}}\left({\widehat{\mathrm{x}}}_{\mathrm{k}},{\tilde{\mathrm{u}}}_{\mathrm{k}},{\widehat{\theta }}_{\mathrm{k}},{\mathrm{y}}_{\mathrm{k}}\right)$$

(3)

Subjects and Clinical Study Experiments

The subjects involved in this study were recruited by the Kovler Diabetes Center, University of Chicago Medical Center, Chicago, Illinois, and were scheduled for a visit at the University of Chicago General Clinical Research Center. The subjects included adults 19-39 years with T1DM. All subjects used CSII pump therapy by adjusting the insulin infusion flow rate based on suggestions from a generalized predictive controller.^6,30,31,38 Subjects wore their personal insulin pump, a continuous glucose monitor (CGM; Medtronic Guardian-Real Time Continuous Glucose Monitors [Medtronic, Northridge, CA] or Dexcom G4 Platinum [Dexcom, San Diego, CA]) and a BodyMedia SenseWear Pro3 (BodyMedia, Pittsburgh, PA) armband reporting physiological signals. Each patient’s visit was approximately 60 hours long during the closed-loop experiment using a multimodule multivariable adaptive artificial pancreas system.^24,29,39-42 The subjects’ own insulin type and pumps were used during the experiments. Subjects participated in two exercise bouts of 20- to 30-minute sessions before and after lunch. Overall 10 subjects participated in the 15 clinical closed-loop experiments. Table 1 shows the characteristics of the participants of the closed-loop studies.

Table 1.

Demographic Information of the Clinical Subjects.

Demographic criteria	Mean ± standard deviation
Age (years)	26.13 ± 5.76
Body weight (kg)	79.29 ± 15.40
Height (m)	1.73 ± 0.08
Body mass index (kg/m2)	26.33 ± 4.18
Total daily insulin dose [basal] (U)	45.75 [26.79] ± 12.73 [11.66]
Duration of time with diabetes (years)	18.30 ± 8.10
Waist (cm)	90.72 ± 13.57

Results

The efficacy of the proposed adaptive system identification algorithm for identifying high fidelity glycemic models is demonstrated using 15 clinical datasets. The output of the model is the predicted CGM measurements and the inputs are the estimates of the PIC, meal effect ( $U_{\mathrm{G}}$ ), and energy expenditure values in units of metabolic equivalents (MET). The two performance indices used to evaluate the performance of the proposed method are the mean absolute error (MAE) and root mean square error (RMSE).

The MAE and RMSE are calculated based on the following equations:

$$MAE=\frac{{{\displaystyle \sum }}_{i=1}^n\left|y_i-{\widehat{y}}_i\right|}{n}$$

(4)

$$RMSE=\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^n{\left(y_i-{\widehat{y}}_i\right)}^2}{n}}$$

(5)

The MAE and RMSE are two of the most common metrics used to measure the prediction accuracy. MAE measures the average magnitude of the errors in a set of predictions, without considering their direction. It’s the average over all samples of the absolute differences between prediction and actual observation where all individual differences have equal weight. RMSE is a quadratic scoring rule that also measures the average magnitude of the error. It is the square root of the average of squared differences between prediction and actual observation. Both the MAE and RMSE express average model prediction error in units of the variable of interest. Both metrics can range from 0 to ∞ and are indifferent to the direction of errors. They are negatively oriented scores, which means lower values are better. Taking the square root of the average squared errors has some interesting implications for the RMSE. Since the errors are squared before they are averaged, the RMSE gives a relatively high weight to large errors. This means the RMSE should be more useful when large errors are particularly undesirable.

The quantitative results for all 15 clinical experiments are presented in Tables 2 and 3. These performance indices are computed based on one- and six-step-ahead predictions. The one- and six-step-ahead prediction are used to show the performance of the proposed technique in predicting the future outputs. The average RMSE (MAE) for one- and six-step-ahead CGM prediction are 4.93 (3.17) and 25.50 mg/dL (18.18 mg/dL), respectively. The RMSE (MAE) values for the one-step-ahead predicted CGM are less than 8.44 mg/dL (4.25 mg/dl), and for the six-step-ahead predicted CGM is less than 32.96 mg/dL (22.02 mg/dL), which demonstrates the effective performance of the adaptive system identification technique to model the glycemic dynamics.

Table 2.

Performance Indices Values for the Modeling Results Using the Proposed Personalized and Adaptive Modeling Algorithm for All Clinical Experiments Based on the One-Step-Ahead Prediction.

No.	MAE (mg/dL)	RMSE (mg/dL)	Median (mg/dL)	First quartile (mg/dL)	Third quartile (mg/dL)
1	3.16	4.50	1.96	0.83	3.80
2	2.86	3.80	2.66	1.28	4.24
3	3.22	5.78	1.66	0.62	3.67
4	3.91	5.22	3.03	1.48	5.65
5	2.97	4.15	2.02	0.93	3.65
6	2.56	4.10	1.86	0.78	3.49
7	2.90	4.50	2.03	0.89	3.45
8	3.16	4.27	2.58	1.15	4.67
9	2.71	4.29	1.72	0.80	3.50
10	3.79	6.58	2.37	1.05	4.37
11	3.89	5.79	2.61	1.19	5.02
12	2.88	4.84	1.85	0.72	3.57
13	2.29	3.30	1.61	0.67	3.15
14	4.25	8.44	2.50	1.21	4.48
15	3.06	4.45	2.13	0.97	4.25
Avg	3.17	4.93	2.17	0.97	4.07

Table 3.

Performance Indices Values for the Modeling Results Using the Proposed Personalized and Adaptive Modeling Algorithm for all Clinical Experiments Based on the Six-Step-Ahead Prediction.

No.	MAE (mg/dL)	RMSE (mg/dL)	Median (mg/dL)	First quartile (mg/dL)	Third quartile (mg/dL)
1	18.28	26.01	10.47	4.80	24.86
2	16.45	20.45	14.23	8.44	22.47
3	17.31	25.94	10.58	4.13	21.85
4	21.45	27.29	16.28	7.30	30.03
5	21.66	30.26	13.90	5.76	27.60
6	16.12	21.91	11.70	5.80	22.53
7	17.71	25.22	11.56	5.21	24.62
8	18.36	24.80	13.54	6.59	24.77
9	17.51	26.92	10.65	4.45	22.27
10	19.35	27.09	13.95	7.18	24.56
11	19.08	25.11	14.08	6.84	26.71
12	16.82	25.43	9.89	4.77	18.93
13	14.71	20.96	10.00	4.65	19.32
14	22.02	32.96	14.87	6.69	27.15
15	15.89	22.13	9.95	4.73	19.74
Avg	18.18	25.50	12.38	5.82	23.83

The proposed adaptive system identification method is able to provide stable and adaptive state-space models at each sampling time. Figure 2 shows the predicted and actual CGM values along with the inputs of the model (the estimates of the PIC and meal effect parameter, and the energy expenditure values) and maximum norm of model eigenvalues for a selected experiment. The results for the clinical experiments (Figure 2) show that the maximum norm of eigenvalues of the models is maintained within the unit circle and evolve over time to best describe the available input and output data, thus ensuring stability of the identified models. Furthermore, the estimated and predicted CGM values closely track the actual CGM measurements. The state-space model is identified such that the available input-output data is well characterized by the system realization while the stability and fidelity of the identified models is guaranteed. Compared to only using the PIC as a model input, additional inputs capturing the effects of meals and physical activity, such as the estimated meal effect parameter and energy expenditure, significantly improve the prediction ability of the modeling approach (Figure 3). Although the model using PIC and meal effect as inputs does not show statistical significance compared to the proposed model using PIC, meal effect and energy expenditures as input variables, the standard deviation is biased by the outliers, which influences the statistical hypothesis test. To elucidate this bias caused by the outliers, a one-sided test for statistical significance is used to compare the model (b) with the proposed model (d), which shows statistical significance (P value = 0.04 for MAE and 0.029 for RMSE). The results demonstrate that the use of the PIC as a filtered insulin input, the incorporation of estimates of the meal effect obtained through the nonlinear observer algorithm, and the inclusion of energy expenditure as an auxiliary input from the BodyMedia SenseWear armband improve the identified glucose-insulin dynamic models.

Figure 2.

Modeling results for a select experiment (subject 7) showing (a) the six-step (30 minutes ahead) predicted and actual CGM values, (b) the estimated plasma insulin concentration and actual PIC values (not used in model development), (c) the estimated meal effect parameter for capturing unannounced meals, (d) the energy expenditure in units of metabolic equivalent (MET) for incorporating unannounced exercise, and (e) the maximum norm of eigenvalues of the stable time-varying identified models.

Figure 3.

Comparison of six-step (30 minutes ahead) prediction accuracies across all experiments using the mean absolute error (MAE) and root mean square error (RMSE) for prediction with different model inputs: (a) PIC; (b) PIC and energy expenditure; (c) PIC and meal effect; and (d) PIC, energy expenditure, and meal effect. All hypothesis tests computed using two-sided paired t-tests at 5% significance level and are relative to the proposed model (d), which considers PIC, energy expenditure, and meal effect as input variables.

The Clarke error grid approach is used to assess the clinical significance of differences between glucose predictions generated by the proposed technique and measured values.^43-46 The method uses a Cartesian diagram, in which the six-step-ahead CGM predictions obtained by the presented technique are displayed on the y-axis, whereas the real (measured) CGM values are displayed on the x-axis. The diagonal represents the perfect agreement between the two, whereas the points below and above the line indicate, respectively, overestimation and underestimation of the actual values. Zone A (acceptable) represents the glucose values that deviate from the real values by ±20% or are in the hypoglycemic range (<70 mg/dl), when the reference is also within the hypoglycemic range. The values within this range are clinically exact and are thus characterized by correct clinical treatment. Zone B (benign errors) is located above and below zone A; this zone represents those values that deviate from the reference values, which are incremented by 20%. The values that fall within zones A and B are clinically acceptable, whereas the values included in areas C-E are potentially dangerous, and there is a possibility of making clinically significant mistakes. Based on the results presented in the Table 4, we can see that 98.87% of values fall within zones A and B. There is no cases fall within zone E that is the most dangerous case and just 1.13% of values fall within zones C and D. Figure 4 shows of the Clarke error grid analysis for a select experiment (subject 7).

Table 4.

Percentages in Each Zone Based on Clarke Error Grid Analysis for the Modeling Results Using the Proposed Personalized and Adaptive Modeling Algorithm for All Clinical Experiments Based on the Six-Step-Ahead Prediction.

No.	A	B	C	D	E
1	77.30	21.44	0.00	1.26	0.00
2	79.67	19.42	0.91	0.00	0.00
3	85.43	14.37	0.20	0.00	0.00
4	66.48	31.25	0.00	2.27	0.00
5	77.08	21.48	0.18	1.26	0.00
6	78.46	21.18	0.00	0.36	0.00
7	82.64	17.15	0.00	0.21	0.00
8	77.84	21.26	0.00	0.90	0.00
9	88.60	10.40	0.00	1.00	0.00
10	76.65	22.37	0.19	0.78	0.00
11	70.48	27.24	0.19	2.10	0.00
12	83.62	16.19	0.00	0.19	0.00
13	81.95	17.13	0.00	0.92	0.00
14	76.04	21.80	0.90	1.26	0.00
15	78.87	19.25	0.00	1.89	0.00
Avg	78.74	20.13	0.17	0.96	0.00

Figure 4.

Figure of the Clarke error grid analysis for a select experiment.

In addition to the MAE and RMSE performance indices and Clarke error grid analysis, the low risk of missing hypoglycemia (below 70 mg/dL) or hyperglycemia (above 200 mg/dL) is a key feature of an efficient and reliable identified glucose-insulin dynamic model. We have also calculated the results for the percentage of the cases which the hypoglycemia or hyperglycemia are missed. Nine out of 15 experiments have CGM values below 70 mg/dl, and the average and standard deviation of missed hypoglycemia for these 9 experiments are 24.09% ± 16.50%. Thirteen out of 15 experiments have CGM values above 200 mg/dl, and the average and standard deviation of missed hyperglycemia for these 13 experiments are 07.83% ± 05.29%.

Discussion

A physiological model derived from the insulin compartment subsystem of Hovorka’s model is incorporated in this work with adaptive data-driven models developed using a recursive subspace identification technique. A computationally efficient nonlinear observer is used to estimate the time-varying parameters and the meal effect, thus facilitating the explicit incorporation of unannounced meals. The results for the proposed integrated compartment model with adaptive subspace identification techniques are promising. The proposed adaptive models may be used in the design of a model-based predictive control algorithm.

There are several benefits to the integration of compartment models with adaptive models. First, the PIC estimates are adaptive and individualized to particular patients, thus providing accurate real-time estimates of the amount of active insulin present in the bloodstream. The plasma insulin information can be readily used to impose constraints on the insulin dosing computation algorithm. Second, in contrast to the discrete basal and acute bolus insulin variations, the estimated PIC from the compartment model readily provides a more appropriate and filtered input variable for model identification. Third, the adaptive subspace identification technique renders the recursively identified glycemic models valid over a diverse range of daily activities without requiring onerous and obscure information such as the amount of carbohydrate consumption or quantification of physical activity levels.

Although meal and exercise announcements may improve the prediction ability of the proposed approach, unannounced meals and exercise are considered in this work to refrain from the inconvenience of user inputs and avoid the potential hazards of erroneous entries. The proposed adaptive modeling approach with an integrated compartment model utilizes additional auxiliary input variables such as the energy expenditure bio-signals to consider the effects of physical activity on glycemic dynamics.

Any insulin model capable to translate the abrupt bolus and discrete basal changes into estimates of PIC can be used. In this work, insulin subsystem of Hovorka’s model has been used since the PIC estimator is designed based on the Hovorka’s glucose-insulin dynamic model and we have developed personalization techniques based on the PLS models to individualize the parameters of this insulin subsystem for each specific patient.

PIC measurements are available from the clinical experiments, however, this information is not used in the development of the PIC estimation model and the on-line and real-time prediction of future CGM values. We use the clinical measurements of PIC to validate the estimates of the PIC obtained through the nonlinear simultaneous state and parameter estimation technique. As the PIC estimates are accurate, we employ the PIC estimates in the identification of adaptive state-space models. This ensures that only CGM measurements, infused insulin information, and physiological measurements from an nonintrusive wearable activity tracker are used in the on-line modeling and predictions of future glucose values.

For any data-driven modeling approaches where the model is obtained based only on input and output data measured form the underlying system, the quality and accuracy of the data play an important role on the performance of identified model. The time and frequency of calibrations have a significant effect on accuracy of the CGM measurements and consequently on the predictions accuracy. Subjects calibrated the CGM at least twice each day according to manufacturer instructions. The participants used their personal glucometers, and the CGM calibration may have occurred in nonideal glucose trends. CGM accuracy is greatest when the CGM is calibrated during a stable glucose period (flat arrows on CGM), but this was not always the case when the subjects chose to calibrate.

The steady-state value of the gut appearance rate (U_G) parameter should be zero during fasting time, provided all model parameters are personalized for each subject and the basal insulin rate stabilizes the BGC around a steady-state value. The purpose of considering the U_G in the PIC estimator is to handle any unknown disturbances present in the system, including the effects of unannounced meals (such as main meals or fast-acting carbohydrates for treatment of hypoglycemia) or any other uncertainties in the model parameters (such as EGP₀, $F_{01}^c$ and F_R parameters). To design the PIC estimator, other than the two model parameters ( $t_{\max ,I}$ and $k_{\mathrm{e}}$ ) that are individualized using the PLS models, the other parameters of Hovorka’s model are not personalized because personalization of all parameters requires sufficient data (more detailed information) on the subject’s glucose-insulin dynamics. Since the steady-state value of U_G is a function of several factors including the parameters of Hovorka’s model, the steady state value of the U_G can be a nonzero value even during fasting period.

The most challenging periods for the CGM predictions are the times that exercise and meal affect the glucose-insulin dynamics. Since, any information about these unknown disturbances such as their time, intensity, amount and duration are not known in advance, the accuracy of multistep-ahead predictions that is just based on current and past data decreases during these periods until the model parameters are updated based on measurements in those specific periods. One of the reasons that a recursively updated model is necessary for a fully automated artificial pancreas system (unannounced meal and exercise) is to update the model parameters at each sampling time using new measurements to adapt the model for that specific period. As shown in Figure 3, including the MET values and meal effects estimates improves the performance of the identified model for predicting the CGM values. In addition to analyzing the global error metrics, investigating the local error variabilities around the unannounced meals and exercises can be useful.

To design a fully automated artificial pancreas system that can operate without any user input, in addition to CGM data, information that could help the artificial pancreas control system to detect exercise and act properly for that period is needed. This information can be physical activity signals measured by wearable devices. To achieve this, different bio-signals have been used in artificial pancreas systems. The metabolic equivalent task (MET), or simply metabolic equivalent, is a physiological measure expressing the energy cost of physical activity (PA) and is defined as the ratio of metabolic rate (and therefore the rate of energy consumption) during a specific PA to a reference metabolic rate. MET is used as a means of expressing the intensity and energy expenditure of activities in a way that is comparable among people. We have developed a real-time MET estimation algorithm by using noninvasive measurements of physiological variables by our research group at Illinois Institute of Technology. In this approach, these MET values are derived using heart rate, galvanic skin response, skin temperature, blood volume pulse, and accelerometer.

The performance of proposed technique has been verified using different indices based on one-step-ahead predictions and six-step-ahead predictions. Different indices are considered to show the ability of the identified model in tracking the real measurements as well as predicting the future trend in CGM measurements. As we extend the prediction horizon, the accuracy of the predicted values decreases (as it can be seen from 1 step to 6 steps ahead of prediction) due to several reasons including unknown phenomena like exercise, stress or meal that may take place. So choosing the prediction horizon is a trade-off between knowing long-term information about the future measurements and their accuracy.

For few clinical experiments (e.g. No. 14 in Table 3) the relative high value of performance indices (MAE and RMSE values) is observed in comparison to the average values. The proposed technique is a data-driven technique and the quality of the measured data (eg, CGM and physiological variables that may be corrupted with noise and outliers) have a significant effect of the performance, as well as the fact that any information about the meal and exercises are not used. Predicting such abnormal measurements is very challenging and updating the glucose-insulin dynamic models using poor quality of data may result in bad performance in predicting the future trend of measurements.

One of the limitations of our clinical experiments is the exercise times are close to the meal times. Although people with type 1 diabetes usually consume carbohydrates before, during, and after exercise to avoid hypoglycemia and the proposed approach demonstrates that the effects of meals and exercise have well been captured by the modeling technique, in our future clinical work, the clinical experiments will be conducted to test the algorithm under different conditions such as. no meals very close to the exercise time.

Conclusion

An adaptive system identification algorithm is proposed to determine reliable dynamic glycemic models for use in AP systems. The proposed modifications to the recursive system identification algorithm assist in handling uncertain disturbances such as unannounced meals and exercise, stochastic measurement noise, and other adverse phenomena. The efficacy of the proposed approach is demonstrated using case studies involving the modeling of time-varying glucose-insulin dynamics in clinical experiments, and the results demonstrate the approach is able to identify reliable, time-varying, and personalized glycemic models.