效果一览
代码功能
代码功能简述
目标:实现贝叶斯多项式曲线拟合,动态展示随着数据点逐步增加,模型后验分布的更新过程。
核心步骤:
数据生成:在区间[0,1]生成带噪声的正弦曲线作为训练数据。
参数设置:定义先验精度(alpha、beta)和多项式阶数(D=8)。
贝叶斯更新:每次迭代使用前j个数据点计算后验分布的均值和方差。
预测分布:对密集采样的X值计算预测均值和置信区间(±0.5σ)。
可视化:
实时绘制数据点、真实曲线、预测均值曲线和置信区间。
支持生成GIF动画,展示拟合过程的动态变化。
数学基础:基于贝叶斯线性回归公式,通过设计矩阵和协方差更新实现高效计算。
关键说明
bsxfun(@power, x2’, 0:D)':生成多项式特征矩阵(例如,x^0, x^1, …, x^D)。
Sinv:权重后验分布的协方差矩阵的逆,结合了先验(aI)和似然项(bφφ^T)。
m(k)和s(k):分别对应预测分布的均值和方差,体现贝叶斯模型对不确定性的量化。
置信区间填充:使用品红色半透明区域表示预测分布的范围,直观展示模型的不确定性随数据增加而降低的过程。
部分代码
clear; % 清空工作区变量
a = 5e-3; % alpha: 权重w的先验精度(控制模型复杂度)
b = 11.1; % beta: 输入噪声的精度(控制数据拟合程度)
D = 8; % 多项式阶数
M = D+1; % 模型参数数量(从0次到D次)
W = 0.5*[1 -1]; % 预测分布的标准差倍数(用于置信区间)
% 随机打乱训练数据顺序(展示算法对数据顺序的鲁棒性)
ind = randperm(N);
t = t(ind); x = x(ind);
% 用户输入是否生成GIF动画
z = input('Make a GIF? (y/n): ','s');
if z == 'y', filename = 'bayesiancurvefit.gif'; end
clc;
%% 主循环:逐步添加数据点并更新贝叶斯模型
for j = 1:N
%% 贝叶斯模型更新
t2 = t(1:j); x2 = x(1:j); % 取前j个数据点
phi = bsxfun(@power, x2', 0:D)'; % 设计矩阵(多项式特征,维度[D+1 x j])
Sinv = a*eye(M) + b*(phi*phi'); % 权重后验协方差的逆(式1.72)
PHI = bsxfun(@power, X', 0:D)'; % 密集采样的X对应的设计矩阵
% 计算预测分布的均值和方差
s = zeros(size(X)); m = s;
for k = 1:length(X)
% 均值计算(式1.70)
m(k) = b * PHI(:,k)' * (Sinv \ (phi * t2'));
% 方差计算(式1.71,包含噪声项1/b)
s(k) = 1/b + PHI(:,k)' * (Sinv \ PHI(:,k));
end
%% 可视化结果
f = figure(1);
clf; scatter(x2, t2, 'b', 'filled'); % 绘制当前数据点
hold on; grid on; box on;
axis([-0.1 1.1 -1.5 1.5]); % 固定坐标轴范围
% 绘制真实曲线(绿色)、预测均值(红色)
plot(X, T, 'g', X, m, 'r-', 'LineWidth', 2);