Simular un valor Poisson

La distribución de Poisson es la otra cara de la moneda de la distribución Exponencial.  Mientras que la distribución exponencial nos sirve para modelizar el tiempo que transcurre entre dos sucesos independientes, la distribución de Poisson modeliza el número de sucesos que tienen lugar en una unidad de tiempo.

Por ejemplo:

• Si queremos saber cuántas llamadas entrarán en un Call Center cada hora, modelizaremos en base a Poisson.
• Si lo que queremos saber es cuánto tiempo transcurrirá entre dos llamadas a ese Call Center, modelizaremos con una distribución Exponencial.

Partiendo de esta idea, y disponiendo ya de la simulación de la distribución exponencial, la simulación de valores Poisson (por ejemplo, llamadas por minuto) puede ser como sigue:

• Se simula un primer valor E utilizando el generador de valores exponenciales. Esto nos indicará el tiempo que ha tardado en entrar la primera llamada.
• Si la llamada ha tardado más de un minuto en entrar, significa que durante el primer minuto han entrado cero llamadas. Por lo tanto, ese valor de Poisson vale cero.
• En caso contrario, ya ha entrado una llamada, y quizá entre una segunda. Simularemos valores exponenciales mientras las sumas de estos valores no superen la unidad de medida (en este caso, un minuto).
• Cada vez que simulamos un valor exponencial, sumamos 1 a P.

Suena un poco complicado, pero en el fondo lo puedes pensar así: arrancamos el cronómetro y lo paramos al minuto. Contamos el total de llamadas que han entrado en ese minuto por una línea. En el momento que entre una llamada después de ese minuto, ya está fuera del contador.

La distribución de Poisson, como la Exponencial, toma un solo parámetro λ, pero ojo porque el valor es el inverso. Es decir, que si en la Exponencial consideramos que cada llamada tarda 30 segundos en entrar, la Poisson indicará que hay 2 llamadas por minuto. En resumen: el parámetro tiene el mismo nombre pero toma el valor inverso. En la Exponencial sería 0,5 y en la Poisson sería 1/0,5 = 2. En el fondo es más intuitivo, aunque de buenas a primeras suena a chino, hay que reconocerlo.

[d3-link]

[/d3-link]

Pues bien, hecha la explicación, vamos a por la simulación de los datos con D3:

[d3-source canvas="chart"]
var
svg  = null,
bars = null,

generator = function() {
var 
lambda = 5,
simulations = 500,
barWidth    = 16.3,
barsScale = d3.scale.linear()
.domain([0, 3 * lambda])
.range([0,29]);
barXScale = d3.scale.linear()
.domain([0, 3 * lambda])
.range([0,270]);
var dataset = poissonGen(lambda, simulations, barsScale),

barYScale = d3.scale.linear()
.domain([0, d3.max(dataset, function(d) { return d;})])
.range([270, 0]),

yAxis = d3.svg.axis()
.scale(barYScale)
.orient("left")
.ticks(5),

xAxis = d3.svg.axis()
.scale(barXScale)
.orient("bottom")
.ticks(4);

if (svg !== null) {
bars = svg.selectAll("rect")
.data(dataset);

bars.enter();
} else {
svg = d3.select(".chart")
.append("svg")
.attr("width", 400)
.attr("height", 300);

bars = svg.selectAll("rect")
.data(dataset)
.enter()
.append("rect");
}

svg.selectAll("g.axis").remove();

svg.append("g")
.attr("class", "axis")
.attr("transform", "translate(30,0)")
.call(yAxis);

svg.append("g")
.attr("class", "axis")
.attr("transform", "translate(30,270)")
.call(xAxis);

bars.attr("x", function(d,i) {
  return 32 + i * (barWidth + 1);
 })
 .attr("y", function(d) {
  return 270-d;
 })
 .attr("width", barWidth)
 .attr("height", function(d) {
  return d;
 })
 .attr("fill", "blue");
},

poissonGen = function(lambda, simulations, barXScale) {
var datum = [];
  for(var row = 0; row < 30; row++) {
  datum.push(0);
  }

  for(var i = 0; i < simulations; i++) {
  var total = 0, poisson = -1;
  do {
  var u = Math.random();
  total += (-Math.log(u)/lambda);
  poisson++;
  } while(total < 1);

  datum[poisson]++;
  }
  return datum;
};
generator();
[/d3-source]

Se puede ver que en general su forma es acampanada. En realidad a partir de valores λ superiores a 30, las formas de la distribución de Poisson y de la distribución normal son casi indistinguibles (salvo, claro, porque la primera toma valores discretos y la segunda, continuos).