# Trajectory phase transitions in non-interacting systems: all-to-all dynamics
and the random energy model

Juan P. Garrahan  School of Physics and Astronomy, University of Nottingham,
Nottingham, NG7 2RD, UK Centre for the Mathematics and Theoretical Physics of
Quantum Non-Equilibrium Systems, University of Nottingham, Nottingham, NG7
2RD, UK  Chokri Manai  Department of Mathematics, TU Munich, Germany  Munich
Center for Quantum Science and Technology, Munich, Germany  Simone Warzel
Department of Mathematics, TU Munich, Germany  Munich Center for Quantum
Science and Technology, Munich, Germany  Department of Physics, TU Munich,
Germany

###### Abstract

We study the fluctuations of time-additive random observables in the
stochastic dynamics of a system of $N$ non-interacting Ising spins. We mainly
consider the case of all-to-all dynamics where transitions are possible
between any two spin configurations with uniform rates. We show that the
cumulant generating function of the time-integral of a normally distributed
quenched random function of configurations, i.e., the energy function of the
random energy model (REM), has a phase transition in the large $N$ limit for
trajectories of any time extent. We prove this by determining the exact limit
of the scaled cumulant generating function. This is accomplished by connecting
the dynamical problem to a spectral analysis of the all-to-all quantum REM. We
also discuss finite $N$ corrections as observed in numerical simulations.

##  1 Introduction

In statistical mechanics we are used to studying static phase transitions from
singularities in partition sums [12]: the value of a control parameter at
which the free-energy becomes non-analytic (in the infinite-size limit)
indicates that the equilibrium ensemble of configurations undergoes a phase
change. The standard equilibrium ensemble method can be generalised
straightforwardly to stochastic dynamics by replacing configurations with
trajectories, static observables with (time-extensive) functions of
trajectories, and the partition sum with the corresponding moment generating
function of the trajectory observable [36, 25]. The “thermodynamics of
trajectories” approach [33, 16, 21] allows to study dynamical or “trajectory”
phase transitions, that is, singular changes in the nature of dynamical
fluctuations that often are not reflected in (thermo)static properties or
occur at different parameters of the model. The singularities of the relevant
large deviation (LD) functions [39] reveal phase transitions in, for example,
the dynamical activity of glassy systems [17, 20, 38], in time-integrated
currents in exclusion processes [15, 5, 22], and in (active) work in active
matter [34].

An interesting question is what occurs in a system of many degrees of freedom
whose dynamics is non-interacting when one considers the fluctuations of a
(quenched) random trajectory-observable that couples them. Our main object of
interest will be a system of $N$ Ising spins which all flip independently from
each other. For the case of non-random local observables and independent spins
with single spin-flip dynamics recent results [40] show that in certain cases
there is a phase transition in the LD function. While, naively, one might
expect nothing interesting to occur due to the non-interacting nature of the
dynamics, these results indicate that the optimal way to generate large
fluctuations is by means of effectively highly correlated dynamics which is
singularly different from the typical dynamics [35].

Here we start addressing the problem of random and long-ranged trajectory
observables by considering the time integral of a function of configurations
whose values are normally distributed with zero mean and variance $N$, that
is, the energy function of the simplest mean-field spin glass: the random
energy model (REM) [14, 10]. For simplicity we will consider dynamics which is
all-to-all, that is, allowed configuration changes are those where any number
of spins can flip simultaneously and independently. We also comment on the
case of single-spin flips, which corresponds to the quantum random energy
model (QREM).

The general problem we consider here has relevance in several areas. One is
the minimisation via trajectory sampling of (quasi) random cost functions
[28], which arises for example when training neural networks. A second one is
in connection to measurement induced phase transitions in quantum systems [27,
37], where the calculation of Renyi entropies reduces to computing the optimal
dynamics of a random coupling function [1, 4] in a system of classical
replicas which evolve independently.

##  2 Unbiased dynamics

Any continuous-time Markov process with trajectories
$\boldsymbol{\omega}:[0,\infty)\to\mathcal{Q}_{N}$ on the configuration space
$\mathcal{Q}_{N}:=\\{-1,1\\}^{N}$of $N$ Ising spins is uniquely characterised
in terms of the transition rates $w_{\boldsymbol{\sigma}\to\boldsymbol{\tau}}$
of spin configurations
$\boldsymbol{\sigma}=(\sigma_{1},\dots,\sigma_{N})\in\mathcal{Q}_{N}$ to any
other configuration $\boldsymbol{\tau}$, and the associated escape rates
$r_{\boldsymbol{\sigma}}:=\sum_{\boldsymbol{\tau}\neq\boldsymbol{\sigma}}w_{\boldsymbol{\sigma}\to\boldsymbol{\tau}}$.
The latter governs the law,
$r_{\boldsymbol{\sigma}}e^{-r_{\boldsymbol{\sigma}}\Delta t}$, of the sojourn
time $\Delta t$ until the next jump out of $\boldsymbol{\sigma}$. In the
following, we choose $w_{\boldsymbol{\sigma}\to\boldsymbol{\tau}}:=N2^{-N}$
independent of the configuration. Since the connectivity of this jump process
is then described by the complete graph on $2^{N}$ vertices (i.e. spin
configurations), this dynamics is called the completely connected or all-to-
all stochastic dynamics on Ising configurations. Using Dirac’s notation, in
which $\left\\{|\boldsymbol{\sigma}\rangle\ |\
\boldsymbol{\sigma}\in\mathcal{Q}_{N}\right\\}$ stands for the canonical
orthonormal basis in the Hilbert space
$\ell^{2}(\mathcal{Q}_{N})\equiv\otimes_{j=1}^{N}\mathbb{C}^{2}$, the
generator of this Markov process is given by

$W:=\sum_{\begin{subarray}{c}\boldsymbol{\sigma},\boldsymbol{\tau}\in\mathcal{Q}_{N}\\\
\boldsymbol{\sigma}\neq\boldsymbol{\tau}\end{subarray}}w_{\boldsymbol{\sigma}\to\boldsymbol{\tau}}\,|\boldsymbol{\tau}\rangle\langle\boldsymbol{\sigma}|-\sum_{\boldsymbol{\sigma}\in\mathcal{Q}_{N}}r_{\boldsymbol{\sigma}}\,|\boldsymbol{\sigma}\rangle\langle\boldsymbol{\sigma}|=N\left(|-\rangle\langle-|-\mathbbm{1}\right)$

in terms of the orthogonal projection $|-\rangle\langle-|$ onto the ’flat
state’ defined by $\langle\boldsymbol{\sigma}|-\rangle=2^{-N/2}$. In its
probabilistic interpretation, $W$ is considered an operator on
$\ell^{1}(\mathcal{Q}_{N})$ and acts on probability distributions
$|p_{t}\rangle$, i.e.
$p_{t}(\boldsymbol{\sigma})\equiv\langle\boldsymbol{\sigma}|p_{t}\rangle\geq
0$ and
$\sum_{\boldsymbol{\sigma}\in\mathcal{Q}_{N}}p_{t}(\boldsymbol{\sigma})=1$.
The dynamics of any initial distribution is governed by the master equation

$\partial_{t}|p_{t}\rangle=W|p_{t}\rangle.$

The completely connected stochastic dynamics can be regarded as a further
simplification of the dynamics of independent spin flips at infinite
temperature. The latter is generated by
$\widehat{W}:=\sum_{j=1}^{N}\left(X_{j}-\mathbbm{1}\right)$, in terms of the
Pauli-$X$ matrices, which flip the $j$th spin, i.e.
$X_{j}|\boldsymbol{\sigma}\rangle=|\sigma_{1},\dots,-\sigma_{j},\dots,\sigma_{N}\rangle$.
Both Markov processes are irreducible and share the equidistribution $p_{\rm
ss}(\boldsymbol{\sigma}):=2^{-N}$ as its unique invariant measure. One
difference is their spectral gap, which governs the rate of approach to the
equidistribution. While the spectral gap is $N$ in the case of $W$, it is $2$
in the case of $\widehat{W}$. In this paper we focus on the completely
connected dynamics $W$ and only comment on the single spin-flip dynamics
$\widehat{W}$.

The dynamics generated by $W$ (and $\widehat{W}$) is “infinite temperature” in
the sense that transitions are completely independent of the initial and final
states. The operator $W$ is therefore bi-stochastic,
$\langle-|W=0,\,W|-\rangle=0$, with the first equality indicating conservation
of probability, and the second that the stationary state is also the flat
state (the stationary probability vector being $2^{-N/2}|-\rangle$). Since the
dynamics of all spins is independent, all correlation functions are
unconnected.

##  3 Trajectory observable and REM

We study the statistics under the above defined all-to-all independent
dynamics of a trajectory observable chosen to explore the energy landscape of
the REM [14, 10]. The REM, $U:\mathcal{Q}_{N}\to\mathbb{R}$, is a Gaussian
random field (with randomness independent of the Markov process) in which the
values $U(\boldsymbol{\sigma})$ are distributed independently for all
$\boldsymbol{\sigma}\in\mathcal{Q}_{N}$ with identical normal law uniquely
characterised by zero mean and covariance $N$. The units are chosen so that
the REM’s large deviations occur on order $N$ which agrees with the norm of
$W$. In this context, we recall [24, 10] that the REM’s minimum (and similarly
for its maximum) satisfies the extremal value statistics:

$\mathbb{P}\left(\min U\geq
u_{N}(x)\right)=\left(1-2^{-N}e^{-x+o(1)}\right)^{2^{N}}$ (1)

