###### proof of Theorem C.16.

From Theorem C.15 and the assumption (131), we have

$\displaystyle\widetilde{P}$ $\displaystyle=$ $\displaystyle
P(\mathbf{Y}|do(\mathbf{X}),\mathbf{Z})\times P(\mathbf{X_{1}})\times
P(\mathbf{X_{2}})\times\cdots\times P(\mathbf{X_{n}})\times P(\mathbf{Z})$
$\displaystyle=$ $\displaystyle P(\mathbf{Y}|\mathbf{X},\mathbf{Z})\times
P(\mathbf{X_{1}})\times P(\mathbf{X_{2}})\times\cdots\times
P(\mathbf{X_{n}})\times P(\mathbf{Z}).$

Thus, from Lemma C.2, we obtain

$e^{-T^{*}(\mathbf{X}_{1},\mathbf{X}_{2},\dots,\mathbf{X}_{n},\mathbf{Z})}=\frac{d\mathbb{Q}}{d\mathbb{P}}=\frac{d\widetilde{P}}{dP}.$
(134)

This completes the proof. ∎

###  C.3 Proofs for Section 7

In this section, we first present a proposition for obtaining the density
ratio between empirical distributions of the source and target distributions.
Next, we present a proposition and lemmas for the early stopping method
proposed in this study.

######  Proposition C.17.

It holds that

$\frac{d\hat{Q}^{(N)}}{d\hat{P}^{(N)}}(\mathbf{x})=\left\\{\begin{array}[]{ll}dQ/dP(\mathbf{x})&\text{if}\
\ 1\leq\exists i\leq N\ \text{s.t.}\ \mathbf{X}_{i}=\mathbf{x},\\\
0&\text{otherwise.}\end{array}\right.$ (135)

###### proof of Proposition C.24.

Let $\nu$ be the countable measure on
$\left\\{\mathbf{X}_{1},\mathbf{X}_{2},\ldots,\mathbf{X}_{N}\right\\}$:

$\nu(\mathbf{x})=\left\\{\begin{array}[]{ll}1&\text{if}\ \ 1\leq\exists i\leq
N\ \text{s.t.}\ \mathbf{X}_{i}=\mathbf{x},\\\
0&\text{otherwise.}\end{array}\right.$ (136)

Then, $\hat{P}^{(N)}\ll\nu$ and $\hat{Q}^{(N)}\ll\nu$.

Note that, from the definitions of $\hat{P}^{(N)}(\mathbf{x})$ and
$\hat{Q}^{(N)}(\mathbf{x})$, we have

$\hat{P}^{(N)}(\mathbf{x})=\frac{1}{N}\sum_{i}\mathbf{1}(\mathbf{X}^{i}_{\sim
P}=\mathbf{x})=\frac{1}{N}\sum_{i}\mathbf{1}(\mathbf{X}^{i}=\mathbf{x})\cdot\frac{dP}{d\mu}(\mathbf{x}),$
(137)

and

$\hat{Q}^{(N)}(\mathbf{x})=\frac{1}{N}\sum_{i}\mathbf{1}(\mathbf{X}^{i}_{\sim
Q}=\mathbf{x})=\frac{1}{N}\sum_{i}\mathbf{1}(\mathbf{X}^{i}=\mathbf{x})\cdot\frac{dQ}{d\mu}(\mathbf{x}),$
(138)

where $\mathbf{1}(\cdot)$ equals one if the statement in parentheses is true
and zero otherwise.

From (137) and (138), if $\mathbf{X}_{i}=\mathbf{x}$, we see

$\frac{d\hat{P}^{(N)}}{d\nu}(\mathbf{x})=\hat{P}^{(N)}(\mathbf{x})=\frac{1}{N}\frac{dP}{d\mu}(\mathbf{x}),$
(139)

and

$\frac{d\hat{Q}^{(N)}}{d\nu}(\mathbf{x})=\hat{Q}^{(N)}(\mathbf{x})=\frac{1}{N}\frac{dQ}{d\mu}(\mathbf{x}).$
(140)

Then, we have

$\frac{d\hat{Q}^{(N)}}{d\hat{P}^{(N)}}(\mathbf{x})=\frac{\frac{d\hat{P}^{(N)}}{d\nu}(\mathbf{x})}{\frac{d\hat{Q}^{(N)}}{d\nu}(\mathbf{x})}=\frac{dQ}{dP}(\mathbf{x}).$
(141)

For
$\mathbf{x}\notin\left\\{\mathbf{X}_{1},\mathbf{X}_{2},\ldots,\mathbf{X}_{N}\right\\}$,
we observe $d\hat{Q}^{(N)}/d\nu(\mathbf{x})=0$. Note that,
$d\hat{Q}^{(N)}/d\hat{P}^{(N)}(\mathbf{x})$ is defined as zero for
$\mathbf{x}\in\Omega$ such that $d\hat{Q}^{(N)}/d\nu(\mathbf{x})=0$.
Subsequently, we see $d\hat{Q}^{(N)}/d\hat{P}^{(N)}(\mathbf{x})=0$. ∎

Next, we present a proposition for the early stopping method proposed in
Section 7.1. We obtain an early stopping step as the step that minimizes the
$W_{1}$ distance of the balanced distribution and target distribution,
$\hat{Q}^{(N)}_{k}$ and $Q$ in (22). To obtain the early stopping step, we
assume that the two distributions differ the worst outside the neighborhood of
the observations because we cannot know the closeness of the two
distributions, $\hat{Q}^{(N)}_{k}$ and $Q$ in (22) except in the neighborhood
of the observations.

We now provide a note on the convergence rate for optimizing the loss function
(16). Let

$f^{(N)}(\mathbf{t})=f^{(N)}(t_{1},t_{2},\ldots,t_{N})=\frac{1}{\alpha}\frac{1}{N}\sum_{i=1}^{N}e^{\alpha\cdot
t_{i}}\cdot\frac{dQ}{d\mu}(\mathbf{x}_{i})+\frac{1}{1-\alpha}\frac{1}{N}\sum_{i=1}^{N}e^{(\alpha-1)\cdot
t_{i}}\cdot\frac{dP}{d\mu}(\mathbf{x}_{i}).$ (142)

Subsequently, let $\mathbf{t}_{K}$ denote a model at step $K$ when optimizing
(142) with a Stochastic Gradient Desent (SGD) algorithm. Because, from
Corollary C.6, $f^{(N)}(\mathbf{t})$ is strongly convex with $\left\|\nabla
f^{(N)}(\mathbf{t})\right\|^{2}\ \leq D^{2}$ ($\exists D\in\mathbb{R}$) around
the optimal point $\mathbf{t}_{*}=\allowbreak(t_{*}^{1},\allowbreak
t_{*}^{2},\allowbreak\ldots,\allowbreak
t_{*}^{N})=\allowbreak(-\log\frac{dQ}{dP}(\mathbf{x}_{1}),\allowbreak-\log\frac{dQ}{dP}(\mathbf{x}_{2}),\allowbreak\ldots,\allowbreak-\log\frac{dQ}{dP}(\mathbf{x}_{N}))$,
an $O(1/K)$ convergence rate can be achieved at step $K$ when optimizing (142)
with SGD algorithms under regular conditions for $\mathbf{t}$:

$E\left[f_{N}(\bar{\mathbf{t}}_{K})\right]-f_{N}(\mathbf{t}_{*})\leq\frac{C}{K+1},$
(143)

where $\bar{\mathbf{t}}_{K}$ is a weighted averaging such that
$\bar{\mathbf{t}}_{K}=\frac{1}{(K+1)\cdot(K+2)}\sum_{k}(k+1)\cdot\mathbf{t}_{k}$
and $C>0$ is constant. Here $E\left[\cdot\right]$ denotes the expectation for
the randomness of batch sampling of SGD.777For the convergence rate of SGD
algorithms, for example, readers can refer to [12]. As assumptions close to
(143), we briefly assume (144) in Assumption E1 and (145) in Assumption E2 to
obtain an early stopping step, which are simpler and more relaxed than (143).

Herein, we make the following assumptions for the early stopping method
presented in Section 7.

  * •

Assumption E1. Let $\\{T_{k}^{(N)}\\}_{k=1}^{\infty}$ be a sequence of
functions in $\mathcal{T}^{\alpha}$ such that
$\lim_{k\rightarrow\infty}\allowbreak
T_{k}^{(N)}(\mathbf{X}_{i})\allowbreak=\allowbreak-\log
dQ/dP(\mathbf{X}_{i})$, for $1\allowbreak\leq\allowbreak\forall\allowbreak
i\allowbreak\leq\allowbreak N$. Suppose that

$\hat{L}_{\alpha}^{(N)}(Q,P;T_{k}^{(N)})-\hat{L}_{\alpha}^{(N)}(Q,P;T_{*})\leq\frac{C_{0}}{K},$
(144)

where $\hat{L}_{\alpha}^{(N)}(Q,P;\cdot)$ is defined as (60) in Lemma C.7 and
$C_{0}>0$ is constant.

  * •

Assumption E2. Let $\\{T_{k}\\}_{k=1}^{\infty}$ be a sequence of functions in
$\mathcal{T}^{\alpha}$ with $E_{P}[e^{-T_{k}(\mathbf{X})}]=1$ such that
$\lim_{k\rightarrow\infty}T_{k}=-\log dQ/dP$, $P$-almost everywhere. Suppose
that

$L_{\alpha}(Q,P;T_{k})-L_{\alpha}(Q,P;T_{*})\leq\frac{C_{1}}{K},$ (145)

where $L_{\alpha}(Q,P;\cdot)$ is defined as (57) in Lemma C.7 and $C_{1}>0$ is
constant.

In addition, we make the following assumptions to simplify the discussion in
the proofs.

  * •

Assumption E3. Let $\Omega$ be a compact set in $\mathbb{R}^{d}$ with
$\mathrm{diam}(\Omega)=1$, where $\mathrm{diam}(\Omega)$ denotes the diameter
of $\Omega$. Then $\lambda$ denotes the Lebesgue measure on $\mathbb{R}^{d}$.

  * •

Assumption E4. Let $Q$ and $P$ be two probabilities on $\Omega$ with
continuous probability densities $p(\mathbf{x})$ and $q(\mathbf{x})$,
respectively. Assume $0<p_{min}\leq p(\mathbf{x})\leq p_{max}$ and
$0<q_{min}\leq q(\mathbf{x})\leq q_{max}$ for all $\mathbf{x}\in\Omega$.

  * •

Assumption E5. For $\\{T_{k}^{(N)}\\}_{k=1}^{\infty}$ in Assumption E1, assume
that each function of $T_{k}^{(N)}(\mathbf{X})$ is Lipschitz continuous: for
$1\leq k\leq\infty$,

$|T_{k}^{(N)}(\mathbf{x})-T_{k}^{(N)}(\mathbf{y})|\leq\rho_{k}\cdot\|\mathbf{x}-\mathbf{y}\|.$
(146)

  * •

Assumption E6. For $\\{T_{k}^{(N)}\\}_{k=1}^{\infty}$ in Assumption E2, assume
that each function of $T_{k}(\mathbf{X})$ is Lipschitz continuous: for $1\leq
k\leq\infty$,

$|T_{k}(\mathbf{x})-T_{k}(\mathbf{y})|\leq\widetilde{\rho}_{k}\cdot\|\mathbf{x}-\mathbf{y}\|.$
(147)

Note that, the Lipschitz coefficient in Assumption E5 does not depend on the
sample size $N$.

######  Lemma C.18.

For $\\{T_{k}^{(N)}\\}_{k=1}^{\infty}$ in Assumption E5, it holds that for
$\mathbf{x}\in\Omega$ and $\|\mathbf{y}-\mathbf{x}\|<D$,

$e^{-T_{k}^{(N)}(\mathbf{y})}=e^{-T_{k}^{(N)}(\mathbf{x})}+e^{-T_{k}^{(N)}(\mathbf{x})}\cdot\left\\{O\left(D\right)+O\left(D^{2}\right)\right\\}+O_{\mathbf{x}}\left(D\right).$
(148)

###### proof of Lemma C.18.

From the intermediate value theorem for the second derivative of $e^{-x}$, we
have

$e^{-y}=e^{-x}-e^{-x}\cdot(y-x)+\frac{e^{-x+\theta\cdot(x-y)}}{2}\cdot(y-x)^{2},$
(149)

where $0<\theta<1$.

By substituting $y=T_{k}^{(N)}(\mathbf{y})$ and $x=T_{k}^{(N)}(\mathbf{x})$
into the above formula, we obtain

$\displaystyle e^{-T_{k}^{(N)}(\mathbf{y})}$ $\displaystyle=$ $\displaystyle
e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{x})}\left(T_{k}^{(N)}(\mathbf{y})-T_{k}^{(N)}(\mathbf{x})\right)$
(150) $\displaystyle\ +\
\frac{e^{-T_{k}^{(N)}(\mathbf{x})+\theta(\mathbf{x},\mathbf{y})\cdot\left(T_{k}^{(N)}(\mathbf{y})-T_{k}^{(N)}(\mathbf{x})\right)}}{2}\left(T_{k}^{(N)}(\mathbf{y})-T_{k}^{(N)}(\mathbf{x})\right)^{2},$

