`output="complex"` for Welch method should only allow `average=None`

Abierto Apto para principiantes
  • ¿Algo no está bien?
#14,495 1 comentario 0 reacciones 0 asignados Ver en GitHub

Los mantenedores suelen responder en 1 día

@iis-MarkKuang ya está trabajando en esto.

Desde el 5/11/2026.

  • #14498 de @iis-MarkKuang — abierto

Evaluación

Dificultad
2/5
Tiempo estimado
1-3 horas
Aptitud para principiantes
74/100
Tipo de issue
Error
Claridad
Bien especificado
Estado de actividad
Activo
Stack tecnológico
python
Área
data

Línea de trabajo

Empieza siguiendo el manejo de Welch de output="complex" y average en epochs.compute_psd(), y luego comprueba cómo Spectrum.plot() convierte los coeficientes en potencia. Añade una comprobación que rechace los valores de average distintos de None para la salida compleja, y busca pruebas existentes que cubran los modos de salida de Welch y la representación gráfica. Se considera terminado cuando la combinación no válida genera un error claro y average=None y la salida de potencia conservan el comportamiento esperado.

Escrito por el modelo de indexación a partir del texto del issue.

Descripción

BUG
Description of the problem

For Welch method, if average is requested with output="complex", segments are aggregated for complex coefficients, which can mess up the phase information in the resulting coefficients. Subsequently computing power from these coefficients gives very different results to if output="power" was used directly. This also affects the conversion to power from coefficients in Spectrum.plot().

Steps to reproduce
import numpy as np

import mne
from mne.time_frequency.psd import _median_biases

data_path = mne.datasets.sample.data_path()
raw = mne.io.read_raw_fif(
    data_path / "MEG" / "sample" / "sample_audvis_filt-0-40_raw.fif", verbose=False
)
raw.pick("eeg", exclude="bads").crop(tmax=60).load_data(verbose=False)
epochs = mne.make_fixed_length_epochs(raw, duration=4.0, preload=True, verbose=False)
kwargs = dict(method="welch", n_fft=128, verbose=False)  # 4 segments per epoch


def coefs_to_power(spectrum):
    """Convert complex Welch coefficients to one-sided power."""
    coefs = spectrum.get_data()
    freqs, sfreq = spectrum.freqs, spectrum.info["sfreq"]
    # double all bins except DC and Nyquist (scipy.signal.spectrogram does this for
    # mode="psd", but not for mode="complex")
    scale = np.where((freqs == 0) | np.isclose(freqs, sfreq / 2), 1.0, 2.0)
    if coefs.ndim == 4:  # unaggregated: (epoch, channel, freq, segment)
        scale = scale[:, np.newaxis]
    return (coefs * coefs.conj()).real * scale


for average in ("mean", "median"):
    power = epochs.compute_psd(output="power", average=average, **kwargs).get_data()

    # output="complex" with average!=None aggregates the complex coefficients across
    # segments, so power computed afterwards is wrong
    spec = epochs.compute_psd(output="complex", average=average, **kwargs)
    power_agg_first = coefs_to_power(spec)

    # correct: average=None, compute power per segment, then aggregate
    spec = epochs.compute_psd(output="complex", average=None, **kwargs)
    power_seg = coefs_to_power(spec)
    if average == "mean":
        power_agg_last = power_seg.mean(axis=-1)
    else:
        n_segments = power_seg.shape[-1]
        power_agg_last = (
            np.median(power_seg, axis=-1) / _median_biases(n_segments)[n_segments]
        )

    print(f"average={average!r}")
    for name, psd in (
        ("complex -> aggregate -> power", power_agg_first),
        ("complex -> power -> aggregate", power_agg_last),
    ):
        ratio = psd / power
        print(
            f"  {name}: ratio to output='power' median={np.median(ratio):.3f} "
            f"min={ratio.min():.3g} max={ratio.max():.3g} "
            f"allclose={np.allclose(psd, power, rtol=1e-6, atol=0)}"
        )
Link to data

No response

Expected results

average!=None with output="complex" is not allowed, as it gives invalid coefficients.

Actual results

Computing output="complex" with average!=None is allowed, and gives very different results when trying to compute power from these coefficients, compared to the proper approach of unaggregated coefficients → unaggregated power → aggregated power (which is what happens with output="power").

average='mean'
  complex -> aggregate -> power: ratio to output='power' median=0.131 min=2.06e-08 max=0.992 allclose=False
  complex -> power -> aggregate: ratio to output='power' median=1.000 min=1 max=1 allclose=True
average='median'
  complex -> aggregate -> power: ratio to output='power' median=0.326 min=9.42e-08 max=36.2 allclose=False
  complex -> power -> aggregate: ratio to output='power' median=1.000 min=1 max=1 allclose=True
Additional information

AI disclosure: Claude Opus 5.5 initially identified the problem, which I further investigated.

Lenguaje dominante
Python
Estrellas
3.6k
Forks
1.6k
Merge medio
1 d 17 h
PR fusionados (30 d)
74

Preparar el entorno

Primeros pasos

  1. Lee el issue completo y luego la guía de contribución del proyecto.
  2. Comenta en el issue que vas a ocuparte — evita que dos personas hagan lo mismo.
  3. Haz un fork del repositorio y trabaja en una rama.
  4. Abre un pull request que haga referencia al número del issue.

Más de mne-tools/mne-python

Todos los issues de mne-tools/mne-python

Issues similares

Más issues de Python

Recibe los nuevos issues en tu correo

Un resumen breve de issues de GitHub para principiantes.