mirror of
https://github.com/mihakralj/QuanTAlib.git
synced 2026-07-27 17:27:43 +00:00
128 lines
4.6 KiB
Plaintext
128 lines
4.6 KiB
Plaintext
// Licensed under the Apache License, Version 2.0
|
|
// © mihakralj
|
|
//@version=6
|
|
indicator("F-Distribution CDF (FDIST)", "FDIST", overlay=false, precision=6)
|
|
|
|
//@function Log-gamma via Lanczos approximation (g=7, 9 coefficients)
|
|
//@param z Input value (z > 0)
|
|
//@returns ln(Gamma(z))
|
|
lnGamma(simple float z) =>
|
|
float g = 7.0
|
|
array<float> c = array.from(
|
|
0.99999999999980993,
|
|
676.5203681218851,
|
|
-1259.1392167224028,
|
|
771.32342877765313,
|
|
-176.61502916214059,
|
|
12.507343278686905,
|
|
-0.13857109526572012,
|
|
9.9843695780195716e-6,
|
|
1.5056327351493116e-7)
|
|
float zz = z - 1.0
|
|
float x = array.get(c, 0)
|
|
for i = 1 to 8
|
|
x += array.get(c, i) / (zz + i)
|
|
float t = zz + g + 0.5
|
|
0.5 * math.log(2.0 * math.pi) + (zz + 0.5) * math.log(t) - t + math.log(x)
|
|
|
|
//@function Regularized incomplete beta function I_x(a,b) via continued fraction (Lentz)
|
|
//@param x Evaluation point (0 <= x <= 1)
|
|
//@param a Shape parameter alpha (a > 0)
|
|
//@param b Shape parameter beta (b > 0)
|
|
//@returns CDF value P(X <= x) for Beta(a,b)
|
|
betaReg(series float x, simple float a, simple float b) =>
|
|
if x <= 0.0
|
|
0.0
|
|
else if x >= 1.0
|
|
1.0
|
|
else
|
|
float lnPfx = a * math.log(x) + b * math.log(1.0 - x) - math.log(a) - lnGamma(a) - lnGamma(b) + lnGamma(a + b)
|
|
float front = math.exp(lnPfx)
|
|
bool flip = x > (a + 1.0) / (a + b + 2.0)
|
|
float xx = flip ? 1.0 - x : x
|
|
float aa = flip ? b : a
|
|
float bb = flip ? a : b
|
|
float lnPfx2 = aa * math.log(xx) + bb * math.log(1.0 - xx) - math.log(aa) - lnGamma(aa) - lnGamma(bb) + lnGamma(aa + bb)
|
|
float front2 = math.exp(lnPfx2)
|
|
float TINY = 1e-30
|
|
float EPS = 1e-10
|
|
int MAXITER = 200
|
|
float f = TINY
|
|
float C = TINY
|
|
float D = 0.0
|
|
float delta = 0.0
|
|
for m = 0 to MAXITER - 1
|
|
float d_val = 0.0
|
|
int mm = m / 2
|
|
if m == 0
|
|
d_val := 1.0
|
|
else if m % 2 == 0
|
|
float mf = float(mm)
|
|
d_val := mf * (bb - mf) * xx / ((aa + 2.0 * mf - 1.0) * (aa + 2.0 * mf))
|
|
else
|
|
float mf = float(mm) + 1.0
|
|
d_val := -(aa + mf - 1.0) * (aa + bb + mf - 1.0) * xx / ((aa + 2.0 * mf - 2.0) * (aa + 2.0 * mf - 1.0))
|
|
D := 1.0 + d_val * D
|
|
if math.abs(D) < TINY
|
|
D := TINY
|
|
D := 1.0 / D
|
|
C := 1.0 + d_val / C
|
|
if math.abs(C) < TINY
|
|
C := TINY
|
|
delta := C * D
|
|
f *= delta
|
|
if math.abs(delta - 1.0) < EPS
|
|
break
|
|
float result = front2 * f
|
|
flip ? 1.0 - result : result
|
|
|
|
//@function Computes F-Distribution CDF for a normalized price series
|
|
//@param source Series to transform
|
|
//@param period Lookback period for min-max normalization
|
|
//@param d1 Numerator degrees of freedom (d1 > 0)
|
|
//@param d2 Denominator degrees of freedom (d2 > 0)
|
|
//@returns CDF value in [0,1]: P(F <= x) via regularized incomplete beta
|
|
//@optimized Lentz continued fraction converges in ~10-20 iterations
|
|
fdist(series float source, simple int period, simple float d1, simple float d2) =>
|
|
if period <= 0
|
|
runtime.error("Period must be greater than 0")
|
|
if d1 <= 0.0
|
|
runtime.error("d1 must be greater than 0")
|
|
if d2 <= 0.0
|
|
runtime.error("d2 must be greater than 0")
|
|
|
|
float minVal = source
|
|
float maxVal = source
|
|
for i = 1 to period - 1
|
|
float v = source[i]
|
|
if not na(v)
|
|
if v < minVal
|
|
minVal := v
|
|
if v > maxVal
|
|
maxVal := v
|
|
|
|
float range = maxVal - minVal
|
|
float x = range > 0.0 ? (source - minVal) / range : 0.5
|
|
|
|
float safeX = math.max(0.0, x)
|
|
float t = d1 * safeX / (d1 * safeX + d2)
|
|
|
|
betaReg(t, d1 / 2.0, d2 / 2.0)
|
|
|
|
// ---------- Main loop ----------
|
|
|
|
// Inputs
|
|
i_source = input.source(close, "Source")
|
|
i_period = input.int(50, "Lookback Period", minval=2, maxval=5000, tooltip="Min-max normalization window")
|
|
i_d1 = input.float(5.0, "d1 (numerator df)", minval=0.1, step=1.0, tooltip="Numerator degrees of freedom")
|
|
i_d2 = input.float(5.0, "d2 (denominator df)", minval=0.1, step=1.0, tooltip="Denominator degrees of freedom")
|
|
|
|
// Calculation
|
|
float result = fdist(i_source, i_period, i_d1, i_d2)
|
|
|
|
// Plot
|
|
plot(result, "FDIST", color=color.yellow, linewidth=2)
|
|
hline(0.5, "Midline", color=color.gray, linestyle=hline.style_dotted)
|
|
hline(0.95, "Upper", color=color.red, linestyle=hline.style_dashed)
|
|
hline(0.05, "Lower", color=color.green, linestyle=hline.style_dashed)
|