¿Qué resultado importa?
Aquí hay un escenario común: se realizó una prueba A/B, donde se seleccionó una muestra aleatoria de unidades (por ejemplo) para una campaña y recibieron tratamiento A. Se seleccionó otra muestra para recibir tratamiento B. “A” podría ser un Comunicación u oferta y “B” no podría ser comunicación o ninguna oferta. “A” podría tener un 10% de descuento y “B” podría ser un 20% de descuento. Dos grupos, dos tratamientos diferentes, donde A y B son dos tratamientos discretos, pero sin pérdida de generalidad a mayores que 2 tratamientos y tratamientos continuos.
Por lo tanto, la campaña se ejecuta y los resultados están disponibles. Con nuestro sistema de backend, podemos rastrear cuáles de estas unidades tomaron la acción de interés (por ejemplo, hicieron una compra) y cuál no. Además, para los que lo hicieron, registramos la intensidad de esa acción. Un escenario común es que podemos rastrear las cantidades de compra para las que compraron. Esto a menudo se llama un monto promedio de pedido o ingresos por métrica del comprador. O cien nombres diferentes que significan lo mismo: para los que compraron, ¿cuánto gastaron, en promedio?
Para algunos casos de uso, el vendedor está interesado en la métrica anterior: la tasa de compra. Por ejemplo, ¿manejamos más compradores (potencialmente por primera vez) en nuestra campaña de adquisición con el tratamiento A o B? A veces, estamos interesados en conducir los ingresos por comprador más alto, por lo que ponemos énfasis en este último.
Sin embargo, más a menudo, estamos interesados en impulsar los ingresos de manera rentable y lo que realmente nos importa es los ingresos que produjo la campaña. en general. ¿El tratamiento A o B generó más ingresos? No siempre tenemos tamaños de muestra equilibrados (tal vez debido al costo o la evitación de riesgos), por lo que dividimos los ingresos medidos por el número de candidatos que fueron tratados en cada grupo (llame a estos recuentos N_A y N_B). Queremos comparar esta medida entre los dos grupos, por lo que el contraste estándar es simplemente:
Este es solo el ingreso promedio para el tratamiento, un ingreso medio menos medio para el tratamiento B, donde esa media se toma sobre todo el conjunto de unidades específicas, independientemente de si respondieron o no. Su interpretación es también sencilla: ¿cuál es el aumento promedio de ingresos por unidad promovida que aumenta del tratamiento A versus el tratamiento B?
Por supuesto, esta última medida representa los dos anteriores: la tasa de respuesta multiplicada por los ingresos medios por respondedor.
¿Incertidumbre?
Cuánto gasta un comprador es muy variable y un par de compras grandes en un grupo de tratamiento u otro pueden sesgar la media significativamente. Del mismo modo, la variación de la muestra puede ser significativa. Por lo tanto, queremos entender cuán seguros estamos en esta comparación de medios y cuantificar la “importancia” de la diferencia observada.
Entonces, arrojas los datos en una prueba t y miras el valor p. ¡Pero espera! Desafortunadamente para el vendedor, la gran mayoría de las veces, la tasa de compra es relativamente baja (a veces muy baja) y, por lo tanto, hay muchos valores de ingresos cero, a menudo la gran mayoría. Los supuestos de la prueba t pueden ser violadas gravemente. Los tamaños de muestra muy grandes pueden llegar al rescate, pero hay una forma más de principios de analizar estos datos que son útiles de múltiples maneras, que se explicará.
Conjunto de datos de ejemplo
Comencemos con el conjunto de datos de muestra para hacer las cosas prácticas. Uno de mis conjuntos de datos de marketing directo favoritos es de la KDD Cup 98.
url="https://kdd.ics.uci.edu/databases/kddcup98/epsilon_mirror/cup98lrn.zip"
filename="cup98LRN.txt"
r = requests.get(url)
z = zipfile.ZipFile(io.BytesIO(r.content))
z.extractall()
pdf_data = pd.read_csv(filename, sep=',')
pdf_data = pdf_data.query('TARGET_D >=0')
pdf_data['TREATMENT'] = np.where(pdf_data.RFA_2F >1,'A','B')
pdf_data['TREATED'] = np.where(pdf_data.RFA_2F >1,1,0)
pdf_data['GT_0'] = np.where(pdf_data.TARGET_D >0,1,0)
pdf_data = pdf_data[['TREATMENT', 'TREATED', 'GT_0', 'TARGET_D']]
En el fragmento de código anterior, estamos descargando un archivo zip (el conjunto de datos de aprendizaje específicamente), extrayendo y leyéndolo en un marco de datos de pandas. La naturaleza de este conjunto de datos es el historial de campaña de una organización sin fines de lucro que buscaba donaciones a través de correos directos. No hay variantes de tratamiento dentro de este conjunto de datos, por lo que estamos fingiendo en su lugar y segmentamos el conjunto de datos en función de la frecuencia de donaciones pasadas. Llamamos a este indicador TRATAMIENTO (como el categórico y crear Tratado como el indicador binario para ‘a’). Considere esto los resultados de un ensayo de control aleatorio donde una porción de la población de la muestra fue tratada con una oferta y el resto no lo fueron. Hacemos un seguimiento de cada individuo y acumulamos la cantidad de su donación.
Entonces, si examinamos este conjunto de datos, vemos que hay alrededor de 95,000 individuos promovidos, generalmente distribuidos por igual en los dos tratamientos:

El tratamiento A tiene una tasa de respuesta mayor, pero en general la tasa de respuesta en el conjunto de datos es de solo alrededor del 5%. Entonces, tenemos 95% ceros.

Para aquellos que donaron, el tratamiento A parece estar asociado con una cantidad de donación promedio más baja.

Combinando juntos a todos los que fueron atacados, el tratamiento A parece estar asociado con una cantidad de donación promedio más alta (la tasa de respuesta más alta supera la cantidad de donación más baja para los respondedores, pero no mucho.

Finalmente, el histograma de la cantidad de donación se muestra aquí, se agrupa sobre ambos tratamientos, lo que ilustra la masa en cero y un sesgo derecho.

Un resumen numérico de los dos grupos de tratamiento cuantifica el fenómeno observado anteriormente, mientras que el tratamiento A parece haber impulsado una respuesta significativamente mayor, los que fueron tratados con un promedio donado menos en promedio. La red de estas dos medidas, la que finalmente buscamos, la donación media general por unidad dirigida, parece ser más alta para el tratamiento A. cuán seguros estamos en ese hallazgo es el tema de este análisis.

Obstáculo gamma
Una forma de modelar estos datos y responder nuestra pregunta de investigación en términos de la diferencia entre los dos tratamientos para generar la donación promedio por unidad dirigida es con la distribución de obstáculos gamma. Similar a la distribución más conocida de Poisson inflado cero (ZIP) o NB (ZinB), esta es una distribución de la mezcla donde una parte se refiere a la masa en cero y la otra, en los casos en que la variable aleatoria es positiva, la densidad gamma función.

Aquí π representa la probabilidad de que la variable aleatoria Y sea> 0. En otras palabras, es la probabilidad del proceso gamma. Del mismo modo, (1- π) es la probabilidad de que la variable aleatoria sea cero. En términos de nuestro problema, esto se refiere a la probabilidad de que se realice una donación y, de ser así, es valor.
Comencemos con las partes componentes del uso de esta distribución en una regresión de regresión: logística y regresión gamma.
Regresión logística
La función logit es la función de enlace aquí, relacionada con las probabilidades de registro con la combinación lineal de nuestras variables predictoras, que con una sola variable como nuestro indicador de tratamiento binario, parece:
Donde π representa la probabilidad de que el resultado sea un evento “positivo” (denotado como 1), como una compra y (1-π) representa la probabilidad de que el resultado sea un evento “negativo” (denotado como 0). Además, π, que es la cantidad de interés anterior, se define por la función de logit inverso:
Ajustar este modelo es muy simple, necesitamos encontrar los valores de las dos betas que maximizan la probabilidad de los datos (el resultado y), que supone que las observaciones de NIID son:

Podríamos usar cualquiera de las múltiples bibliotecas para adaptarse rápidamente a este modelo, pero demostraremos PYMC como el medio para construir una regresión logística bayesiana simple.
Sin ninguno de los pasos normales del flujo de trabajo bayesiano, se ajustamos a este modelo simple usando MCMC.
import pymc as pm
import arviz as az
from scipy.special import expit
with pm.Model() as logistic_model:
# noninformative priors
intercept = pm.Normal('intercept', 0, sigma=10)
beta_treat = pm.Normal('beta_treat', 0, sigma=10)
# linear combination of the treated variable
# through the inverse logit to squish the linear predictor between 0 and 1
p = pm.invlogit(intercept + beta_treat * pdf_data.TREATED)
# Individual level binary variable (respond or not)
pm.Bernoulli(name="logit", p=p, observed=pdf_data.GT_0)
idata = pm.sample(nuts_sampler="numpyro")
az.summary(idata, var_names=['intercept', 'beta_treat'])

Si construimos un contraste de las dos tasas de respuesta media del tratamiento, encontramos que, como se esperaba, el aumento de la tasa de respuesta media para el tratamiento A es 0.026 más grande que el tratamiento B con un intervalo creíble del 94% de (0.024, 0.029).
# create a new column in the posterior which contrasts Treatment A - B
idata.posterior['TREATMENT A - TREATMENT B'] = expit(idata.posterior.intercept + idata.posterior.beta_treat) - expit(idata.posterior.intercept)
az.plot_posterior(
idata,
var_names=['TREATMENT A - TREATMENT B']
)

Regresión gamma
El siguiente componente es la distribución gamma con una de sus parametrizaciones de su función de densidad de probabilidad, como se muestra arriba:
Esta distribución se define para variables aleatorias estrictamente positivas y si se usa en negocios para valores como costos, gastos de demanda de clientes y montos de reclamos de seguro.
Dado que la media y la varianza de gamma se definen en términos de α y β de acuerdo con las fórmulas:
Para la regresión gamma, podemos parametrizar por α y β o por μ y σ. Si hacemos μ definido como una combinación lineal de variables predictoras, entonces podemos definir gamma en términos de α y β usando μ:
El modelo de regresión gamma supone (en este caso, el enlace inverso es otra opción común) el enlace de registro que está destinado a “linealizar” la relación entre el predictor y el resultado:
Después de casi exactamente la misma metodología que para la tasa de respuesta, limitamos el conjunto de datos a solo respondedores y se ajustamos a la regresión gamma usando PYMC.
with pm.Model() as gamma_model:
# noninformative priors
intercept = pm.Normal('intercept', 0, sigma=10)
beta_treat = pm.Normal('beta_treat', 0, sigma=10)
shape = pm.HalfNormal('shape', 5)
# linear combination of the treated variable
# through the exp to ensure the linear predictor is positive
mu = pm.Deterministic('mu',pm.math.exp(intercept + beta_treat * pdf_responders.TREATED))
# Individual level binary variable (respond or not)
pm.Gamma(name="gamma", alpha = shape, beta = shape/mu, observed=pdf_responders.TARGET_D)
idata = pm.sample(nuts_sampler="numpyro")
az.summary(idata, var_names=['intercept', 'beta_treat'])

# create a new column in the posterior which contrasts Treatment A - B
idata.posterior['TREATMENT A - TREATMENT B'] = np.exp(idata.posterior.intercept + idata.posterior.beta_treat) - np.exp(idata.posterior.intercept)
az.plot_posterior(
idata,
var_names=['TREATMENT A - TREATMENT B']
)

Una vez más, como se esperaba, vemos que la elevación media para el tratamiento A tiene un valor esperado igual al valor de muestra de -7.8. El intervalo creíble del 94% es (-8.3, -7.3).
Los componentes, la tasa de respuesta y la cantidad promedio por respondedor que se muestra arriba son lo más simple posible. Pero, es una extensión directa para agregar predictores adicionales para 1) estimar los efectos de tratamiento promedio condicional (CATE) cuando esperamos que el efecto del tratamiento difiera en segmento o 2) reducir la varianza de la estimación del efecto de tratamiento promedio mediante el acondicionamiento Variables previas al tratamiento.
Regresión del modelo de obstáculo (gamma)
En este punto, debería ser bastante sencillo ver dónde estamos progresando. Para el modelo de obstáculo, tenemos una probabilidad condicional, dependiendo de si la observación específica es 0 o mayor que cero, como se muestra arriba para la distribución de obstáculos gamma. Podemos ajustar los dos modelos de componentes (regresión logística y gamma) simultáneamente. Obtenemos gratis, su producto, que en nuestro ejemplo es una estimación del monto de la donación por unidad específica.
No sería difícil ajustar este modelo con el uso de una función de probabilidad con una declaración de conmutación dependiendo del valor de la variable de resultado, pero PYMC tiene esta distribución ya codificada para nosotros.
import pymc as pm
import arviz as az
with pm.Model() as hurdle_model:
## noninformative priors ##
# logistic
intercept_lr = pm.Normal('intercept_lr', 0, sigma=5)
beta_treat_lr = pm.Normal('beta_treat_lr', 0, sigma=1)
# gamma
intercept_gr = pm.Normal('intercept_gr', 0, sigma=5)
beta_treat_gr = pm.Normal('beta_treat_gr', 0, sigma=1)
# alpha
shape = pm.HalfNormal('shape', 1)
## mean functions of predictors ##
p = pm.Deterministic('p', pm.invlogit(intercept_lr + beta_treat_lr * pdf_data.TREATED))
mu = pm.Deterministic('mu',pm.math.exp(intercept_gr + beta_treat_gr * pdf_data.TREATED))
## likliehood ##
# psi is pi
pm.HurdleGamma(name="hurdlegamma", psi=p, alpha = shape, beta = shape/mu, observed=pdf_data.TARGET_D)
idata = pm.sample(cores = 10)
Si examinamos el resumen de rastreo, vemos que los resultados son exactamente los mismos para los dos modelos de componentes.

Como se señaló, la media de la distribución de obstáculos gamma es π * μ para que podamos crear un contraste:
# create a new column in the posterior which contrasts Treatment A - B
idata.posterior['TREATMENT A - TREATMENT B'] = ((expit(idata.posterior.intercept_lr + idata.posterior.beta_treat_lr))* np.exp(idata.posterior.intercept_gr + idata.posterior.beta_treat_gr)) - \
((expit(idata.posterior.intercept_lr))* np.exp(idata.posterior.intercept_gr))
az.plot_posterior(
idata,
var_names=['TREATMENT A - TREATMENT B']
El valor medio esperado de este modelo es 0.043 con un intervalo creíble del 94% de (-0.0069, 0.092). Podríamos interrogar el posterior para ver qué proporción de veces se prevé que la donación por comprador sea mayor para el tratamiento y cualquier otra función de decisión que tuviera sentido para nuestro caso, incluida la adición de un P&L más completo a la estimación (es decir, incluidos los márgenes y el costo) .

Notas: Algunas implementaciones parametrizan el modelo de obstáculo gamma de manera diferente cuando la probabilidad de cero es π y, por lo tanto, la media del obstáculo gamma implica (1-π) en su lugar. También tenga en cuenta que en el momento de este escrito parece haber un asunto Con las muestras de nueces en PYMC y tuvimos que recurrir a la implementación predeterminada de Python para ejecutar el código anterior.
Resumen
Con este enfoque, obtenemos la misma inferencia para ambos modelos por separado y el beneficio adicional de la tercera métrica. El ajuste de estos modelos con PYMC nos permite todos los beneficios del análisis bayesiano, incluida la inyección del conocimiento del dominio anterior y un completo posterior para responder preguntas y cuantificar la incertidumbre.
Créditos:
- Todas las imágenes son los autores, a menos que se indique lo contrario.
- El conjunto de datos utilizado es de la Copa KDD 98 patrocinada por Epsilon. https://kdd.ics.uci.edu/databases/kddcup98/kddcup98.html (CC por 4.0)