# -*- coding: utf-8 -*-
"""
Created on Fri Aug 12 13:20:49 2022

@author: Rene
"""
import numpy as np
from scipy import optimize
import matplotlib.pyplot as plt
import pylab as p

def CDF_T_k(t,k,n):
    num_values = np.arange(t-k+1,t+1,dtype=np.float32)
    denom_values = np.arange(n-k+1,n+1,dtype=np.float32)
    numerator = np.product(num_values)
    denominator = np.product(denom_values)
    value = numerator/denominator
    return(value)

def solve_critical_value(alpha,k,n):
    
    def f(t):
        value = alpha - CDF_T_k(t,k,n)
        return(value)
    
    root = optimize.brentq(f, 1, n)
    return(root)

n = 350
k = 5
t = 61
cdf = CDF_T_k(t,k,n)
alpha = 0.05
t_alpha = solve_critical_value(alpha,k,n)

cdf_t = np.zeros(n+1)
for i in np.arange(1,n+1):
  cdf_t[i] = CDF_T_k(i,k,n)

data = [61,19,56,24,16]
t_5 = max(data)
p_5 = CDF_T_k(t_5,k,n)

# Plotting the analysis results    
LB = 0
UB = 350
off_set_x = 5
off_set_y = 0.1
max_y = 1
    
# Setting the attributes for the Figure
Figure = plt.figure()     
plt.rc('font', family='serif', size='10')
plt.figure(figsize=(9, 4.5))
Figure.subplots_adjust(hspace=0.4, wspace=0.4)
x_lims = (LB-off_set_x,UB+off_set_x)
y_lims = (0-off_set_y,max_y+off_set_y)
plt.axvline(0, lw=2, ls = '-',color = 'lightgray',alpha=0.5)
plt.axhline(0, lw=2, ls = '-',color = 'lightgray',alpha=0.5)
plt.xlim(x_lims)
plt.ylim(y_lims)
plt.xlabel('t')
plt.ylabel('$Pr(T_5 \leq t|N=350$')
plt.plot(np.arange(0,n+1),cdf_t,ls='-',lw=0.5,color='red')

plt.axhline(1, lw=1, ls = ':',color = 'black',alpha=0.5)
plt.axvline(350, lw=1, ls = ':',color = 'black',alpha=0.5)

plt.axvline(t_5, lw=1, ls = '-.',color = 'red',alpha=0.5)
plt.axvline(t_alpha, lw=1, ls = '--',color = 'blue',alpha=0.5)
plt.axhline(alpha, lw=1, ls = '--',color = 'blue',alpha=0.5)


p.arrow( 10, 0.5, t_alpha-20, 0, head_width=0.05, head_length=10,color='green')
p.arrow( t_alpha-10, 0.5, -(t_alpha-20), 0, head_width=0.05, head_length=10,color='green')
plt.text(t_alpha/2,0.525,'Critical Region', horizontalalignment = 'center', verticalalignment = 'bottom', color = 'green')

text_str = '$t_5 = $' +  f'{t_5:2.0f}' 
plt.text(t_5+5,1,text_str, horizontalalignment = 'left', verticalalignment = 'top', color = 'red')

text_str = r'$t_{\alpha} = $' +  f'{t_alpha:2.2f}' 
plt.text(t_alpha+5,-0.05,text_str, horizontalalignment = 'left', verticalalignment = 'center', color = 'blue')

text_str = r'Type 1 Error $ = \alpha = $'+f'{alpha:2.2f}' 
plt.text(10,alpha,text_str, horizontalalignment = 'left', verticalalignment = 'bottom', color = 'blue')

text_str = f'{alpha:2.2f}' 
plt.text(-5,alpha,text_str, horizontalalignment = 'right', verticalalignment = 'center', color = 'blue')

text_str = r"$H_0: N = $" + f'{n:2.0f}'
plt.text(0,0.9,text_str, horizontalalignment = 'left', verticalalignment = 'center', color = 'green')

text_str = r"$H_1: N < $" + f'{n:2.0f}'
plt.text(0,0.8,text_str, horizontalalignment = 'left', verticalalignment = 'center', color = 'green')

text_str = r"German Tank Example: $T_5 = max(X_1,\ldots,X_5)$"
plt.suptitle(text_str,size='14')

plt.savefig('German_Tank_Example_N_350.png', dpi=300)