A one-size-fits-all artificial pancreas for people with type 1 diabetes based on physiological insight and feedback control

Tobias K. S. Ritschel, Asbjørn Thode Reenberg, Emilie B. Lindkvist, Christian Laugesen, Jannet Svensson,
Ajenthen G. Ranjan, Kirsten Nørgaard, Bernd Dammann, John Bagterp Jørgensen This work was partially funded by the IFD Grand Solution project ADAPT-T2D (9068-00056B). A. T. Reenberg, T. K. S. Ritschel, B. Dammann, and J. B. Jørgensen are with the Department of Applied Mathematics and Computer Science, Technical University of Denmark, DK-2800 Kgs. Lyngby, Denmark. E. B. Lindkvist, C. Laugesen, J. Svensson, A. G. Ranjan, and K. Nørgaard are with Steno Diabetes Center Copenhagen, Clinical Research, DK-2730 Herlev, Denmark. Corresponding author: J. B. Jørgensen (E-mail: jbjo@dtu.dk).

Abstract

We propose a model-free artificial pancreas (AP) for people with type 1 diabetes. The algorithmic parameters are tuned to a virtual population of 1,000,000 individuals, and the AP repeatedly estimates the basal and bolus insulin requirements necessary for maintaining normal blood glucose levels. Therefore, the AP can be used without healthcare personnel or engineers customizing the algorithm to each user. The estimates are based on bodyweight, measurements from a continuous glucose monitor (CGM), and estimates of the meal carbohydrate contents. In a virtual clinical trial with all 1,000,000 individuals (i.e., a Monte Carlo closed-loop simulation), the AP achieves a mean time in range of more than 87% and almost 89% of the participants satisfy several glycemic targets.

I Introduction

Diabetes is a chronic disease where the pancreas produces insufficient amounts of insulin or the body is resistant to insulin. More than 10% of the world’s adult population suffers from this disease, and in 2021, USD 966 billion dollars were spent on diabetes (which corresponds to 9% of the global health expenditure) [1]. Type 1 diabetes (T1D) accounts for between 5–10% of all cases, and it is caused by autoimmune destruction of the pancreatic insulin-producing cells. Consequently, the pancreas does not produce any insulin. Therefore, people with T1D require daily insulin treatment in order to prevent high blood glucose concentrations (referred to as hyperglycemia). Long periods of hyperglycemia can cause damage to the nerves and eyes and lead to chronic kidney disease and cardiovascular disease. Additionally, incorrect insulin treatment can lead to low blood glucose concentrations (referred to as hypoglycemia). Severe hypoglycemia can lead to a number of acute complications, e.g., loss of consciousness and seizures.

People with T1D spend significant amounts of time on self-treatment. Therefore, there is considerable interest in developing automated insulin delivery systems which can assist them. Such systems are referred to as artificial pancreases (APs) [2], and they typically consist of 1) a sensor, often a continuous glucose monitor (CGM), 2) a control system, usually a control algorithm implemented on a smartphone or a dedicated device, and 3) an actuator, e.g., an insulin pump.

