# Multimode Mamyshev Oscillator

Henry Haig  School of Applied and Engineering Physics, Cornell University,
Ithaca, New York 14853, USA  Corresponding author: tsh67cornell.edu  Pavel
Sidorenko  School of Applied and Engineering Physics, Cornell University,
Ithaca, New York 14853, USA  Anirban Dhar  Fiber Optics and Photonics
Division, CSIR-Central Glass and Ceramic Research Institute, Kolkata-700032,
India  Nilotpal Choudhury  Fiber Optics and Photonics Division, CSIR-Central
Glass and Ceramic Research Institute, Kolkata-700032, India  Academy of
Scientific and Innovative Research (AcSIR), Gaziabad 21002, India  Ranjan Sen
Fiber Optics and Photonics Division, CSIR-Central Glass and Ceramic Research
Institute, Kolkata-700032, India  Demetrios Christodoulides  CREOL/College of
Optics and Photonics, University of Central Florida, Orlando, Florida 32816,
USA  Frank Wise  School of Applied and Engineering Physics, Cornell
University, Ithaca, New York 14853, USA

###### Abstract

We present a spatiotemporally mode-locked Mamyshev oscillator. A wide variety
of multimode mode-locked states, with varying degrees of spatiotemporal
coupling, are observed. We find that some control of the modal content of the
output beam is possible through the cavity design. Comparison of simulations
to experiments indicates that spatiotemporal mode-locking is enabled by
nonlinear intermodal interactions and spatial filtering, along with the
Mamyshev mechanism. This work represents a first exploration of spatiotemporal
mode-locking in an oscillator with the Mamyshev saturable absorber.

††journal: ol

##  1 Introduction

Mode-locked fiber lasers attract considerable interest in science, industry,
and healthcare due to their compactness, efficiency, and inherent alignment-
free design. However, fiber lasers lag behind their solid-state counterparts
in pulse energy owing mainly to strong nonlinear processes in the waveguide
medium. While nonlinearity enables mode-locking, it eventually destabilizes
pulse evolutions and constrains performance. Consequently, researchers explore
methods to benefit from nonlinearity as well as ways to reduce nonlinear
effects.

One recent innovation in exploiting nonlinearity is use of the Mamyshev
mechanism [1] in fiber oscillators. Mamyshev oscillators contain two offset
spectral filters [2, 3, 4] that extinguish continuous-wave (CW) lasing, but
allow pulses with sufficient peak power to traverse the cavity through
nonlinear spectral broadening between filters. This enacts strong effective
saturable absorbtion [5] which can stabilize highly nonlinear evolutions.
Mamyshev oscillators achieve pulse energies (> 200 nJ in single-mode fiber
[6], > 1 µJ in photonic crystal fiber [7]) and durations (< 40 fs) [6] far
beyond those of fiber lasers with other saturable absorbers.

Increasing fiber core size is another technique to scale energy, by reducing
nonlinearity. This method is limited in principle for single-mode (SM) lasers
by the higher-order modes that exist in large-core fiber. However, the
demonstration of spatiotemporal mode-locking (STML) in multimode (MM) fiber
[8] showed that it is possible to incorporate higher-order modes in mode-
locking. The limits of STML are unknown, and the broad range of nonlinear
multimode propagation phenomena [9] may provide means to generate light with
unprecedented degrees of spatiotemporal control. Some STML lasers build on
well-known SM phenomena, such as dissipative soliton [8] and self-similar
evolutions [10]. Others exploit spatiotemporal degrees of freedom in entirely
new ways, such as beam self-cleaning [11].

Much remains to be explored in MM oscillators. Previous work on STML has
highlighted the importance of linear filtering and nonlinear intermodal
interactions, but the effect of saturable absorption on MM pulse evolution is
not fully understood. To date, STML lasers have used saturable absorbers based
on material absorption [12] or nonlinear polarization evolution [13, 10, 8,
11]. The Mamyshev saturable absorber has not been applied to STML. Most
generally, one would expect that a saturable absorber with deep modulation
would be attractive for STML. However, it is unknown how the Mamyshev
mechanism will affect MM evolutions in general, and whether it can provide
similar or alternative benefits to MM oscillators as it has to SM oscillators.

Here we present a first demonstration of spatiotemporal mode-locking in a
Mamyshev oscillator. The laser exhibits a wide variety of STML states, ranging
from compressible single pulses to spatiotemporally-complex pulse bursts. In
some states, the output profile can be partially controlled by the cavity
alignment. Numerical simulations that account for the main features of the
experimental observations will also be presented.

Figure 1: Cavity Schematic of the Multimode Mamyshev Oscillator. DM, dichroic
mirror; MM, multimode fiber; L, lens; BS, beam sampler; HWP, half-wave plate;
PBS, polarizing beam splitter; OC, output coupling; SF, spectral filter; COL,
collimating lens; ISO, isolator; C, combiner; SM, single mode fiber. Green
indicates Yb doped gain fiber, and black indicates passive fiber.

Fig. 1 is a schematic of the laser cavity. The ring cavity comprises one SM
arm and one MM arm separated by offset spectral filters (4 nm full-width-at-
half-maximum (FWHM) bandwidth and 5 - 10 nm center wavelength separation). The
filters, along with nonlinear spectral broadening in either arm, create an
effective saturable absorber. The SM arm includes 2 meters of passive fiber
and 1.5 meters of active fiber. Both segments are polarization-maintaining
(PM) to suppress polarization dynamics. 8.5 µm core diameter fiber is used to
deliver sufficient average power to the MM arm to saturate the gain in that
segment. The multimode arm consists of 1.5 - 2 m of 30 µm core diameter Yb-
doped graded-index fiber. GRIN MM fiber is chosen for its small intermodal
group velocity differences. Both gain fibers are pumped by wavelength-
stabilized MM pump diodes at 976 nm. Waveplates and a beam splitter provide MM
output coupling and align the polarization with the PM fiber’s axis. These
design choices were motivated by the goal of limiting the complexity of the
evolution, as observed in initial simulations of the cavity.

The GRIN MM fiber was fabricated through the vapor phase chelate delivery
technique [14] by optimization of different fabrication parameters such as
vapor-phase composition, sublimator temperature, deposition temperature, and
number of deposition passes. The vapor phase composition was specifically
adjusted in different core layers to achieve a graded profile along the radial
direction which finally yielded a fiber with numerical aperture of 0.14 and
maximum Yb-ion-concentration of 4000 ppm. The fiber supports 20 transverse
modes of a single linear polarization at 1030 nm according to our mode
calculations [15].

For the experiments described here, the oscillator is started by seeding with
pulses from an ANDi fiber oscillator [16] with 20 nm bandwidth centered at
1040 nm and 1.3 ps (chirped) duration. The seed pulses are coupled to the
cavity by reflection from the spectral filter preceeding the MM arm. While
aligning the laser, we measure the spatially-integrated MM beam spectrum with
a diffuser and highly-MM fiber coupled spectrometer. The temporal behavior is
monitored with a photodiode and oscilloscope, and the near-field MM beam
profile is measured with a camera. When mode-locking occurs, the mode-locked
spectrum and pulse train dominate over the seed spectrum and pulse train. The
seed beam is then blocked for self-sustaining pulsed operation.

Figure 2: Typical RF spectra of an STML state. (a) shows overlapping RF
spectra signal spikes at the cavity fundamental repetition rate (23 MHz)
within a 5 KHz window for different parts of the beam. Colored data correspond
to the crosshairs in (b). Black curve is the background noise floor. (b) shows
the near-field beam profile. The white circle is the core-cladding boundary.
The intensity colormap and scale here are used for all other images presented
in this letter. (c) shows the RF spectrum of the center section of the beam
with harmonics up to the response limit of the photodetector (1 GHz).

To characterize the laser, we use a pinhole to isolate spatial sections of the
beam (< 10% of the beam) for measurements. The laser produces a wide variety
of pulses with various spectral bandwidths (7 – 30 nm), spatial profiles,
energies (7 – 20 nJ), and spatiotemporal complexity. In every mode-locked
state we observe, the RF spectrum is close to uniform across the beam, as
shown in Fig. 2. Contrast between the fundamental frequency and secondary
modulations over 70 dB indicates stable mode locking at the fundamental cavity
repetition rate. We also observe some states for which the RF spectrum and
photodiode trace show multi-pulsing, albeit always uniformly across the beam.
We use autocorrelations based on second-harmonic generation (SHG) to determine
the pulse duration. SHG autocorrelations generally cannot accurately measure
spatiotemporally-coupled beams [17] because of the limited angular acceptance
bandwidth. However, we find that sufficiently small sections of the beam give
autocorrelations that are insensitive to the phase-matching condition, and so
can be measured with this technique. The output pulses are chirped, with
spatially-integrated duration around 1 ps. All mode-locked states that we
observe so far exhibit spatio-spectral complexity, meaning the optical
spectrum varies across the beam. Below, we present two representative STML
states, with minimal and maximal spatiotemporal coupling, respectively.

Measurements of a single-pulse mode-locked state with minimal variation of the
temporal profile across the beam are summarized in Fig. 3. Despite spatio-
spectral variation, the autocorrelation is nearly identical for any section of
the beam. These pulses compress reasonably well by a standard grating
compressor to a FWHM duration of 300 fs, with a transform-limited FWHM
duration of 200 fs. We estimate that under 20% of the total pulse energy is
contained in pedestal structure. As a secondary point, the pulses from the SM
arm (20 nJ with 60 nm FWHM bandwidth) are similar to those from SM Mamyshev
oscillators [6] [18].

Figure 3: Measurements of 15-nJ pulses with weak spatiotemporal coupling. (a)
Near-field beam profile with core boundary in white and fundamental mode cross
section inset at bottom. Measurements in (b), (c), and (d) are taken by
isolating sections of the beam indicated by corresponding colored crosshairs
in (a). (b) normalized spectra, with the spatially-integrated spectrum in
purple. (c) Autocorrelations of pulses directly from laser. (d)
Autocorrelation of compressed pulses. The optimal dispersion for compression
is the same for each section. Figure 4: Measurements of a highly-
spatiotemporally-coupled pulse burst of 20 nJ total energy. (a), (b), and (c)
correspond to the same measurements as in Fig. 3.

In other cases, such as the state shown in Fig. 4, the autocorrelation varies
considerably with the chosen portion of the beam, which is evidence of strong
spatiotemporal coupling. Typically, these states exhibit multiple pulses with
picosecond separation, and the effect of a grating compressor depends strongly
on the chosen section of the beam. These states may be undesirable from an
applications standpoint, since they would require complex spatiotemporal phase
compensation to compress.

Attempts to increase the pulse energy by increasing pump power usually result
in formation of multiple pulses with picosecond separation. Increasing filter
separation is an established technique to prevent multipulsing or increase the
single-pulse energy of Mamyshev oscillators [19] [20] by increasing the
saturation power of the saturable absorber. In this laser, however, we find
that increasing the filter separation past 10 nm disrupts mode-locking before
it prevents multi-pulsing.

The beam profile prior to mode-locking is determined by the coupling of
filtered ASE from the SM arm to the MM arm. We observe qualitatively that this
coupling usually determines the spatial profile of the mode-locked state. This
trend is shown in Fig. 5(a) and (b), where the spatial profile before seeding
and mode-locking strongly resembles that after seeding and mode-locking. Less
often, we observe states for which the profiles are distinct (Fig. 5(c) and
(d)). This phenomenon is not fully understood, but we speculate that the modal
energy distribution in the MM arm is primarily determined by the SM to MM
coupling.

Figure 5: Each pair of images shows the multimode near field beam profile
before (left) and after (right) mode-locking.

Numerical simulations were performed to help understand the MM Mamyshev
oscillator. Drawing precise comparisons between experiment and simulation is
difficult in STML lasers, but some features of a simulated cavity help explain
the MM pulse evolution in this oscillator and align qualitatively with our
experimental observations. To limit the time it takes to run simulations, we
include up to 6 linearly polarized modes. The simulation [21] includes the
effects of Kerr and Raman nonlinearities, self-steepening, and second- and
third- order dispersion. We model the Yb gain by solving the rate equations in
the steady state [22] [23]. To model coupling into and out of the MM arm, we
assume a 1:1 imaging system between the end faces of the SM and MM fibers, and
calculate overlap integrals of each MM fiber mode with the SM fiber mode. This
includes spatial offset between fiber centers in order to excite or couple
modes that are not radially-symmetric. Simulations of the same cavity with
different input and output offsets show that the input offset primarily
determines the modal energy content of the amplified MM pulse despite energy
exchange between modes during propagation. This is consistent with our
experimental observation that the SM to MM coupling controls the modal energy
distribution. However, the general features of the evolution are similar for
all input and output offsets up to 8.5 µm.

Fig. 6 shows various parameters of a typical simulated MM evolution for this
cavity. Nonlinear spectral broadening occurs due to self-phase modulation
(SPM) and cross-phase modulation (XPM) among multiple modes. XPM contributes
more significantly to nonlinear phase than SPM [24][25], so the spectral
broadening per mode is greater than it would be for SM propagation of pulses
with the same energy. As shown in Fig. 6(e) and (d), pulses in several modes
spectrally broaden sufficiently to reach the offset spectral filter, despite
having little of the overall MM pulse energy.

