[FieldTrip] Padding with mtmfft and mtmconvol

Vladimir Litvak litvak.vladimir at gmail.com
Thu Sep 21 12:34:23 CEST 2017


Another thing that I noticed is that in the mtmconvol case padding is added
to the entire trial, not to the short window over which FFT is actually
computed. This might be because you actually use a wavelet which moves
along the data (I didn't check that). Anyhow right now this doesn't make
much difference because padding doesn't affect mtmconvol in such a dramatic
way as mtmfft. However, if you do allow specifying the units as power
rather than density then the way things are now mtmconvol and mtmfft with
the same padding would not be equivalent.

Vladimir

On Thu, Sep 21, 2017 at 11:29 AM, Vladimir Litvak <litvak.vladimir at gmail.com
> wrote:

> Hi Jan-Mathijs,
>
> Yes, you are right about Robert's example. But if you do:
>
> pad = 10;
>
> fsample = 1000;
> time = (1:1000)/fsample;
> dat = randn(size(time));
>
> [spectrum1,ntaper1,freqoi1] = ft_specest_mtmfft(dat, time, 'taper',
> 'hanning', 'pad', pad);
>
> power1 = abs(spectrum1).^2;
> power1 = squeeze(power1);
>
> [spectrum2,ntaper2,freqoi2,timeoi2] = ft_specest_mtmconvol(dat, time,
> 'taper', 'hanning', 'timeoi', mean(time), 'timwin', 0.99, 'pad', pad);
>
> power2 = abs(spectrum2).^2;
> power2 = squeeze(power2);
>
> figure
> plot(freqoi1, power1);
> hold on
> plot(freqoi2, power2, 'r');
>
>
> You will see the problem that I'm talking about. We discussed with Robert
> yesterday and this is indeed 'a feature' which has to do with the fact that
> the outputs of mtmfft and mtmconvol have different units. The former is
> spectral density whereas the latter is spectral power.
>
> Here is what Robert wrote me:
>
>
> the units of computations (also here) are a known and long-standing issue.
> I know for a long time that the two have different scaling, but did not
> think about it for a long time. I recall something like this: To compare
> TFRs over frequencies, you don't want the bandwidth to affect the estimate.
> Shorter wavelets have a larger bandwidth, hence the 1/Hz would affect
> those. E.g. imagine a 10Hz and a 20Hz sine wave, and do a TFR with
> conventional wavelets: at 20Hz the wavelet is 2x shorter, so the spectral
> resolution over which the signal(and noise) spreads is different. If you
> were to compute the TFR in V^2/Hz, the same V at 20Hz would have a
> different value, because the length of the wavelet affects the 1/Hz.
> something related (but nevertheless different) applies to the mtmfft: if
> you want to estimate broadband activity in a window of 1 second or a window
> of 2 seconds, you would get different spectral resolutions. The nyquist is
> the same, but the power gets distributed over more bins between 0 and
> Fnyquist/2. That would cause the values to appear smaller in the 2-s case.
> Hence we compute spectral density, which somehow normalizes for this. I
> never found a really clear explanation, but google got me this
> https://dsp.stackexchange.com/questions/33957/what-is-the-
> difference-between-the-psd-and-the-power-spectrum
> what confuses me is that power (or variance) is already normalized, i.e.
> sum of squared values divided by N. So we have energy (which increases with
> length), power (which does not increase with length), and power density
>
>
> So one issue is that most people don't know about this including me and
> possibly you. I think a good solution would be to add an option to specify
> the output units for all the methods as there might be quite subtle
> considerations for choosing one over the other as Robert suggests.
>
> Vladimir
>
> On Thu, Sep 21, 2017 at 10:53 AM, Schoffelen, J.M. (Jan Mathijs) <
> jan.schoffelen at donders.ru.nl> wrote:
>
>> Hi to all who’s reading along,
>>
>> Perhaps the two cases will become more similar once the ‘timwin’ is
>> increased in length for the mtmconvol case….
>>
>> Best wishes,
>>
>> JM
>>
>> On 20 Sep 2017, at 16:26, Robert Oostenveld <r.oostenveld at donders.ru.nl>
>> wrote:
>>
>> Hi Vladimir,
>>
>> I suggest that you first start with a simpler case, like this
>>
>> fsample = 1000;
>> time = (1:1000)/fsample;
>> dat = randn(size(time));
>>
>> [spectrum1,ntaper1,freqoi1] = ft_specest_mtmfft(dat, time, 'taper',
>> 'hanning');
>>
>> power1 = abs(spectrum1).^2;
>> power1 = squeeze(power1);
>>
>> [spectrum2,ntaper2,freqoi2,timeoi2] = ft_specest_mtmconvol(dat, time,
>> 'taper', 'hanning', 'timeoi', mean(time), 'timwin', 0.5);
>>
>> power2 = abs(spectrum2).^2;
>> power2 = squeeze(power2);
>>
>> figure
>> plot(freqoi1, power1);
>> hold on
>> plot(freqoi2, power2, 'r');
>>
>> Note that these are not the same (albeit similar), which I had expected…
>>
>> best
>> Robert
>>
>>
>>
>> On 20 Sep 2017, at 12:56, Vladimir Litvak <litvak.vladimir at gmail.com>
>> wrote:
>>
>> Dear Fieldtrippers,
>>
>> I'm looking into an issue of one of SPM users who gets different results
>> when doing TF decomposition compared to computing a spectrum for the same
>> time window. I'm not sure I got to the bottom of it yet but one thing I
>> found is that ft_specest_mtmfft and ft_specest_mtmconvol are affected
>> differently by increasing padding. For short padding the results are
>> similar but with increasing padding there are differences both in offset of
>> the spectrum and its overall shape. See attached images where the top one
>> shows original spectra and the bottom one aligns the lowermost bin to zero.
>>
>> Is this a bug or a feature?
>>
>> Below is the script that produces these plots. I could provide the data
>> as well but this could probably be reproduced with any data.
>>
>> Thanks,
>>
>> Vladimir
>>
>>
>> -------------------------------------
>>
>> pad = 0.5;%1%10
>>
>>
>> freqoi = 5:45;
>> timwin = 0.4+0*freqoi;
>>
>> [spectrum,ntaper,freqoi,timeoi] = ft_specest_mtmconvol(data, time,
>> 'taper', 'hanning', 'timeoi', 1.2, 'freqoi', freqoi,...
>>     'timwin', timwin, 'tapsmofrq', 1, 'verbose', 0, 'pad', pad);
>>
>> s1 = squeeze(mean(mean(abs(spectrum), 4), 2));
>>
>> figure;
>> subplot(2,1,1)
>> plot(freqoi, s1);
>> subplot(2,1,2);
>> plot(freqoi, s1-s1(1));
>> %%
>> [spectrum,ntaper,freqoi] = ft_specest_mtmfft(data, time, 'taper',
>> 'hanning', 'freqoi', freqoi,...
>>     'tapsmofrq', 1, 'verbose', 0, 'pad', pad);
>>
>> s2 = squeeze(mean(mean(abs(spectrum), 2), 1));
>>
>> subplot(2,1,1)
>> hold on
>> plot(freqoi, s2, 'r');
>> subplot(2,1,2)
>> hold on
>> plot(freqoi, s2-s2(1), 'r');
>> <pad0_5.png><pad1.png><pad10.png>___________________________
>> ____________________
>> fieldtrip mailing list
>> fieldtrip at donders.ru.nl
>> https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>>
>>
>> _______________________________________________
>> fieldtrip mailing list
>> fieldtrip at donders.ru.nl
>> https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>>
>>
>>
>> _______________________________________________
>> fieldtrip mailing list
>> fieldtrip at donders.ru.nl
>> https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20170921/b05b91a3/attachment.html>


More information about the fieldtrip mailing list