Predictive Control of the Blood Glucose Level in Type I Diabetic Patient Using Delay Differential Equation Wang Model.

Because of increasing risk of diabetes, the measurement along with control of blood sugar has been of great importance in recent decades. In type I diabetes, because of the lack of insulin secretion, the cells cannot absorb glucose leading to low level of glucose. To control blood glucose (BG), the insulin must be injected to the body. This paper proposes a method for BG level regulation in type I diabetes. The control strategy is based on nonlinear model predictive control. The aim of the proposed controller optimized with genetics algorithms is to measure BG level each time and predict it for the next time interval. This merit causes a less amount of control effort, which is the rate of insulin delivered to the patient body. Consequently, this method can decrease the risk of hypoglycemia, a lethal phenomenon in regulating BG level in diabetes caused by a low BG level. Two delay differential equation models, namely Wang model and Enhanced Wang model, are applied as controller model and plant, respectively. The simulation results exhibit an acceptable performance of the proposed controller in meal disturbance rejection and robustness against parameter changes. As a result, if the nutrition of the person decreases instantly, the hypoglycemia will not happen. Furthermore, comparing this method with other works, it was shown that the new method outperforms previous studies.

Introduction

Diabetes is one of the most common endocrine diseases in which insulin secretion is not enough to regulate blood glucose (BG) due to destruction of pancreatic β cells. On the other hand, in case of high BG, glucagon secretion also stops, and thus, BG level exceeds the normal range of 80–140 mg/dl.[1]

The most important goal in the treatment of diabetes is to maintain BG in the normal range. In fact, as depicted in Figure 1, the main objective is to find the optimal control signal for insulin injection rate. According to block diagram in Figure 1, the rate of insulin injections is applied by the pump in diabetic patients as a control signal.

Figure 1

Closed loop controller for a diabetic individual[12]

Some efforts to capture the glucose–insulin mechanism have led to the formulation of various glucose insulin kinetic models. So far, several models have been suggested to predict the dynamic behavior of glucose–insulin system.[2345]

Negative feedback ordinary differential equation (ODE) model of Sturis and Tolic (2000) and delay differential equation (DDE) models of Engelborghs (2001), Bennett and Gourley (2004), and Kuang, Li and Mason (2006), together with the Wang and Li (2007) and the extended Wang (2009) models are among the current and valid ones based on ODE and DDE.[236789]

Primitive models for diabetes could not model time delay from the moment BG level increases till insulin secretion time. In some models to get insulin secretion fluctuations, insulin is divided into two components of plasma and intracellular – this considered a disadvantage for the proposed model. However, in this study, the delayed nonlinear model of Wang and Li is used considering the nonlinear behavior of insulin–glucose interaction for type I diabetic patient.

Closed loop control methods to regulate BG in type I diabetic patients are mainly based on the model and also the experimental data. Methods based on the experimental data identify model parameters using glucose–insulin data and results are then obtained by applying proportional–integral–derivative (PID) controllers.[10] In addition, various control methods including H∞,[1112] adaptive,[13] proportional–derivative (PD),[14] PID,[15] fuzzy,[1617] etc. have been proposed to regulate BG levels in type I diabetic patients. These methods differ in terms of employing control strategies, use of constraints, mathematical models, and ease of implementation, each with its own advantages and disadvantages.

In 2008, Gianni Marchetti proposed an improved PID control strategy for BG control and critically evaluated in silico using a physiologic model of Hovorka.[18] An artificial pancreas strategy using constrained model predictive control is developed to achieve closed-loop glucose control for type I diabetes in 2009. A system of meal detection and meal size estimation is also developed to automatically administer meal insulin boluses as feed-forward action to unmeasured meals.[19] In 2010, a model-based predictive control scheme was applied to a newly developed diabetic patient model. The controller was provided with a feed-forward loop to improve meal compensation, a gain-scheduling scheme to account for different BG levels, and an asymmetric cost function to reduce hypoglycemic risk.[20] In 2011, a system based on a nonlinear model-predictive controller was developed which used a personalized glucose–insulin metabolism model, consisting of two compartmental models and a recurrent neural network. The model took as input the patient’s information regarding meal intake, glucose measurements, and insulin infusion rates, and provided glucose predictions.[21] A novel automatic adaptive control strategy based on frequent glucose measurements and a self-tuning control technique was validated based on a simulation study for 200 virtual patients in 2013. The adaptive control strategy was shown to be highly effective in controlling BG concentration.[22] Control methods based on nonlinear models have been introduced which employ physiological behavior of patients to provide different control methods in regulating BG. It is obvious that if the model behavior is more similar to the patient’s body, the resulting control law is more accurate too. Hence, providing a suitable approach to control delayed nonlinear models for diabetic patients is of great importance.

