如何基于给定二维分布生成伪随机数？逆变换法遇阻求方案

一、二维逆变换法（适配你已有的离散联合CDF） #

你已经计算得到离散化的联合累积分布函数funzione_cumulata，可通过以下步骤实现二维逆变换抽样：

1. 生成M个[0,1]区间的均匀随机数u = rand(M,1)；
2. 对每个u，在二维CDF网格中找到满足funzione_cumulata(i,j) >= u的最小索引(i,j)；
3. 对应网格点的(x(j), y(i))即为样本，若需要更精细结果，可在对应网格单元内做均匀插值。

实现代码 #

% 基于已有联合CDF生成样本
u = rand(M, 1);
samples = zeros(M, 2);

% 展平CDF并记录对应坐标，加速查找
cdf_flat = funzione_cumulata(:);
[x_grid, y_grid] = meshgrid(x, y);
x_flat = x_grid(:);
y_flat = y_grid(:);

for k = 1:M
    % 找到第一个大于等于u(k)的索引
    idx = find(cdf_flat >= u(k), 1, 'first');
    samples(k, :) = [x_flat(idx), y_flat(idx)];
end

二、接受-拒绝抽样法（更高效的替代方案） #

该方法无需依赖联合CDF，针对离散网格形式的PDF实现直观：

1. 找到PDF的最大值max_pdf = max(densita_di_prob(:))；
2. 循环生成候选点：
   • 在[-a,a]×[-a,a]生成均匀随机点x_cand = -a + 2*a*rand()，y_cand = -a + 2*a*rand()；
   • 定位候选点所在网格索引，获取对应PDF值pdf_val；
   • 生成[0, max_pdf]的均匀随机数v = max_pdf*rand()；
   • 若v <= pdf_val则接受该点，否则丢弃；
3. 重复直到收集到M个样本。

实现代码 #

max_pdf = max(densita_di_prob(:));
samples = zeros(M, 2);
count = 0;

while count < M
    % 生成候选点
    x_cand = -a + 2*a*rand();
    y_cand = -a + 2*a*rand();
    
    % 定位候选点对应的网格索引
    idx_x = find(x >= x_cand, 1, 'first');
    idx_y = find(y >= y_cand, 1, 'first');
    idx_x = max(1, min(idx_x, N));
    idx_y = max(1, min(idx_y, N));
    pdf_val = densita_di_prob(idx_y, idx_x);
    
    % 接受-拒绝判断
    if max_pdf*rand() <= pdf_val
        count = count + 1;
        samples(count, :) = [x_cand, y_cand];
    end
end

三、边缘分布+条件分布抽样法 #

利用概率链式法则，分步抽样x和y：

1. 计算x的边缘分布：对PDF按y轴求和得到pdf_x = sum(densita_di_prob, 1)，归一化后计算其CDFcdf_x = cumsum(pdf_x)/sum(pdf_x)；
2. 用一维逆变换生成x样本；
3. 对每个x样本，计算给定x下y的条件分布pdf_y_given_x = densita_di_prob(:, idx_x)/sum(densita_di_prob(:, idx_x))，再用一维逆变换生成y样本。

实现代码 #

% 计算x的边缘分布及CDF
pdf_x = sum(densita_di_prob, 1);
pdf_x = pdf_x / sum(pdf_x);
cdf_x = cumsum(pdf_x);

samples = zeros(M, 2);
for k = 1:M
    % 生成x样本
    u_x = rand();
    idx_x = find(cdf_x >= u_x, 1, 'first');
    x_sample = x(idx_x);
    
    % 计算给定x的y条件分布及CDF
    pdf_y_given_x = densita_di_prob(:, idx_x);
    pdf_y_given_x = pdf_y_given_x / sum(pdf_y_given_x);
    cdf_y_given_x = cumsum(pdf_y_given_x);
    
    % 生成y样本
    u_y = rand();
    idx_y = find(cdf_y_given_x >= u_y, 1, 'first');
    y_sample = y(idx_y);
    
    samples(k, :) = [x_sample, y_sample];
end

额外优化提示 #

你的PDF本质是权重为2:1的两个二维高斯混合分布，可直接用高斯抽样函数生成，效率远高于离散网格抽样：

% 直接生成混合高斯样本（最高效）
samples = zeros(M, 2);
for k = 1:M
    % 按权重选择高斯分量
    if rand() < 2/3
        samples(k, :) = mu1 + sigma*randn(1,2);
    else
        samples(k, :) = mu2 + sigma*randn(1,2);
    end
end
% 可选：截断到[-a,a]范围（与原网格范围一致）
samples = max(min(samples, a), -a);

内容的提问来源于stack exchange，提问作者ETTORE BIOLI