Many control algorithms have been considered for this purpose. They can be divided into model-free and model-based algorithms. Model-based algorithms typically use model predictive control (MPC) [3, 4, 5, 6] where a model is used to predict the body’s response to, e.g., meal carbohydrates and insulin. MPC is a well-proven control methodology that has been applied to many different types of processes [7]. However, it requires an accurate model. Automatic generation of such a model based on historical data (e.g., CGM measurements, administered insulin, and meal carbohydrates) is an ongoing field of research. Consequently, it remains difficult to make model-based APs widely available because a model must be developed for each individual. In contrast, model-free controllers only rely on a few pieces of information about the body for which estimates are readily available, e.g., the bodyweight, the basal insulin requirement, the insulin-to-carb ratio (ICR), and the insulin sensitivity factor (ISF). These controllers are often based on heuristics [8], fuzzy logic [9], or proportional-integral-derivative (PID) control. PID controllers have been successfully applied in many different industrial applications [10], and several researchers have proposed APs based on concepts from PID control. Marchetti et al. [11] proposed a PID controller which is switched off when a meal is announced and switched back on based on heuristical rules. Huyett et al. [12] described a PID control algorithm for intraperitoneal insulin delivery (whereas most APs deliver insulin subcutaneously which results in a more delayed insulin effect). Finally, Jørgensen et al. [13] and Sejersen et al. [14] present APs that 1) use concepts from PID control to estimate the basal rate and 2) compute the insulin bolus as a linear function of the meal carbohydrate content using the ICR (which is assumed to be known). However, these AP algorithms have not been tested on large numbers of real or virtual people. Therefore, it is currently unknown whether they can be adopted without healthcare personnel or engineers tuning the algorithms specifically for each individual (or group of individuals).

In this work, we present a one-size-fits-all AP algorithm with a single set of parameters tuned to a population of 1,000,000 virtual individuals with T1D. It simultaneously estimates the basal insulin and the meal insulin bolus curve (as a function of the meal carbohydrate content normalized with body weight). Therefore, it is straightforward for a user to start using the system. Furthermore, we demonstrate that, for a given objective function, the optimal meal insulin bolus is a nonlinear function of the meal carbohydrate content, and we argue that it can be approximated well by a piecewise linear function. We test the proposed AP using a virtual clinical trial with 1,000,000 participants over 52 weeks, i.e., we perform a Monte Carlo closed-loop simulation. All virtual participants meet the target on time in level 2 hypoglycemia, the mean time in range (TIR) is 87.2%, and almost 89% of the virtual participants meet all targets on TIR, time above range (TAR), time below range (TBR), average glucose, and glucose management indicator (GMI) described in [15].

The remainder of this paper is organized as follows. In Section II, we analyze the optimal insulin bolus as a function of the meal carbohydrate content, and in Section III, we present the AP. We present the results of the virtual clinical trial in Section IV, and conclusions are given in Section V.

II Analysis

In this section, we present a dynamic optimization problem for determining the optimal meal insulin bolus as a function of the meal size. We only use this optimization problem to analyze the meal insulin bolus curve. It is not used in the AP algorithm presented in Section III because it would require a model of each person using the AP.

II-A The dynamic optimization problem

The dynamic optimization problem determining the optimal meal insulin bolus flow rate is in the form


	$\displaystyle\min_{u_{0}}\quad\phi$	$\displaystyle=\int_{t_{0}}^{t_{f}}\rho(z(t))\,dt,$			(1a)
