You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
48 KiB
48 KiB
None
<html lang="en">
<head>
</head>
![No description has been provided for this image]()
![No description has been provided for this image]()
![No description has been provided for this image]()
</html>
In [1]:
# Probability basics 10/02/26
import numpy as np
import matplotlib.pyplot as plt
import math
rng = np.random.default_rng(7)
$$ P(A) = \frac{\text{Success}}{\text{Trial}}$$
In [2]:
def bernoulli(p, rng=rng):
return rng.random() < p
def estimate(trials, fn):
"""Return running estimate of P(A) over time"""
x = np.fromiter((fn() for _ in range(trials)), dtype=np.int8, count=trials)
return np.cumsum(x)/np.arange(1,trials+1)
In [3]:
# Example 1: coin flip as Bernoulli trails
p = 0.5 # fair coin
N = 5000
# running estimate of P(A) over time; A=success
estimation = estimate(N, lambda: bernoulli(p,rng))
plt.plot(estimation,'.k')
plt.axhline(p,color='red')
plt.show()
In [5]:
x = np.fromiter((bernoulli(p,rng) for _ in range(N)), dtype=np.int8, count=N)
proba = np.cumsum(x)/np.arange(1, N+1)
plt.plot(proba, 'ok')
Out[5]:
In [12]:
# Example 2: Law of total probability via simulation
# Objective: Estimate the overall probability of drawing a red ball, when the ball is
# drawn from a randomly selected urn
#-=-=-=-=-=-=-=
# Parameters:
pA = 0.3 # probability that the experiment selects urn A
pRedA = 0.7 # probability of red given urn A was selected P(R|A)
pRedB = 0.2 # P(R|B)
N = 200_000
def bernoulli(p, rng=rng):
return rng.random() < p
def drawRed():
chooseA = bernoulli(pA)
if chooseA:
return bernoulli(pRedA)
else:
return bernoulli(pRedB)
def estimate(trials, fn):
"""Return running estimate of P(A) over time"""
x = np.fromiter((fn() for _ in range(trials)), dtype=np.int8, count=trials)
return np.cumsum(x)/np.arange(1,trials+1)
pHat = estimate(N,drawRed)
pHat[-1]
pExact = pA*pRedA+(1-pA)*pRedB
pHat[-1], pExact
Out[12]:
In [16]:
plt.plot(pHat,'.k')
plt.xlim([0,50])
plt.show()