Prior Distribution Design forMusic Bleeding-Sound Reduction Based on
Nonnegative Matrix Factorization
Yusaku Mizobuchi∗ , Daichi Kitamura∗ , Tomohiko Nakamura† , Hiroshi Saruwatari† ,Yu Takahashi‡ , and Kazunobu Kondo‡∗ National Institute of Technology, Kagawa College, Kagawa, Japan† The University of Tokyo, Tokyo, Japan‡ Yamaha Corporation, Shizuoka, Japan

Abstract—When we place microphones close to a sound sourcenear other sources in audio recording, the obtained audio signal
includes undesired sound from the other sources, which is often
called cross-talk or bleeding sound. For many audio applications
including onstage sound reinforcement and sound editing after alive performance, it is important to reduce the bleeding sound ineach recorded signal. However, since microphones are spatiallyapart from each other in this situation, typical phase-aware
blind source separation (BSS) methods cannot be used. We
propose a phase-insensitive method for blind bleeding-soundreduction. This method is based on time-channel nonnegative
matrix factorization, which is a BSS method using only amplitude
spectrograms. With the proposed method, we introduce the
gamma-distribution-based prior for leakage levels of bleeding
sounds. Its optimization can be interpreted as maximum aposteriori estimation. The experimental results of music bleedingsound reduction indicate that the proposed method is more
effective for bleeding-sound reduction of music signals comparedwith other BSS methods.

I. I NTRODUCTION
When we record a live musical performance, many microphones are usually arranged among the players. Some arelocated very close to each of the audio sources, such as musicalinstruments, vocals, and amplifiers. These close microphones
are placed to pick up the sound from the source other than thatwhich is intended. However, undesirable audio leakage fromthe non-target audio sources is also captured, which is often
called “cross-talk” or “bleeding sound,” as shown in Fig. 1.
In onstage mixing, sound engineers control the balance ofsound levels of individual sources, and the processed sounds
are provided to the audience through loudspeakers and performers through monitor speakers. Bleeding sound makes suchsound reinforcement difficult, degrading musical performance
quality. It is also necessary to avoid sound bleeding for highquality audio editing (remixing) of the recorded signals after alive performance. For these reasons, sound engineers carefully
place close microphones so that the as much bleeding soundis reduced as possible. Putting acoustic barriers between thesound sources and reducing sound reflection in the recording
room are also effective. However, completely avoiding bleeding sound is almost impossible. In other words, sound bleeding
essentially occurs in a live-recording situation.

Fig. 1. Spatial arrangement of sources and close microphones, where M =N = 4. Target sound is contaminated with bleeding sound from other nontarget sources.

Bleeding-sound reduction is similar to the well-investigatedproblem called multichannel audio source separation(MASS) [1], [2], [3], [4], but some conditions are different
from those in MASS, which are listed as follows.
(a) The signal-to-noise ratio (SNR) of the observed signal
is relatively high because of a close miking setup, where
the “signal” is a target source for the close microphone
and the “noise” is the leakage from the other sources.
(b) The observed multichannel signals are already “labeled,”
namely, the target source for each microphone is known
because each microphone is located close to each soundsource.
(c) The microphones are spatially apart from each other(e.g., more than 2 m), resulting in serious spatial aliasing.
(d) The requirement of separation quality is relatively highso as not to degrade the artistic value of the music signal.
Conditions (a) and (b) are advantages of bleeding-sound reduction, which make resolving bleeding sound easier than MASS.
However, conditions (c) and (d) are difficult. In particular,
condition (c) is critical because typical high-quality MASS,
including beamformers [1], [2] and independence-based blindsource separation (BSS) [5], [6], [7], [8], [9], [10], [11], usesphase differences between microphones, which are unreliablein bleeding-sound reduction due to spatial aliasing. To tackle
this problem, phase-insensitive (amplitude- or power-based)
MASS [12], [13], [14], [15] can be applied. Togami et al. [12]applied nonnegative matrix factorization (NMF) [17], [18] tothe time-channel domain in each frequency (hereafter, timechannel NMF: TCNMF), where both the nonnegative mixing
matrix and amplitude activation of each source are estimatedin each frequency bin. TCNMF performs well even undercondition (c) or an asynchronous recording condition [13],[14], although its effectiveness regarding music bleedingsound reduction has not been investigated. A BSS-based
method that ignores the phase information was proposed [15],which is called linear demixed domain multichannel NMF(DMNMF). Similar to TCNMF, this method also estimates
the frequency-wise nonnegative mixing matrix. Das et al. [16]introduced supervised information to accurately reduce thebleeding sound, where the frequency-wise nonnegative mixing
matrix (i.e., leakage levels of the non-target sources for eachclose microphone) is measured before the musical performance or calculated using the solo-played time segments ofeach source. However, to reduce the onsite recording cost forsound reinforcement, such supervision should not be used.
Also, a mismatch between the obtained mixing matrix and actual condition may markedly degrade reduction performance.
We aimed to reduce bleeding sound in a fully blind manner,namely, the spatial locations of sources and microphones are
unknown. We also did not use supervision of sources, suchas solo-played music datasets, to avoid the mismatch between
training and test data; thus, supervised deep-neural-networkbased approaches [19], [20], [21], [22], [23] are out of thescope of this paper. We propose a phase-insensitive method
for blind bleeding-sound reduction, which is a modificationof TCNMF: we introduce an a priori generative model for
diagonal and off-diagonal elements of the frequency-wise
mixing matrix to model relative leakage levels of bleeding
sounds. This method is based on NMF with maximum aposteriori (MAP) estimation, which was originally proposed
by Cemgil [24], and we demonstrate that the proposed method
is suitable for reducing music bleeding sound.
II. C ONVENTIONAL M ETHODS
A. Mixture Model
Let M and N be the numbers of microphones (channels)
and sources, respectively. The source, observed, and estimatedsignals are respectively denoted asTNs̃(t) = [s̃1 (t), s̃2 (t), · · · , s̃n (t), · · · , s̃N (t)] ∈ R ,T(1)

