2.1 Meta label correction

High-quality labels are crucial for enhancing the training performance of neural networks. Following the assumptions of [Reference Charikar, Steinhardt and Valiant12, Reference Veit, Alldrin, Chechik, Krasin, Gupta and Belongie49, Reference Wei, Gao, Chen and Duan51], we assume that labels with noise can be categorized into two distinct groups: one dataset with clean or trusted labels and other datasets with noisy or weak labels. In many cases, the dataset with clean labels is smaller than those with noisy labels, primarily due to the limited availability of clean labels and the associated high cost of labelling. Direct training on these smaller clean datasets can result in suboptimal performance, with a pronounced risk of overfitting. Similarly, training exclusively on noisy datasets—or a mix of both clean and noisy datasets—may also yield unsatisfactory outcomes, as large-scale models can inadvertently fit and memorize the noise [Reference Zhang, Bengio, Hardt and Recht56].

Our objective is to identify precursor features that contribute to seizures, enabling us to distinguish data oscillations caused by other behaviours. For epileptic patients, seizures often last only 30 seconds to 2 minutes, which makes seizure data scarce compared with nonseizure data and leads to data imbalance issues. Moreover, accurate labelling of seizure and nonseizure data serves as the foundation for model training, a process typically requiring the expertise of a physician and demanding substantial time and effort. Yang et al. propose the MMLC method to handle financial data with partial feature information and poor label quality [Reference Yang, Gao, Wei, Dai, Fang and Duan54]. Here, we have far more clean data than noisy data. Our dataset is partitioned into two subsets: a smaller dataset with potentially noisy or no labels during the period when seizures are possible, and a larger dataset with clean and reliable labels during nonseizure period or post-ictal period when the seizures have concluded.

Our aim is to employ clean data to facilitate the labelling of noisy data. Unlike conventional label generation methods, our approach relies on making label judgments for the data based on their predicted time domain, obviating the need to directly obtain labels for the noisy data. This concept operates with a two-layer optimization paradigm. Initially, we employ a meta model to predict labels for the noisy data and subsequently refine our meta model using clean data. The main network uses these corrected labels to classify noisy data.

2.2 Early warning indicators

Under long-term noise disturbances, a dynamical system may experience enormous fluctuations, which could lead to transitions between metastable states [Reference Gao, Duan, Kan and Cheng17, Reference Guo, Gao, Zhang, Han and Duan21]. Sudden transitions between alternative states are a ubiquitous phenomenon within complex systems, spanning a wide spectrum of domains ranging from cellular regulatory networks to neural processes and ecological systems [Reference Ashwin and Zaikin3, Reference Scheffer43]. Numerous instances of abrupt regime shifts have been observed in diverse research fields, including collapses of ecosystems [Reference Hirota, Holmgren, Van Nes and Scheffer24, Reference Wang, Dearing, Langdon, Zhang, Yang, Dakos and Scheffer50], abrupt climate shifts [Reference Drijfhout, Bathiany, Beaulieu, Brovkin, Claussen, Huntingford, Scheffer, Sgubin and Swingedouw14, Reference Lenton, Livina, Dakos, van Nes and Scheffer30] and onset of certain disease states such as atrial fibrillation [Reference Quail, Shrier and Glass37] or epileptic seizures [Reference Meisel and Kuehn31]. Warning signals for impending critical transitions are highly desirable, owing to the formidable challenges associated with restoring a system to its prior state subsequent to the occurrence of such transitions [Reference Folke, Carpenter, Walker, Scheffer, Elmqvist, Gunderson and Holling16, Reference Scheffer, Carpenter, Foley, Folke and Walker42]. Lade and Gross [Reference Lade and Gross28] introduced a method for the identification of early warning signals, offering a comprehensive integration of diverse information sources and data within the context of a generalized model framework. Proverbio et al. [Reference Proverbio, Skupin and Gonçalves36] undertook a systematic exploration of the characteristics and effectiveness of dynamical early warning signals across various conditions of system deterioration.

