概率论基础

import random

# 模拟掷骰子实验
def simulate_dice(n: int = 10000) -> dict[int, float]:
    counts = {i: 0 for i in range(1, 7)}
    for _ in range(n):
        counts[random.randint(1, 6)] += 1
    return {k: v / n for k, v in counts.items()}

print(simulate_dice(10000))

# 蒙特卡洛方法估算圆周率 pi
import random

def estimate_pi(n: int = 100000) -> float:
    inside = 0
    for _ in range(n):
        x, y = random.random(), random.random()
        if x2 + y2 <= 1:
            inside += 1
    return 4 * inside / n

print(f"Pi ≈ {estimate_pi()}")

# 条件概率计算：医学检测示例
# 已知：疾病患病率 1%，检测灵敏度 99%，假阳性率 5%
def bayes_medical():
    p_disease = 0.01          # 先验概率
    p_positive_given_disease = 0.99   # 灵敏度
    p_positive_given_healthy = 0.05   # 假阳性率
    p_healthy = 1 - p_disease
    
    # 全概率公式求 P(阳性)
    p_positive = (p_positive_given_disease * p_disease +
                  p_positive_given_healthy * p_healthy)
    
    # 贝叶斯公式
    p_disease_given_positive = (
        p_positive_given_disease * p_disease / p_positive
    )
    return p_disease_given_positive

print(f"P(患病|阳性) = {bayes_medical():.4f}")

# 验证独立性：掷两枚骰子
import random

def check_independence(n: int = 100000):
    a_and_b = 0  # 两枚都是6
    a_count = 0   # 第一枚是6
    b_count = 0   # 第二枚是6
    
    for _ in range(n):
        d1 = random.randint(1, 6)
        d2 = random.randint(1, 6)
        if d1 == 6: a_count += 1
        if d2 == 6: b_count += 1
        if d1 == 6 and d2 == 6: a_and_b += 1
    
    p_a = a_count / n
    p_b = b_count / n
    p_ab = a_and_b / n
    print(f"P(A) = {p_a:.4f}")
    print(f"P(B) = {p_b:.4f}")
    print(f"P(A∩B) = {p_ab:.4f}")
    print(f"P(A)*P(B) = {p_a * p_b:.4f}")
    print(f"独立? {abs(p_ab - p_a * p_b) < 0.01}")

check_independence()

# 离散随机变量：二项分布模拟
import numpy as np

def binomial_simulation(n: int = 10, p: float = 0.5, trials: int = 10000):
    samples = np.random.binomial(n, p, trials)
    pmf = {}
    for k in samples:
        pmf[k] = pmf.get(k, 0) + 1
    pmf = {k: v / trials for k, v in sorted(pmf.items())}
    
    print(f"理论均值: {n * p}")
    print(f"样本均值: {np.mean(samples):.4f}")
    print(f"理论方差: {n * p * (1 - p)}")
    print(f"样本方差: {np.var(samples):.4f}")
    return pmf

binomial_simulation()

# 连续随机变量：正态分布 PDF 和 CDF
import numpy as np
from scipy.stats import norm

def normal_distribution_analysis():
    mu, sigma = 0, 1  # 标准正态分布
    x = np.linspace(-4, 4, 200)
    
    pdf = norm.pdf(x, mu, sigma)
    cdf = norm.cdf(x, mu, sigma)
    
    # 68-95-99.7 法则验证
    print(f"P(|X| < 1σ) = {norm.cdf(1) - norm.cdf(-1):.4f}")
    print(f"P(|X| < 2σ) = {norm.cdf(2) - norm.cdf(-2):.4f}")
    print(f"P(|X| < 3σ) = {norm.cdf(3) - norm.cdf(-3):.4f}")
    
    # 分位数计算
    print(f"95% 分位数: {norm.ppf(0.95):.4f}")

normal_distribution_analysis()

# 期望和方差的计算与验证
import numpy as np

