Создание калькулятора размера тестовой выборки A/B из блога простого последовательного A/B тестирования Эвана Миллера

Question

Чтобы узнать больше о выборе размеров выборки A/B-тестирования, я пытаюсь использовать популярный блог Эвана Миллера для воссоздания калькулятора размеров выборки ( https://www.evanmiller.org/sequential-ab-testing.html). Однако, похоже, есть ошибка, которая не позволяет мне воссоздать размеры выборки, приведенные в статье. Что бы вы предложили мне пересмотреть, чтобы найти решение моей проблемы?

Эта ошибка должна быть в моих расчетах или в моем прочтении проблемы. Решение уравнений ограничения предполагает небольшой размер выборки ~6. В статье предполагается, что обрезание статистики теста, которое показывает, имеет ли курс лечения более высокую степень конверсии, чем контроль, является функцией размера выборки. Затем перечисляются два неравенства, которые нужно решить для одной переменной, N размера выборки. Как воспроизвести расчеты размера выборки?

Большинство из этих проблем исчезло бы, если бы обрезание статистики теста, d_star, не было функцией N.

d_star = z*Sqrt(N), N - размер выборки, а z- нормальное значение z

Однако в таблицах, представленных в последней половине статьи, d_star = z*Sqrt(N) очень точно связывает N и d_star, что предполагает, что d_star изменяется в зависимости от N.

Ограничения с учетом альфа и бета: сумма R(p=1/(1+ дельта)) > 1 - бета сумма R(p=1/2) <альфа

Я добавлю свой код Python 2.7 и график для каждого уравнения ограничения.

#### Begin Python Code to Calculate Sample Size ####
import random
import scipy.stats
import math
import sys
print sys.version

# Functions and helper functions to 
# calculate the sample size.
def calcStatPowerSum(N, script_delta):
    z = scipy.stats.norm.isf(alpha/2.0) #1.96 for alpha = 0.05
    d_star = z*(N**0.5)
    d_star = int(math.ceil(d_star))
    statPowerSum = 0
    for i in range(1, N+1):
        statPowerSum += (float(d_star)/i)*scipy.stats.binom.pmf((i+d_star)//2, i, 1.0-float(1)/(2+script_delta))
        # p and the (1-p) terms are reversed relative to the binomial distribution

    return statPowerSum

def calcCritValueSum(N, script_delta):
    z = scipy.stats.norm.isf(alpha/2.0) #1.96 for alpha = 0.05
    d_star = z*(N**0.5)
    d_star = int(math.ceil(d_star))
    critValueSum = 0
    for i in range(d_star, N+1, 2):
        critValueSum += (float(d_star)/i)*scipy.stats.binom.pmf((i+d_star)//2, i, 0.5)

    return critValueSum

def determineSampleSize(alpha, beta, script_delta):
    z = scipy.stats.norm.isf(alpha/2.0) #1.96 for alpha = 0.05
    d=1
    N=int(math.ceil(z*z))-1
    statPowerSum = 0
    critValueSum = 1
    while (statPowerSum <= 1 - beta or critValueSum >= alpha) and N<3000:
        d+=1
        N=int(math.floor(d*d/z/z))
        statPowerSum = calcStatPowerSum(N,  script_delta)
        critValueSum = calcCritValueSum(N, script_delta)

    return N


alpha = 0.05
beta = 0.8
lift = script_delta = 0.10
sampleSize = determineSampleSize(alpha, beta, script_delta)
print("beta:   ", beta, ", alpha:   ", alpha, ", sampleSize:   ", sampleSize)

##  The article suggests that N=2922 satisfies the constraint equations.
##  But the calculation suggests otherwise.
print("calcCritValueSum: ", calcCritValueSum(2922, script_delta))
print("statPowerSum:     ", calcStatPowerSum(2922, script_delta))

#### End Python Code ####

Следующий код mathematica воспроизводит графики ограничений как функцию d_star. Скопируйте и вставьте код в https://sandbox.open.wolframcloud.com/ чтобы сгенерировать эти графики. Создание более одного сюжета за один раз превысит лимит бесплатного использования.

(*   alpha constraint Plot   *)
z=1.96    (*alpha=0.05*)
Table[Sum[(d/n)*PDF[BinomialDistribution[Floor[d^2/z^2], 1/2], (d+n)/2], {n, d, Floor[d^2/z^2], 2}],{d,1,130}]
ListLinePlot[%, PlotRange->All, AxesLabel->{d_star,prob},PlotLabel->alpha Constraint Plot]


(*   1-beta constraint Plot   *)
delta=0.10
z=1.96    (*alpha=0.05*)
Table[Sum[(d/n)*PDF[BinomialDistribution[Floor[d^2/z^2], (1+delta)/(2+delta)], (d+n)/2], {n, d, Floor[d^2/z^2], 2}],{d,1,130}]
ListLinePlot[%, PlotRange->All, AxesLabel->{d_star,prob},PlotLabel->1- beta Constraint Plot]