In this study, the proposed control system is capable of predicting BG levels and injecting insulin so that the BG level always lies within normal range, showing suitable performance against meal disturbances and uncertainties in the model.

Delay Differential Equation Model

These models range from simple expressions relating glucose and insulin, to very complicated mathematical models. To simulate the glucose dynamics of the patient’s body, two well-known mathematical models are considered. The first model, proposed by Wang et al. in 2007[2] for glucose–insulin interactions in the body, consists of 2 DDEs describing various sections in the body. They performed both qualitative and quantitative studies of the dynamics of the model. The analytical results showed the existence and uniqueness of a stable periodic solution corresponding to ultradian insulin secretion oscillations. Numerical simulation results of insulin administration based on their model matched with the findings of the clinical studies.

The second model, known as the Enhanced Wang model,[3] consists of 2 DDEs too. The Enhanced Wang model is a nonlinear compartmental model for insulin therapy for both type I and type II diabetes mellitus, in which the insulin degradation rate assumes Michaelis–Menten kinetics.

It is well known that Michaelis–Menten kinetics is suitable for the response function in chemical reaction, when the reaction rate does not increase indefinitely when an excess of resource is available.

However, the existing models for insulin therapies take it for granted that the response function of insulin clearance is proportional to the insulin concentration. Their analysis shows that it is possible to simulate pancreatic insulin secretion by exogenous insulin infusions, and their numerical simulations provide clinical strategies for insulin–administration practices.

The Wang and the Enhanced Wang models have similar dynamical behaviors.

Equations related to Wang and Enhanced Wang models are shown in (1) and (2), respectively:

α > 0, β ∈ [0,1], for type I diabetes: β = 0. (no insulin is secreted from the pancreas)

Here, the uptake of food glucose, modeled by (3),[2] is denoted by Gin:

In these models, G(t) is the output and Gin and Iin(t) are the system inputs. Insulin reduction in body differs from person to person. In Wang model, insulin reduction is defined as d and for Enhanced Wang model with Michaelis–Menten kinetics:

where d1 is the maximum insulin clearance rate, and d2 is the half-saturation value.

Functions of the models are as follows:

f1(G) = R/(1 + exp((C1 − G/V)/a))

f2(G) = U/(1 − exp(−G/(C2V)))

f3(G) = G/(C3 V)

f3(G) = U + (U − U)/(1 + exp(−βln(I/C4 + 1/(Et))))