where $0<\theta(\mathbf{x},\mathbf{y})<1$.

Now, note that, from Assumption E5,

$T_{k}^{(N)}(\mathbf{y})-T_{k}^{(N)}(\mathbf{x})\leq\rho_{k}\cdot\|\mathbf{x}-\mathbf{y}\|\leq\rho_{k}\cdot
D.$ (151)

From (150) and (151), we see

$\displaystyle\left|e^{-T_{k}^{(N)}(\mathbf{y})}-e^{-T_{k}^{(N)}(\mathbf{x})}\right|$
$\displaystyle=$
$\displaystyle\left|e^{-T_{k}^{(N)}(\mathbf{x})}\cdot\left(T_{k}^{(N)}(\mathbf{y})-T_{k}^{(N)}(\mathbf{x})\right)\right.$
(152) $\displaystyle\left.\ \ +\
\frac{e^{-T_{k}^{(N)}(\mathbf{x})+\theta(\mathbf{x},\mathbf{y})\cdot(T_{k}^{(N)}\left(\mathbf{y})-T_{k}^{(N)}(\mathbf{x})\right)}}{2}\cdot\left(T_{k}^{(N)}(\mathbf{y})-T_{k}^{(N)}(\mathbf{x})\right)^{2}\right|$
$\displaystyle\leq$ $\displaystyle
e^{-T_{k}^{(N)}(\mathbf{x})}\cdot\rho_{k}\cdot
D+\frac{e^{-T_{k}^{(N)}(\mathbf{x})+D}}{2}\cdot D^{2}$ $\displaystyle=$
$\displaystyle
e^{-T_{k}^{(N)}(\mathbf{x})}\cdot\left\\{O\left(D\right)+O\left(D^{2}\right)\right\\}+O_{\mathbf{x}}\left(D\right).$

Therefore, we have

$e^{-T_{k}^{(N)}(\mathbf{y})}=e^{-T_{k}^{(N)}(\mathbf{x})}+e^{-T_{k}^{(N)}(\mathbf{x})}\cdot\left\\{O\left(D\right)+O\left(D^{2}\right)\right\\}+O_{\mathbf{x}}\left(D\right).$
(153)

This completes the proof. ∎

######  Lemma C.19.

For $\\{T_{k}\\}_{k=1}^{\infty}$ in Assumption E6, it holds that for
$\mathbf{x}\in\Omega$ and $\|\mathbf{y}-\mathbf{x}\|<D$,

$e^{-T_{k}(\mathbf{y})}=e^{-T_{k}}(\mathbf{x})+e^{-T_{k}(\mathbf{x})}\cdot\left\\{O\left(D\right)+O\left(D^{2}\right)\right\\}+O_{\mathbf{x}}\left(D\right).$
(154)

###### proof of Lemma C.19.

Note that, we use only the Lipschitz continuity of $T_{k}^{(N)}(\mathbf{x})$
to prove Lemma C.18. From Assumption E6, $T_{k}(\mathbf{x})$ is Lipschitz
continuous. Then, (154) can be proven in a manner similar to Lemma C.18.

This completes the proof. ∎

######  Lemma C.20.

Let $T_{*}=-\log dQ/dP$. Under Assumption E3 and E4, it holds that for
$\mathbf{x}\in\Omega$ and $\|\mathbf{y}-\mathbf{x}\|<D$,

$e^{-T_{*}}(\mathbf{y})=e^{-T_{*}(\mathbf{x})}-e^{-T_{*}(\mathbf{x})}\cdot\left\\{O\left(D\right)+O\left(D^{2}\right)\right\\}+O_{\mathbf{x}}\left(D\right).$
(155)

###### proof of Lemma C.20.

Note that, we use only the Lipschitz continuity of $T_{k}^{(N)}(\mathbf{x})$
to prove Lemma C.18. From Assumption E3 and Assumption E4, $T_{*}(\mathbf{x})$
is a bounded continuous function on $\Omega$. Since bounded continuous
functions are Lipschitz continuous, $T_{*}(\mathbf{x})$ is Lipschitz
continuous. Thus, (155) can be proven in a manner similar to Lemma C.18.

This completes the proof. ∎

######  Lemma C.21.

Let $B(\mathbf{x}_{0},D)=\\{\mathbf{y}:\|\mathbf{y}-\mathbf{x}_{0}\|<D\\}$.
Then,

$\displaystyle p_{min}\cdot D^{d}\leq P(B(\mathbf{x}_{0},D))\leq p_{max}\cdot
D^{d}.$ (156)

###### proof of Lemma C.21.

From Assumption E4, $p_{min}\leq p(\mathbf{x})\leq p_{max}$ holds, and by
integrating over $B(\mathbf{x}_{0},D)$ with $\lambda$, we obtain

$\displaystyle\int_{B(\mathbf{x}_{0},D)}p_{min}\ d\lambda\leq$ $\displaystyle
P(B(\mathbf{x}_{0},D))$ $\displaystyle\leq\int_{B(\mathbf{x}_{0},D)}p_{max},\
d\lambda$ $\displaystyle\therefore\ p_{min}\cdot D^{d}\leq$ $\displaystyle
P(B(\mathbf{x}_{0},D))$ $\displaystyle\leq p_{max}\cdot D^{d}.$ (157)

This completes the proof. ∎

######  Lemma C.22.

Let $B(\mathbf{x}_{0},D)=\\{\mathbf{y}:\|\mathbf{y}-\mathbf{x}_{0}\|<D\\}$.
Then,

$\displaystyle P(B(\mathbf{x}_{0},D))=C\cdot p(\mathbf{x}_{0})\cdot D^{d},$
(158)

where $C$ is constant.

###### proof of Lemma C.22.

From Assumption E4, $p$ is a bounded continuous function on $\Omega$. Since
bounded continuous functions are Lipschitz continuous, $p(\mathbf{x})$ is
Lipschitz continuous.

Then, there exist a constant $C$ such that

$\displaystyle p(\mathbf{x})\leq
p(\mathbf{\mathbf{x}_{0}})+C\cdot\|\mathbf{x}-\mathbf{x}_{0}\|,$ (159)

and by integrating over $B(\mathbf{x}_{0},D)$ with $\lambda$, we obtain

$\displaystyle P(B(\mathbf{x}_{0},D))\leq C\cdot p({\mathbf{x}_{0}})\cdot
D^{d}.$ (160)

This completes the proof. ∎

######  Proposition C.23.

For $\\{T_{k}^{(N)}\\}_{k=1}^{\infty}$ in Assumption E1, let
$\hat{Q}^{(N)}_{k}$ be a probability defined as

$d\hat{Q}^{(N)}_{k}=e^{-T_{k}^{(N)}}\cdot dP.$

Then, under Assumpution E1-E6, for a sufficiently large $K>0$, it holds that

$E_{\mathbf{X}_{P}^{(N)}}[W_{1}(Q,\hat{Q}^{(N)}_{K})]\leq 2-N\cdot
K^{-\frac{d}{2}}+K^{-\frac{1}{2}}.$ (161)

######  Corollary C.24.

Let $K_{0}=\allowbreak N^{\frac{2}{d+\delta}}$ with $\delta>0$. Then, under
Assumpution E1-E6, for a sufficiently large $N$, it holds that

$E_{\mathbf{X}_{P}^{(N)}}[W_{1}(Q,\hat{Q}^{(N)}_{K_{0}})]\leq
2-K_{0}^{\frac{\delta}{2}}+K_{0}^{-\frac{1}{2}}.$ (162)

######  Corollary C.25.

In Corollary C.24, let $\delta^{\prime}>0$ such that
$N^{\frac{\delta^{\prime}}{d+\delta^{\prime}}}\allowbreak=\allowbreak 2$, and
let $K_{0}\allowbreak=\allowbreak N^{\frac{2}{d+\delta^{\prime}}}$. Then,
under Assumpution E1-E6, for a sufficiently large $N$, it holds that

$E_{\mathbf{X}_{P}^{(N)}}[W_{1}(Q,\hat{Q}^{(N)}_{K_{0}})]\leq
K_{0}^{-\frac{1}{2}}.$ (163)

Thus, if $N>\left(\frac{1}{\varepsilon}\right)^{d+\delta^{\prime}}$ then
$E_{\mathbf{X}_{P}^{(N)}}[W_{1}({Q},\hat{{Q}}^{(N)}_{K_{0}})]<\varepsilon$.

###### proof of Proposition C.23.

Let $Q_{K}$ be a probability defined as

$dQ_{K}=e^{-T_{K}}\cdot dP.$ (164)

Intuitively, $Q_{K}$ is the true balanced probability distribution at a step
$K$.

First, from the triangle inequality for the $L_{1}$ norm, we have

$E_{\mu}\left|\hat{Q}^{(N)}_{K}-Q\right|\leq
E_{\mu}\Big{|}\hat{Q}^{(N)}_{K}-Q_{K}\Big{|}+E_{\mu}\Big{|}Q_{K}-Q\Big{|}.$
(165)

Considering the expectation $E_{\mathbf{X}_{P}^{(N)}}[\cdot]$ for the both
sides of the above equation, we see

$E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\left|\hat{Q}^{(N)}_{K}-Q\right|\right]\leq
E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\Big{|}\hat{Q}^{(N)}_{K}-Q_{K}\Big{|}\right]+E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\Big{|}Q_{K}-Q\Big{|}\right].$
(166)

Next, we obtain the upper bound of the first term in (166).

Let $\Delta_{i}=B(\mathbf{X}_{i},1/\sqrt{K})$. Subsequently, let
$\Delta=\bigcup_{i=1}^{N}\Delta_{i}$. Then, we have

$\displaystyle E_{\mu}\left|\hat{Q}^{(N)}_{K}-Q_{K}\right|$ $\displaystyle=$
$\displaystyle\int\left|e^{-T_{k}^{(N)}}\cdot\frac{dP}{d\mu}-e^{-T_{k}}\cdot\frac{dP}{d\mu}\right|d\mu$
$\displaystyle=$
$\displaystyle\int\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\frac{dP}{d\mu}d\mu$
$\displaystyle=$ $\displaystyle E_{P}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|$
$\displaystyle=$ $\displaystyle
E_{P}\left[id_{\Delta}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]+E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]$
$\displaystyle\leq$ $\displaystyle
E_{P}\left[id_{\Delta}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]+E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}^{(N)}}\right|\right]+E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}}\right|\right].$

Considering the expectation $E_{\mathbf{X}_{P}^{(N)}}[\cdot]$ for the both
sides of the above equation, we see

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\left|\hat{Q}^{(N)}_{K}-Q_{K}\right|\right]$
(167) $\displaystyle\leq$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Delta}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]\right]$
$\displaystyle\qquad+\
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}^{(N)}}\right|\right]\right]+E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}}\right|\right]\right].$

To obtain the upper bound of the first term in (167), we see

$\displaystyle
E_{P}\left[id_{\Delta}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]$ (168)
$\displaystyle=$ $\displaystyle
E_{P}\left[\sum_{i=1}^{N}id_{\Delta_{i}}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}E_{P}\left[id_{\Delta_{i}}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]$
$\displaystyle=$ $\displaystyle
E_{P}\left[\left|\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}\cdot
e^{-T_{k}^{(N)}}-\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}\cdot
e^{-T_{k}(\mathbf{X}_{i})}\right|\right].$

Subsequently, we have

$\displaystyle
E_{P}\left[\left|\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}\cdot
e^{-T_{k}^{(N)}}-\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}\cdot
e^{-T_{k}(\mathbf{X}_{i})}\right|\right]$ $\displaystyle=$ $\displaystyle
E_{P}\left[\left|\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}^{(N)}(\mathbf{x})}-\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right.\right.$ $\displaystyle\qquad\qquad+\
\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{*}(\mathbf{X}_{i})}$ $\displaystyle\qquad\qquad+\
\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{*}(\mathbf{X}_{i})}-\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{*}(\mathbf{x})}$
$\displaystyle\qquad\qquad+\left.\left.\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{*}(\mathbf{x})}-\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}(\mathbf{x})}\right|\right]$ $\displaystyle=$ $\displaystyle
E_{P}\left[\left|\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left\\{e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right\\}\right.\right.$
$\displaystyle\qquad\qquad+\
\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right\\}$
$\displaystyle\qquad\qquad+\
\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}\left\\{e^{-T_{*}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{x})}\right\\}$
$\displaystyle\qquad\qquad+\
\left.\left.\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left\\{e^{-T_{*}(\mathbf{x})}-e^{-T_{k}(\mathbf{x})}\right\\}\right|\right]$
$\displaystyle\leq$ $\displaystyle
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right|\right]$
$\displaystyle\qquad\qquad+\
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|\right]$
$\displaystyle\qquad\qquad+\
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{*}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{x})}\right|\right]$
$\displaystyle\qquad\qquad+E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{*}(\mathbf{x})}-e^{-T_{k}(\mathbf{x})}\right|\right].$

