Difference in the measured and predicted evoked of Dipole

When I computed the forward solution of my fitted dipole and compared the predicted evoked with the measured evoked, I saw a huge difference in topography. What could be the reason for it? Is it faulty head model or something else?

the left plot (labelled ā€œmeasured fieldā€) and the right plot (labelled ā€œdifferenceā€) look identical to me. Are you sure youā€™re actually plotting the difference correctly? Looking at the left and middle plots, Iā€™m not sure Iā€™d call that a ā€œhugeā€ difference, assuming the colormaps are on the same scale.

Hi @ drammock,

I used the code in the MNE tutorial. As youā€™ve mentioned, I am not sure why the measured field and the one labelled difference look nearly the same.


best_idx = np.argmax(dip.gof)
best_time = dip.times[best_idx]
print('Highest GOF %0.1f%% at t=%0.1f ms with confidence volume %0.1f cm^3'
      % (dip.gof[best_idx], best_time * 1000,
         dip.conf['vol'][best_idx] * 100 ** 3))
# remember to create a subplot for the colorbar
fig, axes = plt.subplots(nrows=1, ncols=4, figsize=[10., 3.4],
                         gridspec_kw=dict(width_ratios=[1, 1, 1, 0.1],
                                          top=0.85))
vmin, vmax = -400, 400  # make sure each plot has same colour range

# first plot the topography at the time of the best fitting (single) dipole
plot_params = dict(times=best_time, ch_type='mag', outlines='skirt',
                   colorbar=False, time_unit='s')
evoked.plot_topomap(time_format='Measured field', axes=axes[0], **plot_params)

# compare this to the predicted field
pred_evoked.plot_topomap(time_format='Predicted field', axes=axes[1],
                         **plot_params)

# Subtract predicted from measured data (apply equal weights)
diff = combine_evoked([evoked, pred_evoked], weights=[1, -1])
plot_params['colorbar'] = True
diff.plot_topomap(time_format='Difference', axes=axes[2:], **plot_params)
fig.suptitle('Comparison of measured and predicted fields '
             'at {:.0f} ms'.format(best_time * 1000.), fontsize=16)
fig.tight_layout()

This is what my BEM looks like. It is not perfect, so I was wondering if it could be the reason for the poor Goodness of fit (only 5%)

@agramfort , I am sharing the files with you here. If possible, would you please take a look at what I could have possibly done wrong? The GOF for my fitted dipole is too low

import matplotlib
import pathlib
import mne
import mne_bids
matplotlib.use('Qt5Agg')
import numpy as np
cd  E:\
fname='phantom_rebuilt_MEG_4232022_L6.fif'
raw=mne.io.read_raw_fif(fname)
raw.load_data()
raw.filter(l_freq=1,h_freq=70)
picks = mne.pick_types(raw.info, meg= True, eog=False,
                       stim=False, exclude='bads')
raw.notch_filter(np.arange(60, 241, 60), picks=picks, filter_length='auto',
                 phase='zero')

event=mne.read_events('E:\L6M.eve')

for i in range(1286):
 event[i][1]=0

for i in range(1286):
    event[i][2]=1 # making all spikes as event 1

for i in range(1285):
    if event[i+1][0]- event[i][0]< 10:
        event[i+1][0]=event[i][0]
        event[i+1][2]=event[i][2] # if two spikes are really close to each other, making them   same event

 event_dict = {'1':1}

tmin, tmax = -.2, 0.3

epochs = mne.Epochs(raw, events=event, event_id=event_dict,tmin=-0.2, tmax=0.3,event_repeated='drop', baseline=None)

epochs.plot()

epochs.save('revised.fif',overwrite=True)

subjects_dir = 'E:\myMRI-20220510T004802Z-001\myMRI'

mne.viz.plot_bem(subject='E:\myMRI-20220510T004802Z-001\myMRI\_my_subject', subjects_dir=subjects_dir,
                 orientation='coronal')
trans='C.fif'
subject = 'E:\myMRI-20220510T004802Z-001\myMRI\_my_subject'
src = mne.setup_source_space(subject='E:\myMRI-20220510T004802Z-001\myMRI\_my_subject',
                             spacing='oct6',  
                             subjects_dir=subjects_dir,
                             add_dist=False)

conductivity = (0.3, 0.006, 0.3)  
model = mne.make_bem_model(subject='E:\myMRI-20220510T004802Z-001\myMRI\_my_subject', ico=4,
                           conductivity=conductivity,
                           subjects_dir=subjects_dir)

bem_sol = mne.make_bem_solution(model)
bem_sol

bem_fname = 'E:\myMRI-20220510T004802Z-001\myMRI\sampleBEM.fif'
mne.bem.write_bem_solution(bem_fname, bem_sol,overwrite=True)
cov=mne.compute_covariance(epochs,tmin=-0.2,tmax=-0.05)
evoked=epochs.average()
evoked.plot()
evoked_full = evoked.copy()
evoked.crop(0.00, 0.01)
dip = mne.fit_dipole(evoked, cov, bem=bem_sol, trans=trans)[0]
dip.plot_locations(trans=trans, subject=subject, subjects_dir=subjects_dir, mode='orthoview',coord_frame='mri',idx='gof')

Here is the google drive link to all the files

https://drive.google.com/drive/folders/1rZL6Zj-nZfzUlLNJKF2ThfTp0oY_93Et?usp=sharing

@Bamford - Not entirely sure of the discrepancy between the diff and original. I looked at your data and noticed a few items that might help. Your spikes (Iā€™m assuming epilepsy) are defined with a maximum at time t=0. Therefore your baseline shouldnā€™t be -0.2-0. If you change it to (-0.2,-0.1), you will get better results. Also at the end of your data you have a lot of noise - see image below. This needs to be removed from the data or it will dominate the signal. Also - I am not sure why, but the EEG doesnā€™t seem to be baseline correcting properly - I donā€™t know if this is a reference issue (??) - but the data near the end is floating off of baseline.

epo = mne.read_epoch('C.fif', preload=True)
epo.apply_baseline(baseline=(-0.2, -0.1))
epo.plot()

Produces this image - these are the last trials and the baseline does not appear to be set (??).

If I just take the first 100 epochs - you get a pretty strong activity vs baseline. The data also appears to be more dipolar (I work with MEG - but I think this is appropriate for an EEG dipole). Fitting this data will give you a better dipole fit / GOF.

epo2 = epo[0:100]
epo2.filter(1.0, None)
epo2.plot()
ave=epo2.average()

One other major component that is going to affect your data: your 3 layer BEM does not resolve the inner and outer-skull boundary well. In EEG this is pretty important since the skull is ~40X the resistance of the other layers. With a typical T1w MRI, it is not uncommon to have this issue.

And the final problem with this type of analysis is that the underlying neurophysiology may not be represented by a single dipole. With these large discharges, it is possible that there is a large synchronous network firing together. With that amount of tissue active - you might not be able to represent all of that activity with a single dipole (typically these results will send the dipole deeper). Its possible that distributed modelling (MNE, Loretaā€¦) would be preferred.

Those are some initial thoughts - but I canā€™t comment too much on the dipole simulation.

-Jeff

3 Likes

@Bamford - I just looked over what you had posted again. It appears that you are doing phantom MEG (ā€˜phantom_rebuilt_MEGā€™) and in your picks, you have meg=True. But in your data it is listed as EEG

Out[74]: 
<Info | 13 non-empty values
 bads: 21 items (E67, E68, E73, E74, E92, E83, E82, E91, E94, E95, ...)
 ch_names: E1, E2, E3, E4, E5, E6, E7, E8, E9, E10, E11, E12, E13, E14, ...
 chs: 256 EEG, 1 misc
 custom_ref_applied: True
 device_info: 1 item (dict)
 dig: 261 items (3 Cardinal, 257 EEG, 1 Extra)

Not sure if this is the expected data. And could result in the poor dipole fit if the sensor type is different from the actual data.

-Jeff

3 Likes

hI @jstout211 Thank you for taking your time out to go through the files. I seem to have attached a wrong epoch file ( the one with EEG here). I have updated the google drive with the new epoch file [ newmeg.fif].

I was able to get a GOF of 97% with the way suggested by you for EEG file, however the GOF for MEG files is still low. Could you give it a check?

Thanks

@Bamford - I looked at your MEG data. The sensor quality is a little noisy and inconsistent between channels. I marked about 10 as bad before making the epochs. I donā€™t know how much this really contributes to your issue - but just wanted to point it out.

I averaged all of the epochs and got the same topography map as you posted. It looks like your sensor data go into a steady state rythmn where the filter rolloff bleeds into the next trial (you can see the persistent activity between the two trials below). I thought that using the initial stimuli may get around this problem, since it wonā€™t be affected by prior trials. I set the baseline as 5 seconds prior to the initial event and everything looks relatively flat until the dipole stimuli. Unfortunately you still get the same topography.

The topography maps (at least the mag data) does not appear to show great phantom dipoles. It should appear to look like a butterfly pattern with a positive and negative lobe (similar to your dipole map), but the topography still looks a little noisy / non-dipolar for having such a high signal to noise level. I think your gradiometer data looks more appropriate, but there is some distant noise there as well.

Another comment on your data - you have an MRI that you are using to fit the dipole - but the MEG data does not have any neurophysiological activity. I would assume you are better off using a spherical head model that represents your phantom head. I dontā€™ really know how to do this in MNE.

I also am not sure what is causing the topography to be distorted, which would lead to a poor diple fit. There are many factors - possibly channel tuning (??). The magnetometers are going to be more sensitive to artifact compared with the gradiometers. Just as a quick test, I would select just the gradiometers (raw.pick_types(meg=ā€˜gradā€™)) and redo the forward model and the dipole fit using the MRI. See if this helps. It could be that the mags are just driven outside of their linearity range.

Did you perform SSS / tSSS on the data? Its possible that this might be able to help. Are you running MEG and EEG concurrently on this phantom? Or are the EEG leads still attached to the phantom, but not to the electronics - they can act like antennae and distort the fields.

Sorry I couldnā€™t help any more with this - I am kind of stumped as to why you are getting this kind of topography.

-Jeff

2 Likes