A Continuous Clustering Algorithm for Detection of Local Sleep in Humans

The likelihood that a certain cortical column enters into a sleep-like state is essentially related to the degree of previous activity (1–3). Therefore, the sleep-like states are distributed spatially in a patchwork structure. Localized sleep-like states can be observed when the study subject is awake (4–6), while localized wake-like states can be observed when the study subject is sleep (1).

In a recent study, we presented evidence based on intracranial subdural grid recordings supporting the concept of local sleep while awake occurring in humans, and we postulated that this underlies the condition commonly described as drowsiness (7). At each moment in time, the electrical activity over each subdural grid electrode was classified as either “active” or “inactive.” The method for determining whether activity at a given subdural electrode was in the active or inactive state relied on the dual observations of an increase in the instantaneous amplitude of the delta band (δA), with a decreased instantaneous frequency in the combined theta-alpha band (θαF), while the study subject was in non-REM (NREM) sleep. The classifier utilized a k-means clustering algorithm that produced the binary outputs (either active or inactive).

One potential flaw in this approach lies in the use of the k-means algorithm. Forcing a binary classification might result in the appearance of a progression toward an increasing number of electrodes with activity appearing consistent with sleep while the study subject is still awake, while the truth might be that the actual underlying values have only shifted slightly on average, moving from one edge of the classifier to the other. A classifier that produces a continuous range of results rather than a binary classification would address this possibility. Accordingly, we have reanalyzed the human electrocorticography (ECoG) data using a fuzzy c-means (FCM) clustering algorithm (8). FCM clustering results in a continuous series of values rather than pure binary classifications, eliminating the potential biases created by the latter method.

One additional benefit to the use of FCM clustering is that it allows for the comparison of the relative contribution of each of the two variables (instantaneous θαF and instantaneous δA) to the overall state classification. While an increase in delta-band spectral power has been demonstrated to be a strong signature of the sleep state, we hypothesized that the instantaneous θαF might serve as a similar signature of the waking state, and the analysis of the individual contributions to the classification could test this proposition.

Methods

Details of the study subjects and the method of ECoG acquisition, along with a very detailed description of the instantaneous amplitude and frequency computation through Hilbert transform, is published elsewhere (7). Briefly, recordings were obtained from adult (age >17 years) patients with medically refractory epilepsy who were undergoing intracranial monitoring before resective epilepsy surgery. We used 12 short recordings (1–1.5 hours) from six study subjects, with six recordings having only the wake state and six recordings having both the wake and sleep states. We also used six long recordings of 24 hours each from four study subjects. Recording days during which seizures occurred were excluded from the analysis. A human observer reviewed the video portion of each sample to identify a 1-hour section during the day when the patient was behaviorally clearly awake (morning circadian peak) and 1 hour at night after sleep onset in which the patient was clearly asleep and remained asleep. The sleep samples were chosen at time points before 80 minutes after sleep onset in an effort to limit the recordings to periods of NREM sleep. When available, concurrent scalp EEG recordings were also reviewed for consistency with NREM sleep. These 1-hour samples served as the baseline exemplar recordings (individual to each patient), from which the salient features associated with the sleep and wake states were derived.

Computation of the Instantaneous Amplitude and Frequency

The instantaneous amplitude of a signal x(n) is the modulus of the analytical signal z(n)=x(n) + jHT[x(n)], where HT[.] is the Hilbert transform (9). The instantaneous frequency is the time derivative of the phase of analytical signal z(n), and it approximates the change with time of the signal x(n) frequency. First, we filtered the z-score of the raw ECoG signal in the delta band (0–4Hz) and theta-alpha band (4–12 Hz). Then, we computed the instantaneous amplitude in the δ band, denoted as δA, and the instantaneous frequency in the θα band, denoted as θαF, from the filtered signals. Finally, we used the median and moving average filtering to smooth and down sample both δA and θαF to 0.5 Hz. To process the 24-hour records, we down sampled δA and θαF to 1 Hz. We denoted as mA the average over all channels of the amplitude δA, , where δA_k is the amplitude of channel k, M is the number of channels, and n is the time moment. We denoted as mF the average over all channels of the frequency θαF, , where θαF_k is the frequency of channel k.

