Nonnegative Matrix Factorization with Transform Learning

TL;DR

This paper introduces a method to jointly learn an optimal short-time orthogonal transform along with nonnegative matrix factorization for improved audio signal decomposition, moving beyond fixed transforms.

Contribution

It proposes a novel joint optimization framework that learns both the transform and the factorization, enhancing traditional NMF techniques.

Findings

Improved audio signal decomposition results.

Effective joint transform learning and factorization.

Practical benefits demonstrated in experiments.

Abstract

Traditional NMF-based signal decomposition relies on the factorization of spectral data, which is typically computed by means of short-time frequency transform. In this paper we propose to relax the choice of a pre-fixed transform and learn a short-time orthogonal transform together with the factorization. To this end, we formulate a regularized optimization problem reminiscent of conventional NMF, yet with the transform as additional unknown parameters, and design a novel block-descent algorithm enabling to find stationary points of this objective function. The proposed joint transform learning and factorization approach is tested for two audio signal processing experiments, illustrating its conceptual and practical benefits.

Figures

Tables

Table 1. Table 1 : Source separation performance.

Method	SDR (dB)		SIR (dB)		SAR (dB)
SNR = -10 dB	${\hat{y}}_{\text{sp}}$	${\hat{y}}_{\text{no}}$	${\hat{y}}_{\text{sp}}$	${\hat{y}}_{\text{no}}$	${\hat{y}}_{\text{sp}}$	${\hat{y}}_{\text{no}}$
Baseline	-9.50	10.00	-9.50	10.00	$\infty$	$\infty$
IS-NMF	-6.75	6.82	-5.00	13.95	4.12	7.93
TL-NMF	1.73	12.29	13.44	13.33	2.22	19.20
SNR = 0 dB	${\hat{y}}_{\text{sp}}$	${\hat{y}}_{\text{no}}$	${\hat{y}}_{\text{sp}}$	${\hat{y}}_{\text{no}}$	${\hat{y}}_{\text{sp}}$	${\hat{y}}_{\text{no}}$
Baseline	0.10	0.08	0.10	0.08	$\infty$	$\infty$
IS-NMF	1.73	0.69	3.06	5.32	9.30	3.65
TL-NMF	6.50	5.81	12.11	9.16	8.16	9.00

Equations

$$\mathbf{V}\approx \mathbf{W}\mathbf{H},$$

$$\underset{\varphi ,\mathbf{W},\mathbf{H}}{\min }D(|\varphi (y){|}^{\circ 2}|\mathbf{W}\mathbf{H})$$

$$\underset{\mathbf{W},\mathbf{H}}{\min }D(|{\mathbf{\Phi }}_{\text{DCT}}\mathbf{Y}{|}^{\circ 2}|\mathbf{W}\mathbf{H})+\lambda \frac{M}{K}||\mathbf{H}|{|}_1$$

$$\text{s.t.}\mathbf{W}\ge 0,\mathbf{H}\ge 0,\forall k,||{\mathbf{w}}_k|{|}_1=1$$

$$C_{\lambda }(\mathbf{\Phi },\mathbf{W},\mathbf{H})\stackrel{\text{mbox}def}{=}D(|\mathbf{\Phi }\mathbf{Y}{|}^{\circ 2}|\mathbf{W}\mathbf{H})+\lambda \frac{M}{K}||\mathbf{H}|{|}_1$$

$$\underset{\mathbf{\Phi },\mathbf{W},\mathbf{H}}{\min }{C}_{\lambda }(\mathbf{\Phi },\mathbf{W},\mathbf{H})$$

$$\text{s.t.}\mathbf{W}\ge 0,\mathbf{H}\ge 0,\forall k,||{\mathbf{w}}_k|{|}_1=1,{\mathbf{\Phi }}^T\mathbf{\Phi }={\mathbf{I}}_M.$$

$$\nabla \stackrel{\text{mbox}def}{=}{\nabla }_{\mathbf{\Phi }}{C}_{\lambda }\left(\mathbf{\Phi },\mathbf{W},\mathbf{H}\right)=2\left(\mathbf{\Delta }\circ \mathbf{X}\right){\mathbf{Y}}^T$$