for any $x$ with the scaling function
$u_{N}(x)\coloneqq-\beta_{c}N+\frac{\ln(N\ln
2)-\ln(4\pi)}{2\beta_{c}}-\frac{x}{\beta_{c}}$, where $\beta_{c}=\sqrt{2\ln
2}$ and $\mathbb{P}$ denotes the joint law of the REM. In particular, the
minimal energy of the REM is roughly at $-\beta_{c}N$.

The trajectory observable we consider is (up to a factor of $t$) the empirical
average of the REM energy along a trajectory $\boldsymbol{\omega}$ of the
Markov process

$U_{t}[\boldsymbol{\omega}]:=\int_{0}^{t}U\left(\boldsymbol{\omega}(s)\right)\
ds,\quad t>0.$

We will be interested in the probability distribution of this quantity under
the law $\mathbb{P}_{t}$ on trajectories associated with $W$ up to time $t$
with the initial spin configurations equally distributed. The main result of
this short note is a proof of a large deviation principle for this
distribution in the limit of large system size $N$ (for trajectories of any
time extent $t$). This large deviation principle is described in terms of the
moment generating function

$Z(t,\lambda):=\int e^{-\lambda U_{t}[\boldsymbol{\omega}]}\
\mathbb{P}_{t}(d\boldsymbol{\omega})=\sum_{\boldsymbol{\sigma},\boldsymbol{\tau}\in\mathcal{Q}_{N}}2^{-N}\langle\boldsymbol{\sigma}|e^{t\left(W-\lambda
U\right)}|\boldsymbol{\tau}\rangle=\langle-|e^{t\left(W-\lambda
U\right)}|-\rangle.$ (2)

Here the second equality is due to the Feynman-Kac formula for the Markov
process under consideration (cf. [23, 26]). Crucially, this formula connects
the question concerning the (a)typical behavior of $U_{t}$ to properties of
the tilted generator

$W_{\lambda}:=W-\lambda U,$

which is a random matrix on $\ell^{2}(\mathcal{Q}_{N})$. Note that by
substituting $W$ by $\widehat{W}$, this random matrix coincides, up to a
constant shift and change of sign, with the Hamiltonian of the QREM – one of
the simplest quantum spin glass models [19, 30, 31, 29]. In our case, the
operator $H_{\lambda}=-W_{\lambda}+N\mathbbm{1}$ instead corresponds to the
Hamiltonian with an all-to-all kinetic energy studied in [2]. Due to the
symmetry of the REM’s distribution the parameter $\lambda$ can be taken non-
negative without loss of generality, and the large deviation function also
known as scaled cumulant generating function (SCGF) is then given by

$\theta(t,\lambda):=\lim_{N\to\infty}\frac{1}{Nt}\ln Z(t,\lambda).$

The SCGF plays the role of a free energy for trajectory ensembles.

It is important to emphasise that what we are considering here is very
different from the study of classical thermal dynamics of the REM under
Glauber or Metropolis schemes, as in e.g. [6, 7, 11, 18]. In those cases the
dynamical Markov generator is interacting (as transitions depend on changes in
$U$) and what is studied are the typical trajectories under that interacting
dynamics. In contrast we study rare trajectories under the non-interacting
dynamics generated by $W$ with large fluctuations of $U_{t}$.

##  4 Trajectory phase diagram

Our main result is the following:

######  Theorem 1.

For any $t>0,\lambda\geq 0$ and almost all realisations of the REM:

$\theta(t,\lambda)=\max\left\\{0,t^{-1}p_{0}(t\lambda)-1\right\\},$ (3)

with

$p_{0}(\beta):=\begin{cases}\frac{\beta^{2}}{2},&\beta\leq\beta_{c}:=\sqrt{2\ln
2}\\\ \beta\beta_{c}-\ln 2,&\beta>\beta_{c}.\end{cases}$ (4)

Before spelling out the short proof of Theorem 1 in Section 7 below, let us
put this result in some context and discuss some consequences. The quantity
defined in (4) is the pressure corresponding to the REM’s static (normalised)
partition function at inverse temperature $\beta$:

$p_{0}(\beta)=\lim_{N\to\infty}\frac{1}{N}\ln\frac{1}{2^{N}}\sum_{\boldsymbol{\sigma}}e^{-\beta
U(\boldsymbol{\sigma})}.$ (5)

The critical value $\beta_{c}=\sqrt{2\ln 2}$ corresponds to the inverse of the
REM’s freezing temperature into a spin glass phase with 1-step replica
symmetry breaking, cf. [10].

The phase diagram resulting from Theorem 1 is thus composed of three regimes
depicted in Fig. 1(a):

  1. 1.

An Active dynamical phase in which the Markov generator $W$ dominates over the
tilting, and which is characterised by $\theta(t,\lambda)=0$ and the specific
activity being unity (see below). It is separated from the remaining regimes
by a first-order transition line. This regime persists for all
$|\lambda|<\beta_{c}(2t)^{-1}+\beta_{c}^{-1}$ in case $t^{-1}<2\beta_{c}^{-2}$
and $|\lambda|<\sqrt{2t^{-1}}$ in case $t^{-1}\geq 2\beta_{c}^{-2}$.

  2. 2.

A regime of vanishing activity which occurs for
$t^{-1}<2|\lambda|\beta_{c}^{-1}-2\beta_{c}^{-2}$ and which is dominated by
the REM’s extreme values where the system localises. This regime is related to
the spin-glass phase of the REM. We call this the Inactive-1 dynamical phase.

  3. 3.

The remaining parameter regime corresponds to a second inactive regime which
we term Inactive-2 dynamical phase. It occurs only if $t^{-1}>2/\beta_{c}^{2}$
and is related to the classical paramagnetic phase of the REM.

In particular, in the long-time limit, $t\to\infty$, the value
$\lambda=\beta_{c}^{-1}$ separates the Active and Inactive-1 phases, the
latter dominating at at large $\lambda$. Not surprisingly, this transition in
the largest eigenvalue of the tilted generator $W_{\lambda}=W-\lambda U$
reflects the known location of its quantum-phase transition. As we will recall
in Section 7 below, the eigenvector corresponding to the largest eigenvalue
changes near $\lambda=\beta_{c}^{-1}$ from a delocalised state resembling
$|-\rangle$ (indicating that trajectories visit all states equally giving rise
to large activity) to a state localised at the REM’s maximising spin
configuration (corresponding to trajectories that are inactive as they do not
move away from this configuration).

The classification above of the trajectory phases in terms of their activity
is obtained as follows. The dynamical activity is the total number of
configuration changes in a trajectory. It can be calculated through the same
tilting method used above for the time-integrated REM energy. Specifically, if
we define the doubly tilted partition sum
$Z(t,\lambda,s):=\langle-|e^{tW_{\lambda,s}}|-\rangle$ with
$W_{\lambda,s}=Ne^{-s}|-\rangle\langle-|-N(1-2^{-N}(1-e^{-s}))\mathbbm{1}-\lambda
U$ (where the additional tilting by $e^{-s}$ of the off-diagonal part of $W$
allows to count jumps in trajectories), we get the activity from
$-\partial_{s}\log Z(t,\lambda,s)|_{s=0}$. Using the results above it is easy
to see that the average activity per unit space and time is unity in the
active phase and zero in the two inactive phases.

Figure 1:  (a) Dynamical phase diagram in the limit of $N\to\infty$. The
abscissa is the counting field $\lambda$ conjugate to the time-integrated REM
energy and the ordinate the inverse of the trajectory length $t$. The full
lines indicate first-order transitions between Active and Inactive-1 or
Inactive-2 trajectory phases, while the double line indicates a 1-step RSB
transition between the Inactive-1 and Inactive-2 phases. (b) Large deviation
function $\varphi(t,u)$ at two different values of $t$. For
$t^{-1}=1<2\beta_{c}^{-2}$ (blue lines) the rate function is one of
coexistence between Active (which has $u=0$) and Inactive-1 (which has
$u=\pm\beta_{c}$, by symmetry). The linear portion of the rate function is the
Maxwell construction indicative of phase coexistence (in time). For
$t^{-1}=10>2\beta_{c}^{-2}$ the rate function describes the coexistence
between the three phases (black). The linear portion between $0$ and
$\sqrt{2t}$ is now the first-order coexistence between Active and Inactive-2.
In Inactive-2 $|u|$ can take values with decreasing probability between
$\sqrt{2t}$ and $\beta_{c}$. The rate function is infinite for any $|u|$
beyond $\beta_{c}$ (indicative of zero probability for such trajectories). (c)
Thermal phase diagram of the all-to-all QREM for comparison to (a).

Via the Gärtner-Ellis theorem [13], the rate function of the large deviation
principle obeyed by $U_{t}$ is given by the Legendre-Fenchel transformation

$\varphi(t,u):=\sup_{\lambda}\left(u\lambda-\theta(t,\lambda)\right)=\begin{cases}|u|\sqrt{\frac{2}{t}},&|u|\leq\min\left\\{\sqrt{2t},\beta_{c}\right\\},\\\
1+\frac{u^{2}}{2t},&\mbox{else},\\\ \infty,&|u|>\beta_{c}.\end{cases}$

Note that, although in Theorem 1 initially defined only for $\lambda\geq 0$ ,
the function $\lambda\to\theta(t,\lambda)$ extends to all real values by
symmetry. The rate function $u\to\varphi(t,u)$ is then symmetric as well. For
times $t>\beta_{c}^{2}/2=\ln 2$, the second case in the above equation is
absent. As a corollary to Theorem 1 and [13, Thm 2.3.6], we thus obtain the
promised large deviation principle

