Muestreo de yacimientos

Algoritmo aleatorio

El muestreo de reservorio es una familia de algoritmos aleatorios para elegir una muestra aleatoria simple , sin reemplazo, de k elementos de una población de tamaño desconocido n en una sola pasada sobre los elementos. El tamaño de la población n no es conocido por el algoritmo y normalmente es demasiado grande para que todos los n elementos quepan en la memoria principal . La población se revela al algoritmo con el tiempo, y el algoritmo no puede mirar hacia atrás a los elementos anteriores. En cualquier punto, el estado actual del algoritmo debe permitir la extracción de una muestra aleatoria simple sin reemplazo de tamaño k sobre la parte de la población vista hasta el momento.

Motivación

Supongamos que vemos una secuencia de elementos, uno a la vez. Queremos conservar 10 elementos en la memoria y queremos que se seleccionen al azar de la secuencia. Si conocemos el número total de elementos n y podemos acceder a ellos de forma arbitraria, la solución es sencilla: seleccionar 10 índices i distintos entre 1 y n con la misma probabilidad y conservar los elementos i -ésimos. El problema es que no siempre conocemos el número exacto de n de antemano.

Sencillo: Algoritmo R

Jeffrey Vitter creó un algoritmo simple y popular pero lento, el algoritmo R. ^[1]

Inicializa una matriz indexada de a , que contiene los primeros k elementos de la entrada . Este es el depósito . ${\displaystyle R}$ ${\displaystyle 1}$ ${\displaystyle k}$ ${\displaystyle x_{1},...,x_{k}}$

Para cada nueva entrada , genere un número aleatorio j de manera uniforme en . Si , entonces establezca . De lo contrario, descarte . ${\displaystyle x_{i}}$ ${\displaystyle \{1,...,i\}}$ ${\displaystyle j\in \{1,...,k\}}$ ${\displaystyle R[j]:=x_{i}.}$ ${\displaystyle x_{i}}$

Regresar después de que se hayan procesado todas las entradas. ${\displaystyle R}$

Este algoritmo funciona por inducción en . ${\displaystyle i\geq k}$

Prueba

Cuando , el algoritmo R devuelve todas las entradas, proporcionando así la base para una prueba por inducción matemática. ${\displaystyle i=k}$

Aquí, la hipótesis de inducción es que la probabilidad de que una entrada particular se incluya en el depósito justo antes de que se procese la -ésima entrada es y debemos demostrar que la probabilidad de que una entrada particular se incluya en el depósito es justo después de que se procese la -ésima entrada. ${\displaystyle (i+1)}$ ${\displaystyle {\frac {k}{i}}}$ ${\displaystyle {\frac {k}{i+1}}}$ ${\displaystyle (i+1)}$

Aplique el algoritmo R a la entrada -ésima. La entrada se incluye con probabilidad por definición del algoritmo. Para cualquier otra entrada , por la hipótesis de inducción, la probabilidad de que se incluya en el depósito justo antes de que se procese la entrada -ésima es . La probabilidad de que todavía esté incluida en el depósito después de que se procese (en otras palabras, que no se reemplace por ) es . Esto último se desprende del supuesto de que el entero se genera de manera uniforme al azar; una vez que queda claro que de hecho ocurrirá un reemplazo, la probabilidad de que en particular se reemplace por es . ${\displaystyle (i+1)}$ ${\displaystyle x_{i+1}}$ ${\displaystyle {\frac {k}{i+1}}}$ ${\displaystyle x_{r}\in \{x_{1},...,x_{i}\}}$ ${\displaystyle (i+1)}$ ${\displaystyle {\frac {k}{i}}}$ ${\displaystyle x_{i+1}}$ ${\displaystyle x_{r}}$ ${\displaystyle x_{i+1}}$ ${\displaystyle {\frac {k}{i}}\times \left(1-{\frac {1}{i+1}}\right)={\frac {k}{i+1}}}$ ${\displaystyle j}$ ${\displaystyle x_{r}}$ ${\displaystyle x_{i+1}}$ ${\displaystyle {\frac {k}{i+1}}{\frac {1}{k}}={\frac {1}{i+1}}}$

