如何用倒数法生成数

How to generate number by inversion method

我有一个计算范围内数字的算法,例如

  1. 计算 F(x),x 从 0 到 n
  2. 生成统一数z If(F(x)<=z<F(x+1)) 然后 k=x

其中,F(x)是由Binomial_distribution计算得到的累积分布函数。

例如,F(x)=0.12, 0.2, 0.5, ..1.

z=0.3 然后 return k=1

我是用matlab实现的。然而,它 returned 错误的结果。你能看到我的实施对我有帮助吗?如有错误,请给出解决方案

clear all;
number_generated = 500000;
p=0.3;
n=20
%% Compute F(x): x=0..20
sum_F=0;
for x=0:n
        sum_F=sum_F+nchoosek(n,x)*p^x*(1-p)^(n-x);
        F(x+1)=sum_F;% matlab index from 1
end    
k_arr=[]; %% Store k
for num=1:number_generated
    %% Generate z_i uniform in the interval (0,1)
    z=rand();
    %% Find k such that F(x)<=z<F(x+1)
    for i=1:length(F)-1 
        if(F(i)<=z & z<F(i+1))
            k=i-1; % matlab index from 1
            break;
        end
    end    
    k_arr=[k_arr k]; %% Record k in array
end

您的代码基本上是正确的,但是您忘记了 0F(1) 之间的第一个间隔。

对您的代码的其他评论:

  1. 总是预分配数组

  2. 在逻辑条件下使用&&代替&

更新代码:

clear all;
close all;
number_generated = 500000;
p=0.3;
n=20;
% Preallocate F array.
F = NaN(n + 2, 1);
% Set first value to 0
F(1) = 0;
% Save binomial distribution for plotting.
bin = zeros(n + 1, 1);
% Compute F(x): x=0..20
sum_F=0;
for x=0:n
    % Save binomial distribution for plotting.
    bin(x + 1) = nchoosek(n,x)*p^x*(1-p)^(n-x);
    sum_F=sum_F+bin(x + 1);
    % This is now x+2 because F(1) is 0.
    F(x+2)=sum_F;% matlab index from 1
end
% Preallocate k_arr array.
k_arr=NaN(number_generated, 1); %% Store k
for num=1:number_generated
    % Generate z_i uniform in the interval (0,1)
    z=rand();
    % Find k such that F(x)<=z<F(x+1)
    for i=1:length(F)-1 
        if(z >= F(i) && z < F(i+1))
            k=i-1; % matlab index from 1
            break;
        end
    end    
    k_arr(num) = k; % Record k in array
end

% Plot expected result
figure
stairs((0:n) - 0.5, bin);
xlim([-1 20]);

% Plot sampled result
figure
histogram(k_arr);
xlim([-1 20]);

预期分布:

样本分布: