Detecting Change-Point via Saddlepoint Approximations

1 Introduction

Consider a sequence of random variables X₁,X₂,⋯,X_n, then the sequence is said to have a change-point at m if X_t for t = 1,2,⋯, m have a common distribution function F₁(x) and X_t for t = m + 1, m + 2,⋯,n have a common distribution function F₂(x), and F₁(x)≠ F₂(x). Change-point analysis is very important and applied in many fields, including macroeconomics[2], environment[3,4], quality control[5–7], medicine[8], geology[9], and etc.

Theoretically, various authors concentrate on the change-point problems. Basically, studies can be divided into two focuses, “change-point in distribution” and “change-point in regression function”. “Change-point in distribution” refers that a pre-change distribution and a post-change distribution exist in the sequence[10]. Studies on this issue include Pettitt[1, 11], Hinkley[12, 13], Talwar and Gentle[14], Carlstein[15], Inclan and Tiao[16], Hawkins and Zamba[35], Kawahara and Sugiyama[18] and etc. “Change-point in regression function”, on the other side, refers the changed regression models, like parameter instability and structural change[19]. Previous studies include Muller[20] and Loader[21], who adapt kernel smoothers, and Davis, et al.[22], Huskova, et al.[23] and Gombay, et al.[24], considering autoregressive models.

The existing methods have several limitations. 1) Large-sample properties have been deduced in most of previous works, while few developments are made on small-sample properties. Thus, these methods usually fail when sample that could be obtained or observed is very small. 2) The computations for existing methods are very complicated and time-consuming. 3) In the case that change-point occurs in the lower and upper tails of the sequence, most of the approaches become ineffective. 4) Even in the common conditions without change-point in the tails, existing methods can only give an estimated interval where the change-point may be located. The limitations mentioned above hinder the applications of these change-point methods in practice.

To settle these problems, we first propose a quantile test of single τ based on the theory of saddlepoint approximations. Further, a “composite quantile test” is constructed to collect the information on each τ of the sequence and give the probability of each location to be a change-point. This detection scheme could effectively solve the problems before, and has several good properties. 1) No assumptions are made about the functional form of the sequence distribution of independent random variables X₁, X₂, ⋯, X_n.2) The computation of the “composite quantile test” is quick and simple. In addition, 3) the detection scheme is very sensitive to the varying information of the sequence. Thus, not only can it give the probable interval to be a change-point, but also can sensitively pinpoint the location. 4) This approach can solve the change-point problem with much smaller sample size than other existing methods. 5) The scheme works well when the change-point occurs in the tails of the sequence. 6) Multiple change-point problems can also be successfully detected.

“Mean shift” problem is an important issue in studies of no matter change-point in distribution or in regression function. Pre-shift and post-shift means are to be estimated concurrently with the change point. Many former studies focus on this problem, like Kokoszka and Leipus[26], Ferger[27] and Shi, et al.[28, 29]. This paper presents the performances of “composite quantile test” in the “mean shift” problem as a representative attempt. Other change-point models are possible to be derived based on this idea.

Suppose X₁, X₂, ⋯, X_m are independent with mean μ, and X_m+1, X_m+2,⋯, X_n are independent with mean μ + δ. These two sequences are also independent with each other, where m is an unknown positive integer. m is the unknown location of a change point X_m. Mean change-point model is formulated as:

where {Z_i, i = 1,2, ⋯, n} is a sequence of random variables with constant mean μ. δ is the unknown amount of change in mean. If δ = 0, there is no change-point at the point of m, otherwise, change-point exists. Under the null hypothesis X₁, X₂,⋯,X_n is an i.i.d. sequence, whereas under the alternative two sequences exist with δ < 0, indicating the mean of X₁, X₂,⋯, X_m greater than the mean of X_m+1, X_m+2,⋯, X_n. In practice, this problem is usually considered as the detection that wether the quality of articles decline[7, 30], and so on. Similarly, δ > 0 is meant that the mean of quantity or quality increases by δ. The objective is to determine wether δ ≠ 0 at m and pinpoint the location of change-point.

2 Quantile Test for Single Quantile

The hypothesis test used here is

• • Null hypothesis H₀: δ = 0, that is, mean shift does not happened in the location of m, m = 2, 3,⋯,n – 1, denoted as η ∈ Θ₀;

• • Alternative hypothesis H₁: δ < 0, and change-point occurs and the mean of the sequence decreases, denoted as η∈Θ0c.

For the case of mean of a sequence increasing, it’s to inverse the order to be the alternative hypothesis above. For the null hypothesis of “no-change” H₀, the test statistic T_τ is constructed as the number of observations among X₁, X₂,⋯, X_m which are greater than the τ quantile of the whole sequence X₁, X₂,⋯, X_n, that is,

where I(·) is the indicative function, which is 1 for X_i > Q_τ(X)and 0 for X_i ≤ Q_τ(X). Q_τ(X) is the τ quantile of the whole sequence.

Proposition 1

Assume that the sequence X₁, X₂,⋯,X_n are i.i.d. with a positive definite variance σ²(X), and X ∈ ℝ. The probability space is (Ω, 𝓐,P). Applying saddlepoint approximation method, the first-order saddlepoint approximation of formula (2.1) is given by

(2.2)P(Tτ≤t)≈1−{n−τ(n+1)}n−τ(n+1)+12{τ(n+1)}τ(n+1)+12(n−m)n−m+12mm+12×{(t+1)(τ(n+1)+t−m+1)}B0(λ^)×{nn+12{n−τ(n+1)−t−1}n−τ(n+1)−t−12{τ(n+1)+t−m+1}τ(n+1)+t−m+32(t+1)t+32(m−t−1)m−t−12{n(t+1)−nm+τ(n+1)m}}−1,

where

λ^=log(t+1)(τ(n+1)+t−m+1)(m−t−1)(n−τ(n+1)−t−1)×(t+1)(n−τ(n+1)−t−1)(m−t−1)(τ(n+1)+t−m+1)(t+1)(n−τ(n+1)−t−1)τ(n+1)+(m−t−1)1/2×{(τ(n+1)+t−m+1)(n−τ(n+1))}−1/2,

and B0(λ^) is Esscher function, given by

B0(λ^)=λ^expλ^22[1−N(λ^)],

where 𝔑(x) is standard normal distribution function,

N(y)=∫−∞yn(s)ds=∫−∞y12πexp−s22ds.

When t≥m[n−τ(n+1)]n, the formula (2.2) is valid. In order to compute the cumulative probability density, P(Tτ ≤ t), let m×τ≥2nn+1.

The procedure of how Proposition 1 is derived is presented in the following. Let ϖ(T_τ) be the statistics relative to T_τ, which helps approximate the distribution of T_τ. The Laplace transform function of ϖ(T_τ) is defined as

φ(θ)=Eexp⁡(θ⋅ϖ(Tτ))=∑ϖ(Tτ)exp⁡(θ⋅ϖ(Tτ))P(ϖ(Tτ)).

A cumulant transform function is defined as

κ(θ)=log⁡φ(θ),

and a characteristic function is defined as

ϕ(t)=Eexp⁡(it⋅ϖ(Tτ)),

where θ ∈ ℝ. The domain of Laplace transform is Θ_θ = {θ ∈ ℝ :φ(θ) < ∞}. The positive definite variance of X₁, X₂, ⋯, X_n enables the cumulant transform κ(θ) strictly convex on Θ_θ[25].

Applications of the saddlepoint to approximate the probability density of a statistic have much more development. Jensen[³²] proposed many methods and ideas to produce a very accurate approximation and we use one of them in our deviation. Let ϖ̅ be the mean of ϖ(T_τ). The saddlepoint approximation to the density of ϖ̅ is

(2.3)fϖ¯(w)≈φ(θ^)exp⁡(−nθ^⋅w)n0(2πΣ(θ^))1/21+1nκ4(θ^)8Σ2(θ^)−5κ3(θ^)224Σ3(θ^),

where κ₃ (·) and κ₄ (·) are the third and forth deviations of κ(·), Σ (θ̂) = σ²(θ̂), and θ̂ is called the saddlepoint satisfying the saddlepoint equation,

κ′(θ^)=w.

The equation is also available in (1.5) of Butler[²⁵], and the density approximation formula (2.3) was derived from Jensen[³²].

For continuous random variables, an saddlepoint approximation to the cumulative density function(CDF) of ϖ̅ is equally straight forward,

(2.4)P(ϖ¯≥w)≈φ(θ^)nexp⁡(−nθ^⋅w)n|θ^|σ(θ^){B0(λ)+sgn(θ^)κ3(θ^)6nΣ3/2(θ^)B3(λ)1nκ4(θ^)24Σ2(θ^)B4(λ)+κ3(θ^)272Σ3(θ^)B6(λ)},

where

Σ(θ^)=σ2(θ^)=κ″(θ^),λ=n|θ^|σ(θ^),B3=−{λ3B0(λ)−(λ3−λ)(2π)−1/2},B4=λ4B0(λ)−(λ4−λ2)(2π)−1/2,B6=λ6B0(λ)−(λ6−λ4+3λ2)(2π)−1/2.

In this case, the saddlepoint θ̂ is obtained by solving κ′(θ̂) = ι(w), i.e., ϵ = 0.

For lattice random variables, use the expanded formula similar to (2.4)

(2.5)P(ϖ¯≥w)≈φ(θ^)nexp⁡(−nθ^⋅w)nσ(θ^){1−exp⁡(−|θ^|)}×{B0(λ)+(n)−1[(σ(θ))−1((|θ^|)−1−exp⁡(−|θ^|)1−exp⁡(−|θ^|))B1(λ)+κ3(θ^)6Σ3/2(θ^)sgn(θ^)B3(λ)]+1n[(exp⁡(−|θ^|)1−exp⁡(−|θ^|)12−(|θ^|)−1+[exp⁡(−|θ^|)1−exp⁡(−|θ^|)]2)B2(λ)Σ(θ^)+1σ(θ^)(|θ^|)−1−exp⁡(−|θ^|)1−exp⁡(−|θ^|)×κ3(θ^)6Σ3/2(θ^)sgn(θ^)B4(λ)+κ4(θ^)24Σ2(θ^)B4(λ)+κ3(θ^)272Σ3(θ^)B6(λ)]},

where

B1=−λ{B0(λ)−(2π)−1/2},B2=λ2{B0(λ)−(2π)−1/2}.

The saddlepoint θ̂ in this case is obtained by solving κ′(θ̂) = ι(w) –ϵ/2, i.e., ϵ= 0. Without loss of generality, let ϵ = 1/n.

The distribution of T_τ can be represented as a conditional distribution. Let ϑ_i = 0 for i ≤ τ (n + 1) and ϑ_i = 1 for i > τ(n + 1). The cumulative density function of T_τ is

(2.6)FTτ(t)=P∑1nϑiUi≤t|∑1nUi=m,t=0,1,⋯,m,

where U₁,U₂,⋯, U_n are i.i.d. with P(Ui=1)=P(Ui=0)=12. Therefore, the statistic ϖ(Tτ)=(∑1nϑiUi,∑1nUi).

In order to approximate the above formula in which θ ∈ ℝ², it is suitable for formula (2.5) to be expanded to the formula (2.2.15) in Jensen[32]

(2.7)P(ϖ1¯(Tτ)≥w1∣ϖ2¯(Tτ)=w2)=exp−r22|Σ22(θ~)|n|Σ(θ^)|1/21|θ^1|B0λ^+OB0(λ^)1+λ^n−1/2,

where r2=2n(θ^−θ~)⋅s−κ(θ^)+κ(θ~) and ϖ̅₁ and ϖ̅₂ are defined as in the Appendix. The procedure of Proposition 1 is given. And details of the deviation from formula (2.6) to formula (2.2) for our case are in the Appendix.

According to the Proposition 1, the probability distribution of T_τ can be obtained via P(T_τ = t) = P(T_τ ≤ t) – P(T_τ ≤ t – 1), which is denoted as fTτ(t)for t = 0,1,⋯,m.

If change-point does not exist at the location m of the sequence, the distribution of T_τ is exactly approximated by saddlepoint approximation method, the formula (2.2). In order to exactly compute its value, formula (2.2) is transformed to

(2.8)P(Tτ≤t)≈1−n−τ(n+1)n−τ(n+1)−t−1n−τ(n+1)−t−12×τ(n+1)τ(n+1)+t−m+1τ(n+1)+12×n−τ(n+1)τ(n+1)+t−m+1t+1×n−mnn+12×τ(n+1)+t−m+1n−mm×mm−t−1m+12×m−t−1t+1t+1×(t+1)−12×{(t+1)[τ(n+1)+t−m+1]B0(λ^)}n(t+1)−nm+τ(n+1)m.

The distribution of T_τ approximated by saddlepoint approximation is under the null hypothesis of “no-change”.

For each location m, the significant level is set to be α = 0.05. The criterion for hypothesis test is the value of T_τ, named as t₀(τ) when the accumulative density distribution of T_τ equals 0.95. That is

t0(τ)=arg⁡{∑t=0t0(τ)fTτ(t)=0.95}.

Here, t₀(τ) is not an integer necessarily, and

∑t=0int(t0(τ))fTτ(t)≤0.95≤∑t=0int(t0(τ))+1fTτ(t),

where int(t) denotes the largest integer smaller than t. The plots of T_τ when m = 10, m = 50 and m = 90 are presented in Figure 1. The curve of T_τ intersects 0.95 at t₀(τ). When k number of τ are considered, for each location m = 4, 5,⋯, n – 1, the threshold matrix of t₀(τ) by {τ_j, j = 1, 2,⋯, k} are plotted in Figure 2.

Figure 1

FTτ: (a) describes FTτ when m = 10 which intersects line 0.95 at t₀(τ); (b) describes FTτ when m = 50 which intersects line 0.95 at t₀(τ); (c) describes FTτwhen m = 90 which intersects line 0.95 at t₀(τ)

Figure 2

The matrix of t₀(τ): (a) describes the surface of t₀(τ) in the large sample case (n = 100); (b) is the surface of t₀(τ) in the small sample case (n = 20)

When “mean shift” happens, the distribution of T_τ of data would largely deviate from the distribution under the null hypothesis. If m is given, the T_τ for target sequence is calculated as t₁(τ) and compared with t₀(τ). If t₁(τ) > t₀(τ), the null hypothesis is rejected and change point exists at m location of the sequence X₁, X₂,⋯, X_n. If t₁(τ) < t₀(τ) instead, null hypothesis cannot be rejected and the change-point cannot be determined. By shifting the value of m from 2 to n – 1, we can detect wether change point exists in the target sequence and pinpoint the location. Therefore, a quantile test that regects H₀ if t1X(τ)>t0(τ),is a test with rejection region R,{tX:t1X(τ)>t0(τ)}.

To compute the power of hypothesis test, under the alternative hypothesis, the distribution of T_τ is obtained by applying the bootstrap method[33]. The basic steps proceed as follows. Assuming that the random sample of size n is drawn from an unspecified probability distribution, Ξ.

• Step 1. Placing a probability of 1/n at each point, x₁, X₂,⋯,x_n. Each sample’s element has the same probability to be drawn.

• Step 2. From the sample x₁, X₂,⋯,x_n, draw a random sample of size d with replacement.

• Step 3. Calculate the statistic of interest, T_τ. For resample, yield Tτ∗=Tτ×nd.

• Step 4. Repeat Steps 2 and 3 N (N > 1000) times.

• Step 5. Construct the relative frequency through the N number of Tτ∗ by placing a probability of 1/N at each point, Tτ∗1,Tτ∗2,⋯,Tτ∗N.The distribution obtained is the bootstrapped estimate of the sampling distribution of T_τ.

The bootstrapped sampling distribution of T_τ is denoted as fTτboot(t). The power of hypothesis test, denoted as β_τ(η), can be calculated by

(2.9)βτ(η)=∑t=[int(t0(τ))+1]mfTτ(t),ifη∈Θ0,1−∑t=0int(t0(τ))fTτboot(t),ifη∈Θ0c,

where, ∑t=[int(t0(τ))+1]mfTτ(t) is the probability of Type I Error the test makes. ∑t=0int(t0(τ))fTτint(t)

is the probability of Type II Error. A good hypothesis test has power near 1 for most η∈Θ0c and near 0 for most η ∈ Θ₀. Due to

supη∈Θ0βτ(η)=α,

quantile test is a size α test.

3 Composite Quantile Test

Instead of adopting one τ quantile (median or other τ) for change-point detection, which will drop much of the valuable information, we try to collect information from several τ of the sequence X₁, X₂,⋯,X_n. A complete detection scheme is constructed to compute the probability of each location m to be a change point. The information across k uniformly spaced quantiles {τ_j, j = 1, 2,⋯, k} is aggregated to maximise the efficiency of making inference, where

The probability of one location m to be a change-point is calculated by

where

and

where

