Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- import mne
- import numpy as np
- import matplotlib.pyplot as plt
- from mne.time_frequency import psd_welch
- # Load the EDF file
- edf_file_path = '*.edf' # Replace with your file path
- raw = mne.io.read_raw_edf(edf_file_path, preload=True)
- # Select channels for analysis
- selected_channels = ['Cz-LE', 'O1-LE', 'O2-LE', 'F3-LE', 'F4-LE', 'Fz-LE', 'Pz-LE']
- raw_selected = raw.copy().pick_channels(selected_channels)
- # Rename the channels to match the standard 10-20 montage naming convention
- rename_mapping = {
- 'F3-LE': 'F3', 'Fz-LE': 'Fz', 'F4-LE': 'F4',
- 'Cz-LE': 'Cz', 'Pz-LE': 'Pz', 'O1-LE': 'O1', 'O2-LE': 'O2'
- }
- raw_selected.rename_channels(rename_mapping)
- # Apply a standard 10-20 EEG montage
- montage = mne.channels.make_standard_montage('standard_1020')
- raw_selected.set_montage(montage)
- # Apply a bandpass filter to retain frequencies between 1 and 45 Hz
- raw_selected.filter(1., 45., fir_design='firwin')
- # Apply Notch filter to remove power line noise at 50 Hz (or 60 Hz depending on your region)
- raw_selected.notch_filter(freqs=[50])
- # Re-reference to the average (common re-referencing step in qEEG analysis)
- raw_selected.set_eeg_reference('average', projection=True)
- # Detect and annotate bad segments (e.g., muscle artifacts, eye blinks)
- raw_selected = mne.preprocessing.annotate_muscle_zscore(raw_selected, ch_type='eeg', threshold=4.0)
- raw_selected = mne.preprocessing.annotate_movement(raw_selected, threshold=5.0)
- raw_selected = mne.preprocessing.annotate_flat(raw_selected, bad_percent=5.0)
- # Interpolate bad channels (if any were marked during artifact detection)
- raw_selected.interpolate_bads()
- # Apply the average reference projection
- raw_selected.apply_proj()
- # Apply ICA for advanced artifact removal
- ica = mne.preprocessing.ICA(n_components=15, random_state=97, max_iter='auto')
- ica.fit(raw_selected)
- # Detect and exclude components related to eye blinks or cardiac artifacts
- ica.detect_artifacts(raw_selected)
- raw_selected = ica.apply(raw_selected)
- # Calculate Power Spectral Density (PSD) for each channel
- psds, freqs = psd_welch(raw_selected, fmin=1, fmax=45, n_fft=2048)
- # Average the PSD values across each frequency band
- band_powers = {
- 'Delta': (0.5, 4),
- 'Theta': (4, 8),
- 'Alpha': (8, 12),
- 'Beta': (12, 30),
- 'Gamma': (30, 45)
- }
- # Create a dictionary to store averaged band power for each channel
- avg_band_powers = {band: [] for band in band_powers}
- for band, (fmin, fmax) in band_powers.items():
- # Find the indices of frequencies within the desired band
- band_idx = np.where((freqs >= fmin) & (freqs <= fmax))[0]
- # Average the power spectral density across the band for each channel
- avg_band_powers[band] = np.mean(psds[:, band_idx], axis=1)
- # Plot topographic maps for each frequency band
- fig, axes = plt.subplots(1, 5, figsize=(20, 4))
- fig.suptitle('Topographic Maps of EEG Frequency Bands', fontsize=16)
- for i, (band, power) in enumerate(avg_band_powers.items()):
- mne.viz.plot_topomap(power, raw_selected.info, axes=axes[i], show=False)
- axes[i].set_title(f'{band} Band')
- plt.show()
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement