Tweak y-axis label (#13298)

Co-authored-by: Daniel McCloy <dan@mccloy.info>

Author: cbrnr

doc/changes/devel/13298.bugfix.rst

@@ -0,0 +1 @@
+Fix y-axis label in PSD plot when ``dB=True`` to show a more conventional label (i.e., dB/Hz or dB/√Hz), by `Clemens Brunner`_.

mne/time_frequency/spectrum.py

@@ -606,8 +606,8 @@ def plot(
             Whether to plot an amplitude spectrum (``True``) or power spectrum
             (``False``).
 
-                .. versionchanged:: 1.8
-                    In version 1.8, the default changed to ``amplitude=False``.
+            .. versionchanged:: 1.8
+                In version 1.8, the default changed to ``amplitude=False``.
         %(xscale_plot_psd)s
         ci : float | 'sd' | 'range' | None
             Type of confidence band drawn around the mean when ``average=True``. If

mne/time_frequency/tfr.py

@@ -1960,7 +1960,7 @@ def plot(
 
             How baseline is computed is determined by the ``mode`` parameter.
         %(mode_tfr_plot)s
-        %(dB_spectrum_plot)s
+        %(dB_tfr_plot)s
         %(combine_tfr_plot)s
 
             .. versionchanged:: 1.3
@@ -2215,7 +2215,7 @@ def plot_joint(
 
             How baseline is computed is determined by the ``mode`` parameter.
         %(mode_tfr_plot)s
-        %(dB_tfr_plot_topo)s
+        %(dB_tfr_plot)s
         %(yscale_tfr_plot)s
         %(vlim_tfr_plot_joint)s
         %(cnorm)s
@@ -2508,7 +2508,7 @@ def plot_topo(
         %(layout_spectrum_plot_topo)s
         %(cmap_tfr_plot_topo)s
         %(title_none)s
-        %(dB_tfr_plot_topo)s
+        %(dB_tfr_plot)s
         %(colorbar)s
         %(layout_scale)s
         %(show)s

mne/utils/docs.py

@@ -1013,32 +1013,36 @@ def _reflow_param_docstring(docstring, has_first_line=True, width=75):
 
 _dB = """
 dB : bool
-    Whether to plot on a decibel-like scale. If ``True``, plots
-    10 × log₁₀({quantity}){caveat}.{extra}
+    Whether to plot on a decibel scale. If ``True``, plots
+    10 × log₁₀({quantity}){ampl}{caveat}.{extra}
 """
-_ignored_if_normalize = " Ignored if ``normalize=True``."
-_psd = "spectral power"
+_psd = "spectral_power/Hz"
 
+# for the legacy func/methods:
 docdict["dB_plot_psd"] = """\
 dB : bool
-    Plot Power Spectral Density (PSD), in units (amplitude**2/Hz (dB)) if
-    ``dB=True``, and ``estimate='power'`` or ``estimate='auto'``. Plot PSD
-    in units (amplitude**2/Hz) if ``dB=False`` and,
-    ``estimate='power'``. Plot Amplitude Spectral Density (ASD), in units
-    (amplitude/sqrt(Hz)), if ``dB=False`` and ``estimate='amplitude'`` or
-    ``estimate='auto'``. Plot ASD, in units (amplitude/sqrt(Hz) (dB)), if
-    ``dB=True`` and ``estimate='amplitude'``.
+    Plot power spectral density (PSD) in units (dB/Hz) if ``dB=True`` and
+    ``estimate='power'``. Plot PSD in units (amplitude**2/Hz) if
+    ``dB=False`` and ``estimate='power'``. Plot amplitude spectral density (ASD) in
+    units (amplitude/sqrt(Hz)) if ``dB=False`` and ``estimate='amplitude'``.
+    Plot ASD in units (dB/sqrt(Hz)) if ``dB=True`` and ``estimate='amplitude'``.
 """
 docdict["dB_plot_topomap"] = _dB.format(
     quantity=_psd,
-    caveat=" following the application of ``agg_fun``",
-    extra=_ignored_if_normalize,
+    ampl="",
+    caveat=", following the application of ``agg_fun``",
+    extra=" Ignored if ``normalize=True``.",
+)
+docdict["dB_spectrum_plot"] = _dB.format(
+    quantity=_psd,
+    ampl=", or 20 × log₁₀(spectral amplitude/√Hz) if ``amplitude=True``",
+    caveat="",
+    extra="",
 )
-docdict["dB_spectrum_plot"] = _dB.format(quantity=_psd, caveat="", extra="")
 docdict["dB_spectrum_plot_topo"] = _dB.format(
-    quantity=_psd, caveat="", extra=_ignored_if_normalize
+    quantity=_psd, ampl="", caveat="", extra=""
 )
-docdict["dB_tfr_plot_topo"] = _dB.format(quantity="data", caveat="", extra="")
+docdict["dB_tfr_plot"] = _dB.format(quantity="data", ampl="", caveat="", extra="")
 
 _data_template = """
 data : ndarray, shape ({})

mne/viz/tests/test_raw.py

@@ -987,6 +987,15 @@ def test_plot_raw_filtered(filtorder, raw, browser_backend):
     RawArray(np.zeros((1, 100)), create_info(1, 20.0, "stim")).plot(lowpass=5)
 
 
+def _check_ylabel_psd(ylabel, amplitude, dB, unit):
+    """Check that the ylabel is correct."""
+    numerator = "dB" if dB else unit
+    denominator = r"\sqrt{\mathrm{Hz}}" if amplitude else r"\mathrm{Hz}"
+    reference_val = rf"\ \mathrm{{re}}\ 1\ \mathrm{{{unit}}}" if dB else ""
+    pattern = rf"\mathrm{{{numerator}}}/{denominator}{reference_val}"
+    assert pattern in ylabel, ylabel
+
+
 def test_plot_raw_psd(raw, raw_orig):
     """Test plotting of raw psds."""
     raw_unchanged = raw.copy()
@@ -1038,24 +1047,20 @@ def test_plot_raw_psd(raw, raw_orig):
     with pytest.warns(UserWarning, match="[Infinite|Zero]"):
         spectrum = raw.compute_psd()
     for dB, amplitude in itertools.product((True, False), (True, False)):
-        with pytest.warns(UserWarning, match="[Infinite|Zero]"):
+        with pytest.warns(UserWarning, match="(Infinite|Zero)"):
             fig = spectrum.plot(average=True, dB=dB, amplitude=amplitude)
         # check grad axes
         title = fig.axes[0].get_title()
         ylabel = fig.axes[0].get_ylabel()
-        unit = r"fT/cm/\sqrt{Hz}" if amplitude else "(fT/cm)²/Hz"
+        unit = "fT/cm" if amplitude else "(fT/cm)^2"
         assert title == "Gradiometers", title
-        assert unit in ylabel, ylabel
-        if dB:
-            assert "dB" in ylabel
-        else:
-            assert "dB" not in ylabel
+        _check_ylabel_psd(ylabel, amplitude, dB, unit)
         # check mag axes
         title = fig.axes[1].get_title()
         ylabel = fig.axes[1].get_ylabel()
-        unit = r"fT/\sqrt{Hz}" if amplitude else "fT²/Hz"
+        unit = "fT" if amplitude else "fT^2"
         assert title == "Magnetometers", title
-        assert unit in ylabel, ylabel
+        _check_ylabel_psd(ylabel, amplitude, dB, unit)
 
     # test xscale value checking
     raw = raw_unchanged

mne/viz/utils.py

@@ -2431,25 +2431,22 @@ def _convert_psds(
         warn(msg, UserWarning)
 
     _check_option("estimate", estimate, ("power", "amplitude"))
+    psds *= scaling * scaling
+    denom = r"\sqrt{\mathrm{Hz}}" if estimate == "amplitude" else r"\mathrm{Hz}"
     if estimate == "amplitude":
         np.sqrt(psds, out=psds)
-        psds *= scaling
-        ylabel = rf"$\mathrm{{{unit}/\sqrt{{Hz}}}}$"
         coef = 20
     else:
-        psds *= scaling * scaling
         if "/" in unit:
             unit = f"({unit})"
-        ylabel = rf"$\mathrm{{{unit}²/Hz}}$"
+        unit = f"{unit}^2"
         coef = 10
+    ylabel = rf"$\mathrm{{{unit}}}/{denom}$"
     if dB:
         np.log10(np.maximum(psds, np.finfo(float).tiny), out=psds)
         psds *= coef
-        ylabel = r"$\mathrm{dB}\ $" + ylabel
-    ylabel = "Power (" + ylabel if estimate == "power" else "Amplitude (" + ylabel
-    ylabel += ")"
-
-    return ylabel
+        ylabel = rf"$\mathrm{{dB}}/{denom}\ \mathrm{{re}}\ 1\ \mathrm{{{unit}}}$"
+    return f"{'Power' if estimate == 'power' else 'Amplitude'} ({ylabel})"
 
 
 def _plot_psd(

tutorials/time-freq/10_spectrum_class.py

@@ -123,7 +123,7 @@
 # (interpolated scalp topography of power, in specific frequency bands). A few
 # plot options are demonstrated below; see the docstrings for full details.
 
-evk_spectrum.plot(picks="data", exclude="bads", amplitude=False)
+evk_spectrum.plot(picks="data", exclude="bads")
 evk_spectrum.plot_topo(color="k", fig_facecolor="w", axis_facecolor="w")
 
 # %%