subject to

	$\displaystyle x(t_{0})$	$\displaystyle=x_{0},$			(1b)
	$\displaystyle\dot{x}(t)$	$\displaystyle=f(x(t),u(t),d(t),\theta),$	$\displaystyle t$	$\displaystyle\in[t_{0},t_{f}],$	(1c)
	$\displaystyle z(t)$	$\displaystyle=g(x(t),\theta),$	$\displaystyle t$	$\displaystyle\in[t_{0},t_{f}],$	(1d)
	$\displaystyle u(t)$	$\displaystyle=u_{k},\quad t\in[t_{k},t_{k+1}[,$	$\displaystyle k$	$\displaystyle=0,\ldots,N-1,$	(1e)
	$\displaystyle d(t)$	$\displaystyle=d_{k},\quad t\in[t_{k},t_{k+1}[,$	$\displaystyle k$	$\displaystyle=0,\ldots,N-1,$	(1f)
	$\displaystyle u_{\min}$	$\displaystyle\leq u_{0}\leq u_{\max}.$			(1g)

The objective function in (1a) is the integral of the penalty function $\rho$ over the time horizon $[t_{0},t_{f}]=[0,12]$ h, and $N$ is the number of control intervals. Furthermore, $t$ is time, $x$ are the state variables, $u$ is a vector of manipulated inputs (i.e., the basal and bolus insulin flow rates), $d$ are disturbance variables (i.e., the meal carbohydrate flow rate), $z$ is the output (i.e., the CGM measurement), and $\theta$ are model parameters. The initial condition (1b) is the steady state corresponding to a blood glucose concentration of 6 mmol/L, and the dynamical constraint (1c) is the model presented by Hovorka et al. [16] extended with a CGM model [17]. Next, (1d) is an output equation, and (1e)–(1f) are zero-order-hold (ZOH) parametrizations of the manipulated inputs and the disturbance variables. Finally, $u_{\min}$ and $u_{\max}$ in (1g) are lower and upper bounds on the manipulated inputs in the first control interval.

In the first control interval ( $k=0$ ), the person consumes a meal with a specified meal carbohydrate flow rate, $d_{0}$ , and an insulin bolus is administered. For the remaining control intervals ( $k>0$ ), the disturbances and the bolus insulin flow rate are zero. The insulin basal rate is equal to its steady state value (corresponding to 6 mmol/L) in all control intervals, i.e., it is not a decision variable. Finally, the lower bound on the insulin bolus flow rate in the first control interval ( $k=0$ ) is 0 and the upper bound is infinity.

The penalty function penalizes the deviation from the setpoint $\bar{z}=6$ mmol/L and the violation of the soft lower bound $z_{\min}=3.9$ mmol/L (see also Table I), i.e.,

\displaystyle\rho(z(t))

\displaystyle=\bar{\rho}(z(t))+\kappa\rho_{\min}(z(t)),

(2)

where


$\displaystyle\bar{\rho}(z(t))$	$\displaystyle=\frac{1}{2}(z(t)-\bar{z})^{2},$	(3a)
$\displaystyle\rho_{\min}(z(t))$	$\displaystyle=\frac{1}{2}\max\{0,z_{\min}-z(t)\}^{2}.$	(3b)

As the main priority is to avoid hypoglycemia, $\kappa=10^{6}$ .

II-B Optimal insulin bolus curves

For 6 virtual people with type 1 diabetes (i.e., 6 different sets of parameters, $\theta$ ), Fig. 1 shows the values of the objective function in (1a) for different combinations of (absolute) meal carbohydrate contents and insulin boluses. The black lines indicate the optimal boluses found by solving the dynamic optimization problem (1) using a single-shooting approach. The optimal insulin bolus curve is more nonlinear for some sets of parameters than others. For instance, a linear insulin bolus curve is a worse approximation for person 6 than for person 4. The first kink in the optimal insulin bolus curve (starting from the left) appears because the total non-insulin-dependent glucose flux decreases when the blood glucose concentration comes below 4.5 mmol/L in the model by Hovorka et al. [16]. Consequently, more insulin is required to decrease the concentration below this value. The second kink arises because of the soft lower bound, $z_{\min}$ , in the penalty function. There are no kinks for person 4 because their blood glucose concentration increases very little after meals, i.e., there would be kinks for meals with more than 150 g carbohydrates. In conclusion, a piecewise linear function is a reasonable approximation of the optimal insulin bolus curve.

Refer to caption — Figure 1: Values of the objective function (1a) for different (absolute) meal carbohydrate contents and insulin boluses. The black lines indicate the optimal insulin boluses as functions of the meal carbohydrate content, i.e., the solutions to the dynamic optimization problem (1).

TABLE I: The five glycemic ranges described by Holt et al. [15].

Category	Range [mmol/L]	Color
Level 2 hyperglycemia	]13.9, $\ \infty\;[$	Orange
Level 1 hyperglycemia	]10.0, 13.9]	Yellow
Normoglycemia	[ 3.9, 10.0]	Green
Level 1 hypoglycemia	[ 3.0, 13.9[	Light red
Level 2 hypoglycemia	[ 0.0, 13.0[	Red

III Algorithm

At time $t_{k}$ [min], the AP receives a CGM measurement, $y_{k}$ [mmol/L], and computes the basal and bolus insulin flow rates, $u_{ba,k}$ [mU/min] and $u_{bo,k}$ [mU/min]. These are clipped and collected in the vector of manipulated inputs,

\displaystyle u_{k}

\displaystyle=\max\left\{0,\min\left\{u_{\max},\begin{bmatrix}u_{ba,k}\\ u_{bo,k}\end{bmatrix}\right\}\right\},

(4)

which is administered over the following control interval, $[t_{k},t_{k+1}[$ . The upper bounds, $u_{\max}$ , on the insulin basal and bolus insulin flow rates are 55 mU/min and 8000 mU/min, respectively, and the control and sampling intervals are identical. Their lengths are $T_{s}=t_{k+1}-t_{k}=5$ min.

If the CGM measurement is above the target value of $\bar{y}=6$ mmol/L, the basal insulin flow rate is the sum of an estimated nominal basal insulin flow rate, $\bar{u}_{ba,k}$ [mU/min], and a microadjustment term, $u_{ma,k}$ [mU/min]. If the measurement is between the target and the safety threshold $y_{s}=3$ mmol/L, the microadjustments are clipped to be non-positive (for safety reasons). Finally, if the measurement is below the safety threshold, the basal rate is zero, i.e.,

\displaystyle u_{ba,k}

\displaystyle=\begin{cases}\bar{u}_{ba,k}+u_{ma,k}&\text{if}~{}\bar{y}\leq y_{k},\\ \bar{u}_{ba,k}+[u_{ma,k}]^{-}&\text{if}~{}\bar{y}>y_{k}>y_{s},\\ 0&\text{otherwise},\end{cases}

(5)

where $[\>\cdot\>]^{-}=\min\{0,\>\cdot\>\}$ . The nominal basal insulin flow rate is estimated using an extended integral (I) controller, and the microadjustments are computed using a proportional-derivative (PD) controller:


$\displaystyle\bar{u}_{ba,k}$	$\displaystyle=I_{ba,k},$	(6a)
$\displaystyle u_{ma,k}$	$\displaystyle=P_{ma,k}+D_{ma,k}.$	(6b)

Here, $I_{ba,k}$ is an integral term (see Section III-A) and $P_{ma,k}$ and $D_{ma,k}$ are proportional and derivative terms (see Section III-B).

Based on the analysis in Section II, we compute the meal bolus insulin flow rate as a continuous piecewise linear function of the estimated normalized meal carbohydrate flow rate, $\hat{d}_{k}$ [g CHO/(kg min)]. That is, $\hat{d}_{k}$ is the amount of meal carbohydrates the user announces they will consume in the $k$ ’th control interval, divided by the product of their bodyweight and the length of the control interval, $T_{s}$ . Specifically, the meal bolus insulin flow rate is given by

\displaystyle u_{bo,k}

\displaystyle=\begin{cases}\alpha_{k}d_{th}+\frac{\alpha_{k}}{\beta}(\hat{d}_{k}-d_{th})&\text{if}~{}\hat{d}_{k}>d_{th},\\ \alpha_{k}\hat{d}_{k}&\text{otherwise}.\end{cases}

(7)

If the normalized meal carbohydrate flow rate is below the threshold $d_{th}=0.1$ g CHO/(kg min), the insulin bolus flow rate is proportional to the meal carbohydrate flow rate and the slope, $\alpha_{k}$ [mU kg/(g CHO)], is essentially the inverse of the ICR. For higher normalized meal carbohydrate contents, the slope is divided by $\beta=2$ (unitless). Both the threshold and $\beta$ were identified by trial-and-error, and we estimate the meal bolus factor $\alpha_{k}$ using another extended I-controller, i.e.,

\displaystyle\alpha_{k}

\displaystyle=I_{bo,k},

(8)

where $I_{bo,k}$ is an integral term described in Section III-C.

III-A Estimation of the basal rate

At time $t_{k}$ , when a CGM measurement becomes available, we update the estimate of the basal rate and ensure that it is non-negative:

\displaystyle I_{ba,k}

\displaystyle=\max\{0,I_{ba,k-1}+\Delta I_{ba,k}\}.

(9)

The increment is

\displaystyle\Delta I_{ba,k}

\displaystyle=w_{ba,k}K_{I,ba}e_{ba,k}T_{s},

(10)

where the unitless binary weight $w_{ba,k}$ is 1 if the last meal was announced more than $\Delta t_{m}=9.5$ h ago and 0 otherwise. Furthermore, the integrator gain is $K_{I,ba}=4\cdot 10^{-4}$ mU L/(mmol min²), and the error is computed using the deadband $[y_{ba}^{\ell},y_{ba}^{u}]=[3.9,8.0]$ mmol/L and a (unitless) hypoglycemia amplification factor, $\gamma=100$ :

\displaystyle e_{ba,k}

\displaystyle=\begin{cases}y_{k}-y_{ba}^{u}&\text{if}~{}y_{k}>y_{ba}^{u},\\ \gamma(y_{k}-y_{ba}^{\ell})&\text{if}~{}y_{k}<y_{ba}^{\ell},\\ 0&\text{otherwise}.\end{cases}

(11)

III-B Microadjustments of the basal rate

The proportional term in the microadjustment of the basal insulin flow rate is

\displaystyle P_{ma,k}

\displaystyle=w_{ma,k}K_{P,ma}e_{k},

(12)

where the gain is $K_{P,ma}=0.3$ mU L/(mmol min), and the error is

\displaystyle e_{k}

\displaystyle=y_{k}-\bar{y}.

(13)

The unitless binary weight $w_{ma,k}$ is 1 if the CGM measurement is below the target or if $w_{ba,k}=1$ . Otherwise, it is zero:

\displaystyle w_{ma,k}

\displaystyle=\begin{cases}1&\text{if}~{}y_{k}<\bar{y}~{}\text{or}~{}w_{ba,k}=1,\\ 0&\text{otherwise}.\end{cases}

(14)

The derivative term is

\displaystyle D_{ma,k}

\displaystyle=w_{ma,k}K_{D,ma}\frac{y_{k}-y_{k-1}}{T_{s}},

(15)

where the gain is $K_{D,ma}=10$ mU L/mmol and we disregard changes in the setpoint (which is also constant in this work).

III-C Estimation of the meal bolus factor

As for the basal rate, we update the estimate of the meal bolus factor whenever a CGM measurement becomes available and ensure that it is non-negative, i.e.,

\displaystyle I_{bo,k}

\displaystyle=\max\{0,I_{bo,k-1}+\Delta I_{bo,k}\},

(16)

where the increment is

\displaystyle\Delta I_{bo,k}

\displaystyle=w_{bo,k}K_{I,bo}e_{bo,k}T_{s}.

(17)

The unitless binary weight $w_{bo,k}$ is 1 for a time period of $\Delta t_{m}$ after every announced meal, i.e., $w_{bo,k}$ and $w_{ba,k}$ never have the same value. Furthermore, the gain is $K_{I,bo}=0.05$ mU kg L/(g CHO mmol min), and we use both a deadband of $[y_{bo}^{\ell},y_{bo}^{u}]=[3.9,10]$ mmol/L, the hypoglycemia amplification factor $\gamma$ , and clipping to compute the error:

\displaystyle e_{bo,k}

\displaystyle=\begin{cases}y_{bo}^{th}-y_{bo}^{u}&\text{if}~{}y_{k}>y_{bo}^{th},\\ y_{k}-y_{bo}^{u}&\text{if}~{}y_{k}\in[y_{bo}^{u},y_{bo}^{th}],\\ \gamma(y_{k}-y_{bo}^{\ell})&\text{if}~{}y_{k}<y_{bo}^{\ell},\\ 0&\text{otherwise}.\end{cases}

(18)

The clipping ensures that all CGM measurements above the threshold $y_{bo}^{th}=13.9$ mmol/L result in the same error.

IV Virtual clinical trial

TABLE II: The compositions of the seasons and the weeks and the meal carbohydrate contents in the protocol described in [18].

Compositions of the seasons
Season	Standard weeks	Active weeks	Vacation weeks
Winter	6	4	3
Spring	6	6	1
Summer	7	3	3
Autumn	9	3	1

Compositions of the weeks
Week type	Standard days	Active days	Movie nights	Late nights
Standard	4	1	1	1
Active	3	3	1	0
Vacation	5	0	0	2

Bodyweight-dependent meal carbohydrate contents
Meal size	Amount of carbohydrates	For a 70 kg person
Large meal	1.29 g CHO/kg	90 g CHO
Medium meal	0.86 g CHO/kg	60 g CHO
Small meal	0.57 g CHO/kg	40 g CHO
Snack	0.29 g CHO/kg	20 g CHO

In this section, we test the AP algorithm described in Section III in a virtual clinical trial with 1,000,000 participants. The trial starts on January 1^st, 2021 and lasts 52 weeks. We use 1) the virtual population and the protocol presented by Reenberg et al. [18] and 2) a previously developed Monte Carlo simulation framework [19, 20]. However, we replace participants for which any time constant is more than 1 order of magnitude smaller or larger than the mean (i.e., we generate new participants). The protocol mimics a Northern European lifestyle, and it consists of 4 seasons lasting 13 weeks each. Each week is categorized as standard, active, or vacation, and all weeks consist of standard days, active days, days with a movie night, and days with a late night (see Table II and Fig. 2). Each participant is represented using the mathematical model presented by Hovorka et al. [16] extended with a CGM model [17] and an exercise model [21]. The initial estimates of the basal rate and the meal bolus factor are zero, i.e., $I_{ba,0}=I_{bo,0}=0$ , and the initial state is the steady state without insulin administration.

Fig. 3 shows the results of the virtual clinical trial for one participant over four different types of days. The basal rate is constant for most parts of the day, and it is decreased when the CGM measurements are below the target of 6 mmol/L (which mostly happens at night). A bolus is administered for each meal, and for this participant, the majority of the CGM measurements are within the normoglycemic range. However, the estimated nominal basal rate is quite low. Therefore, the CGM measurements increase over night.

In the following, we discuss the efficacy of the AP based on the last 48 weeks of the trial as the estimates of the nominal basal insulin flow rate and the meal bolus factor vary significantly during the first 4 weeks. The participant who obtains the lowest CGM measurement during all 52 weeks (specifically, 1.05 mmol/L) is referred to as the worst-case participant, and Fig. 4 shows the mean and worst-case TIRs. The mean of 87.2% TIR exceeds the target of 70%, and the time in level 1 and 2 hypoglycemia is low, even for the worst-case participant. Fig. 5 shows the cumulative distribution of the CGM measurements. The left tail shows that all participants spend less than 1% of the time in level 2 hypoglycemia and less than 8% of the time in level 1 and 2 hypoglycemia. This can also be seen in Fig. 6 which shows box plots of the TIRs. It also shows that most participants do not spend significant amounts of time in level 2 hyperglycemia. Table III shows the percentages of participants satisfying the glycemic targets described by Holt et al. [15]. Almost 82% satisfy all targets, and nearly 89% satisfy all average glucose, GMI, TAR, TIR, and TBR targets. Finally, Fig. 7 shows that most of the participants’ average total daily doses (TDDs) of basal and bolus insulin are in the intervals [7.5, 25] U/day and [5, 20] U/day, respectively. However, the distributions have long tails towards the right indicating that a few participants require high insulin doses.

TABLE III: Satisfaction of the glycemic targets described in [15].

Quantity	Target	Satisfied
Average glucose	$<154$ mg/dL	97.69%
GMI	$<7\%$	97.77%
GV	$\leq 36\%$	82.70%
TAR (level 2 hyperglycemia)	$<5\%$	93.97%
TAR (level 1 and 2 hyperglycemia)	$<25\%$	89.06%
TIR (normoglycemia)	$>70\%$	92.19%
TBR (level 1 and 2 hypoglycemia)	$<4\%$	98.85%
TBR (level 2 hypoglycemia)	$<1\%$	100.00%
All TAR, TIR, and TBR targets		88.92%
All targets except the GV target		88.87%
All targets		81.70%

V Conclusions

We propose a one-size-fits-all AP algorithm for people with T1D, which estimates both the basal insulin flow rate and the meal insulin bolus curve. It is based on physiological insight and concepts from PID control, and it only requires the bodyweight, CGM measurements, and meal carbohydrate estimates. We compute the meal insulin bolus as a piecewise linear function of the meal carbohydrate content normalized with bodyweight, and we test the AP algorithm in a 52 week virtual clinical trial with 1,000,000 participants. The mean TIR is 87.2%, and almost 89% of the participants satisfy targets on average glucose, GMI, TAR, TIR, and TBR.

References

[1] International Diabetes Federation, “IDF diabetes atlas,” 2021, ISBN: 978-2-930229-98-0.
[2] R. A. Lal, L. Ekhlaspour, K. Hood, and B. Buckingham, “Realizing a closed-loop (artificial pancreas) system for the treatment of type 1 diabetes,” Endocrine Reviews, vol. 40, no. 6, pp. 1521–1546, 2019.
[3] D. Boiroux, A. K. Duun-Henriksen, S. Schmidt, K. Nørgaard, S. Madsbad, N. K. Poulsen, H. Madsen, and J. B. Jørgensen, “Overnight glucose control in people with type 1 diabetes,” Biomedical Signal Processing and Control, vol. 39, pp. 503–512, 2018.
[4] D. Boiroux and J. B. Jørgensen, “Nonlinear model predictive control and artificial pancreas technologies,” in 2018 IEEE Conference on Decision and Control (CDC), 2018, pp. 284–290.
[5] A. Chakrabarty, E. Healey, D. Shi, S. Zavitsanou, F. J. Doyle III, and E. Dassau, “Embedded model predictive control for a wearable artificial pancreas,” IEEE Transactions on Control Systems Technology, vol. 28, no. 6, pp. 2600–2607, 2020.
[6] M. Messori, G. P. Incremona, C. Cobelli, and L. Magni, “Individualized model predictive control for the artificial pancreas: In silico evaluation of closed-loop glucose control,” IEEE Control Systems Magazine, vol. 38, no. 1, pp. 86–104, 2018.
[7] M. G. Forbes, R. S. Patwardhan, H. Hamadah, and R. B. Gopaluni, “Model predictive control in industry: Challenges and opportunities,” IFAC-PapersOnLine, vol. 48, no. 8, pp. 531–538, 2015.
[8] I. Capel, M. Rigla, G. Garcia-Saez, A. Rodriguez-Herrero, B. Pons, D. Subias, F. Garcia-Garcia, M. Gallach, M. Aguilar, C. Perez-Gandia, E. J. Gomez, A. Caixas, and M. E. Hernando, “Artificial pancreas using a personalized rule-based controller achieves overnight normoglycemia in patients with type 1 diabetes,” Diabetes Technology and Therapeutics, vol. 16, no. 3, pp. 172–179, 2014.
[9] T. Biester, J. Nir, K. Remus, A. Farfel, I. Muller, S. Biester, E. Atlas, K. Dovc, N. Bratina, O. Kordonouri, T. Battelino, M. Philip, T. Danne, and R. Nimri, “DREAM5: An open-label, randomized, cross-over study to evaluate the safety and efficacy of day and night closed-loop control by comparing the MD-Logic automated insulin delivery system to sensor augmented pump therapy in patients with type 1 diabetes at home,” Diabetes, Obesity and Metabolism, vol. 21, no. 4, pp. 822–828, 2019.
[10] K. J. Åström and R. M. Murray, Feedback Systems: An Introduction for Scientists and Engineers, 2nd ed. Princeton University Press, 2008.
[11] G. Marchetti, M. Barolo, L. Jovanovic, H. Zisser, and D. E. Seborg, “An improved PID switching control strategy for type 1 diabetes,” IEEE Transactions on Biomedical Engineering, vol. 55, no. 3, pp. 857–865, 2008.
[12] L. M. Huyett, E. Dassau, H. C. Zisser, and F. J. Doyle III, “Design and evaluation of a robust PID controller for a fully implantable artificial pancreas,” Industrial & Engineering Chemistry Research, vol. 54, pp. 10 311–10 321, 2015.
[13] J. B. Jørgensen, D. Boiroux, and Z. Mahmoudi, “An artificial pancreas based on simple control algorithms and physiological insight,” IFAC PapersOnLine, vol. 52, no. 1, pp. 1018–1023, 2019.
[14] M. Sejersen, D. Boiroux, S. E. Engell, T. K. S. Ritschel, A. T. Reenberg, and J. B. Jørgensen, “Initial titration for people with type 1 diabetes using an artificial pancreas,” IFAC PapersOnLine, vol. 54, no. 15, pp. 484–489, 2021.
[15] R. I. G. Holt, J. H. DeVries, A. Hess-Fischl, I. B. Hirsch, M. S. Kirkman, T. Klupa, B. Ludwig, K. Nørgaard, J. Pettus, E. Renard, J. S. Skyler, F. J. Snoek, R. S. Weinstock, and A. L. Peters, “The management of type 1 diabetes in adults. A consensus report by the American Diabetes Association (ADA) and the European Association for the Study of Diabetes (EASD),” Diabetologia, vol. 64, pp. 2609–2652, 2021.
[16] 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, and M. E. Wilinska, “Nonlinear model predictive control of glucose concentration in subjects with type 1 diabetes,” Physiological Measurement, vol. 25, no. 4, pp. 905–920, 2004.
[17] A. Facchinetti, S. D. Favero, G. Sparacino, J. R. Castle, W. K. Ward, and C. Cobelli, “Modeling the glucose sensor error,” IEEE Transactions on Biomedical Engineering, vol. 61, no. 3, pp. 620–629, 2014.
[18] A. T. Reenberg, T. K. S. Ritschel, B. Dammann, and J. B. Jørgensen, “High-performance uncertainty quantification in large-scale virtual clinical trials of closed-loop diabetes treatment,” in Proceedings of the 2022 American Control Conference (ACC), 2022.
[19] M. R. Wahlgreen, A. T. Reenberg, M. K. Nielsen, A. Rydahl, T. K. S. Ritschel, B. Dammann, and J. B. Jørgensen, “A high-performance Monte Carlo simulation toolbox for uncertainty quantification of closed-loop systems,” in Proceedings of the 60th Conference on Decision and Control (CDC), 2021, pp. 6755–6761.
[20] DTU Computing Center, “DTU Computing Center resources,” 2021, Technical University of Denmark. [Online]. Available: https://doi.org/10.48714/DTU.HPC.0001
[21] M. Rashid, S. Samadi, M. Sevil, I. Hajizadeh, P. Kolodziej, N. Hobbs, Z. Maloney, R. Brandt, J. Feng, M. Park, L. Quinn, and A. Cinar, “Simulation software for assessment of nonlinear and adaptive multivariable control algorithms: Glucose–insulin dynamics in type 1 diabetes,” Computers & Chemical Engineering, vol. 130, p. 106565, 2019.