The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.
19,757 characters
Detecting Cointegrating Relations in Non-stationary Matrix-Valued Time Series
\newtheorem{remark}{Remark}
\maketitle
\begin{abstract}
This paper proposes a Matrix Error Correction Model to identify cointegration relations in matrix-valued time series.
We hereby allow separate cointegrating relations along the rows and columns of the matrix-valued time series and use information criteria to select the cointegration ranks.
Through Monte Carlo simulations and a macroeconomic application, we demonstrate that our approach provides a reliable estimation of the number of cointegrating relationships.
\end{abstract}
\bigskip
\textsc{Keywords}: Matrix-valued time series, Cointegration rank, Error correction model, Information criteria
\newpage
\section{Introduction}
Understanding the long-run relationships between key macroeconomic variables is a central focus for many economists.
Recently, however, interest has turned to cointegration analysis for matrix-valued time series \citep{li2024coint}.
Modeling matrix-valued time series directly allows researchers to capture long-run relationships across multiple dimensions-- such as between countries (row dimension) and economic indicators (column dimension) --providing a more comprehensive understanding of how different economies interact, co-move, and adjust over time.
This paper introduces the Matrix Error Correction Model and demonstrates that information criteria can be reliably used to determine the cointegrating ranks among the different dimensions of the matrix-valued time series.
To fix ideas, we first review the standard framework for analyzing multivariate cointegrated systems.
Let $\bm y_t$, $t=1, \dots, T$, be an $N$ dimensional time series integrated of order I($1$).
In the presence of cointegration, the vector error correction model (VECM, see e.g., \citealp{johansen1990coint})
\begin{align}
\label{eq:vecm}
\bm \varDelta \bm y_t = \bm d + \bm \alpha \bm \beta' \bm y_{t-1} + \sum_{j=1}^{p} \bm \varPhi_j \bm \varDelta \bm y_{t-j} + \bm e_t,
\end{align}
captures long-run relationships, where $\bm d$ is the vector of deterministic terms,\footnote{For simplicity, the deterministic terms are left unrestricted and are not included within the cointegrating vector.} and $\bm \beta$ and $\bm \alpha$ are the $N \times r$ cointegrating matrix and adjustment coefficients respectively, with $r$ being the cointegrating rank.
Additionally, $\bm \varPhi_j$ represents the $j$th $N \times N$ short-run coefficient matrix, and $\bm e_t$ is the $N$-dimensional error term.
When $N$ is small-- typically ranging from two to five key variables --the canonical correlation approach of \citet{johansen1991metrica} can reliably determine the cointegration rank of the matrix $\bm \varPi = \bm \alpha \bm \beta'$.
However, as $N$ grows, this method becomes less reliable, necessitating alternatives to determine the cointegrating rank (e.g., \citealp{gutierrez2003power}, \citealp{wilms2016forecasting}).
To this end, we exploit the time series' matrix-valued structure, where $N_1$ and $N_2$ denote the number of time series in the rows and columns of the observed matrix over time.
For instance, our empirical analysis examines data from $N_1 = 3$ economic indicators and $N_2 = 4$ countries over $T=116$ observations, resulting in $N_1 N_2 = N = 12$ variables.
This matrix structure allows us to employ a Matrix Error Correction Model (MECM), which jointly accommodates the long-run dependencies across the two dimensions of the matrix, similar to the framework recently proposed by \citet{li2024coint}.
More precisely, we impose a Kronecker structure on the $\bm \alpha$, $\bm \beta$, and $\bm \varPhi_j$ terms of equation \eqref{eq:vecm} (see Section \ref{sec:MECM}).
This structure has three benefits over the traditional VECM.
First, it enables separate analysis of the cointegration dynamics across the rows and columns of the matrix, unlike the vectorized approach in equation \eqref{eq:vecm}.
We thus have a rank associated with the economic indicators ($r_1$, row rank) and a potentially different rank for the countries ($r_2$, column rank).
Second, we allow for a partial full-rank cointegrated system where only one of the two dimensions of the matrix-valued time series is rank-restricted.
Third, the resulting model is typically far more parsimonious, allowing for larger data sets than a traditional VECM.
While \citet{li2024coint} consider the MECM with fixed cointegration ranks mainly from a theoretical point of view, we complement their work by providing practitioners with practical tools, in the form of information criteria, to select the cointegration ranks $r_1$ and $r_2$ (see Section \ref{sec:estimationselection}).
Information criteria have been successfully applied in the cointegration literature \citep{Aznar2002selecting, cheng2009semiparametric} and can flexibly accommodate matrix-valued time series.
We demonstrate the good performance of the information criteria through a Monte Carlo simulation study in Section \ref{sec:simulation}.
Finally, our empirical results in Section \ref{sec:Application} reveal a single (restricted) long-run relationship between three economic indicators for the US, Germany, France, and Great Britain.
Replication material for the simulations and empirical analysis are available at \url{https://github.com/ivanuricardo/MECMrankdetermination}.
\section{The Matrix Error Correction Model}
\label{sec:MECM}
Let $\bm Y_t$ be an $N_1 \times N_2$ matrix-valued time series with I($1$) component series that follows the MECM($p$) model given by
\begin{align}
\label{eq:mecm}
\bm \varDelta \bm Y_t = \bm D + \bm U_1 \bm U_3' \bm Y_{t-1} \bm U_4 \bm U_2' + \sum_{j=1}^{p} \bm \varPhi_{1,j} \bm \varDelta \bm Y_{t-j} \bm \varPhi_{2,j}' + \bm E_t,
\end{align}
where $\bm D$ is the matrix of deterministic terms, $\bm U_3 \in \mathbb{R}^{N_1 \times r_1}$ and $\bm U_4 \in \mathbb{R}^{N_2 \times r_2}$ are the cointegrating matrices for the rows and columns of the matrix-valued time series, $\bm U_1 \in \mathbb{R}^{N_1 \times r_1}$ and $\bm U_2 \in \mathbb{R}^{N_2 \times r_2}$ are the corresponding adjustment coefficients, and $\bm \varPhi_{1,j} \in \mathbb{R}^{N_1 \times N_1}$ and $\bm \varPhi_{2,j} \in \mathbb{R}^{N_2 \times N_2}$ are the matrix autoregressive coefficients \citep{li2024coint}.\footnote{In the case of a stationary matrix AR($p$), different ranks are used for each matrix to distinguish between different right and left null space commonalities in the $\bm U_i$, see \citet{hecq2024reduced}.}
We assume the errors follow a matrix-valued normal distribution \citep{dawid1981matrix}, namely
\begin{align*}
\bm E_t \sim MVN(\mathbf{0}, \bm \varSigma_1, \bm \varSigma_2) \Leftrightarrow \text{vec}(\bm E_t) \sim N (\text{vec}(\mathbf{0}), \bm \varSigma_2 \otimes \bm \varSigma_1),
\end{align*}
where $\bm \varSigma_1 \in \mathbb{R}^{N_1 \times N_1}$ and $\bm \varSigma_2 \in \mathbb{R}^{N_2 \times N_2}$ are positive definite matrices capturing the relations between the rows and columns of the matrix-valued errors, $MVN(\cdot,\cdot,\cdot)$ denotes the matrix-valued normal distribution, and $N(\cdot,\cdot)$ denotes the multivariate normal distribution.
Reorganizing equation \eqref{eq:mecm} to the traditional vector-valued set-up, where $\text{vec}(\bm Y_t) = \bm y_t$, gives the restricted VECM
\begin{align*}
\bm \varDelta \bm y_t = \bm d + \underbrace{(\bm U_2 \otimes \bm U_1)}_{\bm \alpha} \underbrace{(\bm U_4 \otimes \bm U_3)'}_{\bm \beta'} \bm y_{t-1} + \sum_{j=1}^{p} \underbrace{(\bm \varPhi_{2,j} \otimes \bm \varPhi_{1,j})}_{\bm \varPhi_j} \bm \varDelta \bm y_{t-j} + \bm e_t.
\end{align*}
The MECM in equation \eqref{eq:mecm} thus implies a Kronecker structure on the adjustment coefficients $\bm \alpha$, the cointegrating matrix $\bm \beta$, and short-run coefficients $\bm \varPhi_j$ of equation \eqref{eq:vecm}.
The imposed Kronecker structure puts restrictions on the coefficients by separating the cointegrating relations and adjustment coefficients across the two dimensions of the matrix-valued time series.
This separation not only enhances interpretability but also yields a substantial reduction in the number of parameters to be estimated.
The total number of effective parameters, excluding the constant term, is given by
\begin{equation}
\label{eq:numpars}
\psi(r_1, r_2, p) = r_1 (2 N_1 - r_1) + r_2 (2 N_2 - r_2) + p (N_1^2 + N_2^2).
\end{equation}
As an example, consider a scenario where $(N_1, N_2) = (3, 4)$, $(r_1, r_2) = (1,1)$ and $p = 2$, then the MECM requires estimating $62$ parameters, whereas a comparable VECM with $N = 12$, $r = 1$, and $p = 2$ would require estimating $311$ parameters.
\begin{remark}
Model \eqref{eq:mecm} is not uniquely identifiable without additional restrictions on the parameters.
To resolve this, we impose that the top $r_1 \times r_1$ block of $\bm U_3$ and the top $r_2 \times r_2$ block of $\bm U_4$ are the identity matrix and set $\|\bm \varSigma_1\|_F = 1$ such that $\bm \varSigma_1$ is identified up to a sign change.
Identification restrictions on the short-run coefficient matrices $\bm \varPhi_j$ are not required, but in case the dimension-specific matrices $\bm \varPhi_{1,j}$ and $\bm \varPhi_{2,j}$ are of interest, the same restriction may be used, namely, $\| \bm \varPhi_{1,j}\|_F = 1$ for $j = 1, \dots, p$.
\end{remark}
\begin{remark}
We specify the short-run dynamics in equation \eqref{eq:mecm} as a matrix autoregression, following \citet{chen2021mar}.
It is, however, possible to replace this with unrestricted autoregressive components.
This substitution would increase the number of lagged autoregressive parameters from $ N_1^2 + N_2^2$ to $ N_1^2 N_2^2$.
\end{remark}
\begin{remark}
While the Kronecker product structure forms a natural way to reduce the dimensionality in MECMs for $N_1\times N_2$ matrix-valued data, note that it would be interesting to introduce a specification test to investigate whether the data support the presence of such a Kronecker product structure. To this end, an interesting avenue for future research would be to extend the specification test of \cite{chen2021mar} for autoregressive models with matrix-valued time series to the MECM model set-up.
\end{remark}
\begin{remark}
We deliberately focus on MECMs of moderate dimension in this paper.
The current MECM can, however, be extended to high-dimensional settings.
Examples of possible extensions include applying techniques for stationary high-dimensional matrix-valued time series \citep{wang2019factor} to the cointegration case (e.g., nonstationary factor models as in \citealp{trapani2025inference}), as well as adapting methods for vector error correction models (e.g., sparse methods in \citealp{liao2015automated, wilms2016forecasting} or low-rank procedures in \citealp{cubadda2023vecim}) to the high-dimensional MECM set-up.
\end{remark}
\section{Estimation and Selection of Cointegration Ranks}
\label{sec:estimationselection}
The log-likelihood (up to a constant) of the MECM($p$) model in equation \eqref{eq:mecm} for fixed rank $r_1$ and $r_2$ is given by
\begin{align}
\label{eq:loglike}
\mathcal{L}(\bm \varTheta) = -\frac{T N_1}{2} \log |\bm \varSigma_1| - \frac{T N_2}{2} \log |\bm \varSigma_2| - \frac{1}{2}\sum_{t=1}^T \text{tr} ( \bm \varSigma_1^{-1} (\bm \varDelta \bm Y_t - \bm U_1 \bm U_3' \bm Y_{t-1} \bm U_4 \bm U_2' - \sum_{j=1}^p \bm \varPhi_{1,j} \bm \varDelta \bm Y_{t-j} \bm \varPhi_{2,j}' - \bm D) \nonumber \\
\times \bm \varSigma_2^{-1} (\bm \varDelta \bm Y_t - \bm U_1 \bm U_3' \bm Y_{t-1} \bm U_4 \bm U_2' - \sum_{j=1}^p \bm \varPhi_{1,j} \bm \varDelta \bm Y_{t-j} \bm \varPhi_{2,j}' - \bm D)')
\end{align}
where $\bm \varTheta$ collects all parameters.
The objective function is non-convex, but gradient descent can be used to solve problem \eqref{eq:loglike} in a computationally efficient way (see Appendix \ref{sec:gradientdescent}).
In practice, however, the ranks $r_1$ and $r_2$ are unknown and must be selected.
To this end, we use standard information criteria, namely, the Akaike Information Criterion (AIC, \citealp{akaike1974new}) and Bayesian Information Criterion (BIC, \citealp{schwarz1978estimating})
\begin{align*}
\text{AIC}(r_1, r_2, p) &= -2 \mathcal{L}(\widehat{\bm \varTheta}) + 2 \psi(r_1, r_2, p), \\
\text{BIC}(r_1, r_2, p) &= -2 \mathcal{L}(\widehat{\bm \varTheta}) + \ln(T) \psi(r_1, r_2, p),
\end{align*}
where $\mathcal{L}(\widehat{\bm \varTheta})$ is the value of the log-likelihood at the estimated parameters and $\psi(r_1, r_2, p)$ denotes the effective number of parameters (see eq. \ref{eq:numpars}).
\setlength{\tabcolsep}{1em}
\begin{table}[t]
\centering
\caption{MECM(0): Rank selection with AIC or BIC for $T = 100$ and $T = 250$ observations.}
\label{tab:MECM0mat}
\begin{tabular}{
l
>{\centering\arraybackslash}p{2.4cm}
>{\centering\arraybackslash}p{2.4cm}
>{\centering\arraybackslash}p{2.4cm}
>{\centering\arraybackslash}p{2.4cm}
>{\centering\arraybackslash}p{2.4cm}
}
\toprule
\textbf{True Rank} & \textbf{Method} &\textbf{Average Ranks} & \textbf{Standard Deviation} & \textbf{Frequency Correct} \\
\midrule
(1,1) & AIC (100) & (1.02, 1.00) & (0.16, 0.00) & (0.98, 1.00) \\
& BIC (100) & (1.02, 1.00) & (0.13, 0.00) & (0.98, 1.00) \\
& AIC (250) & (1.01, 1.00) & (0.08, 0.00) & (0.99, 1.00) \\
& BIC (250) & (1.01, 1.00) & (0.04, 0.00) & (0.99, 1.00) \\
\hline
(3,1) & AIC (100) & (3.00, 1.07) & (0.00, 0.26) & (1.00, 0.93) \\
& BIC (100) & (3.00, 1.01) & (0.00, 0.04) & (1.00, 0.99) \\
& AIC (250) & (3.00, 1.03) & (0.00, 0.17) & (1.00, 0.97) \\
& BIC (250) & (3.00, 1.00) & (0.00, 0.00) & (1.00, 1.00) \\
\hline
(1,4) & AIC (100) & (1.13, 4.00) & (0.34, 0.03) & (0.87, 1.00) \\
& BIC (100) & (1.06, 4.00) & (0.24, 0.03) & (0.94, 1.00) \\
& AIC (250) & (1.03, 4.00) & (0.18, 0.00) & (0.97, 1.00) \\
& BIC (250) & (1.02, 4.00) & (0.14, 0.00) & (0.98, 1.00) \\
\hline
(3,4) & AIC (100) & (3.00, 4.00) & (0.00, 0.00) & (1.00, 1.00) \\
& BIC (100) & (3.00, 4.00) & (0.00, 0.00) & (1.00, 1.00) \\
& AIC (250) & (3.00, 4.00) & (0.00, 0.00) & (1.00, 1.00) \\
& BIC (250) & (3.00, 4.00) & (0.00, 0.00) & (1.00, 1.00) \\
\bottomrule
\end{tabular}
\end{table}
\section{Simulation Study}
\label{sec:simulation}
We conduct a simulation study to investigate the performance of our rank selection criteria.
The data-generating process (DGP) is examined under two scenarios.
The first is a MECM($0$), which omits the short-run dynamics from the model.
The second is a MECM($1$) with short-run dynamics.
Across all settings, we take $N_1 = 3$ and $N_2 = 4$ in line with our empirical application and generate $T+100$ observations from the MECM detailed in Appendix \ref{sec:signoise}, using the first $100$ as burn-in and taking $T=100$ and $T=250$.
The rank selection criteria then estimates MECMs across all possible combinations of the ranks ($r_1, r_2$), and selects the model with the lowest information criterion value.
We explore four settings of a reduced rank in the MECM: (i) fully reduced ($r_1 = r_2 = 1$), (ii) partially reduced first dimension ($r_1 = 1, r_2 = 4$), (iii) partially reduced second dimension ($r_1 = 3, r_2 = 1$), and (iv) no rank reduction ($r_1 = 3, r_2 = 4$).
To better understand the implications of these ranks in the context of our empirical application with $N_1=3$ economic indicators and $N_2=4$ countries, note that a rank of ($1,1$) indicates one cointegrating relation among the $N=12$ variables, suggesting that all indicators move together across different countries.
A rank of ($1,4$) reflects four distinct cointegrating relations, where all indicators co-move for each country separately.
A rank of ($3,1$) yields three cointegrating relations, where all countries co-move for each indicator separately.
Finally, with no rank reduction, the model captures a stationary process with no co-movements.
Table \ref{tab:MECM0mat} presents the results of the simulation study for the MECM($0$) DGP.
The results for the MECM($1$) DGP are similar and given in Table \ref{tab:MECM1mat} of Appendix \ref{sec:addsim}.
With a true rank of ($1,1$), all information criteria select the correct rank at a rate above $95\%$, with BIC performing best at a rate of over $95\%$.
As the number of observations grows, both AIC and BIC select the correct ranks more often.
Similar conclusions hold for cases with ranks ($3,1$) and ($1,4$).
AIC and BIC select the correct rank at least 85\% of the time with $T=100$, and it goes up to 100\% with $T=250$.
Under the full rank case, we always correctly select the ranks.
\section{Application}
\label{sec:Application}
\begin{figure}
\centering
\includegraphics[width=0.8\linewidth]{globalplots.png}
\caption{Time series plots for macroeconomic economic indicators (rows) and countries (columns).}
\label{fig:globalplot}
\end{figure}
We consider quarterly macroeconomic data from 1991Q1 to 2019Q4 ($T=116$) on $N_1 = 3$ economic indicators across $N_2 = 4$ countries.
This includes the log levels of real gross domestic product (GDP), the log levels of industrial production (PROD), and the levels of long-term interest rates for the United States (USA), Germany (DEU), France (FRA), and Great Britain (GBR).
The twelve time series are shown in Figure \ref{fig:globalplot}.
All variables are found to be integrated of order I($1$) based on Augmented Dickey-Fuller tests.
We estimate a MECM($1$) for the rank selection criteria and both AIC and BIC select ($1,1$) as cointegrating ranks.
This implies one restricted cointegrating relation across all 12 variables, which is visualized in Figure \ref{fig:cointvecs} without adjustment for the short-run dynamics.
The corresponding estimated cointegrating matrices and adjustment coefficients are given in Table \ref{tab:estcoefficients}.
\begin{figure}
\centering
\includegraphics[width=0.65\linewidth]{cointrel.png}
\caption{The cointegrated series $(\bm U_4 \otimes \bm U_3)' \text{vec}(\bm Y_{t-1})$ capturing the long-run dynamics among the 12 variables.}
\label{fig:cointvecs}
\end{figure}
From the estimated indicator-specific cointegrating vector $\widehat{\bm U}_3$, we observe that GDP positively co-moves with industrial production and long-term interest rates.
For the country-specific cointegrating vector $\widehat{\bm U}_4$, the cointegrating relationship between the USA and France is stronger than that between the USA and either Germany or Great Britain.
\begin{table}[t]
\centering
\setlength{\tabcolsep}{12pt}
\begin{tabular}{lcc|lcc}
\hline
\multicolumn{3}{c|}{\textbf{Indicators}} & \multicolumn{3}{c}{\textbf{Countries}} \\
\hline
\rule{0pt}{3ex}
& $\widehat{\bm U}_1$ & $\widehat{\bm U}_3$ & & $\widehat{\bm U}_2$ & $\widehat{\bm U}_4$ \\
\hline
\textbf{GDP} & -0.084 & 1.000 & \textbf{USA} & 0.088 & 1.000 \\
\textbf{PROD} & -0.201 & -0.099 & \textbf{DEU} & 0.774 & 0.055 \\
\textbf{IR} & -7.843 & -0.021 & \textbf{FRA} & 1.103 & -1.092 \\
& & & \textbf{GBR} & 0.674 & -0.196 \\
\hline
\end{tabular}
\caption{Estimated parameter values for indicators and countries.}
\label{tab:estcoefficients}
\end{table}
\bigskip
\noindent
{\bf Acknowledgements.}
We thank the editor and referee for their constructive comments which substantially improved the quality of the manuscript.
The last author was financially supported by the Dutch Research Council (NWO) under grant number VI.Vidi.211.032.
\bibliographystyle{asa}
\bibliography{updatedreferences}
\newpage