贝叶斯网络分析：结构学习、跨尺度稳定性与条件推断 link

本章对应 spatial_structure_analysis/bayesian_network_analysis.py 与 bn_visual_inference.py 的当前实现。输入来自统计单元尺度或控制变量/多目标统计单元尺度生成的 stage_c_quadrat_dataset.csv；当前数据来源列表不直接读取邻域作用尺度批次。

1. 模块定位 link

贝叶斯网络（Bayesian Network, BN）用有向无环图（DAG）表示变量间的条件依赖结构。设节点为 $X_1,\ldots,X_p$，父节点集合为 $Pa(X_i)$，联合分布分解为：

$$ P(X_1,\ldots,X_p)=\prod_{i=1}^{p}P!\left(X_i\mid Pa(X_i)\right). $$

本模块回答：

• 在某一统计单元尺度下，哪些变量之间存在能提升 BIC 的有向依赖；
• 某个目标节点的父节点、子节点和 Markov Blanket 是什么；
• 给定证据状态后，目标节点的后验概率如何变化；
• 同一方向的边在多少个尺度重复出现；
• 哪些网络关系对尺度较稳定，哪些只在特定尺度出现。

它不能仅凭观测数据和自动结构学习证明严格因果关系。网络方向是当前离散化、变量集合、父节点约束、BIC 与搜索路径下的最佳统计方向，应结合理论、时间顺序、实验或准实验设计验证。

[1] 主界面与输入 link

1.1 页面总览 link

图 1-1　贝叶斯网络分析主界面：左侧选择来源、批次、尺度、变量和结构学习参数；右侧浏览数据概览、网络边、Markov Blanket、条件概率表、稳定性、综合报告与网络图推理。

1.2 可用来源 link

可选择上述层级目录或某个具体分析批次。扫描只显示包含可用阶段 C 数据的批次和尺度。

1.3 尺度与变量选择 link

变量候选取自所选尺度样本表的公共可计算字段。以下字段默认不作为网络节点：

• grid_row、grid_col、x、y、lon、lat 等空间索引；
• block_pixels、scale 等尺度/索引字段；
• 常量列、缺失率过高或不能数值/分类编码的字段。

至少选择 2 个变量。做跨尺度稳定性时建议 2–4 个以上尺度；只有 1 个尺度也能学习网络，但“稳定度”只能是 0 或 1，不具有真正的跨尺度比较意义。

变量角色根据字段前缀识别为主结构源、控制变量、响应目标或其他候选，并用于表格、图形颜色和解释；角色不对边方向施加强制约束。

[2] 离散化：连续变量如何变成状态 link

当前结构学习是离散贝叶斯网络。每个变量先变为整数状态 $0,1,\ldots,r_i-1$。

2.1 数值变量 link

缺失填补 link

数值缺失使用本列中位数：

$$ x_j^{fill}=\begin{cases} x_j,&x_j\text{ 有效},\ \operatorname{median}(x),&x_j\text{ 缺失}. \end{cases} $$

若整列无有效数值，所有样本编码为状态 0，并标为 empty_numeric；常量列也编码为状态 0，标为 constant_numeric。

唯一值不多于分箱数 link

若唯一值数量 $u\le B$，将排序后的唯一值逐一映射为状态：

$$ state(x)=\operatorname{index}_{0}\big(x;\operatorname{sort}(unique(x))\big). $$

该方式标为 ordinal_numeric，不会强行合并原有等级。

分位数分箱 link

若 $u>B$，代码先用 rank(method="first") 打破并列，再对秩做等频分箱：

$$ r_j=rank_{first}(x_j),qquad state(x_j)=Q_B(r_j), $$

其中 $Q_B$ 表示最多 $B$ 个分位箱，重复边界会被删除。状态标签保存各箱对应原始数值的最小—最大范围，而不是秩范围。该方式标为 quantile_numeric。

2.2 分类变量 link

分类值按出现频数从高到低排序。若类别数超过 $B$，保留前 $B-1$ 个高频类别，其余合并为 other：

$$ x_j’= \begin{cases} x_j,&x_j\in\mathcal C_{top(B-1)},\ \text{other},&\text{其他}. \end{cases} $$

缺失/文本空值归为 missing 语义状态，再按频次编码。

2.3 分箱数的影响 link

若节点 $X_i$ 有 $r_i$ 个状态，父节点组合数为 $q_i=\prod_{X_j\in Pa_i}r_j$。代码用于 BIC 惩罚的参数量为：

