mne-tools / mne-python

MNE: Magnetoencephalography (MEG) and Electroencephalography (EEG) in Python
https://mne.tools
BSD 3-Clause "New" or "Revised" License
2.72k stars 1.32k forks source link

BUG in DICS beamformer #1549

Closed lukasbreuer closed 10 years ago

lukasbreuer commented 10 years ago

Dear all,

I think in the function "_apply_dics" in the module "_dics.py" (lines 25-140) there is a bug. The spatial filter is computed from the cross-spectral density matrices, i.e. it is a filter in the frequency domain. However, it is applied to the data (line 127), which are in the sensor space, i.e. time domain. I think before applying the filter the data must also be transformed into the frequency domain, e.g. using FFT. As far as I understand the paper from Gross et al. (2001) they also transformed the data in the frequency domain before applying the filter, but maybe I misunderstood the paper.

Best, Lukas

rgoj commented 10 years ago

Hi Lukas,

The cross-spectral density matrix passed to _apply_dics is only calculated for a single frequency (or summed across a range of frequencies), see example line 62. So we're only calculating the spatial filter for that particular frequency [range] and it's only a spatial filter, there's no frequency dimension.

Perhaps what you write about could be achieved by calculating the spatial filters separately for each frequency (e.g. by setting fsum to False for compute_epochs_csd and then using each of the CSDs and resulting spatial filters separately with data at particular frequencies? You'll notice in the example we don't even filter the data using the same fmin and fmax as the CSD, which would be good to do. On the other hand, the frequency range is chosen to contain most of the power in the example, so it's not a big problem for the example itself. But if you wanted to obtain the time courses at different frequencies, you'd have to create a CSD matrix for each frequency and pass it with FFT-filtered data to _apply_dics.

Thanks for bringing this up, the beamformer module and dics code aren't very well tested, so every little bit helps. Does the above address your concerns?

Best, Roman

lukasbreuer commented 10 years ago

Hi Roman,

Thanks for the fast response! I have just reread the paper of Gross and checked your code again. You are right. Your implementation is correct. I think my misinterpretation of your code is maybe based on the name "_apply_dics". What is actually done in the method is to perform something like LCMV beamformer (Gross et al. also used LCMV to show amplitude as function of time courses). Of course you used DICS to create a spatial filter for a particular frequency, but applying it as spatial filter is LCMV beamforming. The main method one has to use for getting the 'real' DICS results as reported in Gross et al. is "dics_source_power". Summarizing I have to thank you a lot! With your explanation you helped me to get a better understand of DICS beamformer and how it is implemented in the mne package.

Best, Lukas

larsoner commented 10 years ago

@lukasbreuer as you work with the code, also let us know if you think of any enhancements / fixes you'd be interested in working on :)

rgoj commented 10 years ago

@lukasbreuer happy to be of help :)