\section{Introduction\label{sec:intro}} 

The vector autoregressive model \citep{Sims1980} is arguably the most fundamental model for multivariate time series \citep{Lutkepohl2005, Tsay2013}.
Applications of the model and its variants can be found in almost any field that involves learning the temporal dependency: economics and finance \citep{Wu_Xia2016}, energy forecasting \citep{Dowell_Pinson2016}, psychopathology \citep{Bringmann2013}, neuroscience \citep{Gorrostieta2012} and reinforcement learning \citep{Recht2018}, etc.

Consider the vector autoregressive model of order one, \textsc{var}(1), in the following form:
\begin{equation*
X_{t+1}= A X_t + \eta_t,
\end{equation*}
where $X_t\in\mathbb{R}^d$ is the observed time series, $A\in\mathbb{R}^{d\times d}$ is the unknown transition matrix, and $\eta_t\in\mathbb{R}^d$ is the innovation. In modern applications, the dimension $d$ is often relatively large. However, since the number of unknown parameters increases as  $d^2$, leading to problems such as over-parametrization, the  model cannot provide reliable estimates or forecasts without further restrictions \citep{Stock_Watson2001}.  A classical approach to dimensionality reduction for vector autoregressive models, which recently has enjoyed a resurgence of interest, advocates the incorporation of prior knowledge into  modelling.  For example, motivated by the fact that in spatio-temporal studies, it is often sufficient to collect information from `neighbors', \cite{Guo_Wang_Yao2016} proposed the banded vector autoregressive model, where the nonzero entries of $A$ are assumed to form a narrow band along the main diagonal, after arranging the $d$ components of $X_t$ by geographic location. To analyse users' time series over large social networks, the network vector autoregressive model of \cite{Zhu_Pan_Li_Liu_Wang2017}  used the follower-followee adjacency matrix to determine the zero-nonzero pattern of $A$, together with equality restrictions to further reduce the dimensionality. In fact, the above models can both be incorporated by the general framework of linear restrictions 
\begin{equation}\label{eq:lr}
\mathcal{C}\text{vec}(A^{ \mathrm{\scriptscriptstyle T} })=\mu,
\end{equation}
where $\mathcal{C}$ is a prespecified restriction matrix, $\mu$ is a known constant vector, and $A^{ \mathrm{\scriptscriptstyle T} }$ is the transpose of $A$. This form of restrictions is traditionally well known by time series modellers; see  books on multivariate time series analysis  such as \cite{Reinsel1993}, \cite{Lutkepohl2005} and \cite{Tsay2013}. 

Meanwhile, drawing inspiration from recent developments in high-dimensional regression, another well studied approach concerns the penalized estimation, where the modeller is agnostic to  the locations of nonzero coordinates in $A$ while assuming certain sparsity \citep{Davis_Zang_Zheng2015, Han_Lu_Liu2015, Basu_Michailidis2015}, or the directions of low-dimensional projections while, e.g., assuming a low rank structure of $A$ \citep{Ahn_Reinsel1988,Negahban_Wainwright2011}. Note that in the former case, once the locations of nonzero coordinates are identified, the model can be formulated as an instance of \eqref{eq:lr}. Although we focus on fully known restrictions, the framework of \eqref{eq:lr} allows us to study the theoretical properties of a much richer variety of restriction patterns in this paper.  

On the other hand, in the literature on large vector autoregressive models, there has been an almost exclusive focus on stable processes. Technically, this means requiring that the spectral radius $\rho(A)<1$, or often more stringently, that the spectral norm $\lVert A \rVert_2<1$.  However, the analysis of stable processes typically cannot be carried over to unstable processes. In this paper, we provide a novel finite-time (non-asymptotic) analysis of the ordinary least squares estimator for stable, unstable and even slightly explosive vector autoregressive models within the general framework of \eqref{eq:lr}. 

Our analysis sheds new light on the phase transition phenomenon of the ordinary least squares estimator across different stability regimes. This is made possible by adopting the non-asymptotic, non-mixing approach of \cite{Simchowitz2018}. Resting upon a generalization of Mendelson's \citeyearpar{Mendelson2014} small-ball method, this approach is particularly attractive because: (i) it unifies stable and unstable cases, whereas in asymptotic theory these two cases would require substantially different techniques; and (ii) in contrast to existing non-asymptotic methods, it can well capture the fundamental  trait that the estimation will be more accurate as $\rho(A)\rightarrow1$. While relaxing the normality assumption in \cite{Simchowitz2018}, we precisely characterize the impact of imposing restrictions on the estimation error.  More importantly, we reveal for the first time that the phase transition from slow to fast error rates depends on the the smallest singular value of $A$, a measure of the least excitable mode of the system.  In addition, we expand the applicable region  $\rho(A)\leq 1$ in the above paper to $\rho(A)\leq 1+c/n$, where  $c>0$ is a universal constant, so  slightly explosive processes are also included.
Compared to \cite{Guo_Wang_Yao2016} which focused on the case with $\lVert A \rVert_2<1$, our assumption on the innovation distribution is stronger, and our error rate for the stable regime is larger by a logarithmic factor which, however, can be dropped under the normality assumption; see  $\S$ \ref{subsec:var_ub} for details.  Note that although \cite{Zhu_Pan_Li_Liu_Wang2017} relied on even milder assumptions on the innovations than \cite{Guo_Wang_Yao2016},  they assumed that the number of unknown parameters is fixed, and hence their theoretical analysis is incomparable to ours; see Example \ref{ex:networkAR} in $\S$ \ref{subsec:var_example}.

Throughout we denote by $\lVert\cdot\rVert$ the Euclidean norm and by $\mathcal{S}^{d-1}=\{\omega\in\mathbb{R}^d: \lVert \omega \rVert=1\}$ the unit sphere in $\mathbb{R}^d$. For a real matrix $A=(a_{ij})$, we let $\lambda_{\max}(A)$ and $\lambda_{\min}(A)$ (or $\sigma_{\max}(A)$ and  $\sigma_{\min}(A)$) be its largest and smallest eigenvalues (or singular values), respectively; additionally, we let $\rho(A)=\vert \lambda_{\max}(A) \vert, \lVert A\rVert_2=\sup_{\lVert\omega\rVert=1} \lVert A\omega\rVert$ and  $\lVert A\rVert_F=(\sum_{i,j}a_{ij}^2)^{1/2}$ be the spectral radius, spectral norm and Frobenius norm of $A$, respectively. For $x\in\mathbb{R}$, let $\lfloor x\rfloor=\max\{k\in\mathbb{Z}: k\leq x\}$ and $\lceil x\rceil=\min\{k\in\mathbb{Z}: k\geq x\}$, where $\mathbb{Z}$ is the set of integers. We write $A\succ0$ (or $A\succeq 0$) if $A$ is a positive definite (or positive semidefinite) matrix. Moreover, for any real symmetric matrices $A$ and $B$, we write $A\prec B$ (or $A \preceq B$) if $B-A\succ0$ (or $B-A\succeq0$), and write $A\nprec B$ (or $A \npreceq B$) if $A\prec B$ (or $A \preceq B$)  does not hold.  For any quantities $X$ and $Y$, we write $X \gtrsim Y$ if there exists a universal constant  $c>0$ independent of $(n,d,m,R,k,\sigma,\delta)$, whose meaning will become clear later, such that $X \geq c Y$. 

\section{Linearly restricted stochastic regression \label{sec:Reg_upperbd}}
\subsection{Problem formulation\label{subsec:formulation}}
Consider a sequence of time-dependent covariate-response pairs $\{(X_t, Y_t)\}_{t=1}^n$ following
\begin{equation}\label{eq:model_regression}
Y_t=A_*X_t+ \eta_t, 
\end{equation}
where $Y_t, \eta_t \in\mathbb{R}^q, X_t\in\mathbb{R}^d$, $A_*\in\mathbb{R}^{q\times d}$, and $E(\eta_t)=0$.  In particular, \eqref{eq:model_regression} becomes the \textsc{var}(1) model  when $Y_t=X_{t+1}$. In model \eqref{eq:model_regression}, the process $\{X_t, t=1,2,\ldots\}$ is adapted to  the filtration \[\mathcal{F}_{t}=\sigma\{\eta_1,\ldots,\eta_{t-1},X_1,\ldots, X_{t}\}.\]

Let $\beta_*=\text{vec}(A_*^{ \mathrm{\scriptscriptstyle T} })\in\mathbb{R}^{N}$, where $N=q d$. The parameter space of the linearly restricted model can be defined as 
\[\mathcal{L} = \{\beta\in\mathbb{R}^N: \mathcal{C}\beta=\mu\},\]
where $\mathcal{C}$ is a known  $(N-m)\times N$ matrix of rank $N-m$, representing $N-m$ independent restrictions, and $\mu\in\mathbb{R}^{N-m}$ is a known constant vector which may simply be set to zero in practice.  Let $\mathcal{C}_{+}$ be an $m\times N$ complement of $\mathcal{C}$ such that $\mathcal{C}_{\text{full}}=(\mathcal{C}_{+}^{{ \mathrm{\scriptscriptstyle T} }}, \mathcal{C}^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }$ is invertible, with its inverse partitioned into two blocks as $\mathcal{C}_{\text{full}}^{-1}=(R, R_{+})$, where $R$ is the  matrix of the first $m$ columns of $\mathcal{C}_{\text{full}}^{-1}$. Additionally, define $\gamma=R_{+}\mu$. Then it holds  $\mathcal{C}\gamma=\mathcal{C}R_{+}\mu=\mu$.  Note that if $\mathcal{C}\beta=\mu$, then $\beta=\mathcal{C}_{\text{full}}^{-1}\mathcal{C}_{\text{full}}\beta=R\mathcal{C}_{+}\beta+R_{+}\mathcal{C}\beta=R\theta+\gamma$, where $\theta=\mathcal{C}_{+}\beta$. Conversely, for any $\theta\in\mathbb{R}^m$, if $\beta=R\theta+\gamma$, then $\mathcal{C}\beta=\mathcal{C}R\theta+\mathcal{C}\gamma=\mu$. Thus, we have 
\[
\mathcal{L} =\{R\theta+\gamma: \theta\in\mathbb{R}^m\},
\]
i.e., the linear space spanned by columns of the restriction matrix $R$, shifted by the constant vector $\gamma$. This immediately implies that, given $(R, \gamma)$, there exists a unique unrestricted parameter $\theta_*\in\mathbb{R}^m$ such that $\beta_*=R\theta_*+\gamma$. Note that  $\gamma=0$ if and only if $\mu=0$. Moreover, the unrestricted model corresponds to the special case where $R=I_N$ and $\gamma=0$. 

The following examples illustrate how the linear restrictions can be encoded by  $(R,\gamma)$ or $(\mathcal{C}, \mu)$, where without loss of generality we set $\mu=\gamma=0$. Let $\beta_{*i}$ denote the $i$th entry of $\beta_*$. 

\begin{example}[Zero restriction]
	The restriction $\beta_{*i}=0$ may be encoded by  setting the $i$th row of $R$ to zero, or  by setting a row of $\mathcal{C}$ to $(0,\ldots, 0, 1, 0,\ldots, 0)\in\mathbb{R}^N$, where the $i$th entry is one.
\end{example}

\begin{example}[Equality restriction]
	Consider the restriction $\beta_{*i}-\beta_{*j}=0$. Suppose that the value of $\beta_{*i}=\beta_{*j}$ is $\theta_{*k}$, the $k$th entry of $\theta_*$. Then this restriction may be encoded by setting both the $i$th and $j$th rows of $R$ to $(0,\ldots, 0, 1, 0,\ldots, 0)\in\mathbb{R}^m$, where the $k$th entry is one. Alternatively, we may set a row of $\mathcal{C}$ to the $1\times N$ vector $c(i,j)$ whose $\ell$th element is defined as $[c(i,j)]_\ell=1(\ell=i)-1(\ell=j)$, where $1(\cdot)$ is the indicator function.
\end{example}

Define $n\times q$ matrices $Y=(Y_1,\ldots, Y_n)^{ \mathrm{\scriptscriptstyle T} }$ and $E=(\eta_1,\ldots, \eta_n)^{ \mathrm{\scriptscriptstyle T} }$, and the $n\times d$ matrix $X=(X_1,\ldots, X_n)^{ \mathrm{\scriptscriptstyle T} }$. Then \eqref{eq:model_regression} has the matrix form $Y=XA_*^{ \mathrm{\scriptscriptstyle T} } +E$.  Let $y=\text{vec}(Y)$,  $\eta=\text{vec}(E)$, $Z=(I_q \otimes X)R$, and $\widetilde{y}=y-(I_q \otimes X)\gamma$.
By vectorization and reparameterization, we have
\[
\widetilde{y}=(I_q \otimes X)(\beta_*-\gamma)+\eta=Z\theta_*+\eta.
\]

As a result, the ordinary least squares estimator of $\beta_*$ for the linearly restricted model is
\begin{equation}\label{eq:ols_theta}
\widehat{\beta}=R\widehat{\theta}+\gamma, \quad \widehat{\theta}=\argmin_{\theta\in\mathbb{R}^m}\lVert \widetilde{y}- Z \theta\rVert^2.
\end{equation}
Notice that $Z\in\mathbb{R}^{q n\times m}$. To ensure the feasibility of \eqref{eq:ols_theta}, we need $q n \geq m$.  Let $R=(R_1^{ \mathrm{\scriptscriptstyle T} }, \ldots, R_q^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }$ and $\gamma=(\gamma_1^{ \mathrm{\scriptscriptstyle T} }, \ldots, \gamma_q^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }$, where $R_i$ are $d\times m$ matrices and $\gamma_i$ are $d\times 1$ vectors. Then, $A_*=(I_q \otimes \theta_*^{ \mathrm{\scriptscriptstyle T} })\widetilde{R}+G$, where
\begin{equation*
\widetilde{R}=(R_1,\ldots, R_q)^{ \mathrm{\scriptscriptstyle T} } \in\mathbb{R}^{m q\times d}, \quad G=(\gamma_1,\ldots, \gamma_q)^{ \mathrm{\scriptscriptstyle T} } \in\mathbb{R}^{q\times d}.
\end{equation*}
Consequently, the ordinary least squares estimator of $A$ is
$\widehat{A}=(I_q \otimes \widehat{\theta}^{ \mathrm{\scriptscriptstyle T} })\widetilde{R}+G$.

\subsection{General upper bounds analysis \label{subsec:sr_ub}}
To derive  upper estimation error bounds for the stochastic regression model in $\S$ \ref{subsec:formulation}, we begin by introducing a key technical ingredient, namely the block martingale small ball condition \citep{Simchowitz2018}. As a generalization of Mendelson's \citeyearpar{Mendelson2014} small-ball method to time-dependent data, this condition can be viewed as a non-asymptotic stability assumption for controlling the lower tail behavior of the Gram matrix $X^{ \mathrm{\scriptscriptstyle T} } X$, or $Z^{ \mathrm{\scriptscriptstyle T} } Z$ in our context. 

\begin{definition}[Block Martingale Small Ball Condition]\label{def:BMSB} (i) 	For a real-valued time series $\{X_t, t=1,2,\ldots\}$ adapted to the filtration $\{\mathcal{F}_t\}$, we say that $\{X_t\}$ satisfies the $(k,\nu,\alpha)$-\textsc{bmsb}\ condition if there exist an integer $k\geq1$ and constants $\nu>0$ and $\alpha\in(0,1)$ such that, for every integer $s\geq 0$,	$k^{-1}\sum_{t=1}^{k} {\mathrm{pr}}(|X_{s+t}|\geq \nu \mid \mathcal{F}_s) \geq \alpha$ with probability one. (ii) For a time series $\{X_t, t=1,2,\ldots\}$ taking values in $\mathbb{R}^d$, we say that  $\{X_t\}$ satisfies the $(k,\Gamma_{\smb},\alpha)$-\textsc{bmsb}\ condition if there exists  $0\prec\Gamma_{\smb}\in\mathbb{R}^{d\times d}$ such that,  for every $\omega\in\mathcal{S}^{d-1}$, the real-valued time series $\{\omega^{ \mathrm{\scriptscriptstyle T} } X_t, t=1,2,\ldots\}$ satisfies the $\{k,  (w^{ \mathrm{\scriptscriptstyle T} }\Gamma_{\smb}w)^{1/2}, \alpha\}$-\textsc{bmsb}\ condition.	
\end{definition}

The value of the probability $\alpha$ is unimportant for our purpose as long as it exists. The thresholding matrix $\Gamma_{\smb}$, or $\nu$ in the univariate case,   captures  the  average cumulative excitability over any size-$k$ block; e.g., if $\{X_t\}$ is a mean-zero vector autoregressive process, then $\Gamma_{\smb}$ will scale proportionally no less than $k^{-1}E(\sum_{t=1}^{k}X_{s+t}X_{s+t}^{ \mathrm{\scriptscriptstyle T} } \mid \mathcal{F}_s)$, which is constant for all $s$. Since every time point is associated with a new shock to the  process,  $E(X_{s+t} X_{s+t}^{ \mathrm{\scriptscriptstyle T} } \mid \mathcal{F}_s)$  will increase as $t$ increases. Consequently, $\Gamma_{\smb}$  will be monotonic increasing in $k$; see Lemma \ref{lem_bmsb} in $\S$ \ref{subsec:verify}. 	
Moreover, by aggregating all the size-$k$ blocks, $\Gamma_{\smb}$  will essentially become the lower bound of $n^{-1}\sum_{t=1}^{n}X_tX_t^{ \mathrm{\scriptscriptstyle T} }$.  Note that for the ordinary least squares estimation, a larger lower bound on the Gram matrix will yield a sharper estimation error bound. Thus, a larger block size $k$ is generally preferred; see Theorem \ref{thm3} in $\S$ \ref{sec:VAR_upperbd} for details.

Let $\Gamma_{\smb}$ and $\overline{\Gamma}$ be  $d\times d$ positive definite matrices, and denote
\begin{equation}\label{notation:Gammaline_R}
\underline{\Gamma}_{R}=R^{ \mathrm{\scriptscriptstyle T} } (I_q \otimes \Gamma_{\smb}) R, \quad \overline{\Gamma}_R=R^{ \mathrm{\scriptscriptstyle T} } (I_q \otimes \overline{\Gamma}) R.
\end{equation}
In our theoretical analysis, properly rescaled matrices $\underline{\Gamma}_{R}$ and $\overline{\Gamma}_R$ will serve as lower and upper bounds of the Gram matrix $Z^{ \mathrm{\scriptscriptstyle T} } Z$, respectively, and the covering numbers derived from them will give rise to the quantity $\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})$ in Theorem \ref{thm1} below.
The regularity conditions  underlying our upper bound analysis are listed as follows:

\begin{assumption}\label{assum:bmsb}
	The covariates process $\{X_t\}_{t=1}^n$ satisfies the $(k,\Gamma_{\smb},\alpha)$-\textsc{bmsb}\ condition.
\end{assumption}

\begin{assumption}\label{assum:upper_matrix}
	For any $\delta\in(0,1)$, there exists $\overline{\Gamma}_R$ defined as in \eqref{notation:Gammaline_R} such that  ${\mathrm{pr}}(Z^{ \mathrm{\scriptscriptstyle T} } Z \npreceq n\overline{\Gamma}_R)\leq \delta$, where $Z=(I_q \otimes X)R$, and  $\overline{\Gamma}_R$ is dependent on $\delta$.
\end{assumption}

\begin{assumption}\label{assum:error_subG}
	For every integer $t\geq 1$, $\eta_t\mid \mathcal{F}_{t}$ is mean-zero and $\sigma^2$-sub-Gaussian.  
\end{assumption}

\begin{theorem}\label{thm1} Let $\{(X_t, Y_t)\}_{t=1}^n$ be generated by the linearly restricted stochastic regression model.
	Fix $\delta\in(0,1)$. Suppose that Assumptions \ref{assum:bmsb}--\ref{assum:error_subG} hold, $0\prec\Gamma_{\smb}\preceq\overline{\Gamma}$, and 
	\begin{equation}\label{thm1_condition}
	n \geq \frac{9k}{\alpha^2} \left \{ m \log\frac{27}{\alpha} + \frac{1}{2}\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})+\log q+\log\frac{1}{\delta}\right \}.
	\end{equation}
	Then, with probability at least $1-3\delta$, we have
	\begin{equation*}
	\lVert \widehat{\beta} - \beta_* \rVert \leq  \frac{9\sigma}{\alpha} \left [\frac{\lambda_{\max}(R \underline{\Gamma}_R^{-1} R^{ \mathrm{\scriptscriptstyle T} })}{n}\left \{12m \log\frac{14}{\alpha}+9\log\det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1}) +6\log\frac{1}{\delta} \right \} \right]^{1/2}.
	\end{equation*}
\end{theorem}

Note that $\lVert \widehat{\beta} - \beta_* \rVert=\lVert \widehat{A} - A_* \rVert_F$. The result in Theorem \ref{thm1} is new even for the unrestricted stochastic regression, where $R=I_N$. For  vector autoregressive processes, we will specify the matrices $\underline{\Gamma}_R$ and $\overline{\Gamma}_R$ in $\S$ \ref{subsec:verify}, where $\underline{\Gamma}_R$ will depend on the block size $k$ through $\Gamma_{\smb}$. As $k$ increases, $\Gamma_{\smb}$ will become larger and hence the factor $\lambda_{\max}(R \underline{\Gamma}_R^{-1} R^{ \mathrm{\scriptscriptstyle T} })$ will become smaller, resulting in a sharper error bound.
However, $k$ cannot be too large due to condition \eqref{thm1_condition}. Specifically, this condition arises from applying  the Chernoff bound technique to lower bound the Gram matrix via aggregation of all the size-$k$ blocks, since the probability guarantee of the Chernoff bound will degrade as the number of blocks decreases; see Lemma \ref{lem_A2} in the supplementary material. 	
Therefore, to apply Theorem \ref{thm1} to vector autoregressive processes, a crucial step  will be to derive a feasible region for $k$ that guarantees condition \eqref{thm1_condition}; see $\S$ \ref{subsec:k_condition} for details. 

Similarly,  we can obtain an analogous upper bound for $\widehat{A}-A_*$ in the spectral norm:

\begin{proposition}\label{prop1} Let $\{(X_t, Y_t)\}_{t=1}^n$ be generated by the linearly restricted  stochastic regression model. Fix $\delta\in(0,1)$. Then,  under the conditions of Theorem \ref{thm1}, with probability at least $1-3\delta$, we have
	\begin{equation*}
	\lVert \widehat{A} - A_* \rVert_2 \leq \frac{9\sigma}{\alpha} \left [\frac{\lambda_{\max}\left (\sum_{i=1}^{q}R_i \underline{\Gamma}_R^{-1} R_i^{ \mathrm{\scriptscriptstyle T} } \right)}{n}\left \{12m \log\frac{14}{\alpha}+9\log\det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1}) +6\log\frac{1}{\delta} \right \} \right]^{1/2}.
	\end{equation*}
\end{proposition}


\section{Linearly restricted vector autoregression\label{sec:VAR_upperbd}}
\subsection{Representative examples\label{subsec:var_example}}
We begin by illustrating how  the formulation in $\S$ \ref{sec:Reg_upperbd} can be used to study vector autoregressive models. Four representative examples will be discussed: \textsc{var}($p$) model, banded vector autoregressive model, network vector autoregressive model, and pure unit root process. 

Consider  the \textsc{var}($1$) model, i.e., model \eqref{eq:model_regression} with $Y_t=X_{t+1}\in\mathbb{R}^d$:
\begin{equation}\label{var1_true}
X_{t+1}=A_*X_t+ \eta_t, 
\end{equation}
subject to $\beta_* = R\theta_*+\gamma$, where  $\beta_*=\text{vec}(A_*^{ \mathrm{\scriptscriptstyle T} })\in\mathbb{R}^{d^2}$, $R=(R_1^{ \mathrm{\scriptscriptstyle T} }, \ldots, R_d^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }\in\mathbb{R}^{d^2\times m}$ with $R_i$  being $d\times m$ matrices, $\theta_*\in\mathbb{R}^m$,  and $\gamma=(\gamma_1^{ \mathrm{\scriptscriptstyle T} }, \ldots, \gamma_d^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }\in\mathbb{R}^{d^2}$ with $\gamma_i\in\mathbb{R}^d$. 

\setcounter{example}{0}
\begin{example}[\textsc{var}($p$) model\label{ex:VARp}] Interestingly, vector autoregressive models of order $p<\infty$ can be viewed as linearly restricted \textsc{var}($1$) models. Consider the \textsc{var}($p$) model
	\begin{equation}\label{model:var_p}
	Z_{t+1}=A_{*1}Z_t+A_{*2}Z_{t-1}+\cdots+A_{*p}Z_{t-p+1}+\varepsilon_t,
	\end{equation}
	where $Z_t, \varepsilon_t \in\mathbb{R}^{d_0}$, and $A_{*i}\in\mathbb{R}^{d_0\times d_0}$ for $i=1,\ldots, p$. Denote $d=d_0 p$, $X_t=(Z_t^{ \mathrm{\scriptscriptstyle T} }, Z_{t-1}^{ \mathrm{\scriptscriptstyle T} }, \ldots, Z_{t-p+1}^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} } \in\mathbb{R}^d$, $\eta_t=(\varepsilon_t^{ \mathrm{\scriptscriptstyle T} }, 0, \ldots, 0)^{ \mathrm{\scriptscriptstyle T} }\in\mathbb{R}^d$, and 
	\begin{equation}\label{eq:var_p_A}
	A_*=\left(\begin{array}{cccc}
	A_{*1} & \cdots & A_{*p-1} & A_{*p}\\ 
	I_{d_0}& \cdots & 0        & 0     \\
	\vdots & \ddots & \vdots   & \vdots\\
	0      & \cdots & I_{d_0}  & 0     
	\end{array}\right) \in\mathbb{R}^{d\times d}.
	\end{equation}
	As a result, \eqref{model:var_p} can be written  into the \textsc{var}($1$) form in \eqref{var1_true}. As shown in  \eqref{eq:var_p_A},  all entries in the last $d-d_0$ rows of  $A_*$ are restricted to either zero or one. The restriction $\beta_{*i}=1$ can be encoded in $(\mathcal{C},\mu)$ in the same way as Example 1 in $\S$ \ref{subsec:formulation} but with the $i$th entry of $\mu$ set to one. Thus,  with or without restrictions, the \textsc{var}($p$) model can be studied by the same method as that for a linearly restricted \textsc{var}(1) model. Note that the special structure of the innovation $\eta_t$ that some entries of $\eta_t$ are fixed at zero will not pose any extra difficulty.
\end{example}

In the following  examples, we consider \textsc{var}(1) models with various structures for $A_*=(a_{*ij})_{d\times d}$, and set $\gamma=0$ so that the restrictions are in the form of $R\theta=0$:

\begin{example}[Banded vector autoregression\label{ex:bandedAR}] \cite{Guo_Wang_Yao2016} proposed the vector autoregressive model with the following zero restrictions:
	\begin{equation}\label{eq:banded}
	a_{*ij}=0,\hspace{5mm} \vert i-j\vert > k_0,
	\end{equation}
	where the integer  $1\leq k_0 \leq \lfloor (d-1)/2 \rfloor$ is called the bandwidth parameter. Let $b_{*i}\in\mathbb{R}^d$ be the transpose of the $i$th row of $A_*$. Hence, $\beta_*=(b_{*1}^{ \mathrm{\scriptscriptstyle T} },\ldots, b_{*d}^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }$. Note that the restrictions  are imposed on each $b_{*i}$ separately.
	As a result, the $b_{*i}$'s are determined by non-overlapping subsets of entries in $\theta_*$; that is, we can write $b_{*i}=R_{(i)}\vartheta_{*i}$, where $R_{(i)}\in\mathbb{R}^{d\times m_i}$, $\vartheta_{*i}\in\mathbb{R}^{m_i}$, $\theta_*=(\vartheta_{*1}^{ \mathrm{\scriptscriptstyle T} }, \ldots, \vartheta_{*d}^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }\in\mathbb{R}^m$, and $m=\sum_{i=1}^{d}m_i$. In this case, $R$ is a block diagonal matrix:
	\begin{equation*
	R=\left(\begin{array}{*{3}c}
	R_{(1)} &  			& 0\\
	& \ddots 	&  \\
	0		& 			& R_{(d)}
	\end{array} \right) \in \mathbb{R}^{d^2\times m},
	\end{equation*}
	and \eqref{eq:banded} can be encoded in $R$ as follows: (1) $m_i=k_0+i$ and $R_{(i)}=(I_{m_i},0)^{ \mathrm{\scriptscriptstyle T} }$ if $1\leq i \leq k_0+1$; (2) $m_i=2k_0+1$ and $R_{(i)}=(0_{m_i\times (i-k_0-1)}, I_{m_i}, 0_{m_i\times (d-i-k_0)})^{ \mathrm{\scriptscriptstyle T} }$ if $k_0+1<i<d-k_0$; and (3) $m_i=k_0+1+d-i$ and $R_{(i)}=(0,I_{m_i})^{ \mathrm{\scriptscriptstyle T} }$ if $d-k_0\leq i \leq d$.
\end{example}

\begin{example}[Network vector autoregression\label{ex:networkAR}] Consider the network model in \cite{Zhu_Pan_Li_Liu_Wang2017}. Let us drop the individual effect and the intercept to ease the notations. This model assumes that all diagonal entries of $A_*$ are equal: $a_{*ii}=\theta_{*1}$ for $1\leq i\leq d$. For the off-diagonal entries,  the zero-nonzero pattern is known and completely determined by the social network: $a_{*ij}\neq 0$ if and only if individual $i$ follows individual $j$. Moreover, all nonzero off-diagonal entries are assumed to be equal: $a_{*ij}=\theta_{*2}$ if $a_{*ij}\neq 0$, for $1\leq i\neq j \leq d$. This model is actually very parsimonious, with only $m=2$, while the network size $d$ can be extremely large. To incorporate the above restrictions, we may define the $d^2\times 2$ matrix $R$ as follows: for  $i=1,\ldots, d^2$, the $i$th row of $R$ is  $(1,0)$ if $\beta_{*i}$ corresponds to a diagonal entry of $A_*$,  $(0,1)$ if $\beta_{*i}$ corresponds to a nonzero off-diagonal entry of $A_*$, and $(0,0)$ if $\beta_{*i}$ corresponds to a zero off-diagonal entry of $A_*$. 
\end{example}

\begin{example}[Pure unit root process\label{ex:pure_unit_root}] Another simple but important case is $A_*=\rho I_d$ with $\rho\in\mathbb{R}$. Then, the smallest true model has $m=1$,  and the corresponding restrictions, $a_{*11}=\cdots=a_{*dd}$ and $a_{*ij}=0$ for $1\leq  i\neq j \leq d$, can be imposed by setting $R=(e_1^{ \mathrm{\scriptscriptstyle T} }, \ldots, e_d^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} } \in\mathbb{R}^{d^2}$, where $e_i$ is the $d\times 1$ vector with all elements zero except the $i$th being one. When $\rho=1$, the underlying model becomes the pure unit root process, a classic example of unstable vector autoregressive processes \citep{Hamilton1994}. In particular, the problem of testing $A_*=I_d$, or unit root testing in panel data, has been extensively studied in the asymptotic literature; see \cite{Chang2004} and \cite{Zhang_Pan_Gao2018} for studies in low and high dimensions, respectively. Note that \cite{Zhang_Pan_Gao2018} focused on asymptotic distributions of the largest eigenvalues of the sample covariance matrix of the pure unit root process under  $\lim_{n,d\rightarrow\infty}d/n=0$. It  does not involve parameter estimation and hence cannot be directly compared to this paper.	
\end{example}

The stochastic regression can also incorporate possibly time-dependent exogenous inputs such as individual effects \citep{Zhu_Pan_Li_Liu_Wang2017} and observable factors \citep{Zhou_Bose_Fan_Liu2018}, leading to the class of \textsc{varx}\ models \citep[see, e.g.,][]{Wilms_Basu_Bien_Matteson2017}. Since \textsc{varx}\ models can be analyzed similarly to vector autoregressive models, we do not pursue the details in this paper.

\subsection{Verification of Assumptions  \ref{assum:bmsb}--\ref{assum:error_subG} in Theorem \ref{thm1}\label{subsec:verify}}

In light of the generalisability to \textsc{var}($p$) models via the \textsc{var}($1$) representation,  to apply the general results in \S\ 2.2 to linearly restricted vector autoregressive models, it suffices to restrict our attention to the \textsc{var}($1$) model in \eqref{var1_true} from now on. 

Following the notations  in $\S$ \ref{sec:Reg_upperbd}, let  $Y=(X_2,\ldots, X_{n+1})^{ \mathrm{\scriptscriptstyle T} }$, $Z=(I_d\otimes X)R$, and $A_*=(I_d\otimes \theta_*^{ \mathrm{\scriptscriptstyle T} })\widetilde{R}+G$, where $\widetilde{R}=(R_1,\ldots, R_d)^{ \mathrm{\scriptscriptstyle T} }$ and $G=(\gamma_1,\ldots, \gamma_d)^{ \mathrm{\scriptscriptstyle T} }$. Note that $q=d$ for the \textsc{var}($1$) model.  In addition, $\{X_t\}$ is adapted to the filtration  
\[\mathcal{F}_{t}=\sigma\{\eta_1,\ldots,\eta_{t-1}\}.\]
The following conditions on $\{X_t\}$ will be invoked in our analysis:

\begin{assumption}\label{assum:var}
	(i) The process $\{X_t\}$ starts at $t=0$ with $X_0=0$; (ii) the innovations $\{\eta_t\}$ are independent and identically distributed with $E(\eta_t)=0$ and $\var(\eta_t)=\Sigma_\eta=\sigma^2 I_d$;  (iii) there is a universal constant $C_0>0$ such that, for every $\nu\in\mathbb{R}^{d}$ with $\nu^{ \mathrm{\scriptscriptstyle T} } \Sigma_\eta \nu\neq0$, the density of $\nu^{ \mathrm{\scriptscriptstyle T} } \eta_t/(\nu^{ \mathrm{\scriptscriptstyle T} } \Sigma_\eta \nu)^{1/2}$ is bounded from above by $C_0$ almost everywhere; and (iv) $\{\eta_t\}$ are $\sigma^2$-sub-Gaussian. 
\end{assumption}

Under Assumption \ref{assum:var}(i),  we can simply write $X_t$ in the finite-order moving average form, $X_t=\sum_{s=0}^{t-1}A_*^{s}\eta_{t-s-1}$ for any $t\geq1$. Then, by Assumption \ref{assum:var}(ii), $\var(X_t)=\sigma^2\Gamma_t$, where 
\begin{equation}\label{eq:Gramian}
\Gamma_t=\sum_{s=0}^{t-1} A_*^s (A_*^{ \mathrm{\scriptscriptstyle T} })^s
\end{equation}
is called the finite-time controllability Gramian \citep{Simchowitz2018}. Note that $\var(X_t)<\infty$ for whatever $A_*$.  By contrast, the typical setup in  asymptotic theory of stable processes assumes that  $\{X_t\}$ starts at $t=-\infty$. In this case, $\{X_t\}_{t\in\mathbb{Z}}$ has the infinite-order moving average form $X_t=\sum_{s=0}^{\infty}A_*^s\eta_{t-s-1}$, so $\var(X_t)=\sigma^2\sum_{s=0}^{\infty} A_*^s (A_*^{ \mathrm{\scriptscriptstyle T} })^s=\sigma^2\lim_{t\rightarrow\infty} \Gamma_{t}<\infty$
if and only if $\rho(A_*)<1$.   Thus, an important benefit of Assumption \ref{assum:var}(i) is that it allows us to  capture the possibly explosive behavior of $\var(X_t)$ over any finite time horizon, and derive upper bounds of the Gram matrix over different stability regimes; see Lemmas \ref{lem_gram1} and \ref{lem_gram} in this subsection.  

The condition $\Sigma_\eta=\sigma^2 I_d$ in Assumption \ref{assum:var}(ii) is imposed for simplicity. However, we can easily extend all the proofs in this paper to the general case with any symmetric matrix $\Sigma_\eta\succeq0$: we only need to rederive all results  with the role of $\Gamma_t$ replaced by $\sum_{s=0}^{t-1} A_*^s \Sigma_\eta (A_*^{ \mathrm{\scriptscriptstyle T} })^s$.
Assumption \ref{assum:var}(iii) is used to establish the block martingale small ball condition, i.e., Assumption \ref{assum:bmsb}, for vector autoregressive processes, and it allows us to  lower bound the small ball probability by leveraging Theorem 1.2 in \cite{Rudelson_Vershynin2015}   on densities of sums of independent random variables; see  also Remark  \ref{remark:normal}. Essentially,  Assumption \ref{assum:var}(iii) only requires that the distribution of any one-dimensional projection of the innovation is well spread on the real line.   Examples of such distributions include multivariate normal and  multivariate $t$ \citep{Kotz_Nadarajah2004} distribution and, more generally, elliptical distributions \citep{Fang_Kotz_Ng1990} with the consistency property in \cite{Kano1994}.  Lastly, it is clear that Assumptions \ref{assum:var}(ii) and (iv) guarantee Assumption \ref{assum:error_subG}.



\begin{remark}
	In asymptotic theory, stable and unstable processes require substantially different techniques, and results derived under $\rho(A_*)<1$ typically cannot be carried over to  unstable processes. For example, the convergence rate of the ordinary least squares estimator for fixed-dimensional unstable \textsc{var}(1) processes  is $n$ instead of $n^{1/2}$, and the limiting distribution is no longer  normal \citep{Hamilton1994}. 
\end{remark}

\begin{remark}
	The controllability Gramian $\Gamma_t$ is interpretable even without Assumption \ref{assum:var}(i). By recursion, we have $X_{s+t}=\sum_{\ell=0}^{t-1}A_*^{\ell}\eta_{s+t-\ell-1}+A_*^tX_{s}$ for any time point $s$ and  duration $t\geq1$. As a result, $\var(X_{s+t}\mid \mathcal{F}_{s})=\sum_{\ell=0}^{t-1} A_*^\ell \Sigma_\eta (A_*^{ \mathrm{\scriptscriptstyle T} })^\ell$, and it simply becomes $\sigma^2\Gamma_t$ if $\Sigma_\eta=\sigma^2 I_d$.  Note that $\var(X_{s+t}\mid \mathcal{F}_{s})$, or equivalently $\Gamma_t$, is a partial sum of a geometric sequence due to the autoregressive structure. Roughly speaking, larger $A_*$ means more persistent impact of $\eta_t$.  
\end{remark}

In the following, we present three lemmas for the linearly restricted vector autoregressive model. Lemma \ref{lem_bmsb} establishes the block martingale small ball condition by specifying  $\Gamma_{\smb}$, while  Lemmas \ref{lem_gram1} and \ref{lem_gram}  verify Assumption \ref{assum:upper_matrix} by providing two possible specifications of $\overline{\Gamma}$. 

Some additional notations to be used in Lemma \ref{lem_gram} are introduced as follows: Denote by $\Sigma_X$  the covariance matrix of the  $d n\times 1$ vector $\vect(X^{ \mathrm{\scriptscriptstyle T} })=(X_1^{ \mathrm{\scriptscriptstyle T} }, \dots, X_n^{ \mathrm{\scriptscriptstyle T} })^{ \mathrm{\scriptscriptstyle T} }$, so the $(t,s)$th $d\times d$ block of $\Sigma_X$ is $E(X_tX_s^{ \mathrm{\scriptscriptstyle T} })$, for $1\leq t,s\leq n$. Then define
\begin{equation}\label{xi}
\xi=\xi(m,d,n,\delta) =2 \left \{\frac{\lambda_{\max}(\Gamma_n)  \psi(m,d,\delta) \lVert \Sigma_X \rVert_2}{\sigma^2 n} \right \}^{1/2} +  \frac{2 \psi(m,d,\delta)\lVert \Sigma_X \rVert_2}{\sigma^2 n},
\end{equation}
where $\psi(m,d,\delta)=C_1 \{m\log 9+\log d+\log(2/\delta)\}$, and $C_1>0$ is a universal constant.


\begin{lemma}\label{lem_bmsb} 
	Suppose that $\{X_t\}_{t=1}^{n+1}$ follows $X_{t+1}=A_*X_t+\eta_t$ for $t=0,1,\dots, n$. Under Assumptions \ref{assum:var}(ii) and (iii), for any $1\leq k \leq \lfloor n/2 \rfloor$, $\{X_t\}_{t=1}^n$ satisfies the $(2k,\Gamma_{\smb},1/10)$-\textsc{bmsb}\ condition, where $\Gamma_{\smb}=\sigma^2\Gamma_{k}/ (4C_0)^2$. 
\end{lemma}

\begin{lemma}\label{lem_gram1}
	Let $\{X_t\}_{t=1}^{n+1}$ be generated by the linearly restricted vector autoregressive model.  Under Assumptions \ref{assum:var}(i) and (ii), for any $\delta\in(0,1)$, it holds  ${\mathrm{pr}}(Z^{ \mathrm{\scriptscriptstyle T} } Z \npreceq n\overline{\Gamma}_R)\leq \delta$, where  $\overline{\Gamma}_R=R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \overline{\Gamma}) R$, with $\overline{\Gamma}=\sigma^2 m\Gamma_{n}/\delta$.
\end{lemma}


\begin{lemma}\label{lem_gram}
	Let $\{X_t\}_{t=1}^{n+1}$ be generated by the linearly restricted vector autoregressive model. Under Assumptions \ref{assum:var}(i) and (ii), if $\{\eta_t\}$ are normal, then for any $\delta\in(0,1)$, it holds  ${\mathrm{pr}}(Z^{ \mathrm{\scriptscriptstyle T} } Z \npreceq n\overline{\Gamma}_R)\leq \delta$, where  $\overline{\Gamma}_R=R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \overline{\Gamma}) R$, with $\overline{\Gamma}=\sigma^2 \Gamma_{n}+ \sigma^2 \xi I_d$, and $\xi=\xi(m,d,n,\delta)$ is defined in \eqref{xi}.
\end{lemma}	

\begin{remark}\label{remark:normal}
	Unlike \cite{Simchowitz2018}, by leveraging \cite{Rudelson_Vershynin2015}, we establish the  block martingale small ball condition without the normality assumption. If Assumption \ref{assum:var}(ii) is relaxed to the general $\var(\eta_t)=\Sigma_\eta\succeq 0$, by a straightforward extension of the proof of Lemma \ref{lem_bmsb}, we can show that Lemma \ref{lem_bmsb} holds with $\Gamma_{\smb}=\sum_{\ell=0}^{k-1} A_*^\ell \Sigma_\eta (A_*^{ \mathrm{\scriptscriptstyle T} })^\ell/ (4C_0)^2$.  

\end{remark}

\begin{remark}\label{remark:hansonw}
	Lemma \ref{lem_gram1} is  a simple consequence of the Markov inequality and the property that $\lambda_{\max}(\cdot)\leq \tr(\cdot)$, so no distributional assumption on $\eta_t$ is required. However,  Lemma \ref{lem_gram} relies on the Hanson-Wright inequality \citep{Vershynin2018}, where the normality assumption is invoked. Although adopting the $\overline{\Gamma}$ in  Lemma \ref{lem_gram}  can eliminate a factor of $\log m$ in the resulting estimation error bounds, the $\overline{\Gamma}$ in Lemma \ref{lem_gram1} actually leads to sharper bounds under certain conditions on $A_*$; see $\S$ \ref{subsec:k_condition} and $\S$ \ref{subsec:var_ub} for details. 
\end{remark}

By Lemma \ref{lem_bmsb}, for any $1\leq k \leq \lfloor n/2 \rfloor$, the matrix  $\underline{\Gamma}_{R}$ in Theorem \ref{thm1} can be specified as
\begin{equation}\label{Gamma_under_R}
\underline{\Gamma}_{R}= \sigma^2 R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \Gamma_k) R/(4C_0)^2.
\end{equation}
By Lemmas \ref{lem_gram1} and \ref{lem_gram}, the matrix $\overline{\Gamma}_R$ in Theorem \ref{thm1} can be chosen as
$\overline{\Gamma}_R=\overline{\Gamma}_{R}^{(1)}$ or $\overline{\Gamma}_{R}^{(2)}$, where
\begin{equation}\label{Gamma_over_R}
\overline{\Gamma}_{R}^{(1)}=\sigma^2 m R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \Gamma_n) R/\delta, \quad 
\overline{\Gamma}_{R}^{(2)}=\sigma^2 R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \Gamma_n) R + \sigma^2 \xi(m,d,n,\delta) R^{ \mathrm{\scriptscriptstyle T} } R;
\end{equation}
recall the definitions of $\underline{\Gamma}_{R}$ and $\overline{\Gamma}_R$  in \eqref{notation:Gammaline_R}, where $q=d$ for the \textsc{var}($1$) model.
Furthermore, observe that $\underline{\Gamma}_{R}$ in \eqref{Gamma_under_R} and  the two $\overline{\Gamma}_R$'s in \eqref{Gamma_over_R}, which serve as lower and upper bounds of $Z^{ \mathrm{\scriptscriptstyle T} } Z$, respectively, are all related to the controllability Gramian $\Gamma_{t}$, and $0\prec I_d \preceq\Gamma_{k}\preceq\Gamma_{n}$. 

\subsection{Feasible region of $k$\label{subsec:k_condition}}

With $\underline{\Gamma}_{R}$ and $\overline{\Gamma}_R$ chosen as in \eqref{Gamma_under_R} and \eqref{Gamma_over_R}, the term $\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})$ in condition \eqref{thm1_condition} in Theorem \ref{thm1} is intricately  dependent on both $k$ and $n$. Thus, in order to apply the theorem to model \eqref{var1_true}, we need to verify the existence of the block size $k$ satisfying  \eqref{thm1_condition}. This boils down to deriving an explicit upper bound of $\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})$ free of $k$. 

By \eqref{eq:Gramian} and \eqref{Gamma_under_R}, it is easy to show that  $\det(\underline{\Gamma}_{R})$ is monotonic increasing in $k$. Then, as $\overline{\Gamma}_R$ is free of $k$, $\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})$ is maximized when $k=1$. As a result, we can first upper bound $\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})$ by its value at $k=1$ to get rid of its dependence on $k$. That is,
\begin{equation}\label{logdet1} 
\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1}) \leq \log \det\{\overline{\Gamma}_R (\sigma^2 R^{ \mathrm{\scriptscriptstyle T} } R)^{-1} (4C_0)^2\}.
\end{equation}

Now it suffices to upper bound the right-hand side of \eqref{logdet1}, where $\overline{\Gamma}_R$ can be chosen from $\overline{\Gamma}_{R}^{(1)}$ and $\overline{\Gamma}_{R}^{(2)}$ in \eqref{Gamma_over_R}.
As shown in the supplementary material,
\begin{equation}\label{logdet} 
\log \det\{\overline{\Gamma}_R (\sigma^2 R^{ \mathrm{\scriptscriptstyle T} } R)^{-1} (4C_0)^2\} \lesssim \begin{cases}
m\log(m/\delta)+\kappa, & \text{if}\quad \overline{\Gamma}_R=\overline{\Gamma}_{R}^{(1)}\\
m\log \{2\max(1,\xi)\}+\kappa, & \text{if}\quad \overline{\Gamma}_R=\overline{\Gamma}_{R}^{(2)}
\end{cases},
\end{equation}
where $\xi=\xi(m,d,n,\delta)$ is defined as in \eqref{xi}, and
\begin{equation}\label{kappa}
\kappa =\log \det\left \{R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \Gamma_n) R (R^{ \mathrm{\scriptscriptstyle T} } R)^{-1} \right \}.
\end{equation}

Obviously, without imposing normality, we can only choose $\overline{\Gamma}_{R}=\overline{\Gamma}_{R}^{(1)}$;  see Lemmas \ref{lem_gram1} and \ref{lem_gram} in $\S$ \ref{subsec:verify}. However, if $\{\eta_t\}$ are normal,  $\overline{\Gamma}_R$ can be set to whichever of $\overline{\Gamma}_{R}^{(1)}$ and $\overline{\Gamma}_{R}^{(2)}$  delivers the sharper upper bound. Note that both $\kappa$ and $\xi$ depend on $n$, and their growth rates with respect to $n$ depend on the magnitude of $A_*$.  This will ultimately affect the choice between $\overline{\Gamma}_{R}^{(1)}$ and $\overline{\Gamma}_{R}^{(2)}$; e.g., if $\kappa$ is the dominating term in both upper bounds in \eqref{logdet}, we will be indifferent between the two.  Assumptions \ref{assum:explosive}--\ref{assum:stable2} below summarize the three cases of $A_*$ we consider:

\begin{assumption}\label{assum:explosive}
	It holds $\rho(A_*)\leq 1+c/n$, where $c>0$ is a universal constant.
\end{assumption}

\begin{assumption}\label{assum:stable1}
	It holds $\rho(A_*)\leq \bar{\rho} <1$ and $\lVert A_*\rVert_2 \leq C$, where $\bar{\rho}, C >0$ are universal constants.
\end{assumption}


\renewcommand{\theassumption}{\arabic{assumption}$^\prime$}
\setcounter{assumption}{5}
\begin{assumption}\label{assum:stable2}
	It holds $\rho(A_*)\leq \bar{\rho} <1$, $\mu_{\min}(\mathcal{A})=\inf_{\lVert z \rVert=1} \lambda_{\min}(\mathcal{A}^*(z)\mathcal{A}(z))\geq \mu_1$, and  $\lVert A_*^t\rVert_2 \leq C \varrho^t$ for any integer $1\leq t\leq n$, where  $\bar{\rho}, \mu_1, C >0$ and $\varrho\in(0,1)$ are universal constants, and $\mathcal{A}(z)=I_d-A_*z$ for any complex number $z$.
\end{assumption}

Assumption  \ref{assum:explosive} is the most general case among the above three, and Assumption \ref{assum:stable1} is weaker than Assumption \ref{assum:stable2}.  Notably, Assumption \ref{assum:stable2} does not require $\lVert A_*\rVert_2 < 1$ because $C$ may be greater than one. \cite{Guo_Wang_Yao2016} assumed  $\lVert A_*^t\rVert_2 \leq \varrho^t$ for $\varrho\in(0,1)$ and any positive integer $t$, while it is unclear if this can be relaxed to $\lVert A_*^t\rVert_2 \leq C \varrho^t$ as in Assumption \ref{assum:stable2}. We need  $\mu_{\min}(\mathcal{A})$ to be bounded away from zero in order to derive a sharp upper bound on $\lVert \Sigma_X \rVert_2$; see Remark \ref{remark:xi} below.  This condition is also necessary for the estimation error rates derived in \cite{Basu_Michailidis2015} for a similar reason. In particular, if $A_*$ is diagonalizable, it is shown in Proposition 2.2 therein that $\mu_{\min}(\mathcal{A})\geq \left[ \{1-\rho(A_*)\}/\cond(S) \right]^2$, where $S$ is defined in \eqref{eq:jordan} below.

\begin{remark}\label{remark:kappa}
	From \eqref{kappa}, $\kappa$ is dependent on $n$ through $\Gamma_n$.  If $\rho(A_*)<1$,  it holds
	$\Gamma_{n}\preceq\Gamma_{\infty}=\lim_{n\rightarrow\infty}\Gamma_{n}<\infty$, and then
	$\kappa \leq \log \det\left \{R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \Gamma_\infty) R ( R^{ \mathrm{\scriptscriptstyle T} } R)^{-1} \right \}$,
	an upper bound free of $n$. By contrast, if $\rho(A_*)\geq 1$, $\Gamma_{\infty}$ no longer exists, and we need to carefully control the growth rate of $\Gamma_{n}$; see Lemma \ref{lem_Gamma_n} in the supplementary material. This is achieved via the Jordan decomposition of $A_*$ in \eqref{eq:jordan} below, and the mildest condition we need is Assumption \ref{assum:explosive}. The upper bound of $\kappa$ under Assumption \ref{assum:explosive} or \ref{assum:stable1} is given by Lemma \ref{lem_kappa} in the supplementary material.
\end{remark}

\begin{remark}\label{remark:xi}
	For $\xi$ defined in \eqref{xi}, it is worth noting that $\lVert \Sigma_X \rVert_2$ depends on $n$, as $\Sigma_X=[E(X_tX_s^{ \mathrm{\scriptscriptstyle T} })]_{1\leq t,s\leq n}$, where $E(X_tX_s^{ \mathrm{\scriptscriptstyle T} })=\sigma^2  A_*^{t-s} \Gamma_s$ for $1\leq s\leq t\leq n$ under Assumptions \ref{assum:var}(i) and (ii). Unlike $\kappa$ discussed in Remark \ref{remark:kappa}, even under Assumption \ref{assum:stable1},  $\lVert \Sigma_X \rVert_2$ is not guaranteed to be bounded by a constant free of $n$; indeed, we need Assumption \ref{assum:stable2} for this purpose, since $\lVert \Sigma_X \rVert_2$ is affected by not only the growing diagonal blocks $\sigma^2\Gamma_1, \dots, \sigma^2\Gamma_n$ but also the growing off-diagonal blocks; see Lemma \ref{lem_Sigma_X} in the supplementary material for details.  The upper bound of $\xi$ under Assumption \ref{assum:explosive} or \ref{assum:stable2} is given by  Lemma \ref{lem_xi} in the supplementary material.
\end{remark}

Let the Jordan decomposition of $A_*$ be 
\begin{equation}\label{eq:jordan}
A_*=SJS^{-1},
\end{equation}
where $J$ has $L$ blocks with sizes $1\leq b_1,\ldots, b_L \leq d$, and both $J$ and $S$ are $d\times d$ complex matrices. Let 
$b_{\max}=\max_{1\leq \ell \leq L} b_\ell,$
and denote the condition number of $S$ by
$\cond(S)=\left \{\lambda_{\max}(S^* S)/\lambda_{\min}(S^*S)\right \}^{1/2}$,
where $S^*$ is the conjugate transpose of $S$.

The following proposition, which follows from \eqref{logdet1}, \eqref{logdet} and upper bounds of $\kappa$ and $\xi$ under Assumptions \ref{assum:explosive}, \ref{assum:stable1} or \ref{assum:stable2}, is proved in $\S$ \ref{sec:S3} of the supplementary material.

\begin{proposition}\label{prop_G}
	For any $A_*\in\mathbb{R}^{d\times d}$, under Assumption \ref{assum:explosive}, we have $$\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1}) \lesssim m \left[\log \{d\cond(S)/\delta\}+ b_{\max} \log n\right]$$ for $\overline{\Gamma}_R=\overline{\Gamma}_R^{(1)}$ or $\overline{\Gamma}_R^{(2)}$. Moreover, if Assumption \ref{assum:stable1} holds, then  $\log \det(\overline{\Gamma}_R^{(1)}\underline{\Gamma}_{R}^{-1}) \lesssim m \log (m/\delta)$. Furthermore,  if Assumption \ref{assum:stable2} holds and  $n  \gtrsim m +\log(d/\delta)$, then  $\log \det(\overline{\Gamma}_R^{(2)}\underline{\Gamma}_{R}^{-1}) \lesssim m$. 
\end{proposition}


The condition $n  \gtrsim m +\log(d/\delta)$ in Proposition \ref{prop_G} is not stringent, because it is necessary for condition \eqref{thm1_condition} in Theorem \ref{thm1}, where $q=d$ for the $\textsc{var}(1)$ model. By Proposition \ref{prop_G},  $\overline{\Gamma}_{R}^{(2)}$ yields a sharper upper bound of $\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})$ than does $\overline{\Gamma}_{R}^{(1)}$  only under Assumption \ref{assum:stable2}. Thus, we shall always set $\overline{\Gamma}_{R}=\overline{\Gamma}_{R}^{(1)}$ unless Assumption \ref{assum:stable2} holds and  $\{\eta_t\}$ are normal.
As a result, under Assumption \ref{assum:var}, the feasible region of $k$ that is sufficient for condition \eqref{thm1_condition} in Theorem \ref{thm1} is
\begin{equation}\label{eq:k_sufficient}
k \lesssim \begin{cases}
\frac{n}{m\left [\log\{d\cond(S)/\delta\}+ b_{\max} \log n \right ]}, & \text{if Assumption \ref{assum:explosive} holds}\\
\frac{n}{m\log(m/\delta)+\log d}, & \text{if Assumption \ref{assum:stable1} holds}\\
\frac{n}{m+\log(d/\delta)}, & \text{if Assumption \ref{assum:stable2} holds and $\{\eta_t\}$ are normal}
\end{cases}.
\end{equation}

Because the upper bound of $\log \det(\overline{\Gamma}_R\underline{\Gamma}_{R}^{-1})$ and the feasible region of $k$ are both dependent on the assumption on $A_*$ and whether  $\{\eta_t\}$ are normal, the resulting estimation error bounds will vary slightly  under different conditions; see Theorems \ref{thm2} and \ref{thm3} in $\S$ \ref{subsec:var_ub}. 


\subsection{Analysis of upper bounds in vector autoregression\label{subsec:var_ub}}

We focus on the upper bound analysis of $\lVert \widehat{\beta} - \beta_* \rVert$; nevertheless, from Proposition \ref{prop1} we can readily obtain analogous results for $\lVert \widehat{A} - A_* \rVert_2$,  which are omitted here.  For simplicity, denote 
\begin{equation*
\Gamma_{R,k}=R \left \{R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \Gamma_k) R\right \}^{-1} R^{ \mathrm{\scriptscriptstyle T} }.
\end{equation*}

The first theorem follows directly from Theorem \ref{thm1}, Lemmas \ref{lem_bmsb}--\ref{lem_gram} and Proposition \ref{prop_G}.
\begin{theorem}\label{thm2}
	Let $\{X_t\}_{t=1}^{n+1}$ be generated by the linearly restricted vector autoregressive model. Fix $\delta\in(0,1)$. For any $1\leq k \leq \lfloor n/2 \rfloor$ satisfying \eqref{eq:k_sufficient}, under Assumption \ref{assum:var}, we have the following results: (i) If Assumption \ref{assum:explosive} holds, with probability at least $1-3\delta$,
	\begin{equation*
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \left (\lambda_{\max}(\Gamma_{R,k})\frac{ m \left[\log \{d\cond(S)/\delta\}+ b_{\max} \log n\right]}{n} \right )^{1/2}.
	\end{equation*}
	(ii) If Assumption \ref{assum:stable1} holds, with probability at least $1-3\delta$,
	\[
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \left \{\lambda_{\max}(\Gamma_{R,k})\frac{m \log(m/\delta)}{n}\right \}^{1/2}.
	\]
	(iii) If Assumption \ref{assum:stable2} holds and $\{\eta_t\}$ are normal,  with probability at least $1-3\delta$,
	\[
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \left \{\lambda_{\max}(\Gamma_{R,k})\frac{m +\log(1/\delta)}{n}\right \}^{1/2}.
	\]		
\end{theorem}

To gain an intuitive understanding of the factor $\lambda_{\max}(\Gamma_{R,k})$ in Theorem \ref{thm2}, consider the asymptotic distribution of $\widehat{\beta}$ under the assumptions that $\rho(A_*)<1$  and that $d, m$ and $A_*$ are all fixed: 
\begin{equation}\label{eq:limit_dist}
n^{1/2}(\widehat{\beta} - \beta_*)\rightarrow N\big(0, \underbrace{R \{R^{ \mathrm{\scriptscriptstyle T} } (I_d\otimes \Gamma_\infty) R\}^{-1} R^{ \mathrm{\scriptscriptstyle T} }}_{\Gamma_{R,\infty}}\big)
\end{equation}
in distribution as $n\rightarrow\infty$, where $\Gamma_{\infty}=\lim_{k\rightarrow\infty}\Gamma_{k}$; see  \cite{Lutkepohl2005}. Thus, $\lambda_{\max}(\Gamma_{R,k})$ in Theorem \ref{thm2} resembles the limiting covariance matrix $\Gamma_{R,\infty}$. However, it is noteworthy that  by adopting a non-asymptotic approach, Theorem \ref{thm2} retains  the dependence of the estimation error on $\Gamma_{R,k}$ across stable, unstable and slightly explosive regimes.  Also note that, similarly to \eqref{eq:limit_dist}, the error bounds in Theorem \ref{thm2} are free of $\sigma^2$, as the scaling effect of $\sigma^2$ on $\eta_t$ is canceled out by that on $X_t$ due to the autoregressive structure. 

As a special case, $b_{\max}=1$ if $A_*$ is diagonalizable. Moreover, if $A_*=\rho I_d$, then $b_{\max}=1$, $\cond(S)=1$, $\Gamma_k=\gamma_k(\rho) I_d$, with $\gamma_k(\rho)=\sum_{s=0}^{k-1}\rho^{2s}$, and thus
\begin{equation}\label{eq:lambda_max_equal}
\lambda_{\max}(\Gamma_{R,k})=\gamma_k^{-1}(\rho) \lambda_{\max}\{R (R^{ \mathrm{\scriptscriptstyle T} } R)^{-1} R^{ \mathrm{\scriptscriptstyle T} }\}=\gamma_k^{-1}(\rho) \lambda_{\max}\{(R^{ \mathrm{\scriptscriptstyle T} } R)^{-1} R^{ \mathrm{\scriptscriptstyle T} } R\}=\gamma_k^{-1}(\rho),
\end{equation}
where the second equality is due to the fact that, for any matrices $A\in\mathbb{R}^{N\times m}$ and $B\in\mathbb{R}^{m\times N}$, $AB$ and $BA$ have the same nonzero eigenvalues \citep[Theorem 1.3.20,][]{Horn_Johnson1985}. 

By Theorem \ref{thm2}, The linear restrictions affect the error bounds through both the factor $\lambda_{\max}(\Gamma_{R,k})$ and the explicit rate function of $m$ and $n$. To further illustrate this, suppose that 
\[\beta_*=R\theta_*+\gamma=R^{(1)}R^{(2)}\theta_*+\gamma,\] 
where $R^{(1)}\in\mathbb{R}^{d^2\times\widetilde{m}}$ has rank $\widetilde{m}$, and $R^{(2)}\in\mathbb{R}^{\widetilde{m}\times m}$ has rank $m$, with $\widetilde{m}\geq m+1$. Then $\mathcal{L}^{(1)}=\{R^{(1)}\theta+\gamma:\theta\in\mathbb{R}^{\widetilde{m}}\}\supseteq\mathcal{L} =\{R\theta+\gamma: \theta\in\mathbb{R}^m\}$. By an argument similar to that in \citet[p.~199]{Lutkepohl2005}, we can show that $\Gamma_{R,k}\preceq \Gamma_{R^{(1)},k}$, so
\begin{equation}\label{eq:lambda_max}
\lambda_{\max}(\Gamma_{R,k}) \leq \lambda_{\max}(\Gamma_{R^{(1)},k}). 
\end{equation}
Note that the parameter space $\mathcal{L}^{(1)}$ has fewer restrictions than $\mathcal{L}$. Therefore, with fewer restrictions,  the effective model size will increase from $m$ to $\widetilde{m}$, and meanwhile $\lambda_{\max}(\Gamma_{R,k})$ will increase to $\lambda_{\max}(\Gamma_{R^{(1)},k})$, both leading to deterioration of the error bound. 

\begin{remark}
	The preservation of the factor $\lambda_{\max}(\Gamma_{R,k})$ in Theorem \ref{thm2} is achieved by bounding $Z^{ \mathrm{\scriptscriptstyle T} } Z$ and $Z^{ \mathrm{\scriptscriptstyle T} } \eta$ simultaneously through the Moore-Penrose pseudoinverse $Z^\dagger$, and note that  $Z^\dagger=(Z^{ \mathrm{\scriptscriptstyle T} } Z)^{-1}Z^{ \mathrm{\scriptscriptstyle T} }$ if $Z^{ \mathrm{\scriptscriptstyle T} } Z\succ 0$; see
	also \cite{Simchowitz2018}.  This key advantage is not enjoyed by the non-asymptotic analyses in \cite{Basu_Michailidis2015} and \cite{Faradonbeh_Tewari_Michailidis2018}.  In their analyses,  $X^{ \mathrm{\scriptscriptstyle T} } X$ and $X^{ \mathrm{\scriptscriptstyle T} } E$, or  $Z^{ \mathrm{\scriptscriptstyle T} } Z$ and $Z^{ \mathrm{\scriptscriptstyle T} } \eta$ in our context, were bounded separately. This would not only break down $\Gamma_{R,k}$, but also cause degradation of the error bound as $\rho(A_*)\rightarrow1$ due to the inevitable involvement of the condition number of $X^{ \mathrm{\scriptscriptstyle T} } X$ in the resulting error bound. 
\end{remark}


\begin{remark}
	If $A_*=\rho I_d$, then \eqref{eq:lambda_max} becomes an equality. However, the equality generally does not hold even for diagonal matrices $A_*$. For example, if $A_*=\diag(\rho_1,\rho_2)\in\mathbb{R}^{2\times 2}$, where $\vert\rho_1\vert>\vert\rho_2\vert$, then $\Gamma_k=\diag\{\gamma_k(\rho_1), \gamma_k(\rho_2)\}$. Let $R=(1,0,0,0)^{ \mathrm{\scriptscriptstyle T} }=R^{(1)}R^{(2)}$, where $R^{(1)}=(I_2,0)^{ \mathrm{\scriptscriptstyle T} }\in\mathbb{R}^{4\times 2}$ and $R^{(2)}=(1,0)^{ \mathrm{\scriptscriptstyle T} }$. Consequently, \eqref{eq:lambda_max} is a strict inequality:
	$\lambda_{\max}(\Gamma_{R,k}) = \gamma_k^{-1}(\rho_1) < \gamma_k^{-1}(\rho_2)=\lambda_{\max}(\Gamma_{R^{(1)},k})$.
\end{remark}


The next theorem sharpens the error bounds in Theorem \ref{thm2} by utilizing the largest possible $k$, since $\lambda_{\max}(\Gamma_{R,k})$ is monotonic decreasing in $k$. The dependence of $\lambda_{\max}(\Gamma_{R,k})$ on $A_*$ will be captured by $\sigma_{\min}(A_*)$, a measure of the least excitable mode of the underlying dynamics.


\begin{theorem}\label{thm3}
	Suppose that the conditions of Theorem \ref{thm2} hold. Fix $\delta\in(0,1)$, and let $c_1>0$ be a universal constant.
	
	(i) Under Assumption \ref{assum:explosive}, if
	\begin{equation}\label{eq:sigmin_slow1}
	\sigma_{\min}(A_*)\leq 1-\frac{c_1m\left [\log\{d\cond(S)/\delta\}+ b_{\max} \log n \right ]}{n},\tag{A.1}
	\end{equation}	
	then, with probability at least $1-3\delta$, 
	\begin{equation}\label{eq:slow1}
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \left (\frac{\{1-\sigma_{\min}^2(A_*)\}  m \left[\log \{d\cond(S)/\delta\}+ b_{\max} \log n\right]  }{n}\right )^{1/2}, \tag{S.1}
	\end{equation}	
	and if inequality \eqref{eq:sigmin_slow1} holds in the reverse direction, then, with probability at least $1-3\delta$, 
	\begin{equation}\label{eq:fast1}
	\lVert \widehat{\beta} - \beta_* \rVert  \lesssim  \frac{m \left[\log \{d\cond(S)/\delta\}+ b_{\max} \log n\right]}{n}. \tag{F.1}
	\end{equation}
	
	(ii) Under Assumption \ref{assum:stable1}, if
	\begin{equation}\label{eq:sigmin_slow2}
	\sigma_{\min}(A_*)\leq 1-\frac{c_1\{m\log(m/\delta)+\log d\}}{n},\tag{A.2}
	\end{equation}
	then, with probability at least $1-3\delta$, 
	\begin{equation}\label{eq:slow2}
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim  \left [\frac{\{1-\sigma_{\min}^2(A_*)\}m \log(m/\delta)}{n}\right ]^{1/2}, \tag{S.2}
	\end{equation}
	and if inequality \eqref{eq:sigmin_slow2} holds in the reverse direction, then, with probability at least $1-3\delta$, 
	\begin{equation}\label{eq:fast2}
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \frac{m \log(m/\delta)+\log d}{n}. \tag{F.2}
	\end{equation}	
	
	(iii) Under Assumption \ref{assum:stable2}, if $\{\eta_t\}$ are normal and
	\begin{equation}\label{eq:sigmin_slow3}
	\sigma_{\min}(A_*)\leq 1-\frac{c_1\{m+\log(d/\delta)\}}{n},\tag{A.3}
	\end{equation}
	then, with probability at least $1-3\delta$, 
	\begin{equation}\label{eq:slow3}
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \left [\frac{ \{1-\sigma_{\min}^2(A_*)\} \{m +\log(1/\delta)\}}{n}\right ]^{1/2}. \tag{S.3}
	\end{equation}	
	and if inequality \eqref{eq:sigmin_slow3} holds in the reverse direction, then, with probability at least $1-3\delta$, 
	\begin{equation}\label{eq:fast3}
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \frac{m +\log(d/\delta)}{n}. \tag{F.3}
	\end{equation}			
\end{theorem}

Theorem \ref{thm3} reveals an interesting phenomenon of phase transition  from slow to fast error rate regimes, i.e., from about $O\{(m/n)^{1/2}\}$ as in \eqref{eq:slow1}--\eqref{eq:slow3} to about $O(m/n)$ as in \eqref{eq:fast1}--\eqref{eq:fast3}, up to logarithmic factors.  Within the slow rate regime, the estimation error decreases as  $\sigma_{\min}(A_*)$ increases. Note that the slow rates in \eqref{eq:slow1}--\eqref{eq:slow3} differ from each other only by logarithmic factors, and so do the fast rates in \eqref{eq:fast1}--\eqref{eq:fast3}.  Moreover, the point at which the transition occurs is dependent on $\sigma_{\min}(A_*)$ instead of $\rho(A_*)$; see conditions \eqref{eq:sigmin_slow1}--\eqref{eq:sigmin_slow3}. Since $\sigma_{\min}(A_*)\leq \rho(A_*)$, conditions \eqref{eq:sigmin_slow1}--\eqref{eq:sigmin_slow3} may be mild as long as $\rho(A_*)$ is not too large. However, the fast rates require the opposite of \eqref{eq:sigmin_slow1}--\eqref{eq:sigmin_slow3},  which cannot be directly inferred from $\rho(A_*)$.


\begin{remark}\label{remark:phases}
	For the special case of $A_*=\rho I_d$ with $\rho\in\mathbb{R}$, Assumptions \ref{assum:stable1} and \ref{assum:stable2} both simply reduce to $\vert\rho\vert<1$, and Assumption \ref{assum:explosive} to $\vert\rho\vert\leq 1+O(1/n)$. Also, $\cond(S)=1$ and $b_{\max}=1$. Thus, under Assumption \ref{assum:var}, by \eqref{eq:fast1}, \eqref{eq:slow2} and \eqref{eq:fast2}, the following holds with high probability:
	\[
	\lVert \widehat{\beta} - \beta_* \rVert \lesssim \begin{cases}
	O[\{(1-\rho^2)m\log m/n\}^{1/2}], & \text{if} \quad \vert\rho\vert\leq 1-O\{(m\log m+\log d)/n\}\\
	O\{(m\log m+\log d)/n\},  & \text{if} \quad 1-O\{(m\log m+\log d)/n\} \leq \vert\rho\vert <1\\
	O\{m\log(d n)/n\}, & \text{if} \quad 1\leq \vert\rho\vert \leq 1+O(1/n)
	\end{cases}.
	\]	
	Moreover, if $\{\eta_t\}$ are normal, by \eqref{eq:slow3} and \eqref{eq:fast3}, we can eliminate all  factors of $\log m$ in the above results; that is, every $m\log m$ will be replaced by $m$.
\end{remark}

\begin{remark}\label{remark:Michailidis}
	\cite{Faradonbeh_Tewari_Michailidis2018} derived error bounds for the ordinary least squares estimator of explosive unrestricted vector autoregressive processes when (1)  $|\lambda_{\min}(A_*)|>1$ or (2)  $A_*$ has no unit eigenvalue. In contrast to case (1), we focus on slightly explosive processes with $\rho(A_*)\leq 1+O(1/n)$.  Moreover, the no unit root requirement of case (2) may be quite restrictive; e.g., it excludes the case of $|\rho|=1$ in Remark \ref{remark:phases}. Thus, the conditions in Theorem \ref{thm3} may be more reasonable in practice.
\end{remark}

\begin{remark}\label{remark:Guo}
	\cite{Guo_Wang_Yao2016} obtained the error rate $O_p\{(m/n)^{1/2}\}$ for the banded vector autoregressive model under weaker conditions on $\{\eta_t\}$ yet stronger conditions on $A_*$ than those in this paper; see Theorem 2 therein. Particularly, they required $\lVert A_*\rVert_2<1$. Note that this rate matches \eqref{eq:slow3}. In view of the lower bounds to be presented in $\S$\ref{sec:lowerbd}, we  conjecture that the rate in \eqref{eq:slow2} is larger than the actual rate by a factor of  $\log m$. 	
\end{remark}

\section{Analysis of lower bounds \label{sec:lowerbd}}

For any $\theta\in\mathbb{R}^m$, let $\beta=R\theta+\gamma$, and the corresponding transition matrix is denoted by $A(\theta)=(I_d\otimes \theta^{ \mathrm{\scriptscriptstyle T} })\widetilde{R}+G$, where $R$, $\gamma$, $\widetilde{R}$ and $G$ are defined as in $\S$ \ref{subsec:var_example}. As $\beta$ is completely determined by $\theta$, it is more convenient to index the probability law of the  model by the unrestricted parameter $\theta$. Thus, we denote by ${\mathrm{pr}}_{\theta}^{(n)}$ the distribution of the sample $(X_1,\ldots, X_{n+1})$ on  $(\mathcal{X}^{n+1}, \mathcal{F}_{n+1})$, where $\mathcal{X}=\mathbb{R}^{d}$ and $\mathcal{F}_{n+1}=\sigma\{\eta_1,\ldots,\eta_{n}\}$.	For any fixed $\bar{\rho}>0$, we write the subspace of $\theta$ such that the spectral radius of $A(\theta)$ is bounded above by $\bar{\rho}$ as
\[
\Theta(\bar{\rho}) = \{\theta\in \mathbb{R}^m:  \rho \{A(\theta)\} \leq \bar{\rho} \}.
\]
The corresponding linearly restricted subspace of $\beta$ is denoted by
$\mathcal{L}(\bar{\rho})=\{R\theta+\gamma: \theta\in\Theta(\bar{\rho})\}$.

The minimax rate of estimation over $\beta\in \mathcal{L}(\bar{\rho})$, or $\theta\in\Theta(\bar{\rho})$, is provided by the next theorem.


\begin{theorem}\label{thm_lower}
	Suppose that $\{X_t\}_{t=1}^{n+1}$ follow the vector autoregressive model $X_{t+1}=AX_t+\eta_t$ with linear restrictions defined as in $\S$ \ref{sec:VAR_upperbd}. In addition, Assumptions \ref{assum:var}(i) and (ii) hold, and  $\{\eta_t\}$ are normally distributed. Fix  $\delta\in(0,1/4)$ and $\bar{\rho}>0$.  Let $\gamma_n(\bar{\rho})=\sum_{s=0}^{n-1}\bar{\rho}^{2s}$. Then, for any $\epsilon\in(0, \bar{\rho}/4]$, we have
	\[
	\inf_{\widehat{\beta}}\sup_{\theta\in \Theta(\bar{\rho})} {\mathrm{pr}}_{\theta}^{(n)} \left \{ \lVert \widehat{\beta}-\beta\rVert \geq \epsilon \right \} \geq \delta,
	\]
	where the infimum is taken over all estimators of $\beta$ subject to $\beta\in\{R\theta+\gamma: \theta\in\mathbb{R}^m\}$,
	for any $n$ such that 
	\begin{equation*
	n \gamma_n(\bar{\rho}) \lesssim \frac{m+\log(1/\delta)}{\epsilon^2}.
	\end{equation*}
\end{theorem}

As a result, we have the following minimax rates of estimation across different values of $\bar{\rho}$.

\begin{corollary}\label{cor_lower}
	For the linearly restricted vector autoregressive model in Theorem \ref{thm_lower}, the minimax rates of estimation over $\beta\in \mathcal{L}(\bar{\rho})$  are given as follows:
	
	(i) $\{{(1-\bar{\rho}^2)m}/{n}\}^{1/2}$,  if $0<\bar{\rho}< (1-1/n)^{1/2}$; 
	
	(ii) $m^{1/2}/n$, if $(1-1/n)^{1/2} \leq \bar{\rho} \leq 1+c/n$ for a fixed $c>0$; and
	
	(iii) $\bar{\rho}^{-n}\{(\bar{\rho}^2-1)m/n\}^{1/2}$, if $\bar{\rho}>1+c/n$.
\end{corollary}

While the phase transition in Corollary \ref{cor_lower} depends on $\rho(A_*)$ instead of  $\sigma_{\min}(A_*)$, we may still compare the lower bounds to the upper bounds in Theorem \ref{thm3}.  For case (i)  in Corollary \ref{cor_lower}, since $\sigma_{\min}(A_*)\leq \rho(A_*)<(1-1/n)^{1/2}<1$, we may expect that condition \eqref{eq:sigmin_slow2} will  hold in most cases, and hence the upper bound $O\{(m \log m/n)^{1/2}\}$ differs from the lower bound by a factor of $\log m$. However, for case (ii), the upper bound may lie in either the slow or fast rate regime, depending on the magnitude of $\sigma_{\min}(A_*)$, whereas the lower bound lies in the fast rate regime. As shown by our first experiment in $\S$ \ref{sec:simul}, the transition from slow to fast error rates actually depends on $\sigma_{\min}(A_*)$ instead of $\rho(A_*)$. This also suggests that the results in Theorem \ref{thm3} are sharp in the sense that they correctly capture the transition behavior.


\begin{remark}\label{remark:phases_all}
	If $A_*=\rho I_d$ with $\rho\in\mathbb{R}$ and $\{\eta_t\}$ are normal, in view of Remark \ref{remark:phases} and Corollary \ref{cor_lower},  a more straightforward comparison of the upper and lower bounds can be made:
	\begin{center}
		\begin{tabular}{l|@{\hskip 1mm}c@{\hskip 1mm}c}
			\multicolumn{1}{c|@{\hskip 1mm}}{Range of $\lvert \rho \rvert$}& Lower bound                            & Upper bound\\\hline
			$(0,1-O\left \{(m+\log d)/n\right \}]$     & $\Omega[\{{(1-\rho^2)m}/{n}\}^{1/2}]$   & $O[\{(1-\rho^2)m/n\}^{1/2}]$\\
			$[1-O\{(m+\log d)/n\},(1-1/n)^{1/2})$       & $\Omega[\{{(1-\rho^2)m}/{n}\}^{1/2}]$   & $O\{(m+\log d)/n\}$\\
			$[(1-1/n)^{1/2},1)$                               & $\Omega(m^{1/2}/n)$               & $O\{(m+\log d)/n\}$\\
			$[1, 1+O(1/n)]$                                  & $\Omega(m^{1/2}/n)$               & $O\{m\log(dn)/n\}$\\
			$(1+O(1/n), \infty)$                             & $\Omega[\lvert\rho\rvert^{-n}\{(\rho^2-1)m/n\}^{1/2}]$& $-$
		\end{tabular}
	\end{center}
	See Figure \ref{fig5} for an illustration of the theoretical bounds and actual rates suggested by simulation results in $\S$ \ref{sec:simul}. Note that the actual rates and the theoretical upper and lower bounds exactly match when $0<\vert \rho \vert\leq 1-O\{(m+\log d)/n\}$. In addition, the suggested actual rate is  $m/n$ for $1-O\{(m+\log d)/n\}< \vert \rho \vert< 1+O(1/n)$ and  even faster for $\vert \rho \vert$ beyond this range. 
\end{remark}

\begin{figure}[h]
	\centering
	\includegraphics[width=\textwidth]{Rplot05_new.eps}
	\caption{Illustration of theoretical upper (dot-dash) and lower (dashes) bounds and the actual rates (solid) suggested by simulation results in $\S$ \ref{sec:simul} for the $\textsc{var}(1)$ model with $A_*=\rho I_d$ and normal innovations.}
	\label{fig5}
\end{figure}

\begin{remark}\label{remark:Han}
	\cite{Han_Xu_Liu2015} considered the estimation  of a class of copula-based stationary vector autoregressive processes which includes the Gaussian \textsc{var}($1$) process as a special case, and extended the theoretical properties to $\textsc{var}(p)$ processes by arguments similar to those in Example \ref{ex:VARp}. Under a strong-mixing condition on the process and  the low-rank assumption $\rank(A_*)\leq r$, the proposed estimator $\widetilde{A}$ was proved to attain the minimax error rate $\lVert \widetilde{A}-A_* \rVert_F=O\{(dr/n)^{1/2}\}$.  While we consider different restrictions and estimation method, the  rate $\{{(1-\bar{\rho}^2)m}/{n}\}^{1/2}$ in Corollary \ref{cor_lower} for the stable regime  resembles that in \cite{Han_Xu_Liu2015} if we regard $dr$  as the effective model size $m$ of the low-rank $\textsc{var}(1)$ model. However, similarly to our upper analysis in $\S$ \ref{sec:VAR_upperbd}, the factor of $(1-\bar{\rho}^2)^{1/2}$ in our lower bound also reveals that the estimation error may decrease as $A_*$ approaches the stability boundary.
\end{remark}

\section{Simulation experiments\label{sec:simul}}
We conduct three simulation experiments to verify the theoretical results in previous sections, including (i) the estimation error rates, (ii) the transition from slow to fast  rates regimes, and (iii) the effect of the ambient dimension $d$ on the estimation. The data are generated from \textsc{var}$(1)$ processes with $\{\eta_t\}$ drawn independently from $N(0, I_d)$ and the following structures of $A_*$:

\DGP1: Banded structure defined by the  zero restrictions $a_{*ij}=0$ if $\vert i-j\vert > k_0$, where $1\leq k_0 \leq \lfloor (d-1)/2 \rfloor$ is the bandwidth parameter; see Example \ref{ex:bandedAR} in $\S$ \ref{subsec:var_example}. As a result, if all restrictions are imposed, the size of the model is $m=d + (2d - 1)k_0 - k_0^2$.

\DGP2: Group structure defined by equality restrictions as follows. Partition the index set $\mathcal{V}=\{1,\ldots, d\}$ of the coordinates of $X_t$  into $K$ groups of size $b=d/K$ as $\mathcal{V}=\bigcup_{k=1}^K \mathcal{G}_k$, where
\[\mathcal{G}_k=\{(k-1)b+1,\ldots, kb\}\quad (k=1,\ldots, K).\] 
In each row of $A_*$, the off-diagonal entries $a_{*ij}$ with $j$ belonging to the same group are assumed to be equal:  for any $1\leq k\leq K$ and $1\leq i\leq d$, all elements of $\{a_{*ij}, j\in\mathcal{G}_k, j\neq i\}$ are equal. Thus, $m=(K+1)d$, as there are $(K+1)$ free parameters in each row of $A_*$.	

\DGP3: $A_*=\rho I_d$, where $\rho \in\mathbb{R}$. Note that the smallest true model with size $m=1$ results from imposing zero restrictions on all off-diagonal entries of $A_*$ and equality restrictions on all diagonal entries; see Example \ref{ex:pure_unit_root} in $\S$ \ref{subsec:var_example}.  

Throughout the experiments, the $\ell_2$ estimation error $\lVert \widehat{\beta} - \beta_* \rVert$ is calculated by averaging over $1000$ replications. Except for \DGP3, nonzero entries of $A_*$ are generated independently from $U[-1,1]$ and then rescaled such that $\rho(A_*)$ is equal to a certain value. 

The first experiment aims to verify the error rates in Theorem \ref{thm3} and the implication of Theorem \ref{thm2} that the restrictions can reduce the estimation error through both the explicit rate $(m/n)^{1/2}$ and the decrease in the factor  $\lambda_{\max}(\Gamma_{R,k})$. Fixing $d=24$ and $\rho(A_*)=0.2,0.8$ or $1$, we generate data from \DGP1 with $k_0=1$, \DGP2 with $K=2$, and \DGP3. 
For \DGP1 and \DGP3, we fit banded vector autoregressive models with $k_0=1,5$ or $7$ such that $m=70, 156$ or $304$, respectively. For \DGP2, we fit the group-structured model with $K=2,8$ or $12$ such that $m=72, 120$ or $312$, respectively. Note that for \DGP1 and \DGP2,  $\sigma_{\min}(A_*)\leq 0.1$ even when the randomly generated matrix $A_*$ has $\rho(A_*)=1$. However, for \DGP3, it holds  $\sigma_{\min}(A_*)=\rho(A_*)$.  The $\ell_2$ estimation error $\lVert \widehat{\beta} - \beta_* \rVert$ is plotted against $(m/n)^{1/2}$ in Figure \ref{fig1}, where we consider $(m/n)^{1/2}\in\{0.15,0.35,0.55,0.75,0.95\}$.  	Our findings are summarized as follows:

\begin{figure}[t!]
	\centering
	\includegraphics[width=\linewidth]{Rplot01_new.eps}
	\caption{Plots of $\lVert \widehat{\beta} - \beta_* \rVert$ against $(m/n)^{1/2}$ for three data generating processes with  $\rho(A_*)=0.2$ (dashed lines, grey-filled symbols), $0.8$ (dotted lines, unfilled symbols), or $1$ (solid lines, black-filled symbols) and different $m$. \DGP1 and \DGP3 were fitted as banded vector autoregressive models with $m=70$ (squares), $156$ (circles), or $304$ (triangles), and \DGP2 was fitted as grouped vector autoregressive models with $m=72$ (squares), $120$ (circles), or $312$ (triangles).}
	\label{fig1}
\end{figure}

(i) When $\rho(A_*)=0.2$, for all data generating processes, the lines for different $m$ coincide completely with each other and scale perfectly linearly with $(m/n)^{1/2}$. This suggests that the actual error rate is $(m/n)^{1/2}$ when $\sigma_{\min}(A_*)$ lies in the slow rate regime of Theorem \ref{thm3}.

(ii) When $\rho(A_*)=0.8$ or 1, for \DGP1 and \DGP2, although $\lVert \widehat{\beta} - \beta_* \rVert$ is still proportional to $(m/n)^{1/2}$, the three lines for the same $\rho(A_*)$  but different $m$ do not coincide: fixing $\rho(A_*)$, the slope increases as $m$ increases, and the variation in slope is greater as $\rho(A_*)$ is larger.  Note that for \DGP1 and \DGP2, $\sigma_{\min}(A_*)$ is very small.
As  finding (i) suggests that the actual error rate is $(m/n)^{1/2}$ for small $\sigma_{\min}(A_*)$, this extra variation in slope  may be partially explained by the factor $\lambda_{\max}(\Gamma_{R,k})$ in the error bound in Theorem \ref{thm2} due to the effect of $R$. However, $\lVert \widehat{\beta} - \beta_* \rVert$ actually depends on the spectrum of $\Gamma_{R,k}$, and its largest eigenvalue merely serves as an upper bound. When $\rho(A_*)$ is smaller,  we will have more control over the spectrum of $A_*$ and hence that of $\Gamma_{R,k}$. This may explain why the variation in slope is smaller when  $\rho(A_*)$ is smaller.

(iii) For \DGP3 with $\rho(A_*)=0.2$ or $0.8$, the three lines corresponding to different $m$ still completely coincide with each other.  This can be explained by the fact that $\lambda_{\max}(\Gamma_{R,k})$ is independent of $R$ when $A_*=\rho I_d$; see \eqref{eq:lambda_max_equal}.

(iv) For \DGP3 with $\rho(A_*)=1$, in sharp contrast to all other cases, the  error rate  appears to be a quadratic function of $(m/n)^{1/2}$. This  matches the implication of  Theorem \ref{thm3} that when $\sigma_{\min}(A_*)=1$, the error rate falls into the fast rate regime.

(v) Fixing both $m$ and $n$, $\lVert \widehat{\beta} - \beta_* \rVert$ always decreases as $\rho(A_*)$ increases. Moreover, when $\sigma_{\min}(A_*)<1$, fixing $m$,  the lines become less steep as $\rho(A_*)$ increases. Note that $\sigma_{\min}(A_*)$ is larger when $\rho(A_*)$ is, due to our method of generating $A_*$. Thus, this finding can be explained by the factor $\{1-\sigma_{\min}^2(A_*)\}^{1/2}$ in the  error bound  for  the slow rate regime in Theorem \ref{thm3}.

\begin{figure}[t!]
	\centering
	\includegraphics[width=0.9\textwidth]{Rplot02_new.eps}
	\caption{Error rates for \DGP3 as $\rho$ is fixed or approaching one at different rates. Left panel: plot of $\lVert \widehat{\beta} - \beta_* \rVert$ against $\{(1-\rho^2)/n\}^{1/2}$ with $\rho=0.2$ (dashed lines, squares), $0.4$ (dotted lines, circles) or $0.6$ (solid lines, triangles), and $m=70$. Right panel: plot of $\lVert \widehat{\beta} - \beta_* \rVert$ against $m/n$  with $\rho=0.99$ (squares), $1-(m+\log d)/n$ (circles), $1+1/n$ (triangles) or $1.01$ (diamonds), and $m=1$ (solid lines, filled symbols) or $70$ (dotted lines, unfilled symbols). The case of $(m,\rho)=(70,1.01)$ is omitted as the process becomes so explosive that the computation is numerically infeasible.}
	\label{fig2a}
\end{figure}

\begin{figure}[t!]
	\centering
	\includegraphics[width=0.9\textwidth]{Rplot03_new.eps}
	\caption{Error rates for \DGP3 when $\rho=1$. Left panel: plot of $\lVert \widehat{\beta} - \beta_* \rVert$ multiplied by $n/8$ (squares), $n/(2\log n)$ (circles) or $n^{1/2}$ (triangles)  against $n$, fixing $m=70$. Right panel: plot of $\lVert \widehat{\beta} - \beta_* \rVert$ divided by $m$  (squares) or $8m^{1/2}$ (triangles) against $m$, fixing $n=400$.}
	\label{fig2b}
\end{figure}

In the second experiment, we focus on \DGP3 to further investigate the error rates and the phase transition. We set $d=24$ and $m=1, 70, 156$ or $304$, where $m=1$ results from fitting the smallest true model, and $m=70, 156$ or $304$ from fitting a banded model with $k_0=1,3$ or $7$, respectively. Figures \ref{fig2a} and \ref{fig2b} display the results, where we have the following findings: 

(i) The combination of results from the first experiment and the left panel of Figure \ref{fig2a} suggests that $\lVert \widehat{\beta} - \beta_* \rVert$ scales as $O[\{(1-\rho^2)m/n\}^{1/2}]$ when $\vert \rho \vert$ is fixed at a level well below one.

(ii) Figure \ref{fig2b} suggests that when $\lvert\rho\rvert=1$ the actual error rate is $m/n$. Specifically, when $m$ is fixed, the left panel shows that $n \lVert \widehat{\beta} - \beta_* \rVert$ becomes stable for $n$ sufficiently large, while $\lVert \widehat{\beta} - \beta_* \rVert$ multiplied by $n^{1/2}$ or $n/\log n$ appears to diminish as $n\rightarrow\infty$. On the other hand, when $n$ is fixed, the right panel shows that $\lVert \widehat{\beta} - \beta_* \rVert/m$ becomes stable for $m$ sufficiently large. 

(iii) The right panel of Figure \ref{fig2a} suggests that the regime of rate $m/n$ is reached as early as $\vert\rho\vert=1-O\{(m+\log d)/n\}$ and maintains even as the process becomes slightly explosive with  $\vert\rho\vert=1+O(1/n)$. This lends support to the boundaries of the fast rate regime suggested by Theorem \ref{thm3}; see Remarks  \ref{remark:phases} and \ref{remark:phases_all}. By contrast, when $\rho=0.99$, the rate appears to be $(m/n)^{1/2}$, similarly to our findings in the first experiment.  On the other hand, when $\rho$ is fixed at a level slightly above one, the rate becomes even faster than $m/n$. This matches the conclusion in Remark \ref{remark:phases_all} that the corresponding lower bound  diminishes at a rate faster than $\lvert\rho\rvert^n$ as $n$ increases. 


\begin{figure}[t!]
	\centering
	\includegraphics[width=0.9\textwidth]{Rplot04_new.eps}
	\caption{Plots of $\lVert \widehat{\beta} - \beta_* \rVert$ against $d$ for \DGP3  with $n=100$ (left panel) and $500$ (right panel), when $\rho=0.2$ (squares), $0.8$ (circles) or $1$ (triangles), and $m=1$ (solid lines, filled symbols) or $20$ (dotted lines, unfilled symbols).}
	\label{fig3}
\end{figure}

The third experiment aims to check if the ambient dimension $d$  directly affects  the estimation error. We generate data from \DGP3 with $\rho=0.2, 0.8$ or $1$, $n=100$ or $500$ and $d\in[25,500]$. For the estimation, we consider $m=1$ or $20$, where $m=1$ corresponds to the smallest true model, and $m=20$ corresponds to a model subject to (1) $a_{*11}=\cdots=a_{*dd}$ and (2) the restriction that all but $m-1$ of the off-diagonal entries of $A_*$ are zero.  To generate the pattern in (2), we sample the $m-1$ positions uniformly without replacement from all off-diagonal positions of $A_*$. 

Figure \ref{fig3} shows that the estimation error is constant in $d$ for almost all cases. This confirms that $d$ does not affect the estimation error when $\lvert \rho\rvert <1-O\left \{(m+\log d)/n\right \}$; see Remark \ref{remark:phases_all}. Moreover, the extra factor of $\log d$ in the theoretical upper bounds for the other regimes might not be necessary. 	For $(n, m, \rho)=(100, 20, 1)$, the estimation error seems less stable when $d\leq 75$; see the left panel of Figure \ref{fig3}. 	
This might be explained  by the indirect effect of $R$ on $\lambda_{\max}(\Gamma_{R,k})$ in Theorem \ref{thm2}. As $m$ is fixed at $20$, different $d$ corresponds to different $R$. The spectrum of $\Gamma_{R,k}$ may be more sensitive to $R$ when $d$ is smaller, and the resulting impact on the estimation error may be more pronounced when $n$ is smaller. However, as $d$ grows, the restrictions will  become relatively more sparse, so eventually the spectrum of $\Gamma_{R,k}$ will  be stable, and the impact of $d$ will be negligible.

\section{Discussion\label{sec:discuss}}
An interesting future direction is dimensionality reduction for vector autoregressive models with data-driven restrictions.
Such a procedure involves first suggesting possible restrictions based on subject knowledge and then selecting the true restrictions by a data-driven approach. Note that the lasso method \citep{Davis_Zang_Zheng2015, Basu_Michailidis2015} can be viewed as a procedure where zero restrictions are initially suggested for all entries of $A$, and then the true zeros are identified by penalized estimation. Adopting a more general point of view, the modeller can initially suggest the general linear restrictions  \eqref{eq:lr}  instead. This will enable a more flexible and data-driven integration of expert knowledge. On the other hand, if it is known that only zero and equality restrictions are true, yet the locations of the restrictions are unknown, we can select the true restrictions efficiently by  the  delete or merge regressors algorithm proposed by \cite{Maj-Kanska2015} based on the Bayesian information criterion. The consistency of this procedure can be easily extended to vector autoregressive models. 


\section*{Acknowledgement}
We thank the editor, associate editor and two referees for their invaluable comments, which have led to substantial improvements of our paper. Cheng's research was partially supported by the U.S. National Science Foundation and the Office of Naval Research, and he wishes to thank the Institute for Advanced Study at Princeton for its hospitality during his visit in Fall 2019.

\section*{Supplementary material}
\label{SM}
Supplementary material available at {\it Biometrika}\ online includes all technical proofs of this paper.


\vspace*{-10pt}


\bibliographystyle{biometrika}
