Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
110,575 characters · 37 sections · 105 citation commands
Near-Optimal, Near-Calibrated, Non-Parametric Sequential Tests and Confidence Sequences with Possibly Dependent Observations
Inference based on randomized experiments forms the basis of important decisions in an incredibly diverse range of domains, from medicine Schulzc332 to development economics banerjee2015miracle to technology business amatriain2012netflix. Toward the aim of making better and faster decisions, it is particularly helpful for statistical procedures to be flexible grunwald. Experimental designs which require a pre-specified sample size can be quite rigid in practice. For example, unless some oracle information is known about the effect size, such designs will ultimately include some over- or under-experimentation; collecting more samples than necessary when the treatment effect is under-estimated, and not collecting enough when over-estimated. Sequential designs, on the other hand, enable faster decision-making as they support the ability to analyze data as and when it arrives, enabling us to stop experimenting when the data strongly supports a conclusion. The ability to continuously monitor experiments in turn leads to better decisions, as it prevents bad practices that arise from attempting to use more rigid statistical procedures in applications where resources and schedules often change dynamically.
The sequential testing problem can in general be modelled by the following framework, which encompasses most settings from the literature we refer to. The analyst receives a stream of random data points $X_1,X_2,\ldots$, all lying in a certain set $\mathcal{X}$, adapted to a filtration $\mathfrak{F}=(\mathcal{F}_t)_{t \geq 1}$. For any $t \geq 1$, we denote $\psi_t = E[X_t \mid \mathcal{F}_{t-1}]$ the conditional mean of $X_t$ given previous observations. We consider the problem of testing the statistical hypothesis $H_0$ that $\psi_t = \psi'_0$ for every $t \geq 1$, for some $\psi'_0 \in \mathbb{R}$, against the composite alternative $\widetilde{H}_{\backslash 0} : \cup_{\psi'_1 \neq \psi'_0} \widetilde{H}(\psi'_1)$, where $\widetilde{H}(\psi'_1) : \lim_{t \to \infty} \psi_t = \psi'_1$. In words, $\widetilde{H}_{\backslash 0}$ is the hypothesis that $\psi_t$ has a limit, and that this limit is not $\psi'_0$. A sequential test is an $\mathfrak{F}$-adapted stopping time $\tau$. If $\tau = t$ for some $t < \infty$, we say that the sequential time rejects $H_0$ at $t$, while if $\tau = \infty$, we say that the test doesn't reject $H_0$.
Analyzing data in this fashion requires special inference that explicitly accounts for the sequential nature of the decision-making process. It is well known that repeated application of classical significance tests to accumulating sets of data results in procedures with drastically inflated type-I error rates armitage. Even under the null hypothesis, the absolute value of the $t$-statistic applied to an independent and identically distributed (i.i.d.) sequence is guaranteed to fall into the rejection region at some point, regardless of the chosen $\alpha$-level strassen1964invariance. An analyst intent on disproving any one hypothesis can keep collecting data until this occurs.
For these reasons, a long line of literature has explored sequential testing and confidence sequences. However, we claim that the existing literature still does not provide all the guarantees necessary for the widespread use of confidence sequences and sequential tests in statistical practice. We believe indeed that practitioners would want to know three things before adopting them.
In our understanding, the prevailing view is that practitioners have essentially two choices at their disposal. One is to use parametric-model-based confidence sequences which tend to under-cover and over-reject. The other is to use concentration-inequalities-based sequences, which tend to over-cover and under-reject. A prominent example of parametric confidence sequences is the mixture SPRT boundaries for Brownian data sequences introduced by robbins1970boundary and further refined in robbins1970statistical. Notable contributions to the field of concentration-inequality-based sequences include howard2021 and howardquantile.
These two types of sequences are not actually the only options but the alternatives may be less well-known as of now. These alternatives are asymptotic in nature and are based on weak and strong invariance principles (WIPs and SIPs). robbins1970boundary propose a sequence of delayed-start mixture SPRT boundaries and show, via Donsker's WIP, that coverage tends to $\alpha$ as the delay before the start of monitoring diverges to infinity. More recently, bibaut2021sequential use mcleish1974dependent's WIP to provide a sequence of delayed-start confidence sequences for dependent data. smith propose a novel definition of asymptotic confidence sequences and demonstrate how SIPs allow to construct such asymptotic sequences. The advantage of these WIP- and SIP-based sequences (or sequences of sequences) is that they hold the promise of tight type-I error control (as opposed to over or anti-conservativity) under minimal assumptions -- typically mere moment assumptions. We discuss more in detail why going for an asymptotic WIP-based or SIP-based solution is the right strategy in the next subsection.
Nevertheless, we claim that, as of the time of writing of this article, no existing work simultaneously proposes confidence sequences and an analysis showing that these satisfy the aforementioned three criteria. In particular, we discuss in (ref) how the current work compares to smith in that respect.
In fixed-sample-size inference, arguably the most used procedure is to construct approximate central-limit-theorem- (CLT) based confidence intervals. The central limit theorem and extensions thereof only require a second-order moment (or a Lindeberg or Lyapunov condition), and allow for an extremely wide range of data-generating processes. The coverage error one makes by using a CLT-based confidence intervals is controlled by Berry-Esseen bounds and decreases rapidly, scaling as the inverse square root of the sample size (provided the observations have a third moment).
Comparatively, parametric-model-based confidence sequences and nonparametric concentration bounds aren't as commonly used in applied statistical practice for the same reasons we mentioned for confidence sequences: the former become anti-conservative as soon as the parametric assumptions do not hold, while the latter are over-conservative in general.
This therefore motivates looking for a CLT-like solution to sequential inference and testing. Fortunately, there exist sequential extensions of CLT results: these are precisely the WIPs and the SIPs we mentioned above. These guarantee convergence (in distribution and almost surely, respectively) of partial sum processes to Gaussian processes, under weak moment assumptions.
We aim to justify the use of existing sequential tests and the corresponding confidence sequences, namely the delayed-start versions robbins_siegmund1974 of the running maximum likelihood SPRT (rmlSPRT) robbins1972class, robbins_siegmund1974, and the normal mixture SPRT (nmSPRT). In other words, we give guarantees supporting that it is safe, efficient, and advisable to use these. We show that, under mild non-parametric assumptions, sequences of these sequential tests indexed by the burn-in time (1) have type-I error converging asymptotically to the prespecified level, and (2) have expected stopping time converging to the lower bounds known in the i.i.d. case both in the regime where $\psi_1' \to \psi_0'$, and in the regime where $\alpha \to 0$.
We show in (ref) that the type-I error of delayed-start nmSPRTs and rmlSPRTs converges to $\alpha$ as burn-in time goes to infinity. We show in (ref) that, in our non-parametric dependent setting, in the asymptotic regime where $\alpha \to 0$, the expected rejection time is asymptotically equivalent, up to a small constant, to the lower bound known for simple-vs-simple SPRTs under parametric i.i.d. data. We also show that, in the asymptotic regime where the effect size $\psi'_1 - \psi'_0$ converges to 0, the expected stopping time converges, up to a constant factor, to the lower bound known in the parametric i.i.d. setting. We propose in (ref) a sample splitting and variance stabilization method to conduct estimating-equations-based inference under sequentially collected data. In (ref), we study heuristically the optimal choice of the hyperparameter in nmSPRTs, and we propose a heuristic procedure that auto-tunes this hyperparameter. In (ref) we conduct extensive simulation studies so as to demonstrate numerically the validity of our theoretical results and the empirical robustness of our heuristic auto-tuned nmSPRTs, and in (ref) we study a real-data application to A/B testing for quality control of the Netflix client application. In (ref) we review the historical and logical progression from wald's simple-vs-simple SPRT to the results we present in this paper. In (ref) we further discuss the related literature.
We work in the sequential setting introduced in (ref). For any $t \geq 1$, we denote $S_t = \sum_{s=1}^t X_s$. Our methods are designed to work with observations for which the conditional variance given the past converges to one. We will make rigorous what are the specific variance convergence requirements in the upcoming sections, but it might be helpful to think throughout of $\mathrm{Var}(X_t \mid \mathcal{F}_{t-1})$ as being approximately one.
For any $t \geq 1$, and $\lambda > 0$, we define the rmlSPRT and nmSPRT statistics as
We refer the reader to (ref) for the motivation behind the expressions of the test statistics and the rationale of their names.
So as to avert any confusion, the reader should keep in mind that while we present the guarantees for both test statistics along one another, the analyst will opt for one of the two and monitor only the chosen one (as opposed for instance to monitoring which of the two crosses a given threshold first, or any variation thereof). Let $t_0 \geq 1$ be the burn-in period, that is the time the analyst waits before starting to monitor whether the test statistic crosses a certain rejection threshold. We will justify that for burn-in period $t_0$ (and parameter $\lambda$ for the nmSPRT), the “right” rejection thresholds for the rmlSPRT and the nmSPRT are $-\log \widetilde{\alpha}_1(\alpha)$ and $-\log \widetilde{\alpha}_{2,\lambda / t_0}(\alpha)$, respectively, where $\widetilde{\alpha}_1(\alpha)$ solves $h_1(-\log \widetilde{\alpha}_1) = \alpha$ and, for any $\eta > 0$, $\widetilde{\alpha}_{2,\eta}(\alpha)$ solves $h_{2,\eta}(- \log \widetilde{\alpha}_{2,\eta}) = \alpha$, where
The combination of a test statistic, a burn-in period, and a rejection level specifies a sequential test. Formally we define the the $t_0$-burn-in rmlSPRT and nmSPRT as the stopping times
For any $t \geq 1$, we define $\mathcal{C}^{\mathrm{rmlSPRT}}_{\alpha, t_0}(t)$ (resp. $\mathcal{C}^{\mathrm{nmSPRT}}_{\alpha, \lambda, t_0}(t)$) as the set of values $\psi'_0 \in \mathbb{R}$ for which the rmlSPRT (resp. the nmSPRT) doesn't reject at $t$ the hypothesis $H_0: \psi_t = \psi'_0 \ \forall t$. It is then straightforward to observe that the sequences $(\mathcal{C}^{\mathrm{rmlSPRT}}_{\alpha, t_0}(t))_{t \geq 1}$ and $(\mathcal{C}^{\mathrm{nmSPRT}}_{\alpha, \lambda, t_0}(t))_{t \geq 1}$ have coverage equal to 1 minus the level of the corresponding tests. Inverting the expressions of the test statistics yields that, for $t\geq t_0$, $\mathcal{C}^{\mathrm{rmlSPRT}}_{\alpha, t_0}(t) = [\pm c^{\mathrm{rmlSPRT}}_{\alpha, t_0}(t)]$ and $\mathcal{C}^{\mathrm{nmSPRT}}_{\alpha, \lambda, t_0}(t) = [\pm c^{\mathrm{nmSPRT}}_{\alpha, \lambda, t_0}(t)]$ with
We now present two applications that exemplify our setting and the use of the delayed-start confidence sequences.
Characterizing the type-I error of the delayed start rmlSPRT and nmSPRT is a priori easier if the data distribution under the null is known. As invariance principles guarantee convergence of partial sum processes of general nonparametric centered data sequences to Wiener (a.k.a. Brownian) processes, we start our type-I error study by characterizing the distributional properties of the test statistics under Brownian data. Let $W = (W(t))_{t \geq 0}$ be a Wiener process adapted to a filtration $\widetilde{\mathfrak{F}} = (\widetilde{\mathcal{F}}(t))_{t \geq 0}$. Let $\widetilde{Y}_t^{\mathrm{rmlSPRT}}$ and $\widetilde{Y}_{\lambda, t}^{\mathrm{nmSPRT}}$ be the Brownian-data counterparts of $Y_t^{\mathrm{rmlSPRT}}$ and $Y_{\lambda, t}^{\mathrm{nmSPRT}}$:
The following result shows that choices $\widetilde{\alpha}_1(\alpha)$ and $\widetilde{\alpha}_{2, \lambda/t_0}(\alpha)$ above ensure that $-\log \widetilde{\alpha}_1(\alpha)$ and $-\log \widetilde{\alpha}_{2, \lambda / t_0}(\alpha)$ are $(1-\alpha)$-boundaries for the Brownian test statistics $\widetilde{Y}_t^{\mathrm{rmlSPRT}}$ and $\widetilde{Y}_{\lambda, t}^{\mathrm{nmSPRT}}$ when monitored continuously in time at all $t \geq t_0$.
The claims on $\widetilde{Y}^{\mathrm{nmSPRT}}_{\lambda, t}$ and $\mathcal{C}^{\mathrm{nmSPRT}}_{\alpha, \lambda, t_0}$ are a direct consequence of Theorem 2 in robbins1970boundary. The claims on $\widetilde{Y}_t^{\mathrm{rmlSPRT}}$ and $\mathcal{C}^{\mathrm{rmlSPRT}}_{\alpha, t_0}$ follow from similar techniques. Version 7 of smith also includes a proof of these. We include a proof of (ref) in the appendix for self-containedness.
We now show that type-I error converges along any joint sequence of rmlSPRTs or nmSPRTs and of values of $\alpha$, such that the burn-in time $t_0$ diverges to $\infty$. The guarantees will hold even with discrete-time, non-normal, dependent observations. The proof hinges on constructing a joint probability space for both these observations and the continuous-time Brownian motion, in which $\sup_{t \geq t_0} |\widetilde{Y}_t^{\mathrm{rmlSPRT}} - Y_t^{\mathrm{rmlSPRT}}|$ and $\sup_{t \geq t_0} |\widetilde{Y}_{\lambda, t}^{\mathrm{nmSPRT}} - Y_{\lambda,t}^{\mathrm{nmSPRT}}|$ are small. Almost-sure approximations of partial sums by Brownian motions under lax conditions are the object of strong invariance principles. We will specifically leverage the SIP result from theorem 1.3 of strassen1967. To do so, we impose the following conditions.
For now, we take (ref) as a primitive assumption. In (ref), we will discuss simple sufficient conditions for it for the case of stabilized estimating equations. Condition (ref) is akin to (but different from) a martingale Lindeberg condition, and in (ref) we satisfy it by assuming a moment higher than 2 exists.
A SIP guarantee is, of course, asymptotic, so we need to consider regimes where the stopping time is large. One such regime is when we set a long burn-in period $t_0$, essentially to wait for the data to look normal. The following theorem provides a type-I error guarantee (equivalently, a coverage guarantee) along sequences $((t_{0,m}, \alpha_m, \lambda_m))_{m\geq 1}$ such that $t_{0,m} \to \infty$.
Before presenting our representation results for the rmlSPRT and nmSPRT statistics, let us first discuss known related facts. Let $\widetilde{S}(t) = \widetilde{\psi} t + W(t)$ be a time-continuous Brownian observation sequence with constant drift $\widetilde{\psi}$. robbins_siegmund1974 show via It\^o's lemma that, under input data process $\widetilde{S}$, the nmSPRT test statistic $\widetilde{Y}^{\mathrm{nmSPRT}}_{\lambda, t} = \frac{1}{2}\{\widetilde{S}^2(t) / (t+\lambda) - \log ((t+\lambda)/\lambda)\}$ can be represented alternatively as $\widetilde{Y}^{\mathrm{nmSPRT}}_{\lambda, t} = \int_{0}^t \widetilde{\psi}_{\lambda}(s) d\widetilde{S}(s) - \frac{1}{2} \widetilde{\psi}_{\lambda}^2(s) ds$, where $\widetilde{\psi}_{\lambda}(t) = \widetilde{S}(t) / (t + \lambda)$ is the posterior mean of the drift $\widetilde{\psi}$ under prior $\mathcal{N}(0, 1/\lambda)$. We refer to this latter representation as a running-posterior-mean SPRT. In the same article, they discuss the stopping time properties of a related test of which the test statistic is $\int_0^t \widetilde{\psi}_0(s) d\widetilde{S}(s) - \frac{1}{2}\widetilde{\psi}^2_0(s) ds$. Observing that $\widetilde{\psi}_0(t)$ is the maximum likelihood estimate of $\widetilde{\psi}$, we refer to this latter test statistic as a running-maximum-likelihood-estimate SPRT. As robbins_siegmund1974 observe, these running-effect-size-estimate SPRT representations are amenable to rejection time analysis. This motivates us to look for such representations of the rmlSPRT and nmSPRT statistics.
While we cannot directly apply It\^o's lemma in our nonparametric discrete-time setting, we derive a finite-differences equivalent of a certain It\^o-derived identity, namely that for $h_{\lambda}(x,t) = \frac{1}{2} x^2 / (t+\lambda)$, $dh_\lambda(\widetilde{S}(t),t) = \widetilde{\psi}_\lambda(t) d\widetilde{S}(t) - \frac{1}{2} \widetilde{\psi}_\lambda^2(t) + \frac{1}{2} (t+\lambda)^{-1}$. The following lemma is our discrete analog of this identity.
As a direct corollary of (ref) we obtain the following representation results for the rmlSPRT and nmSPRT statistics defined in (ref). Let $X^0_t = X_t - E[X_t | \mathcal{F}_{t-1}]$, $S^0_t = \sum_{s=1}^t X^0_s$, $\breve{X}_t = \psi + X^0_t$, $\breve{S}_t = \sum_{s=1}^t \breve{X}_s$. where $\psi = 0$ under $H_0$ and $\psi= \lim_t \psi_t$ under $H_{\backslash 0}$.
(ref) shows that the test statistics $Y^{\mathrm{rmlSPRT}}_t$ and $Y_{\lambda,t}^{\mathrm{nmSPRT}}$ can be represented, up to some remainder terms, as sums of log-likelihood ratios of which the numerator is the likelihood under an non-anticipating (that is, for the $s$-th term, $\mathcal{F}_{s-1}$-measurable) estimate of the parameter. The non-anticipating martingale terminology was introduced by lorden2005nonanticipating to refer to the device introduced by robbins1970statistical and robbins_siegmund1974 to construct their test statistics.
The log non-anticipating martingales can be further expanded to yield the representations in the following theorem.
Observe that $(M^{(0)}_t)_{t \geq 1}$, $(M^{(1)}_{\lambda, t})_{t \geq 1}$ and $(M^\mathrm{skg}_{\lambda,t})_{t\geq 1}$ are martingales with initial expectation 0. This will facilitate their analysis via the optional stopping theorem in the rejection time analysis. The terms with the superscript “skg” are what we refer to as shrinkage terms. Their presence arises from the presence of the shrinkage parameter $\lambda$ in the running mean estimates $(\psi_{\lambda,s})$, and they converge to zero as $\lambda \to 0$. The $\Delta^{\mathrm{qvar}}_{\lambda, t}$ term is the difference between the quadratic variation of a certain discrete time martingale and that of its time-continuous Brownian approximation.
The terms $(\psi^0_{\lambda, s})$ are the centered errors in the estimates $(\psi_{\lambda, s})$. The term $R^\mathrm{adpt}_{\lambda,t}$ behaves roughly as the sum of the squared errors of the estimates $\psi_{\lambda, s}$. As will appear explicitly from (ref), it contributes positively to the rejection time. We interpret it as the cost of adaptively estimating $\psi$ relative to knowing it a priori.
The presence of the burn-in period in the delayed start running-mean estimate SPRTs adds a slight level of complexity to rejection time analysis as compared to when monitoring starts from the beginning of data collection. Fortunately, one can observe that monitoring the crossing of a fixed threshold $A$ by $Y_{\lambda,t}$ after $t_0$ steps is the same as monitoring the crossing of the offset threshold $A - Y_{\lambda, t_0}$ by the offset test statistic $Y_{\lambda, t} - Y_{\lambda, t_0}$. Therefore, conditional on $\mathcal{F}_{t_0}$, we can analyze the rejection time from $t_0$ as we would do in a situation without burn-in, with the difference that the threshold is now the offset threshold $A - Y_{\lambda, t_0}$. We can thus readily obtain a characterization of the expected rejection time given $\mathcal{F}_{t_0}$ given the $\mathcal{F}_{t_0}$-measurable random threshold $A - Y_{\lambda, t_0}$. In what follows, we introduce the shorthand notation $\tau_1 = \tau^\mathrm{rmlSPRT}(\alpha, t_0)$ and $\tau_2 = \tau^\mathrm{nmSPRT}(\alpha, \lambda, t_0)$.
The following lemma connects the test statistics, the stopping times, and the random rejection thresholds.
So as to obtain bounds on the marginal expected (as opposed to conditional on $\mathcal{F}_{t_0}$) rejection times $E[\tau_1]$ and $E[\tau_2]$, we need to characterize the marginal expectations (that is the expectations w.r.t. the distributions of $Y^\mathrm{rmlSPRT}_{t_0}$ and $Y^\mathrm{nmSPRT}_{\lambda, t_0}$) of the random rejection thresholds. We make the following three assumptions.
The following lemma characterizes the behavior of the expected random thresholds as $\alpha \to 0$.
The decompositions of test statistics from (ref), together with (ref) and (ref) imply the following bounds on the expected stopping times.
From (ref) above, we will be able to obtain bounds on the expected rejection times if we can bound the expected remainder terms and the expected quadratic variation difference $E\left[\Delta^{\mathrm{qvar}}_{\lambda,\tau_i} - \Delta^{\mathrm{qvar}}_{\lambda,t_0} \right]$, $i=1,2$. We study these terms in (ref), (ref), and (ref).
In bounding the adaptivity remainder term $E[R^\mathrm{adpt}_{0,\tau_i} - R^\mathrm{adpt}_{0,t_0}]$, $i=1,2$, we use a technique inspired by the proof of lemma 8 in robbins_siegmund1974. As robbins_siegmund1974, we point out that we should expect $E[\tau_i]$, $i=1,2$, to be at least as large as the expected rejection time of the simple-vs-simple SPRT under the alternative $(\psi_t = \psi,\ \forall t\geq 1)$, that is $-2 \psi^{-2} \log \alpha$. Therefore, since $E[(\psi^0_{0,s})^2] \sim s^{-1}$ under (ref), we expect that
as $\psi \to 0$, $t_0 \to \infty$, $\psi \sqrt{t_0} \to 0$, under fixed $\alpha$. The upcoming lemma shows that this is indeed the correct asymptotic equivalent. The result relies on the following assumptions.
Our bound on the bias term relies on the following two assumptions.
Note that (ref) implies (ref), and that (ref) implies (ref).
We obtain upper bounds on the expected rejection times as a direct consequence of (ref), the bounds on the remainder terms, and the fact that under assumptions (ref)-(ref), the rejection time is almost surely finite ((ref) in the appendix). We consider two asymptotic regimes under which these bounds hold. We define these asymptotic regimes below.
We next consider one simple setting where we can apply delayed-start rmlSPRTs and nmSPRTs, and consider interpretable sufficient conditions for our results to hold.
Suppose we observe an i.i.d. sequence $O_1,O_2,\dots$ and we wish to test whether a parameter $\theta$ of the common distribution of the observations is equal to a certain value $\theta_0$. Suppose that $\eta$ is a (potentially infinite dimensional) nuisance parameter lying in a set $\mathcal{T}$ and that we have an estimating function $D(\cdot ; \eta', \theta')$ for $\theta$ defined over the observation space that satisfies the following robustness assumption.
(ref) fits into this framework with $D(z,\theta_0) = z - \theta_0$ to test $\theta = E[Z_1] = \theta_0$, and (ref) with
to test $\theta = E[U_1\mid A_1=1]-E[U_1\mid A_1=0] = \theta_0$.
We now apply our method to this setting. We use a combination of sequential estimation and sample splitting in order to avoid any metric entropy assumptions. Specifically, we will estimate nuisances on half of all the past data and estimate the variance on the other half. Let $\mathcal{I}_{0,t} = \{ t, t-2, t-4,\ldots \}$ and $\mathcal{I}_{1,t} = \{t-1, t-3, t-5, \ldots \}$ index two complementary sample splits. Fix some $\widehat\eta_t$ sequence adapted to $\mathfrak{F}$ such that $\widehat\eta_t$ is independent of $\mathcal{D}_{0,t} = \{O_s : s \in \mathcal{I}_{0,t}\}$ given $\mathcal{D}_{1,t} = \{ O_s : s \in \mathcal{I}_{1,t} \}$, that is, an “estimate" of $\eta$ based only on the ${\lceil t/2\rceil}$ data points from times $\mathcal{I}_{1,t}$. Set
Then, to test $\theta = \theta_0$, set $X_t=\omega_t(\theta_0) D(O_t ; \widehat{\eta}_{t-1}, \theta_0)$, where $\omega_t(\theta')=(\widehat \sigma_{t-1} (\theta') \vee (\chi^{-1} t^{-\iota}))^{-1}$ with some $\chi>0,\iota\in(0,1)$. Then, under (ref), the hypothesis that $X_1,X_2,\dots$ is an MDS holds if and only if $\theta = \theta_0$ holds.
Note that if we have no nuisances, as in (ref), then we can just use all data up to $t-1$ to compute $\hat\sigma^2_t$, not just the past data having the same parity as $t-1$. We can actually do this as long as $\mathcal T$ is sufficiently simple (e.g., a subset of $\mathbb{R}^d$ for some fixed $d$, rather than, say, a space of nonparametric functions). However, to avoid any such assumptions altogether, we focus here on the case where we split the data by parity. Similarly, we here analyze the case where we clip $\omega_t$ by $O(t^\iota)$ to control for the risk of outlying variance estimates, but this is mostly done to make analysis simple. In practice, we do not recommend this, and in our experiments in (ref), we simply recommend using $\omega_t=1$ whenever the variance estimate is zero and not clipping at any other value. This reduces the need to specify hyperparameters.
As mentioned earlier, a confidence sequence $(\mathcal{C}_{\alpha, t_0}(t))_{t \geq 1}$ for $\theta$ may be obtained by setting $\mathcal{C}_{\alpha, t_0}(t)$ to the set of values $\theta_0$ that the test of $\theta = \theta_0$ doesn't reject at $t$. That is, letting $c_{\alpha, t_0}(t)$ to be either $c_{\alpha, t_0}^{\mathrm{rmlSPRT}}(t)$ or $c_{\alpha, \lambda, t_0}^{\mathrm{nmSPRT}}(t)$, $\mathcal{C}_{\alpha, t_0}(t) = \left\lbrace \theta_0 \in \Theta : |S_t| \leq c_{\alpha, t_0}(t) \right\rbrace$. In the case where $D(o; \eta', \theta')=D_1(o; \eta') + D_2(o; \eta') \theta'$ is linear in a scalar parameter $\theta'$, with $\theta =-E[D_2(O; \eta')]^{-1} E[D_1(O; \eta')]$ being the parameter of interest, the resulting confidence sequence simplifies considerably: let
then our confidence sequence is given by $\mathcal{C}_{\alpha, t_0}(t) = [\widehat{\theta}_t \pm \Gamma_t^{-1} c_{\alpha, t_0}(t)]$. Note we do not actually need to require that $\widehat\eta_t\to\eta$, but we may wish this to be the case so as to obtain a smaller confidence sequence. In particular, in the case that $D_2(o; \eta') = d_2$ for a constant $d_2$, the width of the interval may be decomposed as the product of a factor that doesn't depend on the nuisance, and of $(\Gamma_t / t)^{-1} \to d_2^{-1} \sqrt{\mathrm{Var}(D_1(O_1; \eta_1))}$, with $\eta_1$ the limit of $\widehat{\eta}_t$ in an appropriate sense. In various situations, the latter quantity is minimized at $\eta_1 = \eta_0$. For example, it is the case in (ref) that the width of the confidence intervals is minimized at the true regression function $\eta: (a,l) \mapsto E[U_1|A_1 = a, L_1 = l_1]$.
We now verify our assumptions for this simple setting based on simple sufficient conditions. In the following, let $\sigma^2(\eta', \theta')=E[D(O_1;\eta', \theta')^2]-E[D(O_1;\eta', \theta')]^2$. As we make explicit next, it only takes two relatively mild assumptions in addition to (ref) for Assumptions (ref)-(ref) to hold in our estimating equations setting. In fact, for (ref), and therefore the type-I error guarantee to hold, it only takes the following moment condition and variance lower bound condition on the estimating function.
In what follows, denote $\sigma^2(\eta',\theta') = E[D(O_1;\eta',\theta')^2] - E[D(O_1;\eta',\theta')]^2$.
We can now state the type-I error result.
The expected rejection time results take one more assumption, which we now state.
The condition in (ref) formalizes that $\widehat\eta_t$ admits a limit and characterizes the rate of convergence. Usually, this condition would be obtained from a similar bound on $\|\hat\eta_t-\eta_1\|_{\mathcal T}$ in some norm (e.g., Euclidean norm for a vector of nuisances or $L_p$ for nuisance functions) and establishing (or, assuming) that $\sigma(\cdot,\theta_0)$ is Lipschitz (or, H\"older) continuous in this norm. Such guarantees on $\widehat\eta_t$ can be obtained when it is estimated by maximum likelihood or more generally empirical risk minimization van2000empirical. Generally, if (i) $\eta$ is a vector of pathwise differentiable parameters of the distribution of $O_1$, and (ii) the efficient influence function of $\eta$ admits a moment of order $(2+\delta')$, $\delta' > \delta$, then (ref) will hold. If $\eta$ is the best predictor of some $g_1(O_1)$ as a function of $g_2(O_1)$ in $\mathcal T$, then we can generally obtain guarantees in terms of the rate of the critical radii of $\mathcal T$ wainwright2019high.
We can now state the expected rejection time result.
Here we considered just a simple setting with i.i.d. data and a nuisance-invariant estimating equation in order to show how one would verify our assumptions. Our results do also apply to more intricate settings, but further analysis would be needed to verify the assumptions using simple conditions. One example of a possible extension to the simple setting herein where our results still apply is where the estimating equation $D(\cdot ; \eta', \theta')$ is not completely invariant to $\eta'$, but instead we only have Neyman orthogonality chernozukhov2018double in that $\lim_{\epsilon\to0} \epsilon^{-1} E[D(O_1;\eta+\epsilon(\eta'-\eta), \theta')]=0$ for $\eta'\in\mathcal T$. This is, for example, relevant to sequentially observing data from an observational study, where we do not know the propensity score. For example, we may have a sequential trial involving only the intervention of interest (e.g., experimental drug or surgery) and we have an offline pool of controls to compare to, where we assume selection into our trial at random given observed covariates (observed in the trial and in the pool of controls). Another extension may be to parameters that are path differentiable (that is, an influence function for them exists), but may not necessarily be defined in terms of an estimating equation. Yet another example of a possible extension to the simple setting herein where our results still apply is where the data is not i.i.d. but coming instead from an adaptive experiment such as a contextual bandit. In this case, many of our assumptions could be verified in a very similar manner to how theorems 1 and 2 of bibaut2021post are proven, and guarantees for $\widehat\eta_t$ in the form of (ref) can be obtained from bibaut2021risk. These would all simply be applications of our theory.
As the reader might have noticed, we have so far left out the question of tuning $\lambda$ in the delayed-start nmSPRT. We address this question in the present section. To the best of our knowledge, this question hasn't been treated rigorously even in the case of the standard (that is, non-delayed-start) normal-mixture SPRT under normal observations. Our approach in the current section is to first, in the case of normal observations and without burn-in period, aim to uncover experimentally an identity for the optimal value of $\lambda$ as a function of the effect size $\psi$ and the nominal significance level $\alpha$ ((ref)), and then to heuristically justify it mathematically ((ref)). We then propose a heuristic design for a delayed-start sequential test that auto-tunes $\lambda$ as observations are collected.
For the purpose of experimentally calibrating $\lambda$ as a function of $\alpha$ and the effect size, we consider a Brownian observation sequence with drift $\mu$, that is we consider $S_t = \psi t + W(t)$, $\psi \neq 0$, with $W$ a standard Wiener process. We consider various values of $\psi$ and $\alpha$ and scan through values of $\lambda$. We plot the median, first, and third quartile of the rejection time of the burn-in nmSPRT on the left plot of (ref), and we find and represent the optimal value in the $\lambda$ grid for each couple $(\psi, \alpha)$. We work without a burn-in period, that is we set $t_0 = 1$. For the sake of efficiency comparison, also represent the rejection time of the rmlSPRT (which doesn't depend on $\lambda)$ and the oracle simple-vs-simple SPRT where the numerator corresponds to the (a priori unknown to the analyst) true value of $\psi$.
We plot the results of the grid search for the optimal $\lambda$ for each $(\alpha, \mu)$ on the right plot of (ref). The plot seems to indicate, as might be expected from the prior-on-effect-size intuition, that the choice $\lambda = \mu^{-2}$ is close to optimal, at least for small values of $\alpha$.
Under a Brownian observation sequence $\widetilde{S}(t) = \psi t + W(t)$, the rejection time $\tau_2$ of the nmSPRT (without burn-in, which shouldn't matter as the discussion would be the essentially the same under $t_0$ such as $\psi \sqrt{t_0} \to 0$, as per our asymptotic regimes) of level $\alpha$ and parameter $\lambda$ satisfies
As $\alpha \to 0$, $\tau_2 \to \infty$, and therefore, $\psi \tau_2 + W(\tau_2) \sim \mu \tau_2$. Therefore, for small $\alpha$, $\tau_2$ satisfies approximately
Differentiating the above equation with respect to $\lambda$ and using that, at the optimum $\lambda^*$ we have that $(\partial \tau_2 / \partial \lambda)|_{\lambda = \lambda^*} = 0$, yields that
Injecting the last equality in (ref) yields that
which implies that $\tau_2 / \lambda$ must diverge to infinity as $\alpha \to 0$. Therefore, as $\alpha \to 0$, $(\tau_2 + \lambda) / \tau_2 \to 1$, and thus from (ref), we must have that that $\lambda^* / \psi^{-2} \to 1$.
Given the discussion of the previous two sections, it seems natural to use a test statistic that uses an estimate $\widehat{\lambda}_t$ of $\psi^{-2}$ instead of a prespecified value of $\lambda$. In keeping with the design logic of non-anticipating running-estimate SPRT statistics, we propose the following test statistic:
and $\widehat{\psi}_{\lambda,t} = S_t / (t+\lambda)$, for any $t + \lambda >0$, and $\widehat{\psi}_{0,0} = 0$, as defined earlier. So as to fully specify a test of level $\alpha$, it remains to determine the rejection threshold. We speculate that by enforcing a long enough burn-in period, $Y^{\mathrm{adpt}-\lambda}_t$ will behave similarly to an analogous test statistic
obtained from a sequence $(\widetilde{X}_t)_{t \geq 1}$ of standard normal i.i.d. observations, and therefore, that the rejection threshold for the latter should yield approximately the same type-I error for the former. That is, for a burn-in period of duration $t_0$, we are looking for $-\log \alpha^{\mathrm{adpt}-\lambda}(\alpha, t_0)$ such that
that is
Observing that $(\exp(\widetilde{Y}^{\mathrm{adpt}-\lambda}_t - \widetilde{Y}^{\mathrm{adpt}-\lambda}_{t_0})_{ t \geq t_0}$ is a martingale with initial value 1, we should have, from Ville's inequality case of equality, that the above equation is approximately equivalent, for $t_0$ large enough, to
We propose to solve the above equation by performing a grid search over candidate values of $-\log \alpha^{\mathrm{adpt}-\lambda}(\alpha, t_0)- \widetilde{Y}^{\mathrm{adpt}-\lambda}_{t_0}$ and evaluating the above expectation at the grid points by Monte-Carlo simulations. (Note that the Monte-Carlo simulation entails drawing multiple trajectories of an i.i.d. sequence $(\widetilde{X}_t)_{t \geq t_0}$ of standard normal random variables, as opposed to drawing multiple trajectories of the data sequence, which of course is not possible). Our proposed procedure is then the one that starts monitoring the test statistic $Y^{\mathrm{adpt}-\lambda}_t$ after a burn-in period of length $t_0$ and rejects the null hypothesis as soon as it crosses the rejection threshold $-\log \alpha^{\mathrm{adpt}-\lambda}(\alpha, t_0)$ after $t_0$.
We do not analyze formally this procedure in the current version of this work, but we do evaluate it empirically alongside the previously discussed delayed-start nmSPRT and rmlSPRT in the next section.
We now confirm experimentally the type-I error and expected rejection time guarantees for delayed-start rmlSPRT and nmSPRT confidence sequences, and we compare them to alternative confidence sequences.
We start with type-I error experiments. We consider a sequence of i.i.d. observations $O_1,O_2,\ldots \sim \mathrm{Bernoulli}(0.03) - 0.03$. The choice of the relatively small value $0.03$ is to ensure that the sum $S_t$ of the stabilized martingale difference $X_1,X_2,\ldots$ sequence doesn't converge too fast to a normal. It is of course impossible to evaluate the event $\{\tau < \infty \}$ for any sequential test $\tau$, but, since the martingales we consider must follow the law of the iterated logarithm and the sequences we study have $\sqrt{t \log t}$ asymptotics, rejections should happen early on. We evaluate boundary crossings at the points of a time grid $t_1,t_2,\ldots, t_N$ such that $t_i \approx t_{\max} \sum_{j=1}^i j^\beta / \sum_{j=1}^N j^\beta$, with $N=20000$, $\beta = 2$ and $t_{\max} = 10^8$. This specification ensures that the time grid is denser early on the time axis. We refer the reader to our code for implementation details.
(ref) shows convergence of the type-I error with the burn-in period, while the confidence sequences without burn-in period over reject, likely due to erroneous early rejections at time points when $S_t$ is still far from its Wiener process approximations. Note also that for $\lambda = 100$, we don't observe over-rejections even for small values of $t_0$. This is to be expected as setting $\lambda$ trades off early rejections for later tightness. It can be directly observed from the expression of the nmSPRT confidence sequence that the width of the sequence increases early as $\lambda$ increases. is also to be expected from the interpretation of $\lambda^{-1/2}$ as a prior on the effect size: if the effect size is of order $\lambda^{-1/2}$ a confidence sequence optimally targeting that effect size should be loose for $t \ll \lambda$ and tight around $\lambda$, thereby preventing early rejections if $\lambda$ is large.
We now turn to evaluating features of the distribution of the rejection time. As the reader might have noticed, we haven't discussed in detail so far the choice of the parameter $\lambda$ in the nmSPRT expression, beyond the interpretation of $\lambda^{-1/2}$ as an a priori belief on the magnitude of the effect size. We investigate empirically the effect of $\lambda$ on the rejection time of the nmSPRT in the next subsection.
We now plot several measures of sample efficiency for the burn-in nmSPRT sequences and the burn-in rmlSPRT sequences. In particular, we plot the ratio of the median stopping time of our tests over the median stopping time of the oracle (in the sense that it uses the true value of $\mu$) simple-vs-simple SPRT with same burn-in period. (We compute the adjusted $\alpha$ level for the delayed-start simple-vs-simple SPRT in a similar fashion to that of the delayed-start nmSPRT and rmlSPRT). We refer to this ratio as the “relative efficiency” of the sequential tests.
We observe on the left subplot of (ref) that the relative efficiency of the (delayed-start) rmlSPRT and the nmSPRT seem to converge to 1 as $\alpha \to 0$, as implied by (ref). An empirical (as opposed to predicted by any of the theorems of the current article) finding we infer from the right plot of (ref) is that the median stopping time of the nmSPRT at the optimal $\lambda$ value seems to be within a constant factor $c(\alpha) > 1$ of the median stopping time of the oracle simple-vs-simple SPRT, and that $c(\alpha) \to 1$ as $\alpha \to 0$.
We illustrate the use of the delayed-start rmlSPRT and of the delayed-start nmSPRT on a real data case study. We also compare it to an empirical Bernstein sequential boundary based on a gamma-exponential mixture boundary, as proposed in section 4 of howard2021.
Our data comes from an A/B test used to quality-control the release of a new Netflix client application. Users in the treatment and control groups experienced new and existing versions of the Netflix software respectively. The outcome of interest was the delay between requesting playback and the stream starting, which we term “play delay”. Although this dataset was collected according to a pre-determined sample size, we use it to illustrate the performance of our sequential tests. In the following, pre-treatment measurements of play delay are available for all experimental units, which are leveraged for regression adjustment. We work with log-transformed play delay, and use an undisclosed base to protect proprietary information.
We simulate observation trajectories $(X_t)_{t \geq 1}$ as follows: for each $t$, we draw uniformly at random treatment assignment $A_t \in \{0,1\}$ and then we draw outcome $U_t$ and a pre-treatment covariate $L_t$ by sampling with replacement from actual realized pairs $(L, U)$ in arm $A_t$ in the test data. The data structure is an instance of a Bernoulli trial with covariates as described in the second example at the beginning of the article. The outcome $U_t$ is, as we mentioned, the play delay. We use as covariate $L_t$ the pre-treatment outcome of unit $t$, that is the play-delay before treatment assignment. The probability of assignment to cell 1 is $p=0.5$. Our null hypothesis is that there is no difference in play delay on average between the newer and existing Netflix client, that is $H_0: \theta = 0$, where $\theta = E[U_1 \mid A_1 =1] - E[U_1 \mid A_1 =0]$ identifies this mean causal effect. As in our Bernoulli trial with covariates example, we let
where $\eta_{t-1}(a, l)$ is a least-squares estimator of the linear regression of $U$ on $L$, $L \times A$ and an intercept, computed from observations indexed by $\mathcal{I}_{1,t}$, where the index set $\mathcal{I}_{1,t}$ is as introduced in section (ref). We set $X_t = \widehat{\sigma}_{t-1}^{-1} Z_t$ where $\widehat{\sigma}_{t-1}$ is a $\mathcal{D}_{0,t}$-measurable estimator of the standard deviation of $Z_t$.
We set the burn-in period $t_0$ to 1000 observations and the normal mixture SPRT tuning parameter $\lambda$ to $10^2$. Measurements of play delay are bounded from above by $b$, as longer delays are simply abandoned. When using the empirical Bernstein boundary howard2021, $\eta_{t-1}(a, l)$ is the same linear predictor projected onto the interval $[0, b]$, ensuring $Z_t \in [-2b, 2b]$. In their notation, we use a scale $c = 4b$, and a gamma-exponential mixture with parameter $\rho$. Following the recommendation formulated in their section 3.5, we set the gamma-exponential mixture parameter $\rho$ to the same as the Gaussian mixture parameter $\lambda=100$. We make these choices based on domain knowledge of these types of experiments. They correspond to a target effect size of the order of $1$%, and to the fact that we typically get a hundredfold gain in sample size from pre-allocation outcome adjustment. While in the context of our simulation we of course know the effect size (it turns out the point estimate is -0.81%) since we work from an already fully collected data set, the $1$% order of magnitude is what engineers at Netflix expect using domain specific expertise.
So as to make things more concrete, we plot in figure (ref) the trajectory of the test statistics for one arbitrary random draw of the data sequence. We see that if we had received the data in the particular order of the simulated stream, we would have called the test after collecting slightly after the end of the burn-in period $t_0 = 1000$ while the non-asymptotic empirical Bernstein test stops at almost $10^4$ observations.
We now examine the behavior of the test statistics over many draws of the data sequence. We first examine the empirical Type-I error of the nmSPRT by simulating data under the null hypothesis. Specifically, we simulate data for both treatment and control by sampling with replacement from observed control dataset, which we refer to as a simulated A/A test. Figure (ref), left, shows that the Type-I error approaches the nominal 5% the longer the simulation is performed.
We then compare the behavior of the rmlSPRT, the nmSPRT and the non-asymptotic empirical Bernstein test on simulated trajectories of the A/B test (unlike in the A/A test, here observations for treatment and control units are drawn from their respective empirical distributions obtained from the observed dataset, following simulated treatment assignment $A_t$). We plot simulation results in figure (ref). In this example, the rmlSPRT and the nmSPRT perform relatively similary due to the fact that we chose $\lambda$ in the normal mixture to optimize for rejection times of the order of $10^2$, which is earlier than the end of the burn-in time. As expected, the empirical Bernstein test is much more conservative, and we observe it tends to reject a little less than an order of magnitude later.
In this section, we expose the historical and logical progression from Wald's results wald, waldsequential, waldwolfowitz, waldwolfowitzbayes on the optimality of simple-vs-simple parametric SPRTs to our results, that is the type-I error calibration and expected rejection time optimality of the non-parametric delayed-start running-mean-estimate.
\paragraph{Definition and type-I error.} wald, waldsequential introduced the sequential probability ratio test, defined as follows. Suppose that $X_1,X_2,\ldots$ are i.i.d. drawn from a common distribution with density $p$ w.r.t a certain measure. Consider two simple hypotheses $H_0 : p = p_0$ and $H_1 : p = p_1$, with $p_1 / p_0 < \infty$. The SPRT statistic at $t$ is the likelihood ratio
and the $\alpha$-level SPRT is the stopping time $\tau^{\mathrm{svs}} = \inf \{t \geq 1 : Y^{\mathrm{svs}}_t \geq -\log \alpha \}$. That type-I error is at most $\alpha$ follows, via Ville's inequality ville, from the fact that $(Y_t^\mathrm{svs})_{t \geq 1}$ is a martingale under $H_0$ with initial value 1.
\paragraph{Power, expected rejection time.} wald, waldsequential shows that the SPRT has power 1 under the alternative $H_1$ and that its expected rejection time as $d_{KL}(p_1, p_0) \to 0$ is asymptotically equivalent to $- \log\alpha / d_{KL}(p_, p_0)$, where $d_{KL}$ is the Kullback-Leibler divergence. waldwolfowitz further show that, under the present setting, the expected stopping time of the simple-vs-simple SPRT is optimal among all tests of level $\alpha$ under $H_0$, that is, for any sequential test $\tau'$ such that $\mathrm{Pr}_{H_0}[\tau' < \infty] \leq \alpha$, we must have $E_{H_1}[\tau'] \geq E_{H_1}[\tau^\mathrm{svs}]$.
\paragraph{The case for mixture test statistics under composite alternatives.} Consider a family of densities $p_\psi$ indexed by a one-dimensional parameter $\psi$ and the null hypothesis $H_0 : \psi = \psi_0$. It is often the case that the alternative hypothesis is a composite alternative of the form $H_{\backslash 0} = \psi \neq \psi_0$. In clinical trials or A/B tests for instance, experimenters want to test the absence of average treatment effect (ATE) against the hypothesis that the ATE is non-zero, rather than against a specific non-zero value of the ATE.
It can be shown that for $p \in H_{\backslash 0}$ such that $p$ is closer in Kullback-Leibler divergence to $H_0$ than it is to $H_1$, the simple-vs-simple SPRT of $H_0$ against $H_1$ has power strictly smaller than 1. It can further be shown that, even if $d_{\mathrm{KL}}(p, H_1) < d_{\mathrm{KL}}(p, H_0)$, the expected stopping time can get highly suboptimal as $\psi$ gets away from $\psi_1$.
A remedy to this limitation of the simple-vs-simple SPRT that ensures the resulting sequential test has power 1 against any $\psi_1 \neq \psi_0$, is to use a prior $F$ over the possible values of $\psi_1$. This yields so-called mixture SPRTs, introduced by robbins1970boundary, robbins1970statistical, where $F$ is the so-called mixture distribution, in which the test statistic and the test are defined as
Integrating against $F$ preserves the martingale property under $H_0$ and therefore the type-I error guarantee.
\paragraph{Mixture confidence sequences for normal data.} Under an i.i.d. sequence $X_1,X_2 \ldots \sim \mathcal{N}(\psi, 1)$, that is under $p_\psi(x) = (2 \pi)^{-1/2} \exp(-(x-\psi)^2/2)$, $Y_t^{\mathrm{mixt}}$ takes the form
It is common to specify the null hypothesis by setting $\psi_0 = 0$ (think about testing the absence of average treatment effect against non-zero ATE). An analytically convenient mixture distribution $F$ is the normal distribution robbins1970boundary, robbins1970statistical centered around 0 and with variance $\lambda^{-1}$. The variance of the mixture distribution encodes beliefs about the possible values of the effect size if it isn't zero. This yields
Inverting $Y_t^{\mathrm{mixt}}$ yields a confidence sequence for $S_t$. Specifically, for any $\alpha > 0$, that $\Pr_{H_0}[\tau^{\mathrm{mixt}} < \infty] \leq \alpha$ is equivalent to the fact that
\paragraph{Exact calibration for Wiener processes.} robbins1970boundary further show that for a standard Wiener process, the crossing probability of boundary (ref) is exactly $\alpha$, while it is only known that it is at most $\alpha$ in the case of discrete normal i.i.d. data.
\paragraph{Definition of the running-estimate SPRTs.}
A seemingly different strategy to modify the simple-vs-simple SPRT so as to obtain a test of power one against composite alternatives is to harness the sequentiality of data collection by using likelihood ratios of the form $Y^{\mathrm{rngest}}_t=\prod_{s=1}^t p_{\widehat{\psi}_{s-1}}(X_t)/p_{\psi_0}(X_t)$, where $\widehat{\psi}_{s-1}$ is a running $\mathcal{F}_{s-1}$-measurable estimate of $\psi$. This is the approach proposed by robbins1972class to design one-sided tests of $H_{\leq \psi_0} : \psi \leq \psi_0$ against $H_{> \psi_0} : \psi > \psi_0$. A running estimate used by robbins1972class in the case of a sequence of i.i.d data drawn from $\mathcal{N}(\psi,1)$ is the threhsolded maximum likelihood estimator (MLE) $\widehat{\psi}_t = t^{-1} (\sum_{s=1}^t X_s - \psi_0)_+ + \psi_0$, while another one is the posterior mean computed from $X_1,\ldots, X_t$, under a prior $F$ with support $[\psi_0, \infty)$.
\paragraph{Optimality.} It makes sense that the running-estimate SPRTs should have close to optimal expected rejection time, as $\widetilde{\psi}_{s-1}$ is a proxy for $\psi$, and we know from waldwolfowitz that the level-$\alpha$-under-$H_0$ test with optimal rejection time under $\psi = \psi_1$ is the simple-vs-simple SPRT of $H_0$ against $H_1$. robbins_siegmund1974 study the expected rejection time of one-sided tests of $H_{\leq \psi_0}$ against $H_{> \psi_0}$ of the form $\tau^{\mathrm{rngest}} = \inf \{ t \geq 1: Y^{\mathrm{rngest}}_t \geq -\log \alpha \}$ as $\psi \downarrow \psi_0$. They prove that when $\widehat{\psi}_t$ is taken to be the thresholded MLE, $E[\tau^{\mathrm{rngest}}] \sim P_{H_0}[ \tau^{\mathrm{rngest}} = \infty ](\psi - \psi_0)^{-2} \log (\psi - \psi_0)^{-1}$. This isn't too far off the lower bound $2 P_{H_0}[ \tau^{\mathrm{rngest}} = \infty ](\psi - \psi_0)^{-2} \log \log (\psi - \psi_0)^{-1}$ proven by farrell1964asymptotic for such one-sided tests as $\psi \downarrow \psi_0$.
\paragraph{Connection to mixture SPRTs and confidence sequences.} robbins_siegmund1974 show that the running-estimate SPRTs exhibit an interesting connection to mixture SPRTs in the case where the data stream is a time-continuous process of the form $(\widetilde{S}(t))_{t \geq 0}$, with $\widetilde{S}(t) = \psi t + W(t)$ for all $t \geq 0$, where $(W(t))_{t \geq 0}$ is a standard Wiener process. We expose their observations here. Denote $\widetilde{\mathfrak{F}} = (\widetilde{\mathcal{F}}(t))_{t \geq 0}$ the canonical filtration to which $W$ is adapted. The Brownian continuous-time analog of the class of test statistics of the form $Y_t^{\mathrm{rngest}}$ is the class of test statistics of the form:
where $\widetilde{\psi}(s)$ is an $\widetilde{\mathcal{F}}_s$-measurable estimate of $\psi$. Meanwhile, Brownian continuous-time mixture test statistics take the form
where $F$ is the mixing distribution. It\^o's lemma asserts that for any stochastic process $(X(t))_{t \geq 0}$ of the form $X(t) = X(0) + \int_0^t \mu(s) ds + \sigma(s) dW(s)$, where $\sigma(s)$ and $\mu(s)$ are $\widetilde{\mathcal{F}}(s)$-measurable, and any suitably differentiable $(x,t) \mapsto u(x,t)$, it holds that $u(X_t,t) = u(0,0) + \int_0^t \partial_x u(X_s,s)dX(s)+ (\partial_t u(X_s,s) + \frac{1}{2} \partial_{x,x} u(X_s,s))ds$. Applying It\^o's lemma to $u(\widetilde{S}(t),t) = \log f(\widetilde{S}(t),t)$ yields
Notice that
is the posterior mean of $\psi$ given $\widetilde{\mathcal{F}}(s)$ under prior $F$. For $F= \mathcal{N}(0, \lambda^{-1})$, we have that $\widetilde{\psi}'(s) = \widetilde{S}(s) / (s + \lambda)$, that is the shrunken empirical mean with shrinkage parameter $\lambda$. This shows in particular that the normal mixture boundary (ref) is equivalent to $\widetilde{Y}^{\mathrm{rngest}}_t \leq - \log \alpha$ with $\widetilde{\psi}(s) = \widetilde{\psi}'(s) = \widetilde{S}(s) / (s + \lambda)$. robbins_siegmund1974 show in their Theorem 3 that in the normal discrete case, the half-normal prior yields a running-posterior-mean SPRT with expected rejection time behaving as $2 P_{H_0}[\tau^{\mathrm{rngest}} = \infty] (\psi - \psi_0)^{-2} \log (\psi - \psi_0)^{-1}$.
\paragraph{Boundary associated to the running MLE SPRT} As far as we are aware, there doesn't seem to be a mixture distribution corresponding to the MLE or to the thresholded MLE alluded to earlier. However, application of It\^o's lemma to $h(x,t) = x^2 / (2t)$ yields that the running MLE SPRT log test statistic started at 1 can be rewritten as follows:
Under the null $H_0$, $(Y_t^{\mathrm{rmlSPRT}})_{t \geq 1}$ is a martingale with initial value 1. Therefore, from Ville's inequality,
It can readily be checked that putting $x = -\log \widetilde{\alpha}_1(\alpha)$ sets the above quantity to $\alpha$. Therefore, inverting $Y^{\mathrm{rmlSPRT}}_t$ gives that, with probability $1 -\alpha$,
Theorem 2 in robbins1970boundary asserts (we present here a slight two-sided modification of the result) that for a sequence $X_1, X_2, \ldots$ of i.i.d. random variables with mean 0 and variance 1, and a boundary function $(c(u))_{u \geq u_0}$ such that (i) $c(u) u^{-1/2}$ is non-decreasing for $u$ large enough and (ii) $\int_{u_0}^\infty u^{-3/2} c(u) \exp(-c(u)^2 / 2) du < \infty$, it holds that
This result therefore allows, in nonparametric i.i.d. settings, to obtain approximate confidence sequences from a confidence sequence for the Wiener process, by means of a burn-in period and a time-rescaling. Applying this result to (ref) and (ref) yields that
and
The key enabling result in the proof of theorem 2 in robbins1970boundary is Donsker's weak invariance principle for i.i.d. random variables.
\paragraph{Weak invariance principle for sequential testing and confidence sequences.} Theorem 10 in bibaut2021sequential uses mcleish1974dependent's weak invariance principle for martingale-difference triangular arrays to provide a method for constructing non-parametric asymptotic confidence sequences from a confidence sequence for the Wiener process. We present the result here in the case that $X_1,X_2,\ldots$ is an $\mathfrak{F}$-adapted sequence with $\mathrm{Var}(X_t \mid \mathcal{F}_{t-1}) = 1$, $\psi_t = E[X_t \mid \mathcal{F}_{t-1}]$ and $H_0 : \psi_t = 0\ \forall t \geq 1$, that is $(X_t)$ is a martingale difference sequence.
Let $u_0 \in [0,1]$ and et $(c(u))_{u \in [u_0,1]}$ be a symmetric $(1-\alpha)$-confidence sequence for the Wiener process on $[u_0,1]$, that is $\Pr[\forall u \in [u_0,1], |W(u)| \leq c(u)] \geq 1-\alpha$. As a relatively direct corollary of theorem 3.2 in mcleish1974dependent, it holds that
Here $T$ plays the role of the maximum runtime of the experiment, while $u_0$ is the fraction of $T$ the experimenter uses as burn-in time. The need for $T$ is a theoretical limitation, although it might not be a practical one, as experimenters generally have a time budget or sample size budget to spend on a trial. In the case of a martingale data sequence, we conjecture that concentration-inequality-based methods similar to the ones used to prove theorem 2 from robbins1970boundary could be used to show $\lim_{T \to \infty} \Pr[ \forall t \in \mathbb{N} \cap (T, \infty),\ |S_t| \leq \sqrt{T} c(t/T) ] = 0$, and therefore get rid of the need for a maximum experiment runtime. However, doing so might be more complex in the originally intended martingale-difference array setting.
\paragraph{Strong invariance principle for asymptotic time-uniform confidence sequences.} smith introduce the use of strong invariance principles (also known as almost sure invariance principles or strong approximation results) to construct confidence sequences in nonparametric settings. One key contribution of their work is to introduce a definition of asymptotic confidence sequence (AsympCS). They say that a sequence $(c_t)_{t \in \mathbb{N}}$ is a symmetric $(1-\alpha)$-asymptotic confidence sequence for the partial sum process $(S_t)_{t \in \mathbb{N}}$ if there exists a $(1-\alpha)$-exact confidence sequence $(c_t^*)_{t \in \mathbb{N}}$ for $(S_t)_{t \in \mathbb{N}}$ such that $c_t / c^*_t \rightarrow 1$ almost surely.
Successive versions of this article use different strong approximation results. Starting from version 5, they have been using strassen1967's strong invariance principle for martingales and allows for asymptotic confidence sequences for the general setting where $(S_t)_{t \in \mathbb{N}}$ is a martingale. Starting from version 7, they also include type I guarantees for sequences of delayed start confidence sequences. Specifically, they replicate with their techniques (and lift some assumptions for) the type-I error guarantee for sequences of delayed-start normal-mixture sequences, initially proven under martingale data in version 1 of the current paper. They also extend the guarantee (ref) for sequences of the delayed-start running MLE SPRT sequence to martingale data.
In our view, the main focus of smith is proposing a novel definition of AsympCS that is asymptotically close in parameter space to an exact $1-\alpha$ confidence sequence for the parameter of interest. Meanwhile, our work focuses on the testing properties, that is, type-I error and expected rejection time, of sequences of confidence sequences or of tests. In particular, the current version of their paper, version 7 at the time of writing of this work, does not provide a rejection time analysis. One difference between our setting and theirs is that we impose the normalization condition of the conditional variance, which they don't. We use this condition heavily in our rejection time analysis. We leave to future work the discussion of whether this condition is actually necessary for a rejection time analysis.
Classical approaches to hypothesis testing have predominantly dealt with experiments of fixed, predetermined sample sizes, which we refer to as the fixed-n kind. The emphasis on fixed-n tests by early pioneers such as Fisher is presumably a consequence of the motivating applications that drove the development of hypothesis testing procedures in the first half 19th century, in which outcomes of an experiment were only available long after the experiment had been designed, such as in agricultural research armitage2. As tests could only be performed once, fixed-n tests were designed to maximize power subject to a type-I error constraint neymanpearson. Increasingly in modern experiments, however, observations from experimental units become available sequentially instead of simultaneously, providing many opportunities to perform a test instead of just one. The application of fixed-n tests to sequential designs is made difficult because it requires making an undesirable trade-off balancing the competing objectives of detecting large effects early and detecting small effects eventually. Performing the test later risks exposing many experimental units to a potentially large and harmful treatment effect, while performing the test early risks a high type-II error for small effects. These desires have led to bad statistical practices whereby fixed-n procedures are naively applied to accumulating sets of data, see peeking for a discussion pertaining to online A/B tests, which sacrifice type-I error guarantees armitage, permitting the analyst to incorrectly sample to a foregone conclusion anscombe.
For modern sequential designs, sampling until a hypothesis is proven or disproven appears to be a very natural form of scientific inquiry, which requires testing procedures to preserve their type-I/II error guarantees under continuous monitoring. Sequential inference is fundamentally tied to the theory of martingales ramdaskoolen. A test martingale is a statistic that is a nonnegative supermartingale under the null hypothesis. Ville's inequality ville is then used to bound the supremum of the process to provide a time-uniform type-I error guarantee. Research into sequential analysis in the statistics literature began with the introduction of the sequential probability ratio test (SPRT) wald, waldsequential. Although Wald did not reference martingale theory in the exposition of the SPRT, the connection is clear in hindsight by observing that the likelihood ratio is a nonnegative supermartingale under the null. The simple-vs-simple SPRT enjoys the optimality property of being the sequential test that minimizes the average sample number (expected stopping time) among all sequential tests with no larger type-I/II error probabilities waldwolfowitz. This is extended to the continuous-time version in dkw.
The SPRT for simple-vs-simple testing problems and the mixture SPRT (mSPRT) for composite testing problems can be interpreted as Bayes factors jeffreys, kass, forming a bridge between Bayesian, frequentist, and conditional frequentist approaches to sequential testing conditional_frequentist_simple, conditional_frequentist_nested. The SPRT also appears in Bayesian decision-theoretic approaches to sequential hypothesis testing in which there is a constant cost per observation waldwolfowitzbayes, bergerdecision. However, care must be taken when specifying priors in composite testing problems, should one seek to have strict frequentist guarantees deHeide2021. Composite tests in statistical models with group invariances can often be reduced to simple hypothesis tests by constructing invariant SPRTs invariantsprt based on a maximally invariant test statistic lehmann2005testing, LehmCase98. These invariant SPRT test statistics can be obtained as Bayes factors by using the appropriate right-Haar priors on nuisance parameters in group invariant models hendriksen. Such arguments were used by robbins1970statistical to develop sequential tests for location-scale families with unknown scale parameters.
Confidence sequences darling67 can be obtained by inverting a sequential test, and sequential $p$-values can be obtained by tracking the reciprocal of the supremum of the test martingale. Together, these generalize the coverage and type-I guarantees held by fixed-n confidence intervals and $p$-values to hold uniformly through time. Procedures with these guarantees are appropriately referred to as “anytime valid." Relationships between test-martingales, sequential $p$-values and Bayes factors are discussed in shafer. Nonparametric confidence sequences under sub-Gaussian and Bernstein conditions are provided in howard2021. These results are nonasymptotic, yielding valid confidence sequences for all times, but may be conservative. Confidence sequences for quantiles and anytime-valid Kolmogorov-Smirnov tests are provided in howardquantile. smith obtain asymptotic confidence sequences, in the sense that the intervals converge almost surely to a valid confidence sequence with an error that is orders smaller than the width of the latter. They achieve this by approximating the sample average process by a Gaussian process using strong invariance principles strassen1964invariance, strassen1967, Komlos1975,Komlos1976, like us. Their focus is on having approximate confidence sequence width, which need not translate to type-I error guarantees. In particular, there is a risk of rejecting too early when the cumulative sum does not look normal yet. Moreover, they only guarantee that a similar-width confidence sequence has at-least-$\alpha$ coverage, but do not characterize its power, only that the width has the right rate dependence on $t$. Therefore, at the same time, if we do wait, the confidence sequences can be overly conservative.
For certain continuous-time martingales, Ville's inequality is an equality robbins1970boundary. For discretely observed martingales, however, Ville's inequality is generally strict, meaning the type-I-error guarantees it yields for test martingale are conservative. The conservativeness follows from the amount by which the stopped sum process exceeds the rejection boundary (zero in the continuous case) and is often referred to as the “overshoot” problem with the SPRT siegmundbook. Understanding the size of the overshoot is key to understanding how conservative existing bounds are on type-I error and expected stopping times. wald's approximation to the type-I error is obtained by simply ignoring the overshoot. siegmund75 obtains an approximation to the type-I error for the simple-vs-simple SPRT in exponential-family models as complete asymptotic expansions in powers of $\alpha^{-1}$ with exponentially small remainder as $\alpha \rightarrow 0$. With mSPRTs the rejection boundary is curved, and studying the distribution of the overshoot is often tackled via nonlinear renewal theory woodroofe67, woodroofebook, zhang88. As $\alpha\rightarrow 0$, laisiegmund1977, laisiegmund1979 derive asymptotic approximations to the expected value and distribution function of the nmSPRT stopping time under the null so as to study the type-I error resulting from truncated nmSPRT tests. Similar results for the expected stopping times can be found in hagwoodwoodrofe. To our knowledge existing work has focused on asymptotic ($\alpha\rightarrow 0$) approximations to moments of stopping times for parametric SPRTs which yield sharper results than wald's when neglecting the overshoot. While previous authors also use these tools to obtain type-I errors for truncated sequential tests, no attention has been given to calibrating the type-I error for open-ended sequential tests.
As trends in online experimentation shift toward streaming approaches, sequential approaches to A/B testing have seen increased adoption johari2022always, lindon, lindon20. In other applications, particularly in medicine, it may not be possible to test after every new observation. In clinical trials, a small number of interim analyses may be planned, which does not warrant a fully sequential test. Instead, group sequential tests pocock,obrien,demets,jennison1999group can be performed which provide a calibrated sequential test over a fixed and finite number of analyses. Analogous to confidence sequences, repeated confidence intervals provide strict coverage uniformly across all interim analyses rci, jennison84. These procedures are useful when testing on a certain cadence, such as daily, suffices and when a terminal endpoint of the experiment is known. They are, however, not as flexible as fully sequential procedures as they do not allow the experiment to continue past the final analysis, having fully spent their $\alpha$-budget.
Test martingales are closely related to e-processes. An e-variable is a random variable (or statistic) that has expectation at most 1 under the null hypothesis grunwald. An e-process is a nonnegative process, upper bounded by a nonnegative supermartingale, such that the stopped process is an e-variable under any stopping rule ruf22, although it itself may not be a nonnegative supermartingale ramdasexchange. Thanks to this property it is possible to build sequential tests from e-processes. gametheoryav provide a review of test martingales, e-processes, anytime valid inference and game theoretic probability and its applications to sequential testing. See, for example, log-rank tests schure, contingency tables schure and changepoint detection shin. See also the running-MLE sequential likelihood ratio test of wasserman2020universal.