# Robust Multicomponent Tracking of Ultrasonic Vocalizations

Reyhaneh Abbasi, Nicki Holighaus, Peter Balazs, Vincent Lostanlen, Dustin

J. Penn, Sarah M. Zala

► To cite this version:

Reyhaneh Abbasi, Nicki Holighaus, Peter Balazs, Vincent Lostanlen, Dustin J. Penn, et al. Robust Multicomponent Tracking of Ultrasonic Vocalizations. IEEE International Conference on Acoustics, Speech, and Signal Processing (ICASSP), Apr 2025, Hyderabad, India. hal-05002401

HAL Id: hal-05002401

https://hal.science/hal-05002401v1

Submitted on 23 Mar 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Robust Multicomponent Tracking
of Ultrasonic Vocalizations

Reyhaneh Abbasi*, Nicki Holighaus, Peter Balazs

Acoustics Research Institute

Austrian Academy of Sciences

Vienna, Austria

Vincent Lostanlen

Nantes Université, École centrale Nantes

CNRS, LS2N, UMR 6004

F-44000 Nantes, France

Dustin J. Penn, Sarah M. Zala

Konrad Lorenz Institute of Ethology

University of Veterinary Medicine

Vienna, Austria

**Abstract**—Ultrasonic vocalizations (USV) convey information about individual identity and arousal status in mice. We propose to track USV as ridges in the time-frequency domain via a variant of time-frequency reassignment (TFR). The key idea is to perform TFR with empirical Wiener shrinkage and multitapering to improve robustness to noise. Furthermore, we perform TFR over both the short-term Fourier transform and the constant-Q transform so as to detect both the fundamental frequency and its harmonic partial (if any). Experimental results show that our approach effectively estimates multicomponent ridges with high precision and low frequency deviation.

**Index Terms**—Ridge tracking, Multitaper reassignment, Ultrasonic Vocalizations, Bioacoustics

I. INTRODUCTION

Ultrasonic vocalizations (USVs) in mice are considered complex due to their diverse acoustic properties and context-dependent variability. Prior studies [1], [2] have shown that analyzing USVs provide valuable insights about emotional states and identity. Classifying USVs helps these analysis by examining the diversity of USV types and their sequences [3], [4]. Practitioners label USVs manually by inspecting spectrograms [3], [5] or automatically using classification tools such as DeepSqueak [6] and BootSnap [7].

Both manual and automatic classifications rely on spectro-temporal features such as duration, frequency modulation, and complexity. Many of those features are derived from ridges, i.e., the local maxima in the time-frequency domain. Despite the value of ridges in USV classification, the prevalence of low signal-to-noise ratio (SNR) in USV recordings poses a significant challenge for their extraction. Furthermore, current methods for USV ridge extraction only track the fundamental frequency, despite that harmonic vocalizations are closely timed with mice sexual behavior [8], [9]. DeepSqueak, a tool for detection and analysis of USVs, uses a ridge detector based on identifying local amplitude maxima in the time-frequency domain. As an alternative, the authors of [10] propose a method that detects USV ridges through thresholding of smoothed spectrograms. While both methods were shown to be effective at high levels of SNR, they are sensitive to noise and inadequate for multicomponent tracking.

Our work addresses these limitations by introducing a novel method designed to improve the accuracy of USV ridge estimation using more advanced time-frequency processing techniques. Linear time-frequency representations such as short-term Fourier transforms (STFT) or constant-Q transforms (CQT) have limited time-frequency resolution. A data-driven method known as time-frequency reassignment (TFR) has the potential to overcome this limitation as it sharpens

the outline of ridges via local estimates of instantaneous frequency and group delay [11], [12]. Yet, the basic formulation of TFR lacks robustness to noise, preventing its applicability to USV analysis.