Fig. 2. Instantaneous mixture model for bleeding-sound reduction, where
M = N = 4. Color brightness in mixing matrix Ã shows amplitude level ofeach element (brighter is larger). Due to close miking setup, diagonal elementsin Ã have larger amplitudes compared with off-diagonal elements.

on the determined case, which is the most difficult situationin bleeding-sound reduction.
In an instantaneous mixture, the observed and estimatedsignals can respectively be modeled asx̃(t) = Ãs̃(t),(4)ỹ(t) = W̃ x̃(t),(5)

where Ã ∈ RM ×N and W̃ ∈ RN ×M are the time-invariant
mixing and demixing matrices, respectively. The mixture
model (4) is illustrated in Fig. 2. Since the observed signalx̃(t) is “labeled,” as explained in condition (b) in Sect. I,
we define that x̃m (t) is the close-microphone signal forthe mth source s̃m (t) (n = m), as shown in Figs. 1 and2. Thus, x̃m (t) mainly contains the sound from the target
source s̃m (t), although the bleeding sound from the nontarget sources s̃m0 (t) is also included, where m0 6= m. For
this reason, the absolute values of diagonal elements in Ãshould be large enough, and those of off-diagonal elements
become small, which results in high-SNR condition (a) in
Sect. I. Blind bleeding-sound reduction is formulated as an
estimation problem of the mixing matrix Ã or demixing matrix
that satisfies W̃ = Ã−1 from only x̃(t).
In actual recording, the mixing system (4) becomes a
convolutive mixture due to time difference of arrival androom reverberation. To simply model the convolutive mixture,
we assume that the impulse responses (reverberation time)
between microphones and sources are shorter than the window
length used in the short-time Fourier transform (STFT). This
assumption enables us to respectively model the reverberant
observed and estimated signals as1
xij = Ai sij ,(6)

yij = Wi xij ,

(7)

where
sij = [sij1 , sij2 , · · · , sijn , · · · , sijN ]T ∈ CN ,(8)Mx̃(t) = [x̃1 (t), x̃2 (t), · · · , x̃m (t), · · · , x̃M (t)] ∈ R ,(2)

M

xij = [xij1 , xij2 , · · · , xijm , · · · , xijM ] ∈ C ,(9)ỹ(t) = [ỹ1 (t), ỹ2 (t), · · · , ỹn (t), · · · , ỹN (t)]T ∈ RN ,(3)

yij = [yij1 , yij2 , · · · , yijn , · · · , yijN ]T ∈ CN .

(10)

where t = 1, 2, · · · , T , n = 1, 2, · · · , N , and m =1, 2, · · · , M are the indices of discrete time, source, and microphone, respectively. Under the recording condition describedin Sect. I, the mixing system becomes determined (M = N )or overdetermined (M > N ). In this study, we focused onlyT

Here, i = 1, 2, · · · , I and j = 1, 2, · · · , J are the indices of thefrequency bin and time frame, respectively, and Ai ∈ CM ×N
1 Note that roman font signal denotes complex values and italic font signal
denotes real (or nonnegative) values.

