不确定性推理

"""贝叶斯公式：阳性结果下真正患病的概率（基础比率忽略的教训）"""

def bayes(prior, tp_rate, fp_rate):
    """prior=先验患病率；tp_rate=P(阳|病)；fp_rate=P(阳|无病)。"""
    p_positive = tp_rate * prior + fp_rate * (1 - prior)   # 全概率
    return tp_rate * prior / p_positive                    # 贝叶斯

p = bayes(prior=0.001, tp_rate=0.99, fp_rate=0.05)
print(f"检测阳性后真正患病的概率：{p:.4f}，即约 {p*100:.2f}%")

# 观察先验发病率变化对结果的影响
print("\n不同发病率下的阳性预测值：")
for prior in [0.001, 0.01, 0.1, 0.5]:
    post = bayes(prior, 0.99, 0.05)
    print(f"  发病率 {prior:6.3f}  →  P(病|阳性) = {post:6.3f}")

"""贝叶斯网络枚举推理（Burglary-Earthquake-Alarm，AIMA 经典例题）"""
from itertools import product

class BayesNet:
    def __init__(self, nodes):
        # nodes: {变量名: (父节点列表, cpt函数)}
        # cpt(该变量取值, {父变量: 取值}) -> P(变量=取值 | 父取值)
        self.nodes = nodes
        self.vars = list(nodes.keys())

    def parents(self, var):
        return self.nodes[var][0]

    def cpt(self, var, value, parent_values):
        return self.nodes[var][1](value, parent_values)


def enumeration_ask(X, evidence, net):
    """求 P(X | evidence)，返回两个取值的归一化概率分布。"""
    dist = {}
    for xi in (True, False):
        dist[xi] = enumerate_all(net.vars, {**evidence, X: xi}, net)
    total = sum(dist.values())
    return {k: v / total for k, v in dist.items()}

def enumerate_all(variables, assignment, net):
    """按变量列表递归枚举（变量已按拓扑序排列）。"""
    if not variables:
        return 1.0
    first, rest = variables[0], variables[1:]
    if first in assignment:                       # 证据变量：直接取表中概率
        pv = {p: assignment[p] for p in net.parents(first)}
        return net.cpt(first, assignment[first], pv) * enumerate_all(rest, assignment, net)
    total = 0.0
    for val in (True, False):                     # 隐变量：对两种取值求和
        pv = {p: assignment[p] for p in net.parents(first)}
        total += net.cpt(first, val, pv) * enumerate_all(rest, {**assignment, first: val}, net)
    return total


if __name__ == "__main__":
    # 先验
    P_B = lambda v, _: 0.001 if v else 0.999
    P_E = lambda v, _: 0.002 if v else 0.998
    # P(A | B, E)
    def P_A(v, p):
        b, e = p["B"], p["E"]
        table = {(True, True): .95, (True, False): .94,
                 (False, True): .29, (False, False): .001}
        pa = table[(b, e)]
        return pa if v else 1 - pa
    def P_J(v, p):
        pj = .90 if p["A"] else .05
        return pj if v else 1 - pj
    def P_M(v, p):
        pm = .70 if p["A"] else .01
        return pm if v else 1 - pm

    net = BayesNet({
        "B": ([], P_B),
        "E": ([], P_E),
        "A": (["B", "E"], P_A),
        "J": (["A"], P_J),
        "M": (["A"], P_M),
    })

    q1 = enumeration_ask("B", {"J": True, "M": True}, net)
    print(f"John 和 Mary 都来电时，发生盗窃的概率：P(B|J,M) = {q1[True]:.4f}")
    q2 = enumeration_ask("E", {"J": True, "M": True}, net)
    print(f"两人都来电时，发生地震的概率：P(E|J,M) = {q2[True]:.4f}")
    q3 = enumeration_ask("B", {"J": True}, net)
    print(f"仅 John 来电时，发生盗窃的概率：P(B|J)   = {q3[True]:.4f}")

"""主观 Bayes 方法：证据确定出现/不出现时后验几率的更新"""

def p2o(p):
    """概率 -> 几率"""
    return p / (1 - p)

def o2p(o):
    """几率 -> 概率"""
    return o / (1 + o)

def posterior(prior, ls, ln, evidence_state):
    """evidence_state: True=证据出现；False=证据不出现；None=不知道。"""
    if evidence_state is True:
        return o2p(ls * p2o(prior))
    if evidence_state is False:
        return o2p(ln * p2o(prior))
    return prior          # 无证据时维持先验（EH/CP 公式处理部分证据，此处从略）

