Back to the main page.

Bug 2632 - Using simbio to create a volume conduction model and to perform source localisation

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
Depends on:
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 = {'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 = {'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; = {'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; = {'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

- 2014-08-06 14:53:46 +0200

Hi Robert, I am sorry for the delay in following up on this. Unfortunately I was busy with another project. Is it possible that there is a bug in the function determining whether points are inside the brain in a volume created with simbio? I used the following script: cfg = []; cfg.elec = elec_aligned4_128; = 'all'; cfg.grid.resolution = 1; cfg.grid.unit = 'cm'; cfg.vol = vol_download; % or: vol_simbio; grid_pos = ft_prepare_sourcemodel(cfg); Where I only changed the line cfg.vol = ... vol_download is the volume downloaded from the website (created with dipoli) and vol_simbio is the volume I created using simbio. The simbio volume seems to be correct which can be seen on the image attached where I plotted one slice of the volume (taking the first points of the hexahedrons). The other two parts of the image show the grids created from the two volumes. It seems fine for dipoli, but not for simbio. Best regards, Wilhelm