NEJSDS The New England Journal of Statistics in Data Science 2693-71662693-7166 New England Statistical Society NEJSDS51 10.51387/23-NEJSDS51 Methodology Article Statistical Methodology E-detectors: A Nonparametric Framework for Sequential Change Detection ShinJaehyeokshinjaehyeok@google.com∗ RamdasAadityaaramdas@cmu.edu RinaldoAlessandroarinaldo@cmu.edu Google, USA. E-mail address: shinjaehyeok@google.com Carnegie Mellon University, USA. E-mail address: aramdas@cmu.edu Carnegie Mellon University, USA. E-mail address: arinaldo@cmu.edu Corresponding author. 202411220232222926012102023 © 2024 New England Statistical Society2024 Open access article under the CC BY license.

Sequential change detection is a classical problem with a variety of applications. However, the majority of prior work has been parametric, for example, focusing on exponential families. We develop a fundamentally new and general framework for sequential change detection when the pre- and post-change distributions are nonparametrically specified (and thus composite). Our procedures come with clean, nonasymptotic bounds on the average run length (frequency of false alarms). In certain nonparametric cases (like sub-Gaussian or sub-exponential), we also provide near-optimal bounds on the detection delay following a changepoint. The primary technical tool that we introduce is called an e-detector, which is composed of sums of e-processes—a fundamental generalization of nonnegative supermartingales—that are started at consecutive times. We first introduce simple Shiryaev-Roberts and CUSUM-style e-detectors, and then show how to design their mixtures in order to achieve both statistical and computational efficiency. Our e-detector framework can be instantiated to recover classical likelihood-based procedures for parametric problems, as well as yielding the first change detection method for many nonparametric problems. As a running example, we tackle the problem of detecting changes in the mean of a bounded random variable without i.i.d. assumptions, with an application to tracking the performance of a basketball team over multiple seasons.

 Sequential change detection Nonparametric Shiryaev-Roberts CUSUM NSFIIS-2229881DMS-2310718NSFDMS-EPSRC 2015489A. Ramdas was partially supported by NSF grants IIS-2229881 and DMS-2310718. J. Shin and A. Rinaldo were partially supported by NSF grant DMS-EPSRC 2015489. Introduction

Suppose we observe sequentially a stream of random variables X 1 , X 2 , … , whose marginal distributions may change at some unknown time, or changepoint, ν. To take one concrete example that we generalize later, denote the data stream by X 1 , … , X v ∼ P μ 0 and X v + 1 , X v + 2 , … ∼ P μ where P μ 0 and P μ are some probability distributions with parameters μ 0 , μ ∈ R, respectively. Let P ν and E ν denote probability and expectation, respectively, with respect to the distribution of the entire infinite data stream, when the change occurs at time ν. If there is no change, we think of ν as being equal to ∞, and we let P ∞ and E ∞ refer to the corresponding probability and expectation.

We are concerned with designing sequential changepoint detection procedures to determine, at each time, whether a changepoint has occurred in the near past that (i) provide non-asymptotic false alarm guarantees, (ii) allow for non-parametric classes of pre- and post-change distributions and (iii) are computationally efficient. Formally, a sequential changepoint detection algorithm consists of a data-dependent stopping rule N ∗ ≥ 1. If not specified explicitly, the underlying filtration with respect to which N ∗ is defined is assumed to be the natural filtration generated by the data stream X 1 , X 2 , … , but in some cases it is beneficial to coarsen the filtration. If the stopping time N ∗ is finite, then we declare that a changepoint has been detected, in the sense that sufficient evidence has accumulated to support the hypothesis that the data generating distribution has changed. If the algorithm never stops, we set N ∗ = ∞ and no changepoint is proclaimed. To evaluate the performances of detection algorithms, we can check how quickly the algorithm can detect a change in the distribution while controlling the frequency of false alarms.

To control false alarms, we adopt the average run length (ARL) metric [23], defined by ARL : = E ∞ N ∗ . We say that the “ARL is controlled at level α” if E ∞ N ∗ ≥ 1 / α. An equivalent error metric is the False Alarm Rate (FAR), which is the reciprocal of the ARL, and we would like to ensure that the FAR is at most α, and we call a sequential change detection procedure as “valid” if it satisfies the above constraint.

A widely-used measure of the speed of detection after a changepoint is the worst average delays conditioned on the least favorable observations before the change [18], or conditioned on the event that the algorithm stops after the changepoint [25]. These are defined by J L ( N ∗ ) : = sup ν ≥ 0 esssup E ν [ [ N ∗ − ν ] + ∣ F ν ] , J P ( N ∗ ) : = sup ν ≥ 0 E ν [ N ∗ − ν ∣ N ∗ > ν ] , where the subscripts indicate the authors, and it is known that J P ( N ∗ ) ≤ J L ( N ∗ ) [22].

Our implicit objective is to (approximately) minimize J L ( N ∗ ) or J P ( N ∗ ) while guaranteeing that the ARL is controlled at a prespecified level α ∈ ( 0 , 1 ). We differ from other work in that our focus is on nonparametric and composite pre- and post-change distributions, as well as on deriving nonasymptotic guarantees. In several such settings, it is a priori unclear how to define any valid sequential change detection algorithm, let alone an optimal one. Accordingly, we first address the design problem, and then move to questions involving (approximate) optimality for detection delays.

Prior Work and Our Contributions

If both pre- and post-change parameters μ 0 , μ are known and the distributions have densities p μ 0 and p μ with respect to some common reference measure, then the CUSUM procedure [23] has been known to achieve the optimal worst average delay (exactly for J L ( N ∗ ) and asymptotically for J P ( N ∗ ) as α → 0) among all procedures controlling ARL at the same level [18, 21, 31, 14]. Recall that the CUSUM procedure is defined by the stopping time N CU ∗ : = inf { n ≥ 1 : M n CU ≥ c α CU }, where c α CU > 0 is a constant chosen so that E ∞ ( N CU ∗ ) = 1 / α, and the test statistic M n CU is defined by the recursive formula M n CU = p μ ( X n ) p μ 0 ( X n ) · max { M n − 1 CU , 1 } , M 0 CU : = 0 . The Shiryaev–Roberts (SR) procedure [37, 32] is defined by the stopping time N SR ∗ : = inf { n ≥ 1 : M n SR ≥ c α SR }, where c α SR > 0 is a constant chosen so that E ∞ ( N SR ∗ ) = 1 / α, and the test statistic M n SR is obtained recursively as M n SR = p μ ( X n ) p μ 0 ( X n ) · [ M n − 1 SR + 1 ] , M 0 SR : = 0 . Unlike CUSUM, the SR procedure does not achieve exact minimax optimality for J L ( N ∗ ). However, the SR procedure and its generalized versions enjoy strong asymptotic optimality guarantees [26, 27, 41, 39].

The literature sometimes assumes that the pre-change distribution is known or can be approximated with a high precision by using the previous history. However, the post-change distribution is typically unknown and is assumed to belong to a family of distributions P : = { p μ : μ ∈ Θ }. In this case, one natural approach would be to replace the unknown parameter μ with an estimator μ ˆ. If we use the maximum likelihood estimator (MLE), then we obtain the CUSUM procedure based on the generalized likelihood ratio (GLR) rule [e.g., see 1, 49, 38]. For those not familiar with sequential change detection, [15] provides a good overview.

Limitations of Prior Work Though the usage of GLR statistic for the sequential change detection problem has a long history and often yields good empirical performance, the current literature has two main limitations.

First, most existing methodology relies on parametric assumptions on the family of distributions (e.g.: exponential families), both pre-change and post-change. There have been attempts to move away from this setting, and we will discuss these later. However, a general framework for deriving sequential change detection procedures in general nonparametric or composite settings has not been previously presented in the generality that we do here.

Second, the study of statistical properties has typically focused on the asymptotic regime of α → 0, unless the GLR statistic is defined on a well-separated post-change parameter space. In this paper, we guarantee nonasymptotic control on the ARL at a prespecified level (such as α = 0.001). In fact, in many settings considered, we do not know of any existing method to guarantee (even asymptotic) ARL control.

(We think that in theory and practice, the first is a bigger issue than the second. Luckily, our solution for the first automatically handles the second. Indeed, in the composite and nonparametric settings we consider, it is quite unclear how to control the ARL in any asymptotic sense.)

Finally, a direct online implementation of the GLR rule is infeasible since the memory and computation time both increase at least linearly as n → ∞. One natural approach to tackle this online implementation issue is to use window-limited versions of the GLR rule [49, 13]. For instance, a simple form of the window-limited GLR rule can be defined by, at each time n, computing μ ˆ over only times n − W to n for a properly chosen window size W > 0. However, the optimal choice of window size W has been studied only in the asymptotic setting ( α → 0). For a fixed α, the optimal window size depends on the difference between the pre- and post-change distributions, which is unknown.

Our Contributions We present a general framework for sequential change detection, focusing on (parametric and nonparametric) settings with composite pre- and post-change distributions, and nonasymptotic guarantees on the ARL:

We introduce the concept of an e-detector that underlies our construction of sequential change detection procedures. The e-detector utilizes a generalization of the underlying martingale structure of likelihood ratios in classical sequential change detection procedures. This e-detector framework is applicable in nonparametric settings including sub-Gaussian, sub-exponential, and bounded random variables, among many others. In such settings, there is no common reference measure and likelihood ratios cannot be directly defined, thus composite nonnegative supermartingales, or more generally “e-processes”, must be employed in their place.

Despite handling composite pre- and post-change distributions, even without an iid (independent and identically distributed) assumption on the data, our e-CUSUM and e-SR sequential change detection procedures based on e-detectors can always nonasymptotically control the ARL at level α.

Nonasymptotic bounds on the worst average delay are derived in special cases for nonparametric distributions with exponential tail decay (such as sub-Gaussian or sub-exponential), and they match the rate of known lower bounds for exponential families as α → 0.

Computationally feasible algorithms are presented, in order to run our procedures in an online fashion without windowing. Practical strategies to choose hyperparameters are discussed. These are based on an adaptive mixture method, with the number of mixture components growing slowly over time.

Our procedures have natural gambling interpretations, and our work can be viewed as setting the foundations for a game-theoretic approach to sequential change detection. Before discussing the general framework in detail, in the following subsection, we present a motivating real-world example involving bounded random variables to illustrate how our nonparametric framework can be easily used in settings in which it is nontrivial to apply other common methods.

 Example: A Changepoint in Cleveland Cavaliers 2011–2018

The Cleveland Cavaliers are an American professional basketball team. We use the Cavaliers’ game point records over 2010–11 to 2017–18 NBA seasons to illustrate how our proposed nonparametric sequential change detection algorithm can be applied to detect an interesting changepoint in the Cavaliers’ recent history.2

The R code to reproduce all the plots and simulation results of the paper is available at https://github.com/shinjaehyeok/e_detector_paper

The left plot in Figure 1 shows the difference between the scores of the Cavaliers and those of their opposing teams (also known as Plus-Minus) in all the games from the 2010–11 to the 2017–18 regular seasons. Each red line refers to the yearly average difference score in each season. Roughly, this value shows how well the Cavaliers performed against their opponent in terms of scoring. Typically, if a seasonal average is positive (or negative) then we may say that the Cavaliers showed a strong (or poor) performance in the corresponding season. The right plot shows a changepoint detected in early 2015; NBA fans may recall one major cause of the sharp improvement — LeBron James returned to the Cavaliers in 2014. However, how can we detect such a change on the fly by only tracking the Plus-Minus for each game?

This type of question fits well into the sequential change detection framework. Let X 1 , X 2 , … be the sequence of Plus-Minus stats we observe sequentially. After observing a poor performance of the Cavaliers in 2010–11 season, we define the Cavaliers’ pre-change distribution on the Plus-Minus stats as follows: the average Plus-Minus of the team is less than or equal to μ 0 : = − 1. Now, we are interested in detecting a meaningful performance improvement on the fly by defining the post-change distribution as follows: the average Plus-Minus of the team is greater than μ 1 : = 1. Here, the gap | μ 1 − μ 0 | between averages of Plus-Minus in pre- and post-changes refers to the degree of improvement we consider as a significant one.

Left: Plus-Minus of the Cavaliers from 2010–11 to 2017–18 seasons. Each horizontal red line corresponds to the seasonal average. Right: The sample path of (the logarithm of) one of our e-detectors. The horizontal red line is the threshold equal to log ( 1 / α ) controlling the ARL by 1 / α = 1000. In this example, the procedure detects the changepoint in the Plus-Minus of the Cavaliers at the end of 2014–15 season.

Although the formulation of the problem is simple as described above, it is still nontrivial to fit this problem into commonly used sequential change detection procedures for the following reasons. First, it is not easy to choose a proper parametric model to fit observed Plus-Minus stats since they are integer-valued samples with varying mean and variance over seasons as Figure 1 illustrates. Second, even if we can choose a proper model, it is difficult to find a threshold to detect the changepoint since we are interested in detecting any changes larger than | μ 1 − μ 0 | instead of a fixed post-change. Many common methods have been relying on a high-quality simulator or large enough sample history for pre-change observations to compute a valid threshold. For this example, however, it is hard to access such tools since the Cavaliers’ overall Plus-Minus stats are difficult to model directly, and it is tricky to justify using existing records to get a valid threshold as the team’s overall performance varies a lot around the 2010–11 season.

Based on the introduced framework of the sequential change detection procedure using e-detectors, we detour difficulties in the commonly used methods illustrated above as follows. First, due to the nonparametric nature of the new framework, we do not need to choose any parametric model to fit the data. Instead, we simply assume the absolute value of each Plus-Minus stat is bounded by a large number — in this example, we set 80 as the boundary. Though we set a conservatively large boundary, our detection procedure is variance-adaptive so that we can detect the changepoint efficiently without specifying the variance of observations. Second, the nonasymptotic analysis of the new framework makes it possible to choose an explicit detection boundary, which is equal to log ( 1 / α ), to build a sequential change detection procedure controlling the ARL by 1 / α for any given α ∈ ( 0 , 1 ). In this example, we choose α = 10 − 3 to make ARL is larger than at least 10 regular seasons of games. The right plot in Figure 1 shows the log of e-detectors on which we build the sequential change detection procedure. The red horizontal line corresponds to the detection boundary given by log ( 1 / α ). We can check the procedure detects the changepoint of the Plus-Minus stat of the Cavaliers in the middle of the 2014–15 season. See Section 5.2 for the detailed explanation about how we construct the sequential change detection procedure based on the general methodology we introduce in the paper.

Paper Outline The rest of the paper is organized as follows. In Section 2, we introduce a general framework about how to build composite, nonparametric, and nonasymptotic sequential change detection procedures using e-detectors. Section 3 extends the previous framework to the case where we have a set of e-detectors and explain how to use a mixture method to combine multiple e-detectors effectively. In Section 4, we introduce an exponential structure of e-detectors that makes it possible to design a near-optimal detection procedure with an explicit upper bound on worst average delays. Based on the proposed framework, Section 5 presents two canonical examples of Bernoulli (parametric) and bounded random variables (nonparametric) cases with real data applications of the Cavaliers 2011–2018 statistics. We conclude with a discussion, and defer proofs to the supplement.

 Nonparametric Sequential Change Detection Using E-detectors Problem Setup

Let P denote the set of possible pre-change distributions, which could in general be a nonparametric class. We do not assume the observations in the sequence to be independent or identically distributed: each P ∈ P is a distribution over an infinite sequence of random variables.

We will assume throughout that up to the unknown changepoint ν, the observations X 1 , … , X ν follow a distribution P ∈ P. The remaining observations X ν + 1 , X ν + 2 , … are drawn from a distribution Q in a class of post-change distributions Q. In this case, we let P P , ν , Q , E P , ν , Q and V P , ν , Q denote probability, expectation and variance operators over the entire data stream.

If there never is a change, we will use the notation P P , ∞ , E P , ∞ and V P , ∞ . Also, if a change occurs at the beginning ( ν = 0 ) then we use P 0 , Q , E 0 , Q and V 0 , Q . Note that technically P 0 , Q = P Q , ∞ , but we use the former to denote that Q is a post-change distribution and a changepoint has occurred at the very start, and the latter to denote that Q is a pre-change distribution and a changepoint never occurs.

Let F : = { F n } n ≥ 0 be a filtration where we let F 0 : = { ∅ , Ω } for simplicity. Let M : = { M n } n ≥ 0 be a nonnegative adapted process with respect to the filtration F. If required, we define M ∞ : = lim sup n → ∞ M n and F ∞ : = σ ( ⋃ n ≥ 0 F n ) [see, for example, 4]. It is common to consider the natural filtration F n : = σ ( X 1 , … , X n ), but there are situations where restricting the filtration could be advantageous (for example, when there are nuisance parameters). There are yet other situations when enlarging the filtration with external randomness can be useful.

Let T denote the set of all stopping times with respect to F, but we will later see that it typically suffices to consider finite stopping times, or those with finite expectation.

In our paper, the changepoint ν and the post-change observations also do not need to be independent of the pre-change observations, and they do not need to be identically distributed. In other words, ν could be a stopping time, and Q could itself depend on the pre-change data. It may be helpful to imagine an adversary who adaptively decides at each step whether or not to “stop” the pre-change data. If they choose to stop at time ν, there are two options for post-change points: (a) then can pick a distribution Q, draw a sequence Y 1 , Y 2 , … , from Q, and reveal X ν + i = Y i sequentially, or (b) the pick a distribution Q and draw X ν + 1 , X ν + 2 , … from Q ∣ F ν . In situations without a common reference measure, setting (b) may be tricky to formally define (especially if the pre-change data have probability zero under Q), so one may think of setting (a). All of our results on ARL do not require any assumptions on the changepoint ν or post-change distribution of the data. When analyzing detection delay, we will typically assume that X ν + 1 , X ν + 2 , … are independent of the pre-change data, and are drawn from Q as if the time ν was reset to zero; we believe this can be relaxed in future work. But for much of the paper, it is also okay to assume that the post-change data is drawn a data-dependent Q (or drawn from Q ∣ F ν ).

With the appropriate definitions and setup in place, we can now define our central concept, an e-detector.

 What Is an E-detector? (<inline-formula id="j_nejsds51_ineq_087"><alternatives><mml:math> <mml:mi mathvariant="script">P</mml:mi></mml:math><tex-math><![CDATA[$\mathcal{P}$]]></tex-math></alternatives></inline-formula>-e-detector).

The process M is called an e-detector with respect to the class of pre-change distributions P if it satisfies the property E P , ∞ [ M τ ] ≤ E P , ∞ [ τ ] , ∀ τ ∈ T , ∀ P ∈ P .

For brevity, we refer to M as “an e-detector for P”, or if P is understood from context, then simply “an e-detector”. If a stopping time τ has a nonzero probability of being infinite, then inequality (2.1) is trivially satisfied. Thus, the condition is only really required to hold for stopping times with finite expectation under some P. This latter set of stopping times depends on P, and so in order to not complicate notation, we continue to simply consider all stopping times  T. (Also note that the condition can only be satisfied if process M is integrable under any P ∈ P, so this is implicitly assumed to be the case in what follows.)

By the linearity of expectation and Tonelli’s theorem, an average (or “a mixture”) of e-detectors is also an e-detector. More formally, if { M a } a ∈ A is a set of e-detectors (where a is a tuning parameter), then so is ∫ M a d μ ( a ) for any fixed probability distribution μ over A. Later in this paper, we will in particular use finite mixtures of the form ( M 1 + M 2 + ⋯ + M K ) / K in order to adapt to the unknown post-change distribution. (In fact, we will develop more sophisticated mixtures whose number of components grows slowly with time.) For later reference, we state the above as a proposition:

Let { M a } a ∈ A be a set of e-detectors. Then for any probability measure μ on A, the mixture of e-detectors, ∫ M a d μ ( a ) forms a valid e-detector.

An e-detector provides a quantification of evidence for whether a changepoint has occurred or not, and may be continuously monitored, stopped and easily interpreted – e.g, a steep and steady increase of the process in recent times should be taken as an indication that a change has taken place. The following theorem shows how one can immediately obtain a sequential change detection procedure from an e-detector M.

For any α ∈ ( 0 , 1 ) and e-detector M, if we declare a changepoint at the stopping time N ∗ : = inf { n ≥ 1 : M n ≥ 1 / α } , then we have inf P ∈ P E P , ∞ N ∗ ≥ 1 / α . That is, the sequential change detection procedure in (2.2) controls the ARL at level α.

The informal proof is one line long: dropping subscripts, the definition of an e-detector implies that E N ∗ ≥ E M N ∗ , but by definition of N ∗ , we know that M N ∗ ≥ 1 / α if N ∗ is finite. (If N ∗ is not almost surely finite, the claim holds anyway.) The full proof is in Appendix A.1.

 Constructing an E-detector Based on a Sequence of E-processes

The central building block of our e-detector is called an e-process. E-processes are newly-developed tools that have been shown to play a fundamental role in sequential hypothesis testing, especially in composite, nonparametric settings. E-processes are generalizations of nonnegative martingales and supermartingales, and in particular, e-processes are nonparametric and composite generalizations of likelihood ratios. They have strong game-theoretic roots, and have found utility in the meta-analysis, as well as for the purposes of anytime-valid inference in the presence of continuous monitoring [28, 29, 7, 11, 12]. The properties of e-processes have not yet been explored in changepoint analysis, and we undertake this effort here.

To understand their definition, we briefly forget about sequential change detection and consider testing the null hypothesis that X 1 , X 2 , … ∼ P for some P ∈ P. An e-process for P, called a P-e-process, is a sequence of nonnegative random variables ( E t ) t ≥ 1 such that for any P ∈ P and any stopping time τ, we have E P [ E τ ] ≤ 1. As before, the underlying filtration can be that of the data or that of ( E t ), or some enlargement of these. The value of E t measures evidence against the null (larger values, more evidence). A level-α sequential test can be obtained by rejecting the null as soon as E t exceeds 1 / α; this is a consequence of Ville’s inequality [42, 11].

As a result of the optional stopping theorem, nonnegative P-martingales (i.e., the process is a nonnegative P-martingales simultaneously for every P ∈ P) and P-supermartingales are examples of e-processes. However, e-processes are a distinct and more general class of processes. In fact, there exist natural classes P for which the only P-martingales are constants, and the only P-supermartingales are decreasing sequences, but there are e-processes for P that can increase to infinity when the data are not from P. See [29] for one such example, arising from sequentially testing exchangeability of a binary sequence, and [6] for another example arising from testing log-concavity.

In this subsection, we show how to leverage e-processes to build e-detectors.

(<inline-formula id="j_nejsds51_ineq_124"><alternatives><mml:math> <mml:msub> <mml:mrow> <mml:mi mathvariant="normal">e</mml:mi> </mml:mrow> <mml:mrow> <mml:mi mathvariant="italic">j</mml:mi> </mml:mrow> </mml:msub></mml:math><tex-math><![CDATA[${\mathrm{e}_{j}}$]]></tex-math></alternatives></inline-formula>-process).

For any j ≥ 1, Λ ( j ) : = { Λ n ( j ) } n ≥ 1 is called an e j -process for P if it is a nonnegative adapted process such that Λ 1 ( j ) = ⋯ = Λ j − 1 ( j ) = 1, and sup P ∈ P E P , ∞ [ Λ τ ( j ) ∣ F j − 1 ] ≤ 1 , ∀ τ ∈ T .

When j = 1, the e j -process is simply a standard e-process, which tests whether the data distribution is different from the proclaimed null (pre-change) distribution. For j > 1, each e j -process can be viewed as an e-process that begins at time j, and tests whether the data occurring after time j are well explained by the null hypothesis.

(SR and CUSUM e-detectors).

Based on a sequence of e-processes { Λ ( j ) } j ≥ 1 , define SR and CUSUM e-detectors M SR and M CU , respectively by M 0 SR = M 0 CU : = 0 and for each n ≥ 1, M n SR : = ∑ j = 1 n Λ n ( j ) , and M n CU : = max j ∈ [ n ] Λ n ( j ) .

It is not hard to check that the above processes are indeed e-detectors, meaning they satisfy (2.1).

If each e-process starts with an initial weight smaller than 1 such that the sum of all initial weights is less than or equal to 1, then corresponding SR and CUSUM procedures control the probability of the false alarm, sup P ∈ P P P , ∞ ( τ < ∞ ) ≤ α, which is a much more stringent error metric than the ARL. The price to pay is that if there is a change at time ν, then the detection delay will no longer be independent of ν, and will typically increase logarithmically with ν (so that the worst average delay is unbounded).

 Constructing Computationally Efficient E-detectors Using Baseline Increments

In general, it may take O ( n ) time to update the aforementioned e-detectors at time n. In order to construct e-detectors that can be updated online in sublinear time and memory (or even near-constant time and memory), it turns out to be computationally convenient to use a common “baseline” increment in order to build the underlying e j -processes, as we do below. Effectively, this amounts to using e j -processes that are P-supermartingales, which is a special case of particular interest.

(Baseline increment).

A nonnegative, adapted process L : = { L n } n ≥ 1 is called a baseline increment if for each n ≥ 1, we have sup P ∈ P E P , ∞ [ L n ∣ F n − 1 ] ≤ 1 .

It is easy to check that if L 1 and L 2 are baseline increments, and A 1 and A 2 are nonnegative and predictable processes (meaning that A n 1 and A n 2 are both F n − 1 -measurable) such that A 1 + A 2 is strictly positive, then the mixture ( A 1 L 1 + A 2 L 2 ) / ( A 1 + A 2 ) also forms a baseline increment. In short, “predictable mixtures” retain the baseline increment property.

Comparing (2.6) to (2.4), we see that a baseline increment L is not itself an e j -process, because the expectation in (2.6) applies only at fixed times n, with the conditioning being on the previous step n − 1, but (2.4) calculates expectations at any stopping time, and conditions on j − 1. It is best to think about the baseline increment as the multiplicative increment that forms the e j -process, as follows.

(Baseline <inline-formula id="j_nejsds51_ineq_159"><alternatives><mml:math> <mml:msub> <mml:mrow> <mml:mi mathvariant="normal">e</mml:mi> </mml:mrow> <mml:mrow> <mml:mi mathvariant="italic">j</mml:mi> </mml:mrow> </mml:msub></mml:math><tex-math><![CDATA[${\mathrm{e}_{j}}$]]></tex-math></alternatives></inline-formula>-process).

For a given baseline increment L : = { L n } n ≥ 1 , we define the corresponding “baseline e j -process” Λ ( j ) , for each j , n ∈ N, as below: Λ n ( j ) : = 1 if n < j ∏ i = j n L i otherwise .

Under any pre-change distribution P ∈ P, each Λ ( j ) is a nonnegative supermartingale by definition of the baseline increment L i . Therefore, a straightforward application of the optional stopping theorem implies that each Λ ( j ) satisfies condition (2.4), and thus is a valid e j -process.

As an example of a baseline e j -process, consider the case where we have iid observations X 1 , X 2 , … from a distribution p θ parameterized by θ ∈ Θ, and the pre-change distribution is given by θ 0 . Then, for any post-change distribution p θ 1 with θ 1 ≠ θ 0 , the likelihood ratio between two distributions, L n : = p θ 1 ( X n ) / p θ 0 ( X n ) yields a baseline increment process with the inequality in (2.6) being replaced by the equality. Then, each Λ n ( j ) is the likelihood ratio based on X j , … , X n . Further, instead of using a fixed post-change parameter θ 1 , we can also plug-in a running MLE or any other online nonanticipating estimator that is based on the previous history F n − 1 only, say θ ˆ n − 1 , into the likelihood ratio. In this case, although the value L n : = p θ ˆ n − 1 ( X n ) / p θ 0 ( X n ) at time n of the resulting process may depend on the previous history F n − 1 , the inequality in (2.7) will be satisfied as an equality, yelling a valid baseline e j -processes.

While baseline increment processes provide a natural and computationally convenient way to construct e j -process, we emphasize that any e-detector, even one that does not use baseline increments, will automatically control the ARL by Theorem 2.4. To elaborate, baseline e j -processes are composite P-supermartingales (meaning P-supermartingales for every P ∈ P), but there exist other P-e-processes—which are not P-supermartingales—that naturally arise and these can be used to form e-detectors; for example using universal inference [47, Section 8].

 (Baseline SR and CUSUM e-detectors).

When an SR or CUSUM e-detector is constructed using a sequence of baseline e j -processes (2.7), we call it a “baseline SR or CUSUM e-detector”.

Each baseline SR or CUSUM e-detector can be computed recursively like their classical analogs: M n SR = L n · [ M n − 1 SR + 1 ] , M n CU = L n · max { M n − 1 CU , 1 } , with M 0 SR = M 0 CU = 0 for each n ∈ N. The above computational benefit is the primary reason to consider baseline e-detectors, but as mentioned in Remark 2.10 and when introducing e-processes, more general e-detectors are sometimes necessary for certain classes P.

We briefly verify below that the processes M SR and M CU defined above are valid e-detectors. Indeed, for any stopping time τ and pre-change distribution P ∈ P, if P P , ∞ ( τ = ∞ ) > 0 then the condition of the e-detector in (2.1) holds trivially. If not, then τ is finite almost surely, and we have by linearity of expectation and the tower rule: E P , ∞ M τ CU ≤ E P , ∞ M τ SR = E P , ∞ ∑ j = 1 ∞ Λ τ ( j ) 1 ( j ≤ τ ) = ∑ j = 1 ∞ E P , ∞ Λ τ ( j ) 1 ( j ≤ τ ) = ∑ j = 1 ∞ E P , ∞ [ 1 ( j ≤ τ ) E P , ∞ [ Λ τ ( j ) ∣ F j − 1 ] ] ≤ ∑ j = 1 ∞ E P , ∞ 1 ( j ≤ τ ) = E P , ∞ τ , where the first inequality comes from the nonnegativity of e-processes, and the second inequality comes from the definition of the e-process Λ ( j ) for each j ≥ 1. Note that this proof is also applicable to general SR and CUSUM e-detectors.

 Sequential Change Detection Procedures by Thresholding E-detectors

The value of any e-detector process, like M SR or M CU , is directly interpretable without specifying an explicit threshold: a larger value signals an accumulation of evidence of a changepoint. These can be monitored and adaptively stopped. Nevertheless, to explicitly control the ARL at level α, we define SR and CUSUM-style change detection procedures, called “e-SR” and “e-CUSUM” procedures as follows.

(e-SR and e-CUSUM procedures).

Given SR and CUSUM e-detectors M SR and M CU , define e-SR and e-CUSUM procedures by the stopping times N SR ∗ : = inf { n ≥ 1 : M n SR ≥ 1 / α } , N CU ∗ : = inf { n ≥ 1 : M n CU ≥ c α } , where c α is a constant chosen to control the ARL of the e-CUSUM procedure by 1 / α for some α ∈ ( 0 , 1 ). By Theorem 2.4, 1 / α is a valid choice for c α .

We note c α = 1 / α may be a very conservative choice for the e-CUSUM procedure. Indeed, suppose we use the trivial e-processes, that is, we set Λ n ( j ) : = 1 for all j , n. In this setting, the SR e-detector is given by M n SR = n for each n. In contrast, the CUSUM e-detector, M n CU is equal to 1 for all n. Therefore, any valid threshold we can choose for the e-SR procedure must be larger than ⌊ 1 / α ⌋. On the other hand, any threshold above 1 makes N CU ∗ = ∞, which of course controls ARL by 1 / α, but the true ARL is much above the target. Building from this trivial example, it is possible to construct nontrivial examples in which letting α → 0 makes the gap between tight thresholds of e-SR and e-CUSUM procedures arbitrarily different.

Unless we assume the pre-change distribution is time-stationary, known and parametric, or we can access a good sample of the pre-change distribution or large enough historical data, computing a tight or even valid threshold c α can be a challenging task. In the application sections below, we will mainly deal with non-stationary pre-change distributions where pre-change observations may not be identically distributed and thus all observations before the changepoint may follow different distributions. In this case, setting c α = 1 / α seems to be the only reasonable choice, and we recommend using the e-SR procedure rather than e-CUSUM since if we use the same threshold for both procedures, the former always detects the changepoint faster than the latter while provably controlling the ARL at the same level.

 Some Nontrivial Instantiations of E-detectors

E-detectors can be thought of as a general reduction of change detection to sequential testing. Given the recent advances in nonparametric, composite sequential testing using nonnegative supermartingales and e-processes, our e-detectors now make new classes of change detection problems possible. We detail below some interesting nontrivial examples of change detection problems that can now be solved using e-detectors.

(When Likelihood Ratios Are Well-defined).

Consider first the parametric case, when P , Q have a common reference measure and likelihood ratios are well-defined. When the pre- and post-change distributions are known (meaning P , Q are singletons), then the standard likelihood-ratio based CUSUM and SR processes from (1.4) and (1.5) are both e-detectors. If Q is composite, then taking a mixture likelihood ratio yields e-detectors (using either a non-anticipating predictable mixture or using a fixed mixture distribution). If P is also composite, but maximum likelihood estimation is efficient over P, then e-processes can be constructed that take the ratio of mixture likelihoods over the alternative to maximum likelihood under the null, as done in universal inference [47, Section 8] or in [40]. If the “reverse information projection” (RIPr) is computable (analytically or numerically), then one can use the method of [7], though there are some subtleties: by default the method produces e-values over blocks of observations, which can be multiplied across independent blocks to produce a nonnegative supermartingale (and thus e-process) to be used within our framework. But sometimes, the sequence of RIPr’s over increasing sample sizes (nested blocks) automatically produces an e-process, and when this happens, it is more powerful than universal inference (see [7, 30] for details).

 (Change in Distribution).

In this example, P is the set of all iid product distributions over infinite sequences (or its convex closure, the set of all exchangeable distributions), so P = { μ ∞ for some probability distribution μ }. The conformal sequential change detection procedures by [44, 43] are designed to test deviations from exchangeability, meaning that they develop a test martingale (and thus an e-process) for P, meaning that their procedure fits neatly into our framework. Importantly, their filtration is restricted, and is smaller than the natural filtration of the data. This allows nonparametric martingales to exist, but the e-detector property only holds with respect to a smaller class of stopping times. Nevertheless, thresholding our e-detector at 1 / α still controls the ARL at level α, as the latter property is independent of the filtration used to construct the e-detector. As a side note, if one wanted to construct an e-detector that was valid at all stopping times with respect to the natural filtration of the data, e-detectors based on martingales provably do not suffice, but e-detectors based on e-processes can be constructed using the techniques from [29] (at least for categorical distributions).

 (Nonparametric Two-sample Testing).

Suppose we have two streams of (general multivariate) data: X 1 , X 2 , … ∼ P X and Y 1 , Y 2 , … ∼ P Y . For simplicity below, assume that at time t, we observe one point from each stream ( X t , Y t ). Before the changepoint (if one exists), the distributions of X t and Y t are equal, meaning that P = { ( P X , P Y ) ∞ : P X = P Y }. This is a very nonparametric class, since it specifies nothing about the distributions except for the fact that they are equal before the changepoint. After the changepoint, the streams have different distributions (maybe the distribution of X changes, or that of Y changes, or both), thus Q = { ( P X , P Y ) ∞ : P X ≠ P Y }. For this very general nonparametric two-sample testing setup, [34] construct test martingales for P that are provably consistent against Q under minimal assumptions (in particular not requiring any minimum separation between the different distributions after the changepoint, since the tests automatically adapt to the closeness of the unknown alternative). These test martingales fit seamlessly into our e-detector framework, yielding new and practicable e-detectors for detecting a change from homogeneity to non-homogeneity between the streams.

 (Nonparametric Independence Testing).

In this problem setting, we observe a pair of random variables ( X t , Y t ) ∼ P X Y at each step t, where X t , Y t can each lie in a general space. Before the changepoint (if one exists), the data are independent, meaning that P = { P X Y ∞ : P X Y = P X × P Y }. Beyond saying that the joint distribution factorizes into the product of marginals, there is no further structure assumed, making this a rich nonparametric composite class. As before, after the changepoint, ( X , Y ) become dependent, meaning that Q = { P X Y ∞ : P X Y ≠ P X × P Y }. For this general nonparametric independence testing problem, [24] construct test martingales for P that are provably consistent against Q under minimal, weak assumptions (as before, not requiring any separation). Again as before, (in the testing problem) the power of these tests automatically adapts to the difficulty of the unknown alternative. When plugged into our framework, it delivers a novel e-detector for a change from independence to dependence.

We briefly remark that in the two preceding examples (homogeneity and independence), we can move past the iid assumption. The same methods work even when the distribution is allowed to drift within P before the changepoint, and drift within Q after the changepoint. We refer to the original aforementioned papers for details.

(Log-concavity).

Here, the data before the changepoint comes from a log-concave distribution (in a general dimension d ≥ 1), so P = { μ ∞ : μ has a log-concave Lebesgue density }. This is a rich, nonparametric, shape-constrained class. The post-change class of distributions Q consists of, for example, any distribution that has a nonzero KL-divergence and Hellinger distance from every distribution in P. For testing P against Q, [6] show that there exists no nontrivial nonnegative supermartingales, but they design a powerful (and computationally efficient) e-process using universal inference. When plugged into our e-detector, this yields a nontrivial procedure that can detect a deviation from log-concavity.

 (Symmetry).

Suppose P = { μ ∞ : μ is symmetric around 0 } consists of the set of all distributions (in a general dimension d ≥ 1) that are symmetric around the origin, while Q consists its complement (that is, distributions which are not symmetric around the origin). [28] characterize all processes that are nonnegative martingales for P. When used with our e-detector, these provide a clean way to detect a change from symmetry to asymmetry.

 (Change in Mean).

Suppose P C = { μ ∞ : E X ∼ μ [ X ] ≤ 0 , X satisfies C } consists of the set of all univariate distributions with mean less than or equal to zero satisfying come constraint C, while Q C = { μ ∞ : E X ∼ μ [ X ] > 0 , X satisfies C } consists of those with positive mean. [11] provides a large variety of nonnegative supermartingales under various conditions C, such as when X are subGaussian, or bounded from above, or bounded from below, or have only two or three moments; see also [45] for heavy-tailed supermartingales. These can be plugged into the e-detector to yield new nonparametric schemes for changes in mean.

 (Huber-robust Change Detection).

As a final example, suppose we wish to detect a change in mean of heavy-tailed data (as above). But now, suppose that an adversary can also arbitrarily corrupt an ϵ fraction of the data. [46] develop Huber-robust supermartingales for this setting, which can be plugged into e-detectors to yield a valid e-detector in the presence of adversarial corruptions.

Note that in Examples 3, 4, 5 and 6, if P and Q were swapped, the testing problem is much harder, and we are not aware of any powerful test or change detection method. However, if one was simply interested in detecting a change in homogeneity or in dependence, i.e. from some distribution in P ⋃ Q to some other one, there are two possible change detection methods that come to mind. First, an e-detector based on the conformal change detection methods in Example 2 would detect any change from any distribution to any other, though choosing the conformity score may be tricky. As a second and more direct option, one can choose a measure of homogeneity or dependence (like the kernel maximum mean discrepancy or energy distance, or the Hilbert Schmidt independence criterion or distance covariance), construct a confidence sequence for that measure (see [19] for several specific, tight, constructions), and plug it into the recent change detection scheme of [35].

Finally, it is worth noting that it is possible to define e-detectors in cases where there is no common reference measure amongst the pre-change and post-change distributions, and thus no easily-defined likelihood ratio process, and also when there is no nontrivial martingale that can be constructed. This is precisely the utility of e-processes, which are nonparametric and composite generalizations of likelihood ratios. We develop one more interesting nonparametric example (not described above) in the simulations section: detecting change in mean of a bounded random variable.

 Warm-up: Bounds on Worst Average Delays for Baseline E-detectors, when <italic>Q</italic> is Known

Recall that our objective is to minimize worst average delays J L ( N ∗ ) or J P ( N ∗ ), given above in (1.2) and (1.3) respectively, while controlling the ARL E P , ∞ ( N ∗ ) ≥ 1 / α. Note that worst average delays in (1.2) and (1.3) were defined for a fixed pair of pre- and post-change distributions implicitly. In our setting where the pre-change distribution space could be composite, we take an additional supremum over all pre-change distributions when defining both worst average delays for each fixed post-change distribution.

To derive bounds on the worst average delays, we further assume that the post-change observations X ν + 1 , X ν + 2 , … are independent of the pre-change observations and form a strongly stationary process. That is, we assume that, for any finite subset I ⊂ N and any j ∈ N, the joint distributions of { X ν + i } i ∈ I and { X ν + i + j } i ∈ I are equal to each other. About the underlying baseline increment, we further assume that there exist a function f and an integer m ≥ 0 such that L n = f ( X n , X n − 1 , … , X n − m ) for each n. In this warm-up section, assume that we know the post-change distribution Q. Then, as we shall soon see, an optimal choice for L n would set m = 0, but if m is a strictly positive number then we implicitly assume that we can access m observations X 0 , X − 1 , … , X 1 − m from the pre-change distribution in order to build sequential change detection procedures. Under these conditions, the following theorem provides analytically more tractable upper bounds on worst average delays for e-SR and e-CUSUM procedures.

For a given α ∈ ( 0 , 1 ), let N SR ∗ and N CU ∗ be e-SR and e-CUSUM procedures using baseline e-detectors. Under the settings described above, their worst average delays are upper bounded as J P ( N SR ∗ ) ≤ J L ( N SR ∗ ) ≤ E 0 , Q N 1 / α + m , J P ( N CU ∗ ) ≤ J L ( N CU ∗ ) ≤ E 0 , Q N c α + m , respectively, where N c is the stopping time defined for any c > 1 as N c : = inf { n ≥ 1 : ∑ i = 1 n log L i ≥ log ( c ) }, and c α ≤ 1 / α is any threshold that ensures the e-CUSUM procedure (2.11) has an ARL no smaller than 1 / α. Furthermore, if the post-change observations are iid and each L n is a function of X n only (i.e. m = 0) with E 0 , Q log L 1 > 0, then E 0 , Q N c ≤ log ( c ) E 0 , Q log L 1 + V 0 , Q log L 1 [ E 0 , Q log L 1 ] 2 + 1 .

The proof can be found in Appendix A.1.

The stopping time N 1 / α delivers a level-α sequential test for the null hypothesis H 0 : P ∈ P. On the other hand, the stopping time N c α may not necessarily control the type-1 error by α since the threshold c α can be significantly smaller than 1 / α as discussed earlier.

The expected stopping time in the upper bounds (2.12) and (2.13) depends on the parameter m ≥ 0. When Q is known, m = 0 suffices because (2.14) suggests that L i should simply be chosen to maximize E 0 , Q log L 1 , which is identical to the log-optimality criterion used for testing P against Q (as discussed in many recent works, like [33, 7, 48]).

In applications when Q is unknown, the underlying baseline increment may require a long enough sample history (large m) in order to achieve a reasonably small expected stopping time by “learning” Q or using an empirical distribution plug-in for Q. Then, the above results suggest that a reasonable way to choose the window size m is to pick the one minimizing the upper bound on the worst average delays. However, since the optimal choice of the window size should depend on the unknown post-change Q, it remains difficult to minimize the upper bound directly. In our simulations, we often encounter cases where a larger window size is better. In this case, we would choose a window size as large as possible while keeping the procedure computationally tractable. However, the right way to handle unknown Q is dealt with in detail next.

 Combining Baseline E-detectors Using the Method of Mixtures

In the previous section, we discussed how one can construct a valid e-detector and derive upper bounds on worst average delays. However, in most composite and nonparametric sequential change detection scenarios, there is no single optimal e-detector but instead there are often several applicable e-detectors to choose from. In this section, we introduce a practicable and computationally efficient strategy to construct a good e-detector for minimizing the upper bound on worst average delays in Proposition 2.14, especially for the upper bound (2.14) in the m = 0 case.

In detail, suppose we have a set of baseline increments { L λ } λ ∈ Π parametrized by λ ∈ Π. Then, under the additional condition assumed in the second part of Proposition 2.14 (namely, that m = 0 and the post-change observations are from an iid sequence), an ideal choice of the parameter λ op for a post-change distribution Q is given by λ op ( Q ) = arg max λ ∈ Π E 0 , Q log L 1 ( λ ) , minimizing the first term of the upper bound in (2.14), which often becomes a leading term especially for small enough α. In turn, this term is inversely proportional to E 0 , Q log L 1 ( λ op ) : = D ( Q | | P ) , where the second argument P in D ( Q | | P ) explicitly refers to the dependency of the baseline increment L λ op to the class of pre-change distributions P. For the rest of the paper, we will assume that the set of baseline increments, { L λ } λ ∈ Π is rich enough such that D ( Q | | P ) > 0 for all Q ∈ Q. As we observe later, in many canonical cases, we have D ( Q | | P ) = inf P ∈ P KL ( Q | | P ) where KL ( Q | | P ) is the Kullback-Leibler (KL) divergence from Q to P.

Generally, computing λ op is not feasible since it depends on the unknown post-change distribution Q. Next, we show how to build a mixture of baseline e-detectors that can detect the changepoint nearly as quickly as the one with λ op , when known lower and upper bounds λ L and λ U on λ op are available.

Notice that an average of e-detectors is also a valid e-detector, in the sense of satisfying condition (2.1). Therefore, for any mixing distribution W supported on [ λ L , λ U ], we can define mixtures of e-SR and e-CUSUM procedures by following stopping times: N mSR ∗ : = inf { n ≥ 1 : ∫ ∑ j = 1 n ∏ i = j n L i ( λ ) d W ( λ ) ≥ 1 / α } , N mCU ∗ : = inf { n ≥ 1 : ∫ max j ∈ [ n ] ∏ i = j n L i ( λ ) d W ( λ ) ≥ c α } , where c α > 1 is a fixed constant which controls the ARL for some α ∈ ( 0 , 1 ). Since the mixture of e-CUSUM procedure is based on a valid e-detector, we can always set the threshold c α to be equal to 1 / α as same as the threshold of the mixture of e-SR procedures.

Instead of using mixtures, one may be tempted to consider swapping the above integral with a supremum over λ ∈ [ λ L , λ U ]. However, this does not in general yield a valid e-detector.

Computational and Analytical Aspects of Mixtures of Baseline E-detectors

Though any mixing distribution yields a valid e-detector, for computational efficiency, we only consider discrete mixtures where the support of mixing distribution has at most countably many elements. To be specific, let { ω k } k ≥ 1 be a set of nonnegative mixing weights with ∑ k ≥ 1 ω k = 1 and let { λ k } k ≥ 1 be the corresponding supporting set. For ease of notation, we denote L ( λ k ) : = L ( k ) for each k ≥ 1. Based on the set of nonnegative mixing weights and the corresponding set of baseline increments, we define mixtures of SR and CUSUM e-detectors as M 0 mSR = M 0 mCU : = 1, and for each n ∈ N, M n mSR : = ∑ k = 1 ∞ ω k ∑ j = 1 n ∏ i = j n L i ( k ) : = ∑ k = 1 ∞ ω k M n SR ( k ) , M n mCU = ∑ k = 1 ∞ ω k max j ∈ [ n ] ∏ i = j n L i ( k ) : = ∑ k = 1 ∞ ω k M n CU ( k ) . Let K : = | { k : ω k > 0 } | be the number of nonzero mixing weights.

Finite Mixtures If K < ∞, we may for simplicity assume that the first K weights ω 1 , … , ω K are the only nonzero values. In this case, we can compute mixtures of SR and CUSUM e-detectors by M n mSR = ∑ k = 1 K ω k M n SR ( k ) , and M n mCU = ∑ k = 1 K ω k M n CU ( k ) , where M n SR ( k ) and M n CU ( k ) are computed recursively as M n SR ( k ) = L n ( k ) · [ M n − 1 SR ( k ) + 1 ] , M n CU ( k ) = L n ( k ) · max { M n − 1 CU ( k ) , 1 } with M 0 SR ( k ) = M 0 CU ( k ) = 0 for each k ∈ [ K ] and n ∈ N. Therefore, if each computation of L n ( k ) has constant time and space complexities, then the evaluation of mixtures of SR and CUSUM e-detectors at each time n requires O ( K ) time and space complexity.

Infinite Mixtures, Scheduling Functions and Adaptive Re-weighting If K = ∞ or if K is to be chosen adaptively as an increasing function of n we modify our strategy as follows. We first choose an increasing function K : N → N, and let K − 1 : N → N be the generalized inverse function of K defined by K − 1 ( k ) : = inf { j ≥ 1 : K ( j ) ≥ k } for each k ∈ N. Note that K − 1 is also an increasing function. We call such function K as a scheduling function. We intentionally overload notation: in what follows, K ( n ) plays the same role as the constant K in the case of finite support. Note that K − 1 ( k ) ≤ n for any k ≤ K ( n ); we will use this simple fact below when defining nested summations.

Based on a scheduling function K and its generalized inverse K − 1 , we define adaptive SR and CUSUM e-detectors, M n aSR and M n aCU , respectively, as M n aSR = ∑ k = 1 K ( n ) ω k ∑ j = K − 1 ( k ) n γ j ∏ i = j n L i ( k ) : = ∑ k = 1 K ( n ) ω k M n SR ( k ) , M n aCU = ∑ k = 1 K ( n ) ω k max K − 1 ( k ) ≤ j ≤ n γ j ∏ i = j n L i ( k ) : = ∑ k = 1 K ( n ) ω k M n CU ( k ) , where each γ j : = 1 / ∑ k = 1 K ( j ) ω k ≥ 1 is the adaptively re-weighting factor at time j, ensuring that the mixing weights always sum to one at each time. Here, we restrict not only the space over the index k from [ 1 , ∞ ] to [ 1 , K ( n ) ] but also the space over the index j from [ 1 , n ] to [ K − 1 ( k ) , n ]. This choice makes it possible to compute both M n aSR and M n aCU efficiently since each M n SR ( k ) and M n CU ( k ) have following recursive representations M n SR ( k ) = L n ( k ) · [ M n − 1 SR ( k ) + γ n ] , M n CU ( k ) = L n ( k ) · max { M n − 1 CU ( k ) , γ n } , for each n ≥ K − 1 ( k ) and M n aSR ( k ) = M n aCU ( k ) = 0 for all n = 0 , 1 , … , K − 1 ( k ) − 1. Therefore, if each computation of L n ( k ) has constant time and space complexities then the computations of adaptive SR and CUSUM e-detectors at each time n have O ( K ( n ) ) time and space complexities as well. For the purpose of implementing an online algorithm, we are typically interested in the case K ( n ) = O ( log ( n ) ).

Both mixtures of SR and CUSUM e-detectors can be viewed as special cases of their adaptive counterparts where the scheduling function K is understood as a constant function. In this case, we have γ j = ∑ k = 1 K ω k = 1 for each j, and thus M n aSR = M n mSR and M n aCU = M n mCU for each n ∈ N.

Unlike finite mixtures, the mixing distribution deployed in the adaptive SR and CUSUM e-detectors vary over time. Hence, we cannot simply apply Proposition 2.3 to check whether this adaptive scheme yields valid e-detectors. The following proposition formally states the validity of adaptive SR and CUSUM e-detectors. The proof can be found in Appendix A.2.

For any mixing weights { ω k } k ∈ N and a scheduling function K, adaptive SR and CUSUM e-detectors defined in (3.10) and (3.11) form valid e-detectors satisfying the condition (2.1).

Now, based on M n aCU and M n aSR , the adaptive e-SR and e-CUSUM procedures are defined by the stopping times: N aSR ∗ : = inf { n ≥ 1 : M n aSR ≥ 1 / α } , N aCU ∗ : = inf { n ≥ 1 : M n aCU ≥ c α } , where α ∈ ( 0 , 1 ) is a fixed constant and c α is a positive value controlling ARL of the adaptive e-CUSUM procedure by 1 / α. Similar to the usual e-CUSUM procedure case we discussed before, we can always set c α = 1 / α. In this case, from the fact N aSR ∗ ≤ N aCU ∗ , which is implied by M n aSR ≥ M n aCU , we have E P , ∞ N aCU ∗ ≥ E P , ∞ N aSR ∗ ≥ 1 / α , where the last inequality comes from Theorem 2.4 with the fact that M aSR is a valid e-detector. However, the threshold c α for the adaptive e-CUSUM procedure can be chosen to be a significantly smaller value if we have enough knowledge about the pre-change distribution, as discussed in Section 2.5.

 Worst Average Delay Analysis for Adaptive Mixtures of E-detectors

We now derive general upper bounds on worst average delays of the adaptive e-SR and e-CUSUM procedures. As we did before, we further assume that post-change observations X ν + 1 , X ν + 2 , … are independent of the pre-change observations and form a strong stationary process. Also, we further assume that there exist a function f k and an integer m ≥ 0 such that L n ( k ) = f k ( X n , X n − 1 , … , X n − m ) for each k and n. Again, if m is a strictly positive number then we implicitly assume that there exist m observations X 0 , X − 1 , … , X 1 − m from the pre-change distribution we can use to build sequential change detection procedures. Under this additional condition for worst average delay analysis, the following theorem provides analytically more tractable upper bounds on worst average delays for N aSR ∗ and N aCU ∗ .

Under additional conditions described above, worst average delays for N aSR ∗ and N aCU ∗ can be upper bounded as follows: J P ( N aSR ∗ ) ≤ J L ( N aSR ∗ ) ≤ min j > m [ E 0 , Q N 1 / α ( j ) + j − 1 ] , J P ( N aCU ∗ ) ≤ J L ( N aCU ∗ ) ≤ min j > m [ E 0 , Q N c α ( j ) + j − 1 ] , where, for j ∈ N and c > 0, N c ( j ) is the stopping time N c ( j ) : = inf { n ≥ 1 : ∑ k = 1 K ( j ) ω k ∏ i = 1 n L i ( k ) ≥ c } . Here, c α is the same threshold used to build the adaptive e-CUSUM procedure in (3.15). Note that for mixtures of the SR and CUSUM e-detectors where the scheduling function K is a constant function, the stopping times in the upper bounds do not depend on the index j, and thus the upper bounds can be reduced as J P ( N mSR ∗ ) ≤ J L ( N mSR ∗ ) ≤ E 0 , Q N 1 / α + m , J P ( N mCU ∗ ) ≤ J L ( N mCU ∗ ) ≤ E 0 , Q N c α + m , where, for c > 0, N c is the stopping time N c : = inf { n ≥ 1 : ∑ k = 1 K ω k ∏ i = 1 n L i ( k ) ≥ c } .

The proof of upper bounds on worst average delays can be found in Appendix A.2.

Unlike the baseline e-detector case, however, due to mixing weights, it is nontrivial to get further simplified upper bounds on worst average delays as we did in Section 2.7. Next, we present specific adaptive e-SR and e-CUSUM procedures based on exponential baseline increments where we can compute both procedures efficiently and derive upper bounds on worst average delays in explicit forms.

 Exponential Baseline E-detectors and Their Mixtures

Building upon recent advances in time uniform concentration inequalities and sequential testing developed in [11] and [36], below we consider an exponential structure on baseline e-detectors. We show that, in this setting, it is possible to approximate the “oracle” e-SR and e-CUSUM procedures based on the knowledge of the optimal (but unknown) λ op by adaptive procedures built using a mixture of carefully chosen set of baseline increments { L λ k } k ≥ 1 with mixing weights { ω k } k ≥ 1 .

To be specific, assume there exists an extended real-valued convex function ψ on R that is finite and strictly convex on a set Π ⊂ R containing 0 in its interior Π o . Furthermore, assume ψ is continuously differentiable on Π o with ∇ ψ ( 0 ) = 0 = ψ ( 0 ). Then define the “exponential baseline increment” as follows.

(Exponential baseline increment).

For each n ∈ N and λ ∈ Π, define L n λ = exp { λ s ( X n ) − ψ ( λ ) v ( X n ) } , where s is a real-valued function and v is a positive function on the sample space. L λ : = { L λ } n ≥ 1 is called an exponential baseline increment if it satisfies condition (2.6) in Definition 2.8.

Above, s and v are mnemonics for sum and variance. For each Q ∈ Q, define μ ( Q ) : = E 0 , Q s ( X 1 ) , σ 2 : = E 0 , Q v ( X 1 ) and Δ op ( Q ) : = μ ( Q ) σ 2 ( Q ) , where we assume that all expectations are finite. The following proposition provides an explicit expression for D ( Q | | P ) : = max λ ∈ Π E 0 , Q log L 1 ( λ ) and a sufficient condition to have D ( Q | | P ) > 0 when the underlying baseline increments have the form specified in (4.1).

For a fixed Q ∈ Q, suppose there exist λ op ∈ Π o such that Δ op ( Q ) = ∇ ψ ( λ op ). Then, D ( Q | | P ) = E 0 , Q log L 1 ( λ op ) = ψ ∗ ( Δ op ( Q ) ) σ 2 ( Q ) ≥ 0 , where ψ ∗ is the convex conjugate of ψ. Thus, if Δ op ( Q ) ≠ 0, we have D ( Q | | P ) > 0.

The proof of Proposition 4.2 can be found in Appendix B.1. For the rest of the section, we assume that { λ ∈ R : Δ op ( Q ) = ∇ ψ ( λ ) , Q ∈ Q } ⊂ Π , Δ op ( Q ) ≠ 0 , ∀ Q ∈ Q . Then, Proposition 4.2 implies D ( Q | | P ) > 0 for all Q ∈ Q. Also, for ease of notation, we will drop the dependency of Q from related parameters and simply write λ op , μ , σ 2 and Δ.

The exponential structure of the baseline increment in (4.1) results in a simple form of λ op such that λ op = ( ∇ ψ ) − 1 ( Δ op ) = ∇ ψ ∗ ( Δ op ) , where the second equality comes from the fact that λ = ∇ ψ ∗ ∘ ∇ ψ ( λ ) for each λ ∈ Π. Although λ op still depends on the unknown post-change distribution Q via Δ op , in many cases, we can find upper and lower bounds on Δ op . In this section, we explain how to use the knowledge of the range of Δ op to build a mixture of exponential baseline e-detectors that has explicit upper bounds on worst average delays.

Separated Pre- and Post-change Distributions

Suppose we have knowledge of upper and lower bounds on the parameter Δ op given in (4.2), i.e. we know a pair ( Δ L , Δ U ) such that Δ L < Δ op < Δ U . It then follows that λ L < λ op < λ U , where λ L = ∇ ψ ∗ ( Δ L ), λ op = ∇ ψ ∗ ( Δ op ) and λ U = ∇ ψ ∗ ( Δ U ). To simplify presentation, we only consider the one-sided and well-separated case: 0 < λ L < λ U .

Let 1 / α be the target level of the ARL control for a fixed α ∈ ( 0 , 1 ). Let { L ( k ) } k ∈ [ K ] and { ω k } k ∈ [ K ] be K exponential baseline increments and mixing weights whose specific values will be defined later in this subsection. Since each L n ( k ) is a function of the n-th observation X n for each k ∈ [ K ], Theorem 3.4 implies that, if the post-change observations form a strong stationary process then the worst average delays for mixtures of e-SR and e-CUSUM procedures, N mSR ∗ and N mCU ∗ can be upper bounded by E 0 , Q N 1 / α and E 0 , Q N c α , respectively, where N c is the stopping time defined in (3.22). Furthermore, as we can always set the threshold for the e-CUSUM procedure to be c α ≤ 1 / α, we have E 0 , Q N c α ≤ E 0 , Q N 1 / α . Therefore, in this subsection, we construct a set of baseline increments for which we can derive a tight bound on E 0 , Q N 1 / α .

Algorithm 1 describes our methodology for computing mixtures of e-SR procedures in detail. The inputs to the algorithm are the upper and lower bounds Δ U and Δ L on Δ op and the maximal number of baselines processes K max . Mixture of e-CUSUM procedures can be executed similarly by replacing Line 10 by M n CU ( k ) ← exp { λ k s ( X n ) − ψ k v ( X n ) } · max { M n − 1 CU ( k ) , 1 } . Also, for the mixture of e-CUSUM procedures, we can replace the threshold 1 / α with a smaller value c α if we have enough information about the pre-change distribution. For both e-SR and e-CUSUM, at each time n, updates of mixtures of e-detectors have O ( K α ) time and space complexities, which do not depend on n.

Algorithm 1 relies critically on the function computeBaseline in Line 2, which returns a set of parameters and weights to compute a mixture of e-detectors along with a threshold value g α > 0 that will appear in the upper bound on worst average delays given in Theorem 4.3 that will be explained below. The details of computeBaseline are fairly technical and are given in Algorithm 3 in Appendix B.1.

In the main result of this section, we provide bounds on ARL and worst average delays for the mixtures of e-CP procedures obtained with Algorithm 1 that is a function of the parameter λ op and the threshold value g α . The proof can be found in Appendix B.1.

Let N mSR ∗ and N mCU ∗ be the stopping times corresponding to the mixtures of e-SR and e-CUSUM procedures in Algorithm 1 and its variant, respectively. Then, both procedures control the ARL by 1 / α. If we further assume that the post-change observations X ν + 1 , X ν + 2 , … are iid samples from a post-change distribution Q, then the worst average delays for N mSR ∗ and N mCU ∗ can be bounded as max { J L ( N mSR ∗ ) , J L ( N mCU ∗ ) } ≤ g α D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + 1 . The same bound holds also for J P ( N mSR ∗ ) and J P ( N mCU ∗ ).

In Proposition B.2 in Appendix B.1, we show that if the number of baseline processes K max in Algorithm 1 is chosen large enough, then the quantity g α returned by computeBaseline is at most inf η > 1 η [ log ( 1 / α ) + log ( 1 + ⌈ log η ψ ∗ ( Δ U ) ψ ∗ ( Δ L ) ⌉ ) ] , which can be easily evaluated numerically. Expression (B.15) in Appendix B.1 provides a precise formula for how large K max needs to be in order for the above bound to be in effect. In most practical cases, K max = 1000 is a large enough choice. Also, in many canonical examples we will present later, if we choose large enough K max satisfying the condition (B.15) then the first term g α D ( Q | | P ) of the upper bound of worst average delays in Theorem 4.3 become a leading term. In this case, from the inequality (4.7), we can check that this leading term is O ( log ( 1 / α ) / D ( Q | | P ) ) as α → 0.

If there is only one pre-change distribution P and one post-change distribution Q, both from a natural univariate exponential family, then their likelihood ratio forms an exponential baseline increment. In this case, the above upper bound becomes O ( log ( 1 / α ) / KL ( Q | | P ) ) as α → 0, matching the rate of the known lower bounds [18].

The bound on worst average delays in Theorem 4.3 is obtained by analyzing an auxiliary stopping time N ¯ g : = inf { n ≥ 1 : sup λ ∈ ( λ L , λ U ) ∑ i = 1 n log L i ( λ ) ≥ g } , g > 1 . Using the same arguments as in the proof of Proposition 2.14, we immediately have that if the post-change observations are iid from Q, then for any g > 1, E 0 , Q N ¯ g ≤ g D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + 1 . The bound (4.6) is finally established by showing that the stopping time N ¯ g α obtained by using the threshold g α produced by Algorithm 3 is a deterministic upper bound to the stopping times N c α and N 1 / α corresponding to mixtures of SR and CUSUM e-detectors. In detail, it holds that for any stream of observations X 1 , X 2 , … , N c α ≤ N 1 / α ≤ N ¯ g α . This nontrivial result is formally stated in Lemma B.1 in Appendix B.1. Its proof leverages geometric arguments used in [36, Theorem 2] to analyze sequential generalized likelihood ratio tests.

Pseudo-code of the mixture of e-SR procedures.

 Non-separated Pre- and Post-change Distributions

The previous subsection discussed how to build mixtures of e-SR and e-CUSUM procedures with an explicit upper bound on worst average delays when we have known and positive boundary values, λ L and λ U on the unknown λ op via the knowledge of Δ L < Δ op < Δ U . However, in many cases, we may not be fully certain about the boundary values. In this subsection, we discuss how we can generalize the previous argument to the no separation case whereby we only know the sign of λ op ( > 0 ) but do not have specific boundary values.

Recall that, for the well-separated case, we calibrated the mixtures of finitely many exponential baseline e-detectors using the stopping time N ¯ g α in (4.8), which is in turn based on the maximum of underlying baseline increments over the known upper and lower bounds of λ op . Since we no longer have knowledge of the boundary values λ L and λ U , we may use similar stopping times where the range of maximum and the threshold slowly increase over time. In this case, we need an infinite sequence of baseline procedures { L ( k ) } k ∈ N and mixing weights { ω k } k ∈ N to build adaptive e-SR and e-CUSUM procedures.

The bound in Theorem 3.4 along with the fact γ j ≥ 1 for all j ∈ N implies that, for any given scheduling function K : N → N, if the post-change observations form a strong stationary process then worst average delays for adaptive e-SR and e-CUSUM procedures can be upper bounded by min j ≥ 1 [ E 0 , Q N 1 / α ( j ) + j − 1 ] and min j ≥ 1 [ E 0 , Q N c α ( j ) + j − 1 ], respectively, where we recall that N c ( j ) is defined for c > 0 by N c ( j ) : = inf { n ≥ 1 : ∑ k = 1 K ( j ) ω k ∏ i = 1 n L i ( k ) ≥ c } . Again, since we can set the threshold for the e-CUSUM procedure in such a manner that c α ≤ 1 / α (so that E 0 , Q N c α ( j ) ≤ E 0 , Q N 1 / α ( j )), in this subsection, we focus on constructing a set of baseline increments on which we can derive a tight upper bound on min j ≥ 1 [ E 0 , Q N 1 / α ( j ) + j − 1 ].

To derive the set of baseline increments, we use a time-varying boundary function g. Here, we intentionally overload notation: the constant g in the previous subsection for the well-separation case can be viewed as a constant function g in what follows. Let g : [ 1 , ∞ ) → [ 0 , ∞ ) be a nonnegative and nondecreasing continuous function such that the mapping t ↦ g ( t ) / t is nonincreasing and lim t → ∞ g ( t ) / t = 0. For a chosen positive number Δ 0 > 0, let D 0 : = ψ ∗ ( Δ 0 ) and V 0 : = inf { t ≥ 1 : D 0 ≥ g ( t ) / t } . Now, for any fixed η > 1 and j ∈ N, define Δ 1 > Δ 2 > ⋯ as positive solutions of the equations ψ ∗ ( Δ k ) = g ( V 0 η k ) V 0 η k , k = 0 , 1 , 2 , … . Finally, based on the sequence { Δ k } k ≥ 0 , define λ k : = ∇ ψ ∗ ( Δ k ) , and set ω 0 : = α − 1 e − g ( V 0 ) 1 ( g ( V 0 ) > v min D 0 ), ω k : = α − 1 e − g ( V 0 η k ) / η for each k ∈ N where v min : = min x v ( x ), recalling the function v from Definition 4.1.

Based on the quantities defined above, we can construct the stopping time N 1 / α ( j ) for each j. The following lemma shows that we can upper bound the stopping time N 1 / α ( j ) with another stopping time N ¯ g ( j ) from which we can derive an explicit upper bound on its expected stopping time.

For any fixed j ≥ 1, Δ 0 > 0, and tuning parameter η > 1, let N 1 / α ( j ) be the stopping time based on the parameters defined above. Then, we have N 1 / α ( j ) ≤ N ¯ g ( j ) , where N ¯ g ( j ) is a stopping time defined by N ¯ g ( j ) : = inf { n ≥ 1 : sup λ ∈ ( λ K ( j ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V 0 η K ( j ) ) } .

Note that the chosen set of weights { ω k } k ≥ 0 yields valid adaptive e-SR and e-CUSUM procedures if e − g ( V 0 ) 1 ( g ( V 0 ) > v min D 0 ) + ∑ k = 1 ∞ e − g ( V 0 η k ) / η ≤ α . Once the above condition is satisfied, we can use the worst average delay analysis in Section 3.2 with the bound in Lemma 4.5 to get an explicit upper bound on the worst average delay of the adaptive e-SR and e-CUSUM procedures.

In detail, let j op be the smallest integer satisfying λ K ( j op ) < λ op and set K op : = K ( j op ). If we also have λ op < λ 0 then Lemma 4.5 implies E 0 , Q N 1 / α ( j op ) ≤ E 0 , Q N ¯ g ( j op ) ≤ E 0 , Q N op , where the stopping time N op is defined by N op : = inf { n ≥ 1 : ∑ i = 1 n log L i ( λ op ) ≥ g ( V 0 η K op ) } , and the expectation E 0 , Q N op is typically on the order of g ( V 0 η K op ) / D ( Q | | P ). Based on this observation, in the rest of this subsection, we introduce a practical and interpretable way to choose a boundary function g and related tuning parameters which minimize the leading term g ( V 0 η K op ) while satisfying the condition (4.15) on the set of mixing weights.

First note that, although we have no bounds on Δ op in the no separation case, we can still choose Δ L and Δ 0 with Δ L < Δ 0 as tuning parameters that represent our initial guess on the range of the unknown Δ op . Since it is possible that the unknown parameter Δ op of the post-change distribution is outside of the boundary ( Δ L , Δ 0 ), instead of assigning the entire α to the inside of the guessed interval, we split it into two parts by r α and ( 1 − r ) α, respectively where r ∈ ( 0 , 1 ) is another tuning parameter called the importance weight. Roughly speaking, larger r implies we make a higher bet on that the unknown Δ op is inside of our chosen boundaries ( Δ L , Δ 0 ).

Now, given tuning parameters Δ L , Δ 0 and r, we compute the set of { g r α , K L , η } by executing the function computeBaseline, just like in Algorithm 1, except that α is replaced by r α. Then, we can extend the boundary function g to accommodate the case in which the unknown Δ op is not inside the initial interval we had guessed. To be specific, we use the boundary function t ∈ [ 1 , ∞ ) ↦ g ( t ) : = g r α + s η log ( 1 + log η ( t V 0 η K L ∨ 1 ) ) , where V 0 : = g r α / D 0 and s > 1 is a constant obtained as the solution of the equation ζ ( s ) − 1 : = ∑ k = 1 ∞ 1 ( 1 + k ) s = e g r α / η [ α − { e − g r α 1 ( g r α > D 0 ) + K L e − g r α / η } ] . Note that the right hand side of the above equation is approximately equal to ( 1 − r ) α e g r α / η . Therefore, s ≈ ζ − 1 ( 1 + [ 1 − r ] α e g r α / η ) .

Pseudo-code of the adaptive e-SR procedures.

In Algorithm 2, we provide the detailed steps for the adaptive e-SR procedure based on the boundary function in (4.18). The algorithm can be easily modified for the adaptive e-CUSUM procedure by replacing the update in Line 15 with the rule M n CU ( k ) ← exp { λ k s ( X n ) − ψ k v ( X n ) } · max { M n − 1 CU ( k ) , γ } . Also, for the adaptive e-CUSUM procedure, we can replace the threshold 1 / α with a smaller value c α if we have enough information about the pre-change distribution.

In term of computational complexity, in Algorithm 2 we set the scheduling function K : N → N as K ( n ) : = K L + ⌈ m log η n ⌉ , where m ≥ 1 is a tuning parameter. Therefore, for both adaptive e-SR and e-CUSUM procedures, updates of statistics have O ( m log η n ) time and space complexities at each time n. Although it is not a fully online algorithm, logarithm time and space complexities make it feasible to run adaptive e-SR and e-CUSUM procedures in most practical online settings.

From Section 3.1, we know that both procedures control the ARL by 1 / α. The following theorem introduces explicit bounds on the worst average delays for both procedures.

Let N aSR ∗ and N aCU ∗ be stopping times corresponding to the adaptive e-SR procedures in Algorithm 2 and its and e-CUSUM variant, respectively. Then, both procedures control ARL by 1 / α. If we further assume that post-change observations X ν + 1 , X ν + 2 , … are iid samples from a post-change distribution then the worst average delays for N aSR ∗ and N aCU ∗ can be upper bounded as max { J L ( N aSR ∗ ) , J L ( N aCU ∗ ) } ≤ g r α D ( Q | | P ) ψ ∗ ( Δ op ) ψ ∗ ( Δ 0 ) + V 0 , Q [ log L 1 ( λ 0 ) ] [ D ( Q | | P ) ] 2 [ ψ ∗ ( Δ op ) ψ ∗ ( Δ 0 ) ] 2 + 1 if Δ op ≥ Δ 0 g r α D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + 1 if Δ op ∈ ( Δ L , Δ 0 ) g r α + s η log ( 1 + K op − K L ) D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + [ ψ ∗ ( Δ L ) ψ ∗ ( Δ op ) g r α + s η log ( 1 + K op − K L ) g r α ] 1 / m if Δ op ≤ Δ L . Note that η , s > 1 and r ∈ ( 0 , 1 ) do not depend on the unknown Δ op .

 Application to Real Data and Simulation Study Bernoulli Random Variables with Dependent, Time-varying Means

Winning Rates of the Cavaliers To illustrate how sequential change detection procedures based on e-detectors work, we revisit the example of the Cleveland Cavaliers, an American professional basketball team introduced in Section 1.2. Instead of using Plus-Minus stats, in this example, we are monitoring the performance of the Cavaliers by keeping track of wins and losses over all the games. Let X 1 , X 2 , … ∈ { 0 , 1 } be the sequence of win indicators during 2010–11 to 2017–18 regular seasons, where X i = 1 if the Cavaliers won game i. Though Figure 2 presents monthly and seasonal averages for the purpose of visualization, we use the underlying binary sequence to build a sequential change detection procedure.

Modeling Winning Probabilities as a Dependent Sequence of Bernoullis To detect a significant improvement of the performance of the Cavaliers, we assume that before an unknown changepoint ν ∈ N ∪ { ∞ }, the conditional average of winning probability given the sample history is less than or equal to p 0 : = 0.49. That is, under any pre-change distribution P we have p n : = E P , ∞ [ X n ∣ F n − 1 ] ≤ p 0 . (For simplicity, F is taken to be the natural filtration of the data.) Thus, the pre-change class of distributions is P : = { ( p 1 , p 2 , … ) : p i ≤ p 0 , ∀ i ≥ 1 } , where we parameterize each distribution P over binary sequences by the sequence of conditional probabilities.

Our objective is to build mixtures of e-SR and e-CUSUM procedures tuned to quickly detect any significantly improved win rate larger than q 0 : = 0.51 after the changepoint. This can be modeled by assuming that after some changepoint ν, the distribution Q is such that E P , ν , Q [ X n ∣ F n − 1 , n > ν ] : = q n ≥ q 0 . Thus, we may think of the post-change class of distributions as being Q : = { ( q 1 , q 2 , … ) : q i ≥ q 0 , ∀ i ≥ 1 } . In particular, this formalization allows for the winning probabilities to fluctuate over time before and after the changepoint (accounting for factors like form, injuries, etc.).

Deriving Exponential Baseline Processes For each λ > 0, define a baseline increment process L ( λ ) : = { L n ( λ ) } n ≥ 1 as L n ( λ ) : = exp { λ ( X n − p 0 ) − ψ B ( λ ) } , where ψ B ( λ ) : = log ( 1 − p 0 + p 0 e λ ) − λ p 0 is the Bernoulli cumulant generating function. Note that each L ( λ ) is a valid baseline increment as it satisfies the inequality (2.6) in Definition 2.8. That is, under any pre-change distribution P, we have E P , ∞ [ L n ( λ ) ∣ F n − 1 ] = E P , ∞ [ exp { λ X n − log ( 1 − p 0 + p 0 e λ ) } ∣ F n − 1 ] = E P , ∞ [ e λ X n ∣ F n − 1 ] 1 − p 0 + p 0 e λ = 1 − p n + p n e λ 1 − p 0 + p 0 e λ ≤ 1 , ∀ λ ≥ 0 .

To derive exponential baseline processes, we first consider a simplified post-change distribution Q where each post-change observation is identically distributed with E 0 , Q [ X 1 ] : = q ≥ q 0 > p 0 . In this case, the optimal choice of λ ≥ 0 given by λ op : = arg max λ ≥ 0 E 0 , Q exp { λ ( X 1 − p 0 ) − ψ B ( λ ) } . Since the baseline increment has the exponential structure, by Proposition 4.2, we have that D ( Q | | P ) : = E 0 , Q log L 1 ( λ op ) = ψ B ∗ ( q − p 0 ) = KL ( q | | p 0 ) , where KL ( q | | p 0 ) is the Kullback-Leibler (KL) divergence of Bernoulli distributions with parameters q and p 0 written as KL ( q | | p 0 ) : = q log q p 0 + ( 1 − q ) log 1 − q 1 − p 0 , for q , p 0 ∈ ( 0 , 1 ). The appearance of the KL divergence in (5.3) is not a coincidence as the baseline increment can be viewed as a re-parametrized likelihood ratios between two Bernoulli processes. However, the simple geometric structure of the baseline increment make it possible to utilize a prior knowledge about the post-change distribution via Algorithm 1 and 2.

For instance, suppose we know upper and lower bounds of conditional means of the post-change distribution as q n ∈ ( q L , q U ) , ∀ n > ν. Let N mSR ∗ and N mCU ∗ be stopping times of mixtures of e-SR and e-CUSUM procedures in Algorithm 1. In this case, derived sequential change detection procedures do not rely on a specific choice of a post-change distribution Q ∈ Q. However, these procedures can still perform almost as well as the one optimized to a specific choice of the post-change distribution within the same range ( q L , q U ). Typically, if the post-change observations are iid samples from a post-change distribution Q with E 0 , Q [ X 1 ] : = q ∈ ( q L , q U ), then by Theorem 4.3, the worst average delays have the following explicit bound: max { J L ( N mSR ∗ ) , J L ( N mCU ∗ ) } ≤ g α KL ( q | | p 0 ) + q ( 1 − q ) [ log ( 1 − p 0 p 0 q 1 − q ) ] 2 [ KL ( q | | p 0 ) ] 2 + 1 . Typically for small α ≪ 1, from Proposition B.2, we can simplify the above upper bound as max { J L ( N mSR ∗ ) , J L ( N mCU ∗ ) } ⪅ inf η > 1 η [ log ( 1 / α ) + log ( 1 + ⌈ log η KL ( q U | | p 0 ) KL ( q L | | p 0 ) ⌉ ) ] KL ( q | | p 0 ) , which matches the rate of the worst average delays, O ( log ( 1 / α ) / KL ( q | | p 0 ) ) of the oracle sequential change detection procedure as α → 0.

Left: Monthly win rates of the Cavaliers from 2010–11 to 2017–18 seasons (the raw data is Bernoulli, which is harder to visualize). Each red line corresponds to the seasonal average. Right: Paths of log e-detectors (SR: red; CUSUM: green). The horizontal line is the threshold (common to both procedures) equal to log ( 1 / α ), ensuring that the ARL is at least 1 / α = 10 3 , larger than the number of games in 12 seasons (82 per season). The e-SR procedure detects a changepoint during the 2014–15 season.

Implementation of Algorithm 1 and Its Results The lower bound q L can be chosen as q 0 = 0.51 since it is the minimum winning rate we consider as a significant improvement from before the changepoint, when the rates are upper bounded by p 0 = 0.49. We can also safely assume that the win rate cannot be too high given the competitiveness of the NBA, so that the improved win rates cannot be larger than 0.9. In our framework, these considerations can be encoded by setting Δ L : = q 0 − p 0 = 0.02 and Δ U : = 0.41 as input parameters of Algorithm 1. As in Section 1.2, we set α : = 10 − 3 to ensure that the ARL is at least 1 / α : = 10 3 , which is more than the total number of games over 12 years of regular seasons. Finally, we set the maximum number of baselines K max : = 1000. In fact, the computeBaseline function of Algorithm 3 returns only 69 baseline processes, and thus the resulting mixtures of e-SR and e-CUSUM procedures of Algorithm 1 can be computed efficiently in an online fashion.

The right plot in Figure 2 presents the log e-detector values using mixtures of e-SR (red) and e-CUSUM (green) procedures. Although there were a few months in which monthly win rates were higher than p 0 , overall log e-detectors remained at a stable level over the first four seasons. However, after the 2014–15 season starts, the log e-detectors increase rapidly and both procedures detect a changepoint during the 2014–15 season, which is the season that marked the return of LeBron James to the Cavaliers.

 Mean-shift Detection in General Bounded Random Variables

Plus-Minus of the Cavaliers Revisited We return to the Cavaliers 2011–2018 example from Section 1.2. Let X ˜ 1 , X ˜ 2 , … be the sequence of Plus-Minus stats from each game. We assume that the average Plus-Minus of the team is less than or equal to μ < : = − 1 before the changepoint (if any), while after the changepoint it is greater than μ > : = 1. Here, the gap | μ > − μ < | between averages of Plus-Minus in pre- and post-changes refers to the degree of improvement we consider as significant.

For convenience, we first normalize the observed sequence. We assume that the absolute value of each Plus-Minus is bounded by 80, meaning that no team beats another by over 80 points (such an extreme game has never happened in NBA history). Accordingly, define the normalized Plus-Minus, X n : = ( X ˜ n + 80 ) / 160 ∈ [ 0 , 1 ] for each n. Then, the pre-change observations have conditional mean at most m : = ( μ < + 80 ) / 160 = 0.494 and the minimum gap to detect is equal to δ : = | μ > − μ < | / 160 = 0.0125.

Modeling Plus-Minus Stats as a Bounded Sequence with Time-varying, Dependent Means After the normalization above, the Plus-Minus stats form sequence of bounded random variables X 1 , X 2 , … on [ 0 , 1 ]. Each observation may have different distribution (due to seasonal effects, injuries, form, etc.), but we assume that all observations before an unknown changepoint ν have a mean less than or equal to a known boundary m ∈ ( 0 , 1 ), when conditioned on the past sample history. That is, under any pre-change distribution P, we have μ n : = E P , ∞ [ X n ∣ F n − 1 ] ≤ m, ∀ n ≥ 1. In other words, we use P : = { P : μ n ≤ m , ∀ n ≥ 1 } , where other characteristics about P (outside of its sequence of conditional means) are irrelevant. But after the changepoint, all observations have (conditional) mean larger than the boundary m with the minimum gap equal to δ. Thus, Q : = { P : μ n ≥ m + δ , ∀ n ≥ 1 } ,

To build an e-SR or e-CUSUM procedure, we need to choose a baseline increment. To derive it, we first consider a simplified setting where both pre- and post-change observations are independently and identically distributed with E P , ∞ [ X ] ≤ m and E 0 , Q [ X ] ≥ m + δ, respectively. In this simplified case, we simply refer P and Q to marginal pre- and post-change distributions and P and Q to their collections. Then, define KL inf ( Q ; m ) : = inf P ∈ P KL ( Q , P ) to be the smallest KL divergence between Q and P. It is known (see, e.g., [9, 10]) that KL inf has the following variational representation: KL inf ( Q , m ) = sup λ ∈ ( 0 , 1 ) E 0 , Q log ( 1 + λ ( X m − 1 ) ) = : D ( Q | | P ) . Accordingly, for each λ ∈ ( 0 , 1 ), define the baseline increment L λ : = { L n } n ≥ 1 as L n λ : = 1 + λ ( X n m − 1 ) , for each n ∈ N. Though the baseline increment above has been derived in the simplified iid setting, it can be checked that L λ is also a valid baseline increment for the general time-varying, dependent means case since it is nonnegative whenever X n , m ∈ [ 0 , 1 ] as assumed in our setup, and for each pre-change distribution P ∈ P, we have E P , ∞ [ L n λ ∣ F n − 1 ] = 1 + λ ( E P , ∞ [ X n ∣ F n − 1 ] m − 1 ) ≤ 1 , where the inequality comes from the condition μ n ≤ m for any pre-change distribution.

Interestingly, the baseline increments in (5.6) correspond to rescaled increments of the capital process used in [48] to design test martingales for confidence sequences of means of bounded random variables. Though the expressions are essentially identical, ours was obtained via a variational representation of the KL divergence between distributions of bounded random variables, while the derivation presented in [48] is based on a betting interpretation of hypothesis testing.

For any Q ∈ Q, let λ op be the optimal choice of λ ∈ [ 0 , 1 ] given by λ op = arg max λ ∈ [ 0 , 1 ] ) E 0 , Q log ( 1 + λ ( X m − 1 ) ) . Unfortunately, it is typically difficult to compute the optimal λ op since it depends on the unknown post-change distribution Q in a complicated way. In this case, we use a sub-exponential lower bound from [5, 12], given by L ˜ n λ : = exp { λ ( X n m − 1 ) − ψ E ( λ ) ( X n m − 1 ) 2 } ≤ 1 + λ ( X n m − 1 ) = L n λ , where ψ E ( λ ) : = − log ( 1 − λ ) − λ for λ ∈ ( 0 , 1 ). For each λ ∈ ( 0 , 1 ), the process L ˜ λ is itself a valid exponential baseline increment with s ( x ) : = x / m − 1 and v ( x ) : = ( x / m − 1 ) 2 .

The lower bound in (5.8) also implies the lower bound KL inf ( Q , m ) ≥ sup λ ∈ [ 0 , 1 ] { λ μ − ψ E ( λ ) σ 2 } = σ 2 ψ E ∗ ( Δ op ) , where ψ E ∗ ( u ) : = u − log ( 1 + u ) is the convex conjugate of ψ E , while μ , σ 2 and Δ op from Section 4 are: μ : = E 0 , Q s ( X ) = E 0 , Q X − m m , σ 2 : = E 0 , Q v ( X ) = E 0 , Q ( X − m ) 2 m 2 , Δ op : = μ σ 2 = m [ E 0 , Q X − m ] E 0 , Q ( X − m ) 2 . Noting that ψ E ∗ ( u ) ≈ u 2 / 2 for small u, we see that for small Δ op ≪ 1, one has KL inf ( Q , m ) ≳ [ E 0 , Q X − m ] 2 2 E 0 , Q ( X − m ) 2 .

Note that the oracle Δ op depends on the unknown post-change distribution only via first and second moments. Therefore, in contrast to the original set of baseline increments { L λ } λ ∈ ( 0 , 1 ) , the exponential baseline increments { L ˜ λ } λ ∈ ( 0 , 1 ) that lower bound them allow us to more easily set a range ( Δ L , Δ U ) to build mixtures of the e-SR and e-CUSUM procedures. For example, if we assume that the post-change distribution has mean at least m + δ for a positive δ then we can upper and lower bound Δ op by Δ L : = m δ ( 1 − m ) 2 ≤ Δ op ≤ m ( 1 − m ) δ 2 = : Δ U .

Now, given Δ L and Δ U , we can use Algorithm 1 to run the mixture of e-SR or e-CUSUM procedure to detect the changepoint based on the exponential baseline increments { L ˜ λ } λ ∈ ( 0 , 1 ) . It is also straightforward to build the corresponding mixtures of e-SR and e-CUSUM procedures for the original baseline increment { L λ } λ ∈ ( 0 , 1 ) which is always more sample-efficient.

Left: E-detectors (SR: red; CUSUM: green) over eight seasons (82 games per season). Right: Logarithm of e-detectors against date (the sharp rise of e-SR is simply due to the log scale). In both plots, horizontal lines are thresholds equal to 1 / α (left) and log ( 1 / α ) (right) controlling the ARL by 1 / α = 10 3 . The e-SR procedure detects a change during the 2014–15 season, while e-CUSUM takes longer (as expected).

Implementation of Algorithm 1 and Its Results Recall that in the plus-minus stats running example, we use pre-change mean m = 0.494 and the minimum gap δ = 0.0125, which bounds Δ op by Δ L : = 0.024 and Δ U : = m ( 1 − m ) δ 2 = 1600. As before, we choose α = 10 − 3 to make the ARL larger than 12 regular seasons and set the maximum number of baselines K max = 1000. Based on these parameters, we can build mixtures of e-SR and e-CUSUM procedures. Though the difference between Δ L and Δ U may seem to be large, the actual number of baselines returned by the function computeBaseline in Algorithm 1 is 190, which is small enough to update the procedure efficiently on the fly.

Figure 3 shows e-detectors (left) and their logarithms (right). The horizontal line corresponds to the detection boundary given by 1 / α (left) and log α − 1 (right). Similar to the winning rate example, the log e-detectors remained stable during the first four regular seasons, although the difference between SR and CUSUM e-detectors is larger than before. After 2014–15 season started, both e-detectors increased rapidly, and the e-SR procedure detects a changepoint during the 2014–15 season, but e-CUSUM detects the changepoint only in the following season (as expected, since both procedures use the same threshold).

 Simulation-based Comparison with Parametric Methods

In the Bernoulli example of Section 5.1, we showed that, in the simple i.i.d Bernoulli setup, our mixtures of e-SR and e-CUSUM procedures match the rate of the worst average delays O ( log ( 1 / α ) / KL ( q | | p 0 ) ) of the oracle sequential change detection procedure as α → 0. In this subsection, we conduct a simulation study to compare the efficiency of our e-SR procedure with the oracle CUSUM procedure with the exact threshold [21, 31], given by N CU ∗ : = inf { n ≥ 1 : max 0 ≤ j < n ∑ i = j n log f p ∗ ( X i ) f p 0 ( X i ) ≥ c α ∗ } , where p ∗ is the true post-change distribution parameter (hence the oracle designation) and c α ∗ is the value of the threshold so that the ARL is exactly 1 / α. It is well known that the oracle CUSUM procedure (with the appropriate choice of the stopping threshold that controls ARL exactly) minimizes the worst average delay.

We also compare our method with a version of the GLR procedure based on the stopping time N GL ∗ : = inf { n ≥ 1 : max 0 ≤ j < n sup p ∑ i = j n log f p ( X i ) f p 0 ( X i ) ≥ c α } = inf { n ≥ 1 : max 0 ≤ j < n ∑ i = j n log f p ˆ j : n ( X i ) f p 0 ( X i ) ≥ c α } , where each p ˆ j : n is the MLE of the post-change parameter and the exact threshold c α is tuned to control ARL exactly at 1 / α (this is typically only possible in such simple parametric settings, either by analytic derivations or simulations). Unlike the oracle CUSUM procedure, the GLR procedure does not have an iterative update rule, as we need to recompute the MLE of the post-change parameter at each time. As a result, its computational cost at time n is O ( n ), which makes an online implementation very costly. In practice, we may want to use a window-limited GLR procedure to overcome the computational challenge. However, in our study, we deploy the GLR procedure to avoid the additional challenge of picking a window size.

Average detection delay for each changepoint ν = 0 , 100 , … , 500 (each experiment has exactly one changepoint at ν). Three of the methods use an exact threshold calculated via simulation (only possible in this simple, parametric example). Only the Oracle CUSUM method knows the post-change distribution. Even though e-SR uses the conservative log ( 1 / α ) threshold, its detection delay is excellent, often even better than the Oracle CUSUM method (which has optimal average worst-case (across ν) delay).

Pre-change false alarm rates of detection methods for each changepoint ν = 0 , 100 , … , 500 (each experiment has exactly one changepoint at ν). Three of the methods use an exact threshold calculated via simulation (only possible in this simple, parametric example). Only the Oracle CUSUM method knows the post-change distribution exactly. e-SR has the smallest false alert ratio (defined in the text).

Simulation Details Throughout this simulation, we draw pre-change observations as iid Bernoulli random variables with p 0 = 0.5, and post-change observations using p 1 = 0.6. For non-oracle methods, we will only assume that the post-change parameter is known to be in the interval [ 0.51 , 0.99 ]. The e-SR and e-CUSUM procedures in Algorithm 1 will use this range to set Δ L : = 0.01 and Δ U : = 0.49. For the GLR procedure, the MLE of the post-change parameter is p ˆ j : n : = min { max { X ¯ j : n , 0.51 } , 0.99 } , where X ¯ j : n is the sample average over last n − j + 1 observations. Our ARL target is equal to 1 / α : = 500, and each simulation is repeated 5000 times to estimate average delays. The time of the changepoint ν varies in { 0 , 100 , … , 500 }. For simplicity, each run will end no later than time n = 1000.

For the oracle CUSUM, GLR, and e-CUSUM procedures (but not e-SR), we use the same simulation setup to find the exact threshold value that controls the ARL exactly 1 / α for each method. For the e-SR procedure, we simply use the universal threshold of log ( 1 / α ) to demonstrate its efficiency. In particular, practitioners are not required to resort to expensive simulation to identify a good threshold value.

Figure 4 shows the average delays of oracle CUSUM, GLR, e-CUSUM and e-SR procedures for each changepoint ν ∈ { 0 , 100 , … , 500 } (each experiment has only one change at ν). The vertical bar on each point represents 95% confidence interval of the average delay. As the theory guarantees, the oracle CUSUM procedure with an exact threshold results in the smallest worst average delay ( 91.3 ± 1.7), while surprisingly the GLR procedure with an exact threshold shows the largest worst average delay ( 123.7 ± 2.7) despite its high computational complexity. The two e-detectors perform reasonably well, and the e-SR detector in particular performs quite favorably overall despite using its conservative log ( 1 / α ) threshold.

As the value of the changepoint approaches the ARL target of 500, the average delays tend to decrease sharply for e-SR procedure, even falling below the oracle CUSUM method. This is plausibly because e-SR sums the underlying e-processes, in contrast to the other CUSUM-style procedures which take the maximum of the underlying e-processes. We may also intuitively expect the oracle CUSUM delay to be relatively flat across ν because it is known to be minimax optimal, in the sense of minimizing the worst case delay, and minimax procedures often have constant “risk profiles”. Proving these fine-grained behaviors is beyond the scope of the current paper.

Figure 5 illustrates the “pre-change false alarm rate”: the fraction of simulation runs in which the detection procedures stopped before the changepoint at ν (if a change occurs at ν = 0, then it is zero by definition). This is not a common metric, since we provably control ARL at the target level. However, it is an interesting metric, so we plot it. As the time of the true changepoint ν becomes closer to the ARL target of 500, false alarm rates increase across all methods. We notice that the e-SR procedure with the log ( 1 / α ) threshold results in the smallest false alarm rates in most cases.

 Discussion Game-theoretic Interpretation of an E-detector

We briefly mention here a game-theoretic interpretation of an e-detector along the lines of the game-theoretic interpretations of martingales and supermartingales as the wealth of a gambler playing a fair game (well known since the time of [42]). We first summarize the game-theoretic interpretation of a P-e-process, as described in [29].

The standard game-theoretic setup of [33] involves three players: a forecaster, a skeptic, and reality. The forecaster claims at the beginning that P is a plausible model for the yet-to-be-observed data; meaning that the observations are in accordance with (or generated by) some P ∈ P. The skeptic plays (in parallel) a family of games indexed by P ∈ P against nature and begins with one dollar in each game. The objective of the skeptic in the P-th game is to sequentially test whether P is a good explanation for the data by betting against P. At each time step, the skeptic places fair bets (relative to P, in the P-th game) about the next outcome. Then nature reveals the next outcome, and the skeptic’s wealth in every game is updated. The magnitude of the skeptic’s wealth in the P-th game is direct evidence against P being a good explanation; the higher the wealth, the more unlikely the data came from P. Thus in each game, the gambler places different bets, but nature’s moves (the outcomes) are identical across all games. The skeptic’s overall evidence against P is measured by their worst wealth across all the games. If this evidence exceeds 1 / α, it means that the skeptic multiplied their initial capital by at least 1 / α in every game, and if we reject P when this happens, Ville’s inequality implies that we have a valid level-α sequential test.

Since our e-detectors are constructed to be cumulative sums of e-processes started at consecutive times, their game-theoretic interpretation builds on the aforementioned one. Informally, the forecaster not only claims that the data sequence follows P from the start, but that this will not change after some amount of time. The skeptic now wishes to detect a change, if one occurs, as soon as possible. To accomplish this task, the skeptic is provided with one extra dollar every day that they invest (using a P e-process) into testing whether the data from that day onwards is still explained well by P. E-detectors use the wealth in all these games (one against each P ∈ P, starting at each time) as a measure of evidence against the forecaster’s claims. The SR e-detector uses the sum (across time) of the minimum wealth (across P at each time), though it could use the amount that this wealth exceeds n, which is the total dollar amount invested up to time n. The CUSUM e-detector uses the max-min wealth; the maximum (across time) of the minimum wealth (across P). These are only two ways of constructing e-detectors, and we leave other constructions to future work.

 Viewing Lorden’s Reduction to Sequential Testing as an E-detector

Lorden [18] proposed a simple method to construct a change detection method with ARL control via a reduction to sequential testing. We describe this below, first defining a sequential test formally.

A sequential test ϕ is a mapping from increasing amounts of data to a sequence of zeros and ones, where a one represents a rejection of the null hypothesis, and a zero means “continue collecting data”. Formally, define the decision at time t as ϕ t : X t → { 0 , 1 }, and let ϕ : = { ϕ t } t ≥ 0 be the collection of such decisions made one at a time based on the first t datapoints, with ϕ 0 = 0 by default. The sequence of tests ϕ is called a level-α sequential test for P if sup P ∈ P P ( ∃ t ≥ 1 : ϕ t ( X 1 , … , X t ) = 1 ) ≤ α , i.e. if the probability of ever falsely rejecting the null is at most α. By convention, if ϕ t = 1, we set ϕ s = 1 for s > t. This is equivalent to requiring that, for each P ∈ P, P ( ϕ τ ( X 1 , … , X τ ) = 1 ) ≤ α for any stopping time τ .

Let ϕ ( j ) denote a sequential test is started at time j; that is, for ϕ ( j ) , the first observation is actually X j (and not X 1 ), but the test itself can depend on the first j − 1 observations (for example, we can choose our betting strategy based the first j − 1 points, even though our betting score will only be evaluated from time j onwards). Note that ϕ ( 1 ) is simply a standard sequential test as defined above. Ideally, these tests are powerful against alternative Q.

Lorden’s change detection procedure is simple and works as follows. At each time t, start a new sequential test ϕ ( t ) , in addition to the ones that are already running. In other words, consider a sequence of level-α sequential tests ϕ ( 1 ) , ϕ ( 2 ) , … , starting at consecutive times. Lorden declares a change if any of those sequential tests rejects the null P: N Lorden : = inf { n ≥ 1 : max 1 ≤ j ≤ n ϕ n ( j ) ( X n − j + 1 , … , X n ) = 1 } . Lorden proved that this method controls the ARL at 1 / α if the data are iid, and if the same test ϕ ( j ) is deployed at each j (i.e. apart from the delayed start, the tests are identical).

We first observe that Lorden’s method is a special case of an e-detector. Indeed, with each level-α sequential test ϕ ≡ { ϕ t } t ≥ 0 , we can associate an e-process Λ t Lorden : = 1 ( ϕ t = 1 ) α = ϕ t α . Note that Λ only takes on two values: 0 , 1 / α, and when it reaches the latter, it stays there. Furthermore, note that E [ Λ τ Lorden ] ≤ 1 for any stopping time τ, which makes it an e-process as claimed.

Last, note that if we form a “Lorden e-detector” using either the SR or CUSUM methods in (2.5), then both e-detectors start out at 0, and either e-detector jumps to level 1 / α if and only if one of the (delayed start) sequential tests rejects the null, and further our e-detector declares a changepoint at exactly the same instant that Lorden’s does. Thus, Lorden’s procedure can be subsumed within our e-detector framework without any loss of generality or performance.

There are two benefits to viewing Lorden’s method as an e-detector. First, we can dispense with both the aforementioned conditions that Lorden required to prove ARL control: the iid assumption, and the condition that the underlying tests ϕ ( j ) are identical across j. Indeed, our main result guarantees ARL control for any e-detector no matter what the underlying e j -processes are, or whether the data are iid or not.

Second, this viewpoint allows us to see why e-detectors could have much smaller detection delay than Lorden’s method (that is, Lorden’s e-detector). In Lorden’s e-detector, there is no sharing of evidence across different e j -processes: each sequential test acts alone without help from the others, and we need a single e j -process to reach 1 / α before we can declare a change. When a general (say SR-type) e-detector crosses 1 / α, the reason it does so will usually be because of a collaborative effect across various e-processes (caused by the nontrivial sum of e j -processes in the definition of the e-detector), none of which have yet individually reached 1 / α. This will happen much sooner than any individual one reaches that threshold, causing an earlier detection than Lorden’s method. In fact, every level-α sequential test in some sense must be based on threshold an e-process at level 1 / α [28] and using our e-detector with those underlying e-processes will be more statistically efficient (shorter delay) than directly using Lorden’s reduction.

For the sake of future reference, we summarize the above observations below in a “generalized Lorden’s lemma”, whose proof follows immediately from the properties of an e-detector and the discussion above. (Generalized Lorden’s Lemma).

Suppose the data initially come from a distribution in the pre-change class P and, if a change occurs, they later come from a distribution in the post-change class Q (note the lack of any iid assumption). For each j, let ϕ ( j ) denote a (one-sided) level-α sequential test for P against Q that is started at time j (but need not be identical or related in any way to any other ϕ ( k ) , for k ≠ j). If we declare a change at the first time when any one of these sequential tests rejects the null, the resulting change detection procedure has ARL at most 1 / α. Further, this generalization of Lorden’s change detector is a special case of an e-detector.

 Future Directions

There remain a whole host of follow-up directions; we mention only a few below. First, our sequential change detection framework can be straightforwardly generalized to the multi-stream setting where we are monitoring a large number of data streams. In the classical parametric setting, minimum or summation of local CUSUM statistics for multi-stream data were proposed and their asymptotic optimality was studied [8, 20]. Since either minimum or scaled summation (average) of e-detectors also forms a valid e-detector, we can apply the framework in this paper to the multi-stream setting seamlessly. It is interesting to investigate how the framework can be even further generalized to structural multi-stream settings [50, 51, 2].

The kernel sequential change detection is an important class of sequential change detection methods [3, 16]. It is an interesting open direction how to instantiate existing kernel-based methods into our general framework to make it possible to analyze kernel sequential change detection algorithms in a nonasymptotic way. As it currently stands, neither framework is more general than the other, because the kernel methods often assume iid data before the changepoint, while we abstain from such strong assumptions.

Last, throughout this paper, we have only focused on detecting whether a changepoint happened or not but have not dealt with inferential questions surrounding when the change occurred. Future work could study how to perform such inference with e-detectors in our nonparametric settings, either online or post-hoc.

 Summary

We have presented a general framework for sequential change detection based on a new concept called e-detectors. The proposed framework is nonparametric as it does not rely on a parametric assumption on the data-generating distribution (though, when such assumptions are made, we recover well-known parametric methods as special cases). Also, the framework comes with nonasymptotic guarantees, since every component of the framework can be chosen and analyzed explicitly without any asymptotic approximations. By introducing additional structures such as baseline increments and exponential e-detectors on the top of the general framework, we can construct computationally and statistically efficient online algorithms that have explicit upper bounds on worst average delays. Finally, through examples involving Bernoulli and bounded random variables, we explained how one can apply the presented framework in practical settings, with NBA data serving as a running case study.

 Main Proofs Proofs for Statements in Section <xref rid="j_nejsds51_s_004">2</xref>

For any given α ∈ ( 0 , 1 ) let N ∗ be the stopping time defined by N ∗ : = inf { n ≥ 1 : M n ≥ 1 / α } . For any pre-change distribution P ∈ P, we may assume N ∗ < ∞ with probability one under P, without loss of generality. (If not, we have E P , ∞ N ∗ = ∞, which immediately proves the claim.) Then, from the definition of an e-detector, we have E P , ∞ [ M N ∗ ] ≤ E P , ∞ N ∗ , which implies that, for any P ∈ P, E P , ∞ N ∗ ≥ E P , ∞ [ M N ∗ ] = E P , ∞ [ M N ∗ 1 ( N ∗ < ∞ ) ] ≥ E P , ∞ [ 1 α 1 ( N ∗ < ∞ ) ] (by definition of N ∗ ) = 1 α P P , ∞ ( N ∗ < ∞ ) = 1 α , as desired.  □

We first prove (2.12). Note that N SR ∗ = inf { n ≥ 1 : ∑ j = 1 n ∏ i = j n L i ≥ 1 / α } ≤ min j ≥ 1 inf { n ≥ j : ∏ i = j n L i ≥ 1 / α } : = min j ≥ 1 N j . For each fixed P ∈ P and Q ∈ Q, since N j ⫫ F ν for all j > ν + m, we have E P , ν , Q [ [ N SR ∗ − ν ] + ∣ F ν ] ≤ min j ≥ 1 E P , ν , Q [ [ N j − ν ] + ∣ F ν ] ≤ min j > ν + m E P , ν , Q [ N j − ν ∣ F ν ] (since N j ≥ j > ν + m ) = min j > ν + m E P , ν , Q [ N j ] − ν (since N j ⫫ F ν , ∀ j > ν + m ) = min j > ν + m E P , ν , Q [ inf { n ≥ j : ∏ i = j n L i ≥ 1 / α } ] − ν = min j − ν > m E 0 , Q [ inf { n ≥ j − ν : ∏ i = j − ν n L i ≥ 1 / α } ] = min j ′ > m E 0 , Q [ inf { n ≥ j ′ : ∏ i = j ′ n L i ≥ 1 / α } ] (by j ′ : = j − ν ) = min j ′ > m E 0 , Q [ inf { n ≥ 1 : ∏ i = 1 n L i ≥ 1 / α } ] + j ′ − 1 = min j ′ > m E 0 , Q N 1 / α + j ′ − 1 = E 0 , Q N 1 / α + m , where the third equality comes from the fact that the distribution of X j − m , X j − m + 1 , … under P P , ν , Q is equal to the one of X j − ν − m , X j − ν − m + 1 , … under P 0 , Q provided that j − m > ν, and the fifth equality is based on the strong stationarity of the post-change observations. Since the last term does not depend on P, ν or F ν , we obtain the claimed result J P ( N SR ∗ ) ≤ J L ( N SR ∗ ) = sup P ∈ P , ν ≥ 0 esssup E P , ν , Q [ [ N SR ∗ − ν ] + ∣ F ν ] ≤ E 0 , Q N 1 / α + m , as desired.

To prove (2.13), first note that N CU ∗ = inf { n ≥ 1 : max j ∈ [ n ] ∏ i = j n L i ≥ c α } ≤ min j ≥ 1 inf { n ≥ j : ∏ i = j n L i ≥ c α } : = min j ≥ 1 N j . The remaining part of the proof for (2.13) is followed by the same argument for (2.12).

Finally, to prove (2.14), it is enough to show the following inequality holds: E 0 , Q N ( g ) ≤ g E 0 , Q log L 1 + V 0 , Q log L 1 [ E 0 , Q log L 1 ] 2 + 1 , ∀ g > 0 , where N ( g ) : = inf { n ≥ 1 : ∑ i = 1 n log L i ≥ g } for each g > 0. The proof of the above upper bound (A.3) is based on the Lorden’s inequality [17] which can be stated as follows: (Lorden’s inequality [<xref ref-type="bibr" rid="j_nejsds51_ref_017">17</xref>]).

Suppose X 1 , X 2 , … are i.i.d. samples with E X 1 = μ > 0 and E X 1 2 < ∞. For each g > 0, set N ( g ) : = inf { n : S n : = ∑ i = 1 n X i ≥ g } and R g : = S N ( g ) − g. Then, the following inequality holds: sup g > 0 E [ R g ] ≤ E X 1 2 μ = μ + σ 2 μ , where σ 2 = V X 1 .

Now, to prove the upper bound (A.3), fix a constant g > 0. Since E 0 , Q log L 1 > 0, we have E 0 , Q N ( g ) < ∞. Therefore, by Wald’s equation, E 0 , Q log L 1 E 0 , Q N ( g ) = E 0 , Q [ ∑ i = 1 N ( g ) log L i ] . For each g > 0, set R g : = ∑ i = 1 N ( g ) log L i − g. Then, from the Lorden’s inequality, we have E 0 , Q log L 1 E 0 , Q N ( g ) = E 0 , Q [ ∑ i = 1 N ( g ) log L i ] ≤ c + sup g > 0 E 0 , Q [ R g ] ≤ g + V 0 , Q log L 1 E 0 , Q log L 1 + E 0 , Q log L 1 , where the first inequality comes from the definition of N ( g ). By multiplying 1 / E 0 , Q log L 1 on both sides of the inequality (A.6), we have the claimed upper bound, completing the proof.  □

 Proofs for Statements in Section <xref rid="j_nejsds51_s_012">3</xref>

To see adaptive SR and CUSUM e-detectors are actually valid e-detectors, first note that M n aCU ≤ M n aSR for each n ∈ N. Therefore, it is enough to show that E P , ∞ M τ aSR ≤ E P , ∞ τ for any stopping time τ and pre-change distribution P ∈ P. If P P , ∞ ( τ = ∞ ) > 0 then the above inequality holds trivially. Otherwise if P P , ∞ ( τ = ∞ ) = 0, we have that E P , ∞ M τ aSR = E P , ∞ ∑ k = 1 K ( τ ) ω k ∑ j = K − 1 ( k ) τ γ j ∏ i = j τ L i ( k ) : = E P , ∞ ∑ k = 1 K ( τ ) ω k ∑ j = K − 1 ( k ) τ γ j Λ τ ( j ) ( k ) = E P , ∞ ∑ j = 1 τ ∑ k = 1 K ( j ) γ j ω k Λ τ ( j ) ( k ) = ∑ j = 1 ∞ ∑ k = 1 K ( j ) γ j ω k E P , ∞ Λ τ ( j ) ( k ) 1 ( j ≤ τ ) = ∑ j = 1 ∞ ∑ k = 1 K ( j ) γ j ω k E P , ∞ 1 ( j ≤ τ ) E P , ∞ [ Λ τ ( j ) ( k ) ∣ F j − 1 ] ≤ ∑ j = 1 ∞ ∑ k = 1 K ( j ) γ j ω k E P , ∞ 1 ( j ≤ τ ) = ∑ j = 1 ∞ E P , ∞ 1 ( j ≤ τ ) = E P , ∞ τ , as desired. Above, the sole inequality follows since each Λ ( j ) ( k ) is an e j -process, and the following equality invokes the definition of γ j for each j.  □

We first prove the upper bound for the adaptive e-SR procedure in (3.17). Note that N aSR ∗ = inf { n ≥ 1 : ∑ k = 1 K ( n ) ω k ∑ j = K − 1 ( k ) n γ j ∏ i = j n L i ( k ) ≥ 1 / α } ≤ inf { n ≥ 1 : ∑ k = 1 K ( n ) ω k ∑ j = K − 1 ( k ) n ∏ i = j n L i ( k ) ≥ 1 / α } = inf { n ≥ 1 : ∑ j = 1 n ∑ k = 1 K ( j ) ω k ∏ i = j n L i ( k ) ≥ 1 / α } ≤ min j ≥ 1 inf { n ≥ j : ∑ k = 1 K ( j ) ω k ∏ i = j n L i ( k ) ≥ 1 / α } : = min j ≥ 1 N j , where the first inequality follows because γ j ≥ 1 for each j. For each fixed P ∈ P and Q ∈ Q, since N j ⫫ F ν for all j > ν + m, we have E P , ν , Q [ [ N aSR ∗ − ν ] + ∣ F ν ] ≤ min j > ν + m E P , ν , Q [ N j − ν ∣ F ν ] (since N j ≥ j > ν + m ) = min j > ν + m E P , ν , Q [ inf { n ≥ j : ∑ k = 1 K ( j ) ω k ∏ i = j n L i ( k ) ≥ 1 / α } ] − ν (since N j ⫫ F ν , ∀ j > ν + m ) = min j − ν > m E 0 , Q [ inf { n ≥ j − ν : ∑ k = 1 K ( j ) ω k ∏ i = j − ν n L i ( k ) ≥ 1 / α } ] = min j ′ > m E 0 , Q [ inf { n ≥ j ′ : ∑ k = 1 K ( j ′ + ν ) ω k ∏ i = j ′ n L i ( k ) ≥ 1 / α } ] (by setting j ′ : = j − ν ) = min j ′ > m E 0 , Q [ inf { n ≥ 1 : ∑ k = 1 K ( j ′ + ν ) ω k ∏ i = 1 n L i ( k ) ≥ 1 / α } ] + j ′ − 1 ≤ min j ′ > m E 0 , Q [ inf { n ≥ 1 : ∑ k = 1 K ( j ′ ) ω k ∏ i = 1 n L i ( k ) ≥ 1 / α } ] + j ′ − 1 = min j ′ > m E 0 , Q N 1 / α ( j ′ ) + j ′ − 1 , where the second equality comes from the fact that the distribution of X j − m , X j − m + 1 , … under P P , ν , Q is equal to the one of X j − ν − m , X j − ν − m + 1 , … under P 0 , Q provided by j − m > ν and the forth equality is based on the strong stationarity of post-change observations. Since the last term does not depend on P, ν nor F ν , we have the claimed result: J P ( N aSR ∗ ) ≤ J L ( N aSR ∗ ) = sup P ∈ P , ν ≥ 0 esssup E P , ν , Q [ [ N aSR ∗ − ν ] + ∣ F ν ] ≤ min j > m [ E 0 , Q N 1 / α ( j ) + j − 1 ] , as desired.

To prove the adaptive e-CUSUM procedure case in (3.18), first note that N aCU ∗ = inf { n ≥ 1 : ∑ k = 1 K ( n ) ω k max K − 1 ( k ) ≤ j ≤ n γ j ∏ i = j n L i ( k ) ≥ c α } ≤ inf { n ≥ 1 : max j ≤ n ∑ k = 1 K ( j ) ω k ∏ i = j n L i ( k ) ≥ c α } ≤ min j ≥ 1 inf { n ≥ j : ∑ k = 1 K ( j ) ω k ∏ i = j n L i ( k ) ≥ c α } : = min j ≥ 1 N j , where the first inequality comes from γ j ≥ 1 with the fact K − 1 ( k ) ≤ j for any k ≤ K ( j ). The remaining part of the proof of (3.18) follows the same argument used to obtain (3.17).  □

 Remaining Proofs Proofs for Statements in Section <xref rid="j_nejsds51_s_015">4</xref>

To simplify notation, we drop the subscripts { 0 , Q } and 1. The claim reduces to proving that E log L ( λ op ) = ψ ∗ ( Δ op ) σ 2 , where log L ( λ op ) : = λ op s ( X ) − ψ ( λ op ) v ( X ) and Δ op = ∇ ψ ( λ op ). To prove the equality, first note that from the definition of the convex conjugate ψ ∗ of ψ, we have ψ ∗ ( Δ op ) : = sup λ ∈ Π { λ Δ op − ψ ( λ ) } = λ op Δ op − ψ ( λ op ) , which implies E log L ( λ op ) = λ op μ − ψ ( λ op ) σ 2 = σ 2 [ λ op Δ op − ψ ( λ op ) ] = σ 2 ψ ∗ ( Δ op ) , as desired.  □

We first recall the definition of the stopping time in (4.8): N ¯ g : = inf { n ≥ 1 : sup λ ∈ ( λ L , λ U ) ∑ i = 1 n log L i ( λ ) ≥ g } , g > 0 . Then, the same argument used in the proof Proposition 2.14 immediately implies that, i the post-change observations are i.i.d. from Q, then, for any g > 0, E 0 , Q N ¯ g ≤ g D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + 1 .

The claim of the theorem follows from Lemma B.1, whose statement and proof are given below.  □

Let N 1 / α and N c α be stopping times where the underlying mixing weights { ω k } and parameters of baseline increments { λ k } are chosen via Algorithm 3. Let N ¯ g α be the stopping time defined in (4.8) with the threshold given by Algorithm 3. Then, for any stream of observations X 1 , X 2 , … , N c α ≤ N 1 / α ≤ N ¯ g α , deterministically, provided that 1 < c α < 1 / α.

Throughout this proof, we set D L : = ψ ∗ ( Δ L ) < ψ ∗ ( Δ U ) = : D U . The first inequality N c α ≤ N 1 / α follows directly from the definition of the stopping time in (3.22) along with the condition that c α ≤ 1 / α. To prove the second inequality N 1 / α ≤ N ¯ g α , we will exploit on general geometric construction introduced in [36] to analyze the performance of sequential generalized likelihood ratio tests. To that effect, set μ ˆ n : = S n / V n . Then, for each fixed λ > 0 such that Δ = ∇ ψ ( λ ), Proposition 4.2 and the identity λ = ∇ ψ ∗ ( Δ ) imply that the function μ ˆ n ↦ V n − 1 ∑ i = 1 n log L i ( λ ) = λ μ ˆ n − ψ ( λ ) = λ ( μ ˆ n − Δ ) + ψ ∗ ( Δ ) is a mapping of μ ˆ n into the tangent line of the function z ↦ ψ ∗ ( z ) at z = Δ. Next, set V U : = g α / D U and V L : = g α / D L , and define the set R : = { ( z , y ) ∈ [ 0 , ∞ ) 2 : y ≤ ψ ∗ ( z ) } . Then, the stopping event of N ¯ g α can be expressed as { ∃ n ≥ 1 : sup λ ∈ ( λ L , λ U ) ∑ i = 1 n log L i ( λ ) ≥ g α } = { ∃ n ≥ 1 : V n < V U , ( μ ˆ n , g α V n ) ∈ H ( Δ U ) } ∪ { ∃ n ≥ 1 : V n ≥ V U , ( μ ˆ n , g α V n ) ∈ R ∖ H ( Δ L ) } = { ∃ n ≥ 1 : V n < V U , ∑ i = 1 n log L i ( λ U ) ≥ g α } ∪ { ∃ n ≥ 1 : V n ≥ V U , sup λ > λ L ∑ i = 1 n log L i ( λ ) ≥ g α } , H ( Δ U ) and H ( Δ L ) are half spaces contained in and tangent to R at ( Δ U , g α / V U ) and ( Δ L , g α / V L ), respectively. See Figure 6 for an illustration of the stopping event of N ¯ g α .

Illustration of the stopping event of N ¯ g α defined in (4.8), and related regions H ( Δ U ) , H ( Δ L ) and R. The stopping time N ¯ g α is the first time when ( μ ˆ n , g α / V n ) is located in one of the colored areas.

