M-H采样代码 MATLAB实现

M-H采样（Metropolis-Hastings sampling）是一种Markov Chain Monte Carlo（MCMC）方法，用于生成从一个给定分布中抽样的随机数序列。下面是 MATLAB 实现的示例代码：

% 目标分布密度函数
target_pdf = @(x) exp(-x.^2/2) .* (x>=0);

% 提议分布密度函数
proposal_pdf = @(x, sigma) normpdf(x, 0, sigma);

% 初始状态
x = 1;

% M-H采样
N = 10000;  % 采样次数
sigma = 1;  % 提议分布标准差
samples = zeros(N, 1);
for i = 1:N
    % 从提议分布中生成候选样本
    y = x + sigma * randn();
    
    % 接受概率
    alpha = min(target_pdf(y) / target_pdf(x), 1);
    
    % 决定是否接受候选样本
    if rand() < alpha
        x = y;
    end
    
    samples(i) = x;
end

% 绘制采样结果直方图
histogram(samples, 'Normalization', 'pdf');
hold on;
x = -5:0.1:5;
plot(x, target_pdf(x), 'LineWidth', 2);
hold off;

在上面的代码中，`target_pdf` 是目标分布的概率密度函数，`proposal_pdf` 是提议分布的概率密度函数。`x` 是初始状态，`sigma` 是提议分布的标准差。`N` 是采样次数，`samples` 是采样结果。在每一次采样中，从提议分布中生成一个候选样本 `y`，计算接受概率 `alpha`，并根据概率决定是否接受候选样本。最终，将所有采样结果绘制为直方图，并与目标分布比较。