1\documentclass[10pt,a4paper]{article} 2\usepackage{a4wide} 3\usepackage{graphicx,color,alltt} 4\usepackage[authoryear,round,longnamesfirst]{natbib} 5\usepackage{hyperref} 6 7\definecolor{Red}{rgb}{0.7,0,0} 8\definecolor{Blue}{rgb}{0,0,0.8} 9\definecolor{hellgrau}{rgb}{0.55,0.55,0.55} 10 11\newcommand{\E}{\mbox{$\mathsf{E}$}} 12\newcommand{\VAR}{\mbox{$\mathsf{VAR}$}} 13\newcommand{\COV}{\mbox{$\mathsf{COV}$}} 14\newcommand{\p}{\mbox{$\mathsf{P}$}} 15\newcommand{\email}[1]{\href{mailto:#1}{\normalfont\texttt{#1}}} 16\newenvironment{smallexample}{\begin{alltt}\small}{\end{alltt}} 17 18\setlength{\parskip}{0.5ex plus0.1ex minus0.1ex} 19\setlength{\parindent}{0em} 20 21\bibpunct{(}{)}{;}{a}{}{,} 22 23\newcommand{\ui}{\underline{i}} 24\newcommand{\oi}{\overline{\imath}} 25 26\RequirePackage{color} 27\definecolor{Red}{rgb}{0.5,0,0} 28\definecolor{Blue}{rgb}{0,0,0.5} 29\definecolor{hellgrau}{rgb}{0.55,0.55,0.55} 30 31\hypersetup{% 32hyperindex,% 33colorlinks,% 34linktocpage,% 35plainpages=false,% 36linkcolor=Blue,% 37citecolor=Blue,% 38urlcolor=Red,% 39pdfstartview=Fit,% 40pdfview={XYZ null null null}% 41} 42 43 44\begin{document} 45 46\SweaveOpts{engine=R,eps=FALSE} 47%\VignetteIndexEntry{strucchange: An R Package for Testing for Structural Change in Linear Regression Models} 48%\VignetteDepends{strucchange} 49%\VignetteKeywords{structural change, CUSUM, MOSUM, recursive estimates, moving estimates, monitoring, R, S} 50%\VignettePackage{strucchange} 51 52\title{\texttt{strucchange}: An \textsf{R} Package for Testing for Structural 53Change in Linear Regression Models} 54\author{Achim Zeileis \hspace{0.4cm} Friedrich Leisch \hspace{0.4cm} Kurt Hornik \hspace{0.4cm} Christian Kleiber} 55\date{} 56\maketitle 57 58\begin{abstract} 59This introduction to the \textsf{R} package \texttt{strucchange} is a (slightly) 60modified version of 61\cite{Z-papers:Zeileis+Leisch+Hornik:2002}, which 62reviews tests for structural change in linear regression models from 63the generalized fluctuation test framework as well as from the $F$ test (Chow 64test) framework. Since \cite{Z-papers:Zeileis+Leisch+Hornik:2002} various extensions 65were added to the package, in particular related to breakpoint estimation 66\citep[also know as ``dating'', discussed in][]{Z-papers:Zeileis+Kleiber+Kraemer:2003} 67and to structural change tests in other parametric models \citep{Z-papers:Zeileis:2006}. 68A more unifying view of the underlying theory is presented in \cite{Z-papers:Zeileis:2005} 69and \cite{Z-papers:Zeileis+Shah+Patnaik:2010}. 70 71Here, we focus on the linear regression model and 72introduce a unified approach for implementing tests from the fluctuation test 73and $F$~test framework for this model, illustrating how this approach 74has been realized in {\tt strucchange}. 75Enhancing the standard significance test approach the package contains 76methods to fit, plot and test empirical fluctuation processes (like CUSUM, MOSUM 77and estimates-based processes) and to compute, plot and test 78sequences of $F$ statistics with the sup$F$, ave$F$ and exp$F$ test. 79Thus, it makes powerful tools available to display information about 80structural changes in regression relationships and to assess their significance. 81Furthermore, it is described how incoming data can be monitored. 82\end{abstract} 83 84\noindent 85{\bf Keywords:} structural change, CUSUM, MOSUM, recursive estimates, moving 86estimates, monitoring, \textsf{R}, \textsf{S}.\\ 87 88\section{Introduction} \label{sec:intro} 89 90The problem of detecting structural changes in linear regression relationships 91has been an important topic in statistical and econometric research. 92The most important classes of tests on structural change are the tests from the 93generalized fluctuation test framework \citep{Z:Kuan+Hornik:1995} on the 94one hand and tests based on $F$ statistics 95\citep{Z:Hansen:1992,Z:Andrews:1993,Z:Andrews+Ploberger:1994} on the other. 96The first class includes in particular the CUSUM and MOSUM tests and the 97fluctuation test, while the Chow and the sup$F$ test belong to the latter. 98A topic that gained more interest rather recently is to monitor structural 99change, i.e., to start after a history phase (without structural changes) to 100analyze new observations and to be able to detect a structural change as soon 101after its occurrence as possible.\\ 102 103This paper concerns ideas and methods for implementing generalized fluctuation 104tests as well as $F$ tests in a comprehensive and flexible way, that reflects 105the common features of the testing procedures. It also offers facilities to 106display the results in various ways.\\ 107 108This paper is organized as follows: In Section~\ref{sec:model} the standard 109linear regression model upon which all tests are based will be described and the 110testing problem will be specified. Section~\ref{sec:data} introduces a data set 111which is also available in the package and which is used for the 112examples in this paper. The following sections \ref{sec:fluctests}, 113\ref{sec:Ftests} and \ref{sec:monitor} will then explain the tests, how they are 114implemented in {\tt strucchange} and give examples for each. Section 115\ref{sec:fluctests} is concerned with computing empirical fluctuation processes, 116with plotting them and the corresponding boundaries and finally with testing for 117structural change based on these processes. Analogously, 118Section~\ref{sec:Ftests} introduces the $F$ statistics and their plotting and 119testing methods before Section~\ref{sec:monitor} extends the tools from 120Section~\ref{sec:fluctests} for the monitoring case. 121 122 123\section{The model} \label{sec:model} 124Consider the standard linear regression model 125\begin{equation} \label{model1} y_i = x_i^\top \beta_i + u_i 126\qquad (i = 1, \dots, n), \end{equation} 127where at time $i$, $y_i$ is the observation of the dependent 128variable, $x_i = (1, x_{i2}, \dots, x_{ik})^\top$ is a $k \times 1$ vector of 129observations of the independent variables, with the first component equal to 130unity, $u_i$ are iid(0, $\sigma^2$), and $\beta_i$ is the $k \times 1$ vector of 131regression coefficients. Tests on structural change are concerned with testing 132the null hypothesis of ``no structural change'' 133\begin{equation} \label{null-hypothesis} 134H_0: \quad \beta_i = \beta_0 \qquad (i = 1, \dots, n) 135\end{equation} 136against 137the alternative that the coefficient vector varies over time, with 138certain tests being more or less suitable (i.e., having good or poor power) for 139certain patterns of deviation from the null hypothesis.\\ 140 141It is assumed 142that the regressors are nonstochastic with $||x_i|| = O(1)$ and that 143\begin{equation} \label{A5} \frac{1}{n} \sum_{i=1}^n x_i x_i^\top \quad 144\longrightarrow \quad Q\end{equation} 145for some finite regular matrix $Q$. These 146are strict regularity conditions excluding trends in the data which are assumed 147for simplicity. For some tests these assumptions can be 148extended to dynamic models without changing the main properties of the tests; 149but as these details are not part of the focus of this work they are omitted 150here.\\ 151 152In what follows $\hat \beta^{(i, j)}$ is the ordinary least 153squares (OLS) estimate of the regression coefficients based on the 154observations $i+1, \dots, i+j$, and $\hat \beta^{(i)} = \hat \beta^{(0,i)}$ is 155the OLS estimate based on all observations up to $i$. Hence $\hat \beta^{(n)}$ 156is the common OLS estimate in the linear regression model. 157Similarly $X^{(i)}$ is the regressor matrix based on all observations up to $i$. 158The OLS residuals 159are denoted as $\hat u_i = y_i - x_i^\top \hat \beta^{(n)}$ with the 160variance estimate $\hat{\sigma}^2 = \frac{1}{n-k} 161\sum_{i=1}^n \hat u_i^2$. Another type of residuals that are often used in 162tests on structural change are the recursive residuals \begin{equation} 163\label{rr} \tilde u_i \; = \; \frac{y_i - x_i^\top \hat 164\beta^{(i-1)}}{\sqrt{1+x_i^\top \left(X^{(i-1)^\top} X^{(i-1)} \right)^{-1}x_i}} 165\qquad (i = k+1, \dots, n),\end{equation} which have zero mean and variance 166$\sigma^2$ under the null hypothesis. The corresponding variance estimate is 167$\tilde{\sigma}^2 = \frac{1}{n-k} \sum_{i=k+1}^n (\tilde u_i - \bar{ \tilde 168u})^2$. 169 170\section{The data} \label{sec:data} 171The data used for examples throughout this paper are macroeconomic time 172series from the USA. The data set contains the aggregated monthly personal 173income and 174personal consumption expenditures (in billion US dollars) between January 1959 175and February 2001, which are seasonally adjusted at annual rates. 176It was originally taken from \url{http://www.economagic.com/}, a web site for 177economic times series. Both time series are depicted in Figure 178\ref{fig:USIncExp}. 179 180\setkeys{Gin}{width=0.6\textwidth} 181\begin{figure}[htbp] 182\begin{center} 183<<data,fig=TRUE,height=4,width=6,echo=FALSE,results=hide>>= 184library("strucchange") 185data("USIncExp") 186plot(USIncExp, plot.type = "single", col = 1:2, ylab = "billion US$") 187legend(1960, max(USIncExp), c("income", "expenditures"), 188 lty = c(1,1), col = 1:2, bty = "n") 189@ 190\caption{\label{fig:USIncExp} Personal income and personal consumption 191expenditures in the US} 192\end{center} 193\end{figure} 194 195The data is 196available in the {\tt strucchange} package: it can be loaded and a 197suitable subset chosen by 198<<subset>>= 199library("strucchange") 200data("USIncExp") 201USIncExp2 <- window(USIncExp, start = c(1985,12)) 202@ 203 204We use a simple error correction model (ECM) for the consumption function 205similar to \cite{sc:Hansen:1992a}: 206\begin{eqnarray} \label{ecm.model} 207\Delta c_t & = & \beta_1 + \beta_2 \; e_{t-1} + \beta_3 \; \Delta i_t + u_t, \\ 208\label{coint.model} e_t & = & c_t - \alpha_1 - \alpha_2 \; i_t, 209\end{eqnarray} 210where $c_t$ is the consumption expenditure and $i_t$ the income. We estimate 211the cointegration equation (\ref{coint.model}) by OLS and use the residuals 212$\hat e_t$ as regressors in equation (\ref{ecm.model}), in which we will test 213for structural change. Thus, the dependent variable is the increase in 214expenditure and the regressors are the cointegration residuals and the 215increments of income (and a constant). 216To compute the cointegration residuals 217and set up the model equation we need the following steps in \textsf{R}: 218<<ecm-setup>>= 219coint.res <- residuals(lm(expenditure ~ income, data = USIncExp2)) 220coint.res <- lag(ts(coint.res, start = c(1985,12), freq = 12), k = -1) 221USIncExp2 <- cbind(USIncExp2, diff(USIncExp2), coint.res) 222USIncExp2 <- window(USIncExp2, start = c(1986,1), end = c(2001,2)) 223colnames(USIncExp2) <- c("income", "expenditure", "diff.income", 224 "diff.expenditure", "coint.res") 225ecm.model <- diff.expenditure ~ coint.res + diff.income 226@ 227 228Figure~\ref{fig:ts} shows the transformed time series necessary for estimation 229of equation (\ref{ecm.model}). 230 231\begin{figure}[htbp] 232\begin{center} 233<<ts-used,fig=TRUE,echo=FALSE>>= 234plot(USIncExp2[,3:5], main = "") 235@ 236\caption{\label{fig:ts} Time series used -- first differences and 237cointegration residuals} \end{center} 238\end{figure} 239 240 241In the following sections we will apply the methods introduced to test for 242structural change in this model. 243 244 245\section{Generalized fluctuation tests} \label{sec:fluctests} 246 247The generalized fluctuation tests fit a model to the given data and derive an 248empirical process, that captures the fluctuation either in residuals or in 249estimates. For these empirical processes the limiting processes are known, so 250that boundaries can be computed, whose 251crossing probability under the null hypothesis is $\alpha$. If the empirical 252process path crosses these boundaries, the fluctuation is improbably large and 253hence the null hypothesis should be rejected (at significance level $\alpha$). 254 255\subsection{Empirical fluctuation processes: function \texttt{efp}} 256 257Given a formula that describes a linear regression model to be tested the 258function {\tt efp} creates an object of class {\tt "efp"} which contains 259a fitted empirical fluctuation process of a specified type. The types 260available will be described in detail in this section.\\ 261 262{\bf CUSUM processes}: 263The first type of processes that can be computed are CUSUM processes, which 264contain cumulative sums of standardized residuals. 265\cite{Z:Brown+Durbin+Evans:1975} suggested to consider cumulative sums of 266recursive residuals: 267\begin{equation}\label{Rec-CUSUM} W_n(t) \; = \; \frac{1}{\tilde \sigma 268\sqrt{\eta}}\sum_{i=k+1}^{k +\lfloor t\eta \rfloor} \tilde u_i \qquad (0 269\le t \le 1),\end{equation} 270where $\eta = n-k$ is the number of recursive residuals and $\lfloor t\eta 271\rfloor$ is the integer part of $t\eta$.\\ 272 273Under the null hypothesis the limiting process for the empirical fluctuation 274process $W_n(t)$ is the Standard Brownian Motion (or Wiener Process) $W(t)$. 275More precisely the following functional central limit theorem (FCLT) holds: 276\begin{equation}\label{lim(Rec-CUSUM)} W_n \Longrightarrow W, \end{equation} 277as $n \rightarrow \infty$, where $\Rightarrow$ denotes weak convergence 278of the associated probability measures.\\ 279 280Under the alternative, if there is just a single structural change point $t_0$, 281the recursive residuals will only have zero mean up to $t_0$. Hence the path of 282the process should be close to 0 up to $t_0$ and leave its mean afterwards. 283\cite{Z:Kraemer+Ploberger+Alt:1988} 284show that the main properties of the CUSUM quantity remain the same even under 285weaker assumptions, in particular in dynamic models. Therefore {\tt efp} has the 286logical argument {\tt dynamic}; if set to {\tt TRUE} the lagged observations 287$y_{t-1}$ will be included as regressors.\\ 288 289\cite{Z:Ploberger+Kraemer:1992} suggested to 290base a structural change test on cumulative sums of the common OLS residuals. 291Thus, the OLS-CUSUM type empirical fluctuation process is defined by: 292\begin{equation} \label{OLS-CUSUM} W_n^0(t) 293\quad = \quad \frac{1}{\hat \sigma \sqrt{n}} \sum_{i=1}^{\lfloor nt \rfloor} 294\hat u_i \qquad (0 \le t \le 1). \end{equation} 295The limiting process for $W_n^0(t)$ is the standard Brownian bridge $W^0(t) = 296W(t) - t W(1)$. It starts in 0 at $t = 0$ and it also returns to 0 for $t = 2971$. Under a single structural shift alternative the path should have a peak 298around $t_0$.\\ 299 300These processes are available in the function {\tt efp} by specifying the 301argument {\tt type} to be either {\tt "Rec-CUSUM"} or {\tt "OLS-CUSUM"}, 302respectively.\\ 303 304{\bf MOSUM processes}: 305Another possibility to detect a 306structural change is to analyze moving sums of residuals (instead of using 307cumulative sums of the same residuals). 308The resulting empirical fluctuation process does then not contain the sum of all 309residuals up to a certain time~$t$ but the sum of a fixed number of residuals in 310a data window whose size is determined by the bandwidth parameter $h \in (0,1)$ 311and which is moved over the whole sample period. Hence the Recursive MOSUM 312process is defined by 313\begin{eqnarray} \label{Rec-MOSUM} 314M_n(t | h) & = & \frac{1}{\tilde \sigma \sqrt{\eta}} \sum_{i = 315k+\lfloor N_\eta t \rfloor +1}^{k+\lfloor N_\eta t \rfloor +\lfloor 316\eta h \rfloor} \tilde u_i \qquad (0 \le t \le 1-h) \\ 317\label{Rec-MOSUM2} 318& = & W_n \left( \frac{ \lfloor N_\eta t \rfloor + \lfloor \eta 319h\rfloor}{\eta} \right) - W_n \left( \frac{ \lfloor N_\eta t 320\rfloor}{\eta} \right), \end{eqnarray} 321where $N_\eta = (\eta - \lfloor \eta h \rfloor)/(1-h)$. Similarly the OLS-based 322MOSUM process is defined by 323\begin{eqnarray} \label{OLS-MOSUM} 324M_n^0 (t | h) & = & \frac{1}{\hat \sigma \sqrt{n}} \left( \sum_{i = 325\lfloor N_n t \rfloor +1}^{\lfloor N_n t \rfloor +\lfloor nh \rfloor} \hat u_i 326\right) \qquad (0 \le t \le 1 - h) \\ 327\label{OLS-MOSUM2} 328& = & W_n^0 \left( \frac{ \lfloor N_n t \rfloor + \lfloor n 329h\rfloor}{n} \right) - W_n^0 \left( \frac{ \lfloor N_n t \rfloor}{n} 330\right), 331\end{eqnarray} 332where $N_n = (n - \lfloor n h \rfloor)/(1-h)$. As the representations 333(\ref{Rec-MOSUM2}) and (\ref{OLS-MOSUM2}) suggest, the limiting process 334for the empirical MOSUM processes are the increments of a Brownian motion or a 335Brownian bridge respectively. This is shown in detail in 336 \cite{Z:Chu+Hornik+Kuan:1995}.\\ 337 338If again a single structural shift is assumed at $t_0$, then both MOSUM paths 339should also have a strong shift around $t_0$.\\ 340 341The MOSUM processes will be computed if {\tt type} is set to {\tt "Rec-MOSUM"} 342or {\tt "OLS-MOSUM"}, respectively.\\ 343 344{\bf Estimates-based processes}: 345Instead of defining fluctuation processes on the basis of residuals they 346can be equally well based on estimates of the unknown regression coefficients. 347With the same ideas as for the residual-based CUSUM- and MOSUM-type processes 348the $k \times 1$-vector $\beta$ is either estimated recursively with a growing 349number of observations or with a moving data window of constant bandwidth $h$ 350and then compared to the estimates based on the whole sample. 351The former idea leads to the fluctuation process in the spirit of 352\cite{Z:Ploberger+Kraemer+Kontrus:1989} which is defined by 353\begin{equation} 354\label{fluctuation} Y_n \left(t \right) = \frac{\sqrt{i}}{\hat \sigma 355\sqrt{n}} \left({X^{(i)}}^\top X^{(i)} \right)^{\frac{1}{2}} \left( \hat 356\beta^{(i)} - \hat \beta^{(n)} \right), \end{equation} 357where $i = \lfloor k + t(n-k) \rfloor$ with $t \in [0,1]$. And the latter gives 358the moving estimates (ME) process introduced by \cite{Z:Chu+Hornik+Kuan:1995a}: 359\begin{equation} \label{ME} 360Z_n \left( \left. t \right| h \right) = \frac{\sqrt{\lfloor nh 361\rfloor}}{\hat \sigma \sqrt{n}} \left({X^{(\lfloor nt \rfloor, \lfloor nh 362 \rfloor)}}^\top X^{(\lfloor nt \rfloor, 363\lfloor nh \rfloor)} \right)^{\frac{1}{2}} \left( \hat \beta^{(\lfloor nt 364 \rfloor,\lfloor nh 365\rfloor)} - \hat \beta^{(n)} \right), \end{equation} 366where $0 \le t \le 1-h$. 367Both are $k$-dimensional empirical processes. Thus, the limiting processes are a 368$k$-dimensional Brownian Bridge or the increments thereof respectively. Instead 369of rescaling the processes for each $i$ they can also be standardized by 370$\left( {X^{(n)}}^\top X^{(n)} \right)^{\frac{1}{2}}$. This has the 371advantage that it has to be calculated only once, but \cite{Z:Kuan+Chen:1994} 372showed that if there are dependencies between the regressors the rescaling 373improves the empirical size of the resulting test. Heuristically the rescaled 374empirical fluctuation process ``looks'' more like its theoretic counterpart.\\ 375 376Under a single shift alternative the recursive estimates processes should 377have a peak and the moving estimates process 378should again have a shift close to the shift point $t_0$.\\ 379 380For {\tt type="fluctuation"} the function {\tt efp} returns the recursive 381estimates process, whereas for {\tt "ME"} the moving estimates process is 382returned.\\ 383 384All six processes may be fitted using the function {\tt efp}. 385For our example we want to fit an OLS-based CUSUM process, and a moving 386estimates (ME) process with bandwidth $h = 0.2$. The commands are simply 387<<efp>>= 388ocus <- efp(ecm.model, type="OLS-CUSUM", data=USIncExp2) 389me <- efp(ecm.model, type="ME", data=USIncExp2, h=0.2) 390@ 391These return objects of class {\tt "efp"} which contain mainly the 392empirical fluctuation processes and a few further attributes like the process 393type. The process itself is of class {\tt "ts"} (the basic time series class in 394\textsf{R}), which either preserves the time 395properties of the dependent variable if this is a time series (like in our 396example), or which is standardized to the interval $[0,1]$ (or a 397subinterval). For the MOSUM and ME processes the centered interval $[h/2, 3981-h/2]$ is chosen rather than $[0,1-h]$ as in (\ref{Rec-MOSUM}) and 399(\ref{OLS-MOSUM}).\\ 400 401Any other process type introduced in this section can be fitted by setting the 402{\tt type} argument. The fitted process can then be printed, plotted or tested 403with the corresponding test on structural change. For the latter appropriate 404boundaries are needed; the concept of boundaries for fluctuation processes is 405explained in the next section. 406 407\subsection{Boundaries and plotting} 408 409The idea that is common to all generalized fluctuation tests is that the null 410hypothesis of ``no structural change'' should be rejected when the fluctuation 411of the empirical process $\mathit{efp}(t)$ gets improbably large compared to the 412fluctuation of the limiting process. For the one-dimensional residual-based 413processes this comparison is performed by some appropriate boundary $b(t)$, that 414the limiting process just crosses with a given probability $\alpha$. Thus, if 415$\mathit{efp}(t)$ crosses either $b(t)$ or $-b(t)$ for any $t$ then it has to be 416concluded that the fluctuation is improbably large and the null hypothesis can 417be rejected at confidence level $\alpha$. The procedure for the $k$-dimensional 418estimates-based processes is similar, but instead of a boundary for the process 419itself a boundary for $||\mathit{efp}_i(t)||$ is used, 420where $||\cdot||$ is an appropriate functional which is applied component-wise. 421We have implemented the functionals `max' and `range'. The null hypothesis is 422rejected if $||\mathit{efp}_i(t)||$ gets larger than a constant $\lambda$, which 423depends on the confidence level $\alpha$, for any $i = 1, \dots, k$.\\ 424 425The boundaries for the MOSUM processes are also constants, i.e., of form $b(t) = 426\lambda$, which seems natural as the limiting processes are stationary. The 427situation for the CUSUM processes is different though. Both limiting processes, 428the Brownian motion and the Brownian bridge, respectively, are not stationary. 429It would seem natural to use boundaries that are proportional to the standard 430deviation function of the corresponding theoretic process, i.e., 431\begin{eqnarray} 432\label{bound:Rec-CUS1} b(t) &=& \lambda \cdot \sqrt{t}\\ 433\label{bound:OLS-CUS1} b(t) &=& \lambda \cdot \sqrt{t(1-t)} 434\end{eqnarray} 435for the Recursive CUSUM and the OLS-based CUSUM path respectively, 436where $\lambda$ determines the confidence level. But the 437boundaries that are commonly used are linear, because a closed form solution for 438the crossing probability is known. So the standard boundaries for the 439two proccess are of type 440\begin{eqnarray} 441\label{bound:Rec-CUS} b(t) &=& \lambda \cdot (1+2t)\\ 442\label{bound:OLS-CUS} b(t) &=& \lambda. 443\end{eqnarray} 444They were chosen because they are tangential to the boundaries 445(\ref{bound:Rec-CUS1}) and (\ref{bound:OLS-CUS1}) respectively in $t = 0.5$. 446However, \cite{Zo:Zeileis:2000a} examined the properties of the alternative 447boundaries (\ref{bound:Rec-CUS1}) and (\ref{bound:OLS-CUS1}) and showed that 448the resulting OLS-based CUSUM test has better power 449for structural changes early and late in the sample period.\\ 450 451Given a fitted empirical fluctuation process the boundaries can be 452computed very easily using the function {\tt boundary}, which returns a time 453series object with the same time properties as the given fluctuation process: 454<<efp-boundary>>= 455bound.ocus <- boundary(ocus, alpha=0.05) 456@ 457It is also rather convenient to plot the process with its boundaries for some 458confidence level $\alpha$ (by default 0.05) to see whether the path exceeds the 459boundaries or not: This is demonstrated in Figure~\ref{fig:ocus}. 460 461\begin{figure}[hbtp] 462\begin{center} 463<<OLS-CUSUM,fig=TRUE,height=4,width=6>>= 464plot(ocus) 465@ 466\caption{\label{fig:ocus} OLS-based CUSUM process} 467\end{center} 468\end{figure} 469 470It can be seen that the OLS-based CUSUM process exceeds its boundary; hence 471there is evidence for a structural change. Furthermore the process seems to 472indicate two changes: one in the first half of the 1990s and another one at the 473end of 1998.\\ 474 475It is also possible to suppress the boundaries and add them afterwards, e.g. in 476another color 477<<efp-boundary2,eval=FALSE>>= 478plot(ocus, boundary = FALSE) 479lines(bound.ocus, col = 4) 480lines(-bound.ocus, col = 4) 481@ 482For estimates-based processes it is only sensible to use time series plots if 483the functional `max' is used because it is equivalent to rejecting the null 484hypothesis when $\max_{i=1, \dots, k} ||\mathit{efp}(t)||$ gets large or when 485$\max_t \max_{i=1, \dots, k} \mathit{efp}_i(t)$ gets large. This again is 486equivalent to any one of the (one-dimensinal) processes $\mathit{efp}_i(t)$ for 487$i = 1, \dots, k$ exceeding the boundary. The $k$-dimensional process can also 488be plotted by specifying the parameter {\tt functional} (which defaults to {\tt 489"max"}) as {\tt NULL}: 490 491\begin{figure}[h] 492\begin{center} 493<<ME-null,fig=TRUE>>= 494plot(me, functional = NULL) 495@ 496\caption{\label{fig:me-null} 3-dimensional moving estimates process} 497\end{center} 498\end{figure} 499 500The output from \textsf{R} can be seen in Figure~\ref{fig:me-null}, where the 501three parts of the plot 502show the processes that correspond to the estimate of the 503regression coefficients of the intercept, the cointegration residuals and the 504increments of income, respectively. All three paths show two shifts: the first 505shift starts at the beginning of the sample period and ends in about 1991 and 506the second shift occurs at the very end of the sample period. The shift that 507causes the significance seems to be the strong first shift in the process for 508the intercept and the cointegration residuals, because these cross their 509boundaries. Thus, the ME test leads to similar results as the 510OLS-based CUSUM test, but provides a little more information about the nature of 511the structural change. 512 513\subsection{Significance testing with empirical fluctuation processes} 514 515Although calculating and plotting the empiricial fluctuation process with its 516boundaries provides and visualizes most of the information, it might still be 517necessary or desirable to carry out a traditional significance test. This can be 518done easily with the function {\tt sctest} (\underline{s}tructural 519\underline{c}hange \underline{test}) which returns an object of class {\tt 520"htest"} (\textsf{R}'s standard class for statistical test results) containing 521in particular the test statistic and the corresponding $p$ value. The test 522statistics reflect what was described by the crossing of boundaries in the 523previous section. Hence the test statistic is $S_r$ from (\ref{statr}) for the 524residual-based processes and $S_e$ from (\ref{state}) for the estimates-based 525processes: \begin{eqnarray} \label{statr} S_r & = & \max_t 526\frac{\mathit{efp}(t)}{f(t)},\\ \label{state} S_e & = & \max 527||\mathit{efp}(t)||, \end{eqnarray} 528where $f(t)$ depends on the shape of the boundary, i.e., 529$b(t) = \lambda \cdot f(t)$. For most boundaries is $f(t) \equiv 1$, but the 530linear boundary for the Recursive CUSUM test for example has shape $f(t) = 1 + 5312t$.\\ 532 533It is either possible to supply {\tt sctest} with a fitted 534empirical fluctuation process or with a formula describing the model that should 535be tested. Thus, the commands 536<<efp-sctest,eval=FALSE>>= 537sctest(ocus) 538@ 539and 540<<efp-sctest2>>= 541sctest(ecm.model, type="OLS-CUSUM", data=USIncExp2) 542@ 543lead to equivalent results. 544{\tt sctest} is a generic function which has methods not only for 545fluctuation tests, but all 546structural change tests (on historic data) introduced in this paper including 547the $F$ tests described in the next section. 548 549 550\section{$F$ tests} \label{sec:Ftests} 551 552A rather different approach to investigate whether the null hypothesis 553of ``no structural change'' holds, is to use $F$ test statistics. An important 554difference is that the alternative is specified: whereas the generalized 555fluctuation tests are suitable for various patterns of structural changes, the 556$F$ tests are designed to test against a single shift alternative. Thus, the 557alternative can be formulated on the basis of the model (\ref{model1}) 558\begin{equation} \label{H1} 559\beta_i = \left\{ \begin{array}{l} \beta_A \quad (1 \le i \le i_0) \\ 560\beta_B \quad (i_0 < i \le n) \end{array} \right., 561\end{equation} 562where $i_0$ is some change point in the interval $(k, n-k)$. \cite{Z:Chow:1960} 563was the first to suggest such a test on structural change for the case where the 564(potential) change point $i_0$ is known. 565He proposed to fit two separate regressions for the two subsamples defined by 566$i_0$ and to reject whenever 567\begin{equation} \label{Fstat} 568F_{i_0} \; = \; \frac{{\hat u}^\top {\hat u} - {\hat e}^\top {\hat 569e}}{{\hat e}^\top {\hat e}/(n-2k)}. 570\end{equation} 571is too large, where $\hat e = (\hat u_A, \hat u_B)^\top$ are the residuals from 572the full model, where the coefficients in the subsamples are estimated 573separately, and $\hat u$ are the residuals from the restricted model, where the 574parameters are just fitted once for all observations. The test statistic 575$F_{i_0}$ has an asymptotic $\chi^2$ distribution with $k$ degrees of freedom 576and (under the assumption of normality) $F_{i_0}/k$ has an exact $F$ 577distribution with $k$ and $n-2k$ degrees of freedom. The major drawback of this 578``Chow test'' is that the change point has to be known in advance, but there are 579tests based upon $F$ statistics (Chow statistics), that do not require a 580specification of a particular change point and which will be introduced in the 581following sections. 582 583\subsection{$F$ statistics: function \texttt{Fstats}} 584 585A natural idea to extend the ideas from the Chow test is to calculate the $F$ 586statistics for all potential change points or for all potential change points in 587an interval $[\ui, \oi]$ and to reject if any of those statistics get too large. 588Therefore the first step is to compute the $F$ statistics $F_i$ for $k < \ui \le 589i \le \oi < n-k$, which can be easily done using the function {\tt Fstats}. 590Again the model to be tested is specified by a formula interface and the 591parameters $\ui$ and $\oi$ are respresented by {\tt from} and {\tt to}, 592respectively. Alternatively to indices of observations these two parameters can 593also be specified by fractions of the sample; the default is to take {\tt from = 5940.15} and implicitly {\tt to = 0.85}. To compute the $F$ test statistics for all 595potential change points between January 1990 and June 1999 the appropriate 596command would be: 597<<Fstats>>= 598fs <- Fstats(ecm.model, from = c(1990, 1), to = c(1999,6), data = USIncExp2) 599@ 600This returns an object of class {\tt "Fstats"} which mainly contains a time 601series of $F$ statistics. Analogously to the empiricial fluctuation processes 602these objects can be printed, plotted and tested. 603 604\subsection{Boundaries and plotting} 605 606The computation of boundaries and plotting of $F$ statistics is rather similar 607to that of empirical fluctuation processes introduced in the previous section. 608Under the null hypthesis of no structural change boundaries can be computed 609such that the asymptotic probability that the supremum (or the mean) of the 610statistics $F_i$ (for $\ui \le i \le \oi$) exceeds this boundary is $\alpha$. 611So the following command 612\begin{figure}[htbp] 613\begin{center} 614<<Fstats-plot,fig=TRUE,height=4,width=6>>= 615plot(fs) 616@ 617\caption{\label{fig:Fstats} $F$ statistics} 618\end{center} 619\end{figure} 620plots the process of $F$ statistics with its boundary; 621the output can be seen in 622Figure~\ref{fig:Fstats}. As the $F$ statistics cross their boundary, there is 623evidence for a structural change (at the level $\alpha = 0.05$). The process 624has a clear peak in 1998, which mirrors 625the results from the analysis by empirical fluctuation processes and tests, 626respectively, that also indicated a break in the late 1990s.\\ 627 628It is also possible to plot the $p$ values instead of 629the $F$ statistics themselves by 630<<pval-plot,eval=FALSE>>= 631plot(fs, pval=TRUE) 632@ 633which leads to equivalent results. Furthermore it is also possible to set up 634the boundaries for the average instead of the supremum by: 635<<aveF-plot,eval=FALSE>>= 636plot(fs, aveF=TRUE) 637@ 638In this case another dashed line for the observed mean of the $F$ statistics 639will be drawn. 640 641\subsection{Significance testing with $F$ statistics} 642 643As already indicated in the previous section, there is more than one 644possibility to aggregate the series of $F$ statistics into a test statistic. 645\cite{Z:Andrews:1993} and \cite{Z:Andrews+Ploberger:1994} respectively suggested 646three different test statistics and examined their asymptotic distribution: 647\begin{eqnarray} 648\label{supF} \mbox{\textnormal{sup}}F & = & \sup_{\ui\le i \le \oi} F_i, \\ 649\label{aveF} \mbox{\textnormal{ave}}F & = & \frac{1}{\oi - \ui+ 1} 650\sum_{i = \ui}^{\oi} F_i, \\ 651\label{expF} \mbox{\textnormal{exp}}F & = & \log \left( \frac{1}{\oi - \ui+ 6521} \sum_{i = \ui}^{\oi} \exp ( 0.5 \cdot F_i) \right). 653\end{eqnarray} 654The sup$F$ statistic in (\ref{supF}) and the ave$F$ statistic from 655(\ref{aveF}) respectively reflect the testing procedures that have been 656described above. Either the null hypothesis is rejected when the maximal or the 657mean $F$ statistic gets too large. A third possibility is to reject when the 658exp$F$ statistic from (\ref{expF}) gets too large. The ave$F$ and exp$F$ test 659have certain optimality properties \citep{Z:Andrews+Ploberger:1994}. 660The tests can be carried out 661in the same way as the fluctuation tests: either by supplying the fitted {\tt 662Fstats} object or by a formula that describes the model to be tested. Hence the 663commands 664<<Fstats-sctest,eval=FALSE>>= 665sctest(fs, type="expF") 666@ 667and 668<<Fstats-sctest2>>= 669sctest(ecm.model, type = "expF", from = 49, to = 162, data = USIncExp2) 670@ 671lead to equivalent output. 672 673The $p$ values are computed based on \cite{Z:Hansen:1997}.\footnote{The authors 674thank Bruce Hansen, who wrote the original code for computing $p$ values for $F$ 675statistics in \textsf{GAUSS}, for putting his code at disposal for porting to 676\textsf{R}.} 677 678 679\section{Monitoring with the generalized fluctuation test} 680\label{sec:monitor} 681 682In the previous sections we were concerned with the retrospective detection 683of structural changes in \emph{given} data sets. Over the last years 684several structural change tests have been extended to monitoring 685of linear regression models where new data arrive over time 686\citep{Z:Chu+Stinchcombe+White:1996,Z:Leisch+Hornik+Kuan:2000}. Such 687forward looking tests are closely related to sequential tests. When 688new observations arrive, estimates are computed sequentially from all 689available data (historical sample plus newly arrived data) and compared 690to the estimate based only on the historical sample. As in the 691retrospective case, the hypothesis of no structural change is 692rejected if the difference between these two estimates gets too large. 693 694The standard linear regression model~(\ref{model1}) is generalized to 695\begin{equation} 696 y_i = x_i^\top \beta_i + u_i 697 \qquad (i = 1, \dots, n, n+1, \ldots), 698\end{equation} 699i.e., we expect new observations to arrive after time $n$ (when the 700monitoring begins). The sample $\{(x_1,y_1),\ldots,(x_n,y_n)\}$ will 701be called the \emph{historic sample}, the corresponding time period 702$1,\ldots,n$ the \emph{history period}. 703 704Currently monitoring has only been developed for recursive 705\citep{Z:Chu+Stinchcombe+White:1996} and moving 706\citep{Z:Leisch+Hornik+Kuan:2000} estimates tests. The respective 707limiting processes are---as in the retrospective case---the Brownian 708Bridge and increments of the Brownian Bridge. The empirical 709processes are rescaled to map the history period to the interval [0,1] 710of the Brownian Bridge. For recursive estimates there exists a closed 711form solution for boundary functions, such that the limiting Brownian 712Bridge stays within the boundaries on the interval $(1,\infty)$ with 713probability $1-\alpha$. Note that the monitoring period consisting of 714all data arriving after the history period corresponds to the Brownian 715Bridge after time 1. For moving estimates, only the growth rate of the 716boundaries can be derived analytically and critical values have to be 717simulated. 718 719Consider that we want to monitor our ECM during the 7201990s for structural change, using years 1986--1989 as the history 721period. First we cut the historic sample from the complete data set 722and create an object of class \texttt{"mefp"}: 723<<mefp>>= 724USIncExp3 <- window(USIncExp2, start = c(1986, 1), end = c(1989,12)) 725me.mefp <- mefp(ecm.model, type = "ME", data = USIncExp3, alpha = 0.05) 726@ 727Because monitoring is a sequential test procedure, the 728significance level has to be specified \emph{in advance}, i.e., when 729the object of class \texttt{"mefp"} is created. The \texttt{"mefp"} 730object can now be monitored repeatedly for structural changes. 731 732Let us assume we get new observations for the year 1990. Calling function 733\texttt{monitor} on \texttt{me.mefp} automatically updates our 734monitoring object for the new observations and runs a sequential test 735for structural change on each new observation (no structural break is 736detected in 1990): 737<<monitor1>>= 738USIncExp3 <- window(USIncExp2, start = c(1986, 1), end = c(1990,12)) 739me.mefp <- monitor(me.mefp) 740@ 741Then new data for the years 1991--2001 arrive and we repeat the 742monitoring: 743<<monitor2>>= 744USIncExp3 <- window(USIncExp2, start = c(1986, 1)) 745me.mefp <- monitor(me.mefp) 746me.mefp 747@ 748The software informs us that a structural break has been detected at 749observation \#72, which corresponds to December 1991. Boundary and 750plotting methods for \texttt{"mefp"} objects work (almost) exactly as 751their \texttt{"efp"} counterparts, only the significance level 752\texttt{alpha} cannot be specified, because it is specified when the 753\texttt{"mefp"} object is created. The output of 754\texttt{plot(me.mefp)} can be seen in Figure~\ref{fig:monitor}. 755 756\begin{figure}[htbp] 757\begin{center} 758<<monitor-plot,fig=TRUE,echo=FALSE,height=4,width=6>>= 759plot(me.mefp) 760@ 761\caption{\label{fig:monitor} Monitoring structural change with bandwidth $h=1$} 762\end{center} 763\end{figure} 764 765Instead of creating an {\tt "mefp"} object using the formula interface like 766above, it could also be done re-using an existing \texttt{"efp"} object, e.g.: 767<<mefp2>>= 768USIncExp3 <- window(USIncExp2, start = c(1986, 1), end = c(1989,12)) 769me.efp <- efp(ecm.model, type = "ME", data = USIncExp3, h = 0.5) 770me.mefp <- mefp(me.efp, alpha=0.05) 771@ 772If now again the new observations up to February 2001 arrive, we can monitor 773the data 774<<monitor3>>= 775USIncExp3 <- window(USIncExp2, start = c(1986, 1)) 776me.mefp <- monitor(me.mefp) 777@ 778and discover the structural change even two observations earlier as we used 779the bandwidth {\tt h=0.5} instead of {\tt h=1}. Due to this we have not one 780history estimate that is being compared with the new moving estimates, but we 781have a history process, which can be seen on the left in 782Figure~\ref{fig:monitor2}. This plot can simply be generated by 783\texttt{plot(me.mefp)}. 784 785\begin{figure}[htbp] 786\begin{center} 787<<monitor-plot2,fig=TRUE,height=4,width=6,echo=FALSE>>= 788plot(me.mefp) 789@ 790\caption{\label{fig:monitor2} Monitoring structural change with 791bandwidth $h=0.5$} \end{center} 792\end{figure} 793 794The results of the monitoring emphasize the results of the historic tests: the 795moving estimates process has two strong shifts, the first around 1992 and the 796second around 1998. 797 798\section{Conclusions} \label{sec:conclusions} 799 800In this paper, we have described the {\tt strucchange} package that implements 801methods for testing for structural change in linear regression relationships. 802It provides a unified framework for displaying information 803about structural changes flexibly and for assessing their significance 804according to various tests.\\ 805 806Containing tests from the generalized fluctuation test framework as well as 807tests based on $F$ statistics (Chow test ststistics) the package extends 808standard significance testing procedures: There are methods for fitting 809empirical fluctuation processes (CUSUM, MOSUM and estimates-based 810processes), computing an appropriate boundary, plotting these results and 811finally carrying out a formal significance test. Analogously a sequence of $F$ 812statistics with the corresponding boundary can be computed, plotted and tested. 813Finally the methods for estimates-based fluctuation processes have extensions to 814monitor incoming data. 815 816In addition to these methods for the linear regression model, the \texttt{strucchange} 817package contains infrastructure for testing, monitoring, and dating structural 818changes in other parametric models, e.g., estimated by maximum likelihood. Details 819about the underlying theory can be found in \cite{Z-papers:Zeileis:2005}, 820\cite{Z-papers:Zeileis+Hornik:2007}, and \cite{Z-papers:Zeileis+Shah+Patnaik:2010}. 821The corresponding functions in \texttt{strucchange} are presented in 822\cite{Z-papers:Zeileis+Kleiber+Kraemer:2003} and \cite{Z-papers:Zeileis:2006}. 823 824 825\section*{Acknowledgments} 826 827The research of Achim Zeileis, Friedrich Leisch and Kurt Hornik was supported 828by the Austrian Science Foundation (FWF) under grant SFB\#010 (`Adaptive 829Information Systems and Modeling in Economics and Management Science').\\ 830The work of Christian Kleiber was supported by the 831Deutsche Forschungsgemeinschaft, Sonderforschungsbereich 475. 832 833\bibliography{strucchange} 834\bibliographystyle{abbrvnat} 835 836\newpage 837\begin{appendix} 838\section{Implementation details for $p$ values} 839 840An important and useful tool concerning significance tests are $p$ values, 841especially for application in a software package. Their implementation is 842therefore crucial and in this section we will give more detail about the 843implementation in the {\tt strucchange} package.\\ 844 845For the CUSUM tests with linear boundaries there are rather good 846approximations to the asymptotic $p$ value functions given in 847\cite{Zo:Zeileis:2000a}. For the recursive estimates fluctuation test 848there is a series expansion, which is evaluated for the first hundred terms. For 849all other tests from the generalized fluctuation test framework the $p$ values 850are computed by linear interpolation from tabulated critical values. For the 851Recursive CUSUM test with alternative boundaries $p$ values from the interval 852$[0.001, 1]$ and $[0.001, 0.999]$ for the OLS-based version respectively are 853approximated from tables given in \cite{Zo:Zeileis:2000}. The critical values 854for the Recursive MOSUM test for levels in $[0.01, 0.2]$ are taken from 855\cite{Z:Chu+Hornik+Kuan:1995}, while the critical values for the levels 856in $[0.01, 0.1]$ for the OLS-based MOSUM and the ME test are given in 857\cite{Z:Chu+Hornik+Kuan:1995a}; the parameter $h$ is in both cases interpolated 858for values in $[0.05, 0.5]$.\\ 859 860The $p$ values for the sup$F$, ave$F$ and exp$F$ test are approximated based on 861\cite{Z:Hansen:1997}, who also wrote the original code in \textsf{GAUSS}, which 862we merely ported to \textsf{R}. The computation uses tabulated simulated 863regression coefficients. 864 865\end{appendix} 866 867\end{document} 868