[FieldTrip] low-pass filtering

[Diego Lozano-Soldevilla]

Hi Irina,

On the one hand, the smoothing can indeed influence the cluster size. On
the other hand, the "type" of smoothing will influence how sensitive are
you to detect the potential differences too. I would try to illustrate the
reviewer the trade off between the type of differences you may have
(transient vs sustained; regarding the temporal domain) and the filter type
that can be more appropriate to detect the differences to justify your
pipeline. I have no idea about the type of ERF you're studying but I hope
the example is illustrative enough

Look at the attached screenshot. In the first column, I simulated a
transient ERF with added noise (code below) and I ran the non-parametric
cluster-based permutation test that showed some significant and
discontinuous differences. The second column represents the same data but
with a lowpass filter (40 Hz) and I find longer significant clusters. If
you look at the filtered ERFs, the structure is very similar to the noisy,
unfiltered version. The third column is what I understand is requested by
the reviewer. Suddenly, no significant differences are found because the
ground truth simulated difference was too short and weak to survive the
0.25ms averaging. In case of a long and sustained ERF, may be the three
plots would had been more congruent.

I guess that the permutation part of the cluster-based nonparametric test
will take into account the non-independecy in the data resulting from
low-pass filtering. This non-independency will be the same for the original
and the permuted surrogates.

I hope that helps!

Diego



%---------------------------------

fsample = 600;
sigma = 0.5;
ferf = 6.0;
Aerf1 = -0.7;
Aerf2 = -0.2;
tau1 = 0.1;
tau2 = 0.05;
t = -0.8:1/fsample:0.8;
signoise = 0.2;
t0 = 0.05;
t0idx = nearest(t,t0);
ntrials = 100;

data1.label{1}='Oz';
data2.label{1}='Oz';

for k=1:ntrials;
  freq1 = normrnd(ferf,sigma);
  freq2 = normrnd(ferf,sigma);
  for tk = 1:size(t,2);
    phi1 = random('uniform',0,pi,1,1);
    phi2 = random('uniform',0,pi,1,1);
    if tk > t0idx;
      Serf1(1,tk) = Aerf2 * exp(1 - ((t(1,tk)-t0)/tau1)) *
sin(2*pi*freq1*(t(1,tk)-t0)+phi1);
      Serf2(1,tk) = Aerf2 * exp(1 - ((t(1,tk)-t0)/tau2)) *
sin(2*pi*freq2*(t(1,tk)-t0)+phi2);
    end
    Snoise1 = signoise.*randn(size(t));
    Snoise2 = signoise.*randn(size(t));
  end
  data1.trial{k}(1,:) = Serf1 + Snoise1;
  data1.time{k} = t;
  data2.trial{k}(1,:) = Serf2 + Snoise2;
  data2.time{k} = t;
end

cfg = [];
cfg.keeptrials = 'yes';
timelock1 = ft_timelockanalysis(cfg, data1);
timelock2 = ft_timelockanalysis(cfg, data2);

diff = timelock1;
diff.avg = timelock1.avg-timelock2.avg;

