Performance Assessment and Modification of an Adaptive Model

Apr 3, 2019 - Performance Assessment and Modification of an Adaptive Model Predictive Control for Automated Insulin Delivery by a Multivariable Artifi...
1 downloads 0 Views 5MB Size
Subscriber access provided by - Access paid by the | UCSB Libraries

Process Systems Engineering

Performance Assessment and Modication of An Adaptive MPC for Automated Insulin Delivery by A Multivariable Articial Pancreas Iman Hajizadeh, Sediqeh Samadi, Mert Sevil, Mudassir Rashid, and Ali Cinar Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b06202 • Publication Date (Web): 03 Apr 2019 Downloaded from http://pubs.acs.org on April 9, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Performance Assessment and Modication of An Adaptive MPC for Automated Insulin Delivery by A Multivariable Articial Pancreas Iman Hajizadeh,



Sediqeh Samadi,



Mert Sevil,



Mudassir Rashid,



∗,†

and Ali Cinar

†Department of Chemical and Biological Engineering, Illinois Institute of Technology,

Chicago, IL, USA ‡Department of Biomedical Engineering, Illinois Institute of Technology, Chicago, IL, USA E-mail: [email protected]

Phone: +1 (312) 567-3042. Fax: +1 (312) 567-8874

Abstract A controller performance assessment algorithm is developed to analyze the closedloop behavior and modify the parameters of a control system employed in automated insulin delivery. To this end, various performance indices are dened to quantitatively evaluate the controller ecacy in real-time. The controller assessment and modication module also incorporates on-line learning from historical data to anticipate impending disturbances and proactively counteract their eects. A dynamic safety constraint derived from estimates of the physiological states ensures safety of the controlled drug dosing. Using a multivariable simulation platform for type 1 diabetes mellitus, the controller assessment and modication module is applied to the problem of regulating glucose concentrations in people with diabetes by means of automated insulin delivery with an articial pancreas, and the results demonstrate the improvement in controller performance using the performance assessment module. 1

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Introduction Model predictive control (MPC) has gained signicant popularity in diverse applications ranging from chemical and manufacturing process industries to pharmaceutical production and drug delivery. An advantage of MPC formulations is that the nite horizon open-loop optimal control problem solved at each sampling instance can readily handle multivariable systems and systematically deal with state and input constraints 14 . The proliferation of MPC in industrial applications fostered an increased interest in control loop performance assessment and diagnosis 58 . The information gained from the controller assessment and diagnosis routines can be used to retune or redesign the controller for performance improvement, particularly in critical applications like drug delivery. The performance of controllers is typically assessed through comparison with some benchmark. The minimum variance control benchmark, readily computable through simple time series analysis of routine operating data, represents a theoretical lower bound for assessment 913 . A performance index calculated as the ratio of the best achievable variance to the variance of the controlled variable under consideration can indicate the poor performing control loops. Nevertheless, computation of the interactor matrix representing the time delay information of the process is not trivial 1418 . Alternative benchmarks are developed to address the limitations of the minimum variance benchmark. The linear quadratic Gaussian (LQG) or constrained minimum variance benchmark use models of the system and disturbance dynamics to compute the performance curve 1921 . Key performance indexes (KPI) are also developed using the ratio of designed (expected) and achieved (actual) performance value of the MPC cost function 2224 . The expected performance can be a historical period of operation deemed to be satisfactory by control practitioners. A commonality among many control assessment techniques is the introduction of performance benchmarks. MPC performance assessment is challenging because degradation in the closed-loop performance can arise due to model deciencies, poor control design parameters, or inappropriate constraints 2527 . Complex nonlinear dynamical systems such as the metabolic processes 2

ACS Paragon Plus Environment

Page 2 of 44

Page 3 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

in the human body are particularly challenging to control as a variety of disturbances aect the dynamic evolution of the system in uncertain ways 28 . The control of these complex systems requires adaptive models for handling the nonlinear transient dynamics and dynamic constraints that limit input actuation based on the inherent state of the system. Assessing the performance of the resulting adaptive MPC formulations is dicult because of the complications associated with the recursively updated models and dynamic constraints implemented to track the nonlinear transient dynamics and unknown disturbance eects. An example of a system with noteworthy complexities that necessitates adaptive control is the regulation of glucose concentrations in people with type 1 diabetes mellitus (T1DM). T1DM is an autoimmune disease that destroys the insulin-producing beta-cells of the pancreas, and thus people with T1DM must administer exogenous insulin to overcome the insulin deciency and maintain the blood glucose concentration (BGC) within a safe target range 29 . The automation of insulin infusion through technologies that closed the loop between glucose sensing (using continuous glucose monitoring [CGM] sensors) and insulin infusion pumps, termed the articial pancreas (AP) system, is shown to improve glucose control and reduce the likelihood of developing diabetes related complications 3037 . Our novel adaptive and personalized MPC formulation using linear, time-varying glycemic models developed with ecient subspace-based system identication techniques improves the glucose control performance over conventional control strategies. The advantages of dynamically constraining the insulin infusion through the estimated insulin present in the body also help to characterize the protracted eects of previous insulin infusions and avoid dangerous insulin overdosing 38 . Despite the advancements in glycemic control algorithms, the complexity of the control problem in drug delivery has obstructed the reliable assessment of the closed-loop performance. The signicant variability among people with T1DM, and the uctuating temporal dynamics within a subject, have made on-line controller performance assessment a challenging task. Moreover, the uncertain time-varying delays in the glucose-insulin system render the computationally attractive minimum variance benchmarks unreliable. Our previous work 3

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

in assessing the performance of the glycemic control algorithm relied on the LQG performance curve, though the benchmark performance is not practicably achievable 39 . Although KPIs were also developed in our previous work to deduce the potential causes of the performance degradation, the conservative modication of the controller parameters was limited. Motivated by the above considerations, a controller performance assessment system (CPAS) is developed for the adaptive and personalized multivariable AP (mAP) system capable of disturbance rejection without manual user announcements for meals and exercise. A ow chart of the proposed mAP is shown in Fig. 1. Adaptive models identied through system identication techniques are integrated with a physiological compartmental model to characterize the time-varying glucoseinsulin dynamics. The adaptive controller parameters, dynamic insulin constraint, addition of exogenous physiological measurements from wearable devices, feature variables generated from the measurements, and estimates of the uncertain model parameters, including the meal eect, enable the mAP system to eectively compute the optimal insulin infusion over diverse diurnal variations without meal and exercise announcements. A CPAS module is proposed to evaluate the controller ecacy in real-time based on various KPIs and subsequently modify the crucial parameters of the mAP control system to rectify the performance degradation. The CPAS module incorporates on-line learning to identify behaviors and disturbance patterns from the amassed historical data. The identied patterns are then used to proactively mitigate the eects of impending disturbances for improved glucose control. The eectiveness of the proposed CPAS framework to evaluate and modify the adaptive MPC algorithm is demonstrated using a multivariable simulation software platform of virtual subjects with T1DM.