An increase in the variance within the pattern of fluctuations [Reference Carpenter and Brock11] represents a plausible outcome of critical slowing down when a system approaches a critical transition. This phenomenon can be rigorously demonstrated through formal mathematical representations [Reference Biggs, Carpenter and Brock7] and intuitively grasped. As the eigenvalue of the system approaches zero, the influence of perturbations becomes increasingly enduring, resulting in an amplified accumulation effect that subsequently augments the variance of the system’s state variable. In principle, one might anticipate that critical slowing down could impair the system’s capacity to trace these fluctuations, potentially leading to an opposite effect on the variance [Reference Berglund and Gentz5, Reference Berglund and Gentz6]. Nevertheless, investigations into mathematical models consistently reveal that an elevation in variance typically emerges and becomes discernible well in advance of the occurrence of a critical transition [Reference Carpenter and Brock11]. Neubert and Caswell [Reference Neubert and Caswell33] introduced a pair of quantitative metrics for characterizing transient dynamics, namely reactivity and the amplification envelope. These metrics can be computed by using the Jacobian matrix in the vicinity of a stable equilibrium. Reactivity offers an estimation of the system’s susceptibility to immediate perturbation-induced growth, while the maximum amplification quantifies the highest magnitude that a perturbation can transiently achieve relative to the system’s equilibrium.

In this work, one of our goals is to extract fluctuation patterns from data via the meta learning network as an early warning indicator for state transitions in EEG time series signals.

3.1 Automatic label optimization through meta learning

