2017-01-31 27 views
2

Я изменяю силу сигнала для синтетических изображений. Мне нужен сигнал, чтобы он варьировался от 0 до 0,1, но мне нужно сделать это с помощью гамма-распределения, так что больше их приходится на диапазон 0,01/0,02. Проблема в том, что я использую версию Matlab 2010 без инструмента статистики, у которого нет функции gamrnd, являющейся частью ее библиотеки.Как нарисовать случайные числа из гамма-распределения без инструмента статистики?

Любая помощь очень ценится.

ответ

4

Вы можете использовать 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 

Что дает:

gamma dist


Edit:

Без инструментов статистики, вы можете использовать Marsaglia's simple transformation-rejection method для генерации случайных чисел из гамма-распределения с rand и randn:

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 распространяется как гамма (а, б) , Мы можем проверить результат с помощью gamrnd функции:

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

И гистограммы Y и X являются:

gamma dist 2

Однако, имейте в виду, что гамма-распределение может не соответствовать вашим потребностям, потому что он имеет no specific upper bound (т.е. переходит в бесконечность). Таким образом, вы можете использовать другое (ограниченное) распределение, например бета, деленное на 10.

+0

Эта версия Matlab не знает, что означает gaminv. – Kelsey

+0

Странно, в документах он говорит * «Введенный до R2006a» * ... – EBH

+0

@ Kelsey у вас есть ** Statistics Toolbox **? – EBH