Compound Poisson Distribution Model

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.

Answer 1

[Shoaibur Rahman]

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))

[Yasmin Tamimi]

Thank a lot! Can you elaborate on the difference between out using data sample and the out using randperm?

[Shoaibur Rahman]

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.