exponentially distributed random number

7 次查看(过去 30 天)
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.

采纳的回答

John D'Errico
John D'Errico 2020-1-28
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....
  2 个评论
Mos_bad
Mos_bad 2020-1-28
Thank you so much. yeah it is for poisson process. Now I got one more question. I have a while loop whithin a for loop which iterates the lambda values. How I can save the final out put of the nested loops including different generated exponentiall random realization based on different lambda values.
What I am trying to save is "arrivals" for different "m_int" or "lambda_m_int" values.
Here is the code:
M_min=4.5; M_max=8.0;
m=M_min:0.0001:M_max;
T_design=75; %years
a=4.56; b=1.0;
alpha=a*log(10);beta=b*log(10);
syms m
F_M=(1-exp(-beta*(m-M_min))) / (1-(exp(-beta*(M_max-M_min)))); % CDF of Mag.
lambda_m=1-F_M; %annual rat
m_int=4.6:0.6:M_max;
for mm=1:length(m_int)
m_int_=m_int(mm);
lambda_m_int = double(subs(lambda_m, m,m_int_));
arrivals = [];
while true
newarrival = exprnd(1/lambda_m_int,1);
if sum(arrivals) + newarrival < T_design
arrivals(end+1) = newarrival;
else
break
end
end
end
John D'Errico
John D'Errico 2020-1-28
Since the arivals will be a vector of potentially different lengths for every lambda, just save them in a cell array, each set in one cell.

请先登录,再进行评论。

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Error Detection and Correction 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by