Anomaly detection

The bayes_hdc.anomaly module provides calibrated one-class anomaly detection. You fit on “normal” data and get back split-conformal p-values: under exchangeability the false-positive rate is controlled at a target \(\alpha\) with a finite-sample guarantee, no distributional assumptions and no hand-tuned threshold.

The nonconformity score is a distance to the normal region in hypervector space (1 - cosine to the centroid by default, a k-NN mean, or Hamming distance for BSC). The conformal layer turns that score into a p-value via the standard \((1 + \#\{s_\text{cal} \ge s_\text{query}\}) / (n + 1)\) construction (Laxhammar 2014; Bates, Candès, Lei & Romano 2023).

Both classes are registered JAX pytrees, so score / pvalue / predict compose under jit and vmap, and score is differentiable.

Scorer

class bayes_hdc.anomaly.HDCAnomalyScorer(centroid, reference, dimensions, distance_metric='cosine', k_neighbors=1, vsa_model_name='map')[source]

Bases: object

Distance-metric-aware HDC nonconformity scorer.

Given a learned reference set \(\mathcal{R} = \{r_1, \dots, r_n\} \subset \mathbb{R}^d\) (or \(\{0,1\}^d\) for BSC) of “normal” hypervectors, the scorer maps a query \(x \in \mathbb{R}^d\) to a non-negative scalar

\[\begin{split}s(x; \mathcal{R}) \;=\; \begin{cases} 1 - \mathrm{sim}(x, \bar r) & k_{\text{nbrs}} = 1 \\ 1 - \displaystyle \frac{1}{k}\!\!\sum_{r \in \mathrm{topk}(x, k)} \mathrm{sim}(x, r) & k_{\text{nbrs}} = k > 1 \end{cases}\end{split}\]

where \(\bar r\) is the (normalised) centroid of \(\mathcal{R}\) and \(\mathrm{sim}\) is cosine or Hamming similarity depending on the VSA model. Larger \(s\) means more anomalous. The score is bounded in \([0, 2]\) for cosine (typically in \([0, 1]\) after normalisation) and in \([0, 1]\) for Hamming.

This score is the standard HDC nonconformity measure used by Kleyko et al. (2017), Imani et al. (2019), Pandey et al. (2021), Thomas et al. (2021), and Furlong & Eliasmith (2024). The k-nearest variant is a kernel-style alternative that preserves multi-modal structure in the reference set rather than collapsing it into a single centroid.

Variables:
  • centroid (jax.jaxlib._jax.Array) – Bundled / mean reference hypervector of shape (dimensions,). For BSC it is a thresholded majority bundle; for real-valued VSAs it is the unit-normalised mean.

  • reference (jax.jaxlib._jax.Array) – Full reference set of shape (n_reference, dimensions). Stored only when k_neighbors > 1; otherwise filled with zeros so the pytree shape is static.

  • dimensions (int) – Hypervector dimension (static).

  • distance_metric (str) – Either "cosine" or "hamming" (static).

  • k_neighbors (int) – Number of nearest reference vectors to average over. k_neighbors == 1 uses the centroid alone (static).

  • vsa_model_name (str) – Name of the underlying VSA model. Used only to select the default distance metric (static).

Parameters:
  • centroid (Array)

  • reference (Array)

  • dimensions (int)

  • distance_metric (str)

  • k_neighbors (int)

  • vsa_model_name (str)

Example

>>> import jax, jax.numpy as jnp
>>> from bayes_hdc.anomaly import HDCAnomalyScorer
>>> key = jax.random.PRNGKey(0)
>>> normal = jax.random.normal(key, (50, 1024))
>>> scorer = HDCAnomalyScorer.create(dimensions=1024).fit(normal)
>>> in_dist = float(scorer.score(normal[0]))
>>> far_away = float(scorer.score(normal[0] + 100.0))
>>> in_dist < far_away   # in-distribution scores lower than an outlier
True
centroid: Array
reference: Array
dimensions: int
distance_metric: str = 'cosine'
k_neighbors: int = 1
vsa_model_name: str = 'map'
static create(dimensions=10000, vsa_model='map', distance_metric=None, k_neighbors=1, n_reference=0)[source]

Build an empty scorer. Call fit() with normal hypervectors.

Parameters:
  • dimensions (int) – Hypervector dimensionality.

  • vsa_model (str | Any) – VSA model name ("bsc", "map", "hrr", …) or a VSAModel instance. Used only to pick the default distance metric.

  • distance_metric (str | None) – "cosine" or "hamming". Defaults to "hamming" for BSC and "cosine" otherwise.

  • k_neighbors (int) – 1 for centroid mode (cheap, default), or k > 1 to score against the average similarity over the k nearest reference vectors. Larger k is more robust to multi-modal normal distributions.

  • n_reference (int) – Pre-allocate reference buffer to this many rows. Pass the size of the calibration set when k_neighbors > 1 so fit() does not need to reshape the pytree.

Returns:

An unfitted HDCAnomalyScorer.

Return type:

HDCAnomalyScorer

fit(normal_hypervectors)[source]

Learn the reference centroid (and reference set if k_neighbors > 1).

Parameters:

normal_hypervectors (Array) – Array of shape (n, dimensions) of hypervectors drawn from the in-distribution / normal class.

Returns:

A new HDCAnomalyScorer with centroid (and reference) populated.

Raises:
  • ValueError – If normal_hypervectors is empty.

  • ValueError – If the second axis does not match self.dimensions.

Return type:

HDCAnomalyScorer

score(query)[source]

Compute the nonconformity score of a single query.

Returns:

Scalar jax.Array \(s(x; \mathcal{R}) \geq 0\), larger = more anomalous.

Parameters:

query (Array)

Return type:

Array

score_batch(queries)[source]

Vectorised score over a batch of queries.

Parameters:

queries (Array) – Shape (batch, dimensions).

Returns:

Scores of shape (batch,).

Return type:

Array

replace(**updates)[source]

Pytree-friendly functional update.

Parameters:

updates (Any)

Return type:

HDCAnomalyScorer

__init__(centroid, reference, dimensions, distance_metric='cosine', k_neighbors=1, vsa_model_name='map')
Parameters:
  • centroid (Array)

  • reference (Array)

  • dimensions (int)

  • distance_metric (str)

  • k_neighbors (int)

  • vsa_model_name (str)

Return type:

None

Detector

class bayes_hdc.anomaly.ConformalAnomalyDetector(scorer, calibration_scores, n_calibration=0)[source]

Bases: object

Split-conformal one-class anomaly detector with calibrated FPR.

Wraps an HDCAnomalyScorer and a held-out calibration set of normal nonconformity scores \(\{\alpha_1, \dots, \alpha_n\}\). For a query \(x\) with score \(\alpha(x)\), the conformal p-value is

\[p(x) \;=\; \frac{1 + |\{i : \alpha_i \geq \alpha(x)\}|}{n + 1},\]

which is uniformly distributed on \([0, 1]\) under the exchangeability null \(x \sim P_{\text{normal}}\) (Laxhammar, 2014; Lei et al., 2018; Bates et al., 2023). The decision rule

\[\widehat{\text{anomaly}}(x) \;=\; \mathbb{1}\{p(x) \leq \alpha\}\]

has finite-sample false-positive rate at most \(\alpha\) on exchangeable normal test data, with no assumption on the distribution or on the quality of the underlying score. The coverage guarantee is

\[\mathbb{P}_{x \sim P_{\text{normal}}}\!\left(p(x) \leq \alpha\right) \;\leq\; \alpha.\]
Variables:
  • scorer (bayes_hdc.anomaly.HDCAnomalyScorer) – The fitted HDCAnomalyScorer.

  • calibration_scores (jax.jaxlib._jax.Array) – Sorted (ascending) calibration scores of shape (n_calibration,).

  • n_calibration (int) – Calibration-set size (static, for the \((n+1)\) correction).

Parameters:

Example

>>> import jax, jax.numpy as jnp
>>> from bayes_hdc.anomaly import HDCAnomalyScorer, ConformalAnomalyDetector
>>> key = jax.random.PRNGKey(0)
>>> normal = jax.random.normal(key, (200, 1024))
>>> scorer = HDCAnomalyScorer.create(dimensions=1024).fit(normal[:100])
>>> detector = ConformalAnomalyDetector.create(scorer).fit(normal[100:])
>>> p = float(detector.pvalue(normal[0]))    # in-distribution -> not extreme
>>> p_out = float(detector.pvalue(-normal[0]))  # anti-correlated -> anomalous
>>> p > p_out and p_out < 0.05
True
scorer: HDCAnomalyScorer
calibration_scores: Array
n_calibration: int = 0
static create(scorer, n_calibration=0)[source]

Build an unfitted detector around a (possibly unfitted) scorer.

Parameters:
  • scorer (HDCAnomalyScorer) – An HDCAnomalyScorer. Will be re-fit inside fit() if no normal data has been bundled into it yet, otherwise reused as-is.

  • n_calibration (int) – Pre-allocate the calibration buffer to this many entries. Pass the size of the calibration set so the pytree shape is static under JIT.

Returns:

An unfitted ConformalAnomalyDetector.

Return type:

ConformalAnomalyDetector

fit(normal_data_hvs)[source]

Learn the calibration distribution of normal nonconformity scores.

The scorer is left untouched (use HDCAnomalyScorer.fit() on a separate “training” split, then call this fit() on a held-out calibration split — that is the split-conformal protocol of Lei et al. (2018)). If the underlying scorer has not yet been fitted (zero-norm centroid), it is auto-fitted on the same data as a convenience; this collapses into the in-sample variant and the FPR guarantee then holds only asymptotically.

Parameters:

normal_data_hvs (Array) – Calibration hypervectors of shape (n_calibration, dimensions). All assumed to be drawn i.i.d. (or at least exchangeably) from the normal class.

Returns:

A new ConformalAnomalyDetector with calibration scores stored.

Raises:

ValueError – If normal_data_hvs is empty.

Return type:

ConformalAnomalyDetector

score(query_hv)[source]

Raw nonconformity score before conformalisation.

Returns:

Scalar jax.Array; identical to scorer.score(query_hv).

Parameters:

query_hv (Array)

Return type:

Array

pvalue(query_hv)[source]

Conformal p-value \(p(x) \in (0, 1]\).

\(p(x) = \frac{1 + |\{i : \alpha_i \geq \alpha(x)\}|}{n + 1}\). Uniform on \([0, 1]\) under the exchangeability null; smaller values are stronger evidence against normality.

Parameters:

query_hv (Array) – Query hypervector of shape (dimensions,).

Returns:

Scalar p-value as jax.Array.

Return type:

Array

pvalue_batch(queries)[source]

Vectorised p-values over a batch of queries.

Parameters:

queries (Array) – Shape (batch, dimensions).

Returns:

P-values of shape (batch,).

Return type:

Array

predict(query_hv, alpha=0.05)[source]

Boolean anomaly flag at miscoverage level \(\alpha\).

Returns True iff the conformal p-value is at most \(\alpha\). The marginal false-positive rate of this rule is bounded by \(\alpha\) under exchangeability (Laxhammar, 2014; Bates et al., 2023).

Parameters:
  • query_hv (Array) – Query hypervector of shape (dimensions,).

  • alpha (float) – Target false-positive rate in \((0, 1)\). The value is read at call time, so the same fitted detector can be queried at multiple \(\alpha\).

Returns:

Boolean jax.Array scalar.

Raises:

ValueError – If alpha is not strictly in \((0, 1)\).

Return type:

Array

predict_batch(queries, alpha=0.05)[source]

Vectorised predict() over a batch of queries.

Parameters:
  • queries (Array) – Shape (batch, dimensions).

  • alpha (float) – Target false-positive rate in \((0, 1)\).

Returns:

Boolean array of shape (batch,).

Return type:

Array

predict_fdr(queries, q=0.1)[source]

Flag anomalies in a batch with false-discovery-rate control.

Applies the Benjamini-Hochberg (BH) procedure to the conformal p-values of the batch. Because split-conformal p-values are positively dependent (PRDS), BH controls the expected fraction of false discoveries among the flagged points at level \(q\) (Bates, Candès, Lei & Romano, 2023, Testing for outliers with conformal p-values, Annals of Statistics):

\[\mathbb{E}\!\left[ \frac{\#\{\text{true normal points flagged}\}} {\#\{\text{points flagged}\} \vee 1} \right] \le q .\]

This is the right control when you score many points at once and care about the purity of the alarm set, rather than the per-point false-positive rate that predict_batch() controls. The two answer different questions; use FDR control for batch screening (scan logs / a wafer / a patient cohort), per-point \(\alpha\) for an online single-point decision.

Parameters:
  • queries (Array) – Batch of query hypervectors, shape (batch, dimensions).

  • q (float) – Target false-discovery rate in \((0, 1)\).

Returns:

Boolean array of shape (batch,)True where the point is flagged as an anomaly under BH at level q.

Raises:

ValueError – If q is not strictly in \((0, 1)\).

Return type:

Array

replace(**updates)[source]

Pytree-friendly functional update.

Parameters:

updates (Any)

Return type:

ConformalAnomalyDetector

__init__(scorer, calibration_scores, n_calibration=0)
Parameters:
Return type:

None

End-to-end pipeline

bayes_hdc.anomaly.fit_anomaly_pipeline(encoder, normal_data, calibration_data, alpha=0.05, distance_metric=None, k_neighbors=1)[source]

Encode raw features then fit a split-conformal anomaly detector.

Convenience wrapper around the two-step protocol of Lei et al. (2018):

  1. Encode raw features normal_data and calibration_data with encoder.encode_batch to produce hypervectors.

  2. Fit an HDCAnomalyScorer on the encoded normal data (the proper-training split).

  3. Fit a ConformalAnomalyDetector on the encoded calibration data (the calibration split). The scorer is reused — the calibration data only updates the empirical distribution of nonconformity scores, not the centroid.

This protocol gives the strict marginal FPR \(\leq \alpha\) guarantee: normal_data and calibration_data must come from disjoint splits exchangeable with the test distribution. Calling this function with the same data for both arguments is in-sample and gives only an asymptotic guarantee.

Parameters:
  • encoder (Any) – Any object exposing encode_batch(features) -> hypervectors and a dimensions attribute (e.g. ProjectionEncoder).

  • normal_data (Array) – Raw features of shape (n_train, *feature_shape) from the in-distribution / normal class.

  • calibration_data (Array) – Held-out raw features of shape (n_cal, *feature_shape), also from the normal class.

  • alpha (float) – Target false-positive rate (kept on the returned detector for the user’s convenience; the actual decision is made inside ConformalAnomalyDetector.predict()).

  • distance_metric (str | None) – Override the scorer’s distance metric. Defaults to "hamming" for BSC encoders, "cosine" otherwise.

  • k_neighbors (int) – Forwarded to HDCAnomalyScorer.create().

Returns:

A fitted ConformalAnomalyDetector ready for pvalue() and predict().

Raises:

ValueError – If either normal_data or calibration_data is empty, or if alpha is not in \((0, 1)\).

Return type:

ConformalAnomalyDetector

Example

>>> import jax, jax.numpy as jnp
>>> from bayes_hdc.embeddings import ProjectionEncoder
>>> from bayes_hdc.anomaly import fit_anomaly_pipeline
>>> key = jax.random.PRNGKey(0)
>>> enc = ProjectionEncoder.create(input_dim=16, dimensions=512, key=key)
>>> normal = jax.random.normal(jax.random.PRNGKey(1), (200, 16))
>>> det = fit_anomaly_pipeline(enc, normal[:100], normal[100:], alpha=0.05)
>>> bool(det.predict(enc.encode(normal[0]), alpha=0.05))
False

See also

  • Tutorials02_anomaly_detection.py walks through the coverage guarantee, the conformal-vs-naive-threshold comparison, multi-VSA mode, and a tabular fraud-style demo.

  • examples/anomaly_detection_intrusion.py and examples/anomaly_detection_sensors.py — applied demos.