Hemos demostrado que la probabilidad de que un nuevo insumo entre en el yacimiento es igual a la probabilidad de que un insumo existente en el yacimiento se conserve. Por lo tanto, concluimos, por el principio de inducción matemática, que el algoritmo R produce efectivamente una muestra aleatoria uniforme de los insumos.

Si bien es conceptualmente simple y fácil de entender, este algoritmo necesita generar un número aleatorio para cada elemento de la entrada, incluidos los elementos que se descartan. El tiempo de ejecución asintótico del algoritmo es, por lo tanto , . La generación de esta cantidad de aleatoriedad y el tiempo de ejecución lineal hacen que el algoritmo sea innecesariamente lento si la población de entrada es grande. ${\displaystyle O(n)}$

Este es el algoritmo R, implementado de la siguiente manera:

(* S tiene elementos para muestrear, R contendrá el resultado *) ReservoirSample ( S [ 1 .. n ] , R [ 1 .. k ]) // llena la matriz del depósito para i := 1 a k R [ i ] := S [ i ] fin             // reemplaza elementos con probabilidad gradualmente decreciente para i := k + 1 a n (* randomInteger(a, b) genera un entero uniforme del rango inclusivo {a, ..., b} *) j := randomInteger ( 1 , i ) si j <= k R [ j ] := S [ i ] fin fin fin                    

Óptimo: Algoritmo L

Si generamos números aleatorios independientemente, entonces los índices del más pequeño de ellos son una muestra uniforme de los k-subconjuntos de . ${\displaystyle n}$ ${\displaystyle u_{1},...,u_{n}\sim U[0,1]}$ ${\displaystyle k}$ ${\displaystyle \{1,...,n\}}$

El proceso se puede realizar sin saber : ${\displaystyle n}$

Conserve el más pequeño de los que ha visto hasta ahora, así como , el índice del más grande entre ellos. Para cada nuevo , compárelo con . Si , descarte , almacene y establezca que sea el índice del más grande entre ellos. De lo contrario, descarte y establezca . ${\displaystyle k}$ ${\displaystyle u_{1},...,u_{i}}$ ${\displaystyle w_{i}}$ ${\displaystyle u_{i+1}}$ ${\displaystyle u_{w_{i}}}$ ${\displaystyle u_{i+1}<u_{w_{i}}}$ ${\displaystyle u_{w_{i}}}$ ${\displaystyle u_{i+1}}$ ${\displaystyle w_{i+1}}$ ${\displaystyle u_{i+1}}$ ${\displaystyle w_{i+1}=w_{i}}$

Ahora, combine esto con el flujo de entradas . Cada vez que se acepte algo, almacene el correspondiente . Cada vez que se descarte algo, descarte el correspondiente . ${\displaystyle x_{1},...,x_{n}}$ ${\displaystyle u_{i}}$ ${\displaystyle x_{i}}$ ${\displaystyle u_{i}}$ ${\displaystyle x_{i}}$

Este algoritmo aún necesita números aleatorios, por lo que lleva tiempo, pero se puede simplificar. ${\displaystyle O(n)}$ ${\displaystyle O(n)}$

Primera simplificación: no es necesario ir probando uno por uno, ya que la probabilidad de que la siguiente aceptación ocurra en es , es decir, el intervalo de aceptación sigue una distribución geométrica . ${\displaystyle u_{i+1},u_{i+2},...}$ ${\displaystyle u_{i+l}}$ ${\displaystyle (1-u_{w_{i}})^{l-1}\,u_{w_{i}}=(1-u_{w_{i}})^{l-1}-(1-u_{w_{i}})^{l}}$ ${\displaystyle l}$

Segunda simplificación: no es necesario recordar toda la matriz de los más pequeños que se han visto hasta ahora, sino sólo , el más grande entre ellos. Esto se basa en tres observaciones: ${\displaystyle k}$ ${\displaystyle u_{1},...,u_{i}}$ ${\displaystyle w}$