is the complex-valued frequency-wise mixing matrix. Also,
sijn , xijm , and yijn are the complex-valued time-frequencywise elements of the source, observed, and estimated spectrograms Sn ∈ CI×J , Xm ∈ CI×J , and Yn ∈ CI×J ,respectively. In (6), the convolutive mixture is converted to
the frequency-wise instantaneous mixture via STFT.
Typical beamformers [1], [2] and BSS methods [5], [6],[7], [8], [9], [10], [11] are used to estimate the complex-valued
demixing matrix Wi on the basis of a principle of microphonearrays, e.g., time difference of arrival, and these methods
rely on the phase differences between microphones. Whenmicrophones are spatially apart from each other, these methods
cannot precisely estimate Wi because of spatial aliasing. This
problem is salient in bleeding-sound reduction.

Fig. 3. Decomposition model of TCNMF, where M = N = 4 and I = 6.
Note that abs(Xm ) is channel-wise time-frequency matrix, but Ai and Si
are frequency-wise source-channel and time-source matrices, respectively.

The variables Ai and Si can be estimated by solving thefollowing minimization problem [18]:Xmin
DKL (Xi |Ai Si ) s.t. aimn , sinj ≥ 0 ∀i, m, n, j,
A,S

B. DMNMFi(17)

To cope with spatial aliasing, the power-based BSS method
DMNMF was proposed [15]. DMNMF can be interpreted asa phase-insensitive version of independent low-rank matrix
analysis (ILRMA) [4], [10], [11], and the observed signal ismodeled as.2x.2
ij ≈ Ai sij ∀i, j,

(11)

×N
Ai = abs(Ai ) ∈ RM
,
≥0

(12)

xij = abs(xij ) ∈ RM
≥0 ,N
sij = abs(sij ) ∈ R≥0 ,

(13)
(14)

where the dotted exponent ·.q and absolute operation abs(·)
for vectors or matrices return the element-wise qth power and.2absolute, respectively; thus, x.2
ij and sij are the power specM
trogram components of {Xm }m=1 and {Sn }N
n=1 , respectively.
DMNMF approximates (6) by the nonnegative frequency-wise
mixing matrix Ai in the power-spectrogram domain to ignore
the phase information. In addition, the power spectrogram ofeach source is modeled by a low-rank matrix using NMF. After.2estimating Ai and s.2
ij from xij , we can recover the estimatedsignal yij by Wiener filtering.
C. TCNMF
The amplitude-based BSS method TCNMF was proposed [12] and applied [13], [14] to speech enhancement.
Whereas typical NMF is a low-rank decomposition oftime-frequency matrices, TCNMF decomposes frequency-wisetime-channel matrices in the amplitude domain as
Xi ≈ Ai Si ∀i,
Xi = [xi1 xi2 · · ·(15)×JxiJ ] ∈ RM
≥0 ,

(16)

which is illustrated in Fig. 3, where Xi is the frequency-wisetime-channel observed signal in the amplitude domain and×JSi ∈ RN
is a time-source activation matrix: Si involves≥0time-varying gains of each source as the row vectors. By
estimating Ai and Si in the same manner as typical NMF, we
can reconstruct the estimated sources using Wiener filtering.

where
DKL (Xi |Ai Si ) =

X

ximj
n aimn sinj!

ximj log Pm,j− ximj +Xaimn sinj(18)

nis the generalized Kullback-Leibler (KL) divergence that measures the similarity between Xi and Ai Si , A and S are
the sets {Ai }Ii=1 and {Si }Ii=1 , respectively, and ximj , aimn ,and sinj are the elements of Xi , Ai , and Si , respectively.
However, since Ai is an M × N square matrix in thedetermined case, the minimization problem (17) has a trivial
solution, namely, Ai = I for all i, where I is an identity
matrix. To avoid this trivial solution, an L0.5 -norm-based
sparse regularizer was introduced for each time frame [12]as follows:XX
DKL (Xi |Ai Si ) + µ
ksij k0.5
min
A,S

i

i,j

s.t. aimn , sinj ≥ 0 ∀i, m, n, j,(19)

where µ is a weight coefficient for regularization. Notethat sij is a time-frame-wise vector in Si , namely, Si =[si1 si2 · · · siJ ].
III. P ROPOSED M ETHOD
A. Motivation
In bleeding-sound reduction, phase information cannot be
used because of the close miking setup and serious spatial
aliasing. As a phase-insensitive method, DMNMF is a reasonable approach. However, full-blind parameter optimization of DMNMF is difficult and unstable. In fact, a priori
information of steering vectors (column vectors of Ai ) ora phase-aware BSS method is used for pre-estimation [15]to stabilize and improve BSS performance. TCNMF canestimate the source signals without phase information, evenin asynchronous recording [13]. However, its performancefor music BSS or bleeding-sound reduction has
P not been
investigated. In particular, the sparse regularizer i,j ksij k0.5