In this paper, we present a new method for multicomponent tracking of USVs. Our method comprises three preprocessing and postprocessing techniques in addition to TFR: empirical Wiener shrinkage, multitapering, and CQT-based autocorrelation. While the former two techniques are well-known improvements of traditional TFR, multitaper wavelet/CQT reassignment has not been previously reported in the literature, either for USV analysis or broader applications. Indeed, the key idea is to take advantage of the equivariance property in the CQT with respect to the frequency transposition of pure tones. Given an STFT-based estimate of the fundamental ridge and a reassigned CQT, we propose to compute the autocorrelation over adjacent log-frequency bins in the CQT so as to detect multiple components in a way that is both data-driven and informed by prior knowledge about harmonic spectra. To our knowledge, our analysis here demonstrates the first successful application of CQT on mouse USVs, and one of the first in bioacoustics at large.

Experimental results show that the developed approach estimates ridges and harmonic partials with high precision and low mean frequency deviation. Our method will offer a significant improvement in the efficiency of USV analysis, particularly in noisy environments or when harmonic components are present. Another key advantage of the proposed method is that it does not require any hand-annotated data for parametrization.

II. TECHNICAL BACKGROUND

In this section, we recall the central concepts and methods from the time-frequency processing used in the proposed method. Specifically, we compute two different time-frequency transforms, to be used in different processing steps: a STFT [13] and a CQT [14]. Whereas the STFT provides a fixed time and frequency resolution and operates on a linear frequency scale, the time respectively frequency resolution of the CQT are anti-proportional, respectively proportional, to the center frequency of each filter. This is due to the CQT's constant center-frequency-to-bandwidth ratio. It further operates on a logarithmic frequency scale and is therefore particularly suitable for the analysis of harmonic tones.

The coefficients associated with both filter banks can be interpreted as the scalar products of the signal $x \in \mathbb{C}^L$ of interest with a family of time-frequency atoms. Specifically, the STFT is given by

$$V_g x[n, m] = \langle x, g_{n,m} \rangle, \text{ with } g_{n,m}(\bullet) = g[\bullet - na]e^{2\pi i(\bullet)m/M}, \quad (1)$$

where the hop size (or decimation factor) $a$ and the FFT length $M$ are divisors of $L$, $n = 0, \dots, L/a - 1$ and $m = 0, \dots, M - 1$, and $g \in \mathbb{C}^L$ is the *Hann window* function. In this work, we set $a = 150$ (samples) and the FFT length $M = 750$ equals the support length

<footnote>* Corresponding author: Reyhaneh Abbasi (reyhaneh.abbasi@oeaw.ac.at)</footnote>

<footnote>This work was supported by the FWF projects P36446 (R. Abbasi, S. Zala), P34922-N (P. Balazs), P28141-B25 (D. Penn), the WWTF project LS23-014 (N. Holighaus), and the ANR project ANR-23-CE37-0025 (V. Lostanlen). The authors are ordered according to their contribution. We would like to thank Teresa Klaus for providing the mouse recordings.</footnote>
of the window g, corresponding to 2.5 ms. We further denote by $S_g x = |V_g x|^2$ the (STFT) spectrogram of $x$.

The CQT is constructed from bandpass filters $\psi_k \in \mathbb{C}^L$, $k = 0, \dots, K-1$ centered at frequencies $2^{k/B}\xi_{\min}$, where B is the number of frequency bins (or channels) per octave covered by the CQT and $\xi_{\min}$ is the minimum frequency of interest. To ensure the desired constant center-frequency-to-bandwidth ratio, it is common to derive the $\psi_k$ by sampling appropriate dilations of a fixed, continuous bandpass filter $\psi$, revealing the close relation of the CQT to wavelet transforms [15]. Here, we additionally use uniform decimation, such that the CQT is given by

$$C_{\psi}x[n, k] = \langle x, \psi_{n,k} \rangle, \quad \text{with} \quad \psi_{n,k} = \psi_k[\bullet - na]. \quad (2)$$

Specifically, we set $a = 150$, $\xi_{\min} = 28$ kHz, $B = 512$ and $K = 1240$ in order to cover the entire frequency range from $\xi_{\min}$ up to the Nyquist frequency of 150 kHz. Finally, we denote by $S_{\psi}x = |C_{\psi}x|^2$ the CQ spectrogram of $x$.

