Back to the main page.
Bug 3057 - lcmv beamformer source solution mismatch topography
| Status | CLOSED FIXED |
| Reported | 2016-01-31 17:54:00 +0100 |
| Modified | 2019-08-10 12:32:46 +0200 |
| Product: | FieldTrip |
| Component: | inverse |
| Version: | unspecified |
| Hardware: | PC |
| Operating System: | Mac OS |
| Importance: | P5 normal |
| Assigned to: | |
| URL: | |
| Tags: | |
| Depends on: | |
| Blocks: | |
| See also: |
Philipp Ruhnau - 2016-01-31 17:54:24 +0100
Dear ft-developers, I've noticed something strange and I cannot figure out what the issue is (see example screen shots and analysis pipeline aattached). I'm not entirely sure this is a bug, but given the recent problem with localspheres maybe this is related? as can be seen in my screenshot in the attachment, when I beam an erf using the lcmv method with a common filter I get very odd source solutions. I noticed this on the group level (where actually my effect was inverted in polary, and sort of medial). the paradigm is super simple (hit vs. miss of a threshold visual stimulus on the right screen side) with hundreds of trials, thus there should be a quite clear erf for hits, left lateralized. but as i said on group level i get hit < miss and more something medial or, as for the example subject frontal activity. I tried numerous things to name those I still remember: leadfield normalisation, regfactor (0/5/10), creating an average from the mom field (which gives interestingly enough a quite different solution but still an odd one), different grids, fixedori, separate filters, longer covariance time window, solutions based on gradiometers/magnetometers/both, i even calculated the covariance by hand, but still all the same. I think that, if there is a problem and I didn't just kept a typo, then that the problem is somewhere in parts of the code that I do not have the mathematical knowledge to understand.... only observation is that it seems the source reconstructions of both conditions are dominated by the filter. the covariance/data don't seem to do much. I put the data and vol/mri in my dropbox: https://www.dropbox.com/sh/9qi0u1xzbmg6as5/AAAkqOxg2odu7-3wMEKHDs1ua?dl=0 sorry for these large data, but I haven't figured out how to simulate good ERF data for gradiometers.... I was considering writing this to the list first, and apologize if it is a mistake on my side. Best Philipp p.s.: the data in the dropbox contain gradiometers only, to reduce space
Philipp Ruhnau - 2016-01-31 17:55:26 +0100
Created attachment 771 example pipeline and screenshots
Philipp Ruhnau - 2016-01-31 17:59:30 +0100
sorry, I changed my ip during writing this and in result it got submitted twice...
Philipp Ruhnau - 2016-01-31 18:00:59 +0100
*** Bug 3056 has been marked as a duplicate of this bug. ***
Philipp Ruhnau - 2016-02-01 11:40:29 +0100
update: i did also checked different headmodels (singlesphere, localspheres, and i usually use singleshell) - same result more or less though
Philipp Ruhnau - 2016-02-01 19:04:07 +0100
update #2: I can now reconstruct reasonable sources by computing activity from the dipole moments. However, I have to compute the mean of the absolute values (as direction is still present as polarity I assume, which gave me diploe-like solutions for the mom field in the first place), I attach another image with the group level stats based on .avg.pow and .avg.mom field, the pow field still really produces something odd, and I have no clue how that happens, but the mom solution seems to match the topography. I won't do more digging as I found a solution that works for me, but it is a very odd thing. A friend of mine checked my code with some CTF data and there it seemed the final solution was more reasonable, so maybe it has something to do with my sensor type. cheers
Philipp Ruhnau - 2016-02-01 19:04:52 +0100
Created attachment 772 group level stats based on pow and mom field
Jan-Mathijs Schoffelen - 2016-02-01 19:57:01 +0100
Thanks for reporting this Philipp, We'll try and have a look at it soon. I see that Jim has been added (or he added himself) to the CC-list of this bug, so Jim: feel free to chime in :o).
Jan-Mathijs Schoffelen - 2016-02-01 20:18:17 +0100
Hi Philipp, A potential issue I see in your attached script is that line 59 commented out the cfg.removemean = 'no'; for ft_timelockanalysis. It could be that you have tried it with mean non-removal (and that it didn't work), but now I think that this may explain the issue. I quickly reviewed the code in ft_timelockanalysis, which is quite ancient, and I noticed a quirk, in that the cfg.removemean option results in the mean being temporarily subtracted from the data matrix prior to covariance computation, but it is NOT subtracted when computing the average. So the next thing to know would be that the source.avg.pow variable contains number that results from the sensor covariance sandwiched between the spatial filters (i.e. based on a demeaned data), whereas the mom is the sensor data (i.e. tick.avg) left-multiplied with the spatial filter. Taking the abs of the mom (squaring it), and averaging across time should in some occasions yield the same value as the power, but only with a fortuitous combination of cfg-options specified in the call to ft_timelockanalysis that is used for the covariance (and avg) computation. So in your case, since you use a relatively narrow time window for the computation of the covariance (roughly corresponding with the peak in your ERF), in combination with removemean='yes', you are seriously diluting the difference at the level of the covariance matrix, which could lead to a 'polarity change' as observed.
Philipp Ruhnau - 2016-02-01 21:45:49 +0100
Hey Jan-Mathijs, I'm as usual amazed how quickly you find problematic parts, that took me a half day to identify as potential issues... anyway, however, I did check that option (removemean = 'no'), where the default is to remove the mean (isn't that essential for the covariance anyway?), anyhow it didn't change the source solution, which is why I commented it again... what I don't get is why the mean should be removed from the average, but I think that's a different issue (i could check that though, but then again I do trust the filter * data results more these days). oh and I also did this step by hand (averaging and covariance estimation), which yielded the same results to go on: quote: >> Taking the abs of the mom (squaring it), and averaging across time should in some occasions yield the same value as the power, but only with a fortuitous combination of cfg-options specified in the call to ft_timelockanalysis that is used for the covariance (and avg) computation. >> yes I figured, and in the end I am mostly interested in the relationship between conditions not an exact match, what made me believe that it could happen tough was that recently a colleague did that with ctf data and there it matched very well (even numerically)...but still not the issue I guess, what really confuses me is that the solutions change so completely going on : >> So in your case, since you use a relatively narrow time window for the computation of the covariance (roughly corresponding with the peak in your ERF), in combination with removemean='yes', you are seriously diluting the difference at the level of the covariance matrix, which could lead to a 'polarity change' as observed. << maybe, or better: are you sure? we are talking about an absent signal (that is my miss trials) compared to an ERF (hit trials) and in source space the 'no-signal' wins, and that in occipital regions? plus, 100ms in the erp world is a whole age ;) but let's not get into that and just out of curiosity, would this persist even with hundreds of trials? (I have like 300 in each condition) anyway, I did also try the whole thing with a longer covariance window, as I thought that the short window might be the issue (even though I got good results way back when) that is I computed solutions for the first 500ms segments [0 .5] that contain nice ERFs in occipital and also frontal/pre-central areas... still it didn't look anything like the abs mom (again this is a visual ERF, and I get smaller values for the hit (=ERF) in occipital regions compared to the miss (no signal). thank again for looking into this! Philipp
Jim Herring - 2016-02-04 12:40:31 +0100
Hi all, This reminds me a bit of a bug posted by Lin a few months ago (http://bugzilla.fieldtriptoolbox.org/show_bug.cgi?id=2990). There the problem wasn't a bug but rather that calculating the covariance matrix on the single trial level and then averaging made her source reconstruction rather sensetive to induced non-phase locked effects which were overlapping with evoked effects. Although it is not always advised, calculating the covariance matrix on the average gave her a topography that was more similar to her sensor level data.
Philipp Ruhnau - 2016-02-04 13:12:30 +0100
Created attachment 773 lcmv of cov of avg
Philipp Ruhnau - 2016-02-04 13:20:23 +0100
ugnh. ok thanks for the pointer Jim. I tried to check whether there was something similar on the list at some point, but I guess it is a matter of right keywords of course... so I redid the analysis and computed the filters based on the trial-averaged covariance over all trials of both conditions and computed the covariance from the average of trials in each condition that covariance fed into source analysis (and sandwiched...) leads to something very reasonable now indeed! (see recent attachment) important fun facts: - for creation of the filters the average cov over trials is needed (the cov of the avg over all trials delivers something really weird) - the length of my segments has an impact (as suggested by JM). it get's better/more stable with, e.g., a 200 instead of a 100ms segment. one thing I just remembered is that I downsample for space reasons, which for the lcmv pipeline is probably not a good idea (more samples = better cov estimation) ... I still have this feeling that it once worked, but maybe that was just because evoked and induced effects showed in similar directions (which, given the poor variation in my designs over the last few years is not so likely, but still) anyway, thank you guys so much for looking at this. I meanwhile ran some simulations with the code and now that I know the issue I might be able to replicate this on simulated data, which would be interesting. will see. anywho, thank you again Philipp