$${ϵ}^{(i)}=\frac{C_{\lambda }({\mathbf{\Phi }}^{(i)},{\mathbf{W}}^{(i)},{\mathbf{H}}^{(i)})-C_{\lambda }({\mathbf{\Phi }}^{(i-1)},{\mathbf{W}}^{(i-1)},{\mathbf{H}}^{(i-1)})}{C_{\lambda }({\mathbf{\Phi }}^{(i-1)},{\mathbf{W}}^{(i-1)},{\mathbf{H}}^{(i-1)})}$$

$$\mathbf{V}\approx {\mathbf{W}}_{\text{sp}}{\mathbf{H}}_{\text{sp}}+{\mathbf{W}}_{\text{no}}{\mathbf{H}}_{\text{no}},$$

$$C_{\Lambda }(\mathbf{\Phi },{\mathbf{H}}_{\text{sp}},{\mathbf{H}}_{\text{no}})$$

$$+{\lambda }_{\text{sp}}\frac{M}{N_{\text{sp}}}||{\mathbf{H}}_{\text{sp}}|{|}_1+{\lambda }_{\text{no}}\frac{M}{N_{\text{no}}}||{\mathbf{H}}_{\text{no}}|{|}_1$$

$$\text{s.t.}{\mathbf{\Phi }}^T\mathbf{\Phi }=\mathbf{I},{\mathbf{H}}_{\text{sp}}\ge 0,{\mathbf{H}}_{\text{no}}\ge 0,$$

$${\hat{\mathbf{Y}}}_s={\mathbf{\Phi }}^T\left(\frac{{\hat{\mathbf{V}}}_{\text{sp}}}{{\hat{\mathbf{V}}}_{\text{sp}}+{\hat{\mathbf{V}}}_{\text{no}}}\circ (\mathbf{\Phi }\mathbf{Y})\right)$$

$${\nabla }_{\mathbf{\Phi }}{C}_{\Lambda }\left(\mathbf{\Phi },\mathbf{H}\right)=2\left(\mathbf{\Delta }\circ \mathbf{X}\right){\mathbf{Y}}^T+2\left(\mathbf{\Xi }\circ {\mathbf{X}}_{\text{tr}}\right){\mathbf{Y}}_{\text{tr}}^T$$

Full text

Index Terms—  Nonnegative matrix factorization (NMF), transform learning, single-channel source separation

1 Introduction

Nonnegative matrix factorization (NMF) has become a privileged approach to spectral decomposition in several fields such as remote sensing and audio signal processing. In the latter field, it has led to state-of-the-art results in source separation [1] or music transcription [2]. The nonnegative data ${\mathbf{V}}\in\mathbb{R}_{+}^{M\times N}$ is typically the spectrogram $|{\mathbf{X}}|$ or $|{\mathbf{X}}|^{\circ 2}$ of some temporal signal $y\in\mathbb{R}^{T}$, where ${\mathbf{X}}$ is a short-time frequency transform of $y$, $|\cdot|$ denotes the entry-wise absolute value and $\circ$ here denotes entry-wise exponentiation. NMF produces the approximate factorization

[TABLE]

where ${\mathbf{W}}\in\mathbb{R}_{+}^{M\times K}$ is a nonnegative matrix referred to as dictionary that contains spectral patterns characteristic of the data while ${\mathbf{H}}\in\mathbb{R}_{+}^{K\times N}$ is the nonnegative matrix that contains the activation coefficients that approximate the data samples onto the dictionary. The factorization is usually low-rank ($K<\min(M,N)$) but not necessarily so (in which case regularization constraints should apply on ${\mathbf{W}}$ and/or ${\mathbf{H}}$). The decomposition (1) can then be inverted back to the time domain or post-processed in various ways to solve a large panel of audio signal processing problems.

In this traditional setting, a short-time frequency transform such as the short-time Fourier, Cosine or constant-$Q$ transforms acts as a pre-processing of the raw temporal data $y$. This is a potential limitation as any ill-chosen specification of the time-frequency transform may harm the quality of the decomposition. As such, we here propose to learn the transform together with the latent factors ${\mathbf{W}}$ and ${\mathbf{H}}$. We propose to address this task by solving an optimization problem of the form

[TABLE]

subject to structure constraints on $\phi:\mathbb{R}^{T}\rightarrow\mathbb{R}^{M\times N}$ and to nonnegativity of ${\mathbf{W}}$ and ${\mathbf{H}}$, and where $D(\,\cdot\,|\,\cdot\,)$ is a measure of fit. In addition, we study and promote the use of sparsity-inducing penalty on $\mathbf{H}$. We refer to objectives of the form (2) as TL-NMF, which stands for transform-learning NMF.

Connections to other works. TL-NMF is inspired by the work of Ravishankar & Bresler [5] on learning sparsifying transforms. Given a collection of data samples ${\mathbf{Y}}$ (such as a set of images), their work consists in finding an invertible transform $\boldsymbol{\Phi}$ such that the output of $\boldsymbol{\Phi}{\mathbf{Y}}$ is sparse. We are instead looking for a transform $\phi$ such that $|\phi(y)|^{\circ 2}$ can be well approximated by a NMF. TL-NMF can be viewed as finding a one-layer factorizing network, where $y$ acts as the raw data, $\phi$ the linear operator, $|\cdot|^{\circ 2}$ the nonlinearity and ${\mathbf{W}}{\mathbf{H}}$ the output of the network. Recent work has proposed combining deep learning and NMF, but in a different way. For instance, [6] considers a discriminative NMF setting and [7] studies nonnegative auto-encoders. The TL-NMF framework proposed in this work could in addition be extended to fully bridge deep learning and NMF by looking for a cascade of decompositions $f_{L}(\phi_{L}\ldots f_{1}(\phi_{1}(y)))$ such that the output is a NMF. Yet, this is beyond the scope of this paper and left for future work. Moreover, TL-NMF still operates in a transformed domain and is not directly related to synthesis-based NMF models in which the raw data $y(t)$ is modeled as $y(t)=\sum_{k}c_{k}(t)$ where the spectrogram of $c_{k}(t)$ is penalized so as to be closely rank-one [8, 9].

Goals and contributions. The goal of this work is to study the TL-NMF problem of form (2). As a first step, we propose in this paper to gently depart from the traditional short-time Fourier or Cosine transform setting by restricting the transform $\phi(y)$ to be a short-time orthogonal transform (Section 2.1). We consider real-valued transforms for simplicity and use the short-time Cosine transform (STCT) as a baseline. We propose an operational algorithm that returns stationary points of (2) and enables to learn a structured transform $\phi$ together with the NMF factors $\mathbf{W}$ and $\mathbf{H}$ (Section 2.2). The TL-NMF approach is put to test and compared to conventional STCT-based NMF in two benchmark audio signal processing experiments: music decomposition (Section 3) and speech enhancement (Section 4). The results demonstrate that the proposed approach is operational and yields, in both applications, significant benefits (reduced objective function, data-adapted atoms, improved separation accuracy).

2 NMF meets transform learning

2.1 Learning a short-time orthogonal transform

Let us denote by ${\mathbf{Y}}\in\mathbb{R}^{M\times N}$ the matrix that contains adjacent and overlapping short-time frames of size $M$ of $y$ and denote by $\boldsymbol{\Phi}_{\text{DCT}}\in\mathbb{R}^{M\times M}$ the orthogonal real-valued DCT matrix with coefficients $[\boldsymbol{\Phi}_{\text{DCT}}]_{qm}=(2M)^{-1/2}\cos\left(\pi(q+1/2)(m+1/2)/M\right)$. With these notations, the power spectrogram $|\mathbf{X}|^{\circ 2}$ of $y$ is simply given by $|\boldsymbol{\Phi}_{\text{DCT}}{\mathbf{Y}}|^{\circ 2}$. Traditional NMF (with sparsity) may be cast as

[TABLE]

where the notation ${\mathbf{A}}\geq 0$ expresses nonnegativity, $\|{\mathbf{A}}\|_{1}=\sum_{ij}|a_{ij}|$ and ${\mathbf{w}}_{k}$ denotes the $k^{th}$ column of ${\mathbf{W}}$. Regular unpenalized NMF is simply obtained for $\lambda=0$. The normalizing ratio $M/K$ ensures the measure of fit and the penalty term are of the same order of magnitude. We propose in this work to relax the pre-fixed transform $\boldsymbol{\Phi}_{\text{DCT}}$ and learn it jointly with ${\mathbf{W}}$ and ${\mathbf{H}}$. In other words, we consider the general case where $\mathbf{\Phi}\in\mathbb{R}^{M\times M}$ is a parameter of the objective function defined by

[TABLE]

and seek a solution to the TL-NMF problem defined by

[TABLE]

We choose at this stage to impose $\boldsymbol{\Phi}$ to be orthogonal though one could consider relaxing this assumption as well. The orthogonality constraint implicitly keeps $\boldsymbol{\Phi}$ nonsingular and excludes trivial solutions such as $(\boldsymbol{\Phi},{\mathbf{W}},{\mathbf{H}})=({\mathbf{0}},{\mathbf{0}},{\mathbf{0}})$ or $({\mathbf{1}}_{M\times M},{\mathbf{1}}_{M\times 1},|{\mathbf{1}}_{1\times M}{\mathbf{Y}}|^{\circ 2})$, where ${\mathbf{1}}_{M\times N}$ denotes the $M\times N$ matrix filled with ones. In this paper, we also choose the measure of fit $D(\,\cdot\,|\,\cdot\,)$ to be the Itakura-Saito (IS) divergence $D_{\text{IS}}({\mathbf{A}}|{\mathbf{B}})=\sum_{ij}(a_{ij}/b_{ij}-\log(a_{ij}/b_{ij})-1)$. Used with power spectral data, it is known to underlie a variance-structured Gaussian composite model that is relevant to the representation of audio signals [3] and has proven an efficient choice for audio source separation, e.g., [4]. However, the proposed framework can accommodate any other measure of fit.

2.2 Algorithm

We describe a block-coordinate descent algorithm that returns stationary points of problem (2.1). The blocks are the individual variables ${\mathbf{W}}$, ${\mathbf{H}}$ and $\boldsymbol{\Phi}$ that are updated in turn until a convergence criterion is met. As usual of NMF problems, the objective function $C_{\lambda}\left(\mathbf{\Phi},\mathbf{W},\mathbf{H}\right)$ is nonconvex and the returned solution depends on the initialization. We use for ${\mathbf{W}}$ and ${\mathbf{H}}$ the standard multiplicative updates that can be derived from a majorization-minimization procedure [10]. The sum-to-one constraint on the columns of ${\mathbf{W}}$ can be rigorously enforced using a change of variable, like in [11, 12].

Let us now turn our attention towards the update of $\boldsymbol{\Phi}$. We propose to use a gradient-descent procedure with a line-search step selection followed by a projection onto the orthogonal constraint, following the approach of [13]. The main benefit of this approach is that it yields an efficient yet simple algorithm for finding a orthogonal update for $\boldsymbol{\Phi}$. The gradient of the objective function with respect to (w.r.t.) $\mathbf{\Phi}$ can be shown to be given by

[TABLE]

where ${\mathbf{X}}=\boldsymbol{\Phi}{\mathbf{Y}}$, $\mathbf{\Delta}=\hat{\mathbf{V}}^{\circ-1}-\mathbf{V}^{\circ-1}$, ${\mathbf{V}}=|{\mathbf{X}}|^{\circ 2},\hat{{\mathbf{V}}}={\mathbf{W}}{\mathbf{H}}$. The steepest manifold-dependent descent direction is given by the natural gradient $\mathbf{\Omega}=\mathbf{\Phi}\mathbf{\nabla}^{T}\mathbf{\Phi}-\mathbf{\nabla}$. A suitable step-size $\gamma$ is then chosen according to the Armijo rule so that the projection $\pi\left(\mathbf{\Phi}+\gamma\mathbf{\Omega}\right)$ of the updated transform onto the orthogonal constraint induces a significant decrease of the objective function [13]. Our block-coordinate descent algorithm is stopped when the relative variation

[TABLE]

between iteration $i-1$ and $i$ falls below a given threshold $\tau$. The resulting TL-NMF algorithm is summarized in Algorithm 1. The sign ambiguity on $\boldsymbol{\Phi}$ resulting from the squared magnitude operation in $|\mathbf{\Phi}\mathbf{Y}|^{\circ 2}$ is removed by imposing that the entries of the first column of $\boldsymbol{\Phi}$ are positive, likewise the DCT. In the following experiments, we used nonnegative random values for initializing $\mathbf{W}$ and $\mathbf{H}$. Besides, we chose to initialize the transform $\mathbf{\Phi}$ in TL-NMF with a random orthogonal matrix.

3 Experiment 1: Music decomposition