• Cada vez que se selecciona algo nuevo para ingresar al almacenamiento, se descarta una entrada aleatoria uniforme en el almacenamiento. ${\displaystyle x_{i+1}}$
• ${\displaystyle u_{i+1}\sim U[0,w]}$
• ${\displaystyle w_{i+1}}$tiene la misma distribución que , donde todos son independientes. Esto se puede muestrear primero muestreando y luego tomando . ${\displaystyle \max\{u_{1},...,u_{k}\}}$ ${\displaystyle u_{1},...,u_{k}\sim U[0,w]}$ ${\displaystyle u\sim U[0,1]}$ ${\displaystyle w\cdot u^{\frac {1}{k}}}$

Este es el algoritmo L , ^[2] que se implementa de la siguiente manera:

(* S tiene elementos para muestrear, R contendrá el resultado *) ReservoirSample ( S [ 1 .. n ] , R [ 1 .. k ]) // llena la matriz del depósito para i = 1 a k R [ i ] := S [ i ] fin             (* random() genera un número aleatorio uniforme (0,1) *) W := exp ( log ( random ()) / k )    mientras i <= n i := i + floor ( log ( random ()) / log ( 1 - W )) + 1 si i <= n (* reemplazar un elemento aleatorio del depósito con el elemento i *) R [ randomInteger ( 1 , k )] := S [ i ] // índice aleatorio entre 1 y k, inclusive W := W * exp ( log ( random ()) / k ) fin fin fin                          

Este algoritmo calcula tres números aleatorios para cada elemento que pasa a formar parte del reservorio y no dedica tiempo a los elementos que no lo forman. Por lo tanto, su tiempo de ejecución esperado es , ^[2], que es óptimo. ^[1] Al mismo tiempo, es fácil de implementar de manera eficiente y no depende de desviaciones aleatorias de distribuciones exóticas o difíciles de calcular. ${\displaystyle O(k(1+\log(n/k)))}$

Con ordenación aleatoria

Si asociamos a cada elemento de la entrada un número aleatorio generado uniformemente, los k elementos con los valores asociados más grandes (o, equivalentemente, más pequeños) forman una muestra aleatoria simple. ^[3] Un muestreo de reservorio simple mantiene así los k elementos con los valores asociados más grandes actualmente en una cola de prioridad .

(*  S es un flujo de elementos para muestrear  S.Current devuelve el elemento actual en el flujo  S.Next avanza el flujo a la siguiente posición  min-priority-queue admite:  Count -> número de elementos en la cola de prioridad  Minimum -> devuelve el valor de clave mínimo de todos los elementos  Extract-Min() -> Elimina el elemento con la clave mínima  Insert(key, Item) -> Agrega el elemento con la clave especificada *) ReservoirSample ( S [ 1 .. ? ]) H := new min - priority - queue while S tiene datos r := random () // uniformemente aleatorio entre 0 y 1, exclusivo si H . Count < k H . Insert ( r , S . Current ) else // mantiene k elementos con las claves asociadas más grandes si r > H . Minimum H . Extract - Min () H . Insert ( r , S . Current ) end S . Next end devuelve elementos en H end                                  

El tiempo de ejecución esperado de este algoritmo es y es relevante principalmente porque se puede extender fácilmente a elementos con pesos. ${\displaystyle O(n+k\log k\log(n/k))}$

Muestreo aleatorio ponderado

Los métodos presentados en las secciones anteriores no permiten obtener probabilidades de inclusión fijas a priori. Algunas aplicaciones requieren que las probabilidades de muestreo de los elementos sean de acuerdo con los pesos asociados con cada elemento. Por ejemplo, podría requerirse muestrear consultas en un motor de búsqueda con el peso como el número de veces que se realizaron para que la muestra pueda analizarse para determinar el impacto general en la experiencia del usuario. Sea el peso del elemento i y la suma de todos los pesos W . Hay dos formas de interpretar los pesos asignados a cada elemento del conjunto: ^[4] ${\displaystyle w_{i}}$

