MPC-guided Imitation Learning of Neural Network Policies for the Artificial Pancreas

03/03/2020 ∙ by Hongkai Chen, et al. ∙ 13

Even though model predictive control (MPC) is currently the main algorithm for insulin control in the artificial pancreas (AP), it usually requires complex online optimizations, which are infeasible for resource-constrained medical devices. MPC also typically relies on state estimation, an error-prone process. In this paper, we introduce a novel approach to AP control that uses Imitation Learning to synthesize neural-network insulin policies from MPC-computed demonstrations. Such policies are computationally efficient and, by instrumenting MPC at training time with full state information, they can directly map measurements into optimal therapy decisions, thus bypassing state estimation. We apply Bayesian inference via Monte Carlo Dropout to learn policies, which allows us to quantify prediction uncertainty and thereby derive safer therapy decisions. We show that our control policies trained under a specific patient model readily generalize (in terms of model parameters and disturbance distributions) to patient cohorts, consistently outperforming traditional MPC with state estimation.

READ FULL TEXT VIEW PDF
POST COMMENT

Comments

There are no comments yet.

Authors

page 1

page 2

page 3

page 4

Code Repositories

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

The artificial pancreas (AP) is a system for the automated delivery of insulin therapy for Type 1 diabetes (T1D), a disease in which patients produce little or no insulin to regulate their blood glucose (BG) levels and maintain adequate glucose uptake in muscle and adipose tissue. The AP consists of an insulin infusion pump, and a subcutaneous Continuous Glucose Monitor (CGM) for sensing glucose levels. CGM readings are transmitted to a control algorithm that computes the appropriate insulin dose. Such control should maintain a fine balance. Lack of insulin leads to hyperglycemia (i.e., high BG), which, if untreated, can cause complications such as stroke, kidney failure, and blindness. Excessive insulin can lead to hypoglycemia (low BG), a critical state which can result in unconsciousness and death.

Driven by advances in the mathematical modeling of T1D physiology (Hovorka et al., 2004; Man et al., 2014), Model Predictive Control (MPC) (Camacho and Alba, 2013) has become the gold-standard AP algorithm. MPC works by determining the insulin therapy that optimizes the future BG profile, predicted via physiological models. MPC, however, has two important limitations. First, it requires complex (often nonlinear and nonconvex) online optimization, which is infeasible when the algorithm is deployed in resource-constrained medical devices. This is why commercial AP systems either use simplistic linear models within MPC111The algorithm of Typezero® (Forlenza et al., 2019) uses a linearized version of the nonlinear model of (Man et al., 2014). or do not employ MPC at all in favor of simpler control algorithms (e.g., PID control in the Medtronic™ Minimed 670G (Steil, 2013), or fuzzy logic in the Glucositter by DreaMed Diabetes (Atlas et al., 2010)).

Second, and most crucial, MPC requires state estimation (SE) to recover the most recent patient state from CGM measurements (Gondhalekar et al., 2014; Paoletti et al., 2019). Besides its computational cost, SE is error-prone, as it relies strictly on CGM readings, which are an approximate, delayed, and noisy proxy of the target control variable, the BG. Incorrect state estimates might compromise the safety of the insulin therapy.

Our Contributions.

We present a novel method to derive end-to-end insulin control policies, i.e., policies that subsume state estimation and control, directly mapping CGM measurements into patient-optimal insulin dosages. To capture the complex logic of MPC and SE, we consider policies represented as Deep Neural Networks, LSTMs in particular (Hochreiter and Schmidhuber, 1997). Such an approach addresses the main issues surrounding the use of MPC, as it bypasses explicit SE and avoids the cost of MPC’s online optimizations.

Our approach is centered around the use of Imitation Learning (IL) (Ross et al., 2011), where the control policy is trained on examples provided by an MPC-guided teacher policy. Building on the PLATO framework (Kahn et al., 2017), at training time, we instrument the MPC teacher to access the full state of the patient model. As such, the learner policy, using only CGM measurements, will learn to mimic MPC-based decisions, which are based on the true model state and thus highly accurate. In this way, the learned policy incorporates both SE and control.

We employ IL as it alleviates the covariate shift (Storkey, 2009) caused by the discrepancy between the training-time state distribution (induced by the teacher policy) and the test-time one (induced by the learner). Without IL, such a distribution shift can lead to unpredictable behaviour of the learner, and hence jeopardize the patient’s safety.

We learn our insulin policies via approximate Bayesian inference, using Monte Carlo dropout (Gal and Ghahramani, 2016). The resulting stochastic policy can capture model (epistemic) uncertainty, as well as uncertainty in the data, due to, e.g., noisy measurements (Gal, 2016). Uncertainty quantification is key in safety-critical domains like ours, and is not addressed in the PLATO framework. Variations in the patient’s behavior and in physiological dynamics are indeed the norm, and can greatly affect the performance of a deterministic policy, with health consequences for the patient. In contrast, in a Bayesian framework, such variations likely yield high predictive uncertainty, information that we actively use to make robust therapy decisions.

In summary, the main contribution of this paper is an Imitation Learning-based method for deriving Bayesian neural network policies for AP control, a method that overcomes the two main shortcomings of established MPC-based approaches, namely, SE errors and computational cost. We empirically show that: 1) our IL-based approach outperforms Supervised Learning, and with less supervision data; 2) the learned stochastic policies outperform MPC with SE and deterministic policies; 3) our stochastic policies generalize to never-before-seen disturbance distributions and patient parameters that arise in a virtual patient cohort, while in the same setting, MPC with SE show consistent performance degradation. Overall, our best stochastic policy keeps BG in the safe range 8.4%–11.75% longer than MPC with SE and 2.94%–9.07% longer than the deterministic policy.

2 Background on MPC for the AP

We consider an in silico AP system, where the T1D patient is represented by a differential-equation model of BG and insulin metabolism, including absorption, excretion and transport dynamics between different body compartments. In particular, we chose the well-established model of Hovorka et al. (2004), a nonlinear ODE model with 14 state variables. A diagram of the MPC-based AP system is showed in Fig. 3 (a). The state-space description is given in equations (13) below. The notation stands for the indexed sequence .

(1)
(2)
(3)
(4)

Equation (1) describes the T1D patient model, where is the patient state at time , is the insulin input, are the patient parameters with distribution , and is the meal disturbance (amount of ingested carbohydrate) with distribution . Note that both parameters and disturbances are random and time-dependent. Meal disturbances indeed depend on the patient’s eating behaviour, which follows daily and weekly patterns. Similarly, parameters are typically subject to intra-patient variations such as daily oscillations. The parameter distribution can be alternatively defined to capture a population of patients (as we do in the experiment of Section 6).

Variable is the observed output at time , i.e., the CGM measurement, which is subject to Gaussian sensor noise , as done in previous work (Soru et al., 2012; Patek et al., 2007). With we denote the MPC-based control policy, that, given state estimate and future meal disturbances as inputs, with being the MPC prediction horizon, computes the optimal insulin therapy by solving, at each time , the following online optimization problem (Paoletti et al., 2017, 2019)222At each time step, the solution is actually a sequence of insulin values, of which we retain only the first one.:

(5)

subject to

(6)
(7)
(8)
(9)