Computation of the Channel Activity Score

For each ECoG channel, we determined the activity score by using the fuzzy c-means (FCM) clustering (8). In the FCM clustering technique, each data point (in our analysis, the pair δA(n), θαF(n), where n is the time moment) belongs to each cluster to some degree between 0 and 1, which is specified by a membership function. In our analysis, we had two clusters corresponding to the wake (active) and sleep (inactive) states. Thus, for each pair, δA(n), θαF(n), we obtained the degrees of belonging to the active and inactive clusters. These degrees were given by the membership function of the active cluster, MFA(n), and the membership function of the inactive cluster, MFI(n), with MFI(n)=1–MFA(n). We used the MFA as the activity score. All data points with an MFA >0.5 belong to the crisp active set, which is in fact the projection of MFA on the plane of instantaneous amplitude and frequency. Meanwhile, the data points with an MFA <0.5 belong to the crisp inactive set (Figure 1A and 1B). In the fuzzy sets theory, the regular sets, which have the membership function 0 or 1, are called crisp sets (10). For comparison purposes, we also computed the activity scores when only δA (Figure 2A) or θαF (Figure 2B) was used in the FCM clustering. We denoted as mAS the average over all channels of the activity score, , where MFA_k is the activity score of channel k, and n is the time moment (Figure 2C).

Probability of Correct Classification

We computed the probability of correct classification (pcc) of two distributions, X₁ and X₂, in two crisp (regular) sets as follows:

In the above equation, d_p was the discriminability between X₁ and X₂ (11). In our analysis, X₁ was the distribution of (amplitude, frequency) pairs that belonged to the active crisp set, and X₂ was the distribution corresponding to the inactive crisp set (Figure 1A and 1B). For clustering with both amplitude and frequency (X₁ and X₂ were two-dimensional), we computed d_p as follows:

In the above equation, was the average operator, Q was the pooled covariance of matrixes X₁ and X₂, and T was the transpose matrix operator (12). When we used only one feature for classification (e.g., amplitude or frequency), the discriminability d_p was as follows:

In the above equation, m_i was the mean and σ_i was the standard deviation of the one-dimensional distribution X_i (i=1,2) (Figure 3).

Prediction of the Mean Activity Score

For each study subject, we computed the Pearson correlations between mAS (the average over all channels of the activity score), mA (the average over all channels of δA), and mF (the average over all channels of θαF) (Figure 4A). The existence of strong correlations between mAS, mA, and mF allowed us to predict successfully mAS as a function of mA and mF by using a linear regression model as follows:

In the above equation, pAS(n) was the prediction of mAS at the moment n, and b_i coefficients (i=0,1,2,3) were the regression coefficients, computed such that the mean square error between pAS and mAS was minimized (Figure 4). We used the coefficient of determination R² to estimate the prediction performance (13).

Results

Computing the Channel Activity Score

Recently, our group (7) determined the active (wake) and inactive (sleep) states of channels by using the k-means clustering algorithm (14). For each channel, the k-means algorithm partitioned the amplitude-frequency pairs, δA(n), θαF(n), into two clusters (active and inactive) such that each data pair belonged to the cluster with the nearest mean. This binary classification did not provide insightful information about how much either active or inactive was the channel at a certain time point. However, this information could be essential, especially when the study subject is awake and performs a certain cognitive task. We addressed this limitation by using FCM clustering instead of k-means clustering. After applying FCM, for each channel data point, δA(n), θαF(n), the degree of belonging of the data point to the active and inactive class is computed. This degree of belonging, called membership function in the fuzzy sets theory (8), is a number between 0 and 1 and could serve as a measure of how close the data point is to the perfect active and inactive states. We used the MFA class as the activity score of a channel (Figure 1A). The MFI class was equal to 1–MFA (Figure 1A). All data points with an MFA >0.5 were included in the crisp (regular) active set (Figure 1A and 1B), while all data points with an MFI >0.5 were included in the crisp inactive set (Figure 1A and 1B). The active and inactive crisp sets were used to quantify the salience of the amplitude and frequency features for active and inactive state discrimination.