### A. Multitaper Reassignment

Like for any linear filter bank, the local (joint) time-frequency resolution of the STFT and CQT is fundamentally limited by the uncertainty principle [16], which manifests as a smearing in the respective spectrogram that cannot be arbitrarily reduced. Spectrogram reassignment is a technique that sharpens spectrograms by re-allocating time-frequency energy based on its local instantaneous frequency and group delay.

Specifically, every STFT coefficient $V_g x[n, m]$ represents the time-frequency content of $x$ in an area around a specific time-frequency position that we denote by $(t_n, \xi_m)$. Let $t_0[n, m]$ and $\xi_0[n, m]$ be estimates of the local group delay and instantaneous frequency for the coefficient $V_g x[n, m]$ and define

$$n_0[n, m] := \arg \min_{n \in \mathbb{Z}} |t_0[n, m] - t_n| \quad \text{and}$$

$$m_0[n, m] := \arg \min_{m \in \mathbb{Z}} |\xi_0[n, m] - \xi_n|.$$

Then, with $L_{n,m} := \{(l, k) \in \mathbb{Z}^2 : (n_0[l, k], m_0[l, k]) = (n, m)\}$, the reassigned (STFT) spectrogram is defined as

$$S_g^{(R)} x[n, m] = \sum_{(\tilde{n}, \tilde{m}) \in L_{n,m}} S_g x[\tilde{n}, \tilde{m}]. \quad (3)$$

The definition of the reassigned CQ spectrogram $S_{\psi}^{(R)}x[n, m]$ is analogous. Although different techniques for estimating the local group delay and instantaneous frequency exist, this work uses the method proposed in [11] for STFT reassignment and the generalization from [12] for reassignment of the CQT.

To further attenuate noise spectrograms and improve robustness, we use multitaper spectrogram reassignment, which was originally proposed only for the STFT by [17]. Using results from [12], its generalization to the CQT is straightforward. By taking the mean of multiple reassigned spectrograms obtained with orthogonal windows (tapers), the resulting representation emphasizes deterministic signal components while simultaneously reducing variance of background noise. As tapers, we select a fixed number $s \in \mathbb{N}$ of windows $g^{(0)}, \dots, g^{(s-1)}$, resp. prototype filters $\psi^{(0)}, \dots, \psi^{(s-1)}$, derived from orthogonal Hermite functions, as in [17]. In this work, we use $s = 6$ tapers for both the STFT and CQT. The multitaper reassigned spectrograms are then given by

$$V_x^{(R)}[n, m] = s^{-1} \cdot \sum_{j=0}^{s-1} S_{g^{(j)}}^{(R)} x[n, m] \quad (4)$$

and

$$C_x^{(R)}[n, k] = s^{-1} \cdot \sum_{j=0}^{s-1} S_{\psi^{(j)}}^{(R)} x[n, k], \quad (5)$$

respectively. As a final step to further refine the reassigned CQ spectrogram, we implemented a modified version of the pre-denoising step proposed in [17, Section III-B], where we substitute the binary mask by a smooth approximation.

For all calculations involving the CQT, we use the implementations provided in the LTFAT package [18], [19], which are based on [20] and [12].

### III. DATA AND METHOD

Fig. 1 outlines our proposed method for estimating the fundamental frequency ($\xi_0$) and first harmonic partial of USVs. The approach begins by denoising the audio signal with a Wiener filter. Subsequently, multitapered reassigned spectrograms are computed to identify the ridge, through path optimization. Finally, multitapered reassigned CQTs are computed to extract the harmonic partial component using autocorrelation analysis. In the subsections below, we first prescribe the data used (III-A), followed by the different subsections of the proposed method III-B, III-C, III-D, and the evaluation approaches (III-E).

#### A. Data