Previous work on STML lasers [8][17] emphasized the importance of spatial and
spectral filtering to reset temporal walkoff between modes owing to intermodal
group velocity differences. Our simulations indicate that the spectral filter
following the MM arm does not decrease the temporal walkoff between pulses in
different modes. Instead, temporal walk-off that arises from MM propagation
and spectral filtering is corrected by the spatial filtering that occurs on
coupling into the SM fiber (Fig. 6(f)).

In general, MM oscillators and amplifiers are not yet understood well enough
to be fully controlled or employed as practical devices. The results presented
here are a first step toward filling this knowledge gap for the case of
oscillators or regenerative amplifiers that employ the Mamyshev mechanism. The
large space of spatiotemporally-mode-locked states should be a fertile ground
for scientific exploration. Future goals will include demonstration of control
over the modal content of such lasers, scaling of the pulse energy to much-
higher values, and generation of high-energy pulses with high-quality beams.

In conclusion, we have constructed and characterized an STML laser with the
Mamyshev mechanism as saturable absorber. The design features a SM arm in the
cavity that supports the MM evolution by providing spatial filtering,
compensating modal dispersion, and exerting some control over the modal
content in the MM arm. The oscillator produces a wide variety of MM pulses
with varying degrees of spatiotemporal coupling. Simulations and experiments
indicate that nonlinear intermodal interactions, spatial filtering, and the
Mamyshev mechanism together stabilize the evolution.

Figure 6: Converged simulated MM pulse evolution. (a) Intensity of a
spatially-integrated MM pulse (black) at the end of MM gain fiber, along with
pulse in each mode with colors indicated by the colorbar in (d). (b) Power
spectra of the spatially-integrated (black) pulse and individual modes
(colored). (c) time-integrated spatial intensity as would be captured by a
camera. (d) and (e) show energy and RMS bandwidth, respectively, per mode with
evolution coordinate. The blue region captioned SF represents the spectral
filter. (f) shows the time-center coordinate (first moment) of the MM and SM
pulses for a complete cavity round-trip, relative to the simulation frame.
Regions in green, blue, and white represent gain, spectral filtering, and
passive propagation respectively. The time center coordinate is reset during
each round trip. The colors used for the modes are the same as in (a) and (b).

Funding National Science Foundation (ECCS-1912742); Office of Naval Research
(N00014-20-1-2789); Israeli Ministry of Defense (4441069183); Simons
Foundation (733682). Disclosures The authors declare no conflicts of interest.
Data Availability Data underlying the results presented in this paper are not
publicly available at this time but may be obtained from the authors upon
reasonable request.

## References

  * [1] P. Mamyshev, _24th European Conference on Optical Communication. ECOC ’98 (IEEE Cat. No.98TH8398)_ , vol. 1 (1998).
  * [2] K. Sun, M. Rochette, and L. R. Chen, Opt. Express 17, 10419 (2009).
  * [3] M. Piche, _Mode locking through nonlinear frequency broadening and spectral filtering_ , vol. 2041 (SPIE, 1994).
  * [4] K. Regelskis, J. Želudevičius, K. Viskontas, and G. Račiukaitis, Opt. Lett. 40, 5255 (2015).
  * [5] S. Pitois, C. Finot, L. Provost, and D. J. Richardson, J. Opt. Soc. Am. B 25, 1537 (2008).
  * [6] P. Sidorenko, W. Fu, L. G. Wright, M. Olivier, and F. W. Wise, Opt. Lett. 43, 2672 (2018).
  * [7] W. Liu, R. Liao, J. Zhao, J. Cui, Y. Song, C. Wang, and M. Hu, Optica 6, 194 (2019).
  * [8] L. G. Wright, D. N. Christodoulides, and F. W. Wise, Science 358, 94 (2017).
  * [9] K. Krupa, A. Tonello, A. Barthélémy, T. Mansuryan, V. Couderc, G. Millot, P. Grelu, D. Modotto, S. A. Babin, and S. Wabnitz, APL Photonics 4, 110901 (2019).
  * [10] U. Teğin, E. Kakkava, B. Rahmani, D. Psaltis, and C. Moser, Optica 6, 1412 (2019).
  * [11] U. Tegin, B. Rahmani, E. Kakkava, D. Psaltis, and C. Moser, Advanced Photonics 2, 1 (2020).
  * [12] Z.-L. Ma, J.-G. Long, W. Lin, W.-Y. Hong, H. Cui, Z.-C. Luo, W.-C. Xu, M. E. Likhachev, S. S. Aleshkina, V. M. Mashinsky, M. V. Yashkov, and A.-P. Luo, Opt. Express 29, 9465 (2021).
  * [13] Y. Ding, X. Xiao, K. Liu, S. Fan, X. Zhang, and C. Yang, Phys. Rev. Lett. 126, 093901 (2021).
  * [14] N. Choudhury, N. K. Shekhar, A. Dhar, and R. Sen, physica status solidi (a) 216, 1900365 (2019).
  * [15] A. B. Fallahkhair, K. S. Li, and T. E. Murphy, Journal of Lightwave Technology 26, 1423 (2008).
  * [16] A. Chong, J. Buckley, W. Renninger, and F. Wise, Opt. Express 14, 10095 (2006).
  * [17] L. G. Wright, P. Sidorenko, H. Pourbeyram, Z. M. Ziegler, A. Isichenko, B. A. Malomed, C. R. Menyuk, D. N. Christodoulides, and F. W. Wise, Nature Physics 16, 565 (2020).
  * [18] P. Sidorenko, W. Fu, and F. Wise, Optica 6, 1328 (2019).
  * [19] Y.-H. Chen, P. Sidorenko, R. Thorne, and F. Wise, J. Opt. Soc. Am. B 38, 743 (2021).
  * [20] J. Želudevičius, M. Mickus, and K. Regelskis, Opt. Express 26, 27247 (2018).
  * [21] L. G. Wright, Z. M. Ziegler, P. M. Lushnikov, Z. Zhu, M. A. Eftekhar, D. N. Christodoulides, and F. W. Wise, IEEE Journal of Selected Topics in Quantum Electronics 24, 1 (2018).
  * [22] S. K. Turitsyn, A. E. Bednyakova, M. P. Fedoruk, A. I. Latkin, A. A. Fotiadi, A. S. Kurkov, and E. Sholokhov, Opt. Express 19, 8394 (2011).
  * [23] R. Paschotta, J. Nilsson, A. Tropper, and D. Hanna, IEEE Journal of Quantum Electronics 33, 1049 (1997).
  * [24] G. P. Agrawal, _Nonlinear Fiber Optics_ (Academic Press, 1995), 2nd ed.
  * [25] F. Poletti and P. Horak, J. Opt. Soc. Am. B 25, 1645 (2008).

bibliography


# EXoN: EXplainable encoder Network

SeungHwan An1, Hosik Choi2, Jong-June Jeon1

###### Abstract

We propose a new semi-supervised learning method of Variational AutoEncoder
(VAE) which yields a customized and explainable latent space by EXplainable
encoder Network (EXoN). Customization means a manual design of latent space
layout for specific labeled data. To improve the performance of our VAE in a
classification task without the loss of performance as a generative model, we
employ a new semi-supervised classification method called SCI (Soft-label
Consistency Interpolation). The classification loss and the Kullback-Leibler
divergence play a crucial role in constructing explainable latent space. The
variability of generated samples from our proposed model depends on a specific
subspace, called activated latent subspace. Our numerical results with MNIST
and CIFAR-10 datasets show that EXoN produces an explainable latent space and
reduces the cost of investigating representation patterns on the latent space.

##  1 Introduction

Variational AutoEncoder (VAE) (Kingma and Welling 2013; Rezende, Mohamed, and
Wierstra 2014) aims to reconstruct a well-representing latent space and
recover an original observation well. In general, the probabilistic model used
in VAE is parameterized by neural networks defined as two maps from the domain
of observations to the latent space and vice versa. However, since the
probability model of the VAE is not written in a closed form, the maximum
likelihood method is unsuitable for application to the VAE. As an alternative,
the Variational Bayesian method is popularly applied to the model to maximize
the Evidence of Lower Bound (ELBO) (Jordan et al. 1999).

Plenty of semi-supervised learning methods for the VAE (Kingma et al. 2014;
Maaløe et al. 2016; Siddharth et al. 2017; Maaløe, Fraccaro, and Winther 2017;
Li et al. 2019; Feng et al. 2021; Hajimiri, Lotfi, and Baghshah 2021) have
been introduced, and especially (Maaløe, Fraccaro, and Winther 2017; Hajimiri,
Lotfi, and Baghshah 2021) applied the mixture prior distribution such that the
VAE model provides explainable latent space according to labels. However,
existing semi-supervised VAE models still have practical limitations.

(Kingma et al. 2014; Maaløe et al. 2016; Siddharth et al. 2017; Li et al.
2019; Feng et al. 2021) introduced an additional latent space representing
labels and assumed the probabilistic independence of the label and the other
latent variables. While the model is simply formulated, the trained latent
space cannot provide explainable and measurable quantitative information used
to generate a new image by interpolating between two images. So, it is
difficult to impose structural restrictions on the latent space, such as the
proximity of latent features across some labeled observations. Also, the
trained latent space differs according to the training process, so the latent
space is not explained consistently even with the same dataset (Maaløe,
Fraccaro, and Winther 2017; Hajimiri, Lotfi, and Baghshah 2021). For example,
it is difficult to obtain information about interpolated images from latent
space before the trained latent space is entirely investigated by
observations. Thus, a manual design of the latent space is required for an
explainable representation model.

The current study proposes a new semi-supervised VAE model, wherein we have
incorporated the manual design to improve the clarity of the model. The model
employs a Gaussian mixture for prior and posterior latent distributions
(Dilokthanakul et al. 2016; Zheng and Sun 2019; Willetts, Roberts, and Holmes
2019; Mathieu et al. 2019; Guo et al. 2020) and focuses on constructing an
explainable encoder of VAE. In our model, the latent space is customized for
each component’s correspondence to a specific label and proximity. Since the
latent space is shared with the decoder and classifier, each centroid of the
mixture distribution is trained as an identifier representing a specific
label. Therefore, the latent space is manually decomposed by labels, and a
user can obtain a customized explainable latent space.

In addition, we propose a measure of the representation power of latent
variables with which the importance of latent subspace can be investigated. It
is found that the encoder of our VAE model selectively activates only a part
of the latent space, and the latent subspace represents the characteristics of
generated samples. Because the activation is measured by a posterior variance
linked with the encoder, we call the encoder EXoN (EXplainable encoder
Network) by borrowing the term in the field of gene biology.

This paper is organized as follows. Section 2 briefly introduces Variational
AutoEncoder, and Section 3 proposes our VAE model, including its derivation.
Section 4 shows the results of numerical simulations. Concluding remarks
follow in Section 5.

##  2 Preliminary

###  2.1 VAE for Unsupervised Learning

Let ${\bf x}\in\mathcal{X}\subset\mathbb{R}^{m}$ and ${\bf
z}\in\mathcal{Z}\subset\mathbb{R}^{d}$ with $d<m$ be the observed data and the
latent variable. $p({\bf x})$ and $p({\bf z})$ denote the probability density
functions (pdf) of ${\bf x}$ and ${\bf z}$. Let $p({\bf x}|{\bf z})$ be the
conditional pdf of ${\bf x}$ for a given ${\bf z}$ and ${\bf x}|{\bf z}\sim
N(D({\bf z}),\beta\cdot{\bf I}),$ where $D:\mathcal{Z}\mapsto\mathcal{X}$,
${\bf I}$ is the $m\times m$ identity matrix and the observation noise
$\beta>0$. $D({\bf z})$ is a neural network with parameter $\theta$. To
emphasize dependence on parameters, denote $D({\bf z})$ and $p({\bf x}|{\bf
z})$ by $D({\bf z};\theta)$ and $p({\bf x}|{\bf z};\theta,\beta)$,
respectively. $p({\bf x}|{\bf z};\theta,\beta)$ is referred to as the decoder
of VAE. The term of decoder is also used as the map ${\bf z}\mapsto(D({\bf
z};\theta),\beta)$, because $p({\bf x}|{\bf z};\theta,\beta)$ is fully
parameterized by $(D({\bf z};\theta),\beta)$.

The parameter of VAE $(\theta,\beta)$ is estimated by Variational method
(Kingma and Welling 2013; Rezende, Mohamed, and Wierstra 2014) that maximizes
ELBO. The ELBO of $\log p({\bf x};\theta,\beta)$ is obtained from the
inequality

$\displaystyle\log p({\bf x};\theta,\beta)$ $\displaystyle\geq$
$\displaystyle\mathbb{E}_{q({\bf z}|{\bf x})}[\log p({\bf x}|{\bf
z};\theta,\beta)]-\mathcal{KL}\left(q({\bf z}|{\bf x})\|p({\bf z})\right),$

where $\mathcal{KL}(q\|p)$ denotes Kullback-Leibler divergence from $p$ to
$q$. (Kingma and Welling 2013; Rezende, Mohamed, and Wierstra 2014) employ a
neural network with parameter $\phi$ as $q({\bf z}|{\bf x};\phi)$, and obtain
the objective function under finite samples $x_{1},\cdots,x_{n}$,

$\displaystyle\sum_{i=1}^{n}\mathbb{E}_{q}[\log p(x_{i}|{\bf
z};\theta,\beta)]-\mathcal{KL}\left(q({\bf z}|x_{i};\phi)\|p({\bf z})\right).$
(2)

Here, $q({\bf z}|{\bf x};\phi)$ is referred to as the encoder of VAE or the
posterior distribution over the latent variables. Multivariate Gaussian
distribution is a popular choice for $q({\bf z}|{\bf x};\phi)$ in which the
mean and covariance are given by neural network models. Thus, the parameters
of VAE consist of $\phi$ in the encoder and $(\theta,\beta)$ in the decoder.

In practice, the VAE is fitted by maximizing a stochastically approximated
ELBO for (2) with respect to $(\theta,\beta,\phi)$. The approximated ELBO is
given by

$\displaystyle\sum_{i=1}^{n}-\frac{1}{2}\|x_{i}-D(z_{i};\theta)\|^{2}-\beta\cdot\mathcal{KL}(q({\bf
z}|x_{i};\phi)\|p({\bf z}))$ $\displaystyle-\frac{nm}{2}\log 2\pi\beta$ (3)

where $z_{i}$ is a sample from $q({\bf z}|x_{i};\phi)$ and $\|\cdot\|$ is
$l_{2}$-norm. The VAE is fitted by maximizing (2.1) with respect to
$(\theta,\beta,\phi)$ (Lucas et al. 2019). The first term in (2.1) is the
precision of generated samples. Since the KL-divergence is always non-
negative, (Higgins et al. 2016; Mathieu et al. 2019; Li et al. 2019) explain
$\beta$ as a tuning parameter regularizing $\theta$ and $\phi$, and the last
term of (2.1) is considered as constant.

##  3 Proposed Model

We propose a new VAE model with a customized explainable latent space where a
conceptual center of a latent distribution for a specific label can be freely
assigned. Let ${\bf y}\in\mathcal{Y}=\\{1,\cdots,K\\}$ be a discrete random
variable and denote the joint probability density of $({\bf x},{\bf y})$ by
$p({\bf x},{\bf y})$. In this paper, high dimensional observations ${\bf x}$
are partly labeled by ${\bf y}$. Denote the sets of indices for labeled and
unlabeled samples by $I_{L}$ and $I_{U}$.

###  3.1 Model Assumptions

The prior latent distribution for each label is assumed as ${\bf z}|{\bf
y}=k\sim N(m_{k},diag\\{s_{k}^{2}\\})$ and $p({\bf y}=k)=w_{k}$ for
$k\in\mathcal{Y}$ where $diag\\{a\\},a\in\mathbb{R}^{d}$ indicates $d\times d$
matrix whose diagonal elements are $a$. The prior distribution marginalized by
${\bf y}$ is

$\displaystyle p({\bf z})=\sum_{k=1}^{K}w_{k}\cdot\mathcal{N}({\bf
z}|m_{k},diag\\{s_{k}^{2}\\}).$

Here, $m_{k}\in\mathbb{R}^{d}$ and $s_{k}^{2}\in\mathbb{R}_{+}^{d}$ are pre-
determined parameters denoting the conceptual center and dispersion of the
latent variable for each label. The choice of $m_{k}$ and $s_{k}^{2}$ is the
customization of the latent space. Since $m_{k}$ and $w_{k}$ for all $k$ are
fixed, the notation of $m_{k}$ and $diag\\{s_{k}^{2}\\}$ is omitted in $p({\bf
z})$.

The encoder $p({\bf z}|{\bf x})$ is assumed as a mixture distribution, $p({\bf
z}|{\bf x})=\sum_{k=1}^{K}p({\bf y}=k|{\bf x})\cdot p({\bf z}|{\bf y}=k,{\bf
x})$. Since using $p({\bf z}|{\bf x})$ is computationally prohibitive, the
proposal distribution $w({\bf y}=k|{\bf x};\eta)$ for $p({\bf y}=k|{\bf x})$
and $g({\bf z}|{\bf y}=k,{\bf x};\xi)=\mathcal{N}({\bf z}|\mu_{k}({\bf
x};\xi),diag\\{\sigma^{2}_{k}({\bf x};\xi)\\})$ for $p({\bf z}|{\bf y}=k,{\bf
x})$ are introduced. The two proposal distributions are multinomial and
Gaussian, whose parameters are modeled by neural networks. The posterior
distribution is approximated by

$\displaystyle q({\bf z}|{\bf x};\eta,\xi)=\sum_{k=1}^{K}w({\bf y}=k|{\bf
x};\eta)\cdot g({\bf z}|{\bf y}=k,{\bf x};\xi),$ (4)

where the parameters of the neural networks are $\eta$ and $\xi$. Note that
the posterior mixture weight $w({\bf y}=k|{\bf x};\eta)$ decompose the latent
space according to labels. Given image ${\bf x}$ and label ${\bf y}$, the
posterior mixture weight assigns a latent variable to the mixture component
corresponding to the given label ${\bf y}$ and the latent space is separated
by labels.

The decoder is assumed by

$\displaystyle p({\bf x}|{\bf z};\theta,\beta)=\mathcal{N}({\bf x}|D({\bf
z};\theta),\beta\cdot{\bf I}),$

where the true label ${\bf y}$ is not included. Even though a natural decoder
can be considered as a mixture Gaussian $p({\bf x}|{\bf
z};\theta,\beta)=\sum_{k=1}^{K}p({\bf y}=k|{\bf z})\cdot\mathcal{N}({\bf
x}|D_{k}({\bf z};\theta),\beta\cdot{\bf I})$, the mixture distribution is
approximated by a uni-modal Gaussian distribution especially when $p({\bf
y}=k|{\bf z})$ is assumed to be close to 1 for some $k$. The latent space is
separated by $w({\bf y}|{\bf x};\eta)$ in our encoder and we can obtain such
$p({\bf y}|{\bf z})$. Thus, the validity of the assumption mainly depends on
the classification performance of $w({\bf y}|{\bf x};\eta)$. Although the
decoder variance $\beta$ is trainable, we fixed it due to computational
issues.

###  3.2 EXoN: Semi-Supervised VAE

Our proposed VAE model is derived from the joint pdf $p({\bf x},{\bf y})$.
$p({\bf x},{\bf y})$ is decomposed into $p({\bf x})$ and $p({\bf y}|{\bf x})$,
and the derivation of the ELBO only for $\log p({\bf x})$ leads to (5). In the
derivation, $p({\bf y}|{\bf x})$ is approximated by the proposal
classification model $w({\bf y}|{\bf x};\eta)$.

$\displaystyle\log p({\bf x},{\bf y};\theta,\beta)$ (5) $\displaystyle=$
$\displaystyle\log p({\bf x};\theta,\beta)+\log p({\bf y}|{\bf
x};\theta,\beta)$ $\displaystyle\geq$ $\displaystyle\mathbb{E}_{q}[\log p({\bf
x}|{\bf z};\theta,\beta)]-\mathcal{KL}\left(q({\bf z}|{\bf
x};\eta,\xi)\|p({\bf z})\right)$ $\displaystyle+\log p({\bf y}|{\bf
x};\theta,\beta)$ $\displaystyle\simeq$ $\displaystyle\mathbb{E}_{q}[\log
p({\bf x}|{\bf z};\theta,\beta)]-\mathcal{KL}\left(q({\bf z}|{\bf
x};\eta,\xi)\|p({\bf z})\right)$ $\displaystyle+\log w({\bf y}|{\bf x};\eta).$

The first two terms in (5) are the typical ELBO used in the unsupervised VAE
learning, and the remained term is the classification loss. Note that these
terms are coupled with the shared parameter $\eta$. The classification loss
term plays the role of training the posterior mixture weights. Therefore,
$w({\bf y}|{\bf x};\eta)$ with lower classification error can separate the
latent space more clearly by $q({\bf z}|{\bf x};\eta,\xi)$.

Whenever $p({\bf y}|{\bf x};\theta,\beta)=w({\bf y}|{\bf x};\eta)$,

$\displaystyle\log p({\bf x},{\bf y};\theta,\beta)\geq\mathcal{L}({\bf
x};\theta,\beta,\xi,\eta)+\log w({\bf y}|{\bf x};\eta),$

where

$\displaystyle\mathcal{L}({\bf x};\theta,\beta,\xi,\eta)$ $\displaystyle=$
$\displaystyle\mathbb{E}_{q}[\log p({\bf x}|{\bf z};\theta,\beta)]$
$\displaystyle-\mathcal{KL}^{u}\left(q({\bf z}|{\bf x};\eta,\xi)\|p({\bf
z})\right)$

and $\mathcal{KL}^{u}\left(q({\bf z}|{\bf x};\eta,\xi)\|p({\bf z})\right)$ is
the upper bound of $\mathcal{KL}\left(q({\bf z}|{\bf x};\eta,\xi)\|p({\bf
z})\right)$ (Wang et al. 2019; Guo et al. 2020) which is written in closed-
form as

$\displaystyle\mathcal{KL}(w({\bf y}|{\bf x};\eta)\|p({\bf y}))$
$\displaystyle+$ $\displaystyle\sum_{k=1}^{K}w({\bf y}=k|{\bf x};\eta)$
$\displaystyle\times\mathcal{KL}\left(q({\bf z}|{\bf x},{\bf
y}=k;\eta,\xi)\|p({\bf z}|{\bf y}=k)\right).$

Therefore, the lower bound on the joint likelihood for the entire dataset is

$\displaystyle\sum_{i\in I_{L}\cup
I_{U}}\mathcal{L}(x_{i};\theta,\beta,\xi,\eta)+\sum_{i\in I_{L}}\log
w(y_{i}|x_{i};\eta).$ (6)

(Kingma et al. 2014) first employed the classification loss (the second term
in (6)) as a penalty function of VAE and (Maaløe et al. 2016; Li et al. 2019;
Siddharth et al. 2017) use the same penalty function in the subsequent papers.
In our semi-supervised VAE model, the penalty function is applied with a
similar idea. However, the derivation of the regularized objective function
stems from (5) unlike the existing studies.

###  3.3 SCI: Soft-label Consistency Interpolation

The derivation of (6) is mathematically plausible, but it does not guarantee
state-of-the-art semi-supervised classification performance. To improve the
performance of the classification task, we propose a new loss function for the
pseudo-labeling semi-supervised classification method (Iscen et al. 2019;
Berthelot et al. 2019; Arazo et al. 2020) called SCI (Soft-label Consistency
Interpolation) loss.

The SCI loss consists of three parts, 1) an interpolated new image with a pair
of unlabeled images, 2) a pair of pseudo-label of the unlabeled images, and 3)
a convex combination of the cross entropy. Let $({\bf x}^{1},{\bf x}^{2})$ be
a pair of images and let $\tilde{{\bf x}}=\rho\cdot{\bf
x}^{1}+(1-\rho)\cdot{\bf x}^{2}$ be the interpolated image. Denote the
discrete probability $w({\bf y}|{\bf x};\hat{\eta})$ by $f({\bf
x};\hat{\eta})$. The pseudo-labels of ${\bf x}^{1}$ and ${\bf x}^{2}$ are
defined by $\tilde{{\bf y}}^{1}=f({\bf x}_{1};\hat{\eta})$ and $\tilde{{\bf
y}}^{2}=f({\bf x}_{2};\hat{\eta})$ for a given estimate $\hat{\eta}$ (true
label is used for labeled dataset). Note that the pseudo-label is not a
function of $\hat{\eta}$ but a non-trainable quantity determined by
$\hat{\eta}$. Then the SCI loss for $\eta$ with $({\bf x}^{1},{\bf x}^{2})$ is
defined by

$\displaystyle\mbox{SCI}(({\bf x}^{1},\tilde{{\bf y}}^{1}),({\bf
x}^{2},\tilde{{\bf y}}^{2});\eta)$ $\displaystyle=$ $\displaystyle\rho\cdot
H(\tilde{{\bf y}}^{1},f(\tilde{{\bf x}},\eta))+(1-\rho)\cdot H(\tilde{{\bf
y}}^{2},f(\tilde{{\bf x}},\eta)),$

where $H(p,q)$ is the cross entropy of the distribution $q$ relative to a
distribution $p$.

The SCI loss is motivated by the assumption of consistency interpolation.

######  Definition 1.

Given two data point ${\bf x}^{1},{\bf x}^{2}$, if

