From 7ab48cdd7b58c1c6e3b58ca75f48275a58fa60c1 Mon Sep 17 00:00:00 2001 From: Alexander Neubeck Date: Mon, 6 Jul 2026 17:10:15 +0200 Subject: [PATCH 1/3] try a direct function instead of a loop --- crates/geo_filters/src/diff_count/config.rs | 119 ++++++++++++++++++++ crates/geo_filters/src/diff_count/metric.rs | 15 +-- 2 files changed, 125 insertions(+), 9 deletions(-) diff --git a/crates/geo_filters/src/diff_count/config.rs b/crates/geo_filters/src/diff_count/config.rs index b34b2aa..b87e369 100644 --- a/crates/geo_filters/src/diff_count/config.rs +++ b/crates/geo_filters/src/diff_count/config.rs @@ -103,6 +103,79 @@ fn expected_diff_buckets_fast(phi: f64, items: f64) -> (f64, f64) { (sum, derivative) } +/// Euler–Mascheroni constant. +const GAMMA: f64 = 0.577_215_664_901_532_9; +/// `e^gamma`, precomputed. +const E_GAMMA: f64 = 1.781_072_417_990_198; + +// Closed-form approximation of [`expected_diff_buckets`] and its inverse that avoids the +// per-bucket summation (and the Newton iteration on it). It is exact in both limits and stays +// within ~0.35% of the exact model for the predefined configurations (b = 7, 10, 13). +// +// Writing `p_k = (1 - phi) phi^k`, the expected number of one-bits is +// E(N) = 1/2 sum_k (1 - (1 - 2 p_k)^N) ≈ S · Ein(x) + 1/4 (1 - e^-x), +// with `x = 2(1 - phi) N`, `S = 1 / (2 ln(1/phi))`, and `Ein(x) = gamma + ln x + E1(x)` the +// entire exponential-integral bridge (`Ein(x) ≈ x` near 0, `≈ gamma + ln x` for large x). The +// special function `Ein` is replaced by a rational "scaling factor" `rho(x)` inside the logarithm, +// Ein(x) ≈ ln(1 + e^gamma x rho(x)), rho(x) = (e^-gamma + A x + B x^2) / (1 + C x + B x^2), +// where `rho(0) = e^-gamma` reproduces the linear small-N regime and `rho(inf) = 1` reproduces the +// logarithmic asymptote. The coefficients are universal (independent of phi), fitted so the +// rational matches `Ein` to < 0.07%. +const RHO_A: f64 = 0.384_417; +const RHO_B: f64 = 0.130_468; +const RHO_C: f64 = 0.442_202; + +/// Closed-form approximation of [`expected_diff_buckets`]: the expected number of one-bits for +/// `items` items, together with its derivative with respect to `items`. See the comment above for +/// the derivation. Accurate to < 1% for the predefined configurations. +pub(super) fn expected_diff_buckets_approx(phi: f64, items: f64) -> (f64, f64) { + if items <= 0.0 { + return (0.0, 1.0); + } + let s = 0.5 / -phi.ln(); // = 1 / (2 ln(1/phi)) + let a = 2.0 * (1.0 - phi); + let x = a * items; + let e_neg_x = (-x).exp(); + + // Y = e^gamma x rho(x) = (x + e^gamma A x^2 + e^gamma B x^3) / (1 + C x + B x^2) + let num = x * (1.0 + E_GAMMA * RHO_A * x + E_GAMMA * RHO_B * x * x); + let den = 1.0 + RHO_C * x + RHO_B * x * x; + let y = num / den; + let value = s * (1.0 + y).ln() + 0.25 * (1.0 - e_neg_x); + + // derivative w.r.t. items = a * d/dx [ s ln(1 + Y) + 1/4 (1 - e^-x) ] + let num_prime = 1.0 + 2.0 * E_GAMMA * RHO_A * x + 3.0 * E_GAMMA * RHO_B * x * x; + let den_prime = RHO_C + 2.0 * RHO_B * x; + let y_prime = (num_prime * den - num * den_prime) / (den * den); + let derivative = a * (s * y_prime / (1.0 + y) + 0.25 * e_neg_x); + + (value, derivative) +} + +/// Closed-form approximation of the inverse of [`expected_diff_buckets`]: the number of items that +/// produce `ones` one-bits. Uses the analytic large-N seed of the model above followed by three +/// Newton steps on the (cheap) closed-form forward, which converge to full `f64` precision. +pub(super) fn estimate_count_approx(phi: f64, ones: f64) -> f64 { + if ones <= 0.0 { + return 0.0; + } + let s = 0.5 / -phi.ln(); + let a = 2.0 * (1.0 - phi); + // Invert the large-N asymptote `m ≈ S (gamma + ln(a N)) + 1/4` for the initial guess. + let mut items = ((ones - 0.25) / s - GAMMA).exp() / a; + if !items.is_finite() || items <= 0.0 { + items = ones; + } + for _ in 0..3 { + let (m, dm) = expected_diff_buckets_approx(phi, items); + items -= (m - ones) / dm; + if items < f64::MIN_POSITIVE { + items = f64::MIN_POSITIVE; + } + } + items +} + static LOOKUPS: [Lazy; 16] = [ Lazy::new(|| build_lookup(0)), Lazy::new(|| build_lookup(1)), @@ -180,6 +253,52 @@ mod tests { println!("{:?}", estimate_count(phi, 78.450, expected_diff_buckets)); } + #[test] + fn test_expected_diff_buckets_approx() { + // The closed-form forward stays within 1% of the exact model for the predefined configs. + for b in [7u32, 10, 13] { + let phi = 0.5f64.powf(1.0 / (1u64 << b) as f64); + let mut items = 1.0f64; + loop { + let exact = expected_diff_buckets(phi, items).0; + if exact > 2000.0 { + break; + } + if exact > 0.3 { + let approx = expected_diff_buckets_approx(phi, items).0; + let rel = (approx - exact).abs() / exact; + assert!( + rel < 0.01, + "b={b} items={items}: exact={exact} approx={approx} rel={rel}" + ); + } + items *= 1.5; + } + } + } + + #[test] + fn test_estimate_count_approx() { + // The closed-form inverse recovers the item count within 1% of the exact model. + for b in [7u32, 10, 13] { + let phi = 0.5f64.powf(1.0 / (1u64 << b) as f64); + let mut items = 1.0f64; + loop { + let ones = expected_diff_buckets(phi, items).0; + if ones > 2000.0 { + break; + } + let est = estimate_count_approx(phi, ones); + let rel = (est - items).abs() / items; + assert!( + rel < 0.01, + "b={b} items={items} ones={ones}: est={est} rel={rel}" + ); + items *= 1.5; + } + } + } + #[test] fn test_estimation_lut_7() { let c = GeoDiffConfig7::::default(); diff --git a/crates/geo_filters/src/diff_count/metric.rs b/crates/geo_filters/src/diff_count/metric.rs index aa70eb2..a85a2c0 100644 --- a/crates/geo_filters/src/diff_count/metric.rs +++ b/crates/geo_filters/src/diff_count/metric.rs @@ -5,10 +5,8 @@ use std::cmp::Ordering; use std::fmt; use std::ops::Add; -use super::config::expected_diff_buckets; -use crate::config::{ - estimate_count, xor_bit_chunks, BitChunk, GeoConfig, IsBucketType, BITS_PER_BLOCK, -}; +use super::config::{estimate_count_approx, expected_diff_buckets_approx}; +use crate::config::{xor_bit_chunks, BitChunk, GeoConfig, IsBucketType, BITS_PER_BLOCK}; use crate::diff_count::GeoDiffCount; use crate::{Diff, Metric, MetricSpace}; @@ -48,24 +46,23 @@ impl + Default> Metric for OnesMetric { } /// The calibrated size (item count) that this number of one-bits represents, obtained by - /// inverting the expected-buckets function with the same Newton iteration used to build the - /// estimation lookup table. + /// inverting the closed-form approximation of the expected-buckets function. fn to_f32(&self) -> f32 { if *self == Self::infinite() { return f32::INFINITY; } let phi = C::default().phi_f64(); - estimate_count(phi, self.get() as f64, expected_diff_buckets).0 as f32 + estimate_count_approx(phi, self.get() as f64) as f32 } /// The number of one-bits expected for the given calibrated size, obtained by evaluating the - /// forward expected-buckets function. + /// closed-form approximation of the forward expected-buckets function. fn from_f32(value: f32) -> Self { if !value.is_finite() { return Self::infinite(); } let phi = C::default().phi_f64(); - let buckets = expected_diff_buckets(phi, value.max(0.0) as f64).0; + let buckets = expected_diff_buckets_approx(phi, value.max(0.0) as f64).0; Self::new(buckets.round() as usize) } } From 2729729c055a1341e6d045db1cc40f87bb3b38da Mon Sep 17 00:00:00 2001 From: Alexander Neubeck Date: Tue, 7 Jul 2026 15:51:07 +0200 Subject: [PATCH 2/3] document the additional noise introduced by the metric estimate --- crates/geo_filters/README.md | 41 +++++++++++++++++++++ crates/geo_filters/src/diff_count/metric.rs | 6 +++ 2 files changed, 47 insertions(+) diff --git a/crates/geo_filters/README.md b/crates/geo_filters/README.md index 9b18617..7e644c4 100644 --- a/crates/geo_filters/README.md +++ b/crates/geo_filters/README.md @@ -161,6 +161,47 @@ estimate), whereas a `near` candidate is below the bound and therefore scanned i and exact costs match. Since a nearest-neighbor scan rejects far more candidates than it keeps, the early abandon dominates the search cost. Reproduce with `cargo bench --bench nearest_neighbor`. +### Precision and configuration choice + +Ranking candidates by their one-bit distance is a *noisy* estimate of the true symmetric-difference +size (calibrating it to an item count with `Metric::to_f32` is order-preserving, so it does not +change the ranking). The **relative error** below is the standard deviation of `estimate / true − 1` +— i.e. relative to the true number of differing items. It is smallest for small differences and grows +**slowly, like √(ln n)**, as the difference size `n` increases: the dense low bits saturate to a +random ½ (encoding only the *parity* of many items, so they add variance without carrying size +information). Modelling the buckets as independent Bernoulli variables, the count's variance is +exactly `½·expected_diff_buckets(2n)`, which is logarithmic in `n`; the resulting relative error is +about `√((γ + ln(4·(1−φ)·n)) / (2·S))` with `S = 2^b / (2·ln 2)` and `γ` the Euler–Mascheroni +constant. For `GeoDiffConfig7` that formula gives ~8.6% at `n = 100`, ~18% at `10⁴`, ~24% at `10⁶`, +close to the measured ~8% / ~18% / ~27% (real-filter bucket correlations push it a little higher at +the top end). + +| configuration | relative error (1σ), small … 1M-item difference | may reorder candidates within | +| ----------------- | ----------------------------------------------- | ----------------------------- | +| `GeoDiffConfig7` | ~8% … ~27% | ~2.4× | +| `GeoDiffConfig10` | ~2.5% … ~6% | ~1.3× | +| `GeoDiffConfig13` | ~0.8% … ~2% | ~1.1× | + +The one-bit count is *not* exact even for a few items — two items can hash to the same bucket — but +such collisions are rare while the difference is small, so it is a low-variance estimate there. +Compared to the standard windowed estimator `GeoDiffCount::size()` (and `size_with_sketch`), the +one-bit count matches it for small differences but is noisier for large ones: for `GeoDiffConfig7` at +a one-million-item difference it is ~27% here versus ~14% for `size()`, which reads only a fixed +window near the fringe and ignores the saturated bits. Use `size_with_sketch` when you need an +accurate absolute difference size; the one-bit metric exists for fast *ranking* of near neighbors, +where its early-abandon speed matters more than the last few percent of precision. + +Because of the coarse resolution, **the metric is best used with `GeoDiffConfig10` or +`GeoDiffConfig13` (`b ≥ 10`)**. `GeoDiffConfig7` reorders candidates whose true distances lie within +roughly a factor of two of each other, so over a large candidate set it returns an *approximate* +nearest neighbor rather than the exact one — though it never promotes a genuinely far candidate over +a close one, since the noise cannot bridge more than the factor in the last column, even across +hundreds of millions of candidates. + +For an exact result, use the metric as a cheap first pass to shortlist the top *k* candidates, then +re-rank the shortlist with the more robust `Count::size_with_sketch`, which does not accumulate this +variance. + ## Evaluation Accuracy and performance evaluations for the predefined filter configurations are included in the repository: diff --git a/crates/geo_filters/src/diff_count/metric.rs b/crates/geo_filters/src/diff_count/metric.rs index a85a2c0..d618826 100644 --- a/crates/geo_filters/src/diff_count/metric.rs +++ b/crates/geo_filters/src/diff_count/metric.rs @@ -118,6 +118,12 @@ impl> fmt::Debug for OnesMetric { /// search - and, given a bound, abandoned early. The cached count (see [`Self::size`]) also yields /// an O(1) reverse-triangle lower bound on the distance. Once the nearest neighbor is found, /// [`Self::filter`] gives access to the underlying filter for a calibrated size estimate. +/// +/// The one-bit distance is a noisy estimate of the true difference size, so this metric is best +/// used with `GeoDiffConfig10` or `GeoDiffConfig13` (`b >= 10`). The coarser `GeoDiffConfig7` can +/// reorder candidates whose true distances are within roughly a factor of two of each other, +/// returning an approximate rather than exact nearest neighbor; for an exact result, shortlist with +/// this metric and re-rank the top candidates with [`Count::size_with_sketch`](crate::Count). pub struct GeoDiffMetric<'a, C: GeoConfig> { filter: GeoDiffCount<'a, C>, ones: OnesMetric, From 7708b730b0b0f15f738bd19ee8851751fe3aa23e Mon Sep 17 00:00:00 2001 From: Alexander Neubeck Date: Tue, 7 Jul 2026 16:34:16 +0200 Subject: [PATCH 3/3] Update Cargo.toml --- crates/geo_filters/Cargo.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/crates/geo_filters/Cargo.toml b/crates/geo_filters/Cargo.toml index 3c2ba97..38aab10 100644 --- a/crates/geo_filters/Cargo.toml +++ b/crates/geo_filters/Cargo.toml @@ -1,6 +1,6 @@ [package] name = "geo_filters" -version = "0.1.1" +version = "0.2.0" edition = "2021" description = "Geometric filters for set cardinality estimation." repository = "https://github.com/github/rust-gems"