$\displaystyle-\inf_{u\in I^{\circ}}t\varphi(t,u)$
$\displaystyle\leq\liminf_{N\to\infty}\frac{1}{N}\ln\mathbb{P}_{t}\left((Nt)^{-1}U_{t}\in
I\right)$
$\displaystyle\leq\limsup_{N\to\infty}\frac{1}{N}\ln\mathbb{P}_{t}\left((Nt)^{-1}U_{t}\in
I\right)=-\inf_{u\in\overline{I}}t\varphi(t,u)$ (6)

which holds for any Borel set $I\subset\mathbb{R}$ and any $t>0$. The rate
function is shown in Fig. 1(b) for two different times.

Clearly, under the apriori measure $\mathbb{P}_{t}$, which favors rapid
changes of spin configurations at the rate $N(1-2^{-N})$, the typical value of
the REM’s empirical energy density $N^{-1}U_{t}[\boldsymbol{\omega}]$ along
any trajectory $\boldsymbol{\omega}$ is close to zero. The fluctuations about
this typical behavior are described by (4): close to $u=0$, these fluctuations
are linearly suppressed with a rate proportionally to $N\sqrt{2t}$. Tilting
the apriori measure, one encounters one or two phase transitions depending on
whether $t>\ln 2$ or not. If $t<\ln 2$, one enters a regime
$\sqrt{2t}<|u|<\beta_{c}$ with Gaussian fluctuations. Beyond this, i.e., at
energy densities of the order of the REM’s maximum or minimum (1), the energy
density effectively stops fluctuating. Trajectories freeze for long times in
the REM’s extremal values.

##  5 Comparison to the thermal phase diagram

The dynamical partition sum of the stochastic system we are considering is
reminiscent of a quantum (thermal and static) partition sum for the all-to-all
version of the QREM. While the calculation of both is analogous, there are
some important differences. Specifically, if we consider the tilted generator
$W_{\lambda}$ as (minus) a Hamiltonian, the (specific) free energy of the
associated quantum problem at temperature $T$ is

$f(T,\lambda):=\lim_{N\to\infty}\frac{T}{N}\ln\frac{1}{2^{N}}{\operatorname{Tr}\,}e^{(W-\lambda
U)/T}=\max\left\\{-T\ln 2,Tp_{0}(\lambda/T)-1\right\\}.$ (7)

As we will explain in Section 7 below, the last equality follows
straightforwardly from results on the eigenvalues in [2].

Similarly to Theorem 1, from (7) we see that, depending on coupling and
temperature, the all-to-all QREM can be in three different phases, a
delocalised quantum paramagnetic phase (QPM), a localised spin-glass phase
(SG) and a classical paramagnetic phase (PM), see Fig. 1(c). These three
static quantum phases are similar to the dynamical ones of the stochastic
problem. But is worth pointing out that at $T>0$ the (thermo)static phase
transitions described by $f(T,\lambda)$, do not coincide with the dynamic
phase transitions described by $\theta(t,\lambda)$. These differences arise
because of the boundary vectors in the dynamical partition sum versus the
trace in the static quantum one. For a comparison of the phase diagrams, see
Figs. 1(a) and 1(c).

##  6 Numerical illustration of finite size corrections

The exact results above are for the limit $N\to\infty$. At finite $N$ there
are of course finite-size corrections and sample-to-sample fluctuations
between different realisations of the disorder $U$. Using numerics, we now
illustrate some of these finite-size effects. (A comprehensive numerical study
of both the all-to-all and single spin-flip problem will be presented in a
future publication.)

Figure 2:  (a) Dynamical order parameter
$-\partial_{\lambda}\theta(\lambda,t)$ for various times $t$, from exact
diagonalisation (ED) for one disorder realisation with system size $N=12$. (b)
The top panel shows a typical trajectory of time extent $t=10$, corresponding
to unbiased dynamics, $\lambda=0$ (generated from $W$ via standard continuous-
time Monte Carlo from a random initial configuration). The bottom panel shows
a characteristic trajectory of the inactive phase, $\lambda=2$. This was
obtained via transition path sampling (TPS, see main text). The inset to the
lower panel shows the convergence of the trajectory sampling (black): each
iteration is a different trajectory and their time-integrated energy converges
to the ED result (red dashed) with enough TPS iterations. (c) Transition point
$\lambda_{c}(N)$ for $t=\infty$ as a function of system size, averaged over 20
disordered realisations. The dashed line is the large $N$ value of
$\lambda_{c}$.

When the system size is not too large the dynamical partition sum (2) can be
computed numerically using exact diagonalisation (ED). We illustrate results
for one disorder realisation in a system of size $N=12$ of the dynamical order
parameter:

$u(\lambda,t):=\frac{1}{NtZ(t,\lambda)}\int
U_{t}[\boldsymbol{\omega}]e^{-\lambda U_{t}[\boldsymbol{\omega}]}\
\mathbb{P}_{t}(d\boldsymbol{\omega})=-\frac{\partial}{\partial\lambda}\theta_{N}(t,\lambda)$
(8)

Figure 2(a) shows the following: (i) for finite size the phase transitions
turn into crossovers, as expected; (ii) for all $t$ there is a crossover from
$u\approx 0$ at $\lambda=0$ to a large negative $u$ for large $\lambda$,
eventually reaching the minimum of the potential (which changes from sample to
sample); (iii) these crossovers are sharper the longer $t$, also expected due
to the preference of the boundary states in (2) for the delocalised state.

In Fig. 2(b) we show representative trajectories for two values of $\lambda$
for $t=10$. We plot the instantaneous energy as a function of time in the
trajectory. The top panel shows a typical trajectory of the dynamics
corresponding to $\lambda=0$, cf. the red square in Fig. 2(a). This trajectory
generated by $W$ is sampled using standard (continuous-time) Monte Carlo [9].
Since the unbiased dynamics connects all configurations with equal rates the
trajectory jumps between the energy values: it corresponds to the phase which
has high activity and is delocalised. The bottom panel shows a characteristic
trajectory for $\lambda=2$, cf. the red circle in Fig. 2(a). This is a rare
event (exponentially suppressed in $N$ and $t$) of the dynamics, and as such
cannot be easily sampled from running Monte Carlo with $W$ (since
$W_{\lambda}$ is not a stochastic operator). We obtain such rare trajectories
instead by performing importance sampling in trajectory space using transition
path sampling (TPS) [8], essentially a Monte Carlo method in trajectory space
that aims to “equilibrate” to a reweigthed trajectory distribution
$Z(t,\lambda)^{-1}e^{-\lambda U_{t}[\boldsymbol{\omega}]}$ (supplemented with
bridge moves to improve acceptance; we will provide details of this method in
a future publication). The inset to the lower panel shows the convergence of
our TPS approach: it shows the evolution of the sampled trajectories with TPS
iterations by showing their $U_{t}[\boldsymbol{\omega}]$ (per unit time). The
$U_{t}[\boldsymbol{\omega}]$ in the inset converges eventually to the value
expected at $\lambda=2$, showing that TPS converges to the tilted trajectory
ensemble. The trajectory shown in the lower panel is the last trajectory from
TPS. It is very different from the typical one in the upper panel: it has very
low activity and is localised for most of the time in the minium energy
configuration, corresponding to the Inactive-1 dynamical phase. Note that
while we only illustrate the numerics for the size $N=12$, TPS can be used for
larger system sizes in contrast to ED.

Figure 2(c) shows the location of the critical $\lambda$ for $t=\infty$,
averaged over 20 realisations of the disorder, for different systems sizes.
The transition point is inferred from the maximum of the dynamical
susceptibility $\partial_{\lambda}^{2}\theta_{N}(\lambda)$, where
$\theta_{N}(\lambda)$ is the largest eigenvalue of $W_{\lambda}$. This
eigenvalue is calculated using (9) below, which allows to compute it for
larger sizes than those accessible to ED. The figure suggest a convergence to
the limiting value $\beta_{c}^{-1}$ for large $N$, as expected from the
analytics above.

##  7 Proof of the large deviation result

The Feynman-Kac formula (2) reduces the large deviation problem to a spectral
analysis of the random matrix $W_{\lambda}=W-\lambda U$, which – as motivated
in the introduction – may serve as a toy model to the QREM with a simpler all-
to-all kinetic energy term with. Up to a constant shift and rescaling, the
spectrum of $W_{\lambda}$ has been analysed in [2] both on the macro and
microscopic scale of the eigenvalue process. The main technical tool for
studying $W_{\lambda}$ is rank-one perturbation theory according to which $E$
is an eigenvalue of $W_{\lambda}$ if and only if

$\frac{1}{N}=\langle-|(E+N+\lambda
U)^{-1}|-\rangle=\frac{1}{2^{N}}\sum_{\boldsymbol{\sigma}}\frac{1}{E+N+\lambda
U(\boldsymbol{\sigma})}.$ (9)

The corresponding eigenvectors $\psi_{E}$ satisfy for all
$\boldsymbol{\sigma},\boldsymbol{\tau}\in\mathcal{Q}_{N}$:

$\frac{\langle\boldsymbol{\sigma}|\psi_{E}\rangle}{\langle\boldsymbol{\tau}|\psi_{E}\rangle}=\frac{E+N+\lambda
U(\boldsymbol{\tau})}{E+N+\lambda U(\boldsymbol{\sigma})}.$ (10)