The USVs used in this study are part of a dataset from a previous study described in [9]. These USVs were recorded from adult wild-derived house mice (Mus musculus musculus) using the RECORDER USGH-software (Avisoft-RECORDER Version 4.2) at a sampling rate of 300 kHz. The USVs analyzed in this study range from 30 kHz to 120 kHz, with durations between approximately 3 ms and 200 ms. The USVs were classified by biologists experienced with this task into 9 distinct categories based on various acoustic features, as described in [9], i.e., upward ('UP'), down ('D'), 'C2' for 2 component USVs, 'C3' for 3-component USVs, complex ('C'), U-shaped and inverted U-shaped signals labeled as 'U' and 'UI', flat ('F'), and harmonic USVs were denoted as 'H'. For classifying USVs, spectrograms as visual information and A-MUD as automatic detector of USVs [21] were used.

To evaluate our proposed method, described later in III-B,III-C, and III-D, hand-annotated ridges are required. We prepared these annotations by initially estimating ridges for 197 USVs using our proposed method. These samples included a range of low- and high-SNR cases. The initial estimates were then manually refined post-hoc using a custom Python GUI.

#### B. Denoising

To reduce background noise in the USV signals, we apply empirical Wiener shrinkage [22] in the STFT domain, a diagonal approximation of an optimal Wiener filter [23]. An estimate of the power spectral density of the stationary background noise is obtained from the spectrogram of a 2-second recording without the presence of USVs.

#### C. Ridge Estimation

We determine the frequency ridge associated with the fundamental frequency $\xi_0$ using the standard algorithm implemented in the MATLAB function tfridge.

It determines a ridge from the reassigned spectrogram $V_x^{(R)}$ as follows. Define the matrix $\mathcal{L}[n, m]$ recursively by $\mathcal{L}[0, m] = -\ln V_x^{(R)}[0, m]$ and, for $n \ge 1$, set

$$\mathcal{L}[n, m] = -\ln V_x^{(R)}[n, m] + \left( \min_{\tilde{m}=0}^{M-1} L[n-1, \tilde{m}] + \alpha \cdot |m - \tilde{m}|^2 \right).$$
<figure><img src="image_2.png" /><figcaption>Fig. 1. Flowchart of proposed method for multicomponent tracking of USVs.</figcaption></figure>

Here, the parameter $\alpha > 0$ controls how strictly frequency jumps are penalized. Further, define $\mathcal{M}[n, m]$ by $\mathcal{M}[0, m] = 0$ and, for $n \ge 1$,

$$ \mathcal{M}[n, m] = \arg \min_{\tilde{m}=0}^{M-1} L[n-1, \tilde{m}] + \alpha \cdot |m - \tilde{m}|^2. $$

A ridge ($m_0, \dots, m_{N-1}$), where $N-1 := L/a - 1$ is the time index of the final time frame, is then recursively selected by setting $m_{N-1} = \arg \min_{\tilde{m}=0}^{M-1} L[N-1, \tilde{m}]$ and $m_n = M[n+1, m_{n+1}]$ for $n = 0, \dots, N-2$.

In this study, the penalty parameter $\alpha$ is determined heuristically as $\alpha = 6 \times 10^{-4}$, which yielded satisfactory results for all USVs in our data set. After detecting ridges, we refined them by setting a threshold of 10 on the distance between frequency bins. If the distance between adjacent ridge points exceeded this threshold, those points were replaced with NaN values. Subsequently—for those missing values—we applied a linear function to smooth the detected ridges.

## D. Harmonic Partial Component Estimation

The proposed scheme for detection and estimation of harmonic components relies on the simple observation that, due to its logarithmic frequency scale, the CQT is equivariant under transposition. Specifically, if we consider a USV with time-dependent fundamental frequency $\xi_0[n]$ and partials at $\xi_j[n] = c_j\xi_0[n]$, $j = 1, \dots, J$, then the multiplication property of the logarithm yields $\log_2(\xi_j[n]) = \log_2(c_j\xi_0[n]) = \log_2(\xi_0[n]) + \log_2(c_j)$. As a consequence, the distance between partials in the CQ spectrogram (and its reassigned variants) does not depend on the fundamental frequency $\xi_0[n]$. Furthermore, partials are usually expected close to integer multiples of the fundamental frequency, i.e., $c_j \approx j + 1$, such that, e.g., the first partial (if present) is expected to be roughly one octave higher, yielding a ridge at a distance of $\approx B$ bins to the ridge associated with the fundamental.

