Back to the main page.
Bug 2632 - Using simbio to create a volume conduction model and to perform source localisation
| Status | ASSIGNED |
| Reported | 2014-07-03 18:07:00 +0200 |
| Modified | 2014-08-06 14:53:46 +0200 |
| Product: | FieldTrip |
| Component: | inverse |
| Version: | unspecified |
| Hardware: | PC |
| Operating System: | Windows |
| Importance: | P5 normal |
| Assigned to: | Robert Oostenveld |
| URL: | |
| Tags: | |
| Depends on: | |
| Blocks: | |
| See also: |
- 2014-07-03 18:07:15 +0200
Created attachment 649 Detailed bug description with samples Dear all, First of all, thank you very much for this toolbox. I started using it a couple of weeks ago and am really enjoying working with it. I am going to modify head models obtained from MRI scans and it seems to me that volumes created with FEM are more suitable than BEM. That way I ended up using simbio. I found two bugs/suggestions for improvement as outlined below. Since I think it's easier to read when the text is formatted nicely, I also attached the same as a word document. Best wishes, Wilhelm System information: Windows 7, 64 bit Matlab R2014a fieldtrip-20140702 First part: I ran the following script to create a leadfield from a volume conduction model obtained with ‘simbio’: mri = ft_read_mri('Subject01.mri'); %from tutorials cfg = []; cfg.dim = mri.dim; mri_resliced = ft_volumereslice(cfg,mri); cfg = []; cfg.output = {'gray','white','csf','skull','scalp'}; segmentedmri_hexahedral = ft_volumesegment(cfg, mri); cfg = []; cfg.shift = 0.3; cfg.method = 'hexahedral'; mesh = ft_prepare_mesh(cfg,segmentedmri_hexahedral); cfg = []; cfg.method = 'simbio'; cfg.conductivity = [0.33 0.14 1.79 0.01 0.43]; vol_simbio = ft_prepare_headmodel(cfg, mesh); cfg = []; cfg.elec = elec_aligned4_128; % Prepared beforehand cfg.channel = {'AFz', 'T8'}; % Reduced for demonstration cfg.vol = vol_simbio; cfg.grid.resolution = 1; cfg.grid.unit = 'cm'; [grid] = ft_prepare_leadfield(cfg); It’s working fine until the last line: Error using ft_inside_vol (line 120) unrecognized volume conductor model Error in ft_prepare_sourcemodel (line 698) inside = ft_inside_vol(grid.pos, vol, 'grad', sens, 'headshape', cfg.headshape, 'inwardshift', cfg.inwardshift); % this returns a boolean vector Error in ft_prepare_leadfield (line 158) grid = ft_prepare_sourcemodel(tmpcfg); Error in THeadModel_Hexahedral (line 141) [grid] = ft_prepare_leadfield(cfg); If I understand it correctly, only a method to determine the points of the grid that are inside the brain is missing. I think a found a way to work around it by creating a grid from the same MRI image using a BEM method, where vol_download is the volume conduction model from the tutorials created with ‘dipoli’: cfg = []; cfg.elec = elec_aligned4_128; cfg.vol = vol_download; % from tutorials cfg.channel = {'AFz', 'T8'}; cfg.grid.resolution = 1; cfg.grid.unit = 'cm'; [grid_BEM] = ft_prepare_leadfield(cfg); And afterwards run ft_prepare_leadfield(cfg) with cfg.grid = grid_BEM Although it seems to work this way, it would be very nice, if it worked without needing to compute a grid in two ways every time again. ? Second part: Still with the model created above, now doing source analysis. One step which takes a long time to compute when the number of channels is large seems to be done multiple times: In ft_prepare_leadfield(cfg) (see above) and again every time (for contrast plots it’s three times for frequency analysis) in ft_sourceanalysis(cfg,freq_sound) even when vol and sens are supplied to the function as cfg.vol and cfg.sens: cfg = []; cfg.method = 'mtmfft'; cfg.output = 'powandcsd'; cfg.tapsmofrq = 2; cfg.foilim = [36 36]; freq_sound = ft_freqanalysis(cfg, data_sound_m); [vol, sens] = ft_prepare_vol_sens(vol_simbio, elecs); cfg = []; cfg.method = 'dics'; cfg.frequency = 36; cfg.grid = grid; cfg.sens = sens; cfg.vol = vol; cfg.dics.projectnoise = 'yes'; cfg.dics.lambda = '5%'; cfg.dics.keepfilter = 'yes'; cfg.dics.realfilter = 'yes'; cfg.elec = elecs; source_sound_nocon = ft_sourceanalysis(cfg, freq_sound); The output is: converting units from 'mm' to 'cm' converting units from 'mm' to 'cm' Find electrode positions... Calculate transfer matrix... Electrode 2 of 2 Scaling stiffnes matrix... Preconditioning... Finding startvector... Solving equation system... Comes from: sb_transfer (line 32) ft_prepare_vol_sens (line 500) prepare_headmodel (line 94) ft_prepare_leadfield (line 134) And: sb_transfer (line 32) ft_prepare_vol_sens (line 500) prepare_headmodel (line 94) ft_sourceanalysis (line 291) Could FieldTrip check, whether the correct vol and sens are already computed and use those instead?
Robert Oostenveld - 2014-07-04 14:44:06 +0200
Hi Wilhelm, You are right. 1) the inside/outside grid point detection is not yet in place. Your workaround with the BEM model is just fine regarding what it should achieve. I think you could even replace ft_prepare_leadfield with ft_prepare_sourcemodel, which makes it faster. 2) the inefficiency of the code w.r.t. multiple times setting up the same matrix is something I recently discussed with the simbio team. Your suggestion is also what we had identified as the desired solution. Let me first devote some time on your issue #1, as that might be relatively simple to solve. Robert
Robert Oostenveld - 2014-07-06 11:19:09 +0200
Regarding 1: I made an implementation of the inside/outside detection. mac011> svn commit Sending forward/ft_inside_vol.m Adding test/test_tutorial_fem.m Transmitting file data .. Committed revision 9702.
- 2014-07-08 18:52:15 +0200
Hi Robert, Thanks a lot for the quick implementation. I tried the newest version with: cfg = []; cfg.elec = elecs; cfg.vol = vol_simbio; cfg.channel = {'AFz', 'T8'}; cfg.grid.resolution = 1; cfg.grid.unit = 'cm'; cfg.sens = sens; grid = ft_prepare_sourcemodel(cfg); and cfg = []; cfg.elec = elecs; cfg.vol = vol_simbio; cfg.channel = {'AFz', 'T8'}; cfg.grid.resolution = 1; cfg.grid.unit = 'cm'; [grid] = ft_prepare_leadfield(cfg); The first case works fine, I get a grid and the output: 2081 dipoles inside, 5443 dipoles outside brain making tight grid 2081 dipoles inside, 2099 dipoles outside brain In the second case, no points are detected as being inside. grid.inside is empty. However, when I add cfg.sens = sens3 and cfg.vol = vol_simbio3 where sens3 and vol_simbio3 were calculated using [vol_simbio2, sens] = ft_prepare_vol_sens(vol_simbio, elecs) and vol_simbio3 (sens3) being vol_simbio2 (sens) with coordinates in cm and replace in ft_prepare_leadfield (line 134): [vol, sens, cfg] = prepare_headmodel(cfg, data); with: vol = cfg.vol; sens = cfg.sens; I obtain a grid (2069 inside and 2111 outside) and a leadfield. So I guess there is some problem with prepare_headmodel() being used in ft_prepare_leadfield()?
Robert Oostenveld - 2014-07-09 16:47:52 +0200
(In reply to wv12 from comment #3) please do an explicit call to ft_prepare_sourcemodel prior to ft_prepare_leadfield. It will give you better control of the source moddel construction. After cfg = ... grid_pos = ft_prepare_sourcemodel(cfg) you will want to pass it into ft_prepare_leadfield like this ... cfg.grid = grid_pos; grid_lf = ft_prepare_leadfield(cfg) If you only specify cfg.grid.resolution, it will use the span (or range) of the sensor array (electrodes) to make an initial grid with your resolution and will then determine which of those grid points are inside the brain. With 2 electrodes that is likely not to result in a nice grid. You can specify in more detail cfg.grid.xrange etc. Please do check with ft_plot_mesh (or plot3) that the grid points make sense and that inside and outside points (give them different colours) make sense.
- 2014-07-10 14:29:54 +0200
Hi Robert, That's very good to know, thanks! I'll do it from now on. (Normally I use more than two electrodes, only wanted to keep the computation as short as possible.)
- 2014-08-06 14:52:17 +0200
Created attachment 658 grids created from dipoli and simbio