An immediate implication of (9) is the fact that all eigenvalues of
$W_{\lambda}$ aside from the largest one are interlaced with the REM’s
energies and additionally shifted by $-N$ (cf. e.g. [3] and refs. therein for
interlacing and finite-rank perturbation theory). All eigenvalues are almost
surely simple. Moreover, any solution of (9) with $E>(\lambda\beta_{c}-1)N$
(cf. (1)) is independent of the realisation of $U$ up to exponentially small
fluctuations. By the law of large numbers the right-hand side of (9) is then
well approximated as $N\to\infty$ by the integral

$\int_{-\infty}^{\infty}\frac{\exp(-v^{2}/2)}{E+N+\lambda\sqrt{N}v}\frac{dv}{\sqrt{2\pi}}=\frac{1}{E+N}\left(1+\frac{\lambda^{2}}{N(1+E/N)}+\mathcal{O}(N^{-2})\right).$

This explains the following results on the largest eigenvalue
$E_{0}:=\max\sigma(W_{\lambda})$, which are found in [2]:

  1. 1.

In case $\lambda\beta_{c}<1$, on an event with probability exponentially close
to one, the largest eigenvalue is at
$E_{0}:=\max\sigma(W_{\lambda})=\lambda^{2}+\mathcal{O}(N^{-1})$ and the
corresponding eigenvector satisfies
$\langle\boldsymbol{\sigma}|\psi_{E_{0}}\rangle\propto(E_{0}+N+\lambda
U(\boldsymbol{\sigma}))^{-1}$. Since
$2^{-N}\sum_{\boldsymbol{\sigma}}(E_{0}+N+\lambda
U(\boldsymbol{\sigma}))^{-2}$ is of order one up to exponentially small
fluctuations by the law of large numbers, this vector is hence still
delocalised (as in the case $\lambda=0$).

  2. 2.

In case $\lambda\beta_{c}>1$, the largest eigenvalue is at

$E_{0}=\max\sigma(W_{\lambda})=-\lambda\min
U-N-N2^{-N}\left(1-\frac{1}{\lambda\beta_{c}}\right)^{-1}+o(N2^{-N})$ (11)

and the corresponding eigenvector is mostly concentrated on the REM’s
minimising configuration $\boldsymbol{\sigma}_{0}\coloneqq\arg\min U$. This is
specified through the ratios (10). Note that the error term in the above
equation only holds with a probability up to $1-\mathcal{O}(1/N)$, cf. [2, App
A].

In particular, the union of eigenvalues, $\sigma(W_{\lambda})$, when divided
by $N$, converges almost surely to the non-random set
$\\{0\\}\cup[-\lambda\beta_{c}-1,\lambda\beta_{c}-1]$. Together with the
interlacing property, one then also easily arrives at (7) for the free energy
of $W_{\lambda}$.

The proof of Theorem 1 requires slightly more detailed knowledge, since
$\langle-|e^{tW_{\lambda}}|-\rangle$ involves properties of the eigenvectors,
too. The rough picture established in [2] through a more detailed analysis of
the characteristic equation (9) is the following:

  1. 1.

Delocalisation of one eigenstate near energy $0$ is shown to persist up to
$\lambda<\sqrt{2}$. From that value on, $\lambda>\sqrt{2}$, this eigenstates
“melts” into a narrow band of semi-delocalised states near energy $0$.

  2. 2.

The eigenvalue process, when rescaled to order one at some fixed energy
outside $-N$ and $0$, is given by a Poisson process. Correspondingly, outside
those special energies the normalised eigenvectors are localised.

We will need the following result, which is contained in [2, Proof of Thm.
6.3].

######  Proposition 1.

For any $\delta>0$ and any $N$ there is some $a>0$ and an event $\Omega_{N}$
whose complement is summable,
$\sum_{N}\mathbb{P}\left(\Omega_{N}^{c}\right)<\infty$, such that in the event
$\Omega_{N}$ any eigenvalue $E$ of $W_{\lambda}$ with $|E|>\delta N$ and
$|E+N|>\delta N$ has a normalised eigenvector $\psi_{E}$, which satisfies
$|\langle-|\psi_{E}\rangle|^{2}\leq N^{a}\ 2^{-N}$. Moreover, for any such
$E$, there is some $\boldsymbol{\sigma}_{E}\in\\{-1,1\\}^{N}$ such that
$|E+N-\lambda U(\boldsymbol{\sigma}_{E})|\leq\delta N$.

###### Proof of Theorem 1.

The proof proceeds by establishing asymptotically coinciding upper and lower
bounds. For the lower bound, we use Jensen’s inequality to conclude

$\ln\langle-|e^{tW_{\lambda}}|-\rangle\geq
t\langle-|W_{\lambda}|-\rangle=\frac{t\lambda}{2^{N}}\sum_{\boldsymbol{\sigma}}U(\boldsymbol{\sigma}).$

By the law of large numbers, this term converges to zero for almost all
realisations of the REM. For another lower bound, which is sharper in case
$t<p_{0}(t\lambda)$, we estimate

$\displaystyle\langle-|e^{tW_{\lambda}}|-\rangle=\frac{1}{2^{N}}\sum_{\boldsymbol{\sigma},\boldsymbol{\tau}}\langle\boldsymbol{\tau}|e^{tW_{\lambda}}|\boldsymbol{\sigma}\rangle$
$\displaystyle\geq\frac{1}{2^{N}}\sum_{\boldsymbol{\sigma}}\langle\boldsymbol{\sigma}|e^{tW_{\lambda}}|\boldsymbol{\sigma}\rangle$
$\displaystyle\geq\frac{1}{2^{N}}\sum_{\boldsymbol{\sigma}}\exp\left(tN(|\langle\boldsymbol{\sigma}|-\rangle|^{2}-1)-t\lambda
U(\boldsymbol{\sigma})\right),$

where the last step is again by Jensen’s inequality. Using
$|\langle\boldsymbol{\sigma}|-\rangle|^{2}=2^{-N}$ and (5), the combination of
the above estimates yields (3) as a lower bound.

A complementing upper bound is based on Proposition 1. Expanding in
eigenfunctions and splitting the sum over all eigenvalues in three parts
corresponding to energies $E$ with $|E+N|\leq\delta N$, $|E|\leq\delta N$ and
the rest, we write and estimate using Proposition 1:

$\displaystyle\langle-|e^{tW_{\lambda}}|-\rangle=\sum_{E}e^{tE}|\langle-|\psi_{E}\rangle|^{2}\leq
e^{tN(\delta-1)}+e^{tN\delta}+\frac{N^{a}}{2^{N}}e^{tN(\delta-1)}\sum_{\boldsymbol{\sigma}}e^{-t\lambda
U(\boldsymbol{\sigma})}.$ (12)

In the event $\Omega_{N}$ of Proposition 1, we thus conclude

$\limsup_{N\to\infty}\frac{1}{Nt}\ln\langle-|e^{tW_{\lambda}}|-\rangle\leq\max\left\\{0,t^{-1}p_{0}(t\lambda)-1\right\\}+\delta.$

By a Borel-Cantelli argument, this establishes this almost-sure bound on the
upper limit. Since $\delta>0$ is arbitrary, this concludes the proof. ∎

##  8 Outlook: QREM

Let us conclude this note with some conjectures, partial results and
comparison in case $W$ is replaced by the spin-flip dynamics generated by
$\widehat{W}$. In that case, the tilted generator
$\widehat{H}_{\lambda}:=\widehat{W}-\lambda U$ is the QREM. Its low-energy
spectrum as well as the phase transitions in the free energy are well
understood [30, 32]. By the Feynman-Kac formula the dynamical phase transition
is again described in terms of the asymptotic behavior of
$N^{-1}\ln\langle-|e^{t\widehat{H}_{\lambda}}|-\rangle$.

The phase transition in the largest eigenvalue
$\widehat{E}_{0}:=\max\sigma(\widehat{H}_{\lambda})$ occurs on order $N$ at
the same location $\lambda=\beta_{c}^{-1}$ as for $H_{\lambda}$. However, the
finite-volume corrections are different in the localisation regime, i.e. for
all realisations of the REM aside from a set of exponentially small
probability (see [32] for details):

  1. 1.

if $\lambda>\beta_{c}^{-1}$ we have $\widehat{E}_{0}=-\lambda\min
U+(\lambda\beta_{c})^{-1}+\mathcal{O}(N^{-1/4})$,

  2. 2.

if $\lambda<\beta_{c}^{-1}$ we have
$\widehat{E}_{0}=\lambda^{2}+\mathcal{O}(N^{-1/4})$.

Following the steps of the lower bound in the proof of Theorem 1, it is easy
to see that for almost all realizations of the REM one still has:

$\liminf_{N\to\infty}\frac{1}{Nt}\ln\langle-|e^{t\widehat{H}_{\lambda}}|-\rangle\geq\max\left\\{0,t^{-1}p_{0}(t\lambda)-1\right\\}.$
(13)

We conjecture that this bound is sharp. In fact, using the spectral
decomposing as in (12) and decomposing the sum into positive and negative
energies we may again estimate

$\langle-|e^{t\widehat{H}_{\lambda}}|-\rangle\leq\sum_{\begin{subarray}{c}E\sigma(\widehat{H}_{\lambda})\\\
E>0\end{subarray}}e^{tE}\left|\langle-|\psi_{E}\rangle\right|^{2}+1.$

