# Data Segmentation for Time Series Based on a General Moving Sum Approach

Claudia Kirch1,2 and Kerstin Reckruehm1

###### Abstract

In this paper we propose new methodology for the data segmentation, also known
as multiple change point problem, in a general framework including classic
mean change scenarios, changes in linear regression but also changes in the
time series structure such as in the parameters of Poisson-autoregressive time
series. In particular, we derive a general theory based on estimating
equations proving consistency for the number of change points as well as rates
of convergence for the estimators of the locations of the change points. More
precisely, two different types of MOSUM (moving sum) statistics are
considered: A MOSUM-Wald statistic based on differences of local estimators
and a MOSUM-score statistic based on a global estimator. The latter is usually
computationally less involved in particular in non-linear problems where no
closed form of the estimator is known such that numerical methods are
required. Finally, we evaluate the methodology by means of simulated data as
well as using some geophysical well-log data.

11footnotetext: Department of Mathematics, Otto-von-Guericke University;
Magdeburg, Germany. Supported by Deutsche Forschungsgemeinschaft - 314838170,
GRK 2297 MathCoRe.22footnotetext: E-mail<EMAIL_ADDRESS>
##  1 Introduction

Data segmentation or multiple change point estimation is a current topic in
statistics and machine learning and comprises problems in a wide range of
fields like finance, quality control, medicine or climate. For example, Braun
et al., (2000) apply a technique of multiple change point detection on DNA
sequences. Aggarwal et al., (1999) focus on finding structural breaks in the
volatility of stock market returns, while Killick et al., (2010) and Killick
et al., (2012) give an interesting application to oceanography by detecting
changes in the variance of time series for wave heights.

Early literature on change point methodology focused on at-most-one-change
testing and the corresponding estimator for a single change point (see e.g.
Csörgö & Horváth, (1997)). Moving on to more complex data structures this is
still an active field of research, see e.g. Aue & Horváth, (2013), Horváth &
Rice, (2014), Cho & Kirch, (2022+). At the same time, the theoretic
investigation of the multiple change situation with all its additional
challenges has become increasingly popular, where most papers deal with the
detection of multiple changes in the mean; see for example the recent survey
articles Fearnhead & Rigaill, (2020), Cho & Kirch, (2022+).

One approach to the multiple change point in the mean problem with favorable
properties is based on moving sum (MOSUM) statistics first considered in
detail by Eichinger & Kirch, (2018). Corresponding tests had previously been
considered by Bauer & Hackl, (1980), Hušková, (1990), Chu et al., (1995) and
Hušková & Slabỳ, (2001). Cho & Kirch, (2021+) refine the procedure by an
additional post-processing step based on an information criterion. An
efficient implementation of the procedure in an R-package is detailed in Meier
et al., (2021). Yau & Zhao, (2016) use moving sum statistics in combination
with an information criterion for estimating the change points in a linear
autoregressive setting. In the context of renewal processes moving sum
methodology was proposed by Messer et al., (2014) as well as Kirch & Klein,
(2021+).

In this paper, we consider a general framework based on estimating equations
similar to what has been considered in the at-most-one-change situation by
Kirch & Tadjuidje Kamgaing, (2016b) as well as in a sequential change point
context by Kirch & Kamgaing, (2015); Kirch & Weber, (2018). This widely
defined scope contains linear and non-linear (auto-)regressive time series
including such based on neural network approximations. We propose two types of
MOSUM (moving sum) statistics and derive corresponding consistency results for
estimators obtained from one bandwidth. Based on these results post-processing
methods similar to Cho & Kirch, (2021+) can be considered in future work.

The paper is organized as follows: In Section 2 we explain in detail the two
data segmentation algorithms that we propose. The discussion in this paper
occurs in a general framework based on estimating functions where examples are
given in Subsection 2.3. In Section 3 we prove consistency of the data
segmentation algorithms including localization rates for one of the
procedures. Due to the general framework the consistency results are obtained
based on some high-level assumptions. In Section 4 we prove these high-level
assumptions for sufficiently smooth estimating functions under standard moment
conditions. In Section 5 the empirical performance is investigated by means of
a simulation study and using a geophysical well-data set before we give some
conclusions in Section 6. The proofs can be found in an appendix.

##  2 Data Segmentation Based on Moving Sum Statistics

###  2.1 MOSUM-Wald and MOSUM-Score Statistics

We consider the following general setting

$\displaystyle X_{t}=\sum_{j=1}^{q+1}X_{t}^{(j)}\,1_{\\{k_{j-1,n}<t\leq
k_{j,n}\\}},\quad k_{0,n}=0,\;k_{q+1,n}=n.$ (1)

The sequences $\\{X_{t}^{(j)}:{t\geq 1}\\},j=1,\ldots,q+1$, are assumed to be
stationary, where consecutive sequences are distributionally different. Then,
$q$ denotes the number of structural breaks and $k_{1,n},\ldots,k_{q,n}$
denote the change points. For simplicity we assume that $q$ is fixed but all
arguments are valid for a sequence $q_{n}$ as long as $q_{n}$ is bounded.
Under some additional assumptions it is also possible to relax the boundedness
assumption on the number of change points (see e.g. Kirch & Klein, (2021+) or
Section E.3 in Cho & Kirch, (2021+)).

While the construction of the below statistics is based on a given parametric
model, we do not require the observed time series to follow that model.
Instead we explicitly allow for model misspecification in the theoretic
analysis showing that under misspecification changes are detected if
consecutive time series have different best approximating model parameters (in
the sense of Assumption 2.2).

The test statistic is related to estimators of the model parameters based on
estimating functions $\boldsymbol{H}$ also known as M-estimators. These are
obtained as the solution $\boldsymbol{\widehat{\theta}}_{a,b}$ of the
estimating equation system
$\sum_{i=a}^{b}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\theta})\overset{!}{=}\boldsymbol{0}$
for a suitable choice of $\mathbb{X}_{i}$, which can be equal to the
observations, a tuple of response and explanatory variables or (for
autoregressive models) can include lagged observations; see Section 2.3 for
some examples. The estimating function $\boldsymbol{H}$ is vector-valued for
the estimation of multidimensional parameter vectors.

As we adopt a general framework in this paper and explicitly allow for
misspecification, we need to make some high-level assumptions on the
underlying time series. We show their validity under some moment conditions
for sufficiently smooth estimating functions in Section 4.

######  Assumption 2.1.

  1. (a)

For each segment $j$, let $\\{\mathbb{X}_{t}^{(j)}:t\geq 1\\}$ be stationary.

  2. (b)

For a given $\boldsymbol{\theta}$ specified below let
$\boldsymbol{S}_{j}(k,\boldsymbol{{\theta}})=\sum_{i=1}^{k}\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{{\theta}})$
fulfill a strong invariance principle for all $j=1,\ldots,q+1$, i.e. possibly
after changing the probability space there exists some $\nu>0$, a
$p$-dimensional standard Wiener process $\\{\boldsymbol{W}(k):k\geq 0\\}$ and
a symmetric positive definite long-run covariance matrix
$\boldsymbol{\Sigma}_{(j)}(\boldsymbol{\theta})$ such that as $k\to\infty$

$\displaystyle\left\lVert\boldsymbol{S}_{j}(k,\boldsymbol{{\theta}})-E(\boldsymbol{S}_{j}(k,\boldsymbol{{\theta}}))-\boldsymbol{\Sigma}_{(j)}(\boldsymbol{\theta})^{1/2}\,\boldsymbol{W}(k)\right\rVert=O(k^{1/(2+\nu)})\;\;a.s.$

In the literature invariance principles as in (b) have been derived for many
time series – compare also Theorem 4.1 below, where typically the parameter
$\nu$ depends on the number of existing moments.

The following assumption is a Bahadur representation for M-estimators as for
example derived by He et al., (1996) but uniformly in $k$. Effectively, this
is a linearisation of the estimator which is typically used to derive
asymptotic normality. In our case, the representation will also be used to
make a connection between the MOSUM-Wald and the MOSUM-score statistics in the
proof.

######  Assumption 2.2.

There exists a regular matrix $\boldsymbol{V}_{(j)}$ for any $j=1,\ldots,q+1$
such that

$\displaystyle\max_{k_{j-1,n}<k\leq
k_{j,n}-G}\left\lVert\sqrt{\frac{G}{2}}\boldsymbol{V}_{(j)}\left(\boldsymbol{\theta}_{j}-\boldsymbol{\widehat{\theta}}_{k+1,k+G}\right)-\frac{1}{\sqrt{2G}}\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}^{(j)}_{i},\boldsymbol{\theta}_{j})\right\rVert$
$\displaystyle=o_{P}\left((\log(n/G))^{-1/2}\right),$

where $\boldsymbol{\theta}_{j}$ fulfills
$E\left(\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\theta}_{j})\right)=\boldsymbol{0}$
(identifiably unique).

If the parametric assumption underlying the M-estimator is correct, then
$\boldsymbol{\theta}_{j}$ is the true underlying parameter for the $j$th
stationary sequence. Otherwise, this is the best-approximating parameter in
the sense induced by the estimating function. Furthermore, in case of
differentiable estimating functions, it usually holds
$\boldsymbol{V}_{(j)}=E\left(\nabla\boldsymbol{H}(\mathbb{X}^{(j)}_{1},\boldsymbol{\theta}_{j})^{T}\right)$;
see Regularity Condition 4.2 (b).

Moreover, we define

$\displaystyle\boldsymbol{V}_{k}=\boldsymbol{V}_{(j)},\;\boldsymbol{\Sigma}_{k}(\boldsymbol{\theta})=\boldsymbol{\Sigma}_{(j)}(\boldsymbol{\theta}),\;\boldsymbol{\Sigma}_{k}=\boldsymbol{\Sigma}_{(j)}(\boldsymbol{\theta}_{j})\quad\text{for
}k_{j-1,n}<k\leq k_{j,n},$
$\displaystyle\boldsymbol{\Gamma}_{k}=\boldsymbol{V}_{k}^{-1}\boldsymbol{\Sigma}_{k}\left(\boldsymbol{V}_{k}^{-1}\right)^{T},\qquad\boldsymbol{\Gamma}_{(j)}=\boldsymbol{\Gamma}_{k_{j,n}},$
(2)

where $\boldsymbol{\Gamma}_{k}$ is the asymptotic (long-run) covariance matrix
for the estimators from the $j$th stationary segment under the above
assumptions.

A first intuitive approach to the multiple change problem is based on the
following MOSUM-Wald statistic that uses weighted differences of moving
estimators $\boldsymbol{\widehat{\theta}}_{k+1,k+G}$ and
$\boldsymbol{\widehat{\theta}}_{k-G+1,k}$ of the unknown quantity $\theta$
based on the stretch of data $X_{k+1},\ldots,X_{k+G}$ respectively
$X_{k-G+1},\ldots,X_{k}$. Similarly to the Mahalanobis distance the weighting
is done with the asymptotic (long-run) covariance matrix, so that we obtain
for $k=G,\ldots,n-G$

$\displaystyle
T^{(1)}_{k,n}(G)=T^{(1)}_{k,n}(G;\boldsymbol{\Gamma}_{k})=\frac{\sqrt{G}}{\sqrt{2}}\left\lVert\boldsymbol{\Gamma}_{k}^{-1/2}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\widehat{\theta}}_{k-G+1,k}\right)\right\rVert,$
(3)

where the bandwidth $G=G_{n}$ is a tuning parameter that determines the length
of the moving window and $\|\cdot\|$ denotes the Euclidian norm.

In this paper, we consider bandwidths $G$ that are of smaller order than the
sample length $n$. The procedure detects changes that are further apart than
$2G$ such that sublinear changes, whose distance diverges strictly slower than
the sample size, can be detected.

######  Assumption 2.3.

  1. (a)

For $\nu>0$ (specified in Assumption 2.1 (b)) let

$\displaystyle\frac{n}{G}\rightarrow\infty\;\;\text{and}\;\;\frac{n^{\frac{2}{2+\nu}}\log(n)}{G}\rightarrow
0\;\;\text{for}\;\;n\rightarrow\infty.$

  2. (b)

The minimal distance between two neighbouring change points is asymptotically
larger than $2G$ in the sense of

$\lim\inf_{n\to\infty}\min_{j=1,\ldots,q+1}(k_{j,n}-k_{j-1,n})/G>2.$

Close to a change point the above Wald statistic (3) can be expected to be
large such that it can be used to find changes in the parameter vector
$\boldsymbol{\theta}$. However, using the Wald statistic has one major
drawback: Calculating two estimates for each time point, i.e. $2(n-2G)$ in
total, can be computationally challenging in situations where numerical
methods need to be applied as e.g. in non-linear models such as Poisson
autoregressive models. Furthermore, in situations with many (almost) zeroes
(corresponding to many local optima) such as e.g. in a neural-network-
autoregressive situation as in Kirch & Kamgaing, (2012) estimators can be far
apart even though the corresponding stochastic processes are almost identical.
While this is explicitly excluded for our theoretic results by assuming
identifiability, this can easily lead to problems in applications.

In order to avoid these problems and to reduce the computational complexity
and corresponding numerical challenges we consider MOSUM-score statistics
where local parameter estimators are replaced by a global inspection parameter
$\boldsymbol{\widetilde{\theta}}$, which can be fixed or an estimator based on
the same data. The MOSUM-score statistic is based on the following differences
between moving sums of the estimating function at the inspection parameter for
$k=G,\ldots,n-G$

$\displaystyle\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)=\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\widetilde{\theta}})-\sum_{i=k-G+1}^{k}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\widetilde{\theta}}).$

Effectively, this converts a general multiple parameter change problem to a
multiple mean change problem of the transformed sequence
$\\{\boldsymbol{H}(\mathbb{X}_{t},\boldsymbol{\widetilde{\theta}})\\}_{t\geq
1}$. Consequently, a multivariate version of the (univariate) mean-MOSUM
statistic investigated by Eichinger & Kirch, (2018) can be applied, leading to
the following MOSUM-score statistic:

$\displaystyle
T^{(2)}_{k,n}(G)=T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}})=T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}};\boldsymbol{\Sigma}_{k}(\boldsymbol{\widetilde{\theta}}))=\frac{1}{\sqrt{2G}}\left\lVert\boldsymbol{\Sigma}_{k}(\boldsymbol{\widetilde{\theta}})^{-1/2}\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)\right\rVert.$
(4)

If we use data-dependent inspection parameters
${\boldsymbol{\widetilde{\theta}}}_{n}={\boldsymbol{\widetilde{\theta}}}_{n}(\mathbb{X}_{1},\ldots,\mathbb{X}_{n})$,
then we need the following additional assumptions for the MOSUM-score
procedure:

######  Assumption 2.4.

There exists $\boldsymbol{\widetilde{\theta}}$ such that for any
$j=1,\ldots,q$

$\displaystyle(i)\quad\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}\left\lVert\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k)-\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)\right\rVert=o_{P}\left(\sqrt{\frac{G}{\log(n/G)}}\right),$
$\displaystyle(ii)\quad\max_{|k-k_{j,n}|<G}\left\lVert\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k)-\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)\right\rVert=o_{P}\left(\sqrt{G\,\log(n/G)}\right).$

This assumption holds under smoothness and mixing assumptions for any
$\sqrt{n}$-consistent estimator $\boldsymbol{\widetilde{\theta}}_{n}$ (see
Theorem 4.4 below). In particular, $\boldsymbol{\widetilde{\theta}}_{n}$
obtained from the corresponding estimating equations based on the stretch of
data $X_{a},\ldots,X_{b}$ are $\sqrt{n}$-consistent for
$\boldsymbol{\widetilde{\theta}}$, the best-approximating parameter in the
sense of Theorem 4.3. While for the latter approach some theoretical
guarantees towards detectability are given in Section 3.2.1, from a
computational perspective it might be more efficient to use different
estimators (for an empirical example, see Section 5.1).

The latter MOSUM-score approach avoids the numerical drawbacks of the MOSUM-
Wald statistic but can only detect changes for which the given inspection
parameter $\boldsymbol{\widetilde{\theta}}$ results in a change in the
expectation of the transformed series. This problem will be discussed in
detail in Section 3.2.1.

In the mean change model based on moving sample means as considered by
Eichinger & Kirch, (2018) (see also Example 2.3.1 below) the above MOSUM-Wald
and MOSUM-score statistics coincide for any inspection parameter
$\boldsymbol{\widetilde{\theta}}$.

###  2.2 Segmentation Algorithm

MOSUM statistics as introduced in the previous section have peaks close to the
true change points making them particularly suitable for data segmentation.
More precisely, as demonstrated in Figure 1, the MOSUM statistic is a noisy
version of the MOSUM signal, which is a piecewise linear function that is
equal to zero away from the change points and has a single peak at each change
point.

Figure 1: Top panel: Time series with changes in the mean at time points
$100$, $250$ and $700$. Second panel: Signal of the classical MOSUM statistic
for mean changes as investigated by Eichinger & Kirch, (2018) with bandwidth
$G=60$. Third panel: Actual (noisy) MOSUM statistic for the above data set
(the vertical lines indicate the true change points). Fourth panel: Outcome of
the segmentation procedure as described in Section 2.2 in this example.

Clearly, in order to use this noisy signal for segmentation purposes a
threshold is needed to distinguish between significant local maxima that are
due to a close-by change point and local maxima obtained simply from random
fluctuations around zero if no change point is close by. We obtain such
thresholds by globally controlling the random fluctuations of the MOSUM
statistic asymptotically if no change point is present in the time series. As
a consequence the proposed procedure controls the (asymptotic) family-wise
error rate at level $\alpha_{n}$ (see below) for the detected change points.
In Section 2.4 below we detail how to choose such a threshold
$D_{n}(\alpha_{n},G)$ based on asymptotic $\alpha_{n}$-quantiles (of the no-
change situation). These quantiles can also be considered as critical values
for a corresponding uniform (across time) test procedure for the null
hypothesis of no change versus the alternative of one or more changes. An
asymptotic threshold is reasonable given the generality of our procedure, the
(non-parametric) error structure and as it also allows for model
misspecification.

