What it actually does
Force the data through a narrow layer and demand that it come out intact. Because the bottleneck cannot carry everything, the network has to discover what is worth keeping — and whatever survives that squeeze is a representation. The corollary is the useful part: an input the network cannot rebuild is an input unlike anything it learned to compress, which is exactly what an anomaly is.
- Training data is representative of normal, and for anomaly detection genuinely free of the anomalies to be detected
- The data lies near a lower-dimensional manifold, so a bottleneck loses structure rather than signal
- Reconstruction error is a meaningful notion of unusualness for this data — false when the interesting deviation is in a feature the loss barely weights
- Features are on comparable scales, since squared error weights them by variance whether you intended that or not
- Normal behaviour is stable enough between retrains that a threshold calibrated once remains meaningful
Reconstruction objective
The loss says only: rebuild what you were given. What makes that non-trivial is the bottleneck, and the entire art of the model is choosing how to constrain it. Two facts are worth stating together. First, with linear layers and squared error, this provably recovers the principal subspace — an autoencoder is PCA with the linearity removed, and knowing that sets the baseline any nonlinear version must beat. Second, with an unconstrained bottleneck as wide as the input, the optimum is the identity function and the model learns nothing while the loss goes to zero. The regularizer exists precisely because the objective on its own does not prevent that.
- the encoder, mapping input to a latent code of lower dimension
- the decoder, mapping the code back to input space
- an optional constraint on the code — sparsity, contraction, or noise applied to the input instead
- how hard that constraint is enforced; zero gives the plain undercomplete form
- the bottleneck width against the input width — the constraint doing all the work
How it is minimized
Adam on reconstruction error, with the bottleneck or an explicit regularizer supplying the constraint
- the latent code — the thing you actually keep when the model is used for representation
- decoder weights, sometimes tied to the encoder transpose, which halves parameters and regularizes
- the reconstruction residual: the gradient signal AND the anomaly score, which is the model’s central economy
The gradient is the residual times the code, which is worth pausing on because it is the same outer-product form as any linear regression — the network is doing least squares against its own input, layer by layer. Everything interesting is in the constraint rather than the optimizer. An undercomplete bottleneck constrains by width. A denoising variant corrupts the input and asks for the clean version, which forces the model to learn the data manifold rather than a copy operation and is usually the better choice. A sparse variant penalizes mean activation so the code is overcomplete but mostly zero. A contractive one penalizes the Jacobian so the code is insensitive to small input changes. They are different answers to the same question — what stops this collapsing to the identity — and picking one is the actual modelling decision.
bottleneck width
The constraint doing all the work. Too wide and it learns the identity; too narrow and it discards signal (2 to 128, or 1–10% of input width)
depth
Encoder and decoder layers. Deeper buys nonlinearity and costs the clean PCA comparison (1 to 4 each side)
noise level (denoising)
Corruption applied to the input. The most reliable regularizer here, and usually worth using (0.1 to 0.5)
sparsity target
Desired mean activation for a sparse code; enforced by a KL penalty per unit (0.01 to 0.1)
tied weights
Decoder as encoder transpose — halves parameters and regularizes, at some capacity cost (tied / untied)
output activation
Must match the data range. Sigmoid for [0,1], linear for unbounded; a mismatch caps reconstruction permanently (linear / sigmoid)
Easy to train and easy to fool, which is the combination that makes it deceptive. The loss almost always goes down, so it gives no signal about whether the representation is any good — the only honest evaluation is downstream. Three characteristic failures. The identity collapse, where a bottleneck as wide as the input with enough capacity learns to copy: the loss is excellent and the model is useless, and it is invisible without checking the effective rank of the code. Dead latent units, where some dimensions are constantly zero and the effective bottleneck is narrower than declared. And the anomaly-detection failure that catches people repeatedly: train on data that contains the anomalies and the model learns to reconstruct them too, so the detector is silently blind to exactly what it was built for.
O(d · k) per layer per sample for a dense autoencoder, forward and backward. Cheap in absolute terms, fully parallel across the batch, and with no sequential structure at all — the cost is entirely set by input width and depth.
How, where, and why it is used
Time-Series Forecasting
Not applicableThe objective is to reproduce its input, not to predict a future one; there is no mechanism that maps a history to an unseen horizon. It appears in forecasting pipelines as a compressor for a high-dimensional panel before a forecaster runs, but that is preprocessing rather than forecasting, and the model doing the forecast is a different entry.
Anomaly Detection
Primary fitTrain only on normal data, then score each input by how badly it reconstructs. The economy here is unusually clean: the training signal and the anomaly score are literally the same quantity, so there is no separate scoring model to fit, calibrate or maintain. For structured inputs the per-feature residual also localizes the anomaly, which a single scalar score cannot.
- Industrial and equipment monitoring on multivariate sensor data, where faults are rare and unlabelled
- Network intrusion and telemetry monitoring, where the normal profile is abundant and attacks are not
- Manufacturing defect detection on images, where the defect classes are unknown in advance
- Payment and transaction anomaly detection as an unsupervised complement to a supervised classifier
It is the right tool when anomalies are rare, unlabelled and heterogeneous — precisely the case where a supervised classifier cannot be trained and where you would rather model normality than enumerate failure modes. Because the model is nonlinear it also catches anomalies that are unremarkable in every individual feature and only unusual in combination, which is where it beats per-feature thresholding decisively. Against it: the training set must be genuinely clean, and it usually is not; a sufficiently flexible model reconstructs anomalies too, so capacity is a sensitivity knob whether you intended that or not; and on ordinary tabular data an isolation forest is competitive at a fraction of the cost and effort. Reach for this when the input is high-dimensional and structured, not by default.
- Scale features so no single one dominates a squared-error loss, since the residual is the score and an unscaled feature sets the threshold by itself
- Train on a period verified to be normal — the most common failure here is training on data containing the anomalies you want to catch
- Use the denoising variant, which prevents the identity shortcut and generally produces a sharper score gap
- Calibrate the threshold on held-out normal data by percentile rather than picking a round number
- Keep the per-feature residual, not just the norm, so an alert can say which dimension was wrong
Precision and recall at a fixed alert budget on whatever confirmed cases exist, with AUROC as a supporting number rather than the headline — at realistic anomaly rates a high AUROC is compatible with an unusable alert stream. Always compare against an isolation forest and against per-feature thresholding on the same data: those baselines win more often than expected, and skipping them is how an autoencoder pipeline gets deployed without justification.
- Anomalies present in the training set, which teaches the model to reconstruct them and blinds the detector silently
- Too much capacity, so everything reconstructs well and the score gap closes
- An unscaled feature dominating the residual and effectively becoming the entire detector
- Drift in normal behaviour raising every score at once, which reads as a flood of anomalies rather than as staleness
Optimization
AdaptedThe clean result is the PCA equivalence: a linear autoencoder under squared error recovers the principal subspace, so the model is a gradient-descent solution to a problem that has a closed form. That makes it an unusually good case study in when to iterate and when to solve — the closed form is exact and cheap at moderate width, and the iterative version only earns its place once nonlinearity or a data size that will not fit in memory enters the picture. Beyond that, the constraint variants are each a different regularization scheme attacking the same degenerate optimum.
- Linear autoencoder versus SVD as the canonical iterative-versus-closed-form comparison
- Sparsity and contractive penalties as explicit regularizers on a representation rather than on weights
- Denoising as regularization by data corruption, which is a different mechanism from a penalty term
- Nonlinear dimensionality reduction as preprocessing that makes a downstream optimization tractable
Worth studying because the degenerate optimum is so easy to state: with an unconstrained wide bottleneck the identity function is a global minimum with zero loss, and the model learns nothing while reporting perfect success. That is the clearest example in this reference of an objective that is satisfied without solving the problem, and it is why every useful variant is really a story about the constraint rather than about the loss. The PCA equivalence is the other transferable lesson — always know what the closed-form solution to your problem is before reaching for gradient descent, because sometimes it is the same answer for far less effort.
- Compare against PCA at the same latent width first; if the nonlinear model does not beat it, the extra machinery is not earning anything
- Check the effective rank of the learned code, since dead units make the real bottleneck narrower than the declared one
- Prefer denoising to an explicit penalty where possible — corrupting the input is a stronger constraint than penalizing the code
- Whiten or scale inputs before comparing to PCA, or the comparison measures the preprocessing rather than the model
Reconstruction error against PCA at matched latent width, plus downstream task performance on the codes. The second matters more: reconstruction error measures compression, and a representation that compresses better is not automatically one that transfers better.
- Reaching for a nonlinear model when PCA at the same width does as well for a fraction of the cost
- Reading a falling reconstruction loss as evidence the representation is good, which it is not
- Ignoring dead latent units, so the effective and declared bottleneck widths differ silently
Computer Vision
ViableA convolutional encoder-decoder trained to reconstruct images. Used for denoising, for defect detection where the residual map localizes the flaw spatially, and historically as unsupervised pretraining before labelled data was plentiful.
- Image denoising and restoration, where the corrupted-to-clean mapping is the actual task
- Industrial visual inspection, where the pixel-level residual gives a defect map for free
- Compression and latent representation for a downstream generative model, which is how modern latent diffusion is structured
- Unsupervised pretraining on unlabelled imagery, largely superseded by masked-image modelling
On images it earns its place when the residual map itself is the deliverable — a defect localization is more useful than a defect score, and this produces one without any localization supervision. The role in latent diffusion is the other current justification and is a large one: compressing to a latent space is what makes high-resolution diffusion affordable. As a representation learner for classification it has been clearly superseded by contrastive and masked-modelling objectives, which produce features that transfer better, and it is worth saying that plainly rather than presenting it as a live option.
- Use convolutional layers rather than dense ones; a dense autoencoder on pixels discards spatial structure entirely
- Match the output activation to the pixel range, since a mismatch caps reconstruction quality permanently
- Prefer perceptual or adversarial losses over pixel MSE where sharpness matters — squared error produces characteristically blurry output
- Keep the residual as a spatial map rather than a scalar, since localization is most of the value here
PSNR or SSIM for reconstruction quality, but pixel-level localization scores for inspection tasks — those are the numbers that reflect what the model is actually for. Compare against a pretrained-feature detector, which on modern defect benchmarks usually wins.
- Blurry reconstructions from pixel MSE, which is inherent to the loss rather than a capacity problem
- Dense layers on image data, discarding the spatial structure a convolutional encoder would keep
- Using it for representation learning where contrastive or masked-modelling objectives are clearly stronger
Recommendation & Ranking
ViableTreat a user’s interaction vector over the whole catalogue as the input and reconstruct it; the reconstruction over unobserved items is the recommendation score. This is the denoising-autoencoder formulation of collaborative filtering, and it generalizes matrix factorization by allowing a nonlinear encoder.
- Collaborative filtering on implicit feedback, where the denoising formulation is the standard framing
- Session-based recommendation from a bag of recent interactions
- Content-aware recommendation where item features join the input alongside interactions
It generalizes matrix factorization in a defensible way — the same reconstruct-the-interaction-matrix objective with the linearity removed — and the denoising framing matches the problem honestly, since unobserved is not the same as negative and corrupting the input models exactly that ambiguity. The practical limits are real: the input width is the catalogue size, so the first layer is enormous; the model is transductive over items, so a new item needs a retrain; and it optimizes reconstruction rather than ranking, which a model trained on a ranking loss will beat where ranking is what matters.
- Treat missing as missing rather than as zero, or the model learns to predict absence
- Use the denoising variant with input dropout, which is the formulation that actually works here
- Weight observed entries above unobserved ones, since implicit feedback is positive-only and the imbalance is severe
- Normalize by user activity, or heavy users dominate the loss entirely
Recall and NDCG at k on held-out interactions, split strictly by time, with a popularity baseline alongside. Popularity is startlingly hard to beat on recall, and a model that fails to is not worth deploying.
- Zero-filling unobserved entries, which trains the model to predict that items are not interacted with
- A catalogue-width input layer dominating both parameter count and memory
- No path for new items, which is structural rather than a tuning problem
What it costs to run
Minutes to hours on a single GPU for most tabular and sensor applications — among the cheapest models in this category, with no sequential structure and a loss that needs no labels.
One forward pass per input, and for anomaly detection only the encoder-decoder pair is needed. Sub-millisecond on CPU for modest widths, and genuinely streaming — which is why it is practical as an always-on detector rather than a batch one.
Driven entirely by how fast normal behaviour drifts: weekly to monthly for live systems, rarely for stable physical processes. In anomaly detection the retraining question is uncomfortable, because a retrain on recent data absorbs any slow-developing anomaly into the new definition of normal.
- Track the reconstruction-error distribution on known-normal traffic rather than the alert count, since a rising floor means drift and not a genuine increase in anomalies
- Monitor the effective rank of the latent code — collapsing dimensions mean the real bottleneck has narrowed since training
- Watch per-feature residual contributions; one feature quietly dominating means the scaling has drifted
- Alert on the fraction of inputs above threshold, which is the signal that separates model staleness from a real event
- The feature scaler is part of the model. Refitting it at inference changes every score, and because the score is a residual the change is silent and total
- Only the encoder is needed to produce representations, and only the residual for detection — shipping the whole model when half is needed is a common and avoidable serving cost
- Thresholds are not transferable across retrains: the error distribution shifts, so every threshold must be recalibrated with the model
- Anomalies in the training data are the standard failure, and the only defence is a verified-clean training window rather than anything in the model
- The output activation must match the data range; a sigmoid head on unbounded data caps reconstruction and the loss plateaus with no other symptom
Pros and cons, in context
One model, two jobs
The reconstruction residual is simultaneously the training signal and the anomaly score, so there is no separate detector to fit, calibrate or version. An operational simplification that matters more in practice than it sounds.
Nonlinear where PCA is linear
Captures curved manifolds a linear projection cannot, and catches anomalies that are unremarkable in every single feature and only unusual in combination. That combination case is where it beats per-feature thresholding decisively.
Needs no labels at all
Decisive where anomalies are rare, unlabelled and heterogeneous — which is the normal situation in monitoring, and the case where no supervised classifier can be trained.
Per-feature residuals localize the problem
An alert can say which dimension was wrong, or on images which region. That is far more actionable than a scalar score, and it comes for free without any localization supervision.
Cheap and genuinely streaming
One forward pass per input, no sequential structure, sub-millisecond on CPU. Practical as an always-on detector rather than a batch job, which many alternatives are not.
The identity collapse
With a wide enough bottleneck and enough capacity, copying the input is a global optimum with zero loss. The model reports perfect success and has learned nothing, and it is invisible without checking the effective rank of the code.
Training data must be genuinely clean
Anomalies in the training set teach the model to reconstruct them, so the detector is blind to exactly what it was built for. The most common failure in practice, and it has no in-model defence.
Falling loss says nothing about representation quality
Reconstruction error measures compression, not usefulness. The only honest evaluation is downstream, which makes tuning slow and makes it easy to ship something that looked fine.
Capacity is a sensitivity knob in disguise
More capacity reconstructs anomalies too and closes the score gap. Tuning for reconstruction quality and tuning for detection sensitivity pull in opposite directions, which is genuinely awkward to manage.
Superseded for representation learning
Contrastive and masked-modelling objectives produce features that transfer substantially better. This remains excellent for anomaly detection and compression, and is no longer the right choice for learning transferable representations.
Intuitive, idiomatic, optimized
The same algorithm in three languages, each progressing through Kent Beck's three stages. Every optimization names the lever it pulls and the trade-off it accepts — cost profiles are illustrative, not measured benchmarks.
The math, transcribed literally. No dependencies.
"""An autoencoder, transcribed the way the definition reads.
h = sigma(W x + b) encode: squeeze into k dimensions
x_hat = sigma(W' h + b') decode: rebuild d dimensions
loss = ||x - x_hat||^2 the whole objective
Nothing about the loss makes this non-trivial. What does is the BOTTLENECK:
with k as wide as d and enough capacity, the identity function is a global
minimum with zero loss, and the model learns nothing while reporting perfect
success. Every useful variant is a story about what stops that.
Plain loops, no libraries, manual backward pass.
"""
import math
import random
SEED = 11
def sigmoid(value):
if value >= 0.0:
return 1.0 / (1.0 + math.exp(-value))
positive = math.exp(value)
return positive / (1.0 + positive)
def matvec(weight, vector, bias):
"""weight is out_dim rows of in_dim."""
return [
bias[i] + sum(w * v for w, v in zip(row, vector))
for i, row in enumerate(weight)
]
def encode(sample, encoder_weight, encoder_bias):
return [sigmoid(value) for value in matvec(encoder_weight, sample, encoder_bias)]
def decode(code, decoder_weight, decoder_bias):
"""Linear output head: the data here is unbounded, so a sigmoid would cap
reconstruction permanently and the loss would plateau with no other
symptom. The output activation must match the data range."""
return matvec(decoder_weight, code, decoder_bias)
def reconstruction_error(sample, rebuilt):
"""This single quantity is BOTH the training signal and, at inference,
the anomaly score. That economy is the model's central feature."""
return sum((a - b) ** 2 for a, b in zip(sample, rebuilt))
def train_step(sample, params, learning_rate, noise_level, rng):
"""One sample, forward and backward, written out.
The denoising variant: corrupt the INPUT but ask for the CLEAN target.
That is a stronger constraint than any penalty term, because copying the
input is no longer even available as a solution - the network must learn
the data manifold to fill in what was destroyed.
"""
corrupted = [
0.0 if rng.random() < noise_level else value
for value in sample
]
code = encode(corrupted, params['encoder_weight'], params['encoder_bias'])
rebuilt = decode(code, params['decoder_weight'], params['decoder_bias'])
# dL/dx_hat for squared error.
output_delta = [2.0 * (r - s) for r, s in zip(rebuilt, sample)]
# Decoder gradient: residual times code. The same outer-product form as
# any least-squares fit - the network is doing regression against its own
# input, layer by layer.
for i in range(len(params['decoder_weight'])):
for j in range(len(code)):
params['decoder_weight'][i][j] -= learning_rate * output_delta[i] * code[j]
params['decoder_bias'][i] -= learning_rate * output_delta[i]
# Backpropagate into the code, then through the sigmoid.
code_delta = [
sum(
output_delta[i] * params['decoder_weight'][i][j]
for i in range(len(output_delta))
)
for j in range(len(code))
]
hidden_delta = [d * c * (1.0 - c) for d, c in zip(code_delta, code)]
for j in range(len(params['encoder_weight'])):
for k in range(len(corrupted)):
params['encoder_weight'][j][k] -= learning_rate * hidden_delta[j] * corrupted[k]
params['encoder_bias'][j] -= learning_rate * hidden_delta[j]
return reconstruction_error(sample, rebuilt)
def initialise(input_dim, latent_dim, rng):
limit = math.sqrt(6.0 / (input_dim + latent_dim))
return {
'encoder_weight': [
[rng.uniform(-limit, limit) for _ in range(input_dim)]
for _ in range(latent_dim)
],
'encoder_bias': [0.0] * latent_dim,
'decoder_weight': [
[rng.uniform(-limit, limit) for _ in range(latent_dim)]
for _ in range(input_dim)
],
'decoder_bias': [0.0] * input_dim,
}
def train(data, latent_dim, epochs=50, learning_rate=0.01, noise_level=0.2):
"""Note there are no labels anywhere. The target IS the input."""
rng = random.Random(SEED)
params = initialise(len(data[0]), latent_dim, rng)
for _ in range(epochs):
order = list(range(len(data)))
rng.shuffle(order)
for index in order:
train_step(data[index], params, learning_rate, noise_level, rng)
return params
def anomaly_scores(data, params):
"""No separate detector is fitted. The residual IS the score."""
scores = []
for sample in data:
code = encode(sample, params['encoder_weight'], params['encoder_bias'])
rebuilt = decode(code, params['decoder_weight'], params['decoder_bias'])
scores.append(reconstruction_error(sample, rebuilt))
return scores
def per_feature_residual(sample, params):
"""Which dimension was wrong, not just how wrong overall.
A scalar score says something is unusual; this says what. That is most of
the operational value, and it comes free with no localization labels.
"""
code = encode(sample, params['encoder_weight'], params['encoder_bias'])
rebuilt = decode(code, params['decoder_weight'], params['decoder_bias'])
return [(a - b) ** 2 for a, b in zip(sample, rebuilt)]
def effective_rank(codes, tolerance=1e-6):
"""How many latent dimensions are actually being used.
A dead unit makes the real bottleneck narrower than the declared one, and
a code whose rank equals the input width means the model has found the
identity shortcut. Neither is visible in the loss.
"""
latent_dim = len(codes[0])
means = [sum(code[j] for code in codes) / len(codes) for j in range(latent_dim)]
variances = [
sum((code[j] - means[j]) ** 2 for code in codes) / len(codes)
for j in range(latent_dim)
]
return sum(1 for variance in variances if variance > tolerance)
O(d * k) per sample per layer, forward and backward. Illustrative, not a measured benchmark: one sample at a time through nested Python loops means a 200-feature input with a 16-dimensional code is roughly six thousand interpreter iterations per training step, before any batching.
How this language actually wants it written.
"""The same model, with the constraint made explicit and the collapse checked.
Two things change. The variant - what stops this collapsing to the identity -
becomes an explicit type rather than a scattering of flags, because choosing
it IS the modelling decision. And the two failures that the loss cannot see -
the identity collapse and a training set containing anomalies - become
things the code checks rather than things you discover later.
"""
from __future__ import annotations
import math
import random
from dataclasses import dataclass
from enum import Enum
from typing import NamedTuple, Sequence
class Constraint(Enum):
"""What prevents the degenerate solution.
Not a style choice: with an unconstrained wide bottleneck the identity
function is a global minimum with zero loss, so one of these must be in
force or the model learns nothing while reporting perfect success.
"""
UNDERCOMPLETE = 'undercomplete' # the bottleneck alone
DENOISING = 'denoising' # corrupt the input, ask for the clean target
SPARSE = 'sparse' # penalize mean activation
CONTRACTIVE = 'contractive' # penalize the encoder Jacobian
class DegenerateModel(ValueError):
"""Raised when the model can trivially learn the identity.
Its own type because it is the failure with no other symptom: the loss
goes to zero and the representation is worthless.
"""
class ShapeMismatch(ValueError):
"""Raised on a shape violation instead of broadcasting past it."""
class ContaminatedTrainingSet(ValueError):
"""Raised when the training data looks like it contains anomalies.
Its own type because the consequence is specific: the detector learns to
reconstruct the anomalies and is blind to exactly what it was built for,
and nothing downstream will indicate that.
"""
@dataclass(frozen=True)
class AutoencoderConfig:
"""Frozen so the constraint cannot drift after construction."""
input_dim: int
latent_dim: int
constraint: Constraint = Constraint.DENOISING
noise_level: float = 0.2
sparsity_target: float = 0.05
sparsity_weight: float = 0.1
tied_weights: bool = False
def __post_init__(self) -> None:
# Guard clause for the failure the loss cannot see.
if self.constraint is Constraint.UNDERCOMPLETE and self.latent_dim >= self.input_dim:
raise DegenerateModel(
f'a {self.latent_dim}-wide code for {self.input_dim} inputs is not '
'undercomplete; the identity is a zero-loss optimum. Narrow the '
'bottleneck or choose another constraint.'
)
if self.constraint is Constraint.DENOISING and not 0.0 < self.noise_level < 1.0:
raise ShapeMismatch('denoising needs a corruption rate strictly inside (0, 1)')
if self.latent_dim < 1:
raise ShapeMismatch('the code must have at least one dimension')
@property
def compression_ratio(self) -> float:
"""Worth naming: this is the constraint, expressed as a number."""
return self.latent_dim / self.input_dim
class Reconstruction(NamedTuple):
"""Code and rebuilt input together.
Both come back because callers need different halves: representation work
wants the code, detection wants the residual, and returning only one
forces a second forward pass.
"""
code: list[float]
rebuilt: list[float]
def squared_error(self, original: Sequence[float]) -> float:
return math.fsum((a - b) ** 2 for a, b in zip(original, self.rebuilt))
def per_feature_residual(self, original: Sequence[float]) -> list[float]:
"""Which dimension was wrong, not just how wrong overall.
A scalar says something is unusual; this says what. That is most of
the operational value, and it costs nothing extra.
"""
return [(a - b) ** 2 for a, b in zip(original, self.rebuilt)]
def sigmoid(value: float) -> float:
"""Branch on sign so exp cannot overflow on a large negative input."""
if value >= 0.0:
return 1.0 / (1.0 + math.exp(-value))
positive = math.exp(value)
return positive / (1.0 + positive)
class Autoencoder:
"""Encoder and decoder with the constraint applied in the forward pass."""
def __init__(self, config: AutoencoderConfig, seed: int = 11) -> None:
self._config = config
rng = random.Random(seed)
limit = math.sqrt(6.0 / (config.input_dim + config.latent_dim))
self._encoder = [
[rng.uniform(-limit, limit) for _ in range(config.input_dim)]
for _ in range(config.latent_dim)
]
self._encoder_bias = [0.0] * config.latent_dim
# Tied weights halve the parameter count and regularize; the decoder
# is then the encoder transpose and is not stored separately.
self._decoder = (
None
if config.tied_weights
else [
[rng.uniform(-limit, limit) for _ in range(config.latent_dim)]
for _ in range(config.input_dim)
]
)
self._decoder_bias = [0.0] * config.input_dim
def encode(self, sample: Sequence[float]) -> list[float]:
if len(sample) != self._config.input_dim:
raise ShapeMismatch(
f'expected {self._config.input_dim} features, got {len(sample)}'
)
return [
sigmoid(bias + math.fsum(w * v for w, v in zip(row, sample)))
for row, bias in zip(self._encoder, self._encoder_bias)
]
def decode(self, code: Sequence[float]) -> list[float]:
"""Linear output head, because the data range here is unbounded.
Matching the output activation to the data range is not cosmetic: a
sigmoid head on unbounded data caps reconstruction permanently and
the loss plateaus with no other symptom.
"""
if self._decoder is None:
return [
self._decoder_bias[i]
+ math.fsum(self._encoder[j][i] * code[j] for j in range(len(code)))
for i in range(self._config.input_dim)
]
return [
bias + math.fsum(w * c for w, c in zip(row, code))
for row, bias in zip(self._decoder, self._decoder_bias)
]
def forward(self, sample: Sequence[float]) -> Reconstruction:
code = self.encode(sample)
return Reconstruction(code=code, rebuilt=self.decode(code))
def corrupt(self, sample: Sequence[float], rng: random.Random) -> list[float]:
"""Denoising: corrupt the INPUT and ask for the CLEAN target.
A stronger constraint than any penalty term, because copying the
input is no longer available as a solution at all - the network must
learn the data manifold to fill in what was destroyed.
"""
if self._config.constraint is not Constraint.DENOISING:
return list(sample)
return [
0.0 if rng.random() < self._config.noise_level else value
for value in sample
]
def sparsity_penalty(self, codes: Sequence[Sequence[float]]) -> float:
"""KL between the mean activation and the target, per unit.
Applies to an OVERCOMPLETE code: the bottleneck is wide but mostly
zero, so the constraint is on how much is used rather than on how
much is available.
"""
if self._config.constraint is not Constraint.SPARSE:
return 0.0
target = self._config.sparsity_target
penalty = 0.0
for unit in range(self._config.latent_dim):
mean = math.fsum(code[unit] for code in codes) / len(codes)
mean = min(max(mean, 1e-6), 1.0 - 1e-6)
penalty += target * math.log(target / mean) + (1 - target) * math.log(
(1 - target) / (1 - mean)
)
return self._config.sparsity_weight * penalty
def effective_rank(codes: Sequence[Sequence[float]], tolerance: float = 1e-6) -> int:
"""How many latent dimensions are actually being used.
A dead unit makes the real bottleneck narrower than the declared one, and
a rank equal to the input width means the model found the identity
shortcut. Neither is visible in the loss, which is why this is a function
rather than a comment.
"""
if not codes:
raise ShapeMismatch('cannot measure rank with no codes')
latent_dim = len(codes[0])
used = 0
for unit in range(latent_dim):
values = [code[unit] for code in codes]
mean = math.fsum(values) / len(values)
variance = math.fsum((value - mean) ** 2 for value in values) / len(values)
if variance > tolerance:
used += 1
return used
def assert_training_set_clean(
errors: Sequence[float],
tail_fraction: float = 0.01,
ratio_threshold: float = 20.0,
) -> None:
"""Heuristic check that the training set is anomaly-free.
A clean training set produces a reconstruction-error distribution with a
short tail. A long one means the model is being asked to reconstruct
things that do not belong - and if it succeeds, the detector is blind to
exactly what it was built for.
"""
if len(errors) < 100:
raise ShapeMismatch('too few samples to judge the error distribution')
ordered = sorted(errors)
median = ordered[len(ordered) // 2]
tail = ordered[int(len(ordered) * (1.0 - tail_fraction))]
if median > 0.0 and tail / median > ratio_threshold:
raise ContaminatedTrainingSet(
f'the top {tail_fraction:.0%} of reconstruction errors are {tail / median:.0f}x '
'the median; the training set probably contains anomalies'
)
def calibrate_threshold(errors: Sequence[float], false_positive_rate: float) -> float:
"""Threshold as a percentile of held-out normal errors.
Not a round number: the error scale is arbitrary and shifts with every
retrain, so a threshold must be recalibrated alongside the model rather
than carried across.
"""
if not 0.0 < false_positive_rate < 1.0:
raise ShapeMismatch('the false-positive rate must lie inside (0, 1)')
ordered = sorted(errors)
index = min(int((1.0 - false_positive_rate) * len(ordered)), len(ordered) - 1)
return ordered[index]
What changed: Two changes, both aimed at failures the loss cannot see. The constraint — what stops the model collapsing to the identity — becomes an explicit enum rather than a scattering of flags, because choosing it is the actual modelling decision, and the config refuses an undercomplete configuration whose code is as wide as its input: that combination has the identity as a zero-loss global optimum, so the model reports perfect success and has learned nothing. The second is contamination. A training set containing the anomalies it should detect teaches the model to reconstruct them, and the detector is then blind to exactly what it exists for with no downstream signal at all; the heuristic check on the error distribution turns that into something the pipeline can refuse rather than something discovered months later. Alongside those, effective rank becomes a function rather than a comment because dead units make the real bottleneck narrower than the declared one, the threshold is calibrated as a percentile because the error scale is arbitrary and shifts with every retrain, and the reconstruction returns both the code and the rebuilt input since representation work and detection want different halves.
Identical asymptotics; the change is in what fails loudly. Illustrative, not a measured benchmark: fsum and the validation add a small constant and buy two failure modes that would otherwise be invisible.
Faster, and the justification for each change.
"""Batched, fused, and with the PCA baseline computed alongside.
The structural observation: an autoencoder has no sequential structure at
all, so every layer is a dense matrix product and the whole thing batches
trivially. That makes the interesting optimization not the arithmetic but
what you compare against - a linear autoencoder under squared error provably
recovers the principal subspace, so the closed-form SVD answer is available
for free and should be computed before any gradient descent is run.
Three changes:
1. Encode, decode and the backward pass all become GEMMs over the batch.
2. The residual is computed once and reused: it is the loss, the gradient
and the anomaly score, so materializing it three times is pure waste.
3. Tied weights make the decoder a transpose rather than a second matrix,
which halves both parameters and bandwidth.
"""
from __future__ import annotations
import numpy as np
from numpy.typing import NDArray
FLOAT = np.float32
def sigmoid_(values: NDArray[np.float32]) -> NDArray[np.float32]:
"""In-place logistic. The hidden activation is the largest intermediate
in a wide encoder, so avoiding one temporary here is the single biggest
memory saving in the forward pass."""
np.negative(values, out=values)
np.exp(values, out=values)
values += 1.0
np.reciprocal(values, out=values)
return values
class BatchedAutoencoder:
"""Tied-weight autoencoder over a whole batch.
Tied weights mean the decoder IS the encoder transpose: half the
parameters, half the bandwidth, and a regularizer for free. The cost is
capacity, which on most tabular data is not the binding constraint.
"""
def __init__(self, input_dim: int, latent_dim: int, seed: int = 11) -> None:
rng = np.random.default_rng(seed)
limit = np.sqrt(6.0 / (input_dim + latent_dim))
self._weight = rng.uniform(-limit, limit, (input_dim, latent_dim)).astype(FLOAT)
self._encoder_bias = np.zeros(latent_dim, dtype=FLOAT)
self._decoder_bias = np.zeros(input_dim, dtype=FLOAT)
self._latent_dim = latent_dim
def forward(
self, batch: NDArray[np.float32]
) -> tuple[NDArray[np.float32], NDArray[np.float32]]:
"""batch is (n, input_dim). Returns (code, rebuilt)."""
code = batch @ self._weight
code += self._encoder_bias
sigmoid_(code)
# Linear output head: the data range is unbounded, and a sigmoid here
# would cap reconstruction permanently with the loss plateauing and
# no other symptom.
rebuilt = code @ self._weight.T
rebuilt += self._decoder_bias
return code, rebuilt
def train_step(
self,
batch: NDArray[np.float32],
noise_level: float,
learning_rate: float,
rng: np.random.Generator,
) -> np.float32:
"""One batched step. The residual is computed ONCE and reused.
It is simultaneously the loss, the output gradient and the anomaly
score, so materializing it three separate times is pure waste - and
that economy is the model's central feature, not an implementation
detail.
"""
# Denoising: corrupt the INPUT, ask for the CLEAN target. A stronger
# constraint than any penalty, because copying is no longer even
# available as a solution.
mask = (rng.random(batch.shape) >= noise_level).astype(FLOAT)
corrupted = batch * mask
code, rebuilt = self.forward(corrupted)
residual = rebuilt
residual -= batch # in place: `rebuilt` is not needed again
loss = np.einsum('ij,ij->', residual, residual, optimize=True) / batch.shape[0]
scale = FLOAT(2.0 * learning_rate / batch.shape[0])
# Tied weights mean one matrix gets BOTH gradient terms: the decoder
# path and the encoder path. Two GEMMs, one update.
decoder_gradient = code.T @ residual
code_delta = residual @ self._weight
code_delta *= code
code_delta *= 1.0 - code
encoder_gradient = corrupted.T @ code_delta
self._weight -= scale * (decoder_gradient.T + encoder_gradient)
self._decoder_bias -= scale * residual.sum(axis=0)
self._encoder_bias -= scale * code_delta.sum(axis=0)
return loss
def anomaly_scores(self, batch: NDArray[np.float32]) -> NDArray[np.float32]:
"""The residual IS the score. No separate detector is fitted."""
_, rebuilt = self.forward(batch)
rebuilt -= batch
# Row-wise squared norm without materializing the squares: einsum
# contracts the feature axis in one pass.
return np.einsum('ij,ij->i', rebuilt, rebuilt, optimize=True)
def per_feature_residual(self, batch: NDArray[np.float32]) -> NDArray[np.float32]:
"""Which dimension was wrong, not just how wrong overall."""
_, rebuilt = self.forward(batch)
rebuilt -= batch
return np.square(rebuilt, out=rebuilt)
def pca_baseline(
data: NDArray[np.float32], latent_dim: int
) -> tuple[NDArray[np.float32], np.float32]:
"""The closed-form answer, computed before any gradient descent.
A LINEAR autoencoder under squared error provably recovers the principal
subspace, so this is the exact solution to the problem a linear model is
iterating toward. Always know it: if a nonlinear autoencoder does not
beat this at the same latent width, the extra machinery is earning
nothing and the SVD is faster, exact and has no hyperparameters.
"""
centred = data - data.mean(axis=0, dtype=np.float32)
# Economy SVD: the full one computes an (n x n) left factor that is
# immediately discarded, which on a large batch is the dominant cost.
_, singular, components = np.linalg.svd(centred, full_matrices=False)
basis = components[:latent_dim]
projected = centred @ basis.T
reconstruction_error = float(np.square(singular[latent_dim:]).sum() / len(data))
return projected.astype(FLOAT), FLOAT(reconstruction_error)
def effective_rank(codes: NDArray[np.float32], tolerance: float = 1e-6) -> int:
"""How many latent dimensions are actually used.
A dead unit makes the real bottleneck narrower than the declared one, and
a rank matching the input width means the model found the identity
shortcut. Neither is visible in the loss.
"""
return int((codes.var(axis=0, dtype=np.float32) > tolerance).sum())
def calibrate_threshold(errors: NDArray[np.float32], false_positive_rate: float) -> float:
"""Threshold as a percentile of held-out normal errors.
Partition-based, not a sort: a percentile is a selection and ordering the
rest is discarded work. The error scale is arbitrary and shifts with
every retrain, so this must be recomputed alongside the model rather than
carried across.
"""
return float(np.quantile(errors, 1.0 - false_positive_rate))
What changed: The structural point is that an autoencoder has no sequential dependency at all, so every layer is a dense product and the whole model batches trivially — which makes the interesting optimization not the arithmetic but what you compare against. A linear autoencoder under squared error provably recovers the principal subspace, so the closed-form SVD answer is available at a fraction of the cost, and computing it first is the discipline this stage encodes: if the nonlinear model does not beat it at matched latent width, the extra machinery is earning nothing. The SVD itself uses the economy form, since the full factorization computes an n-by-n left factor that is immediately discarded. Within the training step the residual is computed once and reused, because it is simultaneously the loss, the output gradient and the anomaly score — materializing it three times would be waste, and that economy is the model’s defining feature rather than an implementation detail. Tied weights make the decoder a transpose, so one matrix receives both gradient terms and both parameters and bandwidth halve. Squared norms use einsum to contract the feature axis in a single pass rather than squaring into a temporary and summing it.
Replace a closed-form solve with a numerically stabler factorization
The linear case has an exact SVD solution, and the economy factorization gives both the optimal subspace and its reconstruction error without any iteration or hyperparameter.
Trade-off: The SVD needs the whole data matrix resident and is O(n d^2), so past a few hundred thousand rows it stops being the cheap option and the iterative version wins purely on memory.
Vectorize to NumPy/BLAS (@, np.dot, einsum)
Encode, decode and both gradient terms are GEMMs over the batch, and the row-wise squared norms contract the feature axis with einsum instead of squaring into a temporary.
Trade-off: The tied-weight update needs a transpose of the decoder gradient before it can be added, and on a non-square weight that transpose is a genuine copy rather than a view.
Pre-allocate output arrays and use in-place operations
The sigmoid, the residual subtraction and the per-feature squaring all write through buffers that already exist, so a training step allocates only the gradients.
Trade-off: The rebuilt input is destroyed to become the residual, so anything wanting to inspect the reconstruction itself — the natural thing to plot when diagnosing a bad detector — needs a second forward pass.
Fuse adjacent operations so intermediates are never materialized
The sigmoid derivative is computed against the activation already in memory rather than recomputed, and the loss is contracted straight out of the residual.
Trade-off: The pre-activation is never stored, so a saturated encoder — the reason a unit goes dead and the effective rank drops — cannot be diagnosed without an unfused pass.
A training step becomes four GEMMs for the whole batch, against O(n * d * k) interpreter iterations. Illustrative, not a measured benchmark: the PCA baseline on the same data is typically a single factorization and is worth running first precisely because it is so much cheaper than the training loop it might replace.