Preliminaries We rst provide brief descriptions of the nonlinear glucoseinsulin dynamic model and the nonlinear observer developed to estimate the physiological states. The estimates of the phys-

4

ACS Paragon Plus Environment

Page 4 of 44

Page 5 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

iological states, including the estimates of the plasma insulin concentration (PIC), are used in dynamically constraining the insulin dosing control algorithm. Then we review the predictor based subspace identication (PBSID) algorithm for the identication of linear, time-varying state-space models employed in the design of the predictive controller. Finally, we provide details on the performance assessment and modication of the recursively identied models.

Adaptive and Personalized PIC Estimator Hovorka's model, a glucoseinsulin dynamics model detailed in Eq. (I.1) in Appendix I of the supporting information, is used for designing the PIC estimator 4042 , which quanties the available insulin in people with T1DM and is used to impose safety limits for the controller. A nonlinear observer based on the unscented Kalman lter (UKF) algorithm for the simultaneous estimation of the state variables and the time-varying parameters is designed to provide accurate PIC estimates (Ik0 ). The proposed PIC estimator can be individualized by appropriately initializing the time-varying parameters using partial least squares regression. To this end, the readily attainable demographic information of individuals, such as weight, body mass index, and duration with diabetes, is used to estimate the initial values for the time-varying parameters in the model 41,42 . After initialization, the time-varying parameters and the state variables can be estimated on-line using an appropriate observer. The proposed PIC estimator is able to capture the inter- and intra-subject variabilities to provide accurate information on the amount of insulin present in the body. Furthermore, as information about the time and quantity of meals is dicult to ascertain and thus considered as unknown disturbances, the gut absorption rate is also included as an extended state to 0 be estimated (UG,k ).

Recursive Subspace-Based System Identication A stable, reliable, and computationally tractable dynamic model is essential for the design of model-based predictive control algorithms in AP systems. To this end, in this work, 5

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 44

we use the recursive predictor-based subspace identication (PBSID) to identify a linear, time-varying model of the form

x˜k+1 = Ak x˜k + Bk u˜k + Kk ek

(1)

yk = Ck x˜k + ek where x ˜k ∈ Rn˜ x denotes the states of the identied model, u˜k ∈ Rn˜ u is the input variable, in this case the input is the estimates of the PIC, meal eects, and measured physiological 0 METk ]), and yk is the output CGM measurement. Further, Ak , Bk , signals (u ˜k := [Ik0 UG,k

and Ck are the system matrices of appropriate dimensions, Kk is the Kalman gain matrix, and ek is the deviation between the measurement and estimate of the CGM output for feedback to correct the state variables. More details about the identication technique used in this work can be found in the Appendix II of supporting information 4346 .

Performance Assessment and Modication of the Recursively Identied Model Assessing the performance of the recursively identied model is necessary to guarantee the model is able to provide accurate output predictions for use in model-based predictive control algorithms. We propose a performance assessment technique operating with a daily time window for the recursively identied model to check necessary KPIs. The KPIs evaluate the risk of under-predicting undesirably high glucose levels (hyperglycemia), risk of over-predicting dangerously low glucose levels (hypoglycemia), and mean-absolute-error between predicted and measured outputs. The assessment of clinical utility of model-based glucose predictions using Clarke error grid analysis has been presented in the literature for dierent purposes 40,4750 . The Clarke error grid analysis assigns pairs of estimatedtrue glucose concentrations into zones A to E. Zone A implies that correct insulin treatment can be initiated using the estimated glucose concentration. Adjacent zones have progressively lower clinical utility with zone E leading to incorrect (opposite in trend) and potentially dangerous insulin treatment. The predictionerror grid analysis (PREDEGA) is a new version of the 6

ACS Paragon Plus Environment

Page 7 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

continuous glucose-error grid analysis (CGEGA), which is shown to be a reliable and robust method for evaluating the accuracy of the predictions in terms of both accurate nominal values and accurate derivative or rate of change of the predictions by means of clinically accepted metrics 5153 . It includes two interacting components: 1) point-EGA (PEGA) for evaluating the prediction accuracy and 2) rate-EGA (REGA) for assessing the capability of the model to characterize the derivative and rate of change of the measurement values. Although the PEGA and REGA may provide distinct information on the predictive performance of the models, it is meaningful to consider these two analyses simultaneously to accurately detect the prediction accuracy as well as the accuracy in the derivatives of the future predictions 54,55 . In case of decient model performance detected by the KPIs, an optimization problem is solved to obtain new parameters of the model such that the KPIs are minimized for the historical data. The parameters that aect the performance of the recursively identied model are p (the length of the past window of data considered when predicting the future outputs), forgetting factors λ and the initial values of covariance matrices in the RLS problems of Eq. (II.4)-Eq. (II.5) in the Appendix.

Adaptive PIC Cognizant MPC Algorithm In this section, the insulin compartment model of Eq. (I.1) is incorporated with the recursive PBSID approach to build an appropriate form of the glycemic model for use in the AP control system. Then, we describe the adaptive controller setpoint denition, glycemic and PIC risk indexes used in the MPC controller objective function. The safety constraints based on the PIC and a feature extraction method for manipulating these constraints during the meal consumption periods are presented. Furthermore, the performance assessment of the controller during fasting periods and the feature extraction technique are discussed. Subsequently, the adaptive mAP control formulation is presented.

7

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 44

Integrating Insulin Compartment Models with Subspace Identication The manipulated variable of the AP control system is the injected insulin in form of basal/bolus insulin. So, the identied state space model using PBSID described by Eq. (1) cannot be used directly for the model-based control techniques since the PIC is an input of the system 56 . In this subsection, we integrate the data-driven model determined from the PBSID technique with the insulin subsystem from Eq. (I.1) described by Eq. (I.1a)-Eq. (I.1c) yielding

x¯k+1 = Dk x¯k + Ek u¯k

(2)

y¯k = Fk x¯k 

with the new state vector of x ¯k = x˜Tk+3+d S1,k S2,k+1 Ik+2

T

and the output as y¯k =

