Skip to content

BUG: fNIRS interpolation mixes reordered and raw channel indices #13857

@Kallemakela

Description

@Kallemakela

Description of the problem

  1. def _interpolate_bads_nirs has an indexing problem. It constructs a
    distance matrix on faulty indices leading to mixing up channels during
    interpolation. Interpolation produces invalid results if the channels in raw
    are not in lexical order (which often happens since e.g. S10_D10 comes before
    S1_D1 in lexical order).
  2. _validate_nirs_info() reorders successfully validated fnirs channels. This
    leads to downstream issues in e.g. Raw.interpolate_bads() with nirs.
# `mne/channels/interpolation.py:266`
picks_nirs = _validate_nirs_info(inst.info) # This reorders index to [6, 7, 0, 1, 2, 3, 4, 5]
...
# `mne/channels/interpolation.py:275
chs = [inst.info["chs"][i] for i in picks_nirs]
locs3d = np.array([ch["loc"][:3] for ch in chs]) # now loc3d[0] corresponds to channel index 6
dist = pdist(locs3d)
dist = squareform(dist)
for bad in picks_bad:
    # `bad` here is in original raw index, dist is in reordered index
    dists_to_bad = dist[bad]
    ...
    # `mne/channels/interpolation.py:288-289`
    closest_idx = np.argmin(dists_to_bad) + (bad % 2)
    # faulty channel copied to replace bad
    inst._data[bad] = inst._data[closest_idx]

Fixes:

  1. Is trivial to fix and should not cause side-effects.
  2. Should probably not silently reorder valid channels.

PR: #13856

Steps to reproduce

ch_names = [
    "S1_D1 760",
    "S1_D1 850",
    "S2_D2 760",
    "S2_D2 850",
    "S3_D3 760",
    "S3_D3 850",
    "S10_D10 760",
    "S10_D10 850",
]
info = mne.create_info(
    ch_names=ch_names,
    sfreq=1.0,
    ch_types=["fnirs_cw_amplitude"] * len(ch_names),
)
pair_positions = {
    "S1_D1": (0.009, 0.0, 0.0),
    "S2_D2": (0.010, 0.0, 0.0),
    "S3_D3": (0.030, 0.0, 0.0),
    "S10_D10": (0.040, 0.0, 0.0),
}
for idx, ch in enumerate(info["chs"]):
    pair = ch["ch_name"].rsplit(" ", 1)[0]
    ch["loc"][:3] = pair_positions[pair]
    ch["loc"][9] = 760.0 if idx % 2 == 0 else 850.0

data = np.tile(np.arange(len(ch_names))[:, np.newaxis], (1, 5)).astype(float)
raw = mne.io.RawArray(data, info, verbose=False)
raw.info["bads"] = ["S2_D2 760", "S2_D2 850"]

picks_nirs = _validate_nirs_info(raw.info)
picks_bad = pick_channels(raw.ch_names, raw.info["bads"], exclude=[])
bads_mask = np.array([pick in picks_bad for pick in picks_nirs])
locs3d = np.array([raw.info["chs"][pick]["loc"][:3] for pick in picks_nirs])
dist = squareform(pdist(locs3d))

bad = picks_bad[0]
buggy_dists = dist[bad].copy()
buggy_dists[buggy_dists == 0] = np.inf
buggy_dists[bads_mask] = np.inf
buggy_donor = np.argmin(buggy_dists) + (bad % 2)
buggy_dist = buggy_dists[np.argmin(buggy_dists)]

bad_rel = np.where(picks_nirs == bad)[0][0]
correct_dists = dist[bad_rel].copy()
correct_dists[correct_dists == 0] = np.inf
correct_dists[bads_mask] = np.inf
correct_donor = picks_nirs[np.argmin(correct_dists) + (bad % 2)]
correct_dist = correct_dists[np.argmin(correct_dists)]

print("raw channel order:", raw.ch_names)
print("_validate_nirs_info channel order:", [raw.ch_names[pick] for pick in picks_nirs])
print(
    f"expected donor: {raw.ch_names[bad]} > {raw.ch_names[correct_donor]} "
    f"(dist: {correct_dist:g})"
)
print(
    f"buggy donor: {raw.ch_names[bad]} > {raw.ch_names[buggy_donor]} "
    f"(dist: {buggy_dist:g})"
)