Autocorrelation analysis is a simple method for detecting repeating patterns in data, e.g., for pitch estimation in speech signals [24]. Here, we apply it to determine presence of partials in the n-th CQT time frame: For a clean signal, this only involves computing the autocorrelation of the column vector $C_x^{(R)}[n, \bullet]$ and finding maxima in the vicinity of integer multiples of the fundamental frequency. In noisy signals, this process is expected to be unstable, since the background noise might introduce spurious peaks in the autocorrelation. Hence, we propose to use a weighted average over adjacent time frames to ensure that only correlations that persist over some span of time are considered. If $R_n := R(C_x^{(R)}[n, \bullet])$ denotes the autocorrelation of $C_x^{(R)}[n, \bullet]$, then we propose to detect partials via the features

$$ \mathbf{R}_n = C^{-1} \sum_{j=-J}^{J} c_{|j|} R_{n+j}, $$

where $c_0 = 1$, $0 \le c_j \le 1$ and $C = 1 + 2 \sum_{j=1}^{J} c_j$. For this work, we only consider the simplest case $J = 1$ and $c_0 = c_1 = 1$.

Furthermore, since the sampling rate of our USV recordings is not high enough to resolve higher order partials, we restrict the search to the first partial.

From the autocorrelation vector $\mathbf{R}_n$, a candidate bin for the first partial is determined as follows: Set $\mathbf{R}_n[k] = 0$ for $k = 0, \dots, 60, K - 60, \dots K - 1$ and determine $k_1$ as $k_1 = \arg \max_{\tilde{m}=0}^{K-1} \mathbf{R}_n[k]$. If the center frequency $\tilde{\xi} = 2^{k_1/B} \xi_{min}$ associated with the bin $k_1$ is approximately $2\xi_0[n]$, then we consider $\tilde{\xi}$ a candidate for $\xi_1[n]$. If a maximum deviation of 3 kHz is exceeded, we set $\xi_1[n] = NaN$. After performing the above steps for all time frames $n$, we retain only sequences of 5 or more consecutive non-NaN values to represent the ridge of the first partial. The code implementation is available at https://github.com/ReyhanehAbbasi/MultiComponent-USV-Tracking.

## E. Evaluation

We compared the performance of the proposed method with DeepSqueak's USV ridge estimation tool, using both non-denoised audio and reassigned spectrograms of denoised audio. Additionally, we conduct an ablation study to investigate the impact of key modules in the proposed method on performance. The estimated $\xi_0$ and harmonic component $\xi_1$ (if present) for each USV were compared to hand-annotated data. The performance of the methods was evaluated in terms of precision and mean frequency deviation. An estimated ridge was considered a match (true positive) to a manually annotated ridge if the mean difference between the estimated ridge point and the ground truth ridge point did not exceed 500 Hz and a false positive otherwise. For the harmonic component, the same rule applies; however, the threshold for qualifying as true positive is 3 kHz. These thresholds were determined based on preliminary experiments to ensure accurate capture of the relevant features of the USVs. Additionally, false positives may occur if the manual annotations indicates the absence of the first partial and false negatives if a partial is present, but not detected. Hence, we additionally report recall values for harmonic detection. Since DeepSqueak does not provide harmonic estimation, we report metrics for harmonic components only using our proposed method.

# IV. RESULTS