$\displaystyle f(\rho\cdot{\bf x}^{1}+(1-\rho)\cdot{\bf x}^{2};\eta^{*})$
$\displaystyle=$ $\displaystyle\rho\cdot f({\bf x}^{1};\eta^{*})+(1-\rho)\cdot
f({\bf x}^{2};\eta^{*})$

for $\forall\rho\in[0,1]$, then the consistency interpolation (Zhang et al.
2017; Verma et al. 2019) is satisfied for the parameter $\eta^{*}$.

The consistency interpolation assumes the existence of a linear map $f$ from
an image to a pseudo-label. It is well known that such a mix-up strategy can
improve the generalization error (Zhang et al. 2017).

Interestingly, the estimation of $f(\cdot;\eta^{*})$ in our VAE model is also
used in other existing semi-supervised learning methods (Feng et al. 2021;
Arazo et al. 2020) and the following algorithm provides a general framework to
estimate such a $f(\cdot;\eta^{*})$ in the VAE model. Let $\eta^{(t)}$ be an
estimate of $\eta$ obtained by the $t$th step while training the VAE, and
$\eta^{(t+1)}$ be the solution of following optimal interpolation problem
(Feng et al. 2021)

$\displaystyle\min_{\eta}$ $\displaystyle\rho\cdot\mathcal{KL}\left(f({\bf
x}^{1};\eta^{(t)})\|f(\tilde{{\bf x}};\eta)\right)$ (8) $\displaystyle+$
$\displaystyle(1-\rho)\cdot\mathcal{KL}\left(f({\bf
x}^{2};\eta^{(t)})\|f(\tilde{{\bf x}};\eta)\right).$

Then, it is easily shown that

$\displaystyle f(\rho\cdot{\bf x}^{1}+(1-\rho)\cdot{\bf x}^{2};\eta^{(t+1)})$
(9) $\displaystyle=$ $\displaystyle\rho\cdot f({\bf
x}^{1};\eta^{(t)})+(1-\rho)\cdot f({\bf x}^{2};\eta^{(t)}).$

If there exists $\eta^{*}$ such that $\eta^{(t)}\rightarrow\eta^{*}$ as
$t\rightarrow\infty$ given two data point ${\bf x}^{1},{\bf x}^{2}$, then (9)
finally implies the consistency interpolation. Since (8) is equivalent with
(3.3) up to constant, (3.3) is introduced in our proposed objective function.
We also found that using stochastically augmented images for SCI loss helps
our classifier achieve higher accuracy (see Appendix A.4 for detailed
augmentation techniques).

Finally, the objective function is given by

$\displaystyle-$ $\displaystyle\sum_{i\in I_{L}\cup
I_{U}}\mathcal{L}(x_{i};\theta,\beta,\xi,\eta)-\sum_{i\in I_{L}}\log
w(y_{i}|x_{i};\eta)$ (10) $\displaystyle-$
$\displaystyle\frac{\lambda_{1}}{\beta}\sum_{i\in I_{L}}\log
w(y_{i}|x_{i};\eta)$ $\displaystyle+$
$\displaystyle\frac{\lambda_{1}}{\beta}\sum_{i\in
I_{L}}\mbox{SCI}((x^{1}_{i},y^{1}_{i}),(x^{2}_{i},y^{2}_{i});\eta)$
$\displaystyle+$ $\displaystyle\frac{\lambda_{2}(t)}{\beta}\sum_{i\in
I_{U}}\mbox{SCI}((x^{1}_{i},\tilde{y}^{1}_{i}),(x^{2}_{i},\tilde{y}^{2}_{i});\eta),$

where $x_{i}^{1}=x_{i}$, $x_{i}^{2}$ is a randomly chosen sample, and
$\lambda_{2}(t)$ is a ramp-up function. If $\lambda_{1}$ and $\lambda_{2}(t)$
are set to be large enough, the mixture components of $q({\bf z}|{\bf
x};\eta,\xi)$ are shrunk toward those of $p({\bf z})$ according to the true
matching labels, and the latent space is well separated. The shared parameter
$\eta$ leads to this adaption of $q({\bf z}|{\bf x};\eta,\xi)$ to $p({\bf
z})$.

###  3.4 Activated Latent Subspace

In this section, we investigated the meaning of the explainability of the
latent space. Denote the random vector associated with the distribution of the
$k$th component of the mixture distribution (prior or posterior distribution)
by ${\bf z}^{k}=({\bf z}_{1}^{k},\cdots,{\bf z}_{d}^{k})$. Then, what is the
role of the posterior conditional variance $\mbox{Var}({\bf z}_{j}^{k}|{\bf
x};\xi)$ in the encoding process?

To answer above question, consider an extreme case where $\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)$ goes to zero for given $k$. If $\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)\rightarrow 0$, then ${\bf z}_{j}^{k}$ becomes
constant. ${\bf z}_{j}^{k}$ has a deterministic relationship with given ${\bf
x}$, implying that ${\bf z}_{j}^{k}$ can explain the given data point.
Therefore, the value of $\mbox{Var}({\bf z}_{j}^{k}|{\bf x};\xi)$ represents
the explainability of the latent space, and we call $j$ coordinate is
activated.

Interestingly, there is a theorem which shows a connection between our
objective function (10) and the posterior conditional variance
$\mbox{Var}({\bf z}_{j}^{k}|{\bf x};\xi)$.

######  Theorem 1.

Let $p({\bf z})$ the mixture prior distribution defined in Section 3.1 and let
$q({\bf z}|{\bf x};\eta,\xi)$ be (4).

$\displaystyle\frac{d}{2}-\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}}{2}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k}|{\bf
x};\xi)}{\mbox{Var}({\bf z}_{j}^{k})}\right]$ $\displaystyle\leq$
$\displaystyle\mathbb{E}_{{\bf x}}\Big{[}\sum_{k=1}^{K}w({\bf y}=k|{\bf
x};\eta)$ $\displaystyle\times\mathcal{KL}\Big{(}q({\bf z}|{\bf x},{\bf
y}=k;\eta,\xi)\|p({\bf z}|{\bf y}=k)\Big{)}\Big{]}$ $\displaystyle\leq$
$\displaystyle\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}}{2}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k})}{\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)}\right]-\frac{d}{2}.$

###### Proof.

See Appendix A.1. ∎

Theorem 1 means that the KL-divergence upper bound in our objective function
is bounded by the ratio of the prior and posterior variances. The lower bound
of Theorem 1 can be re-written as
$\sum_{k=1}^{K}\sum_{j=1}^{d}w_{k}/2\cdot\mbox{Var}_{{\bf x}|{\bf
y}=k}[\mathbb{E}({\bf z}_{j}^{k}|{\bf x};\xi)]/\mbox{Var}({\bf z}_{j}^{k})$,
so it can be interpreted as the coefficient of determination $R^{2}$ which
measures the proportion of latent variable variation explained by an
observation.

On the other hand, a small $\beta$ increases the weight of the generated
sample precision term and indirectly decreases the relative weight of the KL-
divergence term in (10). It implies that the KL-divergence term is relaxed
(not fully minimized). So, the relaxation of the KL-divergence term induces an
increase in the upper bound of Theorem 1. It means that there exists an
activated coordinate $j$ such that $\mbox{Var}({\bf z}_{j}^{k}|{\bf x};\xi)$
shrinks to zero for some $k$, as the prior mixture weight $w_{k}$ and variance
$\mbox{Var}({\bf z}_{j}^{k})$ are fixed due to pre-design. Thus, Theorem 1
shows that the relaxation of the KL-divergence obtained by tuning $\beta$ is
closely related to controlling the explainability of the latent space for
${\bf x}$.

Based on Theorem 1, we propose a statistics V-nat (VAE-natural unit of
information) which screens such activated latent coordinates,

$\displaystyle\log\mathbb{E}_{{\bf x}|{\bf y}=k}\left[\mbox{Var}({\bf
z}_{j}^{k})/\mbox{Var}({\bf z}_{j}^{k}|{\bf x};\xi)\right].$

Additionally, we define the set of activated latent coordinates as activated
latent subspace for $k$-labeled dataset for some $\delta>0$,

$\displaystyle{\cal A}_{k}(\delta)=\bigg{\\{}j\in[d]:\log\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k})}{\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)}\right]>\delta\bigg{\\}}$

where $[d]=\\{1,\cdots,d\\}$ Our numerical study found that the subspace
represents the informative characteristics of generated samples, and the
subspace can be effectively used to produce a high-quality image (see Section
4.2). These results are consistent with those of VQ-VAE (Van Den Oord, Vinyals
et al. 2017) that an image generation process mainly depends on a nearly
deterministic encoding map, and the refinement of the map can improve the
quality of generated images.

##  4 Experiments

###  4.1 MNIST Dataset

We have used the MNIST dataset (LeCun, Cortes, and Burges 2010) to consider
the 2-dimensional latent space in our VAE model. The values were scaled in the
range of -1 to 1. The encoder returns the 10-component Gaussian mixture
distribution parameters, mixing probability, mean vector, and diagonal
covariance elements. Thus, the encoder maps ${\cal X}$ to
$\left((0,1)\times\mathbb{R}^{2}\times\mathbb{R}_{+}^{2}\right)^{10}$.
Especially, the mixing probabilities are produced by the classifier in the
encoder. Gumbel-Max trick (Gumbel 1954; Jang, Gu, and Poole 2016) is used for
sampling discrete latent variables. Detailed network architectures of the
encoder, decoder and classifier are described in Appendix A.2 Table 4, 5.

Figure 1: Scatter plot of samples from $p({\bf z})$.

The customized $p({\bf z})$ is illustrated in Figure 1. Each label is assigned
counterclockwise from 3 O’clock on the component of the mixture. Note that
Figure 1 illustrates an example of conceptual centers. The $k$th component of
the mixture distribution corresponds to the distribution of the digit $(k-1)$
in the MNIST dataset, for $k=1,\cdots,10$ (see Appendix A.4 for detailed pre-
design settings).

Our model is compared with (Kingma et al. 2014; Maaløe et al. 2016; Li et al.
2019; Hajimiri, Lotfi, and Baghshah 2021; Feng et al. 2021), and the
simulation results show that our fitted model achieves a competitive
classification performance, 3.12% error with 59,900 unlabeled and 100 labeled
images. (see Appendix A.2 Table 6 for comparison results). For implementation
details, see Appendix A.4.

Figure 2: Top row: scatter plot of $z\sim q({\bf z}|x;\eta,\xi)$ given the
test dataset. Bottom row: generated images from grid points on the latent
space.

#### Effect of Regularizations

First, the role of tuning parameter $\beta$ of (10) is investigated in fitted
latent space (Lucas et al. 2019). The top panels of Figure 2 display the
samples from $q({\bf z}|x_{i};\eta,\xi)$ where $x_{i}$ is an observation in
the test dataset, and the bottom panels show images generated from grid points
on the latent space.

The top panels demonstrate that the large $\beta$ regularizes $q({\bf
z}|x;\eta,\xi)$ to $p({\bf z})$ by indirectly increasing the weight of the KL-
divergence term in the objective function. The bottom panels illustrate that
each generated image exactly matches the label defined on the pre-designed
latent space. Also, it is confirmed that generated images are naturally
interpolated on the pre-designed latent space according to our arrangement of
conceptual centers. These results show that the proposed VAE yields an
explainable latent space with the labels. See Appendix A.2 for additional
evaluation results (the negative average single-scale structural similarity
(negative SSIM) (Wang et al. 2004; Zheng and Sun 2019), the classification
error, and the KL-divergence) with various $\beta$ values.

#### Diversity of Generated Images

The bottom panels of Figure 2 show that the images are generated with lower
diversity for larger $\beta$. This result can be explained by the mutual
information $I(\cdot,\cdot)$. The expectation of
$\mathcal{K}\mathcal{L}^{u}(q({\bf z}|{\bf x};\eta,\xi)\|p({\bf z}))$ in (10)
is rewritten in terms of mutual information as

$\displaystyle\mathbb{E}_{p({\bf x})}\left[\mathcal{K}\mathcal{L}^{u}(q({\bf
z}|{\bf x};\eta,\xi)\|p({\bf z}))\right]$ (11) $\displaystyle=$ $\displaystyle
I({\bf x},{\bf y};\eta)+\mathcal{KL}(w({\bf y};\eta)\|p({\bf y}))$
$\displaystyle+\mathbb{E}_{p({\bf y})}[I({\bf x},{\bf z}|{\bf y};\eta,\xi)]$
$\displaystyle+\mathbb{E}_{p({\bf y})}[\mathcal{KL}(q({\bf z}|{\bf
y};\eta,\xi)\|p({\bf z}|{\bf y}))].$

See Appendix A.1 for detailed derivation.

To maximize (10) under a large $\beta$, (11) should be nearly zero and the
mutual information between ${\bf x}|{\bf y}$ and ${\bf z}|{\bf y}$ be close to
zero for all classes ${\bf y}$ as well. The conditional independence of ${\bf
x}$ and ${\bf z}$ implies that $D({\bf z};\theta)$ does not depend on ${\bf
z}$ for each ${\bf y}$, because ${\bf x}|{\bf z}\sim N(D({\bf
z};\theta),\beta\cdot{\bf I})$ is assumed. For this reason, as shown in the
bottom row of Figure 2, the latent mixture component for a specific label is
not able to capture the complex pattern of observations that belong to a
corresponding label when $\beta$ is large.