Based on a MOSUM-Wald or MOSUM-score statistic $T_{k,n}(G)$ with corresponding
threshold $D_{n}(\alpha_{n},G)$ we propose the following MOSUM segmentation
procedure to determine estimators for the number and the locations of the
change points. Below, we suppress the possible dependence on other tuning
parameters such as the inspection parameter $\boldsymbol{\widetilde{\theta}}$
for the MOSUM-score statistic:

Consider all pairs of time points $(v_{j,n},w_{j,n})$,
$v_{1,n}<w_{1,n}<v_{2,n}<w_{2,n}<\ldots$, with

$\displaystyle\;T_{k,n}(G)\geq
D_{n}(\alpha_{n},G)\;\;\text{for}\;\;v_{j,n}\leq k\leq w_{j,n},$
$\displaystyle\;T_{k,n}(G)<D_{n}(\alpha_{n},G)\;\;\text{for}\;\;k=v_{j,n}-1,w_{j,n}+1,$
$\displaystyle\;w_{j,n}-v_{j,n}\geq\varepsilon
G\;\;\;\text{with}\;\;0<\varepsilon<1/2\;\;\text{arbitrary but fixed}.$ (5)

We take the number of these pairs as an estimator for the number of changes:

$\widehat{q}_{n}=\;\;\text{number of pairs}\;(v_{j,n},w_{j,n}).$

Furthermore, we determine the local maxima between $v_{j,n}$ and $w_{j,n}$,
$j=1,\ldots,\widehat{q}_{n}$, and use them as estimators for the locations of
the change points:

$\widehat{k}_{j,n}=\operatorname*{arg\,max}\limits_{v_{j,n}\leq k\leq
w_{j,n}}T_{k,n}(G)$

The pairs $(v_{j,n},w_{j,n})$, $j=1,\ldots,\widehat{q}_{n}$, give start and
end points of intervals on which the statistic exceeds the threshold. For this
reason we call $[v_{j,n},w_{j,n}]$, $j=1,\ldots,\widehat{q}_{n}$, intervals of
exceedings or exceeding intervals.

Condition (5) avoids false positives obtained due to random fluctuation when
the linear signal crosses the threshold. While the linear signal is strictly
monotone when crossing the threshold, the noisy version may not be – possibly
jumping beneath the threshold and then above again simply by chance.

By the below theory (see Section 3), we can be relatively certain (with
asymptotic probability one), that the exceeding intervals contain a change
point. Typically $T_{k,n}(G)$ depends on an estimator for the long-run
covariance matrix $\boldsymbol{\Gamma}_{k}$ resp. $\boldsymbol{\Sigma}_{k}$,
where in many situations a local estimator (depending on the time point $k$)
is required because the long-run covariance matrix depends on the segment.
However, once the exceeding intervals have been determined, different and
possibly better estimators $\widehat{\boldsymbol{\Psi}}_{j,n}$,
$j=1,\ldots,\widehat{q}_{n}$, can be used to obtain the final change point
estimator

$\displaystyle\widehat{k}_{j,n}(\widehat{\boldsymbol{\Psi}}_{j,n})=\operatorname*{arg\,max}\limits_{v_{j,n}\leq
k\leq w_{j,n}}T_{k,n}(G;\widehat{\boldsymbol{\Psi}}_{j,n}),$ (6)

where $\widehat{\boldsymbol{\Psi}}_{j,n}$ is a sequence of symmetric positive
definite matrices fulfilling Assumptions 3.4 below. For a real-valued
estimating function $\mathbf{H}$ this is always equivalent to using the
unscaled maximizer, i.e. $\widehat{\boldsymbol{\Psi}}_{j,n}=1$. For a vector-
valued estimating function it is also possible (e.g. for computational
reasons) to use the identity matrix
$\widehat{\boldsymbol{\Psi}}_{j,n}=\mathbf{\mbox{Id}}$ without jeopardizing
the asymptotic localization rates. Indeed, the theory does not even require
the $\widehat{\boldsymbol{\Psi}}_{j,n}$ to be consistent estimators for the
long-run covariance matrix. While we use the initial estimator in the
simulation study, we obtain the localization rates for the latter. A
corresponding analysis for covariance estimators depending on the time point
$k$ is not possible in our general framework, but can be done for specific
estimators (see e.g. Remark 3.2 and Corollary 3.1 in Eichinger & Kirch, (2018)
for an example).

###  2.3 Examples

The framework that is used in this work covers many different situations (see
e.g. Kirch & Kamgaing, (2015) and Kirch & Weber, (2018) for additional
examples). Here, we concentrate on the following three examples: Uni- and
multivariate changes in expectation where the sample mean but also more robust
estimators can be used, linear regression and Poisson autoregression.

####  2.3.1 Mean Change Model

The univariate mean change model is given by

$\displaystyle X_{i}=\sum_{j=1}^{q+1}\mu_{j}\,1_{\\{k_{j-1,n}<i\leq
k_{j,n}\\}}+\varepsilon_{i},$

where $\\{\varepsilon_{i}:i\geq 1\\}$ is a centered stationary error sequence
with autocovariance function $\gamma(\cdot)$ and long-run variance
$0<\tau^{2}=\sum_{h\in\mathbb{Z}}\gamma(h)<\infty$. Thus, the expected value
of each stretch is given by $\mu_{j}\neq\mu_{j+1}$, $j=1,\ldots,q$. The MOSUM
statistic discussed in Eichinger & Kirch, (2018) is based on the sample mean
$\bar{X}_{1,n}=\frac{1}{n}\sum_{i=1}^{n}X_{i}$ which is the solution of the
estimating equation $\sum_{i=1}^{n}(X_{i}-\mu)\overset{!}{=}0$, hence
$\mathbb{X}_{i}=X_{i}$.

An alternative M-estimator for the expectation for symmetric error
distributions is based on the estimating function
$H(X_{i},\mu)=\frac{2}{\pi}\arctan(\mu-X_{i})$. This estimating function leads
to more robust estimators and can be seen as a smooth approximation of the
sign-function that leads to the median. Henceforth, we will therefore call the
corresponding estimator _median-like estimator_.

A multivariate version $\boldsymbol{X}_{i}=(X_{1,i},\ldots,X_{d,i})^{T}$ is
obtained by replacing $\varepsilon_{j}$ as well as $\mu_{j}$ with multivariate
quantities. The multivariate estimating function is given by
$\boldsymbol{H}(\boldsymbol{X}_{i},\boldsymbol{\mu})=(H(X_{1,i},\mu_{1}),\ldots,H(X_{d,i},\mu_{d}))^{T}$
for $\boldsymbol{\mu}=(\mu_{1},\ldots,\mu_{d})^{T}$, i.e.
$\mathbb{X}_{i}=\boldsymbol{X}_{i}$. Effectively, the procedure uses the
corresponding univariate estimators in each component. Consequently, an
estimator for the long-run covariance matrix is now required which poses a
major difficulty in practice.

####  2.3.2 Linear Regression Model with Structural Breaks

The linear regression model with changes is given by

$\displaystyle
X_{i}=\sum_{j=1}^{q+1}\boldsymbol{Z}_{i}^{T}\boldsymbol{\beta}_{j}\,1_{\\{k_{j-1,n}<i\leq
k_{j,n}\\}}+\varepsilon_{i},$

where $X_{i}$ denotes the response variable, $\boldsymbol{Z}_{i}$ are the
regressors, $\boldsymbol{\beta}_{1},\ldots,\boldsymbol{\beta}_{q+1}$ are the
parameter vectors and $\\{\varepsilon_{i}\\}_{i\geq 1}$ represents an
innovation sequence as in Section 2.3.1. A least squares approach minimising
the sum of squared residuals corresponds to the vector-valued estimating
function

$\displaystyle\boldsymbol{H}\left((X_{i},\boldsymbol{Z}_{i})^{T},\boldsymbol{\beta}\right)=-2\boldsymbol{Z}_{i}\left(X_{i}-\boldsymbol{Z}_{i}^{T}\boldsymbol{\beta}\right),$

i.e. $\mathbb{X}_{i}=(X_{i},\boldsymbol{Z}^{T}_{i})^{T}$. For the mathematical
results some additional assumptions on the regressors such as stationarity and
independence from the innovation sequence are required.

####  2.3.3 Poisson Autoregressive Model with Structural Breaks

In this example we consider the Poisson autoregressive model of order one,
which is a widespread model for integer-valued time series.

A Poisson autoregressive model of order one, also known as INARCH(1) model,
with $q$ change points can be described by (1), where $\\{X_{i}^{(j)}\\}$ are
INARCH(1) time series with parameters
$\boldsymbol{\theta}_{j}=(\theta_{j,1},\theta_{j,2})^{T}$, $j=1,\ldots,q+1$,
i.e.

$\displaystyle X_{i}^{(j)}|\mathcal{F}_{i-1}\sim
Poi(\lambda_{i}),\quad\text{with
}\lambda_{i}=\theta_{j,1}+\theta_{j,2}X_{i-1}.$

In this model the observation $X_{i}$ conditioned on the past is Poisson
distributed with a parameter $\lambda_{i}$ depending on the past observation.
A (partial) likelihood approach leads to the estimating function

$\boldsymbol{H}((X_{i},X_{i-1})^{T},\boldsymbol{\theta})=-2\boldsymbol{X}_{i-1}\left(\frac{X_{i}}{\boldsymbol{X}_{i-1}^{T}\boldsymbol{\theta}}-1\right),$

where $\boldsymbol{X}_{i-1}=\left(1,X_{i-1}\right)^{T}$, i.e. here
$\mathbb{X}_{i}=(X_{i},X_{i-1})^{T}$.

####  2.3.4 Further Examples

In the context of at-most-one-change a posteriori as well as sequential change
point analysis, many procedures proposed in the literature are of this type
(see Kirch & Tadjuidje Kamgaing, (2016a) as well as Kirch & Kamgaing, (2015)
for more examples), such that their extensions to the data segmentation
problem is covered by this framework.

In the context of the mean change model as in Section 2.3.1, in many
applications the variance of the error sequence
$\mathrm{Var}(\varepsilon_{i})=\sigma^{2}$ changes in addition to the mean.
Thus, it is reasonable to assume that the variance is only piecewise constant
such that a change point can be caused by a mean change, a variance change or
a change in both parameters. In order to detect these types of changes one can
use the following vector-valued estimating function
$\boldsymbol{H}(X_{i},\mu,\sigma^{2})=\left(X_{i}-\mu,(X_{i}-\mu)^{2}-\sigma^{2}\right)^{T}$
(which relates to the method of moments but also the maximum likelihood
estimator under Gaussianity assumption). A related MOSUM-Wald statistic has
recently been proposed and analyzed by Messer, (2019).

###  2.4 Threshold Selection

In the following theorem we derive the limit distribution of the maximum of
the MOSUM statistics if no structural break occurs. To this end, let

$\displaystyle a(x)=\sqrt{2\log(x)}\text{ and }$ (7) $\displaystyle
b(x)=2\log(x)+\frac{p}{2}\log(\log(x))-\log\left(\frac{2}{3}\,\Gamma\left(\frac{p}{2}\right)\right),$

where $p$ is the dimension of the parameter space and $\Gamma$ denotes the
gamma function.

######  Theorem 2.1.

Consider model (1) with no change, i.e. $q=0$. Let Assumptions 2.3 on the
bandwidth hold, in addition to Assumptions 2.1, where (b) needs to hold with
$\boldsymbol{\theta}=\boldsymbol{\theta}_{1}$ (as defined in Assumption 2.2)
for the MOSUM-Wald and with
$\boldsymbol{\theta}=\boldsymbol{\widetilde{\theta}}$ for the MOSUM-score
statistic. For the MOSUM-Wald statistic, let also Assumption 2.2 hold.

  1. (a)

Then, for both the MOSUM-Wald ($\ell=1$) as well as the MOSUM-score ($\ell=2$
with $T_{k,n}^{(2)}(G)=T_{k,n}^{(2)}(G,\boldsymbol{\widetilde{\theta}})$)
statistic it holds

$\displaystyle
a(n/G)\,\max_{k=G,\ldots,n-G}T^{(\ell)}_{k,n}(G)-b(n/G)\;\xrightarrow{\mathcal{D}}\;E$

with $E$ a Gumbel distributed random variable fulfilling $P(E\leq
x)=\exp(-2\exp(-x))$ and with $a(x)$ and $b(x)$ as in (7).

  2. (b)

For the MOSUM-score statistic the result remains true for data-dependent
inspection parameters $\boldsymbol{\widetilde{\theta}}_{n}$ if Assumption 2.4
holds. Additionally, the assertions remain true if the long-run covariance
matrix $\boldsymbol{\Sigma}_{1}$ is replaced by an estimator fulfilling (with
$\|\cdot\|$ denoting the spectral norm of a matrix)

$\displaystyle\max_{G\leq k\leq
n-G}\left\lVert\widehat{\boldsymbol{\Sigma}}_{k,n}^{-1/2}-\boldsymbol{\Sigma}_{1}^{-1/2}\right\rVert=o_{P}\biggl{(}\bigl{(}\log(n/G)\bigr{)}^{-1}\biggr{)}.$

  3. (c)

The assertion for the MOSUM-Wald statistic remains true if both
$\boldsymbol{\Sigma}_{1}$ and $\boldsymbol{V}_{1}$ are replaced by estimators
$\widehat{\boldsymbol{\Sigma}}_{k,n}$ respectively
$\widehat{\boldsymbol{V}}_{k,n}$ fulfilling

$\displaystyle\max_{G\leq k\leq
n-G}\left\lVert\boldsymbol{\Sigma}_{1}^{-1/2}\boldsymbol{V}(\boldsymbol{\theta}_{1})-\widehat{\boldsymbol{\Sigma}}_{k,n}^{-1/2}\widehat{\boldsymbol{V}}_{k,n}\right\rVert=o_{P}\left(\left(\log(n/G)\right)^{-1}\right).$
(8)

As threshold in our segmentation algorithm as described in Section 2.2 we use
the asymptotic $(1-\alpha_{n})$-quantiles as given in Theorem 2.1 i.e.

$\displaystyle D_{n}(\alpha_{n},G)=\frac{b(n/G)+c_{\alpha_{n}}}{a(n/G)},\qquad
c_{\alpha_{n}}=-\log\log\frac{1}{\sqrt{1-\alpha_{n}}},$ (9)

where $c_{\alpha_{n}}$ is the $(1-\alpha_{n})$-quantile of the limiting Gumbel
distribution. Furthermore, we require $\alpha_{n}\to 0$ sufficiently slowly in
the sense of:

######  Assumption 2.5.

Let the sequence of significance levels $\alpha_{n}$ fulfill

$\displaystyle\alpha_{n}\rightarrow
0\;\;\;\text{and}\;\;\;c_{\alpha_{n}}=O(b(n/G)).$

######  Remark 2.1.

For fixed $\alpha$ the threshold $D_{n}(\alpha,G)$ can be used as a critical
value in a test procedure based on MOSUM statistics for the null hypothesis of
no change versus the alternative of at least one change. Proposition 3.1 below
shows when such a test procedure has asymptotic power one. While the power is
usually not as good as for the corresponding at-most-one-change statistics
even in the presence of multiple changes, the latter are not as suitable for
the localization of change points.

##  3 Consistency of the MOSUM Segmentation Procedures

###  3.1 MOSUM-Wald Procedure

Part (a) of the following assumption controls the behavior of estimators that
are not contaminated by change points and extends standards results from a
pointwise to a uniform (within a $G$-environment) situation. It follows from
Assumption 2.2 under weak assumptions on the time series and estimating
functions. On the other hand, for contaminated estimators (containing
observations from two segments), the estimator cannot be expected to converge
to either of the best approximating parameters but will be contaminated by
both (see Theorem 4.3 below). Part (b) of the below assumption guarantees that
the estimator differs from the best approximating parameter of the $j$th
segment if at least $\varepsilon G$ of the summands are from a neighboring
regime. For sufficiently smooth estimating functions and mixing time series,
the validity of the assumption will be shown in Section 4 below.

######  Assumption 3.1.

For any $j=1,\ldots,q$ let

$\displaystyle(a)\qquad\max_{\begin{subarray}{c}k_{j-1,n}<k\leq k_{j-1,n}+G\\\
k_{j,n}-2G\leq k\leq
k_{j,n}-G\end{subarray}}\sqrt{G}\,\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\|=O_{P}\left(\sqrt{\log(n/G)}\right).$
$\displaystyle(b)\qquad\sqrt{\frac{G}{\log(n/G)}}\min_{\begin{subarray}{c}k_{j-1,n}-G\leq
k\leq k_{j-1,n}-\varepsilon\,G\\\ k_{j,n}-(1-\varepsilon)\,G\leq k\leq
k_{j,n}\end{subarray}}\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\|\overset{P}{\longrightarrow}\infty.$

Because the long-run covariance matrix
$\boldsymbol{\Gamma}_{k}=\boldsymbol{V}_{k}^{-1}\boldsymbol{\Sigma}_{k}\left(\boldsymbol{V}_{k}^{-1}\right)^{T}$
is usually unknown in applications we need to understand the influence of
estimating this quantity on the procedure. This is reflected by the following
assumption stating that the estimator needs to be consistent away from changes
(related to (8)), while we can allow for a certain amount of misspecification
close to the changes.