where is the control horizon; is the set of admissible insulin values; (7) states that is fixed to the basal insulin rate outside the control horizon. (8) and (9) describe how the predicted state evolves from estimated state at time following the plant model (1). The objective function contains two terms, the former designed to penalize deviations between predicted and target BG (via function ), the latter to avoid abrupt changes of the insulin infusion rate. The coefficient is a hyper-parameter. See (Paoletti et al., 2017, 2019) for more details. Note that in order to predict future BG profiles, MPC requires some information of the future disturbances. In this work, we assume that the true future disturbance values are known. Alternative MPC formulations with weaker assumptions (Hughes et al., 2011; Paoletti et al., 2019) could be considered alternatively.

Finally, (4) describes state estimation (performed by function ), which follows the so-called Moving Horizon Estimation (MHE) method (Rawlings, 2013). MHE is the dual of MPC in that MHE also uses model predictions but, this time, to estimate the most likely state given a sequence of past measurements, control inputs, and disturbances. The estimate is the state that minimizes the discrepancy between observed CGM outputs and corresponding model predictions. As such, the quality of the estimate directly depends on the accuracy of the predictive model and the quality of the measurements. Details of the MHE optimization problem for the AP can be found in (Chen et al., 2019)

. Alternative SE methods include Kalman filters and particle filters 

(Rawlings, 2013), but these are subject to the same kinds of estimation errors as MHE.

3 Overview of the Method

(a)
(b)
Figure 3: a) A typical MPC-based AP system, where the controller requires a state estimate. b) Our learned end-to-end insulin policy instead only requires (noisy) observations.

The main goal of our work is to design an end-to-end insulin control policy for the AP, i.e., policies that directly map noisy system outputs (CGM measurements) into an optimal insulin therapy, without requiring knowledge of the system state. The control loop of such a policy is illustrated in Fig. 3 (b). To this purpose, we follow an imitation learning approach where the learner policy is trained from examples provided by an MPC-based expert, called the teacher policy. Following the PLATO algorithm (Kahn et al., 2017), at training time, the teacher is instrumented with full state information, so that its demonstrations are not affected by SE errors.

Figure 4: Overview of the IL scheme for the AP. Training-time trajectories are generated by applying the adaptive teacher policy, which is a trade-off between the optimal supervision policy and the learner’s policy. For training the learner, visited states are labelled with the optimal action of the supervision policy.

Supervised learning (SL) is not a viable solution to our problem because SL assumes i.i.d. training and test data, while our case is subject to covariate shift: the training data follow the distribution of trajectories visited through the teacher, while at test time we have a different distribution, induced by the learner policy. If the learner cannot imitate perfectly the teacher, in the long-run the learner’s actions will bring the system into a state far from the training distribution, where the behaviour of the learner becomes unpredictable, and with clear safety implications for our AP system. In IL, the teacher should provide demonstrations on trajectories that the learner would explore, but without knowing the learner in advance. To do so, a common solution is to reduce IL into a sequence of SL problems, where at each iteration, the learner is trained on the distribution induced by the previous learners or by a “mixture” of learner and teacher policies (Ross and Bagnell, 2010; Ross et al., 2011).

Our IL method builds on the PLATO algorithm for adaptive trajectory optimization (Kahn et al., 2017). See Section 7 for a summary of our extensions to PLATO. PLATO also reduces IL into a sequence of SL problems, where at each iteration, the teacher’s actions gradually adapt to those of the current learner policy. This adaptive teacher is an MPC-based policy whose objective function is extended to match the learner’s behaviour. As such, it is non-optimal, and is only used to generate trajectories that approach the distribution induced by the learner, thereby alleviating the covariate shift problem. The training data is obtained by labelling the adaptive teacher trajectories with the original (optimal) MPC policy. We call the latter supervision policy. The training process is summarized in Fig. 4. In our approach the learner policy is stochastic, represented as a Bayesian neural network, as explained in Section 5.

4 Imitation Learning Algorithm

The supervision policy is the MPC policy for the AP system of (59). Its control action at time , , is obtained by running the MPC algorithm given the true system state (instead of the estimated state) and future meal disturbances: . The learner policy is represented by an LSTM network with parameters . The choice of a recurrent architecture is natural for our application, because our policies have to subsume both control and SE and, as discussed in Section 2, SE for nonlinear models can be seen as a sequence prediction problem, see (4).

At time , the learner policy takes an input sequence, which includes past observations , past and future disturbances , and past control actions , and predicts the control action at time , , by iterating the following equation:

(10)
(11)
(12)

where is the hidden state of the network at time , and and

are activation functions. The learner results from the composition of

and . Since we use a stochastic policy for (see Section 5), we will denote the learner by the conditional distribution .

The adaptive teacher policy extends the supervision policy with a term to penalize the mismatch between the supervision policy, i.e., the optimal MPC policy, and the learner. The output of is the first control action in the solution of the following MPC problem.

(13)
subject to (69).

The first term of (13) is the cost function of the supervision policy (5). The second term penalizes the distance between the optimal (supervision) control action and the action that the learner policy would perform. The factor determines the relative importance of matching the behavior of the learner policy against minimizing the cost . In our experiments, we set where is the iteration of our IL algorithm, see Algorithm 1. In this way, the relative importance of matching increases as the learner improves (i.e., as increases). By gradually matching the behaviour of the learner, the control actions taken by will lead to exploring a state space similar to the one induced by the learner policy.

The distance term is defined as the -th Wasserstein distance between the output distribution of the learner policy and the optimal control action , or more precisely, the Dirac distribution centered at , . We choose because the two distributions are likely to have disjoint supports, and so other measures like the KL divergence (used in PLATO) are undefined or infinite. As we will see, an analytical form for is not available, but we can sample from it and obtain the empirical approximation , where are the sampled values of . Thus, we have that

(14)

where is a suitable distance. Note that the above equality holds as we are comparing a discrete distribution with a Dirac distribution. We set and .

Algorithm 1 outlines our IL scheme, which consists of a sequence of SL iterations. The learner policy is initialized at random. At each iteration, we start from a random initial state of the plant and generate a trajectory of the system of length . To do so, we first sample a sequence of random meal disturbances, i.e., of ingested carbohydrates (lines 5-7). Then, for each time point of the trajectory, the adaptive teacher computes an insulin infusion rate by solving (13), that is, by optimizing the MPC objective while matching the current learner policy (line 11). The MPC supervision policy is used to compute the optimal control action by solving (5) (line 12). The optimal action is used to label the corresponding training example , which is added to the training set (line 13), while the sub-optimal action by the adaptive teacher is used to evolve the system state (line 14). At the end of each iteration, is trained using the current set of examples . As the teacher gradually matches the learner, the training-time distribution of the system state gradually approximates that at test time.

The time complexity of Algorithm 1 is dominated by the two MPC instances (lines 11 and 12), which are solved times, and the training of the learner , repeated times.

1:  Initialize training dataset .
2:  Randomly initialize learner policy .
3:  Sample random patient parameters
4:  for  to  do
5:     for  to  do
6:        Sample random carb disturbance
7:     end for
8:     Initialize the patient state .
9:     for  to  do
10:        Collect CGM measurements as per (2)
11:        Compute sub-optimal therapy with adaptive teacher: .
12:        Compute optimal therapy with supervision policy: .
13:        Append to .
14:        State evolves as
15:     end for
16:     Train on
17:  end for
Algorithm 1 IL for Artificial Pancreas control policy using adaptive trajectory optimization