#### Customized Latent Space

(a) Parted-VAE (Hajimiri, Lotfi, and Baghshah 2021)

(b) The EXoN ($\beta=5$)

Figure 3: The interpolation path from point A (label 0) to B (label 1) and
reconstructed images.

Figure 3 shows two series of reconstructed images from interpolation on the
prior structure. The latent variable A and B are sampled from mixture
components of labels 0 and 1, respectively. Suppose the interpolation path
from point A to B produces interpolated images with labels 0 and 1. However,
the left panel of Figure 3 shows reconstructed images whose labels are neither
0 nor 1 in the middle of the interpolation path.

It means that the interpolation path is unpredictable before training with
Parted-VAE (Hajimiri, Lotfi, and Baghshah 2021). However, the interpolation
path with our model consists of only interpolated images with labels 0 and 1
since our latent space can be manually designed (see the right panel of Figure
3). Furthermore, we can pre-determine the interpolation path and the
resolution of interpolation by controlling the proximity between mixture
components (see Appendix A.2).

###  4.2 CIFAR-10 Dataset

We apply our VAE model to the CIFAR-10 dataset (Krizhevsky, Hinton et al.
2009) of which images have 10 labels. We evaluated our model with 400 labeled
images per class and 46,000 unlabeled images. All values of the images are
scaled to the range of $(-1,1)$. 256-dimensional latent space and the
10-mixture Gaussian distribution are employed. Each component of the prior
mixture distribution has a separate mean vector, and all components share the
same covariance (see Appendix A.4 for detailed pre-designed priors, hyper-
parameters, and implementation settings). Gumbel-Max trick (Gumbel 1954; Jang,
Gu, and Poole 2016) is used for sampling discrete latent variables. The
network architectures of the encoder, decoder and classifier are shown in
Appendix A.3 Table 8, 9.

#### Comparison

models | error (%) | Inception Score
---|---|---
$\Pi$-model (4.5M) (Laine and Aila 2016) | 17.53 $\pm$ 0.15 | -
VAT (1.5M) (Miyato et al. 2018) | 10.55 $\pm$ 0.05 | -
MixMatch (5.8M) (Berthelot et al. 2019) |  5.91 $\pm$ 0.36 | -
PLCB (4.5M) (Arazo et al. 2020) | 7.82 $\pm$ 0.18 | -
M2 (4.5M) (Kingma et al. 2014) | 27.30 $\pm$ 0.25 | 1.86 $\pm$ 0.02
Parted-VAE (5.8M) (Hajimiri, Lotfi, and Baghshah 2021) | 34.00 $\pm$ 1.67 | 1.73 $\pm$ 0.20
SHOT-VAE (5.8M) (Feng et al. 2021) | 6.45 $\pm$ 0.28 |  3.43 $\pm$ 0.05
Ours (4.5M) ($\beta=0.01$) | 7.45 $\pm$ 0.48 | 3.14 $\pm$ 0.05
Table 1: Mean and standard deviation of test dataset classification error for
4,000 labeled dataset and the Inception Score from 5 replicates of the
experiment. The numbers in parenthesis are the number of classifier
parameters. The Inception Scores are computed using 10,000 generated images. *
means that the result is from the original paper.

(a) SHOT-VAE (between same classes)

(b) EXoN ($\beta=0.05$) (between same classes)

(c) Parted-VAE (between different classes)

(d) EXoN ($\beta=0.05$) (between different classes)

Figure 4:  Series of images obtained by interpolating two points on the latent
space(or subspace), which various models produce. The most left and right are
original images.

Table 1 shows quantitative comparison results of the classification
performance in the CIFAR-10 benchmark setting and image generation
performance. The classification models (Laine and Aila 2016; Miyato et al.
2018; Berthelot et al. 2019; Arazo et al. 2020) are not generative, so the
Inception Score is not computed. Although the EXoN is not the best in the
Inception Score, our model has a remarkable image interpolation quality. It is
shown in Figure 4 which visualizes the reconstructed images from
interpolation. This Figure indicates that our model achieves much better
recovering and interpolation performance than other models (Hajimiri, Lotfi,
and Baghshah 2021; Feng et al. 2021) (see Appendix Figure 11 for more
reconstructed images by our model).

#### Activated Latent Subspace

(a) EXoN ($\beta=0.05$)

(b) EXoN ($\beta=0.05$)

Figure 5:  (a) Plot of V-nat values for automobile class. (b) Top row:
visualization of noises added to the latent variable, bottom row: generated
image from the corresponding perturbed latent variable.

The properties of the activated latent subspace are investigated. For
$\beta=0.05$, the latent space with $|{\cal A}_{k}(1)|=115$ is obtained for
the automobile class. The associated V-nat values are displayed in Figure 5(a)
(see Appendix A.3 for other classes). We generate images on three cases of
latent variables, 1) the original latent variable, 2) the latent variable
perturbed by the uniform noise $U(-2,2)$ on ${\cal A}_{k}(1)$, 3) the latent
variable perturbed by the uniform noise $U(-2,2)$ on ${\cal A}^{c}_{k}(1)$.
Figure 5(b) displays examples of noise values and generated images.

The two images of left and right in Figure 5(b) seem identical even though
their latent variables are far from each other. While, the middle image of the
bottom panel is distorted compared with the other images, which confirms that
only the activated latent subspace determines the features of generated
images. Motivated by these numerical results in Figure 5(b), we manipulated
images using only the activated latent subspace (see results in Appendix A.3).

$\beta$ | 0.01 | 0.05 | 0.1 | 0.5 | 1
---|---|---|---|---|---
Inception Score | 3.14 $\pm$ 0.05 | 3.36 $\pm$ 0.05 | 3.37 $\pm$ 0.10 | 3.17 $\pm$ 0.08 | 2.86 $\pm$ 0.06
$|{\cal A}_{k}(1)|$ | 189.2 $\pm$ 22.7 | 114.6 $\pm$ 4.22 | 82.8 $\pm$ 12.50 | 36.00 $\pm$ 5.87 | 22.00 $\pm$ 1.10
error (%) | 7.45 $\pm$ 0.48 | 7.27 $\pm$ 0.19 | 7.51 $\pm$ 0.25 | 8.18 $\pm$ 0.63 | 9.50 $\pm$ 1.07
Table 2: For various $\beta$ values, mean and standard deviation of the
Inception Score (Salimans et al. 2016), cardinality of activated latent
subspace $|{\cal A}_{k}(1)|$ of automobile class and test dataset
classification error with 4,000 labeled dataset from 5 replicates of the
experiment.

Theorem 1 theoretically supports that relaxation of KL-divergence (shrinkage
of $\beta$) leads to a large set ${\cal A}_{k}(1)$. Meanwhile, (11) implies
that strong regularization of KL-divergence (increasing $\beta$) leads to
decreasing $\mathbb{E}_{p({\bf y})}[I({\bf x},{\bf z}|{\bf y};\eta,\xi)]$,
which results in a small set ${\cal A}_{k}(1)$. These theoretical results are
shown numerically in Table 2.

##  5 Conclusion and Limitations

This paper proposes a new method to construct an explainable latent space in
semi-supervised learning VAE. The explainable latent space can be obtained by
combining classification losses (cross-entropy and SCI loss) and the relaxed
KL-divergence term in our objective function. By the classification, the
latent space can be effectively decomposed according to the mixture components
and the true labels of observations. In addition, it is shown that the
relaxation of the KL-divergence increases the dimension of the activated
latent subspace which determines the characteristics of the images generated
from our proposed model. The activated latent subspace can be discovered
through the V-nat measure.

In short, the explainability of our proposed model is demonstrated by two
points 1) the label-specific conceptual means and variances of the latent
distribution, and 2) the features of the given images can be analyzed in
association with the activated latent subspace. Thus, the practical advantage
of our method is that the manually interpolated images corresponding to the
user’s specific purpose can be achieved by customizing the prior distribution
structure in advance. We guess that the clustering method would replace the
role of the classification loss in our VAE model, and leave the development of
the EXoN in unsupervised learning as our future work.