1. En cada ronda, la probabilidad de que cada elemento no seleccionado sea seleccionado en esa ronda es proporcional a su peso en relación con los pesos de todos los elementos no seleccionados. Si X es la muestra actual, entonces la probabilidad de que un elemento sea seleccionado en la ronda actual es . ${\displaystyle i\notin X}$ ${\displaystyle \textstyle w_{i}/(W-\sum _{j\in X}{w_{j}})}$
2. La probabilidad de que cada elemento se incluya en la muestra aleatoria es proporcional a su peso relativo, es decir, . Tenga en cuenta que esta interpretación podría no ser factible en algunos casos, por ejemplo, . ${\displaystyle w_{i}/W}$ ${\displaystyle k=n}$

Algoritmo A-Res

El siguiente algoritmo fue dado por Efraimidis y Spirakis que utiliza la interpretación 1: ^[5]

(*  S es un flujo de elementos para muestrear  S.Current devuelve el elemento actual en el flujo  S.Weight devuelve el peso del elemento actual en el flujo  S.Next avanza el flujo a la siguiente posición  El operador de potencia está representado por ^  min-priority-queue admite:  Count -> número de elementos en la cola de prioridad  Minimum() -> devuelve el valor de clave mínimo de todos los elementos  Extract-Min() -> Elimina el elemento con la clave mínima  Insert(key, Item) -> Agrega el elemento con la clave especificada *) ReservoirSample ( S [ 1 .. ? ]) H := new min - priority - queue while S tiene datos r := random () ^ ( 1 / S . Weight ) // random() produce un número aleatorio uniforme en (0,1) si H . Count < k H . Insert ( r , S . Current ) de lo contrario // mantiene k elementos con las claves asociadas más grandes si r > H . Minimum H . Extract - Min () H . Insert ( r , S . Current ) end end S . Próximo fin de devolver los artículos en el extremo H                                     

Este algoritmo es idéntico al algoritmo dado en Muestreo de yacimientos con ordenamiento aleatorio, excepto por la generación de las claves de los elementos. El algoritmo es equivalente a asignar a cada elemento una clave donde r es el número aleatorio y luego seleccionar los k elementos con las claves más grandes. De manera equivalente, una formulación numéricamente más estable de este algoritmo calcula las claves como y selecciona los k elementos con las claves más pequeñas . ^{[6] [ verificación fallida ]} ${\displaystyle r^{1/w_{i}}}$ ${\displaystyle -\ln(r)/w_{i}}$

Algoritmo A-ExpJ

El siguiente algoritmo es una versión más eficiente de A-Res, también propuesta por Efraimidis y Spirakis: ^[5]

(*  S es un flujo de elementos para muestrear  S.Current devuelve el elemento actual en el flujo  S.Weight devuelve el peso del elemento actual en el flujo  S.Next avanza el flujo a la siguiente posición  El operador de potencia está representado por ^  min-priority-queue admite:  Count -> número de elementos en la cola de prioridad  Minimum -> clave mínima de cualquier elemento en la cola de prioridad  Extract-Min() -> Elimina el elemento con la clave mínima  Insert(Key, Item) -> Agrega el elemento con la clave especificada *) ReservoirSampleWithJumps ( S [ 1 .. ? ]) H := new min - priority - queue while S tiene datos y H . Count < k r := random () ^ ( 1 / S . Weight ) // random() produce un número uniformemente aleatorio en (0,1) H . Insert ( r , S . Current ) S . Siguiente fin X := log ( random ()) / log ( H . Minimum ) // esta es la cantidad de peso que se debe saltar mientras S tenga datos X := X - S . Weight si X <= 0 t := H . Minimum ^ S . Weight r := random ( t , 1 ) ^ ( 1 / S . Weight ) // random(x, y) produce un número uniformemente aleatorio en (x, y) H . Extract - Min () H . Insert ( r , S . Current )                                                          X := log ( aleatorio ()) / log ( H . Mínimo ) fin S . Siguiente fin devolver elementos en H fin           

