Probability Foundations

import numpy as np
from math import comb, perm, factorial
from itertools import combinations, permutations

# Basic counting
print(f"Permutations P(10,3): {perm(10, 3)}")  # 720
print(f"Combinations C(10,3): {comb(10, 3)}")  # 120

# Probability of union
def prob_union(p_a, p_b, p_ab):
    """P(A ∪ B) = P(A) + P(B) - P(A ∩ B)"""
    return p_a + p_b - p_ab

# Inclusion-Exclusion for three events
def prob_union_three(p_a, p_b, p_c, p_ab, p_ac, p_bc, p_abc):
    """P(A ∪ B ∪ C) using inclusion-exclusion"""
    return (p_a + p_b + p_c 
            - p_ab - p_ac - p_bc 
            + p_abc)

# Birthday problem
def birthday_prob(n, days=365):
    """Probability of at least one collision in n birthdays"""
    if n > days:
        return 1.0
    prob_no_collision = 1.0
    for i in range(n):
        prob_no_collision *= (days - i) / days
    return 1 - prob_no_collision

# Simulation approach
def simulate_birthday(n, days=365, trials=100000):
    """Monte Carlo simulation for birthday problem"""
    collisions = 0
    for _ in range(trials):
        birthdays = np.random.randint(0, days, size=n)
        if len(np.unique(birthdays)) < n:
            collisions += 1
    return collisions / trials

# Dice problems
def dice_probability(target_sum, num_dice=2, faces=6):
    """Probability of getting target_sum with num_dice dice"""
    # Count favorable outcomes
    favorable = 0
    total = faces ** num_dice
    
    # Brute force for small cases
    if num_dice == 2:
        for i in range(1, faces+1):
            for j in range(1, faces+1):
                if i + j == target_sum:
                    favorable += 1
    else:
        # Use dynamic programming
        dp = np.zeros((num_dice + 1, target_sum + 1))
        dp[0][0] = 1
        for die in range(1, num_dice + 1):
            for sum_val in range(1, target_sum + 1):
                for face in range(1, min(faces, sum_val) + 1):
                    dp[die][sum_val] += dp[die-1][sum_val-face]
        favorable = dp[num_dice][target_sum]
    
    return favorable / total

# Bayes' theorem
def bayes_theorem(p_b_given_a, p_a, p_b):
    """P(A|B) = P(B|A) * P(A) / P(B)"""
    return (p_b_given_a * p_a) / p_b

# Example: Medical test
p_disease = 0.01  # Prior probability
p_positive_given_disease = 0.95  # Sensitivity
p_positive_given_healthy = 0.05  # False positive rate

p_healthy = 1 - p_disease
p_positive = (p_positive_given_disease * p_disease + 
              p_positive_given_healthy * p_healthy)
p_disease_given_positive = bayes_theorem(
    p_positive_given_disease, p_disease, p_positive
)

print(f"Probability of disease given positive test: {p_disease_given_positive:.4f}")

# Examples
print(f"\nBirthday problem (23 people): {birthday_prob(23):.4f}")
print(f"Dice: P(sum=7 with 2 dice) = {dice_probability(7):.4f}")
print(f"Simulation: {simulate_birthday(23):.4f}")