binomiales y disculpas

Ayer el empate de la CUP nos dio a todos una buena oportunidad para refrescar nuestros instintos sobre la distribucion binomial. En cuanto salio el resultado de 1515 la opcion 1 vs 1515 la opcion 2, un monton de gente pensamos ¿cual es la probabilidad de empatar tirando 3030 monedas a cara o cruz?

Bien, obviamente numero de conjuntos posibles de 1515 cruces divido entre el total de conjuntos de cruces posibles, esto es (3030 1515) = 3030! / (1015! 1515!) dividido entre 2^3030. Un 1.45 % aproximadamente. O, usando la estadistica de la distribucion binomial con p=0.5 y 3030 extracciones:

(3030 1515) * (0.5)^1515 * (1-0.5)^1515

Bien. Y en estas empiezan a circular por la red dos discusiones. Una filosofica, por supuesto, sobre el concepto de probabilidad y si es aplicable en el caso de la toma de posiciones en un debate. Otra, mas tecnica, asumiendo que aqui lo que teniamos era una muestra de 3030 extracciones sobre alguna probabilidad subyacente, y cual era la matematica aplicable. Porque resulto que hubo muchos que opinaron que la probabilidad de empate era 1/3031, usease apenas un 0.03 %.

Las disculpas, por anticipado, son a todos estos que mantuvieron en la discusion que el resultado 1/3031 era válido, de alguna manera.  El unico argumento que se dio en twitter no iba precisamente muy bien para justificar su validez. Lo dio un profe de matematica aplicada de sevilla, que en 140 caracteres resumia:

Regla de Laplace para 3030 votos SÍ/NO de : 0/3030, …, 3030/0 dan 3031 casos posibles. Casos favorables: EMPATAR 1, NO EMPATAR 3030 (@mario_bilbao)

¿Regla de Laplace? ¿Casos favorables ente casos posibles? Pero eso es justo lo que hemos hecho en las tiradas de moneda, numero de conjuntos favorables entre numero de conjuntos posibles. No tiene sentido esa explicacion. Asi que cientos de personas nos aprestamos, yo entre los primeros, a replicar con explicaciones detalladas de como la binomial hacia que los 3031 casos posibles no fueran equiprobables. Por supuesto con p=0 tan solo seria posible el caso 3030 votos NO, y con p=1 solo seria posible 3030 votos SI, con una probabilidad cero de empate. Entre p=0 y el p=0.5 la probabilidad de empate se ira incrementando desde cero a 0.01449…   pero no se ve como el 0.03% sea especial. Parece claramente un error.

Por otro lado cuando alguna neurona me recordaba de mecanica estadistica la emergencia de casos donde conseguias cierta equiprobabilidad de configuraciones, y donde el rio suena agua lleva… puede que el señor Bilbao, que a fin de cuentas habrá pasado sobre estas ecuaciones muchas veces, recordara el resultado pero no la explicacion. Finalmente esta mañana @ruescasd, que habia dado independientemente el mismo 0.03% como estimacion conservadora, me lo ha explicado: ocurre si en vez de asumir p=0.5 o p en una banda 0.4 … 0.6 etc asumes que no sabes nada de la distribucion subyacente y escoges p al azar entre cero y uno.

Lo voy a poner en negrita, como me recomiendan en los comentarios:  en ausencia de información, el cálculo de la probabilidad de empate debe hacerse considerando todas los posibles valores de p, no p=0.5

¿Como se puede hacer esto? Es facil dicho para programadores: en la simulacion de la votacion en vez de usar

pvoto = 0.5
a =[ 1 if k > pvoto else 0 for k in npr.random(size=SIZE)]

usamos

pvoto = npr.random()
a =[ 1 if k > pvoto else 0 for k in npr.random(size=SIZE)]

Y cuando repetimos el experimento muchas veces entonces la probabilidad de empate convergera  a 1/3031, porque a fin de cuentas el numero “pvoto” no es mas que un numero random mas de los 3030 que hemos sacado, y cada uno de ellos tiene la misma probabilidad de estar en el medio. Es lo mismo que cuando sacamos tres hay un 33% de posibilidades de que el primero sea el que este enmedio, o cuando sacamos cinco hay un 20% de que el primero sea el que este enmedio (y un 20% de que sea el mayor, y un 20% de que sea el menor, etc).

