-
-
Notifications
You must be signed in to change notification settings - Fork 1.5k
[MRG] Improve interpolation of bridged electrodes #11094
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Conversation
| # look for group of bad channels | ||
| nodes = sorted(set(chain(*bridged_idx))) | ||
| G_dense = np.zeros((len(nodes), len(nodes))) | ||
| # fill the edges with a weight of 1 | ||
| for bridge in bridged_idx: | ||
| idx0 = nodes.index(bridge[0]) | ||
| idx1 = nodes.index(bridge[1]) | ||
| G_dense[idx0, idx1] = 1 | ||
| G_dense[idx1, idx0] = 1 | ||
| # look for connected components | ||
| _, labels = connected_components(csr_matrix(G_dense), directed=False) | ||
| groups_idx = [ | ||
| [nodes[j] for j in np.where(labels == k)[0]] for k in set(labels) | ||
| ] | ||
| groups_names = [ | ||
| [inst.info.ch_names[k] for k in group_idx] for group_idx in groups_idx | ||
| ] |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This snippet looks for group of electrodes. It's equivalent to the networkx code:
G = nx.Graph()
for bridge in bridged_idx:
G.add_edge(*bridge)
groups_idx = nx.connected_components(G)
I had no idea scipy had those functionalities, and I never used them. Maybe there is a nicer way to do this?
Also, I used connected_components(csr_matrix(G_dense), directed=False) as per the example of connected_components, but it also works without csr_matrix. Not sure if it's needed..
mne/preprocessing/interpolate.py
Outdated
| def _find_centroid(ch_pos, group_names): | ||
| """Compute the centroid position between N electrodes. | ||
| The centroid should be determined in spherical "head" coordinates which is | ||
| more accurante than cutting through the scalp by averaging in cartesian | ||
| coordinates. | ||
| A simple way is to average the location in cartesian coordinate, convert | ||
| to spehrical coordinate and replace the radius with the average radius of | ||
| the N points in spherical coordinates. | ||
| Parameters | ||
| ---------- | ||
| ch_pos : OrderedDict | ||
| The position of all channels in cartesian coordinates. | ||
| group_names : | ||
| The name of the N electrodes used to determine the centroid. | ||
| Returns | ||
| ------- | ||
| pos_centroid : array of shape (3,) | ||
| The position of the centroid in cartesian coordinates. | ||
| """ | ||
| cartesian_positions = np.zeros((len(group_names), 3)) | ||
| sphere_positions = np.zeros((len(group_names), 3)) | ||
| for i, ch_name in enumerate(group_names): | ||
| cartesian_positions[i, :] = ch_pos[ch_name] | ||
| sphere_positions[i, :] = _cart_to_sph(ch_pos[ch_name]) | ||
| cartesian_pos_centroid = np.average(cartesian_positions, axis=0) | ||
| sphere_pos_centroid = _cart_to_sph(cartesian_pos_centroid) | ||
| # average the radius and overwrite it | ||
| avg_radius = np.average(sphere_positions, axis=0)[0] | ||
| sphere_pos_centroid[0, 0] = avg_radius | ||
| # convert back to cartesian | ||
| pos_centroid = _sph_to_cart(sphere_pos_centroid)[0, :] | ||
| return pos_centroid |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I moved the computation of an average location in a separate function because I wanted to add a test specifically for this. Maybe this function and the associated test could be moved to a different module? Not sure where would be a good fit..
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think this is an okay place for it, we can move if it's useful elsewhere
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I don't mind, but we will definitely forget it 😄
mne/preprocessing/interpolate.py
Outdated
| G_dense = np.zeros((len(nodes), len(nodes))) | ||
| # fill the edges with a weight of 1 | ||
| for bridge in bridged_idx: | ||
| idx0 = nodes.index(bridge[0]) | ||
| idx1 = nodes.index(bridge[1]) | ||
| G_dense[idx0, idx1] = 1 | ||
| G_dense[idx1, idx0] = 1 | ||
| # look for connected components | ||
| _, labels = connected_components(csr_matrix(G_dense), directed=False) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I would do this with https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.coo_matrix.html
Or if you know the matrices will always be small (they probably will) then what you have is fine. I wouldn't bother casting explicitly to csr_matrix, though
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It will always be small, a worst (impossible) case scenario would be a 256 electrodes cap with all electrodes bridged. So 256x256.
So coo_matrix to avoid the for loop and remove the csr_matrix. I'll check it later today.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Up to you if you want to coo_matrix or not. It looks like you might have to loop to get idx0 and idx1. In principle you could maybe use np.searchsorted to avoid even those loops but it starts to be less and less readable...
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'll stick with the current version, I just changed the .index with a better np.searchsorted. I did not find an easy and readable way to convert the information to the format for coo_matrix. I was starting with:
nodes = sorted(set(chain(*bridged_idx)))
rows, cols = zip(*bridged_idx)
# need to convert the idx as with `nodes.index(idx)`
coo_matrix(np.ones(len(rows)), rows, cols, shape=(len(nodes), len(nodes)))
# ...
It doesn't look more readable, and the conversion of the idx requires the same loop.
Co-authored-by: Eric Larson <larson.eric.d@gmail.com>
|
Any way to kill mne-python/.github/workflows/linux_pip.yml Lines 18 to 20 in bdc435d
|
larsoner
left a comment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Just minor comments from me, otherwise LGTM!
|
@larsoner I addressed your last comments, good to merge on my end. |
|
@mscheltienne I disabled the GitHub "notify on commits pushed" so I didn't see that you had addressed my comments (too many emails) so thanks for the ping! I think it's fine to merge as-is, but @alexrockhill let me know if you want to look |
* upstream/main: (25 commits) DOC: Exclude some implicit refs in doc build (mne-tools#10433) MAINT: Better issue forms (mne-tools#11101) [FIX] Typo in example (mne-tools#11118) [MRG] Improve interpolation of bridged electrodes (mne-tools#11094) BUG: Spectrum deprecation cleanup [circle deploy] (mne-tools#11115) Add API entry list and map (mne-tools#10999) Add legacy decorator (mne-tools#11097) [ENH, MRG] Add time-frequency epoch source estimation (mne-tools#11095) Revert "Add error message when conversion of EEG locs to [circle deploy] (mne-tools#11104) MRG: Fixes for mne-tools#11090 (mne-tools#11108) add test for edf units param (mne-tools#11105) BUG: Improve logic for bti (mne-tools#11102) add spectrum class (mne-tools#10184) ENH : add units parameter to read_raw_edf in case units is missing from the file (mne-tools#11099) ENH: Add temperature and galvanic (mne-tools#11090) Add error message when conversion of EEG locs to head space fails (mne-tools#11080) DOC: removed unnecessary line in PSF example (mne-tools#11085) FIX: Update helmet during ICP (mne-tools#11084) Fix various typos (mne-tools#11086) DOC: Rel ...
Closes #11079
Closes #11088
On main,
interpolate_bridged_electrodescreates one virtual channel for each pair of bridged electrodes. On this PR, it creates one virtual channel for each group of bridged electrodes. A warning is emitted if the group includes more than 4 electrodes (not included).On main, the virtual channel location is determined by averaging the 3 spherical coordinates of each pair of electrodes. On this PR, the location is determined by averaging the 3 cartesian coordinates of each electrode in the group. The resulting position is converted to spherical coordinates and the radius is replaced with the average radius of the electrodes in the group before being converted back to cartesian coordinates.
I ran this example locally, just in case: https://mne.tools/stable/auto_examples/preprocessing/eeg_bridging.html#ex-eeg-bridging