Repo contains source code for creating and filtering EEG data from periodic, non-sinusoidal and non-stationary tCS artifacts using weighted comb filters.
Includes also code for artifact removal using adaptive DFT and adaptive PCA, and for simulation of tACS recordings. You might want to read the accompanying pdf for more details about the differences and similarities to SMA and similar template averaging methods, and performance evaluation.
This module is shared under a X11 license. Its development was supported by the BMBF: FKZ 13GW0119.
Upper Limb Bipolar ECG recording during 11 Hz tACS |
|
---|---|
Recover the physiological signal (which is ~120dB weaker than tACS) |
Artifact can be non-stationary and non-sinusoidal, but is required to be periodic. Comb filters natively support only frequencies which are integer divisibles of the sampling frequency. When artacs.kernel.run is used, the signal is automatically resampled, to circumvent this limitation. Note that the method still requires integer frequencies.
By default, the kernel is symmetric and weights are based empirically on the artifacts periodic autocorrelation.
% Add package to path
addpath('.\src\')
% Run a kernel filter, with automatic weights as default
% based on 10 neighbouring periods
% filter a frequency of 11 Hz
% for a signal recorded at 1000 Hz
NumberPeriods = 10;
Freq = 11;
Fs = 1000;
filtered_signal = artacs.kernel.run(Signal,Freq,NumberPeriods,Fs)
% You can also create own kernels
% e.g a symmetric uniform kernel
wfun = 'uniform'; % wfun can be 'uniform', 'linear', 'exp', 'gauss', 'automatic'
symflag = 'symmetric'; % symflag can be 'causal', 'symmetric', 'right', or 'piecewise'.
filtered_signal = artacs.kernel.run(Signal,Freq,NumberPeriods,Fs,symflag,wfun)
Piecewise filtering splits the signal in half (or at a specified latency), filters the left half forward and the right half backwards, and fuses the signal again. Comb filters can cause echos (see causal uniform filter). Piecewise filtering can help to reduce the echo of an ERP in the filtered signal if the ERP latency is known (within roughly half the artifacts period).
filtered_signal = artacs.kernel.run(Signal,Freq,NumberPeriods,Fs,'piecewise','gauss')
Artifact can be non-stationary and non-sinusoidal, but is required to be periodic. Estimates the time-course of artifact amplitude modulation across periods, and removes its prinical components until artifact power is suppressed below the level of neighbouring frequencies. Works only for integer frequencies.
% based on adaptive stepwise removal of prinicipal components
% of the artifact amplitude modulation.
filtered_signal = artacs.template.stepwise(Signal,Freq,Fs)
Artifact is assummed to be sinusoidal, and periodic, but can be non-stationary. The sinusoidal artifacts amplitude is estimated either local (i.e. for the past NumberPeriods periods) or for the complete signal duration and removed. Works natively for any real frequency.
% based on adaptive causal dft
filtered_signal = artacs.dft.causal(Signal,Freq,Fs,NumberPeriods)
% filter harmonics by using a vector for Freq
filtered_signal = artacs.dft.causal(Signal,[1:4]*Freq,Fs,NumberPeriods)
% use fft/ifft on the complete trial duration
filtered_signal = artacs.dft.complete(Signal,Freq,Fs)
Recovery (as R² between filtered and stim-free ECG) for various filtering approaches |
On inspecting kernels in time and frequency domain
- Allow filtering for non-integer frequencies
- Translate the code to Python3