f5(I) = R/(1 + exp(α(I/V − C5))      (4)

f1(G): Glucose-dependent insulin secretion.

f2(G): Insulin-independent glucose consumption by the brain and nerve cells.

f3(G)f4(I): Glucose-dependent insulin consumption by muscle cells and fat.

f5(I): Glucose production controlled by insulin concentration.

Parameters related to functions f2 to f5 are mentioned in Table 1.

Table 1

Parameters of the functions in (4)

Designing Nonlinear Predictive Controller for Insulin Injection System

The method presented in this study is based on predictive nonlinear controller optimized with genetic algorithm.

Figure 2 depicts the basic principle of model predictive control. On the basis of the measurements obtained at time t, the controller forecasts the future.[23] Dynamic manner of the system over a prediction horizon T determines the input such that a predetermined performance objective function is optimized. If there were no model-plant discrepancy, and if the optimization problem were solvable for infinite horizons, then the input function at time t = 0 to the system is applicable for all times t ≥ 0. The resulting manipulated input function will be implemented only until the next measurement is available. The time gap between recalculations can vary; however, it is supposed to be fixed − the measurement will occur every δ sampling time-units. Using the new measurement at time t + δ, the entire process-prediction and optimization will be repeated to find a new input function with the control and prediction horizons moving forward. As shown in Figure 2, the input u is denoted as an arbitrary function of time. The calculation of the applied input based on the predicted system behavior allows the inclusion of constraints on states and inputs as well as the optimization of a given cost function. The stabilization problem for a class of systems is introduced by the following nonlinear set of differential (5):

Figure 2

Principle of model predictive control[18]

X(T) = F(X(T), U(T), X(0)) = X0      (5)

which is subject to input and state constraints of the form:[23]

where x(t) and u(t) represent the states and inputs vector, respectively. U and X constraints are given in (7), where umin, umax and xmin, xmax are constant vectors:[23]

To distinguish between the real system and the system model used to predict the future within the controller, the internal variables in the controller are denoted by a bar (e.g., x̄, ū). The finite horizon optimal control problem commonly described above is mathematically formulated as following:[23]

min J (x(t), u(t); T, T)      (8)

The function F, hereinafter called cost function, specifies the favorable control performance that can arise. The standard quadratic form is the simplest and most often the used one:

F(x, u) = (x − x)T Q(x − x + (u − u)T R(u − u)      (10)

where u and x are the desired input and output, respectively. Q and R are symmetric, definite positive, and weighing matrices. T is the horizon of the predicted output, and T is the control horizon. Eq. (10) gives the error between desired output and model-predicted output. To obtain the output values in the next time interval, the optimal input values in time period T are used, and then, the value of the control variable is set constant in the last calculated value.[23]

The block diagram of nonlinear model predictive control optimized with genetic algorithm is given in Figure 3, composed of a system model and plant.

Figure 3

Block diagram of nonlinear model predictive control optimized with genetic algorithm

The controller is designed based on the data provided by the Wang model. To show the effectiveness and robustness of the designed controller, it is used to regulate BG predicted by both Wang and the Enhanced Wang models.

To design predictive controller for the model, an objective function needs to be designed. Real-time optimization of objective function will lead to design of a control signal that is able to track the suitable reference path by predicting system behavior. According to the studies, the optimal amount of BG is 110 mg/dl; nevertheless, the 80-to-140 interval is also known as green or healthy zone. Therefore, this objective function is proposed:

In (11), Ĝ is the anticipated blood sugar, G is the optimal blood sugar level (110), Gup and Gdown are high (140) and low (80) blood sugar values, respectively, t is present, u(t + j) is future control signal, W(Z−1) is defined to solve a single problem of dynamic matrix that is equal to 1−Z−1, N is the prediction horizon, and M is the control horizon. Most importantly, Q1,2,3 values have to be regulated so that blood sugar levels do not lie in unhealthy conditions.

As a result, parameter tuning is as follows:

Q1 = 1

Q2 = {0 ifĜ < Gup1000ifĜ > Gup

Q3 = {1000 ifĜ < Gdown0if Ĝ > Gdown.      (12)

The R-value in (11), denoting objective function penalty for taking too much insulin, is taken as 100 in simulations. The prediction horizon is set as 45 min.[1]

As the objective function is optimized with genetic algorithm for genetic algorithm parameters, maximum number of repetitions is 30, initial population size is 25, crossover value is 50%, and mutation rate is considered as 40%. In addition, the Rollet Wheel Selection is considered as the method for selecting members of recombination.

Simulations and Results

In this section, the simulation results in the absence of controller (open loop system) are examined, the need to use the controller is checked, and the results of adding controller to the system are discussed. The proposed controller is then tested against disturbances and uncertainties in the system and the obtained results will be finally compared with other studies, all presented in Tables 2 and 3.

Table 2

Numerical simulation results

Table 3

Comparison of the daily infused insulin under the genetics optimal nonlinear model-predictive controller (proposed approach) and that of Lee and Bequette[19]

The proposed simulation method was implemented in MATLAB software.

Open loop response system

Figure 4 depicts the amount of sugar entering the body in standard conditions. The profile provided in this diagram is a standard one defined based on (3).

Figure 4

Glucose intake rate (Gin) in (3)

In the absence of insulin injection, blood sugar is 142 mg/dl and insulin level is 18 U. Please note that the initial system was deliberately put in unhealthy conditions. Results are presented in Figure 5. As seen, blood sugar level always remains in unhealthy situation. In this figure, the amount of injected insulin is zero, but in Figure 6, we set the rate of injection at a nonzero level. Because of lack of a feedback loop in injection system, the amount of BG is decreased until it becomes as low as dangerous situation.

Figure 5

Open loop response system. The rate of insulin (I) and glucose (G) change in 24 h

Figure 6

Open loop response system with constant injection rate. The rate of insulin (I) and glucose (G) change in 24 h

Closed loop response system

As depicted in Figure 7, the amount of BG level is immediately shifted from unauthorized to authorized level when the controller is used.

Figure 7

Closed loop response system. Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

Disturbances and uncertainties in the system

To evaluate the system performance in food disturbances mode, two states of impact noise and increasing blood sugar input were studied. In the first case, the patient is assumed to consume sugar at a time other than the time period defined in sugar intake profile. This process exhibits some changes during the time period 8–12 [Figure 8]. As seen, the controller has an appropriate performance in this case and does not permit above-limit deviation. However, an excessive amount of 1.6% insulin is injected here to compensate for the imposed disturbances.

Figure 8

System behavior for impact noise mode. Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

In this study, another disturbance was studied in which one of the peaks of sugar intake by patient experienced a 20-fold jump [Figure 9]. In this case, to use the controller, an increase of 1.8% in insulin injection is implemented every 45 min.

Figure 9

System behavior while imposing noise in one of the peaks of sugar intake. Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

Insulin reduction parameter varies for different people in Wang model. Therefore, to compare results for various modes, d received 20% increase and also 20% reduction. Results are depicted in Figures 10 and 11. In more critical situations involving d increase, insulin injection increases by 8–10%. However, in d decline mode, injected insulin drops by 20% as shown in Table 2.

Figure 10

System behavior in parametric uncertainty for d increase mode (d (1 + 20%)). Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

Figure 11

System behavior in parametric uncertainty for d reduction mode (d (1−20%)). Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

Another parameter that may vary in different patients is time delay of insulin-dependent glucose consumption τ2 and also time delay related to insulin-dependent glucose secretion τ3. Any time delay of more than 15 min for glucose secretion caused by insulin incline will put the patient at risk for a severe drop in blood sugar. Figure 12 represents system simulation in the event that time delay is 50 min (according to Wang et al.[2]).

Figure 12

System behavior with delay time of 50 instead of 15 min. Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

As seen in above chart, although some drop in blood sugar occurs, the controller is still able to return BG level to permissible range by 10% less insulin infusion. Result of increasing time delay for glucose consumption from 5 to 30 min is shown in Figure 13.

Figure 13

System behavior with delay time of 30 instead of 5 min. Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

Another case studied in this research is the incompatibility between the model used in controller structure for prediction and the real model of glucose in the patient. To do so, the Wang model was used in controller structure, and the incremental Wang model, as an extension to Wang model and suitable for both diabetes types 1 and 2, was used in system structure. In this model, insulin consumption rate is considered more realistically as a function of insulin level in blood. The results of this study verify the appropriate performance of our system even in the absence of agreement between controller structure model and system dynamics model. Simulation results are seen in Figure 14.

Figure 14

Evaluating model uncertainty and the mismatch between controller structure model and system model. Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

Now, we examine the case in which model uncertainty and the noise resulting from impact and increased glucose consumption occur simultaneously. As shown in Figure 15, despite the fact that BG level violates the allowable range in noisy moments, the controller can well bring it back to the range within permissible limit.

Figure 15

Evaluation of model uncertainty along with noise. Glucose intake rate (Gin), insulin injection (Iin), blood glucose level (G), and blood insulin level (I) in 24 h

Comparison of the results

On the basis of Table 3, our proposed controller outperforms in regulating BG levels and also reducing daily insulin dosage compared with other controllers such as fuzzy PD, fuzzy proportional–integral (PI), genetics optimal fuzzy PI, and genetics optimal fuzzy PID that were proposed in Al-Fandi et al.[24] with same glucose–insulin model. Fortunately, Ref.[24] is one of the previous studies in this area that reports the amount of insulin injection per day. In other references, like,[1232526] the main goal is only to protect the BG at safe level. We can refer to inattention to optimization of insulin injection per day as the main drawback of latest studies. In this research, we optimize the amount of injection in addition to maintaining the glucose at the safe level using a nonlinear model predictive control system.

Conclusions

As evident, this study was aimed at providing a predictive control method for improving the performance of a system for automatic injection of insulin to diabetic patients. In this regard, after observing the unstable response of open-loop system (untimely insulin injection), the proposed method was implemented. In this approach, the dynamic model of blood sugar and insulin variations is seen within the structure of the controller, and the controller can predict variations of blood sugar and insulin level using current measurements.

Following the prediction done by the controller, the optimum insulin injection is computed on a real-time basis using genetic algorithm so that the unhealthy blood sugar in a patient is maintained in its lowest possible value in 24 h. To evaluate the performance of the designed controller, various scenarios ranging from normal and noisy conditions with unpredicted factors, to parametric uncertainty, and finally model uncertainty were designed and implemented. The results showed the ability of the controller system to regulate blood sugar levels, ensuring the accuracy of its performance in different conditions. The results of this research can well compete with works done by other researchers.

Financial support and sponsorship

Nil.

Conflicts of interest

There are no conflicts of interest.