In this section, we report results obtained with the proposed algorithm for decomposing real audio data $y(t)$, consisting of a $23$ s excerpt of Mamavatu by Susheela Raman that has been downsampled to $f_{s}=16$ kHz. ${\mathbf{Y}}$ is constructed using $40$ ms-long, 50%-overlapping temporal segments that are windowed with a sine bell. This construction leads to $M=640$ and $N=1191$. The proposed TL-NMF is compared to conventional IS-NMF, which is obtained using Alg. 1 as well, but fixed transform $\boldsymbol{\Phi}=\boldsymbol{\Phi}_{\text{DCT}}$. The two algorithms are run with the same stopping threshold $\tau=10^{-7}$, the same rank $K=10$, and three values of the regularization parameter $\lambda\in\{0,10^{3},10^{6}\}$. The two algorithms are initialized with the same random initializations of $\mathbf{W}$ and $\mathbf{H}$.

**Comparison of objective function values. ** Fig. 1 displays the objective function values w.r.t. iterations for the two approaches for $\lambda=10^{6}$ (the initial objective value is cropped for better readability) and permits the following conclusions. IS-NMF reaches the stopping criterion after fewer iterations than TL-NMF, as expected since the transform $\mathbf{\Phi}$ is fixed for the latter, leading to faster convergence for the remaining variables $\mathbf{W}$ and $\mathbf{H}$. Moreover, TL-NMF starts at a larger objective function value due to the random orthogonal transform initialization for $\mathbf{\Phi}$. Yet, the proposed TL-NMF algorithm monotonously decreases the objective function to a value that, at convergence, is of the order of the final value of IS-NMF.

**Regularization and learnt transform. ** We now examine examples of the atoms returned by TL-NMF (rows ${{\boldsymbol{\phi}}}_{m}$ of $\boldsymbol{\Phi}$) for the three values of $\lambda$. Fig. 2 displays the six atoms which most contribute to the audio signal (i.e., with largest values $||{{\boldsymbol{\phi}}}_{m}{\mathbf{Y}}||_{2}$). Clearly, without regularization ($\lambda=0$), the learnt atoms lack apparent structure. Yet, interestingly, as the value for $\lambda$ is increased, the atoms become oscillatory and smoother. This is a direct consequence of and justifies the use of the sparsity-inducing term in (2.1), which induces a structure-seeking constraint on the transform. Eventually, for strong regularization ($\lambda=10^{6}$), the learnt atoms resemble packet-like, highly regular oscillations, further analyzed in the next paragraph.

**Properties of the learnt atoms. ** In Fig. 3, pairs of atoms $(1,2)$, $(3,6)$, $(4,5)$ (from left to right) for $\lambda=10^{6}$ are plotted, together with the square root of their summed square magnitudes and the sine bell used to window the data. Interestingly, the results demonstrate that the identified pairs of atoms are approximately in quadrature. They hence offer shift invariance properties similar to those of the (complex-valued) short time Fourier transform, which is ubiquitous in audio applications. Yet, the TL-NMF algorithm permits to learn these atoms from random initializations. This provides strong evidence for the relevance of the proposed approach. Further, note that the atoms embrace the sine bell used to taper the observations $\mathbf{Y}$.

4 Experiment 2: Supervised source separation

We now examine whether learning an adaptive transform is actually useful for source separation. To this end, we consider a supervised NMF-based separation setting that follows the approach of [14]. In the following we address the separation of speech from interfering noise, but the method can be applied to any class of sound.

4.1 Principle

We assume that we are given speech and noise training data $y_{\text{sp}}(t)$ and $y_{\text{no}}(t)$ from which we form short-time matrices ${\mathbf{Y}}_{\text{sp}}$ and ${\mathbf{Y}}_{\text{no}}$ of sizes $M\times N_{\text{sp}}$ and $M\times N_{\text{no}}$, as in Section 2.1. Given a noisy speech recording $y(t)$ with short-time matrix ${\mathbf{Y}}$, traditional supervised NMF amounts to estimating activation matrices ${\mathbf{H}}_{\text{sp}}$ and ${\mathbf{H}}_{\text{no}}$ such that

[TABLE]