Este algoritmo sigue las mismas propiedades matemáticas que se utilizan en A-Res, pero en lugar de calcular la clave para cada elemento y comprobar si ese elemento debe insertarse o no, calcula un salto exponencial al siguiente elemento que se insertará. Esto evita tener que crear variables aleatorias para cada elemento, lo que puede resultar costoso. El número de variables aleatorias necesarias se reduce de a en expectativa, donde es el tamaño del depósito y es el número de elementos en el flujo. ^[5] ${\displaystyle O(n)}$ ${\displaystyle O(k\log(n/k))}$ ${\displaystyle k}$ ${\displaystyle n}$

Algoritmo A-Chao

Advertencia: la siguiente descripción es errónea, consulte el artículo original de Chao y la discusión aquí.

El siguiente algoritmo fue proporcionado por MT Chao y utiliza la interpretación 2: ^[7] y Tillé (2006). ^[8]

(*  S tiene elementos para muestrear, R contendrá el resultado  S[i].Weight contiene el peso de cada elemento *) WeightedReservoir - Chao ( S [ 1 .. n ] , R [ 1 .. k ]) WSum := 0 // rellenar la matriz del depósito para i := 1 a k R [ i ] := S [ i ] WSum := WSum + S [ i ] . Peso fin para i := k + 1 a n WSum := WSum + S [ i ] . Peso p := S [ i ] . Peso / WSum // probabilidad para este elemento j := random () ; // uniformemente aleatorio entre 0 y 1 si j <= p // seleccionar elemento según la probabilidad R [ randomInteger ( 1 , k )] := S [ i ] //selección uniforme en el depósito para reemplazo fin fin fin                                                    

Para cada elemento, se calcula su peso relativo y se utiliza para decidir aleatoriamente si el elemento se agregará al depósito. Si se selecciona el elemento, se selecciona de manera uniforme uno de los elementos existentes del depósito y se reemplaza por el nuevo elemento. El truco aquí es que, si las probabilidades de todos los elementos del depósito ya son proporcionales a sus pesos, entonces al seleccionar de manera uniforme qué elemento reemplazar, las probabilidades de todos los elementos siguen siendo proporcionales a su peso después del reemplazo.

Obsérvese que Chao no especifica cómo muestrear los primeros k elementos. Simplemente supone que tenemos otra forma de seleccionarlos en proporción a su peso. Chao: "Supongamos que tenemos un plan de muestreo de tamaño fijo con respecto a S_k en el momento A ; tal que su probabilidad de inclusión de primer orden de X_t es π(k; i) ".

Algoritmo A-Chao con saltos

De manera similar a otros algoritmos, es posible calcular un peso aleatorio jy restar los valores de masa de probabilidad de los elementos, omitiéndolos mientras j > 0se reduce la cantidad de números aleatorios que se deben generar. ^[4]

(*  S tiene elementos para muestrear, R contendrá el resultado  S[i].Weight contiene el peso de cada elemento *) WeightedReservoir - Chao ( S [ 1 .. n ] , R [ 1 .. k ]) WSum := 0 // llena la matriz del reservorio para i := 1 a k R [ i ] := S [ i ] WSum := WSum + S [ i ] . Peso final j := random () // uniformemente aleatorio entre 0 y 1 pNone := 1 // probabilidad de que no se haya seleccionado ningún elemento hasta ahora (en este salto) para i := k + 1 a n WSum := WSum + S [ i ] . Peso p := S [ i ] . Peso / WSum // probabilidad para este artículo j -= p * pNone pNone := pNone * ( 1 - p ) si j <= 0 R [ randomInteger ( 1 , k )] := S [ i ] //selección uniforme en el yacimiento para reemplazo j = random () pNone := 1 fin fin fin                                                                         

Relación con la mezcla Fisher-Yates

Supongamos que uno quisiera sacar k cartas al azar de una baraja de cartas. Un enfoque natural sería barajar la baraja y luego tomar las k cartas superiores. En el caso general, la baraja también debe funcionar incluso si no se conoce de antemano el número de cartas de la baraja, una condición que se cumple con la versión de adentro hacia afuera de la baraja de Fisher-Yates : ^[9]

(* S tiene la entrada, R contendrá la permutación de salida *) Shuffle ( S [ 1 .. n ] , R [ 1 .. n ]) R [ 1 ] := S [ 1 ] para i de 2 a n hacer j := randomInteger ( 1 , i ) // rango inclusivo R [ i ] := R [ j ] R [ j ] := S [ i ] fin fin                       

