mirror of
https://github.com/bitcoinresearchkit/brk.git
synced 2026-07-17 05:58:11 -07:00
global: snapshot
This commit is contained in:
@@ -2,41 +2,36 @@ use brk_types::StoredF32;
|
||||
|
||||
/// Fast expanding percentile tracker using a Fenwick tree (Binary Indexed Tree).
|
||||
///
|
||||
/// Values are discretized to BasisPoints32 precision (×10000) and tracked in
|
||||
/// Values are discretized to 10 BPS (0.1%) resolution and tracked in
|
||||
/// a fixed-size frequency array with Fenwick prefix sums. This gives:
|
||||
/// - O(log N) insert (N = tree size, ~18 ops for 200k buckets)
|
||||
/// - O(log N) insert (N = tree size, ~16 ops for 43k buckets)
|
||||
/// - O(log N) percentile query via prefix-sum walk
|
||||
/// - Exact at BasisPoints32 resolution (no approximation)
|
||||
/// - 0.1% value resolution (10 BPS granularity)
|
||||
#[derive(Clone)]
|
||||
pub(crate) struct ExpandingPercentiles {
|
||||
/// Fenwick tree storing cumulative frequency counts.
|
||||
/// Index 0 is unused (1-indexed). tree[i] covers bucket (i - 1 + offset).
|
||||
tree: Vec<u64>,
|
||||
count: u64,
|
||||
/// Offset so bucket 0 in the tree corresponds to BPS value `offset`.
|
||||
offset: i32,
|
||||
size: usize,
|
||||
/// 1-indexed: tree[0] is unused, tree[1..=TREE_SIZE] hold data.
|
||||
tree: Vec<u32>,
|
||||
count: u32,
|
||||
}
|
||||
|
||||
/// Max BPS value supported. Ratio of 42.0 = 420,000 BPS.
|
||||
/// Bucket granularity in BPS. 10 BPS = 0.1% = 0.001 ratio.
|
||||
const BUCKET_BPS: i32 = 10;
|
||||
/// Max ratio supported: 43.0 = 430,000 BPS.
|
||||
const MAX_BPS: i32 = 430_000;
|
||||
/// Min BPS value supported (0 = ratio of 0.0).
|
||||
const MIN_BPS: i32 = 0;
|
||||
const TREE_SIZE: usize = (MAX_BPS - MIN_BPS) as usize + 1;
|
||||
const TREE_SIZE: usize = (MAX_BPS / BUCKET_BPS) as usize + 1;
|
||||
|
||||
impl Default for ExpandingPercentiles {
|
||||
fn default() -> Self {
|
||||
Self {
|
||||
tree: vec![0u64; TREE_SIZE + 1], // 1-indexed
|
||||
tree: vec![0u32; TREE_SIZE + 1], // 1-indexed
|
||||
count: 0,
|
||||
offset: MIN_BPS,
|
||||
size: TREE_SIZE,
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
impl ExpandingPercentiles {
|
||||
pub fn count(&self) -> u64 {
|
||||
pub fn count(&self) -> u32 {
|
||||
self.count
|
||||
}
|
||||
|
||||
@@ -47,29 +42,27 @@ impl ExpandingPercentiles {
|
||||
|
||||
/// Convert f32 ratio to bucket index (1-indexed for Fenwick).
|
||||
#[inline]
|
||||
fn to_bucket(&self, value: f32) -> usize {
|
||||
fn to_bucket(value: f32) -> usize {
|
||||
let bps = (value as f64 * 10000.0).round() as i32;
|
||||
let clamped = bps.clamp(self.offset, self.offset + self.size as i32 - 1);
|
||||
(clamped - self.offset) as usize + 1 // 1-indexed
|
||||
let bucket = (bps / BUCKET_BPS).clamp(0, TREE_SIZE as i32 - 1);
|
||||
bucket as usize + 1
|
||||
}
|
||||
|
||||
/// Bulk-load values in O(n + N) instead of O(n log N).
|
||||
/// Builds raw frequency counts, then converts to Fenwick in-place.
|
||||
pub fn add_bulk(&mut self, values: &[StoredF32]) {
|
||||
// Build raw frequency counts into tree (treated as flat array)
|
||||
for &v in values {
|
||||
let v = *v;
|
||||
if v.is_nan() {
|
||||
continue;
|
||||
}
|
||||
self.count += 1;
|
||||
let bucket = self.to_bucket(v);
|
||||
self.tree[bucket] += 1;
|
||||
self.tree[Self::to_bucket(v)] += 1;
|
||||
}
|
||||
// Convert flat frequencies to Fenwick tree in O(N)
|
||||
for i in 1..=self.size {
|
||||
for i in 1..=TREE_SIZE {
|
||||
let parent = i + (i & i.wrapping_neg());
|
||||
if parent <= self.size {
|
||||
if parent <= TREE_SIZE {
|
||||
let val = self.tree[i];
|
||||
self.tree[parent] += val;
|
||||
}
|
||||
@@ -83,47 +76,40 @@ impl ExpandingPercentiles {
|
||||
return;
|
||||
}
|
||||
self.count += 1;
|
||||
let mut i = self.to_bucket(value);
|
||||
while i <= self.size {
|
||||
let mut i = Self::to_bucket(value);
|
||||
while i <= TREE_SIZE {
|
||||
self.tree[i] += 1;
|
||||
i += i & i.wrapping_neg(); // i += lowbit(i)
|
||||
i += i & i.wrapping_neg();
|
||||
}
|
||||
}
|
||||
|
||||
/// Find the bucket containing the k-th element (1-indexed k).
|
||||
/// Uses the standard Fenwick tree walk-down in O(log N).
|
||||
#[inline]
|
||||
fn kth(&self, mut k: u64) -> usize {
|
||||
fn kth(&self, mut k: u32) -> usize {
|
||||
let mut pos = 0;
|
||||
let mut bit = 1 << (usize::BITS - 1 - self.size.leading_zeros()); // highest power of 2 <= size
|
||||
let mut bit = 1 << (usize::BITS - 1 - TREE_SIZE.leading_zeros());
|
||||
while bit > 0 {
|
||||
let next = pos + bit;
|
||||
if next <= self.size && self.tree[next] < k {
|
||||
if next <= TREE_SIZE && self.tree[next] < k {
|
||||
k -= self.tree[next];
|
||||
pos = next;
|
||||
}
|
||||
bit >>= 1;
|
||||
}
|
||||
pos + 1 // 1-indexed bucket
|
||||
}
|
||||
|
||||
/// Convert bucket index back to BPS u32 value.
|
||||
#[inline]
|
||||
fn bucket_to_bps(&self, bucket: usize) -> u32 {
|
||||
(bucket as i32 - 1 + self.offset) as u32
|
||||
pos + 1
|
||||
}
|
||||
|
||||
/// Compute 6 percentiles in one call. O(6 × log N).
|
||||
/// Quantiles q must be in (0, 1).
|
||||
/// Quantiles q must be in (0, 1). Output is in BPS.
|
||||
pub fn quantiles(&self, qs: &[f64; 6], out: &mut [u32; 6]) {
|
||||
if self.count == 0 {
|
||||
out.iter_mut().for_each(|o| *o = 0);
|
||||
return;
|
||||
}
|
||||
for (i, &q) in qs.iter().enumerate() {
|
||||
// k = ceil(q * count), clamped to [1, count]
|
||||
let k = ((q * self.count as f64).ceil() as u64).clamp(1, self.count);
|
||||
out[i] = self.bucket_to_bps(self.kth(k));
|
||||
let k = ((q * self.count as f64).ceil() as u32).clamp(1, self.count);
|
||||
out[i] = (self.kth(k) as u32 - 1) * BUCKET_BPS as u32;
|
||||
}
|
||||
}
|
||||
}
|
||||
@@ -141,30 +127,19 @@ mod tests {
|
||||
#[test]
|
||||
fn basic_quantiles() {
|
||||
let mut ep = ExpandingPercentiles::default();
|
||||
// Add ratios 0.01 to 1.0 (BPS 100 to 10000)
|
||||
for i in 1..=1000 {
|
||||
ep.add(i as f32 / 1000.0);
|
||||
}
|
||||
assert_eq!(ep.count(), 1000);
|
||||
|
||||
let median = quantile(&ep, 0.5);
|
||||
// 0.5 ratio = 5000 BPS, median of 1..1000 ratios ≈ 500/1000 = 0.5 = 5000 BPS
|
||||
assert!(
|
||||
(median as i32 - 5000).abs() < 100,
|
||||
"median was {median}"
|
||||
);
|
||||
assert!((median as i32 - 5000).abs() < 100, "median was {median}");
|
||||
|
||||
let p99 = quantile(&ep, 0.99);
|
||||
assert!(
|
||||
(p99 as i32 - 9900).abs() < 100,
|
||||
"p99 was {p99}"
|
||||
);
|
||||
assert!((p99 as i32 - 9900).abs() < 100, "p99 was {p99}");
|
||||
|
||||
let p01 = quantile(&ep, 0.01);
|
||||
assert!(
|
||||
(p01 as i32 - 100).abs() < 100,
|
||||
"p01 was {p01}"
|
||||
);
|
||||
assert!((p01 as i32 - 100).abs() < 100, "p01 was {p01}");
|
||||
}
|
||||
|
||||
#[test]
|
||||
@@ -177,10 +152,9 @@ mod tests {
|
||||
#[test]
|
||||
fn single_value() {
|
||||
let mut ep = ExpandingPercentiles::default();
|
||||
ep.add(0.42); // 4200 BPS
|
||||
assert_eq!(quantile(&ep, 0.0001), 4200);
|
||||
assert_eq!(quantile(&ep, 0.5), 4200);
|
||||
assert_eq!(quantile(&ep, 0.9999), 4200);
|
||||
ep.add(0.42);
|
||||
let v = quantile(&ep, 0.5);
|
||||
assert!((v as i32 - 4200).abs() <= BUCKET_BPS, "got {v}");
|
||||
}
|
||||
|
||||
#[test]
|
||||
|
||||
@@ -61,6 +61,19 @@ where
|
||||
T: Default + SubAssign,
|
||||
{
|
||||
compute_height(&mut self.height)?;
|
||||
self.compute_rest(max_from, windows, exit)
|
||||
}
|
||||
|
||||
/// Compute cumulative + rolling sum from already-populated height data.
|
||||
pub(crate) fn compute_rest(
|
||||
&mut self,
|
||||
max_from: Height,
|
||||
windows: &WindowStarts<'_>,
|
||||
exit: &Exit,
|
||||
) -> Result<()>
|
||||
where
|
||||
T: Default + SubAssign,
|
||||
{
|
||||
self.cumulative
|
||||
.height
|
||||
.compute_cumulative(max_from, &self.height, exit)?;
|
||||
|
||||
@@ -1,6 +1,6 @@
|
||||
use brk_error::Result;
|
||||
use brk_traversable::{Traversable, TreeNode};
|
||||
use brk_types::{BasisPoints16, Cents, Height, Version};
|
||||
use brk_types::{Cents, Height, Version};
|
||||
use vecdb::{AnyExportableVec, Database, ReadOnlyClone, Ro, Rw, StorageMode, WritableVec};
|
||||
|
||||
use crate::indexes;
|
||||
@@ -11,62 +11,6 @@ pub const PERCENTILES: [u8; 19] = [
|
||||
];
|
||||
pub const PERCENTILES_LEN: usize = PERCENTILES.len();
|
||||
|
||||
/// Compute spot percentile rank by interpolating within percentile bands.
|
||||
/// Returns a value between 0 and 100 indicating where spot sits in the distribution.
|
||||
pub(crate) fn compute_spot_percentile_rank(
|
||||
percentile_prices: &[Cents; PERCENTILES_LEN],
|
||||
spot: Cents,
|
||||
) -> BasisPoints16 {
|
||||
if spot == Cents::ZERO && percentile_prices[0] == Cents::ZERO {
|
||||
return BasisPoints16::ZERO;
|
||||
}
|
||||
|
||||
let spot_f64 = f64::from(spot);
|
||||
|
||||
// Below lowest percentile (p5) - extrapolate towards 0
|
||||
let p5 = f64::from(percentile_prices[0]);
|
||||
if spot_f64 <= p5 {
|
||||
if p5 == 0.0 {
|
||||
return BasisPoints16::ZERO;
|
||||
}
|
||||
// Linear extrapolation: rank = 5% * (spot / p5)
|
||||
return BasisPoints16::from((0.05 * spot_f64 / p5).max(0.0));
|
||||
}
|
||||
|
||||
// Above highest percentile (p95) - extrapolate towards 100
|
||||
let p95 = f64::from(percentile_prices[PERCENTILES_LEN - 1]);
|
||||
let p90 = f64::from(percentile_prices[PERCENTILES_LEN - 2]);
|
||||
if spot_f64 >= p95 {
|
||||
if p95 == p90 {
|
||||
return BasisPoints16::ONE;
|
||||
}
|
||||
// Linear extrapolation using p90-p95 slope
|
||||
let slope = 0.05 / (p95 - p90);
|
||||
return BasisPoints16::from((0.95 + (spot_f64 - p95) * slope).min(1.0));
|
||||
}
|
||||
|
||||
// Find the band containing spot and interpolate
|
||||
for i in 0..PERCENTILES_LEN - 1 {
|
||||
let lower = f64::from(percentile_prices[i]);
|
||||
let upper = f64::from(percentile_prices[i + 1]);
|
||||
|
||||
if spot_f64 >= lower && spot_f64 <= upper {
|
||||
let lower_pct = f64::from(PERCENTILES[i]) / 100.0;
|
||||
let upper_pct = f64::from(PERCENTILES[i + 1]) / 100.0;
|
||||
|
||||
if upper == lower {
|
||||
return BasisPoints16::from(lower_pct);
|
||||
}
|
||||
|
||||
// Linear interpolation
|
||||
let ratio = (spot_f64 - lower) / (upper - lower);
|
||||
return BasisPoints16::from(lower_pct + ratio * (upper_pct - lower_pct));
|
||||
}
|
||||
}
|
||||
|
||||
BasisPoints16::ZERO
|
||||
}
|
||||
|
||||
pub struct PercentilesVecs<M: StorageMode = Rw> {
|
||||
pub vecs: [Price<ComputedFromHeight<Cents, M>>; PERCENTILES_LEN],
|
||||
}
|
||||
|
||||
@@ -5,20 +5,20 @@ use derive_more::{Deref, DerefMut};
|
||||
use vecdb::{Database, EagerVec, Exit, PcoVec, Rw, StorageMode};
|
||||
|
||||
use crate::internal::{ComputedFromHeight, Price};
|
||||
use crate::{blocks, indexes, prices};
|
||||
use crate::{indexes, prices};
|
||||
|
||||
use super::ComputedFromHeightRatioExtended;
|
||||
use super::ComputedFromHeightRatio;
|
||||
|
||||
#[derive(Deref, DerefMut, Traversable)]
|
||||
pub struct ComputedFromHeightPriceWithRatioExtended<M: StorageMode = Rw> {
|
||||
pub struct ComputedFromHeightPriceWithRatio<M: StorageMode = Rw> {
|
||||
#[deref]
|
||||
#[deref_mut]
|
||||
#[traversable(flatten)]
|
||||
pub inner: ComputedFromHeightRatioExtended<M>,
|
||||
pub inner: ComputedFromHeightRatio<M>,
|
||||
pub price: Price<ComputedFromHeight<Cents, M>>,
|
||||
}
|
||||
|
||||
impl ComputedFromHeightPriceWithRatioExtended {
|
||||
impl ComputedFromHeightPriceWithRatio {
|
||||
pub(crate) fn forced_import(
|
||||
db: &Database,
|
||||
name: &str,
|
||||
@@ -27,15 +27,14 @@ impl ComputedFromHeightPriceWithRatioExtended {
|
||||
) -> Result<Self> {
|
||||
let v = version + Version::TWO;
|
||||
Ok(Self {
|
||||
inner: ComputedFromHeightRatioExtended::forced_import(db, name, version, indexes)?,
|
||||
inner: ComputedFromHeightRatio::forced_import(db, name, version, indexes)?,
|
||||
price: Price::forced_import(db, name, v, indexes)?,
|
||||
})
|
||||
}
|
||||
|
||||
/// Compute price via closure (in cents), then compute ratio + extended metrics.
|
||||
/// Compute price via closure (in cents), then compute ratio.
|
||||
pub(crate) fn compute_all<F>(
|
||||
&mut self,
|
||||
blocks: &blocks::Vecs,
|
||||
prices: &prices::Vecs,
|
||||
starting_indexes: &Indexes,
|
||||
exit: &Exit,
|
||||
@@ -45,13 +44,9 @@ impl ComputedFromHeightPriceWithRatioExtended {
|
||||
F: FnMut(&mut EagerVec<PcoVec<Height, Cents>>) -> Result<()>,
|
||||
{
|
||||
compute_price(&mut self.price.cents.height)?;
|
||||
self.inner.compute_rest(
|
||||
blocks,
|
||||
prices,
|
||||
starting_indexes,
|
||||
exit,
|
||||
&self.price.cents.height,
|
||||
)?;
|
||||
let close_price = &prices.price.cents.height;
|
||||
self.inner
|
||||
.compute_ratio(starting_indexes, close_price, &self.price.cents.height, exit)?;
|
||||
Ok(())
|
||||
}
|
||||
}
|
||||
|
||||
@@ -23,7 +23,7 @@ pub use derived::{
|
||||
RatioCents64, TimesSqrt,
|
||||
};
|
||||
pub use ratio::{
|
||||
NegRatioDollarsBps32, RatioCentsBp16, RatioCentsBp32, RatioCentsSignedCentsBps32,
|
||||
NegRatioDollarsBps32, RatioCentsBp32, RatioCentsSignedCentsBps32,
|
||||
RatioCentsSignedDollarsBps32, RatioDiffCentsBps32, RatioDiffDollarsBps32, RatioDiffF32Bps32,
|
||||
RatioDollarsBp16, RatioDollarsBp32, RatioDollarsBps32, RatioSatsBp16, RatioU32Bp16,
|
||||
RatioU64Bp16,
|
||||
|
||||
@@ -30,19 +30,6 @@ impl BinaryTransform<Sats, Sats, BasisPoints16> for RatioSatsBp16 {
|
||||
}
|
||||
}
|
||||
|
||||
pub struct RatioCentsBp16;
|
||||
|
||||
impl BinaryTransform<Cents, Cents, BasisPoints16> for RatioCentsBp16 {
|
||||
#[inline(always)]
|
||||
fn apply(numerator: Cents, denominator: Cents) -> BasisPoints16 {
|
||||
if denominator == Cents::ZERO {
|
||||
BasisPoints16::ZERO
|
||||
} else {
|
||||
BasisPoints16::from(numerator.inner() as f64 / denominator.inner() as f64)
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
pub struct RatioCentsBp32;
|
||||
|
||||
impl BinaryTransform<Cents, Cents, BasisPoints32> for RatioCentsBp32 {
|
||||
@@ -143,7 +130,12 @@ pub struct RatioDollarsBp32;
|
||||
impl BinaryTransform<Dollars, Dollars, BasisPoints32> for RatioDollarsBp32 {
|
||||
#[inline(always)]
|
||||
fn apply(numerator: Dollars, denominator: Dollars) -> BasisPoints32 {
|
||||
BasisPoints32::from(f64::from(numerator) / f64::from(denominator))
|
||||
let ratio = f64::from(numerator) / f64::from(denominator);
|
||||
if ratio.is_finite() {
|
||||
BasisPoints32::from(ratio)
|
||||
} else {
|
||||
BasisPoints32::ZERO
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
Reference in New Issue
Block a user