Repository: librosa/librosa
Pull request #1882 (2024-10-25 13:47:35)

Improved CMND implementation for YIN / pYIN
**TL;DR**: We can noticeably improve the accuracy of the YIN implementation by getting rid of the `win_length` parameter entirely and simply calculating the autocorrelation function (ACF) over the entire frame. Mathematical reasoning and some benchmarks provided below.

Apologies in advance for the wall of text. The devil is very much in the details with this one and I wanted to document it thoroughly.

## ACF Calculation

The main point of this pull request is to replace the calculation of the ACF: https://github.com/librosa/librosa/blob/d5aa7e1ae443c917fa53bef31055ed2613d63133/librosa/core/pitch.py#L393-L395 with a call to librosa's existing implementation: 

```python
acf_frames = audio.autocorrelate(y_frames, max_size=max_period+1, axis=-2)
``` 

Looking at this [discussion](https://github.com/librosa/librosa/pull/1063#discussion_r435547835), I believe the intention behind the current implementation is to calculate a type I ACF with a constant window size. Given window length $N$, this is defined as

$$r(k) = \sum_{n=0}^{N-1}x[n]x[n + k]$$

By contrast, the type II ACF uses an integration window that shrinks as the lag increases:

$$r'(k) = \sum_{n=0}^{N-1-k}x[n]x[n + k]$$

I'll go into some detail as to why I believe that 

1. the current implementation does not actually correctly calculate the type I ACF
2. the librosa implementation should use the type II ACF anyway (as implemented in `librosa.core.autocorrelate`)

It's true that the authors of the original YIN paper recommend to use the type I ACF, because the type II ACF implicitly tapers to 0 as $k \rightarrow N$. However, in the formulation of the type I ACF, the $n + k$ index will exceed the window size. This is not an issue discussed in the paper, since the authors consider the signal in its entirety, rather than on a frame-by-frame basis, so they can essentially "look into the future". 

As I understand it, the current implementation addresses this by making a distinction between `frame_length` $N$ and `win_length` $W < N$, where `win_length` essentially places an upper bound on the calculated lags. Mathematically,

$$r(k) = \sum_{n=0}^{W-1}x[n]x[n+k] \text{ for } 0\leq k < N - W$$

However, due to the `win_length:0:-1` indexing, the current code instead calculates:

$$r(k) = \sum_{n=0}^{W-1}x[n]x[n+k+1] \text{ for } 0\leq k < N - W$$

This off-by-one error was introduced [here](https://github.com/librosa/librosa/issues/1425). I think the reason that this change actually improved accuracy at the time is that it might have cancelled out another off-by-one error in the following calculation of the cumulative mean normalized difference (CMND) function, which this PR also fixes (see next section).

But this brings me to the second point of the argument: rather than fixing this calculation, I propose to drop the `win_length` parameter entirely and simply use `librosa.core.audio.autocorrelate` to calculate the type II ACF over the entire frame. This has several benefits:

- it simplifies the function call interface by removing a parameter
- using the type II ACF is standard in frame-by-frame implementations of YIN, in e.g. [aubio](https://github.com/aubio/aubio/blob/master/src/pitch/pitchyin.c) and also in the [original pYIN implementation](https://code.soundsoftware.ac.uk/projects/pyin/repository/entry/YinUtil.cpp#L57), so it aligns librosa with expected behavior.

As an example, the different ACFs are illustrated below for a 2048-sample frame from track 2 of the vocadito dataset. The highest peak of the type II ACF (excluding 0 lag) is at the correct period, while the type I ACF's peak at the correct period is slightly lower than the peak at double the period.

![ACF comparison](https://github.com/user-attachments/assets/eab25b47-8c1b-4927-b188-d31e1ef8054f)


## CMND Calculation

The next step in the YIN algorithm is to calculate the difference function. As mentioned before, the current implementation makes another off-by-one error here: https://github.com/librosa/librosa/blob/d5aa7e1ae443c917fa53bef31055ed2613d63133/librosa/core/pitch.py#L398-L406

Translating line 401 into mathematical notation, we see that the `energy_frame` $E$ is calculated as

$$\begin{align}
E(k) &= \sum_{n=0}^{W + k}x[n]^2 - \sum_{n=0}^{k}x[n]^2\\
&= \sum_{n=k+1}^{W + k}x[n]^2\\
&= \sum_{n=0}^{W - 1}x[n + k + 1]^2
\end{align}$$

Meaning that we end up with the same off-by-one error as seen before in the ACF correlation, which seem to cancel each other out to a degree.

The calculation needs to be adapted in any case to account for switching to the type II ACF $r'(k)$. We'll define 

$$\begin{align}E(k) &= \sum_{n=k}^{N-1}x[n]^2\\
&= \sum_{n=0}^{N-1}x[n]^2 - \sum_{n=0}^{k-1}x[n]^2\\
&= r'(0) - \sum_{n=0}^{k-1}x[n]^2
\end{align}$$

Following the logic of equation (7) from the YIN paper, the difference function can be calculated as

$$\begin{align}
d(k) &= r'(0) + E(k) - 2r'(k)\\
&= r'(0) + \left(r'(0) - \sum_{n=0}^{k-1}x[n]^2\right) - 2r'(k)\\
&= 2\left(r'(0) - r'(k)\right) - \sum_{n=0}^{k-1}x[n]^2
\end{align}$$

And going back from maths to code:

```python
yin_frames = np.cumsum(np.square(y_frames), axis=-2)

# Difference function: d(k) = 2 * (ACF(0) - ACF(k)) - sum_{m=0}^{k-1} y(m)^2
k = slice(1, max_period + 1)
yin_frames[..., 0, :] = 0
yin_frames[..., k, :] = (
    2 * (acf_frames[..., 0, :] - acf_frames[..., k, :]) - yin_frames[..., :k.stop-1, :]
)
``` 

The remaining steps (cumulative mean normalization, peak picking etc) are correct in the current implementation as far as I can tell.

## Benchmarks

I've run some comparisons using `mir_eval` that hopefully show that all this mathematical pedantry was worth it. The full numbers are contained [this Google sheet](https://docs.google.com/spreadsheets/d/15_eYfSo3drQSCM-cWmeYEokGsPsMtqqpg4E5PhnP5bs/edit?gid=0#gid=0), but I'll summarize the highlights here.


Using `fmin=80`, `fmax=800`, `sr=44100` and leaving all other values at their default, the new implementation improves the raw pitch accuracy of YIN (using a tolerance of 50 cents) on the [vocadito](https://zenodo.org/records/5578807) dataset significantly:

- from 95% to 98% using `frame_length=1024`
- from 92% to 98% using `frame_length=2048`

pYIN results also show clear improvements both for raw pitch accuracy and overall accuracy (taking into account V/UV decisions). 

The drawback is an increase in computational cost. The FFT is $O(N\log(N))$, and we've replaced 3 FFT operations with length $N/2$ with 2 operations of length $N$. That means at $N=2048$, we'd expect to take roughly 1.5x longer, which indeed matches the measurements for YIN. However, due to the more costly Viterbi decoding, this does not seem to have the same effect on pYIN running times.

## Additional changes

- I removed this line in both YIN and pYIN which marks $d'(0)$ as a trough if $d'(0) < d'(1)$: https://github.com/librosa/librosa/blob/d5aa7e1ae443c917fa53bef31055ed2613d63133/librosa/core/pitch.py#L612 This check is unnecessary since $d'(0) = 1$ by definition (equation 8 in the YIN paper), and a lag of 0 should never be selected as the period anyway.
- `__check_yin_parameters`
  - Since `win_length` is unused now, I've added a deprecation warning if it is called with any value other than `None`
  - It previously included this check: https://github.com/librosa/librosa/blob/d5aa7e1ae443c917fa53bef31055ed2613d63133/librosa/core/pitch.py#L967 I believe this is an error and should use `fmin`? As I understand it, the intention is an upper bound on the largest possible period, which corresponds to the lowest frequency, not the highest. I've changed this check to throw an error if the largest period `sr / fmin` exceeds `frame_length`.
  - Additionally there's a "softer" check (warning instead of error) for `sr / fmin >= frame_length // 2`, since ideally the frame should fit two periods of the lowest frequency for accurate f0 detection.

bmcfee: Thanks for this!  Dropping a note now to acknowledge that it's on my radar, but I won't have time to dig into it for a week or two.
bmcfee: :wave: thanks again @dsuedholt for the thorough investigation and writeup of your changes here.  Also thanks for your patience while I cleared my stack enough to put some thought into this PR!

Everything you've described makes total sense, and indeed clarifies a few points that popped up in the development of this implementation.

I've enabled the CI to run, and it looks like the changes here do introduce several regressions.  Some of these are down to changes in estimation - I suspect the disabled trough check here  - and others are due to broken handling of multi-channel.  I'll try to sleuth the latter part out to recommend a fix.



> * I believe this is an error and should use `fmin`? As I understand it, the intention is an upper bound on the largest possible period, which corresponds to the lowest frequency, not the highest. I've changed this check to throw an error if the largest period `sr / fmin` exceeds `frame_length`.

:eyes: Huh  We had this right in the issue thread where this was all derived: https://github.com/librosa/librosa/issues/1680#issuecomment-1487241064 ; this was then reverted in the PR at commit 6923fe70a96883a9f1c5b9378ce9dc62d4d1fa25 but I can't at all remember why this was a problem.



dsuedholt: Hi @bmcfee, thanks a lot for the helpful review! Multichannel processing should be fixed, tests are passing now.

I couldn't find the code that generated the ground truth data for the chirps that's being loaded from `tests/data/pitch-yin.npy` 🤔 I changed the tests to use the average instantaneous frequency for every frame as ground truth instead.
bmcfee: > I couldn't find the code that generated the ground truth data for the chirps that's being loaded from `tests/data/pitch-yin.npy`

I believe these outputs came from a prior yin implementation (from the vamp plugin, IIRC) so that we could regression-test for compatibility; see discussion in https://github.com/librosa/librosa/pull/1063 .  I'd prefer to keep that test as it was so that we don't lose track of compatibility.  If we need to relax the test somehow to account for differences introduced in this change, we can do that, but retain documentation (comments) about what changed and why.

I'm happy to have a free-standing unit test against the instantaneous frequency as a separate test though.

Thanks again for all this!
dsuedholt: Makes sense! I got the previous tests against the vamp plugin working now with a few adjustments. I reduced the framesize in the YIN call to 1024; before it was being called with the default `frame_length=2048`, but because of the previous `win_length` implementation, it was effectively only looking at 1024 samples.

For the `pyin` chirp test, I switched to the exact values used in the vamp plugin for fmin, fmax, and resolution, and it seems to pass the test for me.

I did notice one more thing though: the tests assumed that the first frame of the pyin output will always be unvoiced, because the decoding is initialized with a prior of 0 for all voiced bins and `1 / n_pitch_bins` for all unvoiced bins, instead of a uniform `1 / (2 * n_pitch_bins)` for all bins: https://github.com/librosa/librosa/blob/d5aa7e1ae443c917fa53bef31055ed2613d63133/librosa/core/pitch.py#L853-L854

I'm not sure why this decision was made, but it matches neither the original paper (defines uniform V/UV prior in equation 6) nor the Vamp plugin (initializes a uniform prior [here](https://code.soundsoftware.ac.uk/projects/pyin/repository/entry/MonoPitchHMM.cpp#L91)). 

I changed the pyin implementation to a uniform prior, which made the `test_pyin_multi_center` test pass for me when it previously didn't (I'm on a mac, I'm not sure if the mentioned OSX issue still persists).
codecov[bot]: ## [Codecov](https://app.codecov.io/gh/librosa/librosa/pull/1882?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa) Report
Attention: Patch coverage is `77.27273%` with `5 lines` in your changes missing coverage. Please review.
> Project coverage is 98.59%. Comparing base [(`d5aa7e1`)](https://app.codecov.io/gh/librosa/librosa/commit/d5aa7e1ae443c917fa53bef31055ed2613d63133?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa) to head [(`9c6eace`)](https://app.codecov.io/gh/librosa/librosa/commit/9c6eaceb757febc75b8fb427e3c28b7cbf627d7c?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa).

| [Files with missing lines](https://app.codecov.io/gh/librosa/librosa/pull/1882?dropdown=coverage&src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa) | Patch % | Lines |
|---|---|---|
| [librosa/core/pitch.py](https://app.codecov.io/gh/librosa/librosa/pull/1882?src=pr&el=tree&filepath=librosa%2Fcore%2Fpitch.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa#diff-bGlicm9zYS9jb3JlL3BpdGNoLnB5) | 77.27% | [5 Missing :warning: ](https://app.codecov.io/gh/librosa/librosa/pull/1882?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa) |

<details><summary>Additional details and impacted files</summary>


```diff
@@            Coverage Diff             @@
##             main    #1882      +/-   ##
==========================================
- Coverage   98.73%   98.59%   -0.14%     
==========================================
  Files          35       35              
  Lines        4656     4621      -35     
==========================================
- Hits         4597     4556      -41     
- Misses         59       65       +6     
```

| [Flag](https://app.codecov.io/gh/librosa/librosa/pull/1882/flags?src=pr&el=flags&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa) | Coverage Δ | |
|---|---|---|
| [unittests](https://app.codecov.io/gh/librosa/librosa/pull/1882/flags?src=pr&el=flag&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa) | `98.59% <77.27%> (-0.14%)` | :arrow_down: |

Flags with carried forward coverage won't be shown. [Click here](https://docs.codecov.io/docs/carryforward-flags?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa#carryforward-flags-in-the-pull-request-comment) to find out more.

</details>

[:umbrella: View full report in Codecov by Sentry](https://app.codecov.io/gh/librosa/librosa/pull/1882?dropdown=coverage&src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa).   
:loudspeaker: Have feedback on the report? [Share it here](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=librosa).

----
🚨 Try these New Features:

- [Flaky Tests Detection](https://docs.codecov.com/docs/test-result-ingestion-beta) - Detect and resolve failed and flaky tests
dsuedholt: I added tests for the warnings (should fix reduced code coverage) and removed a `numpy>2.0` only feature from one of the new chirp tests (should fix the linter failing).

I'm not sure what to do about the doc workflow failing, I can't get the error to happen on my local machine.
bmcfee: > I added tests for the warnings (should fix reduced code coverage) and removed a `numpy>2.0` only feature from one of the new chirp tests (should fix the linter failing).
> 

Thanks for this - it's all looking great!  I had one minor q about the revised tests, but I'm otherwise happy with where things are at.

> I'm not sure what to do about the doc workflow failing, I can't get the error to happen on my local machine.

The doc CI was failing due to a broken link.  This is almost surely unrelated to the PR, and it happens sporadically when sites we link to are temporarily unreachable.  As long as the tests that execute before link check are passing, I won't worry about it.
bmcfee: > I'm not sure why this decision was made, but it matches neither the original paper (defines uniform V/UV prior in equation 6) nor the Vamp plugin (initializes a uniform prior [here](https://code.soundsoftware.ac.uk/projects/pyin/repository/entry/MonoPitchHMM.cpp#L91)).

That's also a good point.  It must have escaped our attention in the initial CR of #1063 because the test signals generally started on unvoiced frames.  I couldn't find any discussion of it in the comment thread anyway.
bmcfee: I think we're all set with this.  I gave it a spin locally and it's all behaving nicely.

We do take a bit of a performance hit here: speed was about the same, but memory usage is up by about 10%.  I expect this comes from some inefficiencies in either `autocorrelate` or `abs2` (or both), relative to how it was previously implemented.  I can live with this for now - we can always revisit internal optimizations in a later point release.

Unless there are objections, I think this is good to merge.  Thanks again @dsuedholt for the deep dive and scrutiny on this!

Patch:
diff --git a/librosa/core/pitch.py b/librosa/core/pitch.py
index fdc8970b..fc2cd551 100644
--- a/librosa/core/pitch.py
+++ b/librosa/core/pitch.py
@@ -10,9 +10,11 @@ import numba
 
 from .spectrum import _spectrogram
 from . import convert
+from . import audio
 from .._cache import cache
 from .. import util
 from .. import sequence
+from ..util import Deprecated
 from ..util.exceptions import ParameterError
 from numpy.typing import ArrayLike
 from typing import Any, Callable, Optional, Tuple, Union
@@ -360,8 +362,6 @@ def piptrack(
 
 def _cumulative_mean_normalized_difference(
     y_frames: np.ndarray,
-    frame_length: int,
-    win_length: int,
     min_period: int,
     max_period: int,
 ) -> np.ndarray:
@@ -375,10 +375,6 @@ def _cumulative_mean_normalized_difference(
     ----------
     y_frames : np.ndarray [shape=(frame_length, n_frames)]
         framed audio time series.
-    frame_length : int > 0 [scalar]
-        length of the frames in samples.
-    win_length : int > 0 [scalar]
-        length of the window for calculating autocorrelation in samples.
     min_period : int > 0 [scalar]
         minimum period.
     max_period : int > 0 [scalar]
@@ -389,31 +385,26 @@ def _cumulative_mean_normalized_difference(
     yin_frames : np.ndarray [shape=(max_period-min_period+1,n_frames)]
         Cumulative mean normalized difference function for each frame.
     """
-    # Autocorrelation.
-    a = np.fft.rfft(y_frames, frame_length, axis=-2)
-    b = np.fft.rfft(y_frames[..., win_length:0:-1, :], frame_length, axis=-2)
-    acf_frames = np.fft.irfft(a * b, frame_length, axis=-2)[..., win_length:, :]
-    acf_frames[np.abs(acf_frames) < 1e-6] = 0
+    acf_frames = audio.autocorrelate(y_frames, max_size=max_period + 1, axis=-2)
 
     # Energy terms.
-    energy_frames = np.cumsum(y_frames**2, axis=-2)
-    energy_frames = (
-        energy_frames[..., win_length:, :] - energy_frames[..., :-win_length, :]
+    yin_frames = np.square(y_frames)
+    np.cumsum(yin_frames, out=yin_frames, axis=-2)
+
+    # Difference function: d(k) = 2 * (ACF(0) - ACF(k)) - sum_{m=0}^{k-1} y(m)^2
+    k = slice(1, max_period + 1)
+    yin_frames[..., 0, :] = 0
+    yin_frames[..., k, :] = (
+        2 * (acf_frames[..., 0:1, :] - acf_frames[..., k, :]) - yin_frames[..., :k.stop-1, :]
     )
-    energy_frames[np.abs(energy_frames) < 1e-6] = 0
-
-    # Difference function.
-    yin_frames = energy_frames[..., :1, :] + energy_frames - 2 * acf_frames
 
     # Cumulative mean normalized difference function.
     yin_numerator = yin_frames[..., min_period : max_period + 1, :]
     # broadcast this shape to have leading ones
-    tau_range = util.expand_to(
-        np.arange(1, max_period + 1), ndim=yin_frames.ndim, axes=-2
-    )
+    k_range = util.expand_to(np.r_[k], ndim=yin_frames.ndim, axes=-2)
 
     cumulative_mean = (
-        np.cumsum(yin_frames[..., 1 : max_period + 1, :], axis=-2) / tau_range
+        np.cumsum(yin_frames[..., k, :], axis=-2) / k_range
     )
     yin_denominator = cumulative_mean[..., min_period - 1 : max_period, :]
     yin_frames: np.ndarray = yin_numerator / (
@@ -489,7 +480,7 @@ def yin(
     fmax: float,
     sr: float = 22050,
     frame_length: int = 2048,
-    win_length: Optional[int] = None,
+    win_length: Optional[Union[int, Deprecated]] = Deprecated(),
     hop_length: Optional[int] = None,
     trough_threshold: float = 0.1,
     center: bool = True,
@@ -512,34 +503,39 @@ def yin(
     ----------
     y : np.ndarray [shape=(..., n)]
         audio time series. Multi-channel is supported..
+
     fmin : number > 0 [scalar]
         minimum frequency in Hertz.
         The recommended minimum is ``librosa.note_to_hz('C2')`` (~65 Hz)
         though lower values may be feasible.
+
     fmax : number > fmin, <= sr/2 [scalar]
         maximum frequency in Hertz.
         The recommended maximum is ``librosa.note_to_hz('C7')`` (~2093 Hz)
         though higher values may be feasible.
+
     sr : number > 0 [scalar]
         sampling rate of ``y`` in Hertz.
+
     frame_length : int > 0 [scalar]
         length of the frames in samples.
         By default, ``frame_length=2048`` corresponds to a time scale of about 93 ms at
         a sampling rate of 22050 Hz.
-    win_length : None or int > 0 [scalar]
-        length of the window for calculating autocorrelation in samples.
-        If ``None``, defaults to ``frame_length // 2``
+
     hop_length : None or int > 0 [scalar]
         number of audio samples between adjacent YIN predictions.
         If ``None``, defaults to ``frame_length // 4``.
+
     trough_threshold : number > 0 [scalar]
         absolute threshold for peak estimation.
+
     center : boolean
         If ``True``, the signal `y` is padded so that frame
         ``D[:, t]`` is centered at `y[t * hop_length]`.
         If ``False``, then ``D[:, t]`` begins at ``y[t * hop_length]``.
         Defaults to ``True``,  which simplifies the alignment of ``D`` onto a
         time grid by means of ``librosa.core.frames_to_samples``.
+
     pad_mode : string or function
         If ``center=True``, this argument is passed to ``np.pad`` for padding
         the edges of the signal ``y``. By default (``pad_mode="constant"``),
@@ -547,6 +543,12 @@ def yin(
         If ``center=False``,  this argument is ignored.
         .. see also:: `np.pad`
 
+    win_length : Deprecated
+        length of the window for calculating autocorrelation in samples.
+
+        .. warning:: This parameter is deprecated as of 0.11.0 and
+            will be removed in 1.0.
+
     Returns
     -------
     f0: np.ndarray [shape=(..., n_frames)]
@@ -571,12 +573,16 @@ def yin(
     if fmin is None or fmax is None:
         raise ParameterError('both "fmin" and "fmax" must be provided')
 
-    # Set the default window length if it is not already specified.
-    if win_length is None:
-        win_length = frame_length // 2
+    if not isinstance(win_length, Deprecated):
+        warnings.warn(
+            f"The win_length parameter has been deprecated in version 0.11.0 "
+            f"and has no effect. It will be removed in version 1.0.0.",
+            category=FutureWarning,
+            stacklevel=3,
+        )
 
     __check_yin_params(
-        sr=sr, fmax=fmax, fmin=fmin, frame_length=frame_length, win_length=win_length
+        sr=sr, fmax=fmax, fmin=fmin, frame_length=frame_length
     )
 
     # Set the default hop if it is not already specified.
@@ -597,11 +603,11 @@ def yin(
 
     # Calculate minimum and maximum periods
     min_period = int(np.floor(sr / fmax))
-    max_period = min(int(np.ceil(sr / fmin)), frame_length - win_length - 1)
+    max_period = min(int(np.ceil(sr / fmin)), frame_length - 1)
 
     # Calculate cumulative mean normalized difference function.
     yin_frames = _cumulative_mean_normalized_difference(
-        y_frames, frame_length, win_length, min_period, max_period
+        y_frames, min_period, max_period
     )
 
     # Parabolic interpolation.
@@ -650,7 +656,7 @@ def pyin(
     fmax: float,
     sr: float = 22050,
     frame_length: int = 2048,
-    win_length: Optional[int] = None,
+    win_length: Optional[Union[int, Deprecated]] = Deprecated(),
     hop_length: Optional[int] = None,
     n_thresholds: int = 100,
     beta_parameters: Tuple[float, float] = (2, 18),
@@ -681,51 +687,63 @@ def pyin(
     ----------
     y : np.ndarray [shape=(..., n)]
         audio time series. Multi-channel is supported.
+
     fmin : number > 0 [scalar]
         minimum frequency in Hertz.
         The recommended minimum is ``librosa.note_to_hz('C2')`` (~65 Hz)
         though lower values may be feasible.
+
     fmax : number > fmin, <= sr/2 [scalar]
         maximum frequency in Hertz.
         The recommended maximum is ``librosa.note_to_hz('C7')`` (~2093 Hz)
         though higher values may be feasible.
+
     sr : number > 0 [scalar]
         sampling rate of ``y`` in Hertz.
+
     frame_length : int > 0 [scalar]
         length of the frames in samples.
         By default, ``frame_length=2048`` corresponds to a time scale of about 93 ms at
         a sampling rate of 22050 Hz.
-    win_length : None or int > 0 [scalar]
-        length of the window for calculating autocorrelation in samples.
-        If ``None``, defaults to ``frame_length // 2``
+
     hop_length : None or int > 0 [scalar]
         number of audio samples between adjacent pYIN predictions.
         If ``None``, defaults to ``frame_length // 4``.
+
     n_thresholds : int > 0 [scalar]
         number of thresholds for peak estimation.
+
     beta_parameters : tuple
         shape parameters for the beta distribution prior over thresholds.
+
     boltzmann_parameter : number > 0 [scalar]
         shape parameter for the Boltzmann distribution prior over troughs.
         Larger values will assign more mass to smaller periods.
+
     resolution : float in `(0, 1)`
         Resolution of the pitch bins.
         0.01 corresponds to cents.
+
     max_transition_rate : float > 0
         maximum pitch transition rate in octaves per second.
+
     switch_prob : float in ``(0, 1)``
         probability of switching from voiced to unvoiced or vice versa.
+
     no_trough_prob : float in ``(0, 1)``
         maximum probability to add to global minimum if no trough is below threshold.
+
     fill_na : None, float, or ``np.nan``
         default value for unvoiced frames of ``f0``.
         If ``None``, the unvoiced frames will contain a best guess value.
+
     center : boolean
         If ``True``, the signal ``y`` is padded so that frame
         ``D[:, t]`` is centered at ``y[t * hop_length]``.
         If ``False``, then ``D[:, t]`` begins at ``y[t * hop_length]``.
         Defaults to ``True``,  which simplifies the alignment of ``D`` onto a
         time grid by means of ``librosa.core.frames_to_samples``.
+
     pad_mode : string or function
         If ``center=True``, this argument is passed to ``np.pad`` for padding
         the edges of the signal ``y``. By default (``pad_mode="constant"``),
@@ -733,6 +751,12 @@ def pyin(
         If ``center=False``,  this argument is ignored.
         .. see also:: `np.pad`
 
+    win_length : Deprecated
+        length of the window for calculating autocorrelation in samples.
+
+        .. warning:: This parameter is deprecated as of 0.11.0 and
+            will be removed in 1.0.
+
     Returns
     -------
     f0: np.ndarray [shape=(..., n_frames)]
@@ -773,12 +797,16 @@ def pyin(
     if fmin is None or fmax is None:
         raise ParameterError('both "fmin" and "fmax" must be provided')
 
-    # Set the default window length if it is not already specified.
-    if win_length is None:
-        win_length = frame_length // 2
+    if not isinstance(win_length, Deprecated):
+        warnings.warn(
+            f"The win_length parameter has been deprecated in version 0.11.0 "
+            f"and has no effect. It will be removed in version 1.0.0.",
+            category=FutureWarning,
+            stacklevel=3,
+        )
 
     __check_yin_params(
-        sr=sr, fmax=fmax, fmin=fmin, frame_length=frame_length, win_length=win_length
+        sr=sr, fmax=fmax, fmin=fmin, frame_length=frame_length
     )
 
     # Set the default hop if it is not already specified.
@@ -799,11 +827,11 @@ def pyin(
 
     # Calculate minimum and maximum periods
     min_period = int(np.floor(sr / fmax))
-    max_period = min(int(np.ceil(sr / fmin)), frame_length - win_length - 1)
+    max_period = min(int(np.ceil(sr / fmin)), frame_length - 1)
 
     # Calculate cumulative mean normalized difference function.
     yin_frames = _cumulative_mean_normalized_difference(
-        y_frames, frame_length, win_length, min_period, max_period
+        y_frames, min_period, max_period
     )
 
     # Parabolic interpolation.
@@ -850,8 +878,7 @@ def pyin(
     t_switch = sequence.transition_loop(2, 1 - switch_prob)
     transition = np.kron(t_switch, transition)
 
-    p_init = np.zeros(2 * n_pitch_bins)
-    p_init[n_pitch_bins:] = 1 / n_pitch_bins
+    p_init = np.ones(2 * n_pitch_bins) / (2 * n_pitch_bins)
 
     states = sequence.viterbi(observation_probs, transition, p_init=p_init)
 
@@ -944,13 +971,16 @@ def __pyin_helper(
 
 
 def __check_yin_params(
-    *, sr: float, fmax: float, fmin: float, frame_length: int, win_length: int
+    *, sr: float, fmax: float, fmin: float, frame_length: int
 ):
     """Check the feasibility of yin/pyin parameters against
     the following conditions:
 
     1. 0 < fmin < fmax <= sr/2
-    2. frame_length - win_length - 1 > sr/fmax
+    2. sr / fmin < frame_length - 1: At least one period of fmin needs to fit into the frame.
+
+    Give a warning if parameters violate:
+    3. sr / fmin < frame_length / 2: Ideally, at least two periods of fmin should fit into the frame.
     """
     if fmax > sr / 2:
         raise ParameterError(f"fmax={fmax:.3f} cannot exceed Nyquist frequency {sr/2}")
@@ -959,15 +989,21 @@ def __check_yin_params(
     if fmin <= 0:
         raise ParameterError(f"fmin={fmin:.3f} must be strictly positive")
 
-    if win_length >= frame_length:
+    if sr / fmin >= frame_length - 1:
+        fmin_feasible = sr / (frame_length - 1)
+        frame_length_feasible = int(np.ceil(sr / fmin) + 1)
         raise ParameterError(
-            f"win_length={win_length} must be less than frame_length={frame_length}"
+            f"fmin={fmin:.3f} is too small for frame_length={frame_length} and sr={sr}. "
+            f"Either increase to fmin={fmin_feasible:.3f} or frame_length={frame_length_feasible}"
         )
 
-    if frame_length - win_length - 1 <= sr // fmax:
-        fmax_feasible = sr / (frame_length - win_length - 1)
-        frame_length_feasible = int(np.ceil(sr/fmax) + win_length + 1)
-        raise ParameterError(
-            f"fmax={fmax:.3f} is too small for frame_length={frame_length}, win_length={win_length}, and sr={sr}. "
-            f"Either increase to fmax={fmax_feasible:.3f} or frame_length={frame_length_feasible}"
-        )
+    if sr / fmin >= frame_length // 2:
+        fmin_optimal = sr / (frame_length / 2)
+        frame_length_optimal = int(np.ceil(sr / fmin) * 2 + 1)
+
+        warnings.warn(
+            f"With fmin={fmin:.3f}, sr={sr} and frame_length={frame_length}, less than two periods of fmin "
+            f"fit into the frame, which can cause inaccurate pitch detection. "
+            f"Consider increasing to fmin={fmin_optimal:.3f} or frame_length={frame_length_optimal}.",
+            stacklevel=3,
+        )
\ No newline at end of file


Test patch:
diff --git a/tests/test_core.py b/tests/test_core.py
index 53803aad..24a1d328 100644
--- a/tests/test_core.py
+++ b/tests/test_core.py
@@ -1171,39 +1171,70 @@ def test_yin_tone(freq):
 
 
 def test_yin_chirp():
+    # test yin on a chirp, using output from the vamp plugin as ground truth
     y = librosa.chirp(fmin=220, fmax=640, duration=1.0)
-    f0 = librosa.yin(y, fmin=110, fmax=880, center=False)
+    f0 = librosa.yin(y, fmin=110, fmax=880, center=False, frame_length=1024, hop_length=512)
+
+    # adjust frames to the removal of win_length from yin
+    f0 = f0[:-2]
+
     target_f0 = np.load(os.path.join("tests", "data", "pitch-yin.npy"))
     assert np.allclose(np.log2(f0), np.log2(target_f0), rtol=0, atol=1e-2)
 
 
+def test_yin_chirp_instant():
+    # test yin on a chirp, using frame-wise instantaneous frequency as ground truth
+    sr = 22050
+    chirp_min, chirp_max = 220, 640
+
+    t = np.arange(sr) / sr
+    f = chirp_min * (chirp_max / chirp_min) ** t
+
+    fl = 2048
+    hl = 512
+
+    y = librosa.chirp(fmin=chirp_min, fmax=chirp_max, sr=sr, duration=1.0, linear=False)
+    target_f0 = librosa.util.frame(f, frame_length=fl, hop_length=hl).mean(axis=0)
+
+    f0 = librosa.yin(y, fmin=110, fmax=880, sr=sr, frame_length=fl, hop_length=hl, center=False)
+    assert np.allclose(np.log2(f0), np.log2(target_f0), rtol=0, atol=1e-2)
+
+
 @pytest.mark.xfail(raises=librosa.ParameterError)
 @pytest.mark.parametrize(
-    "fmin,fmax,win_length,frame_length",
+    "fmin,fmax,frame_length",
     [
-        (None, None, None, 2048),  # neither
-        (110, None, None, 2048),  # no fmax
-        (None, 880, None, 2048),  # no fmin
-        (110, 880, 2049, 2048),  # win_length >= frame_length
-        (-1, 440, None, 2048),  # Negative fmin
-        (440, 220, None, 2048),  # fmin > fmax
-        (440, 16000, None, 2048),  # fmax > nyquist
-        (10, 21, None, 2048),  # frame_length - win_length - 1 <= sr/fmax
+        (None, None, 2048),  # neither
+        (110, None, 2048),  # no fmax
+        (None, 880, 2048),  # no fmin
+        (-1, 440, 2048),  # Negative fmin
+        (440, 220, 2048),  # fmin > fmax
+        (440, 16000, 2048),  # fmax > nyquist
+        (10, 21, 2048),  # sr / fmin >= frame_length - 1
     ],
 )
-def test_yin_fail(fmin, fmax, win_length, frame_length):
+def test_yin_fail(fmin, fmax, frame_length):
     y = librosa.tone(110, duration=1.0)
     librosa.yin(
-        y, fmin=fmin, fmax=fmax, win_length=win_length, frame_length=frame_length
+        y, fmin=fmin, fmax=fmax, frame_length=frame_length
     )
 
+def test_yin_warn():
+    y = librosa.tone(110, duration=1.0)
+
+    # win_length is deprecated
+    with pytest.warns(FutureWarning, match="deprecated"):
+        librosa.yin(y, fmin=110, fmax=1000, win_length=1024)
+
+    # sr / fmin >= frame_length // 2
+    with pytest.warns(UserWarning, match="two periods"):
+        librosa.yin(y, fmin=20, fmax=1000)
 
 @pytest.mark.parametrize("freq", [110, 220, 440, 880])
 def test_pyin_tone(freq):
     y = librosa.tone(freq, duration=1.0)
     f0, _, _ = librosa.pyin(y, fmin=110, fmax=1000, center=False)
-    # Skip the first frame, since the voicing prior does not allow it
-    assert np.allclose(np.log2(f0[1:]), np.log2(freq), rtol=0, atol=1e-2)
+    assert np.allclose(np.log2(f0), np.log2(freq), rtol=0, atol=1e-2)
 
 
 def test_pyin_multi():
@@ -1256,9 +1287,17 @@ def test_pyin_multi_center():
 
 
 def test_pyin_chirp():
+    # test yin on a chirp, using output from the vamp plugin as ground truth
     y = librosa.chirp(fmin=220, fmax=640, duration=1.0)
     y = np.pad(y, (22050,))
-    f0, voiced_flag, _ = librosa.pyin(y, fmin=110, fmax=880, center=False)
+
+    # default values as set in https://code.soundsoftware.ac.uk/projects/pyin/repository/
+    f0, voiced_flag, _ = librosa.pyin(y, fmin=60, fmax=900, center=False, frame_length=1024, hop_length=512, resolution=0.2)
+
+    # adjust frames to the removal of win_length from yin
+    f0 = f0[:-2]
+    voiced_flag = voiced_flag[:-2]
+
     target_f0 = np.load(os.path.join("tests", "data", "pitch-pyin.npy"))
     # test if correct frames are voiced
     assert np.array_equal(voiced_flag, target_f0 > 0)
@@ -1267,23 +1306,64 @@ def test_pyin_chirp():
         np.log2(f0[voiced_flag]), np.log2(target_f0[target_f0 > 0]), rtol=0, atol=1e-2
     )
 
+def test_pyin_chirp_instant():
+    # test pyin on a chirp, using frame-wise instantaneous frequency as ground truth
+    sr = 22050
+    chirp_min, chirp_max = 220, 640
+
+    t = np.arange(sr) / sr
+    f = chirp_min * (chirp_max / chirp_min) ** t
+    f = np.pad(f, (sr,))
+
+    fl = 2048
+    hl = 512
+
+    target_f0 = librosa.util.frame(f, frame_length=fl, hop_length=hl)
+    target_f0 = target_f0.mean(axis=0, where=target_f0 > 0)
+
+    y = librosa.chirp(fmin=chirp_min, fmax=chirp_max, sr=sr, duration=1.0, linear=False)
+    y = np.pad(y, (sr,))
+
+    f0, voiced_flag, _ = librosa.pyin(y, fmin=110, fmax=880, frame_length=fl, hop_length=hl, center=False)
+
+    # test if correct frames are voiced
+    assert np.array_equal(voiced_flag, target_f0 > 0)
+
+    # test voiced frames are within one cent of the target
+    cents = np.log2(f0[voiced_flag])
+    target_cents = np.log2(target_f0[target_f0 > 0])
+
+    assert np.allclose(
+        np.log2(cents[1:-1]), np.log2(target_cents[1:-1]), rtol=0, atol=1e-2
+    )
+
+    # higher tolerance for the first and last frames, accounting for abrupt start / end
+    assert np.abs(cents[0] - target_cents[0]) <= 1e-1
+    assert np.abs(cents[-1] - target_cents[-1]) <= 1e-1
+
 
 @pytest.mark.xfail(raises=librosa.ParameterError)
 @pytest.mark.parametrize(
-    "fmin,fmax,win_length,frame_length",
+    "fmin,fmax,frame_length",
     [
-        (None, None, None, 2048),
-        (110, None, None, 2048),
-        (None, 880, None, 2048),
-        (110, 880, 2049, 2048),
+        (None, None, 2048),
+        (110, None, 2048),
+        (None, 880, 2048),
     ],
 )
-def test_pyin_fail(fmin, fmax, win_length, frame_length):
+def test_pyin_fail(fmin, fmax, frame_length):
     y = librosa.tone(110, duration=1.0)
     librosa.pyin(
-        y, fmin=fmin, fmax=fmax, win_length=win_length, frame_length=frame_length
+        y, fmin=fmin, fmax=fmax, frame_length=frame_length
     )
 
+def test_pyin_warn():
+    y = librosa.tone(110, duration=1.0)
+
+    # win_length is deprecated
+    with pytest.warns(FutureWarning, match="deprecated"):
+        librosa.pyin(y, fmin=110, fmax=1000, win_length=1024)
+
 
 @pytest.mark.parametrize("freq", [110, 220, 440, 880])
 @pytest.mark.parametrize("n_fft", [1024, 2048, 4096])