基于MCMC算法的MLP与高斯过程贝叶斯推断工具箱

资 源 简 介

本项目是一个基于MATLAB开发的专业贝叶斯推断工具箱，核心依据马尔可夫链蒙特卡尔（MCMC）理论，旨在为多层感知机（MLP）和高斯过程（GP）提供强大的贝叶斯学习框架。该工具箱详细实现了多种高级采样算法，特别是混合蒙特卡尔（Hybrid Monte Carlo, HMC）算法和Metropolis-Hastings算法，用于在高维参数空间中高效地对后验概率分布进行采样。其主要功能包括：1. 贝叶斯神经网络训练，通过MCMC方法对神经网络的权重和偏置进行采样，从而代替传统的点估计，有效解决过拟合问题并提供预测的不确定性量化；2. 高斯过程模型的贝叶斯实现，支持对协方差函数的超参数进行采样，适用于复杂的回归和分类任务；3. 实现了自动相关性确定（ARD）技术，用于输入变量的特征选择和权重重要性分析；4. 提供完整的模型评估与诊断工具，包括采样链的收敛性检查、自相关分析以及预测分布的可视化。该工具适用于机器学习、模式识别及统计建模等领域的研究人员，帮助用户在小样本或含噪数据下构建鲁棒的概率模型。

详 情 说 明

MCMC Methods for MLP and GP V2.1

项目简介

本项目是一个基于 MATLAB 开发的专业贝叶斯推断工具箱，专为多层感知机（MLP）和高斯过程（GP）设计。核心依据马尔可夫链蒙特卡尔（MCMC）理论，旨在为回归任务提供强大的贝叶斯学习框架。该工具箱不依赖深度学习黑盒框架，而是通过底层代码详细实现了混合蒙特卡尔（HMC）算法和 Metropolis-Hastings 算法，用于在高维参数空间中高效地对后验概率分布进行采样，从而实现不确定性量化和鲁棒预测。

功能特性

• 贝叶斯多层感知机 (Bayesian MLP)： 使用混合蒙特卡尔 (HMC) 方法对神经网络权重进行采样，支持 Tanh 激活函数的单隐层结构。
• 贝叶斯高斯过程 (Bayesian GP)： 使用 Metropolis-Hastings算法对 GP 的超参数（长度尺度、信号方差、噪声方差）进行采样。
• 自动相关性确定 (ARD)： 代码框架中集成了 ARD 机制，允许对不同输入维度的重要性进行加权分析。
• 贝叶斯模型平均 (BMA)： 摒弃传统的点估计，通过对后验采样链的模型平均来进行预测，提供 95% 置信区间 (CI)。
• 底层算法实现： 包含自定义的数据标准化、梯度计算（反向传播）、势能计算及自相关分析工具。
• 可视化诊断： 提供完整的模型评估图表，包括预测拟合、能量收敛轨迹、参数后验分布及自相关性检查。

系统要求

• MATLAB: R2018b 或更高版本（代码主要依赖基础矩阵运算，无特定工具箱强依赖）。
• 内存: 建议 4GB 以上，取决于采样数量和网络规模。

使用方法

1. 确保工作路径包含本项目所在文件夹。
2. 直接运行 main 函数。
3. 程序将依次执行：合成数据生成 -> 贝叶斯 MLP 训练 (HMC) -> 贝叶斯 GP 训练 (MH) -> 结果可视化。
4. 运行结束后，控制台将输出迭代进度、接受率及能量值，并弹出一个包含 6 个子图的综合评估窗口。

---

main.m 核心逻辑详解

main.m 是整个项目的入口脚本，其执行流程严格遵循以下逻辑步骤：

1. 数据准备与预处理

程序首先根据设定的随机种子（rng(42)）生成合成回归数据集。数据生成函数为 $y = sin(3x) + 0.1x^2 + epsilon$，其中包含训练集和测试集（测试集范围比训练集更广以测试外推能力）。 关键在于数据标准化：代码手动实现了 Z-score 标准化 (zscore_custom)，分别对输入 $X$ 和目标 $y$ 进行处理。这是保证 MCMC 在参数空间高效收敛的关键步骤。