t₀(τ_j) is the theoretical criterion for hypothesis test at significant level 95% when m and τ_j are given. t₁(τ_j) denotes the number of points in series x₁,x₂,⋯,x_m which are larger than τ_j quantile of the whole sequence x₁, X₂,⋯, X_n. The distance between t₁(τ_j) and t₀(τ_j) indicates the deviation of the distribution of proposed statistics in the target sequence from theoretical criterion in the null hypothesis. Euclidean distance is used here.

The proposed “composite quantile regression” can be considered as a special procedure to combine the results of each single quantile test together by a varying weight. For given j, the distance between t₁(τ_j) and t₀(τ_j) determines the weight of the information from τ_j quantile test. In this regard, “composite quantile regression” is a series of hypothesis tests to the greatest extent collecting the valuable information from the whole sequence, instead of limited data points. This is a more sensitive way to detect the change-point location without restricting the forms of the sequence distribution.

4 Simulation

Monte Carlo simulation is applied to present the performances of the proposed single quantile test and composite quantile test. Three other detecting methods-single tests, CUSUM[5, 34], Student-t statistic[35, 36] and Mann-Whitney statistics[37, 38] are adopted as benchmarks. Three scenarios with different degree of “mean shift” are considered in both large and small sample cases. Multiple change-points situation is also tested.

4.1 Single Quantile Test and Power

We analyze the situations of large size sample n = 100 and small size sample n = 20. Given the simulated sequences with a change-point located at m = 4, 5,⋯, n – 1, the rate of successful detection in the Monte Carlo simulation is defined as S(τ)=N∗N,where N* is the number of successful detection, and N (N = 10⁴) denotes the times of Monte Carlo experiments. The generated sequence X_t for t = 1,2,⋯, m is distributed by Norm(μ₁, σ) and sequence X_t for t = m+1,⋯ , n is distributed by Norm(μ₂, σ), where μ₁ = μ₂+δ ·σ. In the previous empirical studies, the mean shift is usually larger than one σ. Thus, δ is set as —1, —1.5 and —2 in the scenarios respectively indicating the degree of “mean shift”. For simplicity, σ = 1. The results of S (τ) are presented in Figure 3. Most of location m show great success rate S(τ) except for extreme conditions of m ≥n — 3. This is a very good performance in comparison with studies before. Quantile test possesses the advantage that it can find the location of change point exactly.

Figure 3

S (τ): S (τ) in cases of n = 100 and n = 20 with δ = — 1, —1.5, —2, respectively

Sequences without change-points are generated to test the probability of Type I Error. The probability of Type I Error is defined as W(τ)=N^N,where N is the times of Monte Carlo experiments which make Type I Error. As reported in Figure 4, some large W(τ) are at the lower tail of sequence when small τ is adopted. At other location W(τ) are close to zero. Even in the extreme case of small sample (n = 20), for more than half of locations, W(τ) still show very small values. However, small quantile test fail when sample size is extreme small. For large sample size, the Type I Errors of CUSUM, Student-t and Mann-Whitney are 0.9185, 0.0522 and 0.0556, respectively. The corresponding values are 0.9908, 0.0482 and 0.051 for extreme small sample size. CUSUM fails, while Student-t and Mann-Whitney still perform well. The Type I Errors of Student-t and Mann-Whitney are close to the means of Type I Errors of 0.5, 0.6 and 0.7 quantiles obtained from Figure 4, which are 0.0580, 0.0549 and 0.0538 for large sample size, and 0.0639, 0.490 and 0.0573 for extreme small size.

Figure 4

W(τ): W(τ) in cases of n = 100 and n = 20

The power of quantile test is defined as Formula (2.9). Bootstrap resampling is introduced with repeating times N = 10⁴. Figure 5 presents the results of β_τ in three scenarios δ = -1,-1.5,-2 with n = 100 and n = 20. The mean of β_τ is reported in Table 1. β_τ, increasing with higher degree of “mean shift”, shows great performance of the quantile test by single τ. Figure 6 shows the powers of CUSUM, Student-t, Mann-Whitney, and single quantile tests using 0.5, 0.6 and 0.7 quantiles. In the simulation, the sequences are generated with one change point at location m. Single quantile test can exactly find this location, while CUSUM, Student-t and Mann-Whitney can only detect the existence of change point but cannot always give the exact location. In the calculation of power for them in Figure 6, we do not consider whether they find the exact location. The high power of CUSUM is accompanied by extremely high Type I Error. Therefore, it makes little sense to compare the power of CUSUM with other methods. The powers of Student-t and Mann-Whitney are close. For both large and small sample size, the power of single quantile tests is almost always higher than these two methods. Single quantile tests are more powerful than other methods when change point existing in the first half part of sequence. As showed in Table 2, the powers of Student-t and Mann-Whitney are very weak when we consider whether they find the exact location of change point. The power of CUSUM can reach 67 percent and the powers of Student-t and Mann-Whitney are smaller than 20 percent. They are far smaller than that of single quantile test as showed in Table 1.