## References

  * Abadi et al. (2015) Abadi, M.; Agarwal, A.; Barham, P.; Brevdo, E.; Chen, Z.; Citro, C.; Corrado, G. S.; Davis, A.; Dean, J.; Devin, M.; Ghemawat, S.; Goodfellow, I.; Harp, A.; Irving, G.; Isard, M.; Jia, Y.; Jozefowicz, R.; Kaiser, L.; Kudlur, M.; Levenberg, J.; Mané, D.; Monga, R.; Moore, S.; Murray, D.; Olah, C.; Schuster, M.; Shlens, J.; Steiner, B.; Sutskever, I.; Talwar, K.; Tucker, P.; Vanhoucke, V.; Vasudevan, V.; Viégas, F.; Vinyals, O.; Warden, P.; Wattenberg, M.; Wicke, M.; Yu, Y.; and Zheng, X. 2015.  TensorFlow: Large-Scale Machine Learning on Heterogeneous Systems.  Software available from tensorflow.org.
  * Arazo et al. (2020) Arazo, E.; Ortego, D.; Albert, P.; O’Connor, N. E.; and McGuinness, K. 2020.  Pseudo-labeling and confirmation bias in deep semi-supervised learning.  In _2020 International Joint Conference on Neural Networks (IJCNN)_ , 1–8. IEEE.
  * Berthelot et al. (2019) Berthelot, D.; Carlini, N.; Goodfellow, I.; Papernot, N.; Oliver, A.; and Raffel, C. A. 2019.  Mixmatch: A holistic approach to semi-supervised learning.  _Advances in Neural Information Processing Systems_ , 32.
  * Dilokthanakul et al. (2016) Dilokthanakul, N.; Mediano, P. A.; Garnelo, M.; Lee, M. C.; Salimbeni, H.; Arulkumaran, K.; and Shanahan, M. 2016.  Deep unsupervised clustering with gaussian mixture variational autoencoders.  _arXiv preprint arXiv:1611.02648_.
  * Feng et al. (2021) Feng, H.-Z.; Kong, K.; Chen, M.; Zhang, T.; Zhu, M.; and Chen, W. 2021.  SHOT-VAE: semi-supervised deep generative models with label-aware ELBO approximations.  In _Proceedings of the AAAI Conference on Artificial Intelligence_ , volume 35, 7413–7421.
  * Gumbel (1954) Gumbel, E. J. 1954.  _Statistical theory of extreme values and some practical applications: a series of lectures_ , volume 33.  US Government Printing Office.
  * Guo et al. (2020) Guo, C.; Zhou, J.; Chen, H.; Ying, N.; Zhang, J.; and Zhou, D. 2020.  Variational Autoencoder With Optimizing Gaussian Mixture Model Priors.  _IEEE Access_ , 8: 43992–44005.
  * Hajimiri, Lotfi, and Baghshah (2021) Hajimiri, S.; Lotfi, A.; and Baghshah, M. S. 2021.  Semi-Supervised Disentanglement of Class-Related and Class-Independent Factors in VAE.  _arXiv preprint arXiv:2102.00892_.
  * Higgins et al. (2016) Higgins, I.; Matthey, L.; Pal, A.; Burgess, C.; Glorot, X.; Botvinick, M.; Mohamed, S.; and Lerchner, A. 2016.  beta-vae: Learning basic visual concepts with a constrained variational framework.
  * Ioffe and Szegedy (2015) Ioffe, S.; and Szegedy, C. 2015.  Batch normalization: Accelerating deep network training by reducing internal covariate shift.  In _International conference on machine learning_ , 448–456. PMLR.
  * Iscen et al. (2019) Iscen, A.; Tolias, G.; Avrithis, Y.; and Chum, O. 2019.  Label propagation for deep semi-supervised learning.  In _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , 5070–5079.
  * Jang, Gu, and Poole (2016) Jang, E.; Gu, S.; and Poole, B. 2016.  Categorical reparameterization with gumbel-softmax.  _arXiv preprint arXiv:1611.01144_.
  * Jordan et al. (1999) Jordan, M. I.; Ghahramani, Z.; Jaakkola, T. S.; and Saul, L. K. 1999.  An introduction to variational methods for graphical models.  _Machine learning_ , 37(2): 183–233.
  * Kingma and Ba (2014) Kingma, D. P.; and Ba, J. 2014.  Adam: A method for stochastic optimization.  _arXiv preprint arXiv:1412.6980_.
  * Kingma et al. (2014) Kingma, D. P.; Mohamed, S.; Rezende, D. J.; and Welling, M. 2014.  Semi-supervised learning with deep generative models.  In _Advances in neural information processing systems_ , 3581–3589.
  * Kingma and Welling (2013) Kingma, D. P.; and Welling, M. 2013.  Auto-encoding variational bayes.  _arXiv preprint arXiv:1312.6114_.
  * Krizhevsky, Hinton et al. (2009) Krizhevsky, A.; Hinton, G.; et al. 2009.  Learning multiple layers of features from tiny images.
  * Laine and Aila (2016) Laine, S.; and Aila, T. 2016.  Temporal ensembling for semi-supervised learning.  _arXiv preprint arXiv:1610.02242_.
  * LeCun, Cortes, and Burges (2010) LeCun, Y.; Cortes, C.; and Burges, C. 2010.  MNIST handwritten digit database.  _ATT Labs [Online]. Available: http://yann.lecun.com/exdb/mnist_ , 2.
  * Li et al. (2019) Li, Y.; Pan, Q.; Wang, S.; Peng, H.; Yang, T.; and Cambria, E. 2019.  Disentangled variational auto-encoder for semi-supervised learning.  _Information Sciences_ , 482: 73–85.
  * Loshchilov and Hutter (2017) Loshchilov, I.; and Hutter, F. 2017.  Decoupled weight decay regularization.  _arXiv preprint arXiv:1711.05101_.
  * Lucas et al. (2019) Lucas, J.; Tucker, G.; Grosse, R. B.; and Norouzi, M. 2019.  Don’t Blame the ELBO! A Linear VAE Perspective on Posterior Collapse.  In _Advances in Neural Information Processing Systems_ , 9408–9418.
  * Maaløe, Fraccaro, and Winther (2017) Maaløe, L.; Fraccaro, M.; and Winther, O. 2017.  Semi-supervised generation with cluster-aware generative models.  _arXiv preprint arXiv:1704.00637_.
  * Maaløe et al. (2016) Maaløe, L.; Sønderby, C. K.; Sønderby, S. K.; and Winther, O. 2016.  Auxiliary deep generative models.  _arXiv preprint arXiv:1602.05473_.
  * Mathieu et al. (2019) Mathieu, E.; Rainforth, T.; Siddharth, N.; and Teh, Y. W. 2019.  Disentangling disentanglement in variational autoencoders.  In _International Conference on Machine Learning_ , 4402–4412.
  * Miyato et al. (2018) Miyato, T.; Maeda, S.-i.; Koyama, M.; and Ishii, S. 2018.  Virtual adversarial training: a regularization method for supervised and semi-supervised learning.  _IEEE transactions on pattern analysis and machine intelligence_ , 41(8): 1979–1993.
  * Rezende, Mohamed, and Wierstra (2014) Rezende, D. J.; Mohamed, S.; and Wierstra, D. 2014.  Stochastic backpropagation and approximate inference in deep generative models.  _arXiv preprint arXiv:1401.4082_.
  * Salimans et al. (2016) Salimans, T.; Goodfellow, I.; Zaremba, W.; Cheung, V.; Radford, A.; and Chen, X. 2016.  Improved techniques for training gans.  _Advances in neural information processing systems_ , 29: 2234–2242.
  * Siddharth et al. (2017) Siddharth, N.; Paige, B.; Van de Meent, J.-W.; Desmaison, A.; Goodman, N.; Kohli, P.; Wood, F.; and Torr, P. 2017.  Learning disentangled representations with semi-supervised deep generative models.  In _Advances in neural information processing systems_ , 5925–5935.
  * Van Den Oord, Vinyals et al. (2017) Van Den Oord, A.; Vinyals, O.; et al. 2017.  Neural discrete representation learning.  _Advances in neural information processing systems_ , 30.
  * Verma et al. (2019) Verma, V.; Kawaguchi, K.; Lamb, A.; Kannala, J.; Bengio, Y.; and Lopez-Paz, D. 2019\.  Interpolation consistency training for semi-supervised learning.  _arXiv preprint arXiv:1903.03825_.
  * Wang et al. (2019) Wang, W.; Gan, Z.; Xu, H.; Zhang, R.; Wang, G.; Shen, D.; Chen, C.; and Carin, L. 2019.  Topic-guided variational autoencoders for text generation.  _arXiv preprint arXiv:1903.07137_.
  * Wang et al. (2004) Wang, Z.; Bovik, A. C.; Sheikh, H. R.; and Simoncelli, E. P. 2004.  Image quality assessment: from error visibility to structural similarity.  _IEEE transactions on image processing_ , 13(4): 600–612.
  * Willetts, Roberts, and Holmes (2019) Willetts, M.; Roberts, S.; and Holmes, C. 2019.  Disentangling to Cluster: Gaussian Mixture Variational Ladder Autoencoders.  _arXiv preprint arXiv:1909.11501_.
  * Zhang et al. (2017) Zhang, H.; Cisse, M.; Dauphin, Y. N.; and Lopez-Paz, D. 2017.  mixup: Beyond empirical risk minimization.  _arXiv preprint arXiv:1710.09412_.
  * Zheng and Sun (2019) Zheng, Z.; and Sun, L. 2019.  Disentangling latent space for vae by label relevant/irrelevant dimensions.  In _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition_ , 12192–12201.

##  Appendix A Appendix

###  A.1 Theoretical Derivations

#### Upper bound of the KL-divergence

We use the upper bound of the KL-divergence from $p({\bf z})$ to $q({\bf
z}|{\bf x};\eta,\xi)$ (Wang et al. 2019; Guo et al. 2020) by

$\displaystyle\mathcal{KL}^{u}\left(q({\bf z}|{\bf x};\eta,\xi)\|p({\bf
z})\right)$ (12) $\displaystyle\equiv$ $\displaystyle\mathcal{KL}(w({\bf
y}|{\bf x};\eta)\|p({\bf y}))$ $\displaystyle+$
$\displaystyle\sum_{k=1}^{K}w({\bf y}=k|{\bf x};\eta)$
$\displaystyle\times\mathcal{KL}\Big{(}\mathcal{N}({\bf z}|\mu_{k}({\bf
x};\xi),diag\\{(\sigma_{k}({\bf x};\xi))^{2}\\})\|$
$\displaystyle\mathcal{N}({\bf
z}|\mu_{k}^{0},diag\\{(\sigma_{k}^{0})^{2}\\})\Big{)}.$

#### Proof of Theorem 1

Let ${\bf z}=({\bf z}_{1},\cdots,{\bf z}_{d})\sim
N(\mu,diag\\{(\sigma)^{2}\\})$ and ${\bf z}|{\bf x}\sim N(\mu({\bf
x}),diag\\{(\sigma({\bf x}))^{2}\\})$. For understanding following equations,
we denote the $j$th elements of $\mu$, $\mu({\bf x})$, $(\sigma)^{2}$ and
$(\sigma({\bf x}))^{2}$ by $\mathbb{E}({\bf z}_{j})$, $\mathbb{E}({\bf
z}_{j}|{\bf x})$, $\mbox{Var}({\bf z}_{j})$ and $\mbox{Var}({\bf z}_{j}|{\bf
x})$, respectively. The KL-divergence between two univariate normal
distributions leads to a close form of the KL-divergence from
$N(\mu,diag\\{(\sigma)^{2}\\})$ to $N(\mu({\bf x}),diag\\{(\sigma({\bf
x}))^{2}\\})$,

$\displaystyle\mathbb{E}_{{\bf x}}\Big{[}\mathcal{KL}\Big{(}N(\mu({\bf
x}),diag\\{(\sigma({\bf x}))^{2}\\})\|$ $\displaystyle
N(\mu,diag\\{(\sigma)^{2}\\})\Big{)}\Big{]}$ $\displaystyle=$
$\displaystyle\mathbb{E}_{{\bf
x}}\Big{[}\frac{1}{2}\sum_{j=1}^{d}\Big{(}\frac{[\mathbb{E}({\bf z}_{j}|{\bf
x})-\mathbb{E}({\bf z}_{j})]^{2}}{\mbox{Var}({\bf
z}_{j})}+\frac{\mbox{Var}({\bf z}_{j}|{\bf x})}{\mbox{Var}({\bf z}_{j})}$
$\displaystyle+\log\frac{\mbox{Var}({\bf z}_{j})}{\mbox{Var}({\bf z}_{j}|{\bf
x})}-1\Big{)}\Big{]}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{j=1}^{d}\Big{(}\frac{\mbox{Var}_{{\bf
x}}[\mathbb{E}({\bf z}_{j}|{\bf x})]}{\mbox{Var}({\bf
z}_{j})}+\frac{\mathbb{E}_{{\bf x}}[\mbox{Var}({\bf z}_{j}|{\bf
x})]}{\mbox{Var}({\bf z}_{j})}$ $\displaystyle+\mathbb{E}_{{\bf
x}}\left[\log\frac{\mbox{Var}({\bf z}_{j})}{\mbox{Var}({\bf z}_{j}|{\bf
x})}\right]\Big{)}-\frac{d}{2}$ (13) $\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{j=1}^{d}\left(\frac{\mbox{Var}({\bf
z}_{j})}{\mbox{Var}({\bf z}_{j})}+\mathbb{E}_{{\bf
x}}\left[\log\frac{\mbox{Var}({\bf z}_{j})}{\mbox{Var}({\bf z}_{j}|{\bf
x})}\right]\right)-\frac{d}{2}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{j=1}^{d}\mathbb{E}_{{\bf
x}}\left[\log\frac{\mbox{Var}({\bf z}_{j})}{\mbox{Var}({\bf z}_{j}|{\bf
x})}\right].$

Using above equality (A.1), the expectation of the second term in (12) can be
written as:

$\displaystyle\sum_{k=1}^{K}\mathbb{E}_{{\bf x}}\Big{[}w({\bf y}=k|{\bf
x};\eta)$ (14) $\displaystyle\times\mathcal{KL}\Big{(}\mathcal{N}({\bf
z}|\mu_{k}({\bf x};\xi),diag\\{(\sigma_{k}({\bf x};\xi))^{2}\\})\|$
$\displaystyle\mathcal{N}({\bf
z}|\mu_{k}^{0},diag\\{(\sigma_{k}^{0})^{2}\\})\Big{)}\Big{]}$ $\displaystyle=$
$\displaystyle\sum_{k=1}^{K}\int_{{\bf x}}p({\bf x})w({\bf y}=k|{\bf x};\eta)$
$\displaystyle\times\mathcal{KL}\Big{(}\mathcal{N}({\bf z}|\mu_{k}({\bf
x};\xi),diag\\{(\sigma_{k}({\bf x};\xi))^{2}\\})\|$
$\displaystyle\mathcal{N}({\bf
z}|\mu_{k}^{0},diag\\{(\sigma_{k}^{0})^{2}\\})\Big{)}d{\bf x}$
$\displaystyle=$ $\displaystyle\sum_{k=1}^{K}\int_{{\bf x}}p({\bf y}=k)p({\bf
x}|{\bf y}=k;\eta)$ $\displaystyle\times\mathcal{KL}\Big{(}\mathcal{N}({\bf
z}|\mu_{k}({\bf x};\xi),diag\\{(\sigma_{k}({\bf x};\xi))^{2}\\})\|$
$\displaystyle\mathcal{N}({\bf
z}|\mu_{k}^{0},diag\\{(\sigma_{k}^{0})^{2}\\})\Big{)}d{\bf x}$
$\displaystyle=$ $\displaystyle\sum_{k=1}^{K}w_{k}^{0}\cdot\mathbb{E}_{{\bf
x}|{\bf y}=k}\Big{[}$ $\displaystyle\mathcal{KL}\Big{(}\mathcal{N}({\bf
z}|\mu_{k}({\bf x};\xi),diag\\{(\sigma_{k}({\bf x};\xi))^{2}\\})\|$
$\displaystyle\mathcal{N}({\bf
z}|\mu_{k}^{0},diag\\{(\sigma_{k}^{0})^{2}\\})\Big{)}\Big{]}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{k=1}^{K}w_{k}^{0}\sum_{j=1}^{d}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\log\frac{(\sigma^{0}_{kj})^{2}}{(\sigma_{kj}({\bf
x};\xi))^{2}}\right],$

where $\sigma^{0}_{kj}$ and $\sigma_{kj}({\bf x};\xi)$ are the $j$th element
of $\sigma_{k}^{0}$ and $\sigma_{k}({\bf x};\xi)$, respectively. The second
equality above equations holds for $p({\bf y}|{\bf x};\eta)=w({\bf y}|{\bf
x};\eta)$.

Let ${\bf z}^{k}=({\bf z}_{1}^{k},\cdots,{\bf z}_{d}^{k})\sim
N(\mu_{k}^{0},diag\\{(\sigma_{k}^{0})^{2}\\})$, the distribution of the $k$th
component of the mixture distribution $p({\bf z})$. First, the upper bound of
(12) is written as