5 Bayesian inference of control policy

We take a Bayesian approach to learn our control policy, which results in a stochastic policy represented by a neural network with fixed architecture and random parameters . This provides us with a distribution of policy actions, the predictive distribution, from which we can derive uncertainty measures to inform the final therapy decision. Such uncertainty should capture both data uncertainty, due to e.g., noisy measurements, and epistemic uncertainty, i.e., the lack of confidence of the model about a given input (Gal, 2016).

Performing Bayesian inference corresponds, given training data , to computing the posterior from some prior by applying Bayes rule. The distribution of policies is induced by the random parameters . In particular, the predictive distribution , with being the policy inputs, is obtained from as:

(15)

For the non-linearity of the neural network function, precise inference is, however, infeasible and thus one needs to resort to approximate methods (Neal and others, 2011; Blundell et al., 2015). Monte Carlo Dropout (MCD) (Gal and Ghahramani, 2016)

is one of these. Dropout is a well-established regularization technique based on dropping some neurons at random during training with some probability

. Technically, this corresponds to multiplying the weights with a dropout mask

, i.e., a vector of Bernoulli variables with parameter

. Then, at test time, standard dropout derives a deterministic network by scaling back the weights by a factor of . On the other hand, in MCD the random dropout mask is applied at test time too, resulting in a distribution of network parameters. Gal and Ghahramani (2016) show that applying dropout to each weight layer is equivalent to performing approximate Bayesian inference of the neural network. This property is very appealing as it reduces the problem of inferring to drawing Bernoulli samples, which is very efficient.

Sample size.

For the empirical approximation of (15), we choose a sample size of that, by the Dvoretzky-Kiefer-Wolfowitz inequality (Massart, 1990), gives us the following statistical guarantee:

where is the CDF of the predictive distribution and its empirical approximation with observations.

Decision rule.

The output of our policy is the distribution of insulin actions (15). Hence, we need to define a decision rule that produces an individual value out of this distribution, and we want such rule to account for the predictive uncertainty of the policy. Let be the ordered empirical sample of (15), and let be the last measured glucose value. Our rule selects a particular order statistic , i.e., one of the sampled values, depending on the relative distance of w.r.t. the safe BG upper bound and lower bound . We call this adaptive rule because the selected order statistic is adapted on . Formally,

(16)

In this way, if the patient is approaching hypoglycemia ( close to ), we select a conservative insulin value, and we select instead an aggressive therapy if is close to hyperglycemia (). Importantly, as the policy uncertainty increases (and so does the spread of the sample), gets more conservative when is in the lower half of the safe BG range, i.e., we take safer decisions when the policy is less trustworthy333Protecting against hypoglycemia is the primary concern.. For the same principle, when is in the upper half of the range, higher uncertainty yields a more aggressive therapy, but this poses no safety threats because is well away from the hypoglycemia threshold .

In our evaluation, we compare our adaptive rule with the commonly used rule that sets to the sample mean of (15).

6 Experiments

We conducted computational experiments to validate the following claims:

  1. Our IL-based approach outperforms SL, and with less supervision data.

  2. The learned stochastic policies outperform both MPC with SE and deterministic policies.

  3. The stochastic policies generalize well to unseen patient physiological parameters and meal disturbances.

  4. For a stochastic policy, the adaptive decision rule (16) outperforms the mean prediction rule.

  5. The predictive uncertainty of the policy increases with out-of-distribution test inputs.

Runtime cost.

On our workstation444With an Intel i7-8750H CPU and 16GB DDR4 SDRAM the stochastic policy executes in approximately 20 milliseconds, which is well within the CGM measurement period of 5 minutes, and consistently more efficient than MPC-based optimization (on average, 150+ times faster in our experiments)555Optimization for MPC and SE is solved used MATLAB’s implementation of the interior-point algorithm of Byrd et al. (2000).

. Thanks to platforms such as TensorFlow Lite, we believe that a similar runtime can be obtained also after deploying the policy on embedded hardware.

LSTM architecture.

We represent the learner policy as an LSTM network with three layers, of 200 hidden units each, with tanh state activation functions and sigmod

gate activation functions. The network has a fully connected layer, a ReLu activation function layer and a regression layer as final layers. We added a dropout layer before each weight layer with dropout probability

, which empirically gave us the best accuracy among other probability values. We also experimented with feedforward architectures but these performed poorly, confirming the need for a recurrent architecture to adequately represent SE and control.

Performance measures.

To evaluate and compare the insulin policies, we consider three typical performance measures in the AP domain: : the average percentage of time the virtual patient’s BG is in hypoglycemia (i.e., BG 70 mg/dL), hyperglycemia (i.e., BG 180 mg/dL), and in the safe range, or euglycemia (i.e., 70 mg/dL BG 180 mg/dL), respectively. Given the goal of diabetes treatment and BG control, the target is to maximize while minimizing (hypoglycemia leads to more serious consequences than hyperglycemia).

Experimental settings.

During training, we consider a fixed (and time-constant) parameterization of the T1D virtual patient model, and the disturbance distribution of Table 1 (which includes three large meals and three small snacks). The length of each training trajectory is set to . We performed iterations of the IL algorithm, totalling to a training set size of . At test time, we consider the following three configurations of the model parameters:

  1. Fixed Patient. The test-time virtual patient has the same parameters as during training.

  2. Varying Patient. We introduce intra-patient temporal variations to the parameters.

  3. Patient cohort. We introduce both inter-patient and intra-patient variations. The distribution of model parameters represents a cohort of virtual patients, where each patient has time-varying parameters. It might include severe cases when the parameters are sampled at marginal areas of the distribution.

The above distributions for modelling intra- and inter-patient variability are taken from (Wilinska et al., 2010), where the authors derived them from clinical data.

For each configuration and policy, we simulate 90 trajectories of a one-day simulation. For each trajectory, we draw a fresh observation of the random disturbance and patient parameters, but, to ensure a fair comparison, these are kept the same for all control algorithms.

In another set of experiments, we consider (in combination with the above configurations) a meal disturbance distribution different from the one at training time. This meal distribution reflects a late eating habit and a larger probability of larger snacks, which leads to meal profiles with higher carbs. Further details are reported in the supplement.

For each experiment, we evaluate the performance of the following policies:

  • MPC+SI: the supervision policy , i.e., the MPC algorithm with full information on the current state and model parameters. We stress that full observability is an unrealistic assumption, yet useful to establish an ideal performance level. We set control and prediction horizons of MPC to min and min, respectively.

  • MPC+SE: the MPC algorithm with state estimation, as described in (59). This is one of the most common control schemes for the AP. The parameters of the T1D model used within MPC and SE are set to the training-time values. We will show that SE takes a toll on the overall performance, unlike our end-to-end policies.

  • DLP: the deterministic learner policy, where dropout is used at training time but not at test time.

  • SLP-M: the stochastic learner policy with the rule that selects the sample mean of the predictive distribution.

  • SLP-A: the stochastic learner policy with the adaptive rule of (16). As explained above, such rule actively uses uncertainty information by selecting an order statistic of the predictive distribution commensurate to the sensed glucose. We will show that this policy outperforms all the others (but MPC+SI).

IL outperforms SL.

