matlab实现mcmc

MCMC（Markov Chain Monte Carlo）是一种常用的概率统计方法，可以用于模拟复杂的概率分布。在Matlab中，可以使用MCMC工具箱来实现MCMC算法。

MCMC工具箱提供了多种MCMC算法，包括Metropolis-Hastings算法、Gibbs采样算法等。其中，Metropolis-Hastings算法是最常用的MCMC算法之一，它通过构造一个马尔可夫链来模拟目标分布。具体实现过程如下：

1. 定义目标分布函数，即需要模拟的概率分布。
2. 初始化状态，即从目标分布中随机抽取一个样本作为初始状态。
3. 生成新状态，即从一个提议分布中随机抽取一个样本，与当前状态组合生成一个新状态。
4. 计算接受率，即计算接受新状态的概率。
5. 决定是否接受新状态，即根据接受率和一个随机数决定是否接受新状态。
6. 如果接受新状态，则将新状态作为下一次迭代的当前状态；否则，保留当前状态作为下一次迭代的当前状态。

下面是一个简单的Metropolis-Hastings算法的示例代码：

% 定义目标分布函数
target = @(x) exp(-x.^2/2) .* (x>=0);

% 定义提议分布函数
proposal = @(x) normrnd(x,1);

% 初始化状态
x = 0;

% 迭代次数
niter = 10000;

% 存储样本
samples = zeros(niter, 1);

% Metropolis-Hastings算法
for i = 1:niter
    % 生成新状态
    x_new = proposal(x);
    
    % 计算接受率
    alpha = min(1, target(x_new)/target(x));
    
    % 决定是否接受新状态
    if rand < alpha
        x = x_new;
    end
    
    % 存储样本
    samples(i) = x;
end

% 绘制样本分布
histogram(samples, 'Normalization', 'pdf');
hold on;
x = linspace(0, 5, 100);
plot(x, target(x), 'r', 'LineWidth', 2);

该代码实现了对目标分布 $f(x)=e^{-x^2/2}\cdot\mathbb{I}_{x\geq0}$ 的采样，其中 $\mathbb{I}$ 表示指示函数。在代码中，使用正态分布作为提议分布，即 $q(x'|x)=\mathcal{N}(x,1)$。最终，通过绘制样本分布和目标分布的比较，可以看出样本分布已经逼近目标分布。