Initialize NeuroAnalyzer
using NeuroAnalyzer
using DataFrames
using StatisticsWhat is power-line noise?
Removing Power-Line Noise
Power-line noise is interference from the alternating current (AC) of the electrical grid, picked up by the electrodes, cables and the amplifier from nearby power lines, electrical devices and lighting, especially when the electrode impedances are high or unequal or the grounding is poor.
- Frequency: 50 Hz (Europe, most of Asia, Africa, Australia) or 60 Hz (North and Central America, parts of South America and Asia).
- Appearance: a sinusoid at the line frequency, often much larger than the EEG, and sometimes its harmonics (100/150 Hz or 120/180 Hz), caused by non-linear loads (e.g. switching power supplies).
- Identification: a narrow peak at 50/60 Hz in the power spectrum; in the time domain, a regular fast oscillation superimposed on the EEG.
Power-line noise overlaps with the gamma band, so it must be removed before analyzing high-frequency activity, and it inflates the amplitude and power measures of all channels.
Synthetic data
Sixty seconds of a synthetic EEG (pink noise with a 15 μV alpha rhythm at 10 Hz) are contaminated with power-line noise of 20 μV at 50 Hz and its harmonic of 5 μV at 100 Hz.
fs = 256
t = collect(0:(1 / fs):(60 - 1 / fs))
n = length(t)
eeg_signal = generate_noise(n, 10; type = :pink) .+ generate_sine(10, t, 15)
line_noise = generate_sine(50, t, 20, 30) .+ generate_sine(100, t, 5, 70)
syn = create_object(; data_type = "eeg")
syn.header.recording[:sampling_rate] = fs
add_channel!(syn;
data = reshape(eeg_signal .+ line_noise, 1, n, 1),
label = "Cz",
type = "eeg",
unit = "μV")
create_time!(syn; fs = fs)The 50 Hz oscillation dominates the signal:
plot(syn; seg = (0, 1), gui = false)
and appears as narrow peaks at 50 and 100 Hz in the power spectrum:
p_noisy = psd(syn; ch = "Cz", db = true)
plot_psd(p_noisy.f, p_noisy.p[1, :, 1];
xlabel = "Frequency [Hz]",
ylabel = "Power [dB]")
Detecting the power-line frequency
detect_powerline() checks which of the two power-line frequencies (50 or 60 Hz) is present. For each channel (and epoch) the prominence of the power spectrum (Welch’s method, dB) at 50 and at 60 Hz over the median power 2–5 Hz away is calculated; the frequency with the larger prominence is selected if it is at least pr dB (default: 3 dB). The frequency detected in most channels is returned, or 0 if there is no power-line noise.
detect_powerline(syn)50.0detect_powerline!() stores the detected frequency in the object header (obj.header.recording[:line_frequency]), which is the default pl_frq of remove_powerline():
detect_powerline!(syn)Removing power-line noise
remove_powerline() removes the power-line noise and its harmonics with zero-phase IIR notch filters whose bandwidth is optimized for every channel and every peak:
- Peak check: the power spectrum (dB) of the channel must show a peak (prominence ≥
prdB, default: 2 dB) at the power-line frequencypl_frq; channels without a peak are not modified (a warning lists them). - Bandwidth optimization: notch filters of bandwidths
q,2q, … 10 Hz (q= 0.1 Hz) are applied to the signal and the bandwidth that gives the flattest spectrum around the line frequency (the lowest variance of the PSD withinpl_frq ± d,d= 5 Hz) is selected: too narrow a notch leaves a residual peak, too wide a notch creates a dip. - Harmonics: steps 1–2 are repeated for every harmonic (
2 × pl_frq,3 × pl_frq, … below the Nyquist frequency) that shows a peak in this channel (bandwidths up to 5 Hz, window± d / 2).
The filters are applied forward and backward (zero phase), so the EEG waveforms are not shifted in time.
syn_f, pl_bw = remove_powerline(syn;
ch = "Cz",
pl_frq = 50)
pl_bw2×3 DataFrame
| Row | channel | frequency | bandwidth |
|---|---|---|---|
| String | Float64 | Float64 | |
| 1 | Cz | 50.0 | 0.8 |
| 2 | Cz | 100.0 | 0.4 |
The returned DataFrame lists the filtered frequencies and the selected notch bandwidths. Frequency response of the selected notch filter for 50 Hz:
plot_filter(fprototype = :iirnotch,
cutoff = 50,
bw = pl_bw[1, :bandwidth],
fs = fs,
gui = false)
Power spectra of the contaminated signal, the filtered signal and the original EEG (without power-line noise):
syn_clean = deepcopy(syn)
syn_clean.data[1, :, 1] = eeg_signal
p_filtered = psd(syn_f; ch = "Cz", db = true)
p_clean = psd(syn_clean; ch = "Cz", db = true)
p_all = vcat(p_noisy.p[1, :, 1]', p_filtered.p[1, :, 1]', p_clean.p[1, :, 1]')
plot_psd(p_noisy.f, p_all;
clabels = ["with line noise", "remove_powerline()", "without line noise"],
xlabel = "Frequency [Hz]",
ylabel = "Power [dB]")
Around 50 Hz:
plot_psd(p_noisy.f, p_all;
clabels = ["with line noise", "remove_powerline()", "without line noise"],
flim = (40, 60),
xlabel = "Frequency [Hz]",
ylabel = "Power [dB]")
Both peaks are removed and the spectrum of the filtered signal follows the spectrum of the original EEG, except for a small dip at the notch frequencies.
Contaminated (blue) and filtered (black) signal:
plot(syn, syn_f;
seg = (10, 11),
gui = false)
RMS difference between the filtered signal and the original EEG (the power-line noise had an RMS of 14.6 μV):
sqrt(mean(abs2, syn_f.data[1, :, 1] .- eeg_signal))0.44673156995775315Comparison with other methods
Power-line noise can also be attenuated with standard filters, but their parameters have to be chosen manually:
# IIR notch filters at 50 and 100 Hz with a fixed bandwidth of 2 Hz
syn_notch = NeuroAnalyzer.filter(syn; ch = "Cz", fprototype = :iirnotch, cutoff = 50, bw = 2)
NeuroAnalyzer.filter!(syn_notch;
ch = "Cz",
fprototype = :iirnotch,
cutoff = 100,
bw = 2)
# FIR band-stop filter 49–51 Hz
syn_bs = NeuroAnalyzer.filter(syn;
ch = "Cz",
fprototype = :fir,
ftype = :bs,
cutoff = (49, 51),
bw = 1)
# FIR low-pass filter at 45 Hz
syn_lp = NeuroAnalyzer.filter(syn;
ch = "Cz",
fprototype = :fir,
ftype = :lp,
cutoff = 45,
bw = 2)The filtered signals are compared with the original EEG: RMS error, power at 50 and 100 Hz and mean power within 45–55 Hz and 60–95 Hz (dB):
function compare(name, obj)
p = psd(obj; ch = "Cz", db = true)
band(f1, f2) = mean(p.p[1, f1 .<= p.f .<= f2, 1])
return (method = name,
rms_error = sqrt(mean(abs2, obj.data[1, :, 1] .- eeg_signal)),
power_50Hz = p.p[1, p.f .== 50, 1][1],
power_100Hz = p.p[1, p.f .== 100, 1][1],
power_45_55Hz = band(45, 55),
power_60_95Hz = band(60, 95))
end
DataFrame([compare("remove_powerline()", syn_f),
compare("notch, bw = 2 Hz", syn_notch),
compare("band-stop 49–51 Hz", syn_bs),
compare("low-pass 45 Hz", syn_lp),
compare("without line noise", syn_clean)])5×6 DataFrame
| Row | method | rms_error | power_50Hz | power_100Hz | power_45_55Hz | power_60_95Hz |
|---|---|---|---|---|---|---|
| String | Float64 | Float64 | Float64 | Float64 | Float64 | |
| 1 | remove_powerline() | 0.446732 | -17.2976 | -20.2716 | -17.6925 | -19.5095 |
| 2 | notch, bw = 2 Hz | 0.338556 | -27.7894 | -31.7316 | -20.1132 | -19.5819 |
| 3 | band-stop 49–51 Hz | 3.74555 | -23.9624 | 9.12791 | -16.339 | -19.4255 |
| 4 | low-pass 45 Hz | 5.0955 | -38.1942 | -43.0889 | -36.4379 | -41.5198 |
| 5 | without line noise | 0.0 | -17.6715 | -20.9231 | -17.4691 | -19.504 |
remove_powerline()removes both peaks and keeps the spectrum closest to the original EEG.- Notch with a fixed bandwidth: effective, but a too wide notch removes EEG power around the line frequency (the power at 50 Hz and within 45–55 Hz is lower than in the original EEG); a too narrow notch leaves residual noise.
- Band-stop filter: removes only the band given, here not the harmonic at 100 Hz, and its transition bands affect the neighbouring frequencies.
- Low-pass filter: removes the noise, but also all EEG activity above the cutoff (gamma band).
Real EEG
This recording contains strong 50 Hz noise:
eeg = import_edf("files/eeg/eeg.edf")
remove_dc!(eeg;
ch = "eeg")detect_powerline(eeg)50.0eeg_f, pl_bw = remove_powerline(eeg;
ch = "eeg",
pl_frq = 50)
first(pl_bw, 6)6×3 DataFrame
| Row | channel | frequency | bandwidth |
|---|---|---|---|
| String | Float64 | Float64 | |
| 1 | Fp1 | 50.0 | 3.4 |
| 2 | Fp1 | 100.0 | 0.3 |
| 3 | Fp2 | 50.0 | 2.6 |
| 4 | Fp2 | 100.0 | 0.2 |
| 5 | F3 | 50.0 | 3.4 |
| 6 | F3 | 100.0 | 0.3 |
The 50 Hz noise is removed from all channels, together with its harmonic at 100 Hz:
p1 = psd(eeg; ch = "Cz", db = true)
p2 = psd(eeg_f; ch = "Cz", db = true)
plot_psd(p1.f, vcat(p1.p[1, :, 1]', p2.p[1, :, 1]');
clabels = ["original", "remove_powerline()"],
xlabel = "Frequency [Hz]",
ylabel = "Power [dB]")
The broad base of the 50 Hz peak (40–60 Hz) is not removed: it reflects amplitude fluctuations of the power-line noise, which a narrow notch filter cannot remove without distorting the EEG. The peak at 106 Hz is not a harmonic of 50 Hz and is left unchanged.
plot(eeg, eeg_f;
ch = "Cz",
seg = (20, 22),
gui = false)
Options
Narrower search window and lower prominence threshold (e.g. for weak power-line noise):
eeg_f, pl_bw = remove_powerline(eeg;
ch = "eeg",
pl_frq = 50,
d = 3,
pr = 1)Finer bandwidth optimization (slower):
eeg_f, pl_bw = remove_powerline(eeg;
ch = "eeg",
pl_frq = 50,
q = 0.05)In-place filtering (returns the DataFrame of the filtered peaks):
pl_bw = remove_powerline!(eeg;
ch = "eeg",
pl_frq = 50)The selected bandwidths can be used to filter other recordings made in the same conditions with the same notch filters:
eeg2 = import_edf("files/eeg/eeg.edf")
for row in eachrow(pl_bw)
NeuroAnalyzer.filter!(eeg2;
ch = row.channel,
fprototype = :iirnotch,
cutoff = row.frequency,
bw = row.bandwidth)
endremove_powerline() and remove_powerline!() must be applied to a continuous (non-epoched) signal.
Tip: Remove power-line noise before other processing (re-referencing, ICA, epoching): the filters work best on long continuous signals.
Tip: The best way to reduce power-line noise is to prevent it: low and equal electrode impedances, a good ground, short cables away from power cords and switching off unnecessary electrical devices.
Other approaches
Regression-based estimation and subtraction (e.g. CleanLine, Zapline) estimates the amplitude and phase of the power-line sinusoid in short sliding windows and subtracts it, without filtering the background EEG:
- Multi-taper decomposition of the signal in a sliding window.
- Regression: estimation of the amplitude and phase of the sinusoidal line noise in the frequency domain.
- Significance testing (Thomson F-test) of the estimated sinusoid.
- Subtraction of the reconstructed line noise from the signal, repeated until the line noise is no longer significant.
EEG signal ◄─────────────────────────────────────┐
│ │
▼ │
┌── for each sliding window ───────────────┐ │
│ multi-taper spectrum │ │
│ line-noise amplitude and phase │ │
│ Thomson F-test │ │
│ if significant: reconstruct sinusoid │ │
└──────────────────────────────────────────┘ │
│ │
▼ │
subtract reconstructed noise from EEG │
│ │
▼ │
any significant noise? ── yes: repeat ──────────────────┘
│
no
▼
doneThese methods are not available in NeuroAnalyzer yet.