if __name__ == "__main__":
    # PROSPECTOR 风格规则：
    # IF 有石英闪长岩 THEN (LS=300, LN=0.2) 有铜矿脉
    prior_h, ls, ln = 0.01, 300, 0.2

    p_yes = posterior(prior_h, ls, ln, True)
    p_no = posterior(prior_h, ls, ln, False)
    print(f"先验 P(H) = {prior_h}")
    print(f"证据出现   P(H|E)  = {p_yes:.4f}  （LS={ls} 强烈支持）")
    print(f"证据缺席   P(H|~E) = {p_no:.4f}  （LN={ln} 削弱假设）")

    # 多级推理：证据链 E1 -> H1 -> H2
    # 规则1 IF E1 THEN (20, 0.5) H1，H1 先验 0.02
    h1 = posterior(0.02, 20, 0.5, True)
    # 规则2 IF H1 THEN (10, 0.4) H2，H2 先验 0.05，把 H1 的后验当作软证据
    h2 = o2p(20 * p2o(0.05) * (h1 / 0.02))    # 按后验比近似加权（教学演示）
    print(f"\n证据链两级传递后：P(H1|E1)={h1:.4f}, P(H2) 被更新到约 {h2:.4f}")

"""MYCIN 可信度模型：证据合成、规则传递、结论合并"""

def cf_and(*cfs):
    """合取证据：取最小可信度。"""
    return min(cfs)

def cf_or(*cfs):
    """析取证据：取最大可信度。"""
    return max(cfs)

def cf_apply(rule_cf, evidence_cf):
    """不确定性传递：规则强度 × 证据可信度（证据为负则不触发）。"""
    return rule_cf * max(0.0, evidence_cf)

def cf_combine(cf1, cf2):
    """两条独立规则支持同一结论时的可信度合成。"""
    if cf1 >= 0 and cf2 >= 0:
        return cf1 + cf2 - cf1 * cf2
    if cf1 < 0 and cf2 < 0:
        return cf1 + cf2 + cf1 * cf2
    return (cf1 + cf2) / (1 - min(abs(cf1), abs(cf2)))


class CFSolver:
    def __init__(self):
        self.facts = {}      # 证据 -> CF
        self.rules = []      # (前提列表, 结论, 规则CF)

    def set_fact(self, name, cf):
        self.facts[name] = cf

    def add_rule(self, premises, conclusion, rule_cf):
        self.rules.append((premises, conclusion, rule_cf))

    def solve(self, target, trace=True):
        """前向链：反复用规则更新各结论 CF，直到稳定。"""
        cf = dict(self.facts)
        changed = True
        while changed:
            changed = False
            for premises, concl, rule_cf in self.rules:
                if any(p not in cf for p in premises):
                    continue
                evidence_cf = cf_and(*[cf[p] for p in premises])
                new = cf_apply(rule_cf, evidence_cf)
                old = cf.get(concl, 0.0)
                combined = cf_combine(old, new) if old != 0 else new
                if abs(combined - old) > 1e-9:
                    cf[concl] = round(combined, 4)
                    changed = True
                    if trace:
                        print(f"规则 {premises}=>{concl}：证据CF={evidence_cf:.2f}，"
                              f"本次贡献={new:.3f}，合并后 CF({concl})={cf[concl]:.3f}")
        return cf.get(target, 0.0), cf

if __name__ == "__main__":
    # 简化医疗诊断网络（症状可信度来自观察）
    s = CFSolver()
    s.set_fact("发烧", 0.8)
    s.set_fact("咳嗽", 0.6)
    s.set_fact("X光阴影", 0.5)

    s.add_rule(["发烧"], "感染", 0.7)
    s.add_rule(["咳嗽"], "呼吸道感染", 0.75)
    s.add_rule(["感染", "X光阴影"], "肺炎", 0.8)
    s.add_rule(["呼吸道感染", "X光阴影"], "肺炎", 0.5)   # 第二条支持肺炎的规则

    target_cf, all_cf = s.solve("肺炎")
    print("\n最终可信度：")
    for k, v in all_cf.items():
        print(f"  CF({k}) = {v:+.2f}")
    print(f"\n诊断：肺炎的可信度 = {target_cf:+.3f}"
          f"（>0.2 可认为有一定支持，>0.8 为高度可信）")

"""D-S 证据理论：基本概率分配、信任函数、似然函数、Dempster 组合"""
from itertools import combinations

def powerset(iterable):
    s = list(iterable)
    return [set(c) for r in range(len(s) + 1) for c in combinations(s, r)]

def ds_combine(m1, m2):
    """组合两个 mass 函数（键为 frozenset）。返回 (组合mass, 冲突因子K)。"""
    k = 0.0
    raw = {}
    for a, ma in m1.items():
        for b, mb in m2.items():
            inter = a & b
            if not inter:
                k += ma * mb                  # 冲突质量
            else:
                raw[frozenset(inter)] = raw.get(frozenset(inter), 0) + ma * mb
    if k >= 1.0:
        raise ValueError("证据完全冲突（K=1），无法组合")
    norm = 1 - k
    return {a: v / norm for a, v in raw.items()}, k

def belief(mass, subset, theta):
    return sum(v for a, v in mass.items() if a <= subset)

def plausibility(mass, subset, theta):
    complement = set(theta) - subset
    return 1 - belief(mass, frozenset(complement), theta)