We compared the performance of SLP-A, which is trained via the IL algorithm of Section 4, against the corresponding stochastic policy trained using an MPC-guided SL approach, that is, without the adaptive teacher policy that guides the exploration of training trajectories but only using the supervision policy.

(a) Fixed patient
(b) Varying patient
(c) Patient cohort
Figure 8: Mean BG standard deviation for MPC with state estimation vs stochastic policy with adaptive output.
breakfast snack 1 lunch snack 2 dinner snack 3
Occurrence 100 50 100 50 100 50
Probability (%)
CHO Amount 40-60 5-25 70-110 5-25 55-75 5-15
(gram)
Time of 1:00- 5:00- 8:00- 12:00- 15:00- 19:00-
day (hour) 5:00 8:00 12:00 15:00 19:00 21:00
Table 1: Attributes of meal disturbance distribution during training. CHO amounts and starting times are sampled uniformly from the reported intervals.

We found that the IL-based policy has much better control performance than the SL-based one, using the same number of algorithm iterations and training examples. This suggests that the IL-based policy captures more useful state-action information and that the IL algorithm leads the learner policy to efficiently explore the state space. In particular, after 20 iterations of Algorithm 1, SLP-A achieves an average of for (time in range), while the SL-based counterpart only . IL consistently outperforms SL in all other performance metrics too (see Appendix B for more detailed results).

Stochastic policies outperform MPC with SE and deterministic policies.

Results in Table 2 evidence that the stochastic policies (SLP-A and SLP-M) outperform the MPC policy with state estimation in essentially all performance measures and configurations, and in a statistically significant manner. In particular, we observe that, on average, the SLP-A policy stays in range for 8.4%–11.75% longer than MPC+SE, and the SLP-M does so for 5.97%–8.4% longer. This is visible also from the BG profiles of Figure 8. The deterministic policy also performs better than MPC+SE, but lags behind the stochastic ones, with SLP-A achieving a time in range 2.94%–9.07% longer than DLP.

These results suggests that MPC with SE introduces estimation errors that have a detrimental effect on BG control, as also confirmed by the fact that the same control algorithm but with full state information (MPC+SI) is way superior. In realistic settings where the true patient state is not accessible, our analysis shows that an end-to-end policy is to be preferred to explicit SE.

Stochastic policies generalize to unseen patient parameters and disturbances.

From Table 2 and Figure 8, we observe that our stochastic policies are robust to patient parameters unseen during training, with a time in range constantly above 85% for the varying patient and the patient cohort configurations. The superiority of the stochastic policies over the deterministic one under intra-patient and inter-patient variability evidence the role of considering prediction uncertainty when the policy is deployed in environments that deviate from the training one.

Fixed Patient Configuration Varying Patient Configuration Patient Cohort Configuration
(%) (%) (%) (%) (%) (%) (%) (%) (%)
MPC+SI 0.00 0.00 99.84 0.60 0.16 0.60 0.000.00 99.800.83 0.200.83 0.000.00 99.151.81 0.851.81
MPC+SE 0.922.76 80.405.33 18.684.76 1.443.18 79.885.46 18.694.95 0.853.57 77.1715.39 21.9814.30
DLP 0.210.99 85.404.15 14.383.90 0.531.87 82.666.42 16.826.25 0.451.91 82.636.28 16.936.15
SLP-M 0.231.10 86.334.17 13.444.09 0.331.09 86.344.20 13.334.06 0.321.19 85.454.34 14.234.00
SLP-A 0.130.64 91.733.45 8.143.39 0.050.27 91.733.18 8.233.15 0.040.41 85.574.12 14.394.08
Table 2: Performance of MPC with state information (MPC+SI), MPC with state estimation (MPC+SE), deterministic learner policy (DLP), stochastic learner policy with mean value output (SLP-M), and stochastic learner policy with adaptive output (SLP-A). For each column, in bold are the significantly best policies, i.e., with in all pairwise sign tests (one-sided) (Hollander et al., 2013). A complete, detailed table of performance statistics and the processes of hypothesis tests can be found in Appendix A.

We further evaluate the robustness of the SLP-A policy under an unseen meal disturbance distribution characterized by a higher carbs intake. Results, reported in Table 3, evidence that our stochastic policy generalize well also in this case. A hypothesis test can be found in Appendix D, showing that the performance statistics SLP-A achieves under intra-patient and inter-patient variability with two different meal distributions has no significant difference.

Adaptive rule outperforms mean-value rule.

The SLP-A policy obtains a time in range approximately 5% longer than SLP-M in the fixed patient and varying patient configurations, see Table 2. There is an improvement also in the patient cohort experiment, albeit less significant. With the adaptive rule, the policy adopts a more conservative or aggressive therapy depending on the measured glucose and the predictive uncertainty, which can lead to more stable BG trajectories and explain the observed difference.

Policy uncertainty increases at out-of-distribution test inputs.

We verify that inference via Monte Carlo dropout adequately captures epistemic uncertainty. We indeed observe an increased variability of the predictive distribution in the patient cohort configuration, that is, when the policy is applied to out-of-distribution test data resulting from introducing inter- and intra-patient variability. However, uncertainty does not significantly increase from the first to the second configuration, that is, when introducing only inter-patient variability. These results were confirmed via pair-wise two-sample Kolmogorov-Smirnov tests over the coefficients of variation (CoVs), i.e., the mean-normalized standard deviations (differences significant at ). See Appendix C for the elaboration on the Kolmogorov-Smirnov tests and a plot of the CoVs distributions.

Performance Fixed Patient Varying Patient Patient Cohort
Metrics Configuration Configuration Configuration
()
()
()
Table 3: Performance of SLP-A with unseen disturbances.

7 Related Work

Traditional methods for insulin control in the AP mostly rely on MPC (Hovorka et al., 2004; Magni et al., 2007; Lee et al., 2009; Cameron et al., 2011; Paoletti et al., 2019), PID control (Steil et al., 2004; Huyett et al., 2015), and fuzzy rules based on the reasoning of diabetes caregivers (Atlas et al., 2010; Nimri et al., 2012)

. Reinforcement learning approaches have been proposed as well, including policy iteration 

(De Paula et al., 2015; Weng et al., 2017), actor critic methods (Daskalaki et al., 2016), and deep Q-networks for dual-hormone therapy (Zhu et al., 2019)666Dual-hormone therapy allows infusion of both insulin and its antagonist hormone, glucagon, which protects against hypoglycemia. However, it is not as established as insulin-only therapy, and still has safety concerns (Taleb et al., 2017).. Neural networks have been explored in the AP space not just to represent insulin policies (de Canete et al., 2012; Zhu et al., 2019), but also to predict BG concentration based on inputs such as insulin dosage, nutritional intake, and daily activities (Pappada et al., 2008; Pérez-Gandía et al., 2010; Dutta et al., 2018; Li et al., 2019). Our work is different from the above papers as: 1) it takes an imitation learning approach to learn the insulin policy, mitigating potential (and dangerous) test-time distribution drifts; 2) incorporates in the policy uncertainty information obtained via Bayesian inference; 3) it produces end-to-end policies that do not require learning a separate BG prediction model.

