Poisson random number generator

조회 수: 51 (최근 30일)
Ahmed raheem
Ahmed raheem 2012년 2월 6일
댓글: Torsten 2022년 1월 20일
Hi all please i need to know how to generate a Poisson distributed random variable without using the built-in function (poissrnd).

채택된 답변

Andreas Goser
Andreas Goser 2012년 2월 6일
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.

추가 답변 (8개)

Derek O'Connor
Derek O'Connor 2012년 2월 6일
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.
  댓글 수: 1
Ahmed raheem
Ahmed raheem 2012년 2월 7일
thank you for your help....
this one is quite helpful...

댓글을 달려면 로그인하십시오.


Derek O'Connor
Derek O'Connor 2012년 2월 7일
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

Derek O'Connor
Derek O'Connor 2012년 5월 3일
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
PhD Student
PhD Student 2019년 7월 10일
an End is missing
Ian Van Giesen
Ian Van Giesen 2020년 6월 24일
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
Derek O'Connor 2012년 7월 24일
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
John D'Errico
John D'Errico 2019년 7월 10일
Because the function name is spelled poissrnd, not poisrnd.
Binlin Wu
Binlin Wu 2021년 4월 1일
편집: Binlin Wu 2021년 4월 2일
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
Richard Willey 2012년 2월 6일
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.

Ahmed raheem
Ahmed raheem 2012년 2월 7일
First of all, I'd like to thank you all for your cooperation with me. Based on Dirk Kroese, I wrote the following code to generate the Poisson random variable:
The m-file code:
n=1;
lambda=500;
for i=1:10000
x=rand(1);
a=1;
a=a*x;
if a>=exp(-lambda)
n=n+1;
continue
else
X(i)=n-1;
end
end
  • would you please check it if it's correct or not? i feel it is not correct.

Derek O'Connor
Derek O'Connor 2012년 2월 7일
@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;
  댓글 수: 1
Ahmed raheem
Ahmed raheem 2012년 2월 7일
thank you Derek.
it is performing well now...
i putted the code inside a (for loop) and the result was ok.
thanks again.

댓글을 달려면 로그인하십시오.


Ahmed raheem
Ahmed raheem 2012년 2월 7일
This is the final code:
% -------------------------------------------------------------
function X = Poisson(lambda,n); % n represents the number of iterations
% -------------------------------------------------------------
% Generate a random value from the (discrete) Poisson
% distribution with parameter lambda.
% Derek O'Connor, 6 Feb 2012. derekroconnor@eircom.net
%
for i=1:n;
k=1; usave=1;
usave = usave*rand;
while usave >= exp(-lambda)
usave = usave*rand;
k = k+1;
end
X(i) = k;
end
hist(X)
  댓글 수: 2
Abdulramon Adeyiola
Abdulramon Adeyiola 2022년 1월 20일
@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
Torsten 2022년 1월 20일
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

댓글을 달려면 로그인하십시오.

태그

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by