Autoencoders

From deterministic autoencoders to variational and vector-quantised models, and the gradient-estimation tricks (reparameterisation, Gumbel-Softmax, REINFORCE) that make them trainable.

intermediate deep learning unsupervised-learninggenerative-modelsrepresentation-learning
autoencodersvariational-inferencevaevq-vaereparameterisationgumbel-softmaxreinforceelbokl-divergence

Author @Zak-T Updated 2026-06-25

Overview

An autoencoder is a neural network trained to copy its input back to its output, but through a deliberate bottleneck. We compress an input xx with an encoder into a compact latent code zz, and then hand that code to a decoder that tries to reconstruct the input as x′x'. The whole network is trained so that x′x' resembles xx as closely as possible.

The bottleneck is the point. If the code zz is lower-dimensional (or otherwise constrained) than the input, the network cannot simply memorise and pass the input through unchanged. It is forced to discard noise and keep only the structure that matters for reconstruction. A good latent is therefore one that has learned and stored the important features of the input. Because we never need labels — the input is its own target — autoencoders are an unsupervised (or self-supervised) method.

Notation. Throughout the variational sections we follow the convention of the original VAE paper: ϕ\phi are the encoder parameters (the recognition / inference network qϕq_\phi), and θ\theta are the decoder parameters (the generative network pθp_\theta). Some texts swap these letters, so watch for it when reading other sources.


The Standard (Vanilla) Autoencoder

Encoder

Each layer of the encoder is an ordinary neural-network layer. For weights WW, bias bb, and an activation function gg (e.g. ReLU or sigmoid), a layer’s output is:

z=g(Wx+b)z = g(Wx + b)

Stacking several such layers maps the input xx down to the latent code zz.

Decoder

The decoder mirrors this, mapping the code back up to a reconstruction. With its own weights W′W', bias b′b', and activation g′g':

x′=g′(W′z+b′)x' = g'(W'z + b')

Loss function

We want x′x' to match xx, so we penalise their squared distance (the reconstruction error):

