Temporal Specialization in LSTM Models for Enhanced Prediction of Adverse Glycaemic Events in Type 1 Diabetes

Diabetes, as a chronic disease, represents a major public health challenge in the world [1]. There are several types of diabetes, among which the most common are T1D, Type 2 Diabetes (T2D) and gestational diabetes. Each of these types has distinct characteristics, but they all share one thing: a failure in the BG regulation system. This irregularity results from a failure in the production of an essential hormone secreted by the beta cells of the pancreas, called insulin. The latter allows the body's cells to absorb glucose present in the blood. In the absence of insulin, glucose builds up in the blood, leading to serious complications.

T1D is usually diagnosed in children and young adults. In this case, the immune system attacks the beta cells of the pancreas. T2D, which is more prevalent in adults, is distinguished by insulin resistance or insufficient production of this hormone. Thus, cells become less responsive to insulin, making it difficult for cells to absorb glucose, contributing to an increase in BG levels. Gestational diabetes, on the other hand, occurs during pregnancy and can increase the risk of complications for both mother and child. Although this type of diabetes is often temporary, it requires close monitoring to avoid complications during pregnancy and reduce the later risk of developing T2D.

Diabetes management is based on approaches specific to each type of diabetes. Treatment of T2D and gestational diabetes may involve oral medications, insulin injections or other medications, lifestyle changes including eating a balanced diet and regular physical activity. People with T1D rely on insulin, which they give by injection. Two types of insulin are generally used: basal insulin and bolus insulin. Basal insulin aims to maintain a constant level of insulin in the body, working in the background to cover basic needs between meals and during the night. Alternatively, bolus insulin is administered before meals to control the increase in BG caused by food ingestion. Patients often use a set of devices, including insulin pumps and CGM systems. Insulin pumps are portable electronic devices that provide precise and flexible insulin delivery. They eliminate the need for frequent injections. On the other hand, CGM systems are small sensors implanted under the skin that continuously measure the concentration of glucose in interstitial fluid, providing patients with real-time information about their BG levels.

In this section, we present the methodology adopted to carry out our study. This includes the description of the dataset used, the different pre-processing steps applied, defining the input variables and the outcome labels, as well as the modelling and the parameter search procedure. Figure 1 presents an overview of our approach, outlining the fundamental aspects involved in the development and evaluation of the predictive model.

1.jpg

Figure 1. Schematic overview for model prediction development and assessment

4.1 Description of the OHIOTDM dataset

To evaluate the effectiveness of our approach, this study relies on data from OhioT1DM [5]. This dataset was specifically designed to be used for the development of data-driven models for predicting BG levels in T1D patients. The data was collected from 12 patients over a period of 8 weeks (data from 6 patients were made public in two in 2018, followed by 6 others in 2020). Each patient's dataset is divided into two subsets, reserving the last week of data for testing purposes. The features of the dataset show diversity, with variables such as BG measurements, insulin bolus, meal information, sleep information, exercise intensity, time of self-reported stress, illness, and hypoglycaemic events. A second category includes CGM, basal insulin level, galvanic skin response, heart rate, air and skin temperature, and step counter data.

4.2 Pre-processing

The development of prediction models is initiated with rigorous pre-processing steps performed on the dataset. This notably includes handling missing values and normalizing variables. The pre-processing steps are detailed as follows:

Converting XML files to CSV format: we converted the files provided in XML format to CSV format for easier further processing.

CGM data extraction: we collected all the data relating to CGM measurements, serving as inputs for our prediction model.

Filling CGM missing data: we initially replaced missing data points with values calculated by interpolating one step forward and one step backward in time. This method helps to estimate missing values based on the available data before and after the gaps, providing a more accurate representation of the BG levels during those time periods. However, any remaining missing CGM values were subsequently filled with zeros.

To ensure data quality, we implemented a rule regarding missing data rates. If a specific day or night-time period in the dataset exhibited a missing data rate greater than 33%, we opted to remove all data corresponding to that day or night. This strategy was employed to prevent potential bias and inaccuracies that could arise from significant gaps in the data. While this decision may affect the overall size of the dataset, it enhances our confidence in the integrity of the remaining data.

Merging training and testing data: we combined the training and testing datasets for each patient as we used the Leave one patient out cross validation approach by leaving one patient out in each iteration.

Splitting into diurnal and nocturnal sub-datasets: to define the specific time periods chosen for the diurnal and nocturnal models, we focused on leveraging the "sleep" and "basis_sleep" features available in the dataset. These features provided detailed information about when patients were asleep, allowing us to split the dataset into two distinct time periods: diurnal and nocturnal.

