DC offset removal filter

Sometimes slow drifts or DC drifts mask low frequencies that could be of interest, especially when the analysis included fourier-based methods. For example the Contingent Negative Variation (CNV) signal is a very low frequency signal that should be preserved while avoiding to include in the analysis the slow drifts inherent sometimes with EEG recordings. These filters were made for a CNV EEG analysis. Sometimes, in other kinds of recordings an initial drift from baseline exists at the beginning of recordings. For those cases, a specialized DC offset removal filter is necessary for not cutting out parts of data that could be of interest.

Two filters, two solutions are proposed below:

(the examples assume your data are already loaded in EEGLAB)

Solution a: Using an IIR low-pass filter

General idea: Application of a low pass FIR filter (Roll off 0.001 to 0.05 Hz, butter, (IIR), 6 dB attenuation in the stop band), then removal from the signal the filtered lowpass signal

Filtered_low_pass=filterLP(signal)

Final_signal=signal –filterLP(signal).

This is depicted in the following figure for a representative channel.

Figure 1: channel 1 of dataset 101b, raw dataset (blue), the Filtered_low_pass(red) and the final signal Filtered_signal (green). A sinusoid of 0.01 Hz is also plotted for comparison.

If we apply this procedure to 8 selected channels we get the following picture:

Figure 2: filtered (blue) and raw (green).

The lowpass filter characteristics are

h=fdesign.lowpass('Fp,Fst,Ap,Ast',0.001,0.05,1,6,250);

d=design(h,'butter');

fvtool(d)

The program:

DC_offset_removal_21_10_2011_a.m

%% DC offset removal. 21.10.2011 Maria L. Stavrinou

%a)solution 1

%% 1. load the EEG.data and save it as another set

data= EEG.data; % nchan x timeduration_all

%% 2 Low pass filter

h=fdesign.lowpass('Fp,Fst,Ap,Ast',0.001,0.05,1,6,250);

d=design(h,'butter');

fvtool(d)

% check on a single channel

ch1=double(data(1,:));

y=filter(d, ch1);%y=filtfilt(d.Numerator,1,ch1);

figure; plot(ch1); hold on; plot(y, 'r')

hold on; plot((ch1-y), 'g')

%% remove DC from all channels

nchan=size(data,1);

for k=1:nchan

temp=double(data(k,:));

%temp_filt=filtfilt(d.Numerator,1,temp);

temp_filt=filter(d,temp);

data_filt(k,:)=-data(k,:)+temp_filt;

%data_filt(k,:)=(data_filt(k,:)-mean(data_filt(k,:)))

clear temp temp_filt

end

%% 4 check on one file & decide how to cut the file.

% %% see an example file

ch12=data_filt([5 6 7 12 62 72 106 129],:);

figure(100);

nn=size(ch12,1);

for kk=1:nn

a=squeeze(ch12(kk,:));

plot(a);hold on;

clear a

end

hold on;

raw_ch12=data([5 6 7 12 62 72 106 129],:);

for kk=1:nn

a=squeeze(raw_ch12(kk,:));

plot(a, 'g');hold on;legen1=legend('green=raw blue=filtered')

clear a

end

disp('Write down how to cut from the start in datapoints')

%% 5 reload the filtered back to the EEGLAB structure

EEG.data=data_filt;

Solution b: Using an FIR filter without phase shift (time delays)

General idea: Application of a low pass FIR filter (Roll off 0.001 to 0.05 Hz, Equiripple, (FIR), 8 dB attenuation in the stop band), then removal from the signal the filtered lowpass signal

Filtered_low_pass=filterLP(signal)

Final_signal=signal –filterLP(signal).

Filtered_signal=Filtered_signal-mean(Filtered_signal); (for every channel, its own mean is subtracted)

This is depicted in the following figure for a representative channel.

Figure 1: channel 1 of dataset 101b, raw dataset (blue), the Filtered_low_pass(red) and the final signal Filtered_signal (green).

If we apply this procedure to 8 selected channels we get the following picture:

Figure 2: filtered (blue) and raw (green).

The lowpass filter characteristics are

h=fdesign.lowpass('Fp,Fst,Ap,Ast',0.001,0.05,1,8,250);

d=design(h,'equiripple');

fvtool(d)

The program:

DC_offset_removal_21_10_2011_b.m

%% DC offset removal. 21.10.2011 Maria L. Stavrinou

%b)solution 2 %%

%% load eeg data

data= EEG.data; % nchan x timeduration_all

%% 2 Low pass filter

h=fdesign.lowpass('Fp,Fst,Ap,Ast',0.001,0.05,1,8,250);

d=design(h,'equiripple');

fvtool(d)

% check on a single channel

ch1=double(data(1,:));

y=filtfilt(d.Numerator,1,ch1);

figure; plot(ch1); hold on; plot(y, 'r')

hold on; plot((ch1-y), 'g')

%% 3 remove DC from all channels

nchan=size(data,1);

for k=1:nchan

temp=double(data(k,:));

temp_filt=filtfilt(d.Numerator,1,temp);

%temp_filt=filter(d,temp);

data_filt(k,:)=data(k,:)-temp_filt;

data_filt(k,:)=(data_filt(k,:)-mean(squeeze(data_filt(k,:))));

clear temp temp_filt

end

%% 4 check on channels of interest

ch12=data_filt([5 6 7 12 62 72 106 129],:);

figure(100);

nn=size(ch12,1);

for kk=1:nn

a=squeeze(ch12(kk,:));

plot(a);hold on;

clear a

end

hold on;

raw_ch12=data([5 6 7 12 62 72 106 129],:);

for kk=1:nn

a=squeeze(raw_ch12(kk,:));

plot(a, 'g');hold on;legen1=legend('green=raw blue=filtered')

clear a

end

clear nn

disp('Write down how to cut from the start in datapoints')

%% 5 reload the filtered back to the EEGLAB structure

EEG.data=data_filt;

Home