in (19) may degrade the sound quality of estimated signals
in music mixture. This is because the regularizer is based on
a W-disjoint-orthogonality assumption in the time-frequency
domain [25], which is suitable only for speech mixtures. Sincemusic mixtures frequently include both spectral and temporal
overlaps of sources, the sparse regularizer for Si may beinappropriate.
To avoid the trivial solution of Ai in TCNMF, we use
our proposed method to regularize both the diagonal andoff-diagonal elements of the nonnegative mixing matrix Ai
instead of regularizing Si . The proposed method can beinterpreted as a MAP estimation, where the bleeding-soundlevels are assumed to be generated by the gamma distribution
prior.
B. Generative Model of KL-Divergence-Based NMF
Cemgil [24] revealed the generative model of NMF with KL
divergence (KLNMF): the minimization problem in KLNMF
is equivalent to the maximum likelihood (ML) estimation
with the Poisson generative model. For (17), the following
generative model is assumed:
zimnj ∼ P(zimnj ; aimn sinj ),1e−λ λz ,P(z; λ) =Γ(z + 1)

(20)
(21)

where
Pzimnj ∈ N is a random variable that satisfies ximj =e + n zimnj , P(z; λ) is the Poisson distribution with therandom variable z ∈ N and parameter λ > 0, Γ(z + 1) =z! is the gamma function, and e is a random variable thatobeys uniform distribution in the range [0, 1). Also, zimnj
is assumed to be mutually independent w.r.t. i, m, n, and j.
The Poisson random variables have the superposition
property,P
namely, when zn ∼ P(zn ; λn ) and x P= n zn , the marginal
probability is given by p(x) = P(x; n λn ). Therefore, themarginal log-likelihood of Xi is given bylog p(Xi ; Ai , Si )
Y X
= log
p(xinm ; zimnj )p(zimnj ; aimn sinj )m,j zimnj

Y

= log

P(ximj ;P

n aimn sinj )m,j=ximj log

X

(23)
(24)

where δ(a) is the Dirac’s delta distribution and G(a; k, θ) is the
gamma distribution with the random variable a ≥ 0 and shapeand scale parameters k > 0 and θ > 0. Note that the gamma
distribution is a conjugate prior of the Poisson generativemodel (20). In addition, aimn is assumed to be mutuallyindependent w.r.t. i, m, and n; thus, the prior distribution of
Ai becomesYYp(Ai ; k, θ) =p(aimn )p(aimn ; k, θ)
m,n=m

=

Y

m,n6=mδ(aimn − 1)

m,n=m

Y

G(aimn ; k, θ).

m,n6=m(25)
The prior (25) enables us to control the probability of offdiagonal elements of Ai (relative leakage levels of bleeding
sound) by k and θ, while restricting all the diagonal elements
to be unity. As shown in Fig. 4, we can avoid aimn = 0 for
all m 6= n, which is the trivial solution of Ai , by setting theshape parameter to k > 1. Hereafter, we consider k > 1 only.
For the activation matrix Si , we do not assume explicit
structure, but only the nonnegativity prior is used as follows:
sinj ∼ lim

aimn sinjnm,j