Several imitation learning methods have been proposed, such as (Daumé et al., 2009; Ross and Bagnell, 2010; Ross et al., 2010; Ho and Ermon, 2016). Some of these, like our approach, are tailored to work with MPC teachers (Zhang et al., 2016; Kahn et al., 2017; Amos et al., 2018; Lowrey et al., 2019). Our method is also akin to (Cronrath et al., 2018; Menda et al., 2017; Lee et al., 2019) where Bayesian extensions of IL are presented to quantify the learner’s predictive uncertainty and better guide the queries to the teacher policy. Other Bayesian approaches in policy learning include (Deisenroth and Rasmussen, 2011; Gal et al., 2016; Polymenakos et al., 2019). “Frequentist” alternatives to Bayesian uncertainty quantification also exist (Vovk et al., 2005; Wang and Qiao, 2018; Papernot and McDaniel, 2018; Bortolussi et al., 2019; Park et al., 2020).

Our work builds on the PLATO IL algorithm (Kahn et al., 2017) and extends it in three main direction: 1) we consider recurrent architectures, which are more suitable than feedforward ones (used in PLATO) to represent nonlinear state estimation and control; 2) PLATO also derives stochastic policies but, unlike our work, no uncertainty-aware decision-making strategies are considered; 3) PLATO policies do not support systems with external disturbances beyond noise. In our policies instead, random meal disturbances are central.

8 Conclusion

We introduced a method based on MPC-guided Imitation Learning and Bayesian inference to derive stochastic Neural Network policies for insulin control in an artificial pancreas. Our policies are end-to-end in that they directly map CGM values into insulin control actions. By using Bayesian neural networks, we can crucially quantify prediction uncertainty, information that we incorporate in the insulin therapy decisions. We empirically demonstrated that our stochastic insulin policies outperform traditional MPC with explicit state estimation and are more robust than their deterministic counterparts, as they generalize well to unseen T1D patient parameters and meal disturbance distributions.

