¿Por qué esta distribución es uniforme?

12

Estamos investigando las pruebas estadísticas bayesianas, y nos encontramos con un fenómeno extraño (al menos para mí).

Considere el siguiente caso: estamos interesados ​​en medir qué población, A o B, tiene una tasa de conversión más alta. Para una comprobación de cordura, establecemos , es decir, la probabilidad de conversión es igual en ambos grupos. datos artificiales utilizando un modelo binomial, por ejemplo,n ABinomial ( N , p A )pA=pB

nABinomial(N,pA)

Luego tratamos de estimar el pA,pB usando un modelo beta-binomial bayesiano para obtener posteriores para cada tasa de conversión, por ejemplo,

PABeta(1+nA,NnA+1)

Nuestra estadística de prueba se calcula calculando S=P(PA>PB|N,nA,nB) través de Monte Carlo.

Lo que me sorprendió fue que si pA=pB , entonces SUniform(0,1) . Mis pensamientos eran que sería centra alrededor de 0,5, e incluso convergen a 0,5 como el tamaño de la muestra, N , crece.

Mi pregunta es, ¿por qué SUniform(0,1) cuando pA=pB ?


Aquí hay un código de Python para demostrar:

%pylab
from scipy.stats import beta
import numpy as np
import pylab as P

a = b = 0.5
N = 10000
samples = [] #collects the values of S
for i in range(5000):
    assert a==b
    A = np.random.binomial(N, a); B = np.random.binomial(N, b)
    S = (beta.rvs(A+1, N-A+1, size=15000) > beta.rvs(B+1, N-B+1, size=15000)).mean() 
    samples.append(S)

P.hist(samples)
P.show()
Cam.Davidson.Pilon
fuente
Tenga en cuenta que no puede ser exactamente uniforme, porque es una variable discreta. Por lo tanto, está preguntando sobre el comportamiento asintótico. Además, para pequeño (menos de , aproximadamente, con ) la distribución ni siquiera es remotamente uniforme. N 100 / min ( p , 1 - p )SN100/min(p,1p)p=pA=pB
whuber
@whuber S no es discreto, es una probabilidad que puede caer entre 0 y 1. Además, incluso para N bajo, estoy observando un comportamiento uniforme.
Cam.Davidson.Pilon
2
Debo estar malentendiendo su configuración, entonces. Por lo que puedo decir, para cualquier valor dado de el valor de es un número. Por lo tanto, aceptando que y están fijos por el momento (como están en su código), es una función de . Pero este último, al ser la realización de dos distribuciones binomiales, solo puede alcanzar un conjunto discreto de valores. Cuando reproduzco en su código , consigo histogramas decididamente no uniforme por pequeñas . S S ( n A , n B ) NN,nA,nB,Sp BN,pA,pBS(nA,nB)RN
whuber
1
Aunque de hecho su tiene valores entre y , no confunda eso con no discreto: puede tener como máximo valores distintos (y en realidad tiene menos que eso). Es posible que esto no sea perfectamente claro para usted porque su simulación genera estimaciones de lugar de sus valores correctos y las estimaciones esencialmente tienen una distribución continua. 0 1 N 2 SS01N2S
whuber
1
@whuber sí, tienes razón, excelente observación. Todavía estoy atrapado en por qué se ve uniforme entonces.
Cam.Davidson.Pilon

Respuestas:

11

TL; DR: Las mezclas de distribuciones normales pueden verse uniformes cuando los tamaños de los contenedores son grandes.

Esta respuesta toma prestada del código de muestra de @ whuber (que primero pensé que era un error, pero en retrospectiva probablemente fue una pista).

Las proporciones subyacentes en la población son iguales: a = b = 0.5.
Cada grupo, A y B tiene 10000 miembros: N = 10000.
Vamos a realizar 5000 repeticiones de una simulación: for i in range(5000):.

En realidad, lo que estamos haciendo es una de una . En cada una de las 5000 iteraciones haremos . s i m u l a t i o n u n d e r l y i n g s i m u l a t i o n p r i m esimulationprimesimulationunderlyingsimulationprimesimulationunderlying

