# How to verify Bayes theorem output?

I am looking for if I can do some simulation in R for verifying the output of Bayes theorem , if you have some other technique I am open to know that also. Anything that can help me in verification would be useful.

Actually my problem is something like this:-

I have probabilities of dropping a call based on the duration of time customer has waited for before dropping.

Probability of dropping during 0-20sec: $$0.23$$

Probability of dropping during 20-40sec: $$0.11$$

Probability of dropping during 40-50sec: $$0.04$$

Probability of dropping during 50-60sec: $$0.61$$

Now the probability of dropping the call during 50-60sec given that the call has not been dropped during 0-20sec is denoted as P(50-60sec dropped/Not dropped 0-20sec).

And P(50-60sec dropped/Not dropped 0-20sec) is found to be $$\frac{0.61}{1-0.23}=0.79$$.

So now I want to verify if 0.79 is correct or not by some simulation or any other technique.

EDIT: The sum of probabilities adds up to 0.99 because those are calculated from a dataset and rounded off up to 2 digits.

• simulate n realisations of the 4 types of call, keep only those above 20sec, count the proportion among these with duration 50-60sec, done! Oct 31, 2020 at 8:35

Here's some Python that generates 100,000 samples. I'll leave it to you to write it in R.

import numpy as np

n50_60 = n_not0_20 = 0

for i in range(100_000):
u = np.random.uniform()
if u <= 0.23:
pass
elif u <= 0.23 + 0.11:
n_not0_20 += 1
elif u <= 0.23 + 0.11 + 0.04:
n_not0_20 += 1
else:
n_not0_20 += 1
n50_60 += 1

print(n50_60 / n_not0_20)
$$$$


Here is the simulation in R which tries to verify the Probability of dropping during 50-60sec by using P(50-60sec dropped/Not dropped 0-20sec).

#This will have ones with the probability of 0.77
not_20_sec<-rbinom(n=10000, size=1, prob=1-0.23)
#Below code will put 'ones'1' in 'y' with the probabilty of 0.7951
#wherever it was a '1' in  not_20_sec else it would add '0'.
y=c()
for (i in not_20_sec){
if (i==1){
y<-c(y,rbinom(n=10, size=1, prob=0.7951))
}
else{
y<-c(y,replicate(10,0))
}
}
#This is the probability of dropping during 50-60sec.
print(sum(y==1)/length(y)) # 0.609
`