Arbitrary Price Point Probability (APPP)
// This source code is subject to the terms of the Mozilla Public License 2.0 at https://mozilla.org/MPL/2.0/
// © peacefulLizard50262

//@version=5 
indicator("Arbitrary Price Point Probability", "APPP", overlay = true)
import RicardoSantos/MathSpecialFunctionsGamma/1
arb_price = input.price(0.0, "Price", confirm = true)

gamma(float x) =>
    MathSpecialFunctionsGamma.Gamma(x)

t_dist_pdf(float x, int v) =>
    gamma((v + 1) / 2) / (math.sqrt(v * math.pi) * gamma(v / 2)) * math.pow(1 + x * x / v, -(v + 1) / 2)

rising_factorial(float c, int n) =>
    float product = 1.0
    if n > 0
        for i = 0 to n - 1
            product *= (c + i)
    product

hypergeometric(float a, float b, float c, float z) =>
    float sum = 0.0
    for i = 0 to 2
        sum += rising_factorial(a, i) * rising_factorial(b, i) * math.pow(z, i) / (rising_factorial(c, i) * rising_factorial(i, i))
    sum

t_dist_cdf(float x, int v) =>
    // I hate myself for doing this
    1/2.0 + x * gamma((v+1)/2.0) * (hypergeometric(1/2.0, (v+1)/2.0, 3.0/2, x*x/v) / (math.sqrt(v * math.pi)* gamma(v/2)))

mean_length = input.int(20 , "Mean Length", 2)
stdev_length = input.int(20 , "Deviation Length", 2)
// Calculate p-value given the confidence
p = (1 - 0.99) / 2
// We need to find the value of x such that t_dist_cdf(x, length - 1) = p
// We do this by using a binary search
float lower = -100
float upper = 100
float mid = (lower + upper) / 2
while math.abs(t_dist_cdf(mid, stdev_length - 1) - p) > 0.000001
    if t_dist_cdf(mid, stdev_length - 1) > p
        upper := mid
    else
        lower := mid
    mid := (lower + upper) / 2
float t_value = mid
// t-value for the p-value
if t_value < 0
    t_value := -t_value

mean = ta.ema(open, mean_length)
squared_diff = math.pow(open-mean, 2)
summed_squares = math.sum(squared_diff, stdev_length)
// length - 1 because we're using a sample rather than whole population
mean_squared_error = summed_squares / (stdev_length-1)
rmse = math.sqrt(mean_squared_error)
root_standard_error = math.sqrt(1 + 1/stdev_length + squared_diff/summed_squares)
offset = t_value * rmse * root_standard_error
zscore = (arb_price-mean)/offset

probability = t_dist_pdf(zscore, stdev_length - 1)

price = label.new(bar_index + 3, arb_price, text=str.tostring(math.round(probability*100, 2))+"%", style=label.style_label_down, textcolor = arb_price > open ? color.green : color.red, color=color.new(color.white, 100))
lines = line.new(bar_index + 3, arb_price, bar_index + 4, arb_price, color = arb_price > open ? color.green : color.red)
if barstate.isconfirmed
    label.delete(price)
    line.delete(lines)