• "Sleep" feature: this data is reported directly by the patient, indicating the precise times when the patient starts and ends their sleep. It gives insight into the actual sleep habits of each individual, which may vary significantly depending on personal routines and lifestyle factors.
• "Basis_sleep" feature: this feature comes from the dedicated sensor tracking the patient's physiological state, automatically detecting when the patient is asleep. It provides a reliable, sensor-based indication of the patient's sleep schedule, which complements the self-reported "sleep" feature.

To identify the most accurate sleep window, we take the intersection of the two periods by combining the last sleep time with the earliest wake-up time from both features.

In cases where both the "sleep" and "basis_sleep" features were unavailable throughout the night, we standardized the sleep period by defining it from 10:00 p.m. to 6:00 a.m. This choice is based on the typical sleep habits observed among patients in the dataset. Many patients in the dataset follow a routine that aligns with these hours, making it a reasonable approximation for the nocturnal period when specific sleep data is missing.  Accordingly, the diurnal period was defined from the patient's awakening until the onset of sleep, representing the active, wakeful hours of the day.

Normalization of CGM values: to enhance the model's ability to identify patterns and relationships in the data, we applied normalization to the CGM values. Specifically, we used the min-max scaling technique, which transforms features to a fixed range between 0 and 1. This is essential because LSTM models typically perform better and train faster when input features are on a uniform scale. Large differences in feature magnitudes can lead to instability during training. Furthermore, normalization accelerates the convergence of gradient descent-based optimizers, such as the Adam optimizer (Adaptive Moment Estimation) used in this study.

4.3 Establishing input variables and outcome labels

Each data entry in the dataset for a given time point includes the current CGM value at that specific time point x_t, as well as a sequence of previous CGM values within a time window of N timestamps [x_t−1, x_t−2,…, x_t−(N−1)].

Our method consists of evaluating the presence of an adverse event during the PH, rather than looking for the event precisely at the time $t+P H$. Thus, we defined a binary variable y_t+PH indicating the presence (1) or absence (0) of the adverse event in the time interval $[t+1, t+P H]$.

If the adverse event is defined as hypoglycaemia, the threshold is set at 70 mg/dL. The risk label is defined as follows:

$y_{t+P H}=\left\{\begin{array}{lc}1 & { if } \exists \hat{t} \in[t+1, t+P H]: x_{\hat{t}}<70 \\ 0 & { else }\end{array}\right.$    (1)

And similarly, if the adverse event is defined as hyperglycaemia, the threshold is set at 180 mg/dL and the risk label is defined as follows:

$y_{t+P H}=\left\{\begin{array}{lc}1 & {if~} \exists \hat{t} \in[t+1, t+P H]: x_{\hat{t}}>180 \\ 0 & { else }\end{array}\right.$    (2)

The selection of thresholds at 70 mg/dL for hypoglycaemia and 180 mg/dL for hyperglycaemia aligns with the recommendations of the American Diabetes Association (ADA) [32].

The optimal value of N was identified as 6 timestamps (30 minutes) after following experimental tests which aimed to evaluate the performance of the model under different time window configurations.

In this study, we set the PH to 30 minutes. Demonstrations have shown that this duration is effective in initiating corrective measures to deal with the event [33]. Furthermore, it provides the best compromise between minimizing error in prediction outcomes (the longer the PH, the greater the error) and preserving forecast accuracy [3].

To reduce FPs due to the risk of triggering alerts or actions based on temporary fluctuations or unusual single measurements that may not be clinically meaningful; An event or prediction is considered valid if it lasts for at least 10 minutes, corresponding to a minimum of three consecutive timestamps.

4.4 LSTM

In this work, we adopted an approach based on DL, more precisely RNNs with an LSTM architecture in order to predict adverse events, specifically hypoglycaemia and hyperglycaemia. The use of LSTMs is motivated by their remarkable capacity to model the complex temporal dependencies found in data sequences, which makes them particularly suitable for modelling medical time series, where the temporal dynamics of variables is crucial [34]. Furthermore, previous research has demonstrated the effectiveness of LSTM in predicting future BG variability, reinforcing our decision to adopt this architecture to accurately anticipate hypoglycaemic and hyperglycaemic episodes [4, 24]. Formally, an LSTM can be defined as follows:

Let's assume an input sequence $x=\left(x_1, x_2, \ldots, x_t\right)$, where t stands for time. An LSTM processes this sequence step by step, and at each step t, it takes into account the current input $x_t$ as well as the previous hidden state $h_{t-1}$ and the previous memory cell $c_{t-1}$. The equations that describe the dynamics of an LSTM are generally formulated as follows [35]:

Forget gate: It determines what information from the previous memory cell $c_{t-1}$ should be forgotten or retained for the new memory cell $c_t$.

$f_t=\sigma\left(w_f \cdot\left[\begin{array}{ll}h_{t-1}, & x_t\end{array}\right]+b_f\right)$    (3)

Input gate: It decides what new information should be added to the memory cell. $\tilde{c}_t$ represents the candidates for the new memory cell.

$i_t=\sigma\left(w_i \cdot\left[h_{t-1}, x_t\right]+b_i\right)$    (4)

$\tilde{c}_t=tanh \left(w_c \cdot\left[h_{t-1}, x_t\right]+b_c\right)$    (5)

Memory cell update: It is updated based on decisions made by the forget gate and the input gate.

$c_t=f_t \cdot c_{t-1}+i_t \cdot \tilde{c}_t$    (6)

Output gate: It determines the current output $h_t$ based on the updated memory cell and decides what information should be transmitted to the output. In these equations, σ represents the sigmoid function (logistic function), $tanh$ is the hyperbolic tangent function, W and b are the weights and biases of the network, and [$\cdot,\cdot$] indicates the concatenation of vectors.

$o_t=\sigma\left(w_o .\left[h_{t-1}, x_t\right]+b_o\right)$    (7)

$h_t=o_t \cdot \tanh \left(c_t\right)$    (8)

In summary, an LSTM employs gates to control the flow of information within the memory cell, enabling it to capture long-term temporal dependencies in data sequences.

4.5 Modelling

Our work stands out due to the utilization of multi-layer LSTM models, allowing for concurrent predictions across 6 PHs at intervals of 5, 10, 15, 20, 25, and 30 minutes. Thus, we adopted an approach where the initial layers are specifically trained for a PH of 30 minutes, thus facilitating the learning of general characteristics inherent in the data, and the last layers are specific to each PH, allowing the model to capture information finer and specialized contextual information. Each set of specific layers receives information specific to the PH it aims to predict. This multi-output LSTM architecture is illustrated in Figure 2.

Our choice of a multi-output LSTM model represents a strategic divergence from the multi-step approach. Indeed, the multi-step approach involves the propagation of prediction errors from one step to the next, which can lead to an accumulation of errors and a degradation of the model's performance on long-term predictions. By opting for a multi-output model, we circumvented this problem by processing all predictions simultaneously, thus avoiding the transmission of errors from one step to another.

2.jpg

Figure 2. The multi-output LSTM architecture

After the LSTM layers, a dense output layer is added for each specific PHs. Additionally, the SoftMax activation function for binary classification was used to normalize the final outputs, thus providing probabilities associated with each predicted class. Subsequently, a combination approach for the final outputs related to each PH is employed to determine the ultimate prediction, indicating whether an undesirable event will occur. This methodology is detailed in Section 4.7 below.

4.6 Optimization of model parameters

Optimizing a model's predictive performance requires careful selection of hyper-parameters to achieve an optimal balance between model capability and prevention of over fitting.

In this study, our architectural choice was oriented towards the use of an approach where one layer is shared for all PHs, while two specific layers are dedicated to each PH.

The hyper-parameters are adjusted following experiments for a single patient, without resorting to leave-one-out cross-validation. The shared layer was set to 128 LSTM units, while each time horizon-specific LSTM layer is configured with 64 units each. The Adam optimizer [36] was chosen for this study due to its ability to effectively handle the complexities of CGM data, which includes challenges such as BG variability, sensor inaccuracies, noise, and non-stationary patterns. Unlike other optimization algorithms like SGD, RMSProp, and AdaGrad, Adam adapts the learning rates for each parameter based on both the first moment (the mean of the gradients) and the second moment (the variance of the gradients), allowing it to handle noisy and complex datasets efficiently. The learning rate was set to 0.01 initially, but Adam dynamically adjusts it throughout training, ensuring stable convergence and preventing issues like overshooting or local minima. This adaptability makes Adam especially suitable for CGM data. Training spans 50 epochs, with a batch size of 64, and 30% of the data is reserved for validation. The categorical Cross-Entropy loss function is applied to each temporal output. In order to counter over fitting, a dropout layer, with a rate of 20%, is introduced after the shared layer and after the first layer specific to each PH. This approach aims to improve the generalization of the model by randomly deactivating a proportion of the units during the training phase. An early stop mechanism is built in, stopping training if there is no improvement for five consecutive epochs to prevent over fitting.

4.7 Evaluation method

In this work we aim to enhance key performance metrics of the models developed, specifically precision, recall, specificity and F1-score, the definitions of which are provided below.

Precision: measures the proportion of adverse events predicted by the model that are actually predicted. In seeking to increase accuracy, the goal is to minimize FPs, thereby ensuring that the model's positive predictions are highly reliable (see Eq. (9)).

Precision $=\frac{\# \mathrm{TPs}+\text { #FPs }}{\# \mathrm{TPs}}$    (9)

Rappel (Recall): evaluates the model's ability to correctly identify all true adverse events. In seeking to increase recall, emphasis is placed on minimizing False Negatives (FNs), ensuring that the model does not miss important adverse events (see Eq. (10)).

$Rappel =\frac{\text { #TPs }+ \text { #FNs }}{\text { #TPs }}$    (10)

F1-Score: a metric that considers both precision and recall, providing an overall assessment of model performance. In seeking to improve the F1-score, the goal is to achieve an optimal balance between precision and recall to maximize the overall quality of the model's predictions (see Eq. (11)).

F1 - Score $=2 \times\left(\frac{\text { Precision }+ \text { Rappel }}{\text { Precision } \times \text { Rappel }}\right)$   (11) 

In these formulas:

True Positives (TPs): these indicate cases where the model correctly predicted an adverse event.

False Positives (FPs): these denote cases where the model incorrectly predicted an adverse event.

False Negatives (TNs): these represent cases where the model failed to predict an adverse event.

True Negative (TNs): these are cases where the model correctly predicts the absence of an adverse event, and this corresponds to reality.

Our model generates six outputs; to obtain an overall prediction over a 30-minute period, we aggregated the predictions for different PHs using indicator notation ($\mathbb{I}$) defined as follow:

$E=\mathbb{I}\left(\sum_i^n y_i>0\right)$    (12)

where, $y_i$ is the output of the model for the i-th PH. The condition states that an event is considered to have occurred if at least one model output is equal to 1. n is set to 6 corresponding to 6 PHs.

Subsequently, to measure the effectiveness of our approach to predicting adverse events over a 30-minute period, we compared the list of events that actually occurred to that of predicted events. This comparison was based on the use of the algorithm below which is directly inspired by the previous work [3]. This allows quantifying how close or distant the model predictions are from the actual events.

This algorithm is used to evaluate the performance of an event prediction model by comparing predicted events (P) to actual events (V) over a time scale. The idea is to determine TP, FP, and FN based on a tolerance range around the actual time instants.

Here is the explanation of the algorithm:

1. Initialization of sets V and P:

V is the set of temporal instants associated with real events.

P is the set of predicted events.

2. Initializing constants:

k is a positive constant.

ph is the prediction horizon.

3. Initializing result sets:

TP_set, FP_set, and FN_set are empty sets that will be used to store TPs, FPs, and FNs respectively.

4. Iterating through each element tp in P:

For each element tp in the set of predicted events P, the algorithm searches for a real time instant tv in the set of real events V such that the difference between tp et tv is within the tolerance range defined by -k < tp – tv < ph.

5. Searching for a tv in V:

The algorithm goes through all the elements of the set V to find a tv that satisfies the condition.

6. Processing of results:

If a tv is found, tp is considered a TP and tv is removed from the set V to avoid multiple consideration of the same event. If no tv is found, tp is considered a FP.

7. Identification of FN:

The remaining elements in set V after scanning all predicted events are considered FN.

In the context of BG sensing, prediction distances can be interpreted as the temporal difference between when the model predicted an event and the actual time when that event occurred.

$\text{Prediction Distance} =\text{Actual Event Time - Predicted Event Time}$    (13)

A prediction distance close to zero indicates that the model made a prediction very close to the actual time of the event.

A positive prediction distance suggests that the model anticipated the event, while a negative distance suggests a delay from the actual event.

Below are the different conditions that can be associated with a predicted event.

• Δtpv<−k: The event is out of analysis range.
• −k≤Δt_pv<0: The actual event occurs after the predicted time tp, but still within the tolerance range.
• Δt_pv=0: The actual event occurs exactly at the predicted time.
• 0<Δtpv<ph: The actual event tv occurs before the predicted time tp, but still in the future.
• Δtpv≥ph: The actual event has already happened. His prediction failed.

This study assesses the performance of models in predicting hypoglycaemic and hyperglycaemic events, emphasizing the influence of temporal variations on their accuracy. The initial hypothesis suggests that differences in glycaemic behaviour between diurnal and nocturnal periods necessitate a tailored approach, with prediction models developed from distinct datasets specific to each time period.

The performance of models based on the full dataset is slightly better and this can be explained by the fact that they likely benefit from a greater diversity of data, which may contribute to their ability to generalize and predict effectively.

The ability of the specific diurnal and nocturne models to maintain high performance even with smaller datasets suggests some adaptability. This could indicate that specific temporal features can be learned efficiently with smaller datasets, while still maintaining high relevance for prediction.

The relatively lower performance in predicting hypoglycaemia may be attributed to the intrinsic complexity of these events. Hypoglycaemia can occur quickly and requires increased sensitivity of models to detect them early.

The use of RNN, in particular LSTMs, to model the temporal sequences of glycaemic data is consistent with the dynamic and sequential nature of medical data. LSTMs are well suited to capture complex temporal dependencies, making them valuable in the context of BG, where changes can occur quickly and are strongly influenced by past events. This allows the model to take historical information into account, providing a better understanding of BG trends and improving prediction accuracy.

Using multi-layer LSTM models with simultaneous predictions over multiple PHs provides the flexibility to capture both general data characteristics and contextual information specific to each PH, thereby improving the model's ability to anticipate glycaemic events.

In this study, we opted for an architectural approach where one layer is shared across all prediction horizons (PHs), while two distinct layers are allocated to each individual PH. However, it is important to note that this architectural choice is not exhaustive, and other configurations deserve to be explored regarding the number of shared layers and distinct layers. The complexity of glycaemic dynamics can vary significantly across time horizons, and adjusting the number of layers based on each PH output could potentially improve the model's ability to capture these nuances.

The decision to opt for an evaluation which defines a binary variable indicating the presence (1) or absence (0) of the adverse event in the time interval [t+1,t+PH], contributes to maintain model stability and accuracy across PH of 30 minutes.

The distance Δt_pv is a way to quantify the temporal proximity between predicted and actual events, and it is used to make decisions about the accuracy of predictions in a given context. The value of the constants k and ph influences the tolerance granted to this temporal proximity.

The similarity of the results between the models developed based on “Diurnal”, “Nocturnal” and “All day” periods highlight their relative stability throughout the day. This temporal consistency suggests that the models maintain their performance regardless of the specific period considered.

It is worth noting that conducting a thorough comparison of model performance faces challenges due to differences in datasets, variations in dataset sizes, distinctions in data pre-processing methods, and disparities in parameters adjustment. Thus, the comparison of the results displayed in the Table 1 lacks accuracy and substantive significance [3]. Consequently, studies relying on the OHIOT1DM dataset offer the most suitable grounds for comparison. Our approach yields highly comparable results with those presented in Table 1 based on the OHIOT1DM dataset.

From a clinical standpoint, the results of this study carry significant practical clinical implications, particularly due to the use of real-world patient data for both model training and evaluation. The BG measurements were captured by a clinically certified CGM sensor, ensuring that the data reflect real-life conditions.

The performance metrics suggest that the model can be used as a decision-support tool for alerting potential adverse events, especially given the promising results achieved [37]. However, further validation with additional real-world datasets is essential to fully confirm its robustness across diverse populations. While the model demonstrates promising results, areas for improvement remain, including reducing false negatives to minimize the risk of missing critical hypoglycaemic or hyperglycaemic events. Exploring other ML approaches, based on night-time and daytime separation, could improve prediction performances, ensuring that the model is reliable enough for broader real-world clinical use.

The development of predictive models in healthcare, particularly those utilizing CGM data, presents several significant societal concerns that must be addressed. Data privacy and security are paramount, as sensitive CGM information can lead to privacy violations if not properly secured; implementing strong measures like encryption and regulatory compliance is essential. Additionally, algorithmic bias can result in inaccurate predictions for diverse populations if models are trained on unrepresentative datasets, highlighting the need for varied training data and continuous performance monitoring. There is also the risk of over-reliance on AI, which can diminish the essential role of healthcare professionals; thus, AI should be viewed as a supportive tool that requires human oversight. Access to advanced technologies may exacerbate health inequities, particularly in low-income or rural areas, necessitating efforts to enhance access through public healthcare integration and local partnerships. Furthermore, misinterpretation of model predictions can lead to inappropriate medical decisions, making it crucial to provide clear explanations and training for healthcare providers. Lastly, the psychological impacts of AI predictions, such as anxiety from critical alerts, should be mitigated by prioritizing essential notifications and minimizing unnecessary alerts. Proactively tackling these concerns will guarantee the ethical and responsible implementation of BG prediction healthcare applications [37, 38].