######  Assumption 3.2.

Let $\widehat{\boldsymbol{\Gamma}}_{k,n}$ be an estimator sequence of
$\boldsymbol{\Gamma}_{k}$ satisfying

  1. (a)

$\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}\left\lVert\boldsymbol{\widehat{\Gamma}}^{-1/2}_{k,n}-\boldsymbol{\Gamma}_{k}^{-1/2}\right\rVert=o_{P}\left(\log(n/G)^{-1}\right)$
for any $j=1,\ldots,q+1$,

  2. (b)

for any $j=1,\ldots,q$ that

$\displaystyle\sup_{|k-k_{j,n}|\leq(1-\varepsilon)G}\left\lVert\boldsymbol{\widehat{\Gamma}}_{k,n}^{1/2}\right\rVert<\infty,\qquad\sup_{|k-k_{j,n}|\leq(1-\varepsilon)G}\left\lVert\boldsymbol{\widehat{\Gamma}}_{k,n}^{-1/2}\right\rVert<\infty.$

######  Remark 3.1.

In Assumption 3.2 (a) we require a uniformly consistent covariance estimator
at locations well away from change points. This assumption allows to work with
the minimal threshold as discussed in Section 2.4 of exact order
$\sqrt{\log(n/G)}$ which in turns leads to minimal assumptions on the signal
strength in Assumption 3.1(b). In practice, it can be difficult or numerically
expensive to use such estimators (see also Section 5). In such cases, we can
replace this assumption by the weaker one that

$\displaystyle\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}\left\lVert\boldsymbol{\widehat{\Gamma}}^{-1/2}_{k,n}\right\rVert<\infty.$

Then, the below results remain true as long as we use a slightly larger
threshold $\widetilde{D}_{n}(G)$ that fulfills
$\widetilde{D}_{n}(G)/\sqrt{\log(n/G)}\to\infty$ as well as slightly
strengthen Assumption 3.1(b) on the signal strength to

$\displaystyle\frac{\sqrt{G}}{\widetilde{D}_{n}(G)}\min_{\begin{subarray}{c}k_{j-1,n}-G\leq
k\leq k_{j-1,n}-\varepsilon\,G\\\ k_{j,n}-(1-\varepsilon)\,G\leq k\leq
k_{j,n}\end{subarray}}\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\|\overset{P}{\longrightarrow}\infty.$

We are now ready to prove a proposition, from which we will derive a
consistency theorem for the MOSUM segmentation including first rates of
convergence.

######  Proposition 3.1.

Let Assumptions 2.3 on the bandwidth hold, in addition to Assumptions 2.1,
where (b) needs to hold with $\boldsymbol{\theta}=\boldsymbol{\theta}_{j}$,
$j=1,\ldots,q+1$, (as defined in Assumption 2.2) in addition to Assumptions
2.2 and 3.1.

  1. (a)

Then, for $D_{n}(\alpha_{n},G)$ as in (9) fulfilling Assumption 2.5

$\displaystyle(i)\quad$ $\displaystyle
P\left(\max_{j=1,\ldots,q+1}\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}T^{(1)}_{k,n}(G)<D_{n}(\alpha_{n},G)\right)\rightarrow 1,$
$\displaystyle(ii)\quad$ $\displaystyle
P\left(\min_{j=1,\ldots,q}\,\min_{|k-k_{j,n}|\leq(1-\varepsilon)\,G}T^{(1)}_{k,n}(G)\geq
D_{n}(\alpha_{n},G)\right)\rightarrow 1.$

  2. (b)

Furthermore, the statements of (a) remain true if the long-run covariance
matrix $\boldsymbol{\Gamma}_{k}$ is replaced by estimators
$\widehat{\boldsymbol{\Gamma}}_{k,n}$ satisfying Assumption 3.2.

The following theorem shows, that the number of change points are estimated
consistently and that the distance from each change point to the closest
estimator is smaller than $G$.

######  Theorem 3.2.

Let the assumptions of Proposition 3.1 hold (either with known or estimated
long-run covariance matrix). Then, as $n\to\infty$,

$\displaystyle P\left(\widehat{q}_{n}^{(1)}=q,\max_{1\leq j\leq
q}\left|\widehat{k}^{(1)}_{j,n}-k_{j,n}\right|<G\right)\rightarrow 1,$

where $\widehat{q}_{n}^{(1)}$ is the estimated number of change points based
on the MOSUM-Wald statistics and $\widehat{k}^{(1)}_{j,n}$ are the
corresponding (ordered) change point estimators.

As a special case it follows:

######  Remark 3.2.

In the classical multiple change point situation, where
$k_{j,n}=\lfloor\lambda_{j}\,n\rfloor$ for
$0<\lambda_{1}<\ldots<\lambda_{q}<1$ and
$\widehat{\lambda}^{(1)}_{j,n}=\widehat{k}^{(1)}_{j,n}/n$, it holds

$\displaystyle\max_{j=1,\ldots,\min(q,\widehat{q}^{(1)}_{n})}\left|\widehat{\lambda}_{j,n}-\lambda_{j,n}\right|=O_{P}\left(\frac{G}{n}\right)=o_{P}(1).$

###  3.2 MOSUM-Score Procedure

The main disadvantage of the MOSUM-Wald procedure is its large computational
complexity in non-linear problems where $n-G$ estimators need to be calculated
via numerical optimization. Additionally, the use of numerical procedures can
cause statistical problems if distant parameters describe very similar models
(corresponding to many local maxima for the numerical optimization). A
prominent example of this type are neural networks.The MOSUM-score procedure
on the other hand does not suffer from this problem but may (depending on the
underlying model and choice of estimating function) not detect all changes if
only applied with a single inspection parameter.

####  3.2.1 Detectability

A change in the parameter vector $\boldsymbol{\theta}$ can only be detected or
localized by the MOSUM-score statistics if it causes a change in the
expectation of the transformed series
$\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\widetilde{\theta}})$ for the
chosen inspection parameter $\boldsymbol{\widetilde{\theta}}$. Denote the
corresponding set and its cardinality by

$\displaystyle\tilde{Q}=\tilde{Q}(\boldsymbol{\widetilde{\theta}})=\left\\{1\leq
j\leq
q\,\Big{|}\,E\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\widetilde{\theta}})\neq
E\boldsymbol{H}(\mathbb{X}_{1}^{(j+1)},\boldsymbol{\widetilde{\theta}})\right\\},\qquad$
$\displaystyle\tilde{q}=\tilde{q}(\boldsymbol{\widetilde{\theta}})=|\tilde{Q}|.$
(10)

In order to formulate the below results it is helpful to relabel detectable
changes $k_{j,n}$ with $j\in\tilde{Q}$ by
$\tilde{k}_{j,n}=\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})$,
$j=1,\ldots,\tilde{q}$, where
$\tilde{k}_{1,n}<\ldots<\tilde{k}_{\tilde{q},n}$.

Figure 2: The plots give two examples where $\tilde{q}\neq q$, i.e. not all
changes are detectable. The upper panel shows a simulated time series, while
the lower panel shows the corresponding signal of the MOSUM-score statistic
with the global median as inspection parameter and the estimating function of
the median. The change points are indicated by the red vertical lines.

Clearly, the number of detectable changes depends on both the choice of
estimating function and the inspection parameter. Figure 2 gives an example
for $\tilde{q}\neq q$, where the estimating function for the median, i.e.
$H(x,\mu)=\operatorname{sgn}(x-\mu)$, was used with the global median as
inspection parameter. Using a smooth strictly monotone approximation of this
estimating function (compare Section 2.3.1) makes all changes detectable
theoretically but still leads to power loss for some changes, which can be
remedied by using different inspection parameters in addition to some post-
processing (see Remark 3.3 below). This effect is illustrated by means of a
simulation study as well as a data example in Section 5.1 below.

For the classical multiple change point situation, where
$k_{j,n}=\lfloor\lambda_{j}\,n\rfloor$ for
$0=\lambda_{0}<\lambda_{1}<\ldots<\lambda_{q}<\lambda_{q+1}=1$ a canonical
data-driven inspection parameter is
$\boldsymbol{\widetilde{\theta}}_{n}=\widehat{\boldsymbol{\theta}}_{1,n}$
defined as the solution of
$\sum_{i=1}^{n}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\theta})\overset{!}{=}\boldsymbol{0}$.

Theorem 4.3 below shows that under mild conditions this data-dependent
inspection parameter fulfills the required assumptions with
$\boldsymbol{\widetilde{\theta}}_{0,1}$ as the unique solution of
$\sum_{j=1}^{q+1}(\lambda_{j}-\lambda_{j-1})E\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta})\overset{!}{=}\mathbf{0}$.

######  Lemma 3.3.

Let $\boldsymbol{\theta}_{j}$ denote the unique zero of
$E\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta})$, $j=1,\ldots,q+1,$
with $\boldsymbol{\theta}_{j}\neq\boldsymbol{\theta}_{j+1}$ for all
$j=1,\ldots,q$ and $q\geq 1$.

  * (a)

Then, the MOSUM-score procedure with detection parameter
$\boldsymbol{\widetilde{\theta}}_{0,1}$ (or a corresponding data driven
version) detects at least one change, i.e. $\widetilde{q}\geq 1$.

  * (b)

If there are only two distinct regimes then all changes are detectable, i.e.
$\widetilde{q}=q$.

######  Remark 3.3.

In particular, Lemma 3.3 shows that at least one change point is detectable if
the global estimator (computed on the whole sample) is used as inspection
parameter. Consequently, by recursively applying this procedure on detected
segments, all changes are detected with a much smaller computational burden
than for the MOSUM-Wald procedure. An empirical illustration based on
simulated data and a well-log data set is given in Section 5.1 below. This is
of particular interest if the procedure is used for candidate generation to be
combined with a pruning step as suggested by Cho & Kirch, (2021+). Some first
considerations and results in that direction can be found in Chapter 5 of
Reckrühm, (2019).

####  3.2.2 Consistency

The MOSUM-score statistics depends on the long-run covariance matrix
$\boldsymbol{\Sigma}_{k}(\boldsymbol{\widetilde{\theta}})$ as in (2) which is
typically not known. We show that it can be replaced by an estimator as long
as the estimator is consistent away from changes and fulfills some weaker
assumptions close to a change point. The latter is important because local
estimators will typically be contaminated by a close-by change point. For many
models (other than mean changes)
$\boldsymbol{\Sigma}_{k}(\boldsymbol{\widetilde{\theta}})$ differs from one
segment to the next, such that there is is no alternative to local estimation.

######  Assumption 3.3.

Let $\boldsymbol{\widehat{\Sigma}}_{k,n}$ satisfy

  1. (a)

$\max\limits_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}\left\lVert\boldsymbol{\widehat{\Sigma}}^{-1/2}_{k,n}-\boldsymbol{\Sigma}_{k}(\boldsymbol{\widetilde{\theta}})^{-1/2}\right\rVert=o_{P}\left(\log(n/G)^{-1}\right),$
for any $j=1,\ldots,q+1$,

  2. (b)

For any $j=1,\ldots,q$ it holds

$\displaystyle\sup_{|k-k_{j,n}|<G}\left\lVert\boldsymbol{\widehat{\Sigma}}_{k,n}^{1/2}\right\rVert<\infty,\qquad\sup_{|k-k_{j,n}|<G}\left\lVert\boldsymbol{\widehat{\Sigma}}_{k,n}^{-1/2}\right\rVert<\infty.$

######  Remark 3.4.

As in Remark 3.1 it can be difficult or numerically expensive to use an
estimator that is uniformly consistent away from change points. This is even
more important in combination with the MOSUM-score statistics that is mainly
designed to ease the computational burden of the procedure (see also Section
5). Here, too, the below results remain true as long as one uses a threshold
$\widetilde{D}_{n}(T)$ with $\widetilde{D}_{n}(T)/\sqrt{\log(n/G)}\to\infty$.

In analogy to the MOSUM-Wald statistics we prove the following proposition.

######  Proposition 3.4.

Let Assumptions 2.3 on the bandwidth hold, in addition to Assumptions 2.4,
thus allowing for a sequence of inspection parameters
$\boldsymbol{\widetilde{\theta}}_{n}$. Furthermore, let Assumptions 2.1 hold
with $\boldsymbol{\theta}=\boldsymbol{\widetilde{\theta}}$ (as in Assumption
2.4) in (b).

  1. (a)

Then,

$\displaystyle(i)\quad$ $\displaystyle
P\left(\max_{j=1,\ldots,\tilde{q}+1}\max_{\tilde{k}_{j-1,n}+G\leq
k\leq\tilde{k}_{j,n}-G}T^{(2)}_{k,n}(G)<D_{n}(\alpha_{n},G)\right)\rightarrow
1,$ $\displaystyle(ii)\quad$ $\displaystyle
P\left(\min_{j=1,\ldots,\tilde{q}}\,\min_{|k-\tilde{k}_{j,n}|\leq(1-\varepsilon)\,G}T^{(2)}_{k,n}(G)\geq
D_{n}(\alpha_{n},G)\right)\rightarrow 1.$

  2. (b)

Furthermore, the statements of (a) remain true if the long-run covariance
matrix $\boldsymbol{\Sigma}_{k}(\boldsymbol{\widetilde{\theta}})$ is replaced
by estimators $\widehat{\boldsymbol{\Sigma}}_{k,n}$ satisfying Assumption 3.3.

From this we conclude in the following theorem that the number of change
points are estimated consistently and that the distance from each change point
to the closest estimator is smaller than $G$.

######  Theorem 3.5.

Let the assumptions of Proposition 3.4 hold (either with known or estimated
long-run covariance matrix). Then, as $n\to\infty$,

$\displaystyle
P\left(\widehat{q}_{n}^{(2)}(\boldsymbol{\widetilde{\theta}}_{n})=\tilde{q}(\boldsymbol{\widetilde{\theta}}),\max_{1\leq
j\leq\tilde{q}}\left|\widehat{k}^{(2)}_{j,n}(\boldsymbol{\widetilde{\theta}}_{n})-\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})\right|<G\right)\rightarrow
1,$

where $\widehat{q}_{n}^{(2)}(\boldsymbol{\widetilde{\theta}}_{n})$ is the
estimated number of change points based on the MOSUM-score statistics with
inspection parameters $\boldsymbol{\widetilde{\theta}}_{n}$ and
$\widehat{k}^{(2)}_{j,n}(\boldsymbol{\widetilde{\theta}}_{n})$ are the
corresponding (ordered) change point estimators.

######  Remark 3.5.

In the classical multiple change point situation, where
$\tilde{k}_{j,n}=\lfloor\tilde{\lambda}_{j}\,n\rfloor$ for
$0<\tilde{\lambda}_{1}<\ldots<\tilde{\lambda}_{\tilde{q}}<1$ and
$\widehat{\lambda}^{(2)}_{j,n}=\widehat{k}^{(2)}_{j,n}/n$, it holds

$\displaystyle\max_{j=1,\ldots,\min(\tilde{q},\widehat{q}^{(2)}_{n})}|\widehat{\lambda}_{j,n}-\tilde{\lambda}_{j,n}|=O_{P}\left(\frac{G}{n}\right)=o_{P}(1).$

####  3.2.3 Localization Rates

Localization rates quantify the distance between true and estimated change
points. For a bounded number of change points the minimax optimal localization
rates are known to be constant, i.e. it holds
$\max_{j}|\widehat{k}_{j,n}-k_{j}|=O_{P}(1)$; see e.g. Lemma 2 in Wang et al.,
(2020). In this section we show that the MOSUM-score estimators match these
minimax optimal localization rates.

This is important because the MOSUM-score procedure is particularly well
suited as a candidate generating algorithm to be combined with a pruning step
(see Cho & Kirch, (2021+) for a corresponding discussion in the context of the
mean change model). This is due to its computational efficiency and the fact
that combined with several inspection parameters all changes can be detected.
For such pruning methodology to work it is essential that the candidate
generating algorithm produces estimators that are sufficiently close to the
true change points.

In this section, we work with
$\widehat{k}_{j,n}^{(2)}(\boldsymbol{\widetilde{\theta}}_{n};\widehat{\boldsymbol{\Psi}}_{j,n})$
as in (6), i.e. for $j\in\widetilde{Q}$

$\displaystyle\widehat{k}_{j,n}^{(2)}(\boldsymbol{\widetilde{\theta}}_{n};\widehat{\boldsymbol{\Psi}}_{j,n})=\operatorname*{arg\,max}_{v_{j,n}\leq
k\leq
w_{j,n}}\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k)^{T}\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k).$

We need the following assumptions on the matrices
$\boldsymbol{\widehat{\Psi}}_{j,n}$.

######  Assumption 3.4.

$\boldsymbol{\widehat{\Psi}}_{j,n}$ is a sequence (possibly depending on the
data) of symmetric positive definite matrices with

$\displaystyle\|\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\|=O_{P}(1),\qquad\|\boldsymbol{\widehat{\Psi}}_{j,n}^{1/2}\|=O_{P}(1).$

Effectively, the assumption guarantees that the use of
$\boldsymbol{\widehat{\Psi}}_{j,n}$ does not kill the signal nor does it
explode the noise.

Furthermore, we need to slightly strengthen the previous assumptions by
additionally requiring:

######  Assumption 3.5.

  1. (a)

With $\boldsymbol{\widetilde{\theta}}$ as in Assumption 2.4 we require that
for any $\xi_{n}\to\infty$ and any $j=1,\ldots,\tilde{q}_{n}+1$, it holds

