Setting the file. One moment.
Ks2 · Signals Scout Anomaly Detection · PostHog/skills · Skills Docs
Back to the top of the page scripts/ks2.py
scripts/ ks2.py
Python · 130 lines · 5 KB
15
16 The KS statistic D = max|F_a(x) - F_b(x)|; the p-value is the asymptotic Kolmogorov
17 distribution with the Stephens small-sample correction (matches scipy's ks_2samp
18 asymptotic mode closely for n,m >~ 30).
19 """
20
21 from __future__ import annotations
22
23 import sys
24 import json
25 import math
26
27
28 def _ks_pvalue (d: float , n: float , m: float ) -> float :
29 if n <= 0 or m <= 0 or d <= 0 :
30 return 1.0
31 en = math.sqrt(n * m / (n + m))
32 t = (en + 0.12 + 0.11 / en) * d # Stephens correction
33 if t < 1e-12 :
34 return 1.0
35 s = 0.0
36 for k in range ( 1 , 101 ):
37 term = 2.0 * ( - 1 ) ** (k - 1 ) * math.exp( - 2.0 * k * k * t * t)
38 s += term
39 if abs (term) < 1e-10 :
40 break
41 return max ( 0.0 , min ( 1.0 , s))
42
43
44 def ks_2samp (a: list[ float ], b: list[ float ]) -> tuple[ float , float ]:
45 a = sorted (a)
46 b = sorted (b)
47 n, m = len (a), len (b)
48 if n == 0 or m == 0 :
49 return 0.0 , 1.0
50 i = j = 0
51 d = 0.0
52 while i < n and j < m:
53 x = a[i] if a[i] <= b[j] else b[j]
54 while i < n and a[i] <= x:
55 i += 1
56 while j < m and b[j] <= x:
57 j += 1
58 d = max (d, abs (i / n - j / m))
59 return d, _ks_pvalue(d, n, m)
60
61
62 def ks_2samp_binned (a_hist: list[list[ float ]], b_hist: list[list[ float ]]) -> tuple[ float , float ]:
63 """KS on two empirical CDFs given as (value, count) bins — the cheap-payload path."""
64 ca: dict[ float , float ] = {}
65 cb: dict[ float , float ] = {}
66 for v, c in a_hist:
67 ca[v] = ca.get(v, 0.0 ) + c
68 for v, c in b_hist:
69 cb[v] = cb.get(v, 0.0 ) + c
70 na = sum (ca.values())
71 nb = sum (cb.values())
72 if na == 0 or nb == 0 :
73 return 0.0 , 1.0
74 fa = fb = d = 0.0
75 for v in sorted ( set (ca) | set (cb)):
76 fa += ca.get(v, 0.0 ) / na
77 fb += cb.get(v, 0.0 ) / nb
78 d = max (d, abs (fa - fb))
79 return d, _ks_pvalue(d, na, nb)
80
81
82 def changepoint (series: list[ float ], min_seg: int = 8 ) -> dict :
83 """Find the split index whose left/right value distributions differ most (max D).
84
85 The sweep picks the split that maximizes D over many candidates, so the winning
86 `p` is a **scan minimum** — biased low by multiple comparisons and NOT a
87 single-hypothesis p-value. Use `p_adj` (Bonferroni over `tests`) as the calibrated
88 figure, and confirm the chosen split with a direct two-sample KS on
89 seasonality-matched windows before treating it as emit evidence.
90 """
91 n = len (series)
92 if n < 2 * min_seg:
93 return { "changepoint" : None , "reason" : f "need >= { 2 * min_seg } points, got { n } " }
94 best = { "index" : None , "d" : 0.0 , "p" : 1.0 }
95 tests = 0
96 for c in range (min_seg, n - min_seg + 1 ):
97 tests += 1
98 d, p = ks_2samp(series[:c], series[c:])
99 if d > best[ "d" ]:
100 best = { "index" : c, "d" : d, "p" : p}
101 return {
102 "changepoint" : best[ "index" ],
103 "d" : round (best[ "d" ], 4 ),
104 "p" : best[ "p" ], # scan minimum — uncorrected; see p_adj
105 "p_adj" : min ( 1.0 , best[ "p" ] * tests), # Bonferroni over the scan
106 "tests" : tests,
107 "n" : n,
108 }
109
110
111 def main () -> None :
112 req = json.load(sys.stdin)
113 mode = req.get( "mode" , "two_sample" )
114 if mode == "changepoint" :
115 out = changepoint(req[ "series" ], min_seg = req.get( "min_seg" , 8 ))
116 elif "a_hist" in req or "b_hist" in req:
117 if "a_hist" not in req or "b_hist" not in req:
118 out = { "error" : "binned mode needs both a_hist and b_hist" }
119 else :
120 d, p = ks_2samp_binned(req[ "a_hist" ], req[ "b_hist" ])
121 out = { "d" : round (d, 4 ), "p" : p, "na" : sum (c for _, c in req[ "a_hist" ]), "nb" : sum (c for _, c in req[ "b_hist" ])}
122 else :
123 d, p = ks_2samp(req[ "a" ], req[ "b" ])
124 out = { "d" : round (d, 4 ), "p" : p, "n" : len (req[ "a" ]), "m" : len (req[ "b" ])}
125 json.dump(out, sys.stdout)
126 sys.stdout.write( " \n " )
127
128
129 if __name__ == "__main__" :
130 main()