Approach: maintain running sample mean/variance with Welford's algorithm (online, numerically stable) so each new observation updates count, mean, M2 (sum squares). At each look compute the sequential test statistic (e.g., Z = (mean - mu0) / (sigma_hat / sqrt(n)) or t-stat if sigma unknown). Use an alpha-spending function (O'Brien-Fleming or Pocock) to compute cumulative alpha spent α(t) at look t (based on proportion of information/time). Convert α(t) to a critical boundary: for two-sided z-test, z_crit(t) = Phi^{-1}(1 - α(t)/2). Reject at look t if |Z_t| >= z_crit(t). Report adjusted p-value as the smallest cumulative α at which you'd reject (i.e., p_adj = min_t { α(t) : |Z_t| >= z_crit(t)} ) or alternately report sequential p = min_t { two-sided p-value_t / α(t) } capped at 1 (practical reporting: report p_t and compare to α(t)).Example implementation:python
import math
from math import sqrt
from scipy.stats import norm
class RollingSequentialTest:
def __init__(self, mu0=0.0, spending='obrien_fleming'):
self.n = 0
self.mean = 0.0
self.M2 = 0.0
self.mu0 = mu0
self.spending = spending
self.looks = [] # store info fractions
def update(self, x):
# Welford update
self.n += 1
delta = x - self.mean
self.mean += delta / self.n
delta2 = x - self.mean
self.M2 += delta * delta2
def var(self):
return (self.M2 / (self.n - 1)) if self.n > 1 else float('nan')
def info_fraction(self, n_max):
return self.n / n_max
def alpha_spent(self, t):
# t in (0,1] info fraction
if self.spending == 'pocock':
# Pocock approximate: α(t) = α * log(1 + (e-1)*t)
# We'll parameterize later by total alpha when used
return math.log(1 + (math.e - 1) * t) / math.log(math.e)
else: # O'Brien-Fleming-like (very conservative early)
return 1 - norm.cdf(norm.ppf(1 - 0.025) / sqrt(t)) # scaled for two-sided 0.05 baseline
def look(self, n_max, alpha_total=0.05):
if self.n < 2:
return None
s2 = self.var()
se = sqrt(s2 / self.n)
z = (self.mean - self.mu0) / se
t = self.info_fraction(n_max)
alpha_t = alpha_total * self.alpha_spent(t)
z_crit = norm.ppf(1 - alpha_t / 2)
p_two = 2 * (1 - norm.cdf(abs(z)))
rejected = abs(z) >= z_crit
# adjusted p-value: smallest cumulative alpha that would reject; here approximate by p_two
return {'n': self.n, 'mean': self.mean, 'z': z, 'p': p_two, 'alpha_t': alpha_t,
'z_crit': z_crit, 'reject': rejected}