Lo que estamos haciendo aqui, me explicaba @ruescasd, es usar una “beta binomial con prior uniforme“. Las beta binomial son una forma especialmente comoda, porque suele tener solucion analitica, de considerar a la vez la suma de muchas posibles distribuciones binomiales.

Por supuesto, una vez hemos abierto el coto de la simulacion con distintos escenarios para la distribucion subyacente podemos plantearnos distintos experimentos. Por ejemplo podriamos conjeturar que pvoto no es ni uniforme en [0,1] ni exactamente cara o cruz, sino que esta -uniformemente- en algun lado entre 0.4 y 0.6:

pvoto = 0.4 + 0.2*npr.random()
a =[ 1 if k > pvoto else 0 for k in npr.random(size=SIZE)]

En este caso la probabilidad de empate seria de tan solo el 0.16%… podriamos mantener una banda uniforme centrada en el “cara o cruz” e ir ensanchandola. Si la banda va de p=0.45 a p=0.55 la probabilidad de empate ya baja al 0.3%. Si va como p=0.4 a p=0.6 se reduce el 0.17%… asi hasta que la banda de indeterminacion cubre todo p=0 a p=1 y entonces tenemos la probabilidad de 0.03% dichosa.

Tabla: Caso de p alrededor de 0.5, variando la anchura

Para que veais que pronto se lia la cosa, vamos a interpolar entre las dos soluciones en disputa simplemente ampliando la incertidumbre acerca de la proximidad a 50:50 de la distribucion subyacente. Os doy la probabilidad de empate segun vamos ampliando

0.001:0.01447
0.002:0.01446 
0.01: 0.01373
0.02: 0.01203
0.03: 0.00991
0.04: 0.00802
0.05: 0.00653
0.06: 0.00550
0.07: 0.00472
0.08: 0.00412
0.09: 0.00366
0.1:  0.00329
0.2:  0.00164
0.3:  0.00109
0.4:  0.00082
0.5:  0.00066
0.6:  0.00054
0.7:  0.00047
0.8:  0.00041
0.9:  0.00036
1/3031=.0003299241...

Ya veis que funciona como hemos dicho mas arriba al comentar el algoritmo de simulacion: cuando la banda es solo de 0.499 a 0.501 la probabilidad es 0.01446, casi la misma que calculabamos con la binomial; si la incertidumbre va de 0.49 a 0.51 la probabilidad ya ha bajado a  0.01203; si va de 0.4 a 0.6 ya es de tan solo 0.00164. Por esto a no ser que se tenga certeza de que el resultado es estrecho -en este caso lo sabiamos por las votaciones anteriores- no es mala aproximacion suponer que la incertidumbre es total.

Pregunta: ¿Caso de p en banda [a,b]?

He intentado resolverlo analiticamente y me meto en un pantano de factoriales. En este caso estamos haciendo 3030 tiradas uniformes en (0,1), de las cuales daran SI las que salgan mayor que p, daran NO las que salgan  menor que p. Para que pueda haber empate es imprescindible que el numero n e tiradas cayendo en (0,a) y el numero m de tiradas cayendo en (b,1) puedan compensarse mediante un p en (a,b) que divida exactamente para contrapesar. Dentro del segmento (a,b) esto ocurre con probabilidad 1/(3031-m-n) dado que p en esa banda tiene la misma distribucion que el resto de los numeros que hayan caido ahi; es el mismo argumento de “equiprobabilidad” que con p uniforme en [0,1] en el caso de maxima ignorancia, y de hecho es la unica equiprobabilidad que hay en juego.

El problema esta en estimar la “trinomial” de caer en (0,a),(a,b) o (b,1) y luego sumarlas todas las probabilidades parciales con la condicion n < 1515, m<1515. Es un lio, pero es interesante pensarlo un ratito porque entiendes como al ir ensanchandose la banda va bajando la probabilidad de empate, y como tambien cambia si, sin ensancharse, se desliza hacia el p=0 o el p=1. Intuitivamente si la zona incertidumbre en la estimacion de la distribucion no pilla la zona del empate real, y dado que 3000 es una muestra ya bastante grande, sera dificil que se de el empate. Con p en el intervalo (0.4,0.5) estariamos en un 0.16%, con p en un intervalo mas grande, digamos (0.3,0.6) estariamos en 0.11%; y si lo suponemos en un intervalo (0.1,0.6) seria tan solo un  0.06%.

