Exclude parameter for mne.io.raw.compute_psd() not working as expected

💬 Support & Discussions
1
  • MNE version: 1.9.0
  • operating system: e.g. macOS 15.3.2

When calling compute_psd on a raw object, according to the documentation , setting exclude = [] should result in all channels being processed, including bad channels. This does not seem to work for me. In the example below, there are 7 non-data channels :

n_eeg_channels = 26
n_non_data = 7
n_ch_total = 33

Perhaps I have misunderstood this parameter, please let me know. Thanks in advance

Expected behavior (as I undersand it)

When exclude = [], the returned psd data should be of shape n_eeg_channels x n_frequencies.

raw = create_base_raw() # creates a raw object with noise 
print('raw channels : ' , raw.info['ch_names'])
print('raw bads : ' , raw.info['bads'])
raw.info['bads'] = ['Fp1', 'Fp2', 'F7', 'F3', 'Fz'] # manually set bad channels
print('bads after setting : ' , raw.info['bads'])
psd = raw.compute_psd(method='welch',
                fmin=0.5,
                fmax=50,
                exclude=[])
print(" psd channels : " , psd.ch_names , "\nlen : " , len(psd.ch_names))
spectra, freqs = psd.get_data(return_freqs=True)

print(" spectra shape : " , spectra.shape , " freqs shape : " , freqs.shape)

This prints:

Creating RawArray with float64 data, n_channels=33, n_times=1000
    Range : 0 ... 999 =      0.000 ...     9.990 secs
Ready.
raw channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2', 'VPVA', 'VNVB', 'HPHL', 'HNHR', 'Erbs', 'OrbOcc', 'Mass']
raw bads :  []
bads after setting :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz']
Effective window size : 10.000 (s)
 psd channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2'] 
len :  26
 spectra shape :  (21, 496)  freqs shape :  (496,)

If bads = [], then the returned data shape is as expected.

Creating RawArray with float64 data, n_channels=33, n_times=1000
    Range : 0 ... 999 =      0.000 ...     9.990 secs
Ready.
raw channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2', 'VPVA', 'VNVB', 'HPHL', 'HNHR', 'Erbs', 'OrbOcc', 'Mass']
raw bads :  []
bads after setting :  []
Effective window size : 10.000 (s)
 psd channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2'] 
len :  26
 spectra shape :  (26, 496)  freqs shape :  (496,)

Setting picks = ['all'] does not seem to work either.

raw = create_base_raw() # creates a raw object with noise 
print('raw channels : ' , raw.info['ch_names'])
print('raw bads : ' , raw.info['bads'])
raw.info['bads'] = ['Fp1', 'Fp2', 'F7', 'F3', 'Fz'] # manually set bad channels
print('bads after setting : ' , raw.info['bads'])
psd = raw.compute_psd(method='welch',
                fmin=0.5,
                fmax=50,
                picks= ['all'])
print(" psd channels : " , psd.ch_names , "\nlen : " , len(psd.ch_names))
spectra, freqs = psd.get_data(return_freqs=True)

print(" spectra shape : " , spectra.shape , " freqs shape : " , freqs.shape)

prints

Creating RawArray with float64 data, n_channels=33, n_times=1000
    Range : 0 ... 999 =      0.000 ...     9.990 secs
Ready.
raw channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2', 'VPVA', 'VNVB', 'HPHL', 'HNHR', 'Erbs', 'OrbOcc', 'Mass']
raw bads :  []
bads after setting :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz']
Effective window size : 10.000 (s)
 psd channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2', 'VPVA', 'VNVB', 'HPHL', 'HNHR', 'Erbs', 'OrbOcc', 'Mass'] 
len :  33
 spectra shape :  (21, 496)  freqs shape :  (496,)

Helper function used to create raw mne object