$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}^{u}\left(q({\bf z}|{\bf
x};\eta,\xi)\|p({\bf z})\right)$ (15) $\displaystyle=$
$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}(w({\bf y}|{\bf x};\eta)\|p({\bf
y}))$
$\displaystyle+\frac{1}{2}\sum_{k=1}^{K}w_{k}^{0}\sum_{j=1}^{d}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\log\frac{\mbox{Var}({\bf z}_{j}^{k})}{\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)}\right]$ $\displaystyle\leq$
$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}(w({\bf y}|{\bf x};\eta)\|p({\bf
y}))$
$\displaystyle+\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}^{0}}{2}\log\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k})}{\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)}\right]$ $\displaystyle\leq$
$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}(w({\bf y}|{\bf x};\eta)\|p({\bf
y}))$
$\displaystyle+\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}^{0}}{2}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k})}{\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)}\right]-\frac{d}{2},$

by Jensen’s inequality and $\log x\leq x-1$.

Similarly, the lower bound of (12) is derived as

$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}^{u}\left(q({\bf z}|{\bf
x};\eta,\xi)\|p({\bf z})\right)$ (16) $\displaystyle=$
$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}(w({\bf y}|{\bf x};\eta)\|p({\bf
y}))$
$\displaystyle-\frac{1}{2}\sum_{k=1}^{K}w_{k}^{0}\sum_{j=1}^{d}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\log\frac{\mbox{Var}({\bf z}_{j}^{k}|{\bf
x};\xi)}{\mbox{Var}({\bf z}_{j}^{k})}\right]$ $\displaystyle\geq$
$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}(w({\bf y}|{\bf x};\eta)\|p({\bf
y}))$
$\displaystyle-\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}^{0}}{2}\log\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k}|{\bf
x};\xi)}{\mbox{Var}({\bf z}_{j}^{k})}\right]$ $\displaystyle\geq$
$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}(w({\bf y}|{\bf x};\eta)\|p({\bf
y}))$
$\displaystyle-\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}^{0}}{2}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k}|{\bf
x};\xi)}{\mbox{Var}({\bf z}_{j}^{k})}\right]+\frac{d}{2}.$

Therefore, combining (15) and (16),

$\displaystyle\frac{d}{2}-\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}^{0}}{2}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k}|{\bf
x};\xi)}{\mbox{Var}({\bf z}_{j}^{k})}\right]$ $\displaystyle\leq$
$\displaystyle\mathbb{E}_{{\bf x}}\mathcal{KL}^{u}\left(q({\bf z}|{\bf
x};\eta,\xi)\|p({\bf z})\right)-\mathbb{E}_{{\bf x}}\mathcal{KL}(w({\bf
y}|{\bf x};\eta)\|p({\bf y}))$ $\displaystyle\leq$
$\displaystyle\sum_{k=1}^{K}\sum_{j=1}^{d}\frac{w_{k}^{0}}{2}\mathbb{E}_{{\bf
x}|{\bf y}=k}\left[\frac{\mbox{Var}({\bf z}_{j}^{k})}{\mbox{Var}({\bf
z}_{j}^{k}|{\bf x};\xi)}\right]-\frac{d}{2},$ (18)

completes the proof.

#### Mutual Information Derivation

Here, we denote the mutual information as $I(\cdot,\cdot)$ and the expectation
of the KL-divergence upper bound (12) is

$\displaystyle\mathbb{E}_{p({\bf x})}\mathcal{KL}^{u}\left(q({\bf z}|{\bf
x};\eta,\xi)\|p({\bf z})\right)$ (19) $\displaystyle=$
$\displaystyle\mathbb{E}_{p({\bf x})}\left[\mathcal{KL}(w({\bf y}|{\bf
x};\eta)\|p({\bf y}))\right]$ $\displaystyle+\mathbb{E}_{p({\bf
x})}\Big{[}\sum_{k=1}^{K}w({\bf y}=k|{\bf x};\eta)$
$\displaystyle\times\mathcal{KL}\left(q({\bf z}|{\bf x},{\bf y}=k;\xi)\|p({\bf
z}|{\bf y}=k)\right)\Big{]}.$

We assume that $w({\bf y}|{\bf x};\eta)=p({\bf y}|{\bf x};\eta)$.

$\displaystyle\mathbb{E}_{p({\bf x})}[\mathcal{KL}\left(w({\bf y}|{\bf
x};\eta)\|p({\bf y})\right)]$ (20) $\displaystyle=$
$\displaystyle\mathbb{E}_{p({\bf x})}\mathbb{E}_{w({\bf y}|{\bf
x};\eta)}\left[\log\frac{w({\bf y}|{\bf x};\eta)}{w({\bf y})}\right]$
$\displaystyle+\mathbb{E}_{p({\bf x})}\mathbb{E}_{w({\bf y}|{\bf
x};\eta)}[\log w({\bf y};\eta)]$ $\displaystyle-\mathbb{E}_{p({\bf
x})}\mathbb{E}_{w({\bf y}|{\bf x};\eta)}[\log p({\bf y})]$ $\displaystyle=$
$\displaystyle\mathbb{E}_{p({\bf x})}[\mathcal{KL}(w({\bf y}|{\bf
x};\eta)\|w({\bf y};\eta))]$ $\displaystyle+\mathbb{E}_{w({\bf
y};\eta)}\left[\log\frac{w({\bf y};\eta)}{p({\bf y})}\right]$ $\displaystyle=$
$\displaystyle\int p({\bf x})w({\bf y}|{\bf x};\eta)\log\frac{w({\bf y}|{\bf
x};\eta)}{w({\bf y};\eta)}d{\bf y}d{\bf x}$ $\displaystyle+\mathcal{KL}(w({\bf
y};\eta)\|p({\bf y}))$ $\displaystyle=$ $\displaystyle\int p({\bf x},{\bf
y};\eta)\log\frac{p({\bf x},{\bf y};\eta)}{p({\bf x})w({\bf y};\eta)}d{\bf
y}d{\bf x}$ $\displaystyle+\mathcal{KL}(w({\bf y};\eta)\|p({\bf y}))$
$\displaystyle=$ $\displaystyle I({\bf x},{\bf y};\eta)+\mathcal{KL}(w({\bf
y};\eta)\|p({\bf y})),$

where $w({\bf y};\eta)=\mathbb{E}_{p({\bf x})}[w({\bf y}|{\bf x};\eta)]$ and
$p({\bf x},{\bf y};\eta)=p({\bf x})w({\bf y}|{\bf x};\eta)$.

For the second term in (19), we define $q({\bf z}|{\bf y};\eta,\xi)$ as

$\displaystyle\mathbb{E}_{p({\bf x}|{\bf y};\eta)}[q({\bf z}|{\bf x},{\bf
y};\xi)]$ $\displaystyle=$ $\displaystyle\int p({\bf x}|{\bf y};\eta)q({\bf
z}|{\bf x},{\bf y};\xi)d{\bf x}$ $\displaystyle=$
$\displaystyle\int\frac{p({\bf x})w({\bf y}|{\bf x};\eta)}{p({\bf
y})}\frac{q({\bf z},{\bf y}|{\bf x};\eta,\xi)}{w({\bf y}|{\bf x};\eta)}d{\bf
x}$ $\displaystyle=$ $\displaystyle\int\frac{q({\bf z},{\bf y},{\bf
x};\eta,\xi)}{p({\bf y})}d{\bf x}=\frac{q({\bf z},{\bf y};\eta,\xi)}{p({\bf
y})}$ $\displaystyle=$ $\displaystyle q({\bf z}|{\bf y};\eta,\xi)$

where $p({\bf x}|{\bf y};\eta)=p({\bf x})w({\bf y}|{\bf x};\eta)/p({\bf y})$.

Therefore,

$\displaystyle\mathbb{E}_{p({\bf x})}\Big{[}\sum_{k=1}^{K}w({\bf y}=k|{\bf
x};\eta)$ $\displaystyle\times\mathcal{KL}\Big{(}q({\bf z}|{\bf x},{\bf
y}=k;\xi)\|p({\bf z}|{\bf y}=k)\Big{)}\Big{]}$ $\displaystyle=$
$\displaystyle\mathbb{E}_{p({\bf x})w({\bf y}|{\bf x};\eta)q({\bf z}|{\bf
x},{\bf y};\xi)}\left[\log\frac{q({\bf z}|{\bf x},{\bf y};\xi)}{p({\bf z}|{\bf
y})}\right]$ $\displaystyle=$ $\displaystyle\mathbb{E}_{p({\bf x})w({\bf
y}|{\bf x};\eta)q({\bf z}|{\bf x},{\bf y};\xi)}\left[\log\frac{q({\bf z}|{\bf
x},{\bf y};\xi)q({\bf z}|{\bf y};\eta,\xi)}{p({\bf z}|{\bf y})q({\bf z}|{\bf
y};\eta,\xi)}\right]$ $\displaystyle=$ $\displaystyle\mathbb{E}_{p({\bf
x})w({\bf y}|{\bf x};\eta)q({\bf z}|{\bf x},{\bf
y};\xi)}\left[\log\frac{q({\bf z}|{\bf x},{\bf y};\xi)}{q({\bf z}|{\bf
y};\eta,\xi)}\right]$ $\displaystyle+\mathbb{E}_{p({\bf x})w({\bf y}|{\bf
x};\eta)q({\bf z}|{\bf x},{\bf y};\xi)}\left[\log\frac{q({\bf z}|{\bf
y};\eta,\xi)}{p({\bf z}|{\bf y})}\right]$ $\displaystyle=$
$\displaystyle\mathbb{E}_{p({\bf y})p({\bf x}|{\bf y};\eta)q({\bf z}|{\bf
x},{\bf y};\xi)}\left[\log\frac{q({\bf z}|{\bf x},{\bf y};\xi)}{q({\bf z}|{\bf
y};\eta,\xi)}\right]$ $\displaystyle+\mathbb{E}_{p({\bf y})p({\bf x}|{\bf
y};\eta)q({\bf z}|{\bf x},{\bf y};\xi)}\left[\log\frac{q({\bf z}|{\bf
y};\eta,\xi)}{p({\bf z}|{\bf y})}\right]$ $\displaystyle=$
$\displaystyle\mathbb{E}_{p({\bf y})q({\bf x},{\bf z}|{\bf
y};\eta,\xi)}\left[\log\frac{q({\bf z}|{\bf x},{\bf y};\eta,\xi)}{q({\bf
z}|{\bf y};\eta,\xi)}\right]$ (21) $\displaystyle+\mathbb{E}_{p({\bf y})q({\bf
z}|{\bf y};\eta,\xi)}\left[\log\frac{q({\bf z}|{\bf y};\eta,\xi)}{p({\bf
z}|{\bf y})}\right]$ $\displaystyle=$ $\displaystyle\mathbb{E}_{p({\bf
y})}\int q({\bf x},{\bf z}|{\bf y};\eta,\xi)\log\frac{q({\bf x},{\bf z}|{\bf
y};\eta,\xi)}{q({\bf z}|{\bf y};\eta,\xi)p({\bf x}|{\bf y})}d{\bf x}d{\bf z}$
$\displaystyle+\mathbb{E}_{p({\bf y})}[\mathcal{KL}(q({\bf z}|{\bf
y};\eta,\xi)\|p({\bf z}|{\bf y}))]$ $\displaystyle=$
$\displaystyle\mathbb{E}_{p({\bf y})}[I({\bf x},{\bf z}|{\bf y};\eta,\xi)]$
$\displaystyle+\mathbb{E}_{p({\bf y})}[\mathcal{KL}(q({\bf z}|{\bf
y};\eta,\xi)\|p({\bf z}|{\bf y}))],$

where $p({\bf x}|{\bf y};\eta)p({\bf y})=p({\bf x})w({\bf y}|{\bf x};\eta)$
and $q({\bf x},{\bf z}|{\bf y};\eta,\xi)=p({\bf x}|{\bf y};\eta)q({\bf z}|{\bf
x},{\bf y};\xi)$.

In conclusion, the expectation of the KL-divergence upper bound is written in
mutual information terms as

$\displaystyle\mathbb{E}_{p({\bf x})}[\mathcal{KL}^{u}\left(q({\bf z}|{\bf
x};\eta,\xi)\|p({\bf z})\right)]$ (22) $\displaystyle=$ $\displaystyle I({\bf
x},{\bf y};\eta)+\mathcal{KL}(w({\bf y};\eta)\|p({\bf y}))$
$\displaystyle+\mathbb{E}_{p({\bf y})}[I({\bf x},{\bf z}|{\bf y};\eta,\xi)]$
$\displaystyle+\mathbb{E}_{p({\bf y})}[\mathcal{KL}(q({\bf z}|{\bf
y};\eta,\xi)\|p({\bf z}|{\bf y}))].$

###  A.2 MNIST Dataset

#### Evaluation with various tuning parameter $\beta$s

$\beta$ | 0.1 | 0.25 | 0.5 | 0.75 | 1 | 5 | 10 | 50
---|---|---|---|---|---|---|---|---
negative SSIM | 0.43 | 0.438 | 0.44 | 0.445 | 0.435 | 0.439 | 0.418 | 0.316
KL-divergence | 19.654 | 9.654 | 8.481 | 7.684 | 6.75 | 4.064 | 2.729 | 1.855
error(%) | 3.29 | 3.23 | 3.25 | 3.33 | 3.46 | 3.74 | 4.07 | 3.38
Table 3: The negative SSIM, KL-divergence, and the classification error of 100
labeled dataset for various $\beta$ values.

