|Home | About | Journals | Submit | Contact Us | Français|
Control algorithms that regulate blood glucose (BG) levels in individuals with type 1 diabetes mellitus face several fundamental challenges. Two of these are the asymmetric risk of clinical complications associated with low and high glucose levels and the irreversibility of insulin action when using only insulin. Both of these nonlinearities force a controller to be more conservative when uncertainties are high. We developed a novel extended model predictive controller (EMPC) that explicitly addresses these two challenges.
Our extensions to model predictive control (MPC) operate in three ways. First, they explicitly minimize the combined risk of hypoglycemia and hyperglycemia. Second, they integrate the effect of prediction uncertainties into the risk. Third, they understand that future control actions will vary if measurements fall above or below predictions. Using the University of Virginia/Padova Simulator, we compared our novel controller (EMPC) against optimized versions of a proportional-integral-derivative (PID) controller, a traditional MPC, and a basal/bolus (BB) controller, as well as against published results of an independent MPC (IMPC). The BB controller was optimized retrospectively to serve as a bound on the possible performance.
We tuned each controller, where possible, to minimize a published blood glucose risk index (BGRI). The simulated controllers (PID/MPC/EMPC/BB) provided BGRI values of 2.99/3.05/2.51/1.27 as compared to the published IMPC BGRI value of 4.10. These correspond to 73/79/84/92% of BG values lying in the euglycemic range (70–180 mg/dl), respectively, with mean BG levels of 151/156/147/140 mg/dl.
The EMPC strategy extends MPC to explicitly address the issues of asymmetric glycemic risk and irreversible insulin action using estimated prediction uncertainties and an explicit risk function. This controller reduces the avoidable BGRI by 56% (p < .05) relative to a published MPC algorithm studied on a similar population.
Patients with longstanding type 1 diabetes mellitus have little or no endogenous insulin production, leaving the body unable to lower blood glucose (BG) levels without exogenous insulin. Technological advances have enabled the development of a closed-loop artificial pancreas, reducing the burden of glucose management for these patients.1–3 Continuous glucose monitors measure subcutaneous glucose levels while pumps infuse insulin and/or glucagon subcutaneously. A controller then connects these two. While some investigators have proposed and tested controllers using pumps for both insulin and glucagon,4 we focus on the more common approach using only an insulin pump.
In this study, we examined two fundamental challenges faced by any such closed-loop artificial pancreas and developed a novel controller to explicitly address these difficulties.
First, the acute and chronic risks of extreme blood glucose levels are asymmetric.5 Low BG levels (hypoglycemia) are acutely risky as they can result in altered mental state, seizures, and coma. Meanwhile, high BG levels (hyperglycemia) increase the risk of chronic complications such as retinopathy, nephropathy, and cardiovascular disease. In the face of significant uncertainty, this asymmetry should cause a controller to err toward high BG levels. That is, controller behavior should vary with observed uncertainty.
Second, insulin action is irreversible. Using only insulin, the controller cannot actively counter the effects of delivered insulin. This means that a controller should be cautious to avoid slow recoveries from unexpected insulin overdoses. And again, the level of caution should vary with the level of uncertainty.
In this work, we extended model predictive control6 as a basis for a closed-loop artificial pancreas. We explicitly managed the expected future risk to trade off the dangers of hypoglycemia and hyperglycemia. In the process, we explicitly considered uncertainty and the potential effects of future BG measurements and controller actions.
We began by analyzing the challenges of asymmetry and irreversible insulin. We then developed the novel extended model predictive controller (EMPC) that explicitly deals with these challenges. We finally demonstrated the improved performance of this controller on the University of Virginia/Padova (UVa/Padova) Metabolic Simulator7 against a proportional-integral-derivative (PID) controller, a basic model predictive controller (MPC), and a published MPC algorithm.8
The ultimate objective of any BG controller is to minimize the complications resulting from poor BG control. To this end, we proposed incorporating an explicit estimate of risk into our feedback controller. We first discussed known approximations of the risk of diabetic complications versus BG level. Thereafter, we developed a continuous, convex function describing risk that is also suitable to the optimization needed in MPC. We then explored how asymmetric risk should affect control decisions.
The Diabetes Complication and Control Trial (DCCT)9 and the Epidemiology of Diabetes Interventions and Complications (EDIC) trial10 show an increased risk of chronic complications with higher hemoglobin A1c (HbA1c) (glucose levels). They also describe the acute detrimental effects of hypoglycemia. While informative, they do not provide a single numerical metric of the risk of complications.
To provide a single risk metric combining chronic and acute elements, Kovatchev and colleagues5,11 symmetrized the distribution of BG levels. They assumed that less common BG levels are inherently more risky than more common BG levels. The resulting formula is R(g) = 10[1.509ln(g)1.084 – 5.381]2, where g is the glucose concentration in mg/dl. Meanwhile, other research8 proposed a similar function that penalizes low blood glucose levels more: R’(g) = 10[3.5506ln(g).8353 – 3.7932]2. This second metric was “modified with respect to literature values to better suit control performance results.”8
To enable comparisons with published results, we adopted the second approximation in our validations. Both approximations are shown in Figure 1. While these works provide a nice estimate of clinical risk, neither estimate is well suited to estimation of clinical risk as part of a controller, because neither handles negative values and so are not continuous or convex over the set of real numbers. Handling negative values is critical since many proposed closed-loop controllers use simple models of the glucose system that predict negative glucose values when simulating large insulin doses.8,12–15 Continuity and convexity are also important since they make optimizations faster and more robust.
While there are ways to modify the existing risk functions to make them continuous and convex over the real numbers, we chose to generate a new simpler function based directly upon the DCCT and EDIC trials.
The chronic risk of retinopathy correlates roughly linearly with the HbA1c for values between 5.5 and 9.0 (about 111 to 212 mg/dl).9 Because HbA1c correlates with mean BG values,16 we could approximate the risks of chronic complications as a linear function of BG levels.
The acute risks of low BG levels increases sharply with dropping BG level starting with lethargy and mild hypoglycemia and progressing quickly to seizures and sever hypoglycemia. As there is no data to directly support a specific analytic expression, we chose a simple cubic model.
The risk function for control that combines hypoglycemic and hyperglycemic risks is then
and is parameterized by the constants a, b, c, and d. In practice, this expression could be adjusted parametrically or even functionally subject to clinical judgement or at a user’s discretion.
To match the R’(g) function, discussed above, used in Magni and colleagues8 we chose the parameter values
Figure 1 graphs the risk function developed here together with two other published graphs.
The risk of complications is asymmetric about the lowest risk glucose concentration of 140 mg/dl,8 so that a glucose concentration of 190 mg/dl is much less risky than a glucose concentration of 90 mg/dl. Since future BG concentration values are never certain, controllers should always consider the possibility of values that become higher or lower than expected. This means that as the prediction uncertainty increases, less insulin should be delivered to bias BG levels to higher values where risk increases more gently. Figure 2 shows an illustrative example where shifting the range of possible future values upward lowers the average risk. Mathematically, this corresponds to minimizing the expected risk of the prediction distribution instead of minimizing the risk of the mean prediction. This is a significant effect since predictions for several hours into the future can easily have standard deviations of 40 to 100 mg/dl.12
Predictions are inherently uncertain, hence controllers should consider how they will respond when the actual values end up above or below the predicted mean values. With the aid of future glucose measurements, higher than predicted mean values will be rejected by injecting more insulin. Lower than predicted values may be rejected by injecting less insulin provided that injecting less insulin does not mean removing already injected insulin. In that case, the controller can only rely on the endogenous glucose production (EGP) to slowly raise the glucose concentrations back to desirable levels, as shown in Figure 3. Therefore, it is imperative that controllers understand how uncertain predictions are and not inject too much insulin too soon.
This section extends MPC in three steps to explicitly incorporate the risk of complications, prediction uncertainty, and future glucose measurements.
MPC chooses insulin infusions to optimize a cost function.6 We denote u0 as the insulin to be infused at the current time, a vector of future insulin infusions, and a vector of BG predictions based on the current and future insulin. The optimization can be written mathematically as
where all infusions must fall below the pump limit Umax.
Traditionally, cost functions are often quadratic and trade off glucose excursions with infusion rate changes. For example,
where J provides a weighting parameter for tuning controller aggressiveness.
To minimize the risk of complications, we replaced the traditional cost function with the average risk
for all predicted glucose concentrations. We generated the predictions using a published prediction algorithm that incorporates meal detection and estimation.17 The adjusted optimization is
This change allowed us to minimize the risk but ignores both prediction uncertainty and future glucose measurements.
In addition to the predicted future glucose concentrations , the prediction algorithm also supplies the estimated standard deviations of the predicted concentrations. Assuming that the prediction errors are Gaussian, and fully describe the distribution of possible future glucose values.
To approximate the expected risk of the entire distribution of predictions, we calculated the average risk for five parallel predictions that span the full distribution. These five trajectories are , where j [–2.5, –1.25, 0, 1.25, 2.5]. The average risk values for the five trajectories are summed according to their normalized probability Pj. We omitted the constants and standard deviation terms in the calculation of Pj as they either cancel out or do not affect the optimization. The resultant framework is
We chose to approximate the full distribution of predictions by five potential trajectories as this provides minimal error for a reasonable computational cost.
As we were now considering five potential future glucose trajectories, we realized that the controller would have to respond differently in the future based on which of these trajectories will come true. And while future measurements will resolve this multiplicity, until such measurements arrived we had to allow a multiplicity of future control actions. So we defined five parallel future control actions , each associated with a potential glucose trajectory . The final form of the optimization problem solved in the EMPC is thus
We validated the performance of the EMPC on the UVa/Padova Metabolic Simulator.7 The simulator is approved by the U.S. Food and Drug Administration as a substitute for animal trials using a population of 100 patients to model meals, glucose-insulin dynamics, and interpatient variability. The publicly available portion of the simulator includes only 10 adult patients who were chosen to span the full range of the 100 patients. Unfortunately, for such a small a number of patients this implies that these patients are not a typical distribution of patients. In particular, we believe patient 9 encodes an outlier case and is not representative of one tenth of patients. So while we show individual results for all 10 publicly available patients, we exclude patient 9 when averaging over the patient population. This is also detailed in Appendix A.
Each patient was simulated for 36 hours and was provided with six unannounced meals, lasting 20 minutes each and measuring 50 g CHO at 9 a.m., 70 g at 1 p.m., 90 g at 5:30 p.m., 25 g at 8 p.m., 50 g at 9 a.m., and 70 g at 1 p.m. This scenario matches a scenario used for a clinical trial, and further represents a worst-case scenario where the patient announces none of their meals. Research shows that adolescents routinely miss two meal boluses per week while using manual control.18 We feared that the use of closed-loop control will worsen this behavior, making the worst-case scenario an important test.
The controllers used a sample time of 5 minutes, calculating a new insulin infusion every 5 minutes. Each infusion was limited to a Umax of six units and predictions were carried out to a horizon of 60 samples or 300 minutes.
We provided three controllers for direct comparison to the EMPC results.
First, we created a basic PID controller that is fully described by
The setpoint was chosen to match the lowest risk glucose value, 140 mg/dl. The integral term was locked to the basal rate to avoid windup issues associated with saturated inputs and sustained positive disturbances, such as meals. The weights for the proportional (6.375 × 10-5 units per mg/dl) and derivative (0.0046 U per mg/dl) terms were optimized to minimize the average blood glucose risk index (BGRI)8 on the 10 patients.
Second, we created an MPC using the cost function in Equation (6), where J = 300 was chosen to minimize the average BGRI on the 10 patients.
This MPC shares the same prediction algorithm as EMPC but differs in the cost function. Thus, comparing the results of EMPC and MPC controllers should serve to illuminate the benefits of the contributions of this article, isolated from the potential confounding effects of different prediction methods.
Lastly, we used an optimized BB control as a lower bound on the BGRI.8 Specifically, the basal rate and one bolus at the start of each meal were optimized to provide the minimum possible BGRI.8 The basal rate was chosen to cause a steady state value lower than the zero risk value of 140 mg/dl. This reduces the risk effect of meals, and suggests a slow recovery from the last meals in the day. Because so many parameters (basal rates and boluses) were tuned noncausally, on a patient-specific level and with full knowledge of meals, this served as a lower bound on the achievable BGRI.8
These controllers caused different characteristic effects on the controlled BG level. The PID controller reacted slowly to meals allowing larger disturbances due to meals. On the other end, the BB control reacted immediately. Since the EMPC controller was not hampered by the objective weight on large insulin infusion rates, it rejected meal disturbances faster than the MPC controller but slower than the noncausal BB controller.
The EMPC and MPC both show dips in the BG value just before meals as a result of meal anticipation in the prediction between 6 a.m. and 10 p.m. Because the EMPC also uses uncertainty, it avoided the prebreakfast low that occurs for the MPC after 3 hours of anticipated but unrealized meals. Figure 4 shows the average BG levels across all valid patients. Appendix B shows the individual results for all simulated patients.
Tables 1 and 22 provide the value of published performance measures.5 The BGRI, low blood glucose index (LBGI), and high blood glucose index (HBGI) follow Magni and colleagues’ definition used in the published results.8 In Table 1, we include the Student’s t-test significance measure of the changes in the BGRI relative to EMPC. Comparing the EMPC to the PID, MPC, and BB controllers, the t-test assesses whether the two sets of nine BGRI values for the nine valid patients are likely to share the same mean.
We also validated against the published results from an independently created MPC,8 denoted by IMPC. The IMPC was tested on the private set of 100 adult simulator patients. The IMPC controller is advantaged because it uses a linearized version of the simulator model, but disadvantaged because it changes the insulin infusion rate only once every 30 minutes compared to our 5-minute update rate. Consequently, a comparison of the results for the IMPC and EMPC results (Table 3) is informative but not authoritative.
Assuming that the Student’s t-test applies, there is a 98.34% chance that the EMPC provides a lower mean BGRI than the IMPC.
The EMPC controller has shown statistically significant performance benefits on the Metabolic Simulator, as measured by the BGRI. In comparison to the published MPC,8 which differs both in prediction and control algorithm, the EMPC reduced the avoidable risk by 56%. In comparison to our MPC, which shares the same prediction algorithm, the EMPC lowered the avoidable risk by 30%.
These improvements can be traced to two general sources: the prediction algorithm and the control algorithm. With regard to prediction, we know that reduced BGRI correlates with a decreased premeal glucose level. This suggests a benefit from controlling more aggressively in advance of expected meals to help reduce the future meal effect. Indeed, the prediction algorithm used by both the MPC and EMPC anticipates meals, allowing the controllers to implicitly accomplish this. Testing of the EMPC controller with a worst-case missed dinner and evening snack shows a minimal downside for this benefit. In these tests, the glucose concentrations only dropped to a minimum of 84 mg/dl after the missed meals. A concentration of 84 mg/dl has the same BGRI as the hyperglycemic glucose value of 244 mg/dl.
The improvement due to the controller extensions presented in this article can be attributed to the incorporation of uncertainty and thereby the ability to pursue prediction horizons of 5 hours. Other MPC algorithms use maximum prediction horizons between 2 and 4 hours.8,14,19,20 Without evaluating uncertainty, each of those algorithms is forced to treat all predictions with the horizon as equally relevant. Our inclusion of uncertainty can implicitly reduce the weight given to the less certain, long-term predictions. This also allows us to adapt to the increased uncertainty around unannounced meals with unknown sizes.
We also believe the presented MPC extensions will prove robust in practice. The EMPC contains no tuning parameters or indirect weighting parameters. It simply allows the definition and adjustment of the risk function, Equation (3), which has direct clinical relevance.
Of course, the EMPC does leave room for improvement. In particular, we approximated the distribution of potential future glucose values by five distinct glucose trajectories. Other approximations using more trajectories or trajectories that branch further in the future are conceivable. With increasing computational power, it may be feasible and even beneficial to consider greater degrees of multiplicity.
By taking advantage of estimates of prediction uncertainty and understanding the asymmetries associated with the risks of complications and irreversible insulin action, we extended the MPC framework. The resultant, novel MPC-based controller lowers the avoidable risk by 56% relative to a published MPC algorithm when tested on a similar population. We hope and believe that with additional in silico and in vivo testing, these extensions and the controllers and products they enable will help alleviate the burden for type 1 diabetes patients.
Funding from the Juvenile Diabetes Research Foundation Artificial Pancreas Project is gratefully acknowledged.
The UVa/Padova Simulator uses a population of 100 adults to function as an FDA-approved substitute for animal trials. Unfortunately, only 10 adult patients are publicly available. We examined these 10 adult patients to find that patient 9 falls well outside the norm and does not represent a normal patient. While such outliers may be useful in large populations to test worst-case behaviors, they do not represent one tenth of all patients. As such, we were forced to exclude patient 9 from the averaged results.
We first examined insulin action. We see that insulin has a very slow and powerful effect on the BG levels of patient 9. We simulated all 10 patients with a 50 g CHO meal and an optimal meal bolus, as defined by the simulator’s user’s guide. Figure A1 shows the resulting blood glucose values.
We also reviewed EGP suppression following a 50 g CHO meal in Figure A2. Insulin suppressed EGP by 46% more for patient 9 than for any other patient. This effect is compounded by the delayed nature of EGP suppression. Published meal shapes21–24 and insulin time action profiles25–27 are 90% finished after 6 and 6.33 hours or less, respectively. The cumulative effect of the EGP suppression for patient 9 is only 63% complete after 6.33 hours.
This slow, powerful suppression of EGP caused a significant drop in the BG values for patient 9 long after the EGP and meal bolus effects had faded for the other nine patients. While exercise can cause similar symptoms, we were not trying to test for those here.
The following table provides the results individually for all ten patients.