Extracting the time series of activations in a label

Describe the new feature or enhancement

Current API: https://mne.tools/stable/generated/mne.SourceEstimate.html#mne.SourceEstimate.extract_label_time_course
implementation offers 4 modes to extract source time courses in a label.

Modes:
max, mean, mean_flip, pca_flip

Is it possible to extend the mode options to com and min modes?

• com: time course of the vertex that represents spatial center of mass in a label.
• min: minimum value across vertices at each time point within each label.

Reasons:

• For dSPM source time-courses with singed values, sometimes min mode (Minimum value across vertices) is more representative than max mode if the underlying label (dipoles) represent a surface negative source; for example: auditory cortex (Heschl's gyrus: N1 component). With a prior information (if the source/label is well known: surface positive/negative), min/max mode can be more useful than any other strategies to average the times series in a label.
• com: dipole that represent spatial center of mass in a label. com can be an extra option to extract source time-courses in a label i.e., shows similar values as min/max modes.

Example to follow:
https://mne.tools/stable/auto_examples/inverse/label_source_activations.html#sphx-glr-auto-examples-inverse-label-source-activations-py

Current, max mode in MNE is implemented as:

_label_funcs = {
    "mean": lambda flip, data: np.mean(data, axis=0),
    "mean_flip": lambda flip, data: np.mean(flip * data, axis=0),
     "max": lambda flip, data: np.max(np.abs(data), axis=0),
    "pca_flip": _pca_flip,
}

isn't it make more sense not to take the max(absolute(data)) and rather use only max (data) especially, if I consider from the signed prospective? I wrote some code to test this.

Describe your proposed implementation

test code:

import matplotlib.pyplot as plt
import matplotlib.patheffects as path_effects

import mne
from mne.datasets import sample
from mne.minimum_norm import read_inverse_operator, apply_inverse
import numpy as np

data_path = sample.data_path()
label = "Aud-lh"
meg_path = data_path / "MEG" / "sample"
subjects_dir = data_path / "subjects"
label_fname = meg_path / "labels" / f"{label}.label"
fname_inv = meg_path / "sample_audvis-meg-oct-6-meg-inv.fif"
fname_evoked = meg_path / "sample_audvis-ave.fif"

snr = 3.0
lambda2 = 1.0 / snr**2
method = "dSPM"  # use dSPM method (could also be MNE or sLORETA)

# Load data
evoked = mne.read_evokeds(fname_evoked, condition=0, baseline=(None, 0))
inverse_operator = read_inverse_operator(fname_inv)
src = inverse_operator["src"]

# %%
# Compute inverse solution
# ------------------------
pick_ori = "normal"  # Get signed values to see the effect of sign flip
stc = apply_inverse(evoked, inverse_operator, lambda2, method, pick_ori=pick_ori)

label = mne.read_label(label_fname)

stc_label = stc.in_label(label)
modes = ("com", "mne_max", "min", "max",  "mean", "mean_flip", "pca_flip")
tcs = dict()
for mode in modes:
    if mode == "com":
        # compute spatial center of mass (COM)
        label.values.fill(1.0)

        # Restrict the eligible vertices to be those on the surface under
        # consideration and within the label.
        hemi_to_ind = {"lh": 0, "rh": 1}
        surf_vertices = src[hemi_to_ind[label.hemi]]["vertno"]
        restrict_verts = np.intersect1d(surf_vertices, label.vertices)
        com = label.center_of_mass(subjects_dir=subjects_dir, subject='sample', restrict_vertices=restrict_verts)
        com_dipole = np.where(stc_label.vertices[0] == com)[0][0]
        tcs[mode] = [stc_label.data[com_dipole, :]]

    elif mode == "mne_max":  # current implementation in mne
        tcs[mode] = [np.max(abs(stc_label.data), axis=0)]

    elif mode == "min":
        tcs[mode] = [np.min(stc_label.data, axis=0)]

    elif mode == "max":
        tcs[mode] = [np.max(stc_label.data, axis=0)]

    else:
        tcs[mode] = stc.extract_label_time_course(label, src, mode=mode)
print("Number of vertices : %d" % len(stc_label.data))

# %%
# View source activations
# -----------------------

fig, ax = plt.subplots(1)
t = 1e3 * stc_label.times
ax.plot(t, stc_label.data.T, "k", linewidth=0.5, alpha=0.5)
pe = [
    path_effects.Stroke(linewidth=5, foreground="w", alpha=0.5),
    path_effects.Normal(),
]
for mode, tc in tcs.items():
    ax.plot(t, tc[0], linewidth=3, label=str(mode), path_effects=pe)
xlim = t[[0, -1]]
ylim = [-30, 30]
ax.legend(loc="upper right")
ax.set(
    xlabel="Time (ms)",
    ylabel="Source amplitude",
    title="Activations in Label %r" % (label.name),
    xlim=xlim,
    ylim=ylim,
)
mne.viz.tight_layout()
plt.show()

Test result:

Describe possible alternatives

None

Additional context

No response

Start with the SourceEstimate.extract_label_time_course API and the label_source_activations example linked in the issue, then inspect how the existing mean, mean_flip, pca_flip, and max modes are implemented. Confirm the intended signed max, min, and center-of-mass semantics with the issue discussion; done means the accepted modes work for the relevant label/source cases and are covered by tests and example documentation.