$\displaystyle(i)\quad$ $\displaystyle\text{for
}m=\tilde{k}_{j-1,n}(\boldsymbol{\widetilde{\theta}})+G,\;\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})-G,\;\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})\text{
that}$ $\displaystyle\max\limits_{\xi_{n}\leq k\leq
G}\frac{1}{k}\left\lVert\sum\limits_{i=m-k+1}^{m}\left(\boldsymbol{H}\left(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n}\right)-\boldsymbol{H}\left(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}\right)\right)\right\rVert=o_{P}(1),$
$\displaystyle(ii)\quad$ $\displaystyle\text{for
}m=\tilde{k}_{j-1,n}(\boldsymbol{\widetilde{\theta}}),\;\tilde{k}_{j-1,n}(\boldsymbol{\widetilde{\theta}})+G,\;\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})-G\text{
that}$ $\displaystyle\max\limits_{\xi_{n}\leq k\leq
G}\frac{1}{k}\left\lVert\sum\limits_{i=m+1}^{m+k}\left(\boldsymbol{H}\left(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n}\right)-\boldsymbol{H}\left(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}\right)\right)\right\rVert=o_{P}(1).$

  2. (b)

The following backward law of large numbers holds for any $j$ and any sequence
$\xi_{n}\to\infty$

$\displaystyle\max\limits_{\xi_{n}\leq k\leq
G}\frac{1}{k}\left\lVert\sum\limits_{i=G-k+1}^{G}\left(\boldsymbol{H}\left(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}\right)-E\boldsymbol{H}\left(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}\right)\right)\right\rVert=o_{P}(1).$

######  Remark 3.6.

The forward version (as in (a)(ii)) is missing in (b) because it follows
directly from a strong law of large number which (after changing the
probability space) follows directly from the invariance principle as in
Assumption 2.1 (b). However, the backward version does not follow from that
assumption alone - unless the data is i.i.d. Nevertheless, most sequences will
allow for strong ergodic theorems even in backward direction. For example
mixing conditions as in Regularity conditions 4.1(b) are symmetric in the
sense that the backward sequence is also mixing with the same rates.
Consequently, both a forward and backward invariance principle holds and (b)
follows from that.

In the following theorem we derive the minimax optimal localization rates
improving upon Theorem 3.5.

######  Theorem 3.6.

Let the assumptions of Theorem 3.5 hold in addition to Assumptions 3.4 – 3.5.
Then, it holds

$\displaystyle\max_{1\leq
j\leq\min(\tilde{q},\widehat{q}_{n}^{(2)}(\boldsymbol{\widetilde{\theta}}_{n}))}\left|\widehat{k}_{j,n}^{(2)}(\boldsymbol{\widetilde{\theta}}_{n};\widehat{\boldsymbol{\Psi}}_{j,n})-\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})\right|=O_{P}(1).$

##  4 Regularity Conditions for Sufficiently Smooth Estimating Functions

The results of the previous sections were obtained under certain high-level
assumptions, that need to be verified separately for each underlying time
series structure as well as estimating function. On the other hand, many
estimating functions are either sufficiently smooth or can be approximated to
any degree of accuracy by sufficiently smooth functions. Therefore, we will
verify these high-level assumptions exemplarily for sufficiently smooth
estimating function for strongly mixing time series under moment conditions.
This includes the i.i.d. situation as a special case.

The invariance principles of Assumption 2.1 (b) follow from the following
regularity conditions:

######  Regularity Condition 4.1.

For $j=1,\ldots,q+1$

  1. (a)

let there exist a $\tilde{\nu}>0$ such that
$0<E\left\lVert\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta})\right\rVert^{2+\tilde{\nu}}<\infty$
for the $\boldsymbol{\theta}$ of interest,

  2. (b)

let $\\{\mathbb{X}_{t}^{(j)}\\}$ be stationary and strongly mixing sequences
of random vectors with a strong mixing coefficient $\alpha(\cdot)$ satisfying
$\alpha(n)=O(n^{-\beta})$ for some $\beta>1+2/\tilde{\nu}$ as $n\to\infty$,
where $\tilde{\nu}$ is as (a).

######  Theorem 4.1.

  1. (a)

Under Regularity Conditions 4.1 Assumption 2.1 hold for the same
$\boldsymbol{\theta}$ and some $\nu>0$ depending on $\beta$, $\tilde{\nu}$ and
the dimension.

  2. (b)

Additionally, Assumption 3.5 (b) holds if Regularity Conditions 4.1 are
fulfilled with $\boldsymbol{\theta}=\boldsymbol{\widetilde{\theta}}$.

In the following we will always assume that the bandwidth $G$ fulfils
Assumption 2.3 for the above suitable $\nu$.

In order to prove the necessary assumptions for the consistency results of the
above MOSUM statistics the following additional regularity conditions on the
estimating function are sufficient.

######  Regularity Condition 4.2.

  1. (a)

Let $\boldsymbol{\Theta}\subset\mathbb{R}^{p}$ be a compact parameter space
and $\boldsymbol{H}=(H_{1},\ldots,H_{p})^{T}$ twice differentiable such that
for some $\tilde{\nu}>0$

$\displaystyle(i)\quad
E\left\lVert\nabla\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta})\right\rVert^{2+\tilde{\nu}}<\infty\qquad\text{for
all }\boldsymbol{\theta}\in\boldsymbol{\Theta},$ $\displaystyle(ii)\quad
E\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\left\lVert\nabla^{2}H_{l}(\mathbb{X}^{(j)}_{1},\boldsymbol{\theta})\right\rVert^{2+\tilde{\nu}}<\infty,\quad\text{for
all }l=1,\ldots,p,$

for all $j=1,\ldots,q+1$.

  2. (b)

Let
$\boldsymbol{V}_{(j)}(\boldsymbol{\theta})=E\left(\nabla\boldsymbol{H}(\mathbb{X}^{(j)}_{1},\boldsymbol{\theta})\right)^{T}$
be a regular matrix for all $\boldsymbol{\theta}\in\boldsymbol{\Theta}$
fulfilling for all $j=1,\ldots,q+1$

$\displaystyle\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\left\lVert\boldsymbol{V}_{(j)}(\boldsymbol{\theta})^{-1}\right\rVert<\infty.$

  3. (c)

Let $E\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\left\lVert\nabla
H_{l}(\mathbb{X}^{(j)}_{1},\boldsymbol{\theta})\right\rVert<\infty$,
$j=1,\ldots,q+1$.

For specific examples the above assumptions such as the compactness assumption
can be weakened. The important special case of linear regression models has
been discussed in detail in Reckrühm, (2019), Chapter 3.2, under standard
assumptions for linear regression.

######  Theorem 4.2.

Let Regularity Condition 4.1 hold with
$\boldsymbol{\theta}=\boldsymbol{\theta}_{j}$ as defined in Assumption 2.2
identifiably unique, $j=1,\ldots,q+1$, in addition to Regularity Condition 4.2
(a) and (b).

  1. (a)

Then Assumptions 2.2 and Assumption 3.1 (a) hold.

  2. (b)

If additionally Regularity Condition 4.2 (c) hold, we get Assumption 3.1 (b).

The MOSUM-score procedure requires an inspection parameter
$\boldsymbol{\widetilde{\theta}}$ which can also be data-dependent
$\boldsymbol{\widetilde{\theta}}_{n}$. We will show the validity of the
assumptions for situations where $\boldsymbol{\widetilde{\theta}}_{n}$ is a
$\sqrt{n}$-consistent estimator for some $\boldsymbol{\widetilde{\theta}}$ in
Theorem 4.4.

Here, we discuss the choice $\widehat{\boldsymbol{\theta}}_{a,b}$ defined as
the unique zero of
$\sum_{i=a}^{b}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\theta})\overset{!}{=}\boldsymbol{0}$,
which has advantageous properties for detectability (compare Section 3.2.1) in
particular if combined with a subsequent pruning step.

######  Theorem 4.3.

Let the classical multiple change point situation with
$k_{j,n}=\lfloor\lambda_{j}\,n\rfloor$ for
$0=\lambda_{0}<\lambda_{1}<\ldots<\lambda_{q}<\lambda_{q+1}=1$ hold. Consider
$a=\lfloor\gamma_{a}n\rfloor$ and $b=\lfloor\gamma_{b}n\rfloor$ for some
$0\leq\gamma_{a}<\gamma_{b}\leq 1$. Define
$\boldsymbol{\widetilde{\theta}}_{\gamma_{a},\gamma_{b}}$ as the unique
solution of
$\sum_{j=1}^{q+1}(\min(\lambda_{j},\gamma_{b})-\max(\lambda_{j-1},\gamma_{a}))_{+}\,E\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta})\overset{!}{=}\mathbf{0}$,
where $x_{+}=\max(x,0)$. Furthermore, let Regularity Conditions 4.1 hold with
$\boldsymbol{\theta}=\boldsymbol{\widetilde{\theta}}_{\gamma_{a},\gamma_{b}}$
as well as 4.2 (b) and (c). Then,

$\displaystyle\sqrt{n}\left(\widehat{\boldsymbol{\theta}}_{a,b}-\boldsymbol{\widetilde{\theta}}_{\gamma_{a},\gamma_{b}}\right)=O_{P}(1).$

Finally, we show the validity of the necessary assumptions for the MOSUM-score
procedure under suitable regularity conditions on the estimating function.

######  Theorem 4.4.

Let $\boldsymbol{\widetilde{\theta}}_{n}$ and
$\boldsymbol{\widetilde{\theta}}$ with
$\sqrt{n}(\boldsymbol{\widetilde{\theta}}_{n}-\boldsymbol{\widetilde{\theta}})=O_{P}(1)$
in addition to Regularity Condition 4.1 with
$\boldsymbol{\theta}=\boldsymbol{\widetilde{\theta}}$ hold. Additionally, let
Regularity Condition 4.2 (a) and (c) hold. Then Assumption 2.4 as well as 3.5
(a) hold.

##  5 Empirical comparison

For the mean change situation based on the sample mean as in Section 2.3.1
extensive simulation studies can be found in Eichinger & Kirch, (2018); Cho &
Kirch, (2021+); Meier et al., (2021). In this section, we want to give a proof
of concept beyond the mean by including a variety of examples in the
simulation study, where we focus on aspects that differ from the mean
situation. In Section 5.1 we aim at sheding some additional light on the
detectability issue discussed in Section 3.2.1 using the example of a mean
change in combination with the median-like estimator. In Sections 5.2 and 5.3
we focus on differences between the MOSUM-Wald and MOSUM-score statistics in
the empirical performance as well as their computation time. In particular,
this includes the example of linear regression in Section 5.2, where the
estimator is known analytically, as well as the example of Poisson
autogression in Section 5.3, where this is not the case and all estimators
need to be obtained numerically.

All simulations are based on $1000$ repetitions unless otherwise stated,
$\varepsilon$ as in (5) is set to $0.2$ and the threshold from Section 2.4
with $\alpha_{n}=5\%$ is used.

###  5.1 Median-like Estimator and Detectability

$G$ | $[80,120]$ | $[180,220]$ | $[580,620]$ | $[880,920]$
---|---|---|---|---
Inspection Parameter: Sample median $\hat{\mu}_{1,1000}$
$20$ | 0.019 | 1.000 | 0.935 | 0.142
$50$ | 0.343 | 1.000 | 1.000 | 0.665
Inspection Parameter: Sample median $\hat{\mu}_{1,200}$
$20$ | 0.158 | 0.985 | 0.380 | 0.026
$50$ | 0.659 | 1.000 | 0.999 | 0.404
Table 1: Percentage of simulations with a change point estimator in the given
intervals (’detection rate’) based on the MOSUM-score statistic with the
median-like estimating function

In Section 3.2.1 the detectability issue of the MOSUM-score procedure was
discussed in detail. Here, we shed some additional light onto the problem by a
corresponding empirical study and the analysis of well-log data where we use
the median-like estimating function as discussed in Section 2.3.1 with the
sample median (based on different stretches of the original data) as
inspection parameter. If instead the median-like estimator is used as an
inspection parameter, the detection rates as given in Table 1 only differ in
the third digit. However, using the sample median as implemented in R instead
of the median-like estimator (implemented using the uniroot-function from the
R-package stats) saved around 1/3 of computation time in our other-wise
identical simulation.

The procedure is implemented using the mosum function from the R-package mosum
applied to the transformed data sequence
$H(X_{t},\widehat{\boldsymbol{\theta}}_{n})$. In particular, the variance is
estimated empirically with the mosum-window estimator that is the default in
the package (see Meier et al., (2021)). Thus, the signal is effectively given
by the difference in the expectation of the transformed data
$EH(X_{t},\boldsymbol{\widetilde{\theta}})$ with the noise corresponding to
the distribution of $H(X_{t},\boldsymbol{\widetilde{\theta}})$ where
$\boldsymbol{\widetilde{\theta}}$ is the limit of
$\widehat{\boldsymbol{\theta}}_{n}$.

To elaborate we consider a noise sequence of $1000$ independent standard
normally distributed random variables which are then centered around a mean
sequence of $1,2,5,3,4$ with change points at $100,200,600$ and $900$.
Consequently, the global median will be somewhere between $3$ and $4$. As such
the second (jumping from mean $2$ to mean $5$) and to a slightly lesser degree
the third change point (jumping from $5$ to $3$) should be well visible to our
statistic applied with the global median, while the last one (jumping from $3$
to $4$) should be less visible. This is also confirmed by the empirical
detection rates as reported in Table 1. The first change point should be
almost invisible as it jumps from $1$ to $2$ with only few observations
exceeding the global median. While the detection rates in Table 1 are indeed
lower than for the other change points, they are well above the $5\%$
threshold with detection rates of approximately $20\%$ (for $G=20$) and even
above $80\%$ (with $G=50$). The reasons are twofold: First of all, more
observations from the second regime can be expected to cross the threshold
than for the first one. Secondly, with the median-like estimation function,
the data transformation is strictly monotonic such that any mean difference
remains theoretically detectable with any inspection parameter. Furthermore,
the transformation that dampens the signal (if both means are on the same side
of the inspection parameter) will also dampen the noise (i.e. variance) of the
signal. As soon as the sample median of the first $200$ observations is used
as an inspection parameter, which will be between $1$ and $2$, the detection
rates for the first change point increases while the ones for the last change
point decreases illustrating Remark 3.3.

#### Analysis of well-log data

One of the advantages of using the median in combination with the median-like
estimating function is the fact that they are robust with respect to outliers.
Indeed, recently, the distinction between outliers and change points has
received some increased attention. For example, Knoblauch et al., (2018),
Fearnhead & Rigaill, (2019) and Li & Yu, (2021) present robust change point
procedures and illustrate their usefulness using well-log data (see Figure 3
(a)): This geophysical data is collected from a probe being lowered into a
bore-hole, where change points occur when the probe moves from one rock strata
to another. The data notably contains several outliers, which have been
removed before the change point analysis in earlier works such as (Ruanaidh &
Fitzgerald,, 1996, Section 5.7.2), where the data was first discussed, or
Fearnhead & Clifford, (2003), Fearnhead, (2006), Adams & MacKay, (2007),Wyse
et al., (2011), Ruggieri & Antonellis, (2016).

(a) Well-log data with detected change points

(b) MOSUM-Score statistic with global sample median as inspection parameter

(c) MOSUM-Score statistic with median between $1070$ and $2767$ as inspection
parameter

Figure 3: Well-log data with segmentation obtained from the median-like MOSUM-
Score statistic with different inspection parameters. The vertical dashed red
lines give the detected change points using the global median as inspection
parameter, while the vertical solid blue lines give the detected change points
using as an inspection parameter the median of the stretch starting at $1070$
and ending at $2767$.

Using the median-like estimating function in combination with the sample
median yields a robust procedure such that we can use this methodology on the
original data set including the outliers. We use a bandwidth of $G=20$ and
first apply the procedure using the global median as inspection parameter (see
Figure 3 (b)). This results in the detection of 5 change points leaving some
more undetected (as predicted by the discussion in Section 3.2.1). In a second
step, we repeat the procedure using now the sample median between the first
(at $1070$) and last (at $2768$) detected change point as inspection parameter
(see Figure 3 (c)). This results in the detection of 9 additional change
points, which together result in a very reasonable segmentation of the data
(see Figure 3 (a)). Out of the three change points detected by both inspection
parameters in one case the estimator was exactly at the same location (at
$2470$), while the other two were shifted by 2 time points (one in either
direction). This further emphasizes the need for suitable post-processing
methodology that can deal with this type of duplicate estimators obtained from
different inspection parameters or different bandwidths (see also Remark 3.3).
Using additionally inspection parameters obtained from the data up to the
first respectively starting at the last change point (detected with the global
median at $1070$ respectively $2768$) does not yield any additional change
point estimators but some more duplicated ones ($8$ that are exactly equal,
and $6$ slightly shifted).

Concerning the last three detected change points in Figures 3 (a) and (c), the
question may arise as to whether these are truly change points or in fact
outliers. Indeed, the stretch between $3942$ and $3965$ as been removed by
some of the previous works as outliers but not all of them. The distinction
between serial outliers and change points depends on the length of the
’outlier’ stretch. While some authors considered this stretch of $23$ very
small observations as two change points, others have considered them serial
outliers. The MOSUM procedure is robust with respect to classical outliers
that can occur in both directions in an isolated way. Longer (with respect to
the bandwidth $G$) stretches of serial outliers are in some sense
indistinguishable from segments of length of order $G$ and thus may be
detected by the procedure: While the breakdown point of the median is at
$50\%$, smaller stretches of outliers (in the same direction) of length
$\gamma G$, $0<\gamma<0.5$, will lead to the situation that the sample median
of the window of $G$ observations (including the $\gamma G$ serial outliers)
will estimate the $(0.5\pm 0.5\,\gamma/(1-\gamma))$-quantile instead which
typically will differ from the median of the neighboring stretch that does not
contain any serial outliers. This also explains the last change point in
Figures 3 (a) and (c) because there the ’outlier’ stretch (if one decides to
classify it as such) has a length of around $10$ and we use a bandwidth of
only $20$. On the other hand a suitable choice of bandwidth $G$ in combination
with the choice of $\varepsilon$ in (5) will make the procedure robust with
respect to sufficiently short outlier stretches. The latter effect can e.g. be
seen in Figure 3 (c), where some of the other outlier stretches lead to
significant points in the MOSUM score statistic that are sufficiently isolated
to not be considered change points by the procedure due to the choice of
$\varepsilon=0.2$ in (5).

