For a single time series, the distribution is the same distribution as in the two sample Kolmogorov Smirnov Test, namely the distribution of the maximal value of the absolute values of a Brownian bridge. It is computated as follows (van Mulbregt, 2018):
For \(t_n(x) < 1\):
$$
P(t_n(X) \le t_n(x)) =
\sqrt{2 * \pi} / t_n(x) * t (1 + t^8(1 + t^{16}(1 + t^{24}(1 + ...))))$$
up to \(t^{8 * k_{max}},
k_{max} = \lfloor \sqrt{2 - \log(tol)}\rfloor\) where \(t = \exp(-\pi^2 / (8*x^2))\)
else:
$$
P(t_n(X) \le t_n(x)) = 2 * \sum_{k = 1}^{\infty} (-1)^{k - 1} * \exp(-2*k^2*x^2)$$
until
\(|2 * (-1)^{k - 1} * \exp(-2*k^2*x^2) - 2 * (-1)^{(k-1) - 1} * \exp(-2*(k-1)^2*x^2)| \le tol.
\)
In case of multiple time series, the distribution equals that of the maximum of an h dimensional squared Bessel bridge. It can be computed by (Kiefer, 1959):
$$P(t_n(X) \le t_n(x)) =
(4 / (\Gamma(h / 2) 2^{h / 2} t_n^h)) \sum_{i = 1}^{\infty} ( ((\gamma_{(h - 2)/2, n})^{h - 2} \exp(-(\gamma_{(h - 2)/2, n})^2 / (2t_n^2))) / (J_{h/2}(\gamma_{(h - 2)/2, n}))^2 )$$
where \(J_h\) is the Bessel function of first kind and h-th order, \(\Gamma\) is the gamma function and \(\gamma_{h, n}\) denotes the n-th zero of \(J_h\).