Don't trust the source power spectrum results

[rickwassing]

Hi, I have calculated the Power Spectrum in Frequency Bands in Source Space, according to the details below. But I don't trust the results, because the alpha power (8-12 Hz) maximum values are on the bottom of the cortex, but on top of the cortex they are very close to zero (see figure 1). Nonetheless, the power spectral maps in sensor space are very reasonable and look like there should be sources in a frontal/central region and occipital region (see figure 2). I don't know where to start to see where things may have gone wrong.

Thank you very much for your help!

Rick Wassing
University of Sydney

Figure 1

Figure 2

Figure 3

Pipeline

• Imported subject anatomy from FreeSurfer output (15000 vertices), checked the orientation of MRI was correct and marked fiducial points.
• Imported EEGLAB dataset of cleaned sleep-EEG data from stage NREM2 only (~1hr of recording) i.e. bad channels and bad epochs removed.
• Recorded channel locations (Geoscan device) loaded from EEGLAB structure were registered to the head mask using the MRI registration tool.
• Calculated OpenMEEG BEM Forward head model.
• No noise modeling (identity matrix) because this is a "resting-state" recording and no noise estimation can be performed.
• Computed sources using minimum norm imaging, sLORETA (unconstrained)
• I then computed the power spectrum (Welch) on both the sensor-space and source-space data.

[Francois]

Your processing pipeline looks good to me.
To make the import and processing faster and much lighter in terms of memory, you could avoid importing the entire file as a single epoch in the database, by using the menu "Review raw file" instead of "Import MEG/EEG".
https://neuroimage.usc.edu/brainstorm/Tutorials/ChannelFile#Review_vs_Import

Your problem could be due to several things:

• a reference issue (Source Localization Abnormality - EEG reference problem)
• a problem of electrodes positions (have you checked that the electrodes are nicely located on the head? - please post a screen capture that shows the electrodes at the surface of the head)
• a problem in the OpenMEEG computation: try using a spherical model instead
• a failure specific to sLORETA: try using dSPM instead, or a non-normalized map ("Current density map")
• a problem due to the absence of noise covariance information: compute a noise covariance from a short segment of recording where no major event occur

I don't know if any preprocessing step performed by EEGLAB could cause additional trouble...
Have you tried reading and processing the continuous files in Brainstorm?

[rickwassing]

Thank you very much for your response. I have now used the "Review raw file" option, thanks for pointing that out.

Success
The problem resolved by using a Spherical Head model instead of OpenMEEG BEM. I originally used the Freesurfer output from the individual's T1 MRI to construct the BEM head model. To check whether the BEM head model was faulty because of some problem with the Freesurfer surfaces, I ...

• deleted the individual anatomy files
• selected to use Default Anatomy
• re-generated BEM surfaces
• re-registered the sensors to this default anatomy
• re-calculated the OpenMEEG BEM head model
• computed sources using Current Density and sLORETA maps

As can be seen in Figure 1B and C, using the individuals head-model results in unlikely source estimates. However, as can be seen in 1D and 1E, when using the default template head model, both sLORETA and the Current Density map look plausible.

So, it seems that the BEM head model is not right due to a problem with the Freesurfer surfaces. I then recalculated the surfaces using CAT12 and ran through the steps above. However, the results were very similar to using the Freesurfer surfaces, i.e. not good.

To troubleshoot, I tried to project the sources from the individual BEM head model to default anatomy, and then compared the results to those obtained by directly using a default anatomy BEM head model. The results don't overlap at all. If you have another suggestion to troubleshoot this problem, please let me know. Or do you think that my surfaces are not the problem, but something else?

Suggestions that did not matter

• Reference Issue. The EEG data was already referenced to the average. Re-referencing again to the average made no difference.
• The electrode positions are nicely located on the head surface (see Figure 2).
• I could not compute noise covariance, when right-clicking on the channel file > Noise covariance, the options are: no noise modeling, or import from file or from matlab.

Figure 1. Alpha power in (A) sensor-space; individual cortical surface using sLORETA (B) and current density (C), and default template cortical surface using sLORETA (D) and current density (E).

Figure 2. Electrode Positions are ok

[Francois]

re-generated BEM surfaces

Have you checked visually the surfaces?
Open the 3 BEM layers and the low-resolution cortex in the same figure, set a different color for the inner skull for improved readability post the screen capture here.

re-calculated the OpenMEEG BEM head model

• Have you seen any warning in a message box or in the Matlab command window that indicates that some dipoles are not positioned correctly?
• For debugging further the OpenMEEG solution: right-click on the forward model > View lead field vectors, navigate between sensors and see if there are any dipoles that do not follow the trend of their neighbors.
• I suggest you try comparing the OpenMEEG BEM results with a spherical model: compare the source maps, and compare the leadfields (select both the BEM and spherical models, right-click > View lead field vectors)

So, it seems that the BEM head model is not right due to a problem with the Freesurfer surfaces. I then recalculated the surfaces using CAT12 and ran through the steps above. However, the results were very similar to using the Freesurfer surfaces, i.e. not good.

The quality of the BEM solution depends mostly on the BEM surfaces, which are generated by Brainstorm. Whether you use CAT12 or FreeSurfer to generate the cortex surface influences very little these BEM layers.
It is possible that one cortex surface (FreeSurfer or CAT) has more vertices outside the inner skull (computed by Brainstorm) than the other, which could cause instabilities. But this is something you can fix with a right-click on the cortex > Force inside skull.

computed sources using Current Density and sLORETA maps

The minimum norm maps without normalization (current density maps) look similar and OK for both the individual head model and the template head model. The problem is probably not in the surfaces, but rather in the sLORETA normalization (widely used in the EEG literature, but not much promoted by Brainstorm and MNE-Python developers). I suggested you try with dSPM instead, as done in this tutorial:
https://neuroimage.usc.edu/brainstorm/Tutorials/RestingOmega#Source_estimation

Are these recordings only one task (eg. resting or sleeping)? Otherwise, if you have two states (eg. a task vs. rest), compute the noise covariance from the resting periods, or normalize the source maps with a Z-score wrt the rest recordings.

Reference Issue. The EEG data was already referenced to the average. Re-referencing again to the average made no difference.

In the forum thread I pointed out, the problem was related with the type of some data channels (adding or removing the mastoid electrodes to the group of EEG channels), but this is most likely not a problem here, because the recordings were converted to avg ref before being imported in Brainstorm.

I could not compute noise covariance, when right-clicking on the channel file > Noise covariance, the options are: no noise modeling, or import from file or from matlab.

If you want to compute the noise covariance matrix from recordings, you need to right-click on the RECORDINGS, not the channel file.

[Francois]

@Sylvain @John_Mosher
Would you have additional recommendations?

[rickwassing]

First of all, thank you very much for your quick response, that's incredible.

BEM surfaces. Visual inspection does not reveal any problems. I forced to fill any possible holes in the head-surface (right-click > Fill Holes, factor 2; although I understand that BS already does that by default), and I forced all cortical vertices inside skull, but a message appeared this was already the case.

OpenMEEG BEM head model. There was no warning or message to indicate any mis-positioning of dipoles. When inspecting the dipoles, the vector lengths of the majority of vertices were quite long for sensors located on the lower perimeter, i.e. around the side and back of the neck (example channel 82). Would you not expect all vectors to be short(er), just like the vectors far away from channel 9 in the example on the right?

Compare OpenMEEG to Spherical head model. The vectors of the spherical head model are in general quite a bit smaller as compared to those of the OpenMEEG model. In addition, the vector lengths for channel 82 seem to be more uniformly small across all vertices, and for channel 9 only those that are close to the channel are larger.

Compare Alpha Power Results. Ok, based on the results, it seems that the two head models are actually very much equivalent. The Current Density Maps show alpha activations in the Occipital pole (as expected), but do not show any clear activation in central/frontal regions (as shown in the topography). In addition, there seems to be a clear activation in the temporal poles with the OpenMEEG head model, which is not apparent in the topography or with the Spherical Head model. Secondly, the normalised maps (dSPM) show the main activation in the medial cortices, irrespective of the head model.

.

Beta Power Results.* Next, I checked whether the dSPM-normalised Alpha activation in the medial cortices were driven by the central/frontal Alpha activation shown in the topography. The topography of Beta power does not show such a central/frontal activation (except for one left frontal channel), however, the dSPM maps still show a main Beta power activation in the medial cortices, which seems to be not right.