the diagonal and off-diagonal elements of Ai :(δ(aimn − 1)
(m = n)aimn ∼,G(aimn ; k, θ) (m 6= n)1
G(a; k, θ) =ak−1 e−a/θ ,Γ(k)θk1
I[0 ≤ sinj ≤ β]β→∞ β
∝ I[0 ≤ sinj ],"X

Fig. 4. Probabilistic density function of gamma distribution.

(26)

#
−

X

aimn sinj − log Γ(ximj + 1) .

(22)

n

The maximization of (22) w.r.t. aimn and sinj for all i (MLestimation) is equivalent to the minimization of (18).
C. A Priori Generative Model for Bleeding-Sound Levels
With the proposed method, to avoid the trivial solution of
Ai , we introduce the following a priori generative model into

where β is the normalized coefficient and I[·] denotes a binary
distribution that has value one when its argument is true andzero otherwise. Similar to Ai , sinj is mutually independent
w.r.t. i, n, and j, and the prior distribution of Si becomesYp(Si ) =p(sinj )
n,j

∝

Y
n,j

I[0 ≤ sinj ].

(27)

D. Cost Function for MAP Estimation

designed using Jensen’s inequality as follows:

On the basis of the above-mentioned prior distributions, weestimate variables Ai and Si in the MAP sense. The posterior
distribution can be obtained asYYp(Ai , Si ; Xi ) ∝p(Xi ; Ai , Si ) p(Ai ; k, θ)p(Si ) .
{z
}|
{z
}
|

DKL (Xi |Ai Si )
!
X
X
X
c
=
−ximj log
aimn sinj +aimn sinjii

PriorsLikelihood

=

≤

n

aimn sinj X−ximj logξimnj+
aimn sinj
ξimnj
n
n

X

−ximj

aimn sinj X+
aimn sinj
ξimnj log
ξimnj
n
n

J =−

[log p(Xi ; Ai , Si ) + log p(Ai ; k, θ) + log p(Si )] .

i

(29)
From (22), (25), and (27), the cost function J is obtained as
"
X
X
J =
−ximj log
aimn sinj

i,m,j
+

#
+

X

aimn sinj + log Γ(ximj + 1)++

I[aimn = 1]

i,m,n6=m+

X

I[0 ≤ simn ] ,(30)

where I[·] = − log I[·] denotes an indicator function that hasvalue zero when its argument is true and ∞ otherwise. The
MAP estimation of Ai and Si is a minimization problem of(30), and this minimization w.r.t. Ai and Si is equivalent to
the following problem:
A,SX

DKL (Xi |Ai Si ) +ii,m,n6=m

and aimn = 1 ∀m = n.
1−(k − 1) log aimn + aimn
θ

i,n,j

min

i

s.t. aimn , sinj ≥ 0 ∀i, m, n, j, ξimnj > 0 ∀i, m, n, j,

i,m,n=mXcwhere = denotes equality up toPa constant, ξimnj > 0 is an
auxiliary variable that satisfies n ξimnj = 1, and Ξ is a setof ξimnj for all i, m, j, and n. The equality in (33) holds if
and only ifaimn sinj
ξimnj = P∀i, m, j, n.
(34)
00n0 aimn sin j

A,S,Ξ

n

X

(33)

From (33), the MM problem is obtained as
XX
min
D+ (Ai , Si , Ξ) +R(aimn ; k, θ)

n

i,m,j

X

!

X

≡ D (Ai , Si , Ξ),
X

!

X

X
i,m,j

(28)
By taking a negative logarithm of (28), we can decompose theright side of (28) as

n

i,m,j

R(aimn ; k, θ)

i,m,n6=ms.t. aimn , sinj ≥ 0 ∀i, m, n, j and aimn = 1 ∀m = n, (31)(35)

By setting the derivative of the majorization function (35)w.r.t. aimn and sinj to zero and substituting (34) for ξimnj ,we can derive the update rules. Since the regularizer does notaffect sinj , the update rule of sinj is the same as that of simple
KLNMF [18] and expressed as
P
ximj
P
aimn
m
n0 aimn0 sinj 0Psinj ← sinj.
(36)
m aimnFor the off-diagonal elements aimn (m 6= n), we have thefollowing equations from the derivative of the majorization
function:

X
1
1
ξimnj
−ximj
+ sinj − (k − 1)
+ = 0. (37)
a
a
θ
imn
imn
j
Therefore, we have

where


1
R(aimn ; k, θ) = −(k − 1) log aimn + aimn
θaimn =

(32)

is the regularizer that corresponds to the gamma distribution
prior (24) for the off-diagonal elements of Ai .
E. Derivation of Optimization Algorithm
The minimization problem (31) can be solved using amajorization-minimization (MM) algorithm [18], [26], which
is often used in the context of NMF optimization. The majorization function of the fidelity term DKL (Xi |Ai Si ) can be(k − 1) +

P

1
θ +

P

j ximj ξimnj
j sinj.

(38)

The update rule of the off-diagonal elements aimn is derived
by substituting (34) asPx(k − 1) + aimn j P 0 a imj0 s 0 sinj
imn in jnPaimn ←
.
(39)
1
+
s
inj
j
θ
The nonnegativity of aimn and sinj can hold by setting theirinitial values to nonnegative values. Since the value of the
diagonal elements of Ai is restricted, we initialize the diagonal
elements aimn (m = n) with unity and fix them during theiterative optimization of the other variables.

The efficient matrix-form implementation of (36) and (39)is as follows:XiT(k − 1) + Ai
Ai Si Si
∀i,
(40)
Ai ←
1
T
θ + 1Sidiag(Ai ) ← [1, 1, · · · , 1]T ∀i,
Si ← SiXi
AT
i Ai Si
AT
i 1

∀i,

(41)

A. Conditions

F. Balancing Between Fidelity Term and Regularizer
With the proposed method, the diagonal elements of Ai
are restricted to be unity so that the off-diagonal elements
correspond to the relative leakage levels of bleeding sound.
The KL divergence (18) also has a scale-dependent property,
namely,(43)

where α ≥ 0 is an arbitrary coefficient. These facts meanthat an observed gain of Xi , i.e., the signal amplitude ineach
P microphone, affects the balance
P of the fidelity termi DKL (Xi |Ai Si ) and regularizer
i,m,n6=m R(aimn ; k, θ)in (31).
To solve this problem, we also parameterize the observed
gain. The following normalization is carried out for the observed signal x̃(t) before we apply the proposed method:α(44)x̃(t) ← x̃(t) ∀t,v
Tv = max {abs(x̃(t))}t=1 ,(45)where max(·) returns the maximum scalar value of the input set. After the normalization (44), a dynamic range of{x̃(t)}Tt=1 becomes ±α. Similar to µ in (19), we can control
the balance between the fidelity term and regularizer by α. If
we set α to a small value, the regularizer strongly affects the
optimization.
G. Reconstruction of Estimated Signals
Similar to conventional TCNMF, the complex-valued estimated signal Yn can be recovered by applying Wiener filteringto the complex-valued observed signal xijm as follows:(aimm simj )2yijn = Px .
2 ijm
n (aimn sinj )Since aimm = 1, (46) can be implemented as
Si.2
yijn =
xijm ,
.2
A.2
i Si
m,j

IV. E XPERIMENTS(42)

where and the quotient symbol for matrices denote elementwise multiplication and division, respectively, 1 is an M × J
matrix containing only ones, and diag(·) returns a vector thatconsists of the diagonal elements of the input square matrix.
Note that (40) will change the value of the diagonal elementsof Ai , but they are immediately replaced with unity by (41).
It is guaranteed that the iterative calculation of (40)–(42)monotonically decreases the cost function (31).

DKL (αXi |αAi Si ) = αDKL (Xi |Ai Si ) ,where [·]m,j denotes an (m, j) element of the input matrix.
After Wiener filtering, the estimated signal Yn is converted
to the time-domain signal ỹn (t) via the inverse STFT. Then,the signal gain is recovered byvỹ(t) ← ỹ(t) ∀t.
(48)
α

(46)

(47)

To evaluate the performance of the proposed method(proposed TCNMF), we conducted an experiment of blind
bleeding-sound reduction. The observed music mixture signal
was simulated using songKitamura [27], [28], which is an
artificial music dataset. We chose four musical instruments,
clarinet (Cl.), oboe (Ob.), piano (Pf.), and trombone (Tb), as
dry sources Sn and prepared a four-channel observed signal
xij so that M = N = 4. To simulate bleeding sound, we
mixed these instrumental sounds sij using the frequency-wise×N
as follows:
nonnegative random mixing matrix Ai ∈ RM
≥0
xij = Ai sij ,(49)

where the diagonal and off-diagonal elements of Ai are
set to unity and uniformly distributed random values in therange (0, 0.2) for all i, respectively. This mixing system is
an approximation of (6). In this experiment, ten observed
mixtures were prepared using different pseudo-random seeds,
i.e., ten different mixing matrices Ai .
For all signals, we performed STFT using a 4096-point-longhamming window with half-overlap shifting, where a sampling
frequency of the signals was 44.1 kHz. The numbers of
frequency bins and time frames were I = 2049 and J = 109,
respectively. The update rules in the optimization algorithm
were iterated 200 times, and we confirmed the convergence ofthe cost function value.
For DMNMF and the conventional and proposed TCNMFs,
the initial value of Ai was set as follows: the diagonal andoff-diagonal elements were set to unity and the uniformly distributed random value in the range (0, 0.1), respectively. The
other parameters were initialized by the uniformly distributed
random value in the range (0, 1).
As an evaluation criterion, we used the source-to-distortion
ratio (SDR) [29], which indicates total separation quality
including both degree of separation (source-to-interference
ratio: SIR) and absence of artificial distortion (sources-toartificial ratio: SAR). As described in condition (a) in Sect. I,
the SNR and SDR of the observed signals for the bleedingsound reduction are high. In our experiment, the average SDRs
over the ten observed mixture signals of Cl., Ob., Pf., and Tb.
were 18.8, 15.0, 14.7, and 8.6 dB, respectively. We calculated
the improvements from these input SDRs for each source toevaluate the performance of each method.
We compared five methods, i.e., independent vector analysis(IVA) [9], ILRMA [11], DMNMF [15], the conventional
TCNMF [12], and the proposed TCNMF. IVA and ILRMA
estimate the complex-valued demixing matrix Wi , thus aresounds and is based on MAP estimation with the gamma
distribution prior. Experiments using simulated mixture signals
showed that the proposed method could achieve the highest
bleeding-sound-reduction performance. Since the proposed
method has three hyperparameters, an efficient parametertuning method is necessary and is for future work.
ACKNOWLEDGMENT
This work was partly supported by JSPS KAKENHI Grant
Numbers 19K20306 and 19H01116.
Fig. 5. Comparison of SDR improvements, where each bar is average over
10 different observed mixtures and 4 instrumental sources.

phase-aware BSS methods. The other methods are the phaseinsensitive methods that only use amplitude or power spectrograms. The initial value of Wi for IVA and ILRMA was set to
an inverse matrix of the initial mixing matrix used in DMNMF
and the conventional and proposed TCNMFs. We also used
the numerically stable update rule of the demixing matrix
in both IVA and ILRMA, which is called iterative source
steering [30], and the estimated source was recovered using(7). We then applied the projection-back technique [31] to theestimated signal to recover the frequency-wise signal scales.
For DMNMF and the conventional and proposed TCNMFs,
we used Wiener filtering (46) to obtain the estimated source.
For ILRMA and DMNMF, the number of basis vectors in the
NMF source model, L, was set to 10, 30, and 80.
B. Results
Figure 5 shows the average performance comparison among
the five methods, where the hyperparameters of the conventional and proposed TCNMF were experimentally determinedand set to µ = 0.56, k = 1.25, θ = 0.6, and α = 0.006, which
provided the best performance in this experiment. We canconfirm that the phase-aware BSS methods, IVA and ILRMA,
cannot reduce bleeding sound. This is because the observed
mixture signal in this experiment was produced using thenonnegative random mixing matrix Ai as (49), and the phase
information is useless for estimating the demixing matrix. As aresult, many artificial distortions are produced in the estimatedsignals of IVA and ILRMA, degrading their SDR performance.
DMNMF has the potential to reduce bleeding sound, but its
performance did not exceed 0 dB. This result indicates thedifficulty of parameter optimization in DMNMF. For both theconventional and proposed TCNMFs, we can confirm thatthe average SDR improvements exceed 0 dB. In particular,
the proposed TCNMF outperformed the conventional TCNMF
by more than 2.5 dB. This improvement is significant toachieve high-quality post-processing or sound reinforcementof a musical performance.
V. C ONCLUSION
We aimed to reduce the bleeding sound in the observed signal obtained with close microphones. We proposed a TCNMF
method that regularizes the relative leakage levels of bleeding

R EFERENCES[1] M. Brandstein and D. Ward, Microphone Arrays: Signal Processing
Techniques and Applications, Springer-Verlag Berlin Heidelberg, 2001.
[2] H. L. Van Trees, Optimum Array Processing, John Wiley and Sons, New
York, 2002.
[3] X. Yu, D. Hu, J. Xu, Blind Source Separation: Theory and Applications,John Wiley and Sons, New York, 2014.
[4] H. Sawada, N. Ono, H. Kameoka, D. Kitamura, and H. Saruwatari, “A
review of blind source separation methods: Two converging routes to
ILRMA originating from ICA and NMF,” APSIPA Trans. Signal andInfo. Process., vol. 8, no. e12, pp. 1–14, 2019.
[5] P. Smaragdis, “Blind separation of convolved mixtures in the frequency
domain,” Neurocomputing, vol. 22, pp. 21–34, 1998.
[6] H. Saruwatari, T. Kawamura, T. Nishikawa, A. Lee, and K. Shikano,“Blind source separation based on a fast-convergence algorithm combining ICA and beamforming,” IEEE Trans. Audio, Speech, and Lang.
Process., vol. 14, no. 2, pp. 666–678, 2006.
[7] A. Hiroe, “Solution of permutation problem in frequency domain ICA
using multivariate probability density functions,” Proc. Int. Conf. Independent Compon. Anal. Blind Source Separation, pp. 601–608, 2006.
[8] T. Kim, H. T. Attias, S.-Y. Lee, and T.-W. Lee, “Blind source separation
exploiting higher-order frequency dependencies,” IEEE Trans. Audio,
Speech, and Lang. Process., vol. 15, no. 1, pp. 70–79, 2007.
[9] N. Ono, “Stable and fast update rules for independent vector analysis
based on auxiliary function technique,” Proc. IEEE Workshop Appl.
Signal Process. Audio Acoust., pp. 189–192, 2011.
[10] D. Kitamura, N. Ono, H. Sawada, H. Kameoka, and H. Saruwatari, “Determined blind source separation unifying independent vector analysis
and nonnegative matrix factorization,” IEEE/ACM Trans. Audio, Speech,
and Lang. Process., vol. 24, no. 9, pp. 1626–1641, 2016.
[11] D. Kitamura, N. Ono, H. Sawada, H. Kameoka, and H. Saruwatari,“Determined blind source separation with independent low-rank matrix
analysis,” in Audio Source Separation, S. Makino, Ed., pp. 125–155.
Springer, Cham, 2018.
[12] M. Togami, Y. Kawaguch, H. Kokubo, and Y. Obuchi, “Acoustic echo
suppressor with multichannel semi-blind non-negative matrix factorization,” Proc. Asia-Pacific Signal Info. Process. Assoc. Annu. Summit
Conf., pp. 522–525, 2010.
[13] H. Chiba, N. Ono, S. Miyabe, Y. Takahashi, T. Yamada, and S. Makino,“Amplitude-based speech enhancement with nonnegative matrix factorization for asynchronous distributed recording,” Proc. Int. Workshop
Acoustic Signal Enhancement, pp. 203–207, 2014.
[14] Y. Murase, H. Chiba, N. Ono, S. Miyabe, Y. Takahashi, T. Yamada,and S. Makino, “On microphone arrangement for multichannel speechenhancement based on nonnegative matrix factorization in time-channel
domain,” Proc. Asia-Pacific Signal Info. Process. Assoc. Annu. Summit
Conf., 2014.
[15] T. Taniguchi and T. Masuda, “Linear demixed domain multichannel
nonnegative matrix factorization for speech enhancement,” Proc. IEEE
Int. Conf. Acoust., Speech Signal Process., pp. 476–480, 2017.
[16] O. Das, J. O. Smith, and J. S. Abel, “Microphone cross-talk cancellationin ensemble recordings with maximum likelihood estimation,” Proc.
Audio Eng. Soc. Convention, 2021.
[17] D. D. Lee and H. S. Seung, “Learning the parts of objects by nonnegative matrix factorization,” Nature, vol. 401, no. 6755, pp. 788–791,
1999.
[18] D. D. Lee and H. S. Seung, “Algorithms for non-negative matrix
factorization” Proc. Neural Info. Process. Syst., pp. 556–562, 2000.

[19] A. A. Nugraha, A. Liutkus, and E. Vincent, “Multichannel audio source
separation with deep neural networks,” IEEE/ACM Trans. Audio, Speech,
and Lang. Process., vol. 24, no. 9, pp. 1652–1664, 2016.
[20] N. Makishima, S. Mogami, N. Takamune, D. Kitamura, H. Sumino,
S. Takamichi, H. Saruwatari, and N. Ono, “Independent deeply learned
matrix analysis for determined audio source separation,” IEEE/ACM
Trans. Audio, Speech, and Lang. Process., vol. 27, no. 10, pp. 1601–
1615, 2019.
[21] H. Kameoka, L. Li, S. Inoue, and S. Makino, “Supervised determinedsource separation with multichannel variational autoencoder,” NeuralComput., vol. 31, no. 9, pp. 1891–1914, 2019.
[22] N. Makishima, Y. Mitsui, N. Takamune, D. Kitamura, H. Saruwatari,
Y. Takahashi, and K. Kondo, “Independent deeply learned matrix
analysis with automatic selection of stable microphone-wise update andfast sourcewise update of demixing matrix,” Signal Process., vol. 178,
107753, 2021.
[23] T. Nakamura, S. Kozuka, and H. Saruwatari, “Time-domain audio source
separation with neural networks based on multiresolution analysis,”
IEEE/ACM Trans. Audio, Speech, and Lang. Process., vol. 29, pp. 1687–
1701, 2021.
[24] A. T. Cemgil, “Bayesian inference for nonnegative matrix factorisation models,” Computational Intelligence and Neuroscience, vol. 2009,no. 785152, 2009.

[25] O. Yılmaz and S. Rickard, “Blind separation of speech mixtures via
time-frequency masking,” IEEE Trans. Signal Process., vol. 52, no. 7,
pp. 1830–1847, 2004.
[26] Y. Sun, P. Babu, and D. P. Palomar, “Majorization-minimization algorithms in signal processing, communications, and machine learning,”
IEEE Trans. Signal Process., vol. 65, no. 3, 2017.
[27] D. Kitamura, H. Saruwatari, H. Kameoka, Y. Takahashi, K. Kondo,
and S. Nakamura, “Multichannel signal separation combining directional clustering and nonnegative matrix factorization with spectrogram
restoration,” IEEE/ACM Trans. Audio, Speech, and Lang. Process.,
vol. 23, no. 4, pp. 654–669, 2015.
[28] D. Kitamura, “Open dataset: songKitamura,” http://d-kitamura.net/
dataset en.html. Accessed 10 July 2021.
[29] E. Vincent, R. Gribonval, and C. Févotte, “Performance measurement inblind audio source separation,” IEEE Trans. Audio, Speech, and Lang.
Process., vol. 14, no. 4, pp. 1462–1469, 2006.
[30] S. Robin and N. Ono, “Fast and stable blind source separation with
rank-1 updates,” Proc. IEEE Int. Conf. Acoust., Speech Signal Process.,pp.236–240, 2020.
[31] N. Murata, S. Ikeda, and A. Ziehe, “An approach to blind source separation based on temporal structure of speech signals,” Neurocomputing,vol. 41, no. 1–4, pp. 1–24, 2001.