matlab - 如何在没有统计工具箱的情况下从 Gamma 分布中绘制随机数？

Question

我正在改变合成图像的信号强度。我需要信号在 0 到 0.1 之间变化，但我需要使用 Gamma 分布来做到这一点，以便更多信号落在 0.01/.02 范围内。问题是我使用的是 2010 版的 Matlab，没有使用统计工具箱，该工具箱的库中没有 gamrnd 函数。

非常感谢任何和所有帮助。

Answer 1

您可以使用Inverse transform sampling method将均匀分布转换为任何其他分布:

P = rand(1000);
X = gaminv(P(:),2,2); % with k = 2 and theta = 2

这是一个小演示:

for k = [1 3 9]
    for theta = [0.5 1 2]
        X = gaminv(P(:),k,theta);
        histogram(X,50)
        hold on
    end
end

这给出:

---

编辑:

如果没有统计工具箱，可以使用Marsaglia's simple transformation-rejection method使用 rand 和 randn 从 Gamma 分布生成随机数:

N = 10000; % no. of tries
% distribution parameters:
a = 0.5;
b = 0.1;
% Marsaglia's simple transformation-rejection:
d = a - 1/3;
x = randn(N,1);
U = rand(N,1);
v = (1+x./sqrt(9*d)).^3;
accept = log(U)<(0.5*x.^2+d-d*v+d*log(v));
Y = d*(v(accept)).*b;

现在 Y 的分布类似于 gamma(a,b)。我们可以使用 gamrnd 函数测试结果:

n = size(Y,1);
X = gamrnd(a,b,n,1);

Y和X的直方图为:

但是，请记住， Gamma 分布可能不适合您的需求，因为它具有 no specific upper bound (即趋向无穷大)。因此，您可能想要使用另一个(有界)分布，例如 beta 除以 10。