A colleague (@proloyd) and I were recently discussing the issue of EEG referencing in source localization algorithms, which prompted us to look into some old technical notes on the minimum-norm estimate inverse operator. We have a question that may be relevant for MNE-Python's construction of inverse operator, and we hope to raise this question to the developers and seek some clarifications.

In https://mne.tools/stable/overview/implementation.html#the-linear-inverse-operator, the inverse operator M is expressed using the lead field gain matrix G, which can be computed from Maxwell equations. However, we presume the matrix G here is full ranked and a theoretical gain matrix relative to an infinity reference. In real EEG recordings, we can only obtain referenced data, which are typically relative to a reference electrode either on the scalp or elsewhere on the body. The issue is on how to deal with this discrepancy in source localization of real data.

Digging into the MNE-Python inverse operator construction code a bit, we see a function here:

that seems to check whether average referencing is needed. It looks like if source estimates are desired, scalp EEG potentials are forced to be average referenced. Is there a justification to this requirement? Is the intention that once EEG data has been average referenced, it is "as close as it gets" to the theoretical scalp potential that would have been obtained directly from the gain matrix G and therefore the issue of referencing hereafter omitted?

An old "newsletter" between Matti and the author of sLORETA actually proved that for minimum-norm estimates, the unipolar reference applied on the EEG scalp data is not relevant, as long as the lead field gain matrix is average referenced (See equation 7) http://www.uzh.ch/keyinst/NewLORETA/BriefHistory/LORETA-NewsLett2b.pdf

We know that the theoretical gain matrix is whitened and scaled in MNE-Python inverse operator here

But is the gain matrix being average referenced at any point down the road? Perhaps when source estimates are being computed?

This could potentially be a bug if the gain matrix is used as is and never average referenced, because average referenced scalp EEG data is different from the theoretical scalp potential that would have been obtained with infinity reference. If it is indeed considered and taken into account, perhaps a clarification should be mentioned in the implementation details page so that people are not mislead by the expression of linear inverse operator currently in https://mne.tools/stable/overview/implementation.html#the-linear-inverse-operator.

Thank you for your attention!