With larger bandwidths the long stretches at the beginning and the end of the
sequence are also segmented to some extend, which is consistent with
observations made by some previous authors (even after outlier removal)
because these stretches seem not to be mean-stationary which will cause
location segmentation procedures to fit a step-function.

###  5.2 Linear Regression

In this section, we will explore the difference in performance between the
MOSUM-Wald and MOSUM-score methodology when applied to linear regression where
we use the usual least squares methodology that requires no numerical
approximation.

To this end, we consider a time series
$Y_{i}=\boldsymbol{X}_{i}^{T}\boldsymbol{\beta}+\varepsilon_{i}$,
$i=1,\ldots,n,$ of length $n=1000$ with exogenous regressors
$\boldsymbol{X}_{i}=(1,X_{i,1},X_{i,2})^{T}$ with $X_{i,1}\sim N(1,1)$ and
$X_{i,2}\sim N(2,1)$ and i.i.d. standard normal errors. We include three
change points at $200,500$ and $800$ and the regression coefficients
$\boldsymbol{\beta}_{1}=(1,2,2)^{T},\boldsymbol{\beta}_{2}=(1,1,2)^{T},\boldsymbol{\beta}_{3}=(2,1,2)^{T}$
and $\boldsymbol{\beta}_{4}=(2,1,1)^{T}$. We use the global least-squares
regression estimator as inspection parameter
$\widehat{\boldsymbol{\beta}}_{1,n}$ for the MOSUM-score statistic.

As discussed in Remarks 3.1 and 3.4 we can make minimal assumptions on the
covariance estimators theoretically, but the small-sample performance
crucially depends on this estimator.

| Estimated number $\widehat{q}$ | Detection rate
---|---|---
$G$ | $\leq 1$ | $2$ | $3$ | $4$ | $\geq 5$ | $200$ | $500$ | $800$
MOSUM-score with _S-global_ covariance estimator
$50$ | 0.484 | 0.489 | 0.027 | 0 | 0 | 0.494 | 0.027 | 0.993
$100$ | 0.003 | 0.468 | 0.518 | 0.011 | 0 | 0.969 | 0.515 | 0.999
MOSUM-score with _S-local_ covariance estimator
$50$ | 0.110 | 0.502 | 0.353 | 0.034 | 0.001 | 0.804 | 0.430 | 1.000
$100$ | 0 | 0.049 | 0.918 | 0.033 | 0 | 0.985 | 0.917 | 1.000
MOSUM-Wald with _W-local_ covariance estimator
$50$ | 0.018 | 0.445 | 0.501 | 0.035 | 0.001 | 0.963 | 0.539 | 1.000
$100$ | 0 | 0.030 | 0.945 | 0.025 | 0 | 0.998 | 0.938 | 1.000
Table 2: Number of estimated change points and detection rate for all three
change points (i.e. percentage of simulations with a change point estimator in
the interval $[k_{j,n}-20,k_{j,n}+20]$) for the various scenarios in the
linear regression example.

As covariance estimators we use
$\boldsymbol{\widehat{\Sigma}}^{(j)}=\widehat{v}_{n}^{(j)}\,\frac{1}{n}\sum_{i=1}^{n}\boldsymbol{X}_{i}\boldsymbol{X}_{i}^{T}$
with

$\displaystyle\widehat{v}_{n}^{(1)}=\frac{1}{n-1}\sum_{i=1}^{n}(Y_{i}-\boldsymbol{X}_{i}^{T}\widehat{\boldsymbol{\beta}}_{1,n})^{2},$
(_S-global_)
$\displaystyle\widehat{v}_{n}^{(2)}=\frac{1}{2G}\left(\sum_{i=k-G+1}^{k}(\hat{\epsilon}_{i}-\bar{\epsilon}_{k-G+1,k})^{2}+\sum_{i=k+1}^{k+G}(\hat{\epsilon}_{i}-\bar{\epsilon}_{k+1,k+G})^{2}\right),$
$\displaystyle\text{with
}\hat{\epsilon}_{i}:=Y_{i}-\boldsymbol{X}_{i}^{T}\widehat{\boldsymbol{\beta}}_{1,n}\text{
and }\bar{\epsilon}_{l,u}:=\frac{1}{u-l+1}\sum_{i=l}^{u}\hat{\epsilon}_{i},$
(_S-local_)
$\displaystyle\widehat{v}_{n}^{(3)}=\frac{1}{2G}\left(\sum_{i=k-G+1}^{k}(Y_{i}-\boldsymbol{X}_{i}^{T}\widehat{\boldsymbol{\beta}}_{k-G+1,k})^{2}+\sum_{i=k+1}^{k+G}(Y_{i}-\boldsymbol{X}_{i}^{T}\widehat{\boldsymbol{\beta}}_{k+1,k+G})^{2}\right).$
(_W-local_)

We use the first two estimators with the MOSUM-score statistics while the last
one is only used with the Wald statistics where the estimators
$\widehat{\boldsymbol{\beta}}_{t+1,t+G}$ are already available. We do not use
this estimator with the MOSUM-score statistic because their usage cancels the
computational advantage of the MOSUM-score statistic over the MOSUM-Wald
statistic.

The first two estimators $\boldsymbol{\widehat{\Sigma}}^{(j)}$, $j=1,2$, are
only consistent in the no-change situation but do not fulfill assumption
3.3(a) in the presence of change points. Instead of using a threshold as in
Remark 3.4, in this simulation study, we stick to the threshold as in (9) for
all methods.

Table 2 gives the estimated number of change points as well as the detection
rates for all three change points in the various settings. For a bandwidth of
$50$ the MOSUM-Wald procedure (with _W-local_) outperforms the MOSUM-score
procedure (with _S-local_), while both procedures achieve a similar
performance for bandwidth $100$. The MOSUM-score procedure with the global
estimator $(\emph{S-global})$ is clearly inferior to both competing methods
(in particular for smaller bandwidth and the second change point) emphasizing
again the importance of the choice of the covariance estimator for the small
sample performance.

In terms of computation time (for $G=n^{2/3}$) the MOSUM-score clearly
outperforms the MOSUM-Wald statistics: Even with the local estimator (_S-
local_), the MOSUM-score statistic was roughly $22$ times faster in our
simulations than the MOSUM-Wald statistic (with _W-local_) for a time series
of length $1000$ (an average (out of $100$ runs) of $0.03$ seconds as compared
to $0.66$ seconds). For a length of $8000$ it was already more than $31$ times
faster ($0.22$ versus $6.91$ seconds). The numbers only give a qualitative
idea as we did not optimize any of the procedures with respect to computation
time but merely used the R-function rollsum (from the R-package zoo Zeileis &
Grothendieck, (2005)) to calculate the MOSUM-statistics, where the local
estimators for the covariance matrices are implemented naively with a loop and
the local regression parameter for the Wald-statistics is calculated with the
lm-function.

More simulation results including the false alarm rate in the no-change
situation, the results for other bandwidths and covariance estimators and more
information on computing times can be found in Reckrühm, (2019), Section 4.1.

###  5.3 Poisson Autoregressive Model

In this section we consider the Poisson autoregressive model with the
estimation function corresponding to the partial likelihood as in Section
2.3.3. Compared to the previous section this includes two additional
difficulties: First, due to the serial dependence of the data the true scaling
of the procedures depends on the long-run covariance rather than the
covariance matrix. This time-dependency is also the reason behind the larger
bandwidths compared to the previous section. Secondly, there is no analytical
solution to the estimating procedure such that numerical methods are required
which as expected will greatly increase computation time for the MOSUM-Wald
procedure.

| Estimated number $\widehat{q}$ | Detection rate
---|---|---
$G$ | $\leq 1$ | $2$ | $3$ | $4$ | $\geq 5$ | $250$ | $500$ | $750$
MOSUM-score procedure with $\boldsymbol{\widehat{\theta}}_{1,1000}$
$80$ | 0.619 | 0.288 | 0.063 | 0.028 | 0.002 | 0.713 | 0.135 | 0.242
$150$ | 0.056 | 0.321 | 0.449 | 0.137 | 0.037 | 0.921 | 0.583 | 0.623
MOSUM-score procedure with $\boldsymbol{\widehat{\theta}}_{300,700}$
$80$ | 0.100 | 0.397 | 0.300 | 0.143 | 0.060 | 0.936 | 0.199 | 0.734
$150$ | 0.018 | 0.162 | 0.596 | 0.194 | 0.030 | 0.919 | 0.724 | 0.742
MOSUM-Wald procedure
$80$ | 0.069 | 0.295 | 0.373 | 0.199 | 0.064 | 0.890 | 0.603 | 0.645
$150$ | 0.001 | 0.040 | 0.629 | 0.261 | 0.069 | 0.896 | 0.809 | 0.803
Table 3: Number of estimated change points and detection rate for all three
change points (i.e. percentage of simulations with a change point estimator in
the interval $[k_{j,n}-20,k_{j,n}+20]$) for the various scenarios in the
Poisson autoregressive example.

We consider a time series of length $n=1000$ with three change points at times
$250$, $500$ as well as $750$ with the paramaters
$\boldsymbol{\theta}_{1}=(1,0.5)^{T}$,
$\boldsymbol{\theta}_{2}=(2.5,0.5)^{T}$,
$\boldsymbol{\theta}_{3}=(2.5,0.2)^{T}$ as well as
$\boldsymbol{\theta}_{4}=(1,0.5)^{T}$.

For the MOSUM-score statistics we use the global (partial) maximum likelihood
estimator as inspection parameter as well as the one based on the observations
between time point $300$ and $700$ (compare also Remark 3.3). To estimate the
covariance matrix we use the following local estimator

$\displaystyle\widehat{\boldsymbol{\Sigma}}_{k,n}$
$\displaystyle=\frac{1}{2G}\sum_{i=k-G+1}^{k}\left(\boldsymbol{H}(\mathbb{Y}_{i},\boldsymbol{\widehat{\theta}}_{1,n})-\overline{\boldsymbol{H}}_{k-G+1,k}\right)\left(\boldsymbol{H}(\mathbb{Y}_{i},\boldsymbol{\widehat{\theta}}_{1,n})-\overline{\boldsymbol{H}}_{k-G+1,k}\right)^{T}$
$\displaystyle\quad+\frac{1}{2G}\sum_{i=k+1}^{k+G}\left(\boldsymbol{H}(\mathbb{Y}_{i},\boldsymbol{\widehat{\theta}}_{1,n})-\overline{\boldsymbol{H}}_{k+1,k+G}\right)\left(\boldsymbol{H}(\mathbb{Y}_{i},\boldsymbol{\widehat{\theta}}_{1,n})-\overline{\boldsymbol{H}}_{k+1,k+G}\right)^{T},$

where $\overline{\boldsymbol{H}}_{l,u}$ denotes the sample mean of
$\boldsymbol{H}(\mathbb{Y}_{l},\boldsymbol{\widehat{\theta}}_{1,n}),\ldots,\boldsymbol{H}(\mathbb{Y}_{u},\boldsymbol{\widehat{\theta}}_{1,n})$.
However, this estimator does not take the dependence into account and as such
estimates the covariance rather than the long-run covariance matrix, such that
Assumption 3.3(a) is not fulfilled. As before we use a threshold as in (9)
despite Remark 3.4.

Motivated by Weiß, (2010), equations (7) and (8), we use the estimator
$\widetilde{\boldsymbol{\Gamma}}_{k,n}^{-1}=\frac{1}{2}\left(\widetilde{\boldsymbol{\Gamma}}_{k-G+1,k}^{-1}+\widetilde{\boldsymbol{\Gamma}}_{k+1,k+G}^{-1}\right)$
in the MOSUM-Wald procedure, where

$\displaystyle\widetilde{\boldsymbol{\Gamma}}_{l,u}^{-1}=\frac{1}{G}\sum_{i=l}^{u}\frac{1}{(\boldsymbol{Y}_{i-1}^{T}\boldsymbol{\widehat{\theta}}_{l,u}^{ML})^{2}}\left(\begin{array}[]{cc}Y_{i}&Y_{i}Y_{i-1}\\\
Y_{i}Y_{i-1}&Y_{i}Y_{i-1}^{2}\end{array}\right).$

Table 3 gives the estimated number of change points as well as the detection
rates for all three change points in the various setting. In this case, the
MOSUM-Wald statistics outperforms the score procedures in terms of detection
rates for the second and third change point but indeed the MOSUM-score with
the restricted estimator (rather than the global one) slightly outperforms the
MOSUM-Wald for the first one.

The MOSUM-score is computationally much cheaper and scales better with longer
series. Indeed, the median computation time (for $100$ runs and $G=n^{2/3}$)
for the MOSUM-Wald statistics for a length of $n=1000$ (running several
minutes) was more than 272 times slower than the MOSUM-score statistic
(running less than a second). For $n=8000$ the MOSUM-Wald statistics ran for
more than half an hour which was more than 362 times slower than the MOSUM-
score (which ran for less than 6 seconds). In particular, the MOSUM-Wald
statistic is significantly slowed down by the need of numerical optimization
in comparison to the linear regression where no numerical methods are
required. The numbers only give a qualitative idea where a naive loop was used
for the calculation of the statistics and the global and local parameter
estimators for the Poisson regression were calculated with the R-package
tscount Liboschik et al., (2017).

In conclusion, the MOSUM-score statistic can be a good way to generate change
point candidates by means of using different inspection parameters (and
bandwidths) even in combination with a sub-optimal covariance estimation
procedure.

More simulation results including the false alarm rate in the no-change
situation, the results for other bandwidths and covariance estimators, results
for the least-squares estimators and more information on computing times can
be found in Reckrühm, (2019), Section 4.2.

##  6 Conclusion

In this work we propose and analyse a data segmentation procedure in a unified
framework including likelihood-based as well as robust methodology for complex
non-linear time series such as Poisson regression or neural-network based
(auto-)regression. To account for the diversity of models included in the
framework we work with high-level assumptions in Sections 2 and 3, that need
to be verified on a case-by-case basis in general. For sufficiently smooth
estimating functions this high-level assumptions are shown to be valid under
standard moment conditions in Section 4. Assumptions 3.2 for the MOSUM-Wald as
well as 3.3 for the MOSUM-score procedures and their relaxations in Remarks
3.1 and 3.4 on the covariance estimators are very general. Consequently, the
choice of covariance estimators has only a minimal influence on the asymptotic
properties of the corresponding change point estimators, however, it is
critical for the small sample behaviour: On the one hand (away from change
points) the choice of the threshold in combination with the covariance
estimator greatly influences the rate of spuriously detected change points. On
the other hand (close to change points) this combination determines the
detectability rate for a given signal strength. Estimators for the covariance
depend crucially on the given model and estimating function, such that they
need to be chosen separately for every situation.

Another contribution of this paper is the analysis of two different
procedures: The MOSUM-Wald methodology yields better results with a single
bandwidth but at the cost of a higher computation time in particular in
combination with non-linear models requiring numerical optimization
methodology. Furthermore, the MOSUM-Wald procedure is prone to errors in a
situation where many local maxima are present in the optimization procedure
such as is well known for neural network approximations. The reason is that
two parameters that are far away in the parameter space can describe almost
the same model but will be mistaken for different models by the MOSUM-Wald
procedure. The MOSUM-score procedure on the other hand is robust to a multiple
mode situations and in addition is computationally much lighter even for non-
linear models. On the other hand, with a single inspection parameter it might
not be able to detect all change points present (see Section 3.2.1 as well as
5.1). Future work will aim at solving this problem e.g. by combining the
methodology with a pruning step as has been discussed for the mean change
situation based on the sample mean in Cho & Kirch, (2021+). Such a step is
also necessary to deal with truly multiscale heterogeneous situations
containing both frequent large jumps as well as small jumps over long
stretches of stationarity (for a mathematical definition we refer to Cho &
Kirch, (2021+), Definition 2.1). To deal with such multiscale situations
several MOSUM procedures with different bandwidths can be used to generate
change point candidates further emphasizing the need for fast candidate
generating methods. For the MOSUM-score procedure only the computation of a
variety of data-driven inspection parameters requires numerical optimization
while the same inspection parameters can be used with a variety of bandwidths.
The MOSUM-Wald procedure on the other hand requires the (numerical)
computation of such estimators for each bandwidth and each point in time
making the MOSUM-Wald procedure computationally potentially highly
problematic. For a mathematical analysis of the the pruning step in the mean
change situation based on the sample mean, see Cho & Kirch, (2021+), some
first considerations in that direction in the framework of this paper can be
found in Chapter 5 of Reckrühm, (2019).

## Acknowledgements

This work was supported by the Deutsche Forschungsgemeinschaft (DFG, German
Research Foundation) - 314838170, GRK 2297 MathCoRe. Claudia Kirch would also
like to thank the Isaac Newton Institute for Mathematical Sciences for support
and hospitality during the programme ’Statistical scalability’ when work on
this paper was undertaken.