As shown in Tab. I, the proposed method achieves the highest overall precision (89.6%) and the second lowest mean frequency deviation (1.5±5.1 kHz). The incremental improvements observed from rows 3 to 6 highlight the benefits of key modules incorporated into our method: denoising (ridge STFT denoised) improves ridge detection compared to the original audio (ridge STFT), and reassigned STFT (ridge RS STFT denoised) further enhances performance. The proposed approach, which combines denoising with multitapered reassigned STFT, achieves the best results. Also, applying reassigned STFT to denoised audio improved DeepSqueak's baseline performance.
<table><caption>TABLE I
OVERALL AND CLASS-WISE EVALUATION METRIC FOR RIDGE ESTIMATION METHODS.</caption><thead><tr><th rowspan="2">Method</th><th colspan="2">Overall Evaluation Metrics</th><th colspan="9">Class-wise Precision (%)</th></tr><tr><th>Precision (%)</th><th>Deviation (kHz)</th><th>C</th><th>C2</th><th>C3</th><th>D</th><th>F</th><th>H</th><th>U</th><th>UI</th><th>UP</th></tr></thead><tbody><tr><td>DeepSqueak</td><td>80.3</td><td>6.8 ± 11.2</td><td>82.9</td><td>84.3</td><td>80.1</td><td>67.0</td><td>48.7</td><td>86.5</td><td>68.9</td><td>84.0</td><td>63.2</td></tr><tr><td>DeepSqueak RS denoised</td><td>82.8</td><td>3.8 ± 6.4</td><td>85.2</td><td>84.5</td><td>80.2</td><td>72.8</td><td>63.6</td><td>88.5</td><td>74.4</td><td>87.2</td><td>70.2</td></tr><tr><td>ridge STFT</td><td>84.5</td><td>3 ± 8.1</td><td>86.5</td><td>86.6</td><td>83.5</td><td>70.7</td><td>62.8</td><td>92.8</td><td>65.1</td><td>88.4</td><td>72</td></tr><tr><td>ridge STFT denoised</td><td>85.7</td><td>1.7 ± 4.1</td><td>87.9</td><td>86.7</td><td>82.7</td><td>80.2</td><td>76.5</td><td>92.1</td><td>72.7</td><td>89.1</td><td>71.8</td></tr><tr><td>ridge RS STFT denoised</td><td>86</td><td>1.5 ± 3.7</td><td>87.1</td><td>87.4</td><td>83.6</td><td>78.3</td><td>75.5</td><td>92.9</td><td>84.9</td><td>89.1</td><td>76.2</td></tr><tr><td>Proposed</td><td>89.6</td><td>1.5 ± 5.1</td><td>91.5</td><td>91.0</td><td>87.4</td><td>80.6</td><td>78.5</td><td>96</td><td>82.8</td><td>92.3</td><td>79.2</td></tr></tbody></table>

<figure><img src="image_3.png" /><figcaption>Fig. 2. Comparison of ridge estimation methods for USVs under varying SNR conditions.</figcaption></figure>

Examining the class-wise precision in Tab. I, the proposed method shows substantial improvements across all classes compared to all methods. We note particularly large improvements compared to DeepSqueak in the classes D, F, U and UP. Although an in-depth analysis is beyond the scope of this work, a first inspection indicates that in our data, USV recordings associated with these classes have a lower average SNR, providing a further indication that our method is more robust to noise than the baseline. Fig. 2 shows a comparison of the ridge estimation methods for USV under varying SNR conditions. The proposed method consistently outperforms the DeepSqueak and DeepSqueak RS approaches across all SNR conditions. DeepSqueak has more difficulties to accurately estimate ridges for decreasing SNRs. Incorporating a reassigned spectrogram (DeepSqueak RS) enhances its performance, demonstrating the benefits of this technique for ridge estimation in noisy environments. Despite these improvements, in extremely low-SNR scenarios (third column), only the proposed method estimates the ridges correctly.

The results for the partial harmonic is shown in Tab. II. Similar to fundamental frequency detection, harmonic partial detection shows incremental improvements with the successive signal processing steps in our method. Fig. 3 shows harmonic estimation results for USVs under varying precisions and recalls. The first column displays a USV with moderate precision (62%, primarily due to high false positives) and high recall (81%). In the second column, the spectrogram shows a USV with high precision (97%) but lower recall (63%) as a result of high false negatives. The third column shows a USV with both high precision (100%) and high recall (86%).

