Test A/B para el Bandido Multibrazo (Multi-Armed Bandit)

Recientemente hemos visto el problema del Bandido Multibrazo (Multi-Armed Bandit). Una de las posibles soluciones que tenemos en nuestra mano para resolver este problema es utilizar un Test A/B. Esto es, evaluar durante un periodo de tiempo todos los bandidos por igual y decidir una vez finalizado este periodo de prueba cuál es el óptimo. O, si los datos no son concluyentes, volver a realizar otra prueba con más muestras. Vamos a ver cómo simular la utilización de un Test A/B para el Bandido Multibrazo en Python.

Creación de un objeto para simular el bandido

En primer lugar, para llevar a cabo la simulación es necesario crear una clase que represente a cada uno de los bandidos. Clase en la que solamente es necesario crear un método para simular el resultado de una tirada. Para lo que se puede usar una función de distribución como puede ser la binomial. Siendo los parámetros de la función de distribución diferentes en cada uno bandidos. Así se puede crear la siguiente clase llamada Bandit.

En la que el constructor cuenta con una propiedad llamada mean que es la probabilidad de éxito para el bandido. Debido a que se ultima una distribución binomial con un único ensayo este valor coincidirá la recompensa media.

Nótese que en el objeto también se guarda el histórico de recompensas, aunque esto es para facilitar el posterior análisis.

Comprobar cuál es el mejor bandido

Ahora deberíamos simular el comportamiento de dos o más bandidos durante un número de casos. Una vez obtenidos los resultados debemos comprobar cuál es el que ofrece la mayor recompensa promedio y comprobar que el resultado no es debido al azar, esto es, la diferencia entre las dos muestras es estadísticamente significativa. Para lo que se puede utilizar por ejemplo un test Z. Un test que podemos encontrar implementado en la función proportions_ztest de la librería Statsmodels de Python.

La función proportions_ztest requiere dos tuplas como parámetros de entrada. La primera debe contener el número de casos favorables y el segundo el número de observaciones. Como resultado se obtiene una tupla con el resultado del z-test y el p-valor de este. Indicando este último la posibilidad de rechazar la hipótesis nula, que en este caso es que los dos resultados provienen de la misma distribución.

Aplicación a múltiples bandidos

La función proportions_ztest solamente funciona con pares, por lo que si hay más de dos bandidos es necesario utilizar una estrategia para probar todas las combinaciones. Como lo que se quiere es identificar el mejor bandido de todos solamente hay que comprobar que los valores de este pasan el test con el resto de los bandidos. Si esto es así, se puede afirmar que el resultado no es debido al azar.

Implementación en Python

Con todo lo que hemos visto se puede realizar la siguiente implementación en Python.

import numpy as np
from itertools import compress
from statsmodels.stats.proportion import proportions_ztest 

np.random.seed(0)

class Bandit:
    def __init__(self, mean):
        self.mean = mean
        self.rewards = []
    
    def pull(self):
        reward = np.random.binomial(1, self.mean)
        self.rewards.append(reward)
        return reward

original = [Bandit(0.02), Bandit(0.04), Bandit(0.06), Bandit(0.08), Bandit(0.10)]

bandits = original
simulations = 500

num_bandits = len(bandits)
num_iter = 0

while num_bandits > 1 and num_iter < 10:
    num_iter = num_iter + 1
    
    for sim in range(simulations):
        for bandit in range(num_bandits):
            bandits[bandit].pull()

    n_reward = [np.sum(bandit.rewards) for bandit in bandits]
    max_reward = max(n_reward)
    nobs = len(bandits[0].rewards)

    p_value = list(map(lambda a : proportions_ztest((a, max_reward), (nobs, nobs))[1], n_reward))
    bandits = list(compress(bandits, np.array(p_value) > 0.05))
    num_bandits = len(bandits)
    
evaluations = len(bandits[0].rewards)
total_evaluations = np.sum([len(bandit.rewards) for bandit in original])
total_reward = np.sum([np.sum(bandit.rewards) for bandit in original])
avg_reward = total_reward / total_evaluations

En este ejemplo se han creado cinco bandidos con diferente probabilidad de éxito. En un bucle while se simula 500 veces cada uno de los bandidos y posteriormente se extrae el número de casos existimos en cada uno de ellos. Lo que se almacena en el vector n_reward. A partir del cual se puede extraer fácilmente la recompensa máxima que se almacena en max_reward.

Con esto y el número de casos ( nobs) se puede obtener el p-valor ( pvalue) de cada una de las muestras respecto al máximo. Seleccionado aquellos bandidos para los cuales no se puede rechazar la hipótesis nula. Usando para ello una probabilidad de 0.05.

Si esta prueba la pasa solamente un bandido será el mejor. En caso contrario es necesario aumentar el número de muestras hasta que solamente quede uno o se alcance el máximo de iteraciones permitidas. Lo que se tienen que incluir para evitar bucles infinitos.

Resultados

Si se ejecuta el código se puede ver que es necesario 7 iteraciones, 3500, simulaciones para identificar el mejor bandido. Aunque solamente para los dos mejores, ya que los tres peores se han podido descartar después de una sola iteración. Resultados que pueden cambiar ligeramente en caso de que se modifique el valor de la semilla.

Con esto, después de 8500 simulaciones se han obtenido 692 casos positivos, lo que corresponde a un 8,1%. Por debajo del 10% máximo que se esperaría del mejor bandido. Aunque una vez obtenido este resultado ya se puede jugar solo con este bandido y mejorar la recompensa media.

Conclusiones

En esta entrada hemos visto cómo se puede emplear un Test A/B para el Bandido Multibrazo. Pudiendo identificar rápidamente el bandido que ofrece la mayor recompensa. Este es un método que requiere muchas pruebas antes de obtener unos resultados concluyentes, por eso podemos estar mucho tiempo jugando con bandidos que no son óptimos. En futuras entregas veamos como hacer esto con otras estrategias como pueden ser Epsilon-greedy, valores iniciales optimistas o UCB1.

Imagen de PIRO4D en Pixabay