Compound Poisson Distribution Model

8 visualizzazioni (ultimi 30 giorni)
Hey Everyone,
I have the probability density function of a Negative-binomial Distribution (Compound Poisson Distribution) and I would like to generate random numbers based on probability similar to how binornd function works. The following is the compound model written in matlab:
fun = @(lambda) (lambda.^k).*(exp(-1.*lambda)).* gampdf(lambda,alpha,beta)./factorial(k);
P(k)= integral(fun,0,Inf);
I assume the number of trials (k) from 0:7 and the output is the probability of each point P(k). So now I only need to choose random numbers from this distribution based on their probability values.
Thanks in advance.

Risposta accettata

Shoaibur Rahman
Shoaibur Rahman il 14 Gen 2015
I am not sure if I understand your question correctly. I guess once you get P, then you would like to choose its elements randomly, right? If so, you can use the following code:
alpha = 1; beta = 1; % replace by your desired values
for k = 0:7
fun = @(lambda) (lambda.^k).*(exp(-1.*lambda)).* gampdf(lambda,alpha,beta)./factorial(k);
P(k+1) = integral(fun,0,Inf);
end
out = datasample(P,1) % 1 observation, replace 1 by higher value for more observations
This will be uniformly random. If you are interested then you may use:
out = P(randperm(numel(P),1))
  2 Commenti
Yasmin Tamimi
Yasmin Tamimi il 14 Gen 2015
Thank a lot! Can you elaborate on the difference between out using data sample and the out using randperm?
Shoaibur Rahman
Shoaibur Rahman il 14 Gen 2015
Both datasample and randperm take random numbers from uniform random distribution. Further, datasample uses rng and randperm uses rand function to generate random numbers, indices in this case. However, rng is more 'predictable' than rand, like you can set a previous generation stream by using rng. Lets do an example:
gen_setting = rng;
r1 = rand(1,5)
rng(gen_setting)
r2 = rand(1,5)
You will see that both r1 and r2 are same. You can also see the docs for datasample, randperm, rng and rand to explore this a bit more.

Accedi per commentare.

Più risposte (0)

Community Treasure Hunt

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

Start Hunting!

Translated by