Bayesian Loss-based Approach to Change Point Analysis
Laurentiu Hinoveanu
\noindent 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.
Arxiv_R1_DraftV2.tex
\documentclass[a4paper,12pt,authoryear]{article}
\usepackage{xcolor,amsmath}
\usepackage{amsthm,bm}
\usepackage{hyperref}
\usepackage{amssymb}
\usepackage{dsfont}
\usepackage{comment}
\usepackage{float}
\usepackage{graphicx}
\usepackage{slashbox}
\usepackage{pbox}
\usepackage{colortbl}
\usepackage{epstopdf}
\usepackage{mathtools}
\usepackage{multirow}
\usepackage{relsize}
\usepackage{bigints}
\usepackage{cases}
\usepackage[font={footnotesize}]{caption}
\usepackage{hhline}
\usepackage{booktabs}



%\usepackage[justification=centering,margin=2cm]{caption}
%\usepackage{datetime}
%\delimitershortfall=-1pt

\usepackage{natbib}

\newtheorem{example}{\textit{Example}}
\newtheorem{definition}{Definition}
\newtheorem{theorem}{Theorem}
\newtheorem{lemma}{Lemma}
\newtheorem{proposition}{Proposition}
\newtheorem{axiom}{Axiom}
\newtheorem{remark}{\textit{Remark}}

\theoremstyle{remark}
\newtheorem{myremark}{Remark}

\theoremstyle{plain}
\newtheorem{mylemma}{Lemma}


\newcommand\myiid{\mathrel{\overset{\makebox[0pt]{\mbox{\normalfont\tiny\sffamily i.i.d.}}}{\sim}}}
\newcommand{\ie}{\textit{i}.\textit{e}. }
\newcommand*\diff{\mathop{}\!\mathrm{d}}

\title{Bayesian Loss-based Approach to Change Point Analysis}

\usepackage{authblk}

\author[1]{Laurentiu Hinoveanu\thanks{lch36@kent.ac.uk}}
\author[1]{Fabrizio Leisen\thanks{F.Leisen@kent.ac.uk}}
\author[1]{Cristiano Villa\thanks{C.Villa-88@kent.ac.uk}}
\affil[1]{School of Mathematics, Statistics and Actuarial Science, University of Kent}

%\author{Laurentiu Hinoveanu}

%from http://tex.stackexchange.com/questions/180238/custom-date-on-latex
\usepackage[english]{isodate}

%from http://www.howtotex.com/packages/customize-the-date-format-in-your-latex-documents/
%\newdateformat{mydate}{\monthname[\THEMONTH] \THEYEAR} 
%\date{\mydate\today}

\usepackage{parskip}
\setlength{\parindent}{0pt} 

%for define equal sign from http://tex.stackexchange.com/questions/74125/how-do-i-put-text-over-symbols
\newcommand\mydefeq{\mathrel{\overset{\makebox[0pt]{\mbox{\normalfont\tiny\sffamily def}}}{=}}}
\date{}



\begin{document}
\maketitle

\begin{abstract}
\noindent
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.
\end{abstract}
\textbf{Keywords:} Change point; Discrete parameter space; Loss-based prior; Model selection.


\section{Introduction}\label{sc_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 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. \citet{Chernoff1964} focused on the change in the means of normally distributed variables. \citet{Smith1975} 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. \citet{Smith1975} focuses on the binomial and normal distributions. In \citet{Muliere1985} the problem is tackled from a Bayesian nonparametric perspective. The authors consider Dirichlet processes with independent base measures as underlying distributions. In this framework, \citet{Petrone1997} 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. \citet{Raftery} have approached the single change point problem in the context of a Poisson likelihood under both proper and improper priors for the model parameters. \citet{Carlin} build on the work of \citet{Raftery} 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, \citet{Loschi2005} have provided a fully Bayesian treatment for the product partitions model of \citet{Barry1992}. Their application focused on stock exchange data. \citet{Stephens1994} has extended the Gibbs sampler introduced by \citet{Carlin} in the change point literature to handle multiple change points. \citet{Hannart2009} have used Bayesian decision theory, in particular 0-1 cost functions, to estimate multiple changes in homoskedastic normally distributed observations. \citet{Schwaller2017} extend the product partition model of \citet{Barry1992} by adding a graphical structure which could capture the dependencies between multivariate observations. \citet{Fearnhead2007} proposed a filtering algorithm for the sequential  multiple change points detection problem in the case of piecewise regression models.
\citet{Henderson1993} 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 \citet{Tian2009}, who proposed a method which treats the change points as latent variables. 
\citet{ko2015} have proposed an extension to the hidden Markov model of \citet{Chib} 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. \citet{Heard2017} 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 \citet{Moreno2005} and \cite{Giron}. 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 $X_1, X_2,\dots,X_n$ 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 \citet{VillaMass}. 

To illustrate the idea, consider a probability distribution $f(x|m)$, where $m\in \mathbb{M}$ is a discrete parameter. Then, the prior $\pi(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 \citet{berk1966}, 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, $D_{KL}(f(\cdot|m)\|f(\cdot|m^{\prime}))$, where $m^{\prime}$ 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:
\begin{align}
\pi(m)& \propto  \exp\left\lbrace\displaystyle\min_{m^{\prime} \neq m }D_{KL}(f(\cdot|m)\|f(\cdot|m^{\prime}))\right\rbrace-1,
\label{VillaMass}
\end{align} 
where the Kullback--Leibler divergence  \citep{kullback1951} from the sampling distribution with density $f(x|m)$ to the one with density $f(x|m^{\prime})$ is defined as:
\begin{equation*}
D_{KL}(f(\cdot|m)\|f(\cdot|m^{\prime}))=
      \displaystyle\mathlarger{\int}\displaylimits_{\mathcal{X}} f(x|m) \cdot \log\left[\dfrac{f(x|m)}{f(x|m^{\prime})}\right]\diff x.
   \label{KLdef}
\end{equation*}
Throughout the paper, the objective prior defined in equation \eqref{VillaMass} will be referenced as the \textit{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 \citet{VillaModel}, where a method to define a prior on the space of models is proposed.
To illustrate, let us consider $k$ Bayesian models:
\begin{equation}
M_{j}=\{f_{j}(x|\theta_{j}), \pi_{j}(\theta_{j})\}\qquad   j\in\{1,2,\ldots, k\},
\label{modeldef}
\end{equation}

where $f_{j}(x|\theta_{j})$ is the sampling density characterised by $\theta_{j}$ and  $\pi_{j}(\theta_{j})$ represents the prior on the model parameter. 


Assuming the prior on the model parameter, $\pi_{j}(\theta_{j})$, is proper, the model prior probability $\Pr(M_{j})$ is proportional to the expected minimum Kullback--Leibler divergence from $M_j$, where the expectation is considered with respect to  $\pi_{j}(\theta_{j})$. That is: 
\begin{align}
\Pr(M_{j})\propto &\exp\left\lbrace  \mathbb{E}_{\pi_{j}}\left[\inf_{\theta_{i},i\neq j} D_{KL}(f_{j}(x|\theta_{j})\|f_{i}(x|\theta_{i}))\right]\right\rbrace \qquad    j=1, \ldots, k. \label{VillaModelCompactOrg}
\end{align}


%Across this paper, we use the following shorthand notation for equation \eqref{VillaModelCompactOrg}:
%\begin{align*}
%\Pr(M_{j})\propto &\exp\left\lbrace  \mathbb{E}_{\pi_{j}}\left[\inf_{\theta_{i},i\neq j} D_{KL}(M_j\|M_i)\right]\right\rbrace\qquad    j=1, \ldots, k. 
%\end{align*}

The model prior probabilities defined in equation \eqref{VillaModelCompactOrg} can be employed to derive the model posterior probabilities through:
\begin{equation}
\Pr(M_{i}|x)= \left[\sum_{j=1}^{k}\dfrac{\Pr(M_{j})}{\Pr(M_{i})} B_{ji}\right]^{-1},
\label{modelposterior}
\end{equation}
where $B_{ji}$ is the Bayes factor between model $M_{j}$ and model $M_{i}$, defined as 
\begin{equation*}
B_{ji}=\dfrac{\int f_{j}(x|\theta_{j})\pi_{j}(\theta_{j})\diff \theta_{j}}{\int f_{i}(x|\theta_{i})\pi_{i}(\theta_{i})\diff \theta_{i}},
\end{equation*}
with $i\neq j \in \{1, 2, \ldots, k\}$. 

This paper is structured as follows: in Section \ref{sc_locations} we establish the way we set objective priors on both single and multiple change point locations. Section \ref{sc_numberCP} 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 \ref{sc_simulation} and \ref{sc_realdata}, where we work with simulated and real data, respectively. Section \ref{Conclusion} is dedicated to final remarks.

\section{Objective Prior on the Change Point Locations}\label{sc_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 $m_1<m_2<\ldots<m_k$ 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 \ref{sub_sc_OneChgPoint}). Then, we extend the results to the more general case where multiple change points are assumed (see Section \ref{sub_sc_MultipleChgPoints}). 

A well-known objective prior for finite parameter spaces, in cases where there is no structure, is the uniform prior \citep{BBS2012}. As such, a natural choice for the prior on the change points location is the uniform \citep{Koop}. 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.

\subsection{Single Change Point}
\label{sub_sc_OneChgPoint}
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, \ldots, n-1\}$.
 

Let $\mathbf{X}^{(n)}=(X_{1},\ldots, X_{n})$ denote an \textit{n}-dimensional vector of random variables, representing the random sample, and $m$ be our single change point location, that is $m \in \{1, 2, \ldots, n-1\}$, such that
\begin{align}
X_{1}, \ldots, X_{m}| \tilde\theta_{1} &\myiid f_{1}(\cdot|\tilde\theta_{1}) \nonumber \\
X_{m+1}, \ldots, X_{n}| \tilde\theta_{2} &\myiid f_{2}(\cdot|\tilde\theta_{2}).
\label{SingleChangePointCase}
\end{align}
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 $\tilde\theta_1 \neq \tilde\theta_2$ when $f_1=f_2$. The sampling density for the vector of observations  $\mathbf{x}^{(n)}=(x_1,\dots,x_n)$ is:
\begin{equation}
f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})=\prod_{i=1}^{m} f_{1}(x_{i}|\tilde\theta_{1})\prod_{i=m+1}^{n} f_{2}(x_{i}|\tilde\theta_{2}).
\label{likelihood}
\end{equation}
Let $m^{\prime} \neq m$. Then, the Kullback--Leibler divergence between the model parametrised by $m$ and the one parametrised by $m^{\prime}$ is:
\begin{align}
D_{KL}(f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})\|f(\mathbf{x}^{(n)}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))=& \int f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})  \nonumber \\ &\hspace*{1em}\log\left(\dfrac{f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})}{f(\mathbf{x}^{(n)}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2})}\right)\,\mathrm{d}\mathbf{x}^{(n)}.
\label{klremarkOneChg}
\end{align}
Without loss of generality, consider $m<m^{\prime}$. In this case, note that
\begin{align*}
\dfrac{f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})}{f(\mathbf{x}^{(n)}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2})}=\prod_{i=m+1}^{m^{\prime}} \dfrac{f_{2}(x_{i}|\tilde\theta_{2})}{f_{1}(x_{i}|\tilde\theta_{1})},
\end{align*}
leading to
%\begin{align}
%&D_{KL}(f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})\|f(\mathbf{x}^{(n)}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))=\int \left( \prod_{j=1}^{m} f_{1}(x_{j}|\tilde\theta_{1})\right)  \nonumber\\ &\hspace*{1em}\cdot \left( \prod_{i=m+1}^{m^{\prime}} f_{2}(x_{i}|\tilde\theta_{2})\right) \cdot \left( \prod_{j=m^{\prime}+1}^{n} f_{2}(x_{j}|\tilde\theta_{2})\right) \cdot \left[\sum_{i=m+1}^{m^{\prime}} \ln\left(\dfrac{f_{2}(x_{i}|\tilde\theta_{2})}{f_{1}(x_{i}|\tilde\theta_{1})}\right)\right]\,\mathrm{d}\mathbf{x}^{(n)} \nonumber \\
%&=\sum_{i=m+1}^{m^{\prime}} \int \left( \prod_{j=1}^{m} f_{1}(x_{j}|\tilde\theta_{1})\right) \cdot \left( \prod_{j=m^{\prime}+1}^{n} f_{2}(x_{j}|\tilde\theta_{2})\right)\cdot \left(\prod_{\substack{r=m+1\\ r \neq i}}^{m^{\prime}} f_{2}(x_{r}|\tilde\theta_{2})\right) \nonumber  \\ &\hspace*{1em}\cdot \left[ f_{2}(x_{i}|\tilde\theta_{2})\cdot\ln\left(\dfrac{f_{2}(x_{i}|\tilde\theta_{2})}{f_{1}(x_{i}|\tilde\theta_{1})}\right)\right]\,\mathrm{d}\mathbf{x}^{(n)}.
%\label{step1}
%\end{align}
%Integrating with respect to all variables which are not indexed by $i$ in equation \eqref{step1}, we obtain: %\cdot \phantom{\dfrac{1}{2}}\right. \nonumber\\ &\left.
\begin{align}
D_{KL}(f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})&\|f(\mathbf{x}^{(n)}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))=\nonumber \\ &\sum_{i=m+1}^{m^{\prime}}  \int f_{2}(x_{i}|\tilde\theta_{2})\log\left(\dfrac{f_{2}(x_{i}|\tilde\theta_{2})}{f_{1}(x_{i}|\tilde\theta_{1})}\right)\,\mathrm{d}x_{i}.
\label{step2}
\end{align}
On the right hand side of equation \eqref{step2}, we can recognise the Kullback--Leibler divergence from density $f_{2}$ to density $f_{1}$, thus getting:
\begin{align}
D_{KL}(f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})&||f(\mathbf{x}^{(n)}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))=\nonumber \\ 
&\hspace{2cm}(m^{\prime}-m) D_{KL}(f_{2}(\cdot|\tilde\theta_{2})\|f_{1}(\cdot|\tilde\theta_{1})).
\label{step3case1}
\end{align}
%By analogy with the steps taken in Equations \eqref{step1},\eqref{step2} and \eqref{step3case1}, we derive the following relation for the case $m>m^{\prime}$:
In a similar fashion, when $m>m^{\prime}$, we have that:
\begin{align}
D_{KL}(f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})&\|f(\mathbf{x}^{(n)}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))=\nonumber \\ 
&\hspace{2cm}(m-m^{\prime}) D_{KL}(f_{1}(\cdot|\tilde\theta_{1})\|f_{2}(\cdot|\tilde\theta_{2})).
\label{step3case2}
\end{align}

In this single change point scenario, we can consider $m^{\prime}$ as a perturbation of the change point location $m$, that is $m^{\prime}=m\pm l$ where $l \in \mathbb{N}^{*}$, such that $1\leq m^{\prime}<n$. Then, taking into account equations \eqref{step3case1} and \eqref{step3case2}, the Kullback--Leibler divergence becomes:

\begin{align*}
D_{KL}(f(\mathbf{x}^{(n)}|m,\tilde\theta_{1},\tilde\theta_{2})\|f(\mathbf{x}^{(n)}&|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))=\\
&\left\{
\begin{array}{ll}
l \cdot D_{KL}(f_{2}(\cdot|\tilde\theta_{2})\|f_{1}(\cdot|\tilde\theta_{1})), & \mbox{ if } m<m^{\prime}\\
&\\
l \cdot D_{KL}(f_{1}(\cdot|\tilde\theta_{1})\|f_{2}(\cdot|\tilde\theta_{2})), & \mbox{ if } m>m^{\prime},
\end{array}\right.\\
\end{align*}
and
\begin{align}
\min_{m^{\prime}\neq m}&\left[D_{KL}(f(\mathbf{x^{(n)}}|m,\tilde\theta_{1},\tilde\theta_{2})\|f(\mathbf{x^{(n)}}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))\right]=\nonumber \\&=\min_{m^{\prime}\neq m}\lbrace l \cdot D_{KL}(f_{2}(\cdot|\tilde\theta_{2})\|f_{1}(\cdot|\tilde\theta_{1})),l \cdot D_{KL}(f_{1}(\cdot|\tilde\theta_{1})\|f_{2}(\cdot|\tilde\theta_{2}))\rbrace\nonumber \\&= \min_{m^{\prime}\neq m}\lbrace D_{KL}(f_{2}(\cdot|\tilde\theta_{2})\|f_{1}(\cdot|\tilde\theta_{1})),D_{KL}(f_{1}(\cdot|\tilde\theta_{1})\|f_{2}(\cdot|\tilde\theta_{2}))\rbrace  \cdot \underbrace{\min_{m^{\prime}\neq m}\{l\}}_{1}.
\label{KL_no_m}
\end{align}
%From equation \eqref{KL_no_m}, we observe that $\displaystyle\min_{m^{\prime}\neq m}\left[D_{KL}(f(\mathbf{x^{(n)}}|m,\tilde\theta_{1},\tilde\theta_{2})\|\right.\\ \left.f(\mathbf{x^{(n)}}|m^{\prime},\tilde\theta_{1},\tilde\theta_{2}))\right]$ is only a function of $\tilde\theta_{1}$ and $\tilde\theta_{2}$ and does not depend by $m$. Thus, if we apply the Loss-based prior from equation \eqref{VillaMass} to our discrete parameter $m$ and we recall the non-negativity property of the Kullback--Leibler divergence and that $\tilde\theta_{1}\neq \tilde\theta_{2}$, we obtain:
We observe that equation \eqref{KL_no_m} is only a function of $\tilde\theta_{1}$ and $\tilde\theta_{2}$ and does not depend on $m$. Thus, $\pi(m)\propto 1$ and, therefore,
\begin{align}
\pi(m)=\dfrac{1}{n-1} \qquad m\in \{1,\ldots, n-1\}.
\label{SimpleOneChgPriorForm}
\end{align}
This prior was used, for instance, in an econometric context by \citet{Koop} with the rationale of giving equal weight to every possible change point location.

\subsection{Multivariate Change Point Problem}
\label{sub_sc_MultipleChgPoints}

In this section, we address the change point problem in its generality by assuming that there are $1\leq k < n$ change points. In particular, for the data $\mathbf{x}^{(n)}=(x_1,\dots,x_n)$, we consider the following sampling distribution 
\begin{equation}\label{SD_multi}
f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta})=\prod_{i=1}^{m_{1}} f_{1}(x_{i}|\tilde\theta_{1}) \prod_{j=1}^{k-1}\prod_{i=m_{j}+1}^{m_{j+1}} f_{j+1}(x_{i}|\tilde\theta_{j+1})\prod_{i=m_{k}+1}^{n} f_{k+1}(x_{i}|\tilde\theta_{k+1}),
\end{equation}
where $\bm{m}=(m_1,\dots,m_k)$, $1\leq m_1<m_2<\ldots<m_k<n$, is the vector of the change point locations and $\bm{\tilde\theta}=(\tilde\theta_1,\dots,\tilde\theta_k,\tilde\theta_{k+1})$ is the vector of the parameters of the underlying probability distributions. Schematically:
\begin{equation*}
\begin{array}{lcll}
X_{1}&, \ldots, & X_{m_{1}}|\tilde\theta_{1} &\myiid f_{1}(\cdot|\tilde\theta_{1}) \\
X_{m_{1}+1}&, \ldots, &X_{m_2}|\tilde\theta_{2} &\myiid f_{2}(\cdot|\tilde\theta_{2})\\
\vdots&,\ldots, &\vdots &\vdots\ldots\vdots\\
X_{m_{k-1}+1}&, \ldots, &X_{m_k}|\tilde\theta_{k} &\myiid f_{k}(\cdot|\tilde\theta_{k})\\
X_{m_{k}+1}&, \ldots, &X_{n}|\tilde\theta_{k+1} &\myiid f_{k+1}(\cdot|\tilde\theta_{k+1}).\\
\end{array}
\end{equation*}
If $f_1=f_2=\cdots=f_{k+1}$, then it is reasonable to assume that some of the $\theta$'s are different. Without loss of generality, we assume that $\tilde\theta_1\neq\tilde\theta_2\neq\cdots\neq\tilde\theta_k\neq\tilde\theta_{k+1}$. In a similar fashion to the single change point case, we cannot assume $m_k=n$ since we require exactly $k$ change points.

In this case, due to the multivariate nature of the vector $\bm{m}=(m_1,\dots,m_k)$, 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 \ref{lemmaOneChgPoint} (the proof of which is in the Appendix).
In particular, we are able to prove an analogous of equations \eqref{step3case1} and \eqref{step3case2} when only one component is arbitrarily perturbed. 
Let us define the following functions:
\begin{align*}
d^{+1}_j(\bm{\tilde\theta})&=D_{KL}(f_{j+1}(\cdot|\tilde\theta_{j+1})\|f_{j}(\cdot|\tilde\theta_{j}))\\
d^{-1}_j(\bm{\tilde\theta})&=D_{KL}(f_{j}(\cdot|\tilde\theta_{j})\|f_{j+1}(\cdot|\tilde\theta_{j+1})),
\end{align*}
where $j \in \{1, 2, \ldots, k\}$. The following Theorem  is useful to understand the behaviour of the loss-based prior in the general case.\\
\begin{theorem} Let $f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta})$ be the sampling distribution defined in equation \eqref{SD_multi} and consider  $j\in\lbrace 1,\dots,k\rbrace$. Let $\bm{m}'$ be such that $m'_i=m_i$ for $i\neq j$, and let the component $m'_j$ be such that $m'_j\neq m_j$ and $m_{j-1}<m'_j<m_{j+1}$. Therefore,
$$D_{KL}(f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta})\|f(\mathbf{x}^{(n)}|\bm{m}',\bm{\tilde\theta})=|m'_j-m_j| d^S_j(\bm{\tilde\theta}),$$
where $S=\emph{sgn}(m'_j-m_j)$.
%$$m_i'=\left\lbrace\begin{array}{ll}
%m_i&\mbox{ if } i\neq j\\
%p& \mbox{ if } i=j
%\end{array}\right.$$
 
\label{lemmaOneChgPoint}
\end{theorem}

Note that, Theorem \ref{lemmaOneChgPoint} states that the minimum Kullback--Leibler divergence is achieved  when $m'_j=m_j+1$ or $m'_j=m_j-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 \ref{lemmaOneChgPoint} 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 $\bm{m}$ is perturbed by (plus or minus) one unit. As such, the loss-based prior depends on the vector of parameters $\bm{\tilde\theta}$ only, as in the one-dimensional case, yielding the uniform prior for $\bm{m}$.

Therefore, the loss-based prior on the multivariate change point location is  
\begin{align}
\pi(\bm{m})&=\left\lbrace\dbinom{n-1}{k}\right\rbrace^{-1},
\label{MultipleChg}
\end{align}
where $\bm{m}=(m_1,\dots,m_k), 1\leq m_1<m_2<\ldots<m_k<n$.
The denominator in equation \eqref{MultipleChg} 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$ elements, which is $\binom{n-1}{k}$. The same prior was also derived in a different way by \citet{Giron}. 

\section{Loss-based Prior on the Number of Change Points}\label{sc_numberCP}

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 \cite{VillaModel}. 
\begin{figure}[h]
\centering   
\vspace{-6em}
\includegraphics[scale=0.5]{ModelSelectArrowTilde.pdf}
\vspace{-18em}
\caption{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.}
\label{ModelConstruct}
\end{figure}
We proceed as follows. Assume we have to select from $k+1$ possible models. Let $M_{0}$ be the model with no change points, $M_{1}$ the model with one change point and so on. Generalising, model $M_{k}$ 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 $M_{3}$ the first two change point locations will be the same as in the case of model $M_{2}$. To illustrate the way we envision our models, we have provided Figure \ref{ModelConstruct}. It has to be noted that the construction of the possible models from $M_0$ to $M_k$ can be done in a different way to one here described. Obviously, the approach to define the model priors stays unchanged. Consistently with the notation used in Section 1, 
$$\theta_k=\left\{
\begin{array}{ll}
\tilde{\theta}_1,\ldots, \tilde{\theta}_{k+1},m_1, \ldots, m_k &\mbox{ if } k=1,\dots,n-1\\
\tilde{\theta}_1&\mbox{ if }k=0
\end{array}\right.
$$ represents the vector of parameters of model $M_k$, where $\tilde{\theta}_1,\ldots, \tilde{\theta}_{k+1}$ are the model specific parameters and $m_1, \ldots, m_k$ are the change point locations, as in Figure \ref{ModelConstruct}. 

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 \ref{remark_model} (the proof of which is in the Appendix). 
\vspace*{0.5cm}
\begin{theorem}
Let 
$$D_{KL}(M_{i}\|M_{j})=D_{KL}(f(\mathbf{x}^{(n)}|\theta_i)\|f(\mathbf{x}^{(n)}|\theta_j)).$$
For any $0\leq i<j \leq k$ integers, with $k<n$, and the convention $m_{j+1}=n$, we have the following:
\begin{align*}
 &D_{KL}(M_{i}\|M_{j})=\sum_{q=i+1}^{j}\left[(m_{q+1}-m_{q})\cdot D_{KL}(f_{i+1}(\cdot|\tilde{\theta}_{i+1})\|f_{q+1}(\cdot|\tilde{\theta}_{q+1}))\right],
\end{align*}
and
\begin{align*}
 &D_{KL}(M_{j}\|M_{i})=\sum_{q=i+1}^{j}\left[(m_{q+1}-m_{q})\cdot D_{KL}(f_{q+1}(\cdot|\tilde{\theta}_{q+1})\|f_{i+1}(\cdot|\tilde{\theta}_{i+1}))\right].
\end{align*}
\label{remark_model}
\end{theorem}

The result in Theorem \ref{remark_model} is useful when the model selection exercise is implemented. Indeed, the \cite{VillaModel} approach requires the computation of the Kullback--Leibler divergences in Theorem \ref{remark_model}. Recalling equation \eqref{VillaModelCompactOrg}, the objective model prior probabilities are then given by:
\begin{align}
\Pr(M_{j})\propto &\exp\left\lbrace  \mathbb{E}_{\pi_{j}}\left[\inf_{\theta_{i},i\neq j} D_{KL}(M_j\|M_i)\right]\right\rbrace\qquad j=0, 1, \ldots, k. \label{ModelPriors}
\end{align}

For illustrative purposes, in the Appendix we derive the model prior probabilities to perform model selection among $M_0$, $M_1$ and $M_2$. \\

It is easy to infer from equation \eqref{ModelPriors} 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, a sensible choice is the uniform prior which, as shown in Section \ref{sc_locations}, 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 \citep{Intrinsic} may represent a viable solution since they are proper. Nonetheles, the method introduce by \cite{VillaModel} 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. 


{\bf Remark.} In the case where the changes in the underlying sampling distribution are limited to the parameter values, the model prior probabilities defined in \eqref{ModelPriors} follow the uniform distribution. That is, $\Pr(M_j)\propto1$. In the real data example illustrated in Section \ref{sub_sc_British_coal_mine}, we indeed consider a problem where the above case occurs.

\subsection{A special case: selection between $M_0$ and $M_1$}
\label{sub_sc_Model_one_change}
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 $M_0$ (i.e. no change point) and $M_1$ (i.e. one change point). Following our approach, we have:
 
\begin{align}
 \Pr(M_0) \propto \exp\left\{\mathbb{E}_{\pi_{0}}\left[\inf_{\tilde{\theta}_2}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\right]\right\},
 \label{prior_M0_one_chg}
\end{align}
and
\begin{align}
  \Pr(M_1) \propto \exp\left\{\mathbb{E}_{\pi_1}\left[(n-m_{1}) \cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]\right\}
 \label{prior_M1_one_chg}.
\end{align}

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 $\pi_1(m_1,\tilde{\theta}_1, \tilde{\theta}_2)=\pi_1(m_1) \pi_1(\tilde{\theta}_1, \tilde{\theta}_2)$. Let us further recall that, as per equation \eqref{MultipleChg}, $\pi_1(m_1)=1/(n-1)$.  As such, we observe that the model prior probability on $M_1$ becomes:
 \begin{align}
  \Pr(M_1) \propto \exp\left\lbrace\left(\dfrac{n}{2}\right)\mathbb{E}_{\pi_1(\tilde{\theta}_1, \tilde{\theta}_2)}\left[\inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]\right\rbrace.
 \label{prior_M1_trans_one_chg}
\end{align}

We notice that the model prior probability for model $M_1$ is increasing when the sample size increases. This behaviour occurs whether there is or not a change point in the data. We propose to address the above problem by using a non-uniform prior for $m_1$. A reasonable alternative, which works quite well in practice, would be the following shifted binomial as prior:
 \begin{align}
  \pi_1(m_1)=\dbinom{n-2}{m_1-1}\left(\dfrac{n-1}{n}\right)^{m_1-1}\left(\dfrac{1}{n}\right)^{n-m_1-1}, 1\leq m_1 \leq n-1.
  \label{Pseudo_binom}
 \end{align}
To argument the choice of \eqref{Pseudo_binom}, 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 \eqref{prior_M1_one_chg} and \eqref{Pseudo_binom} follows:
 \begin{align}
  \Pr(M_1) \propto \exp\left\lbrace\left(\dfrac{2n-2}{n}\right)\mathbb{E}_{\pi_1(\tilde{\theta}_1, \tilde{\theta}_2)}\left[\inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]\right\rbrace.
 \label{prior_M1_pseudo_binom}
\end{align}
For the more general case where we consider more than two models, the problem highlighted in equation \eqref{prior_M1_trans_one_chg} vanishes. 

\section{Change Point Analysis on Simulated Data}\label{sc_simulation}
In this section, we present the results of several simulation studies based on the methodologies discussed in Sections \ref{sc_locations} and \ref{sc_numberCP}. 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 \citet{VillaModel}.

\subsection{Single sample}

\paragraph*{Scenario 1.} The first scenario concerns the choice between models $M_0$ and $M_1$. Specifically, for $M_0$ we have:
$$X_{1}, X_{2}, \ldots, X_{n}|p \myiid \mbox{Geometric}(p)$$
and for $M_1$ we have:
\begin{align*}
X_{1}, X_{2}, \ldots, X_{m_{1}}&|p \myiid \mbox{Geometric}(p) \\
X_{m_{1}+1}, X_{m_{1}+2}, \ldots, X_{n}&|\lambda \myiid \mbox{Poisson}(\lambda)
\end{align*}

Let us denote with $f_1(\cdot|p)$ and $f_2(\cdot|\lambda)$ the probability mass functions of the Geometric and the Poisson distributions, respectively. The priors for the parameters of $f_1$ and $f_2$ are $
p \sim \mbox{Beta}(a,b)$ and $\lambda \sim \mbox{Gamma}(c,d)$.

In the first simulation, we sample $n=100$ observations from model $M_0$ with $p=0.8$. To perform the change point analysis, we have chosen the following parameters for the priors on $p$ and $\lambda$: $a=2$, $b=2$, $c=3$ and $d=1$. Applying the approach introduced in Section \ref{sc_numberCP}, we obtain $\Pr(M_0)\propto 1.59$ and $\Pr(M_1)\propto 1.81$. These model priors yield the posterior distribution probabilities (refer to equation \eqref{modelposterior}) $\Pr(M_0|\mathbf{x}^{(n)})=0.92$ and $\Pr(M_1|\mathbf{x}^{(n)})=0.08$. As expected, the selection process strongly indicates the true model as $M_0$. Table \ref{Discrete_Results_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 $M_1$, with $p=0.8$ and $\lambda=3$. We have sampled 50 data points from the Geometric distribution and the remaining 50 data points from the Poisson distribution. In Figure \ref{fig_GP}, 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(M_0|\mathbf{x}^{(n)})=0.06$ and $\Pr(M_1|\mathbf{x}^{(n)})=0.94$. Again, the model selection process is assigning heavy posterior mass to the true model $M_1$. These results are further detailed in Table \ref{Discrete_Results_1}.

\begin{figure}[h]
\centering
\includegraphics[scale=0.6]{G_P_model_1.pdf}
\caption{Scatter plot of the data simulated from model $M_1$ in Scenario 1.}
\label{fig_GP}
\end{figure}


\begin{table}[h]
\centering
\begin{tabular}{c|cc}
\hline 
• & \multicolumn{2}{c}{True model} \\ 
%\hline 
• & $M_0$ & $M_1$ \\ 
\hline 
$\Pr(M_0)$ & 0.47 & 0.47 \\ 
%\hline 
$\Pr(M_1)$ & 0.53 & 0.53 \\ 
%\hline 
$B_{01}$ & 12.39 & 0.08 \\ 
%\hline 
$B_{10}$ & 0.08 & 12.80 \\ 
%\hline 
$\Pr(M_0|\mathbf{x}^{(n)})$ & 0.92 & 0.06 \\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.08 & 0.94 \\ 
\hline 
\end{tabular}
\caption{Model prior, Bayes factor and model posterior probabilities for the change point analysis in Scenario 1. We considered samples from, respectively, model $M_0$ and model $M_1$.}
\label{Discrete_Results_1}
\end{table}

\paragraph*{Scenario 2.} In this scenario we consider the case where we have to select among three models, that is model $M_0$:
\begin{align}
X_{1}, X_{2}, \ldots, X_{n}|\lambda, \kappa &\myiid \mbox{Weibull}(\lambda, \kappa) ,\label{Simulation_M_0}
\end{align}
model $M_1$:
\begin{align}
X_{1}, X_{2}, \ldots, X_{m_{1}}|\lambda, \kappa &\myiid \mbox{Weibull}(\lambda, \kappa)\nonumber \\
X_{m_{1}+1}, X_{m_{1}+2}, \ldots, X_{n}|\mu, \tau &\myiid \mbox{Log-normal}(\mu, \tau) \label{Simulation_M_1},
\end{align}
with $1\leq m_{1}\leq n-1$ being the location of the single change point, and model $M_2$:
\begin{align}
X_{1}, X_{2}, \ldots, X_{m_{1}}| \lambda, \kappa &\myiid \mbox{Weibull}(\lambda, \kappa)\nonumber\\
X_{m_{1}+1}, X_{m_{1}+2}, \ldots, X_{m_{2}}| \mu, \tau &\myiid \mbox{Log-normal}(\mu, \tau)\nonumber\\
X_{m_{2}+1}, X_{m_{2}+2}, \ldots, X_{n}|\alpha, \beta &\myiid \mbox{Gamma}(\alpha, \beta), \label{Simulation_M_2}
\end{align}
with $1\leq m_{1}<m_{2}\leq n-1$ representing the locations of the two change points, such that $m_1$ corresponds exactly to the same location as in model $M_1$. Analogously to the previous scenario, we sample from each model in turn and perform the selection to detect the number of change points.

Let $f_{1}(\cdot|\lambda, \kappa), f_{2}(\cdot|\mu, \tau)$ and $f_{3}(\cdot|\alpha, \beta)$ represent the Weibull, Log-normal and Gamma densities, respectively, with $\tilde\theta_1=(\lambda, \kappa)$, $\tilde\theta_2=(\mu, \tau)$ and $\tilde\theta_3=(\alpha, \beta)$. We assume a Normal prior on $\mu$ and Gamma priors on all the other parameters as follows:
$$\lambda\sim\mbox{Gamma}(1.5,1) \qquad \kappa\sim\mbox{Gamma}(5,1) \qquad \mu\sim\mbox{Normal}(0.05,1),$$
$$\quad\tau\sim\mbox{Gamma}(16,1) \qquad \alpha\sim\mbox{Gamma}(10,1) \qquad \beta\sim\mbox{Gamma}(0.2,0.1).$$

In the first exercise, we have simulated $n=100$ observations from model $M_0$, where we have set $\lambda=1.5$ and $\kappa=5$.  We obtain the following model priors: $\Pr(M_0)\propto1.09$, $\Pr(M_1)\propto1.60$ and $\Pr(M_2)\propto1.37$, yielding the posteriors $\Pr(M_0|\mathbf{x}^{(n)})=0.96$, $\Pr(M_1|\mathbf{x}^{(n)})=0.04$ and $\Pr(M_2|\mathbf{x}^{(n)})=0.00$. We then see that the approach assigns high mass to the true model $M_0$. Table \ref{TableResults} reports the above probabilities and the corresponding Bayes factors.
\begin{table}[h]
\centering
\begin{tabular}{c|ccc}
\hline 
• & \multicolumn{3}{c}{True model} \\ 
%\hline 
• & $M_0$ & $M_1$ & $M_2$ \\ 
\hline 
$\Pr(M_0)$ & 0.27 & 0.27 & 0.27 \\ 
%\hline 
$\Pr(M_1)$ & 0.39 & 0.39 & 0.39 \\ 
%\hline 
$\Pr(M_2)$ & 0.34 & 0.34 & 0.34 \\
%hline
$B_{01}$ & 36.55 & $3.24\times10^{-4}$ & $4.65\times10^{-40}$ \\ 
%\hline 
$B_{02}$ & $1.84\times 10^3$ & 0.02 & $1.27\times10^{-45}$ \\ 
%\hline 
$B_{12}$ & 50.44 & 55 & $2.72\times10^{-6}$ \\
%\hline
$\Pr(M_0|\mathbf{x}^{(n)})$ & 0.96 & 0.00 & 0.00\\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.04 & 0.98 & 0.00 \\
%\hline 
$\Pr(M_2|\mathbf{x}^{(n)})$ & 0.00 & 0.02 & 1.00 \\
\hline 
\end{tabular}
\caption{Model prior, Bayes factor and model posterior probabilities for the change point analysis in Scenario 2. We considered samples from, respectively, model $M_0$, model $M_1$ and model $M_2$.}
\label{TableResults}
\end{table}
The second simulation was performed by sampling 50 observations from a Weibull with parameter values as in the previous exercise, and the remaining 50 observations from a Log-normal density with location parameter $\mu=0.05$ and scale parameter $\tau=16$. The data is displayed in Figure \ref{One_change}.
\begin{figure}[H]
\centering
\includegraphics[scale=0.6]{W_L_G_model_1.pdf}
\caption{Scatter plot of the observations simulated from model $M_1$ in Scenario 2.}
\label{One_change}
\end{figure}
The model posterior probabilities are $\Pr(M_0|\mathbf{x}^{(n)})=0.00$, $\Pr(M_1|\mathbf{x}^{(n)})=0.98$ and $\Pr(M_2|\mathbf{x}^{(n)})=0.02$, which are reported in Table \ref{TableResults}. In this case as well, we see that the model selection procedure indicates $M_1$ as the true model, as expected.
\begin{figure}[H]
\centering
\includegraphics[scale=0.6]{W_L_G_model_2.pdf}
\caption{Scatter plot of the observations simulated from model $M_2$ in Scenario 2.}
\label{Two_change}
\end{figure}
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 $\alpha=10$ and $\beta=2$. From Table \ref{TableResults}, we note that the posterior distribution on the model space accumulates on the true model $M_2$.

 \subsection{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. 

%It is not unusual in practice to have relatively large data sets (e.g. financial data). As such, the above is not a limitation of the approach but, actually, a point of strength.

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 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. 

\paragraph*{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 $f_{1}(\cdot|\nu_1), f_{2}(\cdot|\nu_2)$ and $f_{3}(\cdot|\nu_3)$ represent the densities of three standard $t$ distributions, respectively. We assume that $\nu_1, \nu_2$ and $\nu_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:
$$\nu_1\sim 2+\mbox{Poisson}(30) \qquad  \nu_2\sim 2+\mbox{Poisson}(3) \qquad \nu_3\sim 2+\mbox{Poisson}(8).$$
In this experiment, we consider 60 repeated samples, each of size $n=300$ and with the following structure: 
\begin{itemize}
\item $X_1, \dots, X_{100}$ from a Student-$t$ distribution with $\nu_1=30$,
\item $X_{101}, \dots, X_{200}$ from a Student-$t$ distribution with $\nu_2=3$,
\item $X_{201}, \dots, X_{300}$ from a Student-$t$ distribution with $\nu_3=8$.
\end{itemize}
Table \ref{TableResults-T-distrib} reports the frequentist results of the simulation study. First, note that $P(M_1)=P(M_2)=P(M_3)=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 $M_2$. We also note that the above posterior is on average $0.75$  with a variance $0.02$, making the inferential procedure extremely accurate. 
\begin{table}[h]
\centering
\begin{tabular}{c|c|c|c}
\hline 
• & Mean posterior & Variance posterior &  \pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Freq. true model} \\ 
%\hline 
\hline 
%$\Pr(M_0)$ & 0.33 [$3.37\times10^{-9}$]&   \\ 
%%\hline 
%$\Pr(M_1)$ & 0.33 [$3.59\times10^{-9}$]&   \\ 
%%\hline 
%$\Pr(M_2)$ & 0.34 [$1.39\times10^{-8}$] &   \\
%hline
%$B_{01}$ & 36.55 & $3.24\times10^{-4}$ & $4.65\times10^{-40}$ \\ 
%\hline 
%$B_{02}$ & $1.84\times 10^3$ & 0.02 & $1.27\times10^{-45}$ \\ 
%\hline 
%$B_{12}$ & 50.44 & 55 & $2.72\times10^{-6}$ \\
%\hline
$\Pr(M_0|\mathbf{x}^{(n)})$ & 0.01 & $3.84\times10^{-4}$&  0/60\\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.24 & 0.0160  & 0/60 \\
%\hline 
$\Pr(M_2|\mathbf{x}^{(n)})$ & 0.75 & 0.0190 & 60/60 \\
\hline 
%\pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Frequency of identifying \\ the respective model} && \\
%\hline
\end{tabular}
\caption{Average model posterior probabilities, variance and frequency of true model for the Scenario 3 simulation exercise.}
\label{TableResults-T-distrib}
\end{table}


\paragraph*{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 particular, we draw $100$ samples with $n=500$ and $n=1500$. For $n=500$, the loss-based prior probabilities are $P(M_0)=0.18$, $P(M_1)=0.16$ and $P(M_2)=0.66$. For $n=1500$, the loss-based prior probabilities are $P(M_0)=0.015$, $P(M_1)=0.014$ and $P(M_2)=0.971$. The simulation results are reported, respectively, in Table \ref{TableResults-WLG-Conjugate-500} and in Table \ref{TableResults-WLG-Conjugate-1500}.  The two change point locations for $n=500$ are at the $171$st and $341$st observations. For $n=1500$, the first change point is the $501$st observation, while the second is at the $1001$st 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\%$. 


\begin{table}[H]
\centering
\begin{tabular}{c|c|c|c}
\hline 
• & Mean posterior & Variance posterior &  \pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Freq. true model} \\
%\hline 
\hline 
%$\Pr(M_0)$ & 0.18 [$6.66\times10^{-14}$]&   \\ 
%%\hline 
%$\Pr(M_1)$ & 0.16 [$1.53\times10^{-12}$]&   \\ 
%%\hline 
%$\Pr(M_2)$ & 0.66 [$9.54\times10^{-13}$] &   \\
%hline
%$B_{01}$ & 36.55 & $3.24\times10^{-4}$ & $4.65\times10^{-40}$ \\ 
%\hline 
%$B_{02}$ & $1.84\times 10^3$ & 0.02 & $1.27\times10^{-45}$ \\ 
%\hline 
%$B_{12}$ & 50.44 & 55 & $2.72\times10^{-6}$ \\
%\hline
$\Pr(M_0|\mathbf{x}^{(n)})$ & $9.88\times10^{-4}$​ &$2.60\times10^{-5}$&  0/100\\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.63 & 0.0749 & 70/100 \\
%\hline 
$\Pr(M_2|\mathbf{x}^{(n)})$ & 0.37 & 0.0745& 30/100 \\
\hline 
%\pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Frequency of identifying \\ the respective model} && \\
%\hline
\end{tabular}
\caption{Average model posterior probabilities, variance and frequency of true model for the Scenario 4 simulation exercise with $n=500$ and the loss-based prior.}
\label{TableResults-WLG-Conjugate-500}
\end{table}


\begin{table}[h]
\centering
\begin{tabular}{c|c|c|c}
\hline 
• & Mean posterior & Variance posterior &  \pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Freq. true model} \\
%\hline 
\hline 
%$\Pr(M_0)$ & 0.015 [$3.46\times10^{-18}$]&   \\ 
%%\hline 
%$\Pr(M_1)$ & 0.014 [$1.51\times10^{-14}$]&   \\ 
%%\hline 
%$\Pr(M_2)$ & 0.971 [$1.46\times10^{-14}$] &   \\
%hline
%$B_{01}$ & 36.55 & $3.24\times10^{-4}$ & $4.65\times10^{-40}$ \\ 
%\hline 
%$B_{02}$ & $1.84\times 10^3$ & 0.02 & $1.27\times10^{-45}$ \\ 
%\hline 
%$B_{12}$ & 50.44 & 55 & $2.72\times10^{-6}$ \\
%\hline
$\Pr(M_0|\mathbf{x}^{(n)})$ & $1.33\times10^{-13}​$ &$1.76\times10^{-24}$&  0/100\\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.08 &0.0200  & 4/100 \\
%\hline 
$\Pr(M_2|\mathbf{x}^{(n)})$ & 0.92 &0.0200 & 96/100 \\
\hline 
%\pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Frequency of identifying \\ the respective model} && \\
%\hline
\end{tabular}
\caption{Average model posterior probabilities, variance and frequency of true model for the Scenario 4 simulation exercise with $n=1500$ and the loss-based prior. }
\label{TableResults-WLG-Conjugate-1500}
\end{table}

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 \ref{Unif_TableResults-WLG-Conjugate-500} and in Table \ref{Unif_TableResults-WLG-Conjugate-1500}, 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.

\begin{table}[H]
\centering
\begin{tabular}{c|c|c|c}
\hline 
• & Mean posterior & Variance posterior &  \pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Freq. true model} \\
%\hline 
\hline 
%hline
%$B_{01}$ & 36.55 & $3.24\times10^{-4}$ & $4.65\times10^{-40}$ \\ 
%\hline 
%$B_{02}$ & $1.84\times 10^3$ & 0.02 & $1.27\times10^{-45}$ \\ 
%\hline 
%$B_{12}$ & 50.44 & 55 & $2.72\times10^{-6}$ \\
%\hline
$\Pr(M_0|\mathbf{x}^{(n)})$ & $16\times 10^{-4}$& $7.15\times10^{-5}$ &  0/100\\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.82 & 0.0447 & 91/100 \\
%\hline 
$\Pr(M_2|\mathbf{x}^{(n)})$ & 0.18 & 0.0443 & 9/100 \\
\hline 
%\pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Frequency of identifying \\ the respective model} && \\
%\hline
\end{tabular}
\caption{Average model posterior probabilities, variance and frequency of true model for the Scenario 4 simulation exercise with $n=500$ and the uniform prior.}
\label{Unif_TableResults-WLG-Conjugate-500}
\end{table}

\begin{table}[h]
\centering
\begin{tabular}{c|c|c|c}
\hline 
• & Mean posterior & Variance posterior &  \pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Freq. true model} \\
%\hline 
\hline 
%$\Pr(M_0)$ & 0.3333 &   \\ 
%%\hline 
%$\Pr(M_1)$ & 0.3333 &   \\ 
%%\hline 
%$\Pr(M_2)$ & 0.3333 &   \\
%hline
%$B_{01}$ & 36.55 & $3.24\times10^{-4}$ & $4.65\times10^{-40}$ \\ 
%\hline 
%$B_{02}$ & $1.84\times 10^3$ & 0.02 & $1.27\times10^{-45}$ \\ 
%\hline 
%$B_{12}$ & 50.44 & 55 & $2.72\times10^{-6}$ \\
%\hline
$\Pr(M_0|\mathbf{x}^{(n)})$ & $8.64\times10^{-12}$ & $7.45\times10^{-21}$&  0/100\\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.501 &0.1356  & 49/100 \\
%\hline 
$\Pr(M_2|\mathbf{x}^{(n)})$ & 0.499 & 0.1356 & 51/100 \\
\hline 
%\pbox{15cm}{\relax\ifvmode\centering\fi \rule{0pt}{1em} Frequency of identifying \\ the respective model} && \\
%\hline
\end{tabular}
\caption{Average model posterior probabilities, variance and frequency of true model for the Scenario 4 simulation exercise with $n=1500$ and the uniform prior.}
\label{Unif_TableResults-WLG-Conjugate-1500}
\end{table}

\begin{figure}[H]
\centering
\includegraphics[scale=0.7]{WLG.pdf}
\caption{The densities of Weibull($\lambda,\kappa$), Log-normal($\mu,\tau$) and  Gamma($\alpha,\beta$) with the same mean (equal to 5) and the same variance (equal to 2.5).}
%\caption{The density functions of Weibull($\lambda,\kappa$), Log-normal($\mu,\tau$) and  Gamma($\alpha,\beta$) which were plotted using the \texttt{\textbf{R}} function \texttt{density()} where the bandwidth was set to 40 times the default bandwidth}
\label{fig_PDFs_WLG}
\end{figure}

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 \ref{fig_PDFs_WLG}, 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 \ref{Hellinger_TableResults}. In other words, from Figure \ref{fig_PDFs_WLG} 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 densities, because in this case only we would have a sufficient representation of the whole distribution.

\begin{table}[H]
\centering
\begin{tabular}{c|ccc}
\hline 
• & \multicolumn{3}{c}{Hellinger distances} \\ 
%\hline 
• & Weibull($\lambda,\kappa$) & Log-normal($\mu,\tau$) & Gamma($\alpha,\beta$) \\ 
\hline 
Weibull($\lambda,\kappa$) &  & 0.1411996  & 0.09718282 \\ 
%\hline 
Log-normal($\mu,\tau$) &  &  & 0.04899711 \\ 
\hline 
%Gamma($\alpha,\beta$) &  &  & 0.34 \\
%\hline 
\end{tabular}
\caption{Hellinger distances between all the pairs formed from a Weibull($\lambda,\kappa$), Log-normal($\mu,\tau$) and  Gamma($\alpha,\beta$). The six hyperparameters are such that the distributions have the same mean=5 and same variance=2.5.}
\label{Hellinger_TableResults}
\end{table}

\section{Change Point Analysis on Real Data}\label{sc_realdata}
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 \citep{Carlin}. 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 \ref{sub_sc_British_coal_mine}, while the latter in Section \ref{sub_sc_Financial_data}.


\subsection{British Coal-Mining Disaster Data}
\label{sub_sc_British_coal_mine}
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 \ref{British_coalmine} shows the number of deaths per year in the British coal-mining industry from 1851 to 1962.
\begin{figure}[h]
\centering   
\includegraphics[scale=0.6]{Coal_model.pdf}
\caption{Scatter plot of the British coal-mining disaster data.}
\label{British_coalmine}
\end{figure}
As in \citet{Chib}, we assume a Poisson sampling distribution with a possible change in the parameter value. That is
\begin{align}\label{m1_mining}
X_{1}, X_{2}, \ldots, X_{m}|\phi_{1} &\myiid \mbox{Poisson}(\phi_{1})\nonumber \\
X_{m+1}, X_{m+2}, \ldots, X_{n}|\phi_{2} &\myiid \mbox{Poisson}(\phi_{2}),
\end{align}
where $m$ is the unknown location of the single change point, such that $1\leq m \leq n$, and a $\mbox{Gamma}(2,1)$ is assumed for $\phi_1$ and $\phi_2$. The case $m=n$ corresponds to the scenario with no change point, that is model $M_0$. The case $m<n$ assumes one change point, that is model $M_1$.

Let $f_1(\cdot|\phi_1)$ and $f_2(\cdot|\phi_2)$ be the Poisson distributions with parameters $\phi_1$ and $\phi_2$, respectively. Then, the analysis is performed by selecting between model $M_0$, that is when the sampling distribution is $f_1$, and model $M_1$, where the sampling distribution is $f_1$ up to a certain $m<n$ and $f_2$ from $m+1$ to $n$.

As highlighted in the Remark at the end of Section \ref{sc_numberCP}, the prior on the model space is the discrete uniform distribution, that is $\Pr(M_0)=\Pr(M_1)=0.5$. The proposed model selection approach leads to the Bayes factors $B_{01}=1.61 \times 10^{-13}$ and $B_{10}=6.20 \times 10^{12}$, where it is obvious that the odds are strongly in favour of model $M_1$. Indeed, we have $\Pr(M_1|{\bf x}^{(n)})\approx1$.

% Similarly, the infinimum of the Kullback--Leibler divergence between $f_2(\cdot|\phi_2)$ and $f_1(\cdot|\phi_1)$ is attained when $\phi_1=0$. As such, when we apply Theorem \ref{remark_model} and equation \eqref{ModelPriors}, we obtain $\Pr(M_0)=\Pr(M_1)=0.5$. 
%\begin{table}[H]
%\centering
%\resizebox{\columnwidth}{!}{%
%\definecolor{LightCyan}{rgb}{0.88,1,1}
%
%
%    \begin{tabular}{|c|c|c|c\|c|c|c|c|c|}
%        \cline{2-9}
%        \multicolumn{1}{c|}{} & \backslashbox{Results}{No. of Data Points} & Full Data   & $n=20$ & $n=40$ & $n=41$ & $n=42$ & $n=43$ & $n=44$ \\ \hline
%   \multirow{2}{*}{\parbox{3.5cm}{\center{Model posterior probabilities}}}    & \textcolor{red}{$\Pr(M_{0}|\mathbf{X}^{(n)})$} & \textcolor{red}{0} & \textbf{\textcolor{red}{0.73}} & \textbf{\textcolor{red}{ 0.72}}  & \textbf{\textcolor{red}{0.71}} & \textbf{\textcolor{red}{0.62}} & \textbf{\textcolor{red}{0.50}} & \textcolor{red}{0.37}  \\ \hhline{~--------}
%  & \textcolor{red}{$\Pr(M_{1}|\mathbf{X}^{(n)})$} & \textbf{\textcolor{red}{1}}  & \textcolor{red}{0.27} & \textcolor{red}{0.28} & \textcolor{red}{0.29} & \textcolor{red}{0.38} & \textcolor{red}{0.50} & \textbf{\textcolor{red}{0.63}} \\ \hline
%    \multirow{2}{*}{\parbox{3.5cm}{\center{Bayes Factors}}}    &   $B_{01}$ & $1.61 \times 10^{-13}$ &  2.75  &  2.60 & 2.40 & 1.64 &  1.01 & 0.59  \\ \hhline{~--------}
%      &   $B_{10}$ &  $6.20 \times 10^{12}$ &  0.36 & 0.38 & 0.42 & 0.61 & 0.99 & 1.70 \\ \hline
%    \end{tabular}%
%    }
%    \caption{Results for the model posterior probabilities and the Bayes Factors corresponding to the British coal-mining disaster data from \citet{Carlin}}% obtained via the crude Monte Carlo method with 100000 iterations in MATLAB R2015a }
%    \label{ChibTable}
%\end{table}
%
%As shown in Table \ref{ChibTable}, the model posterior probability for model $M_1$ is basically 1. That is, when we consider the full data set, the objective approach based on the Loss-based prior clearly identifies the scenario of one data point. This result is supported by the value of the Bayes factor $B_{10}$. Although we do not specifically estimate the change point location, we can obtain a rough result in a naive fashion. In particular, we perform the analysis on progressively large subsets of the original data. As illustrated in Table \ref{ChibTable}, we consider $n=20$, $n=40$ to $n=44$. We then see that the posterior probability for $M_1$ moves from a value below 1/2 to a value equal or larger than 1/2 when $n\geq43$. As such, it appears that the change point is in the region around the $43^{rd}$ or $44^{th}$ observation. Although a rigorous comparison cannot be made, we see that our naive result is not dissimilar to the $n=41$ provided by \citet{Chib}.

\subsection{Daily S\&P 500 Absolute Log-Return Data}
\label{sub_sc_Financial_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 \ref{fig_sp500}). As underlying sampling distributions we consider the Weibull and the Log-normal \citep{Yu}, and the models among which we select are as follows. $M_0$ is a $\mbox{Weibull}(\lambda,\kappa)$, $M_1$ is formed by a $\mbox{Weibull}(\lambda,\kappa)$ and a $\mbox{Log-normal}(\mu_1,\tau_1)$ and, finally, $M_2$ is formed by a $\mbox{Weibull}(\lambda,\kappa)$, a $\mbox{Log-normal}(\mu_1,\tau_1)$ and a $\mbox{Log-normal}(\mu_2,\tau_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 \cite{KassRaftery}, due to the large sample size of the data set, we could approximate 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.
\begin{figure}[h]
\centering
\includegraphics[scale=0.6]{Financial_model_2.pdf}
\caption{Absolute daily log-returns of the S\&P500 index from 14/01/08 to 30/12/11.}
\label{fig_sp500}
\end{figure}
From the results in Table \ref{Finance_TableResults}, we see that the model indicated by the proposed approach is $M_2$. In other words, there is very strong indication that there are two change points in the data set.
\begin{table}[h]
\centering
\begin{tabular}{c|ccc}
%\hline 
%\hline
\hline 
$\Pr(M_0)$ & 0.36 \\ 
%\hline 
$\Pr(M_1)$ & 0.32 \\ 
%\hline 
$\Pr(M_2)$ & 0.32 \\
%hline
$B_{01}$ & $7.72\times10^{18}$ \\ 
%\hline 
$B_{02}$ & $3.30\times10^{-3}$ \\ 
%\hline 
$B_{12}$ & $4.28\times10^{-22}$ \\
%\hline
$\Pr(M_0|\mathbf{x}^{(n)})$ & 0.00\\ 
%\hline 
$\Pr(M_1|\mathbf{x}^{(n)})$ & 0.00 \\
%\hline 
$\Pr(M_2|\mathbf{x}^{(n)})$ & 1.00 \\
\hline 
\end{tabular}
\caption{Model prior, Bayes factor and model posterior probabilities for the S\&P500 change point analysis.}
\label{Finance_TableResults}
\end{table}
From Table \ref{Finance_TableResults}, we note that the prior on model $M_1$ and $M_2$ 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. 
\section{Conclusion}\label{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.

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 \citep{VillaMass} to make inference. The resulting prior distribution, which is the discrete uniform distribution, it is not new in the literature \citep{Giron}, 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.

\section*{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.

\begin{thebibliography}{xx}

\bibitem[{Barry and Hartigan(1992)}]{Barry1992}
Barry, D., Hartigan, J.~A., 1992. {Product Partition Models for Change Point
  Problems}. {The Annals of Statistics} 20~(1), 260--279.

\bibitem[{Berger et~al.(2012)Berger, Bernardo, and Sun}]{BBS2012}
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.

\bibitem[{Berger and Pericchi(1996)}]{Intrinsic}
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.

\bibitem[{Berk(1966)}]{berk1966}
Berk, R.~H., 1966. {Limiting Behavior of Posterior Distributions when the Model
  is Incorrect}. The Annals of Mathematical Statistics 37~(1), 51--58.

\bibitem[{Carlin et~al.(1992)Carlin, Gelfand, and Smith}]{Carlin}
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.

\bibitem[{Chernoff and Zacks(1964)}]{Chernoff1964}
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.

\bibitem[{Chib(1998)}]{Chib}
Chib, S., 1998. {Estimation and Comparison of Multiple Change-point Models}.
  Journal of Econometrics 86~(2), 221--241.

\bibitem[{Fearnhead and Liu(2007)}]{Fearnhead2007}
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.

\bibitem[{Gir{\'{o}}n et~al.(2007)Gir{\'{o}}n, Moreno, and Casella}]{Giron}
Gir{\'{o}}n, F.~J., Moreno, E., Casella, G., 2007. {Objective Bayesian Analysis
  of Multiple Changepoints for Linear Models (with discussion)}. In: Bernardo,
  {\relax J.M}., Bayarri, {\relax M.J}., Berger, {\relax J.O}., Dawid, {\relax
  A.P}., Heckerman, D., Smith, {\relax A.F.M}., West, M. (Eds.), Bayesian
  Statistics 8. Oxford University Press, London, pp. 227–--252.

\bibitem[{Hannart and Naveau(2009)}]{Hannart2009}
Hannart, A., Naveau, P., 2009. {Bayesian Multiple Change Points and
  Segmentation: Application to Homogenization of Climatic Series}. Water
  Resources Research 45~(10), {W10444}.

\bibitem[{Heard and Turcotte(2017)}]{Heard2017}
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.

\bibitem[{Henderson and Matthews(1993)}]{Henderson1993}
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.

\bibitem[{Kass and Raftery(1995)}]{KassRaftery}
Kass, R.~E., Raftery, A.~E., 1995. {Bayes Factors}. Journal of the American
  Statistical Association 90~(430), 773--795.

\bibitem[{Ko et~al.(2015)Ko, Chong, and Ghosh}]{ko2015}
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.

\bibitem[{Koop and Potter(2009)}]{Koop}
Koop, G., Potter, S.~M., 2009. {Prior Elicitation in Multiple Change-Point
  Models}. International Economic Review 50~(3), 751--772.

\bibitem[{Kullback and Leibler(1951)}]{kullback1951}
Kullback, S., Leibler, R.~A., 1951. {On Information and Sufficiency}. The
  Annals of Mathematical Statistics 22~(1), 79--86.

\bibitem[{Loschi and Cruz(2005)}]{Loschi2005}
Loschi, {\relax R.H}., Cruz, {\relax F.R.B}., 2005. {Extension to the Product
  Partition Model: Computing the Probability of a Change}. Computational
  Statistics \& Data Analysis 48~(2), 255 -- 268.

\bibitem[{Moreno et~al.(2005)Moreno, Casella, and Garcia-Ferrer}]{Moreno2005}
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.

\bibitem[{Muliere and Scarsini(1985)}]{Muliere1985}
Muliere, P., Scarsini, M., 1985. {Change-point Problems: A Bayesian
  Nonparametric Approach}. Aplikace matematiky 30~(6), 397--402.

\bibitem[{Petrone and Raftery(1997)}]{Petrone1997}
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.

\bibitem[{Raftery and Akman(1986)}]{Raftery}
Raftery, A.~E., Akman, V.~E., 1986. {Bayesian Analysis of a Poisson Process
  with a Change-Point}. Biometrika 73~(1), 85--89.

\bibitem[{Schwaller and Robin(2017)}]{Schwaller2017}
Schwaller, L., Robin, S., 2017. {Exact Bayesian Inference for Off-line
  Change-point Detection in Tree-structured Graphical Models}. Statistics and
  Computing 27~(5), 1331--1345.

\bibitem[{Smith(1975)}]{Smith1975}
Smith, A. F.~M., 1975. {A Bayesian Approach to Inference about a Change-Point
  in a Sequence of Random Variables}. Biometrika 62~(2), 407--416.

\bibitem[{Stephens(1994)}]{Stephens1994}
Stephens, D.~A., 1994. {Bayesian Retrospective Multiple-Changepoint
  Identification}. Journal of the Royal Statistical Society. Series C (Applied
  Statistics) 43~(1), 159--178.

\bibitem[{Tian et~al.(2009)Tian, Ng, Li, and Tan}]{Tian2009}
Tian, G.-L., Ng, K.~W., Li, K.-C., Tan, M., 2009. {Non-iterative Sampling-based
  Bayesian Methods for Identifying Changepoints in the Sequence of Cases of
  Haemolytic Uraemic Syndrome}. Computational Statistics and Data Analysis
  53~(9), 3314--3323.

\bibitem[{Villa and Walker(2015{\natexlab{a}})}]{VillaMass}
Villa, C., Walker, S., 2015{\natexlab{a}}. {An Objective Approach to Prior Mass
  Functions for Discrete Parameter Spaces}. Journal of the American Statistical
  Association 110~(511), 1072--1082.

\bibitem[{Villa and Walker(2015{\natexlab{b}})}]{VillaModel}
Villa, C., Walker, S., 2015{\natexlab{b}}. {An Objective Bayesian Criterion to
  Determine Model Prior Probabilities}. Scandinavian Journal of Statistics
  42~(4), 947--966.

\bibitem[{Yu(2001)}]{Yu}
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.

\end{thebibliography}


%\begin{thebibliography}{xx}
%\end{thebibliography}


%Building models such that the number of change points and their locations are appropriately considered, we have seen that our approach based on the method of \citet{VillaModel} is giving satisfactory results in both simulated and real datasets. Moreover, our model selection procedure could be augmented to also identify from an exhaustive list of generating sampling distributions, the ones the observations are generated by, besides detecting the number of changes in the sample. One area we plan to additionally investigate involves the inclusion of change point location estimation in our Loss-based model selection method.

%In conclusion, we have shown that the discrete uniform distribution is the objective prior which should be considered for the case of a single change point location. Furthermore, we have presented an information-theoretical rationale behind the use of the discrete uniform as an objective prior by \citet{Giron} for the multivariate change point problem. These developments allowed us to work on the model selection problem seen from the context of change points. Building models such that the number of change points and their locations are appropriately considered, we have seen that our approach based on the method of \citet{VillaModel} is giving satisfactory results in both simulated and real datasets. Moreover, our model selection procedure could be augmented to also identify the underlying sampling distributions the observations are generated by, besides detecting the number of changes in the sample. One area we plan to additionally investigate involves the inclusion of change point location estimation in our Loss-based model selection method.

%In conclusion, we have shown that the discrete uniform distribution is the objective prior which should be considered for the case of a single change point location. Furthermore, we have conjectured that the same uniform distribution is the prior once we deal with the multivariate change point problem. These developments allowed us to work on the model selection problem seen from the context of change points. Building models such that the number of change points and their locations are appropriately considered, we have seen that our approach based on the method of \citet{VillaModel} is giving satisfactory results in both simulated and real datasets. Moreover, our model selection procedure could be augmented to also identify the underlying sampling distributions the observations are generated by, besides detecting the number of changes in the sample. One area we plan to additionally investigate involves the inclusion of change point location estimation in our Loss-based model selection method.

%\section{Acknowledgement}\label{sc_ackn}
\newpage

%\begin{comment}
%\section*{Research Questions}
%In this paper, we mainly consider two problems:
%\begin{itemize}
%\item[1)]We have $n$ data points, denoted by $X_{1}, X_{2}, \ldots, X_{n}$ and we know that there are $k$ change points (a priori), such that $0<k\leq n$. Let $1\leq m_{1}<\ldots<m_{j}<\ldots<m_{k}\leq n$, so that we have:
%\begin{align*}
%X_{1}, X_{2}, \ldots, X_{m_{1}}&\sim f_{1}(\cdot|\theta_{1})\\
%X_{m_{1}+1}, X_{m_{1}+2}, \ldots, X_{m_{2}}&\sim f_{2}(\cdot|\theta_{2})\\
%\vdots\\
%X_{m_{k}+1}, X_{m_{k}+2}, \ldots, X_{n}&\sim f_{k+1}(\cdot|\theta_{k+1})
%\end{align*}
%Then, we are interested into fixing an objective prior on the positions of the $k$ change points, that is on $m_{1}, m_{2}, \ldots, m_{k}$.
%
%\item[2)] Now, supposing the $k$ is unknown and building the models as seen in Figure \ref{ModelConstruct}, we want to assign a prior probability to the models $\\ M_{0}, M_{1}, \ldots, M_{k}$ in a objective way given $f_{1}, f_{2}, \ldots, f_{k+1}$ and the priors on the parameters.
%
%\begin{figure}[H]
%%\hspace*{-5cm}
%\centering   
%%\vspace{-15em}
%\includegraphics[scale=0.7]{ModelConstructTheta.pdf}
%\vspace{-25em}
%\caption{Diagram showing the way we construct the models under consideration}
%%\vspace{1.5em}
%\label{ModelConstruct}
%\end{figure}
%
%\end{itemize}
%
%\section*{Solution to problem 1}
%By nature, the change position is discrete and this calls for using the Villa-Walker prior. Before we embark into proving some results we need a bit of notation.
%
%
%Let $\mathbf{X}$ denote the \textbf{n}-size vector of random variables, $\mathbf{X^{(n)}}=(X_{1},\ldots, X_{n})$ and 
%\begin{align*}
%\mathbf{m^{(k,n)}}&=(m_{1}, m_{2}, \ldots, m_{k})\\
%\mathbf{\boldsymbol\theta^{(k+1,n)}}&=(\theta_{1}, \theta_{2}, \ldots, \theta_{k}, \theta_{k+1}),
%\end{align*}
%such that $m_{1}, \ldots, m_{k}\in\mathbb{N}^{*}$ with $1<m_{1}<\ldots<m_{k}<\mathbf{n}$ and $\theta_{1}, \ldots, \theta_{k+1}$ are the appropriate statistical parameters. Additionally, let $l\in\mathbb{N^{*}}$.
%
%\begin{myremark}
%In this paper, we use the notation $\mathbf{n}$ and $n$ interchangeably.
%\end{myremark}
%
%Furthermore, let there be $a_{1}, \ldots, a_{k}\in \mathbb{Z}$, such that 
%\begin{align*}
%\mathbf{a^{(k)}}&=(a_{1}, a_{2}, \ldots, a_{k})\\
%\mathbf{m^{(k,n)}+a^{(k)}}&=(m_{1}+a_{1}, m_{2}+a_{2}, \ldots, m_{k}+a_{k}),
%\end{align*}
%where $m_{1}+a_{1}<m_{2}+a_{2}<\ldots<m_{k}+a_{k}$.
%
%Moreover, let us further define:
%\begin{align*}
%\mathbf{a_{j}^{(k)}}&=(a_{1}, a_{2}, \ldots, a_{j}, \ldots, a_{k}) \textrm{ such that } a_{j}=1, a_{i}=0, \forall i\neq j \in\{1, 2, \ldots, k\}
%\end{align*}
%
%Now, let us bring forward the likelihood function of our $k$ change points and $k+1$ statistical parameters, which is denoted with
%\begin{align}
%L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})&=\left( \prod_{i=1}^{m_{1}} f(x_{i}|\theta_{1})\right)\cdot \left( \prod_{i=m_{1}+1}^{m_{2}} f(x_{i}|\theta_{2})\right) \cdot \ldots\nonumber \\ & \ldots \cdot \left( \prod_{i=m_{k}+1}^{n} f(x_{i}|\theta_{k+1})\right).
%\label{likelihood}
%\end{align}
%
%Also, let us use the following notations for the Kullback--Leibler divergence corresponding to a specific change point $m_{j}, \forall j \in \{1, 2, \ldots, k\}$
%\begin{align*}
%d_{+}(\theta_{j}, \theta_{j+1})&=D_{KL}(f(\cdot|\theta_{j+1})\|f(\cdot|\theta_{j}))\\
%d_{-}(\theta_{j}, \theta_{j+1})&=D_{KL}(f(\cdot|\theta_{j})\|f(\cdot|\theta_{j+1}))
%\end{align*}
%
%Now, we denote with
%\begin{align}
%D(\mathbf{m^{(k,n)}}\|\mathbf{m^{(k,n)}\pm \textit{l} \cdot a_{j}^{(k)}}|&\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\
%&\cdot\log\left(\dfrac{L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}{L(\mathbf{m^{(k,n)}\pm \textit{l} \cdot a_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}\right)\,\mathrm{d}\mathbf{X^{(n)}}
%\label{generaldivergence}
%\end{align}
%
%\begin{mylemma}
%$\\D(\mathbf{m^{(k,n)}}\|\mathbf{m^{(k,n)} \pm \textit{l} \cdot a_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})\geq l\cdot d_{\pm}(\theta_{j},\theta_{j+1}) \\ \forall j \in \{1, 2, \ldots, k\}, \forall k<n, \forall n\in \mathbb{N}^{*}$.
%\end{mylemma}
%
%\begin{proof}
%
%Let $\mathbf{a_{j}^{(k)}}$ be an arbitrary, but fixed vector, such that element $a_{j}$ is 1, whereas the other constituent elements are zeroes, that is $ \mathbf{a_{j}^{(k)}}=(0,0, \ldots,\underbrace{1}_{\mathclap{\text{\textit{j}-th element}}}, \ldots,0)$.
%
%We distinguish two cases:
%\begin{itemize}
%\item when the sign is `$+$':
%\begin{align}
%&D(\mathbf{m^{(k,n)}}\|\mathbf{m^{(k,n)}+\textit{l}\cdot a_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})=\nonumber\\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\
%&\cdot\log\left(\dfrac{L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}{L(\mathbf{m^{(k,n)}+\textit{l}\cdot a_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}\right)\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot\left\lbrace\rule{0cm}{1cm}\left[\sum_{i=1}^{m_{1}} \log(f(x_{i}|\theta_{1}))\right]\right.+\nonumber \\
%&+\left.\left[\sum_{i=m_{1}+1}^{m_{2}} \log(f(x_{i}|\theta_{2}))\right]+\ldots+\left[\sum_{i=m_{j-1}+1}^{m_{j}} \log(f(x_{i}|\theta_{j}))\right]+\right.\nonumber\\ 
%&\left.+\left[\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \log(f(x_{i}|\theta_{j+1}))\right]+\left[\sum_{i=m_{j}+\textit{l}+1}^{m_{j+1}} \log(f(x_{i}|\theta_{j+1}))\right]+\right.\nonumber\\
%&+\left.\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \log(f(x_{i}|\theta_{j+2}))\right]+\ldots+\left[\sum_{i=m_{k}+1}^{n} \log(f(x_{i}|\theta_{k+1}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=1}^{m_{1}} \log(f(x_{i}|\theta_{1}))\right]-\left[\sum_{i=m_{1}+1}^{m_{2}} \log(f(x_{i}|\theta_{2}))\right]-\ldots-\right.\nonumber\\
%&\left.-\left[\sum_{i=m_{j-1}+1}^{m_{j}} \log(f(x_{i}|\theta_{j}))\right]-\left[\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \log(f(x_{i}|\theta_{j}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=m_{j}+\textit{l}+1}^{m_{j+1}} \log(f(x_{i}|\theta_{j+1}))\right]-\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \log(f(x_{i}|\theta_{j+2}))\right]-\ldots-\right.\nonumber\\
%&-\left.\left[\sum_{i=m_{k}+1}^{n} \log(f(x_{i}|\theta_{k+1}))\right]\rule{0cm}{1cm}\right\rbrace\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\ 
%\displaybreak \nonumber\\
%&=\int  L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\ &\left[\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \log\left(\dfrac{f(x_{i}|\theta_{j+1})}{f(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber  \\ &\cdot \left[\log\left(\dfrac{f(x_{i}|\theta_{j+1})}{f(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \Bigg\lbrace \left\lbrace \int f(x_{i}|\theta_{j+1}) \cdot  \left[\log\left(\dfrac{f(x_{i}|\theta_{j+1})}{f(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \cdot \nonumber \\
%& \idotsint f(x_{1}|\theta_{1}) \cdot \ldots \cdot f(x_{m_{j}}|\theta_{j}) \cdot \ldots \cdot f(x_{i-1}|\theta_{j+1}) \cdot \nonumber\\& \cdot f(x_{i+1}|\theta_{j+1}) \cdot \ldots \cdot f(x_{n}|\theta_{k+1})\,\mathrm{d}x_{1}\ldots\,\mathrm{d}x_{m_{j}}\ldots\,\mathrm{d}x_{i-1} \nonumber \\ &\,\mathrm{d}x_{i+1}\ldots\,\mathrm{d}x_{n}\Bigg\rbrace \nonumber \\
%&=\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \left\lbrace 1^{n-1} \cdot \int f(x_{i}|\theta_{j+1}) \cdot \left[\log\left(\dfrac{f(x_{i}|\theta_{j+1})}{f(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} D_{KL}(f(x_{i}|\theta_{j+1})\|f(x_{i}|\theta_{j})) \nonumber \\
%&=\textit{l} \cdot D_{KL}(f(\cdot|\theta_{j+1})\|f(\cdot|\theta_{j})) \nonumber \\
%&=\textit{l} \cdot d_{+}(\theta_{j}, \theta_{j+1})
%\label{demopositive}
%\end{align}
%
%\item when the sign is `$-$':
%\begin{align}
%&D(\mathbf{m^{(k,n)}}\|\mathbf{m^{(k,n)}-\textit{l}\cdot a_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})=\nonumber\\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\
%&\cdot\log\left(\dfrac{L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}{L(\mathbf{m^{(k,n)}-\textit{l}\cdot a_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}\right)\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot\left\lbrace\rule{0cm}{1cm}\left[\sum_{i=1}^{m_{1}} \log(f(x_{i}|\theta_{1}))\right]\right.+\nonumber \\
%&+\left.\left[\sum_{i=m_{1}+1}^{m_{2}} \log(f(x_{i}|\theta_{2}))\right]+\ldots+\left[\sum_{i=m_{j-1}+1}^{m_{j}-\textit{l}} \log(f(x_{i}|\theta_{j}))\right]+\right.\nonumber\\ 
%&\left.+\left[\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \log(f(x_{i}|\theta_{j}))\right]+\left[\sum_{i=m_{j}+1}^{m_{j+1}} \log(f(x_{i}|\theta_{j+1}))\right]+\right.\nonumber\\
%&+\left.\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \log(f(x_{i}|\theta_{j+2}))\right]+\ldots+\left[\sum_{i=m_{k}+1}^{n} \log(f(x_{i}|\theta_{k+1}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=1}^{m_{1}} \log(f(x_{i}|\theta_{1}))\right]-\left[\sum_{i=m_{1}+1}^{m_{2}} \log(f(x_{i}|\theta_{2}))\right]-\ldots-\right.\nonumber\\
%&\left.-\left[\sum_{i=m_{j-1}+1}^{m_{j}-\textit{l}} \log(f(x_{i}|\theta_{j}))\right]-\left[\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \log(f(x_{i}|\theta_{j+1}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=m_{j}+1}^{m_{j+1}} \log(f(x_{i}|\theta_{j+1}))\right]-\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \log(f(x_{i}|\theta_{j+2}))\right]-\ldots-\right.\nonumber\\
%&-\left.\left[\sum_{i=m_{k}+1}^{n} \log(f(x_{i}|\theta_{k+1}))\right]\rule{0cm}{1cm}\right\rbrace\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\ 
%\displaybreak \nonumber\\
%&=\int  L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\ &\left[\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \log\left(\dfrac{f(x_{i}|\theta_{j})}{f(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber  \\ &\cdot \left[\log\left(\dfrac{f(x_{i}|\theta_{j})}{f(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \Bigg\lbrace \left\lbrace \int f(x_{i}|\theta_{j}) \cdot  \left[\log\left(\dfrac{f(x_{i}|\theta_{j})}{f(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \cdot \nonumber \\
%& \idotsint f(x_{1}|\theta_{1}) \cdot \ldots \cdot f(x_{m_{j}-\textit{l}}|\theta_{j}) \cdot \ldots \cdot f(x_{i-1}|\theta_{j}) \cdot \nonumber\\& \cdot f(x_{i+1}|\theta_{j}) \cdot \ldots \cdot f(x_{n}|\theta_{k+1})\,\mathrm{d}x_{1}\ldots\,\mathrm{d}x_{m_{j}-\textit{l}}\ldots\,\mathrm{d}x_{i-1} \nonumber \\ &\,\mathrm{d}x_{i+1}\ldots\,\mathrm{d}x_{n}\Bigg\rbrace \nonumber \\
%&=\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \left\lbrace 1^{n-1} \cdot \int f(x_{i}|\theta_{j}) \cdot \left[\log\left(\dfrac{f(x_{i}|\theta_{j})}{f(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} D_{KL}(f(x_{i}|\theta_{j})\|f(x_{i}|\theta_{j+1})) \nonumber \\
%&=\textit{l} \cdot D_{KL}(f(\cdot|\theta_{j})\|f(\cdot|\theta_{j+1})) \nonumber \\
%&=\textit{l} \cdot d_{-}(\theta_{j}, \theta_{j+1})
%\label{demonegative}
%\end{align}
%
%\end{itemize}
%
%From expressions \eqref{demopositive} and \eqref{demonegative}, we see that the lemma holds.
%\end{proof}
%
%\pagebreak
%
%\section*{Solution to problem 2}
%
%Let $\mathbf{X}$ denote the \textbf{n}-size vector of random variables, $\mathbf{X^{(n)}}=(X_{1},\ldots, X_{n})$ and 
%\begin{align*}
%\mathbf{m^{(k,n)}}&=(m_{1}, m_{2}, \ldots,m_{j}, \ldots, m_{k})\\
%\mathbf{\boldsymbol\theta^{(k+1,n)}}&=(\theta_{1}, \theta_{2}, \ldots, \theta_{j}, \ldots, \theta_{k}, \theta_{k+1}),
%\end{align*}
%such that $m_{1}, \ldots,m_{j}, \ldots, m_{k}\in\mathbb{N}^{*}$ with $1\leq m_{1}<\ldots<m_{j}<\ldots<m_{k}<\mathbf{n}$ and $\theta_{1}, \ldots, \theta_{j}, \ldots, \theta_{k+1}$ are the appropriate statistical parameters.
%
%\begin{myremark}
%In this paper, we use the notation $\mathbf{n}$ and $n$ interchangeably.
%\end{myremark}
%
%Let us now construct the models according with the notation seen in \cite[p.~950]{VillaModel} and Figure \ref{ModelConstruct}. The model denoted with $M_{1}$ corresponds to a single change point. Now we expand upon this model, by fixing the position of the change point (denoted with $m_{1}$) and adding a further single change point whose position is to the right of $m_{1}$.  This new model is $M_{2}$. In a similar manner, we construct $M_{3}$, $M_{4}$, \ldots, $M_{k}$. Additionally, let us denote with $M_{0}$, the model with no change points. 
%
%Moreover, let us employ the following notation which we will be using
%for our model parameters:
%\begin{equation*}
%\mathbf{\tilde{\boldsymbol\theta}_{j}}=(\mathbf{m^{(j,n)}}, \mathbf{\boldsymbol\theta^{(j+1,n)}}) \mathrm{ ,with}\: j=1,\ldots, k
%\end{equation*}
%
%Now, let us bring forward the likelihood function of our $k$ change points and $k+1$ statistical parameters, which is denoted with
%\begin{align}
%L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})&=\left( \prod_{i=1}^{m_{1}} f_{1}(x_{i}|\theta_{1})\right)\cdot \left( \prod_{i=m_{1}+1}^{m_{2}} f_{2}(x_{i}|\theta_{2})\right) \cdot \ldots\nonumber \\ & \ldots \cdot \left( \prod_{i=m_{k}+1}^{n} f_{k+1}(x_{i}|\theta_{k+1})\right).
%\label{likelihood}
%\end{align}
%
%
%
%Using the notation introduced previously, as well as the way we denote likelihood corresponding to model $M_{j}$, that is $L(x|\mathbf{\tilde{\boldsymbol\theta}_{j}})=L(\mathbf{m^{(j,n)}}|\mathbf{\boldsymbol\theta^{(j+1,n)}}, \mathbf{X^{(n)}})$.
%
%\pagebreak
%
%\subsection*{Motivating Example}
%In this section, we will compare three models, namely $M_{0}$, $M_{1}$ and $M_{2}$. Furthermore, we will employ the same distributions as studied by \cite[p.~955]{VillaModel}, that is $f_{1}$ is the pdf of a $Weibull(\lambda,\kappa)$, $f_{2}$ is the pdf of a $Log-normal(\mu,\tau)$ and $f_{3}$ is the pdf of a $Gamma(\alpha,\beta)$, with the following parametrizations:
%\begin{align*}
%f_{1}(x|\lambda, \kappa)&=\dfrac{\kappa}{\lambda}\cdot\left(\dfrac{x}{\lambda}\right)^{\kappa-1}\cdot e^{-\left(\dfrac{x}{\lambda}\right)^{\kappa}}\\
%f_{2}(x|\mu, \tau)&=\dfrac{1}{x}\cdot\left(\dfrac{\tau}{2\pi}\right)^{1/2}\cdot e^{-\dfrac{\tau}{2}\cdot(\log(x)-\mu)^2}\\
%f_{3}(x|\alpha, \beta)&=\dfrac{\beta^{\alpha}}{\Gamma(\alpha)}\cdot x^{\alpha-1}\cdot e^{-\beta\cdot x}
%\end{align*}
%
%Following the way we have worked in the context of the first problem and keeping in mind that $\theta_{1}=(\lambda, \kappa)$, $\theta_{2}=(\mu, \tau)$ and $\theta_{3}=(\alpha, \beta)$, we have that:
%\begin{itemize}
%\item the prior probability for model $M_{0}$ is connected to the following quantities:
%\begin{align*}
% D_{KL}(\prod_{i=1}^{n}f_{1}(x|\theta_{1})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}}))&=\ldots=(n-m_{1})\cdot D_{KL}(f_{1}(\cdot|\theta_{1})\|f_{2}(\cdot|\theta_{2}))\\
% D_{KL}(\prod_{i=1}^{n}f_{1}(x|\theta_{1})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}}))&=\ldots=(m_{2}-m_{1})\cdot D_{KL}(f_{1}(\cdot|\theta_{1})\|f_{2}(\cdot|\theta_{2}))+\\
%&+(n-m_{2})\cdot D_{KL}(f_{1}(\cdot|\theta_{1})\|f_{3}(\cdot|\theta_{3}))
%\end{align*}
%
%\item the prior probability for model $M_{1}$ is connected to the following quantities:
%\begin{align*}
% D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}}))&=\ldots=(n-m_{2})\cdot D_{KL}(f_{2}(\cdot|\theta_{2})\|f_{3}(\cdot|\theta_{3}))\\
% D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}})\|\prod_{i=1}^{n}f_{1}(x|\theta_{1}))&=\ldots=(n-m_{1})\cdot D_{KL}(f_{2}(\cdot|\theta_{2})\|f_{1}(\cdot|\theta_{1}))
%\end{align*}
%
%\item the prior probability for model $M_{2}$ is connected to the following quantities:
%\begin{align*}
% D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}}))&=\ldots=(n-m_{2})\cdot D_{KL}(f_{3}(\cdot|\theta_{3})\|f_{2}(\cdot|\theta_{2}))\\
% D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}})\|\prod_{i=1}^{n}f_{1}(x|\theta_{1}))&=\ldots=(m_{2}-m_{1})\cdot D_{KL}(f_{2}(\cdot|\theta_{2})\|f_{1}(\cdot|\theta_{1}))+\\
%&+(n-m_{2})\cdot D_{KL}(f_{3}(\cdot|\theta_{3})\|f_{1}(\cdot|\theta_{1}))
%\end{align*}
%
%
%\end{itemize}
%
%For the sake of being concise, we will use the following notation:
%\begin{align*}
%D_{KL}(\prod_{i=1}^{n}f_{1}(x|\theta_{1})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}}))&=D_{KL}(M_{0}\|M_{1})\\
%D_{KL}(\prod_{i=1}^{n}f_{1}(x|\theta_{1})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}}))&=D_{KL}(M_{0}\|M_{2})\\
%D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}}))&=D_{KL}(M_{1}\|M_{2})\\
%D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}})\|\prod_{i=1}^{n}f_{1}(x|\theta_{1}))&=D_{KL}(M_{1}\|M_{0})\\
%D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{1}}))&=D_{KL}(M_{2}\|M_{1})\\
%D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{2}})\|\prod_{i=1}^{n}f_{1}(x|\theta_{1}))&=D_{KL}(M_{2}\|M_{0})\\
%\end{align*}
%
%\subsection*{Generalization}
%
%\begin{mylemma}
%For any $0\leq i<j \leq k, \forall k\leq n$, with all the models $M_{i}$ and $M_{j}$ constructed appropriately, we have the following:
%\begin{align*}
% &D_{KL}(M_{i}\|M_{j})=(m_{i+2}-m_{i+1})\cdot D_{KL}(f_{i+1}(\cdot|\theta_{i+1})\|f_{i+2}(\cdot|\theta_{i+2}))+ \\
%&+(m_{i+3}-m_{i+2})\cdot D_{KL}(f_{i+1}(\cdot|\theta_{i+1})\|f_{i+3}(\cdot|\theta_{i+3}))+\ldots+\\
%&+(n-m_{j})\cdot D_{KL}(f_{i+1}(\cdot|\theta_{i+1})\|f_{j+1}(\cdot|\theta_{j+1}))
%\end{align*}
%and
%\begin{align*}
% &D_{KL}(M_{j}\|M_{i})=(m_{i+2}-m_{i+1})\cdot D_{KL}(f_{i+2}(\cdot|\theta_{i+2})\|f_{i+1}(\cdot|\theta_{i+1}))+ \\
%&+(m_{i+3}-m_{i+2})\cdot D_{KL}(f_{i+3}(\cdot|\theta_{i+3})\|f_{i+1}(\cdot|\theta_{i+1}))+\ldots+\\
%&+(n-m_{j})\cdot D_{KL}(f_{j+1}(\cdot|\theta_{j+1})\|f_{i+1}(\cdot|\theta_{i+1}))
%\end{align*}
%
%For any $0\leq i<j \leq k, \forall k\leq n$, with all the models $M_{i}$ and $M_{j}$ constructed appropriately, we have the following:
%\begin{align*}
% &D_{KL}(M_{i}\|M_{j})=(m_{i+2}-m_{i+1})\cdot D_{KL}(f_{i+1}(\cdot|\theta_{i+1})\|f_{i+2}(\cdot|\theta_{i+2}))+ \\
%&+(m_{i+3}-m_{i+2})\cdot D_{KL}(f_{i+1}(\cdot|\theta_{i+1})\|f_{i+3}(\cdot|\theta_{i+3}))+\ldots+\\
%&+(m_{j}-m_{j-1})\cdot D_{KL}(f_{i+1}(\cdot|\theta_{i+1})\|f_{j}(\cdot|\theta_{j}))+\\
%&+(n-m_{j})\cdot D_{KL}(f_{i+1}(\cdot|\theta_{i+1})\|f_{j+1}(\cdot|\theta_{j+1}))
%\end{align*}
%and
%\begin{align*}
% &D_{KL}(M_{j}\|M_{i})=(m_{i+2}-m_{i+1})\cdot D_{KL}(f_{i+2}(\cdot|\theta_{i+2})\|f_{i+1}(\cdot|\theta_{i+1}))+ \\
%&+(m_{i+3}-m_{i+2})\cdot D_{KL}(f_{i+3}(\cdot|\theta_{i+3})\|f_{i+1}(\cdot|\theta_{i+1}))+\ldots+\\
%&+(m_{j}-m_{j-1})\cdot D_{KL}(f_{j}(\cdot|\theta_{j})\|f_{i+1}(\cdot|\theta_{i+1}))+\\
%&+(n-m_{j})\cdot D_{KL}(f_{j+1}(\cdot|\theta_{j+1})\|f_{i+1}(\cdot|\theta_{i+1}))
%\end{align*}
%\end{mylemma}
%
%
%\begin{proof}
%Because $M_{j}$ is just an appropriately constructed extension to $M_{i}$ and taking into account how we can link the quantity of interest ($D_{KL}(M_{i}\|M_{j})$ or $D_{KL}(M_{j}\|M_{i})$) with the individual Kullback--Leibler divergences as seen in the motivating example, then the result of this lemma is quite straight-forward.
%\end{proof}
%
%Following the formula shown in \citet[p.~950]{VillaModel}, we have
%
%\begin{equation}
%\Pr(M_{j})\propto \exp\left\lbrace \min_{i\neq j}\left\langle E_{\pi_{j}}\left[\inf_{\mathbf{\tilde{\boldsymbol\theta}_{i}},i\neq j} D_{KL}(L(x|\mathbf{\tilde{\boldsymbol\theta}_{j}})\|L(x|\mathbf{\tilde{\boldsymbol\theta}_{i}}))\right]\right\rangle\right\rbrace,
%\label{VillaModel}
%\end{equation}
%
%with $j=1, \ldots, k$. Alternatively, the statement \ref{VillaModel} can be stated in a ore compact form as
%\begin{equation}
%\Pr(M_{j})\propto \exp\left\lbrace \min_{i\neq j}\left\langle E_{\pi_{j}}\left[\inf_{\mathbf{\tilde{\boldsymbol\theta}_{i}},i\neq j} D_{KL}(M_{j}\|M_{i})\right]\right\rangle\right\rbrace
%\label{VillaModelComapact}
%\end{equation}
%\pagebreak
%
%\section*{Simulation Results}
%In this section we will present the results of the simulations we have run in MATLAB 2015a. For the sake of brevity, we have decided to denote with $
%\mathbf{\underline{X}}$ our simulated data. Keeping in line with the ideas from section 3 of \citet[pp.~950--958]{VillaModel} and the motivating example regarding the solution to problem 2, we have selected the following underlying distributions for our models to simulate data from:
%\begin{align*}
%X_{1}, X_{2}, \ldots, X_{m_{1}}&\sim \textrm{\textit{Weibull}}(\lambda, \kappa)\\
%X_{m_{1}+1}, X_{m_{1}+2}, \ldots, X_{m_{2}} &\sim \textrm{\textit{Log-normal}}(\mu, \tau)\\
%X_{m_{2}+1}, X_{m_{2}+2}, \ldots, X_{n}&\sim \textrm{\textit{Gamma}}(\alpha, \beta),
%\end{align*}
%where $1\leq m_{1}<m_{2}\leq n-1$ represents the positions of the two change points. These distributions have the respective probability density functions:
%\begin{align*}
%f_{1}(x|\lambda, \kappa)&=\dfrac{\kappa}{\lambda}\cdot\left(\dfrac{x}{\lambda}\right)^{\kappa-1}\cdot e^{-\left(\dfrac{x}{\lambda}\right)^{\kappa}}\\
%f_{2}(x|\mu, \tau)&=\dfrac{1}{x}\cdot\left(\dfrac{\tau}{2\pi}\right)^{1/2}\cdot e^{-\dfrac{\tau}{2}\cdot(\log(x)-\mu)^2}\\
%f_{3}(x|\alpha, \beta)&=\dfrac{\beta^{\alpha}}{\Gamma(\alpha)}\cdot x^{\alpha-1}\cdot e^{-\beta\cdot x},
%\end{align*}
%where for the simulated data we fix the parameters to $\lambda=1.5, \kappa=5, \mu=0.05, \tau=16, \alpha=10, \beta=2$, accordingly. Furthermore, for the sake of consistency with the paper of \citet{VillaModel}, we take each of the six parameters of the distributions as independent of each other and we use the following priors on them:
%\begin{align*}
%\lambda &\sim \textrm{\textit{Gamma}}(\alpha_{\lambda}, \beta_{\lambda}),\textrm{where }\alpha_{\lambda}=1.5,  \beta_{\lambda}=1\\
%\kappa &\sim \textrm{\textit{Gamma}}(\alpha_{\kappa}, \beta_{\kappa}),\textrm{where }\alpha_{\kappa}=5,  \beta_{\kappa}=1\\
%\mu &\sim \textrm{\textit{Normal}}(\mu_{1}, \sigma_{1}^{2}),\textrm{where }\mu_{1}=0.05,  \sigma_{1}^{2}=1\\
%\tau &\sim \textrm{\textit{Gamma}}(\alpha_{\tau}, \beta_{\tau}),\textrm{where }\alpha_{\tau}=16,  \beta_{\tau}=1\\
%\alpha &\sim \textrm{\textit{Gamma}}(\alpha_{1}, \beta_{1}),\textrm{where }\alpha_{1}=10,  \beta_{1}=1\\
%\beta &\sim \textrm{\textit{Gamma}}(\alpha_{2}, \beta_{2}),\textrm{where }\alpha_{2}=0.2,  \beta_{2}=0.1
%\end{align*}
%
%We have selected the twelve aforementioned hyperparameters($\alpha_{\lambda}, \beta_{\lambda}, \alpha_{\kappa}, \beta_{\kappa},\\ \mu_{1}, \sigma_{1}^{2}, \alpha_{\tau}, \beta_{\tau}, \alpha_{1}, \beta_{1}, \alpha_{1}$ and $\beta_{1}$) in such a way that the prior on each of our three models (the no change point, the one change point and the two change points models) are reasonably close to one another, which will be much more clearer from the results of Table \ref{TableResults}.
%
%The next figures represent the simulated data, where Figure \ref{No_change} shows data which contains no change points, whereas in Figure \ref{One_change} and Figure \ref{Two_change} we have data with one and two change points, respectively. 
%\begin{figure}[H]
%%\hspace*{-5cm}
%\centering   
%%\vspace{-15em}
%\includegraphics[scale=0.6]{No_change_shuffle.pdf}
%%\vspace{-25em}
%\caption{Here the data was simulated from a \textit{Weibull}(1.5, 5)}
%%\vspace{1.5em}
%\label{No_change}
%\end{figure}
%
%\begin{figure}[H]
%%\hspace*{-5cm}
%\centering   
%%\vspace{-15em}
%\includegraphics[scale=0.6]{One_change_shuffle.pdf}
%%\vspace{-25em}
%\caption{Here the data was simulated from a \textit{Weibull}(1.5, 5) until position 50 and afterwards we simulated from a \textit{Log-normal}(0.05, 16)}
%%\vspace{1.5em}
%\label{One_change}
%\end{figure}
%
%\begin{figure}[H]
%%\hspace*{-5cm}
%\centering   
%%\vspace{-15em}
%\includegraphics[scale=0.6]{Two_change_shuffle.pdf}
%%\vspace{-25em}
%\caption{Here the data was simulated from a \textit{Weibull}(1.5, 5) until position 50 and afterwards we simulated from a \textit{Log-normal}(0.05, 16) until position 70, while after position 70 we simulated from a \textit{Gamma}(10, 2)}
%%\vspace{1.5em}
%\label{Two_change}
%\end{figure}
% 
%
%Therefore, when we simulate 100 data points from the underlying distributions, we obtain the following results:
%
%\begin{table}[H]
%\centering
%\definecolor{LightCyan}{rgb}{0.88,1,1}
%
%
%    \begin{tabular}{|c|c|c|c|}
%        \hline
%        \backslashbox{Results}{Model}        & \pbox{10cm}{No change point} & \pbox{10cm}{One change point \\ ($m_{1}=50$)} & \pbox{10cm}{Two change points \\ ($m_{1}=50$ and $m_{2}=70$)}\\ \hline
%     Unnormalised $\Pr(M_{0})$  & 1.0942  & 1.0942  & 1.0942 \\ \hline
%     Unnormalised $\Pr(M_{1})$ & 1.6036 & 1.6037 & 1.6045\\ \hline
%     Unnormalised $\Pr(M_{2})$  & 1.3738  & 1.3738  &  1.3738\\ \hline
%                  $\Pr(M_{0})$ & 0.2687 & 0.2687 & 0.2687 \\ \hline
%                  $\Pr(M_{1})$ & 0.3938 & 0.3939 & 0.3940 \\ \hline
%                  $\Pr(M_{2})$ & 0.3374 & 0.3374 & 0.3373 \\ \hline
%   \textcolor{red}{$\Pr(M_{0}|\mathbf{\underline{X}})$} & \textcolor{red}{0.9099} & \textcolor{red}{$9.8391 \times 10^{-4}$}  & \textcolor{red}{$1.2307 \times 10^{-38}$} \\ \hline
%   \textcolor{red}{$\Pr(M_{1}|\mathbf{\underline{X}})$} & \textcolor{red}{0.0893} & \textcolor{red}{0.9660} & \textcolor{red}{0.1826} \\ \hline
%   \textcolor{red}{$\Pr(M_{2}|\mathbf{\underline{X}})$} & \textcolor{red}{$7.7895 \times 10^{-4}$} & \textcolor{red}{0.0330} & \textcolor{red}{0.8174} \\ \hline
%        $B_{01}$ &  14.9366  & 0.0015 & $9.8825 \times 10^{-38}$\\ \hline
%        $B_{02}$ & 1.4667$\times 10^{3}$ & 0.0374  & $1.8903 \times 10^{-38}$ \\ \hline
%        $B_{12}$ &  98.1973  & 25.0733 & 0.1913 \\ \hline
%    \end{tabular}
%\caption{The results for the normalised posterior model probabilities are shown in red for a dataset of 100 points where we have used parallelised loops in a cluster of 10 threads to compute the integrals via the crude Monte Carlo method where the number of repetitions was set to 100000. Note that $B_{01}$ is the Bayes Factor between model $M_{0}$ (no change point) and model $M_{1}$ (one change point) and so on. }
%\label{TableResults}
%\end{table}
%
%
%\begin{table}[H]
%\centering
%\definecolor{LightCyan}{rgb}{0.88,1,1}
%
%
%    \begin{tabular}{|c|c|c|c|}
%        \hline
%        \backslashbox{Results}{Model}        & \pbox{10cm}{No change point} & \pbox{10cm}{One change point \\ ($m_{1}=50$)} & \pbox{10cm}{Two change points \\ ($m_{1}=50$ and $m_{2}=70$)}\\ \hline
%     Unnormalised $\Pr(M_{0})$  & 1.0942  & 1.0942  & 1.0942 \\ \hline
%     Unnormalised $\Pr(M_{1})$ & 1.6029 & 1.6036 & 1.6011\\ \hline
%     Unnormalised $\Pr(M_{2})$  & 1.3743  & 1.3732  &  1.3737\\ \hline
%                  $\Pr(M_{0})$ & 0.2687 & 0.2688 & 0.2689 \\ \hline
%                  $\Pr(M_{1})$ & 0.3937 & 0.3939 & 0.3935 \\ \hline
%                  $\Pr(M_{2})$ & 0.3375 & 0.3373 & 0.3376 \\ \hline
%   \textcolor{red}{$\Pr(M_{0}|\mathbf{\underline{X}})$} & \textcolor{red}{0.9608} & \textcolor{red}{$2.1756 \times 10^{-4}$}  & \textcolor{red}{$1.0083 \times 10^{-45}$} \\ \hline
%   \textcolor{red}{$\Pr(M_{1}|\mathbf{\underline{X}})$} & \textcolor{red}{0.0385} & \textcolor{red}{0.9845} & \textcolor{red}{$3.1704 \times 10^{-6}$} \\ \hline
%   \textcolor{red}{$\Pr(M_{2}|\mathbf{\underline{X}})$} & \textcolor{red}{$6.5463 \times 10^{-4}$} & \textcolor{red}{0.0153} & \textcolor{red}{1.0000} \\ \hline
%        $B_{01}$ &  36.5471  & $3.2388 \times 10^{-4}$ & $4.6538 \times 10^{-40}$\\ \hline
%        $B_{02}$ & $1.8435 \times 10^{3}$ & 0.0178 & $1.2659 \times 10^{-45}$ \\ \hline
%        $B_{12}$ &  50.4419  & 55.0020 & $2.7202 \times 10^{-6}$ \\ \hline
%    \end{tabular}
%\caption{The results for the normalised posterior model probabilities are shown in red for a dataset of 100 points where we have used parallelised loops in a cluster of 10 threads to compute the integrals via the crude Monte Carlo method where the number of repetitions was set to 100000. Note that $B_{01}$ is the Bayes Factor between model $M_{0}$ (no change point) and model $M_{1}$ (one change point) and so on. }
%\label{TableResults}
%\end{table}
%
%
%In Table \ref{TableResults}, the model labelled as `No change point' corresponds to the whole 100 data points simulated from a \textit{Weibull}(1.5, 5), whereas model `One change point' has the first 50 data points follow a \textit{Weibull}(1.5, 5), while the next 50 data points come from a \textit{Log-normal}(0.05,16) and `Two change points' model's data is split as follows: the first 50 data points come from a \textit{Weibull}(1.5, 5), the next 20 data points are simulated from a \textit{Log-normal}(0.05, 16) and the last 30 points come from a \textit{Gamma}(10, 2).
%\end{comment}

\section*{Appendix}
\setcounter{subsection}{0}
\renewcommand{\thesubsection}{\Alph{subsection}}
\subsection{Model prior probabilities to select among models $M_0$, $M_1$ and $M_2$}
Here, we show how model prior probabilities can be derived for the relatively simple case of selecting among scenarios with no change points ($M_0$), one change point ($M_1$) or two change points ($M_2$). First, by applying the result in Theorem \ref{remark_model}, we derive the Kullback--Leibler divergences between any two models. That is:
\begin{itemize}
\item the prior probability for model $M_{0}$ depends on the following quantities:
\begin{align*}
 D_{KL}(M_0\|M_1)=&(n-m_{1})\cdot D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\\
  D_{KL}(M_0\|M_2)=&(m_{2}-m_{1})\cdot D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\\
&+(n-m_{2})\cdot D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{3}(\cdot|\tilde{\theta}_{3}))
\end{align*}

\item the prior probability for model $M_{1}$ depends on the following quantities:
\begin{align*}
 D_{KL}(M_{1}\|M_{2})=&(n-m_{2})\cdot D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{3}(\cdot|\tilde{\theta}_{3}))\\
 D_{KL}(M_1\|M_0)=&(n-m_{1})\cdot D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))
\end{align*}

\item the prior probability for model $M_{2}$ depends on the following quantities:
\begin{align*}
 D_{KL}(M_2\|M_1)=&(n-m_{2})\cdot D_{KL}(f_{3}(\cdot|\tilde{\theta}_{3})\|f_{2}(\cdot|\tilde{\theta}_{2}))\\
 D_{KL}(M_2\|M_0)=&(m_{2}-m_{1})\cdot D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\\
&+(n-m_{2})\cdot D_{KL}(f_{3}(\cdot|\tilde{\theta}_{3})\|f_{1}(\cdot|\tilde{\theta}_{1}))
\end{align*}
\end{itemize}
%Recalling the need for $\displaystyle\inf_{\theta_{i},i\neq j} D_{KL}(M_j\|M_i)$, we obtain:
The next step is to derive the minimum Kullback--Leibler divergence computed at each model:
\begin{itemize}
\item for model $M_{0}$:
\begin{align*}
 \inf_{\theta_1}D_{KL}(M_0\|M_1)=&\underbrace{\left[\inf_{m_{1}\neq n}(n-m_{1})\right]}_{1}\cdot \left[\inf_{\tilde{\theta}_2}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\right]\\=& \inf_{\tilde{\theta}_2}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\\
   \inf_{\theta_2}D_{KL}(M_0\|M_2)=&\underbrace{\left[\inf_{m_{1}\neq m_{2}}(m_{2}-m_{1})\right]}_{1}\cdot \left[\inf_{\tilde{\theta}_{2}}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\right]\\
&+\underbrace{\left[\inf_{m_{2} \neq n}(n-m_{2})\right]}_{1}\cdot \left[\inf_{\tilde{\theta}_3}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{3}(\cdot|\tilde{\theta}_{3}))\right]\\ =& \inf_{\tilde{\theta}_{2}}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))+\inf_{\tilde{\theta}_3}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{3}(\cdot|\tilde{\theta}_{3}))
\end{align*}

\item for model $M_{1}$:
\begin{align*}
 \inf_{\theta_2}D_{KL}(M_{1}\|M_{2})=&\underbrace{\left[\inf_{m_2 \neq n}(n-m_{2})\right]}_{1}\cdot \left[\inf_{\tilde{\theta}_3}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{3}(\cdot|\tilde{\theta}_{3}))\right]\\=&\inf_{\tilde{\theta}_3}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{3}(\cdot|\tilde{\theta}_{3}))\\
 \inf_{\theta_0=\tilde{\theta}_1}D_{KL}(M_1\|M_0)=&(n-m_{1})\cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))
\end{align*}

\item for model $M_{2}$:
\begin{align*}
 \inf_{\theta_1}D_{KL}(M_2\|M_1)=&(n-m_{2})\cdot \inf_{\tilde{\theta}_2}D_{KL}(f_{3}(\cdot|\tilde{\theta}_{3})\|f_{2}(\cdot|\tilde{\theta}_{2}))\\
 \inf_{\theta_0=\tilde{\theta}_1}D_{KL}(M_2\|M_0)=&(m_{2}-m_{1})\cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\\
&+(n-m_{2})\cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{3}(\cdot|\tilde{\theta}_{3})\|f_{1}(\cdot|\tilde{\theta}_{1}))
\end{align*}
\end{itemize}
Therefore, the model prior probabilities can be computed through equation \eqref{ModelPriors}, so that:
\begin{itemize}
\item the model prior probability $\Pr(M_{0})$ is proportional to the exponential of the minimum between:
\begin{align*}
 &\left\lbrace\mathbb{E}_{\pi_{0}}\left[\inf_{\tilde{\theta}_2}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\right],\mathbb{E}_{\pi_{0}}\left[\inf_{\tilde{\theta}_{2}}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{2}(\cdot|\tilde{\theta}_{2}))\right.\right. \nonumber \\ &\left.\left. \hspace*{0.3cm}+\inf_{\tilde{\theta}_3}D_{KL}(f_{1}(\cdot|\tilde{\theta}_{1})\|f_{3}(\cdot|\tilde{\theta}_{3}))\right]\right\rbrace
\end{align*}

\item the model prior probability $\Pr(M_{1})$ is proportional to the exponential of the minimum between:
\begin{align*}
 &\left\lbrace\mathbb{E}_{\pi_{1}}\left[\inf_{\tilde{\theta}_3}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{3}(\cdot|\tilde{\theta}_{3}))\right], \right. \nonumber \\ &\left. \hspace*{0.3cm} \mathbb{E}_{\pi_1}\left[(n-m_{1}) \cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]\right\rbrace
\end{align*}

\item the model prior probability $\Pr(M_{2})$ is proportional to the exponential of the minimum between:
\begin{align*}
  &\left\lbrace\mathbb{E}_{\pi_{2}}\left[(n-m_{2})\cdot \inf_{\tilde{\theta}_2}D_{KL}(f_{3}(\cdot|\tilde{\theta}_{3})\|f_{2}(\cdot|\tilde{\theta}_{2}))\right],\right.\nonumber\\&\left.\hspace*{0.3cm}\mathbb{E}_{\pi_{2}}\left[(m_{2}-m_{1})\cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))+(n-m_{2})\right.\right.\nonumber\\&\left.\left.\hspace*{0.3cm}\cdot\inf_{\tilde{\theta}_1}D_{KL}(f_{3}(\cdot|\tilde{\theta}_{3})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]\right\rbrace
\end{align*}
\end{itemize}
\newpage
\subsection{Proofs}
\subsubsection*{Proof of Theorem \ref{lemmaOneChgPoint}}
%\begin{proof}
We distinguish two cases: $S=+1$ and $S=-1$. When $S=+1$, equivalent to $m_j<m_j^{\prime}$:
\begin{align}
&D_{KL}(f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta})\|f(\mathbf{x}^{(n)}|\bm{m}',\bm{\tilde\theta}))=\int f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta}) \cdot \ln\left(\dfrac{f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta})}{f(\mathbf{x}^{(n)}|\bm{m}',\bm{\tilde\theta})}\right)\,\mathrm{d}\mathbf{x}^{(n)} \nonumber \\
=&\int  f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta}) \cdot \left[\sum_{i=m_{j}+1}^{m_{j}^{\prime}} \ln\left(\dfrac{f_{j+1}(x_{i}|\tilde\theta_{j+1})}{f_{j}(x_{i}|\tilde\theta_{j})}\right)\right]\,\mathrm{d}\mathbf{x}^{(n)} \nonumber \\
=&\mathlarger{\sum_{i=m_{j}+1}^{m_{j}^{\prime}}} \int f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta}) \cdot \left[\ln\left(\dfrac{f_{j+1}(x_{i}|\tilde\theta_{j+1})}{f_{j}(x_{i}|\tilde\theta_{j})}\right)\right]\,\mathrm{d}\mathbf{x}^{(n)} \nonumber \\
=&\mathlarger{\sum_{i=m_{j}+1}^{m_{j}^{\prime}}} \left\lbrace 1^{n-1} \cdot \int f_{j+1}(x_{i}|\tilde\theta_{j+1}) \cdot \left[\ln\left(\dfrac{f_{j+1}(x_{i}|\tilde\theta_{j+1})}{f_{j}(x_{i}|\tilde\theta_{j})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
=&\sum_{i=m_{j}+1}^{m_{j}^{\prime}} D_{KL}(f_{j+1}(x_{i}|\tilde\theta_{j+1})\|f_{j}(x_{i}|\tilde\theta_{j})) \nonumber \\
=&(m_j^{\prime}-m_{j}) \cdot D_{KL}(f_{j+1}(\cdot|\tilde\theta_{j+1})\|f_{j}(\cdot|\tilde\theta_{j})) \nonumber \\
=&(m_j^{\prime}-m_{j}) \cdot d_j^{+1}(\bm{\tilde\theta}).
\label{demopositive}
\end{align}
When $S=-1$, equivalent to $m_j>m_j^{\prime}$, in a similar fashion, we get 
\begin{align}
&D_{KL}(f(\mathbf{x}^{(n)}|\bm{m},\bm{\tilde\theta})\|f(\mathbf{x}^{(n)}|\bm{m}',\bm{\tilde\theta}))
=(m_{j}-m_{j}^{\prime}) \cdot d_j^{-1}(\bm{\tilde\theta})
\label{demonegative}
\end{align}
From equations \eqref{demopositive} and \eqref{demonegative}, we get the result in Theorem \ref{lemmaOneChgPoint}.
%\end{proof}

\subsubsection*{Proof of Theorem \ref{remark_model}}
%\begin{proof}
We recall that the model parameter $\theta_i$ is the vector $(m_1, m_2, \ldots, m_i, \tilde{\theta}_1, \tilde{\theta}_2, \ldots, \tilde{\theta}_{i+1})$, where $i=0, 1, \ldots, k$. Here, $\tilde{\theta}_1, \tilde{\theta}_2,\ldots, \tilde{\theta}_{i+1}$ represent the parameters of the underlying sampling distributions considered under model $M_i$ and $m_1, m_2, \ldots, m_i$ are the respective $i$ change point locations. In this setting, 
\begin{align}
f(\mathbf{x}^{(n)}|\theta_i)=\prod_{r=1}^{m_{1}} f_{1}(x_{r}|\tilde{\theta}_{1})\prod_{t=1}^{i-1}\prod_{r=m_{t}+1}^{m_{t+1}} f_{t+1}(x_{r}|\tilde{\theta}_{t+1})\prod_{r=m_{i}+1}^{n} f_{i+1}(x_{r}|\tilde{\theta}_{i+1})
\label{proof2_0}
\end{align}

We proceed to the computation of $D_{KL}(M_{i}\|M_{j})$, that is the Kullback--Leibler divergence  introduced in Section \ref{sc_numberCP}. Similarly to the proof of Theorem \ref{lemmaOneChgPoint}, we obtain the following result.  
\begin{align*}
D_{KL}(M_{i}\|M_{j})=&\sum_{r=m_{i+1}+1}^{m_{i+2}} \int f(\mathbf{x}^{(n)}|\theta_i)\ln\left(\dfrac{f_{i+1}(x_{r}|\tilde\theta_{i+1})}{f_{i+2}(x_{r}|\tilde\theta_{i+2})}\right)\,\mathrm{d}\mathbf{x}^{(n)}\nonumber \\
&+ \sum_{r=m_{i+2}+1}^{m_{i+3}} \int f(\mathbf{x}^{(n)}|\theta_i)\ln\left(\dfrac{f_{i+1}(x_{r}|\tilde\theta_{i+1})}{f_{i+3}(x_{r}|\tilde\theta_{i+3})}\right)\,\mathrm{d}\mathbf{x}^{(n)}+\nonumber \\
&\ldots+ \sum_{r=m_{j}+1}^{n} \int f(\mathbf{x}^{(n)}|\theta_i) \ln\left(\dfrac{f_{i+1}(x_{r}|\tilde\theta_{i+1})}{f_{j+1}(x_{r}|\tilde\theta_{j+1})}\right)\,\mathrm{d}\mathbf{x}^{(n)}.
\end{align*}
Given equation \eqref{proof2_0}, if we integrate out the variables not involved in the logarithms, we obtain
\begin{align*}
 D_{KL}(M_{i}\|M_{j})=&(m_{i+2}-m_{i+1})\cdot D_{KL}(f_{i+1}(\cdot|\tilde{\theta}_{i+1})\|f_{i+2}(\cdot|\tilde{\theta}_{i+2})) \\
&+(m_{i+3}-m_{i+2})\cdot D_{KL}(f_{i+1}(\cdot|\tilde{\theta}_{i+1})\|f_{i+3}(\cdot|\tilde{\theta}_{i+3}))+\\
&\ldots+(n-m_{j})\cdot D_{KL}(f_{i+1}(\cdot|\tilde{\theta}_{i+1})\|f_{j+1}(\cdot|\tilde{\theta}_{j+1})).
\end{align*}
In a similar fashion, it can be shown that
\begin{align*}
 D_{KL}(M_{j}\|M_{i})=&(m_{i+2}-m_{i+1})\cdot D_{KL}(f_{i+2}(\cdot|\tilde{\theta}_{i+2})\|f_{i+1}(\cdot|\tilde{\theta}_{i+1}))\\
&+(m_{i+3}-m_{i+2})\cdot D_{KL}(f_{i+3}(\cdot|\tilde{\theta}_{i+3})\|f_{i+1}(\cdot|\tilde{\theta}_{i+1}))+\\
&\ldots+(n-m_{j})\cdot D_{KL}(f_{j+1}(\cdot|\tilde{\theta}_{j+1})\|f_{i+1}(\cdot|\tilde{\theta}_{i+1}))
\end{align*}
%\end{proof}

%\subsubsection*{Proof of Equation \eqref{prior_M1_trans_one_chg}}
%%\begin{proof}
%First, let us recall that $\pi_1$ is a shorthand notation for  $\pi_1(m_1,\tilde{\theta}_1, \tilde{\theta}_2)$. Also, let the prior on $m_1$ be independent of the prior on the parameters of the underlying distributions,  $\pi_1(m_1,\tilde{\theta}_1, \tilde{\theta}_2)=\pi_1(m_1)\cdot \pi_1(\tilde{\theta}_1, \tilde{\theta}_2)$. Furthermore, let us denote by
%\begin{align}
%u(\tilde{\theta}_1, \tilde{\theta}_2)&\equiv\mathbb{E}_{\pi_1(\tilde{\theta}_1, \tilde{\theta}_2)}\left[\inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right] \label{Proof_in_M1_Unif_Case}\\
%&=\int\left[\inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]\cdot \pi_1(\tilde{\theta}_1, \tilde{\theta}_2)\,\mathrm{d}\tilde{\theta}_1\,\mathrm{d}\tilde{\theta}_2 \nonumber
%\end{align}
%
%We observe that $u(\tilde{\theta}_1, \tilde{\theta}_2)$ does not depend on $m_1$.
%In this proof, we make use of the following prior:
%\begin{align*}
%  &\pi_1(m_1)=\dfrac{1}{n-1}, \qquad 1\leq m_1\leq n-1, 
%\end{align*}
%and we employ the expectation of a discrete uniform random variable on the support $\{1, 2, \ldots, n-1\}$, that is $\mathbb{E}_{\pi_1(m_1)}(m_1)=\dfrac{n}{2}$.
% 
% Therefore, by expanding Equation \eqref{prior_M1_one_chg} and by using Equation \eqref{Proof_in_M1_Unif_Case}, as well as the independence between the change point location prior and the prior on the parameters of the underlying distributions, we get:
%\begin{footnotesize}
%\begin{align}
%  \mathbb{E}_{\pi_1}\left[(n-m_{1}) \cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]=& \sum_{m_1=1}^{n-1}(n-m_1)\cdot \pi_1(m_1)\cdot u(\tilde{\theta}_1, \tilde{\theta}_2)\nonumber\\
%  =&u(\tilde{\theta}_1, \tilde{\theta}_2) \cdot \sum_{m_1=1}^{n-1}(n-m_1)\cdot \pi_1(m_1) \nonumber\\
%  =&u(\tilde{\theta}_1, \tilde{\theta}_2) \cdot \left\lbrace \left[n\cdot\underbrace{\sum_{m_1=1}^{n-1}\pi_1(m_1)}_{1}\right]\right. \nonumber\\
%  &\left.-\underbrace{\left[\sum_{m_1=1}^{n-1}m_1\cdot\pi_1(m_1)\right]}_{\mathbb{E}_{\pi_1(m_1)}(m_1)}\right\rbrace \nonumber\\
%  =&u(\tilde{\theta}_1, \tilde{\theta}_2) \cdot \left(n-\dfrac{n}{2}\right)=\dfrac{n}{2} \cdot u(\tilde{\theta}_1, \tilde{\theta}_2)
%  \label{Proof_in_M1_Unif_Case_2}
%\end{align} 
%\end{footnotesize}
%
%From Equations \eqref{Proof_in_M1_Unif_Case} and \eqref{Proof_in_M1_Unif_Case_2}, we obtain  
%\begin{footnotesize}
%\begin{align*}
%  \mathbb{E}_{\pi_1}\left[(n-m_{1}) \cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]=&\left(\dfrac{n}{2}\right) \\
%  &\cdot \mathbb{E}_{\pi_1(\tilde{\theta}_1, \tilde{\theta}_2)}\left[\inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right],
%\end{align*}
%\end{footnotesize}
%thus arriving at Equation \eqref{prior_M1_trans_one_chg}.
%%\end{proof}
%
%\subsubsection*{Proof of Equation \eqref{prior_M1_pseudo_binom}}
%%\begin{proof}
%It is similar to the proof for Equation \eqref{prior_M1_trans_one_chg}. In this situation, we make use of the following prior:
%\begin{align*}
%  \pi_1(m_1)=\dbinom{n-2}{m_1-1}\cdot\left(\dfrac{n-1}{n}\right)^{m_1-1}\cdot\left(\dfrac{1}{n}\right)^{n-m_1-1}, 1\leq m_1 \leq n-1,
% \end{align*}
%and we employ the expectation of a shifted binomial random variable on the support $\{1, 2, \ldots, n-1\}$, where the probability of success in each trial is $\dfrac{n-1}{n}$, that is $\mathbb{E}_{\pi_1(m_1)}(m_1-1)=(n-2)\cdot \left(\dfrac{n-1}{n}\right)$.
%
% 
% Therefore, by expanding Equation \eqref{prior_M1_one_chg} and using Equation \eqref{Proof_in_M1_Unif_Case}, as well as the independence between the change point location prior and the prior on the parameters of the underlying distributions, we get:
%\begin{footnotesize}
%\begin{align}
%  \mathbb{E}_{\pi_1}\left[(n-m_{1}) \cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]=& \sum_{m_1=1}^{n-1}(n-m_1)\cdot \pi_1(m_1)\cdot u(\tilde{\theta}_1, \tilde{\theta}_2)\nonumber\\
%  =&u(\tilde{\theta}_1, \tilde{\theta}_2) \cdot \sum_{m_1=1}^{n-1}[n-(m_1-1)-1]\cdot \pi_1(m_1) \nonumber\\
%  =&u(\tilde{\theta}_1, \tilde{\theta}_2) \cdot \left\lbrace \left[(n-1) \cdot \underbrace{\sum_{m_1=1}^{n-1}\pi_1(m_1)}_{1}\right]\right. \nonumber\\
%  &\left.-\underbrace{\left[\sum_{m_1=1}^{n-1}(m_1-1)\cdot\pi_1(m_1)\right]}_{\mathbb{E}_{\pi_1(m_1)}(m_1-1)}\right\rbrace \nonumber\\
%  =&u(\tilde{\theta}_1, \tilde{\theta}_2) \cdot \left(n-1-\dfrac{(n-2)\cdot(n-1)}{n}\right)\nonumber\\
%  =&\left(\dfrac{2n-2}{n}\right) \cdot u(\tilde{\theta}_1, \tilde{\theta}_2)
%  \label{Proof_in_M2_Unif_Case_2}
%\end{align} 
%\end{footnotesize}
%
%From Equations \eqref{Proof_in_M1_Unif_Case} and \eqref{Proof_in_M2_Unif_Case_2} we obtain  
%\begin{footnotesize}
%\begin{align*}
%  \mathbb{E}_{\pi_1}\left[(n-m_{1}) \cdot \inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right]=&\left(\dfrac{2n-2}{n}\right) \\
%  &\cdot \mathbb{E}_{\pi_1(\tilde{\theta}_1, \tilde{\theta}_2)}\left[\inf_{\tilde{\theta}_1}D_{KL}(f_{2}(\cdot|\tilde{\theta}_{2})\|f_{1}(\cdot|\tilde{\theta}_{1}))\right],
%\end{align*}
%\end{footnotesize}
%thus arriving at Equation \eqref{prior_M1_pseudo_binom}.
%\end{proof}

%\subsubsection*{Proof that the Infimum of the KL divergence between two Log-normal Distributions is 0}
%%\begin{proof}
%Let us denote with $f_2(\cdot|\mu_2,\tau_2)$ and $f_3(\cdot|\mu_3,\tau_3)$, the probability densities of two Log-normal distributions which are parametrized as following:
%\begin{align*}
%f_{2}(x|\mu_2, \tau_2)&=\dfrac{1}{x}\cdot\left(\dfrac{\tau_2}{2\pi}\right)^{1/2}\cdot e^{-\dfrac{\tau_2}{2}\cdot[\ln(x)-\mu_2]^2}\\
%f_{3}(x|\mu_3, \tau_3)&=\dfrac{1}{x}\cdot\left(\dfrac{\tau_3}{2\pi}\right)^{1/2}\cdot e^{-\dfrac{\tau_3}{2}\cdot[\ln(x)-\mu_3]^2}.
%\end{align*}
%Then, the Kullback--Leibler divergence from $f_2(\cdot|\mu_2,\tau_2)$ to $f_3(\cdot|\mu_3,\tau_3)$ is:
%\begin{align}
%D_{KL}(f_2(x|\mu_2,\tau_2)\|f_3(x|\mu_3,\tau_3))=&\int\displaylimits_{0}^{\infty} f_2(x|\mu_2,\tau_2) \cdot \ln\left[\dfrac{f_2(x|\mu_2,\tau_2)}{f_3(x|\mu_3,\tau_3)}\right]\,\mathrm{d}x \nonumber \\
%=&\int\displaylimits_{0}^{\infty} f_2(x|\mu_2,\tau_2) \cdot \ln\left\lbrace\dfrac{\dfrac{1}{x}\cdot\left(\dfrac{\tau_2}{2\pi}\right)^{1/2}}{\dfrac{1}{x}\cdot\left(\dfrac{\tau_3}{2\pi}\right)^{1/2}} \right. \nonumber \\ 
%&\left.\cdot \dfrac{e^{-\dfrac{\tau_2}{2}\cdot[\ln(x)-\mu_2]^2}}{e^{-\dfrac{\tau_3}{2}\cdot[\ln(x)-\mu_3]^2}}\right\rbrace\,\mathrm{d}x \nonumber \\
%=&\int\displaylimits_{0}^{\infty} f_2(x|\mu_2,\tau_2) \cdot \left\lbrace \dfrac{1}{2} \cdot\ln\left(\dfrac{\tau_2}{\tau_3}\right) \right. \nonumber \\
%&\left.-\left(\dfrac{\tau_2}{2}\right) \cdot [\ln^2(x)-2\mu_2\ln(x)+\mu_2^2] \right. \nonumber \\
%&\left.+\left(\dfrac{\tau_3}{2}\right) \cdot [\ln^2(x)-2\mu_3\ln(x)+\mu_3^2]\vphantom{\dfrac{1}{2}}\right\rbrace\,\mathrm{d}x  \nonumber \\
%=&\int\displaylimits_{0}^{\infty} f_2(x|\mu_2,\tau_2) \cdot \left\lbrace \dfrac{1}{2} \cdot\ln\left(\dfrac{\tau_2}{\tau_3}\right) \right. \nonumber \\
%&\left.+\left[\left(\dfrac{\tau_3}{2}-\dfrac{\tau_2}{2}\right)\cdot \ln^2(x)\right]+\left[\left(\tau_2\mu_2\right.\right.\right. \nonumber \\ 
%&\left.\left.\left.-\tau_3\mu_3\right)\cdot \ln(x)\right]+\left(\dfrac{\tau_3\mu_3^2}{2}-\dfrac{\tau_2\mu_2^2}{2}\right)\vphantom{\dfrac{1}{2}}\right\rbrace\,\mathrm{d}x  \nonumber \\
%\displaybreak \nonumber \\
%=&\dfrac{1}{2} \cdot\ln\left(\dfrac{\tau_2}{\tau_3}\right)+\left(\dfrac{\tau_3}{2}-\dfrac{\tau_2}{2}\right)\cdot \mathbb{E}_2[\ln^2(x)] \nonumber \\
%&+\left(\tau_2\mu_2-\tau_3\mu_3\right)\cdot  \mathbb{E}_2[\ln(x)] \nonumber \\
%&+\left(\dfrac{\tau_3\mu_3^2}{2}-\dfrac{\tau_2\mu_2^2}{2}\right)
%\label{KL_lognormals_1}
%\end{align}
%According to \citet[pp.~207--208]{johnson1994}, the logarithm of a log-normally distributed random variable with parameters $\mu_2$ and $\tau_2$ has a normal distribution with location $\mu_2$ and scale $1/\tau_2$. So, we have that $\mathbb{E}_2[\ln(x)]=\mu_2$ and $\mathbb{E}_2[\ln^2(x)]=\mu_2^{2}+\dfrac{1}{\tau_2}$. As such, Equation \eqref{KL_lognormals_1} becomes:
%\begin{align}
%D_{KL}(f_2(x|\mu_2,\tau_2)\|f_3(x|\mu_3,\tau_3))=&\dfrac{1}{2} \cdot\ln\left(\dfrac{\tau_2}{\tau_3}\right)+\left(\dfrac{\tau_3}{2}-\dfrac{\tau_2}{2}\right)\cdot \left(\mu_2^2+\dfrac{1}{\tau_2}\right) \nonumber \\
%&+\left(\tau_2\mu_2-\tau_3\mu_3\right)\cdot \mu_2 +\left(\dfrac{\tau_3\mu_3^2}{2}-\dfrac{\tau_2\mu_2^2}{2}\right)
%\label{KL_lognormals_2}
%\end{align}
%To find the infinimum with respect to $\mu_3$ and $\tau_3$ of the Kullback--Leibler divergence from $f_2(\cdot|\mu_2,\tau_2)$ to $f_3(\cdot|\mu_3,\tau_3)$, we need to differentiate once Equation \eqref{KL_lognormals_2} with respect to $\mu_3$ and $\tau_3$ and set the resulting equations to 0. Therefore, we obtain:
%\begin{align}
%\dfrac{\partial D_{KL}(f_2(x|\mu_2,\tau_2)\|f_3(x|\mu_3,\tau_3))}{\partial \mu_3}&= -\tau_3\mu_2+\tau_3\mu_3=0 \label{KL_lognormals_3} \\
%\dfrac{\partial D_{KL}(f_2(x|\mu_2,\tau_2)\|f_3(x|\mu_3,\tau_3))}{\partial \tau_3}&= -\dfrac{1}{2\tau_3}+\dfrac{\mu_2^2}{2}+\dfrac{1}{2\tau_2}-\mu_2\mu_3+\dfrac{\mu_3^2}{2}=0 \label{KL_lognormals_4}
%\end{align}
%From Equations \eqref{KL_lognormals_3} and \eqref{KL_lognormals_4}, as well as using the constraint that $\tau_2,\tau_3>0$, we have that $\mu_2=\mu_3$ and $\tau_2=\tau_3$, respectively. As such, when we plug back those values into Equation \eqref{KL_lognormals_2}, we end up with:
%\begin{align*}
% \inf_{\mu_3,\tau_3}D_{KL}(f_2(x|\mu_2,\tau_2)\|f_3(x|\mu_3,\tau_3))=0.
%\end{align*}
%In a similar fashion, it can be shown that:
%\begin{align*}
% \inf_{\mu_2,\tau_2}D_{KL}(f_3(x|\mu_3,\tau_3)\|f_2(x|\mu_2,\tau_2))=0.
%\end{align*}
%\end{proof}


%\subsection{Infima of the KL divergences needed in Section \ref{sc_simulation}}
%From \citet{VillaModel}, we have the following results necessary for obtaining the model prior probabilities (needed in equations \eqref{priorM0}, \eqref{priorM1} and \eqref{priorM2}), where $f_{1}(\cdot|\lambda, \kappa), f_{2}(\cdot|\mu, \tau)$ and $f_{3}(\cdot|\alpha, \beta)$ are the Weibull, Log-normal and Gamma distributions, respectively:
%\begin{align*}
%\inf_{\mu,\tau} D_{KL}(f_1(\cdot|\lambda,\kappa)\|f_2(\cdot|\mu, \tau))=&\dfrac{\ln(2\pi)}{2}+\ln(\pi)-\gamma-\dfrac{\ln(6)}{2}-\dfrac{1}{2},\\
%\inf_{\lambda,\kappa} D_{KL}(f_2(\cdot|\mu,\tau)\|f_1(\cdot|\lambda, \kappa))=&1-\dfrac{\ln(2\pi)}{2},\\
%\inf_{\alpha,\beta} D_{KL}(f_1(\cdot|\lambda,\kappa)\|f_3(\cdot|\alpha, \beta))=&\ln(\kappa)-\gamma-1-\alpha\ln(\beta)+\ln[\Gamma(\alpha)]\\&-\alpha\ln(\lambda)+\alpha\cdot \dfrac{\gamma}{\kappa}+\beta\lambda\Gamma\left(1+\dfrac{1}{\kappa}\right),\\
%\inf_{\alpha,\beta} D_{KL}(f_2(\cdot|\mu,\tau)\|f_3(\cdot|\alpha, \beta))=&\dfrac{\ln(\tau)}{2}-\dfrac{\ln(2\pi)}{2}-\dfrac{1}{2}-\alpha\ln(\beta)+\ln[\Gamma(\alpha)]\\&-\alpha\mu+\beta\cdot e^{\mu+\dfrac{1}{2\tau}},\\
%\inf_{\lambda,\kappa} D_{KL}(f_3(\cdot|\alpha,\beta)\|f_1(\cdot|\lambda, \kappa))=&-\ln[\Gamma(\alpha)]+\alpha\Psi(\alpha)-\alpha-\ln(\kappa)\\&-\kappa\Psi(\alpha)+\kappa\ln(\beta)+\kappa\ln(\lambda)+1,\\
%\inf_{\mu,\tau} D_{KL}(f_3(\cdot|\alpha,\beta)\|f_2(\cdot|\mu, \tau))=&-\ln[\Gamma(\alpha)]+\alpha\Psi(\alpha)-\alpha\\&+\dfrac{\ln[\Psi^{\prime}(\alpha)]}{2}+\dfrac{\ln(2\pi)}{2}+\dfrac{1}{2},
%\end{align*}
%where $\gamma$ is the Euler-Mascheroni constant and $\Gamma$, $\Psi$ and $\Psi^{\prime}$ are the Gamma, Digamma and Trigamma functions, respectively.
%
%%Now, let us denote with $f_1(\cdot|\textrm{p}_1)$ and $f_2(\cdot|\textrm{p}_2)$, the probability mass functions of two Geometric distributions, and with $f_1(\cdot|\lambda_1)$ and $f_2(\cdot|\lambda_2)$, the probability mass functions of two Poisson distributions. Then, following \citet{VillaModel}, we have:
%Now, let us denote with $f_1(\cdot|\textrm{p})$ and $f_2(\cdot|\lambda)$ , the probability mass functions of a Geometric distribution and a Poisson distribution, respectively. Then, following \citet{VillaModel}, we have:
%\begin{align*}
%\inf_{\lambda} D_{KL}(f_1(x|\textrm{p})\|f_2(x|\lambda))=&\ln(\textrm{p})+\dfrac{1-\textrm{p}}{\textrm{p}}\cdot \ln(\textrm{p})+\dfrac{1-\textrm{p}}{\textrm{p}}\\
%\displaybreak
%&+\left[\sum_{x=0}^{\infty} \textrm{p} \cdot (1-\textrm{p})^x \cdot \ln(x!)\right],\\
%\inf_{\textrm{p}} D_{KL}(f_2(x|\lambda)\|f_1(x|\textrm{p}))=&-\lambda+\lambda\cdot\ln(1+\lambda)+\ln(1+\lambda)\\
%&-\left[\sum_{x=0}^{\infty} \dfrac{\lambda^{x}}{x!}\cdot \ln(x!) \cdot e^{-\lambda} \right].\\
%\end{align*}

%\begin{comment}
%\begin{align*}
%&=\sum_{i=m_{j}^{\prime}+1}^{m_{j}} \int \left( \prod_{r=1}^{m_{j}^{\prime}} f_{j}(x_{r}|\theta_{j})\right) \cdot \left( \prod_{r=m_{j}+1}^{n} f_{j+1}(x_{r}|\theta_{j+1})\right)\cdot \nonumber  \\ &\cdot \left[\left( \prod_{i=m_{j}^{\prime}+1}^{m_{j}} f_{j}(x_{i}|\theta_{j})\right)\cdot \ln\left(\dfrac{f_{j}(x_{i}|\theta_{j})}{f_{j+1}(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}\mathbf{x}^{(n)} \nonumber \\
%\end{align*}
%
%\subsection*{Proof of Lemma \ref{SimpleOneChgPoint}}
%\begin{proof}
%When we expand equation \ref{klremarkOneChg}, we need to take into account two cases in this context:
%\begin{itemize}
%\item when $m_{1}<m_{1}^{\prime}$, we have 
%\begin{align}
%&D_{KL}(f(\mathbf{X^{(n)}}|m_{1},\theta_{1},\theta_{2})\|f(\mathbf{X^{(n)}}|m_{1}^{\prime},\theta_{1},\theta_{2}))=\int f(\mathbf{X^{(n)}}|m_{1},\theta_{1},\theta_{2}) \cdot \nonumber\\ &\left\lbrace\left[\sum_{i=1}^{m_{1}} \ln(f_{1}(x_{i}|\theta_{1}))\right]+\left[\sum_{i=m_{1}+1}^{m_{1}^{\prime}} \ln(f_{2}(x_{i}|\theta_{2}))\right]+\right.\nonumber\\&\left.+\left[\sum_{i=m_{1}^{\prime}+1}^{\textbf{n}} \ln(f_{2}(x_{i}|\theta_{2}))\right]-\left[\sum_{i=1}^{m_{1}} \ln(f_{1}(x_{i}|\theta_{1}))\right]-\right. \nonumber\\ &-\left.\left[\sum_{i=m_{1}+1}^{m_{1}^{\prime}} \ln(f_{1}(x_{i}|\theta_{1}))\right]-\left[\sum_{i=m_{1}^{\prime}+1}^{\textbf{n}} \ln(f_{2}(x_{i}|\theta_{2}))\right]\right\rbrace \,\mathrm{d}\mathbf{X^{(n)}}\nonumber\\
%&=\int \left( \prod_{j=1}^{m_{1}} f_{1}(x_{j}|\theta_{1})\right)\cdot \left( \prod_{i=m_{1}+1}^{m_{1}^{\prime}} f_{2}(x_{i}|\theta_{2})\right) \cdot \left( \prod_{j=m_{1}^{\prime}+1}^{\textbf{n}} f_{2}(x_{j}|\theta_{2})\right)\cdot \nonumber  \\ &\cdot \left[\sum_{i=m_{1}+1}^{m_{1}^{\prime}} \ln\left(\dfrac{f_{2}(x_{i}|\theta_{2})}{f_{1}(x_{i}|\theta_{1})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\sum_{i=m_{1}+1}^{m_{1}^{\prime}} \int \left( \prod_{j=1}^{m_{1}} f_{1}(x_{j}|\theta_{1})\right) \cdot \left( \prod_{j=m_{1}^{\prime}+1}^{\textbf{n}} f_{2}(x_{j}|\theta_{2})\right)\cdot \nonumber  \\ &\cdot \left[\left( \prod_{i=m_{1}+1}^{m_{1}^{\prime}} f_{2}(x_{i}|\theta_{2})\right)\cdot\ln\left(\dfrac{f_{2}(x_{i}|\theta_{2})}{f_{1}(x_{i}|\theta_{1})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%\displaybreak \nonumber\\
%&=\sum_{i=m_{1}+1}^{m_{1}^{\prime}} \left\lbrace \int f_{2}(x_{i}|\theta_{2}) \cdot  \left[\ln\left(\dfrac{f_{2}(x_{i}|\theta_{2})}{f_{1}(x_{i}|\theta_{1})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \cdot \nonumber \\
%& \idotsint f_{1}(x_{1}|\theta_{1}) \cdot \ldots \cdot f_{1}(x_{m_{1}-1}|\theta_{1}) \cdot \ldots \cdot f_{2}(x_{i-1}|\theta_{2}) \cdot f_{2}(x_{i+1}|\theta_{2}) \cdot  \nonumber \\
%&\ldots \cdot f_{2}(x_{\textbf{n}}|\theta_{2})\,\mathrm{d}x_{1}\ldots\,\mathrm{d}x_{m_{1}-1}\ldots\,\mathrm{d}x_{i-1}\,\mathrm{d}x_{i+1}\ldots\,\mathrm{d}x_{\textbf{n}} \nonumber \\
%&=\sum_{i=m_{1}+1}^{m_{1}^{\prime}} \left\lbrace 1^{n-1} \cdot \int f_{2}(x_{i}|\theta_{2}) \cdot  \left[\ln\left(\dfrac{f_{2}(x_{i}|\theta_{2})}{f_{1}(x_{i}|\theta_{1})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{1}+1}^{m_{1}^{\prime}} \left\lbrace \int f_{2}(x_{i}|\theta_{2}) \cdot  \left[\ln\left(\dfrac{f_{2}(x_{i}|\theta_{2})}{f_{1}(x_{i}|\theta_{1})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{1}+1}^{m_{1}^{\prime}} D_{KL}(f_{2}(x_{i}|\theta_{2})\|f_{1}(x_{i}|\theta_{1}))\nonumber\\
%&=(m_{1}^{\prime}-m_{1}) \cdot D_{KL}(f_{2}(x_{i}|\theta_{2})\|f_{1}(x_{i}|\theta_{1}))
%\end{align}
%
%\item when $m_{1}>m_{1}^{\prime}$, we have 
%\begin{align}
%&D_{KL}(f(\mathbf{X^{(n)}}|m_{1},\theta_{1},\theta_{2})\|f(\mathbf{X^{(n)}}|m_{1}^{\prime},\theta_{1},\theta_{2}))=\int f(\mathbf{X^{(n)}}|m_{1},\theta_{1},\theta_{2}) \cdot \nonumber\\ &\left\lbrace\left[\sum_{i=1}^{m_{1}^{\prime}} \ln(f_{1}(x_{i}|\theta_{1}))\right]+\left[\sum_{i=m_{1}^{\prime}+1}^{m_{1}} \ln(f_{1}(x_{i}|\theta_{1}))\right]+\right.\nonumber\\&\left.+\left[\sum_{i=m_{1}+1}^{\textbf{n}} \ln(f_{2}(x_{i}|\theta_{2}))\right]-\left[\sum_{i=1}^{m_{1}^{\prime}} \ln(f_{1}(x_{i}|\theta_{1}))\right]-\right. \nonumber\\ &-\left.\left[\sum_{i=m_{1}^{\prime}+1}^{m_{1}} \ln(f_{2}(x_{i}|\theta_{2}))\right]-\left[\sum_{i=m_{1}+1}^{\textbf{n}} \ln(f_{2}(x_{i}|\theta_{2}))\right]\right\rbrace \,\mathrm{d}\mathbf{X^{(n)}}\nonumber\\
%&=\int \left( \prod_{j=1}^{m_{1}^{\prime}} f_{1}(x_{j}|\theta_{1})\right)\cdot \left( \prod_{i=m_{1}^{\prime}+1}^{m_{1}} f_{1}(x_{i}|\theta_{1})\right) \cdot \left( \prod_{j=m_{1}+1}^{\textbf{n}} f_{2}(x_{j}|\theta_{2})\right)\cdot \nonumber  \\ &\cdot \left[\sum_{i=m_{1}^{\prime}+1}^{m_{1}} \ln\left(\dfrac{f_{1}(x_{i}|\theta_{1})}{f_{2}(x_{i}|\theta_{2})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%\displaybreak \nonumber\\
%&=\sum_{i=m_{1}^{\prime}+1}^{m_{1}} \int \left( \prod_{j=1}^{m_{1}^{\prime}} f_{1}(x_{j}|\theta_{1})\right) \cdot \left( \prod_{j=m_{1}+1}^{\textbf{n}} f_{2}(x_{j}|\theta_{2})\right)\cdot \nonumber  \\ &\cdot \left[\left( \prod_{i=m_{1}^{\prime}+1}^{m_{1}} f_{1}(x_{i}|\theta_{1})\right)\cdot \ln\left(\dfrac{f_{1}(x_{i}|\theta_{1})}{f_{2}(x_{i}|\theta_{2})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\sum_{i=m_{1}^{\prime}+1}^{m_{1}} \left\lbrace \int f_{1}(x_{i}|\theta_{1}) \cdot  \left[\ln\left(\dfrac{f_{1}(x_{i}|\theta_{1})}{f_{2}(x_{2}|\theta_{2})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \cdot \nonumber \\
%& \idotsint f_{1}(x_{1}|\theta_{1}) \cdot \ldots \cdot f_{1}(x_{m_{1}^{\prime}-1}|\theta_{1}) \cdot \ldots \cdot f_{1}(x_{i-1}|\theta_{1}) \cdot f_{1}(x_{i+1}|\theta_{1}) \cdot  \nonumber \\
%&\ldots \cdot f_{2}(x_{\textbf{n}}|\theta_{2})\,\mathrm{d}x_{1}\ldots\,\mathrm{d}x_{m_{1}^{\prime}-1}\ldots\,\mathrm{d}x_{i-1}\,\mathrm{d}x_{i+1}\ldots\,\mathrm{d}x_{\textbf{n}} \nonumber \\
%&=\sum_{i=m_{1}^{\prime}+1}^{m_{1}} \left\lbrace 1^{n-1} \cdot \int f_{1}(x_{i}|\theta_{1}) \cdot  \left[\ln\left(\dfrac{f_{1}(x_{i}|\theta_{1})}{f_{2}(x_{i}|\theta_{2})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{1}^{\prime}+1}^{m_{1}} \left\lbrace \int f_{1}(x_{i}|\theta_{1}) \cdot  \left[\ln\left(\dfrac{f_{1}(x_{i}|\theta_{1})}{f_{2}(x_{i}|\theta_{2})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{1}^{\prime}+1}^{m_{1}} D_{KL}(f_{1}(x_{i}|\theta_{1})\|f_{2}(x_{i}|\theta_{2}))\nonumber\\
%&=(m_{1}-m_{1}^{\prime}) \cdot D_{KL}(f_{1}(x_{i}|\theta_{1})\|f_{2}(x_{i}|\theta_{2}))
%\end{align}
%\end{itemize}
%\end{proof}
%\subsubsection*{Proof of Lemma \ref{lemmaOneChgPoint}}
%\begin{proof}
%
%Let $\mathbf{e_{j}^{(k)}}$ be an arbitrary, but fixed vector, such that the $j$-th element is 1, whereas the other constituent elements are zeros, that is $ \mathbf{e_{j}^{(k)}}=(0,0, \ldots,\underbrace{1}_{\mathclap{\text{\textit{j}-th element}}}, \ldots,0)$.
%
%We distinguish two cases:
%\begin{itemize}
%\item when the sign is `$+$':
%\begin{align}
%&D(\mathbf{m^{(k,n)}}\|\mathbf{m^{(k,n)}+\textit{l}\cdot e_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})=\nonumber\\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\
%&\cdot\ln\left(\dfrac{L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}{L(\mathbf{m^{(k,n)}+\textit{l}\cdot e_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}\right)\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot\left\lbrace\rule{0cm}{1cm}\left[\sum_{i=1}^{m_{1}} \ln(f_{1}(x_{i}|\theta_{1}))\right]\right.+\nonumber \\
%&+\left.\left[\sum_{i=m_{1}+1}^{m_{2}} \ln(f_{2}(x_{i}|\theta_{2}))\right]+\ldots+\left[\sum_{i=m_{j-1}+1}^{m_{j}} \ln(f_{j}(x_{i}|\theta_{j}))\right]+\right.\nonumber\\ 
%&\left.+\left[\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \ln(f_{j+1}(x_{i}|\theta_{j+1}))\right]+\left[\sum_{i=m_{j}+\textit{l}+1}^{m_{j+1}} \ln(f_{j+1}(x_{i}|\theta_{j+1}))\right]+\right.\nonumber\\
%&+\left.\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \ln(f_{j+2}(x_{i}|\theta_{j+2}))\right]+\ldots+\left[\sum_{i=m_{k}+1}^{\textbf{n}} \ln(f_{k+1}(x_{i}|\theta_{k+1}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=1}^{m_{1}} \ln(f_{1}(x_{i}|\theta_{1}))\right]-\left[\sum_{i=m_{1}+1}^{m_{2}} \ln(f_{2}(x_{i}|\theta_{2}))\right]-\ldots-\right.\nonumber\\
%\displaybreak \nonumber\\
%&\left.-\left[\sum_{i=m_{j-1}+1}^{m_{j}} \ln(f_{j}(x_{i}|\theta_{j}))\right]-\left[\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \ln(f_{j}(x_{i}|\theta_{j}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=m_{j}+\textit{l}+1}^{m_{j+1}} \ln(f_{j+1}(x_{i}|\theta_{j+1}))\right]-\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \ln(f_{j+2}(x_{i}|\theta_{j+2}))\right]-\ldots-\right.\nonumber\\
%&-\left.\left[\sum_{i=m_{k}+1}^{\textbf{n}} \ln_{k+1}(f(x_{i}|\theta_{k+1}))\right]\rule{0cm}{1cm}\right\rbrace\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\ 
%&=\int  L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\ &\left[\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} \ln\left(\dfrac{f_{j+1}(x_{i}|\theta_{j+1})}{f_{j}(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\mathlarger{\sum_{i=m_{j}+1}^{m_{j}+\textit{l}}} \int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber  \\ &\cdot \left[\ln\left(\dfrac{f_{j+1}(x_{i}|\theta_{j+1})}{f_{j}(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\mathlarger{\sum_{i=m_{j}+1}^{m_{j}+\textit{l}}} \Bigg\lbrace \left\lbrace \int f_{j+1}(x_{i}|\theta_{j+1}) \cdot  \left[\ln\left(\dfrac{f_{j+1}(x_{i}|\theta_{j+1})}{f_{j}(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \cdot \nonumber \\
%& \idotsint f_{1}(x_{1}|\theta_{1}) \cdot \ldots \cdot f_{j}(x_{m_{j}-1}|\theta_{j}) \cdot \ldots \cdot f_{j+1}(x_{i-1}|\theta_{j+1}) \cdot \nonumber\\& \cdot f_{j+1}(x_{i+1}|\theta_{j+1}) \cdot \ldots \cdot f_{k+1}(x_{\mathbf{n}}|\theta_{k+1})\,\mathrm{d}x_{1}\ldots\,\mathrm{d}x_{m_{j}-1}\ldots\,\mathrm{d}x_{i-1} \nonumber \\ &\,\mathrm{d}x_{i+1}\ldots\,\mathrm{d}x_{\mathbf{n}}\Bigg\rbrace \nonumber \\
%&=\mathlarger{\sum_{i=m_{j}+1}^{m_{j}+\textit{l}}} \left\lbrace 1^{n-1} \cdot \int f_{j+1}(x_{i}|\theta_{j+1}) \cdot \left[\ln\left(\dfrac{f_{j+1}(x_{i}|\theta_{j+1})}{f_{j}(x_{i}|\theta_{j})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{j}+1}^{m_{j}+\textit{l}} D_{KL}(f_{j+1}(x_{i}|\theta_{j+1})\|f_{j}(x_{i}|\theta_{j})) \nonumber \\
%\displaybreak \nonumber\\
%&=\textit{l} \cdot D_{KL}(f_{j+1}(\cdot|\theta_{j+1})\|f_{j}(\cdot|\theta_{j})) \nonumber \\
%&=\textit{l} \cdot d_{+}(\theta_{j}, \theta_{j+1})
%\label{demopositive}
%\end{align}
%
%\item when the sign is `$-$':
%\begin{align}
%&D(\mathbf{m^{(k,n)}}\|\mathbf{m^{(k,n)}-\textit{l}\cdot e_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})=\nonumber\\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\
%&\cdot\ln\left(\dfrac{L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}{L(\mathbf{m^{(k,n)}-\textit{l}\cdot e_{j}^{(k)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}})}\right)\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot\left\lbrace\rule{0cm}{1cm}\left[\sum_{i=1}^{m_{1}} \ln(f_{1}(x_{i}|\theta_{1}))\right]\right.+\nonumber \\
%&+\left.\left[\sum_{i=m_{1}+1}^{m_{2}} \ln(f_{2}(x_{i}|\theta_{2}))\right]+\ldots+\left[\sum_{i=m_{j-1}+1}^{m_{j}-\textit{l}} \ln(f_{j}(x_{i}|\theta_{j}))\right]+\right.\nonumber\\ 
%&\left.+\left[\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \ln(f_{j}(x_{i}|\theta_{j}))\right]+\left[\sum_{i=m_{j}+1}^{m_{j+1}} \ln(f_{j+1}(x_{i}|\theta_{j+1}))\right]+\right.\nonumber\\
%&+\left.\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \ln(f_{j+2}(x_{i}|\theta_{j+2}))\right]+\ldots+\left[\sum_{i=m_{k}+1}^{\textbf{n}} \ln(f_{k+1}(x_{i}|\theta_{k+1}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=1}^{m_{1}} \ln(f_{1}(x_{i}|\theta_{1}))\right]-\left[\sum_{i=m_{1}+1}^{m_{2}} \ln(f_{2}(x_{i}|\theta_{2}))\right]-\ldots-\right.\nonumber\\
%&\left.-\left[\sum_{i=m_{j-1}+1}^{m_{j}-\textit{l}} \ln(f_{j}(x_{i}|\theta_{j}))\right]-\left[\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \ln(f_{j+1}(x_{i}|\theta_{j+1}))\right]-\right.\nonumber\\
%&\left.-\left[\sum_{i=m_{j}+1}^{m_{j+1}} \ln(f_{j+1}(x_{i}|\theta_{j+1}))\right]-\left[\sum_{i=m_{j+1}+1}^{m_{j+2}} \ln(f_{j+2}(x_{i}|\theta_{j+2}))\right]-\ldots-\right.\nonumber\\
%\displaybreak \nonumber\\
%&-\left.\left[\sum_{i=m_{k}+1}^{\textbf{n}} \ln(f_{k+1}(x_{i}|\theta_{k+1}))\right]\rule{0cm}{1cm}\right\rbrace\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\ 
%&=\int  L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber\\ &\left[\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} \ln\left(\dfrac{f_{j}(x_{i}|\theta_{j})}{f_{j+1}(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\mathlarger{\sum_{i=m_{j}-\textit{l}+1}^{m_{j}}} \int L(\mathbf{m^{(k,n)}}|\mathbf{\boldsymbol\theta^{(k+1,n)}}, \mathbf{X^{(n)}}) \cdot \nonumber  \\ &\cdot \left[\ln\left(\dfrac{f_{j}(x_{i}|\theta_{j})}{f_{j+1}(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}\mathbf{X^{(n)}} \nonumber \\
%&=\mathlarger{\sum_{i=m_{j}-\textit{l}+1}^{m_{j}}} \Bigg\lbrace \left\lbrace \int f_{j}(x_{i}|\theta_{j}) \cdot  \left[\ln\left(\dfrac{f_{j}(x_{i}|\theta_{j})}{f_{j+1}(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \cdot \nonumber \\
%& \idotsint f_{1}(x_{1}|\theta_{1}) \cdot \ldots \cdot f_{j}(x_{m_{j}-\textit{l}-1}|\theta_{j}) \cdot \ldots \cdot f_{j}(x_{i-1}|\theta_{j}) \cdot \nonumber\\& \cdot f_{j}(x_{i+1}|\theta_{j})  \cdot \ldots \cdot f_{j+1}(x_{m_{j}+2}|\theta_{j+1}) \cdot \ldots \cdot f_{k+1}(x_{\mathbf{n}}|\theta_{k+1}) \,\mathrm{d}x_{1} \ldots\nonumber \\ &\,\mathrm{d}x_{m_{j}-\textit{l}-1}\ldots\,\mathrm{d}x_{i-1}  \,\mathrm{d}x_{i+1}\ldots\,\mathrm{d}x_{m_{j}+2}\ldots\,\mathrm{d}x_{\mathbf{n}}\Bigg\rbrace \nonumber \\
%&=\mathlarger{\sum_{i=m_{j}-\textit{l}+1}^{m_{j}}}\left\lbrace 1^{n-1} \cdot \int f_{j}(x_{i}|\theta_{j}) \cdot \left[\ln\left(\dfrac{f_{j}(x_{i}|\theta_{j})}{f_{j+1}(x_{i}|\theta_{j+1})}\right)\right]\,\mathrm{d}x_{i}\right\rbrace \nonumber \\
%&=\sum_{i=m_{j}-\textit{l}+1}^{m_{j}} D_{KL}(f_{j}(x_{i}|\theta_{j})\|f_{j+1}(x_{i}|\theta_{j+1})) \nonumber \\
%&=\textit{l} \cdot D_{KL}(f_{j}(\cdot|\theta_{j})\|f_{j+1}(\cdot|\theta_{j+1})) \nonumber \\
%&=\textit{l} \cdot d_{-}(\theta_{j}, \theta_{j+1})
%\label{demonegative}
%\end{align}
%
%\end{itemize}
%
%From expressions \eqref{demopositive} and \eqref{demonegative}, we see that the lemma holds.
%\end{proof}
%\end{comment}
\appendix





\end{document}

1702.05462.extracted_bib
[{Barry and Hartigan(1992)]{Barry1992 Barry, D., Hartigan, J.~A., 1992. {Product Partition Models for Change Point Problems. {The Annals of Statistics 20~(1), 260--279.
1702.05462.extracted_bib
[{Berger et~al.(2012)Berger, Bernardo, and Sun]{BBS2012 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.
1702.05462.extracted_bib
[{Berger and Pericchi(1996)]{Intrinsic 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.
1702.05462.extracted_bib
[{Berk(1966)]{berk1966 Berk, R.~H., 1966. {Limiting Behavior of Posterior Distributions when the Model is Incorrect. The Annals of Mathematical Statistics 37~(1), 51--58.
1702.05462.extracted_bib
[{Carlin et~al.(1992)Carlin, Gelfand, and Smith]{Carlin 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.
1702.05462.extracted_bib
[{Chernoff and Zacks(1964)]{Chernoff1964 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.
1702.05462.extracted_bib
[{Chib(1998)]{Chib Chib, S., 1998. {Estimation and Comparison of Multiple Change-point Models. Journal of Econometrics 86~(2), 221--241.
1702.05462.extracted_bib
[{Fearnhead and Liu(2007)]{Fearnhead2007 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.
1702.05462.extracted_bib
[{Gir{\'{on et~al.(2007)Gir{\'{on, Moreno, and Casella]{Giron Gir{\'{on, F.~J., Moreno, E., Casella, G., 2007. {Objective Bayesian Analysis of Multiple Changepoints for Linear Models (with discussion). In: Bernardo, {\relax J.M., Bayarri, {\relax M.J., Berger, {\relax J.O., Dawid, {\relax A.P., Heckerman, D., Smith, {\relax A.F.M., West, M. (Eds.), Bayesian Statistics 8. Oxford University Press, London, pp. 227–--252.
1702.05462.extracted_bib
[{Hannart and Naveau(2009)]{Hannart2009 Hannart, A., Naveau, P., 2009. {Bayesian Multiple Change Points and Segmentation: Application to Homogenization of Climatic Series. Water Resources Research 45~(10), {W10444.
1702.05462.extracted_bib
[{Heard and Turcotte(2017)]{Heard2017 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.
1702.05462.extracted_bib
[{Henderson and Matthews(1993)]{Henderson1993 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.
1702.05462.extracted_bib
[{Kass and Raftery(1995)]{KassRaftery Kass, R.~E., Raftery, A.~E., 1995. {Bayes Factors. Journal of the American Statistical Association 90~(430), 773--795.
1702.05462.extracted_bib
[{Ko et~al.(2015)Ko, Chong, and Ghosh]{ko2015 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.
1702.05462.extracted_bib
[{Koop and Potter(2009)]{Koop Koop, G., Potter, S.~M., 2009. {Prior Elicitation in Multiple Change-Point Models. International Economic Review 50~(3), 751--772.
1702.05462.extracted_bib
[{Kullback and Leibler(1951)]{kullback1951 Kullback, S., Leibler, R.~A., 1951. {On Information and Sufficiency. The Annals of Mathematical Statistics 22~(1), 79--86.
1702.05462.extracted_bib
[{Loschi and Cruz(2005)]{Loschi2005 Loschi, {\relax R.H., Cruz, {\relax F.R.B., 2005. {Extension to the Product Partition Model: Computing the Probability of a Change. Computational Statistics \& Data Analysis 48~(2), 255 -- 268.
1702.05462.extracted_bib
[{Moreno et~al.(2005)Moreno, Casella, and Garcia-Ferrer]{Moreno2005 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.
1702.05462.extracted_bib
[{Muliere and Scarsini(1985)]{Muliere1985 Muliere, P., Scarsini, M., 1985. {Change-point Problems: A Bayesian Nonparametric Approach. Aplikace matematiky 30~(6), 397--402.
1702.05462.extracted_bib
[{Petrone and Raftery(1997)]{Petrone1997 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.
1702.05462.extracted_bib
[{Raftery and Akman(1986)]{Raftery Raftery, A.~E., Akman, V.~E., 1986. {Bayesian Analysis of a Poisson Process with a Change-Point. Biometrika 73~(1), 85--89.
1702.05462.extracted_bib
[{Schwaller and Robin(2017)]{Schwaller2017 Schwaller, L., Robin, S., 2017. {Exact Bayesian Inference for Off-line Change-point Detection in Tree-structured Graphical Models. Statistics and Computing 27~(5), 1331--1345.
1702.05462.extracted_bib
[{Smith(1975)]{Smith1975 Smith, A. F.~M., 1975. {A Bayesian Approach to Inference about a Change-Point in a Sequence of Random Variables. Biometrika 62~(2), 407--416.
1702.05462.extracted_bib
[{Stephens(1994)]{Stephens1994 Stephens, D.~A., 1994. {Bayesian Retrospective Multiple-Changepoint Identification. Journal of the Royal Statistical Society. Series C (Applied Statistics) 43~(1), 159--178.
1702.05462.extracted_bib
[{Tian et~al.(2009)Tian, Ng, Li, and Tan]{Tian2009 Tian, G.-L., Ng, K.~W., Li, K.-C., Tan, M., 2009. {Non-iterative Sampling-based Bayesian Methods for Identifying Changepoints in the Sequence of Cases of Haemolytic Uraemic Syndrome. Computational Statistics and Data Analysis 53~(9), 3314--3323.
1702.05462.extracted_bib
[{Villa and Walker(2015{\natexlab{a)]{VillaMass Villa, C., Walker, S., 2015{\natexlab{a. {An Objective Approach to Prior Mass Functions for Discrete Parameter Spaces. Journal of the American Statistical Association 110~(511), 1072--1082.
1702.05462.extracted_bib
[{Villa and Walker(2015{\natexlab{b)]{VillaModel Villa, C., Walker, S., 2015{\natexlab{b. {An Objective Bayesian Criterion to Determine Model Prior Probabilities. Scandinavian Journal of Statistics 42~(4), 947--966.
1702.05462.extracted_bib
[{Yu(2001)]{Yu 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.