References

  • B. Amos, I. D. J. Rodriguez, J. Sacks, B. Boots, and J. Z. Kolter (2018) Differentiable MPC for end-to-end planning and control. In Proceedings of the 32nd International Conference on Neural Information Processing Systems, NIPS’18. Cited by: §7.
  • E. Atlas, R. Nimri, S. Miller, E. A. Grunberg, and M. Phillip (2010) MD-logic artificial pancreas system: a pilot study in adults with type 1 diabetes. Diabetes care 33 (5), pp. 1072–1076. Cited by: §1, §7.
  • C. Blundell, J. Cornebise, K. Kavukcuoglu, and D. Wierstra (2015) Weight uncertainty in neural networks. arXiv preprint arXiv:1505.05424. Cited by: §5.
  • L. Bortolussi, F. Cairoli, N. Paoletti, S. A. Smolka, and S. D. Stoller (2019) Neural predictive monitoring. In International Conference on Runtime Verification (RV 2019), pp. 129–147. Cited by: §7.
  • R. H. Byrd, J. C. Gilbert, and J. Nocedal (2000) A trust region method based on interior point techniques for nonlinear programming. Mathematical programming 89 (1), pp. 149–185. Cited by: footnote 5.
  • E. F. Camacho and C. B. Alba (2013) Model predictive control. Springer Science & Business Media. Cited by: §1.
  • F. Cameron, B. W. Bequette, D. M. Wilson, B. A. Buckingham, H. Lee, and G. Niemeyer (2011) A closed-loop artificial pancreas based on risk management. Journal of diabetes science and technology 5 (2), pp. 368–379. Cited by: §7.
  • H. Chen, N. Paoletti, S. Smolka, and S. Lin (2019) Committed moving horizon estimation for meal detection and estimation in type 1 diabetes. In Proc. of the 2019 American Control Conference (ACC), Cited by: §2.
  • C. Cronrath, E. Jorge, J. Moberg, M. Jirstrand, and B. Lennartson (2018) BAgger: a Bayesian algorithm for safe and query-efficient imitation learning. In Machine Learning in Robot Motion Planning – IROS 2018 Workshop, Cited by: §7.
  • E. Daskalaki, P. Diem, and S. Mougiakakou (2016) Model-free machine learning in biomedicine: feasibility study in type 1 diabetes. PloS one. Cited by: §7.
  • H. Daumé, J. Langford, and D. Marcu (2009) Search-based structured prediction. Machine learning 75 (3), pp. 297–325. Cited by: §7.
  • J. F. de Canete, S. Gonzalez-Perez, and J. Ramos-Diaz (2012) Artificial neural networks for closed loop control of in silico and ad hoc type 1 diabetes. Computer methods and programs in biomedicine 106 (1), pp. 55–66. Cited by: §7.
  • M. De Paula, L. O. Ávila, and E. C. Martínez (2015) Controlling blood glucose variability under uncertainty using reinforcement learning and gaussian processes. Applied Soft Computing 35, pp. 310–332. Cited by: §7.
  • M. Deisenroth and C. E. Rasmussen (2011) PILCO: a model-based and data-efficient approach to policy search. In Proceedings of the 28th International Conference on machine learning (ICML-11), pp. 465–472. Cited by: §7.
  • S. Dutta, T. Kushner, and S. Sankaranarayanan (2018) Robust data-driven control of artificial pancreas systems using neural networks. In International Conference on Computational Methods in Systems Biology, pp. 183–202. Cited by: §7.
  • G. P. Forlenza, L. Ekhlaspour, M. Breton, D. M. Maahs, R. P. Wadwa, M. DeBoer, L. H. Messer, M. Town, J. Pinnata, G. Kruse, et al. (2019) Successful at-home use of the tandem control-iq artificial pancreas system in young children during a randomized controlled trial. Diabetes technology & therapeutics 21 (4), pp. 159–169. Cited by: footnote 1.
  • Y. Gal and Z. Ghahramani (2016)

    Dropout as a bayesian approximation: representing model uncertainty in deep learning

    .
    In international conference on machine learning, pp. 1050–1059. Cited by: §1, §5.
  • Y. Gal, R. McAllister, and C. E. Rasmussen (2016) Improving pilco with bayesian neural network dynamics models. In Data-Efficient Machine Learning workshop, ICML, Vol. 4, pp. 34. Cited by: §7.
  • Y. Gal (2016) Uncertainty in deep learning. Ph.D. Thesis, PhD thesis, University of Cambridge. Cited by: §1, §5.
  • R. Gondhalekar, E. Dassau, and F. J. Doyle (2014) Moving-horizon-like state estimation via continuous glucose monitor feedback in MPC of an artificial pancreas for type 1 diabetes. In Decision and Control (CDC), 2014 IEEE 53rd Annual Conference on, pp. 310–315. Cited by: §1.
  • J. Ho and S. Ermon (2016) Generative adversarial imitation learning. In Advances in neural information processing systems, pp. 4565–4573. Cited by: §7.
  • S. Hochreiter and J. Schmidhuber (1997) Long short-term memory. Neural computation 9 (8), pp. 1735–1780. Cited by: §1.
  • M. Hollander, D. A. Wolfe, and E. Chicken (2013) Nonparametric statistical methods. Vol. 751, John Wiley & Sons. Cited by: Table 2.
  • R. Hovorka, V. Canonico, L. J. Chassin, U. Haueter, M. Massi-Benedetti, M. O. Federici, T. R. Pieber, H. C. Schaller, L. Schaupp, T. Vering, et al. (2004) Nonlinear model predictive control of glucose concentration in subjects with type 1 diabetes. Physiological measurement 25 (4), pp. 905. Cited by: §1, §2, §7.
  • C. Hughes, S. D. Patek, M. Breton, and B. P. Kovatchev (2011) Anticipating the next meal using meal behavioral profiles: a hybrid model-based stochastic predictive control algorithm for t1dm. Computer methods and programs in biomedicine 102 (2), pp. 138–148. Cited by: §2.
  • L. M. Huyett, E. Dassau, H. C. Zisser, and F. J. Doyle III (2015) Design and evaluation of a robust PID controller for a fully implantable artificial pancreas. Industrial & engineering chemistry research 54 (42), pp. 10311–10321. Cited by: §7.
  • G. Kahn, T. Zhang, S. Levine, and P. Abbeel (2017) PLATO: policy learning using adaptive trajectory optimization. In Robotics and Automation (ICRA), 2017 IEEE International Conference on, pp. 3342–3349. Cited by: §1, §3, §3, §7, §7.
  • H. Lee, B. A. Buckingham, D. M. Wilson, and B. W. Bequette (2009) A closed-loop artificial pancreas using model predictive control and a sliding meal size estimator. Journal of Diabetes Science and Technology 3 (5), pp. 1082–1090. Cited by: §7.
  • K. Lee, K. Saigol, and E. A. Theodorou (2019) Safe end-to-end imitation learning for model predictive control. In International conference on robotics and automation (ICRA), Cited by: §7.
  • K. Li, C. Liu, T. Zhu, P. Herrero, and P. Georgiou (2019) GluNet: a deep learning framework for accurate glucose forecasting. IEEE journal of biomedical and health informatics. Cited by: §7.
  • K. Lowrey, A. Rajeswaran, S. Kakade, E. Todorov, and I. Mordatch (2019) Plan online, learn offline: efficient learning and exploration via model-based control. In International Conference on Learning Representations, Cited by: §7.
  • L. Magni, D. M. Raimondo, L. Bossi, C. Dalla Man, G. De Nicolao, B. Kovatchev, and C. Cobelli (2007) Model predictive control of type 1 diabetes: an in silico trial. Journal of Diabetes Science and Technology 1 (6), pp. 804–812. Cited by: §7.
  • C. D. Man, F. Micheletto, D. Lv, M. Breton, B. Kovatchev, and C. Cobelli (2014) The UVA/PADOVA type 1 diabetes simulator: new features. Journal of diabetes science and technology 8 (1), pp. 26–34. Cited by: §1, footnote 1.
  • P. Massart (1990) The tight constant in the Dvoretzky-Kiefer-Wolfowitz inequality. The annals of Probability, pp. 1269–1283. Cited by: §5.
  • K. Menda, K. Driggs-Campbell, and M. J. Kochenderfer (2017) Dropoutdagger: a bayesian approach to safe imitation learning. arXiv preprint arXiv:1709.06166. Cited by: §7.
  • R. M. Neal et al. (2011) MCMC using hamiltonian dynamics.

    Handbook of markov chain monte carlo

    2 (11), pp. 2.
    Cited by: §5.
  • R. Nimri, E. Atlas, M. Ajzensztejn, S. Miller, T. Oron, and M. Phillip (2012) Feasibility study of automated overnight closed-loop glucose control under MD-logic artificial pancreas in patients with type 1 diabetes: the DREAM project. Diabetes technology & therapeutics 14 (8), pp. 728–735. Cited by: §7.
  • N. Paoletti, K. S. Liu, H. Chen, S. Smolka, and S. Lin (2019) Data-driven robust control for a closed-loop artificial pancreas. IEEE/ACM transactions on computational biology and bioinformatics. Cited by: §1, §2, §2, §7.
  • N. Paoletti, K. S. Liu, S. A. Smolka, and S. Lin (2017) Data-driven robust control for type 1 diabetes under meal and exercise uncertainties. In International Conference on Computational Methods in Systems Biology, pp. 214–232. Cited by: §2, §2.
  • N. Papernot and P. McDaniel (2018) Deep k-nearest neighbors: towards confident, interpretable and robust deep learning. arXiv preprint arXiv:1803.04765. Cited by: §7.
  • S. M. Pappada, B. D. Cameron, and P. M. Rosman (2008) Development of a neural network for prediction of glucose concentration in type 1 diabetes patients. Journal of diabetes science and technology 2 (5), pp. 792–801. Cited by: §7.
  • S. Park, O. Bastani, N. Matni, and I. Lee (2020) PAC confidence sets for deep neural networks via calibrated prediction. In International Conference on Learning Representations, Cited by: §7.
  • S. D. Patek, M. D. Breton, Y. Chen, C. Solomon, and B. Kovatchev (2007) Linear quadratic gaussian-based closed-loop control of type 1 diabetes. SAGE Publications. Cited by: §2.
  • C. Pérez-Gandía, A. Facchinetti, G. Sparacino, C. Cobelli, E. Gómez, M. Rigla, A. de Leiva, and M. Hernando (2010) Artificial neural network algorithm for online glucose prediction from continuous glucose monitoring. Diabetes technology & therapeutics 12 (1), pp. 81–88. Cited by: §7.
  • K. Polymenakos, A. Abate, and S. Roberts (2019) Safe policy search using gaussian process models. In Proceedings of the 18th International Conference on Autonomous Agents and MultiAgent Systems, pp. 1565–1573. Cited by: §7.
  • J. B. Rawlings (2013) Moving horizon estimation. Encyclopedia of Systems and Control, pp. 1–7. Cited by: §2.
  • S. Ross and D. Bagnell (2010) Efficient reductions for imitation learning. In

    Proceedings of the thirteenth international conference on artificial intelligence and statistics

    ,
    pp. 661–668. Cited by: §3, §7.
  • S. Ross, G. Gordon, and D. Bagnell (2011) A reduction of imitation learning and structured prediction to no-regret online learning. In Proceedings of the fourteenth international conference on artificial intelligence and statistics, pp. 627–635. Cited by: §1, §3.
  • S. Ross, G. Gordon, and J. Bagnell (2010) A reduction of imitation learning and structured prediction to no-regret online learning. Journal of Machine Learning Research - Proceedings Track. Cited by: §7.
  • P. Soru, G. De Nicolao, C. Toffanin, C. Dalla Man, C. Cobelli, L. Magni, A. H. Consortium, et al. (2012) MPC based artificial pancreas: strategies for individualization and meal compensation. Annual Reviews in Control 36 (1), pp. 118–128. Cited by: §2.
  • G. M. Steil, A. E. Panteleon, and K. Rebrin (2004) Closed-loop insulin delivery—the path to physiological glucose control. Advanced drug delivery reviews 56 (2), pp. 125–144. Cited by: §7.
  • G. M. Steil (2013) Algorithms for a closed-loop artificial pancreas: the case for proportional-integral-derivative control. Journal of diabetes science and technology 7 (6), pp. 1621–1631. Cited by: §1.
  • A. Storkey (2009) When training and test sets are different: characterizing learning transfer. Dataset shift in machine learning, pp. 3–28. Cited by: §1.
  • N. Taleb, A. Haidar, V. Messier, V. Gingras, L. Legault, and R. Rabasa-Lhoret (2017) Glucagon in artificial pancreas systems: potential benefits and safety profile of future chronic use. Diabetes, Obesity and Metabolism 19 (1), pp. 13–23. Cited by: footnote 6.
  • V. Vovk, A. Gammerman, and G. Shafer (2005) Algorithmic learning in a random world. Springer Science & Business Media. Cited by: §7.
  • W. Wang and X. Qiao (2018)

    Learning confidence sets using support vector machines

    .
    In Advances in Neural Information Processing Systems, pp. 4929–4938. Cited by: §7.
  • W. Weng, M. Gao, Z. He, S. Yan, and P. Szolovits (2017) Representation and reinforcement learning for personalized glycemic control in septic patients. arXiv preprint arXiv:1712.00654. Cited by: §7.
  • M. E. Wilinska, L. J. Chassin, C. L. Acerini, J. M. Allen, D. B. Dunger, and R. Hovorka (2010) Simulation environment to evaluate closed-loop insulin delivery systems in type 1 diabetes. Journal of diabetes science and technology 4 (1), pp. 132–144. Cited by: §6.
  • T. Zhang, G. Kahn, S. Levine, and P. Abbeel (2016) Learning deep control policies for autonomous aerial vehicles with MPC-guided policy search. In 2016 IEEE international conference on robotics and automation (ICRA), pp. 528–535. Cited by: §7.
  • T. Zhu, K. Li, and P. Georgiou (2019) A dual-hormone closed-loop delivery system for type 1 diabetes using deep reinforcement learning. arXiv preprint arXiv:1910.04059. Cited by: §7.