V. CONCLUSION

Ultrasonic vocalizations (USVs) provide valuable insights into animal communication and their classification is key to analysis. USV classification is based on ridge estimation, which is challenging, especially at realistic signal-to-noise ratios (SNR). We propose a new method based on time-frequency processing of USV signals

<table><caption>TABLE II
EVALUATION METRICS FOR HARMONIC PARTIAL ESTIMATION.
PRECISION AND RECALL ARE IN %, AND DEVIATION IS IN KHZ. AC STANDS FOR AUTOCORRELATION.</caption><thead><tr><th rowspan="2">Input</th><th colspan="3">Evaluation Metric</th></tr><tr><th>Precision</th><th>Recall</th><th>Deviation</th></tr></thead><tbody><tr><td>AC CQT</td><td>73</td><td>51</td><td>1 ± 0.27</td></tr><tr><td>AC CQT denoised</td><td>80.7</td><td>54</td><td>0.8 ± 0.26</td></tr><tr><td>AC RS CQT denoised</td><td>82.3</td><td>60</td><td>0.8 ± 0.24</td></tr><tr><td>Proposed</td><td>89</td><td>59</td><td>0.8 ± 0.15</td></tr></tbody></table>

<figure><img src="image_4.png" /><figcaption>Fig. 3. harmonic partial estimated by the proposed method for three harmonic USVs under varying precisions and recalls.</figcaption></figure>

that significantly improves precision and frequency deviation for detecting $\xi_0$ and, furthermore, detects the first harmonic partial$\xi_1$ with high precision and low frequency deviation. Harmonic USVs present a particular challenge for automatic classification, especially during sexual and social interactions when high levels of background noise and competing sounds, such as female broadband vocalizations, are present. In contrast to prior methods, which usually neglected harmonics entirely, our method facilitates the correct classification of USVs of harmonic type.

For our hand-annotated test data, the harmonic detection scheme achieves a moderate recall. Future work should incorporate advanced noise reduction methods to improve harmonic detection. On the other hand, annotating harmonics is challenging and prone to inaccuracies due to low SNR at high frequencies. In future work, we will verify harmonic detection using synthetic data with known ground truth. We will further integrate the proposed method into USV classification algorithms for enhanced performance. In order to improve computational efficiency in that context, all components of the method will be transitioned to a CQT-based framework. Although the current method works effectively for any USV type without requiring parameter adjustments, it still relies on predefined thresholds for harmonic detection. By automatically learning these parameters directly from the data, future approaches will become fully automated.
REFERENCES

[1] Y. K. Matsumoto and K. Okanoya, "Phase-specific vocalizations of male mice at the initial encounter during the courtship sequence," *PloS one*, vol. 11, no. 2, p. e0147102, 2016.

[2] G. Ehret, "Characteristics of vocalization in adult mice," in *Handbook of behavioral neuroscience*. Elsevier, 2018, vol. 25, pp. 187-195.

[3] D. Nicolakis, M. A. Marconi, S. M. Zala, and D. J. Penn, "Ultrasonic vocalizations in house mice depend upon genetic relatedness of mating partners and correlate with subsequent reproductive success," *Frontiers in zoology*, vol. 17, pp. 1-19, 2020.

[4] J. Chabout, A. Sarkar, D. B. Dunson, and E. D. Jarvis, "Male mice song syntax depends on social contexts and influences female preferences," *Frontiers in behavioral neuroscience*, vol. 9, p. 76, 2015.

[5] M. A. Marconi, D. Nicolakis, R. Abbasi, D. J. Penn, and S. M. Zala, "Ultrasonic courtship vocalizations of male house mice contain distinct individual signatures," *Animal Behaviour*, vol. 169, pp. 169-197, 2020.

[6] K. R. Coffey, R. E. Marx, and J. F. Neumaier, "Deepsqueak: a deep learning-based system for detection and analysis of ultrasonic vocalizations," *Neuropsychopharmacology*, vol. 44, no. 5, pp. 859-868, 2019.

