# [FieldTrip] Frequency shift after source reconstruction

Luca La Fisca 534843 at umons.ac.be
Wed May 27 10:44:42 CEST 2020

```Dear Jan-Mathijs,

It is exactly the problem!

Thank you very much :)

Have a nice day,

*Luca La Fisca*
*PhD Student*
Service d'Information, Signal et Intelligence Artificielle

Boulevard Dolez, 31

7000 Mons

+32 (0)65/37.40.83

*Luca.LAFISCA at umons.ac.be <Luca.LAFISCA at umons.ac.be>*

Garanti
sans virus. www.avast.com

On Wed, May 27, 2020 at 10:15 AM Schoffelen, J.M. (Jan Mathijs) <
jan.schoffelen at donders.ru.nl> wrote:

> Caro Luca,
>
> I haven’t read your code in detail, but I suspect that you perform a
> spectral analysis on the source.avg.pow, which is the magnitude(-squared)
> version of the signal’s time course.
> As a consequence of this, your frequency ‘doubles’. Consider a sine wave
> with a given frequency (i.e. a period of 1/f), say 1 Hz. If you take that
> absolute(-squared) of that sine wave, you get a signal with ’two bumps’ per
> second, which reflects a periodicity of 2 Hz.
>
> You should work with the ‘mom’ field at the source level.
>
> Best wishes,
>
> Jan-Mathijs
>
>
> On 27 May 2020, at 09:33, Luca La Fisca <534843 at umons.ac.be> wrote:
>
> Dear community,
>
> I'm Luca La Fisca, PhD student in Belgium and I'm validating the source
> reconstruction algorithm using "ft_dipolesimulation".
>
> After having performed a frequency analysis, I noticed the reconstructed
> sources have twice the initial frequencies as shown in the following
> figures (initial dipole frequencies = 10 and 85 Hz) :
>
> Timelock:
> <image.png>
>
> Reconstructed sources:
> <image.png>
>
> Here is the code used to get the timelock signal:
> fs = 2048;
> fs_down = 512;
>
> % generate EEG pseudo-signal
> cfg      = [];
> cfg.headmodel = vol;
> cfg.elec = elec_aligned;
> cfg.channel = elec_aligned.label(1:64);
>
> cfg.dip.pos = [
>    -1 65 36           % dipole left
>    1 -62 43          % dipole right
>    ];
>
> unif_dip = 3^(1/2)/3;
> cfg.dip.mom = ...
>   [ unif_dip unif_dip unif_dip 0 0 0 ]' + ...
>   [ 0 0 0 unif_dip unif_dip unif_dip ]';
>
> event = ft_read_event(EEG_FILE);
> event = event(3:end);
> signal1 = zeros(length(eeg_data.trial));
> signal2 = zeros(length(eeg_data.trial));
> val = [event.value];
> for s = drange([event(val==10).sample, event(val==11).sample,
> event(val==12).sample])
>     time = (0:ceil(fs*0.5))/fs;
>     signal1(s:s+ceil(fs*0.5)) = 10*sin(10*time*2*pi);
>     signal2(s:s+ceil(fs*0.5)) = 30*cos(85*time*2*pi);
> end
>
> cfg.dip.signal = {[signal1; signal2]};
> cfg.fsample = fs;
> raw1 = ft_dipolesimulation(cfg);
>
> % Timelock analysis
> % define trials
> cfg = [];
> cfg.dataset = EEG_FILE;
> cfg.trialdef.eventtype  = 'STATUS';
> cfg.trialdef.eventvalue = [10 11 12];
> cfg.trialdef.prestim    = 0.5;
> cfg.trialdef.poststim   = 1;
> cfg = ft_definetrial(cfg);
> data_trial = ft_redefinetrial(cfg, raw1);
>
> % downsample
> cfg = [];
> cfg.resamplefs = fs_down;
> data_trial = ft_resampledata(cfg,data_trial);
>
> cfg = [];
> timelock = ft_timelockanalysis(cfg, data_trial);
>
>
> Here the code related to the source reconstruction:
> % Estimating the leadfield
> cfg = [];
> cfg.sourcemodel = sourcemodel;    %% where are the sources?
> cfg.headmodel   = vol;      %% how do currents spread?
> cfg.elec        = elec_aligned; %% where are the sensors?
> cfg.channel = elec_aligned.label(1:64);
>
> % Source reconstruction
> cfg                     = [];
> cfg.method              = 'mne';
> cfg.sourcemodel         = sourcemodel_and_leadfield;        %the
> cfg.headmodel           = vol;                      %the head model
> cfg.elec                = elec_aligned;             %the electrodes
> cfg.channel             = elec_aligned.label(1:64); %the useful channels
> cfg.mne.prewhiten       = 'yes';                    %prewhiten data
> cfg.mne.lambda          = 0.01;                     %regularisation
> parameter
> cfg.mne.scalesourcecov  = 'yes';                %scaling the source
> covariance matrix
> minimum_norm_eeg        = ft_sourceanalysis(cfg,timelock);
>
> Here the code to perform the frequency analysis:
> tmp_pow = minimum_norm_eeg.avg.pow;
> pow_fft = fft(tmp_pow'); %freq analysis of reconstructed sources
> % pow_fft = fft(timelock.avg'); %freq analysis of timelock signal
> % (un)comment depending on the desired one
>
> L = length(timelock.avg);
> f = fs_down*(0:(L/2))/L;
>
> P2 = abs(pow_fft/L);
> P1 = P2(1:L/2+1,:);
> P1 = P1.*f';
>
> figure()
> plot(f,P1)
> title('Single-Sided Amplitude Spectrum of Reconstructed Sources')
> xlabel('f (Hz)')
> ylabel('|P1(f)|')
>
> Do you have any idea of what is responsible to this double frequency?
> Do you think it is a general issue?
> Is it a major issue leading to wrong analysis of the reconstructed signal
> (e.g. connectivity analysis)?
>
> I thank you in advance.
>
> Best regards,
>
> *Luca La Fisca*
> *PhD Student*
> Service d'Information, Signal et Intelligence Artificielle
>
> Boulevard Dolez, 31
>
> 7000 Mons
>
> +32 (0)65/37.40.83
>
> *Luca.LAFISCA at umons.ac.be <Luca.LAFISCA at umons.ac.be>*
>
>
>
>
>
> sans virus. www.avast.com
> _______________________________________________
> fieldtrip mailing list
> https://mailman.science.ru.nl/mailman/listinfo/fieldtrip
> https://doi.org/10.1371/journal.pcbi.1002202
>
>
> _______________________________________________
> 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/20200527/96972713/attachment.htm>
```