It's noted that the number of folks in a stationary system will maintain an average equal to the rate of arrival multiplied by the mean of the service distribution.

The formula $L = \lambda w$ is valid for any queueing model in steady-state, where $L$ and $w$ are long-term steady-state average values respectively, and $\lambda$ denotes the arrival rate to the system.

We can add up the total service time in the system as follows:

$$\sum w_{j} = \sum i T_i$$

where $T_i$ represents the time units in which $i$ entities were in the system.

But things are always more interesting with sample sizing and simulating results, so say in $5000$ iterations we estimate a state of a system in 1-min intervals, between the first minute and the last arrival at a determined time.

So suppose we use a random interarrival rate of $\lambda = 2$ per minute and the service distribution is $N(8,1)$ minutes for the system.

How can I simulate this model in R, using rexp() and rnorm()? I also want to display in ts() to show in time plot.

  • 2
    $\begingroup$ While I think the topic does fit into our site, asking for a model in R will probably get more attention on SO (see stackoverflow.com/search?q=%5BR%5D+lambda for some questions there). But let's see what our other users will come up to answer your question. $\endgroup$
    – JakobS
    Jul 12, 2019 at 22:56
  • $\begingroup$ I came up with a portion, perhaps there can be come feedback see below $\endgroup$ Jul 13, 2019 at 0:17
  • $\begingroup$ Is your question about OR concepts in the queuing model, e.g., some aspect of Little’s Law, or is it about the implementation of some functions in R? If the latter, I agree with @JakobS that the question might get more response on SO. If the former, I think the question could use some clarification about what exactly you are asking. $\endgroup$ Jul 14, 2019 at 13:23
  • $\begingroup$ It was more oh do we set this up programmatically, but also a question could be why would a distribution of outcomes (ala simulation) provide and more rigor than any other methods $\endgroup$ Jul 14, 2019 at 20:33
  • $\begingroup$ I am happy to expand my answer if you can refine your question appropriately. We are all seeking to keep the answered rate high. $\endgroup$
    – CMichael
    Jul 15, 2019 at 22:05

2 Answers 2


Not directly answering your question of how to code it manually but for discrete simulation of queues in R I would strongly recommend the simmer package. The minimal code for your example would look like this (adapted from the tutorial).

lambda <- 2

queue <- trajectory() %>%
  seize("server", amount=1) %>%
  timeout(function() {rnorm(1, 8, 1)}) %>%
  release("server", amount=1)

env <- simmer() %>%
  add_resource("server", capacity=1) %>%
  add_generator("arrival", queue, function() {rnorm(1, lambda)}) %>%

resources <- get_mon_resources(env)
arrivals <- get_mon_arrivals(env)

plot(resources, metric = "usage", items = "server")
plot(arrivals, metric = "flow_time")

However, looking at the plot of flow times I am unsure if you conceived a stable queuing system - for mean = 1 in the service time it looks better:

Service Time Mean 8: mean = 8

Service Time Mean 1: mean = 1

  • $\begingroup$ Thanks for the suggestion! I’ll share what I came up with! $\endgroup$ Jul 14, 2019 at 21:14
  • 1
    $\begingroup$ Have fun - it is a powerful yet neat package. I updated my answer to be a bit more illustrative. $\endgroup$
    – CMichael
    Jul 14, 2019 at 21:24

So far this is what I have come up with

lambda <- 2
interarrivals <- rexp(5000,lambda)  ## (2 items per minute)

Provided the $\mu$ we expect that the interarrivals is about half a minute

mean(interarrivals) <- 0.516
service.times <- rnorm(5000, mean=8,sd=1)

where the service distribution is $N(8,1)$

arrival.times <- cumsum(interarrivals)
departure.times <- arrival.times+service.times

Where I am having issues is determining how many individuals remain in the system at time 1-min to the last.

I am having issues conceptualizing the model - this is what you are supposed to estimate (Little's law), which states that rate of arrival = 2 units per minute × mean of the service distribution = 8.

  • 1
    $\begingroup$ my intent is to convert the system state to a time.series ts() so I can plot. $\endgroup$ Jul 13, 2019 at 0:28

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.