def expectation_and_variance():
    # 模拟：掷两个骰子之和
    n = 1000000
    d1 = np.random.randint(1, 7, n)
    d2 = np.random.randint(1, 7, n)
    S = d1 + d2  # 两枚骰子之和
    
    # 理论值
    e_d = 3.5
    var_d = 35/12  # 单个骰子方差
    e_s_theory = 2 * e_d  # 期望的线性性质
    var_s_theory = 2 * var_d  # 独立时方差可加
    
    print(f"E(S) 理论={e_s_theory:.4f}, 样本={np.mean(S):.4f}")
    print(f"Var(S) 理论={var_s_theory:.4f}, 样本={np.var(S):.4f}")
    
    # 验证线性性：E(X+Y) = E(X) + E(Y)
    print(f"E(X)+E(Y) = {np.mean(d1) + np.mean(d2):.4f}")

expectation_and_variance()

# 协方差矩阵与 PCA 的关系
import numpy as np

def covariance_demo():
    np.random.seed(42)
    # 生成相关数据
    x = np.random.randn(200)
    y = 0.5 * x + np.random.randn(200) * 0.3
    
    # 计算协方差矩阵
    cov_matrix = np.cov(x, y)
    corr_matrix = np.corrcoef(x, y)
    
    print(f"协方差矩阵:
{cov_matrix}")
    print(f"相关系数矩阵:
{corr_matrix}")
    
    # 特征值分解
    eigenvalues = np.linalg.eigvals(cov_matrix)
    print(f"特征值: {np.sort(eigenvalues)[::-1]}")
    print(f"第一主成分解释方差比: {max(eigenvalues)/sum(eigenvalues):.4f}")

covariance_demo()

# 三大离散分布对比
import numpy as np
from scipy.stats import bernoulli, binom, poisson

def discrete_distributions():
    # 伯努利分布
    p = 0.3
    print(f"伯努利: P(X=1) = {bernoulli.pmf(1, p):.3f}")
    
    # 二项分布
    n, p_bin = 20, 0.3
    k_values = list(range(n + 1))
    pmf_bin = binom.pmf(k_values, n, p_bin)
    mode_bin = np.argmax(pmf_bin)
    print(f"二项 B({n},{p_bin}): 众数={mode_bin}, 均值={n*p_bin}")
    
    # 泊松分布
    lam = 6
    k_pois = list(range(15))
    pmf_pois = poisson.pmf(k_pois, lam)
    print(f"泊松 λ={lam}: P(X=6) = {pmf_pois[6]:.3f}")
    print(f"泊松均值={lam}, 方差={lam}")

discrete_distributions()

# 正态分布的性质可视化
import numpy as np
from scipy.stats import norm

def normal_properties():
    # 不同参数的正态分布
    params = [(0, 1), (0, 2), (2, 1)]
    x = np.linspace(-6, 8, 300)
    
    for mu, sigma in params:
        pdf = norm.pdf(x, mu, sigma)
        peak = norm.pdf(mu, mu, sigma)
        print(f"N({mu}, {sigma}²): 峰值 f(μ) = {peak:.4f}")
        print(f"  P(μ-σ < X < μ+σ) = {norm.cdf(mu+sigma, mu, sigma) - norm.cdf(mu-sigma, mu, sigma):.4f}")
    
    # 标准正态分布的尾部概率
    print(f"P(Z > 1.96) = {1 - norm.cdf(1.96):.4f}")
    print(f"P(Z > 2.58) = {1 - norm.cdf(2.58):.4f}")

normal_properties()

# 大数定律的可视化验证
import numpy as np

def law_of_large_numbers():
    n_max = 100000
    true_mean = 0.5  # 均匀分布 U(0,1) 的期望
    
    # 生成样本并计算累积均值
    samples = np.random.uniform(0, 1, n_max)
    cumulative_mean = np.cumsum(samples) / np.arange(1, n_max + 1)
    
    # 验证收敛性
    for n in [100, 1000, 10000, 100000]:
        err = abs(cumulative_mean[n-1] - true_mean)
        print(f"n={n:>7d}: 样本均值={cumulative_mean[n-1]:.6f}, 误差={err:.6f}")
    
    print(f"
真实期望: {true_mean}")
    print(f"大数定律: 样本均值 → 理论期望")

law_of_large_numbers()

# 中心极限定理演示：任意分布的样本均值趋向正态
import numpy as np
from scipy.stats import norm, kstest

def central_limit_theorem():
    np.random.seed(42)
    sample_size = 50  # 每组样本量
    n_groups = 10000   # 组数
    
    # 用指数分布（明显非正态）
    group_means = []
    for _ in range(n_groups):
        sample = np.random.exponential(scale=2, size=sample_size)
        group_means.append(np.mean(sample))
    
    sample_means = np.array(group_means)
    theoretical_mean = 2  # 指数分布均值
    theoretical_std = 2 / np.sqrt(sample_size)
    
    # KS 检验正态性
    z_scores = (sample_means - theoretical_mean) / theoretical_std
    ks_stat, p_value = kstest(z_scores, "norm")
    
    print(f"样本均值分布: mean={np.mean(sample_means):.4f}, std={np.std(sample_means):.4f}")
    print(f"理论值: mean={theoretical_mean}, std={theoretical_std:.4f}")
    print(f"KS 检验 p={p_value:.4f} (>0.05 则接受正态假设)")

central_limit_theorem()

# 贝叶斯定理：医疗检测的完整分析
from typing import NamedTuple

class TestResult(NamedTuple):
    sensitivity: float    # 真阳性率
    specificity: float    # 真阴性率
    prevalence: float     # 患病率

def bayes_analysis(test: TestResult) -> dict:
    p_d = test.prevalence
    p_h = 1 - p_d
    sens = test.sensitivity
    spec = test.specificity
    
    # 全概率
    p_pos = sens * p_d + (1 - spec) * p_h
    p_neg = (1 - sens) * p_d + spec * p_h
    
    # 贝叶斯
    ppv = sens * p_d / p_pos      # 阳性预测值
    npv = spec * p_h / p_neg      # 阴性预测值
    
    return {"PPV": ppv, "NPV": npv, "P(pos)": p_pos}

# 案例1：罕见病检测
rare = TestResult(0.99, 0.95, 0.001)
print(f"罕见病: PPV={bayes_analysis(rare)['PPV']:.4f}")

# 案例2：常见病检测
common = TestResult(0.99, 0.95, 0.3)
print(f"常见病: PPV={bayes_analysis(common)['PPV']:.4f}")

# 贝叶斯拼写纠正器
import re
from collections import Counter

class SpellCorrector:
    def __init__(self, corpus: str):
        words = re.findall(r"\w+", corpus.lower())
        self.word_freq = Counter(words)
        self.total = len(words)
        self.alphabet = "abcdefghijklmnopqrstuvwxyz"
    
    def P(self, word: str) -> float:
        """先验概率：词在语料中的频率"""
        return self.word_freq.get(word, 1) / self.total
    
    def edits(self, word: str) -> set[str]:
        """编辑距离为1的候选词"""
        splits = [(word[:i], word[i:]) for i in range(len(word)+1)]
        deletes = [L+R[1:] for L,R in splits if R]
        inserts = [L+c+R for L,R in splits for c in self.alphabet]
        return set(deletes + inserts)
    
    def correct(self, word: str) -> str:
        """贝叶斯：argmax P(c|w) ∝ P(w|c)·P(c)"""
        word = word.lower()
        if word in self.word_freq:
            return word
        candidates = self.edits(word) & set(self.word_freq)
        if not candidates:
            return word
        return max(candidates, key=self.P)

corpus = "this is a simple test corpus with some common english words and the word spelling is very important"
sc = SpellCorrector(corpus)
print(f"corect → {sc.correct('corect')}")
print(f"importnt → {sc.correct('importnt')}")