Bayesian Queueing Networks?

I am interested in using PyMC to estimate parameters in queueing network models. To date I have found the following tools outside PyMC:

  • The queue standard library has queue.Queue, queue.LifoQueue, and queue.PriorityQueue. Combined with a tool to handle a network of queues, like Networkx, this
  • "CiW is a discrete event simulation library for open queueing networks. Its core features include the capability to simulate networks of queues, multiple customer classes, and implementation of Type I blocking for restricted networks. A number of other features are also implemented, including priorities, baulking, reneging, schedules, batch arrivals, service disciplines, dynamic customer classes, and deadlock detection.
  • Queueing-tool is a ready-to-use queueing network library.
  • SimPy is a discrete event simulation tool with an emphasis on asynchronous even handling.

None of the above tools are Bayesian per se. My guess is that using queue.PriorityQueue might be the easiest path for integration with PyMC b/c it does not have any 3rd party dependencies. SimPy is likewise built on CPython with few dependencies. However, CiW and Queuing-tool are already ready-made implementations of queue networks.

Is there a way to implement PyMC queueing networks or patch existing tools to work with PyMC?

This blog post provides a minimal example of a Markov chain in PyMC3. However, modern queueing networks are often more detailed than just the Markovian state transitions between states. See queueing theory for an overview.

Objective: Implementing Markov Chains model
Author: Venkatesh Nagilla
Date: 2020-08-09

import pymc3 as pm
model = pm.Model()

import pymc3.distributions.continuous as pmc
import pymc3.distributions.discrete as pmd
import pymc3.math as pmm

with model:
    passenger_onboarding = pmc.Wald('Passenger Onboarding', mu=0.5, lam=0.2)
    refueling = pmc.Wald('Refueling', mu=0.25, lam=0.5)
    departure_traffic_delay = pmc.Wald('Departure Traffic Delay', mu=0.1, lam=0.2)

    departure_time = pm.Deterministic('Departure Time',
                                      12.0 + departure_traffic_delay +
                                      pmm.switch(passenger_onboarding >= refueling,

    rough_weather = pmd.Bernoulli('Rough Weather', p=0.35)

    flight_time = pmc.Exponential('Flight Time', lam=0.5 - (0.1 * rough_weather))
    arrival_traffic_delay = pmc.Wald('Arrival Traffic Delay', mu=0.1, lam=0.2)

    arrival_time = pm.Deterministic('Arrival time',
                                    departure_time +
                                    flight_time +

nb_samples = 500
with model:
    samples = pm.sample(draws=nb_samples, random_seed=1000)

Some related links:

PyMC will need the logp (and even better the gradient with respect to the parameters) of the observed events of your stochastic process.

To make it Bayesian you just add priors to the parameters of the process.

If you can write the logp using PyTensor operations you get the gradient for free (if you can write it using numpy you should be able to translate it to PyTensor pretty much 1-to-1).

Otherwise, if any library can calculate it for you, you can wrap arbitrary code in PyTensor like described in these examples:

It’s also easy to integrate with JAX so you can see if there’s a good JAX library that gives you what you need

1 Like

Quite more experimental you may also try approximate inference with a Simulator: Simulator — PyMC 5.10.3 documentation

1 Like