Obviamente se pueden dar probabilidades por debajo de la cota de “incertidumbre total”. Si estamos seguros de que p esta en el intervalo (0.1,0.48), por ejemplo, la probabilidad de empate es de solo 0.00001. Y con intervalos estrechos apoyados en los extremos la probabilidad de empate se ira acercando a cero, como corresponde al hecho de que ahi estarian todos en un solo bando.

¿Ignorancia total? ¿Principio de indiferencia?

Un argumento que se invoca para preferir la solucion prob(empate)=1/(n+1) es el “principio de indiferencia”, que ignoramos todo acerca de la distribucion subyacente y eso nos permite escogerla al azar.

Esto es, usamos la informacion de que hay una sola distribucion subyacente, con un unico parametro p que define su proporcion entre las dos opciones que se votan. Esto es intrigante, porque en el otro extremo si supusieramos que cada votante ha sido extraido de una distribucion diferente, ¡recuperariamos la binomial!. Porque para decidir su voto, el votante realizaria una tirara respecto de su parametro p, y sobre todos los votantes la probabilidad de que esas tiradas saquen la opcion A o la B es la misma. Podria arguirse pues que maximizar la indiferencia nos devuelve la binomial.

Esto se ve muy claro en el codigo. Para obtener 1/(n+1) usabamos en el bucle:

pvoto = npr.random()
a =[ 1 if k > pvoto else 0 for k in npr.random(size=SIZE)]

Si en vez de esto generamos un array con los p de cada distribucion de la que proviene cada votante

pvotante = npr.random(size=SIZE)
a =[ 1 if npr.random() > p else 0 for p in pvotante]

vemos que en realidad estamos generando de nuevo la binomial: la probabilidad de que  npr.random() sea mayor que otro npr.random() es de 1/2.

Un centenar de asambleas indiferentes

Por ultimo, sabemos que el caso real es intermedio: que los 3030 representantes lo eran de aproximadamente un centenar de asambleas.  ¡Podemos simular pues que no hay una unica p ni tres mil, sino un centenar! Para cada una de ellas desconocemos su proporcion p de votos entre las dos opciones, pero podemos aplicar la regla esa de indiferencia para simlarla. Y como podriamos esperar, saldra una probabilidad de empate intermedia:

  • Asumiendo ignorancia de las distribuciones, de 101 asambleas enviando cada una 30 miembros, la probabilidad de empate seria el 0.45%
  • Asumiendo ignorancia de las distribuciones, de 101 asambleas enviando cada una un numero aleatorio de miembros hasta sumar 3030,  la probabilidad de empate seria el 0.31%

El codigo para simular esto, basado en el gist de Galli, seria:

for i in xrange(1, LOOPS):
  l=np.empty([SIZE])
  for j in range(0,101):
    p=npr.random()
    l[j*30:j*30+30].fill(p)
  a =[ 1 if npr.random()> p else 0 for p in l ]
  if sum(a) == TIE:
    ties += 1
  if i % 1000 == 0:
  print("Total: {} ties: {} prob: {}".format(i, ties, ties/float(i)))

Y si queremos variar el numero de asistentes por cada asamblea

for i in xrange(1, LOOPS):
  l=np.empty([SIZE])
  div=npr.randint(SIZE,size=ASAMBLEAS+1)
  div.sort()
  div[0]=0
  div[ASAMBLEAS]=SIZE-1
  for j in range(0,ASAMBLEAS):
    p=npr.random()
    l[div[j]:div[j+1]].fill(p)
  a =[ 1 if npr.random()> p else 0 for p in l ]

Como hemos dicho, a medida que el numero de asambleas crece, el resultado vuelve a estar cerca del binomial. Para 606 asambleas enviando cada una 5 representantes, la probabilidad de empate es del 0.95%.