In conclusion, for each channel, after applying FCM, we obtained an activity score sequence, MFA(n), which measured in the continuous interval ([0, 1]) the degree of activity of each amplitude-frequency data pair (Figure 1C). By putting together the activity of all the channels, we obtained the study subjects’ restgram results, a matrix of numbers in the interval. The restgrams were used in this study to illustrate the degree of local sleep in the study subjects.

Comparing Frequency and Amplitude for Activity Classification

We computed the channel activity score by clustering only the channel amplitude (Figure 2A), or only the channel frequency (Figure 2B), and compared the obtained average activity scores (mAS) with the average score for clustering with both amplitude and frequency (Figure 2C). During the wake state, we observed that the amplitude and frequency clustering resulted in a lower activity score than amplitude clustering (Figure 2C) but a higher activity score than frequency clustering (Figure 2C). For a conservative estimate of the activity score during the wake state, the frequency would be the best feature. The amplitude would be a better feature for a conservative estimate of the channel activity during the inactive state.

By clustering both the amplitude and frequency, we computed their statistics for the active (mAS >0.5) and inactive (mAS <0.5) periods. During the active periods, the standard deviation of θαF was significantly higher than the standard deviation of δA. However, the opposite was true for the inactive periods (p<0.001, rank-sum test) (Figure 2D), suggesting that the frequency was the most salient feature for classification during active periods, while the amplitude was the most salient feature during inactive periods. We did not find a significant difference between the modulus of correlation between mA and mAS, R(mA,mAS), and the correlation between mF and mAS, R(mA,mAS), for either the active or inactive period (t-paired test, p>0.05) (Figure 2E), indicating that for the prediction of mean activity, mAS, the mA and mF were equally important predictors.

After FCM clustering, for every channel, a distribution corresponding to the features included in the crisp active set (MFA >0.5) and another distribution corresponding to the crisp inactive set (MFA <0.5) can be derived (Figure 1A and 1B). The mean and the standard deviation of the amplitudes and frequencies included in the active and inactive crisp sets, made for all of the channels of subject SC539, are shown in Figure 3. The separation of features in two crisp sets allowed us to compute the pcc in these sets. Importantly, the pcc offered a quantitative way of estimating the importance of amplitude and frequency in the computation of the activity score. For example, Figure 3 also shows the pcc computed from active and inactive distributions when we used both the amplitude and frequency for clustering as well as only the amplitude and only the frequency. For subject SC539, the pcc computed after frequency classification (pccF) was significantly better than the pcc computed after amplitude classification (pccA). However, this result did not hold for all study subjects, and in general we did not observe that the θαF (the feature introduced by our group previously [7]) was significantly a more salient feature than the δA. Yet interestingly enough, we did find that for all study subjects, the pcc computed from both the amplitude and the frequency (pccAF) (Figure 4), was always better or significantly better than the pccA or pccF (Figure 3).

Linear Regression Successfully Predicts the Mean Activity Score From the Mean Amplitude and Frequency

Similar to our previous analysis completed with a binary classification scheme (7), the present continuous classification scheme showed strong correlations between the mean instantaneous frequency mF and mean instantaneous amplitude mA, both with each other (Figure 4) and between each separately and the mean activity score (Figure 4A). Considering these results, we predicted the mean activity score mAS with a linear regression scheme using as predictors the mA and mF (Figure 4). For all records, we obtained a good prediction performance of the mean activity score mAS (average coefficient of determination, R²=0.91) (Figure 4A). This demonstrates that after the regression coefficients for a study subject were determined, the iterative FCM classification process could be successfully replaced by the much-simple computations of linear regression, without affecting the accuracy of the mean activity score.