Appendix

Appendix A Significance Tests Details

Detailed statistics for all control policies and all configurations are available in Table 4. We show extra performance measures in this table: and , average of maximum and minimum of BG values (in mg/dL); , average insulin infusion rate (in mU/min). We performed significance tests out of the collected data to validate our hypotheses regarding the performance of the learned stochastic policies. Details are given below.

Configuration (%) (%) (%) (mg/dL) (mg/dL) u (mU/min)
MPC+SI Fixed Patient 0.00 + 0.00 99.84+ 0.60 0.16+ 0.60 181.28+13.64 123.40+3.85 28.11+0.95
Varying Patient 0.00+0.00 99.80+0.83 0.20+0.83 181.55+12.72 123.40+4.68 28.17+1.05
Patient Cohort 0.00+0.00 99.15+1.81 0.85+1.81 184.56+18.06 115.82+9.68 20.74+4.62
MPC+SE Fixed Patient 0.92+2.76 80.40+5.33 18.68+4.76 242.68+12.84 84.87+14.68 27.98+1.17
Varying Patient 1.44+3.18 79.88+5.46 18.69+4.95 239.79+12.22 81.44+14.80 28.16+1.05
Patient Cohort 0.85+3.57 77.17+15.39 21.98+14.30 241.17+32.43 95.45+21.73 20.55+6.52
DLP Fixed Patient 0.32+1.19 85.45+4.34 14.23+4.00 232.63+10.49 80.98+9.22 27.79+0.78
Varying Patient 0.53+1.87 82.66+6.42 16.82+6.25 225.66+10.40 82.27+10.91 28.59+0.76
Patient Cohort 0.45+1.91 82.63+6.28 16.93+6.15 225.58+10.41 83.14+10.95 28.59+0.79
SLP-M Fixed Patient 0.28+0.99 86.33+4.17 13.44+4.09 230.41+10.27 81.20+9.31 28.39+0.75
Varying Patient 0.33+1.09 86.34+4.20 13.33+4.06 230.16+10.37 81.01+9.77 28.38+0.73
Patient Cohort 0.21+0.99 85.40+4.15 14.38+3.90 232.72+10.27 81.13+8.98 28.39+0.72
SLP-A Fixed Patient 0.13+0.64 91.73+3.45 8.14+3.39 216.65+8.05 82.42+7.86 29.01+0.67
Varying Patient 0.05+0.27 91.73+3.18 8.23+3.15 216.98+7.68 82.50+8.35 28.97+0.69
Patient Cohort 0.04+0.41 85.57+4.12 14.39+4.08 232.02+8.92 86.84+9.58 28.61+0.73
Table 4: Full statistics of the experiments in which virtual T1D patient(s) are regulated by the different controllers/control policies

a.1 Hypotheses Formulation

We define our null and alternative hypotheses in the following for each of the performance measures: , , , , and .

  • Null hypothesis . The median of the distribution of the difference of a specific performance metric from two sets of experiments regulated by two control policies is zero. Note the hypothesis of zero median of is not equivalent to a hypothesis of equal median of and .

  • Alternative hypothesis . The median of distribution of the differences of a specific performance metric from two sets of experiments regulated by two control policies is larger (smaller) than zero.

This formulation is for a paired significance test, which is suitable when comparing, like we do, experimental results about same configuration (and hence, same sampled parameters and disturbances) but different policies.

a.2 Calculating statistic tests

For each combination of the compared policies, we use sign tests to test the hypothesis that medians of distribution of the difference between their , and are larger than zero, and medians of distribution of the difference between their , and are smaller than zero. In the main paper, we set the significance level to (i.e., confidence).

Configuration MPC+SE vs SLP-A DLP vs SLP-A SLP-M vs SLP-A MPC+SE vs DLP MPC+SE vs SLP-M DLP vs SLP-M
Fixed Patient 5.9090e-03 5.0000e-01 5.0000e-01 2.0695e-02 1.3302e-02 5.0000e-01
3.7664e-20 7.1194e-20 3.8869e-14 6.2463e-07 2.1929e-10 4.2922e-11
3.7664e-20 1.3449e-19 2.6012e-13 1.8975e-06 9.7412e-10 1.3568e-12
9.8208e-23 6.5759e-18 6.3970e-16 2.0826e-03 5.4510e-06 5.6818e-08
8.2858e-01 5.4194e-01 7.0079e-01 9.1488e-01 9.7770e-01 2.9921e-01
1.0000e+00 1.0000e+00 1.0000e+00 9.9897e-01 9.9897e-01 6.2398e-01
Varying Patient 2.4617e-05 5.8594e-03 4.6143e-02 1.4725e-02 6.2705e-03 5.8810e-01
2.1623e-21 6.3970e-16 5.4061e-19 5.6672e-02 3.6187e-10 3.8285e-05
5.4061e-19 5.2608e-15 5.4061e-19 8.5121e-02 2.9933e-07 9.3878e-05
3.7664e-20 4.0354e-09 3.8869e-14 9.7412e-10 6.2463e-07 9.9598e-01
5.4194e-01 5.6672e-02 1.7142e-01 3.7602e-01 6.2398e-01 5.4194e-01
1.0000e+00 9.9991e-01 1.0000e+00 9.9257e-01 9.9598e-01 2.2299e-02
Patient Cohort 1.9531e-02 3.5156e-02 1.0742e-02 2.9053e-01 6.8547e-01 7.7275e-01
1.0301e-03 6.6763e-04 4.5806e-01 1.2305e-01 2.1896e-04 4.0230e-03
1.0301e-03 2.7725e-03 6.2398e-01 1.2305e-01 2.1896e-04 1.0301e-03
1.2305e-01 9.9792e-01 3.7602e-01 1.3151e-02 3.7602e-01 9.9999e-01
9.9999e-01 7.4334e-03 2.0826e-03 9.9996e-01 9.9999e-01 7.0079e-01
1.0000e+00 7.6960e-01 9.1488e-01 1.0000e+00 1.0000e+00 4.0230e-03
Table 5: P values of the significance tests for all combinations of compared policies. In each column, the alternative hypothesis is the medians of the observations is larger than zero of the difference in the metrics (i.e., observations from the first policy minus those from the second) in the categories 1) , 2) and 3) , and is smaller than zero in the rest.

