I am very new to Bayesian modelling and MCMC - I would like to know if the problem I describe below can be solved. It seems to be there is too much missing information but I wanted to get your thoughts. Consider the following:

I have a road intersection with one entrance and two exits, A and B. My goal is to estimate the number of cars that pass through this intersection in a given day, which equals to the number of cars that pass through the entrance. I post two people, one at exit A and the other at exit B, to count the number of cars coming out of their exits. They are both not very good so they only capture C percent of the cars that actually pass through their respective exits. I do not know what C is. To make matters worse, the person at exit B lost his record so I only have the numbers from A. I think this situation should be described by the following graph:

enter image description here

From historical data, I know on average p% of ppl go through A and (1-p)% people go through B, but on this given day, I have no information. In this example I only have two exits, but in general I may have more (e.g. 3-5).

Is is possible to estimate the distribution for "# cars thru entrance" with the data I have? If so, what distributions would you assign to each random variable?

  • $\begingroup$ By definition, $100*p$ percent of people must use A and $100*(1-p)$ percent of people must use B, unless there is a third exit, so that change in notation will simplify your graph. The same goes for number of cars thru A and B: they must add up to the total number. $\endgroup$ – Sycorax Jun 1 '14 at 15:02

I'll record my modelling thoughts as I create it:

import pymc as pm

How can I recreate this data? Well, N_T cars enter. Do I know N_T? No, so it's a random variable. For simplicity, I'll say its a discrete uniform with max 1000 (you can change this to, say, a Poisson)

N_T = pm.DiscreteUniform('N_T', 0, 1000)

Of these N_T cars that enter, I know a car will enter A with probability p. Thus the number of cars through A is a binomial rv. But do I know p for certain? Probably not, you mentioned it was a historical average, so we should model this stochastically (using the beta-binomial model):

p_A = pm.Beta('p_A', 50, 50) #centered at 0.5, but with some variance. This depends on your historical data.
N_A = pm.Binomial('N_A', N_T, p_A)

Next, only C percent of cars are recorded by A's observer. Let's just assume that C is completely accurate, i.e. there is no uncertainty about what C is:


Then the number of cars recorded is again a binomial, with N_A trials. Furthermore, suppose observer A record 20 cars.

obs_A = pm.Binomial('obs', N_A, C, observed=True, value = 20)  

Let's fire the mcmc canon:

mcmc = pm.MCMC( [obs_A, N_A, p_A, N_T] )
mcmc.sample(10000, 5000)

from pymc.Matplot import plot as mcplot

And I get this posterior for N_T:

enter image description here

  • $\begingroup$ Fun question, btw =) $\endgroup$ – Cam.Davidson.Pilon Jun 1 '14 at 16:34
  • $\begingroup$ Thanks for the reply and glad you found it interesting. All I know about this topic is from your EXCELLENT Bayesian Methods for Hackers book. I definitely need to re-read it again soon as I made a very quick pass (mainly because it was so engaging and readable) the first time. Any plans on completing chapters X1 and X2? $\endgroup$ – mchangun Jun 1 '14 at 16:44
  • $\begingroup$ Yes, but they will be in a different medium. Details soon ;) $\endgroup$ – Cam.Davidson.Pilon Jun 1 '14 at 16:47
  • $\begingroup$ What if I don't know the value of C and need to model that as a random variable as well? It's actually quite crucial to this problem that I don't assume I know it's value. Would I still be able to get reasonable answers? $\endgroup$ – mchangun Jun 1 '14 at 16:50
  • $\begingroup$ Of course: you can make C a pymc stochastic variable: C = pm.Beta('C', some_alpha, some_beta). The implication of doing this would be wider posteriors for N_T. $\endgroup$ – Cam.Davidson.Pilon Jun 1 '14 at 17:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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