拒绝采样MATLAB代码

拒绝采样是蒙特卡洛方法中的一种采样技术，用于从不易直接采样的目标分布中抽取样本。在拒绝采样中，我们首先选择一个容易采样的建议分布，这个分布需要比目标分布更加“宽松”，即在目标分布的所有值上都有较高的概率值。然后，从建议分布中抽取样本，并根据一定的拒绝率来决定是否接受这个样本，以确保最终样本符合目标分布。

以下是一个简单的拒绝采样MATLAB代码示例，目标分布是指数分布，而建议分布是均匀分布：

% 参数设置
target_density = @(x) exp(-x); % 目标分布的概率密度函数
proposal_density = @(x) ones(size(x)); % 建议分布的概率密度函数，这里用均匀分布
proposal_support = [0, 10]; % 建议分布的支持范围
n = 1000; % 希望抽取的样本数量

% 拒绝采样
accepted_samples = [];
while length(accepted_samples) < n
    % 从建议分布中抽取样本
    proposal_samples = proposal_support(1) + (proposal_support(2) - proposal_support(1)) * rand(1, n);
    
    % 计算目标分布与建议分布的概率密度比值
    ratio = target_density(proposal_samples) ./ proposal_density(proposal_samples);
    
    % 生成接受概率
    accept_prob = min(ratio, 1);
    
    % 根据接受概率决定是否接受样本
    accept = rand(1, n) < accept_prob;
    accepted_samples = [accepted_samples; proposal_samples(accept)];
end

% 移除空行
accepted_samples = accepted_samples';
accepted_samples = accepted_samples(~isempty(accepted_samples), :);

% 绘制结果
histogram(accepted_samples, 'Normalization', 'pdf');
hold on;
f = @(x) target_density(x);
x_values = linspace(min(accepted_samples), max(accepted_samples), 100);
plot(x_values, f(x_values), 'r', 'LineWidth', 2);
title('拒绝采样结果');
xlabel('样本值');
ylabel('概率密度');
legend('样本分布', '目标分布');

在这段代码中，我们首先定义了目标分布和建议分布的概率密度函数，然后通过循环从建议分布中抽取样本，并计算每个样本被接受的概率。如果样本被接受，就将其加入到最终样本集中。最后，我们绘制了样本的直方图和目标分布的曲线图，以可视化结果。