The first sum is estimated trivially by $e^{t\widehat{E}_{0}}$. In case
$\lambda<\beta_{c}^{-1}$ this yields the upper bound
$\limsup_{N\to\infty}\frac{1}{Nt}\ln\langle-|e^{t\widehat{H}_{\lambda}}|-\rangle\leq
0$, which coincides with the lower bound. In case $\lambda>\beta_{c}^{-1}$, we
know from [32] that eigenvalues with energies $E>0$ are in one-to-one
correspondence with values $U(\boldsymbol{\sigma}_{E})=E+\mathcal{O}(1)$. We
conjecture that the local density of states at these energies satisfies
$\lim_{N\to\infty}N^{-1}\ln\langle-|1_{(E-\delta_{N},E+\delta
N)}(\widehat{H}_{\lambda})|-\rangle=-\ln 2$ for all sufficiently small
$\delta>0$. This would prove that (13) is indeed sharp.

## Acknowledgements

JPG acknowledges financial support from EPSRC Grant no. EP/R04421X/1 and the
Leverhulme Trust Grant No. RPG-2018-181. SW thanks the DFG for support under
grant EXC-2111 – 390814868. Numerical simulations were performed using the
Sulis Tier 2 HPC platform funded by EPSRC Grant EP/T022108/1 and the HPC
Midlands+ consortium.

Data: Research data are available from the Nottingham Research Data Management
Repository at http://doi.org/10.17639/nott.7196.

## References

  * [1] U. Agrawal, A. Zabalo, K. Chen, J. H. Wilson, A. C. Potter, J. Pixley, S. Gopalakrishnan, and R. Vasseur.  Entanglement and charge-sharpening transitions in u (1) symmetric monitored quantum circuits.  Preprint arXiv:2107.10279, 2021.
  * [2] M. Aizenman, M. Shamis, and S. Warzel.  Resonances and partial delocalization on the complete graph.  Annales Henri Poincaré, 16(9):1969–2003, 2015.
  * [3] M. Aizenman and S. Warzel.  Random Operators: Disorder Effects on Quantum Spectra and Dynamics.  Graduate Studies in Mathematics Volume 168. AMS, 2015.
  * [4] A. Altland, M. Buchhold, S. Diehl, and T. Micklitz.  Dynamics of measured many-body quantum chaotic systems, 2021.
  * [5] C. Appert-Rolland, B. Derrida, V. Lecomte, and F. van Wijland.  Universal cumulants of the current in diffusive systems on a ring.  Phys. Rev. E, 78(2):021122, aug 2008.
  * [6] G. B. Arous, A. Bovier, and V. Gayrard.  Glauber dynamics of the random energy model.  Communications in Mathematical Physics, 236(1):1–54, 2003.
  * [7] G. B. Arous, A. Bovier, and V. Gayrard.  Glauber dynamics of the random energy model.  Communications in Mathematical Physics, 235(3):379–425, 2003.
  * [8] P. G. Bolhuis, D. Chandler, C. Dellago, and P. L. Geissler.  Transition path sampling: throwing ropes over rough mountain passes, in the dark.  Annu. Rev. Phys. Chem., 53:291, 2002.
  * [9] A. B. Bortz, M. H. Kalos, and J. L. Lebowitz.  New algorithm for Monte-Carlo simulation of Ising spin systems.  J. Comp. Phys., 17(1):10–18, 1975.
  * [10] A. Bovier.  Statistical Mechanics of Disordered Systems: A Mathematical Perspective.  Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 2006.
  * [11] J. Černý and T. Wassmer.  Aging of the metropolis dynamics on the random energy model.  Probability Theory and Related Fields, 167(1):253–303, 2017.
  * [12] D. Chandler.  Introduction to Modern Statistical Mechanics, 1987.
  * [13] A. Dembo and O. Zeitouni.  Large Deviations Techniques and Applications.  Springer, 2nd edition edition, 1998.
  * [14] B. Derrida.  Random-energy model: Limit of a family of disordered models.  Physical Review Letters, 45(2):79–82, 07 1980.
  * [15] B. Derrida.  Non-equilibrium steady states: fluctuations and large deviations of the density and of the current.  J. Stat. Mech., 2007(07):P07023, 2007.
  * [16] J. P. Garrahan.  Aspects of non-equilibrium in classical and quantum systems: Slow relaxation and glasses, dynamical large deviations, quantum non-ergodicity, and open quantum dynamics.  Physica A, 504:130–154, 2018.
  * [17] J. P. Garrahan, R. L. Jack, V. Lecomte, E. Pitard, K. van Duijvendijk, and F. van Wijland.  Dynamical First-Order Phase Transition in Kinetically Constrained Models of Glasses.  Phys. Rev. Lett., 98(19):195702, 2007.
  * [18] V. Gayrard, and L. Hartung.  Dynamic Phase Diagram of the REM,  in: V. Gayrard, L.-P. Arguin, N. Kistler, and I. Kourkova, editors.  Statistical Mechanics of Classical and Disordered Systems, Springer, 2019.
  * [19] Y. Y. Goldschmidt.  Solvable model of the quantum spin glass in a transverse field.  Physical Review B, 41(7):4858–4861, 03 1990.
  * [20] L. O. Hedges, R. L. Jack, J. P. Garrahan, and D. Chandler.  Dynamic order-disorder in atomistic models of structural glass formers.  Science, 323(5919):1309, 2009.
  * [21] R. L. Jack.  Ergodicity and large deviations in physical systems with stochastic dynamics.  Eur. Phys. J. B, 93(4):74, 2020.
  * [22] R. L. Jack, I. R. Thompson, and P. Sollich.  Hyperuniformity and Phase Separation in Biased Ensembles of Trajectories for Diffusive Systems.  Phys. Rev. Lett., 114(6):060601, 2015.
  * [23] M. Keller, D. Lenz, and R. K. Wojciechowski.  Graphs and Discrete Dirichlet Spaces.  Grundlehren der mathematischen Wissenschaften Volume 358. Springer, 2021\.
  * [24] M. R. Leadbetter, G. Lindgren, and H. Rootzén.  Extremes and Related Properties of Random Sequences and Processes.  Springer, 1983.
  * [25] V. Lecomte, C. Appert-Rolland, and F. van Wijland.  Thermodynamic formalism for systems with markov dynamics.  J. Stat. Phys., 127(1):51, 2007.
  * [26] H. Leschke, S. Rothlauf, R. Ruder, and W. Spitzer.  The free energy of a quantum Sherrington–Kirkpatrick spin-glass model for weak disorder.  Journal of Statistical Physics, 182(3):55, 2021.
  * [27] Y. Li, X. Chen, and M. P. A. Fisher.  Quantum zeno effect and the many-body entanglement transition.  Phys. Rev. B, 98:205136, 2018.
  * [28] J. Mair, D. Rose, and J. Garrahan.  Training of neural network ensembles via trajectory sampling.  in preparation, 2022.
  * [29] C. Manai and S. Warzel.  Generalized random energy models in a transversal magnetic field: Free energy and phase diagrams.  Preprint arXive:2007.03290. To appear in: Probab. Math. Phys., 2020.
  * [30] C. Manai and S. Warzel.  Phase diagram of the quantum random energy model.  Journal of Statistical Physics, 180(1):654–664, 2020.
  * [31] C. Manai and S. Warzel.  The de Almeida–Thouless line in hierarchical quantum spin glasses.  Journal of Statistical Physics, 186(1):14, 2021.
  * [32] C. Manai and S. Warzel.  Spectral analysis of the quantum random energy model.  Preprint arXive:2202.00334, 2022.
  * [33] M. Merolle, J. P. Garrahan, and D. Chandler.  Space-time thermodynamics of the glass transition.  Proc. Natl. Acad. Sci. USA, 102(31):10837, 2005.
  * [34] T. Nemoto, E. Fodor, M. E. Cates, R. L. Jack, and J. Tailleur.  Optimizing active work: Dynamical phase transitions, collective motion, and jamming.  Phys. Rev. E, 99:022605, 2019.
  * [35] P. T. Nyawo and H. Touchette.  A minimal model of dynamical phase transition.  EPL, 116(5):50009, 2016.
  * [36] D. Ruelle.  Thermodynamic formalism.  Cambridge University Press, 2004.
  * [37] B. Skinner, J. Ruhman, and A. Nahum.  Measurement-induced phase transitions in the dynamics of entanglement.  Phys. Rev. X, 9:031009, 2019.
  * [38] T. Speck, A. Malins, and C. P. Royall.  First-order phase transition in a model glass former: Coupling of local structure and dynamics.  Phys. Rev. Lett., 109:195703, 2012.
  * [39] H. Touchette.  The large deviation approach to statistical mechanics.  Phys. Rep., 478(1-3):1, 2009.
  * [40] L. M. Vasiloiu, T. H. E. Oakes, F. Carollo, and J. P. Garrahan.  Trajectory phase transitions in noninteracting spin systems.  Phys. Rev. E, 101:042115, 2020.


# Multi-Channel Multi-Speaker ASR Using 3D Spatial Feature

###### Abstract

Automatic speech recognition (ASR) of multi-channel multi-speaker overlapped
speech remains one of the most challenging tasks to the speech community. In
this paper, we look into this challenge by utilizing the location information
of target speakers in the 3D space for the first time. To explore the strength
of proposed the 3D spatial feature, two paradigms are investigated. 1) a
pipelined system with a multi-channel speech separation module followed by the
state-of-the-art single-channel ASR module; 2) a “All-In-One” model where the
3D spatial feature is directly used as an input to ASR system without explicit
separation modules. Both of them are fully differentiable and can be back-
propagated end-to-end. We test them on simulated overlapped speech and real
recordings. Experimental results show that 1) the proposed ALL-In-One model
achieved a comparable error rate to the pipelined system while reducing the
inference time by half; 2) the proposed 3D spatial feature significantly
outperformed (31% CERR) all previous works of using the 1D directional
information in both paradigms.

