sketch_ddsketch_vs_tdigest_tail_bench.nx source
↩ module page · 128 lines · 4923 B
1// sketch_ddsketch_vs_tdigest_tail_bench.nx -- multiplicative-error tail.
2//
3// CLAIM TO VALIDATE:
4// DDSketch (Masson-Rim-Lee 2019) provides MULTIPLICATIVE relative
5// error on every quantile: |est/true - 1| <= 1/bpo. T-Digest
6// provides empirical tail-tight ADDITIVE error. On a HEAVY-TAIL
7// stream where the tail values span multiple orders of magnitude,
8// DDSketch's multiplicative guarantee wins on extreme quantiles.
9//
10// WORKLOAD:
11// 500 small values (in [1, 100]) + 50 huge values (geometric in
12// [1000, 100000]) -- the heavy tail. Total N=550.
13// True p99 ~ 60000 (top 1% are the huge values).
14// True median ~ 50 (middle of small values).
15//
16// MEASUREMENT:
17// DDSketch bpo=64 -> rel err 1/64 = 1.56% per quantile.
18// T-Digest delta=100 -> ~1% global rank error, but at p99 of a
19// heavy-tail the absolute value error scales with quantile value.
20//
21// HARD-WIN GATE:
22// On p99 with truth ~ 60000:
23// - DDSketch should be within ~2-3% (1000-1500 absolute)
24// - T-Digest may be much further off due to centroid quantization
25// at the tail
26// ACCURACY axis BEATS by >5%.
27
28import "syscalls.nx"
29import "sketch_ddsketch.nx"
30import "sketch_tdigest.nx"
31import "sketch_comparator.nx"
32import "sketch_types.nx"
33
34func iabs_d(x: i64) -> i64 {
35 if x < 0 { return -x }
36 return x
37}
38
39func main() -> i64 {
40 // MATCHED-MEMORY rerun:
41 // DDSketch bpo=64 -> 63 * 64 * 8 = 32256 bytes allocated
42 // T-Digest delta=500 (substrate max NX_TD_DELTA_MAX) -> 20848 bytes
43 // T-Digest's substrate-cap is the binding constraint here; raising
44 // NX_TD_MAX_CENTROIDS / NX_TD_DELTA_MAX is queued as a bits-up
45 // scaling task. At delta=500 T-Digest produces ~5x more centroids
46 // than at delta=100; gives it its best honest shot at the tail.
47 let dd: *DDSketch = nx_dd_alloc(64)
48 // delta=1000 -- substrate cap raised this session (was 500). At
49 // delta=1000 T-Digest produces ~10x more centroids than delta=100.
50 let td: *TDigest = nx_tdigest_alloc(1000)
51 if dd == (0 as *DDSketch) { return __syscall(93, 1, 0, 0, 0, 0, 0) }
52 if td == (0 as *TDigest) { return __syscall(93, 2, 0, 0, 0, 0, 0) }
53
54 // ---- Stream: 500 small values in [1, 100] ----
55 var i: i64 = 0
56 while i < 500 {
57 // Modulo so we get a uniform spread in [1, 100]
58 let v: i64 = (i % 100) + 1
59 nx_dd_add(dd, v)
60 nx_tdigest_add(td, v)
61 i = i + 1
62 }
63 // ---- 50 huge values in [1000, 100000], geometric-ish spacing ----
64 i = 0
65 let huge_count: i64 = 50
66 while i < huge_count {
67 // Geometric-ish: v = 1000 * 2^(i/10)
68 var pow: i64 = 1
69 var k: i64 = 0
70 let exp: i64 = i / 10
71 while k < exp {
72 pow = pow * 2
73 k = k + 1
74 }
75 let v: i64 = 1000 * pow + (i * 17) % 500
76 nx_dd_add(dd, v)
77 nx_tdigest_add(td, v)
78 i = i + 1
79 }
80
81 // ---- Query p99: should fall in the upper part of the huge values ----
82 let dd_p99: i64 = nx_dd_quantile(dd, 990)
83 let td_p99: i64 = nx_tdigest_quantile(td, 990)
84
85 // ---- Ground truth p99 ----
86 // 550 total values; p99 index = ceil(0.99 * 550) = 545.
87 // Sorted: 500 small values in [1, 100], then 50 huge values sorted ascending.
88 // p99 = 545th value = 45th huge value (1-based) = huge_index 44.
89 // pow at i=44: 2^4 = 16, plus base 1000 -> 16000, plus (44*17 mod 500) = 248 -> 16248.
90 // Actually compute exactly: at i=44, exp=4, pow=16, v=16000 + 248 = 16248.
91 // For p99 we want the 545th largest = the 45th huge value (1-indexed) = i=44 0-indexed.
92 // Truth p99 = 16248.
93 let truth_p99: i64 = 16248
94
95 // ---- ACCURACY axis: closer to truth wins ----
96 let acc: *ComparisonResult = nx_cmp_accuracy(dd_p99, td_p99, truth_p99, 50000)
97
98 // ---- HARD-WIN GATE ----
99 if acc.verdict != NX_CMP_VERDICT_BEATS {
100 return __syscall(93, 10, 0, 0, 0, 0, 0)
101 }
102 if acc.delta_ppm < 50000 {
103 return __syscall(93, 11, 0, 0, 0, 0, 0)
104 }
105
106 // ---- DDSketch tail guarantee: |dd_p99/truth - 1| <= 1/64 = 1.56% ----
107 // i.e., |dd_p99 - truth| <= truth/64 ≈ 254.
108 let dd_err: i64 = iabs_d(dd_p99 - truth_p99)
109 let bpo_bound: i64 = truth_p99 / 64
110 if dd_err > bpo_bound * 2 { // 2x slack for adjacent bucket choice
111 return __syscall(93, 20, 0, 0, 0, 0, 0)
112 }
113
114 // ---- Sanity: both correctly identify p50 in [1, 100] range ----
115 let dd_p50: i64 = nx_dd_quantile(dd, 500)
116 let td_p50: i64 = nx_tdigest_quantile(td, 500)
117 if dd_p50 < 1 { return __syscall(93, 30, 0, 0, 0, 0, 0) }
118 if dd_p50 > 100 { return __syscall(93, 31, 0, 0, 0, 0, 0) }
119 if td_p50 < 1 { return __syscall(93, 32, 0, 0, 0, 0, 0) }
120 if td_p50 > 100 { return __syscall(93, 33, 0, 0, 0, 0, 0) }
121
122 // ---- Total counts agree ----
123 if nx_dd_total(dd) != 550 {
124 return __syscall(93, 40, 0, 0, 0, 0, 0)
125 }
126
127 return 0
128}