Note that the first decomposition part { ∃ n ≥ 1 : V n < V U , ∑ i = 1 n log L i ( λ U ) ≥ g α } is nonempty only if V U > min x v ( x ) : = v min , which is equivalent to g α > v min D U . For the second part, a straightforward extension of Lemma 1 in the appendix of [36] implies that, for any fixed η > 1, the second part can be further decomposed by sets of simple events as follows: { ∃ n ≥ 1 : V n ≥ V U , sup λ > λ L ∑ i = 1 n log L i ( λ ) ≥ g α } ⊂ ⋃ k = 1 K ( η ) − 1 { ∃ n ≥ 1 : V n ∈ [ V U η k − 1 , V U η k ) , ∑ i = 1 n log L i ( λ k ) ≥ g α / η } ∪ { ∃ n ≥ 1 : V n ≥ V U η K ( η ) − 1 , ∑ i = 1 n log L i ( λ K ( η ) ) ≥ g α / η } , where K ( η ) is a positive integer defined by K ( η ) : = ⌈ log η ( D U D L ) ⌉ , and, for k = 1 , … , K ( η ) − 1, λ k is given by λ k : = ψ ∗ ( Δ k ), with Δ k the solution with respect to z > 0 of the equation ψ ∗ ( z ) = D U η k , while Δ K ( η ) : = Δ L . It can be checked that λ U : = λ 0 > λ 1 > λ 2 > ⋯ > λ K ( η ) = λ L . Decomposing the stopping event of N ¯ g α , we can lower bound the stopping time N ¯ g α for any η > 1 as: N ¯ g α = inf { n ≥ 1 : sup λ ∈ ( λ L , λ U ) ∑ i = 1 n log L i ( λ ) ≥ g α } ≥ inf { n ≥ 1 : e − g α ∏ i = 1 n L i ( λ 0 ) 1 ( g α > v min D U ) + ∑ k = 1 K ( η ) e − g α / η ∏ i = 1 n L i ( λ k ) ≥ 1 } = inf { n ≥ 1 : α − 1 e − g α 1 ( g α > v min D U ) × exp { λ 0 S n − ψ ( λ 0 ) V n } + ∑ k = 1 K ( η ) α − 1 e − g α / η exp { λ k S n − ψ ( λ k ) V n } ≥ 1 / α : = inf { n ≥ 1 : ∑ k = 0 K ( η ) ω k ( η ) exp { λ k S n − ψ ( λ k ) V n } ≥ 1 / α } : = N ( η ) , where ω 0 ( η ) : = α − 1 e − g α 1 ( g α > v min D U ) and ω k ( η ) : = α − 1 e − g α / η for each η > 1 and k = 1 , … , K ( η ). As a quick remark, note that ω 0 ( η ) does not depend on η and each ω k ( η ) in fact does not depend on the index k but we use this notation just for consistency.

Pseudo-code of computeBaseline function.

Finally, set η α : = ( D U D L ) 1 / K α where K α is the integer defined in (B.11) of Algorithm 3. It can be easily checked that K ( η α ) = K α . Therefore, once we choose η = η α then, from the definition of the mixing weights in Algorithm 3, we conclude that N ( η α ) ≥ N 1 / α provided that W : = e − g α 1 ( g α > v min D U ) + K α e − g α / η α ≤ α , where the constant g α is the constant used in Algorithm 3 and given by g α : = inf { g > log ( 1 / α ) : e − g 1 ( g > v min D U ) + min k ∈ [ K max ] k exp { − g ( D U D L ) − 1 / k } ≤ α } . By the definition of K α , we can immediately check that the inequality (B.13) holds, which proves the claimed inequality N 1 / α ≤ N ( η α ) ≤ N ¯ g α , as desired.  □

We conclude this section with a formal proof of the validity of the upper bound in (4.7).

Let K max be a large enough integer such that K α = arg min k ∈ { 1 , … , K max } k exp { − g α ( D U D L ) − 1 / k } = arg min k ∈ N k exp { − g α ( D U D L ) − 1 / k } . Then, the quantity g α Algorithm 1 specified in (B.10) is such that g α < inf η > 1 η [ log ( 1 / α ) + log ( 1 + ⌈ log η ψ ∗ ( Δ U ) ψ ∗ ( Δ L ) ⌉ ) ] .

Once K max is large enough to satisfy (B.15) then the constant g α can be written as g α = inf { g > log ( 1 / α ) : e − g 1 ( g > v min D U ) + inf η > 1 ⌈ log η ( D U D L ) ⌉ e − g / η ≤ α } , because the following equality holds for each g > 0 and D U > D L . min k ∈ N k exp { − g ( D U D L ) − 1 / k } = inf η > 1 ⌈ log η ( D U D L ) ⌉ e − g / η . Finally, to prove the claimed bound in (B.16), first note that from (B.17), we have g α ≤ g ( η ) where g ( η ) > 0 is given by g ( η ) : = inf { g > log ( 1 / α ) : e − g + ⌈ log η ( D U D L ) ⌉ e − g / η ≤ α } ≤ η [ log ( 1 / α ) + log ( 1 + ⌈ log η D U D L ⌉ ) ] , for each η > 1. By taking infimum over η > 1, we get the upper bound in (B.16), as desired.  □

The proof of Lemma 4.5 is similar to the one of Lemma B.1 except the previous threshold g α being replaced with g ( V 0 η K ( j ) ). Also note that, in this proof, the terms Δ 0 and V 0 play a similar role of Δ U and V U in the previous proof of Lemma B.1.

As same as the Lemma B.1 case, set S n : = ∑ i = 1 n s ( X i ), V n : = ∑ i = 1 n v ( X i ), and μ ˆ n : = S n / V n . Then, for each fixed λ > 0 with Δ = ∇ ψ ( λ ), the function μ ˆ n ↦ V n − 1 ∑ i = 1 n log L i ( λ ) = λ μ ˆ n − ψ ( λ ) = λ ( μ ˆ n − Δ ) + ψ ∗ ( Δ ) is a mapping of μ ˆ n to the tangent line of the function z ↦ ψ ∗ ( z ) at z = Δ. Now, since the boundary function g is non-decreasing, the stopping event of N ¯ g ( j ) can be bounded as { ∃ n ≥ 1 : sup λ ∈ ( λ K ( j ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V 0 η K ( j ) ) } ⊂ { ∃ n ≥ 1 : V n < V 0 , ∑ i = 1 n log L i ( λ 0 ) ≥ g ( V 0 ) } ∪ { ∃ n ≥ 1 : V n ∈ [ V 0 , V 0 η K ( j ) ) , sup λ ∈ ( λ K ( j ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V n ) } ∪ { ∃ n ≥ 1 : V n ≥ V 0 η K ( j ) , ∑ i = 1 n log L i ( λ K ( j ) ) ≥ g ( V 0 η K ( j ) ) } = { ∃ n ≥ 1 : V n < V 0 , ( μ ˆ n , g ( V 0 ) V n ) ∈ H ( Δ 0 ) } ∪ { ∃ n ≥ 1 : V n ∈ [ V 0 , V 0 η K ( j ) ) , ( μ ˆ n , g ( V n ) V n ) ∈ R } ∪ { ∃ n ≥ 1 : V n ≥ V 0 η K ( j ) , ( μ ˆ n , g ( V 0 η K ( j ) ) V n ) ∈ H ( Δ K ( j ) ) } = { ∃ n ≥ 1 : V n < V 0 , ∑ i = 1 n log L i ( λ 0 ) ≥ g ( V 0 ) } ∪ { ∃ n ≥ 1 : V n ∈ [ V 0 , V 0 η K ( j ) ) , sup λ ∈ ( λ K ( j ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V n ) } ∪ { ∃ n ≥ 1 : V n ≥ V 0 η K ( j ) , ∑ i = 1 n log L i ( λ K ( j ) ) ≥ g ( V 0 η K ( j ) ) } , where V 0 : = inf { t ≥ 1 : D 0 ≥ g ( t ) / t } and the set R is defined by R : = { ( z , y ) ∈ [ 0 , ∞ ) 2 : y ≤ ψ ∗ ( z ) } , and H ( Δ 0 ) and H ( Δ K ( j ) ) are half spaces contained in and tangent to R at ( Δ 0 , g ( V 0 ) V 0 ) and ( Δ L , g ( V 0 η K ( j ) ) V 0 η K ( j ) ), respectively. See Figure 7 for an illustration of the upper bound of the stopping event of N ¯ g ( j ).

Illustration of the stopping event of N ¯ g ( j ) defined in (4.14), and related regions H ( Δ 0 ) , H ( Δ K ( j ) ) and R. The stopping time N ¯ g ( j ) is the first time when ( μ ˆ n , g ( V ‾ n ( j ) ) / V n ) is located in one of the colored areas, where V ‾ n ( j ) : = max { V 0 , min { V n , V 0 η K ( j ) } }.

Note that the first decomposition part { ∃ n ≥ 1 : V n < V 0 , ∑ i = 1 n log L i ( λ 0 ) ≥ g ( V 0 ) } is nonempty only if V 0 > v min , which is equivalent to g ( V 0 ) > v min D 0 . For the second part, a straightforward extension of Lemma 1 in the appendix of [36] implies that, for any fixed η > 1, the second part can be further decomposed by sets of simple events as follows: { ∃ n ≥ 1 : V n ∈ [ V 0 , V 0 η K ( j ) ) , sup λ ∈ ( λ K ( j ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V n ) } ⊂ ⋃ k = 1 K ( j ) { ∃ n ≥ 1 : V n ∈ [ V 0 η k − 1 , V 0 η k ) , ∑ i = 1 n log L i ( λ k ) ≥ g ( V 0 η k ) / η } , where each λ k is given by λ k : = ψ ∗ ( Δ k ) and each Δ k is the solution of the equation (4.12) for k = 1 , … , K ( j ).

From this decomposition of the stopping event of N ¯ g ( j ), for any fixed η > 1, we can lower bound the stopping time N ¯ g ( j ) as follows: N ¯ g ( j ) = inf { n ≥ 1 : sup λ ∈ ( λ K ( j ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V 0 η K ( j ) ) } ≥ inf { n ≥ 1 : e − g ( V 0 ) ∏ i = 1 n L i ( λ 0 ) 1 ( g ( V 0 ) > v min D 0 ) + ∑ k = 1 K ( j ) e − g ( V 0 η k ) / η ∏ i = 1 n L i ( λ k ) ≥ 1 } = inf { n ≥ 1 : ∑ k = 0 K ( j ) ω k ∏ i = 1 n L i ( λ k ) ≥ 1 / α } = N 1 / α ( j ) , which proves the claimed inequality.  □

From Theorem 3.4 and Lemma 4.5, worst average delays of adaptive e-SR and e-CUSUM procedures are upper bounded by max ν ≥ 0 min j ≥ 1 [ E 0 , Q N G ( j + ν ) + j − 1 ] where N ¯ g ( j ) is a stopping time defined by N ¯ g ( j ) : = inf { n ≥ 1 : sup λ ∈ ( λ K ( j ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V 0 η K ( j ) ) } , for each j ≥ 1. Now, let us first consider the case λ op ≥ λ 0 . In this case, we use the following simple upper bound: min j ≥ 1 [ E 0 , Q N G ( j ) + j − 1 ] ≤ E 0 , Q N G ( 1 ) = E 0 , Q inf { n ≥ 1 : sup λ ∈ ( λ K ( 1 ) , λ 0 ) ∑ i = 1 n log L i ( λ ) ≥ g ( V 0 η K ( 1 ) ) } ≤ E 0 , Q inf { n ≥ 1 : ∑ i = 1 n log L i ( λ 0 ) ≥ g r α } . Since E 0 , Q log L 1 ( λ 0 ) = σ 2 ( λ 0 Δ op − ψ ( λ 0 ) ) ≥ σ 2 ( λ 0 Δ 0 − ψ ( λ 0 ) ) = σ 2 ψ ∗ ( Δ 0 ), by the same argument of Proposition 2.14, the last term above can be further upper bounded by g r α D ( Q | | P ) ψ ∗ ( Δ op ) ψ ∗ ( Δ 0 ) + V 0 , Q [ log L 1 ( λ 0 ) ] [ D ( Q | | P ) ] 2 [ ψ ∗ ( Δ op ) ψ ∗ ( Δ 0 ) ] 2 + 1 . Since we are in the case λ op ≥ λ 0 , we have ψ ∗ ( Δ op ) ψ ∗ ( Δ 0 ) ≥ 1, which can be understood as a measure of inefficiency due to the misspecified upper bound of the oracle λ op .

Now, consider the case λ op < λ 0 where we correctly specified the upper bound. In this case, let j op be the smallest integer satisfying λ K ( j op ) : = λ K op < λ op < λ 0 . Then, we can further upper bound the worst average delays by min j ≥ 1 [ E 0 , Q N G ( j ) + j − 1 ] ≤ E 0 , Q inf { n ≥ 1 : ∑ i = 1 n log L i ( λ op ) ≥ g ( V 0 η K op ) } + j op − 1 : = E 0 , Q N op + j op − 1 . By Equation (4.9), we have the following intermediate upper bound on the worst average delays, E 0 , Q N op + j op − 1 ≤ g ( V 0 η K op ) D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + j op . Note that if j op = 1 then λ K ( 1 ) = λ L < λ op . Thus, in this case, we also correctly specified the lower bound, and the above bound is reduced to the same upper bound on the worst average delays in Theorem 4.3 of the well-separation case except the ARL parameter α being replaced by r α.

Pseudo-code of computeThreshold function.

Finally, to get an explicit upper bound on j op for the case j op > 1, fist note that, from the definition of λ K ( j op − 1 ) with the fact λ op < λ 1 ( j op − 1 ) ⇔ Δ op < Δ K ( j op − 1 ) , we have g ( V 0 η K ( j op − 1 ) ) V 0 η K ( j op − 1 ) = ψ ∗ ( Δ K ( j op − 1 ) ) > ψ ∗ ( Δ op ) . Also, the condition K ( j ) ≥ K L + m log η j implies j ≤ [ η − K L V 0 V 0 η K ( j ) ] 1 / m , for each j ≥ 1. By combining two inequalities above, we have j op − 1 ≤ [ η − K L V 0 V 0 η K ( j op − 1 ) ] 1 / m < [ 1 V 0 η K L g ( V 0 η K ( j op − 1 ) ) ψ ∗ ( Δ op ) ] 1 / m ≤ [ ψ ∗ ( Δ L ) ψ ∗ ( Δ op ) g ( V 0 η K op ) g r α ] 1 / m . In sum, by combining all bounds above, we have min j ≥ 1 [ E 0 , Q N G ( j ) + j − 1 ] ≤ g r α D ( Q | | P ) ψ ∗ ( Δ op ) ψ ∗ ( Δ 0 ) + V 0 , Q [ log L 1 ( λ 0 ) ] [ D ( Q | | P ) ] 2 [ ψ ∗ ( Δ op ) ψ ∗ ( Δ 0 ) ] 2 + 1 if λ op ≥ λ 0 g r α D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + 1 if λ op ∈ ( λ L , λ 0 ) g r α + s η log ( 1 + K op − K L ) D ( Q | | P ) + V 0 , Q [ log L 1 ( λ op ) ] [ D ( Q | | P ) ] 2 + [ ψ ∗ ( Δ L ) ψ ∗ ( Δ op ) g r α + s η log ( 1 + K op − K L ) g r α ] 1 / m if λ op ≤ λ L , as desired.  □

 An Explicit Way to Compute the Threshold in Algorithm <xref rid="j_nejsds51_fig_009">3</xref>

The following pseudo-code describes how to compute the threshold g α defined by g α : = inf { g > log ( 1 / α ) : e − g 1 ( g > v min D U ) + min k ∈ [ K max ] k exp { − g ( D U D L ) − 1 / k } ≤ α } .

References Barnard, G. A. (1959). Control charts and stochastic processes. Journal of the Royal Statistical Society: Series B (Methodological) 21(2) 239–257. Chen, Y., Wang, T. and Samworth, R. J. (2022). High-dimensional, multiscale online changepoint detection. Journal of the Royal Statistical Society Series B: Statistical Methodology 84(1) 234–266. https://doi.org/10.1111/rssb.12447. MR4400396 Desobry, F., Davy, M. and Doncarli, C. (2005). An online kernel change detection algorithm. IEEE Transactions on Signal Processing 53(8) 2961–2974. https://doi.org/10.1109/TSP.2005.851098. MR2169647 Durrett, R. (2019) Probability: Theory and Examples 49. Cambridge University Press. https://doi.org/10.1017/9781108591034. MR3930614 Fan, X., Grama, I. and Liu, Q. (2015). Exponential inequalities for martingales with applications. Electronic Journal of Probability 20, 1. https://doi.org/10.1214/EJP.v20-3496. MR3311214 Gangrade, A., Rinaldo, A. and Ramdas, A. (2023). A Sequential Test for Log-Concavity. arXiv preprint arXiv:2301.03542. Grünwald, P., de Heide, R. and Koolen, W. (2023). Safe testing. Journal of the Royal Statistical Society, Series B. 195 47–63. https://doi.org/10.1016/j.jspi.2017.09.014. MR3760837 Hadjiliadis, O., Zhang, H. and Poor, H. V. (2009). One shot schemes for decentralized quickest change detection. IEEE Transactions on Information Theory 55(7) 3346–3359. https://doi.org/10.1109/TIT.2009.2021311. MR2598025 Honda, J. and Takemura, A. (2010). An Asymptotically Optimal Bandit Algorithm for Bounded Support Models. In Conference on Learning Theory 67–79. Citeseer. Honda, J. and Takemura, A. (2015). Non-asymptotic analysis of a new bandit algorithm for semi-bounded rewards. Journal of Machine Learning Research 16 3721–3756. MR3450551 Howard, S. R., Ramdas, A., McAuliffe, J. and Sekhon, J. (2020). Time-uniform Chernoff bounds via nonnegative supermartingales. Probability Surveys 17 257–317. https://doi.org/10.1214/18-PS321. MR4100718 Howard, S. R., Ramdas, A., McAuliffe, J. and Sekhon, J. (2021). Time-uniform, nonparametric, nonasymptotic confidence sequences. The Annals of Statistics 49(2) 1055–1080. https://doi.org/10.1214/20-aos1991. MR4255119 Lai, T. L. (1995). Sequential changepoint detection in quality control and dynamical systems. Journal of the Royal Statistical Society: Series B (Methodological) 57(4) 613–644. MR1354072 Lai, T. L. (1998). Information bounds and quick detection of parameter changes in stochastic systems. IEEE Transactions on Information Theory 44(7) 2917–2929. https://doi.org/10.1109/18.737522. MR1672051 Lai, T. L. (2001). Sequential analysis: some classical problems and new challenges. Statistica Sinica 11(2) 303–351. MR1844531 Li, S., Xie, Y., Dai, H. and Song, L. (2015). M-Statistic for Kernel Change-Point Detection. In Advances in Neural Information Processing Systems 28. Lorden, G. (1970). On excess over the boundary. The Annals of Mathematical Statistics 41 520–527. https://doi.org/10.1214/aoms/1177697092. MR0254981 Lorden, G. (1971). Procedures for reacting to a change in distribution. The Annals of Mathematical Statistics 42 1897–1908. https://doi.org/10.1214/aoms/1177693055. MR0309251 Manole, T. and Ramdas, A. (2023). Martingale methods for sequential estimation of convex functionals and divergences. IEEE Transactions on Information Theory 69(7) 4641–4658. MR4613565 Mei, Y. (2010). Efficient scalable schemes for monitoring a large number of data streams. Biometrika 97(2) 419–433. https://doi.org/10.1093/biomet/asq010. MR2650748 Moustakides, G. V. (1986). Optimal stopping times for detecting changes in distributions. The Annals of Statistics 14(4) 1379–1387. https://doi.org/10.1214/aos/1176350164. MR0868306 Moustakides, G. V. (2008). Sequential change detection revisited. The Annals of Statistics 36(2) 787–807. https://doi.org/10.1214/009053607000000938. MR2396815 Page, E. S. (1954). Continuous inspection schemes. Biometrika 41(1/2) 100–115. https://doi.org/10.1093/biomet/41.1-2.100. MR0088850 Podkopaev, A., Blöbaum, P., Kasiviswanathan, S. P. and Ramdas, A. (2023). Sequential Kernelized Independence Testing. In: International Conference on Machine Learning. Pollak, M. and Siegmund, D. (1975). Approximations to the expected sample size of certain sequential tests. The Annals of Statistics 3(6) 1267–1282. MR0403114 Pollak, M. (1985). Optimal detection of a change in distribution. The Annals of Statistics 13(1) 206–227. https://doi.org/10.1214/aos/1176346587. MR0773162 Polunchenko, A. S. and Tartakovsky, A. G. (2010). On optimality of the Shiryaev-Roberts procedure for detecting a change in distribution. The Annals of Statistics 38(6) 3445–3457. https://doi.org/10.1214/09-AOS775. MR2766858 Ramdas, A., Ruf, J., Larsson, M. and Koolen, W. (2020). Admissible anytime-valid sequential inference must rely on nonnegative martingales. arXiv preprint arXiv:2009.03167. Ramdas, A., Ruf, J., Larsson, M. and Koolen, W. M. (2022). Testing exchangeability: Fork-convexity, supermartingales and e-processes. International Journal of Approximate Reasoning. 141 83–109. https://doi.org/10.1016/j.ijar.2021.06.017. MR4364897 Ramdas, A., Grünwald, P., Vovk, V. and Shafer, G. (2023). Game-theoretic statistics and safe anytime-valid inference. Statistical Science 38(4) 576–601. https://doi.org/10.1214/23-sts894. MR4665027 Ritov, Y. (1990). Decision theoretic optimality of the CUSUM procedure. The Annals of Statistics 18(3) 1464–1469. https://doi.org/10.1214/aos/1176347761. MR1062720 Roberts, S. (1966). A comparison of some control chart procedures. Technometrics 8(3) 411–430. https://doi.org/10.2307/1266688. MR0196887 Shafer, G. and Vovk, V. (2019) Game-Theoretic Foundations for Probability and Finance 455. John Wiley & Sons. https://doi.org/10.1002/0471249696. MR1852450 Shekhar, S. and Ramdas, A. (2023). Nonparametric two sample testing by betting. IEEE Transactions on Information Theory. https://doi.org/10.1109/TIT.2023.3305867 Shekhar, S. and Ramdas, A. (2023). Sequential change detection via backward confidence sequences. In International Conference on Machine Learning. Shin, J., Ramdas, A. and Rinaldo, A. (2021). Nonparametric iterated-logarithm extensions of the sequential generalized likelihood ratio test. IEEE Journal on Selected Areas in Information Theory. Shiryaev, A. N. (1963). On optimum methods in quickest detection problems. Theory of Probability & Its Applications 8(1) 22–46. MR0155708 Siegmund, D. and Venkatraman, E. (1995). Using the generalized likelihood ratio statistic for sequential detection of a change-point. The Annals of Statistics 23 255–271. https://doi.org/10.1214/aos/1176324466. MR1331667 Tartakovsky, A., Nikiforov, I. and Basseville, M. (2014) Sequential Analysis: Hypothesis Testing and Changepoint Detection. CRC Press. MR3241619 Tartakovsky, A. G. (2014). Nearly optimal sequential tests of composite hypotheses revisited. Proceedings of the Steklov Institute of Mathematics 287(1) 268–288. https://doi.org/10.1134/S0081543814080161. MR3484334 Tartakovsky, A. G., Pollak, M. and Polunchenko, A. S. (2012). Third-order Asymptotic Optimality of the Generalized Shiryaev–Roberts Changepoint Detection Procedures. Theory of Probability & Its Applications 56(3) 457–484. https://doi.org/10.4213/tvp4406. MR3136464 Ville, J. (1939). Étude critique de la notion de collectif. Monographies des Probabilités 3 144. MR3533075 Vovk, V. (2021). Testing randomness online. Statistical Science 36(4) 595–611. https://doi.org/10.1214/20-sts817. MR4323055 Vovk, V., Gammerman, A. and Shafer, G. (2005) Algorithmic Learning in a Random World. Springer Science & Business Media. MR2161220 Wang, H. and Ramdas, A. (2023). Catoni-style confidence sequences for heavy-tailed mean estimation. Stochastic Processes and Applications 163 168–202. https://doi.org/10.1016/j.spa.2023.05.007. MR4610125 Wang, H. and Ramdas, A. (2023). Huber-robust confidence sequences. In International Conference on Artificial Intelligence and Statistics 9662–9679. PMLR. Wasserman, L., Ramdas, A. and Balakrishnan, S. (2020). Universal inference. Proceedings of the National Academy of Sciences 117(29) 16880–16890. https://doi.org/10.1073/pnas.1922664117. MR4242731 Waudby-Smith, I. and Ramdas, A. (2023). Estimating means of bounded random variables by betting. Journal of the Royal Statistical Society Series B: Statistical Methodology. Willsky, A. and Jones, H. (1976). A generalized likelihood ratio approach to the detection and estimation of jumps in linear systems. IEEE Transactions on Automatic Control 21(1) 108–112. https://doi.org/10.1109/tac.1976.1101146. MR0396019 Xie, Y. and Siegmund, D. (2013). Sequential multi-sensor change-point detection. In 2013 Information Theory and Applications Workshop (ITA) 1–20. IEEE. https://doi.org/10.1214/13-AOS1094. MR3099117 Zou, S., Veeravalli, V. V., Li, J. and Towsley, D. (2019). Quickest detection of dynamic events in networks. IEEE Transactions on Information Theory 66(4) 2280–2295. https://doi.org/10.1109/TIT.2019.2948350. MR4087685