[PYTHON] Analyse nach Bayes'schem Denken (1) ... Was ist besser, A oder B?

Beispiel einer Bayes'schen Argumentationsanalyse: Was ist besser, A oder B?

Aus Bayes'schen Überlegungen in Python

"In Python erlebte Bayes-Inferenz"

Was ist besser, Standort A oder Standort B, in Bayes'scher Argumentation? Ist das Problem.

Frage: Was ist die wahre Conversion für Site A?

$ P_A $ sei die Wahrscheinlichkeit, dass ein Benutzer, der Site A sieht, schließlich zu einer profitablen Aktion (dies wird als Konvertierung bezeichnet) führt, z. B. zum Anfordern von Materialien oder zum Kaufen. Angenommen, $ N $ Personen sehen sich Site A an und $ n $ Personen führen zu Conversions. Dann

p_A=n/N

Sie könnten denken, aber es ist nicht. Das heißt, wir wissen nicht, ob $ p_A $ und $ n / N $ gleich sind. Es gibt einen Unterschied zwischen der beobachteten Frequenz und der wahren Frequenz, und die wahre Frequenz kann als die Wahrscheinlichkeit interpretiert werden, dass ein Ereignis eintritt. Mit anderen Worten, seit ich Site A gesehen habe, weiß ich nicht, ob dies zu der Aktion geführt hat. Am Beispiel eines Würfels beträgt die wahre Häufigkeit des Würfelns und des Erhaltens einer 1 $ 1/6 , aber selbst wenn Sie die Würfel tatsächlich 6 Mal fallen lassen, erhalten Sie möglicherweise nicht einmal eine 1. Das ist die Häufigkeit der Beobachtung. Leider ist die Realität zu kompliziert und Rauschen ist unvermeidbar, sodass wir die wahre Frequenz nicht kennen und keine andere Wahl haben, als die wahre Frequenz aus der beobachteten Frequenz zu schätzen. Aus der Bayes'schen Statistik wird die wahre Häufigkeit ( p_A ) aus den vorherigen Verteilungs- und Beobachtungsdaten ( n / N $) geschätzt.

Zunächst müssen wir die vorherige Verteilung für die unbekannte Zahl $ p_A $ berücksichtigen. Was halten Sie von $ p_A $, wenn Sie die Daten nicht haben? Ich bin mir noch nicht sicher, also nehmen wir an, dass $ p_A $ einer gleichmäßigen Verteilung von [0,1] folgt. Geben Sie der Variablen "p" eine gleichmäßige Verteilung von 0 bis 1 mit "@ pm.Uniform".

import pymc3 as pm

# The parameters are the bounds of the Uniform.
with pm.Model() as model:
    p = pm.Uniform('p', lower=0, upper=1)

Beispiel: $ p_A = 0,05 $ und $ N = 1500 $ Simulieren wir, wie viele Benutzer konvertiert haben (gekaufte oder angeforderte Materialien, indem sie sich die Site angesehen haben), als der Benutzer Site A angesehen hat. Eigentlich kenne ich $ p_A $ nicht, aber ich gehe davon aus, dass ich es vorerst kenne und simuliere es, um Beobachtungsdaten zu erstellen.

Da es N Personen gibt, um N Versuche zu simulieren, [Bernouy-Verteilung](https://ja.wikipedia.org/wiki/%E3%83%99%E3%83%AB%E3%83% 8C% E3% 83% BC% E3% 82% A4% E5% 88% 86% E5% B8% 83) wird verwendet. Da die Bernoulli-Verteilung eine Wahrscheinlichkeitsverteilung für binäre Variablen ist (Variablen, die entweder 0 oder 1 annehmen), ist sie eine gute Wahrscheinlichkeitsverteilung für die Berechnung, ob konvertiert werden soll oder nicht.

 X\sim \text{Ber}(p) 

Zu diesem Zeitpunkt ist $ X $ $ 1 $ mit einer Wahrscheinlichkeit von $ p $ und $ 0 $ mit einer Wahrscheinlichkeit von $ 1-p $.

Das Python-Skript, das dies simuliert, ist

import scipy.stats as stats
import numpy as np

#set constants
p_true = 0.05  # remember, this is unknown.
N = 1500

# sample N Bernoulli random variables from Ber(0.05).
# each random variable has a 0.05 chance of being a 1.
# this is the data-generation step
occurrences = stats.bernoulli.rvs(p_true, size=N)

print(occurrences) # Remember: Python treats True == 1, and False == 0
print(np.sum(occurrences))

#[0 0 0 ... 0 0 0]
#67

# Occurrences.mean is equal to n/N.
print("What is the observed frequency in Group A? %.4f" % np.mean(occurrences))
print("Does this equal the true frequency? %s" % (np.mean(occurrences) == p_true))

#What is the observed frequency in Group A? 0.0447
#Does this equal the true frequency? False

Nachdem wir die beobachteten Daten (Vorkommen) haben, setzen wir die beobachteten Daten in die obs-Variable von "PyMC" und führen den Inferenzalgorithmus aus, um ein Diagramm der posterioren Verteilung von unbekanntem $ p_A $ zu erstellen.

#include the observations, which are Bernoulli
with model:
    obs = pm.Bernoulli("obs", p, observed=occurrences)
    # To be explained in chapter 3
    step = pm.Metropolis()
    trace = pm.sample(18000, step=step,cores =1)
    burned_trace = trace[1000:]

import matplotlib.pyplot as plt
plt.title("Posterior distribution of $p_A$, the true effectiveness of site A")
plt.vlines(p_true, 0, 90, linestyle="--",color="r", label="true $p_A$ (unknown)")
plt.hist(burned_trace["p"], bins=25, histtype="bar",density=True)
plt.legend()
plt.show()

Das resultierende Diagramm sieht folgendermaßen aus. 0811_pagoodFigure 2020-08-11 161933.png

Diesmal ist der Graph der Nachwahrscheinlichkeit so, während die Vorwahrscheinlichkeit des wahren Wertes 0,5 beträgt. Wenn Sie dies wiederholen, erhalten Sie jedes Mal ein etwas anderes Diagramm.

Frage: Vergleich von Standort A und Standort B, welcher ist besser?

Wenn Stelle B auf die gleiche Weise berechnet wird, kann die hintere Wahrscheinlichkeit von $ p_B $ erhalten werden. Also habe ich beschlossen, $ delta = p_A-p_B $ zu berechnen.

Wie $ p_A $ kennt $ p_B $ den wahren Wert nicht. Nehmen wir an, dass $ p_B = 0,4 $ ist. Dann ist $ delta = 0,1 $. Berechnen Sie $ p_B $ und $ delta $ auf die gleiche Weise wie $ p_A $. Ich habe alles noch einmal berechnet.

# -*- coding: utf-8 -*-

import pymc3 as pm
import scipy.stats as stats
import numpy as np
import matplotlib.pyplot as plt

#these two quantities are unknown to us.
true_p_A = 0.05
true_p_B = 0.04

#notice the unequal sample sizes -- no problem in Bayesian analysis.
N_A = 1500
N_B = 750

#generate some observations
observations_A = stats.bernoulli.rvs(true_p_A, size=N_A)
observations_B = stats.bernoulli.rvs(true_p_B, size=N_B)
print("Obs from Site A: ", observations_A[:30], "...")
print("Obs from Site B: ", observations_B[:30], "...")

#Obs from Site A:  [0 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] ...
#Obs from Site B:  [0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] ...


print(np.mean(observations_A))
#0.04666666666666667
print(np.mean(observations_B))
#0.03866666666666667

# Set up the pymc3 model. Again assume Uniform priors for p_A and p_B.
with pm.Model() as model:
    p_A = pm.Uniform("p_A", 0, 1)
    p_B = pm.Uniform("p_B", 0, 1)
    
    # Define the deterministic delta function. This is our unknown of interest.
    delta = pm.Deterministic("delta", p_A - p_B)

    
    # Set of observations, in this case we have two observation datasets.
    obs_A = pm.Bernoulli("obs_A", p_A, observed=observations_A)
    obs_B = pm.Bernoulli("obs_B", p_B, observed=observations_B)

    # To be explained in chapter 3.
    step = pm.Metropolis()
    trace = pm.sample(20000, step=step,cores=1)
    burned_trace=trace[1000:]


p_A_samples = burned_trace["p_A"]
p_B_samples = burned_trace["p_B"]
delta_samples = burned_trace["delta"]

#histogram of posteriors

ax = plt.subplot(311)

plt.xlim(0, .1)
plt.hist(p_A_samples, histtype='stepfilled', bins=25, alpha=0.85,
         label="posterior of $p_A$", color="#A60628", density=True)
plt.vlines(true_p_A, 0, 80, linestyle="--", label="true $p_A$ (unknown)")
plt.legend(loc="upper right")
plt.title("Posterior distributions of $p_A$, $p_B$, and delta unknowns")

ax = plt.subplot(312)

plt.xlim(0, .1)
plt.hist(p_B_samples, histtype='stepfilled', bins=25, alpha=0.85,
         label="posterior of $p_B$", color="#467821", density=True)
plt.vlines(true_p_B, 0, 80, linestyle="--", label="true $p_B$ (unknown)")
plt.legend(loc="upper right")

ax = plt.subplot(313)
plt.hist(delta_samples, histtype='stepfilled', bins=30, alpha=0.85,
         label="posterior of delta", color="#7A68A6", density=True)
plt.vlines(true_p_A - true_p_B, 0, 60, linestyle="--",
           label="true delta (unknown)")
plt.vlines(0, 0, 60, color="black", alpha=0.2)
plt.legend(loc="upper right");

plt.show()

# Count the number of samples less than 0, i.e. the area under the curve
# before 0, represent the probability that site A is worse than site B.
print("Probability site A is WORSE than site B: %.3f" % \
    np.mean(delta_samples < 0))

print("Probability site A is BETTER than site B: %.3f" % \
    np.mean(delta_samples > 0))

#Probability site A is WORSE than site B: 0.201
#Probability site A is BETTER than site B: 0.799

Das resultierende Diagramm sieht folgendermaßen aus. 0811_web_abFigure_1.png

Es ist ersichtlich, dass die Basis der posterioren Verteilung von $ p_B $ etwas breiter ist. Dies liegt daran, dass Site B weniger Daten enthält. Dies weist darauf hin, dass wir aufgrund des Mangels an Daten auf Site B nicht sicher sind, welchen Wert $ p_B $ hat.

Der dritte Graph, $ delta = p_A-p_B $, zeigt einen Wert größer als 0. Dies bedeutet, dass Site A als Conversion besser auf Benutzer reagiert.

Sie können ein besseres Verständnis erhalten, indem Sie die Werte von $ p_A $, $ p_B $ und die Anzahl der Daten ändern.

Recommended Posts

Analyse nach Bayes'schem Denken (1) ... Was ist besser, A oder B?
Was ist besser, PyPy oder Python?
Analyse durch Bayes'sches Denken (2) ... Testen Sie den Algorithmus zur Erkennung von Konserven
[Linux] Ende des Prozesses oder Jobs, was ist besser?
Was ist besser, wenn Pythons Standardeingabe input () oder sys.stdin empfängt?
Dies ist schneller zu verwenden, wenn Anweisung oder Wörterbuchtyp beim Konvertieren einer Zeichenfolge (a-> b) in Python
[Anfänger sind besorgt] Was ist besser, Ruby, PHP oder Python?
[Python] gibt A [oder / und] B zurück
Wird der Leerzeichen bei der Verarbeitung der prozentualen Codierung durch ein Pluszeichen oder% 20 ersetzt?
Die Netzwerkanalyse ist eine Weblinkstruktur ①
Die Netzwerkanalyse ist eine Weblinkstruktur ②
Was ist schneller, Python Shuffle oder Sample?