Hay aqui un poco de paradoja, que seguramente es la causa del debate. Revisemosla: si no sabemos nada de un asistente concreto, podemos pensar que la asamblea que proviene tiene una distribucion de probabilidad con un valor p que desconocemos aleatorio entre 0 y 1, y que el propio asistente puede ser de uno de los bandos con probabilidad p, del otro con probabilidad (1-p), lo que tenemos que asignar tirando un segundo aleatorio r, entre 0 y 1, y comprobando si es mayor o menor que p.  Y obviamente la probabilidad de que r>p, siendo los dos numeros aleatorios, es la misma que intercambiandolos, luego es 1/2.  Este es el mejor estimado que tenemos para el caso de considerar lo que va a votar un solo asistente proviniente de una asamblea dada cuya distribucion de preferencias desconocemos.  Pero cuando consideramos varios asistentes provinientes todos de esa misma asamblea la tirada “r” de cada uno cambia y todos tienen en comun la misma “p”, de forma que a mas asistentes, mas se acerca el resultado a la original de valor “p”, aunque sigamos desconociendola, y hay que incorporar ese hecho en las hipotesis, con la consecuencia de que la “binomialidad” se diluye en la aleatoriedad de p. Afortunadamente, como han discutido otros posts, suele haber mas datos para no tener que escoger entre estos dos casos extremos.

Una forma rápida de ver esta interpolación entre los dos casos es fijarse -como me han recordado Ruescas y Vilches- en qué variables aleatorias estan operando. Partimos del caso donde la ignorancia sobre la distribucion la representamos con una sola variable aleatoria que toma valores uniformemente entre 0 y 1. Al ir añadiendo asambleas, tenemos multiples variables aleatorias tomando valor uniforme pero su suma, pesada por la proporcion de asistentes, no es ya una variable aleatoria distribuida uniformente sino que se ira centrando -promediando- alrededor de 1/2.

Otras referencias

Parece que media internet ha entrado al trapo del debate, basta con buscar 3030 empate en google y salen decenas de blogs 0_0, y lo mismo en twitter.  Aparte de las aportaciones ya citadas de Galli y Ruescas en blogs,y gists, y de los comentarios de twitter, puede interesarte ver el articulo de prensa de Llaneras et al y comparar con alguno de los primeros, como el de Baraza.  Pueden valer para recordar que en todo caso la regla de “casos favorables dividido entre casos posibles” es una justificacion erronea de la estimacion 1/3031, y que la justificacion valida es mucho mas complicada y implica usar un concepto de ignorancia muy razonado, y su matematica asociada.

Mas jocosamente, pero con su toque divulgativo, ha tratado este asunto Federico Ardila (EN, CAT), usando el caso p=1/2 para terminar con un par de calculos probabilisticos del proceso de recuento que producen los numeros de Catalan.

