NeuroAnalyzer tutorials: Remove power-line noise

Load data:

using NeuroAnalyzer
using ContinuousWavelets
eeg = import_edf("files/eeg.edf");
[ Info: Precompiling NeuroAnalyzer [b40dcafa-b65c-47ad-a230-3174ba5cadd1] (cache misses: include_dependency fsize change (2), wrong dep version loaded (2), incompatible header (4), mismatched flags (10))
[ Info: NeuroAnalyzer v0.25.4-dev
[ Info: NeuroAnalyzer path: /home/eb/Documents/Code/NeuroAnalyzer.jl
[ Info:  Preferences:
[ Info:     Use CUDA: false
[ Info: Progress bar: true
[ Info:      Verbose: true
[ Info: Exclude bads: false
[ Info: Preparing resources
[ Info: Loading plugins:
[ Info:  Loaded: na_test_plugin.jl
[ Info:  Loaded: plot_env.jl
[ Info:  Loaded: plot_ispc.jl
[ Info:  Loaded: plot_itpc.jl
[ Info:  Loaded: plot_pli.jl
[ Info: Imported: EEG (24 × 308480 × 1; 1205.0 s)

Power-line noise (at 50 or 60 Hz, depending on the country) may be automatically attenuated using remove_powerline(). Currently the following methods are implemented:

  • IIR notch filtering + automated detection and removal of power-line peaks; for best results the notch filter bandwidth is calculated separately for each peak for each channel (this may take approx. 2 minutes per object)

Detect power-line noise frequency using linear regression method

The detection algorithm is based a series of linear regression models that estimates the relationships between the input signal and the sum of sine and cosine waveforms of consecutive frequencies (from 1 Hz to the input signal Nyquist frequency). For each model the sum of squared coefficients for the sine and cosine components of the model is used to calculate the noise power. Noise frequency is determined from the frequency that produces the highest power in these models. This experimental modeling is based on GLM.lm() and requires further validations.

The output is in Hz and shaped as channels × epochs.

noise_frq = detect_powerline(eeg)
24×1 Matrix{Float64}:
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0
 50.0

Remove power-line noise using IIR notch filtering method

Plot PSD before filtering:

plot_psd(eeg, ch="Fp1")

Plot spectrogram before filtering:

e10 = epoch(eeg, ep_len=10)
plot_spectrogram(e10, ch="Fp1", ep=2)

Detect and remove power-line noise at 50 Hz:

remove_powerline!(eeg, ch="all", method=:iir, pl_frq=50)
[ Info: Removing power line noise at 50 Hz and its harmonics
24×6 DataFrame
Row channel power line bandwidth peak 1 frequency peak 2 frequency peak 1 bandwidth peak 2 bandwidth
String Float64 Float64 Float64 Float64 Float64
1 Fp1 3.4 100.0 106.0 0.3 0.2
2 Fp2 2.6 100.0 106.0 0.2 0.1
3 F3 3.4 100.0 106.0 0.3 0.2
4 F4 3.3 100.0 106.0 0.3 0.3
5 C3 3.4 100.0 106.0 0.3 0.2
6 C4 3.3 100.0 106.0 0.2 0.2
7 P3 3.4 100.0 106.0 0.3 0.2
8 P4 3.3 100.0 106.0 0.3 0.1
9 O1 3.4 100.0 106.0 0.3 0.2
10 O2 3.4 100.0 106.0 0.3 0.2
11 F7 3.4 100.0 106.0 0.4 0.2
12 F8 3.3 100.0 106.0 0.2 0.2
13 T3 3.3 100.0 106.0 0.3 0.2
14 T4 3.3 100.0 106.0 0.2 0.1
15 T5 3.3 100.0 106.0 0.3 0.2
16 T6 3.3 100.0 106.0 0.3 0.2
17 A1 3.4 100.0 106.0 0.3 0.2
18 A2 3.3 100.0 106.0 0.3 0.2
19 Fz 2.8 100.0 106.0 0.2 0.1
20 Cz 3.3 100.0 106.0 0.2 0.1
21 Pz 3.4 100.0 106.0 0.3 0.2
22 EOG1 3.1 100.0 106.0 0.2 0.1
23 EOG2 3.3 100.0 106.0 0.3 0.2
24 ECG 2.7 100.0 106.0 0.3 0.2

Plot PSD after filtering:

plot_psd(eeg, ch="Fp1")

Plot spectrogram after filtering:

e10 = epoch(eeg, ep_len=10)
plot_spectrogram(e10, ch="Fp1", ep=2)

(!) This method is based on analyzing a series of IIR notch filters of various bandwidths. The optimal bandwidth is selected using the filtered signal. The variance within the window of peak power-line noise frequency ±5.0 Hz (default value) is calculated for each peak and the IIR notch filter bandwidth that gives the lowest variance is used for noise filtering. Next, the procedure is repeated for other detected peaks of defined prominence (default is 2.0 dB). This is and experimental modeling and requires further validations.

The output data frame may be used to remove peaks manually:

NeuroAnalyzer.filter!(eeg, ch=1, fprototype=:iirnotch, cutoff=50, bw=1.6)
NeuroAnalyzer.filter!(eeg, ch=1, fprototype=:iirnotch, cutoff=100, bw=0.2)
NeuroAnalyzer.filter!(eeg, ch=1, fprototype=:iirnotch, cutoff=106, bw=0.2)