yk+3+d . Further, u¯k includes delayed infused basal/bolus insulin information by a delay of order d, estimates of meal eects, and measured physiological signals. All other parameters of Eq. (2) are detailed in the Appendix III of supporting information.

Adaptive Controller Set-point Denition In this work, an adaptive learning technique is proposed to modify the controller set-point based on historical information during meal and exercise periods. Historical data, including the CGM measurements and physiological variables, are used to deduce the probable times of meals consumption and physical activities. This valuable information on daily life behaviors and habits can be used in the mAP system to anticipate and proactively mitigate the eects of disturbances. Therefore, the performance of the control system can be improved if the controller set-point is appropriately modied in advance for the anticipated periods of the disturbance eects. The nominal value of the reference glucose set-point target in the controller is r = 110 mg/dL. Using historical data, we can assign dierent values to 8

ACS Paragon Plus Environment

Page 9 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

the controller set-point around specic times that disturbances are expected to aect the glucose values. We propose a feature extraction method that is used for meal detection, while the physiological variables are employed to detect the presence of physical activity. At the end of each day, the time windows detected for meal disturbances are assigned a lower set-point value of r = 80 mg/dL to make the controller more aggressive around those periods. Analogously, during exercise times, a higher value set-point value of r = 160 mg/dL makes the controller more conservative to reduce insulin infusion in response to possible physical activities. The controller set-point values for the anticipated meal and exercise disturbance time windows of the new day are calculated by averaging the previously assigned values. For example, for one subject who usually eats breakfast between [06:00, 08:00], lunch in the time interval [12:00, 14:00], dinner around [18:00, 20:00], and does two exercise sessions around [09:30, 10:30] and [15:30, 16:30], the new modied controller set-points are obtained as shown in Fig. 2.

Adaptive Glycemic Risk Index An adaptive glycemic risk index (GRI) is used to determine the weighting matrix for penalizing the deviations of the outputs from their nominal set-point 38,57 . To this end, the time-varying positive semi-denite weighting matrix Qk := Q (¯ yk ) is dened as Q (¯ yk ) :=

α (¯ yk ) · Qnom where Qnom denotes a nominal weight and     phypo × y¯k2 + phypo × y¯k + phypo  I II III     0 α (¯ yk ) :=   phyper × y¯k2 + phyper × y¯k + phyper  I II III      1

if y¯k

  ∈ 50 + rsp , rk

if y¯k

= rk

if y¯k

∈ rk , 300 + rsp



(3)

otherwise

where rk is the controller set-point at sampling time k and rsp = rk − 110. The glycemic risk index asymmetrically increases the set-point tracking weight when y¯k diverges from the

9

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 44

target range. Since hypoglycemic events have serious short-term implications, the set-point penalty increases rapidly in response to hypoglycemic excursions and more gradually in hyperglycemic excursions. A plot of the glycemic risk index is given in Fig. 3.

Plasma Insulin Risk Index A plasma insulin risk index (PIRI) is dened to manipulate the weighting matrix for penalizing the amount of input actuation (aggressiveness of insulin dosing) depending on the estimated PIC, thus suppressing the infusion rate if sucient insulin is present in the bloodstream 38,57 . To this end, the time-varying positive denite weighting matrix Rk := R (γk ) is developed from the PIRI as R (γk ) := γk · Rnom , with Rnom as a nominal weight and γk dened as

γk :=

Ik0

!2

Ibasal,k

(4)

where

Ibasal,k :=

udb,k VI · ke,k

(5)

and udb is the patient specic (possibly time-varying) basal insulin rate that is known in practice, and VI and ke are parameters of Hovorka's model. Furthermore, the parameter ke is estimated on-line using the UKF and the CGM output measurements. As it is impractical to directly consider the estimates of the PIC to dene parameters of the MPC due to the variability among subjects, the normalized value of the PIC is employed, which eliminates the dependency of the PIC estimates to a particular subject by standardizing with the known patient-specic basal PIC value. A plot of the plasma insulin risk index (Fig. 4) indicates that as the penalty weight on the input action increases, and dosing becomes less aggressive, if the estimated PIC in high.

10

ACS Paragon Plus Environment

Page 11 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Plasma Insulin Concentration Bounds In the proposed MPC, the estimated future PIC is dynamically bounded depending on the value of the CGM measurements. For instance, if the CGM values are elevated, the bounds on the PIC are increased to ensure sucient insulin is administered to regulate the glucose concentration. Furthermore, the PIC bounds also constrain the search space in the optimization problem, thus improving the computational tractability of the proposed MPC. The PIC bounds are determined based on the CGM measurements as XPIC := (Pfasting + Pmeal )X (¯ yk ), where XPIC denes the lower and upper bounds and a desired target for the normalized PIC through the predicted CGM y¯k . Pmeal is a parameter which modies the PIC bounds when there is a rapid increase in the CGM measurements. Pfasting is a patient specic parameter that denes the controller aggressiveness/conservativeness during the fasting period. Fig. 5 depicts the bounds and the reference target for the normalized PIC as a function of the CGM measurement, and the MPC solution should satisfy the PIC constraints while maintaining the PIC close to the desired value. The nominal PIC bounds can be determined by multiplying the normalized PIC bounds with the basal PIC value (Ibasal ). Therefore, appropriate PIC bounds can be determined based on each subject's basal PIC value and the CGM measurement.

Assessment of the Controller During Fasting Period The Pfasting parameter plays an important role in the aggressiveness/conservativeness of the controller. The initial value of the Pfasting is set to a safe initial value of unity. The Pfasting parameter is modied for each subject by analyzing the past CGM values and injected insulin information during fasting periods. For example, if during fasting periods the BGC remains below the controller set-point due to over-delivery of insulin, it indicates that the controller is overly aggressive in insulin dosing and the value of Pfasting needs to be reduced. Finding an appropriate value for this parameter can also be helpful for the case where the patient specic basal insulin rate is high, thus increasing the chance of hypoglycemia. Moreover, this 11

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 44

parameter should be increased for a patient that has low patientspecic basal insulin rates such that the CGM values are typically higher than the desired set-point value. Therefore, by nding an appropriate value for the Pfasting , the safety and eciency of the mAP can be improved. In this work, we consider a 3 hours of postprandial period for each detected meal and 2 hours of recovery period after each detected exercise. The rest of the time is considered as fasting period. The collected data during this period is used to nd the appropriate value of the Pfasting parameter. At the end of each day, the value of αf (¯ yk ) parameter is obtained based on the CGM values in the fasting period:

   −1       −α(¯ yk )    αf (¯ yk ) := 0      α(¯ yk )       1

  if y¯k ∈ 40, 50 + rsp   if y¯k ∈ 50 + rsp , rk (6)

if y¯k = rk if y¯k ∈ rk , 300 + rsp



 if y¯k ∈ 300 + rsp , 400

where the α (¯ yk ) is a parameter used in the adaptive glycemic risk index (GRI). At the end new ) using the following of each day, the value of Pfasting is modied for the next day (Pfasting

equation: new old Pfasting = Pfasting +α ¯ f (¯ yk )

(7)

where α ¯ f (¯ yk ) is the average of all αf (¯ yk ) values calculated using CGM data during the fasting period. For example, if the average value gets a negative value, the controller needs to be more conservative by reducing the PIC bounds for the upper, lower and desired values. Similar to the adaptive determination of controller set-point and glycemic risk index, the

X (¯ yk ) features (slope and intercept of the lines for the upper, lower, and desired PIC bounds) are modied over time by learning from historical data to be more aggressive during meal times and more conservative during exercise periods.

12

ACS Paragon Plus Environment

Page 13 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Feature Extraction for Manipulating Constraints Meal consumption can be automatically detected using qualitative descriptions of glucose time-series data, which is useful in modifying the aggressiveness of the adaptive MPC 58,59 . In this work, features are generated from the data to describe the recent trajectory of the  glycemic measurements. To this end, a p-order polynomial yi = f ti , θkM with parameters

θkM is tted to the most recent l glucose measurements yi:i−l := [yi yi−1 . . . yi−l ] at each sampling time using ordinary least squares where ti denotes the sampling index of the recent measurements. Then the derivatives of the polynomial are obtained and the rst- and secondorder derivatives, denoted f (1) and f (2) , are analyzed to derive parameter Pmeal for detecting carbohydrate consumption as

Pmeal

 (1)  f   m   1  c(1) = f   cm 2     0

(2) ≥0 if f (1) ≥ cm 0 and f (2) if f (1) ≥ cm 0, Rk := R (γk ) is a strictly positive denite symmetric matrix to penalize the manipulated input variables. At each iteration, the quadratic programming problem in Eq. (10) is solved, and uk := v0 is the optimal solution implemented to infuse insulin over the current sampling interval with the MPC computation repeated at subsequent sampling instances using new CGM measurements, updated states, and newly calculated penalty weights of the objective function.

Results The ecacy of the proposed mAP is demonstrated by using a multivariable simulator developed by our research group at Illinois Institute of Technology based on a modied Hovorka's glucose-insulin dynamics model that takes into account the eects of dierent physical activities 60 . In addition to the CGM values, the multivariable simulator generates physiological variable signals reported by noninvasive wearable devices. For this purpose, aerobic exercises with treadmill and bicycle are considered for testing the mAP system. Twenty virtual subjects are simulated for multiple days with varying times and quantities of meals consumed on each day and physical activities with dierent types, intensities, and durations as detailed in Table 1 and Table 2. The meal and physical activity information are not entered manually to the AP system as the AP controller is designed to regulate the BGC in the presence of signicant unknown disturbances such as unannounced meals and exercises. The Metabolic Equivalent (MET) values computed by the simulator are used as physiological signals in the mAP system. In the MPC formulation, the umin is zero and umax is dened as

umax

103 = udb + uBolus,max · ∆t

(11)

where ∆t is 5 minutes sampling time and the maximum allowable bolus in units of U is

16

ACS Paragon Plus Environment

Page 17 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Table 1: Meal scenario for 3 weeks (21 days) closed-loop experiment using the multivariable metabolic simulator Range for values Time Amount (g)

Meal Breakfast Lunch Dinner Snack

[06:00, 07:00] [12:00, 13:00] [18:00, 19:00] [22:00, 23:00]

[40, 80] [40, 80] [40, 80] [20, 30]

Table 2: Exercise scenario for 3 weeks (21 days) closed-loop experiment using the multivariable metabolic simulator Exercise Treadmill Bicycling

Range for values Duration (min) Speed

Time

[09:00, 10:00] or [15:00, 16:00] [09:00, 10:00] or [15:00, 16:00]

[50, 70] [50, 70]

[4, 6] -

Incline

Power

[0, 2] -

[70, 110]

deed as

 uBolus,max := max

 y¯k − rk Ik0 · VI − ,0 . CFk 103

(12)

where CFk is the patient-specic correction factor, Ik0 is the estimated PIC in units of mU/L,

VI is the insulin distribution volume in units of L, and therefore, Ik0 · VI provides an estimate of the total insulin present in the body in units of mU . The safety insulin constraints dened for the desired, upper and lower bounds are based on the normalized PIC values that are calculated by dividing the PIC by the basal PIC. The patient-specic basal PIC is the value of the PIC if only basal insulin is administered. The MPC prediction horizon is nP = 24 sampling times or 2 hours. In order to avoid any potential hypoglycemia, a predictive hypoglycemia alarm module is designed to suggest carbohydrates to be consumed for the safety of the subject. The details of this module can be found in the Appendix IV of supporting information. The purpose of the CPAS is to minimize the deviation of outputs from controller setpoint and cost of insulin injection. To show the signicance of the proposed CPAS, the value of the following terms are calculated based on the real outputs (measured CGM values) and 17

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 44

injected insulin for each day using: Yterm

nd 1 X := Qi (yi − ri ) nd i=1

Uterm

nd 1 X := Ri (vi ) nd i=1

(13)

where yi , ri , and vi are the measured CGM values, set-point, and injected insulin, respectively. nd is the total number of samples in each day. We consider the Qi and Ri parameters in calculating the Yterm and Uterm terms to give meaningful weights based on the situation of the patient. As an example, for 20 mg/dl deviation from the controller set-point, if the CGM value is lower than the controller set-point, the weight is higher due to the risk of hypoglycemia.

Case 1 In this case, the quantitative evaluation of the closed-loop results are shown based on the proposed algorithms in Table V.2 and Table V.3 of Appendix V, where it is assumed that the meal and exercise specications are same every day. The purpose of these simulation results is to show the capability of the mAP in learning patients' habits and activities over time and also signicant improvement in the closed-loop performance due to an ecient CPAS. It is readily observed that the average percentage of time spent in the target ranges (BGC ∈ [70, 140] mg/dL and BGC ∈ [70, 180] mg/dL) improves signicantly for all subjects over time. There is no hypoglycemia event as the BGC is never below 70 mg/dL. The minimum and maximum observed BGC values across all experiments are 75 and 231 ml/dL, respectively. There is also a signicant improvement in mean, standard deviation (SD) and coecient of variation (CV) of CGM values. Overall, the results demonstrate that the proposed mAP is able to regulate BGC eectively without requiring meal and exercise announcement. The closed-loop results for all subjects are shown in Fig. 6 and Fig. 7. We can see that the performance of the mAP improves over time due to the CPAS and by learning from historical data.

Case 2 In this case, more realistic scenarios for meal consumption and exercise are con18

ACS Paragon Plus Environment

Page 19 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

sidered. The meal and exercise specications are randomly selected from the specied ranges mentioned in Table 1 and Table 2. The quantitative evaluation of the closed-loop results based on the proposed algorithms in Table V.4 and Table V.5, where the meal and exercise specications are dierent every day. The purpose of these simulation results is to show that the mAP integrated with the CPAS and adaptive learning technique is robust and reliable in learning patients' habits and activities that vary over time, and also handling them. In this case, the average percentage of time spent in the target ranges (BGC ∈ [70, 140] mg/dL and BGC ∈ [70, 180] mg/dL) improves signicantly for all subjects over time as well. There is no hypoglycemia event as the BGC is never below 70 mg/dL. The minimum and maximum observed BGC values across all experiments are 84 and 224 ml/dL, respectively. There is also a signicant improvement in mean, standard deviation (SD) and coecient of variation (CV) of CGM values. Overall, the results demonstrate that the proposed mAP is able to regulate BGC eectively in presence of signicant disturbances caused by the diverse timing and amounts of meals and exercise specications, while mitigating severe hypo- and hyperglycemic excursions. The closed-loop results for all subjects in Fig. 8 and Fig. 9 indicate that the performance of the mAP improves over time due to the CPAS and by learning from historical data.

Case 3 In this case, more challenging scenarios for meal consumption and exercise are considered. The meal and exercise specications are randomly selected from the specied ranges mentioned in the Table 1 and Table 2 but to challenge the mAP control system, in the last day of simulation (day 21), the meals consumption and physical activities are dened according to the Table 3 and Table 4. As listed in these tables, the times of meal consumption and exercise are delayed which is in contrast to the historical data. The quantitative evaluation of the closed-loop results based on the proposed algorithms for this case are presented in Table V.6 and Table V.7. The purpose of these simulation results is to show that the mAP integrated with the CPAS and adaptive learning technique is robust and reliable in facing very dierent conditions from what the patient does regularly 19

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 44

Table 3: Meal scenario for day 21 Meal

values Time Amount (g)

Breakfast Lunch Dinner Snack

[09:00] [14:00] [19:00] [23:00]

[60] [50] [55] [20]

Table 4: Exercise scenario for day 21 Exercise Treadmill Bicycling

values Duration (min) Speed

Time

[17:00] [12:00]

[61] [50]

[5] -

Incline

Power

[1] -

[90]

every day. In this case, there is no hypoglycemia event as well, as the BGC is never below

70 mg/dL. The minimum and maximum observed BGC values across all experiments are 73 and 278 ml/dL, respectively. The closed-loop results for all subjects are shown in Fig. 10. The performance of the mAP in presence of signicant changes in patient's daily behavior is safe and reliable without causing any hypo or hyperglycemia events.

Case 4

In this case, the performance of the controller is tested in presence of more

atypical days. To this end, We dened a longer scenario of 60 days with a signicant number of atypical days to evaluate the performance of controller against more frequent non-regular days. In this case, 28 out of 60 days are considered to be atypical. These atypical days have been dened randomly to challenge thoroughly the mAP system. For example, the times of meal consumption and exercise are delayed, the order of meal consumption and exercise is changed, the meal and exercise specications are randomly selected, and the frequency of exercise sessions and meal consumptions are varied in each day. According to the results shown in the Fig. 11 and Fig. 12, although the performance of the mAP controller deviates from its normal operation during atypical days, it still ensures reliability and safety against unexpected behaviors. As the number of atypical days increases, the learning process and performance of the controller would be adversely aected. In our future work, we will enhance 20

ACS Paragon Plus Environment

Page 21 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

the proposed learning methods to automatically detect atypical days and consequently avoid modication of key controller parameters using non-regular historical data. As an example, the modication of the Pfasting parameter by the CPAS and using historical data is shown in Fig. 13. The Pfasting starts from an initial value and changes to its optimum one. The adaptive modied controller set-point is also shown (Fig. 14). By using historical data over time, the proposed algorithm learns the habits and behaviors of the patients and consequently denes more realistic set-points for the controller. The controller set-point modication is done automatically using only collected historical data without requiring any information from a patient. To show the performance of the recursively identied models, based on one-step-ahead (5 minutes) and 4-steps-ahead (20 minutes) prediction of CGM, the prediction-error grid analysis (PRED-EGA) is done to quantify the accuracy of glucose predictions generated by the proposed technique as compared to the measured values. The result of this test for all 20 subjects is shown in Fig. 15. The predictions are classied as accurate if the predictions fall within the A or B zones of both the P-EGA and the R-EGA. The benign errors refer to predictions with acceptable point accuracy (A or B zones in P-EGA) yet with more signicant errors in the rate accuracy (C, D, or E zones in R-EGA). The benign errors are unlikely to lead to negative clinical consequences as the blood glucose predictions are relatively close to the measured values with only the derivative or rate of change of the blood glucose values not predicted accurately. The erroneous predictions (C, D, or E zones of the P-EGA regardless of the corresponding R-EGA zone) refers to the samples that have signicant prediction errors in the nominal blood glucose values, which could potentially lead to negative clinical actions and outcomes. As you can see in the Fig. 15, for foursteps (20 min) aheadCGMprediction, most predicted CGM values are classied as accurate (99.93 %) and just a few of them are in the benign error region (0.07 %). These results show the reliability and ecacy of the identied models to be used in the mAP system. 21

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Due to the very complex nature of BGC dynamics in people with T1DM, the performance of the MPC techniques developed for the AP system may not be always be as expected. Performance degradation may occur due to the time-varying characteristics of a system and its nonlinear behavior, presence of stochastic and unknown disturbances, and uncertain timevarying delays. As the human body is a system that has such characteristics, it is challenging to automate the monitoring of the AP performance, diagnosing the root causes of poor performance, and modifying the relevant AP modules. In the literature, various theoretical techniques have been proposed for designing control performance assessment and modication systems. In the majority of these techniques, it is assumed that the dynamic model of underlying system including the disturbance model and time delays are perfectly known. These methods also often utilize idealized criterion with which the current control performance is compared. For the AP systems, neither the available physiological models nor the datadriven techniques can describe the dynamic behavior of BGC variations in human body perfectly. The eects of major disturbances to the BGC like physical activities, meals, and stress are not well established yet. The infused insulin as the manipulated variable aects the BGC with an uncertain timevarying delay. Considering all these limitations, the answer to the question of "Who will control the AP controller itself?" is not easy. In this work, to design the CPAS, we were inspired by the theoretical concepts developed in the context of CPAS in the literature, and we also took advantage of historical clinical data of people with T1DM to dene realistic performance indices to assess the AP performance in a medically-acceptable way under dierent conditions such as postprandial periods or exercise times.

Conclusions In this work, a CPAS module is proposed for a multivariable AP controlled by MPC. The mAP eciency is evaluated online, based on dierent key performance indexes (KPI). The

22

ACS Paragon Plus Environment

Page 22 of 44

Page 23 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

mAP modules are capable to learn over time the patient's habits and activities using historical data such as CGM data and biometric variables. The recursively identied models provide accurate prediction of the output as well as any potential danger for hypoglycemia and hyperglycemia in advance. A feature extraction method automatically detects meal consumption using qualitative descriptions of glucose time-series data, which provides useful information for modifying the aggressiveness of the adaptive MPC. Safety constraints for the controller are dened based on plasma insulin concentration to guarantee the ecacy of the mAP. In order to avoid any potential hypoglycemia, a predictive hypoglycemia alarm module is designed to suggest carbohydrates to be consumed by patient for the safety. The proposed adaptive learning mAP system integrated with the CPAS is evaluated by using a multivariable simulator including 20 virtual subjects with T1DM. Dierent scenarios for meal consumption and exercise were dened to signicantly challenge the AP system and assess its performance. These validations are necessary before going into full practical application of an AP controller with real subjects. The next step is to test oline the proposed techniques in clinical experiments to further evaluate the insulin suggestions. After the successful performance of this step, the proposed mAP will be utilized in clinical experiments to automate insulin delivery in real-time.

Acknowledgement This work is supported by the National Institute of Diabetes and Digestive and Kidney Diseases grants DP3 DK101075-01 and DP3 DK101077-01 and Juvenile Diabetes Research Foundation grant A18-0036-001 made possible through collaboration between the JDRF and The Leona M. and Harry B. Helmsley Charitable Trust.

Supporting Information Available • Appendix I: Hovorka's glucoseinsulin dynamics model. 23

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

• Appendix II: Predictor-based system identication technique. • Appendix III: Derivations for calculating the parameters of integrated glycemic model. • Appendix IV: Hypoglycemia detection and carbohydrate suggestion. • Appendix V: CGM glucose metrics for closed-loop results.

References (1) Garcia, C. E.; Prett, D. M.; Morari, M. Model predictive control: Theory and practice A survey. Automatica

1989, 25, 335  348.

(2) Rawlings, J. B.; Mayne, D. Q. Model predictive control: Theory and design ; Nob Hill Pub. Madison, Wisconsin, 2009. (3) Zavala, V. M.; Biegler, L. T. The advanced-step NMPC controller: Optimality, stability and robustness. Automatica

2009, 45, 86  93.

(4) Mayne, D. Q. Model predictive control: Recent developments and future promise. Au-

tomatica 2014, 50, 2967  2986. (5) Ko, B. S.; Edgar, T. F. Performance assessment of constrained model predictive control systems. AIChE Journal

2001, 47, 13631371.

(6) Argawal, N.; Huang, B.; Tamayo, E. Assessing model predictive control (MPC) performance. Ind. Eng. Chem. Res

2007, 46, 81018111.

(7) Yu, J.; Qin, S. J. Statistical MIMO controller performance monitoring. Part I: Datadriven covariance benchmark. J. of Process Control

2008, 18, 277296.

(8) Harrison, C. A.; Qin, S. J. Minimum variance performance map for constrained model predictive control. J. of Process Control

2009, 19, 11991204.

24

ACS Paragon Plus Environment

Page 24 of 44

Page 25 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(9) Årström, K. Assessment of achievable performance of simple feedback loops. Int. J.

Adapt. Control Signal Process. 1991, 5, 319. (10) Åström, K. J. Introduction to stochastic control theory ; Courier Corporation, 2012. (11) Harris, T. Assessment of control loop performance. Can. J. Chem. Eng.

1989, 67,

856861. (12) Harris, T.; Seppala, C.; Desborough, L. A review of performance monitoring and assessment techniques for univariate and multivariate control systems. J. Process Control

1999, 9, 1  17. (13) Desborough, L.; Harris, T. Performance assessment measures for univariate feedback control. Can. J. Chem. Eng.

1992, 70, 11861197.

(14) Huang, B.; Shah, S.; Kwok, E. Good, bad or optimal? Performance assessment of multivariable processes. Automatica

1997, 33, 11751183.

(15) Huang, B.; Shah, S. L.; Fujii, H. The unitary interactor matrix and its estimation using closed-loop data. J. Process Control

1997, 7, 195  207.

(16) Harris, T. J.; Boudreau, F.; MacGregor, J. F. Performance assessment of multivariable feedback controllers. Automatica

1996, 32, 15051518.

(17) McNabb, C. A.; Qin, S. J. Projection based MIMO control performance monitoring: I-covariance monitoring in state space. J. Process Control

2003, 13, 739757.

(18) Huang, B.; Ding, S. X.; Thornhill, N. Practical solutions to multivariate feedback control performance assessment problem: reduced a priori knowledge of interactor matrices.

J. Process Control 2005, 15, 573583. (19) Patwardhan, R. S.; Shah, S. L. Issues in performance diagnostics of model-based controllers. J. Process Control

2002, 12, 413427. 25

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 44

(20) Huang, B.; Shah, S. L. Performance assessment of control loops: theory and applica-

tions ; Springer Science & Business Media, 2012. (21) Zagrobelny, M.; Ji, L.; Rawlings, J. B. Quis custodiet ipsos custodes? Annu. Rev.

Control 2013, 37, 260270. (22) Schäfer, J.; Cinar, A. Multivariable MPC system performance assessment, monitoring, and diagnosis. J. Process Control

2004, 14, 113129.

(23) Patwardhan, R. S.; Shah, S. L.; Qi, K. Z. Assessing the performance of model predictive controllers. Can. J. Chem. Eng.

2002, 80, 954966.

(24) Julien, R. H.; Foley, M. W.; Cluett, W. R. Performance assessment using a model predictive control benchmark. J. Process Control

2004, 14, 441  456.

(25) Sun, Z.; Qin, S. J.; Singhal, A.; Megan, L. Performance monitoring of model-predictive controllers via model residual assessment. J. Process Control

2013, 23, 473482.

(26) Zhao, Y.; Chu, J.; Su, H.; Huang, B. Multi-step prediction error approach for controller performance monitoring. Control Eng. Pract.

2010, 18, 112.

(27) Pannocchia, G.; De Luca, A.; Bottai, M. Prediction Error Based Performance Monitoring, Degradation Diagnosis and Remedies in Oset-Free MPC: Theory and Applications. Asian J. Control

2014, 16, 9951005.

(28) Doyle III, F. J.; Bequette, B. W.; Middleton, R.; Ogunnaike, B.; Paden, B.; Parker, R. S.; Vidyasagar, M. Control in biological systems. In T. Samad, & A. M.

Annaswamy (Eds.), The impact of control technology, part 1. IEEE Control Systems Society 2011, (29) Cinar, A.; Turksoy, K. Advances in Articial Pancreas Systems: Adaptive and Multi-

variable Predictive Control ; Springer, 2018.

26

ACS Paragon Plus Environment

Page 27 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(30) Gondhalekar, R.; Dassau, E.; Doyle, F. J. Periodic zone-MPC with asymmetric costs for outpatient-ready safety of an articial pancreas to treat type 1 diabetes. Automatica

2016, 71, 237  246. (31) Chakrabarty, A.; Zavitsanou, S.; Doyle, F. J.; Dassau, E. Event-triggered model predictive control for embedded articial pancreas systems. IEEE Trans. Biomed. Eng.

2018, 65, 575586. (32) Cameron, F.; Niemeyer, G.; Bequette, B. W. Extended multiple model prediction with application to blood glucose regulation. J. of Process Control

2012, 22, 1422  1432.

(33) Boiroux, D.; Duun-Henriksen, A. K.; Schmidt, S.; Nørgaard, K.; Madsbad, S.; Poulsen, N. K.; Madsen, H.; Jørgensen, J. B. Overnight glucose control in people with type 1 diabetes. Biomed. Signal Process. Control

2018, 39, 503512.

(34) El Fathi, A.; Smaoui, M. R.; Gingras, V.; Boulet, B.; Haidar, A. The Articial Pancreas and Meal Control: An Overview of Postprandial Glucose Regulation in Type 1 Diabetes.

IEEE Control Syst. Mag. 2018, 38, 6785. (35) Messori, M.; Incremona, G. P.; Cobelli, C.; Magni, L. Individualized model predictive control for the articial pancreas: In silico evaluation of closed-loop glucose control.

IEEE Control Syst. Mag. 2018, 38, 86104. (36) Shi, D.; Dassau, E.; Doyle III, F. J. Adaptive Zone Model Predictive Control of Articial Pancreas Based on Glucose-and Velocity-Dependent Control Penalties. IEEE Trans.

Biomed. Eng. 2018, (37) Incremona, G. P.; Messori, M.; Toanin, C.; Cobelli, C.; Magni, L. Model predictive control with integral action for articial pancreas. Control. Eng. Pract. 94.

27

ACS Paragon Plus Environment

2018, 77, 86 

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(38) Hajizadeh, I.; Rashid, M.; Sevil, M.; Brandt, R.; Samadi, S.; Hobbs, N.; Cinar, A. Adaptive Model Predictive Control for Nonlinearity in Biomedical Applications. IFAC-

PapersOnLine 2018, 51, 368  373. (39) Feng, J.; Hajizadeh, I.; Yu, X.; Rashid, M.; Turksoy, K.; Samadi, S.; Sevil, M.; Hobbs, N.; Brandt, R.; Lazaro, C.; Maloney, Z.; Littlejohn, E.; Philipson, L. H.; Cinar, A. Multi-level supervision and modication of articial pancreas control system. Comput. Chem. Eng.

2018, 112, 57  69.

(40) Hovorka, R.; Canonico, V.; Chassin, L. J.; Haueter, U.; Massi-Benedetti, M.; Federici, M. O.; Pieber, T. R.; Schaller, H. C.; Schaupp, L.; Vering, T.; Wilinska, M. E. Nonlinear model predictive control of glucose concentration in subjects with type 1 diabetes. Physiol. Meas.

2004, 25, 905.

(41) Hajizadeh, I.; Rashid, M.; Turksoy, K.; Samadi, S.; Feng, J.; Frantz, N.; Sevil, M.; Cengiz, E.; Cinar, A. Plasma Insulin Estimation in People with Type 1 Diabetes Mellitus.

Ind. Eng. Chem. Res. 2017, 56, 98469857. (42) Hajizadeh, I.; Rashid, M.; Samadi, S.; Feng, J.; Sevil, M.; Hobbs, N.; Lazaro, C.; Maloney, Z.; Brandt, R.; Yu, X.; Turksoy, K.; Littlejohn, E.; Cengiz, E.; Cinar, A. Adaptive and Personalized Plasma Insulin Concentration Estimation for Articial Pancreas Systems. J. Diabetes Sci. Technol.

2018, 12, 639649.

(43) Houtzager, I.; van Wingerden, J.-W.; Verhaegen, M. Recursive Predictor-Based Subspace Identication With Application to the Real-Time Closed-Loop Tracking of Flutter. IEEE Trans. Control Syst. Technol.

2012, 20, 934949.

(44) Hajizadeh, I.; Rashid, M.; Turksoy, K.; Samadi, S.; Feng, J.; Sevili, M.; Frantz, N.; Lazaro, C.; Maloney, Z.; Littlejohn, E.; Cinar, A. Multivariable Recursive Subspace Identication with Application to Articial Pancreas Systems. IFAC-PapersOnLine

2017, 909  914. 28

ACS Paragon Plus Environment

Page 28 of 44

Page 29 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(45) Hajizadeh, I.; Rashid, M.; Cinar, A. Ensuring Stability and Fidelity of Recursively Identied Control-Relevant Models. The 18th IFAC Symposium on System Identication (SYSID). 2018; pp 927932. (46) Hajizadeh, I.; Rashid, M.; Turksoy, K.; Samadi, S.; Feng, J.; Sevil, M.; Hobbs, N.; Lazaro, C.; Maloney, Z.; Littlejohn, E.; Cinar, A. Incorporating Unannounced Meals and Exercise in Adaptive Learning of Personalized Models for Multivariable Articial Pancreas Systems. J. Diabetes Sci. Technol.

2018, 12, 953966.

(47) Del Favero, S.; Facchinetti, A.; Cobelli, C. A glucose-specic metric to assess predictors and identify models. IEEE Trans. Biomed. Eng.

2012, 59, 12811290.

(48) Reifman, J.; Rajaraman, S.; Gribok, A.; Ward, W. K. Predictive monitoring for improved management of glucose levels. J. Diabetes Sci. Technol.

2007, 1, 478486.

(49) Contreras, I.; Oviedo, S.; Vettoretti, M.; Visentin, R.; Vehí, J. Personalized blood glucose prediction: A hybrid approach using grammatical evolution and physiological models. PloS one

2017, 12, e0187754.

(50) Pappada, S. M.; Cameron, B. D.; Tulman, D. B.; Bourey, R. E.; Borst, M. J.; Olorunto, W.; Bergese, S. D.; Evans, D. C.; Stawicki, S. P.; Papadimos, T. J. Evaluation of a model for glycemic prediction in critically ill surgical patients. PloS one

2013,

8, e69475. (51) Kovatchev, B. P.; Gonder-Frederick, L. A.; Cox, D. J.; Clarke, W. L. Evaluating the accuracy of continuous glucose-monitoring sensors: continuous glucoseerror grid analysis illustrated by TheraSense Freestyle Navigator data. Diabetes Care

2004, 27, 19221928.

(52) Clarke, W. L.; Anderson, S.; Farhy, L.; Breton, M.; Gonder-Frederick, L.; Cox, D.; Kovatchev, B. Evaluating the clinical accuracy of two continuous glucose sensors using Continuous glucoseerror grid analysis. Diabetes care

29

ACS Paragon Plus Environment

2005, 28, 24122417.

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 44

(53) Sivananthan, S.; Naumova, V.; Man, C. D.; Facchinetti, A.; Renard, E.; Cobelli, C.; Pereverzyev, S. V. Assessment of blood glucose predictors: the prediction-error grid analysis. Diabetes Technol. Ther.

2011, 13, 787796.

(54) Yu, X.; Rashid, M.; Feng, J.; Hobbs, N.; Hajizadeh, I.; Samadi, S.; Sevil, M.; Lazaro, C.; Maloney, Z.; Littlejohn, E.; Quinn, L.; Cinar, A. Online Glucose Prediction Using Computationally Ecient Sparse Kernel Filtering Algorithms in Type-1 Diabetes. IEEE

Trans. Control Syst. Technol. 2018, 113. (55) Yu, X.; Turksoy, K.; Rashid, M.; Feng, J.; Hobbs, N.; Hajizadeh, I.; Samadi, S.; Sevil, M.; Lazaro, C.; Maloney, Z.; Littlejohn, E.; Quinn, L.; Cinar, A. Model-fusionbased online glucose concentration predictions in people with type 1 diabetes. Control

Eng Pract. 2018, 71, 129  141. (56) Hajizadeh, I.; Rashid, M.; Cinar, A. Integrating Compartment Models with Recursive System Identication. American Control Conference (ACC). 2018; pp 35833588. (57) Rashid, M.; Hajizadeh, I.; Cinar, A. Plasma Insulin Cognizant Predictive Control for Articial Pancreas. American Control Conference (ACC). 2018; pp 35893594. (58) Samadi, S.; Turksoy, K.; Hajizadeh, I.; Feng, J.; Sevil, M.; Cinar, A. Meal detection and carbohydrate estimation using continuous glucose sensor data. IEEE J. Biomed.

Health Inform. 2017, 21, 619627. (59) Samadi, S.; Rashid, M.; Turksoy, K.; Feng, J.; Hajizadeh, I.; Hobbs, N.; Lazaro, C.; Sevil, M.; Littlejohn, E.; Cinar, A. Automatic detection and estimation of unannounced meals for multivariable articial pancreas system. Diabetes Technol. Ther.

2018, 20,

235246. (60) Samadi, S.; Rashid, M.; Sevil, M.; Hajizadeh, I.; Hobbs, N.; Kolodziej, P.; Feng, J.; Park, M.; Quinn, L.; Cinar, A. Multivariable Simulation Platform For Type 1 Diabetes Mellitus. The 18th Annual Diabetes Technology Meeting. 2018. 30

ACS Paragon Plus Environment

Page 31 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 1: Flow chart of the proposed mAP integrated with the CPAS

31

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2: Adaptive controller set-point denition

Figure 3: Glycemic risk index 32

ACS Paragon Plus Environment

Page 32 of 44

15 10 5 0 0

1

2 Normalized PIC (I'k/Ibasal,k)

3

4

Figure 4: Plasma insulin risk index

3.5 3 Normalized PIC

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Plasma insulin risk index (PIRI)

Page 33 of 44

2.5 2 1.5 1 0.5 0

100

200 CGM (mg/dL)

Lower Bound for PIC

Upper Bound for PIC

300

Desired Value for PIC

Figure 5: Plasma insulin concentration bounds

33

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 6: Closed-loop results for all subjects based on the dened scenario in case 1

34

ACS Paragon Plus Environment

Page 34 of 44

Page 35 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 7: Plot of changes in Yterm and Uterm over time for case 1

35

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 8: Closed-loop results for all subjects based on the dened scenario in case 2

36

ACS Paragon Plus Environment

Page 36 of 44

Page 37 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 9: Plot of changes in Yterm and Uterm over time for case 2

37

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 10: Closed-loop results for all subjects based on the dened scenario in case 3

38

ACS Paragon Plus Environment

Page 38 of 44

Page 39 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 11: Closed-loop results for all 20 subjects for case 4. The bottom and top of the boxes are the rst and third quartiles and the line inside the box is the median. The whiskers ends represent minimum and maximum values and + indicates mean values.

39

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 12: Plot of changes in Yterm and Uterm over time for case 4

40

ACS Paragon Plus Environment

Page 40 of 44

Page 41 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 13: Plot of modication in Pfasting parameter over time

41

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 14: Adaptive modication of the controller set-point using historical data

42

ACS Paragon Plus Environment

Page 42 of 44

Page 43 of 44 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 15: Plot of prediction-error grid analysis (PRED-EGA) based on one-step-ahead (5 minutes) and 4-steps-ahead (20 minutes) prediction of CGM. The P-EGA compares the predicted CGM to the actual measured CGM. The R-EGA compares the rate of change of the predicted CGM to the actual rate of change of the measured CGM.

43

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

For Table of Contents Only

44

ACS Paragon Plus Environment

Page 44 of 44