a.3 Conclusion about

We can see the tests reject when comparing and between SLP-A and all other policies, indicating they are statistical significant larger (smaller) when a fixed or varying patient is regulated by SLP-A than those by other policies. The p-values of the sign test on also support that SLP-A can significantly reduce hyperglycemia for a fixed or varying patient compared to all others.

In a patient cohort, the null hypothesis cannot be rejected for metrics and when comparing SLP-A VS SLP-M, indicating there is no significant difference in using these two stochastic policies. However, when comparing these two policies with any other policies, the null is rejected for both and . We can conclude that SLP-A and SLP-M significantly reduce hyperglycemia while maintaining short time in hypoglycemia.

Appendix B Comparison Between Supervised Learning and Imitation Learning

b.1 Performance improvement

Figure 11 (a) show the proportion of time blood glucose is in every state (eu-, hyper- and hypo-glycemia) at each training iterations using supervised learning, which improves much slower than Figure 11 (b).

We can see the learner policy trained by IL outperforms it by SL in reducing hypo-/hyper-glycemia and keeping patients’ blood glucose in healthy range (euglycemia). The mean of maximum blood glucose decreases faster while the mean of minimum blood glucose increases to prevent dangerous hypoglycemia when SLP-A is trained by IL methods.

(a) Supervised Learning
(b) Imitation Learning
Figure 11: Mean BG standard deviation on hypoglycemia, euglycemia, hyperglycemia for SLP-A trained by supervised learning v.s. imitation learning. Imitation learning methods explore the state space whose distribution that has a small distance to that in test time, so the SLP-A improves faster with IL-based methods.

Figure 12 is the root mean squared error (RMSE) of both learning methods at each training iteration. As we can see from the curves, the training RMSE for SL and IL goes down as the number of training episodes increases. The RMSE of IL starts larger but decrease dramatically especially after the first few episodes. Although the RMSE is smaller, stochastic learner policy trained by SL methods has a limited training state space and the out-of-distribution disturbance leads to slow learning process from supervision (i.e., MPC controller with perfect state information).

Figure 12: RSME of SL and IL at each training iteration.

Appendix C Uncertainty Analysis of Stochastic Learner Policy

In our machine learning model for the stochastic learner policy, we quantify the uncertainty of our predictive action. Quantifying the uncertainty is critical in real-world settings such as artificial pancreas (AP) system, in which often involve out-of-distribution disturbance and patient variations. We analyze the uncertainty quantifying ability of our Monte-Carlo Dropout Bayesian method to investigate the effect of covariate shift and virtual patient model variation.

c.1 Coefficients of variation

The coefficients of variation (CV), aslo known as relative standard deviation, is defined as the division of the standard deviation and the mean.

Figure 13

is the cumulative distribution functions (CDF) of the CV of the stochastic learner policy outputs when the learner policy is the control algorithm of an AP system under three patient configuration as illustrated in Section 6 of our paper. The curve of patient cohort experiments is the slowest to reach a probability of 1, and the curve of the fixed patient experiments is the fastest but the curve of varying patient experiment is close. We remark that, in this case, a faster growing CDF implies a higher probability of smaller CV values.

c.2 Kolmogorov-Smirnov test

In statistics, the Kolmogrov-Smirnov test, also known as K-S test, tests null hypothesis for the equality of two probability distributions, or a sample distribution with a reference probability distribution.

We use K-S test to quantify the distance between the empirical CDF of the CV of SLP-A outputs under three patient configurations with significance level . The null hypotheses are rejected when compare the empirical CDF of CV from patient cohort experiment and others in 5% confidence level. Although the CDf curve of fixed patient experiment is faster but the K-S test between it and the CDF of varying patient experiment returns a decision that the null hypothesis is accepted with p value 0.2453, thus they are likely from the same CDF. The p values are believed to be accurate since the sample sizes of the CDFs are large enough, such that , with and are the sample sizes of compared CDFs.

Figure 13: CDFs of Coefficients of variations of the stochastic learner policy outputs under three different patient configuration. The uncertainty of the patient model caused by intra-/inter-patient parameters variations is the largest in a patient cohort and smallest in a fixed virtual patient.

Appendix D Generalization to different diet behaviors

The stochastic learner policy can not only generalize its control ability onto different virtual patients but also to different meal distributions.

d.1 Meal distribution

The attributes of the distinct meal distribution we used to evaluate the generalization ability of our stochastic learner policy is summarized in Table 6. The new distribution has a higher probability of occurrence of snacks and larger snacks with larger carbohydrate amounts. This describes a larger system disturbance to the virtual patient models, and is more challenging for our stochastic learner policy to regulate the patients’ blood glucose level.

breakfast snack 1 lunch snack 2 dinner snack 3
Occurrence 100 80 100 80 100 80
Probability (%)
CHO Amount 40-60 15-30 70-110 15-30 55-75 15-30
(gram)
Time of 3:00- 7:00- 10:00 14:00 17:00 21:00
day (hour) 7:00 10:00 14:00 17:00 21:00 23:00
Table 6: Attributes of a different meal disturbance distribution during testing to that during training. CHO amounts and starting times are sampled uniformly from the reported intervals.

d.2 Performance

We also test our stochastic learner policy in three patient configurations when the patients’ meals are subject to the distribution in Table 6. The results of the performance are shown in Table 7.

Performance
Measures
Fixed Patient
Configuration
Varying Patient
Configuration
Patient Cohort
Configuration
(%) 0.351.08 0.671.83 0.000.00
(%) 91.413.67 90.774.44 85.923.96
(%) 8.243.43 8.563.73 14.083.96
(mg/dL) 218.258.05 217.036.70 233.2210.48
(mg/dL) 76.526.02 76.15+6.70 87.699.30
u (mU/min) 29.260.75 29.320.74 28.700.67
Table 7: The full performance of SLP-A with unseen disturbances.

d.3 Significance tests

We use above described sign test to evaluate, for each configuration, if there are significant performance differences when training-time VS never-before-seen disturbances are used. The two sets of observations are from Table 7 and those from SLP-A in Table 4. The p values of the statistical tests are listed in Table 8.

Performance
Metrics
Fixed Patient
Configuration
Varying Patient
Configuration
Patient Cohort
Configuration
(%) 5.9235e-02 2.2125e-03 1.0000e+00
(%) 2.9921e-01 3.6275e-02 7.6960e-01
(%) 4.5806e-01 2.3040e-01 7.0079e-01
(mg/dL) 3.7602e-01 5.4194e-01 4.5806e-01
(mg/dL) 5.6818e-08 1.4832e-05 8.2858e-01
u (mU/min) 2.2299e-02 2.2299e-02 2.9921e-01
Table 8: P values of sign tests for experiments with SLP-A under trained meal distribution VS unseen meal distribution.

We can observe that all p-values in the patient cohort experiments are very high (the smallest is ), indicating that, as expected, no significant performance difference is observed in patient cohort: SLP-A is able to generalize its control ability to an unseen meal distribution during training phase. In the other two patient configurations, the null hypotheses are also accepted except in the second column and both .