<jupyter_start><jupyter_text>Estimating the weight of a biased coin- There is an unknown $\mathsf{Bernoulli}(p)$ object.- How many samples are needed in order to obtain a **good** estimate $\widehat{p}$?<jupyter_code># unknown distribution
# p = np.array([0.1, 0.6, 0.005, -1])
# p[-1] = 1 - sum(p[:-1])

p = np.random.uniform(0, 1, size=100)
p = p/p.sum()
assert np.isclose(sum(p), 1)

# print(p)

# item whose weight we want to estimate, i.e. p[σ]
σ = 0

def get_samples(n):
    """Get n random samples from `p`."""
    return np.random.choice(np.arange(len(p)), size=n, p=p, replace=True)
ϵ = 0.05
δ = 0.05

n = math.ceil((3/ϵ**2) * math.log(2/δ))
# n = math.ceil((0.5/ϵ**2) * math.log(2/δ))
print(f"{n=:,d}")

xs = get_samples(n)

emp = np.array([
    (xs == i).sum() / n for i in range(len(p))
])

l1 = np.abs(emp - p).sum()
l1
σ = 2
mle = (xs == σ).sum() / n
print(mle, p[σ])

print(f"{np.abs(mle - p[σ]) <=        ϵ = }")
print(f"{np.abs(mle - p[σ]) <= p[σ] * ϵ = }")<jupyter_output><empty_output><jupyter_text>Additive Chernoff bounds for distinguishing between two differently biased coins$\begin{cases}\texttt{coin}_1: \Pr(H) = \frac12 + \epsilon_1 \\\texttt{coin}_2: \Pr(H) = \frac12 - \epsilon_2\end{cases}$<jupyter_code>class Coin:
    def __init__(self, p):
        self.p = p
        
    def sample(self, n=1):
        return np.random.binomial(n=1, p=self.p, size=n)
    
    def __str__(self):
        return f'Coin({p})'
    
    def __repr__(self):
        return str(self)
    

def err(coin, n, m):
    x = []
    for _ in tqdm(range(m)):
        x += [coin.sample(n).sum()/n - coin.p]

    x = np.array(x)

    sns.histplot(x, bins=2**5)
    plt.show()
coin = Coin(p=1/2)
err(coin, n=1000, m=100_000)
ϵ1 = 0.1
ϵ2 = 0
δ  = 0.01
λ  = (1 + ϵ1 - ϵ2) / 2 # best threshold
assert 0 <= ϵ1 < 1 and 0 <= ϵ2 < 1 and 0 < δ < 1

coin1 = Coin(p=1/2 + ϵ1)
coin2 = Coin(p=1/2 - ϵ2)

n = math.ceil(2 * math.log(1/δ) / (ϵ1 + ϵ2)**2)

print(f"best threshold {λ=}")

# make the bound a function of n
P_wrong_bound = lambda n: 0.5 * (math.exp(-2 * n * (λ-0.5+ϵ2)**2) + math.exp(-2 * n * (λ-0.5-ϵ1)**2))
print(f'P(wrong | {n=}) <= {P_wrong_bound(n):.5f}')

# one can actuall
ys = np.array([P_wrong_bound(x) for x in np.arange(1, 10000)])
print(f'smallest n that will get P(wrong) < δ = {np.argwhere(ys <= δ)[0].item()}')
def experiment(n):
    # choose coin uniformly at random
    coin = np.random.choice([coin1, coin2])
    
    # count number of heads in n trials
    X = coin.sample(n).sum()
    
    if X/n > λ:
        return coin is coin1
    else:
        return coin is coin2

N = 100_000
correct = 0
diff = []

for _ in tqdm(range(N), desc='Experiments'):
    correct += experiment(n)
    
P_wrong = 1 - (correct / N)

assert P_wrong <= P_wrong_bound(n)
print(f'experimental {P_wrong=:.7f} <= {P_wrong_bound(n)=:.7f}')
xs = [coin1.sample(921).sum() / n > λ for _ in tqdm(range(1_000_000))]<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.19890822244967588, 'fraction_numerical': 0.04674172637325145, 'mean_word_length': 3.0553250345781464, 'pattern_counts': {'":': 0, '<': 22, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '482641', 'license': 'MIT', 'n_tokens_mistral': 1367, 'n_tokens_neox': 1283, 'n_words': 307, 'path': 'src/chernoff.ipynb', 'repo_name': 'alexandru-dinu/notebooks', 'size': '2931'}