## References

  * Adams & MacKay, (2007) Adams, Ryan Prescott, & MacKay, David JC. 2007.  Bayesian online changepoint detection.  arXiv preprint arXiv:0710.3742.
  * Aggarwal et al., (1999) Aggarwal, Reena, Inclan, Carla, & Leal, Ricardo. 1999.  Volatility in emerging stock markets.  Journal of Financial and Quantitative Analysis, 34(1), 33–55.
  * Aue & Horváth, (2013) Aue, A., & Horváth, L. 2013.  Structural breaks in time series.  Journal of Time Series Analysis, 34(1), 1–16.
  * Bauer & Hackl, (1980) Bauer, Peter, & Hackl, Peter. 1980.  An extension of the MOSUM technique for quality control.  Technometrics, 22(1), 1–7.
  * Braun et al., (2000) Braun, J. V., Braun, R. K., & Müller, H.-G. 2000.  Multiple changepoint fitting via quasilikelihood, with application to DNA sequence segmentation.  Biometrika, 87(2), 301–314.
  * Cho & Kirch, (2021+) Cho, Haeran, & Kirch, Claudia. 2021+.  Two-stage data segmentation permitting multiscale changepoints, heavy tails and dependence.  Ann. Inst. Stat. Math. To appear.
  * Cho & Kirch, (2022+) Cho, Haeran, & Kirch, Claudia. 2022+.  Data segmentation algorithms: Univariate mean change and beyond.  Econometrics and Statistics.  To appear.
  * Chu et al., (1995) Chu, Chia-Shang J., Hornik, Kurt, & Kaun, Chung-Ming. 1995.  MOSUM tests for parameter constancy.  Biometrika, 82(3), 603–617.
  * Csörgö & Horváth, (1997) Csörgö, Miklós, & Horváth, Lajos. 1997.  Limit theorems in change-point analysis.  John Wiley & Sons Inc.
  * Eichinger & Kirch, (2018) Eichinger, Birte, & Kirch, Claudia. 2018.  A MOSUM procedure for the estimation of multiple random change points.  Bernoulli, 24(1), 526–564.
  * Fearnhead, (2006) Fearnhead, Paul. 2006.  Exact and efficient Bayesian inference for multiple changepoint problems.  Statistics and computing, 16(2), 203–213.
  * Fearnhead & Clifford, (2003) Fearnhead, Paul, & Clifford, Peter. 2003.  On-line inference for hidden Markov models via particle filters.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 65(4), 887–899.
  * Fearnhead & Rigaill, (2019) Fearnhead, Paul, & Rigaill, Guillem. 2019.  Changepoint detection in the presence of outliers.  Journal of the American Statistical Association, 114(525), 169–183.
  * Fearnhead & Rigaill, (2020) Fearnhead, Paul, & Rigaill, Guillem. 2020.  Relating and Comparing Methods for Detecting Changes in Mean.  Stat, e291.
  * He et al., (1996) He, Xuming, Shao, Qi-Man, et al. 1996\.  A general Bahadur representation of M-estimators and its application to linear regression with nonstochastic designs.  The Annals of Statistics, 24(6), 2608–2630.
  * Horváth & Rice, (2014) Horváth, L., & Rice, G. 2014.  Extensions of some classical methods in change point analysis.  Test, 23(2), 219–255.
  * Hušková, (1990) Hušková, Marie. 1990.  Asymptotics for robust MOSUM.  Commentationes Mathematicae Universitatis Carolinae, 31(2), 345–356.
  * Hušková & Slabỳ, (2001) Hušková, Marie, & Slabỳ, Aleš. 2001.  Permutation tests for multiple changes.  Kybernetika, 37(5), 605–622.
  * Killick et al., (2010) Killick, Rebecca, Eckley, Idris A., Ewans, Kevin, & Jonathan, Philip. 2010.  Detection of changes in variance of oceanographic time-series using changepoint analysis.  Ocean Engineering, 37, 1120–1126.
  * Killick et al., (2012) Killick, Rebecca, Fearnhead, Paul, & Eckley, Idris A. 2012.  Optimal detection of changepoints with a linear computational cost.  Journal of the American Statistical Association, 107(500), 1590–1598.
  * Kirch & Kamgaing, (2012) Kirch, Claudia, & Kamgaing, Joseph Tadjuidje. 2012.  Testing for parameter stability in nonlinear autoregressive models.  Journal of Time Series Analysis, 33(3), 365–385.
  * Kirch & Kamgaing, (2015) Kirch, Claudia, & Kamgaing, Joseph Tadjuidje. 2015.  On the use of estimating functions in monitoring time series for change points.  Journal of Statistical Planning and Inference, 161, 25–49.
  * Kirch & Klein, (2021+) Kirch, Claudia, & Klein, Philipp. 2021+.  Moving sum data segmentation for stochastics processes based on invariance.  Statistica Sinica.  To appear.
  * Kirch & Tadjuidje Kamgaing, (2016a) Kirch, Claudia, & Tadjuidje Kamgaing, Josef. 2016a.  Detection of change points in discrete valued time series.  Pages 219–244 of: Davis, Richard A., Holan, Scott H., Lund, Robert, & Ravishanker, Nalini (eds), Handbook of discrete valued time series.  CRC Press.
  * Kirch & Tadjuidje Kamgaing, (2016b) Kirch, Claudia, & Tadjuidje Kamgaing, Joseph. 2016b.  Detection of Change Points in Discrete-Valued Time Series.  Pages 239–264 of: Handbook of Discrete-Valued Time Series.  Chapman and Hall/CRC.
  * Kirch & Weber, (2018) Kirch, Claudia, & Weber, Silke. 2018.  Modified sequential change point procedures based on estimating functions.  Electronic Journal of Statistics, 12(1), 1579–1613.
  * Knoblauch et al., (2018) Knoblauch, Jeremias, Jewson, Jack E, & Damoulas, Theodoros. 2018.  Doubly Robust Bayesian Inference for Non-Stationary Streaming Data with $\beta$-Divergences.  Advances in Neural Information Processing Systems, 31.  Software, simulation code and data available at https://github.com/alan-turing-institute/rbocpdms/.
  * Kuelbs & Philipp, (1980) Kuelbs, James, & Philipp, Walter. 1980.  Almost sure invariance principles for partial sums of mixing B-valued random variables.  The Annals of Probability, 8(6), 1003–1036.
  * Li & Yu, (2021) Li, Mengchu, & Yu, Yi. 2021.  Adversarially robust change point detection.  Advances in Neural Information Processing Systems, 34, 22955–22967.
  * Liboschik et al., (2017) Liboschik, Tobias, Fokianos, Konstantinos, & Fried, Roland. 2017.  tscount: An R Package for Analysis of Count Time Series Following Generalized Linear Models.  Journal of Statistical Software, 82(5), 1–51.
  * Meier et al., (2021) Meier, A., Cho, H., & Kirch, C. 2021.  mosum: A package for moving sums in change point analysis.  Journal of Statistical Software, 97, 1–42.
  * Messer, (2019) Messer, Michael. 2019.  Bivariate change point detection: joint detection of changes in expectation and variance.  arXiv preprint arXiv:1904.01320.
  * Messer et al., (2014) Messer, Michael, Kirchner, Marietta, Schiemann, Julia, Roeper, Jochen, Neininger, Ralph, & Schneider, Gaby. 2014.  A multiple filter test for the detection of rate changes in renewal processes with varying variance.  The Annals of Applied Statistics, 8(4), 2027–2067.
  * Rao, (1962) Rao, R. Ranga. 1962.  Relations between weak and uniform convergence of measures with applications.  The Annals of Mathematical Statistics, 33(2), 659–680.
  * Reckrühm, (2019) Reckrühm, Kerstin. 2019.  Estimating Multiple Structural Breaks in Time Series: A Generalized MOSUM Approach Based on Estimating Functions.  Ph.D. thesis, Otto-von-Guericke-Universität Magdeburg.
  * Ruanaidh & Fitzgerald, (1996) Ruanaidh, Joseph JK O, & Fitzgerald, William J. 1996.  Numerical Bayesian methods applied to signal processing.  Springer Science & Business Media.
  * Ruggieri & Antonellis, (2016) Ruggieri, Eric, & Antonellis, Marcus. 2016.  An exact approach to Bayesian sequential change point detection.  Computational Statistics & Data Analysis, 97, 71–86.
  * Steinebach & Eastwood, (1996) Steinebach, Josef, & Eastwood, Vera R. 1996.  Extreme value asymptotics for multivariate renewal processes.  Journal of multivariate analysis, 56(2), 284–302.
  * Wang et al., (2020) Wang, Daren, Yu, Yi, & Rinaldo, Alessandro. 2020.  Univariate mean change point detection: Penalization, cusum and optimality.  Electronic Journal of Statistics, 14(1), 1917–1961.
  * Weiß, (2010) Weiß, Christian H. 2010.  The INARCH (1) model for overdispersed time series of counts.  Communications in Statistics-Simulation and Computation, 39(6), 1269–1291.
  * Wyse et al., (2011) Wyse, Jason, Friel, Nial, & Rue, Håvard. 2011.  Approximate simulation-free Bayesian inference for multiple changepoint models with dependence within segments.  Bayesian Analysis, 6(4), 501–528.
  * Yau & Zhao, (2016) Yau, Chun Yip, & Zhao, Zifeng. 2016.  Inference for multiple change points in time series via likelihood ratio scan statistics.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 78(4), 895–916.
  * Zeileis & Grothendieck, (2005) Zeileis, Achim, & Grothendieck, Gabor. 2005.  zoo: S3 Infrastructure for Regular and Irregular Time Series.  Journal of Statistical Software, 14(6), 1–27.

## Appendix

##  Appendix A Proofs of Section 2.4

###### Proof of Theorem 2.1.

We first prove the assertions for the MOSUM-score statistics ($\ell=2$): By
the invariance principle in Assumption 2.1 (b) as well as Assumption 2.3 on
the bandwidth we get (with
$\boldsymbol{\Sigma}=\boldsymbol{\Sigma}_{1}(\widetilde{\theta})$) that

$\displaystyle\max_{G\leq k\leq
n-G}\frac{1}{\sqrt{2G}}\left\lVert\boldsymbol{\Sigma}^{-1/2}\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)\right\rVert$
$\displaystyle=\max_{G\leq k\leq
n-G}\frac{1}{\sqrt{2G}}\left\lVert\boldsymbol{W}(k+G)-2\boldsymbol{W}(k)+\boldsymbol{W}(k-G)\right\rVert+o_{P}\left(a(n/G)^{-1}\right)$
$\displaystyle=\sup_{r\in[G,n-G]}\frac{1}{\sqrt{2G}}\left\lVert\boldsymbol{W}(r+G)-2\boldsymbol{W}(r)+\boldsymbol{W}(r-G)\right\rVert+o_{P}\left(a(n/G)^{-1}\right)$
$\displaystyle\overset{D}{=}\sup_{t\in[1,n/G-1]}\frac{1}{\sqrt{2}}\left\lVert\boldsymbol{W}(t+1)-2\boldsymbol{W}(t)+\boldsymbol{W}(t-1)\right\rVert+o_{P}\left(a(n/G)^{-1}\right),$

where we used the self-similarity of the Wiener process in the last step. By
an application of Lemma 3.1 in combination with Remark 3.1 of Steinebach &
Eastwood, (1996) with $\alpha=1$ and $C_{1}=\ldots=C_{p}=3/2$, assertion (a)
for the MOSUM-score statistics follows. The assertion for the MOSUM-Wald
statistics ($\ell=1$) follows from this and Assumption 2.2. The assertions in
(b) follow successively by an application of the triangular inequality in
combination with the consistency of the spectral matrix norm with the
Euclidean vector norm (as the corresponding induced norm). The assertions in
(c) are obtained analogously. ∎

##  Appendix B Proofs of Section 3.1

###### Proof of Proposition 3.1.

First, by Theorem 2.1,

$\displaystyle P\left(\max_{j=1,\ldots,q+1}\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}T^{(1)}_{k,n}(G)\geq D_{n}(\alpha_{n},G)\right)$
$\displaystyle\leq\sum_{j=1}^{q+1}P\left(a(n/G)\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}T^{(1)}_{k,n}(G)-b(n/G)\geq c_{\alpha_{n}}\right)$
$\displaystyle\leq\sum_{j=1}^{q+1}(\alpha_{n}+o(1))\to 0,$

showing (a) (i). Furthermore, because the spectral matrix norm is induced by
the Euclidean vector norm, it holds
$\|\mathbf{A}\mathbf{x}\|\leq\|\mathbf{A}\|\|\mathbf{x}\|$ as well as
$\|\mathbf{x}\|\leq\|\mathbf{A}^{-1}\|\|\mathbf{A}\mathbf{x}\|$. Then, by
Assumption 3.1 (a) applied to the second term (which is dominated by the
maximum in the assumption as can be seen by an index shift),

$\displaystyle\min_{k_{j,n}-(1-\varepsilon)\,G\leq k\leq
k_{j,n}}T^{(1)}_{k,n}(G)$
$\displaystyle\geq\min_{k_{j,n}-(1-\varepsilon)\,G\leq k\leq
k_{j,n}}\sqrt{\frac{G}{2}}\left\|\boldsymbol{\Gamma}_{k}^{-1/2}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)\right\|$
$\displaystyle\qquad-\max_{k_{j,n}-(1-\varepsilon)\,G\leq k\leq
k_{j,n}}\sqrt{\frac{G}{2}}\left\|\boldsymbol{\Gamma}_{k}^{-1/2}\left(\boldsymbol{\widehat{\theta}}_{k-G+1,k}-\boldsymbol{\theta}_{j}\right)\right\|$
$\displaystyle\geq\left\|\boldsymbol{\Gamma}_{(j)}^{1/2}\right\|^{-1}\,\min_{k_{j,n}-(1-\varepsilon)\,G\leq
k\leq
k_{j,n}}\sqrt{\frac{G}{2}}\left\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right\|$
$\displaystyle\qquad-\left\|\boldsymbol{\Gamma}_{(j)}^{-1/2}\right\|\,\max_{k_{j,n}-(1-\varepsilon)\,G\leq
k\leq
k_{j,n}}\sqrt{\frac{G}{2}}\left\|\boldsymbol{\widehat{\theta}}_{k-G+1,k}-\boldsymbol{\theta}_{j}\right\|$
$\displaystyle=\left\|\boldsymbol{\Gamma}_{(j)}^{1/2}\right\|^{-1}\,\min_{k_{j,n}-(1-\varepsilon)\,G\leq
k\leq
k_{j,n}}\sqrt{\frac{G}{2}}\left\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right\|+O_{P}\left(\sqrt{\log(n/G)}\right).$

We get by Assumption 2.5 and 3.1 (b)

$\displaystyle P\left(\min_{j=1,\ldots,q}\min_{k_{j,n}-(1-\varepsilon)\,G\leq
k\leq k_{j,n}}T^{(1)}_{k,n}(G)<D_{n}(\alpha_{n},G)\right)$
$\displaystyle\leq\sum_{j=1}^{q}P\left(\min_{k_{j,n}-(1-\varepsilon)\,G\leq
k\leq
k_{j,n}}\sqrt{\frac{G}{2}}\left\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right\|<O_{P}\left(\sqrt{\log(n/G)}\right)\right)$
$\displaystyle\to 0.$

Similar arguments deal with the minimum over $k_{j,n}<k\leq
k_{j,n}+(1-\varepsilon)\,G$, such that assertion (a)(ii) follows.

The proof of (b) follows along the same lines taking Assumption 3.2 into
account. ∎

###### Proof of Remark 3.1.

The proof is analogous to the above proofs where we use that by Theorem 2.1

$\displaystyle\max_{j=1,\ldots,q+1}\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}\left\|T^{(1)}_{k,n}(G)\right\|=O_{P}\left(\sqrt{\log(n/G)}\right).$

∎

###### Proof of Theorem 3.2.

The assertion follows immediately from Proposition 3.1 on noting that

$\displaystyle\left\\{\max_{j=1,\ldots,q+1}\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}T^{(1)}_{k,n}(G)<D_{n}(\alpha_{n},G)\right\\}$
$\displaystyle\quad\cap\left\\{\min_{k_{j,n}-(1-\varepsilon)\,G\leq k\leq
k_{j,n}+(1-\varepsilon)\,G}T^{(1)}_{k,n}(G)\geq D_{n}(\alpha_{n},G)\right\\}$
$\displaystyle\subset\left\\{\widehat{q}_{n}^{(1)}=q\right\\}\cap\left\\{\max_{1\leq
j\leq q}\left|\widehat{k}^{(1)}_{j,n}-k_{j,n}\right|<G\right\\}.$

∎

##  Appendix C Proofs of Section 3.2

###### Proof of Lemma 3.3.

We will show by contradiction that
$E\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\widetilde{\theta}}_{0,1})\neq
E\boldsymbol{H}(\mathbb{X}_{1}^{(j+1)},\boldsymbol{\widetilde{\theta}}_{0,1})$
holds for at least one $j\in\\{1,\ldots,q\\}$. Assume that all these
expectations are equal, then by definition of
$\boldsymbol{\widetilde{\theta}}_{0,1}$ we get for all $j=1,\ldots,q+1$

$\displaystyle\boldsymbol{0}=\sum_{l=1}^{q+1}\left(\lambda_{l}-\lambda_{l-1}\right)E\boldsymbol{H}(\mathbb{X}_{1}^{(l)},\boldsymbol{\widetilde{\theta}}_{0,1})=E\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\widetilde{\theta}}_{0,1}),$

