How do i divide the compressed spectral array (CSA) of an EEG into REM and NREM stages?

I calculated the CSA of a nine hours EEG (that I previously filtered). Frequency of sampling = 512 Hz.

To do this, I divided the signal into windows of 300 seconds, with a window overlap of 50%, and estimated the power spectral density (PSD) for each window, with the autoregressive Yule-Walker method. The CSA is saved as a matrix, where the n-th row represents the PSD for the n-th time window (n=1: 0 - 5 // n=2: 2.5 - 7.5 etc. [min]).

At the very least, I'd want to highlight the windows where there's a clear predominance of NREM sleep, that is, a higher content of delta, theta waves, in general with high amplitude and thus corresponding to higher PSD values. If i can distinguish the intervals of consecutive windows where there's NREM sleep, I can assume the rest as REM, ideally. How can i do this?

Answer 1

[Shubham]

Hey @Davide Vezzari,

I understand that you want to highlight intervals of NREM sleep based on the power spectral density (PSD) values in the CSA matrix.

For identifying the intervals having high PSD values, you can try comparing it with a threshold value. The intervals with PSD values greater than the threshold could be considered as NREM sleep while the others would be considered as REM sleep.

Please refer to the following code snippet:

% Compute average PSD for each time window 
average_psd = mean(CSA,2); 


% Define threshold value for NREM sleep 
threshold = mean(CSA,"all") % Adjust this value based on your data 
 
% Identify windows of NREM sleep 
nrem_windows = zeros([length(average_psd),1]); 
start_index = 1; 
for i = 1:length(average_psd) 
    if average_psd(i) > threshold 
        nrem_windows(i)=1; 
    end 
end 
 
% Extracting the intervals 
nrem_intervals=[]; 
interval_start = false; 
startidx=1; 
for i=1:length(nrem_windows) 
    if nrem_windows(i)==1 
        if ~interval_start 
            interval_start=true; 
            startidx=i; 
        end 
    else 
        if interval_start 
            nrem_intervals=[nrem_intervals;startidx,i]; 
            interval_start=false; 
        end 
    end 
end 

% Display results 
disp(nrem_windows) 
disp(nrem_intervals) 

In the above code snippet, I used the overall average of PSD values as the threshold. I calculated the mean of PSD values for each window and compared it with the threshold to identify the window containing NREM sleep and then calculated the intervals.

The above code produces the following intervals as output:

The intervals represent the row indices of the matrix.

You can also try using the number of intervals above threshold to identify NREM sleep.

Hope this helps!!