subject to sparsity of ${\mathbf{H}}_{\text{sp}}$ and ${\mathbf{H}}_{\text{no}}$, where ${\mathbf{V}}=|\boldsymbol{\Phi}_{\text{DCT}}{\mathbf{Y}}|^{\circ 2}$, ${\mathbf{W}}_{\text{sp}}=|\boldsymbol{\Phi}_{\text{DCT}}{\mathbf{Y}}_{\text{sp}}|^{\circ 2}$, ${\mathbf{W}}_{\text{no}}=|\boldsymbol{\Phi}_{\text{DCT}}{\mathbf{Y}}_{\text{no}}|^{\circ 2}$ [14]. Temporal source and noise estimates are then reconstructed in a second step by so-called Wiener filtering [3], based on the spectrogram estimates $\hat{{\mathbf{V}}}_{\text{sp}}={\mathbf{W}}_{\text{sp}}{\mathbf{H}}_{\text{sp}}$ and $\hat{{\mathbf{V}}}_{\text{no}}={\mathbf{W}}_{\text{no}}{\mathbf{H}}_{\text{no}}$.

In this section, we generalize this procedure by again learning an optimal transform within the separation procedure. To this end, we propose to build an approximation like (8) but where the fixed transform $\boldsymbol{\Phi}=\boldsymbol{\Phi}_{\text{DCT}}$ is now relaxed and learnt together with ${\mathbf{H}}_{\text{sp}}$ and ${\mathbf{H}}_{\text{no}}$. This means we propose to minimize

[TABLE]

with $\Lambda=\left(\lambda_{\text{sp}},\lambda_{\text{no}}\right)$ defining the possibly different weights in front of the sparsity terms. Note that $\boldsymbol{\Phi}$ now appears in both sides of the data-fitting term $D_{\text{IS}}(\cdot|\cdot)$, as the same transform is applied to the mixed data ${\mathbf{Y}}$ and the training data ${\mathbf{Y}}_{\text{sp}}$ and ${\mathbf{Y}}_{\text{no}}$. This requires to slightly modify the gradient of $C_{\Lambda}$ w.r.t. $\boldsymbol{\Phi}$ as compared to Section 2 and as described in next section. Given a solution to (9) and ${\mathbf{V}}=|\boldsymbol{\Phi}{\mathbf{Y}}|^{\circ 2}$ along with speech and noise spectrogram estimates $\hat{{\mathbf{V}}}_{\text{sp}}=|\boldsymbol{\Phi}{\mathbf{Y}}_{\text{sp}}|^{\circ 2}{\mathbf{H}}_{\text{sp}}$ and $\hat{{\mathbf{V}}}_{\text{no}}=|\boldsymbol{\Phi}{\mathbf{Y}}_{\text{no}}|^{\circ 2}{\mathbf{H}}_{\text{no}}$, temporal estimates may still be produced with Wiener filtering, i.e.,

[TABLE]

followed by standard overlap-adding of the columns of $\hat{{\mathbf{Y}}}_{\text{sp}}$ to return $\hat{y}_{\text{sp}}(t)$, and likewise for the noise. This is exactly the same procedure than in traditional NMF-based separation except that $\boldsymbol{\Phi}_{\text{DCT}}$ and $\boldsymbol{\Phi}_{\text{DCT}}^{T}$ are replaced by $\boldsymbol{\Phi}$ and $\boldsymbol{\Phi}^{T}$.

4.2 Algorithm

Denote ${\mathbf{Y}}_{\text{tr}}=\left[{\mathbf{Y}}_{\text{sp}},{\mathbf{Y}}_{\text{no}}\right]$, ${\mathbf{X}}_{\text{tr}}=\boldsymbol{\Phi}{\mathbf{Y}}_{\text{tr}}$, ${\mathbf{W}}=|{\mathbf{X}}_{\text{tr}}|^{\circ 2}$, ${\mathbf{H}}=\left[{\mathbf{H}}_{\text{sp}}^{T},{\mathbf{H}}_{\text{no}}^{T}\right]^{T}$ and $\hat{{\mathbf{V}}}={\mathbf{W}}{\mathbf{H}}$. Given ${\mathbf{W}}$, ${\mathbf{H}}$ can be updated with multiplicative rules derived from majorization-minimization as in [10]. We use again a gradient-descent approach for the update of $\boldsymbol{\Phi}$. The gradient of the objective function (9) can be expressed as

[TABLE]

where $\mathbf{\Delta}=\hat{\mathbf{V}}^{\circ-1}-\mathbf{V}^{\circ-1}$ and $\mathbf{\Xi}=\mathbf{\Delta}^{\prime}{\mathbf{H}}^{T}$ with $\mathbf{\Delta}^{\prime}=\frac{\hat{\mathbf{V}}-\mathbf{V}}{\hat{\mathbf{V}}^{\circ 2}}$.

