feature.mfcc: expose `power_to_db` parameters, as proposed in #1734

Summary

This is a follow-up to #1734 rather than a new report. In that thread @bmcfee diagnosed the cause correctly and proposed the fix:

This is basically another iteration on an issue identified in the multichannel PR, wherein the noise floor calculation became sensitive to information leaking across different channels.

Probably what we ought to do here is expose some of the parameters to power_to_db, ie, so that you can explicitly set this in the call to mfcc rather than having to split it up into multiple steps.

The issue was closed three days later once the reporter's own problem was solved by the two-step workaround. The proposed change does not appear to have been made: in 1.0.0, feature.mfcc still takes no top_db, and the call is

S = power_to_db(melspectrogram(y=y, sr=sr, norm=mel_norm, **kwargs))

power_to_db receives no arguments, so top_db stays at its default of 80.0, and **kwargs is documented as going to melspectrogram. In 1.0.0 power_to_db's axes='auto' resolves to (-2, -1) for 2-D input, so the reduction on line

log_spec = np.maximum(log_spec, log_spec.max(axis=axes, keepdims=True) - top_db)

is over both mel and time — i.e. the whole clip for mono input.

I am opening this because I ran into it while building conformance vectors for MFCC, and the severity turned out to be larger than #1734 shows. Not a new diagnosis, just a measurement that might help prioritise the change that was already agreed.

What the clamp costs

Self-contained, no downloads, librosa 1.0.0:

import numpy as np, librosa

SR, HOP, WIN, NFFT = 16000, 160, 400, 512
rng = np.random.default_rng(0)

quiet = 0.001 * rng.standard_normal(HOP * 200 + WIN)   # 200 frames, ~2 s
kw = dict(sr=SR, n_mfcc=13, n_fft=NFFT, win_length=WIN, hop_length=HOP,
          center=False, n_mels=26)

alone = librosa.feature.mfcc(y=quiet, **kw)

print(f"{'appended burst':>16}{'max|diff|':>12}{'energy in c1..c12':>20}{'frames flat':>13}")
for amp in [0.0, 0.05, 0.9, 5.0, 50.0]:
    y = quiet if amp == 0 else np.concatenate([quiet, amp * rng.standard_normal(HOP * 100)])
    got = librosa.feature.mfcc(y=y, **kw)[:, : alone.shape[1]]
    flat = int((np.abs(got[1:]).max(axis=0) < 1e-9).sum())
    label = "none" if amp == 0 else f"+{20*np.log10(amp/0.001):.0f} dB"
    print(f"{label:>16}{np.abs(alone-got).max():>12.4f}"
          f"{np.abs(got[1:]).sum():>20.2f}{flat:>9d}/{alone.shape[1]}")

center=False and a length that is an exact multiple of hop_length mean the 200 compared frames are byte-identical in both computations; only what follows them differs.

  appended burst   max|diff|   energy in c1..c12  frames flat
            none      0.0000             5199.94        0/200
          +34 dB      0.0000             5199.94        0/200
          +59 dB      0.0000             5199.94        0/200
          +74 dB     15.1276             1720.86        0/200
          +94 dB    116.8845                0.00      200/200

The last row is the part I did not expect. Once the burst is loud enough that the quiet frames sit more than top_db below the file peak, every mel bin in those frames clamps to the same floor, the frame becomes constant across mel, and the DCT of a constant is zero everywhere except c0. All 200 frames go flat and c1..c12 carry exactly nothing. A door slam in an otherwise quiet room recording is roughly that far above the noise floor.

The threshold is sharp rather than gradual: nothing at +59 dB, partial at +74 dB, total at +94 dB.

Confirming the clamp is the sole cause, using the two-step route from the mfcc docstring:

def unclamped(y):
    M = librosa.feature.melspectrogram(y=y, sr=SR, n_fft=NFFT, win_length=WIN,
                                       hop_length=HOP, center=False, n_mels=26)
    return librosa.feature.mfcc(S=librosa.power_to_db(M, top_db=None), n_mfcc=13)

top_db=None via the two-step route: max|diff| = 0.000000000

Same frames, same everything, clamp removed, bitwise identical.

Why it might be worth doing now

• Streaming and batch inference compute different features from identical audio by default, so a model trained on whole files and served on chunks sees a distribution shift with no error and no warning.
• A single loud transient degrades features everywhere else in the same file, which is common in field recordings, meeting audio, and anything with a door, a cough, or clipping.
• The workaround exists and is in the docstring, but it is only reachable by knowing to avoid feature.mfcc, and the failure is silent for anyone who does not.

I appreciate this is a default-behaviour change and the compatibility question is real. Adding top_db (and possibly ref/amin) as explicit feature.mfcc arguments keeping the present defaults would make it opt-out without changing anything for existing callers, which I think matches what was proposed in #1734.

Happy to open a PR for that if it would be useful — including melspectrogram and any other feature that routes through power_to_db, if you would prefer it done consistently.

Versions: librosa 1.0.0, numpy 2.5.2, Python 3.12, macOS.