Bayesian Loss-based Approach to Change
Point Analysis
arXiv:1702.05462v2 [stat.ME] 7 Jan 2018

Laurentiu Hinoveanu∗1 , Fabrizio Leisen†1 , and Cristiano Villa‡1
1

School of Mathematics, Statistics and Actuarial Science,
University of Kent

Abstract
In this paper we present a loss-based approach to change point analysis. In particular, we look at the problem from two perspectives. The
first focuses on the definition of a prior when the number of change
points is known a priori. The second contribution aims to estimate
the number of change points by using a loss-based approach recently
introduced in the literature. The latter considers change point estimation as a model selection exercise. We show the performance of the
proposed approach on simulated data and real data sets.

Keywords: Change point; Discrete parameter space; Loss-based prior; Model
selection.

1

Introduction

There are several practical scenarios where it is inappropriate to assume
that the distribution of the observations does not change. For example, financial data sets can exhibit alternate behaviours due to crisis periods. In
this case it is sensible to assume changes in the underlying distribution. The
change in the distribution can be either in the value of one or more of the
parameters or, more in general, on the family of the distribution. In the
∗

lch36@kent.ac.uk
F.Leisen@kent.ac.uk
‡
C.Villa-88@kent.ac.uk
†

1

latter case, for example, one may deem appropriate to consider a normal
density for the stagnation periods, while a Student t, with relatively heavy
tails, may be more suitable to represent observations in the more turbulent
stages of a crisis. The task of identifying if, and when, one or more changes
have occurred is not trivial and requires appropriate methods to avoid detection of a large number of changes or, at the opposite extreme, seeing no
changes at all. The change point problem has been deeply studied from a
Bayesian point of view. Chernoff and Zacks (1964) focused on the change in
the means of normally distributed variables. Smith (1975) looked into the
single change point problem when different knowledge of the parameters of
the underlying distributions is available: all known, some of them known
or none of them known. Smith (1975) focuses on the binomial and normal
distributions. In Muliere and Scarsini (1985) the problem is tackled from
a Bayesian nonparametric perspective. The authors consider Dirichlet processes with independent base measures as underlying distributions. In this
framework, Petrone and Raftery (1997) have showed that the Dirichlet process prior could have a strong effect on the inference and may lead to wrong
conclusions in the case of a single change point. Raftery and Akman (1986)
have approached the single change point problem in the context of a Poisson
likelihood under both proper and improper priors for the model parameters.
Carlin et al. (1992) build on the work of Raftery and Akman (1986) by considering a two level hierarchical model. Both papers illustrate the respective
approaches by studying the well-known British coal-mining disaster data set.
In the context of multiple change points detection, Loschi and Cruz (2005)
have provided a fully Bayesian treatment for the product partitions model
of Barry and Hartigan (1992). Their application focused on stock exchange
data. Stephens (1994) has extended the Gibbs sampler introduced by Carlin et al. (1992) in the change point literature to handle multiple change
points. Hannart and Naveau (2009) have used Bayesian decision theory, in
particular 0-1 cost functions, to estimate multiple changes in homoskedastic
normally distributed observations. Schwaller and Robin (2017) extend the
product partition model of Barry and Hartigan (1992) by adding a graphical
structure which could capture the dependencies between multivariate observations. Fearnhead and Liu (2007) proposed a filtering algorithm for the
sequential multiple change points detection problem in the case of piecewise
regression models. Henderson and Matthews (1993) introduced a partial
Bayesian approach which involves the use of a profile likelihood, where the
aim is to detect multiple changes in the mean of Poisson distributions with an
application to haemolytic uraemic syndrome (HUS) data. The same data set
was studied by Tian et al. (2009), who proposed a method which treats the
change points as latent variables. Ko et al. (2015) have proposed an exten2

sion to the hidden Markov model of Chib (1998) by using a Dirichlet process
prior on each row of the regime matrix. Their model is semiparametric, as
the number of states is not specified in advance, but it grows according to
the data size. Heard and Turcotte (2017) have proposed a new sequential
Monte Carlo algorithm to infer multiple change points.
Whilst the literature covering change point analysis from a Bayesian perspective is vast when prior distributions are elicited, the documentation referring
to analysis under minimal prior information is limited, see Moreno et al.
(2005) and Girón et al. (2007). The former paper discusses the single change
point problem in a model selection setting, whilst the latter paper, which is
an extension of the former, tackles the multivariate change point problem
in the context of linear regression models. Our work aims to contribute to
the methodology for change point analysis under the assumption that the information about the number of change points and their location is minimal.
First, we discuss the definition of an objective prior for change point location, both for single and multiple changes, assuming the number of changes
is known a priori. Then, we define a prior on the number of change points via
a model selection approach. Here, we assume that the change point coincides
with one of the observations. As such, given X1 , X2 , . . . , Xn data points, the
change point location is discrete. To the best of our knowledge, the sole
general objective approach to define prior distributions on discrete spaces is
the one introduced by Villa and Walker (2015a).
To illustrate the idea, consider a probability distribution f (x|m), where
m ∈ M is a discrete parameter. Then, the prior π(m) is obtained by objectively measuring what is lost if the value m is removed from the parameter space, and it is the true value. According to Berk (1966), if a
model is misspecified, the posterior distribution asymptotically accumulates
on the model which is the most similar to the true one, where the similarity
is measured in terms of the Kullback–Leibler (KL) divergence. Therefore,
DKL (f (·|m)kf (·|m0 )), where m0 is the parameter characterising the nearest
model to f (x|m), represents the utility of keeping m. The objective prior is
then obtained by linking the aforementioned utility via the self-information
loss:


0
π(m) ∝ exp min
DKL (f (·|m)kf (·|m )) − 1,
(1)
0
m 6=m

where the Kullback–Leibler divergence (Kullback and Leibler, 1951) from the
sampling distribution with density f (x|m) to the one with density f (x|m0 )

3

is defined as:

Z
0

DKL (f (·|m)kf (·|m )) =




f (x|m)
f (x|m) · log
dx.
f (x|m0 )

X

Throughout the paper, the objective prior defined in equation (1) will be
referenced as the loss-based prior. This approach is used to define an objective prior distribution when the number of change points is known a priori.
To obtain a prior distribution for the number of change points, we adopt a
model selection approach based on the results in Villa and Walker (2015b),
where a method to define a prior on the space of models is proposed. To
illustrate, let us consider k Bayesian models:
Mj = {fj (x|θj ), πj (θj )}

j ∈ {1, 2, . . . , k},

(2)

where fj (x|θj ) is the sampling density characterised by θj and πj (θj ) represents the prior on the model parameter.
Assuming the prior on the model parameter, πj (θj ), is proper, the model
prior probability Pr(Mj ) is proportional to the expected minimum Kullback–
Leibler divergence from Mj , where the expectation is considered with respect
to πj (θj ). That is:



Pr(Mj ) ∝ exp Eπj inf DKL (fj (x|θj )kfi (x|θi ))
j = 1, . . . , k. (3)
θi ,i6=j

The model prior probabilities defined in equation (3) can be employed to
derive the model posterior probabilities through:
#−1
" k
X Pr(Mj )
Bji
,
(4)
Pr(Mi |x) =
Pr(M
)
i
j=1
where Bji is the Bayes factor between model Mj and model Mi , defined as
R
fj (x|θj )πj (θj ) dθj
Bji = R
,
fi (x|θi )πi (θi ) dθi
with i 6= j ∈ {1, 2, . . . , k}.
This paper is structured as follows: in Section 2 we establish the way we set
objective priors on both single and multiple change point locations. Section
3 shows how we define the model prior probabilities for the number of change
point locations. Illustrations of the model selection exercise are provided in
Sections 4 and 5, where we work with simulated and real data, respectively.
Section 6 is dedicated to final remarks.
4

2

Objective Prior on the Change Point Locations

This section is devoted to the derivation of the loss-based prior when the
number of change points is known a priori. Specifically, let k be the number
of change points and m1 < m2 < . . . < mk their locations. We introduce the
idea in the simple case where we assume that there is only one change point
in the data set (see Section 2.1). Then, we extend the results to the more
general case where multiple change points are assumed (see Section 2.2).
A well-known objective prior for finite parameter spaces, in cases where there
is no structure, is the uniform prior (Berger et al., 2012). As such, a natural
choice for the prior on the change points location is the uniform (Koop and
Potter, 2009). The corresponding loss-based prior is indeed the uniform, as
shown below, which is a reassuring result as the objective prior for a specific
parameter space, if exists, should be unique.

2.1

Single Change Point

As mentioned above, we show that the loss-based prior for the single change
point case coincides with the discrete uniform distribution over the set {1, 2, . . . , n−
1}.
Let X(n) = (X1 , . . . , Xn ) denote an n-dimensional vector of random variables,
representing the random sample, and m be our single change point location,
that is m ∈ {1, 2, . . . , n − 1}, such that
i.i.d.

X1 , . . . , Xm |θ̃1 ∼ f1 (·|θ̃1 )
i.i.d.

Xm+1 , . . . , Xn |θ̃2 ∼ f2 (·|θ̃2 ).

(5)

Note that we assume that there is a change point in the series, as such the
space of m does not include the case m = n. In addition, we assume that
θ̃1 6= θ̃2 when f1 = f2 . The sampling density for the vector of observations
x(n) = (x1 , . . . , xn ) is:
f (x

(n)

|m, θ̃1 , θ̃2 ) =

m
Y

f1 (xi |θ̃1 )

i=1

n
Y

f2 (xi |θ̃2 ).

(6)

i=m+1

Let m0 6= m. Then, the Kullback–Leibler divergence between the model

5

parametrised by m and the one parametrised by m0 is:
Z
(n)
(n)
0
DKL (f (x |m, θ̃1 , θ̃2 )kf (x |m , θ̃1 , θ̃2 )) = f (x(n) |m, θ̃1 , θ̃2 )
log

f (x(n) |m, θ̃1 , θ̃2 )
f (x(n) |m0 , θ̃1 , θ̃2 )

!
dx(n) .
(7)

Without loss of generality, consider m < m0 . In this case, note that
0

m
Y
f (x(n) |m, θ̃1 , θ̃2 )
f2 (xi |θ̃2 )
=
,
f (x(n) |m0 , θ̃1 , θ̃2 ) i=m+1 f1 (xi |θ̃1 )

leading to
DKL (f (x(n) |m, θ̃1 , θ̃2 )kf (x(n) |m0 , θ̃1 , θ̃2 )) =
m0 Z
X
f2 (xi |θ̃2 ) log
i=m+1

f2 (xi |θ̃2 )
f1 (xi |θ̃1 )

!
dxi .

(8)

On the right hand side of equation (8), we can recognise the Kullback–Leibler
divergence from density f2 to density f1 , thus getting:
DKL (f (x(n) |m, θ̃1 , θ̃2 )||f (x(n) |m0 , θ̃1 , θ̃2 )) =
(m0 − m)DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )).

(9)

In a similar fashion, when m > m0 , we have that:
DKL (f (x(n) |m, θ̃1 , θ̃2 )kf (x(n) |m0 , θ̃1 , θ̃2 )) =
(m − m0 )DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 )).

(10)

In this single change point scenario, we can consider m0 as a perturbation
of the change point location m, that is m0 = m ± l where l ∈ N∗ , such
that 1 ≤ m0 < n. Then, taking into account equations (9) and (10), the
Kullback–Leibler divergence becomes:

DKL (f (x(n) |m, θ̃1 , θ̃2 )kf (x(n) |m0 , θ̃1 , θ̃2 )) =

 l · DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )),


l · DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 )),

6

if m < m0
if m > m0 ,

and
h
i
(n)
(n)
0
min
DKL (f (x |m, θ̃1 , θ̃2 )kf (x |m , θ̃1 , θ̃2 )) =
0

m 6=m

= min
{l · DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )), l · DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))}
0
m 6=m

= min
{DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )), DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))} · min
{l} .
m0 6=m
m0 6=m
| {z }
1

(11)
We observe that equation (11) is only a function of θ̃1 and θ̃2 and does not
depend on m. Thus, π(m) ∝ 1 and, therefore,
π(m) =

1
n−1

m ∈ {1, . . . , n − 1}.

(12)

This prior was used, for instance, in an econometric context by Koop and
Potter (2009) with the rationale of giving equal weight to every possible
change point location.

2.2

Multivariate Change Point Problem

In this section, we address the change point problem in its generality by
assuming that there are 1 ≤ k < n change points. In particular, for the data
x(n) = (x1 , . . . , xn ), we consider the following sampling distribution
f (x

(n)

|m, θ̃) =

m1
Y

f1 (xi |θ̃1 )

i=1

j+1
k−1
Y mY

n
Y

fj+1 (xi |θ̃j+1 )

j=1 i=mj +1

fk+1 (xi |θ̃k+1 ),

i=mk +1

(13)
where m = (m1 , . . . , mk ), 1 ≤ m1 < m2 < . . . < mk < n, is the vector
of the change point locations and θ̃ = (θ̃1 , . . . , θ̃k , θ̃k+1 ) is the vector of the
parameters of the underlying probability distributions. Schematically:
X1
Xm1 +1
..
.

, . . . , Xm1 |θ̃1
, . . . , Xm2 |θ̃2
.
, . . . , ..

i.i.d.

∼ f1 (·|θ̃1 )
∼ f2 (·|θ̃2 )
.. ..
.....

i.i.d.

i.i.d.

Xmk−1 +1 , . . . , Xmk |θ̃k ∼ fk (·|θ̃k )
i.i.d.
Xmk +1
, . . . , Xn |θ̃k+1 ∼ fk+1 (·|θ̃k+1 ).
If f1 = f2 = · · · = fk+1 , then it is reasonable to assume that some of the θ’s
are different. Without loss of generality, we assume that θ̃1 6= θ̃2 6= · · · 6=
7

θ̃k 6= θ̃k+1 . In a similar fashion to the single change point case, we cannot
assume mk = n since we require exactly k change points.
In this case, due to the multivariate nature of the vector m = (m1 , . . . , mk ),
the derivation of the loss-based prior is not as straightforward as in the one
dimensional case. In fact, the derivation of the prior is based on heuristic
considerations supported by the below Theorem 1 (the proof of which is in
the Appendix). In particular, we are able to prove an analogous of equations
(9) and (10) when only one component is arbitrarily perturbed. Let us define
the following functions:
d+1
j (θ̃) = DKL (fj+1 (·|θ̃j+1 )kfj (·|θ̃j ))
d−1
j (θ̃) = DKL (fj (·|θ̃j )kfj+1 (·|θ̃j+1 )),
where j ∈ {1, 2, . . . , k}. The following Theorem is useful to understand the
behaviour of the loss-based prior in the general case.
Theorem 1. Let f (x(n) |m, θ̃) be the sampling distribution defined in equation (13) and consider j ∈ {1, . . . , k}. Let m0 be such that m0i = mi for i 6= j,
and let the component m0j be such that m0j 6= mj and mj−1 < m0j < mj+1 .
Therefore,
DKL (f (x(n) |m, θ̃)kf (x(n) |m0 , θ̃) = |m0j − mj |dSj (θ̃),
where S = sgn(m0j − mj ).
Note that, Theorem 1 states that the minimum Kullback–Leibler divergence
is achieved when m0j = mj + 1 or m0j = mj − 1. This result is not surprising
since the Kullback–Leibler divergence measures the degree of similarity between two distributions. The smaller the perturbation caused by changes in
one of the parameters is, the smaller the Kullback–Leibler divergence between
the two distributions is. Although Theorem 1 makes a partial statement
about the multiple change points scenario, it provides a strong argument for
supporting the uniform prior. Indeed, if now we consider the general case of
having k change points, it is straightforward to see that the Kullback–Leibler
divergence is minimised when only one of the components of the vector m is
perturbed by (plus or minus) one unit. As such, the loss-based prior depends
on the vector of parameters θ̃ only, as in the one-dimensional case, yielding
the uniform prior for m.
Therefore, the loss-based prior on the multivariate change point location is

−1
n−1
π(m) =
,
(14)
k
8

where m = (m1 , . . . , mk ), 1 ≤ m1 < m2 < . . . < mk < n. The denominator
in equation (14) has the above form because, for every number of k change
points, we are interested
 in the number of k-subsets from a set of n − 1
n−1
elements, which is k . The same prior was also derived in a different way
by Girón et al. (2007).

3

Loss-based Prior on the Number of Change
Points

Here, we approach the change point analysis as a model selection problem.
In particular, we define a prior on the space of models, where each model
represents a certain number of change points (including the case of no change
points). The method adopted to define the prior on the space of models is the
one introduced in Villa and Walker (2015b). We proceed as follows. Assume
1

1
M1

n

𝜃̃1

M0
𝜃̃1

m1

n

𝜃̃2

1

1
M2

θ1

𝜃̃1

m1

1
M3

𝜃̃1

𝜃̃2

m2

n

𝜃̃3

2

θ1

m1

θ1

𝜃̃2

m2

𝜃̃3

m3

n

𝜃̃4

⋮
1
Mk

𝜃̃1

m1

𝜃̃2

m2

𝜃̃3

m3

𝜃̃4

⋯ ⋯ ⋯

𝜃̃𝑘

mk

n
𝜃̃𝑘+1

Figure 1: Diagram showing the way we specify our models. The arrows indicate that the
respective change point locations remain fixed from the previous model to the current one.

we have to select from k + 1 possible models. Let M0 be the model with no
change points, M1 the model with one change point and so on. Generalising,
model Mk corresponds to the model with k change points. The idea is that
the current model encompasses the change point locations of the previous
model. As an example, in model M3 the first two change point locations will
be the same as in the case of model M2 . To illustrate the way we envision our
models, we have provided Figure 1. It has to be noted that the construction
of the possible models from M0 to Mk can be done in a different way to
one here described. Obviously, the approach to define the model priors stays

9

unchanged. Consistently with the notation used in Section 1,

θ̃1 , . . . , θ̃k+1 , m1 , . . . , mk if k = 1, . . . , n − 1
θk =
θ̃1
if k = 0
represents the vector of parameters of model Mk , where θ̃1 , . . . , θ̃k+1 are the
model specific parameters and m1 , . . . , mk are the change point locations, as
in Figure 1.
Based on the way we have specified our models, which are in direct correspondence with the number of change points and their locations, we state
Theorem 2 (the proof of which is in the Appendix).
Theorem 2. Let
DKL (Mi kMj ) = DKL (f (x(n) |θi )kf (x(n) |θj )).
For any 0 ≤ i < j ≤ k integers, with k < n, and the convention mj+1 = n,
we have the following:
DKL (Mi kMj ) =

j
h
X

i
(mq+1 − mq ) · DKL (fi+1 (·|θ̃i+1 )kfq+1 (·|θ̃q+1 )) ,

q=i+1

and
DKL (Mj kMi ) =

j
h
X

i
(mq+1 − mq ) · DKL (fq+1 (·|θ̃q+1 )kfi+1 (·|θ̃i+1 )) .

q=i+1

The result in Theorem 2 is useful when the model selection exercise is implemented. Indeed, the Villa and Walker (2015b) approach requires the computation of the Kullback–Leibler divergences in Theorem 2. Recalling equation
(3), the objective model prior probabilities are then given by:



Pr(Mj ) ∝ exp Eπj inf DKL (Mj kMi )
j = 0, 1, . . . , k.
(15)
θi ,i6=j

For illustrative purposes, in the Appendix we derive the model prior probabilities to perform model selection among M0 , M1 and M2 .
It is easy to infer from equation (15) that model priors depend on the prior
distribution assigned to the model parameters, that is on the level of uncertainty that we have about their true values. For the change point location,
10

a sensible choice is the uniform prior which, as shown in Section 2, corresponds to the loss-based prior. For the model specific parameters, we have
several options. If one wishes to pursue an objective analysis, intrinsic priors
(Berger and Pericchi, 1996) may represent a viable solution since they are
proper. Nonetheles, the method introduce by Villa and Walker (2015b) does
not require, in principle, an objective choice as long as the priors are proper.
Given that we use the latter approach, here we consider subjective priors for
the model specific parameters.
Remark. In the case where the changes in the underlying sampling distribution are limited to the parameter values, the model prior probabilities
defined in (15) follow the uniform distribution. That is, Pr(Mj ) ∝ 1. In the
real data example illustrated in Section 5.1, we indeed consider a problem
where the above case occurs.

3.1

A special case: selection between M0 and M1

Let us consider the case where we have to estimate whether there is or not
a change point in a set of observations. This implies that we have to choose
between model M0 (i.e. no change point) and M1 (i.e. one change point).
Following our approach, we have:



Pr(M0 ) ∝ exp Eπ0 inf DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 )) ,

(16)

θ̃2

and



Pr(M1 ) ∝ exp Eπ1 (n − m1 ) · inf DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )) .

(17)

θ̃1

Now, let us assume independence between the prior on the change point
location and the prior on the parameters of the underlying sampling distributions, that is π1 (m1 , θ̃1 , θ̃2 ) = π1 (m1 )π1 (θ̃1 , θ̃2 ). Let us further recall that,
as per equation (14), π1 (m1 ) = 1/(n−1). As such, we observe that the model
prior probability on M1 becomes:


 
n
Eπ1 (θ̃1 ,θ̃2 ) inf DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )) .
(18)
Pr(M1 ) ∝ exp
2
θ̃1
We notice that the model prior probability for model M1 is increasing when
the sample size increases. This behaviour occurs whether there is or not a
11

change point in the data. We propose to address the above problem by using
a non-uniform prior for m1 . A reasonable alternative, which works quite well
in practice, would be the following shifted binomial as prior:


m −1  n−m1 −1
n−2
n−1 1
1
π1 (m1 ) =
, 1 ≤ m1 ≤ n − 1.
(19)
m1 − 1
n
n
To argument the choice of (19), we note that, as n increases, the probability
mass will be more and more concentrated towards the upper end of the
support. Therefore, from equations (17) and (19) follows:




2n − 2
Pr(M1 ) ∝ exp
Eπ1 (θ̃1 ,θ̃2 ) inf DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )) . (20)
n
θ̃1
For the more general case where we consider more than two models, the
problem highlighted in equation (18) vanishes.

4

Change Point Analysis on Simulated Data

In this section, we present the results of several simulation studies based
on the methodologies discussed in Sections 2 and 3. We start with a scenario involving discrete distributions in the context of the one change point
problem. We then show the results obtained when we consider continuous
distributions for the case of two change points. The choice of the underlying
sampling distributions is in line with Villa and Walker (2015b).

4.1

Single sample

Scenario 1. The first scenario concerns the choice between models M0 and
M1 . Specifically, for M0 we have:
i.i.d.

X1 , X2 , . . . , Xn |p ∼ Geometric(p)
and for M1 we have:
i.i.d.

X1 , X2 , . . . , Xm1 |p ∼ Geometric(p)
i.i.d.

Xm1 +1 , Xm1 +2 , . . . , Xn |λ ∼ Poisson(λ)
Let us denote with f1 (·|p) and f2 (·|λ) the probability mass functions of the
Geometric and the Poisson distributions, respectively. The priors for the
parameters of f1 and f2 are p ∼ Beta(a, b) and λ ∼ Gamma(c, d).
12

In the first simulation, we sample n = 100 observations from model M0 with
p = 0.8. To perform the change point analysis, we have chosen the following
parameters for the priors on p and λ: a = 2, b = 2, c = 3 and d = 1.
Applying the approach introduced in Section 3, we obtain Pr(M0 ) ∝ 1.59
and Pr(M1 ) ∝ 1.81. These model priors yield the posterior distribution
probabilities (refer to equation (4)) Pr(M0 |x(n) ) = 0.92 and Pr(M1 |x(n) ) =
0.08. As expected, the selection process strongly indicates the true model
as M0 . Table 1 reports the above probabilities including other information,
such as the appropriate Bayes factors.
The second simulation looked at the opposite setup, that is we sample n =
100 observations from M1 , with p = 0.8 and λ = 3. We have sampled 50 data
points from the Geometric distribution and the remaining 50 data points from
the Poisson distribution. In Figure 2, we have plotted the simulated sample,
where it is legitimate to assume a change in the underlying distribution.
Using the same prior parameters as above, we obtain Pr(M0 |x(n) ) = 0.06 and
Pr(M1 |x(n) ) = 0.94. Again, the model selection process is assigning heavy
posterior mass to the true model M1 . These results are further detailed in
Table 1.
8

6

4

2

0
1

20

40

60

80

100

Figure 2: Scatter plot of the data simulated from model M1 in Scenario 1.

Scenario 2. In this scenario we consider the case where we have to select
among three models, that is model M0 :
i.i.d.

X1 , X2 , . . . , Xn |λ, κ ∼ Weibull(λ, κ),

(21)

model M1 :
i.i.d.

X1 , X2 , . . . , Xm1 |λ, κ ∼ Weibull(λ, κ)
i.i.d.

Xm1 +1 , Xm1 +2 , . . . , Xn |µ, τ ∼ Log-normal(µ, τ ),
13

(22)

Pr(M0 )
Pr(M1 )
B01
B10
Pr(M0 |x(n) )
Pr(M1 |x(n) )

True model
M0
M1
0.47 0.47
0.53 0.53
12.39 0.08
0.08 12.80
0.92 0.06
0.08 0.94

Table 1: Model prior, Bayes factor and model posterior probabilities for the change point
analysis in Scenario 1. We considered samples from, respectively, model M0 and model
M1 .

with 1 ≤ m1 ≤ n−1 being the location of the single change point, and model
M2 :
i.i.d.

X1 , X2 , . . . , Xm1 |λ, κ ∼ Weibull(λ, κ)
i.i.d.

Xm1 +1 , Xm1 +2 , . . . , Xm2 |µ, τ ∼ Log-normal(µ, τ )
i.i.d.

Xm2 +1 , Xm2 +2 , . . . , Xn |α, β ∼ Gamma(α, β),

(23)

with 1 ≤ m1 < m2 ≤ n − 1 representing the locations of the two change
points, such that m1 corresponds exactly to the same location as in model
M1 . Analogously to the previous scenario, we sample from each model in
turn and perform the selection to detect the number of change points.
Let f1 (·|λ, κ), f2 (·|µ, τ ) and f3 (·|α, β) represent the Weibull, Log-normal and
Gamma densities, respectively, with θ̃1 = (λ, κ), θ̃2 = (µ, τ ) and θ̃3 = (α, β).
We assume a Normal prior on µ and Gamma priors on all the other parameters as follows:
λ ∼ Gamma(1.5, 1)

κ ∼ Gamma(5, 1)

τ ∼ Gamma(16, 1)

α ∼ Gamma(10, 1)

µ ∼ Normal(0.05, 1),
β ∼ Gamma(0.2, 0.1).

In the first exercise, we have simulated n = 100 observations from model
M0 , where we have set λ = 1.5 and κ = 5. We obtain the following model
priors: Pr(M0 ) ∝ 1.09, Pr(M1 ) ∝ 1.60 and Pr(M2 ) ∝ 1.37, yielding the
posteriors Pr(M0 |x(n) ) = 0.96, Pr(M1 |x(n) ) = 0.04 and Pr(M2 |x(n) ) = 0.00.
We then see that the approach assigns high mass to the true model M0 .
Table 2 reports the above probabilities and the corresponding Bayes factors.
The second simulation was performed by sampling 50 observations from a
Weibull with parameter values as in the previous exercise, and the remaining
14

Pr(M0 )
Pr(M1 )
Pr(M2 )
B01
B02
B12
Pr(M0 |x(n) )
Pr(M1 |x(n) )
Pr(M2 |x(n) )

M0
0.27
0.39
0.34
36.55
1.84 × 103
50.44
0.96
0.04
0.00

True model
M1
M2
0.27
0.27
0.39
0.39
0.34
0.34
3.24 × 10−4 4.65 × 10−40
0.02
1.27 × 10−45
55
2.72 × 10−6
0.00
0.00
0.98
0.00
0.02
1.00

Table 2: Model prior, Bayes factor and model posterior probabilities for the change point
analysis in Scenario 2. We considered samples from, respectively, model M0 , model M1
and model M2 .

50 observations from a Log-normal density with location parameter µ = 0.05
and scale parameter τ = 16. The data is displayed in Figure 3.

2.5
2
1.5
1
0.5
0
1

20

40

60

80

100

Figure 3: Scatter plot of the observations simulated from model M1 in Scenario 2.

The model posterior probabilities are Pr(M0 |x(n) ) = 0.00, Pr(M1 |x(n) ) = 0.98
and Pr(M2 |x(n) ) = 0.02, which are reported in Table 2. In this case as well,
we see that the model selection procedure indicates M1 as the true model,
as expected.

15

12
10
8
6
4
2
0
1

20

40

60

80

100

Figure 4: Scatter plot of the observations simulated from model M2 in Scenario 2.

Finally, for the third simulation exercise we sample 50 and 20 data points
from, respectively, a Weibull and a Log-normal with parameter values as
defined above, and the last 30 observations are sampled from a Gamma
distribution with parameters α = 10 and β = 2. From Table 2, we note that
the posterior distribution on the model space accumulates on the true model
M2 .

4.2

Frequentist Analysis

In this section, we perform a frequentist analysis of the performance of the
proposed prior by drawing repeated samples from different scenarios. In
particular, we look at a two change points problem where the sampling distributions are Student-t with different degrees of freedom. In this scenario,
we perform the analysis with 60 repeated samples generated by different
densities with the same mean values.
Then, we repeat the analysis of Scenario 2 by selecting 100 samples for
n = 500 and n = 1500. We consider different sampling distributions with
the same mean and variance. In this scenario, where we added the further
constraint of the equal variance, it is interesting to note that the change in
distribution is captured when we increase the sample size, meaning that we
learn more about the true sampling distributions.
We also compare the performances of the loss-based prior with the uniform
prior when we analyse the scenario with different sampling distributions.
Namely, Weibull/Log-normal/Gamma. It is interesting to note that the uniform prior is unable to capture the change in distribution even for a large
sample size. On the contrary, the loss-based prior is able to detect the number of change points when n = 1500. Furthermore, for n = 500, even though
16

both priors are not able to detect the change points most of the times, the
loss-based prior has a higher frequency of success when compared to the
uniform prior.
Scenario 3. In this scenario, we consider the case where the sampling
distributions belong to the same family, that is Student-t, where the true
model has two change points. In particular, let f1 (·|ν1 ), f2 (·|ν2 ) and f3 (·|ν3 )
represent the densities of three standard t distributions, respectively. We
assume that ν1 , ν2 and ν3 are positive integers strictly greater than one so
to have defined mean for each density. Note that this allows us to compare
distributions of the same family with equal mean. The priors assigned to the
number of degrees of freedom assume a parameter space of positive integers
strictly larger than 1. As such, we define them as follows:
ν1 ∼ 2 + Poisson(30)

ν2 ∼ 2 + Poisson(3)

ν3 ∼ 2 + Poisson(8).

In this experiment, we consider 60 repeated samples, each of size n = 300
and with the following structure:
• X1 , . . . , X100 from a Student-t distribution with ν1 = 30,
• X101 , . . . , X200 from a Student-t distribution with ν2 = 3,
• X201 , . . . , X300 from a Student-t distribution with ν3 = 8.
Table 3 reports the frequentist results of the simulation study. First, note
that P (M1 ) = P (M2 ) = P (M3 ) = 1/3 as per the Remark in Section 3. For
all the simulated samples, the loss-based prior yields a posterior with the
highest probability assigned to the true model M2 . We also note that the
above posterior is on average 0.75 with a variance 0.02, making the inferential
procedure extremely accurate.
Pr(M0 |x(n) )
Pr(M1 |x(n) )
Pr(M2 |x(n) )

Mean posterior
0.01
0.24
0.75

Variance posterior
3.84 × 10−4
0.0160
0.0190

Freq. true model
0/60
0/60
60/60

Table 3: Average model posterior probabilities, variance and frequency of true model for
the Scenario 3 simulation exercise.

Scenario 4. In this scenario, we perform repeated sampling from the setup
described in scenario 2 above, where the true model has two change points. In
17

particular, we draw 100 samples with n = 500 and n = 1500. For n = 500, the
loss-based prior probabilities are P (M0 ) = 0.18, P (M1 ) = 0.16 and P (M2 ) =
0.66. For n = 1500, the loss-based prior probabilities are P (M0 ) = 0.015,
P (M1 ) = 0.014 and P (M2 ) = 0.971. The simulation results are reported,
respectively, in Table 4 and in Table 5. The two change point locations
for n = 500 are at the 171st and 341st observations. For n = 1500, the
first change point is the 501st observation, while the second is at the 1001st
observation. We note that there is a sensible improvement in detecting the
true model, using the loss-based prior, when the sample size increases. In
particular, we move from 30% to 96%.

(n)

Pr(M0 |x )
Pr(M1 |x(n) )
Pr(M2 |x(n) )

Mean posterior
9.88 × 10−4
0.63
0.37

Variance posterior
2.60 × 10−5
0.0749
0.0745

Freq. true model
0/100
70/100
30/100

Table 4: Average model posterior probabilities, variance and frequency of true model for
the Scenario 4 simulation exercise with n = 500 and the loss-based prior.

(n)

Pr(M0 |x )
Pr(M1 |x(n) )
Pr(M2 |x(n) )

Mean posterior
1.33 × 10−13
0.08
0.92

Variance posterior
1.76 × 10−24
0.0200
0.0200

Freq. true model
0/100
4/100
96/100

Table 5: Average model posterior probabilities, variance and frequency of true model for
the Scenario 4 simulation exercise with n = 1500 and the loss-based prior.

To compare the loss-based prior with the uniform prior we have run the
simulation on the same data samples used above. The results for n = 500
and n = 1500 are in Table 6 and in Table 7, respectively. Although we can
observe an improvement when the sample size increases, the uniform prior
does not lead to a clear detection of the true model for both sample sizes.

Pr(M0 |x(n) )
Pr(M1 |x(n) )
Pr(M2 |x(n) )

Mean posterior
16 × 10−4
0.82
0.18

Variance posterior
7.15 × 10−5
0.0447
0.0443

Freq. true model
0/100
91/100
9/100

Table 6: Average model posterior probabilities, variance and frequency of true model for
the Scenario 4 simulation exercise with n = 500 and the uniform prior.

18

Mean posterior
8.64 × 10−12
0.501
0.499

(n)

Pr(M0 |x )
Pr(M1 |x(n) )
Pr(M2 |x(n) )

Variance posterior
7.45 × 10−21
0.1356
0.1356

Freq. true model
0/100
49/100
51/100

Table 7: Average model posterior probabilities, variance and frequency of true model for
the Scenario 4 simulation exercise with n = 1500 and the uniform prior.

0.3
Weibull
Log-normal
Gamma

0.25

Pdfs

0.2

0.15

0.1

0.05

0

0

1

2

3

4

5

6

7

8

9

10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25

Figure 5: The densities of Weibull(λ, κ), Log-normal(µ, τ ) and Gamma(α, β) with the
same mean (equal to 5) and the same variance (equal to 2.5).

Finally, we conclude this section with a remark. One may wonder why the
change point detection requires an increasing in the sample size, and the
reply can be inferred from Figure 5, which displays the density functions
of the distributions employed in this scenario. As it can be observed, the
densities are quite similar, which is not surprising since these distributions
have the same means and the same variances. The above similarity can
also be appreciated in terms of Hellinger distance, see Table 8. In other
words, from Figure 5 we can see that the main differences in the underlying
distributions are in the tail areas. It is therefore necessary to have a relatively
large number of observations in order to be able to discern differences in the
19

densities, because in this case only we would have a sufficient representation
of the whole distribution.

Weibull(λ, κ)
Log-normal(µ, τ )

Hellinger distances
Weibull(λ, κ) Log-normal(µ, τ ) Gamma(α, β)
0.1411996
0.09718282
0.04899711

Table 8: Hellinger distances between all the pairs formed from a Weibull(λ, κ), Lognormal(µ, τ ) and Gamma(α, β). The six hyperparameters are such that the distributions
have the same mean=5 and same variance=2.5.

5

Change Point Analysis on Real Data

In this section, we illustrate the proposed approach applied to real data. We
first consider a well known data set which has been extensively studied in
the literature of the change point analysis, that is the British coal-mining
disaster data (Carlin et al., 1992). The second set of data we consider refers
to the daily returns of the S&P 500 index observed over a period of four years.
The former data set will be investigated in Section 5.1, while the latter in
Section 5.2.

5.1

British Coal-Mining Disaster Data

The British coal-mining disaster data consists of the yearly number of deaths
for the British coal miners over the period 1851-1962. It is believed that the
change in the working conditions, and in particular, the enhancement of the
security measures, led to a decrease in the number of deaths. This calls for
a model which can take into account a change in the underlying distribution
around a certain observed year. With the proposed methodology we wish
to detect if the assumption is appropriate. In particular, if a model with
one change point is more suitable to represent the data than a model where
no changes in the sampling distribution are assumed. Figure 6 shows the
number of deaths per year in the British coal-mining industry from 1851 to
1962. As in Chib (1998), we assume a Poisson sampling distribution with a
possible change in the parameter value. That is
i.i.d.

X1 , X2 , . . . , Xm |φ1 ∼ Poisson(φ1 )
i.i.d.

Xm+1 , Xm+2 , . . . , Xn |φ2 ∼ Poisson(φ2 ),
20

(24)

7
6
5
4
3
2
1
0
1851

1871

1891

1911

1931

1951

Year

Figure 6: Scatter plot of the British coal-mining disaster data.

where m is the unknown location of the single change point, such that 1 ≤
m ≤ n, and a Gamma(2, 1) is assumed for φ1 and φ2 . The case m = n
corresponds to the scenario with no change point, that is model M0 . The
case m < n assumes one change point, that is model M1 .
Let f1 (·|φ1 ) and f2 (·|φ2 ) be the Poisson distributions with parameters φ1 and
φ2 , respectively. Then, the analysis is performed by selecting between model
M0 , that is when the sampling distribution is f1 , and model M1 , where the
sampling distribution is f1 up to a certain m < n and f2 from m + 1 to n.
As highlighted in the Remark at the end of Section 3, the prior on the model
space is the discrete uniform distribution, that is Pr(M0 ) = Pr(M1 ) = 0.5.
The proposed model selection approach leads to the Bayes factors B01 =
1.61 × 10−13 and B10 = 6.20 × 1012 , where it is obvious that the odds are
strongly in favour of model M1 . Indeed, we have Pr(M1 |x(n) ) ≈ 1.

5.2

Daily S&P 500 Absolute Log-Return Data

The second real data analysis aims to detect change points in the absolute
value of the daily logarithmic returns of the S&P500 index observed from
the 14/01/2008 to the 31/12/2011 (see Figure 7). As underlying sampling
distributions we consider the Weibull and the Log-normal (Yu, 2001), and
the models among which we select are as follows. M0 is a Weibull(λ, κ), M1 is
formed by a Weibull(λ, κ) and a Log-normal(µ1 , τ1 ) and, finally, M2 is formed
by a Weibull(λ, κ), a Log-normal(µ1 , τ1 ) and a Log-normal(µ2 , τ2 ). An interesting particularity of this problem is that we will consider a scenario where
the changes are in the underlying distribution or in the parameter values of
the same distribution. As suggested in Section 4.1.3 of Kass and Raftery
(1995), due to the large sample size of the data set, we could approximate
21

the Bayes factor by using the Schwarz criterion. Therefore, in this case the
specification of the priors for the parameters of the underlying distributions
is not necessary. From the results in Table 9, we see that the model indicated
0.12
0.1
0.08
0.06
0.04
0.02
0
14/01/08

29/10/08

14/08/09

31/05/10

16/03/11

30/12/11

Figure 7: Absolute daily log-returns of the S&P500 index from 14/01/08 to 30/12/11.

by the proposed approach is M2 . In other words, there is very strong indication that there are two change points in the data set. From Table 9, we note
Pr(M0 )
Pr(M1 )
Pr(M2 )
B01
B02
B12
Pr(M0 |x(n) )
Pr(M1 |x(n) )
Pr(M2 |x(n) )

0.36
0.32
0.32
7.72 × 1018
3.30 × 10−3
4.28 × 10−22
0.00
0.00
1.00

Table 9: Model prior, Bayes factor and model posterior probabilities for the S&P500
change point analysis.

that the prior on model M1 and M2 assigned by the proposed method are the
same. This is not surprising as the only difference between the two models
is an additional Log-normal distribution with different parameter values.

6

Conclusion

Bayesian inference in change point problems under the assumption of minimal
prior information has not been deeply explored in the past, as the limited
literature on the matter shows.
22

We contribute to the area by deriving an objective prior distribution to detect change point locations, when the number of change points is known a
priori. As a change point location can be interpreted as a discrete parameter,
we apply recent results in the literature (Villa and Walker, 2015a) to make
inference. The resulting prior distribution, which is the discrete uniform distribution, it is not new in the literature (Girón et al., 2007), and therefore
can be considered as a validation of the proposed approach.
A second major contribution is in defining an objective prior on the number
of change points, which has been approached by considering the problem
as a model selection exercise. The results of the proposed method on both
simulated and real data, show the strength of the approach in estimating the
number of change points in a series of observations. A point to note is the
generality of the scenarios considered. Indeed, we consider situations where
the change is in the value of the parameter(s) of the underlying sampling
distribution, or in the distribution itself. Of particular interest is the last
real data analysis (S&P 500 index), where we consider a scenario where we
have both types of changes, that is the distribution for the first change point
and on the parameters of the distribution for the second.
The aim of this work was to set up a novel approach to address change point
problems. In particular, we have selected prior densities for the parameters
of the models to reflect a scenario of equal knowledge, in the sense that
model priors are close to represent a uniform distribution. Two remarks are
necessary here. First, in the case prior information about the true value of
the parameters is available, and one wishes to exploit it, the prior densities
will need to reflect it and, obviously, the model prior will be impacted by
the choice. Second, in applications it is recommended that some sensitivity
analysis is performed, so to investigate if and how the choice of the parameter
densities affects the selection process.

Acknowledgements
Fabrizio Leisen was supported by the European Community’s Seventh Framework Programme [FP7/2007-2013] under grant agreement no: 630677. Cristiano Villa was supported by the Royal Society Research Grant no: RG150786.

23

References
Barry, D., Hartigan, J. A., 1992. Product Partition Models for Change Point
Problems. The Annals of Statistics 20 (1), 260–279.
Berger, J. O., Bernardo, J. M., Sun, D., 2012. Objective Priors for Discrete Parameter Spaces. Journal of the American Statistical Association
107 (498), 636–648.
Berger, J. O., Pericchi, L. R., 1996. The Intrinsic Bayes Factor for Model
Selection and Prediction. Journal of the American Statistical Association
91 (433), 109–122.
Berk, R. H., 1966. Limiting Behavior of Posterior Distributions when the
Model is Incorrect. The Annals of Mathematical Statistics 37 (1), 51–58.
Carlin, B. P., Gelfand, A. E., Smith, A. F. M., 1992. Hierarchical Bayesian
Analysis of Changepoint Problems. Journal of the Royal Statistical Society.
Series C (Applied Statistics) 41 (2), 389–405.
Chernoff, H., Zacks, S., 1964. Estimating the Current Mean of a Normal
Distribution which is Subjected to Changes in Time. The Annals of Mathematical Statistics 35 (3), 999–1018.
Chib, S., 1998. Estimation and Comparison of Multiple Change-point Models. Journal of Econometrics 86 (2), 221–241.
Fearnhead, P., Liu, Z., 2007. On-line Inference for Multiple Changepoint
Problems. Journal of the Royal Statistical Society: Series B (Statistical
Methodology) 69 (4), 589–605.
Girón, F. J., Moreno, E., Casella, G., 2007. Objective Bayesian Analysis of
Multiple Changepoints for Linear Models (with discussion). In: Bernardo,
J.M., Bayarri, M.J., Berger, J.O., Dawid, A.P., Heckerman, D., Smith,
A.F.M., West, M. (Eds.), Bayesian Statistics 8. Oxford University Press,
London, pp. 227–252.
Hannart, A., Naveau, P., 2009. Bayesian Multiple Change Points and Segmentation: Application to Homogenization of Climatic Series. Water Resources Research 45 (10), W10444.
Heard, N. A., Turcotte, M. J. M., 2017. Adaptive Sequential Monte Carlo for
Multiple Changepoint Analysis. Journal of Computational and Graphical
Statistics 26 (2), 414–423.
24

Henderson, R., Matthews, J. N. S., 1993. An Investigation of Changepoints in
the Annual Number of Cases of Haemolytic Uraemic Syndrome. Journal of
the Royal Statistical Society. Series C (Applied Statistics) 42 (3), 461–471.
Kass, R. E., Raftery, A. E., 1995. Bayes Factors. Journal of the American
Statistical Association 90 (430), 773–795.
Ko, S. I. M., Chong, T. T. L., Ghosh, P., 2015. Dirichlet Process Hidden
Markov Multiple Change-point Model. Bayesian Analysis 10 (2), 275–296.
Koop, G., Potter, S. M., 2009. Prior Elicitation in Multiple Change-Point
Models. International Economic Review 50 (3), 751–772.
Kullback, S., Leibler, R. A., 1951. On Information and Sufficiency. The Annals of Mathematical Statistics 22 (1), 79–86.
Loschi, R.H., Cruz, F.R.B., 2005. Extension to the Product Partition Model:
Computing the Probability of a Change. Computational Statistics & Data
Analysis 48 (2), 255 – 268.
Moreno, E., Casella, G., Garcia-Ferrer, A., 2005. An Objective Bayesian
Analysis of the Change Point Problem. Stochastic Environmental Research
and Risk Assessment 19 (3), 191–204.
Muliere, P., Scarsini, M., 1985. Change-point Problems: A Bayesian Nonparametric Approach. Aplikace matematiky 30 (6), 397–402.
Petrone, S., Raftery, A. E., 1997. A Note on the Dirichlet Process Prior in
Bayesian Nonparametric Inference with Partial Exchangeability. Statistics
& Probability Letters 36 (1), 69 – 83.
Raftery, A. E., Akman, V. E., 1986. Bayesian Analysis of a Poisson Process
with a Change-Point. Biometrika 73 (1), 85–89.
Schwaller, L., Robin, S., 2017. Exact Bayesian Inference for Off-line Changepoint Detection in Tree-structured Graphical Models. Statistics and Computing 27 (5), 1331–1345.
Smith, A. F. M., 1975. A Bayesian Approach to Inference about a ChangePoint in a Sequence of Random Variables. Biometrika 62 (2), 407–416.
Stephens, D. A., 1994. Bayesian Retrospective Multiple-Changepoint Identification. Journal of the Royal Statistical Society. Series C (Applied Statistics) 43 (1), 159–178.
25

Tian, G.-L., Ng, K. W., Li, K.-C., Tan, M., 2009. Non-iterative Samplingbased Bayesian Methods for Identifying Changepoints in the Sequence of
Cases of Haemolytic Uraemic Syndrome. Computational Statistics and
Data Analysis 53 (9), 3314–3323.
Villa, C., Walker, S., 2015a. An Objective Approach to Prior Mass Functions for Discrete Parameter Spaces. Journal of the American Statistical
Association 110 (511), 1072–1082.
Villa, C., Walker, S., 2015b. An Objective Bayesian Criterion to Determine
Model Prior Probabilities. Scandinavian Journal of Statistics 42 (4), 947–
966.
Yu, J., 2001. Chapter 6 - Testing for a Finite Variance in Stock Return
Distributions. In: Knight, J., Satchell, S. E. (Eds.), Return Distributions
in Finance (Quantitative Finance). Butterworth-Heinemann, Oxford, pp.
143 – 164.

26

Appendix
A

Model prior probabilities to select among models
M0 , M1 and M2

Here, we show how model prior probabilities can be derived for the relatively
simple case of selecting among scenarios with no change points (M0 ), one
change point (M1 ) or two change points (M2 ). First, by applying the result
in Theorem 2, we derive the Kullback–Leibler divergences between any two
models. That is:
• the prior probability for model M0 depends on the following quantities:
DKL (M0 kM1 ) =(n − m1 ) · DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))
DKL (M0 kM2 ) =(m2 − m1 ) · DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))
+ (n − m2 ) · DKL (f1 (·|θ̃1 )kf3 (·|θ̃3 ))
• the prior probability for model M1 depends on the following quantities:
DKL (M1 kM2 ) =(n − m2 ) · DKL (f2 (·|θ̃2 )kf3 (·|θ̃3 ))
DKL (M1 kM0 ) =(n − m1 ) · DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 ))
• the prior probability for model M2 depends on the following quantities:
DKL (M2 kM1 ) =(n − m2 ) · DKL (f3 (·|θ̃3 )kf2 (·|θ̃2 ))
DKL (M2 kM0 ) =(m2 − m1 ) · DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 ))
+ (n − m2 ) · DKL (f3 (·|θ̃3 )kf1 (·|θ̃1 ))
The next step is to derive the minimum Kullback–Leibler divergence computed at each model:

27

• for model M0 :

 
inf DKL (M0 kM1 ) = inf (n − m1 ) · inf DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))
θ1
m1 6=n
θ̃2
|
{z
}


1

= inf DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))
θ̃

2
 
inf DKL (M0 kM2 ) = inf (m2 − m1 ) · inf DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))
θ2
m1 6=m2
θ̃2
{z
}
|
1

 

+ inf (n − m2 ) · inf DKL (f1 (·|θ̃1 )kf3 (·|θ̃3 ))
m2 6=n
θ̃3
|
{z
}
1

= inf DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 )) + inf DKL (f1 (·|θ̃1 )kf3 (·|θ̃3 ))
θ̃2

θ̃3

• for model M1 :
 

inf DKL (M1 kM2 ) = inf (n − m2 ) · inf DKL (f2 (·|θ̃2 )kf3 (·|θ̃3 ))
θ2
m2 6=n
θ̃3
|
{z
}


1

= inf DKL (f2 (·|θ̃2 )kf3 (·|θ̃3 ))
θ̃3

inf DKL (M1 kM0 ) =(n − m1 ) · inf DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 ))
θ0 =θ̃1

θ̃1

• for model M2 :
inf DKL (M2 kM1 ) =(n − m2 ) · inf DKL (f3 (·|θ̃3 )kf2 (·|θ̃2 ))
θ1

θ̃2

inf DKL (M2 kM0 ) =(m2 − m1 ) · inf DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 ))
θ0 =θ̃1

θ̃1

+ (n − m2 ) · inf DKL (f3 (·|θ̃3 )kf1 (·|θ̃1 ))
θ̃1

Therefore, the model prior probabilities can be computed through equation
(15), so that:
• the model prior probability Pr(M0 ) is proportional to the exponential
of the minimum between:




Eπ0 inf DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 )) , Eπ0 inf DKL (f1 (·|θ̃1 )kf2 (·|θ̃2 ))
θ̃2
θ̃2

+ inf DKL (f1 (·|θ̃1 )kf3 (·|θ̃3 ))
θ̃3

28

• the model prior probability Pr(M1 ) is proportional to the exponential
of the minimum between:



Eπ1 inf DKL (f2 (·|θ̃2 )kf3 (·|θ̃3 )) ,
θ̃

 3
Eπ1 (n − m1 ) · inf DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 ))
θ̃1

• the model prior probability Pr(M2 ) is proportional to the exponential
of the minimum between:



Eπ2 (n − m2 ) · inf DKL (f3 (·|θ̃3 )kf2 (·|θ̃2 )) ,
θ̃2

Eπ2 (m2 − m1 ) · inf DKL (f2 (·|θ̃2 )kf1 (·|θ̃1 )) + (n − m2 )
θ̃1

· inf DKL (f3 (·|θ̃3 )kf1 (·|θ̃1 ))
θ̃1

29

B

Proofs

Proof of Theorem 1
We distinguish two cases: S = +1 and S = −1. When S = +1, equivalent
to mj < m0j :
!
Z
(n)
f
(x
|m,
θ̃)
DKL (f (x(n) |m, θ̃)kf (x(n) |m0 , θ̃)) = f (x(n) |m, θ̃) · ln
dx(n)
f (x(n) |m0 , θ̃)
 0
!
mj
Z
X
fj+1 (xi |θ̃j+1 )  (n)
= f (x(n) |m, θ̃) · 
dx
ln
fj (xi |θ̃j )
i=mj +1
0

=

mj Z
X

"
f (x(n) |m, θ̃) · ln

i=mj +1

=

m0j (
X

1n−1 ·

fj+1 (xi |θ̃j+1 )
fj (xi |θ̃j )
"

Z

fj+1 (xi |θ̃j+1 ) · ln

i=mj +1

!#
dx(n)

fj+1 (xi |θ̃j+1 )
fj (xi |θ̃j )

!#

)
dxi

m0

=

j
X

DKL (fj+1 (xi |θ̃j+1 )kfj (xi |θ̃j ))

i=mj +1

=(m0j − mj ) · DKL (fj+1 (·|θ̃j+1 )kfj (·|θ̃j ))
=(m0j − mj ) · d+1
j (θ̃).

(25)

When S = −1, equivalent to mj > m0j , in a similar fashion, we get
DKL (f (x(n) |m, θ̃)kf (x(n) |m0 , θ̃)) = (mj − m0j ) · d−1
j (θ̃)

(26)

From equations (25) and (26), we get the result in Theorem 1.
Proof of Theorem 2
We recall that the model parameter θi is the vector (m1 , m2 , . . . , mi , θ̃1 , θ̃2 , . . . , θ̃i+1 ),
where i = 0, 1, . . . , k. Here, θ̃1 , θ̃2 , . . . , θ̃i+1 represent the parameters of the
underlying sampling distributions considered under model Mi and m1 , m2 , . . . , mi
are the respective i change point locations. In this setting,
f (x

(n)

|θi ) =

m1
Y
r=1

f1 (xr |θ̃1 )

i−1 m
t+1
Y
Y

ft+1 (xr |θ̃t+1 )

t=1 r=mt +1

30

n
Y
r=mi +1

fi+1 (xr |θ̃i+1 )

(27)

We proceed to the computation of DKL (Mi kMj ), that is the Kullback–Leibler
divergence introduced in Section 3. Similarly to the proof of Theorem 1, we
obtain the following result.
!
mi+2
X Z
f
(x
|
θ̃
)
i+1
r
i+1
dx(n)
DKL (Mi kMj ) =
f (x(n) |θi ) ln
f
(x
|
θ̃
)
i+2 r i+2
r=mi+1 +1
!
Z
mi+3
X
f
(x
|
θ̃
)
i+1 r i+1
+
f (x(n) |θi ) ln
dx(n) +
fi+3 (xr |θ̃i+3 )
r=mi+2 +1
!
Z
n
X
f
(x
|
θ̃
)
i+1
r
i+1
... +
f (x(n) |θi ) ln
dx(n) .
f
(x
|
θ̃
)
j+1 r j+1
r=mj +1
Given equation (27), if we integrate out the variables not involved in the
logarithms, we obtain
DKL (Mi kMj ) =(mi+2 − mi+1 ) · DKL (fi+1 (·|θ̃i+1 )kfi+2 (·|θ̃i+2 ))
+ (mi+3 − mi+2 ) · DKL (fi+1 (·|θ̃i+1 )kfi+3 (·|θ̃i+3 ))+
. . . + (n − mj ) · DKL (fi+1 (·|θ̃i+1 )kfj+1 (·|θ̃j+1 )).
In a similar fashion, it can be shown that
DKL (Mj kMi ) =(mi+2 − mi+1 ) · DKL (fi+2 (·|θ̃i+2 )kfi+1 (·|θ̃i+1 ))
+ (mi+3 − mi+2 ) · DKL (fi+3 (·|θ̃i+3 )kfi+1 (·|θ̃i+1 ))+
. . . + (n − mj ) · DKL (fj+1 (·|θ̃j+1 )kfi+1 (·|θ̃i+1 ))

31