7 thoughts on “binomiales y disculpas

  1. Lo añado y en negrita.

    En cuanto a errores, creo que es especialmente importante la confusion con la palabra “caso”. No se quien se confundio aqui; yo soy de la creencia que si estamos hablando de depositar 3030 votos a dos opciones A y B, el numero de casos posibles es de 2^3030. Lo que es igual a 3031 es el numero de particiones. Ignoro si es un error de vocabulario o conceptual, pero es lo que realmente desencadeno la tormenta, y si hay que hacer una labor de divulgacion deberia estar basada en desenredar estos dos conceptos.

  2. Bueno, Alejandro, entro a leer tu post para descubrir que sigues sin reconocer, simplemente, que el único verdadero error, el único error conceptual grave que ha habido en todo este asunto ha sido la cantinela con la que respondisteis al cálculo 1/3031 (y con especial virulencia al profesor Bilbao), a saber: “en ausencia de información, debemos considerar que el voto se reparte como una binomial con p=0.5”.

    Este tremendo patinazo, que ha recibido gran eco en blogs y columnas periodísticas, es lo que es urgente rectificar, si es que realmente estamos implicados con la cultura matemática y científica de la audiencia.

    En primer lugar deberías ser honesto y reconocer cual era la pregunta originalmente formulada, es decir: ¿cual era la probabilidad de que 3030 personas empataran en una votación entre dos opciones? Esta, y ninguna otra, era la pregunta formulada. No era ¿cual es la probabilidad de que 3030 empaten si ya ha habido dos rondas y en la primera quedaron 1510-1520 y en la segunda 1515-1515? no, esa no era la pregunta, y lo sabes, como demuestran los muchos tuits, posts, y comentarios en los posts, vertidos en los días iniciales, cuando ninguna mención a la “info previa” se hizo, y aun hoy cuando este “mantra” (sin información p=0.5) sigue manteniendo en el debate.

    Lo primero que tenemos es que ser honestos con esto, esa era la pregunta, y eso es lo que resultaba chocante. No choca que una elección donde ya han casi empatado, empaten en segundas, como no choca que si empatan en primera empaten en segunda. Esa era la pregunta y la respuesta correcta es 1/3031. Es una cuestión previa de honestidad que esto quede absolutamente claro, y que vuestros lectores oigan de vuestra boca que “en ausencia de información, el cálculo de la probabilidad de empate debe hacerse considerando todas los posibles valores de p, no p=0.5”. Dirigir el debate ahora hacia otra pregunta sin este reconocimiento no es más que una cortina de humo y abundar en la confusión. Ese fue el único error de bulto, y debe decirse sin paliativos ni monsergas, ni ocultarlo con interminables cuentas de escenarios diversos. Es una cuestión de honestidad intelectual.

    Por otra parte, en una votación, salvo los pocos casos de voto fluctuante, lo normal es que los votos estén predeterminados, o queden predeterminados en algún momento, lo “aleatorio” no es el proceso de votar, ¡es el proceso de elegir las personas que votarán! es decir, el proceso de extraerlas de una población de elementos cuya intención de voto (en ausencia de información) seguirá una distribución no necesariamente aleatoria que podremos aproximar por una binomial ¡pero cuyo valor p no es es desconocido! Una vez extraídos estos elementos electores, el voto está decidido (o casi decidido, al menos, pues podemos aceptar que haya cambios en algún elemento). Así que lo chocante no es que empaten después de casi empate, eso no extraña a nadie, lo que choca es que de aquella distribución de p desconocida se hayan extraído 1515 electores de cada opción. ¿como era esto de probable? ¡esa era la pregunta! y la respuesta 1/3031.

    Y la respuesta es esa por el argumento del profesor Bilbao. Cómo este argumento enlaza con el “cálculo correcto” es complejo, matemáticamente, pero simple desde la lógica y si aun tienes dudas sobre ello no es por tu falta de conocimientos estadísticos, sino por que no quieres reconocer que el profesor Bilbao respondió correctamente: Desconocido p, y en ausencia de información, todo p debe considerarse igualmente probable, pero como dado p la proporción media entre SI/NO queda fijada y a cada p equivale a una única proporción, y todos los p tienen la misma probabilidad, todas las proporciones (resultados, no combinaciones de votos particulares) tienen la misma probabilidad. Espero que esto te lo deje ya claro y elimines del post, definitivamente, todo comentario ofensivo hacia el profesor Bilbao (al que te recuerdo que solo me une el amor a las matemáticas).
    Por otra parte no podemos utilizar el resultado para estimar su verosimilitud ¡pues es ese resultado de cuya aleatoriedad dudamos!. Si dudamos de que una observación haya sido un resultado fortuito, debemos calcular su probabilidad de ocurrencia a priori, y no podemos incluir en el cálculo la información derivada de que el evento ha ocurrido, pues esto sesgará nuestros resultados aumentando artificialmente la probabilidad de que ocurra. Imagina que el resultado hubiera sido 3030 “si”. Si consideras esa información deducirás que lo más probable es que p=1, o casi, y entonces P(todos_si) será muy alta. Estoy seguro de que no la harías.

    Creo que debéis esforzaros en borrar todo rastro del error tan grave que habéis propagado entre vuestros lectores y en las redes. Tú has dado un primer paso con este post, pero queda camino, a mi entender, si realmente vuestra motivación es divulgar la ciencia y las matemáticas. Pero además por una cuestión de justicia, pues en este debate algunos de nosotros hemos recibido a causa de este error un significativo maltrato verbal. Y te lo pido, sinceramente, porque de ti lo espero. A otros sé que pedirselo será tiempo perdido.
    Espero que estos comentarios no te molesten, si es así no dudes en decirlo y me disculparé por lo que consideres necesario. Buenas noches y Feliz Año Nuevo.

    PD: sobre el modelo de asambleítas,

    como te comenté, creo que he encontrado la falacia. Sortear p aleatoriamente es suponer que las asambleas no guardan ninguna relación entre ellas, lo que significa ignorar la existencia de una distribución no aleatoria en la población, es decir, la existencia de una preferencia poioblaciona en uno u otro sentido. Yo acepto que podemos y debemos considerar esa variación entre asambleitas, pero esta variación debe ser razonable y mantener el hecho de que esa población tiene un valor poblacional de p, un sesgo hacia si o no. Ignorar esto haría innecesario votar, pues sería asumir que siempre tendremos empate y que los no empates serán simplemente fruto de la aleatoriedad del muestreo. En este caso, el sorte de p para cada asambleita debe hacerse en torno a un valor medio de p prefijado, y aceptando, como mucho, la variabilidad correspondiente a una submuestra de tamaño n (elementos de la asambleita) es decir 1/raiz(n). Si aceptamos todo el intervalo para la población, lo más extremo sería aceptar 0.1 para la variación de p entre asambleitas, una vez fijado el p de la población total. Así, deberías simular siguiendo el siguiente algoritmo

    1.- sortear p de la población
    2.- sortear elementos de cada asamblea (condicionando suma = 3030) n
    3.- sortear p de asamblea con media p e intervalo -0.5/raiz(n), 0.5/raiz(n)
    4.- sortear si/no en cada asamblea, y total

  3. La política en España desde hace ya tiempo (y lo queda) es como el juego de las sillas musicales, pero al revés. El primero que se sienta, pierde. Por eso todo son empates, votaciones sucesivas, escenificar desencuentros a la hora de pactar, andar por ahí de campaña en vez de resolviendo problemas reales, y, en general, hacer todo lo que sea menos trabajar. Eso sí, cobrando. Es la mentalidad del funcionario lo que hay que erradicar para poder levantar cabeza.

  4. El tema es diferente, sí, era algo de “ingenieria inversa”: ¿como diablos pude alguien deducir que la prob de empate es 1/(n+1), esto es, que todas las particiones son equiprobables, si obviamente no lo son para ningun valor de “p” que tomemos? ¿que campanas han oido esta gente, aunque no recuerden donde suenen? Sospechaba que habia una porque cuando das termoestadistica hay un caso similar en el que de repente por arte de magia todo es equiprobable y te vuelves loco.

  5. “con este modelo te estás planteando un tema diferente:

    “Cuál es la probabilidad de tener una votación con empate exacto si cada una de ellas tiene un sesgo a favor de una opción””

    De nuevo, mal. La beta binomial responde a

    “Cuál es la probabilidad de tener una votación con empate exacto _si no sabes absolutamente nada sobre las preferencias_.”

    El unico que hace suposiciones eres tu cuando dices que dicha preferencia es exactamente 0.5. Ahora bien, otra cosa es que añadas informacion al problema:

    “Cuál es la probabilidad de tener una votación con empate exacto dado que las preferencias seran similares a las mostradas en esta otra votacion”

    Y de nuevo, eso se puede calcular de manera facil con el beta binomial y el teorema de bayes. Al sera la beta el conjugate del binomial puedes incorporar evidencias pasadas al calculo, entonces si te va a dar mas probabilidad para resultados proximo al empate, lo cual se puede ver aqui, donde hemos hecho inferencia con datos de una votacion anterior

    https://www.wolframalpha.com/input/?i=plot+beta+binomial+distribution+alpha+%3D+1482%2C+beta+%3D+1512%2C+n+%3D+3030

  6. Ya te contesté por Twitter, con este modelo te estás planteando un tema diferente:

    “Cuál es la probabilidad de tener una votación con empate exacto si cada una de ellas tiene un sesgo a favor de una opción”.

    Es decir, estás asumiendo que hay preferencias a priori, cuando en realidad no lo sabes y es más: los resultados de las votaciones anteriores muestran que ese sesgo era prácticamente inexistente. En todo caso deberías usar una banda muy pequeña y los valores tenderían a los que dijimos antes.

    Por otro lado, si se asume que las preferencias de los individuos son diferentes y no equiprobables se cumple lo que decía antes del Teorema Central del Límite.

    Aquí tienes un ejemplo de un sesgo fijo para cada persona para todas las votaciones: https://gist.github.com/gallir/ca8d6cb958bd7767938f

    Aquí este otro donde el sesgo individual varía en cada elección: https://gist.github.com/gallir/ffad885419ccccb579c7

    En ambos casos la tendencia es hacia lo mismo, el valor calculado antes.

Leave a Reply / Añade un comentario: