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.
63,714 characters · 12 sections · 74 citation commands
Recently, \citet*{Kalliovirta+Meitz+Saikkonen:2015} introduced a mixture autoregressive model based on Gaussian distribution with very attractive features. The Gaussian mixture autoregressive (GMAR) model has linear Gaussian autoregressions as its component models and mixing weights that, for a $p$th order model, depend on the full distribution of the $p$ past observations. The specific formulation of the mixing weights leads to ergodicity and full knowledge of the stationary distribution of $p+1$ consecutive observations. Moreover, it allows regime switches to depend on the level, variability, and temporal dependence of the past observations.
Meitz+Preve+Saikkonen:2018 proposed a mixture autoregressive model closely related to the GMAR model but based on Student's $t$-distribution. The Student's $t$ mixture autoregressive (StMAR) model has linear Student's $t$ autoregressions as its component models and mixing weights constructed analogously to the GMAR model, leading to similar theoretical and practical properties. The linear Student's $t$ autoregressions have the same form for the conditional mean as the Gaussian autoregressions (a linear function of the past observations) but different conditional variance. In particular, the conditional variances of the Student’s $t$ autoregressions depend on quadratic forms of past observations, whereas in the Gaussian case the conditional variances of the component models are constants. Utilization of the $t$-distribution does hence not only allow the StMAR model to account for larger kurtosis than the GMAR model but also stronger forms of conditional heteroskedasticity.
In this paper, we propose a generalization of the GMAR and StMAR models. The G-StMAR model accommodates both Gaussian autoregressions and Student's $t$ autoregressions as its component models, and its mixing weights are constructed analogously to the GMAR and StMAR models, leading to similar attractive features. It thus enables to model series which consist of regimes with time varying conditional variance and excess kurtosis as well as regimes with constant conditional variance and zero excess kurtosis. It turns out that the G-StMAR model is a limiting case of a StMAR model with the $t$-distributions of some regimes tending to normal distributions as the degrees of freedom parameters tend to infinity. As opposed to the limiting StMAR model, the advantage of the G-StMAR model is that it removes the redundant degrees of freedom parameters from the model and is free from numerical problems induced by weak identification of very large degrees of freedom parameters.
We demonstrate the usefulness of the G-StMAR model in an empirical application to the monthly U.S. interest rate spread between the 3-month Treasury bill (TB) rate and the effective federal funds (FF) rate. Our G-StMAR model identifies three regimes for the spread, with a GMAR type regime mainly appearing after the financial crisis in 2008 when the zero lower bound limits movements of the spread. The remaining regimes are of the StMAR type, one accommodating eras of low mean and high variability and the other high mean and moderate variability. The former StMAR type regime dominates often when the market possibly anticipates decreases in the FF rate or has increased preference for safety, whereas the latter one mostly prevails when the Fed is arguably not expected to significantly decrease the FF rate target. Our findings are consistent with \citet*{Sarno+Thornton:2003} who found that the FF rate seems to adjust to the TB rate, supporting the hypothesis that the market anticipates movements of the FF rate, moving the TB rate, and hence the spread, in advance.
The rest of this article is organized as follows. Section (ref) first introduces the component processes of the G-StMAR model and then proceeds to define the G-StMAR model and discusses its theoretical properties. Section (ref) discusses maximum likelihood (ML) estimation of the model parameters and establishes the asymptotic properties of the ML estimator. It is, in particular, discussed how the accompanying R package "uGMAR" uGMAR estimates the model parameters in practice with a two-phase procedure. Section (ref) describes a simple model selection procedure and discusses numerical consequences of very large degrees of freedom parameter estimates. Section (ref) presents the empirical application to the interest rate spread and Section (ref) concludes. Details of the estimation procedure employed by uGMAR, as well as proofs for the stated theorems are given in an Appendix.
Throughout this paper, we use the following notation. We write $\boldsymbol{x}=(x_1,...,x_n)$ for the column vector $\boldsymbol{x}$ where the components $x_i$ may be either scalars or (column) vectors. The notation $\boldsymbol{x}\sim n_d(\boldsymbol{\mu},\boldsymbol{\Gamma})$ signifies that the random vector $\boldsymbol{x}$ has a $d$-dimensional Gaussian distribution with mean $\boldsymbol{\mu}$ and (positive definite) covariance matrix $\boldsymbol{\Gamma}$. Similarly, $\boldsymbol{x}\sim t_d(\boldsymbol{\mu},\boldsymbol{\Gamma},\nu)$ signifies that $\boldsymbol{x}$ has a $d$-dimensional $t$-distribution with mean $\boldsymbol{\mu}$, (positive definite) covariance matrix $\boldsymbol{\Gamma}$, and degrees of freedom $\nu$ (assumed to satisfy $\nu>2$). The density functions and some properties of the multivariate Gaussian and Student's $t$-distributions are given in an Appendix. The vectorization operator $vec$ stacks columns of a matrix on top of each other and, $\iota_{d}$ is the $d$ dimensional vector $(1,0,...,0)$, $I_d$ signifies the identity matrix of dimension $d$, and $\otimes$ denotes the Kronecker product. Moreover, $\boldsymbol{1}_d$ and $\boldsymbol{0}_d$ denote $d$ dimensional vectors of ones and zeros, respectively.
We consider mixture autoregressive models in which each observation is generated by a mixture component that is randomly selected according to the probabilities pointed by the mixing weights. The mixture components are either (linear) conditionally homoscedastic Gaussian autoregressions as in the GMAR model Kalliovirta+Meitz+Saikkonen:2015 or conditionally heteroscedastic Student's $t$ autoregressions as in the StMAR model Meitz+Preve+Saikkonen:2018. The mixing weights are functions of the past observations constructed in a way that, for a $p$th order model, leads to ergodicity and full knowledge of the stationary distribution $p+1$ consecutive observations. Moreover, as the mixing weights depend on the full distribution of the past $p$ observations, they allow regime switches to depend on the level, variability, and temporal dependence of the past observations. In this section, we first introduce the component processes of the G-StMAR model and then proceed define of the G-StMAR model and discuss its properties.
To develop theory and notation, we first consider the component processes of the G-StMAR model. For a linear $p$th order Gaussian or Student $t$ autoregression $z_t$, we have
where $\sigma_t>0$, $\varphi_{0}\in\mathbb{R}$, and the autoregressive (AR) parameter $\boldsymbol{\varphi}=(\varphi_1,...,\varphi_p)$ satisfies the stationarity condition $\boldsymbol{\varphi}\in\mathbb{S}^p$ where
In the case of Gaussian autoregression, the distribution of the errors terms $\varepsilon_t$ is standard normal and $\sigma_t$ is a constant $\sigma$ for all $t$. Denoting $\boldsymbol{z}_t=(z_t,...,z_{t-p+1})$ and $\mu=\text{E}[z_t]$, $\gamma_j=\text{Cov}(z_t,z_{t-j})$, and $\boldsymbol{\gamma}_p=(\gamma_1,...,\gamma_p)$, it is well know that the stationary solution to ((ref)) for the Gaussian autoregression satisfies
where $\mu=\varphi_0/(1-\boldsymbol{\varphi}'\boldsymbol{1}_p)$ $\boldsymbol{\gamma}_p=\boldsymbol{\Gamma}_p\boldsymbol{\varphi}$, and the covariance matrices $\boldsymbol{\Gamma}_{p}$ and $\boldsymbol{\Gamma}_{p+1}$ are Toeplitz matrices given as (see, e.g., L\"utkepohl (2005), eq. (2.1.39))
Using the same notation as in ((ref))-((ref)) for $\boldsymbol{z}_{t-1}$, $\mu$, and $\boldsymbol{\Gamma}_p$, the Student's $t$ autoregressions utilized by Meitz+Preve+Saikkonen:2018 (which have also appeared at least in Spanos:1994 and Heracleus+Spanos:2006) are obtained by letting $\varepsilon_t \sim t_1(0,1,\nu+p)$ with $\nu>2$ in ((ref)) and defining
This definition (which requires the stationarity condition of the AR parameter) guarantees stationarity of the Student's $t$ autoregressions. Distributional properties of such stationary Student's $t$ autoregressions are similar to the Gaussian case, in particular Meitz+Preve+Saikkonen:2018,
The aforementioned properties of the component processes are essential in the following discussions and will be exploited implicitly. Gaussian component processes of the G-StMAR model are referred to as GMAR type and Student's $t$ component processes as StMAR type since they are identical to the component processes of the GMAR model Kalliovirta+Meitz+Saikkonen:2015 and the StMAR model Meitz+Preve+Saikkonen:2018, respectively.
Let $y_t$ ($t=1,2,...$) be the real valued time series of interest, and let $\mathcal{F}_{t-1}$ denote the $\sigma$-algebra generated by the random variables $\left\lbrace y_{t-j},j>0\right\rbrace$. For a G-StMAR model with $M$ mixture components and autoregressive order $p$, we have
where $\sigma_{m,t}>0$ are $\mathcal{F}_{t-1}$-measurable, $\varepsilon_{m,t}$ are independent of $\mathcal{F}_{t-1}$, $\varphi_{m,0}\in\mathbb{R}$, $\boldsymbol{\varphi}_m\in\mathbb{S}^p$ (the set $\mathbb{S}^p$ is defined in ((ref))), and $s_{1,t},...,s_{M,t}$ are unobservable regime variables such that for each $t$, exactly one of them takes the value one and the others take the value zero. Given the past of $y_t$, $s_{m,t}$ and $\varepsilon_{m,t}$ are assumed to be conditionally independent, and the conditional probability for regime $m$ occurring at the time $t$ is expressed in terms of the mixing weights $\alpha_{m,t}\equiv \text{Pr}\left(s_{m,t}=1\right|\mathcal{F}_{t-1})$ that satisfy $\sum_{m=1}^{M}\alpha_{m,t}=1$ (for all $t=1,2,...$). Each observation is thus generated by a linear autoregression corresponding to some (unobserved) mixture component $m$ which is selected randomly according to the probabilities determined by the mixing weights.
The first $M_1$ mixture components are (linear) Gaussian autoregressions and the rest $M_2\equiv M-M_1$ are Student's $t$ autoregressions. Regarding equation ((ref)), this means that for $m=1,...,M_1$, the terms $\varepsilon_{m,t}$ have standard normal distributions and the variances $\sigma_{m,t}^2$ are constants $\sigma_m^2$. For $m=M_1+1,...,M$, the terms $\varepsilon_{m,t}$ follow the $t$-distribution $t_1\left(0,1,\nu_m+p\right)$ and the variances $\sigma_{m,t}^2$ are as in equation ((ref)) except that $\boldsymbol{z}_{t-1}$ is replaced with $\boldsymbol{y}_{t-1}=(y_{t-1},...,y_{t-p})$ and the regime specific parameters $\varphi_{m,0},\boldsymbol{\varphi}_m,\sigma_m^2,\nu_m$ are used to define $\mu$ and $\bold{\Gamma}_p$ therein. The component specific conditional means $\mu_{m,t}$ are defined by equation ((ref)) for all the components.
Based on the above specifications, the conditional density function of a G-StMAR model with autoregressive order $p$ is given as
where the conditional densities $n_1(y_t;\mu_{m,t},\sigma_m^2)$ and $t_1\left(y_t;\mu_{m,t},\sigma_{m,t}^2,\nu_m+p\right)$ are obtained from the properties of the component processes (using the regime specific parameters). The form of the Student's $t$ density function in ((ref)) is given in online Appendix. The G-StMAR model adds to the class of mixture models introduced by Le+Martin+Raftery:1996 and further developed by Wong+Li:2000, Wong+Li:2001, Wong+Li2:2001, Glasbey:2001, Lanne+Saikkonen:2003, and Wong+Chan+Kam:2009, to name a few.
In order to specify the mixing weights $\alpha_{m,t}$ in ((ref)), we first define the following function for notational convenience. Let
where the $p$-dimensional densities $n_p(\boldsymbol{y};\mu_m\mathbf{1}_p,\boldsymbol{\Gamma}_m)$ and $t_p(\boldsymbol{y};\mu_m\mathbf{1}_p,\boldsymbol{\Gamma}_m,\nu_m)$ correspond to the stationary distribution of the $m$th component process (given in the equations ((ref)) and ((ref))). Denoting $\boldsymbol{y}_{t-1}=(y_{t-1},...,y_{t-p})$, the mixing weights of the G-StMAR model are defined as
where the parameters $\alpha_1,...,\alpha_M$ satisfy $\sum_{m=1}^M\alpha_m=1$. The mixing weights are thus weighted ratios of densities of the component processes corresponding to the $p$ previous observations. This specific definition of the mixing weights is appealing as it states that an observation is more likely to be generated from a regime with higher relative weighted likelihood. Moreover, it allows the probabilities of each regime occurring to depend on the level, variability, and temporal dependence of the past observations. This is not only convenient for forecasting but it also allows the researcher to associate specific characteristics to different regimes. It turns out that this formulation of the mixing weights also leads to attractive theoretical properties such as fully known stationary distribution of realizations $(y_t,...,y_{t-h})$, $h=0,1,...,p$, and ergodicity of the process. These theoretical properties are formally stated in Theorem (ref) below.
Before stating the theorem, a few notational conventions are provided. We collect the parameters of the G-StMAR model to a $(M(p+3)+M_2-1)\times1$ vector $\boldsymbol{\theta} \equiv (\boldsymbol{\theta}^-,\boldsymbol{\nu})$ where $\boldsymbol{\theta}^-=(\boldsymbol{\vartheta}_1,...,\boldsymbol{\vartheta}_M,\alpha_1,...,\alpha_{M-1})$, $\boldsymbol{\vartheta}_m=(\varphi_{m,0},\boldsymbol{\varphi}_m,\sigma^2_m)$, $\boldsymbol{\varphi}_m=(\varphi_{m,1},...,\varphi_{m,p})$, $m=1,...,M$, and $\boldsymbol{\nu}=(\nu_{M_1+1},...,\nu_M)$. The parameter $\alpha_M$ is omitted because it is obtained from the restriction $\sum_{m=1}^M\alpha_m=1$. The parameter space for the G-StMAR model is
where the restriction $\nu_m>2$ ($m=M_1+1,...,M$) is made to ensure existence of finite second moments and the set $\mathbb{S}^p$ is as in ((ref)). A G-StMAR model with autoregressive order $p$, $M_1$ GMAR type regimes, and $M_2$ StMAR type regimes is referred to as the G-StMAR($p,M_1,M_2$) model, whenever clarity of the presentation requires.
The stationary distribution of $\boldsymbol{y}_t$ is a mixture of $p$-dimensional normal and $t$-distributions with constant mixing weights $\alpha_m$. By the well known properties of the normal and the $t$-distribution, all its moments lower than $\min\lbrace\nu_{M_1+1},...,\nu_{M}\rbrace$ exist and are finite. Moreover, as shown in the proof of Theorem (ref), for any $h=0,1,...,p$, the marginal stationary distribution of the vector $(y_t,..,y_{t-h})$ is also a mixture of normal and $t$-distributions. This gives the parameters $\alpha_m$ an interpretation as the unconditional probabilities for the observation $y_t$ being generated from the $m$th component process. Similarly to the GMAR and the StMAR process, the mean, variance, and first $p$ autocovariances of $y_t$ are thus
where $\gamma_{m,j}$ is the $j$:th autocovariance of the $m$:th component process.
The conditional mean and variance of the G-StMAR process are obtained from the definition of the model as $ \text{E}[y_t|\mathcal{F}_{t-1}]=\sum_{m=1}^{M}\alpha_{m,t}\mu_{m,t} $ and
The conditional mean shares a common form with the GMAR model and StMAR model but differs from them in the definition of the mixing weights. The conditional variance includes three components; the first one is related to the conditional variances of the GMAR type components and the second one to the StMAR type components, whereas the third term encapsulates heteroskedasticity caused by variations in the conditional mean.
Notice that the GMAR model Kalliovirta+Meitz+Saikkonen:2015 can be obtained as a special of the G-StMAR model by setting $M_1=M$ and $M_2=0$, and similarly the StMAR model Meitz+Preve+Saikkonen:2018 is obtained by setting $M_1=0$ and $M_2=M$. We simply need to drop the corresponding terms from the formulas, and all the definitions and results stated in this and in the next section also hold for to the GMAR and StMAR models individually. However, some theory developed for the GMAR model, such as geometric ergodicity Kalliovirta+Meitz+Saikkonen:2015, has not been established for the StMAR and G-StMAR models. The GMAR model also requires less (currently) unverified assumptions than the StMAR and G-StMAR models for concluding asymptotic normality of the maximum likelihood estimator (see Kalliovirta+Meitz+Saikkonen:2015, Kalliovirta+Meitz+Saikkonen:2015, Section 2, Meitz+Preve+Saikkonen:2018, Meitz+Preve+Saikkonen:2018, Theorem 3, and Theorem (ref) of this paper)
Parameters of the G-StMAR model can be estimated with the method of maximum likelihood (ML). Because the stationary distribution of the process is known, the exact log-likelihood function can be used. Suppose the observed time series is $y_{-p+1},...,y_0,y_1,...,y_T$ and that the initial values are stationary. Then the log-likelihood function of the G-StMAR model takes the form
where
and the density functions $n_d(\cdot;\cdot)$ and $t_d\left(\cdot;\cdot\right)$ follow the notation described in Section (ref). If stationarity of the initial values seems unreasonable, one can condition on the initial values by dropping the first term on the right hand side of ((ref)) and base the estimation on the resulting conditional log-likelihood function.
In what follows, we assume estimation based on the conditional log-likelihood function $L_T^{(c)}(\boldsymbol{\theta})=T^{-1}\sum_{t=1}^Tl_t(\boldsymbol{\theta})$, i.e., that the ML estimator $\hat{\boldsymbol{\theta}}_T$ maximizes $L_T^{(c)}(\boldsymbol{\theta})$. We have scaled the conditional log-likelihood function with the sample size $T$ so that the notation is consistent with the referred literature.
To investigate the asymptotic properties of the ML estimator $\hat{\boldsymbol{\theta}}_T$, the parameter space $\boldsymbol{\Theta}$ given in ((ref)) needs to be restricted in a way that guarantees identification of the parameters. This amounts to requiring that components of the G-StMAR model cannot be "relabelled" so that one ends up with the same model with different parameter vector; that is,
The restrictions required to establish asymptotic properties of the ML estimator are summarized in the following assumption.
Asymptotic properties of the ML estimator under the conventional high-level conditions are stated in the following theorem (which is similar to Theorem 3 in Meitz+Preve+Saikkonen:2018 on the ML estimator of the StMAR model). Denote $\mathcal{I}(\boldsymbol{\theta}) = \text{E}\big[\frac{\partial l_t(\boldsymbol{\theta})}{\partial\boldsymbol{\theta}}\frac{\partial l_t(\boldsymbol{\theta})}{\partial\boldsymbol{\theta}'} \big]$ and $\mathcal{J}(\boldsymbol{\theta})=\text{E}\big[\frac{\partial^2 l_t(\boldsymbol{\theta})}{\partial\boldsymbol{\theta}\partial\boldsymbol{\theta}'}\big]$.
If one is willing to assume validity of the conditions (i)-(iii) of Theorem (ref), the ML estimator $\hat{\boldsymbol{\theta}}_T$ has the conventional limiting distribution, implying that approximative standard errors for the estimates are obtained as usual. Moreover, standard likelihood based tests are applicable as long as the orders $M_1$ and $M_2$ are correctly specified. If $M_1$ or $M_2$ is chosen too large, some of the parameters are not identified causing the result of Theorem (ref) to break down. This particularly happens when one tests for the number of regimes as the null hypothesis would imply that some regime is reduced from the model\footnote{Meitz+Saikkonen:2017 have, however, recently developed such tests for mixture models with Gaussian conditional densities.} Kalliovirta+Meitz+Saikkonen:2015. Similar caution also applies for testing whether a regime is of the GMAR type against the alternative that it is of the StMAR type, as under the null hypothesis $\nu_m=\infty$ for the StMAR type regime $m$ being tested, violating Assumption (ref). Numerical consequences of the weak identification of very large degrees of freedom parameters are briefly discussed in Section (ref).
Finding the ML estimates amounts to maximizing the log-likelihood function ((ref)) over the high dimensional parameter space ((ref)) satisfying several constraints. Due to the complexity of the log-likelihood function, finding an analytical solution is infeasible, so numerical optimization methods are required. The EM algorithm Redner+Walker:1984 has been a popular choice for estimating mixture models Wong+Li:2000, Wong+Li:2001, Wong+Li2:2001 as it is suitable for problems where all the data relevant to estimation is not observed (for mixture models that is the origin of each observation; in our case, the random variables $s_{1,t},...,s_{M,t}$ in ((ref))). For the G-StMAR model the EM algorithm is not, however, particularly useful because in each maximization step one faces a new optimization problem that is not much simpler than the original one. This is because in the G-StMAR model the mixing weights also depend on the AR parameters (in a complex way). Conventional gradient based algorithms, on the other hand, tend to converge to some local maximum near the starting point, making them generally insufficient for maximizing multimodal objective functions such as ((ref)) that require thorough exploration of the parameter space.
Several optimization algorithms capable of escaping from local maxima have been proposed for maximization of complicated multimodal objective functions. Such robust methods, which include simulated annealing and the genetic algorithm (see, e.g., \citeauthor*{Goffe+Ferrier+Rogers:1994}, Goffe+Ferrier+Rogers:1994 and Dorsey+Mayer:1995, Dorsey+Mayer:1995), often perform well but they are computationally heavy and tend to converge slowly when near the global maximum point Dorsey+Mayer:1995. Following Dorsey+Mayer:1995 Meitz+Preve+Saikkonen:2018,Meitz+Preve+Saikkonen2:2018, we hence suggest employing a hybrid estimation procedure where a genetic algorithm is used to find starting values for a gradient based method which then accurately converges to a nearby local maximum or saddle point.
Even with the two-phase estimation procedure, parameters of the G-StMAR model can be challenging to estimate. We have therefore accompanied this paper with the CRAN distributed R package "uGMAR" uGMAR in which the genetic algorithm has been modified to improve its performance.\footnote{In addition to the G-StMAR model, uGMAR also accomodates the GMAR and StMAR models.} Brief descriptions of the employed genetic algorithm and its modifications are given in an Appendix. After running the genetic algorithm, the estimation is finalized with a variable metric algorithm Nash:1990 using central difference approximation for the gradient of the log-likelihood function. Because of the presence of multiple local maxima, a (sometimes large) number of estimation rounds should be performed to obtain reliable results, for which uGMAR makes use of parallel computing to shorten the estimation time.
In empirical applications, building a G-StMAR model amounts to finding a suitable autoregressive order $p$, the number of GMAR type regimes $M_1$, and the number of StMAR type regimes $M_2$. Different strategies for choosing the number of each type of regimes may be considered depending on the application. We propose a simple model selection procedure which takes advantage of the observation that the G-StMAR model is a limiting case of the StMAR model\footnote{The definition of the StMAR($p,M$) model is technically the same as of the G-StMAR($p,0,M$) model.}.
It is easy to see that the linear Gaussian autoregression defined in Section (ref) is obtained as a limiting case of the Student's $t$ autoregression with the degrees of freedom parameter tending to infinity. As the mixing weights ((ref)) are weighted ratios of the component process densities, it then follows that the G-StMAR($p,M_1,M_2$) model is obtained as a limiting case of a StMAR($p,M$) model with the parameters $\nu_1,...,\nu_{M_1}$ limiting to infinity. Consequently, if a StMAR($p,M$) model is fitted to data generated by a G-StMAR($p,M_1,M-M_1$) process, then asymptotically, the $M_1$ regimes of the fitted StMAR model are expected to get large degrees of freedom estimates. We therefore suggest building a G-StMAR model by first finding a suitable StMAR model, and then estimating the appropriate G-StMAR model if the fitted StMAR model contains large degrees of freedom estimates. A StMAR model can be specified, for example, by using information criteria together with quantile residual diagnostics Kalliovirta:2012.
Overly large degrees of freedom estimates in a StMAR model are redundant but their weak identification also causes several inconveniences in numerical analysis of the model. They lead to nearly numerically singular Hessian matrix of the log-likelihood function when evaluated at the estimate, making the approximate standard errors often unavailable. Weakly identified degrees of freedom parameters also cause inconvenience in quantile residual based model diagnostics. In particular, the quantile residual tests proposed by Kalliovirta:2012 require a positive definite approximation of the Hessian matrix (evaluated the ML estimate). The tests are thus not applicable for StMAR models with too large degrees of freedom estimates, whereas they are for the corresponding G-StMAR models. Applicability of Kalliovirta:2012's (Kalliovirta:2012) tests, which take into account the uncertainty caused by estimation of the parameters, might have consequences in model selection when sheer graphical analysis of the quantile residuals fails to reveal inadequacies. We demonstrate such a case in the empirical application.
We consider the monthly U.S. interest rate spread between the 3-month Treasury bill (TB) secondary market rate and the effective federal funds (FF) rate, covering the period from 1954VII to 2019VII (781 observations). The series is plotted in Figure (ref) (top left) along with the 3-month TB and FF rates, and with the shaded areas indicating the periods of (NBER based) U.S. recessions. All the data were taken from the Federal Reserve Bank of St. Louis database.
Treasury bills are short-term pure discount bonds which are backed by the U.S. government and therefore generally considered to be almost free from default-risk. The effective federal funds rate is the averaged rate at which depository institutes loan federal funds to each other overnight. The overnight FF lending agreements are one of the most liquid financial asset, but unlike TBs, they are subject to a notable default-risk. The relationship between TB and FF rates has been studied, among others, by Simon:1990 and Sarno+Thornton:2003, while Kishor+Marfatia:2013 examine the relationship between TB and FF futures rate.
According to term structure theory, a long-term interest rate should reflect the current and expected future short-term rates, and also perceptions of risk and liquidity in the form of (possibly time varying) premium. Simon:1990 studied the predictive power of the weekly spread between the 3-month TB and FF rates on the future levels of the FF rate in 1972-1987. He argued that the current and expected future FF rates affect the spread between the TB and FF rates through the repurchase agreement (repo) market\footnote{In a repo, the borrower sells a security to the lender and agrees to repurchase it in the future (often in the next day). Effectively, repos function similarly to collateralized loans. See Baklanova+Copeland+McCaughrin:2015 for an overview of the U.S. repo market.} because repos are closely linked to the FF rate, and corporations with funds to invest can buy TBs alternatively to investing in consecutive overnight repos. TB rates are linked to the FF rates also because security dealers finance the bulk of their TB inventories in the repo market, which is closely tied to the FF market. Furthermore, when trust in solidity of the banking system weakens, the increased demand for safety lowers TB rates relatively to FF rates. Simon:1990 accounted for this by employing the spread between the 3-month Eurodollar time deposit\footnote{Eurodollar time deposit is a U.S. dollar-denominated deposit at a bank outside the U.S. with a fixed maturity.} and TB rates as a risk premium for bank safety. He found that the spread between the 3-month TB and FF rate had significant predictive power on future levels of the FF rate in the volatile nonborrowed reserves operating period (late 1979 - late 1982) but less or none in the other subperiods.
Sarno+Thornton:2003 identified an error correction model (ECM) between the daily 3-month TB and effective FF rate (covering the period from 1974 to 1999) and showed that their ECM, which allows for asymmetries and nonlinearities, outperforms the alternative of a linear ECM. One of their main findings was that the FF rate (which is controlled by the Fed) seems to adjust to the TB rate and not vice versa, supporting the hypothesis that the market anticipates changes in the FF rate, moving the TB rate in advance. Moreover, it appears that the adjustment speed depends on the sign and size of the deviation from the long-run equilibrium. Sarno+Thornton:2003 argued that although there has been a number of procedural changes affecting predictability of the FF rate, their results implicate that the changes have been statistically unimportant. Furthermore, their robustness checks indicate that their findings on the adjustments from disequilibria also hold for monthly data. Variations and asymmetries in the adjustment speed, on the other hand, indicate that the dynamics of the spread between the TB and FF rates might fluctuate along with the level of the spread. This suggests that a mixture model, such as the G-StMAR model, which is able encapsulate such behaviour could be an appropriate choice of model.
Kishor+Marfatia:2013 argued that the results in Sarno+Thornton:2003 are not very surprising since the effective FF rate always tends to revert back to the FF target rate, and it does not incorporate markets expectations of the changes in the future FF rate. To get around that, they studied the relationship between the 3-month TB rate and the 1-month FF futures rate which does incorporate information about market's anticipations on the future FF rate. They fitted a linear ECM to a daily series from 1989 to 2008, and found that the TB rate and the FF futures rate both seem to move to correct a short-run disequilibrium.
Interestingly, the spread between the 3-month TB rate and the effective FF rate is most of the time (covered in our sample period) negative. Sarno+Thornton:2003 made a similar observation for their daily series and suggested that only a small fraction of the negative difference could be attributed to the low default-risk of TBs, but that a more plausible explanation is that the interest on TBs is exempt from some local and state taxes. As the smaller taxes have larger effect on paid net interest (relative to interest paid on federal funds) when the interest rates are higher, some movements of the spread could be partially caused by the differences in taxation.
We employ the method of maximum likelihood based on the exact log-likelihood function for estimating the parameters of the considered models. Adequacy of the estimated models is examined using quantile residual diagnostics in the framework presented in Kalliovirta:2012. The quantile residuals of a correctly specified G-StMAR model are asymptotically independent with standard normal distributions Kalliovirta:2012, so they can be used for graphical analysis in a similar fashion to conventional Pearson's residuals. In addition to graphical analysis of the quantile residuals, we perform Kalliovirta:2012's (Kalliovirta:2012) asymptotic tests (which take into account the uncertainty caused by estimation of the parameters) for testing normality, autocorrelation, and conditional heteroskedasticity of the quantile residuals. The estimation, quantile residual diagnostics, and other numerical analysis of the models is conducted using the R package uGMAR uGMAR which is available through the CRAN repository.\footnote{There is also Matlab code available for the StMAR model in the form of StMAR MATLAB Toolbox by Meitz+Preve+Saikkonen2:2018.} uGMAR estimates the model parameters using the two-phase procedure described in Section (ref).
Following the model selection procedure described in Section (ref), we started by finding a suitable StMAR model. First, we estimated the StMAR($p, M$) model with one mixture component, $M = 1$, and autoregressive orders $p=1,...,24$ and found that the order $p = 6$ yields the largest likelihood. Adequacy of the StMAR($6, 1$) model was clearly rejected by the quantile residual tests (see Table 2), so we estimated the StMAR($p, M$) models with orders $p = 1,...,6$ and $M = 2,3$. The order $(p,M)=(5,2)$ minimized the Schwarz-Bayesian (BIC) and the Hannan-Quinn (HQIC) information criteria, whereas the Akaike's information criterion (AIC) was minimized by the order $(p,M)=(5,3)$. Inappropriate estimates extremely near the border of the stationarity region were discarded as they are not solutions of interest (but maximize the likelihood for rather a technical reason), so in such cases the next-best local maximum of the log-likelihood function was considered instead. In both the StMAR($5,2$) and the StMAR($5,3$) model, a very large degrees of freedom estimate for one regime was obtained (approximately $99000$ and $95000$, respectively), so we estimated the corresponding G-StMAR($5,1,1$) and G-StMAR($5,1,2$) models. Removing the weakly identified degrees of freedom parameters by switching to the G-StMAR models enabled us to compute approximate standard errors of the estimates and to calculate Kalliovirta:2012's (Kalliovirta:2012) test statistics (see Section (ref)). The values of the information criteria are reported in Table (ref) and the parameter estimates of the G-StMAR models are reported in Table (ref) with the approximate standard errors for the estimates in brackets.
Estimates regarding the GMAR type regime are quite similar for the two G-StMAR models, and their standard errors are relatively large. This is because for both of the models the GMAR type regime mainly occurs in the period of near-zero interest rates after 2008 and there are hence only few observations from that regime (regime 1 in Figure (ref), bottom left, which displays the mixing weights of the G-StMAR($5,1,2$) model; the mixing weights of the G-StMAR($5,1,1$) model are not shown). The three zeros in the variance parameter estimates (and in their standard errors) signify that the estimates (and their standard errors) round to zero in three digits accuracy\footnote{More accurate values for the ML estimate of $\sigma^2_1$ and its standard error are $3.237\times 10^{-4}$ and $6.884\times 10^{-5}$ for the G-StMAR($5,1,1$) model, $3.070\times 10^{-4}$ and $6.092\times 10^{-5}$ for the G-StMAR($5,1,2$) model, and $3.593\times 10^{-4}$ and $5.552\times 10^{-5}$ for the G-StMAR($5,1,2)^r$ model, respectively.}, implying that the GMAR type regime exhibits very low variability (conditionally and unconditionally). The small mixing weight parameter estimates, interpreted as the unconditional probability for the GMAR type regime occurring, reflect the observation that eras of such a low variability have been rare in the sample period. Also, a remarkably large standard error for the second regime's variance parameter sticks out for both of the models. Examination of the profile log-likelihood functions (not shown) does not, however, reveal anything notable.
Since the AR parameter estimates for the G-StMAR($5,1,2$) model are somewhat similar in all regimes, we estimated a StMAR($5,3$) model with the AR parameters restricted to be the same in all regimes, allowing for changes in the level, variability, and kurtosis only. The degrees of freedom estimate for one regime was very large (approximately $97000$), so we estimated the corresponding restricted G-StMAR model which we refer to as the G-StMAR($5,1,2)^r$ model. The parameter estimates of this model are also presented in Table (ref) with the related statistics, and the values of the information criteria in Table (ref). The standard errors of the AR parameters are notably smaller than in the non-restricted models because the AR parameters are common for all the regimes.
Figure (ref) presents the time series, normal quantile plot, and the sample autocorrelation function of the quantile residuals, and the sample autocorrelation function of the squared quantile residuals for the G-StMAR models presented in Table (ref). Graphical analysis of the quantile residuals does not show significant signs of inadequacy for any of the models. A slightly too fat lower tail in the quantile residuals' distributions and somewhat large, approximately $0.1$, sample autocorrelation at lag $12$ sticks out for each of the three models, however.
In order to further study adequacy of the models, we employed Kalliovirta:2012's (Kalliovirta:2012) tests, and tested for normality, autocorrelation, and conditional heteroskedasticity of the quantile residuals, taking into account $1,3,6$, and $12$ lags in the autocorrelation and heteroskedasticity tests. The $p$-values obtained from the tests are reported in Table (ref). The normality test rejects for all the three models at $1\%$ level of significance, possibly because of the fat lower tails in the quantile residuals' distributions. More interestingly, despite the similarities in the graphical analysis, the autocorrelation tests unambiguously reject adequacy of the G-StMAR($5,1,1$) model, whereas the $p$-values are reasonable for the G-StMAR($5,1,2$) model which also passes the heteroskedasticity tests. The $p$-values for the autocorrelation tests are rather small also for the restricted G-StMAR($5,1,2)^r$ model, which is preferred by the information criteria, showing some evidence of inadequacy. We therefore prefer the unrestricted G-StMAR($5,1,2$) model whose overall adequacy seems quite satisfactory. The fact that the restricted model has information criteria values superior to the unrestricted models, however, suggests that imposing the autocorrelation structure to be the same for all regimes would also be a reasonable modelling choice.\footnote{For comparison, we also estimated the GMAR($p,M$) model with orders $p=1,...,6$ and $M=1,...,4$. The values of the information criteria were, however, found inferior to our G-StMAR models, with the GMAR($3,4$) model minimizing BIC ($-432$) and the GMAR($5,4$) model minimizing HQIC ($-517$) and AIC ($-572$).}
Our model selection procedure led to the (unrestricted) G-StMAR($5,1,2$) model which identifies three statistical regimes for the spread between the 3-month TB secondary market rate and the effective FF rate. The mixing weights of the model are presented in Figure (ref) (bottom left) along with the interest rate spread series (top left). The GMAR type regime (red) dominates the period of near-zero interest rates occurring after 2008, where also the spread stays close to zero and exhibits very low variability. The second regime (green) identifies periods of high variability and low mean, spanning through most of the recessions, whereas the third regime (blue) often occurs\footnote{By a regime occurring at a point of time we mean that according to the estimated mixing weights, the process generated an observation from that regime with a probability close to one.} after the recessions when the spread moderately varies around zero. These characteristics of the regimes are also highlighted in Figure (ref) (right) where a kernel density estimate of the spread (black solid line) is presented with the model implied density (grey dashed line) and the regime densities (red, green, and blue dotted lines; regime densities are multiplied by the mixing weight parameter estimates $\alpha_m$, $m=1,2,3$). The model implied density matches fairly well to the skewed distribution of the observations, but peakiness of the distribution seems a bit exaggerated and the lower tail is not fat enough.
Based on our G-StMAR($5,1,2$) model, the regime specific unconditional mean of the spread varies from the $-0.06$ %-units of the first (GMAR type) and third regime to the $-0.52$ %-units of the second regime, with each regime regularly occurring for several consecutive months. As the second regime dominates during most of the recessions, and also often occurs before the recessions when the interest rates are relatively high, it seems plausible that part of the larger negative mean is explained by expectations of a decrease in the near-future FF rate. The third regime, on the other hand, mostly occurs after the recessions when the interest rates seem relatively low, possibly indicating that the larger mean of the regime could be related to the lack of expected decreases in the FF rate. These findings are consistent with Sarno+Thornton:2003 who found that the FF rate corrects disequilibriums from the long-run relationship, supporting the hypothesis that market's anticipations in the future movements of the FF rate are reflected in the TB rate.
Sarno+Thornton:2003 also found that the adjustment speed of FF rate towards the long-run equilibrium depends on the sign and size of the deviation. Specifically, FF rate below the long-run trend or larger deviation implies faster adjustment, suggesting that too high values of the spread would be corrected faster than too low values. This might partially explain why the low mean second regime usually occurs when the interest rates are declining, but a rise in the FF rate is not always accompanied with a switch to the higher mean third regime. Another possibility is that market's predictions on the future movements of the FF rate are sometimes rather poor or a premium has an increased effect on the opposite direction. During the savings and loan crisis in 80's and 90's, increased preference for the safety of TBs would seem like a plausible partial explanation for the moderately negative spread despite of the mainly increasing FF rate from late 1986 to early 1989.
Overall, the three statistical regimes of our G-StMAR model identify three economic regimes, with the first regime dominating the period in which the movements of the interest rates are limited by the zero lower bound. The second regime arguably occurs often when the market anticipates decreases in the FF rate or possibly has increased preferences for the safety of the almost default-risk free TBs. The third regime seems to mostly occur at times when the Fed is arguably not expected to significantly decrease the FF rate target (because the recession has already passed and the interest rates are relatively low).
This article introduced a mixture autoregressive model which is a combination of the Gaussian mixture autoregressive (GMAR) model Kalliovirta+Meitz+Saikkonen:2015 and the Student's $t$ mixture autoregressive (StMAR) model Meitz+Preve+Saikkonen:2018. This model, referred to as the G-StMAR model, has several attractive theoretical and practical properties that are analogous to those of the GMAR and StMAR model. In addition to discussing the properties, it was noted that estimating the parameters of the G-StMAR model can be challenging in practice. Following Dorsey+Mayer:1995 Meitz+Preve+Saikkonen:2018,Meitz+Preve+Saikkonen2:2018, we suggested using a two-phase estimation procedure where a genetic algorithm is used to find starting values for a gradient based method and accompied the paper with the R package uGMAR uGMAR which implements the two-phase estimation procedure with a modified version of a genetic algorithm.
We stated that the G-StMAR model is a limiting case of a StMAR model with some degrees of freedom parameters tending to infinity, and found that large degrees of freedom estimates in a StMAR model are not only redundant but also cause several inconveniences in numerical analysis of the model. In particular, weak identification of large degrees of freedom parameters was found to lead to numerically nearly singular approximation of the observed information matrix when evaluated at the estimate, making the approximate standard errors for the estimates and Kalliovirta:2012's (Kalliovirta:2012) diagnostic tests often unavailable. Removing the redundant degrees of freedom parameters by switching to a G-StMAR model was concluded to obviate the problems.
As an empirical application, we considered the monthly U.S. interest rate spread between the 3-month Treasury bill rate and the effective federal funds rate. Our G-StMAR model identified three regimes for the spread, with a switch from a StMAR type regime to a GMAR type regime arising from a switch in the economic regime, namely, to a regime where the zero lower bound limits the movements of the interest rates. The two StMAR type regimes accommodate eras of low mean and high variability and high mean and moderate variability. The first StMAR type regime arguably occurs often when the market anticipates decreases in the FF rate or possibly has increased preferences for safety, whereas the second one mostly occurs when the Fed is arguably not expected to significantly decrease the FF rate target. As opposed to modelling the series with a StMAR model containing an overly large degrees of freedom estimate, switching to the more parsimonious G-StMAR model allowed us to numerically compute approximate standard errors for the estimates, and moreover, to perform the Kalliovirta:2012's (Kalliovirta:2012) quantile residual tests which turned out to have significance in the model selection.
The author thanks Markku Lanne, Mika Meitz, and Pentti Saikkonen who commented the work and gave insightful suggestions that helped to improve the paper substantially. The author also thanks the Academy of Finland for financing the project.