Considering the expectation $E_{\mathbf{X}_{P}^{(N)}}[\cdot]$ for the both
sides of the above equation, we obtain

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Delta}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]\right]$
(169) $\displaystyle\leq$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right|\right]\right]$
$\displaystyle\qquad\qquad+\ E_{\mathbf{X}_{P}^{(N)}}\left[\
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|\right]\right]$
$\displaystyle\qquad\qquad+\ E_{\mathbf{X}_{P}^{(N)}}\left[\
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{*}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{x})}\right|\right]\right]$
$\displaystyle\qquad\qquad+\
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{*}(\mathbf{x})}-e^{-T_{k}(\mathbf{x})}\right|\right]\right].$

Now, from Lemma C.19, we have, for $\mathbf{x}\in
id_{B(\mathbf{X}_{i},1/\sqrt{K})}$

$\left|e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right|=e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{X}_{i}}\left(\frac{1}{K}\right).\\\
$

Then, we see

$\displaystyle
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right|\right]$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{X}_{i}}\left(\frac{1}{\sqrt{K}}\right)\right\\}\right]$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\right]\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{X}_{i}}\left(\frac{1}{\sqrt{K}}\right)\right\\}$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}P\left(B(\mathbf{X}_{i},1/\sqrt{K})\right)\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{X}_{i}}\left(\frac{1}{\sqrt{K}}\right)\right\\}$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}P\left(B(\mathbf{X}_{i},1/\sqrt{K})\right)\cdot
O_{\mathbf{X}_{i}}\left(\frac{1}{\sqrt{K}}\right)$ $\displaystyle\leq$
$\displaystyle\sum_{i=1}^{N}p_{max}\cdot
O_{\mathbf{X}_{i}}\left(\frac{1}{\left(\sqrt{K}\right)^{d}}\right)\cdot
O_{\mathbf{X}_{i}}\left(\frac{1}{\sqrt{K}}\right)$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}O_{\mathbf{X}_{i}}\left(K^{-\frac{d+1}{2}}\right).$

Here, we obtain (C.3) by using Lemma C.21.

Considering the expectation $E_{\mathbf{X}_{P}^{(N)}}[\cdot]$ for the both
sides of the above equation, we have

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right|\right]\right]$
(171) $\displaystyle=$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[\sum_{i=1}^{N}O_{\mathbf{X}_{i}}\left(K^{-\frac{d+1}{2}}\right)\right]$
$\displaystyle=$ $\displaystyle\sum_{i=1}^{N}O\left(K^{-\frac{d+1}{2}}\right)$
$\displaystyle=$ $\displaystyle N\cdot O\left(K^{-\frac{d+1}{2}}\right).$

In addition, from Lemma C.19 and C.20, it holds that, for $\mathbf{x}\in
id_{B(\mathbf{X}_{i},1/\sqrt{K})}$,

$\left|e^{-T_{k}(\mathbf{x})}-e^{-T_{k}(\mathbf{X}_{i})}\right|=e^{-T_{k}(\mathbf{X}_{i})}\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{X}_{i}}\left(\frac{1}{K}\right),$

and

$\left|e^{-T_{*}(\mathbf{x})}-e^{-T_{*}(\mathbf{X}_{i})}\right|=e^{-T_{*}(\mathbf{X}_{i})}\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{X}_{i}}\left(\frac{1}{K}\right).$

In the similar manner to obtain (171), it holds that, for $\mathbf{x}\in
id_{B(\mathbf{X}_{i},1/\sqrt{K})}$, we have

$E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}(\mathbf{x})}-e^{-T_{k}(\mathbf{X}_{i})}\right|\right]\right]=N\cdot
O\left(K^{-\frac{d+1}{2}}\right),$ (172)

and

$E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{*}(\mathbf{x})}-e^{-T_{*}(\mathbf{X}_{i})}\right|\right]\right]=N\cdot
O\left(K^{-\frac{d+1}{2}}\right).$ (173)

Here, we obtain the upper bounds of the first, third, and fourth terms in
(169).

We now have the upper bounds of the second term in (169).

First, we obtain

$\displaystyle
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|\right]$
(174) $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\right]\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}P\left(B(\mathbf{X}_{i},1/\sqrt{K}\right)\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|.$

Then, from Lemma C.22, we have

$\displaystyle\sum_{i=1}^{N}P\left(B(\mathbf{X}_{i},1/\sqrt{K}\right)\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|$
(175) $\displaystyle=$ $\displaystyle
C\,\sum_{i=1}^{N}p\left(\mathbf{X}_{i}\right)\cdot
O\left(\frac{1}{\left(\sqrt{K}\right)^{d}}\right)\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|.$

Next, note that,

$\sum_{i=1}^{N}p\left(\mathbf{X}_{i}\right)\cdot\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|=N\cdot
E_{\nu}\left|\hat{Q}_{k}^{(N)}-Q^{(N)}\right|.$ (176)

Here, $\nu$ is the countable measure on
$\left\\{\mathbf{X}_{1},\mathbf{X}_{2},\ldots,\mathbf{X}_{N}\right\\}$ defined
as (136).

In addition, since

$\hat{L}_{\alpha}^{(N)}(Q,P;T)=L_{\alpha}(\hat{Q}^{(N)},\hat{P}^{(N)};T)\\\ $

holds, we obtain, from Propsition C.13 and Assumption E1,

$\displaystyle\sqrt{\frac{\alpha}{2}}\cdot
E_{\nu}\Big{|}\hat{Q}_{K}^{(N)}-Q^{(N)}\Big{|}$ $\displaystyle\leq$
$\displaystyle\sqrt{L_{\alpha}(Q^{(N)},P^{(N)};T_{K}^{(N)})-L_{\alpha}(Q^{(N)},P^{(N)};T_{*})}$
(177) $\displaystyle=$
$\displaystyle\sqrt{\hat{L}_{\alpha}^{(N)}(Q,P;T_{K}^{(N)})-\hat{L}_{\alpha}^{(N)}(Q,P;T_{*})}$
$\displaystyle=$ $\displaystyle O\left(\frac{1}{\sqrt{K}}\right).$

Finally, from (174), (175), (176), and (177), we have

$\displaystyle
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|\right]$
(178) $\displaystyle\leq$ $\displaystyle
O\left(\frac{1}{\left(\sqrt{K}\right)^{d}}\right)\sum_{i=1}^{N}p\left(\mathbf{X}_{i}\right)\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|.$
$\displaystyle\leq$ $\displaystyle
O\left(\frac{1}{\left(\sqrt{K}\right)^{d}}\right)\cdot N\cdot
E_{\nu}\left|\hat{Q}_{k}^{(N)}-Q\right|$ $\displaystyle\leq$ $\displaystyle
O\left(\frac{1}{\left(\sqrt{K}\right)^{d}}\right)\cdot N\cdot
O\left(\frac{1}{\sqrt{K}}\right)$ $\displaystyle=$ $\displaystyle N\cdot
O\left(K^{-\frac{d+1}{2}}\right).$

Considering the expectation $E_{\mathbf{X}_{P}^{(N)}}[\cdot]$ for the both
sides of the above equation, we have

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|\right]\right]$
(179) $\displaystyle=$ $\displaystyle E_{\mathbf{X}_{P}^{(N)}}\left[N\cdot
O\left(K^{-\frac{d+1}{2}}\right)\right]$ $\displaystyle=$ $\displaystyle
N\cdot O\left(K^{-\frac{d+1}{2}}\right).$

Summarizing (169), (171), (172), (173), and (179), we obtain

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Delta}\left|e^{-T_{k}^{(N)}}-e^{-T_{k}}\right|\right]\right]$
(180) $\displaystyle\leq$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{x})}-e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\right|\right]\right]$
$\displaystyle\qquad\qquad+E_{\mathbf{X}_{P}^{(N)}}\left[\
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{k}^{(N)}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{X}_{i})}\right|\right]\right]$
$\displaystyle\qquad\qquad+E_{\mathbf{X}_{P}^{(N)}}\left[\
E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{*}(\mathbf{X}_{i})}-e^{-T_{*}(\mathbf{x})}\right|\right]\right]$
$\displaystyle\qquad\qquad+E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\left|e^{-T_{*}(\mathbf{x})}-e^{-T_{k}(\mathbf{x})}\right|\right]\right]$
$\displaystyle=$ $\displaystyle N\cdot O\left(K^{-\frac{d+1}{2}}\right)+N\cdot
O\left(K^{-\frac{d+1}{2}}\right)+N\cdot
O\left(K^{-\frac{d+1}{2}}\right)+N\cdot O\left(K^{-\frac{d+1}{2}}\right)$
$\displaystyle=$ $\displaystyle N\cdot O\left(K^{-\frac{d+1}{2}}\right).$

Here, we see the upper bound of the first term in (167).

Next, we obtain the upper bound of the second and third term in (167).

First, we obtain the upper bound of the second term in (167). Now, we have

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}}\right|\right]\right]$
(181) $\displaystyle=$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Omega\setminus\bigcup_{i=1}^{N}\Delta_{i}}(\mathbf{x})\cdot\left|e^{-T_{k}^{(N)}}\right|\right]\right]$
$\displaystyle=$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[\left\\{1-id_{\bigcup_{i=1}^{N}\Delta_{i}}(\mathbf{x})\right\\}\cdot
e^{-T_{k}^{(N)}}\right]\right]$ $\displaystyle=$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[1-E_{P}\left[\sum_{i=1}^{N}id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}^{(N)}}\right]\right]$ $\displaystyle=$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[1-\sum_{i=1}^{N}E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}^{(N)}}\right]\right].$

Then, from Lemma C.18 and Lemma C.21, we have

$\displaystyle E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}^{(N)}(\mathbf{x})}\right]$ (182) $\displaystyle=$ $\displaystyle
E_{P}\left[id_{B\left(\mathbf{X}_{i},1/\sqrt{K}\right)}(\mathbf{x})\cdot\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}+e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\cdot\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{x}}\left(\frac{1}{\sqrt{K}}\right)\right\\}\right]$
$\displaystyle=$ $\displaystyle
E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\right]\cdot\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}+e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\cdot\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{x}}\left(\frac{1}{\sqrt{K}}\right)\right\\}$
$\displaystyle=$ $\displaystyle
P\left(B(\mathbf{X}_{i},1/\sqrt{K})\right)\cdot\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}+e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\cdot\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{x}}\left(\frac{1}{\sqrt{K}}\right)\right\\}$
$\displaystyle\geq$ $\displaystyle p_{min}\cdot
O_{\mathbf{X}_{i}}\left(\frac{1}{\left(\sqrt{K}\right)^{d}}\right)\cdot\left\\{e^{-T_{k}^{(N)}(\mathbf{X}_{i})}+e^{-T_{k}^{(N)}(\mathbf{X}_{i})}\cdot\left\\{O\left(\frac{1}{\sqrt{K}}\right)+O\left(\frac{1}{K}\right)\right\\}+O_{\mathbf{x}}\left(\frac{1}{\sqrt{K}}\right)\right\\}$
$\displaystyle=$ $\displaystyle
O_{\mathbf{X}_{i}}\left(K^{-\frac{d}{2}}\right)+O_{\mathbf{X}_{i}}\left(K^{-\frac{d+1}{2}}\right).$

From (181) and (182), we see

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}^{(N)}}\right|\right]\right]$
(183) $\displaystyle=$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[1-\sum_{i=1}^{N}E_{P}\left[id_{B(\mathbf{X}_{i},1/\sqrt{K})}(\mathbf{x})\cdot
e^{-T_{k}^{(N)}}\right]\right]$ $\displaystyle\leq$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[1-\sum_{i=1}^{N}O_{\mathbf{X}_{i}}\left(K^{-\frac{d}{2}}\right)+O_{\mathbf{X}_{i}}\left(K^{-\frac{d+1}{2}}\right)\right]$
$\displaystyle=$ $\displaystyle
1-N\cdot\left\\{O\left(K^{-\frac{d}{2}}\right)+O\left(K^{-\frac{d+1}{2}}\right)\right\\}$
$\displaystyle=$ $\displaystyle 1-N\cdot O\left(K^{-\frac{d}{2}}\right)$

In the similar manner to obtain (183), we obtain the upper bound of the third
term in (167):

$E_{\mathbf{X}_{P}^{(N)}}\left[E_{P}\left[id_{\Omega\setminus\Delta}\left|e^{-T_{k}}\right|\right]\right]=1-N\cdot
O\left(K^{-\frac{d}{2}}\right).$ (184)

Summarizing (167), (180), (183) and (184), we have