def create_base_raw(n_times: int = 1000, sfreq: float = 100.0, init_zeros: bool = False) -> mne.io.RawArray:
    """
    Create a base MNE Raw object with synthetic EEG data using the standard 32-channel setup.
    
    Args:
        n_times: Number of time points
        sfreq: Sampling frequency in Hz
        neg_proportion: Proportion of values that should be negative (0.0-1.0)
        
    Returns:
        mne.io.RawArray: Base MNE Raw object
    """
    Config.set_global_seed()
    # Define channel types and names
    ch_types = ['eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg',
                'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg', 'eeg',
                'eog', 'eog', 'eog', 'eog', 'ecg', 'eog', 'emg']

    ch_names = ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3',
                'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2', 'VPVA', 'VNVB', 'HPHL', 'HNHR', 'Erbs',
                'OrbOcc', 'Mass']

    # Define channel positions
    dict_ch_pos_array = {'Fp1' : np.array([-0.02681, 0.08406, -0.01056]),
            'Fp2' : np.array([0.02941, 0.08374, -0.01004]),
            'F7'  : np.array([-0.06699, 0.04169, -0.01596]),
            'F3'  : np.array([-0.04805, 0.05187, 0.03987]),
            'Fz'  : np.array([0.00090, 0.05701, 0.06636]),
            'F4'  : np.array([0.05038, 0.05184, 0.04133]),
            'F8'  : np.array([0.06871, 0.04116, -0.01531]),
            'FC3' : np.array([-0.05883, 0.02102, 0.05482]),
            'FCz' : np.array([0.00057, 0.02463, 0.08763]),
            'FC4' : np.array([0.06029, 0.02116, 0.05558]), 
            'T7'  : np.array([-0.08336, -0.01652, -0.01265]), 
            'C3'  : np.array([-0.06557, -0.01325, 0.06498]),
            'Cz'  : np.array([0.000023, -0.01128, 0.09981]),
            'C4'  : np.array([0.06650, -0.01280, 0.06511]),
            'T8'  : np.array([0.08444, -0.01665, -0.01179]), 
            'CP3' : np.array([-0.06551, -0.04848, 0.06857]),
            'CPz' : np.array([-0.0042, -0.04877, 0.09837]), 
            'CP4' : np.array([0.06503, -0.04835, 0.06857]), 
            'P7'  : np.array([-0.07146, -0.07517, -0.00370]), 
            'P3'  : np.array([-0.05507, -0.08011, 0.05944]), 
            'Pz'  : np.array([-0.00087, -0.08223, 0.08243]),
            'P4'  : np.array([0.05351, -0.08013, 0.05940]), 
            'P8'  : np.array([0.07110, -0.07517, -0.00369]), 
            'O1'  : np.array([-0.02898, -0.11452, 0.00967]),  
            'Oz'  : np.array([-0.00141, -0.11779, 0.01584]),
            'O2'  : np.array([0.02689, -0.11468, 0.00945])
            }

    # Generate synthetic data
    n_channels = len(ch_names)
    if init_zeros:
        data = np.zeros_like(np.random.randn(n_channels, n_times))
    else:
        data = np.random.randn(n_channels, n_times)
    
    # Create montage
    montage = mne.channels.make_dig_montage(ch_pos=dict_ch_pos_array, coord_frame='head')
    
    # Create info object
    info = mne.create_info(ch_names=ch_names, ch_types=ch_types, sfreq=sfreq)
    info.set_montage(montage=montage, on_missing='raise')
    
    # Create Raw object
    raw = mne.io.RawArray(data, info)
    return raw
2

Hi,

After calling psd = raw.compute_psd(method='welch', fmin=0.5, fmax=50, picks=["all"]) you have 33 channels.

I think you need to do the same when calling: spectra, freqs = psd.get_data(return_freqs=True, picks=["all"]).

Cheers,

3

Ah indeed. For the behavior I want I need to do this :

from eeglearn.utils.augmentations import create_base_raw

raw = create_base_raw() # creates a raw object with noise 
print('raw channels : ' , raw.info['ch_names'])
print('raw bads : ' , raw.info['bads'])
raw.info['bads'] = ['Fp1', 'Fp2'] # manually set bad channels
print('bads after setting : ' , raw.info['bads'])
psd = raw.compute_psd(method='welch',
                fmin=0.5,
                fmax=50,
                picks= ['all'])
print(" psd channels : " , psd.ch_names , "\nlen : " , len(psd.ch_names))
spectra, freqs = psd.get_data(return_freqs=True, picks = ['eeg'], exclude = [])

print(" spectra shape : " , spectra.shape , " freqs shape : " , freqs.shape)

which prints the following. I have asked for only the eeg-channels (26), without excluding any bad channels. This can also be achieved with picks = ['data'] to pick only the data channels.

Creating RawArray with float64 data, n_channels=33, n_times=1000
    Range : 0 ... 999 =      0.000 ...     9.990 secs
Ready.
raw channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2', 'VPVA', 'VNVB', 'HPHL', 'HNHR', 'Erbs', 'OrbOcc', 'Mass']
raw bads :  []
bads after setting :  ['Fp1', 'Fp2']
Effective window size : 10.000 (s)
 psd channels :  ['Fp1', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8', 'FC3', 'FCz', 'FC4', 'T7', 'C3', 'Cz', 'C4', 'T8', 'CP3', 'CPz', 'CP4', 'P7', 'P3', 'Pz', 'P4', 'P8', 'O1', 'Oz', 'O2', 'VPVA', 'VNVB', 'HPHL', 'HNHR', 'Erbs', 'OrbOcc', 'Mass'] 
len :  33
 spectra shape :  (26, 496)  freqs shape :  (496,)

settings picks = ['all'] returns all 33 channels as expected :slight_smile:

Thanks @sotpapad !

1 Like