Performing the fuzzy clustering of the amplitude and frequency for 24-hour records, we observed that the mean frequency mF had large variations during the wake periods, while the mean amplitude mA had large variations during the sleep period (Figure 5A and 5B). This finding underscores again that the frequency carries more information than the amplitude during the wake periods, while the amplitude is more indicative during the sleep periods (Figure 2D). In addition, we observed strong correlations between mAS and mF, between mAS and mA, and between mF and mA (e.g., –0.89 in Figure 5B). These strong correlations allowed us to predict successfully (R²=0.94) (Figure 5C) the mean activity mAS (Figure 5C) from mA and mF, extending to 24-hour records the prediction results for short records of 1–2 hours.

Discussion

By substituting a continuous FCM clustering for the binary k-means classifier, we addressed the potential bias of binary classification, and the findings of our previous study still hold (7). Utilizing a continuous score for estimating the active state of a channel has the important advantage of increasing the accuracy in determining the drowsiness of the study subject. In our previous work, the drowsiness was measured as the fraction of inactive channels during the wake state. This drowsiness measure is acceptable when the number of channels is large (as in ECoG recordings) but less acceptable when the number of channels is reduced (as in the EEG recordings). Thus, for a small number of electrodes, a continuous score for the electrodes’ state is preferable to a binary score.

As an extension of our previous work, we proposed in the present study a method to quantify the importance of features (amplitude and frequency) for classification in two crisp sets (the active and inactive sets obtained after the FCM clustering). This method, based on computing pcc in the two sets, proved numerically that using together the amplitude and frequency for FCM clustering allowed us to obtain a higher pcc than by using only amplitude or only frequency. In addition, using only amplitude or frequency did not ensure significantly different pccs, indicating that amplitude and frequency were equally salient for FCM clustering. Importantly, the pccs computed for single features had high values, suggesting that even when using only one feature for classification it is possible to achieve an acceptable separability of the active and inactive sets. The fact that the instantaneous frequency and amplitude contribute equally to the classification value suggests that hypotheses attempting to explain the changes occurring in brain electrical activity with increasing time during the wake state must take both into account. The conclusion from our previous work—that the study “provides additional support for the concept of local sleep while awake, the first based on intracranial subdural grid recordings in humans”—remains intact (7).

Another essential, practical advantage offered by using the continuous activity scores is the direct connection that can be made between the mF or mA and the mAS. As a result of the high correlations between mF and mAS, as well as between mA and mAS, one could predict mAS by using only one predictor (mF or mA), and, consequently, the direct computations of mF or mA during the wake state could provide an estimate of the drowsiness of the study subject. This possibility raises the very important question, What feature (amplitude or frequency) is the best to estimate the drowsiness? We found that the frequency shows more variability than the amplitude during the active periods (mAS >0.5) and, consequently, could offer a more sensitive measure of drowsiness. Our previous work demonstrated a link between specific frequency in the alpha band and the variability of reaction time during an attention task, a study that was inspired by the clinical observation of slowing of the background alpha frequency accompanying drowsiness in routine EEGs (15). The present study is consistent with, and further expands upon, this finding. Given this perspective, future experimental studies relating the instantaneous θα frequency to the study subject’s performance on cognitive tasks will be instrumental in finding a suitable measure of drowsiness.

Conclusions

The detection of local sleep, based on fuzzy clustering, offers not only a continuous estimate of the sleep-like states during whole wake or wake-like states during whole sleep, but it also allows for direct connection between drowsiness and measurable characteristics of the ECoG signal, such as instantaneous frequency or amplitude, in certain EEG bands. This continuous approach could be instrumental in computing the drowsiness score for records with a small number of electrodes, such as those from the scalp.