Index Terms—  multi-channel ASR, multi-speaker ASR, audio-visual ASR, speech
separation, 3D feature

##  1 Introduction

With the development of speech techniques and deep neural networks, dramatic
improvement has been achieved on multiple automatic speech recognition (ASR)
benchmarks [1, 2, 3, 4]. However, it remains a challenging task for multi-
channel multi-speaker overlapped speech recognition due to the interfering
speakers or background noise [5, 6, 7]. To solve this problem, [8, 9, 10]
adopt a pipeline-based paradigm where a preceding speech separation module
separates clean target speech for the back-end standard ASR system to
transcribe. These modules can be jointly trained with ASR criteria to achieve
a boost on the final ASR performance.

To handle the multi-talker speech, the common approach is to use the
spectrogram masking based separation [11, 12, 13], where the weight (mask) of
the target speaker at each time-frequency (T-F) bin is estimated. To identify
the target speech from the mixture, additional information about the target
speaker is required. Common clues including speaker’s voiceprint based
features [14, 15], vision clues [6, 16, 17] and location based features [18,
19, 20] have been proven to be beneficial. This paper focuses on location
based features due to its robustness and accuracy of visual face detection
[6].

However, there is one main issue with these location-based methods. To our
best knowledge, all existing systems using spatial features only consider
speech in a 1D scenario [6, 10, 18]. They assume that the sound sources
traveled like a 2D planar wave and only the 1D directional of arrival (DoA)
information (i.e. azimuth) is considered. As a result, when the coming
directions of two simultaneous speech are close, even if they are far apart in
distance, these systems will fail to separate them (see Fig. 2 (a)).

Contributions: Previous works in this area are extended in two important
directions. First, the traditional 1D directional features are generalized to
3D spatial features, which can leverage all the 3D location information of
target speakers, including the azimuth, elevation, and distance respectively.
Second, to further explore the strength of 3D spatial features and simplify
the pipeline, we also design an All-In-One system that totally removes the
separation module which can reduce training and inference time by more than a
half, while still maintaining comparable results to the pipelined system.
Experiments on both simulated and real data are reported.

Fig. 1: Two paradigms for multi-channel multi-speaker ASR. a) pipelined system
with a multi-channel speech separation module followed by a single-channel ASR
module; 2) “All-In-One” model where 3D spatial features are directly used
without separation module. Note the face and body tracking algorithm provides
the location of each speaker, including the distance which is obtained from
depth or binocular cameras.

##  2 3D Spatial Feature

###  2.1 Spectral Feature

The first step is to transform the $M$-channel noisy speech $y^{1:M}$ from
time domain to frequency domain by short time Fourier transform (STFT). Given
a window function $w$ with length $N$, the complex spectrogram ${\bf Y}^{m}$
of the $m$-th channel can be calculated as:

${\bf Y}^{m}_{t,f}=\sum_{n=0}^{N-1}y^{m}[n]w[n-t]\cdot\text{exp}(-i\frac{2\pi
n}{N}f)$ (1)

where $t$ and $f$ denote the frame index and frequency band index,
respectively. The logarithm power spectrum (LPS) of the reference channel (the
first channel in this work) is served as the spectral feature, calculated by
$\text{LPS}=\text{log}|{\bf Y}^{m=1}|^{2}$.

###  2.2 Interchannel Phase Difference (IPD)

IPD has been shown to be a strong cue as a spatial feature in many previous
works [18, 21] for speech separation. We follow the practice in [6] of
computing the phase difference between the selected $P$ pairs of microphones
for both efficiency and better performance:

$\text{IPD}^{(m)}_{t,f}=\angle{\bf Y}^{m_{1}}_{t,f}-\angle{\bf
Y}^{m_{2}}_{t,f}$ (2)

where $\angle{\bf Y}^{m_{1}}_{t,f}$ is the angle of this complex number and
$(m)$ is the index for mic-pair $(m_{1},m_{2})$. These are concatenated to
form the IPD features:
$\textbf{IPD}=[...,\text{IPD}^{(m)}_{t,f},...,]_{P\times T\times F}$ .

###  2.3 3D Spatial Feature

Spatial Feature (SF) or angle feature (AF) previously proposed in [20], is
used to indicate the dominance of the target source in the T-F bins. We define
a target-dependent phase difference (TPD) as a phase delay caused by a
theoretical wave (with frequency $f$) that travels from target location
$\mathbb{L}_{t}$, measured at the $m$-th microphone pair at time $t$ [6].
Generally, the spatial features of a target speaker are formulated as the
cosine distance between the target-dependent phase difference (TPD) and the
interchannel phase differences (IPD):

$\mathbb{SF}_{t,f}=\sum_{p=1}^{P}\langle{\bm{e}}^{\text{TPD}_{t,f}^{(m)}(\mathbb{L})},{\bm{e}}^{\text{IPD}_{t,f}^{(m)}}\rangle\vspace{-0.1cm}$
(3)

where ${\bm{e}}^{(\cdot)}=\begin{bmatrix}\cos(\cdot)\\\
\sin(\cdot)\end{bmatrix}$ and $\text{TPD}_{t,f}^{(m)}(\mathbb{L})$ depends on
target speaker’s location $\mathbb{L}$ obtained from the vision. It measures
the phase difference at mic-pair $(m)$ if there is a wave spread from location
$\mathbb{L}$. Intuitively, if the speaker at location $\mathbb{L}$ is
speaking, this _theoretical_ phase difference
$\text{TPD}_{t,f}^{(m)}(\mathbb{L})$ from vision will be very similar to the
_observed_ phase difference $\text{IPD}_{t,f}^{(m)}$. Thus $\mathbb{SF}_{t,f}$
will be close to 1, otherwise close to 0.

Fig. 2: The illustration of 1D and 3D scenarios. Fig 2 a) shows inseparable
issue when azimuths of two simultaneous speech are close.

####  2.3.1 Target Phase Difference in 1D

In 1D cases, the planar wave propagation model was assumed and only the
directional information–azimuth $\theta_{a}$ of the target speaker is known.
Thus, the location $\mathbb{L}=[\theta_{a}]$ and the target phase difference
$\text{TPD}_{t,f}^{(m)}(\mathbb{L})$ can be obtained by:

$\text{TPD}_{t,f}^{(m)}(\theta_{a})=\frac{2\pi f}{c(F-1)}\cdot f_{s}\cdot
d_{(m)}\text{ cos }\theta_{a}\vspace{-0.1cm}$ (4)

where $d^{(m)}$ is the distance between the microphone-pair $(m)$, $f_{s}$ is
the sampling rate, $c$ is the sound velocity and $F$ is the total number of
frequency bands.

####  2.3.2 Target Phase Difference in 3D

One problem of the $\mathbb{SF}_{t,f}$ in 1D is its limited discriminative
power when any of the interference speakers shares close azimuth with the
target speaker, even they are far apart. As we observed in the empirical study
[18, 20], when the speakers are close, e.g., the azimuth difference (AD)
between speakers is less than 15 degrees, the spatial information calculated
from azimuth is not sufficient to distinguish between speakers. Thus, it is
necessary to generalize the TPD to 3D space with full location information
$\mathbb{L}=[\theta_{a},\theta_{e},d_{o}]$. We illustrate it in Fig. 2 and the
generalized $\text{TPD}_{t,f}^{(m)}(\mathbb{L})$ in 3D becomes

$\begin{split}&\text{TPD}_{t,f}^{(m)}(\theta_{a},\theta_{e},d_{o})=\frac{2\pi
f}{c(F-1)}\cdot f_{s}\cdot(d_{m_{1}}-d_{m_{2}})\\\
&d_{m_{i}}=\sqrt{d^{2}_{om_{i}}+d^{2}_{o}-2d_{om_{i}}d\cos\theta_{a}\cos\theta_{e}}\leavevmode\nobreak\
\leavevmode\nobreak\ _{\forall i\in\\{1,2\\}}\end{split}\vspace{-0.1cm}$ (5)

where $\theta_{a}$ is azimuth, $\theta_{e}$ is elevation. $d_{o}$, $d_{m_{1}}$
and $d_{m_{2}}$ are distances between the target speaker and the camera, the
$m_{1}$-th and the $m_{2}$-th microphone, respectively. $d_{om_{i}}$ is the
distance between the camera and the $m_{i}$-th microphone (as shown in Fig.
2). Here the camera is at the center of microphone array.

Fig. 3 is an example of a mixed speech of two. Comparing Fig. 3 (d) and (e)
shows that the proposed 3D $\mathbb{SF}_{t,f}$ not only exhibits the dominance
of speakers in T-F bins precisely, but also possesses stronger contrast than
its 1D counterpart.

(a) Overlap speech

(b) Target

(c) Interference

(d) 1D $\mathbb{SF}_{t,f}$

(e) 3D $\mathbb{SF}_{t,f}$

Fig. 3: An illustration of spatial features in 1D and 3D. Black solid
rectangle marks the T-F bins dominated by target speaker while the dashed
rectangle marks T-F bins dominated by interfering speaker.

##  3 Systems for multi-speaker ASR

###  3.1 Pipelined System

There are two paradigms for multi-speaker ASR — Pipelined approaches and All-
In-One approaches. As shown in Fig. 1, the pipelined system takes multi-
channel $Y^{1:M}$ along with target’s location information
$\mathbb{L}=(\theta_{a},\theta_{e},d_{o})$ as input, followed by a single-
channel ASR module.

####  3.1.1 Separation Module