$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\left|\hat{Q}^{(N)}_{K}-Q_{K}\right|\right]$
(185) $\displaystyle\leq$ $\displaystyle N\cdot
O\left(K^{-\frac{d+1}{2}}\right)+\left\\{1-N\cdot
O\left(K^{-\frac{d}{2}}\right)\right\\}+\left\\{1-N\cdot
O\left(K^{-\frac{d}{2}}\right)\right\\}$ $\displaystyle=$ $\displaystyle
2-N\cdot O\left(K^{-\frac{d}{2}}\right).$

For the upper bound of the second term in (166), from Propsition C.13 we
obtain

$\sqrt{\frac{\alpha}{2}}E_{\mu}\Big{|}Q_{K}-Q\Big{|}\leq\sqrt{L_{\alpha}(Q,P;T_{K})-L_{\alpha}(Q,P;T_{*})}.$
(186)

Thus, under Assumption E2, we see

$E_{\mu}\Big{|}Q_{K}-Q\Big{|}\leq\frac{C_{0}^{\prime}}{\sqrt{K}},$ (187)

where $C_{0}^{\prime}=\sqrt{(2\cdot C_{0})/\alpha}$.

Considering the expectation $E_{\mathbf{X}_{P}^{(N)}}[\cdot]$ for the both
sides of the above equation, we have

$E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\Big{|}Q_{K}-Q\Big{|}\right]=O\left(K^{-\frac{1}{2}}\right).$
(188)

Finally, (165), (185), and (188), we have

$\displaystyle E_{\mathbf{X}_{P}^{(N)}}[W_{1}(Q,\hat{Q}^{(N)}_{K_{0}})]$
$\displaystyle\leq$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[\mathrm{diam}(\Omega)\cdot
E_{\mu}\left|\hat{Q}^{(N)}_{K}-Q\right|\right]$ $\displaystyle=$
$\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\left|\hat{Q}^{(N)}_{K}-Q\right|\right]$
$\displaystyle\leq$ $\displaystyle
E_{\mathbf{X}_{P}^{(N)}}\left[E_{\mu}\Big{|}\hat{Q}^{(N)}_{K}-Q_{K}\Big{|}\right]+E_{\mathbf{X}_{P}^{(N)}}\Big{[}E_{\mu}\left|Q_{K}-Q\right|\Big{]}$
$\displaystyle=$ $\displaystyle 2-N\cdot
O\left(K^{-\frac{d}{2}}\right)+O\left(K^{-\frac{1}{2}}\right).$

From this, for sufficiently large $K>0$, we see

$E_{\mathbf{X}_{P}^{(N)}}[W_{1}(Q,\hat{Q}^{(N)}_{K_{0}})]\leq 2-N\cdot
K^{-\frac{d}{2}}+K^{-\frac{1}{2}}.$

Here, we show (161).

This completes the proof. ∎

###### proof of Corollary C.24.

For 161, substituting $K_{0}$ for $K$, and $N=K_{0}^{\frac{d+\delta}{2}}$, we
have

$\displaystyle E_{\mathbf{X}_{P}^{(N)}}[W_{1}(Q,\hat{Q}^{(N)}_{K_{0}})]$
$\displaystyle\leq$ $\displaystyle 2-K_{0}^{\frac{d+\delta}{2}}\cdot
K_{0}^{-\frac{d}{2}}+K_{0}^{-\frac{1}{2}}$ (190) $\displaystyle=$
$\displaystyle 2-K_{0}^{\frac{\delta}{2}}+K_{0}^{-\frac{1}{2}}.$

This completes the proof. ∎

###### proof of Corollary C.25.

For the setting of the proposition, we have
$K_{0}^{\frac{\delta^{\prime}}{2}}=N^{\frac{\delta^{\prime}}{d+\delta^{\prime}}}=2$.
Thus, for 190, we see

$\displaystyle E_{\mathbf{X}_{P}^{(N)}}[W_{1}(Q,\hat{Q}^{(N)}_{K_{0}})]$
$\displaystyle\leq$ $\displaystyle 2-K_{0}^{\frac{d+\delta}{2}}\cdot
K_{0}^{-\frac{d}{2}}+K_{0}^{-\frac{1}{2}}$ $\displaystyle=$ $\displaystyle
2-2+K_{0}^{-\frac{1}{2}}$ $\displaystyle=$ $\displaystyle
K_{0}^{-\frac{1}{2}}.$

This completes the proof. ∎

##  Appendix D Neumerical Experiments

In this section, we report the results of numerical experiments conducted in
this study.

###  D.1 Experiments on convergence for different values of $\alpha$

In this section, we report the results of the numerical experiments related to
the discussion in Section 5: the results of the numerical experiments on the
convergence of learning for different values of $\alpha$ are presented.

#### Experimental Setup.

For $\alpha=-3,-2,-1,0.2,0.5,0,8,2.0,3.0$, and $4.0$, we generated training
and test dataset, and then trained an NGB model with the training dataset
while estimating the $\alpha$ divergence at each learning step with the test
dataset. One hundred numerical simulations were performed for each $\alpha$.
As a result of the experiment, the median of the estimated value and ranges
between the 45th and 55th percentile quartiles and between the 5th and 95th
percentile quartiles at each learning step are reported.

#### Synthetic Data.

We generated synthetic data of size 5000 from 5-dimantional normal
distribution $\\{X_{1},X_{2},\ldots,X_{5}\\}$ such that $E[X_{i}]=0$,
$\mathrm{Var}[X_{i}]=1$ and $E[X_{i}\cdot X_{j}]=0.8$ ($i\neq j$), for each of
the training and test datasets.

#### Estimating the $\alpha$ divergence.

The $\alpha$ divergence was estimated in the following way

$\displaystyle\hat{D}_{\alpha}(Q||P)(t)$ $\displaystyle=$
$\displaystyle\frac{1}{\alpha\cdot(1-\alpha)}-\frac{1}{\alpha}\hat{E}_{Q}\left[e^{\alpha\cdot
T_{\theta_{t}}(\mathbf{x}^{te})}\right]-\frac{1}{1-\alpha}\hat{E}_{P}\left[e^{(\alpha-1)\cdot
T_{\theta_{t}}(\mathbf{x}^{te})}\right]$ (191) $\displaystyle=$
$\displaystyle\frac{1}{\alpha\cdot(1-\alpha)}-\mathcal{L}_{\alpha}(\theta_{t}),$
(192)

where $T_{\theta_{t}}$ is a model at learning step $t$ in Algorithm 1 and
$\mathbf{x}^{te}$ denotes the test dataset. Note that, decreasing of the
estimated divergence $\hat{D}_{\alpha}(Q||P)(t)$ in (191) implies increasing
of the loss $\mathcal{L}_{\alpha}(\theta_{t})$ in (192).

#### Implementation and Training Details.

We used a neural network which has 3 hidden layers of 100 units in each layer.
The Adam algorithm in PyTorch was used. For the hyperparameters in the
training, the learning rate was 0.001, BathSize was 2500, and the number of
epochs was 500. A NVDIA Tesla K80 GPU was used. It took approximately four
hours to conduct all simulations for each value of $\alpha$.

(a) $\alpha=-3$ |  (b) $\alpha=-2$
---|---
(c) $\alpha=-1$ |  (d) $\alpha=2.0$
(e) $\alpha=3.0$ |  (f) $\alpha=4.0$
Figure 1: Results of estimating the $\alpha$ divergence for
$\alpha=-3,-2,-1,2,3$ and $4$, over the number of learning steps during the
optimization. The $y$-axis of each graph represents the estimated value of the
$\alpha$ divergence, and the $x$-axis of each graph represents the learning
step. The solid blue line shows the median of the estimates of the $\alpha$
divergence. The dark blue area shows the ranges of the estimates between the
45th and 55th percentiles, and the light blue area shows the range of the
estimates between the 5th and 95th percentile quartiles.

(a) *

(g) $\alpha=0.2$ |  (h) $\alpha=0.5$
---|---
(i) $\alpha=0.8$ |

Figure 2: Results of estimating the $\alpha$ divergence for $\alpha=0.2,0.5$
and $0.8$, over the number of learning steps during the optimization. The
$y$-axis of each graph represents the estimated value of the $\alpha$
divergence, and the $x$-axis of each graph represents the learning step. The
solid blue line shows the median of the estimates of the $\alpha$ divergence.
The dark blue area shows the ranges of the estimates between the 45th and 55th
percentiles, and the light blue area shows the range of the estimates between
the 5th and 95th percentile quartiles.

#### Results.

Figure 1 and 2 show the results of estimating the $\alpha$ divergence over the
number of learning steps during the optimization. Figure 1 is for
$\alpha=-3,-2,-1,2,3$ and $4$, and Figure 2 is for $\alpha=0.2,0.5$, and
$0.8$. The $y$-axis of each graph represents the estimated value of the
$\alpha$ divergence, and the $x$-axis of each graph represents the learning
step. The solid blue line shows the median of the estimates of the $\alpha$
divergence. The dark blue area shows the ranges of the estimates between the
45th and 55th percentiles, and the light blue area shows the range of the
estimates between the 5th and 95th percentile quartiles.

As shown in Figure 1, the estimates of the $\alpha$ divergence diverged. This
corresponds to a negative divergence of the loss function
$\mathcal{L}_{\alpha}(\theta_{t})$ in (192), and then implies that
$E_{Q}[e^{T_{\theta_{t}}}]\rightarrow 0$ for $\alpha>1$, and
$E_{Q}[e^{T_{\theta_{t}}}]\rightarrow\infty$ for $\alpha<0$ in (191). The
discussion in Section 5 suggests that
$E[\nabla_{\theta}\mathcal{L}_{\alpha}(\theta)]\rightarrow\vec{0}$. That is,
the gradients of the neural networks in this case vanished for
$\alpha=-3,-2,-1,2,3$ and $4$. However, as shown in Figure 2, the estimates of
the $\alpha$ divergence converge stably for $\alpha=0.2,0.5$, and $0.8$.

###  D.2 Experiments to confirm the relationship between dimensions of dataset
and steps in training

In this section, we report the results of numerical experiments related to the
discussion in Section 7: the results of numerical experiments to confirm the
relationship between dimensions of dataset and steps in training are
presented.

#### Experimental Setup.

We generated training and test datasets of dimensions $d=2$, $3$, $4$, $5$,
$6$, and $7$, and then trained an NGB model with the training dataset while
estimating the $\alpha$ divergence at each learning step with the test
dataset. One hundred numerical simulations were performed for each dimension
$d$. As a result of the experiment, the median of the estimated value and
ranges between the 5th and 95th percentile quartiles at each learning step are
reported.

#### Synthetic Data.

For each $d=2$, $3$, $4$, $5$, $6$, and $7$, we generated the training and
test datasets of size 5000 from $d$-dimantional normal distribution
$\\{X_{1},X_{2},\ldots,X_{d}\\}$, such that $E[X_{i}]=0$,
$\mathrm{Var}[X_{i}]=1$ and $E[X_{i}\cdot X_{j}]=0.8$ ($i\neq j$).

#### Estimating the $\alpha$ divergence.

The $\alpha$ divergence was estimated in the following way:

$\displaystyle\hat{D}_{\alpha}(Q||P)(t)$ $\displaystyle=$
$\displaystyle\frac{1}{\alpha\cdot(1-\alpha)}-\frac{1}{\alpha}\hat{E}_{Q}\left[e^{\alpha\cdot
T_{\theta_{t}}(\mathbf{x}^{te})}\right]-\frac{1}{1-\alpha}\hat{E}_{P}\left[e^{(\alpha-1)\cdot
T_{\theta_{t}}(\mathbf{x}^{te})}\right]$ (193) $\displaystyle=$
$\displaystyle\frac{1}{\alpha\cdot(1-\alpha)}-\mathcal{L}_{\alpha}(\theta_{t}),$
(194)

where $T_{\theta_{t}}$ is a model at learning step $t$ in Algorithm 1 and
$\mathbf{x}^{te}$ denotes the test dataset. Note that, decreasing of the
estimated divergence $\hat{D}_{\alpha}(Q||P)(t)$ in (193) implies increasing
of the loss $\mathcal{L}_{\alpha}(\theta_{t})$ in (194).

#### Implementation and Training Details.

We used a neural network which has 3 hidden layers of 100 units in each layer.
The Adam algorithm in PyTorch was used. For the hyperparameters in the
training, the learning rate was 0.001, BathSize was 2500, and the number of
epochs was 500. A NVDIA Tesla K80 GPU was used. It took approximately four
hours to conduct all simulations for each $d$.

