Approximations for the boundary crossing probabilities of moving sums of random variables

TL;DR

This paper develops highly accurate approximation methods for boundary crossing probabilities of moving sums of i.i.d. normal variables, applicable to various scenarios including non-normal data and unequal weights, with practical implications for ARL estimation.

Contribution

It introduces a novel approximation approach transforming a discrete problem into a continuous one, improving accuracy for small window lengths and non-normal data.

Findings

Approximations are highly accurate even for small window lengths.

Method extends well to non-normal and unequal weight scenarios.

Provides simple, accurate formulas for average run length (ARL).

Abstract

In this paper we study approximations for the boundary crossing probabilities of moving sums of i.i.d. normal r.v. We approximate a discrete time problem with a continuous time problem allowing us to apply established theory for stationary Gaussian processes. By then subsequently correcting approximations for discrete time, we show that the developed approximations are very accurate even for small window length. Also, they have high accuracy when the original r.v. are not exactly normal and when the weights in the moving window are not all equal. We then provide accurate and simple approximations for ARL, the average run length until crossing the boundary.

Tables

Table 1. Table 1 : Values of λ L , 1 ​ ( h ) subscript 𝜆 𝐿 1 ℎ {\lambda}_{L,1}(h) , λ L , 2 ​ ( h ) subscript 𝜆 𝐿 2 ℎ {\lambda}_{L,2}(h) and μ L ​ ( h ) subscript 𝜇 𝐿 ℎ {\mu}_{L}(h) with L = 20 𝐿 20 L=20 for different h ℎ h .

	$h$=0	$h$=0.5	$h$=1	$h$=1.5	$h$=2	$h$=2.5	$h$=3	$h$=3.5	$h$=4
​​​​ ${\lambda }_{L,1}\left(h\right)$ ​​​​	0.28494	0.46443	0.65331	0.81186	0.91687	0.97090	0.99209	0.99835	0.99974
​​​​ ${\lambda }_{L,2}\left(h\right)$ ​​​​	0.25744	0.43811	0.63472	0.80239	0.91348	0.97005	0.99195	0.99833	0.99974
​​​​ ${\mu }_L\left(h\right)$ ​​​​	0.25527	0.43677	0.63432	0.80241	0.91353	0.97007	0.99195	0.99833	0.99974

Table 2. Table 2 : Average values from 100 evaluations of Approximation 1 for different h ℎ h along with maximum and minimum with L = 5 𝐿 5 L=5 and T = 100 𝑇 100 T=100 .

	$h$=2.5	$h$=2.75	$h$=3	$h$=3.25	$h$=3.5	$h$=3.75	$h$=4
​​​​ $\overline{x}$ ​​	0.855957	0.627299	0.376337	0.191122	0.086253	0.033769	0.013156
​​​​ $s$ ​​	0.004127	0.008588	0.013805	0.015181	0.012826	0.008510	0.005131
​​​​ $Max\left(x_i\right)-\overline{x}$ ​​	0.010665	0.023748	0.029819	0.027066	0.025629	0.016208	0.011609
​​​​ $\overline{x}-Min\left(x_i\right)$ ​​	0.012176	0.021268	0.033211	0.041322	0.041350	0.022650	0.018146
​​​​ Approximation 3​​	0.854844	0.625113	0.373863	0.188933	0.083981	0.033833	0.012551
​​​​ $P_L(T,h)$ ​​	0.855429	0.627463	0.376681	0.191625	0.085697	0.034675	0.013116

Table 3. Table 3 : Average values from 100 evaluations of Approximation 1 for different h ℎ h along with maximum and minimum with L = 20 𝐿 20 L=20 and T = 100 𝑇 100 T=100 .

	$h$=2.5	$h$=2.75	$h$=3	$h$=3.25	$h$=3.5	$h$=3.75	$h$=4
​​​​ $\overline{x}$ ​​	0.952007	0.802073	0.554613	0.315085	0.155331	0.066113	0.025608
​​​​ $s$ ​​	0.001479	0.004856	0.012540	0.015050	0.015160	0.011647	0.008129
​​​​ $Max\left(x_i\right)-\overline{x}$ ​​	0.004746	0.013360	0.027078	0.030940	0.033991	0.024111	0.030014
​​​​ $\overline{x}-Min\left(x_i\right)$ ​​	0.003662	0.010894	0.031463	0.037715	0.041021	0.043283	0.016997
​​​​ Approximation 3​​	0.952475	0.802100	0.555109	0.316076	0.153803	0.066438	0.026143
​​​​ $P_L(T,h)$ ​​	0.952818	0.803078	0.555530	0.315784	0.153446	0.066642	0.026244

Table 4. Table 4 : Average values from 100 evaluations of Approximation 1 for different h ℎ h along with maximum and minimum with L = 100 𝐿 100 L=100 and T = 100 𝑇 100 T=100 .

	$h$=2.5	$h$=2.75	$h$=3	$h$=3.25	$h$=3.5	$h$=3.75	$h$=4
​​​​ $\overline{x}$ ​​	0.979027	0.878031	0.661247	0.402887	0.211894	0.093329	0.039110
​​​​ $s$ ​​	0.000884	0.005502	0.014418	0.021283	0.018493	0.020459	0.015536
​​​​ $Max\left(x_i\right)-\overline{x}$ ​​	0.001995	0.009243	0.039695	0.040615	0.063578	0.064306	0.037958
​​​​ $\overline{x}-Min\left(x_i\right)$ ​​	0.002414	0.020613	0.025530	0.093876	0.038484	0.05694	0.033748
​​​​ Approximation 3​​	0.979119	0.878481	0.660662	0.405674	0.209313	0.094517	0.038529

Table 5. Table 5 : Approximations for ARL h ​ ( 𝕏 ) subscript ARL ℎ 𝕏 {\rm ARL}_{h}(\mathbb{X}) and S ​ D ​ ( τ h ​ ( 𝕏 ) ) 𝑆 𝐷 subscript 𝜏 ℎ 𝕏 SD(\tau_{h}(\mathbb{X})) with L = 10 𝐿 10 L=10 .

$h$	2	2.25	2.5	2.75	3	3.25	3.5
(7.4)	126	217	395	759	1551	3375	7837
(7.6)	126 ​$\pm$ ​​1	218 ​$\pm$ ​​2	394 ​$\pm$ ​​5	756 ​$\pm$ ​​17	1545 ​$\pm$ ​​65	3388 ​$\pm$ ​​300	7791 ​$\pm$ ​​1100
${\mathrm{ARL}}_h\left(𝕏\right)$	127	218	396	757	1550	3344	7721

Table 6. Table 6 : Approximations for ARL h ​ ( 𝕏 ) subscript ARL ℎ 𝕏 {\rm ARL}_{h}(\mathbb{X}) and S ​ D ​ ( τ h ​ ( 𝕏 ) ) 𝑆 𝐷 subscript 𝜏 ℎ 𝕏 SD(\tau_{h}(\mathbb{X})) with L = 50 𝐿 50 L=50 .

$h$	2	2.25	2.5	2.75	3	3.25	3.5
(7.4)	471	791	1392	2587	5099	10695	23918
(7.6)	471 ​$\pm$ ​​3	791 ​$\pm$ ​​7	1393 ​$\pm$ ​​25	2597 ​$\pm$ ​​75	​​​ 5101 ​$\pm$ ​​270 ​​​	​​​​ 10708 ​$\pm$ ​​1250 ​​​​	​​​​ 24639 ​$\pm$ ​​5800 ​​​​
${\mathrm{ARL}}_h\left(𝕏\right)$	472	792	1397	2588	5085	10749	24131

Equations

$$S_{n,L}:=\sum _{j=n+1}^{n+L}{\epsilon }_j(n=0,1,\dots ).$$

$${\mathcal{P}}_{\mathbb{S}}(M,H,L):=\mathrm{Pr}\left(\underset{n=0,1,\dots ,M}{\max }{S}_{n,L}\ge H\right),$$

$$\displaystyle\xi_{n,L}:=\frac{S_{n,L}-\operatorname{\mathbb{E}}S_{n,L}}{\sqrt{{\rm var}(S_{n,L})}}=\frac{S_{n,L}-\theta L}{\sigma\sqrt{L}}\,,\;\;\mbox{$n=0,1,\ldots\,,$}$$

$$\phi (x):=\frac{1}{\sqrt{2\pi }}{e}^{-x^2/2},\Phi (t):={\int }_{-\infty }^t\phi (x)dx.$$

$$\displaystyle{\rm Corr}(\xi_{n,L},\xi_{n+k,L})=\max\{0,1-{k}/{L}\}=\left\{\begin{array}[]{cl}1-{k}/{L}&\;\;\;\;{\rm for}\;0\leq k\leq L\\ 0&\;\;\;\;{\rm for}\;k>L\,.\end{array}\right.$$

$$h=\frac{H-\theta L}{\sigma \sqrt{L}}\text{mbox}sothatH=\theta L+\sigma h\sqrt{L}.$$

$$P_L(T,h):=\text{Pr}\left(\underset{n=0,1,\dots ,TL}{\max }{\xi }_{n,L}\ge h\right).$$

$${\mathcal{P}}_{\mathbb{S}}(M,H,L)=P_L(T,h)\text{mbox}foranyH,L\text{mbox}andT=M/L.$$

$$F_L(T,h)=\text{Pr}\left(\underset{n=0,1,\dots ,TL}{\max }{\xi }_{n,L}<h\right).$$

$$P_L(T,h)\simeq 1-F_L(2,h){\left[\frac{F_L(2,h)}{F_L(1,h)}\right]}^{T-2},$$

$$S_L(t)=\frac{1}{\Delta }\left[(t_n-t){\xi }_{n-1,L}+(t-t_{n-1}){\xi }_{n,L}\right]\mathrm{for}t\in [t_{n-1},t_n],n=1,\dots ,M.$$

$$P(T,h):=\mathrm{Pr}\left\{\underset{0\le t\le T}{\max }S(t)\ge h\right\}.$$

$$F(T,h)=\mathrm{Pr}\left\{\underset{0\le t\le T}{\max }S(t)<h\right\}=1-P(T,h).$$

$$\displaystyle{F}(T,h\,|\,x):={\rm Pr}\Big{\{}S(t)\!<\!h\;{\rm for\;all\;}t\!\in\![0,\!{T}]\,|\,S(0)=x\Big{\}}\,.\;\;$$

$$F(n,h|x)=\frac{1}{\phi (x)}{\int }_{D_x}\det [\phi (y_i-y_{j+1}+h){]}_{i,j=0}^{n}d{y}_2\dots d{y}_{n+1}$$

$$F(1,h)$$

$$F(2,h)$$

$$F(n,h|x)=\frac{1}{\phi (x)}{\int }_{-\infty }^h\dots {\int }_{-\infty }^h\det [\phi (s_i+a_{i,j}){]}_{i,j=0}^{n}d{s}_1\dots d{s}_n,$$

$$\displaystyle a_{i,j}=y_{i+1}\!-\!y_{j+1}=\left\{\begin{array}[]{cl}0&\text{ for }i=j\\ (i-j)h\!-\!s_{j+1}\!-\!\ldots\!-\!s_{i+1}&\text{ for }i>j\\ (i-j)h+s_{i+1}+\ldots+s_{j}&\text{ for }i<j\,.\par\end{array}\right.$$

$$p(s_0,s_1,\dots s_n)=\frac{1}{\phi (s_0)F(n,h|s_0)}\det [\phi (s_i+a_{i,j}){]}_{i,j=0}^n.$$

$$q_h^{(0,n)}(x\to s_n)=\frac{1}{\phi (x)}{\int }_{-\infty }^h\dots {\int }_{-\infty }^h\det [\phi (s_i+a_{i,j}){]}_{i,j=0}^{n}d{s}_1\dots d{s}_{n-1}.$$

$$p_h^{(0,n)}(z):={\int }_{-\infty }^h{q}_h^{(0,n)}(x\to z)\phi (x)dx$$

$$q^{(0,1)}_{h}(x\to z)=\frac{1}{\varphi(x)}\det\left(\begin{array}[]{cc}\varphi(x)&\varphi(x\!-\!h\!+\!z)\\ \varphi(h)&\varphi(z)\\ \end{array}\right)=\varphi(z)\left[1-e^{-(h-z)(h-x)}\right]\,,\,\,z=s_{1}<h.$$

$$p_h^{(0,1)}(z)={\int }_{-\infty }^h{q}_h^{(0,1)}(x\to z)\phi (x)dx=\Phi (h)\phi (z)-\Phi (z)\phi (h)$$

$$q_h^{(j,n)}(x\to z)=\frac{1}{p_h^{(0,j)}(z)}{\int }_{-\infty }^h\dots {\int }_{-\infty }^h\det [\phi (s_i+a_{i,j}){]}_{i,j=0}^{n}d{s}_{0}d{s}_1\dots d{s}_{j-1}d{s}_{j+1}\dots d{s}_{n-1},$$

$$q_h^{(1,2)}(x\to z)$$

$$=$$

$$\mathrm{\Omega }$$

$$\displaystyle{F}(T,h\,|\,x)=\int\cdots\int{\rm Pr}\{{\rm\,\Omega}\,\big{|}\,W(i)=x_{i},\,\,i=0,1,2,\ldots,T+1,\,\,\,W(0)=0,\,\,W(0)-W(1)=x\}$$

$$\displaystyle\times\,{\rm Pr}\{W(i)\in dx_{i},\,\,i=0,1,2,\ldots,T+1,\,\,\big{|}\,\,W(0)=0,\,\,W(0)-W(1)=x\}.$$

$$\displaystyle{\rm Pr}\{{\rm\Omega}\,\,\big{|}\,\,W(i)=x_{i},\,\,i=0,1,2,\ldots,T+1,\,\,\,W(0)=0,\,\,W(0)-W(1)=x\}$$

Full text

1 Introduction: Statement of the problem

Let $\varepsilon_{1},\varepsilon_{2},\ldots$ be a sequence of i.i.d. normal random variables (r.v.) with mean $\theta$ and variance $\sigma^{2}>0$. For a fixed positive integer $L$, the moving sums are defined by

[TABLE]

The sequence of the moving sums (1.1) will be denoted by $\mathbb{S}$ so that $\;\mathbb{S}=\{S_{0,L},S_{1,L},\ldots\}$.

The main aim of this paper is development of accurate approximations for the boundary crossing probability (BCP) for the maximum of the moving sums:

[TABLE]

where $M$ is a given positive integer and $H$ is a fixed threshold. Note that the total number of r.v. $\varepsilon_{i}$ used in (1.2) is $M+L$ and ${\cal P}_{\mathbb{S}}(M,H,L)\to 1$ as $M\to\infty$, for all $H$ and $L$. We will mostly be interested in deriving accurate approximations when $M\geq L$. The case of $M\leq L$ is much simpler and is comprehensively covered in (AandZ2019, , Section 3), see Section 4.6 for discussion.

Developing accurate approximations for the BCP ${\cal P}_{\mathbb{S}}(M,H,L)$ for generic parameters $H$, $M$ and $L$ is very important in various areas of statistics, predominantly in applications related to change-point detection; see, for example, papers Bau2 ; Chu ; Glaz2012 ; MZ2003 ; Xia and especially books glaz2001scan ; glaz2009scan2 . Engineering applications of MOSUM (moving sums charts) are extremely important and have been widely discussed in literature; see e.g. Chu ; eiauer1978use ; glaz2001scan ; glaz2009scan2 ; waldmann1986bounds . The BCP ${\cal P}_{\mathbb{S}}(M,H,L)$ is an ($M+1$)-dimensional integral and therefore direct evaluation of this BCP is hardly possible even with modern software.

To derive approximations for the BCP (1.2) one can use standard tools and approximate the sequence of moving sums with a continuous-time process and then use some continuous-time approximations, see e.g. Haiman ; these approximations, however, are not accurate especially for small window length $L$; see discussion in Section 4.7. There is, therefore, a need for derivation of specific approximations for the BCP (1.2). Such a need was well understood in the statistical community and indeed very accurate approximations for the BCP and the Average Run Length (ARL) have been developed in a series of quality papers by J. Glaz and coauthors, see for example Glaz_old ; Glaz2012 ; wang2014variable ; wang2014multiple (the methodology was also extended to the case when $\varepsilon_{j}$ are integer-valued r.v., see glaz1991tight ). We will call these approximations ‘Glaz approximations’ by the name of the main author of these papers; they will be formally written down in Sections 2.2 and 7.

The accuracy of the approximations developed in the present paper is very high and similar to the Glaz approximations; this is discussed in Sections 6 and 7. The methodologies of derivation of Glaz approximations and the approximations of this paper are very different. The practical advantage of our approximations (they require approximating either a one-dimensional integral or an eigenvalue of an integral operator) is their relative simplicity as to compute the Glaz approximations one needs to numerically approximate $L+1$ and $2L+1$ dimensional integrals. This is not an easy task even taking into account the fact of existence of a sophisticated software; see references in Section 2.2.

The paper is structured as follows. In Section 2 we reformulate the problem, state the Glaz approximation and discuss how to approximate our discrete-time problem with a continuous-time problem. In Section 3 we provide exact formulas for the first-passage probabilities (in the continuous-time setup) due to L. Shepp Shepp71 and give their alternative representation which will be crucial for deriving some of our approximations. In Section 4 we adapt the methodology of D. Siegmund to correct Shepp’s formulas for discrete time and define a version of the Glaz approximation which we will call Glaz-Shepp-Siegmund approximation. In Section 5 we develop continuous-time approximations based on approximating eigenvalues of integral operators and subsequently correct them for discrete time. In Sections 4.7 and 6 we present results of large-scale simulation studies evaluating the performance of the considered approximations (also, in the cases when the original r.v. $\varepsilon_{j}$ are not normal and the weights in the moving window are not equal). In Section 7, we develop an approximation for ARL and compare its accuracy to the one developed in Glaz2012 .

2 Boundary crossing probabilities: discrete and continuous time

2.1 Standardisation of the moving sums

The first two moments of $S_{n,L}$ are $\operatorname{\mathbb{E}}S_{n,L}=\theta L$ and ${\rm var}(S_{n,L})=\displaystyle\sigma^{2}L.$ Define

[TABLE]

which are the standardized versions of $S_{n,L}$. All r.v. $\xi_{n,L}$ are $N(0,1)$; that is, they have the probability density function and c.d.f.

[TABLE]

Unlike the original r.v. $\varepsilon_{i}$, the r.v. $\xi_{0,L},\xi_{1,L},\ldots$ are correlated so that for all $k=0,1,\ldots$ we have ${\rm Corr}(\xi_{0,L},\xi_{k,L})={\rm Corr}(\xi_{n,L},\xi_{n+k,L})$ and

[TABLE]

Proof of (2.5) is straightforward, see (AandZ2019, , Lemma 1).

Set $T=M/L$ and

[TABLE]

Define the BCP for the sequence of r.v. $\xi_{0,L},\xi_{1,L},\ldots$:

[TABLE]

From (2.1) and (2.6), the BCPs ${\cal P}_{\mathbb{S}}(M,H,L)$ and ${P}_{L}(T,h)$ are equal:

[TABLE]

Note also that ${P}_{L}(T,h)=1-F_{L}(T,h)$, where

[TABLE]

In accordance with the terminology of Shepp71 and slepian1961first we shall call $F_{L}(T,h)$ ‘first-passage probability’. In the following sections, we derive approximations for (2.7). These approximations will be based on approximating the sequence of r.v. $\{\xi_{0,L},\xi_{1,L},\ldots,\xi_{M,L}\}$ by a continuous-time random process and subsequently correcting the obtained approximations for discreteness. Before doing this, we formulate the approximation which is currently the state-of-the-art.

2.2 Glaz approximation for ${P}_{L}(T,h)$

The approximation for the BCP ${P}_{L}(T,h)$ developed in Glaz_old ; Glaz2012 ; wang2014variable ; wang2014multiple and discussed in the introduction is as follows.

Approximation 1. (Glaz approximation) For $T\geq 2$,

[TABLE]

*where to approximate the first-passage probabilities $F_{L}(1,h)$ and $F_{L}(2,h)$, which are $L+1$ and $2L+1$ dimensional integrals respectively, it is advised to use the so-called ‘GenzBretz’ algorithm for numerical evaluation of multivariate normal probabilities; see genz2009computation ; GenzR .

Unless $h$ is large (say, $h>3$), Approximation 1 is very accurate. However, its computational cost is also high, especially for large $L$. Moreover, the main option in the ‘GenzBretz’ package requires the use of Monte-Carlo simulations so that for reliable estimation of high-dimensional integrals one needs to make a lot of averaging; see Section 6.1 and 7 for more discussion on these issues.

2.3 Continuous-time (diffusion) approximation

For the purpose of approximating the BCP ${P}_{L}(T,h)$, we replace the discrete-time process $\xi_{0,L},\ldots,$ $\xi_{M,L}$ with a continuous process $S(t)$, $t\in[0,T]$, where $T=M/L$ (we will then correct the corresponding first-passage probabilities for discreteness). We do this as follows.

Set $\Delta=1/L$ and define $t_{n}=n\Delta\in[0,{T}]\;n=0,1,\ldots,{M}.$ Define a piece-wise linear continuous-time process ${S_{L}(t)},$ $t\in[0,T]:$

[TABLE]

By construction, the process ${S_{L}(t)}$ is such that ${{S_{L}(t_{n})}}=\xi_{n,L}\;{\rm for}\;n=0,\ldots,{M}$. Also we have that ${S_{L}(t)}$ is a second-order stationary process in the sense that $\operatorname{\mathbb{E}}{S_{L}(t)},\,$ ${\rm var}({S_{L}(t)})$ and the autocorrelation function $R^{(L)}(t,t+k\Delta)={\rm Corr}({S_{L}(t)},{S_{L}(t+k\Delta)})$ do not depend on $t$.

Lemma 1

Assume $L\to\infty$. The limiting process $S(t)$ = $\lim_{L\rightarrow\infty}{S_{L}(t)}$, where $t\in[0,T]$, is a Gaussian second-order stationary process with marginal distribution $S(t)\sim N(0,1)$ for all $t\in[0,T]$ and autocorrelation function $R(t,t+s)=R(s)=\max\{0,\;1\!-\!|s|\}\,$.

This lemma is a simple consequence of (2.5).

2.4 Diffusion approximations: definition and their role in this study

The above approximation of a discrete-time process $\{\xi_{0,L},\xi_{1,L},\ldots,\xi_{M,L}\}$ with a continuous process $S(t),\,t\in[0,T]$, allows us to approximate the BCP ${P}_{L}(T,h)$ by a continuous-time analogue as follows.

By the definition of a diffusion approximation, the BCP ${P}_{L}(T,h)$ is approximated by

[TABLE]

Note that approximating the discrete process of moving sums by a continuous-time process $S(t)$ and subsequent approximation of the BCP ${P}_{L}(T,h)$ by ${P}(T,h)$ is by no means new. This has been done, in particular, in Haiman .

We will call (2.10) and any approximation to (2.10), which does not involve the knowledge of $L$, ‘diffusion approximation’. These approximations can be greatly improved with the help of the methodology developed by D.Siegmund and adapted to our setup in Section 4. The importance of the discrete-time correction is illustrated by Figures 1 and 2, where for a fixed $h$ and $T$ we can see a significant difference in values of the BCPs ${P}_{L}(T,h)$ for different values of $L$. As seen from Figure 2, even for very large $L=1000$, the discrete-time correction is still needed. Hence we are not recommending to use any approximation for ${P}(T,h)$ (including rather sophisticated ones like the one developed in Haiman ) as an approximation for ${P}_{L}(T,h)$. In the next section we will discuss a diffusion approximation that, after correcting for discrete time, will be a cornerstone for all approximations developed in this paper.

In what follows, it will also be convenient to use the first-passage probability

[TABLE]

Since $\xi_{0,L}=S(0)\sim N(0,1)$, we have ${F}(0,h)=1-{P}(0,h)=\Phi(h)$.

3 Exact formulas for the first-passage probabilities in the continuous-time case

3.1 Shepp’s formulas

Define the conditional first-passage probability

[TABLE]

Since ${F}(T,h\,|\,x)=0$ for $x>h$, for the unconditional first-passage probability ${F}(T,h)$ we have ${F}(T,h)=\int_{-\infty}^{h}{F}(T,h\,|\,x)\varphi(x)dx$.

The result of (Shepp71, , p.949) states than if $T=n$ is a positive integer then

[TABLE]

where $y_{0}=0,y_{1}=h-x,$ $D_{x}=\{y_{2},\dots,y_{n+1}\>|\>h-x<y_{2}<y_{3}<\ldots<y_{n+1}\}$. For non-integer $T\geq 1$, the exact formula for ${F}(T,h\,|\,x)$ is even more complex (the integral has the dimension $\lceil 2T\rceil$) and completely impractical for computing ${P}(T,h)$ with $T>2$, see (Shepp71, , p.950).

For $n=1$, we obtain

[TABLE]

For $n=2$, (3.2) yields

[TABLE]

The three-dimensional integral in (3.4) can be reduced to a one-dimensional, see (4.11) below with $h_{L}=h$.

3.2 An alternative representation of the Shepp’s formula (3.2)

Set $s_{i}=h+y_{i}-y_{i+1}$ ($i=0,1,\ldots,n$) with $s_{0}=x$, $y_{0}=0,y_{1}=h-x$. It follows from Shepp’s proof of (3.2) that $s_{0},s_{1},\ldots,s_{n}$ have the meaning of the values of the process $S(t)$ at the times $t=0,1,\ldots,n$: $S(i)=s_{i}$ ($i=0,1,\ldots,n$). The range of the variables $s_{i}$ is $(-\infty,h)$. Changing the variables in (3.2), we obtain

[TABLE]

where

[TABLE]

3.3 Joint density for the values $\{S(i)\}$ and associated transition densities

From (3.5), we obtain the following expression for the joint probability density function for the values $S(0),S(1),\ldots,S(n)$ under the condition $S(t)<h$ for all $t\in[0,n]$:

[TABLE]

From this formula, we can derive the transition density from $s_{0}\!=\!x$ to $s_{n}$ conditionally $S(t)\!<\!h,$ $\forall t\in[0,n]$:

[TABLE]

For this transition density, $\int_{-\infty}^{h}q^{(0,n)}_{h}(x\to z)dz={F}(n,h\,|\,x)$. Moreover, since $S(0)\sim N(0,1)$, the non-normalized density of $S(n)$ under the condition $S(t)<h$ for all $t\in[0,n]$ is

[TABLE]

with $z<h$ and $\int_{-\infty}^{h}{p}^{(0,n)}(z)dz={F}(n,h)$. In the case $n=1$, (3.8) gives

[TABLE]

From this and (3.9) we get

[TABLE]

with $z<h$ and $\int_{-\infty}^{h}{p}^{(0,n)}(z)dz={F}(1,h)$.

Rather than just recovering the transition density from $s_{0}=x$ to $s_{n}$, we can also use (3.7) and (3.9) to obtain the transition density from $x=s_{j}$ to $z=s_{n}$, $0<j<n$, under the condition $S(t)<h$ for all $t\in[0,n]$:

[TABLE]

where $s_{j}=x$ and $s_{n}=z$. For $j=1$ and $n=2$ we obtain the transition density from $x=s_{1}$ to $z=s_{2}$ under the condition $S(t)<h$ for all $t\in[0,2]$:

[TABLE]

4 Correcting Shepp’s formula (3.2) for discrete time

4.1 Rewriting (3.2) in terms of the Brownian motion

Let $W(t)$ be the standard Brownian Motion process on $[0,\infty)$ with $W(0)=0$ and $\operatorname{\mathbb{E}}W(t)W(s)=\min(t,s).$ Recall the conditional probability ${F}(T,h\,|\,x)$ defined in (3.1). Suppose $T\geq 1$ is an integer and define the event

[TABLE]

If $W(i)=x_{i}$, $i=0,1,\ldots,T+1$, we obtain from (Shepp71, , p.948)

[TABLE]

It follows from the proof of (3.2) that to correct (4.1) for discrete time, one must correct the following probability for discrete time

[TABLE]

where $W_{i}(t)=\frac{\sqrt{2}}{2}[W(t+i-1)-W(t+i)]$, $i=1,2,\ldots,T$. Due to the conditioning on the rhs of (4.2), the processes $W_{i}(t)$ can be treated as independent Brownian motion processes. Therefore, the independent increments of the Brownian motion means correcting formula (3.2) for discrete time is equivalent to correcting the probability ${\rm Pr}(\sqrt{2}\,W(t)<h,\,\,\forall\,\,\,0\leq t\leq 1\,\,)$ for discrete time.

4.2 Discrete-time correction for the BCP of cumulative sums.

Let $X_{1},X_{2},\ldots$ be i.i.d. $N(0,1)$ r.v’s and set $Y_{n}=X_{1}+X_{2}+\ldots+X_{n}$. Consider the sequence of cumulative sums $\{Y_{n}\}$ and define the stopping time $\tau_{Y,a,b}=\inf\{n\geq 1:Y_{n}\geq a+bn\}$ for $a>0$ and $b\in\mathbb{R}$. Consider the problem of evaluating

[TABLE]

Exact evaluation of (4.3) is difficult even if $N$ is not very large but it was accurately approximated by D.Siegmund see e.g. (Sieg_paper, , p.19). Let $W(t)$ be the standard Brownian Motion process on $[0,\infty)$. For $a>0$ and $b\in\mathbb{R}$, define $\tau_{W,a,b}=\inf\{t:W(t)\geq a+bt\}$ so that

[TABLE]

In Sieg_paper , (4.4) was used to approximate (4.3) after translating the barrier $a+bt$ by a suitable scalar $\rho\geq 0$. Specifically, the following approximation has been constructed:

[TABLE]

where the constant $\rho$ approximates the expected excess of the process $\{Y_{n}\}$ over the barrier $a+bt$. From (Sieg_book, , p. 225)

[TABLE]

4.3 Discretised Brownian motion

Define $\epsilon$ = $1/L$ and let $t^{\prime}_{n}=n\epsilon\in[0,{1}],$ $n=0,1,\ldots,{L}.$ Let $X_{1},X_{2},\ldots$ be i.i.d. $N(0,1)$ r.v’s and set $W(t^{\prime}_{n})=\sqrt{\epsilon}\sum_{i=1}^{{n}}X_{i}.$ For $a>0$ define the stopping time

[TABLE]

and consider the problem of approximating

[TABLE]

As $L\rightarrow\infty$, the piecewise linear continuous-time process $W^{\epsilon}(t)$, $t\in[0,1]$, defined by:

[TABLE]

converges to $W(t)$ on $[0,1]$ as so we can refer to $W(t^{\prime}_{n})$ as discretised Brownian motion. We make the following connection between $\sqrt{2}W(t^{\prime}_{n})$ and the random walk $Y_{n}$:

[TABLE]

Then by using (4.5), we approximate the expected excess over the boundary for the process $\sqrt{2}W(t^{\prime}_{n})$ by

[TABLE]

We have deliberately rounded the value $\sqrt{2}\rho\simeq 0.8239...$ to $0.82$ as for small $h$ and small $L$ it provides marginally better approximation (4.9).

4.4 Corrected version of (3.2)

Set ${h}_{L}=h+\omega_{L}$. To correct (3.2) for discrete time we substitute the barrier $h$ with $h_{L}$. From this and the relation ${F}(T,h)=\int_{-\infty}^{h}{F}(T,h\,|\,x)\varphi(x)dx$, the discrete-time corrected form of ${F}(T,h)$ is

[TABLE]

where $y_{0}=0,y_{1}={h_{L}}-x,$ and $D_{x}=\{y_{2},\dots,y_{T+1}\>|\>{h_{L}}-x<y_{2}<y_{3}<\ldots<y_{T+1}\}.$

4.5 A generic approximation involving corrected Shepp’s formula

Approximation 2. For integral $T\geq 1$, the discrete-time correction for the BCP (2.7) is

[TABLE]

*where ${F}_{L}(T,{h},h_{L})$ is given in (4.4).

Whilst Approximation 2 is very accurate (see the next subsection), computation of ${P}(T,{h},h_{L})$ requires numerical evaluation of a $T+1$ dimensional integral which is impractical for large $T$. To overcome this, in Section 5.2 we develop approximations that can be easily used for any $T>0$ (which is not necessarily integer).

4.6 Particular cases: $T=1$ and $T=2$

For $T=1$, evaluation of (4.4) yields

[TABLE]

In our previous work AandZ2019 we have derived approximations $\hat{P}_{L}(T,h)$ for the BCP ${P}_{L}(T,h)$ with $T\leq 1$. The approximations $\hat{P}_{L}(T,h)$ developed in AandZ2019 are also discrete-time corrections of the continuous-time probabilities ${P}(T,h)$ but they are based almost exclusively on the fact that the process $S(t)$ is conditionally Markov on the interval $t\in[0,1]$; hence the technique of AandZ2019 cannot be extended for intervals $t\in[0,T]$ with $T>1$. The approximation $\hat{P}_{L}(1,h)$ of AandZ2019 is different from ${P}(1,{h},h_{L})=1-{F}(1,{h},h_{L})$ of (4.10). It appears that $\hat{P}_{L}(1,h)$ is more complicated and less accurate approximation than ${P}(1,{h},h_{L})$.

For $T=2$, (4.4) can be expressed (after some manipulations) as follows:

[TABLE]

Only a one-dimensional integral has to be numerically evaluated for computing ${F}(2,{h},h_{L})$.

4.7 Simulation study

In this section, we assess the quality of the approximations (4.10) and (4.11) as well as the sensitivity of the BCP ${P}_{L}(T,h)$ to the value of $L$. In Figures 1 and 2, the black dashed line corresponds to the empirical values of the BCP ${P}_{L}(T,h)$ (for $T=M/L=1,2$) computed from 100 000 simulations with different values of $L$ and $M$ (for given $L$ and $M$, we simulate $L+M$ normal random variables 100 000 times). The solid red line corresponds to Approximation 2. The axis are: the $x$-axis shows the value of the barrier $h$ in Figure 1 and value of $L$ in Figure 2; the $y$-axis denotes the probabilities of reaching the barrier. The graphs, therefore, show the empirical probabilities of reaching the barrier $h$ (for the dashed line) and values of considered approximations for these probabilities. From these graphs we can conclude that Approximation 2 is very accurate, at least for $T=1,2$. We can also conclude that the BCP ${P}_{L}(T,h)$ is very sensitive to the value of $L$. From Figure 2 we can observe a counter-intuitive fact that even for very high value $L=1000$, the BCP ${P}_{L}(T,{h})$ is not even close to ${P}_{\infty}(T,{h})={P}(T,{h})$ from (2.10). This may be explained by the fact that for any fixed $T$ and $h$, the inaccuracy $|{P}_{L}(T,h)-{P}(T,{h})|$ decreases with the rate const$/\sqrt{L}$ as $L\to\infty$.

4.8 The Glaz-Shepp-Siegmund approximation

Combining (2.9) and the approximation (4.9) for Shepp’s formula (3.2), we arrive at the following approximation to which we suggest the name ‘Glaz-Shepp-Siegmund approximation’.

Approximation 3. *For all $T>0,$ *

[TABLE]

*where ${F}(1,{h},h_{L})$ and ${F}(2,{h},h_{L})$ are defined in (4.10) and (4.11) respectively.

*Approximations 1 and 3 look similar but computing Approximation 1 is very hard and Approximation 3 is very easy (only a one-dimensional integral should be numerically computed).

5 Approximations for the BCP ${P}_{L}(T,h)$ through eigenvalues of integral operators

5.1 Continuous time: approximations for $F(T,h)$

Let $m$ be a positive integer, and $q(x\to z)$ be the transition density ${q}_{h}^{(m-1,m)}(x\to z)$ defined by (3.10) for $m=1$ (3.15) for $m=2$ and (3.11) for $m>2$.

Let us approximate the distributions of the values $s_{i}=S(i)$ for integral $i>m$ in the following way. Let $p_{i}(x)$ be the density of $S(i)$ under the condition that $S(t)$ does not reach $h$ for $t\in[0,i]$. By ignoring the past values of $S(t)$ in $[0,i)$, the non-normalized density of $S(i+1)$ under the conditions that $S(i)\sim p_{i}(x)$ and $S(t)$ does not reach $h$ for $t\in[i,i+1]$ is

[TABLE]

We can then define ${p}_{i+1}(x)=\tilde{p}_{i}(x)/c_{i},\;x<h,$ where $c_{i}=\int_{-\infty}^{h}\tilde{p}_{i}(x)dx$. We then replace formula (5.1) with

[TABLE]

where $p(x)$ is an eigenfunction of the integral operator with kernel (3.10) corresponding to the maximum eigenvalue $\lambda_{m}(h)$:

[TABLE]

This eigenfunction $p(x)$ is a probability density on $(-\infty,h]$ with $p(x)>0$ for all $x\in(-\infty,h)$ and $\int_{-\infty}^{h}p(x)dx=1\,.$ Moreover, the maximum eigenvalue $\lambda_{m}(h)$ of the operator with kernel $K(x,y)={q}_{h}^{(m-1,m)}(x\to z)$ is simple and positive. The fact that such maximum eigenvalue $\lambda_{m}(h)$ is simple and real (and hence positive) and the eigenfunction $p(x)$ can be chosen as a probability density follows from the Ruelle-Krasnoselskii-Perron-Frobenius theory of bounded linear positive operators, see e.g. Theorem XIII.43 in ReedSimon .

Using (5.2) and (5.3), we derive recursively: ${F}(i+1,h)\simeq F(i,h)\lambda_{m}(h)$ ($i=m,m+1,\ldots$). By induction, for any integer $T\geq m$ we then have

[TABLE]

The approximation (5.4) can be used for any $T>0$ which is not necessarily an integer. The most important particular cases of (5.4) are with $m=1$ and $m=2$. In these two cases, the kernel ${q}_{h}^{(m-1,m)}(x\to z)$ and hence the approximation (5.4) will be corrected for discrete time in the next section.

5.2 Correcting approximation (5.4) for discrete time

To correct the approximation (5.4) for discrete time we need to correct: (a) the first-passage probability ${F}(m,h)$ and (b) the kernel ${q}_{h}^{(m-1,m)}(x\to z)$. The discrete-time correction of ${F}(m,h)$ can be done using ${F}_{L}(m,h,h_{L})$ from (4.4) so that what is left is to correct the kernel ${q}_{h}^{(m-1,m)}(x\to z)$ and hence $\lambda_{m}(h)$.

5.2.1 Correcting the transition kernels for discrete time

As explained in Section 4, to make a discrete-time correction in the Shepp’s formula (3.2) we need to replace the barrier $h$ with ${h_{L}}=h+\omega_{L}$ in all places except for the upper bound for the initial value $S(0)$. Therefore, using the notation of Section 3.2, the joint probability density function for the values $S(0),S(1),\ldots,S(m)$ under the condition $S(t)<h$ for all $t\in[0,m]$ corrected for discrete time is:

[TABLE]

with $-\infty<s_{0}<h$, $-\infty<s_{j}<h_{L}$ $(j=1,\ldots,m)$,

[TABLE]

This gives us the discrete-time corrected transition density from $s_{0}\!=\!x$ to $s_{m}$ conditionally $S(t)\!<\!h,$ $\forall t\in[0,m]$:

[TABLE]

which is exactly (3.8) with $h_{L}$ is substituted for $h$. In a particular case $m=1$, the corrected transition density is

[TABLE]

with $s_{1}<{h_{L}}$.

Let us now make the discrete-time correction of the transition density ${q}_{h}^{(1,2)}(x\to z)$. Denote by ${p}_{h,L}^{(0,1)}(z)$, $z<h,$ the non-normalized density of $S(1)$ under the condition $S(t)<h$ for all $t\in[0,1]$ corrected for discrete time; it satisfies $\int_{-\infty}^{h}{p}_{h,L}^{(0,1)}(z)dz={F}(1,h,h_{L})$. Using (5.10), we obtain

[TABLE]

From (5.5) and (5.10), the transition density from $x=s_{1}$ to $z=s_{2}$ under the condition $S(t)<h$ for all $t\in[0,2]$ corrected for discrete time (the corrected form of (3.15)) is given by

[TABLE]

Unlike the transition density (5.7) (and (5.10) in the particular case $m=1$), which only depends on $h_{L}$ and not on $h$, the transition density ${q}_{h,L}^{(1,2)}(x\to z)$ depends on both $h$ and $h_{L}$ and hence the notation. The dependence on $h$ has appeared from integration over the $s_{0}\in(-\infty,h)$.

5.2.2 Approximations for the BCP ${P}_{L}(T,h)$

With discrete-time corrected transition densities ${q}_{h}^{(0,1)}(x\to z)$ and ${q}_{h}^{(1,2)}(x\to z)$, we obtain the corrected versions of the approximations (5.4).

Approximation 4: ${P}_{L}(T,h)\simeq 1-{F}(1,h,h_{L})\cdot\left[{\lambda}_{L,1}(h)\right]^{T-1},$* where ${F}(1,h,h_{L})$ is given in (4.10) and ${\lambda}_{L,1}(h)$ is the maximal eigenvalue of the integral operator with kernel $K(x,z)=q^{(0,1)}_{h_{L}}(x\to z)$ defined in (5.10).*

Approximation 5: ${P}_{L}(T,h)\simeq 1-{F}(2,h,h_{L})\cdot\left[{\lambda}_{L,2}(h)\right]^{T-2},$* where ${F}(2,h,h_{L})$ is given in (4.11) and ${\lambda}_{L,2}(h)$ is the maximal eigenvalue of the integral operator with kernel $K(x,z)=q^{(1,2)}_{h,L}(x\to z)$ defined in (5.14).

Similarly to $\lambda_{m}(h)$ from (5.3), the maximum eigenvalues ${\lambda}_{L,1}(h)$ and ${\lambda}_{L,2}(h)$ of the operators with kernels $K(x,z)=q^{(0,1)}_{h_{L}}(x\to z)$ and $K(x,z)=q^{(1,2)}_{h,L}(x\to z)$ are simple and positive; the corresponding eigenfunctions $p(x)$ can be chosen as probability densities. Both approximations can be used for any $T>0$.

In numerical examples below we approximate the eigenvalues ${\lambda}_{L,k}(h)$ ($k=1,2$) using the methodology described in Quadrature , p.154. This methodology is based on the Gauss-Legendre discretization of the interval $[-c,h]$, with some large $c>0$, into an $N$-point set $x_{1},\ldots,x_{N}$ (the $x_{i}$’s are the roots of the $N$-th Legendre polynomial on $[-c,h]$), and the use of the Gauss-Legendre weights $w_{i}$ associated with points $x_{i}$; $\lambda_{L,k}(h)$ and $p(x)$ are then approximated by the largest eigenvalue and associated eigenvector of the matrix $D^{1/2}AD^{1/2},$ where $D=\text{diag}({w}_{i})$ and $A_{i,j}=K(x_{i},x_{j})$ with the respective kernel $K(x,z)$. If $N$ is large enough then the resulting approximation to $\lambda_{L,k}(h)$ is arbitrarily accurate. With modern software, computing Approximations 4 and 5 (as well as Approximation 3) with high accuracy takes only milliseconds on a regular laptop.

As discussed in the next section, Approximation 5 is more accurate than Approximation 4, especially for small $h$; the accuracies of Approximations 3 and 5 are very similar. Note also that a version of Approximation 4 has been developed in our previous work AandZ2019 ; this version was based on a different discrete-time approximation (discussed in Section 4.6) of the continuous-time BCP probability $P(T,h)$.

6 Simulation study

6.1 Accuracy of approximations for the BCP ${P}_{L}(T,h)$

In this section we study the quality of Approximations 4 and 5 for the BCP ${P}_{L}(T,h)$ defined in (2.7). Approximation 3 is visually indistinguishable from Approximation 5 and is therefore not plotted (see Table 1). Without loss of generality, $\varepsilon_{j}$ in (1.1) are normal r.v.’s with mean [math] and variance $1$. The style of Fig. 3 is exactly the same as of Fig. 1 and is described in the beginning of Section 4.7. In Fig. 3, the dashed green line corresponds to Approximation 4 and the solid red line corresponds to Approximation 5.

From Figure 3 we see that the performance of Approximations 4 and 5 is very strong even for small $L$. For small $h$, Approximation 5 is more precise than Approximation 4 in view of its better accommodation to the non-Markovian nature of the process $S(t)$.

In Table 1, we display the values of ${\lambda}_{L,1}(h)$, ${\lambda}_{L,2}(h)$ and ${\mu}_{L}(h)$ with $L=20$ for a number of different $h$. From this table, we see only a small difference between ${\lambda}_{L,2}(h)$ and ${\mu}_{L}(h)$; this difference is too small to visually differentiate between Approximations 3 and 5 in Fig. 3.

In Tables 2, 3 and 4 we numerically compare the performance of Approximations 1 and 3 for approximating ${P}_{L}(T,h)$ across different values of $L$ and $h$. Since Approximation 1 relies on Monte-Carlo methods, we present the average over 100 evaluations and denote this by $\bar{x}$. We have also provided values for the standard deviation and maximum and minimum of the 100 runs to illustrate the randomised nature of this approximation. These are denoted by $s$, $Max(x_{i})$ and $Min(x_{i})$ respectively. The values of ${P}_{L}(T,h)$ presented in the tables below are the empirical probabilities of reaching the barrier $h$ obtained by $10^{6}$ simulations. We have not included Approximation 5 in these tables as results are identical to Approximation 3 up to four decimal places.

From Tables 2, 3 and 4 we see that with this choice of $T=100$, the errors of approximating $F_{L}(2,h)$ and $F_{L}(1,h)$ via the ’GenzBretz’ algorithm can accumulate and lead to a fairly significant variation of Approximation 1. This demonstrates the need to average the outcomes of Approximation 1 over a significant number of runs, should one desire an accurate approximation. This may require rather high computational cost and run time, especially if $L$ is large. On the other hand, evaluation of Approximation 3 is practically instantaneous for all $L$. Even for a very small choice of $L=5$, Table 2 shows that Approximation 3 still remains very accurate. As $L$ increases from 5 to $20$, Table 3 shows that the accuracy of Approximation 3 increases. The averaged Approximation 1 is also very accurate but a larger $L$ appears to produce a larger range for $Max(x_{i})$ and $Min(x_{i})$ when $h$ is large; this is seen in Table 4. Note we have not included empirical values of ${P}_{L}(T,h)$ in Table 4 due to the large computational cost.

6.2 Approximation for the BCP in the case of non-normal moving sums

Approximations 3, 4 and 5 remain very accurate when then the original $\varepsilon_{i}$ in (1.1) are not exactly normal. We consider two cases: (a) $\varepsilon_{i}$ are uniform r.v’s on [0,1] and (b) $\varepsilon_{i}$ are Laplace r.v’s with mean zero and scale parameter 1. Simulation results are shown in Figure 4; this figure has the same style as figures in Sections 4.7 and 6.1.

Some selected values used for plots in Figure 4 are:

$h=2,L=20\!:\;\;$ Emp: $0.6045\pm 0.0030\;(0.6123\pm 0.0030)\;[0.5894\pm 0.003];\;\,$ Ap. 4(5): 0.5921(0.6054);

$h=2,L=\!100\!:\;$ Emp: $0.6771\pm 0.0029\;(0.6801\pm 0.0029)\;[0.6722\pm 0.003];\;\,$ Ap. 4(5): 0.6633(0.6775);

$h=3,L=20\!:\;\;$ Emp: $0.0788\pm 0.0017\;(0.0710\pm 0.0016)\;[0.0915\pm 0.002];\;\,$ Ap. 4(5): 0.0777(0.0789);

$h=3,L=\!100\!:\;$ Emp: $0.1039\pm 0.0019\;(0.1033\pm 0.0019)\;[0.1048\pm 0.002];\;\,$ Ap. 4(5): 0.1022(0.1034).

Here we provided means and 95% confidence intervals for the empirical (Emp) values of the BCP ${P}_{L}(T,h)$ (with $T=M/L=10$) computed from 100 000 Monte-Carlo runs of the sequences of the moving sums (1.1) with normal (no brackets), uniform (regular brackets) and Laplace (square brackets) distributions for $\varepsilon_{i}$ in (1.1). Values of Approximations (Ap.) 4 and 5 are also given.

From Figure 4 and associated numbers we can make the following conclusions: (a) the BCP ${P}_{L}(T,h)$ for the case where $\varepsilon_{i}$ in (1.1) are uniform is closer to the case where $\varepsilon_{i}$ are normal, than for the case where $\varepsilon_{i}$ have Laplace distribution; (b) as $L$ increases, the probabilities ${P}_{L}(T,h)$ in the cases of uniform and Laplace distributions of $\varepsilon_{i}$ become closer to the BCP for the case of normal $\varepsilon_{i}$ and hence the approximations to the BCP become more precise; (c) accuracy of Approximation 5 is excellent for the case of normal $\varepsilon_{i}$ and remains very good in the case of uniform $\varepsilon_{i}$; it is also rather good in the case when $\varepsilon_{i}$ have Laplace distribution; (d) Approximation 4 is slightly less accurate than Approximation 5 (and Approximation 3) for the case of normal and uniform $\varepsilon_{i}$ (this is in a full agreement with discussions in Sections 5.2.2 and 6.1); however, Approximation 4 is very simple and can still be considered as rather accurate.

6.3 Approximation for the BCP in the case of moving weighted sums

We have also investigated the performance of Approximation 5 (and 3) after introducing particular weights into (1.1). We explored the following two ways of incorporating weights:

• (i)

$L$ random weights $w_{1},w_{2},\ldots,w_{L}$, with $w_{i}$ i.i.d. uniform on $[0,2]$, are associated with a position in the moving window; this results in the moving weighted sum

[TABLE]

• (ii)

$M+L$ random weights $w_{1},\ldots,w_{M+L}$ are associated with r.v. $\varepsilon_{1},\ldots,\varepsilon_{M+L}$; here $w_{j}$ are i.i.d. uniform r.v’s on [0,2]; this gives the moving weighted sum

[TABLE]

Simulations results are shown in Fig. 5. In both cases, we have repeated simulations 1,000 times and plotted all the curves representing the BCP as functions of $h$ in grey colour and Approximation 5 for the BCP for the non-weighted case (when all weights $w_{j}=1$) as red dashed line. We can see that for both scenarios the Approximation 5 for the BCP in the non-weighted case gives fairly accurate approximation for the weighted BCP. Similar results have been observed for other values of $L$ and $T$.

7 Approximating Average Run Length (ARL)

In this section, we provide approximations to the probability distribution of the moment of time $\tau_{H}(\mathbb{S})\!:=\!\min\{n\!\geq\!0\!:S_{n,L}\!\geq\!H\}$ when the sequence $\mathbb{S}=\{S_{0,L},S_{1,L},\ldots\}$ reaches the threshold $H$ for the first time. Note that $\tau_{H}(\mathbb{S})=\tau_{h}(\mathbb{X})$, where $\tau_{h}(\mathbb{X})\!:=\!\min\{n\!\geq\!0\!:\xi_{n,L}\!\geq\!h\}$ and $\mathbb{X}=\{\xi_{0,L},\xi_{1,L},\ldots\}$. The BCP ${\cal P}_{\mathbb{S}}(M,H,L)$, considered as a function of $M$, is the c.d.f. of this probability distribution: ${\cal P}_{\mathbb{S}}(M,H,L)={\rm Pr}\left(\tau_{H}(\mathbb{S})\leq M\right)$. The average run length (ARL) until $\mathbb{S}$ reaches $H$ for the first time is

[TABLE]

Note that ${\rm ARL}_{H}(\mathbb{S})={\rm ARL}_{h}(\mathbb{X})$. The diffusion approximation to the time moment $\tau_{h}(\mathbb{X})$ is $\tau_{h}(S(t)):=\min\{{t\geq 0}:\;S(t)\geq h\}$, which is the time moment when the process $S(t)$ reaches $h$. The distribution of $\tau_{h}(S(t))$ has the form:

[TABLE]

where $\delta_{0}(ds)$ is the delta-measure concentrated at 0 and

[TABLE]

The function $q(s,h,S(t))/\Phi(h)$, considered as a function of $s$, is a probability density function on $(0,\infty)$ since

[TABLE]

From this, the diffusion approximation for ${\rm ARL}_{H}(\mathbb{X})/L$ is

[TABLE]

The diffusion approximation (7.3) should be corrected for discrete time; otherwise it is poor, especially for small $L$. As shown in Section 6, Approximations 3 and 5 are very accurate approximations for ${P}_{L}(T,h)$ and can be used for all $T>0$. We shall use Approximation 3 to formulate our approximations but note that the use of Approximation 5 would give very similar results.

We define the approximation $\hat{q}(s,h)$ for the probability density function of \tau_{h}(\mathbb{X})$$/L by

[TABLE]

The corresponding approximation for ${\rm ARL}_{h}(\mathbb{X})$ is

[TABLE]

The standard deviation of $\tau_{h}(\mathbb{X})$, denoted $SD(\tau_{h}(\mathbb{X}))$, is approximated by:

[TABLE]

In this paper, we define ARL in terms of the number of random variables $\xi_{n,L}$ rather than number of random variables $\varepsilon_{j}$. This means we have to modify the approximation for ARL of Glaz2012 by subtracting $L$. The standard deviation approximation in Glaz2012 is not altered.

The Glaz approximations for ${\rm ARL}_{h}(\mathbb{X})$ and $SD(\tau_{h}(\mathbb{X}))$ are as follows:

[TABLE]

[TABLE]

where $x={F}_{L}(2,h)/{F}_{L}(1,h)$.

In Tables 5 and 6 we assess the accuracy of the approximations (7.4) and (7.5) and also Glaz approximations (7.6) and (7.7). In these tables, the values of ${\rm ARL}_{h}(\mathbb{X})$ and $SD(\tau_{h}(\mathbb{X}))$ have been calculated using $100,000$ simulations. Since the Glaz approximations rely on Monte Carlo methods, in the tables we have reported value $2s$-confidence intervals computed from 150 evaluations.

Tables 5 and 6 show that the approximations developed in this paper perform strongly and are similar, for small or moderate $h$, to the Glaz approximations. For $h\!\geq\!3$, the Glaz approximation produces rather large uncertainty intervals and the uncertainty quickly deteriorates with the increase of $h$. This is due to the fairly large uncertainty intervals formed by Approximation 1 when approximating ${P}_{L}(T,h)$ with large $h$ and hence small ${P}_{L}(T,h)$, as discussed in Section 6.1. The approximations developed in this paper are deterministic and are much simpler in comparison to the Glaz approximations. Moreover, they do not deteriorate for large $h$.

Acknowledgment

The authors are grateful to the referees for careful reading of the manuscript and useful comments.