L(x,x′)=∥x−x′∥2\mathcal{L}(x, x') = \lVert x - x' \rVert^2

Substituting the decoder, and then the encoder, makes the full composition explicit:

L(x,x′)=∥x−g′(W′z+b′)∥2=∥x−g′(W′ g(Wx+b)+b′)∥2\mathcal{L}(x, x') = \lVert x - g'(W'z + b') \rVert^2 = \lVert x - g'\big(W'\, g(Wx + b) + b'\big) \rVert^2

This is fully deterministic and differentiable, so it trains by ordinary gradient descent.

Why the latent space is a problem

A standard autoencoder is trained on a finite, discrete set of data points, so it only ever learns to place those points somewhere in latent space. This creates two issues:

Both problems stem from the same root cause: the encoder maps each input to a single point, with no pressure to organise the space smoothly. This is exactly what the variational autoencoder fixes.


From a Point to a Distribution: the Variational Autoencoder

Deterministic vs variational encoders

A deterministic encoder takes an input xx and returns one latent vector zz — same xx in, same zz out, every time.

A variational encoder does not output a code directly. Instead it outputs the parameters of a probability distribution over zz. For a Gaussian latent, the encoder reads xx and produces a mean vector μ(x)\mu(x) and a variance σ2(x)\sigma^2(x). The actual code is then a sample from the distribution those parameters define:

z∼N ⁣(μ(x), σ2(x))z \sim \mathcal{N}\!\big(\mu(x),\, \sigma^2(x)\big)

Why is zz a random variable?

Because we built it to be one. The encoder produces a distribution, and zz is a draw from it.

Why design zz to be random?

There are two complementary reasons.

The Bayesian reason. Given an observed xx, which latent codes zz could have produced it? Rarely a single one. Many codes could plausibly have generated the same image, and some are more likely than others. The spread of plausible explanations is a distribution — the posterior p(z∣x)p(z \mid x). The posterior is inherently a distribution because it encodes our uncertainty about the hidden cause. The encoder qϕ(z∣x)q_\phi(z \mid x) is our learned approximation to that posterior, so it has to be a distribution too, and a sample zz from it has to be random.

The practical reason. Because zz is sampled with genuine noise around μ(x)\mu(x), and because a regularising term (the KL term, below) continually pulls qϕ(z∣x)q_\phi(z\mid x) toward a smooth prior p(z)p(z), the latent space gets organised into a continuous, well-behaved region. Nearby points decode to similar outputs, and there are no “dead” regions that decode to nonsense. This is what lets us sample a fresh zz from the prior and decode a brand-new plausible data point, or smoothly interpolate between two encodings. A deterministic autoencoder gives no such guarantee.


Bayesian Foundations

The encoder is learning to infer the latent zz from an observation xx — exactly the situation Bayes’ theorem describes:

p(z∣x)=p(x∣z) p(z)p(x)p(z \mid x) = \frac{p(x \mid z)\, p(z)}{p(x)}

Each piece has a name and a role:

Writing the denominator out shows where the trouble starts:

p(z∣x)=p(x∣z) p(z)∑zp(x∣z) p(z)p(z \mid x) = \frac{p(x \mid z)\, p(z)}{\sum_{z} p(x \mid z)\, p(z)}

To compute p(x)p(x) we would have to sum (or integrate) over every possible value of zz. For any realistic latent space this is intractable. Since we cannot compute the exact posterior, we approximate it — and the tool for that is variational inference, leading to the ELBO.


Variational Inference and the ELBO

A family of approximate posteriors

We pick a family of tractable distributions Q\mathcal{Q} and search within it for the member that best approximates the true posterior. We index a member by a variational parameter λ\lambda:

qλ(z∣x)q_\lambda(z \mid x)

If the family is Gaussian, then λ\lambda is just the mean and variance: knowing those pins down the exact distribution. For the ii-th data point xix_i:

λxi=(μxi, σxi2)\lambda_{x_i} = \big(\mu_{x_i},\, \sigma^2_{x_i}\big)

Because λ\lambda is chosen per data point, qλ(z)=qλ(z∣x)q_\lambda(z) = q_\lambda(z \mid x).

Measuring the approximation with KL divergence

The Kullback–Leibler (KL) divergence measures how far one distribution is from another. It is the expected log-ratio of the two densities; it is always ≥0\ge 0 and equals 00 only when the distributions are identical. We use it to ask: how well does qλ(z)q_\lambda(z) approximate p(z∣x)p(z\mid x)?

DKL(qλ(z) ∥ p(z∣x))=∑zqλ(z)log⁡qλ(z)p(z∣x)=Ez∼q ⁣[log⁡qλ(z)p(z∣x)]D_{\mathrm{KL}}\big(q_\lambda(z) \,\|\, p(z \mid x)\big) = \sum_z q_\lambda(z) \log \frac{q_\lambda(z)}{p(z \mid x)} = \mathbb{E}_{z \sim q}\!\left[ \log \frac{q_\lambda(z)}{p(z \mid x)} \right]

Now we expand it, one step at a time. Split the log of a ratio into a difference of logs:

=Ez∼q ⁣[log⁡qλ(z)]−Ez∼q ⁣[log⁡p(z∣x)]= \mathbb{E}_{z \sim q}\!\left[\log q_\lambda(z)\right] - \mathbb{E}_{z \sim q}\!\left[\log p(z \mid x)\right]

Rewrite the posterior using p(z∣x)=p(z,x)p(x)p(z \mid x) = \tfrac{p(z, x)}{p(x)}:

=Ez∼q ⁣[log⁡qλ(z)]−Ez∼q ⁣[log⁡p(z,x)p(x)]= \mathbb{E}_{z \sim q}\!\left[\log q_\lambda(z)\right] - \mathbb{E}_{z \sim q}\!\left[\log \frac{p(z, x)}{p(x)}\right]

Split that second log too, and note log⁡p(x)\log p(x) does not depend on zz, so its expectation is just log⁡p(x)\log p(x):

=Ez∼q ⁣[log⁡qλ(z)]−Ez∼q ⁣[log⁡p(z,x)]+log⁡p(x)= \mathbb{E}_{z \sim q}\!\left[\log q_\lambda(z)\right] - \mathbb{E}_{z \sim q}\!\left[\log p(z, x)\right] + \log p(x)

Defining the ELBO

The first two terms are exactly the negative of what we call the Evidence Lower BOund (ELBO):

−ELBO(λ)=Ez∼q ⁣[log⁡qλ(z)]−Ez∼q ⁣[log⁡p(z,x)]-\mathrm{ELBO}(\lambda) = \mathbb{E}_{z \sim q}\!\left[\log q_\lambda(z)\right] - \mathbb{E}_{z \sim q}\!\left[\log p(z, x)\right]

Substituting that name back in gives the central identity of variational inference:

log⁡p(x)=ELBO(λ)+DKL(qλ(z) ∥ p(z∣x))\log p(x) = \mathrm{ELBO}(\lambda) + D_{\mathrm{KL}}\big(q_\lambda(z) \,\|\, p(z \mid x)\big)

Rearranged (this is the correct direction of the rearrangement):

DKL(qλ(z) ∥ p(z∣x))=log⁡p(x)−ELBO(λ)D_{\mathrm{KL}}\big(q_\lambda(z) \,\|\, p(z \mid x)\big) = \log p(x) - \mathrm{ELBO}(\lambda)

Why this is the whole point

Look at the identity carefully:

So ELBO(λ)≤log⁡p(x)\mathrm{ELBO}(\lambda) \le \log p(x) — it is genuinely a lower bound on the evidence, hence the name. And because log⁡p(x)\log p(x) is constant, maximising the ELBO is exactly the same as minimising the KL divergence between our approximation and the true posterior. The optimal approximation is:

qλ∗(z)=arg min⁡q∈QDKL(qλ(z) ∥ p(z∣x))=arg max⁡q∈QELBO(λ)q^*_\lambda(z) = \operatorname*{arg\,min}_{q \in \mathcal{Q}} D_{\mathrm{KL}}\big(q_\lambda(z) \,\|\, p(z \mid x)\big) = \operatorname*{arg\,max}_{q \in \mathcal{Q}} \mathrm{ELBO}(\lambda)

The intractable p(x)p(x) has dropped out of anything we need to optimise. We can differentiate the ELBO and run gradient descent. That is the breakthrough that makes the VAE trainable.

The per-data-point ELBO: reconstruction + regulariser

We assume local latent variables: each data point xx has its own latent zz, shared with no other point. The ELBO then splits into a sum over data points, and using p(z,x)=p(x∣z) p(z)p(z, x) = p(x \mid z)\,p(z) we can rewrite each term as:

ELBOi(λ)=Ez∼qλ(z∣xi) ⁣[log⁡p(xi∣z)]−DKL(qλ(z∣xi) ∥ p(z))\mathrm{ELBO}_i(\lambda) = \mathbb{E}_{z \sim q_\lambda(z \mid x_i)}\!\big[\log p(x_i \mid z)\big] - D_{\mathrm{KL}}\big(q_\lambda(z \mid x_i) \,\|\, p(z)\big)

(The derivation: Eq[log⁡p(x∣z)+log⁡p(z)]−Eq[log⁡q]=Eq[log⁡p(x∣z)]−DKL(q∥p(z))\mathbb{E}_q[\log p(x\mid z) + \log p(z)] - \mathbb{E}_q[\log q] = \mathbb{E}_q[\log p(x\mid z)] - D_{\mathrm{KL}}(q\|p(z)).)

Since the encoder is the network qϕq_\phi and the decoder is pθp_\theta, training minimises the negative ELBO, which is the VAE loss:

Li(θ,ϕ)=−Ez∼qϕ(z∣xi) ⁣[log⁡pθ(xi∣z)]⏟reconstruction term+DKL(qϕ(z∣xi) ∥ p(z))⏟regulariser\mathcal{L}_i(\theta, \phi) = -\underbrace{\mathbb{E}_{z \sim q_\phi(z \mid x_i)}\!\big[\log p_\theta(x_i \mid z)\big]}_{\text{reconstruction term}} + \underbrace{D_{\mathrm{KL}}\big(q_\phi(z \mid x_i) \,\|\, p(z)\big)}_{\text{regulariser}}

Both θ\theta and ϕ\phi appear, so a single gradient step on this loss updates both networks at once.

Computing the Gaussian loss in closed form

For diagonal-Gaussian qϕ=N(μ,σ2)q_\phi = \mathcal{N}(\mu, \sigma^2) (dimension JJ) against the standard-normal prior p(z)=N(0,I)p(z) = \mathcal{N}(0, I), the KL term has a clean closed form — no sampling needed:

DKL(N(μ,σ2) ∥ N(0,I))=12∑j=1J(μj2+σj2−log⁡σj2−1)D_{\mathrm{KL}}\big(\mathcal{N}(\mu, \sigma^2) \,\|\, \mathcal{N}(0, I)\big) = \frac{1}{2}\sum_{j=1}^{J}\Big(\mu_j^2 + \sigma_j^2 - \log \sigma_j^2 - 1\Big)

Reading the terms: μj2\mu_j^2 pulls the mean toward 00; σj2−log⁡σj2−1\sigma_j^2 - \log\sigma_j^2 - 1 is minimised at σj2=1\sigma_j^2 = 1, pulling the variance toward 11. Together they pull each latent dimension toward the unit Gaussian.

The reconstruction term depends on the decoder’s chosen output distribution:

[TODO] Full derivation of the reconstruction term from −log⁡pθ(x∣z)-\log p_\theta(x\mid z) for the Gaussian and Bernoulli cases, including the role of the (often fixed) decoder variance.

Connection to Expectation–Maximisation

The ELBO is not unique to VAEs — it is the same object that powers the EM algorithm. EM maximises the marginal likelihood log⁡pθ(X)\log p_\theta(X), where the latent variables ZZ have been marginalised out (pθ(X)=∫pθ(X,Z) dZp_\theta(X) = \int p_\theta(X, Z)\,dZ; note it is ZZ, not θ\theta, that is integrated away). For any distribution q(Z)q(Z) over the latents:

log⁡pθ(X)=Eq ⁣[log⁡pθ(X,Z)]+H(q)⏟F(q,θ) = ELBO+DKL(q(Z) ∥ pθ(Z∣X))\log p_\theta(X) = \underbrace{\mathbb{E}_q\!\big[\log p_\theta(X, Z)\big] + H(q)}_{F(q, \theta)\ =\ \text{ELBO}} + D_{\mathrm{KL}}\big(q(Z) \,\|\, p_\theta(Z \mid X)\big)

where H(q)=−Eq[log⁡q(Z)]H(q) = -\mathbb{E}_q[\log q(Z)] is the entropy of qq. The two parts are the ELBO F(q,θ)F(q,\theta) and the KL gap between qq and the true posterior. EM alternates:

A VAE is essentially amortised variational EM: instead of solving for the optimal qq exactly for every data point (a per-point E-step), an encoder network qϕq_\phi predicts the variational parameters in one forward pass, and both ϕ\phi and θ\theta are improved together by gradient ascent on the ELBO.


Computing the Gradients: Why Backprop Breaks

Training is gradient ascent on the ELBO (equivalently, descent on the negative-ELBO loss above). The decoder θ\theta is no trouble: once a value of zz is given, the gradient flows straight through log⁡pθ(x∣z)\log p_\theta(x\mid z). The encoder ϕ\phi is the problem, because ϕ\phi lives inside the distribution we sample from, z∼qϕ(z∣x)z \sim q_\phi(z\mid x).

Two conditions backprop requires

For gradients to flow, every node from inputs to loss must be:

  1. a deterministic function of its inputs, and
  2. differentiable.

A deterministic node, e.g. y=3μy = 3\mu. Ask: if I nudge μ\mu up by δ\delta, how does yy change? Answer: yy goes up by exactly 3δ3\delta, every single time. The map μ↦y\mu \mapsto y is a genuine function with a well-defined derivative ∂y∂μ=3\tfrac{\partial y}{\partial \mu} = 3.

A sampling node, e.g. z∼N(μ,σ2)z \sim \mathcal{N}(\mu, \sigma^2). Ask the same question: if I nudge μ\mu up by δ\delta, how does zz change? There is no answer. For a fixed μ\mu, zz is not a single number — it is a freshly drawn random value that could land anywhere under the bell curve. The map μ↦z\mu \mapsto z is one-to-many: one input, infinitely many possible outputs. No derivative exists.

Why not just differentiate the expectation directly?

The loss is an expectation, which we can write as an integral. (Using a generic ϕ\phi for the parameters of the sampling distribution and f(z)=log⁡pθ(x∣z)f(z) = \log p_\theta(x\mid z) for the thing inside.)

L(ϕ)=Ez∼qϕ(z∣x) ⁣[f(z)]=∫qϕ(z∣x) f(z) dzL(\phi) = \mathbb{E}_{z \sim q_\phi(z\mid x)}\!\big[f(z)\big] = \int q_\phi(z \mid x)\, f(z)\, dz

Push the gradient inside and apply the product rule:

∇ϕL=∫∇ϕ[qϕ(z∣x) f(z)] dz=∫([∇ϕqϕ(z∣x)]f(z)+qϕ(z∣x)[∇ϕf(z)]) dz\nabla_\phi L = \int \nabla_\phi\big[q_\phi(z \mid x)\, f(z)\big]\, dz = \int \Big(\big[\nabla_\phi q_\phi(z\mid x)\big] f(z) + q_\phi(z\mid x)\big[\nabla_\phi f(z)\big]\Big)\, dz

Here is the key observation: f(z)=log⁡pθ(x∣z)f(z) = \log p_\theta(x \mid z) depends only on the decoder θ\theta and the value zz. The encoder ϕ\phi does not appear in ff at all, so:

∇ϕf(z)=0\nabla_\phi f(z) = 0

and the whole expression collapses to:

∇ϕL=∫[∇ϕqϕ(z∣x)] f(z) dz\nabla_\phi L = \int \big[\nabla_\phi q_\phi(z \mid x)\big]\, f(z)\, dz

💡 This is why backprop breaks. The encoder ϕ\phi exerts its influence entirely through the shape of the distribution qϕq_\phi that we sample from. The gradient therefore lives inside the qϕq_\phi object — a density that measures how the parameters reshape the latent space — while f(z)f(z) only answers “how well do we reconstruct if the latent happens to be zz?”. If we replace “the distribution” with “one concrete sample zz”, we throw away every trace of ϕ\phi: the sampled zz is just a number, and ϕ\phi has vanished from the computation graph. The autograd engine, seeing only the graph after zz was sampled, faithfully reports a zero gradient — because in the graph it can see, there is no path from ϕ\phi to the loss.

The differentiation graph can only begin after zz has been sampled, and from there it differentiates ff, which doesn’t depend on ϕ\phi. ϕ\phi’s influence is locked inside the sampling step, where backprop can’t reach.


The Reparameterisation Trick

The fix is to relocate ϕ\phi‘s influence out of the sampling step and into a deterministic node. Instead of drawing z∼N(μϕ(x),σϕ2(x))z \sim \mathcal{N}(\mu_\phi(x), \sigma_\phi^2(x)) directly, we write:

z=μϕ(x)+σϕ(x)⊙ϵ,ϵ∼N(0,I)z = \mu_\phi(x) + \sigma_\phi(x) \odot \epsilon, \qquad \epsilon \sim \mathcal{N}(0, I)

Now the randomness lives entirely in ϵ\epsilon, an external input with no parameters. Given ϵ\epsilon, the latent zz is a deterministic, differentiable function of ϕ\phi (through μϕ\mu_\phi and σϕ\sigma_\phi). The graph now has a clean path ϕ→z→f→loss\phi \to z \to f \to \text{loss}, and ∇ϕ\nabla_\phi flows through it normally.

Why this is valid: the location–scale property

This works because of a specific property of the Gaussian, the location–scale property. The distribution

z∼N(μ,σ2)z \sim \mathcal{N}(\mu, \sigma^2)

is identical to

z=μ+σ⋅ϵ,ϵ∼N(0,1)z = \mu + \sigma \cdot \epsilon, \qquad \epsilon \sim \mathcal{N}(0, 1)

A standard normal can be shifted by μ\mu (location) and scaled by σ\sigma (scale) to produce any Gaussian. So the reparameterised samples follow exactly the same distribution as before — we have changed how we draw the sample, not what we are drawing — but the new form is differentiable in the parameters.

The resulting gradient is itself an expectation we can Monte-Carlo estimate by drawing a few ϵ\epsilon:

∇ϕ Ez∼qϕ ⁣[f(z)]=Eϵ∼N(0,I) ⁣[∇ϕf(μϕ+σϕ⊙ϵ)]\nabla_\phi\, \mathbb{E}_{z \sim q_\phi}\!\big[f(z)\big] = \mathbb{E}_{\epsilon \sim \mathcal{N}(0, I)}\!\big[\nabla_\phi f(\mu_\phi + \sigma_\phi \odot \epsilon)\big]

This estimator tends to have low variance, which is one reason it works so well in practice (contrast with the score-function estimator below).


Discrete Latents: the VQ-VAE

A vanilla VAE assumes the prior and posterior are continuous Gaussians. But sometimes a discrete latent is the better model.

Motivation. Images may belong to categories like “cat” and “car”, and it makes little sense to interpolate between such categories — a half-cat-half-car is not a meaningful image. Discrete codes are also easier to model: each category is a single value. The VQ-VAE (Vector-Quantised VAE) therefore replaces the continuous normal latent with a discrete one:

  1. the encoder produces continuous vectors which are matched to a categorical codebook,
  2. sampling returns an integer index into a dictionary of embedding vectors,
  3. the indexed embedding is passed to the decoder.

Dimensions

nn = batch size, hh = image height, ww = image width, cc = input channels, dd = hidden channels, kk = number of codebook entries.

The Vector-Quantisation layer

This layer takes the encoder’s continuous embeddings ze(x)z_e(x), snaps each to its nearest codebook entry, and outputs the quantised tensor zqz_q. The decoder then consumes zqz_q and outputs x′x'.

  1. Reshape. Combine all dimensions except the last, giving n⋅h⋅wn\cdot h\cdot w vectors, each of dimensionality dd.
  2. Compute distances. For each of the n⋅h⋅wn\cdot h\cdot w vectors, measure its distance to each of the kk codebook embeddings, producing a matrix of shape (n⋅h⋅w, k)(n\cdot h\cdot w,\, k).
  3. Argmin. For each vector, find the index of the nearest of the kk codebook vectors.
  4. Index. Replace each vector with its nearest codebook vector, giving zqz_q.
  5. Reshape back to (n,h,w,d)(n, h, w, d).
  6. Copy gradients. The arg min⁡\operatorname*{arg\,min} step has zero gradient almost everywhere, so backprop cannot flow through it. We approximate by copying the gradients arriving at zqz_q back onto zez_e, so some training signal still reaches the encoder.

The straight-through estimator

How do we “copy the gradients” and pretend the VQ layer isn’t there during backprop? Note that zez_e and zqz_q are vectors of the same dimensionality dd, and zqz_q is the codebook entry nearest to zez_e, so geometrically they sit close together in the same space. The straight-through estimator says: only in the backward pass, pretend the quantisation step was the identity function, and pass whatever gradient arrives at zqz_q straight through to zez_e.

The justification is that the copied gradient points in a sensible direction — “produce outputs that make the decoder happier”. The forward pass uses the true (discrete) quantisation; the backward pass uses the identity. It is an approximation, but it works well in practice.

The VQ-VAE loss (three components)

Here sg[⋅]\text{sg}[\cdot] is the stop-gradient operator: it is the identity in the forward pass but lets no gradient flow through whatever it wraps in the backward pass.

L=−log⁡p(x∣zq)⏟reconstruction  +  ∥sg[ze(x)]−e∥2⏟codebook (VQ) loss  +  β∥ze(x)−sg[e]∥2⏟commitment loss\mathcal{L} = \underbrace{-\log p(x \mid z_q)}_{\text{reconstruction}} \;+\; \underbrace{\big\lVert \text{sg}[z_e(x)] - e \big\rVert^2}_{\text{codebook (VQ) loss}} \;+\; \underbrace{\beta \big\lVert z_e(x) - \text{sg}[e] \big\rVert^2}_{\text{commitment loss}}

Reparameterising a Discrete Distribution

We just saw how the reparameterisation trick relocates randomness so we can differentiate through continuous latents. The Gumbel-Max trick is the analogue for a categorical distribution.

The core difficulty: the map from continuous parameters (class probabilities / logits) to a discrete outcome is a step function. Nudging the class probabilities by a tiny amount produces a discontinuous jump in the (one-hot) output — there is no smooth response to differentiate. For a continuous distribution the output can move continuously with the parameters; a discrete random variable cannot, because its output lives in a finite set.

The Gumbel-Max trick

Suppose the categorical is specified by logits α1,…,αk\alpha_1, \dots, \alpha_k (unnormalised log-probabilities). The trick gives a routine to draw a sample:

  1. Draw independent Gumbel noise g1,…,gkg_1, \dots, g_k, where each gi=−log⁡(−log⁡ui)g_i = -\log(-\log u_i) with ui∼Uniform(0,1)u_i \sim \text{Uniform}(0, 1).
  2. Output the index arg max⁡i(αi+gi)\operatorname*{arg\,max}_i (\alpha_i + g_i).

This produces an exact draw from the softmax categorical. But arg max⁡\operatorname*{arg\,max} is still a non-differentiable step function — so this alone doesn’t solve our gradient problem; it just isolates the non-differentiability into a single clean operation.

What is the Gumbel distribution?

It comes from extreme value theory. If you take the maximum of many independent random variables and rescale it appropriately, the result converges to one of three limiting shapes. When the underlying distribution has a light tail (decaying exponentially or faster) — like the normal, the exponential, or the logistic — the maximum converges to the Gumbel distribution. Since the trick is built around an arg max⁡\operatorname*{arg\,max} over (noised) logits, the Gumbel is the natural noise to inject.

Why are there two logs?

This is inverse-transform sampling: a way to turn plain uniform randomness into a sample from any target distribution. If U∼Uniform(0,1)U \sim \text{Uniform}(0,1) and FF is the CDF (cumulative distribution function) of the target, then X=F−1(U)X = F^{-1}(U) has exactly that distribution. The proof is short:

P(X≤x)=P(F−1(U)≤x)=P(U≤F(x))=F(x)P(X \le x) = P\big(F^{-1}(U) \le x\big) = P\big(U \le F(x)\big) = F(x)

The last step holds because UU is uniform, so P(U≤t)=tP(U \le t) = t. The Gumbel’s CDF is F(x)=exp⁡(−exp⁡(−x))F(x) = \exp(-\exp(-x)). Inverting it (set u=exp⁡(−exp⁡(−x))u = \exp(-\exp(-x)) and solve for xx) gives:

F−1(u)=−log⁡(−log⁡u)F^{-1}(u) = -\log(-\log u)

— hence the two nested logs.

[TODO] Proof of the Gumbel-Max identity: that   P(arg max⁡i(αi+gi)=j)=eαj∑keαk  \;P\big(\operatorname*{arg\,max}_i (\alpha_i + g_i) = j\big) = \dfrac{e^{\alpha_j}}{\sum_k e^{\alpha_k}}\;, i.e. the routine really does sample from the intended categorical.

The Gumbel-Softmax (Concrete) trick

To make the sample differentiable, replace the hard arg max⁡\operatorname*{arg\,max} with a temperature-controlled softmax:

yi=softmax⁡ ⁣(αi+giτ)=exp⁡ ⁣((αi+gi)/τ)∑jexp⁡ ⁣((αj+gj)/τ)y_i = \operatorname{softmax}\!\left(\frac{\alpha_i + g_i}{\tau}\right) = \frac{\exp\!\big((\alpha_i + g_i)/\tau\big)}{\sum_{j} \exp\!\big((\alpha_j + g_j)/\tau\big)}

This produces a soft vector — a point in the probability simplex (the set of non-negative vectors summing to 11), rather than a hard one-hot. Crucially, the softmax is smooth, so gradients flow back to the logits α\alpha.

The temperature τ\tau controls how peaky the output is, and sets a trade-off:

In practice τ\tau is often annealed (started high and decreased during training) to balance faithful samples against usable gradients.


The Score-Function Estimator (REINFORCE)

Motivation

The recurring problem is optimising an objective that is an expectation over a distribution that itself depends on the parameters:

J(θ)=Ex∼pθ(x) ⁣[f(x)]J(\theta) = \mathbb{E}_{x \sim p_\theta(x)}\!\big[f(x)\big]

where f(x)f(x) scores a sample xx, and xx is drawn through the very parameters θ\theta we want to optimise. If we naively move the gradient inside the integral:

∇θJ(θ)=∫∇θpθ(x) f(x) dx\nabla_\theta J(\theta) = \int \nabla_\theta p_\theta(x)\, f(x)\, dx

we hit a snag: ∇θpθ(x)\nabla_\theta p_\theta(x) is not a probability distribution — it doesn’t integrate to 11 and can be negative. So this integral is not an expectation under pθp_\theta, and we cannot directly approximate it by Monte-Carlo sampling (Monte Carlo requires a genuine expectation). Both REINFORCE and the reparameterisation trick are ways of rewriting this quantity as a proper expectation we can sample.

The log-derivative (score-function) trick

Start from the derivative of a log:

∇θlog⁡pθ(x)=1pθ(x)∇θpθ(x)\nabla_\theta \log p_\theta(x) = \frac{1}{p_\theta(x)} \nabla_\theta p_\theta(x)

Rearranging gives an identity that reintroduces the density pθ(x)p_\theta(x) as a multiplicative factor:

∇θpθ(x)=pθ(x) ∇θlog⁡pθ(x)\nabla_\theta p_\theta(x) = p_\theta(x)\, \nabla_\theta \log p_\theta(x)

The term ∇θlog⁡pθ(x)\nabla_\theta \log p_\theta(x) is called the score.

Now the derivation. Write the objective as an integral, push in the gradient, and use that ff does not depend on θ\theta (so ∇θf=0\nabla_\theta f = 0, eliminating one product-rule term):

∇θEx∼pθ(x) ⁣[f(x)]=∇θ∫f(x) pθ(x) dx=∫f(x) ∇θpθ(x) dx\nabla_\theta \mathbb{E}_{x \sim p_\theta(x)}\!\big[f(x)\big] = \nabla_\theta \int f(x)\, p_\theta(x)\, dx = \int f(x)\, \nabla_\theta p_\theta(x)\, dx

Substitute the log-derivative identity:

=∫f(x) pθ(x) ∇θlog⁡pθ(x) dx=Ex∼pθ(x) ⁣[f(x) ∇θlog⁡pθ(x)]= \int f(x)\, p_\theta(x)\, \nabla_\theta \log p_\theta(x)\, dx = \mathbb{E}_{x \sim p_\theta(x)}\!\big[f(x)\, \nabla_\theta \log p_\theta(x)\big]

This is a proper expectation under pθp_\theta, so we can estimate it by Monte-Carlo sampling:

∇θJ(θ)≈1N∑i=1Nf(xi) ∇θlog⁡pθ(xi)\nabla_\theta J(\theta) \approx \frac{1}{N} \sum_{i=1}^{N} f(x_i)\, \nabla_\theta \log p_\theta(x_i)

Notice the estimator only ever evaluates ff — it never differentiates it. So ff and xx need not be differentiable in θ\theta, which is exactly why this estimator works for discrete distributions where reparameterisation does not apply.

REINFORCE vs the Reparameterisation Trick

Same objective J(θ)=Ex∼pθ[f(x)]J(\theta) = \mathbb{E}_{x \sim p_\theta}[f(x)], two different rewrites.

The reparameterisation trick applies a transformation x=gθ(ϵ)=μ+σϵx = g_\theta(\epsilon) = \mu + \sigma\epsilon with ϵ∼p(ϵ)\epsilon \sim p(\epsilon), turning the objective into an expectation over the parameter-free noise ϵ\epsilon:

J(θ)=Eϵ∼p(ϵ) ⁣[f(gθ(ϵ))]J(\theta) = \mathbb{E}_{\epsilon \sim p(\epsilon)}\!\big[f(g_\theta(\epsilon))\big]

Now the gradient can move inside and be computed by the chain rule:

∇θJ=Eϵ ⁣[∇θf(gθ(ϵ))]=Eϵ ⁣[f′(gθ(ϵ)) ∇θgθ(ϵ)]\nabla_\theta J = \mathbb{E}_{\epsilon}\!\big[\nabla_\theta f(g_\theta(\epsilon))\big] = \mathbb{E}_{\epsilon}\!\big[f'(g_\theta(\epsilon))\, \nabla_\theta g_\theta(\epsilon)\big]

The contrast:

This is precisely why continuous Gaussian latents (the VAE) use reparameterisation, while discrete latents (the VQ-VAE) must fall back on the straight-through estimator, Gumbel-Softmax, or REINFORCE.

[TODO] Variance reduction for REINFORCE — baselines and control variates. The estimator’s high variance is its main practical drawback, and baselines are the standard remedy.

Learning Path

What should I learn before Autoencoders?

Start with gradient-descent before studying Autoencoders.

How does Autoencoders connect to other AI topics?

Autoencoders belongs to Deep Learning and connects to unsupervised-learning, generative-models, representation-learning, autoencoders, variational-inference, vae, vq-vae, reparameterisation, gumbel-softmax, reinforce, elbo, kl-divergence.

Further Reading

  1. Auto-Encoding Variational Bayes (Kingma & Welling, 2013)
  2. Neural Discrete Representation Learning (van den Oord et al., 2017)
  3. Categorical Reparameterization with Gumbel-Softmax (Jang et al., 2016)