Table 2: The early stop step (“$N^{2/d}$”) and the median of the steps at which the estimated divergence reaches its maximum (“$\operatorname{median}(K_{\operatorname{max}})$”), for each dimension $d=2$, $3$, $4$, $5$, $6$, and $7$. | $d=2$ | $d=3$ | $d=4$ | $d=5$ | $d=6$ | $d=7$ |
---|---|---|---|---|---|---|---
$N^{2/d}$ | 5000 | 292 | 71 | 30 | 17 | 11 |
$\operatorname{median}(K_{\operatorname{max}})$ | 130 | 112 | 130 | 136 | 50 | 50 |
(j) $d=2$ |  (k) $d=3$
---|---
(l) $d=4$ |  (m) $d=5$
(n) $d=6$ |  (o) $d=7$
Figure 3: Results of estimating the $\alpha$ divergence for
$\alpha=-3,-2,-1,2,3$ and $4$, over the number of learning steps during the
optimization. The $y$-axis of each graph represents the estimated value of the
$\alpha$ divergence divided by the true value of the divergence, and the
$x$-axis of each graph represents the learning step. The dashed red line
indicates Y=1, which corresponds to the theoretical value of the estimate for
each $d$. The solid blue line shows the median of the estimates of the
$\alpha$ divergence. The dark blue area shows the ranges of the estimates
between the 45th and 55th percentiles, and the light blue area shows the range
of the estimates between the 5th and 95th percentile quartiles.

#### Results.

Let $K_{\operatorname{max}}$ denote the step at which the estimated divergence
reaches its maximum:

$K_{\operatorname{max}}=\underset{t}{\operatorname{argmax}}\hat{D}_{\alpha}(Q||P)(t).$
(195)

Table 3(c) lists $N^{2/d}$, the early stop step obtained from (23), and the
median of $K_{\operatorname{max}}$, for each dimension $d=2$, $3$, $4$, $5$,
$6$, and $7$. In Figure 3, we show the results of estimating the $\alpha$
divergence over the number of learning steps during the optimization. Since
the value of the $\alpha$ divergence changes as the dimension of the dataset
changes, we divided by the the estimated value of the divergence by the true
value of the divergence to normalize the results of each dimension. The
$y$-axis of each graph represents the estimated value of the $\alpha$
divergence divided by the true value of the divergence, and the $x$-axis of
each graph represents the learning step. The solid blue line shows the median
of the estimates of the $\alpha$ divergence. The light blue area shows the
range of the estimates between the 5th and 95th percentile quartiles. The
dashed red line indicates Y=1, which corresponds to the theoretical value of
the estimate for each $d$.

As shown in Table 2, the steps from the early stop method and those at which
the estimates decreased were approximately consistent, except in the case of
$d=2$. However, the estimates of the data of the low dimensions, particularly
$d=2$, decreased earlier than the early stop method suggests. This may be
because $C$ in (23) for the data of low dimensions can be small because the
neural network learns quickly when the dimensions of the data are low.
However, Figure 3 shows that the estimates of the divergence decreased slowly
when the dimensions of the data are low, and they decreased more quickly when
the dimensions of the data were higher. These results suggest that the curse
of dimensionality of balancing is easier to observe when dimensions of data
are higher.

###  D.3 Experiments for esitimating causal effects of joint and
multidimensional interventions with different sample sizes

In this section, we report the results of numerical experiments related to the
discussion in Section 8: the results of numerical experiments for esitimating
causal effects of joint and multidimensional interventions with different
sample sizes are presented.

#### Experimental Setup.

The following two experiments were conducted, in which synthetic data of size
$N=\allowbreak 1000,\allowbreak 10000$, and $100000$ were generated using the
method developed by Vegetabile et al.(2021).

  * •

Experiment 1. An experiment on estimating the causal effect of a single
intervention, especially for continuous intervention,
$E_{\widetilde{P}_{1}}[Y|A,\mathbf{X}]$ where
$\widetilde{P}_{1}=P(Y,\overline{do}(A),\mathbf{X})$.

  * •

Experiment 2. An experiment to estimate the causal effect of a mixture of both
arbitrary discrete and continuous interventions,
$E_{\widetilde{P}_{2}}[Y|A,\mathbf{X}]$ where
$\widetilde{P}_{2}=P(Y,\overline{do}(A),\overline{do}(\mathbf{X}_{1}),\overline{do}(X_{2}),\overline{do}(\mathbf{X}_{3})]$.

#### Experimental Details.

Experiments 1 and 2 were conducted using the following steps.

Step 1: We created training dataset of size $N=1000,\allowbreak 10000$, and
$100000$, and test dataset with size $N=1000$. The training dataset were
generated using the method developed by Vegetabile et al.(2021). The test
dataset were generated from the following distribution:

  * •

Experiment 1. $\widetilde{P}_{1}=P(Y,\overline{do}(A),\mathbf{X})$,

  * •

Experiment 2.
$\widetilde{P}_{2}=P(Y,\overline{do}(A),\overline{do}(\mathbf{X}_{1}),\overline{do}(X_{2}),\overline{do}(\mathbf{X}_{3}))$,

where $P$ denotes the distribution of the training dataset. To create the test
dataset, we shuffled the dataset generated from the same distribution as the
training dataset.

Step 2: The balancing weights were estimated for each experiment. We estimated
$\mathit{BW}(A,\mathbf{X}:T_{\theta_{*}})$ for Experiment 1, and
$\mathit{BW}(A,\mathbf{X}_{1},X_{2},\mathbf{X}_{3}:T_{\theta_{*}})$ for
Experiment 2.

Step 3: We created models for each experiment using the linear regression (LR)
or the gradient boosting tree (GBT) algorithm with our weights from the
previous step. The hyperparameters were tuned to create models of GBT.

Step 4: We estimate the average causal effects
$E_{\widetilde{P}_{1}}[\allowbreak
Y\allowbreak|A,\allowbreak\mathbf{X}\allowbreak]$ and
$E_{\widetilde{P}_{2}}[\allowbreak
Y\allowbreak|A,\allowbreak\mathbf{X}\allowbreak]$ using the predictions of the
models from Step 2 with the test dataset. Finally, we report the mean squared
error (RMSE) between the true and estimated values.

#### Baseline Method.

The main baseline method used in our experiments is entropy balancing [31]. We
compared our method with the method for balancing $\mathbf{X}$ with $A$ for
each of the moments from 1 to 4. For Experiment 1, both our method and the
baseline method estimated the same target:
$E_{\widetilde{P}_{1}}[Y|A,\mathbf{X}]$ where
$\widetilde{P}_{1}(Y,A,\mathbf{X})=P(Y,\overline{do}(A),\mathbf{X})$. However,
no existing method can fully deal with the target of Experiment 2:
$E_{\widetilde{P}_{2}}[Y|A,\mathbf{X}]$ where
$\widetilde{P}_{2}=P(Y,\overline{do}(A),\overline{do}(\mathbf{X}_{1}),\overline{do}(X_{2}),\overline{do}(\mathbf{X}_{3})]$.
Therefore, the same entropy balancing as in Experiment 1 was used in
Experiment 2. This may be an unfair comparison to the baseline method. In
addition, we included a “naïve” estimation, using algorithms with no sample
weights, as a baseline. For the calculation of entropy balancing weights,
WeightIt library in R was used.
888https://cran.r-project.org/web/packages/WeightIt/index.html.

#### Training Data Set.

Specifically, we used the following steps to generate the dataset. First,
$\mathbf{W}=(W_{1},W_{2},W_{3},W_{4},W_{5})$ were generated independently,
such that $W_{1}\sim\mathcal{N}(-0.5,1)$, $W_{2}\sim\mathcal{N}(1,1)$,
$W_{3}\sim\mathcal{N}(0,1)$, $W_{4}\sim\mathcal{N}(1,1)$, and
$\mathcal{X}_{W_{5}}=\\{0,1,2\\}$ with $P(W_{5}=0)=0.70$ and
$P(W_{5}=1)=P(W_{5}=2)=0.15$. Second, $A$ and $Y$ were generated as follows:

$\displaystyle A$ $\displaystyle\sim$
$\displaystyle\mathscr{X}^{2}(df=3,\,\mu_{A}(W_{1},W_{2},W_{4},W_{5})),$
$\displaystyle Y$ $\displaystyle=$
$\displaystyle\frac{1}{50}\left[\left(-0.15A^{2}+A(W_{1}^{2}+W_{2}^{2})-15\right)\right.$
(196) $\displaystyle\quad+\left((W_{1}+3)^{2}+2(W_{2}-25)^{2}+W_{3}\right)$
$\displaystyle\quad\left.-C+\varepsilon\right],$

where $\mu_{A}(W_{1},W_{2},W_{4},W_{5})=5|W_{1}|+6|W_{2}|+|W_{4}|+a$, and
$a=0$ if $W_{5}=0$, and $a=1$ if $W_{5}=1$, and $a=5$ if $W_{5}=2$, and
$C=E[(W_{1}+3)^{2}]+2E[(W_{2}-25)^{2}]+E[W_{3}]$, and
$\varepsilon\sim\mathcal{N}(0,1)$. Here, $\mathscr{X}^{2}(df=n,\mu)$ is the
noncentral $\chi^{2}$ distribution with n degrees of freedom and a noncentral
parameter $\mu$. Finally, we create new variables
$\mathbf{X}=(\mathbf{X}_{1},X_{2},\mathbf{X}_{3})$, as observed values of
$\mathbf{W}$ using the following transformation:

$\displaystyle\mathbf{X}_{1}$ $\displaystyle=$
$\displaystyle(X_{(1,1)},X_{(1,2)},X_{(1,3)}),$ where $\displaystyle
X_{(1,1)}=\exp\left(W_{1}/2\right),$ $\displaystyle
X_{(1,2)}=W_{2}/(1+\exp(W_{1}))+10,$ $\displaystyle
X_{(1,3)}=W_{1}W_{3}/25+0.6,$ $\displaystyle X_{2}$ $\displaystyle=$
$\displaystyle(W_{4}-1)^{2},$ (198) $\displaystyle\mathbf{X_{3}}$
$\displaystyle=$ $\displaystyle\begin{cases}(1,0)\quad\text{if $W_{5}=0$},\\\
(0,1)\quad\text{if $W_{5}=1$},\\\ (0,0)\quad\text{if $W_{5}=2$}.\end{cases}$
(199)

#### Test Data Set.

We first generated dataset from the same distribution as the training dataset.
Second, the dataset were shuffled by the index, with the following divided
parts treated as a single piece of data: for Experiment 1, $A$ and
$\mathbf{X}$ were shuffled by the index, and for Experiment 2, each of $A$,
$\mathbf{X}_{1}$, $X_{2}$ and $\mathbf{X}_{3}$ were shuffled by the index.
Third, using the inverse transformation of Eq. (D.3)-(199), we calculated
$(W_{1},W_{2},W_{3},W_{4},W_{5})$ from
$\mathbf{X_{1}}=(X_{(1,1)},X_{(1,2)},X_{(1,3)})$, $X_{2}$, and
$\mathbf{X_{3}}$ of the shuffled dataset:

$\displaystyle W_{1}$ $\displaystyle=$ $\displaystyle 2\log X_{(1,1)},\quad
W_{2}=X_{(1,2)}\cdot(1+X_{(1,1)}^{2}),$ $\displaystyle W_{3}$ $\displaystyle=$
$\displaystyle\frac{25(X_{(1,3)}-0.6)}{2\log X_{(1,1)}},\quad
W_{4}=\sqrt{X_{2}}+1,$ $\displaystyle W_{5}$ $\displaystyle=$
$\displaystyle\begin{cases}0\quad\text{if $\mathbf{X_{3}}=(1,0)$},\\\
1\quad\text{if $\mathbf{X_{3}}=(0,1)$},\\\ 2\quad\text{if
$\mathbf{X_{3}}=(0,0)$}.\end{cases}$

Finally, the true values of $Y$ for causal effects were calculated using the
terms in Eq. (196) without the term $\varepsilon$.

#### Implementation and Training Details.

$N=1000$: For experiments with the dataset of size $N=1000$, we used a neural
network which has 10 hidden layers of 100 units in each layer. $\alpha=0.5$
was used to estimate the divergence. The Adam algorithm in PyTorch was used.
For the hyperparameters in the training, the learning rate was 0.0001,
BathSize was 1000, and the number of epochs was 70. A NVDIA Tesla K80 GPU was
used. It took approximately 40 min to conduct all the simulations for each
experiment.

$N=10000$: For experiments with the dataset of size $N=10000$, We used a
neural network which has 10 hidden layers of 100 units in each layer.
$\alpha=0.5$ was used to estimate the divergence. The Adam algorithm in
PyTorch was used. For the hyperparameters in the training, the learning rate
was 0.0001, BathSize was 2500, and the number of epochs was 200. A NVDIA Tesla
K80 GPU was used. It took approximately 7 h to conduct all the simulations for
each experiment.

$N=100000$: For experiments with the dataset of size $N=100000$, We used a
neural network which has 10 hidden layers of 100 units in each layer.
$\alpha=0.5$ was used to estimate the divergence. The Adam algorithm in
PyTorch was used. For the hyperparameters in the training, the learning rate
was 0.0001, BathSize was 2500, and the number of epochs was 200. A NVDIA Tesla
K80 GPU was used. It took approximately 78 h to conduct all the simulations
for each experiment.