Features: For short time Fourier transform (STFT) setting, we use 32ms sqrt
hann window and 16ms hop size at a sampling rate of 16kHz, resulting in a
frame length of 512 and hop size of 256. 512-point FFT is used to extract
257-dimensional LPS. 6 pairs of microphones are selected for IPD and TPD
computation, which are $(0,7),(1,6),(2,5),(3,4),(4,7),(3,4)$. The total
dimension of the input feature after concatenation is $257\times(1+6+1)=2056$.

Architecture: A Conv-TasNet is used for speech separation with the same
settings as [6]. It consists of multiple stacked dilated 1-dimensional
convolution blocks to increase the receptive field over the input sequence. We
use batch normalization to reduce the training and evaluation duration. An
iSTFT convolution 1d layer is placed at the end of the model to convert the
estimated target speaker complex spectrogram back to the single-channel
waveform.

Pretraining: We pretrain the separation module on the simulated data with
scale-invariant signal to noise ratio (Si-SNR) as a criterion and the
reverberant clean speech as supervision. The signal-to-distortion ratio (SDR)
and perceptual evaluation of speech quality (PESQ) scores computed with the
reverberant clean speech as the reference are used as the metric for speech
separation. Our best pretrained separation module achieves 12.8 dB of SDR and
3.35 of PESQ on the test set.

####  3.1.2 ASR Module

Features: Log Mel-Filterbank (LFB) acoustic features with 40 bins are used,
which are extracted with a 25ms window, 10ms hop-length at a sample rate of
16kHz.

Architecture: For better performance on overlapped speech as shown in [8, 9],
Lattice-free MMI (LFMMI) loss in k2 [22] is used as the training criterion. As
to model structure, two convolution sub-sampling layers with kernel size 3*3
and stride 2 are used in the front of the 12-layer 4-head Conformer encoder
[23]. Each conformer layer uses 384 attention dimensions and 2048 feed-forward
dimensions.

Pretraining: We pretrain the ASR module on either dry clean (without
reverberation) or reverberant clean data. Adam optimizer is used with a
learning rate schedule with 5,000 warm-up steps. Table 1 shows the character
error rates (CERs) we get from different training set. Note that ”separated”
represents the output from the best pretrained separation module we trained.
Since we are not doing any dereverberation in this system, separated speech is
also reverberant. In this case, the best CER on reverberant clean test data
(8.66%) can be treated as the upper bound of the ASR module when given perfect
separated speech. And the best CER on separated speech (16.41%) can be
interpreted as the lower bound, which shows how the pipelined system performs
when we simply connect two pretrained modules.

####  3.1.3 Joint Training

It is exhibited in previous works that end-to-end joint training of separation
module and ASR module using ASR criterion (LFMMI loss in our case) can bring
further improvement to the pipelined system [8, 9, 24]. We compare several
different training strategies as shown in Table 1. We notice that when there
is no pretraining (System 1), the pipelined system fails to converge. The best
practice is to start from pretrained models, and do joint finetuning on both
modules using mixed speech (System $4$). It gives the best CER of 10.50% on
simulated test and 5.03% on the real recordings, which is close to the upper
bound of 8.66% in Table 2.

sys |  | Pipelined System
---
Separation module (SISNR loss)$\rightarrow|\rightarrow$ASR module (LFMMI loss)
| All-in-One System
---
(LFMMI loss)
|  |  | SDR of | CER(%)
ID | Inputs |  | pretrain
---
(overlap)
| joint
---
finetune
Inputs |  | pretrain
---
(nonoverlap)
| joint
---
finetune
| Inputs
---
# param (M) |  Training time/iter (s) |  Decoding RTF |  Separation module | simulated | real
1 | LPS+IPD+3D | ✗ | ✓ | LFB | ✗ | ✓ |  | 45.4 | 2.35 | 0.0264 | – | not converge
2 | LPS+IPD+3D | ✓ | ✗ | LFB | ✓ | ✗ |  | 45.4 | 2.35 | 0.0264 | 12.8 dB | 16.41 | 8.03
3 | LPS+IPD+3D | ✓ | ✗ | LFB | ✓ | ✓ |  | 45.4 | 2.35 | 0.0264 | 12.8 dB | 14.65 | 7.08
4 | LPS+IPD+3D | ✓ | ✓ | LFB | ✓ | ✓ | N/A | 45.4 | 2.35 | 0.0264 | 11.0 dB | 10.50 | 5.03
5 | LPS+IPD+1D | ✓ | ✓ | LFB | ✓ | ✓ |  | 45.4 | 2.31 | 0.0261 |  9.7 dB | 17.16 | 9.15
6 | LFB+3D | ✓ | ✓ | LFB | ✓ | ✓ |  | 44.3 | 2.29 | 0.0252 | 10.0 dB | 15.03 | 7.74
7 |  |  |  |  |  |  | LPS+IPD+3D | 37.5 | 0.92 | 0.0157 |  | 18.24 | 8.16
8 |  |  |  |  |  |  | LFB+IPD+3D | 37.2 | 0.91 | 0.0157 |  | 17.51 | 7.85
9 | N/A | LFB+3D | 36.4 | 0.78 | 0.0121 | N/A | 14.37 | 6.35
10 |  |  |  |  |  |  | LFB+1D | 36.4 | 0.78 | 0.0121 |  | 22.75 | 9.22

Table 1: Comparing pipelined systems and All-In-One systems, in terms of
decoding speed and the CER (%) on both Aishell-1-mixed simulated test and real
recordings. The joint finetuning of the pipelined system uses LFMMI criteria
for both separation and ASR modules. Decoding RTF was tested on the GPU
(Nvidia V100, batch size=1) with perfetched input and no model quantization.
Note, the all-in-one system 9 performs better than the pipelined system 3 in
CER and with only half the decoding time. Another advantage of the system 9 is
it uses LFB and 3D spatial features only, which makes the all-in-one model
independent to microphone arrays typologies.

###  3.2 All-In-One System

Motivation: Although pipelined systems achieve great results with 3D spatial
feature, there still exists several limitations: 1) Requirement of parallel
reverberant clean data, 2) mismatch of training scale between the separation
module (chunk) and ASR module (utterance) and c) lack of usage of contextual
and semantic information. For these reasons, we believe it is beneficial to
simplify the pipelined system into a ”All-In-One” model that removes the
separation module.

Architecture: In the All-In-One system, to focus on the 3D spatial feature we
propose, we further reduce the input features to only contain $\mathbb{SF}$
and 40-dimension LFB. We unify the window length and hop size to be 25ms and
10ms, resulting in a 400-point FFT and 201 dimensions for $\mathbb{SF}$. We
simply concatenate 2 features and the final input feature dimension is
$40+201=241$. One advantage of this setup is all the spatial and multi-channel
information have been baked into $\mathbb{SF}$, thus the model becomes
independent to the microphone array typologies. The model structure is the
same 12-layer 4-head conformer in the ASR module described above.

Train Test CER(%) |  | dry
---
clean
| reverberant
---
clean
| overlapped
---
speech
| separated
---
speech
dry clean |  7.00111We achieved CER=6.3% with model averaging [25] during decoding which is comparable to the state-of-the-art results in [26]. For simplicity we did not apply model averages for the rest of experiments. | 13.05 | 101.03 | 26.92
reverberant clean | 7.06 | 8.66 | 100.54 | 16.41
Table 2: CERs (%) on different Aishell-1 test set trained with dry/reverberant
clean speech.

##  4 Experiments

As there is no existing dataset containing overlapped speech with 3D location
information, we prepare the training data by simulation. The open-source
Chinese Mandarin speech corpus AISHELL-1 [27] is selected as the clean data,
which contains a 150-hr training set, a 10-hr development set, and a 5-hr test
set. For simulated data, 100,000 sets of room impose response (RIRs) are
generated. The room size is ranging from [3, 3, 3] to [10, 8, 5] meters. An
8-element non-uniform linear array with spacing 15-10-5-20-5-10-15 cm is used
as the recorded sensors. The multi-channel signals are generated using the
image-source method (ISM) [28]. The reverberation time T60 is ranging 0.05s to
0.7 seconds. The signal-to-interference rate (SIR) is ranging from -6 to 6 dB.
The sampling rate is 16kHz. Each utterance is simulated 3 times, mixed with
random interfering utterances from the same dataset split.

To evaluate the proposed 3D spatial features in real-life applications, the
systems were also tested on our in-car data which contains 2000 utterances
recorded using our device in daily car scenarios and $10000$-hr simulated data
is used to train the in-car models. The target speech of in-car tests is the
general voice commands of controlling the car, such as ”turn on the air
conditioner” or ”open the back sit window”. The interference speech is
replayed from loudspeakers. There are 4 potential speakers in the car: the
main driver, the co-driver and two passengers sitting in the back. The main
driver’s voice is taken as the target. Note the azimuths of the main driver
and the main-back passenger (loudspeaker) are very close, though their
elevations and distances are far apart. In our experiments on simulated data,
we use the oracle location information in both training and test time. On car
datasets, this information can be obtained by face detection and tracking
techniques with depth cameras [29].

A complete comparison between the pipelined systems and the All-In-One system
is shown in Table 2. Systems 1-6 are pipelined systems and Systems 7-10 are
All-In-One systems. When compared with the pipelined systems, although the
All-In-One systems still fall behind in CER, the training time and decoding
time is speedup by $3\times$ and $2\times$ respectively. In addition, All-In-
One systems are trained from scratch which further reduces the effort on
pretraining. Thus we believe when empowered by well-designed features like 3D
spatial features, the All-In-One system has good potential and is a promising
direction for further study.

##  5 Conclusions