Table 1

β̅_τ in three scenarios for both n = 100 and n = 20

δ	−1	−1.5	−2
n = 100	0.9088	0.9763	0.9895
n = 20	0.7941	0.9575	0.9577

Table 2

Summary of statistics about the powers of three methods under the requirement of finding exact location of change point

n = 100	δ = −1	Max	Min	Median	Mean
	CUSUM	0.2983	0	0.2245	0.2006
	Student-t	0.1113	0	0.1014	0.0895
	Mann-Whitney	0.1160	0	0.1013	0.0884
	δ = −1.5	Max	Min	Median	Mean
	CUSUM	0.4908	0	0.3737	0.3353
	Student-t	0.1319	0	0.1246	0.1174
	Mann-Whitney	0.1496	0	0.1218	0.1143
	δ = −2	Max	Min	Median	Mean
	CUSUM	0.6538	0.0087	0.5085	0.4559
	Student-t	0.1205	0	0.1109	0.1090
	Mann-Whitney	0.1616	0	0.1129	0.1071
n = 20	δ = −1	Max	Min	Median	Mean
	CUSUM	0.3281	0	0.2783	0.2468
	Student-t	0.0451	0	0.0374	0.0325
	Mann-Whitney	0.0526	0	0.0414	0.0337
	δ = −1.5	Max	Min	Median	Mean
	CUSUM	0.5153	0	0.4413	0.3906
	Student-t	0.0906	0	0.0776	0.0667
	Mann-Whitney	0.1051	0	0.0793	0.0655
	δ = −2	Max	Min	Median	Mean
	CUSUM	0.6726	0.2430	0.5747	0.5373
	Student-t	0.1079	0	0.0948	0.0869
	Mann-Whitney	0.1365	0	0.0941	0.0845

Figure 5

β_τ : β_τ in cases of n = 100 and n = 20 under three scenarios δ = –1, –1.5, –2, respectively

Figure 6

Power comparison of CUSUM, Student-t, Mann-Whitney and single quantile tests in cases of n = 100 and n = 20 under three scenarios δ = −1, -1.5, -2, respectively

4.2 Composite Quantile Test and Multiple Change-Points Situation

Based on quantile test of single τ, “composite quantile test” is constructed as Formula (3.1). In this part, single change-point is set at location m̂ of the generated sequence. m̂ shifts from 4 to n - 1. For each given m̂ , the probability P(m) of each location m to be a change-point is calculated by collecting the information of each single quantile test. The results are presented in Figure 7 both for n = 100 and n = 20. A series of peaks can be observed on the diagonal of the surface, where the simulated change-points are actually located. That means, by applying “composite quantile test”, the highest probabilities pinpoint the actual change-point locations. Similar patterns can also be observed for small sample cases. The ridge-like surfaces of Figure 7 reflect the sensitive performance of “composite quantile test”. Even in the extreme conditions that generated change-point occurs in the tails of the sequence, the proposed method still responds sensitively and gives the highest probability to the right locations.

Figure 7

P(m): P(m) in cases of n = 100 and n = 20 under three scenarios δ = –1,–1.5,–2, respectively

Multiple change-points cases are further investigated to confirm the advantage of the proposed method. For two change-points condition, the sequences which have two change-points are generated. The two change-points at first are on the two tails of the sequence, one on the lower and the other on the upper. The two points approach to each other and eventually join together. The degree of “mean shift” is set as δ = –1. The length of sequence is set as n = 150. It is the same with three change-points condition, except for a third change-point which is constantly in the middle of a n = 200 length sequence. The results of two and three change-points situations are plotted in Figure 8.

Two ridges are observed in Figure 7(a) and (b), exactly positioning the locations of change-points. Similar patterns can be found for three change-points situation, except that the change-point in the middle is more emphasized by the method when two other points are on the tails. Multiple change-points situation indicates that the “composite quantile test” works robustly in dealing with unfavorable data conditions.

