code wiki / (root) / sketch_stream_stats.nx

sketch_stream_stats.nx source

↩ module page · 207 lines · 6613 B

1// sketch_stream_stats.nx -- streaming mean / variance / stddev / min / max. 2// 3// Fundamental single-pass streaming primitive. Maintains six i64 fields 4// across the stream: 5// count = N 6// sum = Σ x_i 7// sum_sq = Σ x_i² 8// min = min over stream 9// max = max over stream 10// 11// Algebraic-identity-based estimators (NOT Welford, which is f64-stable but 12// loses precision when arithmetic is exact-integer): 13// mean = sum / count 14// variance = sum_sq / count - mean² (population) 15// stddev = isqrt(variance) 16// 17// OVERFLOW BUDGET: 18// For sample values bounded by |x| <= V and count N: 19// sum: N * V — must fit i64. Safe to N * V < 2^62. 20// sum_sq: N * V² — same constraint applied to V² instead of V. 21// For V = 2^30 (1 billion): N up to ~2^32 (4 billion) safe for sum, 22// but sum_sq overflows at N * 2^60 -- safe to N=4. 23// For V = 2^20 (~1M): N up to 2^22 (~4M) safe for sum_sq. 24// For V = 2^16 (65K): N up to 2^30 (1B) safe. 25// Caller responsibility to size domain; nx_stats_safe_p returns 1 iff 26// adding `value` would NOT overflow. 27// 28// COMPLEMENTS the sketch suite: 29// - sketch_reservoir: arbitrary-statistic estimation via sample-based math 30// - sketch_kll / sketch_tdigest: quantiles (median, p99) 31// - sketch_stream_stats (here): exact-integer mean/variance over the FULL 32// stream (no sampling, no probabilistic bounds) 33// 34// LOSSLESS-LANGUAGE DISCIPLINE: nx_stats_query_mean and friends return 35// ApproxI64 with envelope_kind = NX_ENV_ABS, param_a = 0 (mean and 36// variance are EXACT under the overflow budget), conf_ppb = 1e9. 37 38import "syscalls.nx" 39import "sketch_types.nx" 40 41struct StreamStats { 42 count: i64, 43 sum: i64, 44 sum_sq: i64, 45 min_val: i64, 46 max_val: i64, 47 has_data: i64, // 0 if empty, 1 if any add happened 48} 49 50// === construction ================================================= 51 52func nx_stats_alloc() -> *StreamStats { 53 let raw: *u8 = sys_mmap(56) 54 let s: *StreamStats = raw as *StreamStats 55 s.count = 0 56 s.sum = 0 57 s.sum_sq = 0 58 s.min_val = 0 59 s.max_val = 0 60 s.has_data = 0 61 return s 62} 63 64// === overflow guard ================================================= 65// 66// Returns 1 iff adding `value` is overflow-safe given the current 67// accumulator state. Conservative: declares unsafe if sum_sq + v² would 68// hit upper 2 bits of i64 (margin for next-step arithmetic). 69 70const NX_STATS_SAFE_HI: i64 = 0x2000000000000000 // 2^61 71 72func nx_stats_safe_p(s: *StreamStats, value: i64) -> i64 { 73 var v: i64 = value 74 if v < 0 { v = -v } 75 if v >= 0x40000000 { return 0 } // |v| >= 2^30 -> v² overflows i64 76 let v_sq: i64 = v * v 77 if s.sum > NX_STATS_SAFE_HI - v { return 0 } 78 if s.sum_sq > NX_STATS_SAFE_HI - v_sq { return 0 } 79 return 1 80} 81 82// === add ========================================================== 83 84func nx_stats_add(s: *StreamStats, value: i64) -> i64 { 85 if nx_stats_safe_p(s, value) == 0 { return -1 } 86 s.count = s.count + 1 87 s.sum = s.sum + value 88 s.sum_sq = s.sum_sq + value * value 89 if s.has_data == 0 { 90 s.min_val = value 91 s.max_val = value 92 s.has_data = 1 93 } 94 if s.has_data == 1 { 95 if value < s.min_val { s.min_val = value } 96 if value > s.max_val { s.max_val = value } 97 } 98 return 0 99} 100 101// === queries ====================================================== 102 103func nx_stats_mean(s: *StreamStats) -> i64 { 104 if s.count == 0 { return 0 } 105 return s.sum / s.count 106} 107 108func nx_stats_variance(s: *StreamStats) -> i64 { 109 if s.count == 0 { return 0 } 110 let m: i64 = nx_stats_mean(s) 111 let e_sq: i64 = s.sum_sq / s.count 112 let m_sq: i64 = m * m 113 if e_sq < m_sq { return 0 } // shouldn't happen but guards against rounding 114 return e_sq - m_sq 115} 116 117// Integer square root via Newton's method. For x >= 0, returns floor(sqrt(x)). 118func nx_stats_isqrt(x: i64) -> i64 { 119 if x < 0 { return 0 } 120 if x == 0 { return 0 } 121 if x < 4 { return 1 } 122 var g: i64 = x 123 // Initial overestimate: x itself (will converge fast). 124 // Better: g = x >> 1 + 1 for rough first guess. 125 g = (x >> 1) + 1 126 var iter: i64 = 0 127 while iter < 64 { 128 let next_g: i64 = (g + x / g) / 2 129 if next_g >= g { iter = 64 } 130 if next_g < g { 131 g = next_g 132 iter = iter + 1 133 } 134 } 135 return g 136} 137 138func nx_stats_stddev(s: *StreamStats) -> i64 { 139 return nx_stats_isqrt(nx_stats_variance(s)) 140} 141 142func nx_stats_min(s: *StreamStats) -> i64 { 143 return s.min_val 144} 145 146func nx_stats_max(s: *StreamStats) -> i64 { 147 return s.max_val 148} 149 150func nx_stats_count(s: *StreamStats) -> i64 { 151 return s.count 152} 153 154// === typed queries ================================================ 155// 156// All exact under the overflow budget -- envelope conf_ppb = 1e9. 157 158func nx_stats_query_mean(s: *StreamStats) -> *ApproxI64 { 159 let m: i64 = nx_stats_mean(s) 160 return nx_approx_new(m, NX_ENV_ABS, 0, 1000000000, 161 NX_MATURITY_PRODUCTION, 162 NX_ADV_HONEST) 163} 164 165func nx_stats_query_variance(s: *StreamStats) -> *ApproxI64 { 166 let v: i64 = nx_stats_variance(s) 167 return nx_approx_new(v, NX_ENV_ABS, 0, 1000000000, 168 NX_MATURITY_PRODUCTION, 169 NX_ADV_HONEST) 170} 171 172// === merge (Chan 1979 parallel combine, integer variant) ========== 173// 174// Combine two independent streams. All sums add exactly. 175 176func nx_stats_merge(a: *StreamStats, b: *StreamStats) -> *StreamStats { 177 let out: *StreamStats = nx_stats_alloc() 178 out.count = a.count + b.count 179 out.sum = a.sum + b.sum 180 out.sum_sq = a.sum_sq + b.sum_sq 181 if a.has_data == 1 { 182 if b.has_data == 1 { 183 out.min_val = a.min_val 184 if b.min_val < out.min_val { out.min_val = b.min_val } 185 out.max_val = a.max_val 186 if b.max_val > out.max_val { out.max_val = b.max_val } 187 out.has_data = 1 188 } 189 if b.has_data == 0 { 190 out.min_val = a.min_val 191 out.max_val = a.max_val 192 out.has_data = 1 193 } 194 } 195 if a.has_data == 0 { 196 if b.has_data == 1 { 197 out.min_val = b.min_val 198 out.max_val = b.max_val 199 out.has_data = 1 200 } 201 } 202 return out 203} 204 205func nx_stats_memory_bytes(s: *StreamStats) -> i64 { 206 return 56 207}