raw.interpolate_bads(origin=(0.0, 0.0, 0.0), verbose=False)
correct_pair = raw.ch_names[correct_donor].rsplit(" ", 1)[0]
buggy_pair = raw.ch_names[buggy_donor].rsplit(" ", 1)[0]
bad_pair = raw.ch_names[bad].rsplit(" ", 1)[0]
label_width = max(len(correct_pair), len(buggy_pair), len(bad_pair))

print("interpolated value:")
print(f"{bad_pair:<{label_width}} values:", raw.get_data(picks=picks_bad)[:, 0].tolist())
print("donor values:")
print(f"{correct_pair:<{label_width}} values:", raw.get_data(picks=[correct_donor, correct_donor + 1])[:, 0].tolist())
print(f"{buggy_pair:<{label_width}} values:", raw.get_data(picks=[buggy_donor, buggy_donor + 1])[:, 0].tolist())

Link to data

None.

Expected results

Bad S2_D2 copies nearest good S1_D1.

Actual results

Bad S2_D2 copies S10_D10.

Output:

raw channel order: ['S1_D1 760', 'S1_D1 850', 'S2_D2 760', 'S2_D2 850', 'S3_D3 760', 'S3_D3 850',
 'S10_D10 760', 'S10_D10 850']
_validate_nirs_info channel order: ['S10_D10 760', 'S10_D10 850', 'S1_D1 760', 'S1_D1 850', 'S2_D
2 760', 'S2_D2 850', 'S3_D3 760', 'S3_D3 850']
expected donor: S2_D2 760 > S1_D1 760 (dist: 0.001)
buggy donor: S2_D2 760 > S10_D10 760 (dist: 0.021)
interpolated value:
S2_D2   values: [6.0, 7.0]
donor values:
S1_D1   values: [0.0, 1.0]
S10_D10 values: [6.0, 7.0]

Additional information

Platform             macOS-26.2-arm64-arm-64bit-Mach-O
Python               3.14.4 | packaged by conda-forge | (main, Apr  8 2026, 02:33:53) [Clang 20.
1.8 ]
Executable           /Users/kallemakela/mambaforge/envs/mnedev/bin/python
CPU                  Apple M1 Pro (8 cores)
Memory               32.0 GiB

Core
├☑ mne               1.13.0.dev22+gf85db9d20.d20260421 (development, latest release is 1.12.1)
├☑ numpy             2.4.3 (OpenBLAS 0.3.32 with 8 threads)
├☑ scipy             1.17.1
└☑ matplotlib        3.10.8 (backend=macosx)

Numerical (optional)
├☑ sklearn           1.8.0
├☑ numba             0.65.0
├☑ nibabel           5.4.0
├☑ nilearn           0.13.1
├☑ dipy              1.12.0
├☑ openmeeg          2.5.16
├☑ pandas            3.0.2
├☑ h5io              0.2.5
├☑ h5py              3.16.0
└☐ unavailable       cupy

Visualization (optional)
├☑ pyvista           0.47.3 (OpenGL 4.1 Metal - 90.5 via Apple M1 Pro)
├☑ pyvistaqt         0.11.4
├☑ vtk               9.6.1
├☑ qtpy              2.4.3 (PySide6=6.11.0)
├☑ ipympl            0.10.0
├☑ pyqtgraph         0.14.0
├☑ mne-qt-browser    0.7.4
├☑ ipywidgets        8.1.8
├☑ trame_client      3.11.4
├☑ trame_server      3.10.0
├☑ trame_vtk         2.11.7
└☑ trame_vuetify     3.2.1

Ecosystem (optional)
├☑ mne-bids          0.18.0
├☑ mne-connectivity  0.8.1
├☑ neo               0.14.4
├☑ eeglabio          0.1.3
├☑ edfio             0.4.13
├☑ curryreader       0.1.2
├☑ mffpy             0.10.0
├☑ pybv              0.7.6
├☑ pymef             1.4.8
├☑ antio             0.6.1
├☑ defusedxml        0.7.1
└☐ unavailable       mne-nirs, mne-features, mne-icalabel, mne-bids-pipeline

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions