sketch_kmeans1d.nx source
↩ module page · 185 lines · 5652 B
1// sketch_kmeans1d.nx -- streaming 1D k-means clustering.
2//
3// Online single-pass clustering via Lloyd's-style centroid updates.
4// For each new value x:
5// 1. Find nearest centroid i = argmin |x - centroid[j]|
6// 2. count[i] += 1
7// 3. centroid[i] += (x - centroid[i]) / count[i] (running mean)
8//
9// Centroids initialized via simple priming: first k distinct values
10// observed seed the centroids; subsequent values cluster against them.
11//
12// USE CASES:
13// - log-message latency clustering (fast / slow / outlier)
14// - sensor value bucketing
15// - online categorical discretization
16// - simple anomaly via distance-to-nearest-centroid
17//
18// COMPLEMENTS:
19// - sketch_histogram: pre-defined uniform bins
20// - sketch_kmeans1d: adaptive bins driven by data density
21//
22// LOSSLESS-LANGUAGE DISCIPLINE: nx_kmeans_query_centroid returns the
23// running-mean centroid value with NX_ENV_REL_STDDEV ~ 1/sqrt(count).
24// Production tier (exact integer running means, no probabilistic error).
25
26import "syscalls.nx"
27import "sketch_types.nx"
28
29const NX_KM_MIN_K: i64 = 2
30const NX_KM_MAX_K: i64 = 1024
31
32struct KMeans1D {
33 centroids: *i64, // k values
34 counts: *i64, // k per-cluster counts
35 k: i64,
36 n_seeded: i64, // # centroids primed so far
37 total_obs: i64,
38}
39
40// === construction =================================================
41
42func nx_kmeans_alloc(k: i64) -> *KMeans1D {
43 if k < NX_KM_MIN_K { return 0 as *KMeans1D }
44 if k > NX_KM_MAX_K { return 0 as *KMeans1D }
45 let raw: *u8 = sys_mmap(40)
46 let m: *KMeans1D = raw as *KMeans1D
47 m.centroids = sys_mmap(k * 8) as *i64
48 m.counts = sys_mmap(k * 8) as *i64
49 var i: i64 = 0
50 while i < k {
51 m.centroids[i] = 0
52 m.counts[i] = 0
53 i = i + 1
54 }
55 m.k = k
56 m.n_seeded = 0
57 m.total_obs = 0
58 return m
59}
60
61// === abs ==========================================================
62
63func nx_kmeans_iabs(x: i64) -> i64 {
64 if x < 0 { return -x }
65 return x
66}
67
68// === nearest centroid =============================================
69
70func nx_kmeans_nearest(m: *KMeans1D, x: i64) -> i64 {
71 if m.n_seeded == 0 { return -1 }
72 var best: i64 = 0
73 var best_dist: i64 = nx_kmeans_iabs(x - m.centroids[0])
74 var i: i64 = 1
75 while i < m.n_seeded {
76 let d: i64 = nx_kmeans_iabs(x - m.centroids[i])
77 if d < best_dist {
78 best = i
79 best_dist = d
80 }
81 i = i + 1
82 }
83 return best
84}
85
86// === observe =====================================================
87
88func nx_kmeans_observe(m: *KMeans1D, x: i64) -> i64 {
89 m.total_obs = m.total_obs + 1
90 if m.n_seeded < m.k {
91 // Seed phase: each new value seeds a new centroid if it's
92 // distinct from existing ones (or if we still have capacity).
93 // Simple policy: every distinct value up to k seeds.
94 var found: i64 = 0
95 var i: i64 = 0
96 while i < m.n_seeded {
97 if m.centroids[i] == x { found = 1 }
98 i = i + 1
99 }
100 if found == 0 {
101 m.centroids[m.n_seeded] = x
102 m.counts[m.n_seeded] = 1
103 m.n_seeded = m.n_seeded + 1
104 return 0
105 }
106 // x matches an existing seed -- update that cluster.
107 }
108 // Standard Lloyd's update.
109 let idx: i64 = nx_kmeans_nearest(m, x)
110 if idx < 0 { return -1 }
111 m.counts[idx] = m.counts[idx] + 1
112 // centroid += (x - centroid) / count (integer running mean)
113 let cnt: i64 = m.counts[idx]
114 let diff: i64 = x - m.centroids[idx]
115 m.centroids[idx] = m.centroids[idx] + diff / cnt
116 return 0
117}
118
119// === queries =====================================================
120
121func nx_kmeans_predict(m: *KMeans1D, x: i64) -> i64 {
122 return nx_kmeans_nearest(m, x)
123}
124
125func nx_kmeans_centroid(m: *KMeans1D, i: i64) -> i64 {
126 if i < 0 { return 0 }
127 if i >= m.n_seeded { return 0 }
128 return m.centroids[i]
129}
130
131func nx_kmeans_count(m: *KMeans1D, i: i64) -> i64 {
132 if i < 0 { return 0 }
133 if i >= m.n_seeded { return 0 }
134 return m.counts[i]
135}
136
137func nx_kmeans_n_clusters(m: *KMeans1D) -> i64 {
138 return m.n_seeded
139}
140
141func nx_kmeans_total(m: *KMeans1D) -> i64 {
142 return m.total_obs
143}
144
145// === isqrt ========================================================
146
147func nx_kmeans_isqrt(x: i64) -> i64 {
148 if x < 0 { return 0 }
149 if x == 0 { return 0 }
150 if x < 4 { return 1 }
151 var g: i64 = (x >> 1) + 1
152 var iter: i64 = 0
153 while iter < 64 {
154 let next_g: i64 = (g + x / g) / 2
155 if next_g >= g { iter = 64 }
156 if next_g < g {
157 g = next_g
158 iter = iter + 1
159 }
160 }
161 return g
162}
163
164// === typed envelope ==============================================
165//
166// Centroid stddev ~ within-cluster-std / sqrt(count_i). We declare
167// 1/sqrt(count) as the rel_stddev bound (conservative for narrow clusters).
168
169func nx_kmeans_query_centroid(m: *KMeans1D, i: i64) -> *ApproxI64 {
170 let c: i64 = nx_kmeans_centroid(m, i)
171 let cnt: i64 = nx_kmeans_count(m, i)
172 var stderr_ppb: i64 = 1000000000
173 if cnt > 0 {
174 let isq: i64 = nx_kmeans_isqrt(cnt)
175 if isq > 0 { stderr_ppb = 1000000000 / isq }
176 }
177 return nx_approx_new(c, NX_ENV_REL_STDDEV, stderr_ppb,
178 682700000,
179 NX_MATURITY_PRODUCTION,
180 NX_ADV_HONEST)
181}
182
183func nx_kmeans_memory_bytes(m: *KMeans1D) -> i64 {
184 return 40 + m.k * 16
185}