The meta-model comprises two neural networks: a meta network and a main network. The meta network’s objective is to enhance label accuracy for noisy data, while the main network’s role is to achieve correct classification. We denote the parameters of the meta network, considered as meta-knowledge, as $\boldsymbol {\alpha }$ , denoting the meta network as $g_{\boldsymbol {\alpha }}(X', Y')$ . Conversely, the main network’s parameters are denoted as $\boldsymbol {\omega }$ , and it computes the function $f_{\boldsymbol {\omega }}(X', y_c)$ , where $y_c$ represents the outcome obtained from the meta network.

The training process is illustrated in Figure 1. The data we feed into the meta network are divided into two parts. One part is the historical data $X'$ , input as a time window of length $h\in \mathbb {R}^+$ , that is, $X'$ represents the EEG time series from time t to $t+h$ . The other part is the known data $Y'$ , and its time window size is denoted as $m\in \mathbb {R}^+$ , that is, $Y'$ represents the EEG time series from time $t+h$ to time $t+h+m$ ; both $X'$ and $Y'$ are sequential in time. The meta network learns the label $y_c$ of data $Y'$ , where $y_c\in \{0,1\}$ . If the pattern of $Y'$ indicates that the patient had seizures during this time period, the meta network should output $y_c=1$ ; otherwise, it should output $y_c=0$ . As for the main network, its objective is to learn the label of the next (unknown) EEG signal window Y of the current time series data X, that is, the main network does not know the sequential data Y and is supposed to predict the label of unknown Y with known data X. This approach enables us to extract additional information about the prediction horizon and historical data without being influenced by noisy labels. We have already mentioned that labels given manually by experts can be subjective and are certainly expensive; our approach can avoid this problem by generating noisy data labels.

Figure 1 The framework of meta label generation model.

In every batch, we obtain both noisy and clean sample data. We initially input the noisy samples, along with their corresponding $(X', Y')$ pairs, into the meta network to derive accurate labels. Subsequently, the noisy data $X'$ are used as input for the current task’s classifier to generate the corresponding predicted labels. We update the main network classifier’s parameters based on the loss between the corrected and predicted labels. Following this, the clean data pair $(X, y)$ is introduced to the new classifier, and the classification loss is computed to facilitate the meta network’s update.

The bilevel optimization of our model can be formulated as follows:

$$ \begin{align*} \min_{\boldsymbol{\alpha}}\ \ & \mathbb{E}_{(X,y)\in D}\mathcal{L}_D(f_{\boldsymbol{\omega}^*_{\boldsymbol{\alpha}}}(X),y)\quad\text{where}\\ \boldsymbol{\omega}^*_{\boldsymbol{\alpha}} &=\arg\,\min\limits_{\boldsymbol{\omega}}\mathbb{E}_{(X',Y')\in D'}[\mathcal{L}_{D'}(f_{\boldsymbol{\omega}}(X'),g_{\boldsymbol{\alpha}}(X',Y'))], \end{align*} $$

and where the loss function is defined as

$$ \begin{align*}\mathcal{L}(a,b)=\frac{1}{N}\sum_{i=1}^N (b_i\ln{(a_i)}+(1-b_i)\ln{(1-a_i)}),\end{align*} $$

that is, the cross entropy loss for $a, b\in \mathbb {R}^{N\times 2}$ . Our algorithm is shown in Algorithm 1.

Updating the meta model presents a formidable challenge, as the resulting loss function does not explicitly feature $\alpha $ . Consequently, we have to compute $\partial \boldsymbol {\omega }^*_{\boldsymbol {\alpha }}/\partial \boldsymbol {\alpha }$ , which involves second-order gradient computations. Several different methods exist for this purpose. To circumvent resource-intensive demands, we employ the Taylor expansion method, as detailed by Zheng et al. [Reference Zheng, Awadallah and Dumais57], to compute the Hessian matrix $({\partial ^2}/{\partial \boldsymbol {\omega }^2}) \mathcal {L}_{D^{\prime }}(\boldsymbol {\alpha }, \boldsymbol {\omega })$ .

Denote the meta loss as

$$ \begin{align*}\mathcal{L}_{D}(\boldsymbol{\omega}^*_{\boldsymbol{\alpha}}) \triangleq \mathbb{E}_{(\mathbf{X}, y) \in D} \ell(y, f_{\mathbf{w}^*_{\boldsymbol{\alpha}}}(\mathbf{X})).\end{align*} $$

Then, we have

$$ \begin{align*} \min _{\boldsymbol{\alpha}} \mathcal{L}_{D}(\mathbf{w}^*_{\boldsymbol{\alpha}}) & \approx \mathcal{L}_{D}(\mathbf{w}^{\prime}_{\boldsymbol{\alpha}}) \\ & =\mathcal{L}_{D}(\boldsymbol{\omega}-\eta \nabla_{\boldsymbol{\omega}} \mathcal{L}_{D^{\prime}}(\boldsymbol{\alpha}, \boldsymbol{\omega})). \end{align*} $$

Subsequently, the update rule for the meta-knowledge $\boldsymbol {\alpha }$ can be written as follows:

(3.1) $$ \begin{align} \boldsymbol{\alpha}=\boldsymbol{\alpha}-\mu \nabla_{\boldsymbol{\alpha}} \mathcal{L}_{D}(\boldsymbol{\omega}^{\prime}_{\boldsymbol{\alpha}}). \end{align} $$

Then, the meta-parameter gradient is shown as

$$ \begin{align*} \frac{\partial \mathcal{L}_{D}(\boldsymbol{\omega}^{\prime}_{\boldsymbol{\alpha}})}{\partial \alpha}&= g_{\boldsymbol{\omega}^{\prime}}(I-\eta \nabla_{\boldsymbol{\omega}, \boldsymbol{\omega}} \mathcal{L}_{D^{\prime}}(\boldsymbol{\alpha}, \boldsymbol{\omega})) \frac{g_{\boldsymbol{\omega}^{\top}}}{\|g_{\boldsymbol{\omega}}\|^2} \frac{\partial \mathcal{L}_{D}(\boldsymbol{\omega})}{\partial \boldsymbol{\alpha}} \\ &\quad -\eta \nabla_{\boldsymbol{\alpha}}(\nabla_{\boldsymbol{\omega}}^{\top} \mathcal{L}_{D^{\prime}}(\boldsymbol{\alpha}, \boldsymbol{\omega}) \nabla_{\boldsymbol{\omega}^{\prime}} \mathcal{L}_{D}(\boldsymbol{\omega}^{\prime})). \end{align*} $$

3.2 Early warning probability indicator

The capacity to discern transitions between meta-stable states assumes a crucial role in the prediction and regulation of brain function. In our model, the objective of the main network is to learn the label of the next sequential data Y (unknown) based on the current time series data X, that is, the main network does not know the sequential data Y and should predict the label of unknown Y with known data X. Namely, the output of the main network is the probability of epileptic seizures occurring within the subsequent sequential time period following $[t,t+h]$ .

An escalation in the variance within the pattern of fluctuations, as established by Carpenter and Brock [Reference Carpenter and Brock11], can serve as a credible indicator when a system is near a tipping point and is experiencing a critical slowing down. We adopt the value of the main network’s output as an early warning indicator in the context of epilepsy, that is,

(3.2) $$ \begin{align} I_{P}=f_{\boldsymbol{\omega}}(X), \end{align} $$

which could be seen as the probability of ictal seizure occurring in the next time window. We conduct experiments in Section 4.2 to compare our early warning indicator $I_{P}$ with the variance of time series data.

4.1 Data processing

The scalp EEG data are recorded via the bipolar montage method [Reference Zaveri, Duckrow and Spencer55], where the voltage of each electrode is linked and compared with neighbouring electrodes to form a chain of electrodes. Figure 2(a) shows a sample of the data, and the red line represents time intervals identified as seizure occurrences by experienced physicians in this case. The horizontal axis represents the time of the test, while the vertical axis corresponds to specific brain regions where electrodes were inserted. Each segment of EEG data was recorded with a 23-channel sampling rate of 256 Hz and recording duration of 1 hour, during which five seizures were recorded for the subject. Our dataset is derived from this EEG data, collected at a frequency of 32 HzFootnote ¹ . Due to the high cost of acquiring such data, the size of the sample under examination remains modest. Consequently, the development of an early warning system based on this limited dataset becomes especially challenging. We collected data from five epileptic episodes of an epileptic patient and used four epileptic episodes to generate the training set and the remaining one episode as the test set. Figure 2(a) shows one epileptic episode.

Figure 2 A sample of our EEG data: (a) time series for the electrodes; the red part shows seizure states and the grey part represents rest states; (b) variance of the EEG data.

To establish the benchmark, we employ a sliding window technique, shifting the window left and right with a 0.5-second step. We set the time window of X to encompass 10 seconds, 15 seconds and 20 seconds, aligning precisely with durations of historical data. For Y, however, we have considered time windows spanning 5, 10, 15 and 20 seconds; these windows indicate the length of our prediction horizon.

Our dataset is divided into clean data and noisy data. Noisy data refer to data instances that are especially challenging for labelling due to their inherent complexity and lack of clear, unequivocal categorization. To be precise, we define an indeterminate zone around the clinically established onset time of epilepsy as the demarcation point, that is, we define a time window spanning $\pm 10$ seconds around this moment, amounting to 40 samples in every epileptic episode, as noisy data. Then, we expand the time window both forward and backward, collecting 40 clean data samples with labels 0 and 1, respectively, on every side. In instances where the right boundary of the Y-window intersects with this indeterminate, we designate the corresponding window X and Y as noisy data; we denote noisy data as $X'$ and $Y'$ , and obtain the resulting noisy dataset $D'=\{X',Y'\}$ . Noisy data do not have fixed determined labels, and our objective is to assign accurate labels to this data category without relying on expert annotations. In essence, we aim to generate precise labels for this noisy dataset autonomously. The rest of the data, which fall outside the indeterminate zones around onset times, is categorized as clean data. The clean dataset is partitioned into two distinct segments: one segment corresponds to the data obtained during the resting state, with label $0$ , and the other segment comprises data recorded during the epileptic state, with label $1$ . The clean dataset is denoted as $D=\{X,y\}$ , where $y\in \{0,1\}$ is the label.

Different electrodes have data with different orders of magnitude, so all data were normalized, and then noisy and clean data were divided according to different tasks. Figure 2(b) shows a sample of the covariance of the EEG data time series, with the red dotted line showing an epileptic seizure in the sample. The figure illustrates how growing variance could be seen as an indicator for an upcoming abrupt transition, which is consistent with [Reference Carpenter and Brock11].

4.2 Experimental results

4.2.1 Performance comparison

To test the performance of our approach, we have conducted experiments with a number of different models. During testing, we rely on the seizure onset moments delineated by physicians. Specifically, if the sample with X-window or Y-window intersects with the ictal moment demarcated by physicians, this sample is assigned with a label $y=1$ . Conversely, if neither the X-window nor the Y-window intersect with seizure onset moment identified by the physician, it is labelled by $y=0$ . Then accuracy is computed by comparing the model’s output labels with the corresponding ground truth labels y, that is, the number of true positives and true negatives divided by total sample size.

Table 1 provides a comprehensive overview of our classification results, comparing the performance of the baseline LSTM model, the ResNet model and our proposed model. Here, we use LSTM for the meta network and ResNet for the main network. Columns in Table 1 show the length of the input (historical) time window $|X|$ , and the rows correspond to different sizes of the prediction horizon $|Y|$ . It is clear that as the prediction horizon $|Y|$ (shown in the first column) grows, it becomes increasingly difficult to make accurate predictions for all three models. The accuracy of the LSTM model rises significantly as the length of the history duration becomes longer, because longer sequences in the input data contain more information. ResNet, however, has no particularly pronounced dependence on the input time length. Additionally, LSTM shows a very promising predictive ability over short periods of time, but the accuracy decreases rapidly as the prediction horizon length increases. The ResNet model, however, still shows a good degree of adaptation for prediction lengths greater than 10 seconds.

Table 1 The accuracy of LSTM, ResNet and our method. Columns correspond to the duration of the input data time X (shown in row 2) and rows correspond to different lengths of the prediction window Y.

Throughout all settings, our model shows substantial improvements in prediction accuracy. Naturally, all three models achieve higher accuracy values at shorter prediction ranges, but our model demonstrates excellent results for longer noisy data as well. Our approach reduces the error rate by at least $2\times $ in most cases: for 5 seconds prediction horizons our model reduces the error rate from 5–7% to 2–2.5%, for 10 seconds, from 10–15% to 5–8%, and for 15 seconds and 20 seconds the error goes down from 20–30% or higher to 15–25%. These results suggest that our approach based on meta learning can very significantly improve the performance of baseline models in epileptic seizure prediction tasks.

To show the efficiency of our predictions more intuitively, we compare predicted ictal probabilities with the true ictal moment identified by experienced physicians. Figure 3 illustrates the results with input duration $|X|=20$ and prediction horizons $|Y|=5$ (Figures 3(a), 3(b), 3(c)) and $|Y|=10$ (Figures 3(d), 3(e), 3(f)). The horizontal axis shows 120 consecutive time samples with an interval of 0.5 seconds, and the vertical axis represents the output probability of the model. Solid blue lines show the ictal probabilities learned by different models, and red dashed lines represent the true ictal moment. It is clear that for the LSTM model and ResNet model, the predicted onset time is approximately 5 seconds late compared with the real onset time. Both models almost always output probability values of $0$ or $1$ , with very few moments of output probability between $0$ and $1$ . This means that these two models can only learn features of clean data, that is, they cannot capture the features of noisy data and make a timely prediction. In contrast, our model has a clear fluctuation that occurs just before and after the onset of the seizure (see the position of the red dashed line in Figures 3(c) and 3(f)), which suggests that our model captures the changes in noisy data adequately.

Figure 3 Predicted ictal probabilities (blue lines) and true ictal onset moments (red dashed lines): (a–c) for $|X|=20$ and $|Y|=5$ ; (d–f) for $|X|=20$ and $|Y|=10$ .

4.2.2 Early warning prediction

Early warning of epileptic seizures is of paramount importance for individuals with epilepsy. In addition to the accuracy associated with predicted label generation, our research also places a significant emphasis on evaluating the early warning capacity inherent to the model.

Inspired by [Reference Carpenter and Brock11], we compare the tipping signal in terms of our indicator $I_{TP}$ and the variance indicator. Figure 4 shows the variance from Figure 2(b) averaged over a window of $0.5$ seconds (dashed blue line) compared with the predicted probability averaged over five experiments (solid red line), which we call the $I_P$ indicator. Figure 4 illustrates that there is an earlier tipping point in our indicator $I_{P}$ than in the variance indicator on the same data. It is clear that for any given reasonable threshold, our indicator $I_{P}$ always gives an earlier warning than the variance indicator, which suggests that our indicator $I_{P}$ is an effective early warning signal.

Figure 4 A comparison of the proposed indicator $I_{P}$ (red solid line) with the variance indicator (blue dashed line) for EEG data: (a) $|X|=20$ , $|Y|=5$ ; (b) $|X|=20$ , $Y=10$ .