Table 3: Average RMSE for estimation in Experiments 1 and 2 for dataset of
size $N=1000$, $10000$, and $100000$. For entropy balancing, the number to the
right side of the method name, “($\mathit{m}$),” denotes the number of moments
that are balanced. The results from 100 simulations are in the form of “mean
(std. err.)” . (a) $N=1000$

| Experiment 1 | Experiment 2
---|---|---
Method | LR | GBT | LR | GBT
Unweighted | 1.347(0.039) | 0.739(0.066) | 1.347(0.033) | 0.741(0.068)
Entropy Balancing(1) | 1.303(0.056) | 0.724(0.058) | 1.303(0.052) | 0.726(0.060)
Entropy Balancing(2) | 1.206(0.029) | 0.693(0.056) | 1.206(0.026) | 0.698(0.055)
Entropy Balancing(3) | 1.201(0.026) | 0.690(0.054) | 1.201(0.024) | 0.698(0.061)
Entropy Balancing(4) | 1.203(0.027) | 0.699(0.057) | 1.203(0.025) | 0.699(0.061)
NBW | 1.347(0.039) | 0.745(0.065) | 1.347(0.034) | 0.738(0.063)

(b) $N=10000$

| Experiment 1 | Experiment 2
---|---|---
Method | LR | GBT | LR | GBT
Unweighted | 1.342(0.030) | 0.489(0.035) | 1.342(0.026) | 0.489(0.039)
Entropy Balancing(1) | 1.295(0.033) | 0.486(0.026) | 1.295(0.030) | 0.487(0.035)
Entropy Balancing(2) | 1.194(0.025) | 0.466(0.036) | 1.194(0.025) | 0.468(0.041)
Entropy Balancing(3) | 1.187(0.025) | 0.459(0.032) | 1.187(0.024) | 0.457(0.036)
Entropy Balancing(4) | 1.189(0.024) | 0.457(0.035) | 1.189(0.023) | 0.452(0.034)
NBW | 1.274(0.038) | 0.488(0.035) | 1.273(0.031) | 0.485(0.032)

(c) $N=100000$

| Experiment 1 | Experiment 2
---|---|---
Method | LR | GBT | LR | GBT
Unweighted | 1.342(0.027) | 0.457(0.048) | 1.342(0.023) | 0.459(0.044)
Entropy Balancing(1) | 1.299(0.029) | 0.453(0.037) | 1.298(0.027) | 0.455(0.036)
Entropy Balancing(2) | 1.195(0.025) | 0.391(0.034) | 1.194(0.023) | 0.386(0.039)
Entropy Balancing(3) | 1.186(0.024) | 0.361(0.025) | 1.186(0.023) | 0.360(0.023)
Entropy Balancing(4) | 1.188(0.024) | 0.353(0.022) | 1.187(0.023) | 0.356(0.020)
NBW | 1.239(0.095) | 0.376(0.033) | 1.252(0.080) | 0.388(0.030)

#### Results.

We report the average and standard errors of the root mean squared error
(RMSE) between the estimated and true values of the average causal effects for
synthetic data of size $N=1000$, $10000$, and $100000$. Table 3(c) lists the
results of Experiments 1 and 2 for each $N$. Each result is in the form of
“mean (std. err.)” from 100 simulations.

As shown in all the results, the results of NBW were less accurate than those
of the entropy-balancing method. Moreover, the results for $N=1000$ shows that
NBW were less accurate than the unweighted estimation. However, as seen in all
results for $N=100000$, the accuracy of NBW was superior to that of the
unweighted estimation, which was close to the accuracy of the entropy-
balancing method. These results imply that the sample size requirements of the
proposed method are larger than those of the entropy balancing method.

##  Appendix E Back-Propagation Algorithm using Neural Balancing Weights

We show a back-propagation algorithm using NBW for MSE loss in Algorithm 3.
The MSE loss here is calculated by the mean of the element wise product of
both the original squared errors and the balancing weights.

Algorithm 3 Back-Propagation Algorithm using Neural Balancing Weights

1:Data
$(y,\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{n},\mathbf{z})=\\{(y^{i},\mathbf{x}_{1}^{i},\mathbf{x}_{2}^{i},\ldots,\mathbf{x}_{n}^{i},\mathbf{z}^{i})|i=1,2,\ldots,N\\}$

2: A Neural Balancing Weight Model $T$

3:A Neural Network Model $f_{\phi}$ for Estimating
$E_{\widetilde{P}}[Y|\mathbf{X},\mathbf{Z}]$

4:repeat

5: $\hat{y}\leftarrow
f_{\phi}(\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{n},\mathbf{z})$ //
Forward Propagation

6:
$\mathit{BW}(\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{n},\mathbf{z})\leftarrow\frac{e^{-T(\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{n},\mathbf{z})}}{\mathit{MEAN}(e^{-T(\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{n},\mathbf{z})})}$

7: $\mathit{Err}_{\phi}\leftarrow y-\hat{y}$ // Obtaining Errors for
$f_{\phi}$

8:
$\mathcal{L}({\phi})\leftarrow\mathit{MEAN}(\mathit{Err}_{\phi}\otimes\mathit{Err}_{\phi}\otimes\mathit{BW}(\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{n},\mathbf{z}))$
// Calculating Loss $\mathcal{L}({\phi})$

9: $\phi\leftarrow\phi-\nabla\mathcal{L}({\phi})$

10:until convergence


[affiliation=1,2]ThemosStafylakis [affiliation=3]AnnaSilnova
[affiliation=3]JohanRohdin [affiliation=3]OldrichPlchot
[affiliation=3]LukasBurget

# Challenging margin-based speaker embedding extractors by using the
variational information bottleneck

###### Abstract

Speaker embedding extractors are typically trained using a classification loss
over the training speakers. During the last few years, the standard
softmax/cross-entropy loss has been replaced by the margin-based losses,
yielding significant improvements in speaker recognition accuracy. Motivated
by the fact that the margin merely reduces the logit of the target speaker
during training, we consider a probabilistic framework that has a similar
effect. The variational information bottleneck provides a principled mechanism
for making deterministic nodes stochastic, resulting in an implicit reduction
of the posterior of the target speaker. We experiment with a wide range of
speaker recognition benchmarks and scoring methods and report competitive
results to those obtained with the state-of-the-art Additive Angular Margin
loss.

###### keywords:

speaker recognition, variational information bottleneck

##  1 Introduction

During the last several years, speaker embeddings extracted with speaker-
discriminatively trained deep neural networks have attained impressive
performance on several datasets. During this period, the field has witnessed
numerous architectures, such as LSTM-based models ([1]), 1D (TDNN, ECAPA-TDNN,
xi-vectors, TitaNet [2, 3, 4, 5]) and 2D CNNs (ResNet, ResNeXt, Res2Net [6, 7,
8]), and more recently self-supervised pretrained Transformer models (Wav2Vec
2.0, HuBERT, WavLM [9, 10]), fined-tuned for the particular downstream task.

One of the most effective and architecture-agnostic methods to improve speaker
embedding discriminability has been the margin-based loss, especially the
Additive-Angular Margin (AAM [11, 12]) variant. Originating from face
recognition, AAM combines directional statistics with a margin penalty on the
angle between the embedding and the prototype of the target class. It comes as
a drop-in replacement of regular softmax/cross-entropy loss and achieves
substantial improvements in many speaker recognition benchmarks.

On the other hand, the effectiveness of margin in other, more challenging
datasets is questionable. Margin-based losses do not retain much of their
strength when a backend model (e.g., probabilistic linear discriminant
analysis, PLDA) needs to be introduced in the pipeline [13, 14]. This is
crucial since most of the successful systems in NIST-SRE do make use of such a
backend. Furthermore, in industry-level voice biometrics, where a single
extractor typically serves several different deployments (often in different
languages), a backend model trained on in-domain data is usually compulsory
for attaining state-of-the-art performance. Finally, margin-based losses are
not based on a probabilistic framework, which is a desired property for tasks
and methods like calibration, self-supervised training, and multitask
learning.

In this paper, we propose an alternative to margin-based loss, based on the
variational information bottleneck (VIB [15]). The method has been introduced
in deep learning for some years, and although it has been applied to many
fields (e.g., speech antispoofing, computer vision, NLP, explainability, a.o.
[16, 17, 18, 19, 20]), it has not been sufficiently explored for training
speaker embedding extractors. In [21], an approach similar to ours is
presented, however, the experiments are conducted only with VoxCeleb1, which
is considered low-resource, with limited test-set, and having low speaker
variability compared to, e.g., NIST-SRE or CNCeleb. Furthermore, to the best
of our knowledge, we are the first to underscore the similarities between VIB
and margin-based losses, conduct experiments on challenging benchmarks such as
NIST-SRE, and combine it with backend models (e.g., PLDA).

In the rest of the paper, we examine whether making the embedding stochastic
and sampling from its conditional distribution during training may have a
similar impact on the speaker-discriminability of the embeddings. We
experiment with a wide collection of speaker recognition benchmarks, namely
VoxCeleb [22, 23], CNCeleb [24], and the latest NIST-SREs [25, 26, 27].
Moreover, we examine different backend methods, ranging from pure cosine-
scoring to PLDA and Toroidal Probabilistic Spherical Discriminant Analysis
(T-PSDA) [28]. Finally, we provide arguments in favor of further exploration
of VIB and propose future research directions.

##  2 Training the network with VIB

In this section, the method for training the embedding extractor with VIB is
demonstrated. We begin by providing some intuition of the Information
Bottleneck principle and the derivation of VIB. We then discuss similarities
of VIB with other methods, while providing implementation details of our
method.

###  2.1 The Information Bottleneck principle

The Information Bottleneck [29] builds on the desired properties of a model
utilizing some internal representation: such representation has to be
effective at performing the task the model is used for; at the same time, the
information about the input data contained in such representation should be as
compressed as possible. In IB, these general considerations are formalized by
introducing the following maximization problem:

$R_{IB}=I(Z,Y;\bm{\theta})-\beta I(Z,X;\bm{\theta}),$ (1)

where the second term corresponds to the mutual information between the input
$X$ and its internal representation $Z$ that should be minimized, and the
first one is the mutual information between $Z$ and the output $Y$ to be
maximized. The scalar $\beta>0$ controls the amount of compression of $X$ in
$Z$, while $\bm{\theta}$ is a vector of model parameters.

###  2.2 Derivation of VIB

In practice, for an arbitrary model, estimating mutual information terms in
eq. (1) can be challenging. Deep Variational Information Bottleneck (VIB) [15]
addresses this issue by introducing two variational approximations to the IB
by parameterizing it as a neural network. The first is $q(y|z;\bm{\psi})$,
which is an approximation to the true conditional distribution $p(y|z)$, where
the former is parametrized by the decoder's parameters $\bm{\psi}$, which in
our case corresponds to the linear layer that maps $Z$ to the logits (i.e., it
is a collection of speaker prototypes, one for each training speaker). The
second variational approximation is $r(z)$, with which we approximate the
marginal distribution of $Z$, i.e., $r(z)\approx p(z)=\int
p(z|x;\bm{\theta})p(x)dx$. Although the parameters of $r(z)$ can be learnable,
it is common to keep them fixed and equal to those of a standardized
multivariate normal (MVN), i.e., $r(z)={\cal N}(z|0,I)$ – the same is used in
the experiments of this paper.

Using the above variational approximations, we obtain

$I(Z,Y;\bm{\theta})\geq\int p(x)p(y|x)p(z|x;\bm{\theta})\log
q(y|z;\bm{\psi})dxdydz$ (2)

which is equal to the negative expected cross-entropy, where the expectation
is with respect to the training data $p(x)p(y|x)$ but also with respect to
$p(z|x;\bm{\theta})$.

Similarly, $I(Z,X;\bm{\theta})$ is upper-bounded by the following expression

$I(Z,X;\bm{\theta})\leq\int
p(x)p(z|x;\bm{\theta})\log\frac{p(z|x;\bm{\theta})}{r(z)}dxdz$ (3)

which is equal to the expected Kullback-Leibler (KL) divergence between
$p(z|x;\bm{\theta})$ and $r(z)$, with respect to the training data
distribution $p(x)$. Therefore, $R_{IB}$ is lower-bounded by the following
expression

$R_{IB}\geq\mathbb{E}_{p(x)p(y|x)p(z|x;\bm{\theta})}[\log q(y|z;\bm{\psi})]\\\
-\beta\mathbb{E}_{p(x)p(z|x;\bm{\theta})}\left[\log\frac{p(z|x;\bm{\theta})}{r(z)}\right]\\\
=-\mathbb{E}_{p(x)p(z|x;\bm{\theta})}[\text{CE}(p(y|x),q(y|z;\bm{\psi})]\\\
-\beta\mathbb{E}_{p(x)}[\text{KL}(p(z|x;\bm{\theta}),r(z))]\\\
\approx\frac{1}{n}\sum_{i=1}^{n}\bigg{[}\frac{1}{m}\sum_{j=1}^{m}\left(\log
q(y_{i}|z^{(j)};\bm{\psi})\right)\\\
-\beta\text{KL}(p(z|x_{i};\bm{\theta}),r(z))\bigg{]}$ (4)

###  2.3 Sampling from the conditional distribution via the reparametrization
trick

Similarly to variational autoencoders (VAEs [30]), VIB models
$p(z|x;\bm{\theta})$ by a MVN with diagonal covariance matrix ${\cal
N}(z|f_{\mu}(x),\text{diag}(f_{\sigma}^{2}(x)))$. Therefore, the KL term in
eq. (4) for a given $x$ is simply the KL divergence between two MVN
distributions. For calculating the CE term in eq. (4) for each given $x$,
samples from $p(z|x;\bm{\theta})$ are required. The reparametrization trick of
VAEs is utilized, where samples $e^{(j)}$ are generated from a standardized
MVN and transformed to samples $z^{(j)}$ from the target distribution
$p(z|x;\theta)$ as follows

$z^{(j)}=f_{\sigma}(x)\odot e^{(j)}+f_{\mu}(x),\,e^{(j)}\sim{\cal N}(0,I).$
(5)

where $j=1,\ldots,J$. We implement $f_{\mu}(x)$ and $f_{\sigma}(x)$ with two
linear layers having the statistics pooling layer as common input. Note that
statistics pooling is defined as mean and std pooling for each frequency-
channel pair of the feature maps of the last convolutional layer. These mean
and std pooled features should not be confused with $f_{\mu}(x)$ and
$f_{\sigma}(x)$, as the latter correspond to the mean and std parametrizing
the distribution of the stochastic embedding $Z$. The non-negativity of the
diagonal elements of $f_{\sigma}(x)$ is enforced by the softplus function,
although the exponential function may also be applied. Finally, the encoder is
a ResNet-34 with statistics pooling layer, although other architectures may be
considered.

###  2.4 VIB and VAE

As we observe, VIB and VAE have several similarities, since (a) they both make
use of a stochastic internal representation $Z$, (b) they prevent the collapse
of $p(z|x;\bm{\theta})$ to a point mass by penalizing the KL divergence
between it and $p(z)$, and (c) they make use of the reparametrization trick to
sample from $p(z|x)$. From this perspective, one may consider VIB as a
supervised learning analog of VAEs. On the other hand, there are certain
differences, the main of which are the following.

  * •

The classifier-decoder of the VIB, as opposed to the (deep) reconstruction-
based decoder of VAE. Note that this allows generating many samples $z^{(j)}$
(i.e., $J\gg 1$) per training example to reduce the variance of the estimator
with a minor increase in computation since the classifier is merely a linear
layer (we set $J=10$).

  * •

The information-theoretic formulation of VIB, as opposed to the variational-
Bayes (VB) formulation of VAEs, which among others provides a better
justification of the tunable trade-off parameter $\beta>0$ (The VB formulation
implies $\beta=1$).

###  2.5 Similarities with margin-based losses

Apart from VAEs, VIB has certain similarities with margin-based losses. The
AAM loss is a variant of Softmax/CE, where the logit of the $i$th training
example and $j$th training speaker is modified as
$l_{i,j}=s\cos(\theta_{z_{i},y_{j}}+m\delta_{j,y_{i}})$, where $s$ is the
scale (typically $s=30$), $m$ is the margin, $\delta_{\cdot,\cdot}$ is the
Dirac function, and $\theta_{z_{i},y_{j}}$ is the (positive) angle between the
embedding $z_{i}$ and the $j$th speaker prototype. It therefore makes the
classification task artificially harder, by intervening to the logits and
reducing the one corresponding to the target class. This intervention creates
a safety zone around the class boundaries, which increases the
discriminability of the speaker embeddings, especially for speakers not
appearing in the training set.

In the VIB framework, $Z$ is assumed stochastic, its conditional distribution
$p(z|x;\bm{\theta})$ is estimated and samples are generated from it, which
should be classified using regular Softmax/CE. Having to classify samples from
$p(z|x;\bm{\theta})$ instead of the expected value
$\mathbb{E}_{p(z|x;\bm{\theta})}[Z]$ makes the task harder, as it increases
the within-class variability of the embeddings. This results from the fact
that $\mathbb{E}_{p(z|x;\bm{\theta})}\log q(y|z;\bm{\psi})\leq\log
q(y|\mathbb{E}_{p(z|y;\bm{\psi})}[Z];\bm{\psi})$ where $y$ is the ground-truth
speaker. As sampling is applied only during training, it creates a safety zone
around the class boundaries, similarly to the AAM loss.

Note that other regularization techniques are employed for making the task
harder during training, such as data augmentation (by adding noise and
reverberation or applying SpecAugment), short training utterances, or dropout.
Nonetheless, the similarity between VIB and margin stems from the fact that
they both act on the classifier of the architecture, while the other methods
are either input-level or act on the intermediate layers. Therefore, as with
margin-based losses, the other regularization techniques are orthogonal to
VIB.

###  2.6 Angular margin, unit-length normalization, and VIB

Margin-based losses are typically combined with angular distances, e.g., by
unit-length normalizing the embedding and the class/speaker prototypes. For
example, penalizing the target logit with a fixed additive margin without
applying such unit-length normalizations may not be effective, as the network
would be free to overcome the penalty by increasing the average magnitude of
the embeddings, and therefore increasing the average scale of the logits,
rendering the margin ineffective. As shown in Sect. 2.5, the AAM unit-length
normalizes both the embeddings and speaker prototypes and bounds the target
logit to $[-s,s\cos(m)]$.

On the other hand, VIB does not necessitate unit-length normalization of
neither the embedding nor the prototypes. As discussed in Sect 2.2, the KL
term enforces the aggregated conditional distribution $p(z)=\int
p(z|x;\bm{\theta})p(x)dx$ of the training set to be as close as possible to
$r(z)$ preventing the magnitude (length) of $Z$ from increasing beyond the
(soft) boundaries of the support of $r(z)$. This inherent mechanism of VIB
provides us with flexibility in choosing whether or not to length-normalize
embeddings and/or prototypes, e.g., depending on the benchmark and the scoring
method to be used.

##  3 Experiments

We perform three sets of experiments, NIST-SRE, VoxCeleb, and CNCeleb, to
analyze the performance of VIB regularization. All of them were conducted
using WeSpeaker ([31]) toolkit111https://github.com/wenet-e2e/wespeaker and
closely followed experimental setups of the corresponding recipes including
training set, augmentations, and training hyperparameters like optimizer,
learning rate, etc. The only exception is that we used longer training
examples of 300 frames (instead of the default 200) in all experiments. Also,
we use our custom implementation of scoring backends for NIST-SRE experiments.
When implementing the training with VIB we reuse the same scheduler as used
for the margin in WeSpeaker: for the first 20 epochs, the margin (or $\beta$
in VIB) is set to zero and then exponentially increased to its final value in
the course of the next 20 epochs after which it is kept constant for the rest
of the training.

###  3.1 NIST-SRE

In these experiments, we train all embedding extractors on NIST CTS superset
[32] and test their performance on evaluation sets of three editions of NIST
SRE: SRE2016 [25], telephone condition from SRE2018 [26], and audio-only part
of SRE 2021 [27]. In all cases, the embeddings were centered and length-
normalized. Then, we applied linear discriminant analysis (LDA) reducing the
dimensionality of the embeddings from 256 to 100. Finally, we trained two
scoring backends on preprocessed embeddings: PLDA and T-PSDA. When PLDA is
used for scoring, both speaker and channel subspaces have a dimensionality of
100 (i.e., we use the two-covariance version of the PLDA model); when T-PSDA
is used, we do not attempt to tune its hyperparameters but rather adopt the
same values as found optimal in [28]: we use one 60-dimensional speaker
variable and two 5-dimensional channel variables. The parameters for
centering, LDA, and the backends were estimated on the embedding extractor
training set. The performance is reported in terms of equal error rate (EER)
and minimum cost ($\mathrm{min\\_C}$) as defined by the respective evaluation
plan and computed by the scoring tools provided by the evaluation organizers.

Table 1 shows the comparison of the baseline embedding extractors trained with
regular Softmax/CE or AAM objectives versus the VIB version of the same
objective. In the first column of the table, the number in the brackets
corresponds to the value of the margin $m$ for AAM and to $\beta$ for VIB. The
table is separated into two parts for two different scoring backends allowing
not only comparison between the backends themselves but also showing the
effectiveness of the VIB approach across different scoring methods. By
analyzing the results, we see that VIB consistently outperforms the
corresponding loss not utilizing the margin (VIB vs. CE and VIB_LN vs.
AAM(0.0)) while in most of the cases being competitive to the margin-based
AAM(0.2).

Table 1: Results on NIST SRE evaluation sets. |  | SRE16 yue | SRE16 tgl | SRE18 CMN2 | SRE21 audio
---|---|---|---|---|---
loss | backend | $\mathrm{\mathbf{min\\_C}}$ | EER(%) | $\mathrm{\mathbf{min\\_C}}$ | EER(%) | $\mathrm{\mathbf{min\\_C}}$ | EER(%) | $\mathrm{\mathbf{min\\_C}}$ | EER(%)
CE | PLDA | .470 | 4.50 | .998 | 20.17 | .569 | 7.81 | .564 | 10.54
VIB(0.002) | PLDA | .326 | 3.36 | .951 | 14.79 | .517 | 6.56 | .555 | 10.50
AAM(0.0) | PLDA | .336 | 3.39 | .988 | 15.19 | .502 | 6.55 | .552 | 10.83
AAM(0.2) | PLDA | .336 | 3.34 | .849 | 13.14 | .487 | 6.19 | .562 | 10.29
VIB_LN(0.002) | PLDA | .315 | 3.16 | .984 | 14.91 | .491 | 6.19 | .538 | 9.42
CE | T-PSDA | .441 | 4.80 | .872 | 16.16 | .569 | 7.96 | .573 | 9.98
VIB(0.002) | T-PSDA | .394 | 4.00 | .807 | 12.65 | .541 | 6.55 | .565 | 10.44
AAM(0.0) | T-PSDA | .355 | 3.86 | .742 | 11.37 | .518 | 6.48 | .569 | 10.26
AAM(0.2) | T-PSDA | .399 | 3.96 | .758 | 10.87 | .518 | 6.30 | .599 | 10.37
VIB_LN(0.002) | T-PSDA | .360 | 3.64 | .779 | 12.12 | .525 | 6.29 | .588 | 9.51

###  3.2 VoxCeleb

When running the experiments on VoxCeleb dataset, we follow a commonly adopted
setup: we use the development part of VoxCeleb2 [23] to train the embedding
extractor and the whole VoxCeleb1 [22] for testing. As a backend, we use
simple cosine scoring with only centering and length-normalization of the
embeddings as preprocessing. The centering is done with the mean computed on
VoxCeleb2 development set. The same set was used as a cohort for adaptive
score normalization, where we selected 300 highest scores from the cohort to
estimate normalization parameters.

Table 2 displays the results achieved with the baseline extractors trained
with AAM with margins 0.2 and 0.0 along with the one that uses VIB. Similar to
NIST-SRE case, we observe that VIB consistently outperforms AAM loss with the
margin set to zero, although when the margin is used AAM is clearly superior.
These observations are valid for both cases: whether we use score-
normalization or not. Apart from this, the table shows the large-margin fine-
tuned model along with an analogous setting for VIB: in both cases, the length
of the training examples was increased to 6 seconds for the last 10 epochs
with the margin increased to 0.5 and the value of $\beta$ kept the same as for
the rest of the training (0.004 in this case).

Table 2: Results on VoxCeleb with cosine scoring, without and with score normalization (wo/w). |  | Vox1-O | Vox1-E | Vox1-H
---|---|---|---|---
|  | $\mathrm{\mathbf{minDCF_{0.01}}}$ | EER(%) | $\mathrm{\mathbf{minDCF_{0.01}}}$ | EER(%) | $\mathrm{\mathbf{minDCF_{0.01}}}$ | EER(%)
AAM(0.0) |  | .150/.131 | 1.28/1.11 | .154/.132 | 1.33/1.18 | .245/.192 | 2.55/2.13
AAM(0.2) |  | .115/.090 | 0.96/0.83 | .114/.104 | 0.98/0.89 | .182/.160 | 1.86/1.63
VIB_LN(0.004) |  | .109/.094 | 0.99/0.88 | .130/.115 | 1.11/1.02 | .204/.174 | 2.08/1.82
AAM(0.2)+FT(0.5) |  | .074/.056 | 0.88/0.74 | .101/.092 | 0.95/0.89 | .173/.151 | 1.69/1.53
VIB_LN+FT(0.004) |  | .113/.078 | 0.96/0.80 | .121/.108 | 1.06/0.95 | .194/.160 | 1.96/1.71

###  3.3 CNCeleb

Following WeSpeaker CNCeleb recipe, we use a combination of CNCeleb2 and the
development set of CNCeleb1 as the training set and evaluate on the test set
of CNCeleb1. For multi-enrollment trials, the embeddings for all sessions are
extracted and averaged to obtain a single enrollment embedding.

We performed a similar set of experiments to the ones presented for VoxCeleb
dataset: we compared the performance of several embedding extractors trained
with AAM or VIB with and without finetuning on the long training examples. In
all cases, simple cosine scoring was used, with optionally using score
normalization with the network training set used as a cohort. The results are
given in Table 3 and show that VIB is competitive to AAM loss in most of the
metrics.

Table 3:  Results on the CNCeleb evaluation set, with cosine scoring, without and with score normalization (wo/w). |  | $\mathrm{\mathbf{minDCF_{0.01}}}$ | EER(%)
---|---|---|---
AAM(0.0) |  | .406/.367 | 7.97/7.25
AAM(0.2) |  | .394/.360 | 7.22/6.61
VIB_LN(0.004) |  | .406/.372 | 7.24/6.87
AAM(0.2)+FT(0.5) |  | .393/.356 | 7.23/6.56
VIB_LN+FT(0.004) |  | .399/.364 | 7.29/6.78

###  3.4 Discussion

The experiments we conducted on VoxCeleb and CNCeleb show that our
implementation of VIB is comparable to the state-of-the-art AAM loss,
recovering most of its performance gains compared to AAM with zero margin
(i.e., only length-normalization of embeddings and prototypes). The
experiments on NIST-SRE show that our implementation of VIB is competitive to
the best setting for each evaluation test we examine.

We emphasize that there are numerous directions one may consider to boost the
performance of the method, such as deeper and less mutually-coupled branches
for estimating $f_{\mu}(x)$ and $f_{\sigma}(x)$, different $r(z)$ such as
Gaussian Mixture Models or von Mises-Fisher distribution, and several
optimization techniques that can be found in the rich literature of VIB and
VAE. We should also mention that most of the hyperparameters we used for
training are optimized for the AAM loss, leaving room for further improvement
simply by hyperparameter optimization.

Furthermore, the solid probabilistic/information-theoretic framework of VIB
facilitates research in many directions, such as (a) unsupervised domain
adaptation and multi-domain training (by adapting or conditioning the marginal
$p(z)$ to each domain), (b) uncertainty propagation (by propagating the
uncertainly of the conditional distribution $p(z|x;\bm{\theta})$ in the
probabilistic backend [33, 34]), (c) disentangled speaker embedding extractors
(e.g., [35]) and finally, (d) incorporating the embedding extractor to more
general architectures (e.g., diarization, target-speaker extractor, voice
conversion, joint speaker and speech recognition architecture, a.o.). Joint
architectures typically behave better with well-defined distributions and
losses, and VIB is a step towards this direction.

##  4 Conclusions

In this paper, we examined the strength of VIB in training speaker embedding
extractors. Motivated by the wide adoption of margin-based losses, their
strength in boosting the performance of virtually any architecture, but also
their weaknesses, we tried to address the question of whether or not the VIB
can eventually replace margin-based losses. We experimented with a diverse
collection of speaker recognition benchmarks and used different scoring
methods. Our experiments show that our implementation of VIB yields
improvement over margin-based loss in several NIST-SRE test sets, while it can
recover much of the performance gains of the margin-based loss in VoxCeleb and
CNCeleb. Finally, several directions to improve its performance are provided,
together with other speaker-related applications and settings in which VIB may
be beneficial.

##  5 Acknowledgements

The work was partially supported by Czech Ministry of Interior project No.
VB02000060, 2024-2026 "NABOSO" and Horizon 2020 Marie Sklodowska-Curie grant
ESPERANTO, No. 101007666. Computing on IT4I supercomputer was supported by the
Ministry of Education, Youth and Sports of the Czech Republic through e-INFRA
CZ (ID:90254).

## References

  * [1] G. Heigold, I. Moreno, S. Bengio, and N. Shazeer, ``End-to-end text-dependent speaker verification,'' in _2016 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2016, pp. 5115–5119.
  * [2] D. Snyder, D. Garcia-Romero, G. Sell, D. Povey, and S. Khudanpur, ``X-vectors: Robust dnn embeddings for speaker recognition,'' in _2018 IEEE international conference on acoustics, speech and signal processing (ICASSP)_. IEEE, 2018, pp. 5329–5333.
  * [3] B. Desplanques, J. Thienpondt, and K. Demuynck, ``ECAPA-TDNN: Emphasized Channel Attention, Propagation and Aggregation in TDNN Based Speaker Verification,'' in _Interspeech_ , 2020.
  * [4] K. A. Lee, Q. Wang, and T. Koshinaka, ``Xi-vector embedding for speaker recognition,'' _IEEE Signal Processing Letters_ , vol. 28, pp. 1385–1389, 2021.
  * [5] N. R. Koluguri, T. Park, and B. Ginsburg, ``Titanet: Neural model for speaker representation with 1d depth-wise separable convolutions and global context,'' in _ICASSP 2022-2022 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2022, pp. 8102–8106.
  * [6] H. Zeinali, S. Wang, A. Silnova, P. Matějka, and O. Plchot, ``But system description to voxceleb speaker recognition challenge 2019,'' _arXiv preprint arXiv:1910.12592_ , 2019.
  * [7] M. Rybicka, J. Villalba, P. Zelasko, N. Dehak, and K. Kowalczyk, ``Spine2net: Spinenet with res2net and time-squeeze-and-excitation blocks for speaker recognition.'' in _Interspeech_ , 2021, pp. 496–500.
  * [8] T. Zhou, Y. Zhao, and J. Wu, ``Resnext and res2net structures for speaker verification,'' in _2021 IEEE Spoken Language Technology Workshop (SLT)_. IEEE, 2021, pp. 301–307.
  * [9] S. Novoselov, G. Lavrentyeva, A. Avdeeva, V. Volokhov, and A. Gusev, ``Robust speaker recognition with transformers using wav2vec 2.0,'' _arXiv preprint arXiv:2203.15095_ , 2022.
  * [10] J. Peng, T. Stafylakis, R. Gu, O. Plchot, L. Mošner, L. Burget, and J. Černockỳ, ``Parameter-efficient transfer learning of pre-trained transformer models for speaker verification using adapters,'' in _ICASSP 2023-2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2023, pp. 1–5.
  * [11] J. Deng, J. Guo, N. Xue, and S. Zafeiriou, ``ArcFace: Additive Angular Margin Loss for Deep Face Recognition,'' in _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR)_ , June 2019.
  * [12] X. Xiang, S. Wang, H. Huang, Y. Qian, and K. Yu, ``Margin matters: Towards more discriminative deep neural network embeddings for speaker recognition,'' in _2019 Asia-Pacific Signal and Information Processing Association Annual Summit and Conference (APSIPA ASC)_. IEEE, 2019, pp. 1652–1656.
  * [13] A. Silnova, T. Stafylakis, L. Mošner, O. Plchot, J. Rohdin, P. Matějka, L. Burget, O. Glembek, and N. Brümmer, ``Analyzing speaker verification embedding extractors and back-ends under language and channel mismatch,'' in _Odyssey 2022: The speaker and Language Recongnition Workshop, Beijing_ , 2022. [Online]. Available: https://arxiv.org/abs/2203.10300
  * [14] J. Villalba, N. Chen, D. Snyder, D. Garcia-Romero, A. McCree, G. Sell, J. Borgstrom, L. P. García-Perera, F. Richardson, R. Dehak _et al._ , ``State-of-the-art speaker recognition with neural network embeddings in nist sre18 and speakers in the wild evaluations,'' _Computer Speech & Language_, vol. 60, p. 101026, 2020.
  * [15] A. A. Alemi, I. Fischer, J. V. Dillon, and K. Murphy, ``Deep variational information bottleneck,'' in _International Conference on Learning Representations_ , 2016.
  * [16] Y. Eom, Y. Lee, J. Um, and H.-R. Kim, ``Anti-spoofing using transfer learning with variational information bottleneck,'' in _23rd Annual Conference of the International Speech Communication Association, INTERSPEECH 2022_. ISCA, 2022, pp. 3568–3572.
  * [17] Y. Du, J. Xu, H. Xiong, Q. Qiu, X. Zhen, C. G. Snoek, and L. Shao, ``Learning to learn with variational information bottleneck for domain generalization,'' in _Computer Vision–ECCV 2020: 16th European Conference, Glasgow, UK, August 23–28, 2020, Proceedings, Part X 16_. Springer, 2020, pp. 200–216.
  * [18] S. Cui, J. Cao, X. Cong, J. Sheng, Q. Li, T. Liu, and J. Shi, ``Enhancing multimodal entity and relation extraction with variational information bottleneck,'' _IEEE/ACM Transactions on Audio, Speech, and Language Processing_ , 2024.
  * [19] T. Gu, G. Xu, and J. Luo, ``Sentiment analysis via deep multichannel neural networks with variational information bottleneck,'' _IEEE Access_ , vol. 8, pp. 121 014–121 021, 2020.
  * [20] S. Bang, P. Xie, H. Lee, W. Wu, and E. Xing, ``Explaining a black-box by using a deep variational information bottleneck approach,'' in _Proceedings of the AAAI conference on artificial intelligence_ , vol. 35, no. 13, 2021, pp. 11 396–11 404.
  * [21] D. Wang, Y. Dong, Y. Li, Y. Zi, Z. Zhang, X. Li, and S. Xiong, ``Variational information bottleneck based regularization for speaker recognition.'' in _Interspeech_ , 2021, pp. 1054–1058.
  * [22] A. Nagrani, J. S. Chung, W. Xie, and A. Zisserman, ``Voxceleb: Large-scale speaker verification in the wild,'' _Computer Speech and Language_ , 2019.
  * [23] J. Chung, A. Nagrani, and A. Zisserman, ``Voxceleb2: Deep speaker recognition,'' _Interspeech 2018_ , 2018.
  * [24] L. Li, R. Liu, J. Kang, Y. Fan, H. Cui, Y. Cai, R. Vipperla, T. F. Zheng, and D. Wang, ``Cn-celeb: multi-genre speaker recognition,'' _Speech Communication_ , vol. 137, pp. 77–91, 2022.
  * [25] S. O. Sadjadi, T. Kheyrkhah, A. Tong, C. Greenberg, D. Reynolds, E. Singer, L. Mason, and J. Hernandez-Cordero, ``The 2016 nist speaker recognition evaluation,'' _Interspeech 2017_ , 2017.
  * [26] O. Sadjadi, C. Greenberg, E. Singer, D. Reynolds, L. Mason, and J. Hernandez-Cordero, ``The 2018 nist speaker recognition evaluation.'' INTERSPEECH, Graz, AT, 2019-09-15 00:09:00 2019. [Online]. Available: https://tsapps.nist.gov/publication/get_pdf.cfm?pub_id=927673
  * [27] S. O. Sadjadi, C. Greenberg, E. Singer, L. Mason, and D. Reynolds, ``The 2021 nist speaker recognition evaluation,'' _arXiv preprint arXiv:2204.10242_ , 2022.
  * [28] A. Silnova, N. Brümmer, A. Swart, and L. Burget, ``Toroidal probabilistic spherical discriminant analysis,'' in _ICASSP 2023-2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2023, pp. 1–5.
  * [29] N. Tishby, F. C. Pereira, and W. Bialek, ``The information bottleneck method,'' in _The 37th annual Allerton Conference on Communication, Control, and Computing_ , 1999, pp. 368–377.
  * [30] D. P. Kingma and M. Welling, ``Auto-encoding variational bayes,'' in _In Proceedings of the 2nd International Conference on Learning Representations (ICLR)_ , 2014.
  * [31] H. Wang, C. Liang, S. Wang, Z. Chen, B. Zhang, X. Xiang, Y. Deng, and Y. Qian, ``Wespeaker: A research and production oriented speaker embedding learning toolkit,'' in _ICASSP 2023-2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2023, pp. 1–5.
  * [32] O. Sadjadi, ``NIST SRE CTS Superset: A large-scale dataset for telephony speaker recognition,'' 2021-08-16 04:08:00 2021. [Online]. Available: https://tsapps.nist.gov/publication/get_pdf.cfm?pub_id=933116
  * [33] P. Kenny, T. Stafylakis, P. Ouellet, M. J. Alam, and P. Dumouchel, ``Plda for speaker verification with utterances of arbitrary duration,'' in _2013 IEEE International Conference on Acoustics, Speech and Signal Processing_. IEEE, 2013, pp. 7649–7653.
  * [34] Q. Wang, K. A. Lee, and T. Liu, ``Incorporating uncertainty from speaker embedding estimation to speaker verification,'' in _ICASSP 2023-2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_. IEEE, 2023, pp. 1–5.
  * [35] T. Liu, K. A. Lee, Q. Wang, and H. Li, ``Disentangling voice and content with self-supervision for speaker recognition,'' _Advances in Neural Information Processing Systems_ , vol. 36, 2024.