[7] R. Abbasi, P. Balazs, M. A. Marconi, D. Nicolakis, S. M. Zala, and D. J. Penn, "Capturing the songs of mice with an improved detection and classification method for ultrasonic vocalizations (bootsnap)," *PLoS Computational Biology*, vol. 18, no. 5, p. e1010049, 2022.

[8] J. L. Hanson and L. M. Hurley, "Female presence and estrous state influence mouse ultrasonic courtship vocalizations," *PloS one*, vol. 7, no. 7, p. e40782, 2012.

[9] T. Klaus, B. Wernisch, S. M. Zala, and D. J. Penn, "Courtship vocalizations of wild house mice show highly dynamic changes and correlate with male copulatory success," *Animal Behaviour*, p. 123024, 2024.

[10] X. Dou, S. Shirahata, and H. Sugimoto, "Functional clustering of mouse ultrasonic vocalization data," *PloS one*, vol. 13, no. 5, p. e0196834, 2018.

[11] F. Auger and P. Flandrin, "Improving the readability of time-frequency and time-scale representations by the method of re-assignment," *IEEE Trans. Signal Process.*, vol. 43, no. 5, pp. 1068-1089, May 1995.

[12] N. Holighaus, Z. Průša, and P. L. Søndergaard, "Reassignment and synchrosqueezing for general time-frequency filter banks, subsampling and processing," *Signal Processing*, vol. 125, pp. 1-8, 2016.

[13] J. Allen, "Short term spectral analysis, synthesis, and modification by discrete fourier transform," *IEEE transactions on acoustics, speech, and signal processing*, vol. 25, no. 3, pp. 235-238, 1977.

[14] J. C. Brown, "Calculation of a constant q spectral transform," *The Journal of the Acoustical Society of America*, vol. 89, no. 1, pp. 425-434, 1991.

[15] I. Daubechies, "Ten lectures on wavelets," *Society for Industrial and Applied Mathematics*, 1992.

[16] K. Gröchenig, "Uncertainty principles for time-frequency representations," *Advances in Gabor analysis*, pp. 11-30, 2003.

[17] J. Xiao and P. Flandrin, "Multitaper time-frequency reassignment for nonstationary spectrum estimation and chirp enhancement," *IEEE Transactions on Signal Processing*, vol. 55, no. 6, pp. 2851-2860, 2007.

[18] Z. Prša, P. L. Søndergaard, N. Holighaus, C. Wiesmeyr, and P. Balazs, "The large time-frequency analysis toolbox 2.0," in *Sound, Music, and Motion: 10th International Symposium, CMMR 2013, Marseille, France, October 15-18, 2013. Revised Selected Papers 10*. Springer, 2014, pp. 419-442.

[19] P. L. Søndergaard, B. Torrésani, and P. Balazs, "The linear time frequency analysis toolbox," *International Journal of Wavelets, Multiresolution and Information Processing*, vol. 10, no. 04, p. 1250032, 2012.

[20] N. Holighaus, M. Dörfler, G. A. Velasco, and T. Grill, "A framework for invertible, real-time constant-Q transforms," *IEEE Trans. Audio Speech Lang. Process.*, vol. 21, no. 4, pp. 775-785, 2013.

[21] S. M. Zala, D. Reitschmidt, A. Noll, P. Balazs, and D. Penn, "Automatic mouse ultrasound detector (A-MUD): A new tool for processing rodent vocalizations," *PLOS ONE*, 2017.

[22] S. P. Ghael, A. M. Sayeed, and R. G. Baraniuk, "Improved wavelet denoising via empirical wiener filtering," in *Wavelet Applications in Signal and Image Processing V*, vol. 3169. SPIE, 1997, pp. 389-399.

[23] N. Wiener, *Extrapolation, interpolation, and smoothing of stationary time series: with engineering applications*. The MIT press, 1949.

[24] L. Rabiner, "On the use of autocorrelation analysis for pitch detection," *IEEE transactions on acoustics, speech, and signal processing*, vol. 25, no. 1, pp. 24-33, 1977.