Tenga en cuenta que, aunque se barajan el resto de las cartas, solo las primeras k son importantes en el presente contexto. Por lo tanto, la matriz R solo necesita rastrear las cartas en las primeras k posiciones mientras realiza la baraja, lo que reduce la cantidad de memoria necesaria. Al truncar R a una longitud k , el algoritmo se modifica en consecuencia:

(* S tiene elementos para muestrear, R contendrá el resultado *) ReservoirSample ( S [ 1..n ] , R [ 1..k ] ) R [ 1 ] := S [ 1 ] para i de 2 a k hacer j := randomInteger ( 1 , i ) // rango inclusivo R [ i ] := R [ j ] R [ j ] := S [ i ] fin para i de k + 1 a n hacer j := randomInteger ( 1 , i ) // rango inclusivo si ( j <= k ) R [ j ] : = S [ i ] fin fin fin                                              

Dado que el orden de las primeras k tarjetas es irrelevante, se puede eliminar el primer bucle y se puede inicializar R para que sean los primeros k elementos de la entrada. Esto produce el algoritmo R .

Limitaciones

El muestreo de yacimientos parte del supuesto de que la muestra deseada cabe en la memoria principal , lo que a menudo implica que k es una constante independiente de n . En aplicaciones en las que nos gustaría seleccionar un subconjunto grande de la lista de entrada (por ejemplo, un tercero, es decir, ), se deben adoptar otros métodos. Se han propuesto implementaciones distribuidas para este problema. ^[10] ${\displaystyle k=n/3}$

Referencias

1. Vitter, Jeffrey S. (1 de marzo de 1985). "Muestreo aleatorio con un reservorio" (PDF) . ACM Transactions on Mathematical Software . 11 (1): 37–57. CiteSeerX  10.1.1.138.784 . doi :10.1145/3147.3165. S2CID  17881708.
2. Li, Kim-Hung (4 de diciembre de 1994). "Algoritmos de muestreo de yacimientos de complejidad temporal O(n(1+log(N/n)))". ACM Transactions on Mathematical Software . 20 (4): 481–493. doi : 10.1145/198429.198435 . S2CID  15721242.
3. Fan, C.; Muller, ME; Rezucha, I. (1962). "Desarrollo de planes de muestreo mediante técnicas de selección secuencial (ítem por artículo) y computadoras digitales". Revista de la Asociación Estadounidense de Estadística . 57 (298): 387–402. doi :10.1080/01621459.1962.10480667. JSTOR  2281647.
4. Efraimidis, Pavlos S. (2015). "Muestreo aleatorio ponderado sobre flujos de datos". Algoritmos, probabilidad, redes y juegos . Apuntes de clase en informática. Vol. 9295. págs. 183–195. arXiv : 1012.0256 . doi :10.1007/978-3-319-24024-4_12. ISBN . 978-3-319-24023-7. Número de identificación del sujeto  2008731.
5. Efraimidis, Pavlos S.; Spirakis, Paul G. (16 de marzo de 2006). "Muestreo aleatorio ponderado con un reservorio". Cartas de procesamiento de la información . 97 (5): 181–185. doi :10.1016/j.ipl.2005.11.003.
6. Arratia, Richard (2002). Bela Bollobas (ed.). "Sobre la cantidad de dependencia en la factorización prima de un entero aleatorio uniforme". Combinatoria contemporánea . 10 : 29–91. CiteSeerX 10.1.1.745.3975 . ISBN  978-3-642-07660-2.
7. Chao, MT (1982). "Un plan de muestreo de probabilidad desigual de propósito general". Biometrika . 69 (3): 653–656. doi :10.1093/biomet/69.3.653.
8. Tillé, Yves (2006). Algoritmos de muestreo . Springer. ISBN 978-0-387-30814-2.
9. Consejo Nacional de Investigación (2013). Fronteras en el análisis masivo de datos. The National Academies Press. pág. 121. ISBN 978-0-309-28781-4.
10. Muestreo de yacimientos en MapReduce