if __name__ == "__main__":
    theta = ["感冒", "流感", "健康"]

    # 证据1：体温计读数，主要支持{流感}，部分支持{感冒,流感}
    m1 = {frozenset(["流感"]): 0.6, frozenset(["感冒", "流感"]): 0.3,
          frozenset(theta): 0.1}
    # 证据2：问诊，支持{感冒}，把一部分信任留给整个框架（表示不知道）
    m2 = {frozenset(["感冒"]): 0.5, frozenset(theta): 0.5}

    combined, K = ds_combine(m1, m2)
    print(f"冲突因子 K = {K:.2f}\n")
    print("组合后的基本概率分配：")
    for a, v in sorted(combined.items(), key=lambda x: -x[1]):
        print(f"   m({sorted(a)}) = {v:.3f}")

    for name in ["感冒", "流感", "健康"]:
        sub = frozenset([name])
        bel = belief(combined, sub, theta)
        pl = plausibility(combined, sub, theta)
        print(f"{name}：信任区间 [{bel:.3f}, {pl:.3f}]")

"""模糊集合：隶属函数、交并补、模糊关系的 max-min 合成"""
import numpy as np

def trimf(x, a, b, c):
    """三角形隶属函数"""
    x = np.asarray(x, dtype=float)
    left = (x - a) / (b - a + 1e-12)
    right = (c - x) / (c - b + 1e-12)
    return np.clip(np.minimum(left, right), 0, 1)

def trapmf(x, a, b, c, d):
    """梯形隶属函数"""
    x = np.asarray(x, dtype=float)
    return np.clip(np.minimum((x - a) / (b - a + 1e-12),
                              (d - x) / (d - c + 1e-12)), 0, 1)

def gaussmf(x, mean, sigma):
    return np.exp(-((x - mean) / sigma) ** 2)

# ---- 模糊集合运算（Zadeh 算子）----
fuzzy_and  = np.minimum     # 交：隶属度取小
fuzzy_or   = np.maximum     # 并：隶属度取大
fuzzy_not  = lambda a: 1 - a

def fuzzy_compose(R, S):
    """模糊关系合成 R∘S：max-min 合成。"""
    R, S = np.asarray(R), np.asarray(S)
    n, m1 = R.shape
    m2, k = S.shape
    T = np.zeros((n, k))
    for i in range(n):
        for j in range(k):
            T[i, j] = max(min(R[i, t], S[t, j]) for t in range(m1))
    return T

if __name__ == "__main__":
    temp = np.arange(0, 41, 5)
    cold = trapmf(temp, 0, 0, 10, 20)
    hot  = trapmf(temp, 20, 30, 40, 40)
    print("温度论域：", list(temp))
    print("'冷'隶属度：", np.round(cold, 2))
    print("'热'隶属度：", np.round(hot, 2))
    print("20℃ 属于'不冷'的程度：", round(fuzzy_not(cold)[list(temp).index(20)], 2))

    # 模糊关系：R=人长相对身高 印象；S=身高-体重；合成=长相-体重
    R = [[0.8, 0.2], [0.3, 0.7]]
    S = [[0.6, 0.4], [0.1, 0.9]]
    print("\nmax-min 合成 R∘S =\n", fuzzy_compose(R, S))

"""Mamdani 模糊推理：服务/菜品质量 -> 小费比例（重心法去模糊化）"""
import numpy as np

# ---- 输入隶属函数（论域 0~10 分）----
def f_service_poor(x):   return np.clip((3 - x) / 3, 0, 1)
def f_service_good(x):   return np.where(x < 5, x / 5, (10 - x) / 5).clip(0, 1)
def f_service_excellent(x): return np.clip((x - 7) / 3, 0, 1)
def f_food_bad(x):       return np.clip((3 - x) / 3, 0, 1)
def f_food_great(x):     return np.clip((x - 7) / 3, 0, 1)

# ---- 输出论域：小费 0%~30% ----
tip = np.linspace(0, 30, 601)
def tip_cheap(t):  return np.clip((12 - t) / 12, 0, 1)
def tip_average(t):
    return np.where(t < 15, (t - 8) / 7, (22 - t) / 7).clip(0, 1)
def tip_generous(t): return np.clip((t - 18) / 12, 0, 1)

def fuzzy_tip(service, food):
    # 1) 模糊化：求各语言值的隶属度
    w_rules = [
        # (前件强度, 后件隶属函数)
        (min(f_service_poor(service), f_food_bad(food)), tip_cheap),       # 服务差且菜差 -> 小费低
        (f_service_excellent(service), tip_generous),                       # 服务极好 -> 小费高
        (f_service_good(service), tip_average),                             # 服务一般 -> 小费中
    ]
    # 2) Mamdani 推理：每条规则用前件强度 min-截断后件，再对各规则取 max 聚合
    aggregated = np.zeros_like(tip)
    for w, mf in w_rules:
        aggregated = np.maximum(aggregated, np.minimum(w, mf(tip)))
    # 3) 重心法（centroid）去模糊化
    if aggregated.sum() == 0:
        return 15.0
    return np.sum(tip * aggregated) / np.sum(aggregated)

if __name__ == "__main__":
    for s, f in [(2, 2), (5, 5), (8, 3), (10, 10), (3, 9)]:
        print(f"服务 {s:2d} 分、菜品 {f:2d} 分  →  建议小费 {fuzzy_tip(s, f):5.2f}%")