MATLAB: Exponentially distributed random number

exponential distributionexprndrandom number

I want to generate exponentially distributed randm numbers (with mean=lambda) untill the summation of generated numbers is less or equal than T_design=75. thereby, I need to save the random numbers.

Best Answer

Wow. Its been 40 years since I did this sort of stuff. ;-)
I assume you have a stochastic process that you wish to model, as that is what it looks like. So you have some (Poisson) process with exponential inter-arrival time, with rate parameter lambda. You want to generate a set of such arrivals over a span of time T_design.
The simple solution is to just generate those samples, one at a time until you exceed T_design, using a loop.
Be careful here though, because if you use exprnd to generate the samples, exprnd defines the parameter as the mean. Some exponential distribution definitions use the inverse of that. For example, if you read wikipedia, it uses the inverse mean definition of an exponential distribution. (Sigh.)
However, a simple code might look like this:
T_design = 75;
lambda = 3; % I just picked an exponential rate parameter out of a hat

arrivals = [];
while true
  newarrival = exprnd(lambda,1);
  if sum(arrivals) + newarrival < T_design
    arrivals(end+1) = newarrival;
  else
    break
  end
end
So a very simple scheme, that grows a list of arrivals, until you get one that pushes you over the end point.
arrivals
arrivals =
  Columns 1 through 11
       2.6221        9.008      0.30704      0.17039       2.1348       2.1446       3.2566       0.3159       2.9889       6.5892      0.74441
  Columns 12 through 22
       2.8268       4.2603       2.7197       7.0161       6.0755      0.17909      0.13457       1.6591       8.4513       4.3473       3.1225
  Column 23
      0.59099
      
sum(arrivals)
ans =
       71.665
The only flaw to this is it forces you to dynamically grow the array arrivals. That is inefficient. So, perhaps better is to preallocate arrivals to some reasonable size. I expect to see roughly T_design/lambda events over that time span. If I pre-allocate arrivals to be 50% larger than the expected number, I should normally be ok. (We can do better than that estimate, as a function of lambda, of course.) However, this code should be reasonably efficient:
T_design = 75;
lambda = 3; % I just picked an exponential rate parameter out of a hat
arrivals = NaN(1,ceil(1.5*T_design/lambda));
counter = 0;
sumarrivals = 0;
while true
  newarrival = exprnd(lambda,1);
  if sumarrivals + newarrival < T_design
    sumarrivals = sumarrivals + newarrival;
    counter = counter + 1;
    arrivals(counter) = newarrival;
  else
    arrivals = arrivals(1:counter);
    break
  end
end
Can we do better? Of course....
Related Question