Denote the set of indices for the test dataset $I_{test}$. The KL-divergence
measures how different the trained posterior distribution is from the pre-
designed prior distribution and is defined by:

$\displaystyle\frac{1}{|I_{test}|}\sum_{i\in
I_{test}}\mathcal{KL}^{u}\left(q({\bf z}|x_{i};\eta,\xi)\|p({\bf z})\right).$
(23)

The negative average single-scale structural similarity (negative SSIM) of
${\bf X}$ is defined by

$\displaystyle\frac{1}{2}\bigg{(}1-\frac{1}{|{\bf X}|^{2}}\sum_{({\bf x},{\bf
x}^{\prime})\in{\bf X}\times{\bf X}}\mbox{SSIM}({\bf x},{\bf
x}^{\prime})\bigg{)},$ (24)

where $\mbox{SSIM}({\bf x},{\bf x}^{\prime})$ for ${\bf x},{\bf
x}^{\prime}\in{\bf X}$ is the similarity measure between two images ${\bf x}$,
${\bf x}^{\prime}$ and has a value on $[-1,1]$ (Wang et al. 2004; Zheng and
Sun 2019). So, negative SSIM has a value on $[0,1]$ and indicates how many
diverse images ${\bf X}$ consists of. For ${\bf X}$ being the set of images
generated from a fitted VAE model, $\mbox{negative SSIM}({\bf X})$ indicates
how expressive the model is. We consist ${\bf X}$ with the images generated
from our trained decoder, where each image is produced by 31$\times$31 equally
spaced grid points on the 2-dimensional latent space.

The classification error is given by

$\displaystyle\frac{1}{|I_{test}|}\sum_{i\in
I_{test}}\mathbb{I}\left(y_{i}\neq\arg\max_{k}(w({\bf
y}=k|x_{i};\eta))\right),$ (25)

where $\mathbb{I}(\cdot)$ is a indicator function. The classification error
shows the degree of discrepancy between the assigned label by the posterior
probability and the true observation label. Thus, the VAE model with low
classification error can separate the latent space into subsets on which the
labels of observations are well identified.

Table 3 indicates that the KL-divergence mostly depends on $\beta$. Because
$\beta$ indirectly controls the weight of the KL-divergence loss term, the KL-
divergence for a large $\beta$ dominates in our objective function (Lucas et
al. 2019). Also, we can observe that the diversity of generated sample
(negative SSIM) increases as $\beta$ decreases. For all $\beta$s, the latent
space is clearly separated according to the true labels. It means that the
classification performance does not depend on $\beta$ since the weight of
additional classification loss terms is scaled by $1/\beta$.

#### Controlling Proximity and Interpolation

Figure 6: Interpolated images produced by various pre-designs of the prior
distribution. From top to bottom, the distance between the center of the two
components increases as 8, 16, 24, and 32.

We investigate the patterns of interpolated images according to various pre-
designs of the proximity of the prior mixture components. We use a subset of
the MNIST dataset with only 0 and 1 labels, so the 2-component Gaussian
mixture distribution is used. All Gaussian components have diagonal covariance
matrices, and their diagonal elements are all 4. We set the location
parameters as $(-r,0),(r,0)$, which determines the proximity. We set the
distances between the two location parameters as 8, 16, 24, and 32. The images
were generated from 11 points of equally spacing line segment from $(-10,0)$
to $(10,0)$ on the 2-dimensional latent space. All experiment settings are
equally given by those of MNIST dataset analysis in Section 4.1. Figure 6
shows that the interpolated images vary more slowly if the two location
parameters of $p({\bf z})$ are farther from each other and the latent space is
effectively adapted to pre-designed characteristics.

encoder | decoder
---|---
(28, 28, 1) image |  $d$-dimension latent variable
Flatten | Dense(128), ReLU
Dense(256), ReLU | Dense(256), ReLU
Dense(128), ReLU | Dense(784), ReLU
2 $\times$ K $\times$ Dense($d$), Linear | Reshape(28, 28, 1)
Table 4: Model descriptions of the encoder and decoder. Here, $K$ denotes the
number of classes, and $d$ is the dimension of the latent variable. classifier
---
input: 28 $\times$ 28 $\times$ 1 Image
Conv2D(32, 5, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
MaxPool2D(pool size=(2, 2), strides=2)
SpatialDropout2D(rate=0.5)
Conv2D(64, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
MaxPool2D(pool size=(2, 2), strides=2)
SpatialDropout2D(rate=0.5)
Conv2D(128, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
MaxPool2D(pool size=(2, 2), strides=2)
SpatialDropout2D(rate=0.5)
GlobalAveragePooling2D
Dense(64), BN, ReLU
Dense(K), softmax
Table 5: Network structures of the classifier. Here, $K$ denotes the number of classes, and BN is Batch-Normalization layer (Ioffe and Szegedy 2015). Conv2D hyper-parameters are filters, kernel size, strides, and padding. models | error(%)
---|---
M2(2014) | 11.97
M1+M2(2014) | 3.33
ADGM(2016) | 0.96
Disentangled-VAE(2019) | 2.71
Parted-VAE(2021) | 9.79
SHOT-VAE(2021) | 3.12
Ours ($\beta=0.25$) | 3.23
Table 6: Semi-supervised test classification errors with 100 labeled datasets
(except for Parted-VAE, all results are from original papers).

###  A.3 CIFAR-10 Dataset

#### V-nat of the EXoN

Figure 7: Visualization of 256-dimensional V-nat vector of the EXoN for all
classes where $\beta=0.05$. Figure 8: Correlation matrix between V-nat vectors
of all classes with $\beta=0.05$

Figure 7 shows V-nat vectors for all classes in CIFAR-10 datasets, and $j$th
element of V-nat vector is $\log\mathbb{E}_{{\bf x}|{\bf
y}=k}\left[\mbox{Var}({\bf z}_{j}^{k})/\mbox{Var}({\bf z}_{j}^{k}|{\bf
x};\xi)\right]$. We found that the activated latent subspace of each class
(latent dimensions which have a higher V-nat value than $\delta=1$) does not
differ significantly from each other.

To show that the variability of generated images depends on the almost same
latent subspace for all classes, the correlation matrix between V-nat vectors
of all classes with $\beta=0.05$ is shown in Figure 8. The correlation plot
indicates that the V-nat vectors are remarkably correlated (correlation
coefficients are almost 1), implying that the diversity of generated images
depends on almost the same subspace.

#### Manipulation on the activated latent subspace

(a) $\beta=0.01$

(b) $\beta=0.1$

Figure 9: Manipulation on the single activated latent subspace and generated
images.

(a) $\beta=0.01$

(b) $\beta=0.1$

Figure 10: Manipulation on all activated latent subspace and generated images.

In Figure 9, 5 series of images are generated from latent variables where the
value of the activated latent axis is changed from -3 to 3, and values of
other latent dimensions are fixed. From top to bottom, latent axes having the
top 5 most significant V-nat values are used for each $\beta$. As the value of
each activated latent subspace changes, some different features of generated
samples are varied (like brightness, the color and the shape of the object, or
the color of the background).

Note that characteristics of images change more significantly where the
decoder variance parameter $\beta$ is large because a single activated latent
dimension determines relatively more image characteristics due to the small
cardinality of the activated latent subspace when $\beta$ is large.
Disentangling properties represented by the activated latent subspace is left
for our future work.

Furthermore, we manipulated all activated latent subspace, and it is
visualized in Figure 10 for each $\beta$. It shows images generated from
latent variables perturbed by the uniform noise on ${\cal A}_{k}(1)$, where
the uniform noise is sampled from $U(-1.5,1.5)$ and automobile class $k$.

Most tops left images are reconstructed images given the unperturbed latent
variable. As in manipulating the single activated latent subspace, this Figure
confirms that a single activated latent dimension determines relatively more
characteristics of a given image because the cardinality of the activated
latent subspace is small for a large $\beta$.

Figure 11: Recovered images given the train dataset where $\beta=0.05$. It
confirms that the proposed model could recover an original observation without
loss of information. Figure 12: Recovered images given the test dataset where
$\beta=0.05$.

###  A.4 Experiment settings

We run all experiments using Geforce GTX 2080 Ti GPU and 16GB RAM, and our
experimental codes are all available with Tensorflow (Abadi et al. 2015). Pre-
designs of the prior distribution:

  * •

MNIST: The conceptual centers and variances are given by
$\mu_{k}^{0}=r\cdot\big{(}\cos(\frac{\pi}{10}(2k-1)),\sin(\frac{\pi}{10}(2k-1))\big{)}$,
where $r=4/\sin(\frac{\pi}{10})$ and $(\sigma_{k}^{0})^{2}=(4,4)$ for
$k=1,\cdots,10$.

  * •

CIFAR-10: The mean vector of $p({\bf z})$ consists of two subvectors, a
10-dimensional label-relevant mean which is the one-hot vector for the
corresponding label, and a 246-dimensional label-irrelevant mean which has all
zero values (Zheng and Sun 2019). All covariance matrices are commonly given
by $diag\\{({0.1}_{10},{1}_{246})\\}$ where ${a}_{k}$ for ${a}\in\mathbb{R}$
and $k\in\mathbb{N}$ is a $k$-dimensional vector whose elements are all $a$.

Stochastic image augmentations we used:

  * •

MNIST: random rotation, random cropping

  * •

CIFAR-10: random horizontal flip, random cropping

dataset | epochs | batch size($U,L$) | drop rate | $\lambda_{1}$ | $\lambda_{2}(e)$
---|---|---|---|---|---
MNIST | 100 | (128, 32) | 0.5 | 6000 | $\lambda_{1}\exp\left(-5(1-\min\\{1,\frac{c}{10}\\})^{2}\right)$
CIFAR-10 | 600 | (128, 32) | 0.1 | 5000 | $\lambda_{1}\exp\left(-5(1-\min\\{1,\frac{c}{50}\\})^{2}\right)$
Table 7: Detailed experiment settings for the MNIST and the CIFAR-10 datasets.

In implementation, we used Adam(Kingma and Ba 2014) optimizer for both
datasets. The initial learning rate of the MNIST experiments is 0.003, and we
decayed the learning rate exponentially by multiplying
$\exp\left(-5(1-\frac{100-c}{100})^{2}\right)$ after 10 epoch, where $c$ is
current epoch number. And the initial learning rate of CIFAR-10 experiments is
0.001, and we multiplied the learning rate by 0.1 at 250, 350, 450, and 550
epochs. After 550 epoch, we decayed the learning rate exponentially by
multiplying $\exp\left(-5(1-\frac{600-c}{600-550})^{2}\right)$. For both
datasets, we applied decoupled weight decay (Loshchilov and Hutter 2017) with
a factor of 0.0005. Other experiment settings are shown in Table 7.

encoder | decoder
---|---
input: 32 $\times$ 32 $\times$ 3 Image | input: $d$-dimension latent variable
Conv2D(32, 5, 2, ‘same’), BN, ReLU | Dense(4 $\times$ 4 $\times$ 512), BN, ReLU
Conv2D(64, 4, 2, ‘same’), BN, ReLU | Reshape(4, 4, 512)
Conv2D(128, 4, 2, ‘same’), BN, ReLU | Conv2DTranspose(256, 5, 2, ‘same’), BN, ReLU
Conv2D(256, 4, 1, ‘same’), BN, ReLU | Conv2DTranspose(128, 5, 2, ‘same’), BN, ReLU
Flatten | Conv2DTranspose(64, 5, 2, ‘same’), BN, ReLU
Dense(1024), BN, Linear | Conv2DTranspose(32, 5, 1, ‘same’), BN, ReLU
mean: K $\times$ Dense($d$), Linear | Conv2D(3, 4, 1, ‘same’), tanh
log-variance: K $\times$ Dense($d$), Linear | -
Table 8: Model descriptions of the encoder and decoder. Here, $K$ denotes the
number of classes, $d$ is the dimension of the latent variable, and BN is
Batch-Normalization layer (Ioffe and Szegedy 2015). Conv2D and Conv2DTranspose
hyper-parameters are filters, kernel size, strides, and padding. classifier
---
input: 32 $\times$ 32 $\times$ 3 Image
Conv2D(128, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
Conv2D(128, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
Conv2D(128, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
MaxPool2D(pool size=(2, 2), strides=2)
SpatialDropout2D(rate=0.1)
Conv2D(256, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
Conv2D(256, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
Conv2D(256, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
MaxPool2D(pool size=(2, 2), strides=2)
SpatialDropout2D(rate=0.1)
Conv2D(512, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
Conv2D(256, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
Conv2D(128, 3, 1, ‘same’), BN, LeakyReLU($\alpha$=0.1)
GlobalAveragePooling2D
Dense(K), softmax
Table 9: Network structures of the classifier. Here, $K$ denotes the number of
classes, and BN is Batch-Normalization layer (Ioffe and Szegedy 2015). Conv2D
hyper-parameters are filters, kernel size, strides, and padding.