In this work, we study the application of 3D spatial features for multi-
channel multi-speaker overlapped speech recognition. It shows great
improvement over the previous 1D feature and achieves 31+% relative CERR on
both simulated and real data. We also explore the idea of the All-In-One model
with the proposed 3D feature as input and obtain comparable results. We
believe when paid more effort and attention, the All-In-One model has the
potential to surpass the prevalent pipelined systems in terms of both
efficiency and accuracy.

## References

  * [1] Geoffrey Hinton, Li Deng, Dong Yu, George E Dahl, Abdel-rahman Mohamed, Navdeep Jaitly, Andrew Senior, Vincent Vanhoucke, Patrick Nguyen, Tara N Sainath, et al.,  “Deep neural networks for acoustic modeling in speech recognition: The shared views of four research groups,”  IEEE Signal processing magazine, vol. 29, no. 6, pp. 82–97, 2012\.
  * [2] Yiming Wang, Tongfei Chen, Hainan Xu, Shuoyang Ding, Hang Lv, Yiwen Shao, Nanyun Peng, Lei Xie, Shinji Watanabe, and Sanjeev Khudanpur,  “Espresso: A fast end-to-end neural speech recognition toolkit,”  in 2019 IEEE Automatic Speech Recognition and Understanding Workshop (ASRU). IEEE, 2019, pp. 136–143.
  * [3] Yiwen Shao, Yiming Wang, Daniel Povey, and Sanjeev Khudanpur,  “Pychain: A fully parallelized pytorch implementation of lf-mmi for end-to-end asr,”  arXiv preprint arXiv:2005.09824, 2020.
  * [4] Yu Zhang, James Qin, Daniel S Park, Wei Han, Chung-Cheng Chiu, Ruoming Pang, Quoc V Le, and Yonghui Wu,  “Pushing the limits of semi-supervised learning for automatic speech recognition,”  arXiv preprint arXiv:2010.10504, 2020.
  * [5] Aswin Shanmugam Subramanian, Chao Weng, Meng Yu, Shi-Xiong Zhang, Yong Xu, Shinji Watanabe, and Dong Yu,  “Far-field location guided target speech extraction using end-to-end speech recognition objectives,”  in ICASSP 2020-2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2020, pp. 7299–7303.
  * [6] Rongzhi Gu, Shi-Xiong Zhang, Yong Xu, Lianwu Chen, Yuexian Zou, and Dong Yu,  “Multi-modal multi-channel target speech separation,”  IEEE Journal of Selected Topics in Signal Processing, vol. 14, no. 3, pp. 530–541, 2020.
  * [7] Naoyuki Kanda, Xuankai Chang, Yashesh Gaur, Xiaofei Wang, Zhong Meng, Zhuo Chen, and Takuya Yoshioka,  “Investigation of end-to-end speaker-attributed asr for continuous multi-talker recordings,”  in 2021 IEEE Spoken Language Technology Workshop (SLT). IEEE, 2021, pp. 809–816.
  * [8] Jianwei Yu, Shi-Xiong Zhang, Jian Wu, Shahram Ghorbani, Bo Wu, Shiyin Kang, Shansong Liu, Xunying Liu, Helen Meng, and Dong Yu,  “Audio-visual recognition of overlapped speech for the lrs2 dataset,”  in ICASSP 2020-2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2020, pp. 6984–6988.
  * [9] Jianwei Yu, Shi-Xiong Zhang, Bo Wu, Shansong Liu, Shoukang Hu, Xunying Liu, Helen M Meng, and Dong Yu,  “Audio-visual multi-channel integration and recognition of overlapped speech,”  IEEE/ACM Transactions on Audio, Speech, and Language Processing, 2021.
  * [10] Aswin Shanmugam Subramanian, Chao Weng, Shinji Watanabe, Meng Yu, Yong Xu, Shi-Xiong Zhang, and Dong Yu,  “Directional asr: A new paradigm for e2e multi-speaker speech recognition with source localization,”  in ICASSP 2021-2021 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2021, pp. 8433–8437.
  * [11] Yuxuan Wang, Arun Narayanan, and DeLiang Wang,  “On training targets for supervised speech separation,”  IEEE/ACM transactions on audio, speech, and language processing, vol. 22, no. 12, pp. 1849–1858, 2014.
  * [12] DeLiang Wang and Jitong Chen,  “Supervised speech separation based on deep learning: An overview,”  IEEE/ACM Transactions on Audio, Speech, and Language Processing, vol. 26, no. 10, pp. 1702–1726, 2018.
  * [13] John R Hershey, Zhuo Chen, Jonathan Le Roux, and Shinji Watanabe,  “Deep clustering: Discriminative embeddings for segmentation and separation,”  in 2016 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2016, pp. 31–35.
  * [14] Kateřina Žmolíková, Marc Delcroix, Keisuke Kinoshita, Tsubasa Ochiai, Tomohiro Nakatani, Lukáš Burget, and Jan Černockỳ,  “Speakerbeam: Speaker aware neural network for target speaker extraction in speech mixtures,”  IEEE Journal of Selected Topics in Signal Processing, vol. 13, no. 4, pp. 800–814, 2019.
  * [15] Quan Wang, Hannah Muckenhirn, Kevin Wilson, Prashant Sridhar, Zelin Wu, John R Hershey, Rif A Saurous, Ron J Weiss, Ye Jia, and Ignacio Lopez Moreno,  “Voicefilter: Targeted voice separation by speaker-conditioned spectrogram masking,”  Proc. Interspeech, pp. 2728–2732, 2019.
  * [16] Ariel Ephrat, Inbar Mosseri, Oran Lang, Tali Dekel, Kevin Wilson, Avinatan Hassidim, William T Freeman, and Michael Rubinstein,  “Looking to listen at the cocktail party: A speaker-independent audio-visual model for speech separation,”  arXiv preprint arXiv:1804.03619, 2018.
  * [17] Jian Wu, Yong Xu, Shi-Xiong Zhang, Lian-Wu Chen, Meng Yu, Lei Xie, and Dong Yu,  “Time domain audio visual speech separation,”  in IEEE Automatic Speech Recognition and Understanding Workshop (ASRU). IEEE, 2019, pp. 667–673.
  * [18] Zhuo Chen, Xiong Xiao, Takuya Yoshioka, Hakan Erdogan, Jinyu Li, and Yifan Gong,  “Multi-channel overlapped speech recognition with location guided speech extraction network,”  in 2018 IEEE Spoken Language Technology Workshop (SLT). IEEE, 2018, pp. 558–565.
  * [19] Zhong-Qiu Wang and DeLiang Wang,  “Combining spectral and spatial features for deep learning based blind speaker separation,”  IEEE/ACM Transactions on Audio, Speech, and Language Processing, vol. 27, no. 2, pp. 457–468, 2018.
  * [20] Rongzhi Gu, Lianwu Chen, Shi-Xiong Zhang, Jimeng Zheng, Yong Xu, Meng Yu, Dan Su, Yuexian Zou, and Dong Yu,  “Neural spatial filter: Target speaker speech separation assisted with directional information,”  Proc. Interspeech, pp. 4290–4294, 2019.
  * [21] Zhong-Qiu Wang, Jonathan Le Roux, and John R Hershey,  “Multi-channel deep clustering: Discriminative spectral and spatial embeddings for speaker-independent speech separation,”  in 2018 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2018, pp. 1–5.
  * [22] Daniel Povey et al.,  “k2,” [Online] Available: https://github.com/k2-fsa/k2.
  * [23] Anmol Gulati, James Qin, Chung-Cheng Chiu, Niki Parmar, Yu Zhang, Jiahui Yu, Wei Han, Shibo Wang, Zhengdong Zhang, Yonghui Wu, et al.,  “Conformer: Convolution-augmented transformer for speech recognition,”  arXiv preprint arXiv:2005.08100, 2020.
  * [24] Xuankai Chang, Wangyou Zhang, Yanmin Qian, Jonathan Le Roux, and Shinji Watanabe,  “Mimo-speech: End-to-end multi-channel multi-speaker speech recognition,”  in 2019 IEEE Automatic Speech Recognition and Understanding Workshop (ASRU). IEEE, 2019, pp. 237–244.
  * [25] Pavel Izmailov, Dmitrii Podoprikhin, Timur Garipov, Dmitry Vetrov, and Andrew Gordon Wilson,  “Averaging weights leads to wider optima and better generalization,”  arXiv preprint arXiv:1803.05407, 2018.
  * [26] Shinji Watanabe, Takaaki Hori, Shigeki Karita, Tomoki Hayashi, Jiro Nishitoba, Yuya Unno, Nelson Enrique Yalta Soplin, Jahn Heymann, Matthew Wiesner, Nanxin Chen, et al.,  “Espnet: End-to-end speech processing toolkit,”  arXiv preprint arXiv:1804.00015, 2018.
  * [27] Hui Bu, Jiayu Du, Xingyu Na, Bengu Wu, and Hao Zheng,  “Aishell-1: An open-source mandarin speech corpus and a speech recognition baseline,”  in 2017 20th Conference of the Oriental Chapter of the International Coordinating Committee on Speech Databases and Speech I/O Systems and Assessment (O-COCOSDA). IEEE, 2017, pp. 1–5.
  * [28] Eric A Lehmann and Anders M Johansson,  “Prediction of energy decay in room impulse responses simulated with an image-source model,”  The Journal of the Acoustical Society of America, vol. 124, no. 1, pp. 269–277, 2008.
  * [29] Davis E King,  “Dlib-ml: A machine learning toolkit,”  The Journal of Machine Learning Research, vol. 10, pp. 1755–1758, 2009.