En cada iteración del vamos a simular un número aleatorio de A y B que son 'éxitos' (también conocido como convertida) dadas las proporciones subyacentes iguales definidos anteriormente: . Nominalmente, esto producirá A = 5000 y B = 5000, pero A y B varían de una ejecución sim a otra y se distribuyen entre las 5000 ejecuciones de simulación de forma independiente y (aproximadamente) normalmente (volveremos a eso).simulationprimeA = np.random.binomial(N, a); B = np.random.binomial(N, b)

Pasemos ahora a para una única iteración de en la que A y B han tenido el mismo número de éxitos (como será el promedio). En cada iteración de , dados A y B, crearemos variantes aleatorias de la distribución beta para cada grupo. Luego los compararemos y descubriremos si , produciendo un VERDADERO o FALSO (1 o 0). Al final de una ejecución de , hemos completado 15000 iteraciones y tenemos 15000 valores VERDADERO / FALSO. El promedio de estos arrojará un solo valor de la distribución de muestreo (aproximadamente normal) de la proporción de s i m B e tsimulationunderlying s i m u l a t i o n u n d e r l y i n g B e t un a >simulationprimesimulationunderlying simulatio n u n d e r l y i n g B e t a A > B e t a BBetaA>BetaBsimulationunderlyingBetaA>BetaB .

Excepto que ahora seleccionará los valores 5000 A y B. A y B rara vez serán exactamente iguales, pero las diferencias típicas en el número de éxitos de A y B se ven reducidas por el tamaño total de la muestra de A y B. Los As y Bs típicos producirán más tirones de su distribución de muestreo de proporciones de , pero los que están en los bordes de la distribución A / B también serán eliminados.B e t a A > B e t a BsimulationprimeBetaA>BetaB

Entonces, lo que en esencia detenemos muchas ejecuciones de simulación es una combinación de distribuciones de muestreo de para combinaciones de A y B (con más tirones de las distribuciones de muestreo hechas de los valores comunes de A y B que los valores poco comunes de A y B). Esto da como resultado mezclas de distribuciones normales-ish. Cuando los combina en un tamaño de contenedor pequeño (como es el valor predeterminado para la función de histograma que utilizó y se especificó directamente en su código original), termina con algo que parece una distribución uniforme.BetaA>BetaB

Considerar:

a = b = 0.5
N = 10
samples = [] #collects the values of S
for i in range(5000):
    assert a==b
    A = np.random.binomial(N, a); B = np.random.binomial(N, b)
    S = (beta.rvs(A+1, N-A+1, size=15000) > beta.rvs(B+1, N-B+1, size=15000)).mean() 
    samples.append(S)

P.hist(samples,1000)
P.show()
russellpierce
fuente
1
Entonces, hay una diferencia entre el mío y su código. Muestro A y B en cada ciclo, lo muestras una vez y calculas S 5000 veces.
Cam.Davidson.Pilon
1
La discrepancia radica en sus llamadas a rbinom, que devuelve un vector. La llamada posterior al rbetainterior replicatese vectoriza, por lo que el bucle interno (interno) está utilizando un y diferente para cada una de las 15000 variables aleatorias generadas (envolviendo las 5000 finales desde su ). Mira para más. Esto difiere del código de @ Cam con un solo y fijo utilizado en todas las 15000 llamadas de variable aleatoria para cada uno de los 5000 bucles de muestreo ( ). B A BABNSIM = 10000?rbetaABreplicate
cardenal
1
Aquí está la salida para aquellos curiosos: imgur.com/ryvWbJO
Cam.Davidson.Pilon
1
Lo único que sé que son potencialmente pertinentes a nivel conceptual es que a) la distribución esperada de resultados es simétrica, b) un tamaño de bin de 1 siempre es uniforme, c) un tamaño de bin de 2 para una distribución simétrica también siempre parecerá uniforme, d) el número de posibles distribuciones de muestreo que pueden extraerse de los aumentos con N, e) los valores de S no pueden acumularse en 0 o 1 solo porque beta no está definido cuando hay 0 éxitos en cualquiera de los grupos , y f) las muestras están restringidas entre 0 y 1.
russellpierce
1
Solo como observación podemos ver que las distancias entre los centroides de las distribuciones de muestreo se reducen a medida que los centroides de las distribuciones de muestreo se alejan de .5 (probablemente relacionado con el punto f anterior). Este efecto tiende a contrarrestar la tendencia de las altas frecuencias de observaciones para los éxitos más comunes casi iguales en el caso del grupo A y del grupo B. Sin embargo, dar una solución matemática de por qué es eso o por qué debería producir distribuciones normales para ciertos tamaños de contenedores no está cerca de mi territorio.
russellpierce
16

