Bayesian Inference: Modeling checkout times at a store

Question

I am currently learning how to use Bayesian inference. I have been making up problems (by defining some population parameters) and then trying to infer those values from samples. I recently made up a problem that I cannot figure out how to solve.

I would like to model the time that it takes to checkout at a grocery store. In this model, there are two main parts. The first is, the time it takes the cashier to ring up and bag the items. The second is, the time it takes to pay the cashier.

I have chosen to model the time it takes the cashier to scan and bag $k$ items (assume for simplicity that the same number of items are always purchased) as a gamma distribution:

$S \sim \textrm{Gamma}(k, \theta_1)$

The $\theta_1$ parameter is the average time it takes to scan and bag one item. So this is modeling the sum of $k$ exponential random variables (with rate parameter $1/\theta_1$).

I have chosen to model the time it takes to pay the cashier as an exponential variable:

$P \sim \textrm{Exponential}(\theta_2)$

The total time it takes to checkout is the sum of the time it takes to scan and bag plus the time it takes to pay:

$W_i \sim S + P$

If I directly observe the random variables $s_i$ and $p_i$ (and I know the number of items k), then the inference problem for $\theta_1$ and $\theta_2$ is easy. The likelihood function for $s_i$ is the gamma distribution and the likelihood function for $p_i$ is the exponential distribution. Standard Monte Carlo sampling techniques can be used to obtain their posterior distributions.

However, if all that I know are the total waiting times $w_i$ (and $k$), then I do not know what the likelihood function should be. From the research I have done, it appears that in general sums of random variables do not have closed form PDFs.

How can this kind of problem be solved? Are there certain approximations that could be made to obtain a closed form likelihood?

I am currently using the pymc3 library to do Bayesian inference. I realize that it is possible to perform a numerical integration to evaluate the convolution of two PDFs, but this is computationally expensive and I'm not sure it would be possible to obtain the gradient of such a function, which would rule out the use of Hamiltonian Monte Carlo sampling (e.g. NUTS).

Answer 1

Hello, qsc. I am very sympathetic to your quest of learning Bayesian statistics by posing problems to yourself, and then trying to solve them. However, this method also meets many roadblocks at times, such as here.

You are not trying to infer $\theta_1$ and $\theta_2$ individually (which is easy), so I assume you are attempting to infer them from their convolution, and working backward, given your assumptions on the distributions taken by $S$ and $P$.

So the first step is solving this convolution problem. Based on wikipedia's convolution page, the convolution of Exponential's is a Gamma, and the convolution of Gamma's is yet another Gamma, so based on some intuition, I imagine the conv. you are looking for would also produce a Gamma. Here, you are in luck, because as you found, most convolutions do not have closed-form answers as easy as this.

However, where you are out of luck, is that there does not seem to be a known, documented way of performing convolution via Bayesian methods, and then using that to go backwards somehow. This is a simple roadblock that I hope you figure out how to overcome with enough research.