We only have one "task", i.e. sleeping. We could estimate the noise level based on the first few minutes of the recording where the participant has their eyes still open. But we're studying sleep disorders where we're interested in wake-like brain activity during sleep, so I guess we'd best use an identity matrix as noise estimates? And, I see, duh, that I have to right-click on the recording to estimate noise from a recording (face-palm).

All EEG channels are of type "EEG", including the original reference channel (Cz), but indeed, the recording is converted to av-ref.

[Francois]

Would you not expect all vectors to be short(er), just like the vectors far away from channel 9 in the example on the right?

The size of the blue arrows is normalized, you can't make any interpretation about their length. What is meaningful is their relative length between different parts of the brain for one given sensor. All the arrows of similar length and orientation indicates that the all dipoles are influenced in a similar way by the selected electrode, which is expected for a distant electrode.
@tmedani Is this correct?

In addition, there seems to be a clear activation in the temporal poles with the OpenMEEG head model, which is not apparent in the topography or with the Spherical Head model.

This is mostly a matter of colormap setting. In all cases, these maps should be normalized, you shouldn't try to interpret these maps directly. After computing group-level statistics with many subjects, the results obtained with the two forward models might not be very important.

Secondly, the normalised maps (dSPM) show the main activation in the medial cortices, irrespective of the head model.

What clearly fails is the normalization of your minimum norm maps (both dSPM and sLORETA).
@Sylvain @John_Mosher: Could you please share your opinion on this?

But we're studying sleep disorders where we're interested in wake-like brain activity during sleep, so I guess we'd best use an identity matrix as noise estimates?

That's what we wrote in the noise covariance tutorial, indeed (https://neuroimage.usc.edu/brainstorm/Tutorials/NoiseCovariance#Variations_on_how_to_estimate_sample_noise_covariance), but it doesn't seem to work in your case.

Can you please post a screen capture illustrating what you obtain with a noise covariance estimated on a short segment of "calm" segments (early stages of sleep maybe)? Just to make sure that this is the problem for dSPM/sLORETA.

[Sylvain]

Please make sure this is not an issue caused by an EEG reference issue. I recommend you change the reference to "Average". I would also be curious to look at a portion of the time series on which you have computed the PSD. The scalp topography of alpha power is a bit off wrt to the expected typical parieto-occipital scalp patterns.

[tmedani]

Hello,

@Francois

Would you not expect all vectors to be short(er), just like the vectors far away from channel 9 in the example on the right?

The size of the blue arrows is normalized, you can't make any interpretation about their length. What is meaningful is their relative length between different parts of the brain for one given sensor. All the arrows of similar length and orientation indicates that the all dipoles are influenced in a similar way by the selected electrode, which is expected for a distant electrode.
@tmedani Is this correct?

@Francois your explanation is correct, this can be interpreted as the sensitivity of the selected lead/sensors to a source locating within the source space, this sensitivity varies by the localization of the source.

regarding the length of the vectors, I think when the methods are displayed separately the scale is normalized according to the current data. And when displayed together, they should have similar lengths and orientations. A difference may appear according to the used method and the parameters.

[rickwassing]

Hi everyone, again thanks for your responses, and my apologies for my late response. In the last 6 weeks I've used multiple resting-state recordings from different studies to make sure the issue is not related to one particular recording. But to no avail. The issue with the normalisation step is still there.

In order to facilitate ease of interpretation, I have now stepped away from sleep data, and loaded a resting-state recording (~30 seconds eyes open and 60 seconds eyes closed) directly from the raw recording (.mff, Philips-EGI Netstation Recording).

• Corrected for ECG and Blink artefacts using SSP as per your tutorials
• Manually entered bad channels
• Calculated the Average reference SSP (As suggested by Sylvain above)
• Bandpass filter between 0.5 and 48 Hz including applying the SSP-projections

I zipped the protocol directory from the BST database, and uploaded it here
https://cloudstor.aarnet.edu.au/plus/s/nULGnASXSDTQGzL

Here a screenshot of a subset of channels showing ~30 seconds of eyes-open data, which were used to calculate the noise covariance matrix.

The noise covariance matrix looks like this (colormap 0 - 0.1)

This is 30 seconds (out of 60 seconds) of eyes closed resting state data with clear alpha activity.

Next, I calculated the sources with the MNE method and 2 measures: Current Density and dSPM (otherwise same settings: unconstrained, use depth weighting, diagonal covariance, etc.)

And then I computed the power spectrum density (Welch method). My database looks like this:

And finally, for the alpha-band,

• on top row, the topoplot, in 2D and 3D, as well as the sensor locations
• in the middle row, the Current Density Maps (looks good)
• in the bottom row, the dSPM map (does not look good)
  
  

I have tried all other Noise Covariance Regularization options (Regularize, Median Eigenvalue, No Covariance, Shrinkage) but the results are similar.

Your help is again very much appreciated.

Best,

Rick Wassing

[Francois]

• Calculated the Average reference SSP (As suggested by Sylvain above)
• Bandpass filter between 0.5 and 48 Hz including applying the SSP-projections

I recommend you follow the pipeline in this tutorial:
https://neuroimage.usc.edu/brainstorm/Tutorials/Epilepsy#Pre-process_recordings
Filter first and then apply any linear operator (re-referencing, SSP, ICA...)

The noise covariance matrix looks like this (colormap 0 - 0.1)

Some channels look particularly noisy, you should probably remove then before any other preprocessing step.

(otherwise same settings: unconstrained, use depth weighting, diagonal covariance, etc.)

Our recommendation was to use constrained orientations:
https://neuroimage.usc.edu/brainstorm/Tutorials/RestingOmega#Source_estimation
https://neuroimage.usc.edu/brainstorm/Tutorials/SourceEstimation#Source_estimation_options

• in the middle row, the Current Density Maps (looks good)
• in the bottom row, the dSPM map (does not look good)

Then do not use the dSPM normalization.
Why do you need it in the first place?
If you use normalized power maps, it should not change (much) the results if you use one or the other.

I zipped the protocol directory from the BST database, and uploaded it here

The noise recordings you shared are not what you posted in your message, it has a lot of noise in the beginning and at the end:

In the noise covariance, all the diagonal points that are clearly higher than the other should be marked as bad channels. Or use only 4-18s in these recordings.

It looks like you recordings have a lot of eye movements, you could try getting rid of these with ICA.

Doing all these adjustments and computing the relative power, as done in the resting state tutorial, the maps for the alpha band look great.

[Francois]

The min norm / dSPM localize the alpha generators very broadly, possibly they are a bit deep and and extended in space.

Once divided with the total power, it shows exactly what we expect:

[rickwassing]

Fantastic, I could replicate your screenshots. Thank you for your answers and persistence in helping out! I see now that what I assumed was an issue, is not actually an issue at all. Below my thoughts to make sure I understand it now.

It was just related to the normalisation step in dSPM where the superficial vertices have relatively larger noise covariance estimates as compared to the deeper vertices (this helps to reduce the bias towards superficial sources). Therefore, when the absolute power-spectral source estimates were normalized (divided) by these noise covariance estimates, the deeper source estimates have higher values than the superficial ones. Please correct me if this is not correct.

May I please still ask one more question? Because I can't wrap my head around it. Below the comparison of the topography (top row), non-normalised Current Density (middle row) and dSPM (bottom row), between absolute alpha power (first column) and relative alpha power (second column). Why is the dSPM map of relative alpha power exactly equal to the non-normalised current density map of relative alpha power? I calculated the noise covariances over a different segment of data than the segment over which I calculated the Welch power spectrum density. So, it can't be that the noise covariance estimate and total power is equal, right?

[Francois]

Why is the dSPM map of relative alpha power exactly equal to the non-normalised current density map of relative alpha power?

If you apply two successive normalization procedures processing each vertex of the source maps individually, only the second one matters.
If you compute a Z-score (https://neuroimage.usc.edu/brainstorm/Tutorials/SourceEstimation#Standardization_of_source_maps) of the non-normalized current density maps or the dSPM, you get exactly the same result (minus some rounding errors). Same goes for this normalization of the power maps wrt to the total power.