Para tener una idea de lo que está sucediendo, siéntase libre de hacer que muy grande y, al hacerlo, ignorar el comportamiento de y explotar los teoremas asintóticos que establecen que las distribuciones Beta y Binomial se vuelven aproximadamente normales. (Con algunos problemas, todo esto puede hacerse riguroso). Cuando hacemos esto, el resultado surge de una relación específica entre los diversos parámetros.O ( 1 / N )NO(1/N)


Debido a que planeamos usar aproximaciones normales, prestaremos atención a las expectativas y variaciones de las variables:

  • Como Binomial variables aleatorias, y tienen expectativas de y varianzas de . En consecuencia y tienen expectativas de y la varianza .n A n B p N p ( 1 - p ) N α = n A / N β = n B / N p p ( 1 - p ) / N(N,p)nAnBpNp(1p)Nα=nA/Nβ=nB/Npp(1p)/N

  • Como una Beta , tiene una expectativa de y una varianza de . Aproximadamente, encontramos que tiene una expectativa de(nA+1,N+1nA)PA(nA+1)/(N+2)(nA+1)(N+1nA)/[(N+2)2(N+3)]PA

    E(PA)=α+O(1/N)

    y una varianza de

    Var(PA)=α(1α)/N+O(1/N2),

    con resultados similares para .PB

Por lo tanto, aproximaremos las distribuciones de y con Normal y Normal (donde el segundo parámetro designa la varianza ). La distribución de consecuencia es aproximadamente Normal; esto es,PAPB(α,α(1α)/N)(β,β(1β)/N)PAPB

PAPBNormal(αβ,α(1α)+β(1β)N).

Para muy grande , la expresión no variará apreciablemente de excepto con muy baja probabilidad (otro término descuidado ). En consecuencia, dejando que sea ​​el CDF normal estándar,Nα(1α)+β(1β)p(1p)+p(1p)=2p(1p)O(1/N)Φ

Pr(PA>PB)=Pr(PAPB>0)Φ(αβ2p(1p)/N).

Pero dado que tiene media cero y varianza es una Normal estándar variante (al menos aproximadamente). es su probabilidad de transformación integral ; es uniforme .2 p ( 1 - p ) / N , Z = α - βαβ2p(1p)/N, ΦΦ(Z)Z=αβ2p(1p)/NΦΦ(Z)

whuber
fuente
1
No estoy contigo hasta ... luego te vas a otra dirección que no seguí del todo. ¿Se define dos veces, una como el CDF normal estándar y luego como la transformación integral de probabilidad? Espero que pueda expandir su descripción en torno a estos pasos y relacionarlos con el código / problema inicial. Quizás retroceda y repita qué parámetros específicos producen el resultado uniforme. ΦPAPBNormalΦ
russellpierce
1
@rpierce (1) La diferencia es aproximadamente normal porque y son independientes y cada una es aproximadamente normal. La media es la diferencia de las medias y la varianza es la suma de las varianzas. (2) La transformación integral de probabilidad es el CDF: es el caso de cualquier variable aleatoria con distribución continua , que es uniforme. P A P B X F F ( X )PAPBPAPBXFF(X)
whuber
1
Oh, obtuve 1, fueron las cosas después de donde me perdí. Esto será increíblemente tonto, pero ¿por qué lo mismo que el CDF? Pr(PA>PB)
russellpierce
1
@rpierce Eso se deduce directamente de la definición, pero hay un ligero giro en el que se invoca la simetría de la distribución Normal. Estamos hablando de una variable aleatoria normal supone que tienen una expectativa de y la varianza . Estandarizando , es natural reescribir la probabilidad comoX=PAPBμ=αβσ2=2p(1p)/NX
Pr(X>0)=Pr((Xμ)/σ>(0μ)/σ)=1Φ(μ/σ)=Φ(μ/σ).
whuber
3
@whuber esto es bastante sorprendente. Eres un maestro maravilloso Le agradezco tanto la respuesta suya como la de rpierce, todavía le daré crédito porque resolvió nuestro problema, y ​​usted ha demostrado por qué ocurre el comportamiento. Ty!
Cam.Davidson.Pilon