Figure 8

P(m): P(m) for two and three change-point δ = –1

5 Real Data

Two sets of data are analyzed. The first one is a much analysed[39–43] data set of intervals between British coal-mining disasters which resulted in 10 or more fatalities during the 112-year period 1851–1962, which was gathered by Maguire et al.[44], extended and corrected by Jarrett[45]. And Carlin, et al.[46, 47] converted Jarrett’s table into yearly intervals. The second data is a 60-year (1912–1971) annual series of temperatures from Libby, Montana. The documented change point is 1938 (resulting from a change in latitude and elevation)[48]. The plots of these two time series are showed in Figure 9(a) and (b). As presented in Figure 9(c) and (d), the probability of each location to be a change point is computed by “composite quantile test”.

The aim is to find a location m* cutting off the sequence making F₁(x) ≠ F₂(x). F₁(x) and F₂(x) are common distribution function for X_t, t = 1, 2, · · · , m and X_t, t = m+ 1, m + 2, ⋯ , n, respectively. For the first data, two peaks with very near locations are identified at m = 36 and m = 39. The probabilities of these two points are also very close. The standard deviation is σ_L = 1.6475, which is the empirical basis to choose the value of δ in Monte Carlo experiments. In detail, at location m = 36, μ¯136=3.250 is the mean of the sequence X₁, X₂, · · · , X_m and μ¯236=0.9737 is the mean of the sequence X_m+1, X_m+2, · · · , X_n. Thus the mean difference is μ¯136−μ¯236=1.382σL. At location m = 39, on the other hand, μ¯139=3.1538andμ¯239=0.9315. And μ¯139−μ¯239=1.349σL. It should be noticed that the mean difference at m does not necessarily lead to a high P(m). Change-point m* makes the sample points before and after m* coming from different distributions. According to the result of composite quantile test, two special points are pinpointed and more likely to be change-points. In previous studies, a probable interval is usually given without exact location to be a change-point. For example, Carlin, et al.[49] gave the probable interval of change-point 1889–1892, which is m = 39 to m = 42 in our context. The point m = 36 is not detected and the interval m = 39 to m = 42 is also too wide to be accurate. Student-t detect only one change point at m = 36 and also Mann-Whitney detect a change point at m = 41. In this regard, “composite quantile test” responds more sensitively to change “tendency” of the sequence than methods before.

Figure 9

Real data: (a) and (b) describe the time series of dataset 1 and 2. (c) and (d) describe P(m), the probability of being change-point

In the second data, one peak is observed at m = 19, which is the year of 1930. The standard deviation is σ_S = 0.8769. The mean of points before this peak is μ̄₁ = 6.5626 and the mean of points after that point is μ̄₂ = 7.3521. The mean difference is μ̄₂ – μ̄₁ = 0.9σ_S. The proposed method senses a tendency that the distribution of sample is going to change. Student-t and Mann-Whitney detect the same change point location. This data is also studies by Reeves, et al[48]. They applied eight methods, where the standard normal homogeneity (SNH) test and nonparametric SNH procedure detect a change point in 1930, Sawa’s Bayes criteria, two-phase regression with a common trend and modified vincent’s method detect a change point in 1947, and Akaike’s information criteria, two-phase regression (LR) and a modified LR method detect a change point in 1945. Our method detect the same change point location as SNH test and nonparametric SNH procedure. This location differs from the documented change point. This discrepancy maybe caused by other undocumented change points in the series or changed information happening before the location of true change point.

6 Conclusion

“Mean shift” problem in change-point studies is discussed in this paper. Saddlepoint approximation method is applied to approximate the cumulative density function of statistic T_τ for null hypothesis of “no-change”, and a quantile test of single τ is proposed. In order to capture the overall tendency information of the sequence, we construct “composite quantile test” to calculate the probability of each location to be a change-point. As showed in Monte Carlo experiments and real data analysis, the proposed detection scheme performs better than existing methods. The types of functional form of sequence distribution are not restrained. The sample sizes (n = 100, 20) considered in our context are much smaller than previous studies. In the cases that change-points at lower or upper tails of sequence, the proposed test still works well. In addition, the robustness of this scheme is confirmed by multiple change-points conditions. The ideas behind this paper are extendable to other complex issues in change-point studies. In addition, the extension of this testing to two-sided is left as a future research topic.