which by the identifiability of $\boldsymbol{\theta}_{j}$ implies
$\boldsymbol{\theta}_{j}=\boldsymbol{\widetilde{\theta}}_{0,1}$, for all
$j=1,\ldots,q+1$, contradicting the assumption. If there are only two possible
regimes, then clearly if one change is detectable all of them are. ∎

###### Proof of Proposition 3.4.

Analogously to the proof of Proposition 3.1 (a) (i) we get

$\displaystyle P\left(\max_{j=1,\ldots,q+1}\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}})\geq
D_{n}(\alpha_{n},G)\right)\to 0.$

The statement remains true when a sequence of inspection parameters
$\boldsymbol{\widetilde{\theta}}_{n}$ is used because by Assumption 2.4 (i) it
holds for any $j=1,\ldots,q$ that

$\displaystyle a(n/G)\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})$
$\displaystyle=a(n/G)\max_{k_{j-1,n}+G\leq k\leq
k_{j,n}-G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}})+o_{P}(1).$

Additionally, we need the statement for environments of non-detectable change
points $k_{j,n}$ with $j\not\in\tilde{Q}$. Indeed, it holds for
$|k-k_{j,n}|\leq G$ with $j\not\in\tilde{Q}$ that
$E\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)=0$ by definition of
$\tilde{Q}$ as in (3.2.1). Consequently, by Assumptions 2.1 and 2.3 it holds
for $j\not\in\tilde{Q}$

$\displaystyle\max_{k_{j,n}\leq
k<k_{j,n}+G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}})$
$\displaystyle=o_{P}(1)+O_{P}\left(\max_{k_{j,n}\leq
k<k_{j,n}+G}\frac{1}{\sqrt{G}}\left\|\boldsymbol{W}(k+G)-2\boldsymbol{W}(k)+\boldsymbol{W}(k_{j,n})\right\|\right)$
$\displaystyle\qquad+O_{P}\left(\max_{k_{j,n}\leq
k<k_{j,n}+G}\frac{1}{\sqrt{G}}\left\|\boldsymbol{\Sigma}^{-1/2}_{(j+1)}\boldsymbol{\Sigma}^{1/2}_{(j)}\left(\boldsymbol{W}(k_{j,n})-\boldsymbol{W}(k-G))\right)\right\|\right)$
$\displaystyle=O_{P}(1)=o_{P}(D_{n}(\alpha_{n},G)),$ (11)

where the last line follows by the self-similarity of Wiener processes, the
stationarity of its increments and the continuous sample paths. An analogous
assertion holds for $k_{j,n}-G\leq k<k_{j,n}$ showing that

$\displaystyle
P\left(\max_{j\not\in\tilde{Q}}\max_{|k-k_{j,n}|<G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}})\geq
D_{n}(\alpha_{n},G)\right)\to 0,$

completing the proof of (a) (i) for a fixed inspection parameter
$\boldsymbol{\widetilde{\theta}}$. Here, the statement remains true for a
sequence of inspection parameters, because by Assumption 2.4 (ii) it holds for
$j\not\in\tilde{Q}$

$\displaystyle\max_{|k-k_{j,n}|<G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})=\max_{|k-k_{j,n}|<G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}})+o_{P}(\sqrt{\log(n/G)})$
(12) $\displaystyle=o_{P}(D_{n}(\alpha_{n},G)).$

The assertion with estimated long-run covariances as in (b) can be obtained
along the same lines by using the consistency of the spectral matrix norm with
the Euclidean vector norm and Assumptions 3.3 (b).

Concerning (ii) first observe that for
$\tilde{k}_{j,n}<k\leq\tilde{k}_{j,n}+(1-\varepsilon)G$,
$j=1,\ldots,\tilde{q}(\boldsymbol{\widetilde{\theta}})$ it holds

$\displaystyle E\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)$
$\displaystyle=\sum_{i=k+1}^{k+G}E\boldsymbol{H}(\mathbb{X}_{i}^{(j+1)},\boldsymbol{\widetilde{\theta}})-\sum_{i=k-G+1}^{\tilde{k}_{j,n}}E\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}})-\sum_{i=\tilde{k}_{j,n}+1}^{k}E\boldsymbol{H}(\mathbb{X}_{i}^{(j+1)},\boldsymbol{\widetilde{\theta}})$
$\displaystyle=\left(G-|k-\tilde{k}_{j,n}|\right)\,\boldsymbol{d}_{j},$

where we denote the signal by

$\displaystyle\boldsymbol{d}_{j}=E\boldsymbol{H}(\mathbb{X}_{1}^{(j+1)},\boldsymbol{\widetilde{\theta}})-E\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\widetilde{\theta}}),\quad
j=1,\ldots\tilde{q}(\boldsymbol{\widetilde{\theta}}).$ (13)

For $\tilde{k}_{j,n}-(1-\varepsilon)G\leq k\leq\tilde{k}_{j,n}$ we arrive at
the same conclusion. Consequently, it holds for all $j\in\tilde{Q}$ and
$|k-\tilde{k}_{j,n}|\leq(1-\varepsilon)G$ by the consistency of the spectral
matrix with the Euclidean vector norm

$\displaystyle\left\lVert\boldsymbol{\Sigma}_{k}^{-1/2}E\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)\right\rVert\geq\left\|\boldsymbol{\Sigma}_{k}^{1/2}\right\|^{-1}\,\left(G-|k-\tilde{k}_{j,n}|\right)\,\left\lVert\boldsymbol{d}_{j}\right\rVert\geq
c\,G,$

for some $c>0$ (depending on $\varepsilon$, the difference in expectation and
the long-run covariances, noting that $\boldsymbol{\Sigma}_{k}$ is constant on
each segment).

By analogous arguments as in (C) (but involving the necessary centering due to
$\tilde{Q}$) and (12) it holds

$\displaystyle\min_{|k-\tilde{k}_{j,n}|\leq(1-\varepsilon)\,G}T^{(2)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})$
$\displaystyle=\min_{|k-\tilde{k}_{j,n}|\leq(1-\varepsilon)\,G}\frac{1}{\sqrt{2G}}\,\left\lVert\boldsymbol{\Sigma}_{k}^{-1/2}E\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}}(k)\right\rVert+o_{P}(\sqrt{\log(n/G)})$
$\displaystyle\geq c\,\sqrt{\frac{G}{2}}+o_{P}\left(\sqrt{\log(n/G)}\right).$

The proof can now be concluded as in the proof of Proposition 3.1 (a) and (b)
(ii). ∎

###### Proofs of Theorem 3.5 and Remark 3.4.

The proofs are completely analogous to the proofs of Theorem 3.2 respectively
Remark 3.1 and therefore omitted. ∎

The following lemma helps simplify several arguments.

######  Lemma C.1.

For a sequence of real random variables $\\{X_{n}\\}$ it holds for
$n\to\infty$

$\displaystyle X_{n}=O_{P}(1)\quad\iff\quad P\left(|X_{n}|>\xi_{n}\right)\to
0\quad\text{for any }\xi_{n}\to\infty.$

###### Proof.

The proof of the only-if-part is straightforward. We prove the if-part by
contradiction. If $X_{n}$ is not stochastically bounded, then there exists
$\eta>0$ such that for any bound $C>0$ and any $n_{0}\geq 0$, there exists
$n_{1}(n_{0},C)>n_{0}$ such that $P\left(|X_{n_{1}}|>C\right)>\eta.$ Setting
$N_{0}=0$ and recursively $N_{l}=n_{1}(N_{l-1},l)$ as well as
$\xi_{N_{l-1}+1}=\ldots=\xi_{N_{l}}=l$, we get $N_{n}\to\infty$,
$\xi_{n}\to\infty$ as well as by construction

$\displaystyle P\left(|X_{N_{l}}|>\xi_{N_{l}}\right)>\eta,$

which is a contradiction. ∎

The proof technique of the below proof is well known in change point analysis,
for example it has been used in the context of MOSUM statistics for the mean
change problem by Eichinger & Kirch, (2018) (Proof of Theorem 3.2).

###### Proof of Theorem 3.6.

By finiteness of $q$ and Lemma C.1 it is sufficient to prove that for any
sequence $\xi_{n}\to\infty$ (arbitrarily slow) it holds

$\displaystyle
P\left(\widehat{k}^{(2)}_{j,n}(\boldsymbol{\widetilde{\theta}}_{n};\widehat{\boldsymbol{\Psi}}_{j,n})<\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})-\xi_{n}\right)\to
0,$ $\displaystyle
P\left(\widehat{k}^{(2)}_{j,n}(\boldsymbol{\widetilde{\theta}}_{n};\widehat{\boldsymbol{\Psi}}_{j,n})>\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})+\xi_{n}\right)\to
0.$

We will prove the first assertion in detail, the second one follows
analogously. For simplicity of notation denote
$\tilde{k}_{j,n}=\tilde{k}_{j,n}(\boldsymbol{\widetilde{\theta}})$ throughout
this proof.

On the asymptotic 1-set of Theorem 3.5 and where
$\min_{j=1,\ldots,q+1}|k_{j,n}-k_{j-1,n}|>2G$ (which holds for $n$ large
enough by Assumption 2.3 (b)), it holds for any
$j=1,\ldots,\tilde{q}(\boldsymbol{\widetilde{\theta}})$ with
$\boldsymbol{d}_{j}$ as in (13)

$\displaystyle\widehat{k}_{j,n}^{(2)}(\boldsymbol{\widetilde{\theta}}_{n};\widehat{\boldsymbol{\Psi}}_{j,n})=\operatorname*{arg\,max}_{v_{j,n}\leq
k\leq
w_{j,n}}V^{(j)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n}),\quad\text{where
}$ $\displaystyle
V^{(j)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})=\left\lVert\boldsymbol{\widehat{\Psi}}_{j,n}^{-1/2}\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k)\right\rVert^{2}-\left\lVert\boldsymbol{\widehat{\Psi}}_{j,n}^{-1/2}\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k_{j,n})\right\rVert^{2}$
$\displaystyle=-\left(\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(\tilde{k}_{j,n})-\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k)\right)\,\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\,\left(\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(\tilde{k}_{j,n})+\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k)\right)$
$\displaystyle=:-\left(\boldsymbol{E}_{1}(k,G,\boldsymbol{\widetilde{\theta}}_{n})+\boldsymbol{d}_{j}(\tilde{k}_{j,n}-k)\right)\,\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\,\left(\boldsymbol{E}_{2}(k,G,\boldsymbol{\widetilde{\theta}}_{n})+\boldsymbol{d}_{j}(2G+k-\tilde{k}_{j,n})\right).$

Denote
$\Delta\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n})=\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n})-\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}})$
and
$\boldsymbol{H}_{0}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}})=\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}})-E\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}})$.
Then, it holds for $k<\tilde{k}_{j,n}$

$\displaystyle\boldsymbol{E}_{1}(k,G,\boldsymbol{\widetilde{\theta}}_{n})=\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(\tilde{k}_{j,n})-\boldsymbol{M}_{\boldsymbol{\widetilde{\theta}}_{n}}(k)-\boldsymbol{d}_{j}(\tilde{k}_{j,n}-k)$
$\displaystyle=\sum_{i=k+G+1}^{\tilde{k}_{j,n}+G}\Delta\boldsymbol{H}(\mathbb{X}_{i}^{(j+1)},\boldsymbol{\widetilde{\theta}}_{n})+\sum_{i=k-G+1}^{\tilde{k}_{j,n}-G}\Delta\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n})-2\sum_{i=k+1}^{\tilde{k}_{j,n}}\Delta\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n})$
$\displaystyle\qquad+\sum_{i=k+G+1}^{\tilde{k}_{j,n}+G}\boldsymbol{H}_{0}(\mathbb{X}_{i}^{(j+1)},\boldsymbol{\widetilde{\theta}})+\sum_{i=k-G+1}^{\tilde{k}_{j,n}-G}\boldsymbol{H}_{0}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}})-2\sum_{i=k+1}^{\tilde{k}_{j,n}}\boldsymbol{H}_{0}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}),$

where by Assumption 3.5 (see also Remark 3.6 for the situation when
$k>\tilde{k}_{j,n}$) and stationarity of the segments it follows

$\displaystyle\sup_{\xi_{n}<\tilde{k}_{j,n}-k\leq
G}\frac{\|\boldsymbol{E}_{1}(k,G,\boldsymbol{\widetilde{\theta}}_{n})\|}{\tilde{k}_{j,n}-k}=o_{P}(1).$
(14)

Furthermore,

$\displaystyle\boldsymbol{E}_{2}(k,G,\boldsymbol{\widetilde{\theta}}_{n})$
$\displaystyle=-\boldsymbol{E}_{1}(k,G,\boldsymbol{\widetilde{\theta}}_{n})$
$\displaystyle\quad+2\sum_{i=\tilde{k}_{j,n}+1}^{\tilde{k}_{j,n}+G}\Delta\boldsymbol{H}(\mathbb{X}_{i}^{(j+1)},\boldsymbol{\widetilde{\theta}}_{n})-2\sum_{i=\tilde{k}_{j,n}-G+1}^{\tilde{k}_{j,n}}\Delta\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n})$
$\displaystyle\quad+2\sum_{i=\tilde{k}_{j,n}+1}^{\tilde{k}_{j,n}+G}\boldsymbol{H}_{0}(\mathbb{X}_{i}^{(j+1)},\boldsymbol{\widetilde{\theta}})-2\sum_{i=\tilde{k}_{j,n}-G+1}^{\tilde{k}_{j,n}}\boldsymbol{H}_{0}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}),$

such that by (14), Assumption 3.5 (a) and stationarity of the segments in
combination with the law of large number (that follows from Assumption 2.1,
see also Remark 3.6) it holds

$\displaystyle\sup_{\xi_{n}<\tilde{k}_{j,n}-k\leq
G}\frac{\|\boldsymbol{E}_{2}(k,G,\boldsymbol{\widetilde{\theta}}_{n})\|}{G}=o_{P}(1).$

By Assumptions 3.4 and the consistency of the spectral matrix norm with the
Euclidean vector norm we conclude

$\displaystyle\sup_{\xi_{n}<\tilde{k}_{j,n}-k\leq
G}V^{(j)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})$
$\displaystyle\leq(\boldsymbol{d}_{j}\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\boldsymbol{d}_{j}+o_{P}(1))\sup_{\xi_{n}<\tilde{k}_{j,n}-k\leq
G}[-(\tilde{k}_{j,n}-k)(2G+k-\tilde{k}_{j,n})].$

Finally,

$\displaystyle\sup_{\xi_{n}<\tilde{k}_{j,n}-k\leq
G}[-(\tilde{k}_{j,n}-k)(2G+k-\tilde{k}_{j,n})]\leq-\xi_{n}G<0,$

and

$\displaystyle\boldsymbol{d}_{j}\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\boldsymbol{d}_{j}=\|\boldsymbol{\widehat{\Psi}}_{j,n}^{-1/2}\boldsymbol{d}_{j}\|^{2}\geq\|\boldsymbol{\widehat{\Psi}}_{j,n}^{1/2}\|^{-2}\|\boldsymbol{d}_{j}\|^{2}$

Consequently,

$\displaystyle
P\left(\widehat{k}^{(2)}_{j,n}(\boldsymbol{\widetilde{\theta}}_{n};\widehat{\boldsymbol{\Psi}}_{j,n})<\tilde{k}_{j,n}-\xi_{n}\right)$
$\displaystyle\leq P\left(\sup_{\tilde{k}_{j,n}-G\leq
k<\tilde{k}_{j,n}-\xi_{n}}V^{(j)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})\geq\sup_{\tilde{k}_{j,n}-\xi_{n}\leq
k\leq\tilde{k}_{j,n}+G}V^{(j)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})\right)+o(1)$
$\displaystyle\leq P\left(\sup_{\tilde{k}_{j,n}-G\leq
k<\tilde{k}_{j,n}-\xi_{n}}V^{(j)}_{k,n}(G,\boldsymbol{\widetilde{\theta}}_{n})\geq
0\right)+o(1)$ $\displaystyle\leq
P\left((\boldsymbol{d}_{j}\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\boldsymbol{d}_{j}+o_{P}(1))\sup_{\xi_{n}<\tilde{k}_{j,n}-k\leq
G}[-(\tilde{k}_{j,n}-k)(2G+k-\tilde{k}_{j,n})]\geq 0\right)+o(1)$
$\displaystyle\leq
P\left(|o_{P}(1)|\geq\boldsymbol{d}_{j}\boldsymbol{\widehat{\Psi}}_{j,n}^{-1}\boldsymbol{d}_{j}\right)+o(1)\leq
P\left(|o_{P}(1)|\,\|\boldsymbol{\widehat{\Psi}}_{j,n}^{1/2}\|^{2}\geq\|\boldsymbol{d}_{j}\|^{2}\right)+o(1)$
$\displaystyle=o(1),$

where the last line follows from Assumption 3.4, concluding the proof. ∎

##  Appendix D Proofs of Section 4

###### Proof of Theorem 4.1.

The assertion in (a) follows by Theorem 4 of Kuelbs & Philipp, (1980) on
noting that the mixing rate of
$\boldsymbol{H}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta})$ is at least as good
as the one of $\\{\mathbb{X}_{1}^{(j)}\\}$ by definition. Because the time
series in backward time is also mixing with the same rate, Assumption 3.5 (b)
follow from the invariance principle in backward time (see also Remark 3.6). ∎

The remaining assumptions all correspond to well known results in statistics
if a global estimator based on estimating functions is used. However, here, it
is maximized over an increasing number of windows. To this end, we require
versions of uniform laws of large numbers taking these moving windows into
account as given in the following lemma.

######  Lemma D.1.

Let $\\{\mathbb{Y}_{t}\\}$ be $p$-dimensional random vectors fulfilling
Regularity Condition 4.1 (b) with $\tilde{\nu}$ as below (in (a),(b),(d) and
arbitrary in (c)), $\boldsymbol{\Theta}\subset\mathbb{R}^{p}$ be a compact
parameter space and
$\boldsymbol{F}=(F_{1},\ldots,F_{p})^{T}:(\mathbb{R}^{p},\boldsymbol{\Theta})\to\mathbb{R}^{p}$
measurable.

  1. (a)

If
$0<E\left\|\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\right\|^{2+\tilde{\nu}}<\infty$
for some $\tilde{\nu}>0$ and some $\boldsymbol{\theta}$, then for the same
$\boldsymbol{\theta}$

$\displaystyle\sup_{0\leq k\leq
n-G}\left\|\sum_{i=k+1}^{k+G}\left(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\left(\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\right)\right)\right\|=O_{P}\left(\sqrt{G\,\log(n/G)}\right),$

  2. (b)

If for some $\tilde{\nu}>0$ it holds
$0<E\left\|\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\right\|^{2+\tilde{\nu}}<\infty$
for all $\boldsymbol{\theta}$ as well as
$E\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\|\nabla\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\|^{2+\tilde{\nu}}<\infty,$
then

$\displaystyle\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\max_{0\leq
k\leq
n-G}\frac{1}{G}\left\lVert\sum_{i=k+1}^{k+G}(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\left(\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\right))\right\rVert=o_{P}(1).$

  3. (c)

If
$E\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\|\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\|<\infty$,
then for any sequence $G\to\infty$ it holds

$\displaystyle(i)\quad$
$\displaystyle\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\max_{1\leq
k\leq
G}\frac{1}{k}\left\lVert\sum_{i=G-k+1}^{G}\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\right\rVert=O_{P}(1),\quad$
$\displaystyle\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\max_{1\leq
k\leq
G}\frac{1}{k}\left\lVert\sum_{i=1}^{k}\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\right\rVert=O_{P}(1).$
$\displaystyle(ii)\quad$
$\displaystyle\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\max_{1\leq
k\leq
G}\frac{1}{G}\left\lVert\sum_{i=1}^{k}(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\left(\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\right))\right\rVert=o_{P}(1),$
$\displaystyle\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\max_{1\leq
k\leq
G}\frac{1}{G}\left\lVert\sum_{i=k}^{G}(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\left(\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\right))\right\rVert=o_{P}(1).$

  4. (d)

If
$E\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\|\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\|^{2+\tilde{\nu}}<\infty$,
then

$\displaystyle\sup_{0\leq k\leq
n-G}\sum_{i=k+1}^{k+G}\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\left\|\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\right\|=O_{P}\left(G\right).$

###### Proof.

Analogously to Theorem 4.1
$\\{\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\\}$ fulfills an
invariance principle from which assertion (a) follows by similar arguments as
in the proof of Theorem 2.1 (see Reckrühm, (2019), Theorem E.2.12 for
details).

The proof technique for (b) is well known (and we only use a basic version
thereof). Thus, we only sketch the proof. First note that by the compactness
assumption on $\boldsymbol{\Theta}$ for each $\delta>0$ there exist
$M=M(\delta)\geq 1$ and
$\boldsymbol{\xi}_{1},\ldots,\boldsymbol{\xi}_{M}\in\boldsymbol{\Theta}$ such
that for any $\boldsymbol{\theta}\in\boldsymbol{\Theta}$ there is an
$m=1,\ldots,M$ with
$\left\lVert\boldsymbol{\theta}-\boldsymbol{\xi}_{m}\right\rVert<\delta$. We
get for any $\boldsymbol{\xi},\boldsymbol{\theta}$

$\displaystyle\max_{0\leq k\leq
n-G}\frac{1}{G}\sum_{i=k+1}^{k+G}\left\lVert\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-\left(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi})-E\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi})\right)\right\rVert$
$\displaystyle\leq\left(2E\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\|\nabla\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\|+o_{P}(1)\right)\,\|\boldsymbol{\theta}-\boldsymbol{\xi}\|=\|\boldsymbol{\theta}-\boldsymbol{\xi}\|\,O_{P}(1),$

where the last line follows from a first order Taylor expansion in addition to
a moving law of large numbers as in (a) applied to the time series
$\\{\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\|\nabla\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\|\\}$
(see also (d)). For given $\eta_{1},\eta_{2}>0$ we can now choose
$\delta=\delta(\eta_{1},\eta_{2})>0$ such that

$\displaystyle
P\Big{(}\sup_{\|\boldsymbol{\theta}-\boldsymbol{\xi}\|<\delta}\max_{0\leq
k\leq
n-G}\frac{1}{G}\sum_{i=k+1}^{k+G}\left\lVert\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\right.$
$\displaystyle\phantom{P(\max_{\|\boldsymbol{\theta}-\boldsymbol{\xi}\|<\delta}\max_{0\leq
k\leq
n-G}\frac{1}{G}\sum_{i=k+1}^{k+G}}\quad\left.-\left(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi})-E\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi})\right)\right\rVert\geq\eta_{1}\Big{)}\leq\eta_{2}$

for all $n\geq n_{0}(\eta_{1},\eta_{2})$. For the (to $\delta$) corresponding
$\boldsymbol{\xi_{1}},\ldots,\boldsymbol{\xi_{M}}$ it holds by another
application of (a)

$\displaystyle P\left(\max_{m=1,\ldots,M}\max_{0\leq k\leq
n-G}\frac{1}{G}\left\lVert\sum_{i=k+1}^{k+G}(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi}_{m})-E\left(\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\xi}_{m})\right))\right\rVert\geq\eta_{1}\right)\leq\eta_{2}.$

for all $n\geq n_{1}(\eta_{1},\eta_{2})$. Combining these arguments yields (b)
on noting that for any $\delta$ and corresponding $\boldsymbol{\xi}_{m}$,
$m=1,\ldots,M$, it holds

$\displaystyle\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\max_{0\leq
k\leq
n-G}\frac{1}{G}\left\lVert\sum_{i=k+1}^{k+G}(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\left(\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\theta})\right))\right\rVert$
$\displaystyle\leq\max_{m=1,\ldots,M}\max_{0\leq k\leq
n-G}\frac{1}{G}\left\lVert\sum_{i=k+1}^{k+G}(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi}_{m})-E\left(\boldsymbol{F}(\mathbb{Y}_{1},\boldsymbol{\xi}_{m})\right))\right\rVert$
$\displaystyle\qquad+\sup_{\|\boldsymbol{\theta}-\boldsymbol{\xi}\|<\delta}\max_{0\leq
k\leq
n-G}\frac{1}{G}\sum_{i=k+1}^{k+G}\left\lVert\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})-E\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\right.$
$\displaystyle\phantom{\qquad+\sup_{\|\boldsymbol{\theta}-\boldsymbol{\xi}\|<\delta}\max_{0\leq
k\leq
n-G}\frac{1}{G}\sum_{i=k+1}^{k+G}\quad}\left.-\left(\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi})-E\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\xi})\right)\right\rVert.$

By Rao, (1962), Theorem 6.5, a uniform (in $\boldsymbol{\theta}$) strong law
of large numbers holds for
$\\{\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\\}$ because
stationarity and mixing implies ergodicity (both forward and backward). By the
almost sure convergence standard arguments give the assertions in (c).

The proof of (d) follows along the same lines as the proof of (a) but applied
to the function
$\sup_{\boldsymbol{\theta}\in\boldsymbol{\Theta}}\left\|\boldsymbol{F}(\mathbb{Y}_{i},\boldsymbol{\theta})\right\|$.
The necessary centering is of the order $G$. ∎

###### Proof of Theorem 4.2.

For $k_{j-1,n}<k\leq k_{j,n}-G$ a Taylor expansion in
$\boldsymbol{\widehat{\theta}}_{k+1,k+G}$ around $\boldsymbol{\theta}_{j}$
yields that there exists
$\left\lVert\boldsymbol{\xi}^{(j)}_{k,n}-\boldsymbol{\theta}_{j}\right\rVert\leq\left\lVert\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right\rVert$
such that

$\displaystyle-\frac{1}{\sqrt{G}}\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\theta_{j}})$
$\displaystyle=\left(\frac{1}{G}\sum_{i=k+1}^{k+G}\nabla\boldsymbol{H}(\mathbb{X}^{(j)}_{i},\boldsymbol{\xi}_{k,n}^{(j)})\right)^{T}\sqrt{G}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)$
$\displaystyle=\left(o_{P}(1)+\boldsymbol{V}_{(j)}(\boldsymbol{\xi}_{k,n}^{(j)})\right)\sqrt{G}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)\text{
uniformly in }k,$

where the last line follows by Regularity Conditions 4.2 (a) in combination
with Lemma D.1 (b). By Lemma D.1 (a), Regularity Condition 4.1 and the
definition of $\boldsymbol{\theta}_{j}$ we get

$\displaystyle\sup_{1\leq k\leq
n-G}\frac{1}{\sqrt{G}}\left\|\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\theta_{j}})\right\|=O_{P}\left(\sqrt{\log(n/G)}\right).$

In combination with Regularity Conditions 4.2(b) this yields

$\displaystyle\max_{j=1,\ldots,q+1}\max_{k_{j-1,n}<k\leq
k_{j,n}-G}\sqrt{G}\left\lVert\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right\rVert=O_{P}\left(\sqrt{\log(n/G)}\right).$
(15)

In particular, this shows the validity of Assumption 3.1 (a).
Moreover, for each $l=1,\ldots,p$ and $k_{j-1,n}<k\leq k_{j,n}-G$ a second
order Taylor expansion yields the existence of
$\left\lVert\boldsymbol{\xi}^{(j)}_{l,n,k}-\boldsymbol{\theta}_{j}\right\rVert\leq\left\lVert\widehat{\boldsymbol{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right\rVert$
with

$\displaystyle-\sum_{i=k+1}^{k+G}H_{l}(\mathbb{X}_{i}^{(j)},\boldsymbol{\theta}_{j})$
$\displaystyle=\left(\sum_{i=k+1}^{k+G}\left(\nabla
H_{l}(\mathbb{X}_{i}^{(j)},\boldsymbol{\theta}_{j})\right)\right)^{T}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)$
$\displaystyle\quad+\frac{1}{2}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)^{T}\left(\sum_{i=k+1}^{k+G}\nabla^{2}H_{l}(\mathbb{X}^{(j)}_{i},\boldsymbol{\xi}^{(j)}_{l,n,k})\right)\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)$
$\displaystyle=\left(E\left(\nabla
H_{l}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta}_{j})\right)+O_{P}\left(\sqrt{\frac{\log(n/G)}{G}}\right)\right)^{T}\;G\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)$
$\displaystyle\quad+O_{P}\left(G\right)\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\|^{2}\quad\text{uniformly
in }k,$

where the last line follows by an application of Lemma D.1 both (a) and (d).
Thus, an application of (15) yields uniformly in $k$

$\displaystyle-\frac{1}{\sqrt{2G}}\sum_{i=k+1}^{k+G}H_{l}(\mathbb{X}_{i}^{(j)},\boldsymbol{\theta}_{j})$
$\displaystyle=E\left(\nabla
H_{l}(\mathbb{X}_{1}^{(j)},\boldsymbol{\theta}_{j})\right)^{T}\sqrt{\frac{G}{2}}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{\theta}_{j}\right)+o_{P}\left(\left(\log
n/G\right)^{-1/2}\right),$

showing the validity of Assumption 2.2, concluding the proof of (a). By the
strong law of large numbers (and similar arguments as in the proof of Lemma
D.1 (c)) and by definition of $\boldsymbol{\theta}_{j}$ it holds

$\displaystyle\max_{k_{j-1,n}-G\leq k\leq
k_{j-1,n}-\varepsilon\,G}\left\|\frac{1}{G}\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\theta}_{j})-\frac{k_{j-1,n}-k}{G}E\boldsymbol{H}(\mathbb{X}_{i}^{(j-1)},\boldsymbol{\theta}_{j})\right\|$
$\displaystyle=\max_{k_{j-1,n}-G\leq k\leq
k_{j-1,n}-\varepsilon\,G}\left\|\frac{1}{G}\sum_{i=k+1}^{k_{j-1,n}}\left(\boldsymbol{H}(\mathbb{X}_{i}^{(j-1)},\boldsymbol{\theta}_{j})-E\boldsymbol{H}(\mathbb{X}_{i}^{(j-1)},\boldsymbol{\theta}_{j})\right)\right.$
$\displaystyle\phantom{=\max_{k_{j-1,n}-G\leq k\leq
k_{j-1,n}-\varepsilon\,G}}\quad+\left.\frac{1}{G}\sum_{i=k_{j-1,n}+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\theta}_{j})\right\|=o_{P}(1).$

By the identifiable uniqueness of $\boldsymbol{\theta}_{j}$ it holds
$E\boldsymbol{H}(\mathbb{X}_{i}^{(j-1)},\boldsymbol{\theta}_{j})\neq 0$, such
that

$\displaystyle\sqrt{\frac{G}{\log(n/G)}}\,\min_{k_{j-1,n}-G\leq k\leq
k_{j-1,n}-\varepsilon\,G}\left\|\frac{1}{G}\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\theta}_{j})\right\|$
$\displaystyle\geq\sqrt{\frac{G}{\log(n/G)}}\left(\varepsilon\,\left\|E\boldsymbol{H}(\mathbb{X}_{i}^{(j-1)},\boldsymbol{\theta}_{j})\right\|+o_{P}(1)\right)\overset{P}{\longrightarrow}\infty.$

By Lemma D.1 (c) it holds

$\displaystyle\sup_{\boldsymbol{\theta}\in\Theta}\max_{k_{j-1,n}-G\leq k\leq
k_{j-1,n}-\varepsilon\,G}\left\|\frac{1}{G}\sum_{i=k+1}^{k+G}\nabla\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\theta})\right\|_{F}=O_{P}(1),$

where $\|\cdot\|_{F}$ denotes the Frobenius matrix norm.

Furthermore, a Taylor expansion of $\boldsymbol{\widehat{\theta}}_{k+1,k+G}$
around $\boldsymbol{\theta}_{j}$ yields for some $\boldsymbol{\xi}_{k,n}$

$\displaystyle\sqrt{\frac{G}{\log(n/G)}}\,\left\|\frac{1}{G}\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\theta}_{j})\right\|$
$\displaystyle=\left\|\left(\frac{1}{G}\sum_{i=k+1}^{k+G}\nabla\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\xi}_{k,n})\right)^{T}\sqrt{\frac{G}{\log(n/G)}}\left(\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{{\theta}}_{j}\right)\right\|$
$\displaystyle\leq\left\|\frac{1}{G}\sum_{i=k+1}^{k+G}\nabla\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\xi}_{k,n})\right\|_{F}\,\sqrt{\frac{G}{\log(n/G)}}\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{{\theta}}_{j}\|,$

where we used the consistency of the Frobenius matrix norm with the Euclidean
vector norm in the last step. Thus, we have shown that the left hand side
diverges to infinity stochastically, while the first term on the right hand
side is stochastically bounded – both in an appropriate uniform sense. By
standard arguments this shows that indeed

$\displaystyle\sqrt{\frac{G}{\log(n/G)}}\min_{k_{j-1,n}-G\leq k\leq
k_{j-1,n}-\varepsilon\,G}\|\boldsymbol{\widehat{\theta}}_{k+1,k+G}-\boldsymbol{{\theta}}_{j}\|\overset{P}{\longrightarrow}\infty.$

The assertion for $k_{j,n}-(1-\varepsilon)\,G\leq k\leq k_{j,n}$ follows
analogously concluding the proof.

∎

###### Proof of Theorem 4.3.

The proof is analogous to the proof of (15) where the sum

$-\frac{1}{\sqrt{n}}\sum_{i=a}^{b}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\widetilde{\theta}}_{\gamma_{a},\gamma_{b}})$
is considered instead. The better rate compared to (15) is due to the fact
that a piecewise application of the central limit theorem (which follows from
the mixing condition) to each regime yields

$\displaystyle\frac{1}{\sqrt{n}}\left\|\sum_{i=a}^{b}\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\widetilde{\theta}}_{\gamma_{a},\gamma_{b}})\right\|=O_{P}\left(1\right).$

The convergence of
$\frac{1}{n}\sum_{i=a}^{b}\nabla\boldsymbol{H}(\mathbb{X}_{i},\boldsymbol{\xi}_{n}^{(\gamma_{a},\gamma_{b})})$
under these regularity conditions follows also e.g. from a piecewise
application of Lemma D.1 (c)(ii) (with $G$ replaced by $n$). ∎

###### Proof of Theorem 4.4.

A first order Taylor expansion yields

$\displaystyle\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}}_{n})-\sum_{i=k+1}^{k+G}\boldsymbol{H}(\mathbb{X}_{i}^{(j)},\boldsymbol{\widetilde{\theta}})$
$\displaystyle=\left(\sum_{i=k+1}^{k+G}\nabla\boldsymbol{H}(\mathbb{X}^{(j)}_{i},\boldsymbol{\xi}_{k,n}^{(j)})\right)^{T}\left(\boldsymbol{\widetilde{\theta}}_{n}-\boldsymbol{\widetilde{\theta}}\right)=O_{P}(G/\sqrt{n})=o_{P}\left(\sqrt{\frac{G}{\log(n/G)}}\right),$

where the last line follows uniformly in $k$ by Lemma D.1 (b). This shows the
validity of Assumption 2.4 (i). The proof of (ii) and Assumption 3.5 (a) are
analogous by using Lemma D.1(c) instead.

∎