figure;plot(t,[timelock1.avg' timelock2.avg' diff.avg']);
legend('condition 1', 'condition 2', 'difference');

cfg = [];
cfg.design = [ 1*ones(1,ntrials) 2*ones(1,ntrials) ];
cfg.ivar = 1;
cfg.method = 'montecarlo';
cfg.statistic = 'indepsamplesT';
cfg.correctm = 'cluster';
cfg.numrandomization = 5000;
stat = ft_timelockstatistics(cfg, timelock1, timelock2);

figure;
subplot(4,1,1); plot(stat.time, stat.stat); ylabel('t-value');
title('raw data: unfiltered')
subplot(4,1,2); plot(diff.time, diff.avg); ylabel('avg1-avg2 (uV)');
subplot(4,1,3); semilogy(stat.time, stat.prob); ylabel('prob');
subplot(4,1,4); plot(stat.time, stat.mask); ylabel('significant');


%% low pass 40 Hz
cfg          = [];
cfg.lpfilter = 'yes';
cfg.lpfreq   = 40;
data1_lp40   = ft_preprocessing(cfg, data1);
data2_lp40   = ft_preprocessing(cfg, data2);

cfg = [];
cfg.keeptrials = 'yes';
timelock1 = ft_timelockanalysis(cfg, data1_lp40);
timelock2 = ft_timelockanalysis(cfg, data2_lp40);

diff = timelock1;
diff.avg = timelock1.avg-timelock2.avg;

figure;plot(t,[timelock1.avg' timelock2.avg' diff.avg']);


cfg = [];
cfg.design = [ 1*ones(1,ntrials) 2*ones(1,ntrials) ];
cfg.ivar = 1;
cfg.method = 'montecarlo';
cfg.statistic = 'indepsamplesT';
cfg.correctm = 'cluster';
cfg.numrandomization = 5000;
stat = ft_timelockstatistics(cfg, timelock1, timelock2);

figure;
subplot(4,1,1); plot(stat.time, stat.stat); ylabel('t-value');
title('lowpass filtering 40Hz')
subplot(4,1,2); plot(diff.time, diff.avg); ylabel('avg1-avg2 (uV)');
subplot(4,1,3); semilogy(stat.time, stat.prob); ylabel('prob');
subplot(4,1,4); plot(stat.time, stat.mask); ylabel('significant');


%% boxcar 0.25ms
cfg = [];
cfg.boxcar = 0.25;
data1_box   = ft_preprocessing(cfg, data1);
data2_box   = ft_preprocessing(cfg, data2);

cfg = [];
cfg.keeptrials = 'yes';
timelock1 = ft_timelockanalysis(cfg, data1_box);
timelock2 = ft_timelockanalysis(cfg, data2_box);

diff = timelock1;
diff.avg = timelock1.avg-timelock2.avg;

figure;plot(t,[timelock1.avg' timelock2.avg' diff.avg']);

cfg = [];
cfg.design = [ 1*ones(1,ntrials) 2*ones(1,ntrials) ];
cfg.ivar = 1;
cfg.method = 'montecarlo';
cfg.statistic = 'indepsamplesT';
cfg.correctm = 'cluster';
cfg.numrandomization = 5000;
stat = ft_timelockstatistics(cfg, timelock1, timelock2);

figure;
subplot(4,1,1); plot(stat.time, stat.stat); ylabel('t-value');
title('boxcar filter 0.25s')
subplot(4,1,2); plot(diff.time, diff.avg); ylabel('avg1-avg2 (uV)');
subplot(4,1,3); semilogy(stat.time, stat.prob); ylabel('prob');
subplot(4,1,4); plot(stat.time, stat.mask); ylabel('significant');





On 26 June 2018 at 17:16, Simanova, I. (Irina) <i.simanova at donders.ru.nl>
wrote:

> Dear experts,
>
> We recently submitted a paper where we use a cluster permutation analysis
> of ERFs (testing conditions exchangeability). The EEG was sampled at 500
> Hz and low-pass filtered at 40 Hz during preprocessing. One of the
> reviewers has indicated that using this temporal resolution for the
> analysis seems unnecessary given the low-pass filter, which makes each time
> point not independent due to smoothing. He further asks: "Do all key
> significant components replicate if the analysis is performed using a mean
> amplitude that is averaged for larger time bins (25 ms bins) instead of
> individual time points?"
>
> I understand that when doing parametric analysis one might want to reduce
> the number of comparisons by averaging ERF samples over time bins. But here
> we solve the MCP with the cluster analysis. But it also seems like the
> reviewer is concerned about non-independecy in the data resulting from
> low-pass filtering. I wanted to check with you: can smoothing time-series
> indeed compromise the cluster-based analysis (e.g. affect the cluster size)?
>
> Thank you,
> Kind regards,
> Irina
>
>
>
>
> _______________________________________________
> fieldtrip mailing list
> https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
> https://doi.org/10.1371/journal.pcbi.1002202
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20180627/98b3c909/attachment-0002.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: Screenshot from 2018-06-27 11-39-54.png
Type: image/png
Size: 87542 bytes
Desc: not available
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20180627/98b3c909/attachment-0002.png>