everythingfades/Math-Stats-AP
1
1import gradio as gr2import os3os.system("pip install scipy")4def birthdayCheck(people):5 p = 1;6 for i in range(people):7 p = p*(365-i)/3658 return 1-p9import random10import matplotlib.pyplot as plt11def alter(x,y):12 return x, y+313def getSame(birthdays):14 cnt = 015 x = list()16 y = list()17 z = list()18 for i in range(len(birthdays)):19 for j in range(i+1,len(birthdays)):20 if birthdays[i] == birthdays[j]:21 cnt += 122 x.append(i)23 y.append(j)24 z.append(birthdays[i])25 return cnt,x,y,z26 27def simulation(people):28 birthdays = list()29 peoples = list()30 for i in range(people):31 birthdays.append(random.randint(1,365))32 peoples.append(i)33 return birthdays,peoples34def tran(z):35 a = [1,31,59,90,120,151,181,212,243,273,304,334]36 str = ''37 for j in z:38 for i in range(len(a)):39 if j < a[i]:40 str += "{}/{}".format(i,j-a[i-1]+1)41 str += " | "42 break43 return str if str != '' else "None"44#birthdays,peoples = simulation(30)45import matplotlib.pyplot as plt46def alter(x,y):47 return x, y+348def getSame(birthdays):49 cnt = 050 x = list()51 y = list()52 z = list()53 for i in range(len(birthdays)):54 for j in range(i+1,len(birthdays)):55 if birthdays[i] == birthdays[j]:56 cnt += 157 x.append(i)58 y.append(j)59 z.append(birthdays[i])60 return cnt,x,y,z61import gradio as gr62 63def getSame_gradio(people):64 birthdays,peoples = simulation(people)65 getSame(birthdays)66 cnt = 067 x = list()68 y = list()69 z = list()70 for i in range(len(birthdays)):71 for j in range(i+1,len(birthdays)):72 if birthdays[i] == birthdays[j]:73 cnt += 174 x.append(i)75 y.append(j)76 z.append(birthdays[i])77 temp = ""78 for i in range(len(x)):79 temp += ("No.{}and No.{} have same birthdays".format(x[i]+1,y[i]+1))80 import matplotlib.pyplot as plt81 import matplotlib as mpl82 from matplotlib.patches import Ellipse83 fig, ax = plt.subplots(figsize=(16,10), dpi= 80)84 ax.hlines(y=[31,59,90,120,151,181,212,243,273,304,334,365], xmin=0, xmax=len(peoples)+1, color='gray', alpha=0.7, linewidth=1, linestyles='dashdot')85 ax.scatter(y=birthdays, x=peoples, s=75, color='black', alpha=0.7)86 87 mpl.rcParams['font.sans-serif']=['SimHei']88 89 # Title, Label, Ticks and Ylim90 ax.set_title('birthdays:{}'.format(len(peoples)), fontdict={'size':22})91 ax.set_xlabel('people No.')92 ax.set_ylabel('date')93 ax.set_yticklabels = ['Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep','Oct','Nov','Dec']94 ax.set_yticks([31,59,90,120,151,181,212,243,273,304,334,365])95 for i in range(len(birthdays)):96 textx,texty = alter(peoples[i],birthdays[i])97 ax.text(textx,texty,birthdays[i])98 if cnt != 0:99 for i in range(len(x)):100 plt.gcf().gca().add_artist(Ellipse(xy = (x[i],z[i]), width = 1, height = 17, fill = False, color = "red"))101 plt.gcf().gca().add_artist(Ellipse(xy = (y[i],z[i]), width = 1, height = 17, fill = False, color = "red"))102 return tran(birthdays), "we have {} pairs of same birthdays".format(cnt), "the repeated date is" + tran(z), temp, plt.gcf()103def multi_simulation(times,people):104 cnt_list = list()105 cnt_check = list()106 for i in range(times):107 cnt = one_simulation(people)108 cnt_list.append(cnt)109 cnt_check.append(0 if cnt ==0 else 1)110 return cnt_list, cnt_check111#demo.launch()112 113def one_simulation(people):114 birthdays,peoples = simulation(people)115 cnt,x,y,z = getSame(birthdays)116 return cnt117def draw_one_sim_one_zero_gradio(times, people):118 fig, ax = plt.subplots(figsize=(16,10), dpi= 80)119 cnt_list, cnt_check = multi_simulation(times,people)120 all_colors = list(plt.cm.colors.cnames.keys())121 random.seed(100)122 c = random.choices(all_colors, k=len(list(set(cnt_list))))123 # Plot Bars124 ax.bar(list(set(cnt_check)), find_class(cnt_check), color=c, width=.5)125 #for i in range(len(list(set(cnt_check)))):126 # plt.annotate(float(find_class(cnt_check)[i]), 127 # xy=(list(set(cnt_check))[i],find_class(cnt_check)[i] + 4),#箭头末端位置128 # 129 # xytext=(list(set(cnt_check))[i],find_class(cnt_check)[i] + 4),#文本起始位置130 # 131 # #箭头属性设置132 # arrowprops=dict(facecolor='#74C476', 133 # shrink=1,#箭头的收缩比134 # alpha=0.6,135 # width=7,#箭身宽136 # headwidth=40,#箭头宽137 # hatch='--',#填充形状138 # frac=0.8,#身与头比139 # #其它参考matplotlib.patches.Polygon中任何参数140 # ),141 # )142 return fig,float(find_class(cnt_check)[0]),float(find_class(cnt_check)[1]),float(find_class(cnt_check)[1]/times),birthdayCheck(people)143import matplotlib144import scipy.stats as stats145import numpy as np146def multi_multi_simulation_gradio(people, times_for_p, iter_times):147 fig, ax = plt.subplots(figsize=(16,10), dpi= 80)148 results_30,a,b = multi_fixed_simulation(people, times_for_p, iter_times)149 mu = np.mean(results_30)150 std = np.std(results_30)151 Bin = 80152 n, bins, patches = ax.hist(results_30, bins = Bin)153 y = stats.norm.pdf(bins,mu,std)154 #y = matplotlib.mlab.normpdf(bins,mu,std)155 print(mu,std)156 plt.xlim((birthdayCheck(people)-0.08, birthdayCheck(people)+0.08)) 157 ax.plot(bins,y,'r--')158 return fig159 160def find_class(list1):161 classes = list()162 for i in list(set(list1)):163 temp_cnt= 0164 for j in list1:165 if i==j:166 temp_cnt += 1167 classes.append(temp_cnt)168 if len(classes) == 1:169 classes.append(classes[0])170 classes[0] = 0171 return classes172def multi_fixed_simulation(people,times, iter_times):173 results = list()174 for i in range(iter_times):175 results.append(get_ratio(people,times))176 if i % 100 == 0:177 print(i)178 return results,people, times179def get_ratio(people, times):180 cnt_list1, cnt_check1 = multi_simulation(times,people)181 return(find_class(cnt_check1)[1]/len(cnt_check1))182import math183import numpy as np184import scipy185def random_gen_data(population,p):186 results = np.zeros(int(population))187 for i in range(results.size):188 if np.random.rand(1)[0] <= p:189 results[i] = 1190 return results191def get_ratio1(results):192 #print(results.size,results.sum(),"#########")193 return results.sum()/results.size194def simulation(p,SRS_size,population,times):195 results = np.zeros(times)196 for i in range(times):197 sample1 = np.array(random.sample(list(population),int(SRS_size)))198 results[i] = get_ratio1(sample1)199 return results#,population200def get_CI(p,alpha,sample_size):#size = (SRS_size)201 sigma = math.sqrt(p*(1-p)/sample_size)202 z = scipy.stats.norm.ppf(1-(1-alpha)/2)203 return(sigma*z)204def get_CI_new(alpha,sample):#size = (SRS_size)205 p_hat = sample206 sigma = np.sqrt(p_hat*(1-p_hat)/sample.shape[0])207 z = scipy.stats.norm.ppf(1-(1-alpha)/2)208 return(sigma*z)209def draw(p,alpha,SRS_size,population_size,times):210 population = random_gen_data(population_size,p)211 #print(population)212 sample = simulation(p,SRS_size,population,times)213 CI = get_CI_new(alpha,sample)214 #print(CI)215 all_colors = list(plt.cm.colors.cnames.keys())216 random.seed(42)217 c = random.choices(all_colors, k=population_size)218 219 220 fig,ax=plt.subplot_mosaic([['A','A','A'],221 ['C','C','C'],222 ['B','B','B'],223 ['B','B','B'],224 ['B','B','B']],figsize=(16, 10), constrained_layout=True)225 ax['A'].set_xlim(0, 1)226 #ax['B'].set_xlim(sample.min()-CI, sample.max()+CI)227 ax['B'].set_xlim(0,1)228 #ax['C'].set_xlim(sample.min()-CI, sample.max()+CI)229 ax['C'].set_xlim(0,1)230 mu = population.mean()231 std = population.std()232 n, bins, patches = ax['A'].hist(population,bins = 40)233 y = normfun(bins,mu,std)234 ax['A'].vlines(p,ax['A'].get_ylim()[0],ax['A'].get_ylim()[1],linestyles ="solid", colors ="k")235 ax['C'].plot(bins,y,'r--',label = "population dis")236 237 mu = sample.mean()238 std = sample.std()239 n, bins, patches = ax['C'].hist(sample,bins = 40)240 y = normfun(bins,mu,std)241 ax['C'].plot(bins,y,':',label = "sample dis")242 ax['C'].vlines(p,ax['C'].get_ylim()[0],ax['C'].get_ylim()[1],linestyles ="solid", colors ="k")243 ax['C'].vlines(p-CI,ax['C'].get_ylim()[0],ax['C'].get_ylim()[1],linestyles ="dashed", colors ="k")244 ax['C'].vlines(p+CI,ax['C'].get_ylim()[0],ax['C'].get_ylim()[1],linestyles ="dashed", colors ="k")245 ax['C'].legend()246 247 ax['B'].hlines(y = np.arange(0,times),xmin = sample-CI, xmax = sample+CI,color = "gray")248 ax['C'].vlines(p,ax['C'].get_ylim()[0],ax['C'].get_ylim()[1],linestyles ="solid", colors ="k")249 250 ax['B'].scatter(sample,np.arange(0,times),c = "blue")251 out = np.zeros(sample.shape[0])252 ax['B'].vlines(p,ax['B'].get_ylim()[0],ax['B'].get_ylim()[1],linestyles ="solid", colors ="k")253 ax['B'].vlines(p-CI,ax['B'].get_ylim()[0],ax['B'].get_ylim()[1],linestyles ="dashed", colors ="k")254 ax['B'].vlines(p+CI,ax['B'].get_ylim()[0],ax['B'].get_ylim()[1],linestyles ="dashed", colors ="k")255 check = (np.abs(sample - p)>CI)*sample256 ax['B'].scatter(check,np.arange(0,times),c = "red")257 ax['B'].scatter(np.zeros(times),np.arange(0,times),c = "white")258 #print(check)259 return fig260 261import math262import numpy as np263import scipy264def random_gen_data_mean(population,mu,sigma):265 results = np.random.normal(mu, sigma,population) # 均值,标准差,数量266 return results267def get_ratio1_mean(results):268 #print(results.size,results.sum(),"#########")269 return results.sum()/results.size270def simulation_mean(SRS_size,population,times):271 results = np.zeros(times)272 for i in range(times):273 sample1 = np.array(random.sample(population.tolist(),int(SRS_size)))274 results[i] = get_ratio1_mean(sample1)275 return results#,population276def get_CI_mean(sample,alpha):#size = (SRS_size)277 #print(sample)278 sigma = sample.std()279 sigma1 = sigma/math.sqrt(sample.shape[0])280 z = scipy.stats.norm.ppf(1-(1-alpha)/2)281 return(sigma*z)282def get_CI_mean_new(sample,alpha):#size = (SRS_size)283 #print(sample)284 sigma = sample285 sigma1 = sigma/math.sqrt(sample.shape[0])286 z = scipy.stats.norm.ppf(1-(1-alpha)/2)287 return(sigma*z)288def normfun(x,mu,std):289 pdf = np.exp(-((x - mu)**2)/(2*std**2)) / (std * np.sqrt(2*np.pi))290 return pdf291def normfun_d(x):292 u = x-x.mean()293 pdf = x * np.exp(-u*u/(2*x.std()*x.std()))/(x.std()*x.std()*x.std()*math.sqrt(math.pi))294# 295import matplotlib.pyplot as plt296import random297import math298import numpy as np299import scipy300import gradio as gr301from scipy.stats import norm302def draw_mean(mu,sigma,alpha,SRS_size,population_size,times):303 population = random_gen_data_mean(population_size,mu, sigma)304 #print(population)305 sample = simulation_mean(SRS_size,population,times)306 CI = get_CI_mean(sample,alpha)307 #print(CI)308 #all_colors = list(plt.cm.colors.cnames.keys())309 random.seed(42)310 #c = random.choices(all_colors, k=population_size)311 312 313 fig,ax=plt.subplot_mosaic([['A','A','A'],314 ['C','C','C'],315 ['B','B','B'],316 ['B','B','B'],317 ['B','B','B']],figsize=(16, 10), constrained_layout=True)318 #ax['C'].set_ylim(0,0.5)319 ax['A'].set_xlim(population.min()-CI, population.max()+CI)320 ax['B'].set_xlim(population.min()-CI, population.max()+CI)321 ax['C'].set_xlim(population.min()-CI, population.max()+CI)322 mu = population.mean()323 std = population.std()324 n, bins, patches = ax['A'].hist(population,bins = 40)#,density = True)325 y = normfun(bins,mu,std)326 ax['A'].vlines(mu,ax['A'].get_ylim()[0],ax['A'].get_ylim()[1],linestyles ="solid", colors ="k")327 #ax['A'].vlines(mu-CI,ax['A'].get_ylim()[0],ax['A'].get_ylim()[1],linestyles ="dashed", colors ="k")328 #ax['A'].vlines(mu+CI,ax['A'].get_ylim()[0],ax['A'].get_ylim()[1],linestyles ="dashed", colors ="k")329 330 ax['C'].vlines(mu,ax['C'].get_ylim()[0],ax['C'].get_ylim()[1],linestyles ="solid", colors ="k")331 ax['C'].vlines(mu-CI,ax['C'].get_ylim()[0],ax['C'].get_ylim()[1],linestyles ="dashed", colors ="k")332 ax['C'].vlines(mu+CI,ax['C'].get_ylim()[0],ax['C'].get_ylim()[1],linestyles ="dashed", colors ="k")333 334 ax['C'].plot(bins,y,'r--',label = "population dis",color = "red")335 336 mu = sample.mean()337 std = sample.std()338 339 340 n, bins, patches = ax['C'].hist(sample,bins = 40)#,density = True,weights = weights)341 y = normfun(bins,mu,std)342 ax['C'].plot(bins,y,'-',label = "sample dis",color = "red")343 ax['C'].legend()344 345 ax['B'].hlines(y = np.arange(0,times),xmin = sample-CI, xmax = sample+CI,color = "gray")346 ax['B'].vlines(mu,ax['B'].get_ylim()[0],ax['B'].get_ylim()[1],linestyles ="solid", colors ="k")347 ax['B'].vlines(mu-CI,ax['B'].get_ylim()[0],ax['B'].get_ylim()[1],linestyles ="dashed", colors ="k")348 ax['B'].vlines(mu+CI,ax['B'].get_ylim()[0],ax['B'].get_ylim()[1],linestyles ="dashed", colors ="k")349 ax['B'].scatter(sample,np.arange(0,times))350 351 check = (np.abs(sample - mu)>CI)*sample352 ax['B'].scatter(check,np.arange(0,times),c = "red")353 ax['B'].scatter(np.zeros(times),np.arange(0,times),c = "white")354 return fig355 356 357with gr.Blocks() as demo:358 with gr.Tab(label='Birthday paradox', open=False):359 with gr.Accordion(label="intro", open = True):360 gr.Markdown("The Birthday paradox illustrates a scnario to find whether or how many people are there with the same birthdays among the 365 days a year")361 gr.Markdown("Basically, if we have n people, and each of the birthdays were randomly selected, the predicted probability for at least two people with the same birthdays is:")362 gr.Markdown("1 - (365 * 364 * 363 * (3650n+1))/(365^n)")363 gr.Markdown("since for the first person, his birthday could be any day. For the second, to avoid have same birthdays, there are 364 remaining, and so on")364 gr.Markdown("we could use 1 minus the value to get the final probability")365 birthday_intro_input1 = gr.Slider(2,1000,label = "number of people")366 birthday_intro_output1 = gr.Textbox(label = "the probability to find same birthdays among the gien value of people is")367 birthday_intro_button = gr.Button()368 with gr.Accordion(label="single simulation", open = False):369 gr.Markdown("this simulation is for one possible outcome given that the number of people is fixed(the slider)")370 gr.Markdown("because I am to discuss the long term probability of multiple simulations combined, so the tag is single simulation")371 birthday_single_simulation_input1 = gr.Slider(0,100,label = "number of people")372 with gr.Row():373 birthday_single_simulation_output1 = gr.Textbox("birthdays")374 birthday_single_simulation_output2 = gr.Textbox()375 birthday_single_simulation_output3 = gr.Textbox()376 birthday_single_simulation_output4 = gr.Textbox()377 birthday_single_simulation_output5 = gr.Plot()378 birthday_single_simulation_button = gr.Button()379 with gr.Accordion(label = "fix number of people, simulate n times to find approximate probability",open = False):380 with gr.Row():381 one_zero_input2 = gr.Slider(0,100,label = "number of people")382 one_zero_input1 = gr.Slider(10,10000,10, label = "number of experiments")383 one_zero_output = gr.Plot()384 with gr.Row():385 one_zero_output1 = gr.Textbox(label = "experiments with no same birthdays")386 one_zero_output2 = gr.Textbox(label = "experiments with same birthdays")387 one_zero_output3 = gr.Textbox(label = "probability of same birhtdays")388 one_zero_output4 = gr.Textbox(label = "estimated probability")389 one_zero_button = gr.Button()390 with gr.Accordion("normal distribution", open = False):391 gr.Markdown("it may took some time when the numbers are big")392 multi_input1 = gr.Slider(1,100,label = "number of people")393 multi_input2 = gr.Slider(10,1000,label = "number of experiments to get one probability")394 multi_input3 = gr.Slider(1,1000,label = "number of probabilities")395 multi_output = gr.Plot()396 multi_button = gr.Button()397 with gr.Tab(label="Confidence Interval"):398 with gr.Accordion("CI for ratio"):399 gr.Markdown("given that the probability to succeed in each trial is p")400 gr.Markdown("find the Confidence Interval under certain significance level")401 with gr.Row():402 input1 = gr.Slider(0,1,0.001,label = "probability to succeed in each trial")403 input2 = gr.Slider(0,1,0.001,label = "significane level")404 input3 = gr.Slider(0,1000,1,label = "sample size")405 input4 = gr.Slider(1,100000,1,label = "population")406 input5 = gr.Slider(1,1000,1,label = "number of trials")407 output = gr.Plot()408 text_button = gr.Button()409 with gr.Accordion("CI for mean"):410 gr.Markdown("tip: the y axis has been normalized")411 with gr.Row():412 input11 = gr.Slider(0,100,1,label = "mean")413 input21 = gr.Slider(0,100,1,label = "std")414 input31 = gr.Slider(0,1,0.001,label = "signifigance level")415 input41 = gr.Slider(0,1000,1,label = "sample size")416 input51 = gr.Slider(1,100000,1,label = "population")417 input61 = gr.Slider(1,1000,1,label = "sample times")418 output1 = gr.Plot()419 text_button1 = gr.Button()420 text_button1.click(draw_mean, inputs=[input11,input21,input31,input41,input51,input61], outputs=[output1])421 birthday_intro_button.click(birthdayCheck, inputs = [birthday_intro_input1], outputs = [birthday_intro_output1])422 birthday_single_simulation_button.click(getSame_gradio, inputs = [birthday_single_simulation_input1], outputs = [birthday_single_simulation_output1,birthday_single_simulation_output2,birthday_single_simulation_output3,birthday_single_simulation_output4,birthday_single_simulation_output5])423 one_zero_button.click(draw_one_sim_one_zero_gradio, inputs=[one_zero_input1,one_zero_input2], outputs=[one_zero_output,one_zero_output1,one_zero_output2,one_zero_output3,one_zero_output4])424 multi_button.click(multi_multi_simulation_gradio,inputs = [multi_input1,multi_input2,multi_input3],outputs = multi_output)425 text_button.click(draw, inputs=[input1,input2,input3,input4,input5], outputs=[output])426demo.launch()