Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
143,391 characters · 18 sections · 118 citation commands
Supplementary Appendix for 'Selective linear segmentation for detecting relevant parameter changes'
\onehalfspacing \justifying
\pagenumbering{arabic}
Long time series are standard in this period of large publicly available datasets. Care is required when modeling such a time series since many of them span over critical events that may change the series dynamic. At least two statistical solutions exist to take into account these changes. On the one hand, a process with fixed parameters can be used but it needs to exhibit a rich and complex dynamic. This complexity often makes the model difficult to estimate and to interpret (see, for instance, long memory processes such as GPH83). On the other hand, one can rely on time-varying parameter (TVP) models and in particular Markov-switching and change-point (CP) processes since they allow for abrupt changes in the model parameters when a critical event affects the series dynamic Hamilton89,Bauwens2013. This paper deals with CP linear regression models where we allow the mean parameters to change over time. \\ The CP literature dates back to chernoff1964estimating and is nowadays vast. Just focusing on linear regressions, andrews1993tests, BaiPerron1998, killick2012optimal, fryzlewicz2014wild and yau2016inference develop prominent procedures to detect breakpoints. On the Bayesian side, there also exist many ways to estimate structural breaks and important contributions can be found in Stephens1994, Chib98, Fearnhead2007, rigaill2012exact and Maheu13. While all these methods differ in the criterion or in the algorithm used to detect the changes, most of them rely on the assumption that, when a break is detected (that may be triggered by the change in only one model parameter), a new segment is created and a new set of parameters needs to be estimated. Although the assumption seems harmless, it creates two important drawbacks:
In this paper, we propose a method to relax the assumption that a break triggers a change in all the model parameters. To do so, we first estimate the potential break dates exhibited by the series and then we use a penalized likelihood approach to detect which parameters change. Because some segments in the CP regression can be small, we opt for a (nearly) unbiased penalty function, called the seamless-L0 (SELO) penalty function, recently proposed by dicker2013variable. We prove the consistency of the SELO estimator in detecting which parameters indeed vary over time and we suggest using a deterministic annealing expectation-maximisation (DAEM) algorithm to deal with the multimodality of the objective function ueda1998deterministic. Since the SELO penalty function depends on two tuning parameters, we use a criterion (new in this literature) to choose the best tuning parameters and as a result the best model. This new criterion exhibits a Bayesian interpretation which makes possible to assess the parameters' uncertainty as well as the model's uncertainty. This last feature is determinant when predicting a time series since the Bayesian model averaging technique, that typically improves forecast accuracy, is readily applicable raftery2010online,koop2012forecasting.
We are aware of five other papers that also relax the assumption on the number of parameters that changes when a break is detected. In the frequentist literature, the influential paper of BaiPerron1998 proposes a method that also operates when only a subset of parameters can break. However, the number of possibilities grows exponentially with the number of breaks as well as with the number of parameters that can break. From a Bayesian perspective, GiordaniKohn2008, Eo2012, huber2019should and dufays2019relevant propose flexible state-space models to capture which parameters vary over time. However, all these estimation procedures break down when the number of parameters is large (see Supplementary Appendix (ref) for more details).
We believe that our method exhibits several advantages over the existing alternatives. Firstly, it operates for small and large dimensions. Secondly, the estimation is fast compared to the Bayesian alternatives. As a final advantage, we relax the assumption on breakpoints once the structural breaks have been detected which makes our approach operating in combination with any existing CP methods. In this paper, we illustrate our approach with the CP procedure of yau2016inference but any other CP method could have been used.
A final reference close to our framework is chan2014group who propose a penalized regression for segmenting time series in piecewise linear models. The paper uses a group Lasso penalty function yuan2006model to get an overestimated number of segments and in a second phase, an information criterion is used to improve the estimation. Nevertheless, we differ from their methods in many aspects. First, we use an almost unbiased penalty function and from a theoretical perspective, as we use the penalized regression on a potential break date set, our assumptions for a consistent estimator are different and in line with the standard penalized regression literature. We also use a Bayesian criterion to select among the promising models uncovered by the penalty function which allows for model uncertainty and for Bayesian model averaging. Also, our estimation procedure is fast compared to chan2014group since we iterate on closed-form expressions and because our model exhibits fewer parameters. As a final difference, we provide break uncertainty.
Eventually, we apply our method on Hedge funds (HF) returns. As highlighted by Fung2008hedgefunds, by meligkotsidou2008detecting, by Bollen2009hedgefundriskdynamics, and more recently by Patton2015change, the dynamics of HF risk exposures and the nonlinear generating process of HF returns should be associated with market events and structural breaks. For a sample of 14 monthly Credit Suisse HF indices spanning from March 1994 to March 2016, and using the asset based style pricing model introduced by fung2001risk, we show that our modeling is particularly appealing to detect time-varying exposures in HF tradings. In particular, our results report the relative role played by static and dynamic parameters and factors in the decomposition of HF returns. We also investigate the prediction performance of our approach and it turns out that the selective segmentation approach compares favorably in terms of root mean squared forecast errors and cumulative log-predictive densities with respect to other CP processes. In particular, it almost systematically dominates the CP model which assumes that all the parameters vary when a break is detected.
The paper is organized as follows. Section (ref) documents the model specification and the SELO penalty function. Section (ref) explains how the DAEM algorithm is applied to our framework. In Section (ref), we detail the criterion used to select the SELO tuning parameters and we relate it to the Bayesian paradigm. Section (ref) documents the CP method of yau2016inference and discusses how it can be slightly improved. An extensive Monte Carlo study is proposed in Section (ref). We end the paper by applying the method on HF returns in Section (ref). All the proofs are given in the Supplementary Appendix (SA).
We consider a standard linear regression specified as
where $ \epsilon_t \sim MDS(0,\sigma^2)$ (in which $MDS$ stands for the martingale difference sequence), $\mathbf x_t = (1, x_{t,2}, \ldots, x_{t,K})'$ and $\boldsymbol{\beta}_1 = (\beta_1, \beta_2, \ldots, \beta_{K})'$. Typically, if a linear model is estimated over a long period, the parameters are subject to abrupt changes over time. To take this time-varying dynamic into account, we allow for $m-1$ structural breaks in the model parameters as follows,
in which $\boldsymbol{\beta}_i^*$, is the true parameter of the explanatory variables over the regime $i$, $\boldsymbol{\tau}_0 = \{\tau_0,\ldots,\tau_m\} \in \mathbb{N}^{m+1}$ where $\tau_0=0$, $\tau_m = T$ and $\tau_i < \tau_{i+1}$ $\forall i \in [0,m-1]$. In this paper, we are interested in capturing which parameters are subject to breaks and which do not vary over time. To do so, we reframe the model (ref) as follows,
where $\mathbf 1_{\{ {x>a} \}}=1$ if $x>a$ and zero otherwise, $\Delta \boldsymbol{\beta}_j^* = \boldsymbol{\beta}_{j}^* - \boldsymbol{\beta}_{j-1}^*$, for $j\in [2,m]$, stands for the model parameters in first-difference, $\mathbf y = (y_1,\ldots,y_T)'$, $\mathbf X_{\boldsymbol{\tau}} = (\tilde{\mathbf X}_{\tau_0}, \tilde{\mathbf X}_{\tau_1},\ldots,\tilde{\mathbf X}_{\tau_{m-1}})$ with $\tilde{\mathbf X}_{\tau_i} = (\mathbf 0,\mathbf 0,\ldots,\mathbf 0,\mathbf x_{\tau_i+1},\ldots,\mathbf x_{T})'$, $\boldsymbol{\epsilon} = (\epsilon_1,\epsilon_2,\ldots,\epsilon_T)'$ and $\boldsymbol{\beta}^*=(\boldsymbol{\beta}_1^{*\prime},\Delta \boldsymbol{\beta}_2^{*\prime},\ldots,\Delta \boldsymbol{\beta}^{*\prime}_m)^{\prime} \in \Re^{Km\times 1}$. Note that the matrix $\tilde{\mathbf X}_{\tau_0}$ stands for the standard regressors since $\tau_0=0$. Regarding the notations, the first-difference parameter in regime $j$ is a K-dimensional vector $\Delta \boldsymbol{\beta}_{j}^*$ such that $\Delta \boldsymbol{\beta}_{j}^* =(\Delta \beta_{j1}^*,\dots, \Delta \beta_{jK}^*)^{\prime}$. Let us also denote $A = \{(j,k); \Delta \beta_{jk}^*\neq 0, \text{ for } j\in[2,m] \text{ and for } k\in[1,K]\}$, the set of indices defining the true model.
Our strategy to uncover which parameters truly vary over time consists in first finding where are the potential break dates $\boldsymbol{\tau}$, then, in a second phase, in detecting which parameters evolve. Note that even when we know the true break dates $\boldsymbol{\tau}$, the problem of finding which parameters vary when a break occurs is not straightforward as the number of models to consider amounts to $2^{(m-1)K}$. Consequently, it is infeasible to carry out an exhaustive model selection when $K$ or $m$ is large. We propose a penalized likelihood approach to explore this large model space and to select which parameters experience breaks. To focus on our selective segmentation approach, we shall first assume that we have obtained a set of potential break dates $\boldsymbol{\tau}$. We discuss how we estimate this set in Section (ref).
As emphasized by Equation (ref), given a set of break dates $\boldsymbol{\tau}$, the problem of finding which parameters abruptly change when a break occurs boils down to a penalized linear regression problem. Specifically, one can solve the following optimization problem
where $||~ . ~||_p$ denotes the $L_p$ norm and $\text{pen}(\Delta \beta_{jk})$ stands for a penalty function. Popular choices of $\text{pen}(\Delta \beta_{jk})$ are the Lasso penalty function (i.e., $\text{pen}(\Delta \beta_{jk}) = \lambda ||\Delta \boldsymbol{\beta}_{jk}||_1$, see Tibshirani94) or the rigde function (i.e., $\text{pen}(\Delta \beta_{jk}) = \lambda \Delta ||\beta_{jk}||_2^2$, see, for instance, Ishwaran2005).
Following Fan2001, standard desirable properties induced by a penalty function are i) unbiasdness, ii) sparsity and iii) continuity. For instance, the ridge function is only continuous while the Lasso penalty function achieves sparsity and continuity (beside at zero). However one standard issue with these popular penalty functions is that they provide biased (but typically consistent) estimators. In our framework, this drawback is problematic since a segment can sometimes contain a small amount of observations that makes consistency results not sufficient. Recently, dicker2013variable propose a penalty function, called seamless-$L_0$ (SELO), that exhibits all the desirable properties. For a model parameter denoted $\omega$, the penalty function reads as
where the parameter $\zeta$ controls for the concavity of the function and $\lambda$ stands for the penalty imposed when $\omega \neq 0$. We slightly modify their function to end up with parameters that are directly interpretable. In fact, we use the following penalty function,
where $\zeta = \frac{2^y-2}{1-2^y}$ with $y \in (0,1)$ and we set $y=0.99$. In most cases, the parameter $a$ can be interpreted as an interval $\omega \in [-a,a]$ in which $\omega$ will be biased with respect to the OLS estimate since $\mathcal{P}_{\mathbf{SELO}}(a) = \lambda y$. Intuitively, when $|\omega|>a$, we have $\mathcal{P}_{\mathbf{SELO}}(\omega) \approx \lambda$ and $\frac{d\mathcal{P}_{\mathbf{SELO}}(\omega)}{d\omega}|_{|\omega|\geq a} \approx 0$ for large values of $a$. Figure (ref) shows the SELO penalty function with $\{a,\lambda\}=\{1,0.9\}$ and illustrates that the function is almost flat for absolute values greater than $a$. To be more precise about how large $a$ must be, when $|\bar{\beta}|\geq a$ with $a\geq \frac{\zeta}{\ln 2[\zeta^2 + 3\zeta + 2]} = 0.0099$, we have that $\frac{d\mathcal{P}_{\mathbf{SELO}}(\omega)}{d\omega}|_{|\omega|\geq a} \leq \lambda$ which implies that the bias imposed by the SELO penalty function is smaller than the one of the Lasso function (i.e. $\lambda |\omega|)$.
The interval $[-a,a]$ in which a parameter is biased is likely to change with the variable to which it refers. Furthermore, if we assume that this interval is fixed over time, we should set a new parameter $a$ for each variable on the $m$ segments. Unlike dicker2013variable who define a single parameter for all the variables, we use $K$ parameters $a_1,\dots,a_K$, that is, one per explanatory variable. Thus, the objective function to minimize is given by
Before discussing how to maximize the objective function, we present the main results about the modified SELO estimator. As highlighted in dicker2013variable, the SELO estimator is consistent under reasonable conditions. Proposition (ref) shows that this consistency result also applies in our framework. To do so, we consider the following assumptions (in which a sequence $\omega_T \rightarrow \omega $ is understood as $\text{lim}_{T \rightarrow \infty} \omega_T = \omega$):
We first discuss the assumptions before detailing our consistency result. Assumptions (ref) to (ref) are similar to those found in the variable selection literature fan2004nonconcave,dicker2013variable and in the CP literature BaiPerron1998,yau2016inference. Condition (ref) assumes that the estimated CPs are the true locations. However, the SELO estimator maintains the same asymptotic properties with a set of potential breakpoints as long as it contains the true break dates (see the adapted assumption (ref) below). In such case, Proposition (ref) also ensures that the number of breakpoints is consistently estimated. Note that condition (ref) implies that the length of each segment increases linearly with $T$. Although unattractive, this condition is generally made in the CP literature perron2006dealing,yau2016inference. For interested readers, perron2006dealing motivate this assumption in details. Assumption (ref) allows the minimum break size to decrease with the sample size but at a slower rate than $T^{-\frac{1}{2}}$. Conditions (ref) are related to the eigenvalues and are standard in the variable selection literature zhang2010regularization. However, we show in SA (ref) that this condition is not innocuous and that it implies a fixed number of regimes as well as $\min_j\{\delta_{\tau_j}\}>0$; that is $\delta_{\tau_j}$ does not drift to 0 as $T \rightarrow \infty$. Avoiding this assumption would imply stronger conditions on the process $\{y_t,\mathbf x_t\}$ chan2014group. The assumption (ref) refers to ergodicity and stationarity of each segment and imposes the standard exogeneity hypothesis. This assumption ensures that sampled counterparts of the first two moments of $\{\mathbf x_t\epsilon_t\}$ are converging to finite values. Importantly, it does not rule out conditional heteroskedasticity. Eventually, condition (ref) defines restrictions on the tuning parameters rate. The same condition applies in dicker2013variable. The consistency of SELO estimator is given by the following Proposition.
The objective function to minimize is given by
in which $q_{k}(\Delta \beta_{jk}) = \left(\frac{2(\frac{|\Delta \beta_{jk}|}{a_k}) + \zeta}{(\frac{|\Delta \beta_{jk}|}{a_k}) + \zeta}\right)^{\left(\frac{T \lambda}{\ln(2)}\right)}$. Due to the penalty function, we cannot find any analytical expression of the minimizer. In addition to that, the function likely exhibits many local modes which complicates the optimization. We address the problem of finding the global mode by using a deterministic annealing expectation-minimization (DAEM) algorithm ueda1998deterministic. To do so, we first approximate the penalty function by a mixture of two Normal components (to take into account the large tail of the SELO penalty function), the details of it are given in SA (ref). Secondly, since minimizing the sum of squared residuals is identical to maximizing a likelihood function when the error term is normally distributed, we work with the following model, $\mathbf y = \mathbf X_{\boldsymbol{\tau}} \boldsymbol{\beta} + \eta$, where $\eta \sim \mathcal{N}(\mathbf 0, \sigma^2 I_T)$. The modified model implied the following objective function to maximize with respect to $(\boldsymbol{\beta},\sigma^2)$:
where $g_k(\Delta \beta_{jk}) = \sum_{i=1}^{2}\omega_{i}^{(k)} f_N(\Delta \beta_{jk}|\mu_{i}^{(k)},s_{i}^{(k)})$, $f_N(x|\mu,s)$ stands for the normal density function evaluated at $x$ with expectation and variance given by $\mu$ and $s$ respectively and $\omega_{i}^{(k)} \in (0,1)$ with $\sum_{i=1}^{2}\omega_{i}^{(k)} = 1$. Note that the function $f(\mathbf y|\boldsymbol{\beta},\sigma^2)$ in Equation (ref) is proportional to the posterior log-density of the parameter distribution $\boldsymbol{\beta},\sigma^2|\mathbf y$ from a Bayesian perspective with prior distributions given by $f(\sigma^2,\boldsymbol{\beta}_1) \propto 1$ and $f(\Delta \beta_{jk}) = g_k(\Delta \beta_{jk})$ for $j\in[2,m]$ and $k\in [1,K]$. In particular, the distribution $f(\Delta \beta_{jk})$ can be understood as a spike and slab prior George_1993 and our optimization procedure fits into the framework of rovckova2014emvs which proposes tackling the linear variable selection problem with the EM algorithm and its DAEM variant. The optimization is therefore equivalent to finding the mode of $\boldsymbol{\beta},\sigma^2|\mathbf y$. Using a data augmentation approach, we add latent variables $\mathbf z = (z_{21},z_{22},\ldots,z_{mK})'$ such that $f(z_{jk}=i) = \omega_i^{(k)}$, $\forall j \in [2,m]$, $\forall k \in [1,K]$ and $\forall i \in [1,2]$. With these latent variables, we can write the prior distribution of $\Delta \beta_{jk}$ in a convenient hierarchical way as follows,
By fixing $\boldsymbol{\theta} = \{\boldsymbol{\beta},\sigma^2\}$, the EM algorithm (and its DAEM variant) solves the following optimization at iteration $n$,
One can easily show that maximizing $Q(\boldsymbol{\theta}_{n}|\boldsymbol{\theta}_{n-1})$ implies that $f(\boldsymbol{\theta}_n|\mathbf y) \geq f(\boldsymbol{\theta}_{n-1}|\mathbf y)$.
To apply the DAEM algorithm, we need to find an expression of $Q(\boldsymbol{\theta}|\boldsymbol{\theta}_{n-1})$. Given a set of parameter $\boldsymbol{\theta}_{n-1}$, we have that
where
with $p_{jk}^{(i)} = f(z_{jk}=i|\mathbf y,\boldsymbol{\theta}_{n-1})$ $\forall i \in [1,2]$,$\forall j \in [2,m]$ and $\forall k \in [1,K]$ . Importantly, the difference between the EM algorithm and its DA version only appears in the quantities $p_{jk}^{(i)}$. In fact, the DAEM algorithm introduces an increasing function $\phi(r):[1,N]\rightarrow (0,1]$ such that $0<\phi(1) \leq 1$ and $\phi(N) = 1$. For each value $r=1, \ldots, N$, it applies recursively the EM algorithm (that starts with the final estimate of the previous EM algorithm) where the posterior probabilities $p_{jk}^{(i)}$ are denoted $p_{jk}^{(i,\phi(r))}$ and are modified as follows,
When $r=N$, the increasing function $\phi(r)=1$ and the standard EM algorithm is run (but with a promising starting point). To find the maximum of $Q(\boldsymbol{\theta}|\boldsymbol{\theta}_{n-1})$, we sequentially maximize $\boldsymbol{\beta}$ given $\sigma^2$ and then $\sigma^2$ with respect to $\boldsymbol{\beta}$. This approach, called coordinate iterative ascent, operates in two steps:
At the end of the two steps, we necessarily have $Q(\boldsymbol{\beta}_{n-1},\sigma_{n-1}^2|\boldsymbol{\theta}_{n-1}) \leq Q(\boldsymbol{\beta}_{n},\sigma_{n-1}^2|\boldsymbol{\theta}_{n-1}) \leq Q(\boldsymbol{\beta}_{n},\sigma_{n}^2|\boldsymbol{\theta}_{n-1})$. The maximisation of $\boldsymbol{\beta}$ given $\sigma_{n-1}^2$ leads to
The update of $\sigma^2$ conditional to $\boldsymbol{\beta}_n$ is given by
We summarize the DAEM procedure in Algorithm (ref). In practice, the minimum distance $e$ indicating a convergence of the algorithm is set to $10^{-5}$ and the number of DAEM iteration $N$ is fixed to 10.
The EM and the DAEM algorithms are sensitive to starting values. Inspired by zhao2012random, we mitigate this issue by randomly exploring the model space using a swapping approach before applying the DAEM algorithm. To be specific, we generate $N_{\text{init}}$ values as explained in Algorithm (ref) and we initialize the DAEM algorithm with the parameter estimates that minimize the penalized function given in Equation (ref). In practice, we set $N_{\text{init}} = \text{min}(2^{(m-1)K-1},3000)$.
The SELO penalty function exhibits two tuning parameters $\boldsymbol{a}$ and $\lambda$. The standard approach to fix them consists in considering a grid of values of these parameters and in selecting the parameters that maximize a (generally consistent) information criterion zhang2010regularization. Instead of relying on a standard information criterion and select the tuning parameters $\boldsymbol{a}$ and $\lambda$ that maximize it, we consider each pair $(\boldsymbol{a},\lambda)$ as a model to take into account the model uncertainty. For a given value of $(\boldsymbol{a},\lambda)$, the DAEM algorithm exposed in Section (ref) provides an estimate $\hat{\Delta \boldsymbol{\beta}}$ of $\Delta \boldsymbol{\beta}$ which delivers an estimate of $\hat{A}$, i.e., the set of indices with $\hat{\Delta \beta_{jk}} \neq 0$ for $ j\in [2,m]$ and for $k \in [1,K]$. This set tells us which covariates should be included in the linear regression and which should not. Let us denote by $\tilde{\mathbf X}_{\boldsymbol{\tau}}^{\hat{A}}$ the covariates related to the first-difference estimates that are different from zero. We use the following criterion for selecting $\boldsymbol{a}$ and $\lambda$:
where $s_{\tilde{\mathbf X}_{\tau_0}}$ stands for the residual sum of squares (RSS) from the ordinary least squares (OLS) with $\mathbf X=\tilde{\mathbf X}_{\tau_0}$ (i.e., a regression without break), $s_{\tilde{\mathbf X}_{\tau_0},\tilde{\mathbf X}_{\boldsymbol{\tau}}^{\hat{A}}}$ is the RSS from the OLS with $\mathbf X=(\tilde{\mathbf X}_{\tau_0},\tilde{\mathbf X}_{\boldsymbol{\tau}}^{\hat{A}})$, the value $k_{\hat{A}} = |\hat{A}|$ denotes the number of first-difference parameters different from zero in the model and $g_{\hat{A}}$ is a user parameter. We properly derive the criterion in SA (ref). fernandez2001benchmark show that the criterion (ref) is consistent in the sense that it selects asymptotically the true subset of regressors when $g_{\hat{A}} = w(T)^{-1}$ as stated in proposition (ref). \\
In fernandez2001benchmark, they advocate for setting $g_{\hat{A}} = \text{min}(T^{-1},(k_{\hat{A}}+K)^{-2})$ as this prior empirically delivers good results for selecting the true covariates in standard linear regressions. However, we deviate from this benchmark prior by fixing $g_{\hat{A}} = \frac{1}{T^{\alpha}-1}$ with $\alpha=1$ when $k_{\hat{A}}=0$ and $\alpha = \frac{k_{\hat{A}} + \hat{m}_{\hat{A}}-1}{k_{\hat{A}}}>1$ when $k_{\hat{A}}>0$ in which $\hat{m}_{\hat{A}}$ denotes the number of active segments. When $\alpha>1$, we show in SA (ref) that the criterion in Equation (ref) asymptotically converges in probability to
The asymptotic value is equivalent to the Bayesian information criterion (BIC) of a linear regression model exhibiting a number of parameters of $\alpha k_{\hat{A}}$.\footnote{The BIC of a linear regression model with K parameters is given by $-\frac{T}{2} \ln(\frac{s_{\tilde{\mathbf X}_{\tau_0},\tilde{\mathbf X}_{\boldsymbol{\tau}}^{\hat{A}}}}{T}) -\frac{K}{2} \ln T$. So the marginal likelihood criterion of Equation (ref) converges to the BIC up to an additive constant (that is $\frac{T}{2} \ln T$).} Consequently, the model penalty takes additionally into account the number of active breakpoints when $\alpha = \frac{k_{\hat{A}} + \hat{m}_{\hat{A}}-1}{k_{\hat{A}}}$. This stronger penalty works empirically well and is motivated by several CP papers advocating for stronger penalties than the BIC as it tends to overfit the number of regimes in finite sample LiuAl1997,zhang2007modified,kim2016consistent.
Interestingly, criterion (ref) stands for a marginal likelihood in the Bayesian paradigm under $\epsilon \sim \mathcal{N}(0,\sigma^2 I_T)$ and the following prior,
where $\mathbf M_{\tilde{\mathbf X}_{\tau_0}} = I_T - \tilde{\mathbf X}_{\tau_0}((\tilde{\mathbf X}_{\tau_0})'\tilde{\mathbf X}_{\tau_0})^{-1}(\tilde{\mathbf X}_{\tau_0})'$. The prior distributions given by Equations (ref) lead to simple posterior inference. The posterior distribution of the model parameters are given by, see SA (ref) for derivations,
in which $\mathcal{IG}(-,-)$ denotes the Inverse-Gamma distribution. Consequently, we can go beyond selecting the best pair $(\boldsymbol{a}_p,\lambda_p)$ (i.e., the pair that maximizes the criterion (ref)) and can take the uncertainty of this selection into account. Given a set of models $M_z = (\boldsymbol{a}_z,\lambda_z)$, with $z=1,...,Z$, we can directly assess the posterior probability of a specific model as follows
where $f(M_z|\boldsymbol{\tau})$ denotes the prior probability of model $M_z$. In this paper, we assume uninformative prior, so $f(M_z|\boldsymbol{\tau}) = Z^{-1}$. The posterior probability can be used to account for uncertainty on the selected regressors. In fact, we have
It is worth emphasizing that the consistent property of the criterion (ref) does not depend on the normality assumption. Only, the posterior distribution of the model parameters does. We do not see this as a limitation since one can easily extend the model with another distributional assumption and compute the posterior distribution by numerical integrations.
Equation (ref) shows how to take into account the uncertainty of the model parameters with respect to the selection of the SELO parameters. The Bayesian paradigm also provides a simple tool to forecast the series taking this uncertainty into account. In particular, the predictive density $f(y_{T+1:T+h}|\mathbf y)$, for $h\geq 1$, is related to the posterior density as follows
where $\{\boldsymbol{\beta}_1^{(i)},\Delta \boldsymbol{\beta}^{(i)},(\sigma^2)^{(i)},M^{(i)}\}_{i=1}^{N}$ are independent draws from the posterior distribution (i.e., $\boldsymbol{\beta}_1,\Delta \boldsymbol{\beta},\sigma^2,M|\mathbf y,\boldsymbol{\tau}$). From (ref), it is apparent that the predictive density takes the model uncertainty into account.\footnote{Using the full marginal likelihood for weighting the models' predictions could raise concerns as only the last segment matters in CP processes. However, as marginal likelihood is frequently used for selecting the number of regimes in the literature and because it is also informative about the fit of the last regime, this average should give large weights to the models exhibiting a good fit at the end of the sample. Nevertheless, we could also weight the models' predictions using the predictive marginal likelihood $f(y_{t_1+1:T}|y_{1:t_1},\boldsymbol{\tau})$ in which $t_1$ is a user-defined value.} This feature should be contrasted with the standard penalized regression literature in which forecasting is performed using one unique set of parameter estimates; i.e., the estimates given by one penalty parameter selected, for instance, by cross-validation or by an information criterion. \\ In practice, simulations from the posterior distribution are not required for evaluating the predictive density. Assuming that the future covariates $\mathbf x_{T+1:T+h}$ are observed at time $T$, the predictive distribution of $y_{T+1:T+h}$ given a model $M_z$ turns out to be a multivariate student distribution. Supplementary Appendix (ref) documents the analytical expression of $f(y_{T+1:T+h}|\mathbf y,M_z)$. Therefore, we can efficiently take into account model uncertainty in the predictive density since Equation (ref) simplifies into
When the number of models to consider is too large to directly explore the model space using the criterion (ref) (i.e., when $(m-1)K>10$, see remark (ref)), we rely on the SELO penalty function to uncover the promising explanatory variables. While the asymptotic result of Proposition (ref) is reassuring, it only applies if the parameters $\lambda$ and $\boldsymbol{a}$ are adequately chosen. Similar to what is generally done in the penalized regression literature, we propose to explore many values of $\lambda$ and $\boldsymbol{a}$ and consider each couple as a model that would be ultimately discriminated via criterion (ref). For the parameter $\boldsymbol{a}$, we use a value of $a_i = \kappa \times \text{std}(\hat{\beta}_{j1})$ for each $j$ of the $K$ parameters per regime where $\text{std}(\hat{\beta}_{j1})$ stands for the standard deviation of the OLS estimate $\hat{\beta}_{j1}$ when we assume no break in the linear regression (i.e., $\mathbf X = \mathbf X_{\tau_0}$). We test several values for the parameter $\kappa$, namely $\kappa \in \{0.1,1\}$. Regarding the penalty parameter $\lambda$, we test 50 different values uniformly spaced in the interval $(0,\bar{\lambda}]$ in which $\bar{\lambda} = 2\ln T$. The penalty imposes by the upper bound $\bar{\lambda}$ is conservative enough as it is stronger than standard information criteria such as the BIC (that corresponds to a penalty of $\frac{1}{2}\ln T$) and the modified BIC.
In this Section, we present one approach to obtain a set of potential break dates. Before going into details, it is worth emphasizing that our method for detecting which parameters vary when a break occurs is independent of the segmentation detection procedure used in the first phase. To build the break date set, we could, for instance, adapt the dynamic programming method of bai2003computation for the marginal likelihood given by Equation (ref) and therefore propose our own CP detection method. We could also detect the locations of the segments using one of the standard segmentation approaches such as BaiPerron1998, killick2012optimal or korkas2017multiple. Even better, we could apply several CP detection algorithms and discriminate between the sets of breakpoints by comparing their marginal likelihoods once the SELO optimization has been carried out on each set. However, as the emphasis of the paper is not on the break detection, we prefer relying on one break detection procedure, the one documented in yau2016inference, because i) it delivers a set of potential break dates with a computational complexity of $\mathcal{O}(T \left(\log(T)\right)^2)$ (which is faster than $\mathcal{O}(T^2)$, i.e., the complexity of the dynamic programming method of bai2003computation) and because ii) we slightly improve their CP detection procedure. In particular, their estimated breakpoints depend on one tuning parameter, the radius $h$. Instead of fixing it, we use multiple values of $h$ and we also adapt their approach to end up with a potential breakpoint set. \\ It is worth noting that, as the paper combines model selection and CP detection methods, our approach only requires a set of potential break dates that includes the correct break dates. By penalizing the parameter variation between two consecutive regimes, the spurious break dates are consistently deleted (see remarks (ref) and (ref)).
yau2016inference propose a likelihood ratio scan method in three steps for estimating multiple break dates in piecewise stationnary processes. They also establish the consistency of the estimated number and location fractions of the CPs. We apply their three steps to detect the break dates but we modify them to reduce the computational burden and to keep at the end of the procedure a potential break date set (that could overestimate the true number of regimes). We now detail the three steps that we use to segment the data.\\ First step. Fix a window radius $h\in [K+1,T-K]$. For $t=h$ to $T-h$, compute the likelihood ratio scan statistic given by,
where $\hat{L}_{t_1:t_2}(\hat{\boldsymbol{\beta}},\hat{\sigma})$ denotes the maximum value of the log-likelihood of model ((ref)) over the segment $t\in [t_1,t_2]$, assuming that $\epsilon_t \sim \mathcal{N}(0,\sigma^2)$. Then, the set $\Gamma(h)$ of potential break dates is given by,
where $S_h(t) = 0$ for $t<h$ and $t>T-h$. As the window radius $h$ is crucial, we differ from yau2016inference by using a grid of M values uniformly-spaced in the interval $[\frac{h_{\text{YZ}}}{2},2h_{\text{YZ}}]$ in which $h_{\text{YZ}}$ denotes their advocated value that is $h_{\text{YZ}} =\max\left\{25, \left(\log(T)\right)^2\right\}$ when $T < 800$ and $h_{\text{YZ}} = \max\left\{50, 2\left(\log(T)\right)^2\right\}$ otherwise. So, at the end of the first step, we end up with M potential break date sets, i.e., $\Gamma(h_1),\ldots,\Gamma(h_M)$. \\ Second step. For every $z\in [1,M]$ and $i\in [1,m_{h_z}-1]$ where $m_{h_z} = |\Gamma(h_z)|+1$, we re-estimate each break date location $\tau_i^{(z)} \in \Gamma(h_z)$ as follows
in which $\left\lfloor x \right\rceil$ stands for the nearest integer to $x$. Gathering all the new locations in the set $\hat{\Gamma}(h_z) = \{\hat{\tau}_1^{(z)},\ldots,\hat{\tau}_{m_{h_z}}^{(z)}\}$, it is clear from Theorems 1 to 3 in yau2016inference that for any $j \in \{1,\dots,m-1\}$, there exist $\hat{\tau}_i^{(z)} \in \hat{\Gamma}(h_z)$ with $i\in [1,m_{h_z}-1]$ such that $\hat{\tau}_i^{(z)} - \tau_j = \mathcal{O}_p(1)$.
Third step. We select the best breakpoints among the M potential break date sets by minimizing the Minimum Description Length (MDL) defined by, for $z\in [1,M]$,
where $\hat{\tau}_{0}=0$, $\hat{\tau}_{m_{h_z}}=T$ and $\{\hat{\tau}_j\}_{j=2,\dots,m_{h_z}-1}=\hat{\Gamma}(h_z)$. In practice, we fix $M=30$.\\
Given a set of break dates obtained either from the procedure described in Section (ref) or from any other existing break detection method such as the one of BaiPerron1998, our method to uncover the partial structural changes can be undertaken. Let us denote by $M_*=(\boldsymbol{a}_*,\lambda_*)$ the SELO parameters maximizing the marginal likelihood criterion (ref) and their corresponding break dates $J = \{\bar{\tau}_0 =0, \bar{\tau}_1, \ldots, \bar{\tau}_{\hat{m}-1}, \bar{\tau}_{\hat{m}} = T\}$. To provide break uncertainty, we shall infer the posterior distribution of the structural breaks; i.e., $\boldsymbol{\tau} \equiv \tau_1,\ldots,\tau_{\hat{m}-1}|\mathbf y,M_*$. To do so, we first assume uninformative priors for the break dates using the set $J$. For $i=1,...,\hat{m}-1$, the break parameter $\tau_i$ is driven by a Uniform distribution as follows
in which $\left\lfloor x \right\rfloor$ stands for the nearest integer less than or equal to $x$ and $\gamma=(K+1)$ is a minimum duration parameter ensuring that the marginal likelihood criterion (ref) can be computed for any break parameters complying with the prior distributions given by Equation (ref). The posterior density is proportional to
As shown in SA (ref), the marginal likelihood $f(\mathbf y|M_*,\boldsymbol{\tau})=f(\mathbf y|\boldsymbol{a}_*,\lambda_*,\boldsymbol{\tau})$ exhibits a closed form expression. Several solutions exist to sample the break parameters Stephens1994,Liao08. In this paper, we use the D-DREAM algorithm developed in DREAM11. It builds a symmetric proposal distribution inspired by the Differential Evolution optimization literature and draws from this proposal distribution are accepted or rejected through a Metropolis step in a Markov-chain Monte Carlo (MCMC) algorithm. As shown in DREAM11, the D-DREAM algorithm complexity is $\mathcal{O}(T)$ and leads to a rapidly mixing MCMC algorithm since the break parameters are jointly sampled from the proposal distribution. To infer the break parameters, we apply the following steps:
In practice, we set the number of MCMC iterations to $4000$ and start collecting the draws after $\text{round}[\frac{M}{2}]$ MCMC iterations. In addition, we assess the convergence of the MCMC algorithm using the multivariate Potential Scale Reduction Factor test proposed in brooks1998general. For the two in-sample applications below in which credible intervals of the breakpoints are computed, the convergence statistics amount to 1.014 and 1.052, respectively. These values meet the threshold of 1.1 commonly used to validate the convergence of MCMCs.
In this Section, we document a Monte Carlo study to assess the accuracy of the SELO approach. We first rely on nine different data generating processes (DGPs) that are documented in Table (ref). For each DGP, we simulate 1000 series with a sample size equal to $T=1024$ and we investigate i) the performance of detecting the break dates using the approach in Section (ref) and ii) the performance of the SELO method for detecting which parameter truly varies when a break occurs. The nine DGPs differ in their mean parameter specifications. For each of them, we study the SELO performance when the innovation is either homoskedastic or driven by a GARCH process. \\ Regarding the DGPs, the first six DGPs are piecewise stationary AR models directly taken from yau2016inference while the others cover situations with exogenous explanatory variables. DGP A and E do not exhibit any breakpoint. They aim at showing the performance of the SELO approach when only spurious break dates are detected. DGPs B and C are weakly persistent piecewise stationary AR models exhibiting three regimes. Simulated series from DGP D experience a break after 50 observations. This DGP should highlight the performance in a short regime context. DGPs E and F are highly persistent piecewise stationary AR models but DGP F differs by exhibiting breaks in the mean parameters. Eventually, DGPs G, H and I include exogenous variables. While DGPs G only exhibits exogeneous regressors, DGPs H and I stand for ARX processes by mixing the parameters of the DGPs B and G.
Table (ref) documents the percentage of detecting a number of regimes per model parameter over the 1000 simulated series per DGP for the SELO method. Overall, the detection rates of identifying the true number of regimes per parameter are excellent and besides DGP F, they are at least equal to 86.4%. Interestingly, this detection rate does not deteriorate when the innovation is driven by a GARCH process. The worst detection rates arise for the DGP F. Even though this DGP is highly persistent with an autocorrelation structure that barely varies over time, the SELO method correctly identifies that the intercept does not experience abrupt switches 69.7% of the times. Note that the potential breakpoint sets for this DGP poorly identify the true breakpoints since only 25.5% of the sets exhibit at least one potential CP close to every true breakpoints. Therefore, the SELO detection rate could hardly exceed this bound. As exemplified by DGPs G, H and I, the detection rates of the SELO method remain excellent when exogenous variables kick in even in the presence of heteroscedasticity. The Table also documents the rate of detecting the true model (i.e. jointly the correct number of regimes) with a posterior probability of at least 10%.\footnote{For this simulation study, the number of explanatory variables ranges from 2 to 5 and the maximum number of potential regimes observed for each DGP is as follows: DGP A (5), DGP B (5), DGP C (6), DGP D (6), DGP E (6), DGP F (6), DGP G (5), DGP H (9), DGP I (9). It can thus lead to a number of models amounting to $2^{40}$. In such a case, the set of the models exhibiting a probability equal or greater than 10% has a prior probability of containing the true model that is approximately equal to $\frac{1000}{2^{40}}\%$.} For all the DGPs but DGP F, the correct detection amounts to at least 83.1% and 85.1% for the constant and the GARCH innovation dynamics, respectively. These excellent results highlight that model uncertainty should be taken into account since several models often exhibit high posterior probabilities.
DGPs from Table (ref) are frequently used in the CP literature to assess the performance of a new segmentation method (see, e.g., cho2015multiple, yau2016inference and korkas2017multiple). Nevertheless, our empirical exercise implies more explanatory variables and a smaller sample size. To assess the SELO performance in such environment, we also consider fourteen variants of an 'empirical DGP' given by
where $T=256$ as in the application, $\epsilon_t \sim \mathcal{N}(0,1)$ and $\mathbf x_t=(1,x_{t,1},\ldots,x_{t,12})'$. The explanatory variables are close to the risk factors used in our empirical exercise. In particular, they are generated from AR models whose coefficients and AR orders are estimated using the risk factors of the application. The parameter values of $\boldsymbol{\beta}_1$ are equal to the OLS estimates of the Hedge fund Index (HFI) regression without breakpoints (see Table (ref) in the empirical application). We consider 14 variants of the DGP given by Equation (ref) that differ by the number of parameters experiencing a breakpoint. Defining $\boldsymbol{\beta}_{2}=(\beta_{2,1},\ldots,\beta_{2,13})'$ and considering the $ith$ DGP, with $i=1,\ldots,14$, we have $\beta_{2,j} = \beta_{1,j} + 3\omega_j \text{sign}(\beta_{1,j})$ for $j< i$ and $\beta_{2,j} = \beta_{1,j}$ for $j\geq i$. The size of the break given by $\omega_j$ is equal to the standard deviation of the $jth$ OLS estimate of the Hedge fund Index (HFI) regression without breakpoints (see Table (ref)). To summarize, the first DGP does not exhibit a breakpoint while the 14th one exhibits a structural change in all its parameters. As before, we consider homoskedastic errors with $\sigma_t^2\equiv \bar{\omega}^2 = 1.7$ $\forall t$ and heteroskedastic ones with $\sigma_t^2= 0.05\bar{\omega}^2 + 0.05 \epsilon_{t-1}^2 + 0.9\sigma_{t-1}^2$ for $t>1$ in which $\bar{\omega}^2$ stands for the OLS variance estimate of the Hedge fund Index (HFI) regression without breakpoints (see Table 7 of the paper). Figure (ref) displays one simulated series drawn from some of the fourteen variants with heteroskedastic errors.
To carry out the Monte Carlo study, we have drawn 100 series from the 14 variants of the empirical DGP. For each simulated series, we have also generated the explanatory variables using the AR models. Table (ref) provides the percentage of the number of regimes detected by the selective segmentation and by the Lasso methods. The Lasso approach consists in using a Lasso penalty function instead of the SELO penalty function and in choosing the best Lasso penalty value and the corresponding model using the marginal likelihood proposed in Section (ref).\footnote{We replace the SELO with the Lasso penalty function. In particular, conditional on the potential breakpoints $\boldsymbol{\tau}$, we minimize the following objective function:
Using the Matlab Lasso toolbox of mcilhagga2016penalized. We also tested the matlab glmnet toolbox of LassoMatlab which leads to similar results. In particular, we also observe that the Lasso estimates over-estimate the number of regimes.} First, we observe that the selective segmentation approach delivers high detection rates whatever the dynamic of the variance. In addition, the detection rate does not deteriorate when the parameter experiences a break. Note also that all the average detections are above 89%. In contrast, the Lasso method does not provide good results when only a subset of the parameters exhibits a CP. In fact, the average detection rates follow a 'U-shape' function meaning that the Lasso method is good at detecting no breakpoint or when all the parameters are experiencing a break. When partial breaks occur, the Lasso approach typically over-estimates the number of regimes.
To further illustrate the issue with the Lasso method, Table (ref) shows the detailed results based on 100 simulated series when the first five parameters exhibit a CP (equivalent to the variant called '\# CPs = 5' in Table (ref)). While the selective segmentation method accurately detects the number of regimes for each parameter, the Lasso approach finds two regimes for most of the parameters that are constant over the sample.
We end this simulation section with a "big data" example motivated by the fact that when the number of explanatory variables is large, the current Bayesian alternatives do not work GiordaniKohn2008,Eo2012,huber2019should,dufays2019relevant (see SA (ref) for more details). To do so, we propose the DGP J that is specified by 100 explanatory variables and one CP as follows:
DGP J: piecewise linear model with big data\\ $$ Y_{t} =
$$
where $T=1024$, $\forall$ $t\in [1,T]$ and for $i=1,...,100$, $\mathbf{x}_{t,i} \sim \mathcal{N}{\left( 0,1\right)}$ and $\varepsilon_t \sim \mathcal{N}(0,1)$. The parameter values of $\boldsymbol{\beta_1}$ are uniformly and randomly set to $-1$ or $1$. In the second regime, the parameter values of $\boldsymbol{\beta_2}$ are equal to $\boldsymbol{\beta_1}$ except for 10 of them randomly chosen that are set to the opposite value (i.e. $-\boldsymbol{\beta}_1$). Thus, 10 parameters of DGP J does experience a break at observation 500. \\ We simulate 100 series from DGP J to assess the SELO performance in detecting which parameters experience a breakpoint. For every simulation, the selective segmentation approach identifies 10 parameters that experience one breakpoint in the sample while the others remain constant. In addition, the exact model specification was always among the specification exhibiting a posterior probability of at least 10%.
We illustrate the selective segmentation method with 14 monthly Credit Suisse HF indices spanning from March 1994 to March 2016. These indices are the weighted average of HF returns following specific trading strategies. fung2004hedge suggested a risk-based approach to model HF returns and identified seven factors on which HF strategies are generally exposed fung2001risk. Since this seminal work, many other risk factors have been uncovered. So, we include five other risk factors that are also popular in the literature. We add two Fung and Hsieh trend following risk factors, PTFSIR, returns on PTFS short term interest rate lookback straddle, and PTFSSTK, returns on PTFS stock index lookback straddle. Following Agarwal2004Risks, among many others, we also use the Up-minus-Down (UMD) factor Carhart-JF-1997. As suggested by chen1986economic, we include the expected inflation, the log relative of US Consumer Price Index (CPI). Finally, we also take into account a factor for real estate risk relevant to explain HF and stocks returns Ambrose-DLima-JREFE-2016,Carmichael-Coen-REE-2018. Table (ref) documents the fourteen strategies on which we focus as well as the twelve factors. \\
It is well acknowledged in the financial literature that HF strategies (or trading techniques) are time-varying. Their changing risk exposures are directly related to market events and economic fluctuations (see, e.g., Agarwal2004Risks, Fung2008hedgefunds or Patton2015change among others). Hedge fund time-varying risk dynamics has important implications for performance appraisal. As pointed out by Mitchell2001characteristics, the changes can be in response to arbitrage opportunities. The cycles of mergers and acquisitions in the 1990s and the 2000s and the corresponding level of risk arbitrage led by HF are illustrations of these changing dynamics. In standard linear asset pricing models, the intercept and risk factor loadings are not constant but time-varying. Moreover, HF returns exhibit significant non-linearities. Therefore, there is a need of dynamic models able to capture non-linearities and changes in risk exposures.\\ Following meligkotsidou2008detecting, we suggest the use of CP risk factor models. This class of models is suited for studying the changes in risk exposures and their time-varying parameters. However, instead of directly focusing on the twelve factors, we take a slightly different approach since we additionally take into account autocorrelations of the returns.\footnote{As reported by getmansky2004econometric, the analysis of serial dependence of returns is a reasonable way of assessing the liquidity of hedge fund investments.} To do so, we first look at the best autoregressive model that fits the returns. In particular, for each HF returns, we estimate ARX(q) models with $q$ ranging from 0 to 4 and in which the explanatory variables are the twelve factors (and an intercept) and we select the best AR order using the Bayesian information criterion (BIC). Table (ref) documents the best order for each strategy.
As reported by fung2004hedge, composites obtained from the individual funds may be contaminated with severe survivorship, selection and instant history biases. Therefore, to avoid these problems, we use the Credit Suisse indices that provide full transparency about their constituents.
Section (ref) discusses in-sample results of our selective segmentation method and we compare them to those of standard CP models and time-varying parameter models. We then illustrate the difference of our approach with the CP method of meligkotsidou2008detecting in Section (ref). Section (ref) documents a forecasting exercise in which we assess the predictive performance of the selective segmentation approach with respect to flexible alternatives. Importantly, all the subsequent results include the optimal AR order documented in Table (ref) as additional explanatory variables.
fung2004hedge focus on linear models. However, as the period covers critical events such as the Long Term Capital Management (LTCM) collapse, the dot-com crisis and the global financial crisis (GFC), one could argue that CP models are more appropriate. In this Section, we focus on two specific indices, namely the Hedge Fund Index (HFI) and the HF returns that are applying a Fixed-Income Arbitrage (FIA) strategy. Results for all the other returns are available upon request. \\ Tables (ref) and (ref) show how the selective segmentation method can improve the interpretation of CP models. The Tables document how the results evolve from a standard linear risk model to a selective segmentation model passing by a standard CP process. As expected, for the two HF returns, ignoring breakpoints can be misleading as the CP results emphasize that they modify the risk exposition of the returns. Also, although one can study in details the results of the standard CP model, the selective segmentation model offers a straightforward picture of the relevant risk factors and how the risk exposition evolves. It also estimates more accurately the parameters that do not change when a break occurs. As the CP model detects three breakpoints for the HFI and six abrupt changes for the FIA strategy, the number of models to consider amounts to $2^{36}$ and $2^{84}$ respectively. Our selective segmentation strategy explores these large model spaces and find the most promising configurations in several minutes on a standard laptop. Let us now discuss in more details the results of the two returns.\\
Hedge Fund Index\\ As documented in Table (ref), the CP model with breakpoints determined by the approach in Section (ref) finds four regimes (hereafter CP-YZ). The relevant breakpoints occur in April 2000, in December 2001 and in August 2014. Interestingly, these dates coincide with the dot-com crash that spanned from March 2000 to October 2002 and when stocks suffer steepest drop in 2014 (introducing the first significant stock market scares after the GFC). It is well acknowledged in the financial literature that the end of the dot-com bubble had important consequences for financial markets in the early 2000s. While all the parameters change for the CP model, the selective segmentation mainly identifies that the factors related to the breaks are the market factor (PMKT), the default risk factor (DEF), the momentum risk factor (UMD) and the real estate risk factor (NAREIT). Moreover, it discards two spurious breaks occurring in December 2001 and in August 2014 making the model even more parsimonious. We can notice that the market factor decreases from 0.34 during the first period to 0.21 during the second period. HFI is indeed more conservative during the 2000s and less correlated with the financial markets. We observe the same trend for the default risk factor increasing from -9.15 to -2.42. The momentum factor, UMD, sharply declines after the nineties known as very volatile as reported by Shiller-PUP-2015 and documented by Campbell-JF-2000 who both highlight the bullish market of this decade. The momentum is still significant since 2000 but its impact has significantly weakened. It decreases from 0.22 to 0.05. The real estate risk factor also exhibits a breakdown in the early 2000s. It is indeed well acknowledged in the real estate economics literature that the 1990s are considered as the new era of real estate investment trusts (REITs) and the period beginning in the early 2000s as the maturity REITs era (Pagliari-Scherer-Monopoli-REE-2005, Ambrose-Lee-Peek-REE-2007 and Carmichael-Coen-REE-2018 among others). From the new REITs era to the maturity REITs era, the real estate risk factor declines from 0.15 to -0.01. These results are consistent with the important variations of interest rates during these two sub-periods and the important increase of credit risk in the 2000s. As a final note, the credible interval of the breakpoint is narrow which indicates a sharp change in the risk exposition in March 2000 (see also Figure (ref) below). \\
Figure (ref) shows the posterior medians over time and their corresponding credible intervals of the parameters related to the MKT, DEF and UMD factors given by our method (see Section (ref) for the related Bayesian model and how the breakpoints are integrated out) and the time-varying parameter (TVP) model (see SA (ref) for the model specification). As with the CP model, one can easier interpret the time-varying dynamics of the parameters given by the selective segmentation method than those of the TVP model. For instance, while the exposition to the default factor seems fixed over the sample due to the smooth transition of the parameter, it is clear that the exposition is changing when we look at the selective segmentation results. Regarding the market factor, we also observe with the TVP model that the exposition seems different before and after the dot-com crash but the credible intervals are too wide to confirm the statement.
Fixed Income Arbitrage (FIA) strategy\\ The FIA strategy is based on the exploitation of inefficiencies in the pricing of bonds and interest rate derivatives (including futures, options, swaps and also mortgage back securities). It was very appreciated among hedge fund managers until the collapse of the LTCM fund in September 1998. After this incident, a change of behavior among managers has been observed for this strategy on financial markets. \\ The results of the FIA returns from a standard linear regression to the selective segmentation model are documented in Table (ref). Focusing on the latter model, breakpoints are detected in September 1997, in March 1999, August 2000 which are related to the Russia financial crisis and the dot-com crash. In addition, three other breakpoints capture the financial crisis: the beginning of the global financial crisis in late June 2007, the turmoil of September and October 2008 after the collapse of Lehman Brothers and the sovereign debt crisis in the euro zone. The selective segmentation specification highlights the role played by the market factor (i.e., PMKT, before and after the global financial crisis with estimates of 0.02 and 0.08 respectively), the variation of the size effect, the default risk factor, the momentum factor, the inflation factor and four trend following risk factors, especially during the sub-prime crisis from late June 2007 to October 2008. The size effect changes from 0.02 before the crisis to -0.05 afterwards. The bond trend following factor, PTFSBD, the currency trend following factor, PTFSFX, the commodity trend following factor, PTFSCOM and the stock index trend following factor, PTFSSTK are highly significant during the GFC. The credit spread factor, DEF, is highly significant and constant after the dot com crisis (estimate amounting to -1.18). Significant during the first and the second periods related to the bullish market of the 1990s, the momentum effect, UMD, dramatically changes during the GFC, with a negative estimate of -0.34 (as expected). Our results also report the impact of inflation risk, CPI, during the crisis (with an estimate of 2.69) until the sovereign debt crisis (fifth and sixth periods). The real estate risk factor, NAREIT, is significant and positive during the 1990s (first, second and third periods) and during the GFC, with a negative estimate (as expected) of -0.16. As shown by Figures (ref) to (ref) in SA (ref), the selective segmentation method allows easier detection of the relevant factors as compared to the TVP model.
We now compare our results with those of meligkotsidou2008detecting. meligkotsidou2008detecting rely on CP models to capture the risk exposition of HF returns over time. In particular, they consider the $4096$ distinct combinations of the twelve risk factors and for each of them, they estimate a CP model exhibiting several numbers of segments $m$ (from one to ten). Eventually, they use the marginal likelihood to select the best model among the set of $m \times 2^K$ estimated processes (i.e., 40960 models since $m=10$ and $K=12$). Their approach consists therefore in first selecting the relevant factors and then, in investigating if the exposition to them is time-varying. \\ There is a striking difference with our approach since, for each breakpoint, our procedure can detect what are the time-varying factors. In fact, our approach discriminates between $2^{m \times K}$ models; a number of models that exponentially increases with the amount of breaks. Note that we could also search for the best regressors to include by considering all the 4096 distinct combinations of the twelve factors. In such a case, the number of models to consider would reach $2^{(m+1) \times K}$.\\ We reproduce the results of meligkotsidou2008detecting on our data by additionally taking the autocorrelation structure into account. Fixing the AR order $q$ to the value given in Table (ref), for each possible combination of the factors, we estimate CP-ARX(q) models with different numbers of breaks (ranging from 1 to 10) by (globally) minimizing the MDL criterion. Then, we report the combination of factors exhibiting the best MDL value. Hereafter, we denote this model by CP-MV.\footnote{Our approach is slightly different as the one used in meligkotsidou2008detecting since we minimize the MDL criterion instead of maximizing the marginal likelihood of a Bayesian CP model for finding the best combination of the factors and the breakpoints. This is motivated by the fact that the MDL criterion consistently selects the true number of regimes while there is no equivalent proof for the marginal likelihood used in meligkotsidou2008detecting. In addition, Dufays2019 show that the MDL criterion is equal to minus the marginal log-likelihood of a CP Bayesian model with particular g-prior distributions. So, our approach can be understood as the method of meligkotsidou2008detecting with different hyper-parameters. We globally minimize the MDL criterion using the dynamic programming of bai2003computation.} \\ Table (ref) documents the factors of the CP-MV model for the two HF strategies. It also reports the factors which exhibit significant parameter estimates at least at one period for the TVP model and the selective segmentation approach. For the HFI, the selected factors by the selective segmentation methods seem more complete than those of the TVP model in light of the current literature. While the TVP selects eight risk factors, the selective segmentation reports the same factors, except the SMB factor, and adds two trend following risk factors, PTFSBD and PTFSSTK, and the real estate risk factor, NAREIT. For the FIA, the results are much more contrasted. All risk factors are selected by the selective segmentation methods whereas five factors are omitted by the TVP (SMB, TERM, PTFSFX, PTFSCOM and UMD). The CP-MV approach does not select SMB, NAREIT and four trend following risk factors (PTFSBD, PTFSCOM, PTFSSIR and PTFSSTK) for the HFI. The analysis of the single strategy, Fixed Income Arbitrage (FIA), also highlights important and significant differences as far as only two factors are selected, DEF and PTFSFX, by the CP-MV process. It may be very surprising since the market premium is (almost) always used in linear asset pricing models as pointed out by fung2001risk. We may also note that the look-back straddles on commodities, PTFSCOM, on bond, PTFSBD and on stock index and PTFSSTK designed to capture non-linearities especially during changes in international economic policies, are not selected whereas the phenomenon is observed just after the GFC.
Table (ref) does not inform on the dynamic of the selected factors by the CP-MV process. Although the preferred specification of the CP-MV model does not include all the factors, the risk exposure of the HF strategies is still abruptly changing over time. Regarding the HFI, two breakpoints are detected and occur in April 2000 and in March 2005, respectively. Table (ref) shows how the selective segmentation method improves the interpretation of the CP-MV results. It also illustrates the improvement in economic modelling as far as it highlights the relative role played by static and dynamic parameters. First, we observe that the alpha, the currency lookback straddle, PTFSFX, and the consumer prince index, CPI, are quite static during the full period. Interestingly, the selective segmentation also reports the CP in the bullish market in the early 2000s. As mentioned earlier, the financial markets were indeed very volatile during the 1990s. This break is clearly reported by the dynamic risk factors, PMKT, UMD and DEF. PMKT declines from 0.41 to 0.21 during the 2000s and afterwards. The trend is more striking for the momentum factor, UMD, with a decline from 0.18 to 0.05, and for the credit default risk factor, DEF, rising from -11.71 to -2.21. The term structure risk factor, TERM, is not statistically significant after the rise of the housing price index in 2005-2006, announcing the collapse of the financial markets in 2008 (as anticipated by Shiller-PUP-2015 among others) and the following quantitative easing policies with very low inters rates.
The FIA strategy exhibits seven regimes which makes the CP-MV model heavily parametrized (i.e., $K \times m = 35 $ parameters). This large number of regimes is probably related to the fact that more breakpoints are needed to adequately fit the FIA returns since the CP-MV specification includes only two risk factors, DEF and PTFSFX and/or because in this specific economic modelling the breakpoints also capture the variance dynamic. Using the selected factors and the breakpoints of the best CP-MV specification, we estimate the selective segmentation model to uncover what are the static and the dynamic parameters. First, we must acknowledge that the best specification selected by the CP-MV model is doubtful with regards to the financial literature and the practice. Nevertheless, we compare the CP-MV approach with the selected segmentation to highlight the contribution of the latter. In particular, we observe that the 'alpha' is varying and statistically positive during the 1990s until LTCM collapse. As expected, the default risk factor is negative, time-varying and high during crises (-9.49 during the LTCM collapse and -6.17 during the GFC). After the GFC, the default factor is constant, negative and not statistically significant. This result is consistent with the trend observed on financial markets (especially on fixed incomes markets after the GFC). The currency trend following factor, PTFSFX, is very low, time-varying and statistically significant during the global financial crisis (as expected) and before the impact of the quantitative easing policies starting in the late 2010. After this date (11/2010) the factor is not statistically significant. This is an illustration of the impact of quantitative easing on fixed income arbitrage.
As illustrated by these empirical results, our method uncovers which parameters truly vary when a CP is detected. This technical improvement induces financial consequences and especially cost reductions.\footnote{We thank an anonymous referee for this relevant comment.} Moreover, our method should imply more accurate and thus less frequent portfolio rebalancing strategies. Investor could indeed change his timing by using our approach and decide to rebalance parsimoniously (and thus efficiently) his investment when a break is detected, with a special focus on the relevant benchmark.
Sections (ref) and (ref) highlight the in-sample advantages of detecting which parameter truly varies when a break is detected. In addition to that, since the selective segmentation method can more accurately estimate parameters that do not change when a break occurs, we could also expect some prediction gains with respect to the standard CP model. In this Section, we investigate this aspect using the root mean squared forecast errors (RMSFE) and the cumulative log predictive density (CLPD), two standard loss functions specified as,
in which $\hat{y}_t$ is the conditional mean of $y_t$ given the information up to period $t$ (i.e., $\mathbb{E}(y_t|y_{1:t-1},\mathbf x_t)$), $f(y_{t}|y_{1:t-1},\mathbf x_t)$ denotes the predictive density of the model and $\underline{t}+1$ denotes the beginning of the out-of-sample forecasting period. In our prediction exercise, the training set is fixed to 20% of the sample size and the 80% remaining observations are used to assess the forecast performance (i.e., $\underline{t} = 0.2T$). Since our data comprise 265 monthly returns, the out-of-sample set of observations amounts to 212 months. Each time we move forward by one month, all the considered models are re-estimated and a forecast for the next period is produced.\\ As competitors to our model, we consider three other processes: i) a linear regression, ii) a standard CP model with breakpoints determined by the modified method of yau2016inference documented in Section (ref) (hereafter CP-YZ), iii) a CP model with the number and the locations of the breakpoints selected by minimizing the MDL criterion (hereafter CP-MDL).\footnote{We do not compare with the CP model of meligkotsidou2008detecting since the model is computationally too involved due to the number of explanatory variables. When an AR(2) model is selected, the number of models to consider at each iteration of the prediction exercise amounts to $10 \times 2^{15} = 327680$.} The minimization of the MDL criterion is carried out using the dynamic programming of bai2003computation.\footnote{See eckley2011analysis for a discussion on the implementation of the algorithm for the MDL criterion. Minimum regime duration is set to $\frac{3}{2}(K+1)$ to avoid capturing outliers. This choice is in favor of the standard CP model as the parameter estimates of the new regimes are based on at least $\frac{3}{2}(K+1)$ observations.} In addition to the factors and an intercept, we also account for the autocorrelation of the HF returns by fixing the AR order to the value given in Table (ref). Regarding the CLPD metric, we assume a normal distribution for the error term and we also use the prior distributions given in Equation (ref) for the linear and the full CP models.
Table (ref) documents the RMSFE and the CLPD criteria for all the Credit Suisse HF returns. For both metrics, we observe that the linear model dominates at least half of the times. Overall, the selective segmentation method improves the RMSFE and the CLPD for 6 and 5 out of 14 HF returns, respectively. Importantly, Table (ref) highlights that the selective segmentation process provides the most robust predictions. In particular, our approach delivers at least the second best predictive performance for all the HF returns. This is evidence that model averaging stabilizes the forecast by reducing its variance as argued in RapachEtAl2009. Interestingly, our method compares extremely well with respect to the two CP models since it outperforms them 13 out of 14 HF returns for both metrics. Since the CP models are based on the same breakpoints as the selective segmentation processes, it is remarkable that the latter models almost systematically dominate CP models where all the parameters are time-varying. From this small sample of series, we could argue that the selective segmentation approach should replace the CP process as it would likely improve the forecast performance.
Since the seminal work of chernoff1964estimating, many CP detection methods for linear models have been proposed. Most of these CP models have in common to assume, at least in practice, that all the model parameters have to change when a break is detected. In this paper, we propose to go beyond this standard framework by capturing which parameters vary when a structural break occurs. Even when conditioning to the break dates, detecting the parameters that vary from one segment to the next is not straightforward since the number of possibilities grows exponentially with the number of breaks and the number of explanatory variables. To solve this dimensional problem, we propose a penalized regression method to explore the model space and we select the best specification by maximizing a criterion that can be interpreted as a marginal likelihood in the Bayesian paradigm. \\ To carry out the model space exploration, we use an almost unbiased penalty function, a desirable property in CP frameworks that is not exhibited by standard penalty functions (e.g., LASSO and Ridge estimators). Also, we prove the consistency of our estimator and we show how to estimate it using the DAEM algorithm. To apply the DAEM algorithm in our context, we transform the penalty function into a mixture of Normal distributions. This simple transformation greatly speeds the estimation as the DAEM algorithm iterates over closed-form expressions.\\ Once the promising models have been uncovered by the penalized regression approach, selecting the parameters of the penalty function is carried out by maximizing a marginal likelihood. Thanks to the Bayesian interpretation of this consistent criterion, we can take model uncertainty into account and can do Bayesian model averaging, a feature that generally improves forecast performance. A simulation study highlights that our selective segmentation method works well in practice for a range of diversified data generating processes.\\ We illustrate our approach with HF returns. The selective segmentation model has two main advantages. First, as the standard CP models, it detects the breakpoints and the corresponding regimes. Second, it highlights the time-varying dynamics of the changing risk factors. When we compare our model with previous advanced CP models, we observe that it is particularly appealing to capture the time-varying dynamics of risk exposures. Then, we test the predictive performance of the selective segmentation approach with respect to the linear regression and standard CP processes. We note that our method produces the most robust forecasts and almost systematically dominates the CP processes based on the same breakpoints. These encouraging results suggest promising developments and applications in financial economics. \\ Importantly, an R package for estimating the model is available on the corresponding author's web page. This package stands for a building block of our future research that will include a dynamic variance and multivariate models.
{\setstretch{0.9} }