$$ k_i=\max!\left[1,(r_i-1)q_i\right]. $$

因此分箱数和最大父节点数会乘法增加参数量。小样本不应同时使用多分箱和多父节点。

[3] BIC 评分与 Hill-Climb 结构学习 link

3.1 局部对数似然 link

对节点 $X_i$，其状态 $k$、父配置 $j$ 的样本计数为 $N_{ijk}$，父配置总计 $N_{ij}=\sum_kN_{ijk}$。离散条件多项式的最大似然对数似然为：

$$ \ell_i=\sum_{j=1}^{q_i}\sum_{k=1}^{r_i}N_{ijk}\ln\frac{N_{ijk}}{N_{ij}}. $$

根节点没有父配置分组，此时：

$$ \ell_i=\sum_{k=1}^{r_i}N_{ik}\ln\frac{N_{ik}}{n}. $$

3.2 本实现的 BIC link

代码使用“越大越好”的 BIC 记号：

$$ BIC_i=\ell_i-\frac12k_i\ln n, \qquad BIC(G)=\sum_iBIC_i. $$

与某些文献采用 $-2\ell+k\ln n$ 且越小越好的定义只差线性变换；阅读本软件结果时必须按“总 BIC 越大越优”解释。

边操作增益为：

$$ \Delta BIC=BIC(G_{after})-BIC(G_{before}). $$

BIC增益 越大，表示接受该方向调整带来的局部评分改善越大；它不是效应大小或概率。

3.3 贪心搜索 link

每轮枚举可能的父—子对，并评估：

1. 添加边 $X_a\rightarrow X_b$；
2. 删除已有边；
3. 反转已有边 $X_a\rightarrow X_b$ 为 $X_b\rightarrow X_a$。

候选操作必须：

• 不形成有向环；
• 不超过子节点最大父节点数；
• $\Delta BIC$ 大于 min_improvement；
• 是本轮增益最大的操作。

接受后进入下一轮，直到没有满足阈值的操作或达到最大迭代数。伪数学表示为：

$$ G^{(t+1)}=arg\max_{G’\in\mathcal N(G^{(t)})}BIC(G’), $$

仅当 $BIC(G^{(t+1)})-BIC(G^{(t)})>\delta_{min}$ 时更新。

3.4 多重随机重启 link

第 1 次重启从空网络开始；后续重启随机生成满足 DAG 与父节点上限的初始网络。每尺度重启次数为：

$$ R=\max!\left(1,\left\lceil\frac{W_{request}}{L_{valid}}\right\rceil\right), $$

$W_{request}$ 为请求并行任务数，$L_{valid}$ 为有效尺度数。每个尺度从有效重启中选择总 BIC 最大的结果：

$$ G_s^*=\arg\max_{r=1,\ldots,R}BIC(G_{s,r}). $$

并行任务数大于尺度数时，额外资源会用于单尺度多重重启，而不是闲置。随机重启减轻局部最优，但不能保证找到全局最优 DAG。

[4] 结构学习参数 link

4.1 自动配置规则 link

自动配置根据尺度数 $L$、变量数 $p$、各尺度最小样本 $n_{min}$ 和复杂度比 $c=p/n_{min}$ 设置参数。

分箱数 link

$$ B= \begin{cases} 3,&n_{min}<40\text{ 或 }c>0.60,\ 5,&n_{min}\ge250\text{ 且 }p\le24,\ 4,&\text{其他}. \end{cases} $$

当前代码在 5 箱条件之后还写有“$n_{min}\ge500$ 且 $p\le16$ 时设为 6”的 elif，但该样本会先满足前面的 5 箱条件，因此自动配置实际上不会到达 6 箱分支；需要 6 箱时应在自动配置后手工设置。最终值仍会限制在界面 2–12 的范围内。

最大父节点数 link

基础为 2；满足更充足的样本、更少变量和更低复杂度时依次提高到 3、4 或 5，并且不超过 $p-1$。

最大迭代 link

$$ I=30+5p+12L+18P_{max}, $$

若 $n_{min}<50$，先限制到不超过 120，最终截断到 1–500。

最小改进阈值 link

$$ \delta_{min}= \begin{cases} 0.05,&n_{min}<40,\ 0.02,&40\le n_{min}<120,\ 0.01,&120\le n_{min}<300,\ 0.005,&n_{min}\ge300. \end{cases} $$

稳定边阈值 link

自动配置是基于规模的保护规则，不替代研究者对变量语义、状态数和网络复杂度的判断。

[5] 网络边 link

图 5-1　网络边结果：按尺度列出父节点、子节点、双方角色与 BIC 增益；可调用 AI 对当前标签进行辅助解释。

结果列解释：

边方向的正确表述 link

• 可以说：“在该尺度和变量集合下，BIC 搜索选择了 $X\rightarrow Y$ 的条件依赖方向。”
• 不应仅凭该边说：“改变 $X$ 一定导致 $Y$ 改变。”
• 同一统计等价类可能包含多个难以由观测数据区分的 DAG；未纳入变量也会改变方向。

[6] 条件概率表（CPT） link

6.1 最大似然频率估计 link

对父配置 $pa_j$：

$$ \widehat P(X_i=k\mid Pa_i=pa_j)=\frac{N_{ijk}}{N_{ij}}. $$

根节点为：

$$ \widehat P(X_i=k)=\frac{N_{ik}}{n}. $$

当前代码不做 Laplace/Dirichlet 平滑，只导出观测到的状态/父配置组合。因此小样本、多分箱或多父节点会产生稀疏 CPT；未观测组合不能被理解为稳定的真实零概率。

图 6-1　条件概率表：每行给出尺度、节点、父节点状态配置、节点状态、计数和观测概率。

读表示例 link

若一行显示 父节点配置: A=2, B=0、状态: 1、概率: 0.475，表示：

$$ \widehat P(X=1\mid A=2,B=0)=0.475. $$

状态编号必须结合 bn_discretization.csv 的“状态标签/状态范围”翻译回原始数值区间，不能把 state=2 直接解释为原变量值 2。

[7] Markov Blanket link

对节点 $X$，Markov Blanket 由父节点、子节点以及子节点的其他父节点（配偶节点）组成：

$$ MB(X)=Pa(X)\cup Ch(X)\cup\bigcup_{Y\in Ch(X)}Pa(Y)\setminus{X}. $$

在贝叶斯网络成立的前提下：

$$ X\perp V\setminus({X}\cup MB(X))\mid MB(X). $$

即给定 Blanket 后，网络中其他节点不再为 $X$ 提供额外条件信息。

图 7-1　Markov Blanket 结果：列出每个节点的父节点、子节点、完整 Blanket 与 Blanket 大小。

7.1 解释价值 link

• 目标特征筛选：目标的 Blanket 是预测该目标的局部候选集合；
• 控制路径识别：控制变量若进入响应目标 Blanket，说明其处于直接条件依赖邻域；
• 网络复杂度：平均 Blanket 很大可能表示网络过密、变量冗余或分箱/父节点上限过高；
• 跨尺度比较：同一目标 Blanket 在不同尺度变化，提示条件依赖邻域具有尺度特异性。

Blanket 是当前模型中的条件独立结构，不是保证充分无偏的因果调整集。

[8] 跨尺度稳定性 link

对有向边 $e=(X\rightarrow Y)$，若在 $L$ 个成功分析尺度中出现于 $L_e$ 个尺度，方向稳定度为：

$$ Stab(e)=\frac{L_e}{L}. $$

当：

$$ Stab(e)\ge\tau_{stable} $$

时标记为“稳定边”。反向边 $Y\rightarrow X$ 被视为另一条边，不计入 $X\rightarrow Y$ 的支持。

图 8-1　跨尺度稳定性：显示有向边出现尺度数、方向稳定度、稳定边判断及具体尺度目录。

8.1 如何解释稳定度 link

稳定度不包含每个尺度的 BIC 增益大小，也不衡量边方向是否因果正确。最好联合报告稳定度、出现尺度、各尺度 BIC、样本数和变量状态分布。

[9] 网络图与精确条件推断 link

9.1 图模式 link

界面还可选择分层/力导布局、全图/焦点邻域/目标 Blanket/关键节点子图、标签密度和 12–600 的节点显示上限。

9.2 变量消元 link

给定证据 $E=e$ 和查询节点 $Q$：

$$ P(Q\mid e)=\frac{1}{Z}\sum_{\mathbf H}\prod_{i=1}^{p}P(X_i\mid Pa_i), $$

$\mathbf H=V\setminus({Q}\cup E)$ 为需消元的隐藏变量，$Z$ 为归一化常数。

代码执行：

1. 用证据限制每个 CPT 因子；
2. 构造因子图并用 min-fill 启发式选择消元顺序；
3. 相乘所有包含待消元变量的因子；
4. 对该变量求和；
5. 最后保留查询节点并归一化。

因子乘积与求和分别为：

$$ (\phi_a\phi_b)(\mathbf z)=\phi_a(\mathbf z_a)\phi_b(\mathbf z_b), \qquad \phi’(\mathbf z_{-X})=\sum_X\phi(\mathbf z_{-X},X). $$

若稀疏 CPT 和证据导致无法形成正概率分布，当前实现退回查询状态的均匀分布；这应视为数据/模型可支持性警告，而不是实质推断结论。

9.3 先验—后验变化 link

目标状态 $k$ 的变化为：

$$ \Delta P_k=P(Q=k\mid e)-P(Q=k). $$

全网节点影响排序使用最大绝对概率变化：

$$ \Delta_{max}(X)=\max_k|P(X=k\mid e)-P(X=k)|, $$

以及 Jensen–Shannon 散度。令 $P$ 为先验、$Q$ 为后验、$M=(P+Q)/2$：

$$ JS(P,Q)=\frac12KL_2(P|M)+\frac12KL_2(Q|M), $$

$$ KL_2(P|M)=\sum_{k:P_k>0,M_k>0}P_k\log_2\frac{P_k}{M_k}. $$

JS 越大表示整个状态分布改变越明显；$\Delta_{max}$ 强调变化最大的单一状态。

图 9-1　网络图推理：左侧显示单尺度概率网络，右侧设置目标与证据状态、查看目标先验/后验、全网影响排序和推断摘要。

9.4 操作步骤 link

1. 切换到“单尺度概率网络”并选择尺度。
2. 选择目标节点。
3. 选择证据节点和离散状态；状态标签显示原始数值区间。
4. 点击“加入证据”，可加入多个节点—状态组合。
5. 点击“执行条件推断”。
6. 查看目标节点的先验/后验概率与 $\Delta P$。
7. 查看全网影响排序，选中节点可使网络图聚焦。
8. 按需导出自包含 HTML 或推断摘要 Markdown。

证据表示“条件假设该节点处于某离散状态”，不是软件对原始栅格实施干预。

[10] 综合报告 link

图 10-1　贝叶斯网络综合报告：汇总尺度数、变量数、总边数、稳定边、关键关系及可选 AI 详细解释。

报告的程序生成部分包括：

• 分析批次、成功尺度数、变量数、总边数与稳定边数量；
• 离散化、Hill-Climb+BIC、CPT 频率估计和稳定度定义；
• 达到当前阈值的关键稳定边；
• 每尺度样本数、节点数、边数和平均 Blanket 大小。

AI 解读可以分别作用于网络边、Blanket、CPT、跨尺度稳定性、综合报告和网络图推理。AI 只读取已生成结果与方法上下文，不重新估计网络；任何 AI 给出的因果或政策建议都应人工核对。

[11] 输出文件 link

11.1 分析根目录 link

11.2 每尺度子目录 link

网络图 HTML 和推断摘要 Markdown 由用户在“网络图推理”页另行导出。临时重启快照只用于恢复，正式汇总成功后会清理。

[12] 推荐分析流程 link

1. 在前置统计单元模块检查 stage_c_quadrat_dataset.csv，确保变量语义和样本数可靠。
2. 选择正确来源类型和批次，至少勾选 2–4 个具有代表性的尺度。
3. 不要一开始勾选数百个高度冗余分量；优先按研究问题筛选主结构、关键控制和响应变量。
4. 点击“自动配置参数”作为起点，再根据样本/状态数检查 CPT 参数量。
5. 小样本使用 3 箱、2 个父节点和较高最小改进阈值；样本充足再逐步增加复杂度。
6. 运行后先看数据概览和总 BIC，再检查网络边、Blanket 和 CPT 稀疏性。
7. 做跨尺度比较时关注同方向边、目标 Blanket 和网络密度是否稳定。
8. 在单尺度网络中设置有业务意义且样本支持的证据状态，检查先验—后验变化。
9. 用不同分箱数、父节点上限和变量子集做敏感性分析。
10. 报告中明确写出：观测性结构学习、离散化规则、BIC 定义、重启次数和稳定阈值。

[13] 常见问题与解释边界 link

贝叶斯网络最有价值的用途是提出可检验的机制假设、识别目标的局部条件依赖邻域和比较结构是否跨尺度复现，而不是把自动学习的箭头直接当作最终因果结论。