Note that the first term of (11) is the gradient in (6). The second term is nothing but the gradient of the data-fitting term $D_{\text{IS}}$ with its first argument fixed. Based on (11), we again use a line-search step selection in the steepest natural gradient direction followed by a projection, like in Section 2.2 and following [13]. The resulting algorithm is summarized in Alg. 2.

4.3 Speech enhancement experiment

We consider clean speech and noise data from the TIMIT corpus [15] and the CHIME challenge,111http://spandh.dcs.shef.ac.uk/chime_challenge respectively. For speech training data $y_{\text{sp}}(t)$, we use all utterances but the first one in the train/dr1/fcjf0 directory (about $21$ s in total). For noise training data $y_{\text{no}}(t)$, we use $30$ s of the file BGD_150204_010_BUS.CH1.wav, which contains noise recorded in a bus. A simulated mixed signal $y(t)$ of duration 3 s is generated by mixing the remaining speech utterance with another segment of the noise file (as such, the test data is not included in the training data), using signal-to-noise (SNR) ratios of -10 dB and 0 dB. The audio files’ sampling frequency is $f_{s}=16$ kHz, and short-term matrices ${\mathbf{Y}}$, ${\mathbf{Y}}_{\text{sp}}$ and ${\mathbf{Y}}_{\text{no}}$ are constructed using $40$ ms-long, 50%-overlapping windowed segments like in Section 3, leading to $M=640$, $N=149$, $N_{\text{sp}}=1059$ and $N_{\text{no}}=1517$.

Our supervised TL-NMF approach is compared to the traditional supervised NMF procedure (with the IS divergence) described in Section 4.1, based on the same training data and using the same regularization parameters (only the transform $\boldsymbol{\Phi}$ differs between the two approaches). Source separation performance was assessed using the standard BSS_eval criteria [16]. For simplicity we set $\lambda_{\text{sp}}=\lambda_{\text{no}}=\lambda$ and tested various order of magnitudes of $\lambda$. The performance of both IS-NMF and TL-NMF were rather robust to the choice of $\lambda$, and also offered best individual performance for the same range of values. We report results with the optimal values $\lambda=10^{-1}$ and $\lambda=10^{-4}$ for SNR = -10 dB and SNR = 0 dB, respectively. We also compute the performance criteria obtained by $\hat{y}_{\text{sp}}=\hat{y}_{\text{no}}=y/2$ as an indicative baseline. Table 1 reports the comparison results.

The results show that the extra adaptability offered by TL-NMF is clearly beneficial as far as source separation capabilities are concerned. Indeed, TL-NMF dramatically improves the signal to distortion and interference ratios for the speech source by $8.5$, $18.5$ (SNR = -10 dB) and $4.73$, $9.05$ (SNR = 0 dB) as compared to IS-NMF. However, the signal to artifact ratios are slightly better with IS-NMF, with improvements of $-1.9$ (SNR = -10 dB) and $-1.14$ (SNR = 0 dB). Audio samples are available online.222https://www.irit.fr/~Cedric.Fevotte/extras/audio_icassp2018.zip

5 Conclusion and future work

We addressed the task of learning the transform underlying NMF-based signal decomposition jointly with the factorization. Specifically, we have proposed a block-coordinate descent algorithm that enables us to find a orthogonal transform $\boldsymbol{\Phi}$ jointly with the dictionary ${\mathbf{W}}$ and the activation matrix ${\mathbf{H}}$. To our knowledge, the proposed algorithm is the first operational procedure for learning a transform in the context of NMF. Our preliminary experiments with real audio data indicate that automatically adapting the transform to the signal pays off when seeking latent factors that accurately represent the data. In particular, we were able to retrieve a naturally structured transform which embeds some form of invariance in a first musical example. Then, we achieved source separation performance that compares very favorably against the state-of-the-art in a speech enhancement experiment. Note that although our presentation focused on the processing of audio data, the approach can be adapted to many other settings where NMF is applied to preprocessed data.

Future work will include the estimation of complex-valued transforms which can be directly compared to the STFT, the influence of the initialization of $\boldsymbol{\Phi}$, the influence of the value of $K$ on the learnt transform as well as relaxations of the orthogonal constraint on $\boldsymbol{\Phi}$ to mere nonsingular matrices. Also, the use of alternative optimization strategies that lend themselves well to dealing with nonconvex problems in high dimension, including stochastic gradient descent, will be investigated.