Poisson random number generator
7 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
Hi all please i need to know how to generate a Poisson distributed random variable without using the built-in function (poissrnd).
0 commentaires
Réponse acceptée
Andreas Goser
le 6 Fév 2012
If this is an acadamic exercise - you can look at the literature refererence
% References:
% [1] Devroye, L. (1986) Non-Uniform Random Variate Generation,
% Springer-Verlag.
Plus de réponses (8)
Derek O'Connor
le 6 Fév 2012
Dirk Kroese has excellent notes here: http://www.maths.uq.edu.au/~kroese/mccourse.pdf, which are based on his book:
D.P. Kroese, T. Taimre, Z.I. Botev: Handbook of Monte Carlo Methods. John Wiley & Sons, 2011.
His notes and book have lots of Matlab examples.
Derek O'Connor
le 7 Fév 2012
I prefer this:
% -------------------------------------------------------------
function S = PoissonSamp(lambda,ns);
% -------------------------------------------------------------
% Generate a random sample S of size ns from the (discrete)
% Poisson distribution with parameter lambda.
% Derek O'Connor, 6 Feb 2012. derekroconnor@eircom.net
%
S = zeros(ns,1);
for i = 1:ns
k=1; produ = 1;
produ = produ*rand;
while produ >= exp(-lambda)
produ = produ*rand;
k = k+1;
end
S(i) = k;
end
0 commentaires
Derek O'Connor
le 3 Mai 2012
I would like to thank Kang Wook Lee of Berkeley for pointing out an error in the code above. The last line should be S(i) = k-1;
% -------------------------------------------------------------
function S = PoissonSamp(lambda,ns);
% -------------------------------------------------------------
% Generate a random sample S of size ns from the (discrete)
% Poisson distribution with parameter lambda.
% Fixed error: changed S(i) = k; to S(i) = k-1;
% Derek O'Connor, 3 May 2012. derekroconnor@eircom.net
%
S = zeros(ns,1);
for i = 1:ns
k=1; produ = 1;
produ = produ*rand;
while produ >= exp(-lambda)
produ = produ*rand;
k = k+1;
end
S(i) = k-1;
end
2 commentaires
Ian Van Giesen
le 24 Juin 2020
Why was the last line changed? Was it to make the events where probability for 'success' a null event? Sorry in advance for a confusing question, still wrapping my head around this, but many thanks for the code!
Derek O'Connor
le 24 Juil 2012
This is a cleaner fix of PoissonSamp
% -------------------------------------------------------------
function S = PoissonSamp3(lambda,ns);
% -------------------------------------------------------------
% Generate a random sample S of size ns from the (discrete)
% Poisson distribution with parameter lambda.
% Fixed error:
% CHANGED k = 1; produ = 1; produ = produ*rand
% TO k = 0; produ = rand;
% Derek O'Connor, 24 July 2012. derekroconnor@eircom.net
%
S = zeros(ns,1);
for i = 1:ns
k = 0;
produ = rand;
while produ >= exp(-lambda)
produ = produ*rand;
k = k+1;
end
S(i) = k;
end
3 commentaires
Binlin Wu
le 1 Avr 2021
Modifié(e) : Binlin Wu
le 2 Avr 2021
Would just like to make some minor changes to accept any array size:
function S = PoissonSamp3(lambda,varargin);
% -------------------------------------------------------------
% Generate a random sample S of size ns from the (discrete)
% Poisson distribution with parameter lambda.
% Fixed error:
% CHANGED k = 1; produ = 1; produ = produ*rand
% TO k = 0; produ = rand;
% Derek O'Connor, 24 July 2012. derekroconnor@eircom.net
%
nn = [varargin{:}];
S = zeros([nn,1]);
for i = 1:numel(S(:))
k = 0;
produ = rand;
while produ >= exp(-lambda)
produ = produ*rand;
k = k+1;
end
S(i) = k;
end
Richard Willey
le 6 Fév 2012
Mark Steyvers has written a nice book titled "Computational Statistics with MATLAB" which can be downloaded from
The first chapter has a very good section describing inverse transform sampling which provides everything you need to know.
0 commentaires
Derek O'Connor
le 7 Fév 2012
@Ahmed, you're correct, it is not correct.
The function below is a Matlab translation of Kroese's algorithm. It seems to work ok but needs to be thoroughly tested. You should do this and let us know the results. Note that Poisson(L) ~ Norm(L,L), for large L.
The one thing I don't like about Kroese is his awful algorithm and programming style. Why does he use GOTOs instead of proper WHILEs etc?
% -------------------------------------------------------------
function X = Poisson(lambda);
% -------------------------------------------------------------
% Generate a random value from the (discrete) Poisson
% distribution with parameter lambda.
% Derek O'Connor, 6 Feb 2012. derekroconnor@eircom.net
%
k=1; produ = 1;
produ = produ*rand;
while produ >= exp(-lambda)
produ = produ*rand;
k = k+1;
end
X = k;
Ahmed raheem
le 7 Fév 2012
2 commentaires
Abdulramon Adeyiola
le 20 Jan 2022
@Ahmed raheem I tried using your code to generate 100 samples from Poisson distribution with parameter 1, but the mean of the resulting samples was way above 1. This is strange!
Torsten
le 20 Jan 2022
L = exp(-lambda);
for i=1:n
k=0; usave=1;
while usave > L
k = k+1;
usave = usave*rand;
end
X(i) = k-1;
end
Voir également
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!