ENH: labels_near_sensors function

[larsoner]

Another thing I'm having to do multiple times in fNIRS analyses is say which anatomical region a channel corresponds to, basically a label-mapping analogous to stc_near_sensors.

The easiest solution is just to take each channel and do a cdist-like operation to the pial surface of fsaverage to say which region is closest. However, I like the idea better of taking a stc_near_sensors approach, where for the Evoked data we use np.eye(n_chan) to obtain a (n_dense_vertices, n_channels) array. Then the "easiest" solution just gets vertex numbers from np.argmax(..., axis=0) of this array when mode='nearest', but using other modes like sum or weighted is more like a np.argsort(..., axis=0)-type operation where you get a varying number of labels per channel that it maps into, along with the weights in those labels.

That way, the labels_near_sensors can:

1. Have an API that matches stc_near_sensors very closely, with the addition of a labels argument (which would typically be obtained from read_labels_from_annot for example).
2. Have an output that also mirrors whatever stc_near_sensors does to project activation onto the cortical surface.

This is another one that could in theory go in MNE-NIRS (cc @rob-luke), but to me MNE-Python is probably a better scope since stc_near_sensors was actually designed for ECoG originally, and it's nice to have labels_near_sensors mirror stc_near_sensors in terms of namespace. To avoid crowding the mne. namespace where stc_near_sensors lives, we could have this new function in mne.label.labels_near_sensors.

[rob-luke]

The correct way to determine what anatomical structure your fNIRS channels are sensitive too is to use a photon migration simulation. There are several tools available to do this and it's usually done when designing your experiment/montage.

The tool I tend to use is the FOLD toolbox. And there are convenience functions in MNE-nirs to read the fold data and to determine what structures standard positions are sensitive to and the sensitivity to each structure. https://mne.tools/mne-nirs/generated/mne_nirs.io.fold_channel_specificity.html#mne_nirs.io.fold_channel_specificity

Your suggestion would work with non standard (10-05) positions, a large benefit. But I don't think the approach you propose would be easily accepted in a fNIRS journal, as there are established methods for doing this.

But until we have our own photon migration integration that would allow for non standard 10-05 optode positions mne-tools/mne-nirs#405 your suggestion would be practically useful. But I think it's better suited to MNE-Python than MNE-NIRS.

[larsoner]

Nice, I'll try that instead! We can reopen this later if someone is interested. Just for posterity, this is code I wrote to approximately accomplish this task:

    ch_weights = list()
    labels = mne.read_labels_from_annot(
        'fsaverage', parc, subjects_dir=subjects_dir)
    offsets = dict(lh=0, rh=len(stc.vertices[0]))
    label_map = -np.ones(stc.data.shape[0], int)
    for li, label in enumerate(labels):
        assert label.hemi in ('lh', 'rh')  # need to handle "both" eventually
        hi = 0 if label.hemi == 'lh' else 1
        idx = offsets[label.hemi] + np.where(
            np.in1d(stc.vertices[hi], label.vertices))[0]
        assert (label_map[idx] == -1).all()
        label_map[idx] = li
    assert (label_map >= 0).all()  # complete (enough) parcellations
    # go through channels
    for ci, d in enumerate(stc_eye.data.T):
        ch_weights.append(list())
        assert d.shape == label_map.shape
        w = np.bincount(label_map, weights=d, minlength=len(labels))
        assert np.isclose(w.sum(), d.sum())  # all mapped somewhere
        w /= w.sum()  # normalize
        order = np.argsort(w)[::-1]
        for idx in order:
            perc = w[idx] * 100
            if perc == 0:
                break
            ch_weights[-1].append(f'{labels[idx].name}({perc:0.3f}%)')

    with open(op.join(results_path, f'ch_mapping_{parc}.csv'), 'w') as fid:
        fid.write('channel,labels\n')
        for ci in range(info['nchan']):
            fid.write(f'{ci + 1},')
            fid.write(':'.join(ch_weights[ci]) + '\n')

[rob-luke]

The fold data is precomputed for standard 10-20 etc positions only. So I think we will be revisiting this discussion soon 😉

[rob-luke]

Also https://mne.tools/mne-nirs/auto_examples/general/plot_12_group_glm.html#relating-responses-to-brain-landmarks