2. 贝叶斯 MLP 训练 (基于 HMC)

• 网络配置： 定义了一个结构为 Input -> 10 Hidden (Tanh) -> 1 Output 的神经网络。
• 参数初始化： 所有权重和偏置被扁平化为一个向量 theta。
• HMC 采样循环：

* 设定采样数（2000）和预烧期（500）。 * 循环调用 hmc_step（HMC 单步采样），该过程利用 Leapfrog 积分器模拟哈密顿动力学。 * 在计算势能时，结合了负对数似然（MSE 损失）和负对数先验（L2 正则化），其中输入层权重正则化项包含了 ARD 参数 alpha_ard。 * 虽然代码保留了 Gibbs 采样更新超参数的接口，但在当前的循环演示中，ARD 参数和噪声精度保持初始化状态，重点展示 HMC 对权重空间的探索。

• BMA 预测： 采样结束后，去除预烧期的样本，对剩余有效样本进行前向传播，计算所有模型的预测均值和方差，并反归一化还原到原始数据尺度。

3. 贝叶斯 GP 训练 (基于 Metropolis-Hastings)

• 参数空间： 采样对象为 GP 的超参数，构建为对数空间向量 gp_theta，包含：每个输入维度的 Length Scale (ARD)、Signal Std 和 Noise Std。
• MH 采样循环：

* 采用随机游走 Metropolis (Random Walk Metropolis) 策略，使用高斯分布作为提议分布。 * 计算接受率 alpha 并根据 Metropolis 准则决定是否接受新参数。 * 核心目标函数是高斯过程的对数边缘似然 (Log Marginal Likelihood)。

• 预测与不确定性： 通过稀疏提取采样链中的超参数，分别构建 GP 模型进行预测。最终的预测方差不仅包含了数据噪声，还融合了参数后验分布带来的认知不确定性（Total Variance = Expected Variance + Variance of Expectations）。

4. 结果可视化与评估

通过 figure 窗口展示 6 个维度的评估结果：

• 预测拟合图 (MLP & GP)： 展示真实数据点、模型预测均值曲线以及 95% 置信区间区域填充。
• 采样链诊断： 绘制 MLP 的势能 (Potential Energy) 轨迹，用于判断 HMC 是否达到稳态分布。
• 自相关分析 (Autocorrelation)： 对比 MLP 的能量链和 GP 的噪声参数链的自相关系数，评估混合效率。
• 参数分布： 绘制 GP 中 ARD Length Scale 参数的直方图。
• ARD 权重分析： 计算 MLP 输入层权重的均方根 (RMS)，直观展示输入特征的重要性。

---

关键算法与实现细节

1. 混合蒙特卡尔 (HMC) 势能与梯度 代码中的 mlp_grad_energy 函数不仅计算势能 $U(theta)$，还手动实现了反向传播算法来计算势能对参数的梯度 $nabla_theta U$。

• $U(theta) = -log P(D|theta) - log P(theta)$
• 似然项采用高斯噪声假设下的平方误差和 (SSE)。
• 先验项采用高斯先验（对应 L2 正则化），并引入了特征相关的超参数 alpha_ard，这使得梯度更新方向受到这两种力量的共同约束。

2. 自定义统计工具 为了保持代码的独立性，项目内嵌了 zscore_custom (标准化) 和 autocorr_custom (自相关计算) 函数，不依赖 MATLAB 统计工具箱的高级函数，展示了底层的计算逻辑。

3. 参数扁平化与重构 在 MLP 实现中，利用 unpack_weights 函数实现了从一维采样向量 theta 到神经网络矩阵形式 $(W1, b1, W2, b2)$ 的无缝转换，这是通用 MCMC 采样器适配神经网络结构的关键技巧。