A simple modeling framework for prediction in the human glucose-insulin system - PubMed (original) (raw)
A simple modeling framework for prediction in the human glucose-insulin system
Melike Sirlanci et al. Chaos. 2023.
Abstract
Forecasting blood glucose (BG) levels with routinely collected data is useful for glycemic management. BG dynamics are nonlinear, complex, and nonstationary, which can be represented by nonlinear models. However, the sparsity of routinely collected data creates parameter identifiability issues when high-fidelity complex models are used, thereby resulting in inaccurate forecasts. One can use models with reduced physiological fidelity for robust and accurate parameter estimation and forecasting with sparse data. For this purpose, we approximate the nonlinear dynamics of BG regulation by a linear stochastic differential equation: we develop a linear stochastic model, which can be specialized to different settings: type 2 diabetes mellitus (T2DM) and intensive care unit (ICU), with different choices of appropriate model functions. The model includes deterministic terms quantifying glucose removal from the bloodstream through the glycemic regulation system and representing the effect of nutrition and externally delivered insulin. The stochastic term encapsulates the BG oscillations. The model output is in the form of an expected value accompanied by a band around this value. The model parameters are estimated patient-specifically, leading to personalized models. The forecasts consist of values for BG mean and variation, quantifying possible high and low BG levels. Such predictions have potential use for glycemic management as part of control systems. We present experimental results on parameter estimation and forecasting in T2DM and ICU settings. We compare the model's predictive capability with two different nonlinear models built for T2DM and ICU contexts to have a sense of the level of prediction achieved by this model.
© 2023 Author(s). Published under an exclusive license by AIP Publishing.
Conflict of interest statement
The authors have no conflicts to disclose.
Figures
FIG. 1.
This schema shows how we divided T2DM patients’ data into training and test time windows. For each patient, the first week of data is the training time window used to estimate the model. The estimated model represents each patient’s personalized BG behavior and is used to forecast BG values over the test time window, which is of length three weeks and follows the training time window.
FIG. 2.
This schema shows our experimental design for the ICU setting. Each training time window has a length of approximately 24 h, which is used for model estimation. Then, we forecast the first BG value measured after the training time window. We perform this prediction for the whole dataset by moving the time windows.
FIG. 3.
Parameter estimation and uncertainty quantification in the T2DM setting. Panels (a)–(c) show results obtained with optimization and panels (d)–(f) show results obtained with MCMC for patients 1, 2, and 3, respectively. Both approaches are used in a patient-specific manner. We see that the point estimates obtained with two approaches are very close to each other in most cases. Also the width of the 1- and 2-stdev intervals, which are obtained with Laplace approximation (optimization) and directly from the approximate posterior samples (MCMC), are also in agreement with each other. Here, these intervals quantify the uncertainty in the point-estimates of the parameters. The parameter estimates and agree with real physiological values and the non-stationary behavior of the glucose dynamics of T2DM patients is reflected in the time-evolving behavior of the estimated parameters. All these features enforce the reliability of the parameter estimation results.
FIG. 4.
In panel (a), the model output of the estimated linear stochastic model is shown over the week of the training data along with the true BG measurements. Model output is a stochastic process and described by a mean and variance. The red circles show true BG measurements, the blue curve shows the mean of the model output, and the gray area represents the estimated 2-stdev band around the mean. Here 2-stdev band is used to quantify the oscillations of BG levels, which are not aimed to be tracked by the mean of the model output, but rather to be encapsulated by the gray region. The peaks in the model output show the BG response to the nutrition. Since the model aims to track the mean BG behavior (by blue curve) and capture the amplitude of BG oscillations (by the gray region), we plot the model output using a curve and a region. In panels (b)–(d), kernel density estimate (KDE) of 1000 different realizations of the estimated model output and BG measurements are shown for each patient. Comparison of the true BG measurements, which are assumed to be a realization of the model output, with the mean and 2-stdev band of the stochastic process—being the model output—along with the KDE plots in panels (b)–(d) implies that BG measurements could indeed be considered as a realizations of the random process.
FIG. 5.
Panels (a)–(c) show forecasting results in the T2DM setting obtained via models formed by using the estimated parameters with the optimization approach for patients 1, 2, and 3, respectively. In each plot, the red circles show the true BG measurements, the blue curve shows the mean of the model output, and the gray region is the estimated 2-stdev band around the mean of the model output, quantifying possible low and high values of the forecasted BG levels. These forecasting results show that the proposed model mean, when equipped with confidence bands found from standard deviations, estimate the BG levels accurately, and in a patient specific way. This reinforces the claim that the model parameters could be used to provide information about the health condition of individual patients.
FIG. 6.
Panels (a)–(c) show parameter estimation and uncertainty quantification results in the ICU setting obtained with MCMC for patients 4, 5, and 6, respectively. In each plot, the blue stars are the point-estimate of each parameter and the gray area is the 2-stdev band around the point-estimates (both obtained from the resulting random samples). Here, the gray region represents the uncertainty in the point-estimates. The estimated model parameters exhibit biophysically realistic values. Also, 1- and 2-stdev bands enforces the reliability of the estimated parameters, especially,
Gb
and
σ
, which are the most influential parameters in predicting the mean and variance of BG levels.
FIG. 7.
Panels (a)–(c) show forecasting results obtained based on parameters estimated with MCMC in the ICU setting for patients 4, 5, and 6, respectively. In each plot, the red circles show the true BG measurements, the blue stars show the mean of the model output, and the gray region shows the 2-stdev band around this mean, obtained from the model output of the stochastic model, forecasting the magnitude of the BG oscillations. These results are, in general, very close to those obtained using the optimization approach, and the most relevant properties are shared by them both. Obtaining similar results with another numerical solution technique based on the same mechanistic model shows the reliability of the model and estimated model parameters.
FIG. 8.
Smoothing piecewise constant nutrition function that is used for ICU patients.
FIG. 9.
BG simulations are shown with respective to the estimated parameters over respective training time window. In each plot, the light blue curve is the glucose rate in the nutrition delivered to the patient (right y-axis), the red circles show the true BG measurements (left y-axis), the red curve is the mean of the model output (left y-axis), and the gray area is the 2-stdev band around the mean of the model output (left y-axis). These figures show two main cases that could arise as a result of parameter estimation in the ICU setting: Panel (a): the input (nutrition rate) is reflected in the output (BG measurements), panels (b) and (c): the input is not reflected in the output, which makes it impossible to estimate the decay rate
γ
accurately.
FIG. 10.
Panels (a)–(c) show the UQ bands when we use the optimization approach for parameter estimation in the ICU setting for patients 4, 5, and 6, respectively. The blue stars are the point estimates for the respective parameters and the gray region quantifies the uncertainty in the estimation. Because of the parameter identifiability issues, the Hessian matrix is ill-conditioned. When we compute its inverse to obtain the variance in the estimated model parameters, we obtain unreasonably wide UQ bands.
FIG. 11.
Panels (a)–(c) show histograms for the elapsed time after meal time until the first BG measurement for all recorded meal data in the training time window for patients 1, 2, and 3, respectively. Lack of variability in the data of patient 2 causes suboptimal parameter estimation and forecasting results for this patient.
References
- Li J. and Kuang Y., “Analysis of a model of the glucose-insulin regulatory system with two delays,” SIAM J. Appl. Math. 67, 757–776 (2007). 10.1137/050634001 - DOI
MeSH terms
Substances
LinkOut - more resources
Full Text Sources
Medical