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.
84,709 characters
Honest Inference for High-Dimensional Linear Models
\maketitle
\begin{abstract}
\noindent We propose a heteroskedasticity-robust and dimention-adaptive test for multiple linear restrictions in sparse high-dimensional linear regression models. It combines desparsified LASSO estimation with a ridge-regularized quadratic form statistic and critical values that adapt to the restriction dimension. Under stated regularity conditions, our approach is asymptotically valid under a fixed or diverging number of restrictions without specifying the asymptotic regime in advance. Our random-design theory allows independent non-identically distributed covariate rows with marginal sub-Weibull tails and heteroskedastic errors satisfying finite polynomial moments. Monte Carlo simulations demonstrate improved finite-sample size control relative to conventional Wald tests, particularly as the number of restrictions increases. \\
\noindent\textbf{Keywords:}
High-Dimensional Inference, Multiple Linear Restrictions, Desparsified LASSO, Strong Approximation, Heteroskedasticity, Heavy-Tailed Distributions. \\
\noindent\textbf{JEL Classification:} C12, C13, C55.
\end{abstract}
\setcounter{oldtocdepth}{\value{tocdepth}}
\addtocontents{toc}{\protect\setcounter{tocdepth}{-10}}
\section{Introduction}
We consider the linear model
\begin{align}
Y = {\bf X}\beta^0 + \varepsilon,
\label{eq:model}
\end{align}
where $Y \in \mathbb{R}^n$ is the outcome vector,
${\bf X} \in \mathbb{R}^{n \times p}$ is the \textit{design matrix},
which may be either fixed or random, and $\varepsilon \in \mathbb{R}^n$ is
the error vector. We assume that, conditional on ${\bf X}$, the individual
errors $\{\varepsilon_i\}_{i\in[n]}$ are independent but not necessarily
identically distributed, thereby allowing for heteroskedasticity. In the
random-design setting, we also allow the rows
$\{X_i\}_{i\in[n]}$ of ${\bf X}$ to be non-identically distributed.
We are concerned with the high-dimensional setting in which $p$ may be
much larger than $n$ and $\beta^0$ is sparse or approximately
sparse.\footnote{We focus on the case in which $\beta^0$ is exactly sparse,
but our theory can be extended to the approximately sparse case provided
that the approximation error is sufficiently small.} Suppose an analyst is interested in testing a system of linear
restrictions involving the subvector
$\beta_G^0 := (\beta_j^0)_{j\in G}$, where
$G\subseteq\{1,\ldots,p\}$ is an arbitrary index set with $|G|=m$.
Specifically, consider $r$ restrictions of the form
\begin{align*}
R_G\beta_G^0=q,
\qquad
R_G\in\mathbb{R}^{r\times m},
\qquad
\operatorname{rank}(R_G)=r\leq m,
\qquad
q\in\mathbb{R}^r,
\end{align*}
where $R_G$ and $q$ are specified by the analyst. Let the
\textit{restriction matrix} $R\in\mathbb{R}^{r\times p}$ be obtained by
padding $R_G$ with zero columns outside $G$. Then
$R\beta^0=R_G\beta_G^0$, and the testing problem can be written as
\begin{align}
H_0:R\beta^0=q
\qquad \text{vs.} \qquad
H_1:R\beta^0\neq q.
\label{eq:H_0}
\end{align}
By construction, $R$ has the same full row rank as $R_G$ and has nonzero
columns only on the coordinates indexed by $G$. When $m=r=1$, the testing
problem reduces to inference on a single linear restriction involving one
coefficient. More generally, we allow both $m$ and $r$ to diverge with
$n$, subject to the rate conditions imposed below. An important special
case is the joint test that a group of coefficients is zero, for which
$R_G=I_m$, $r=m$, and $q=0_m$.
\subsection{Contribution to the Literature}
In this paper, we propose a heteroskedasticity-robust dimension-adaptive test for a broad class of linear restrictions in sparse high-dimensional regressions. This procedure is asymptotically valid across a fixed or diverging number of restrictions, avoiding the need to classify problems into low or increasing dimensional asymptotic regimes. We also establish the power properties of the proposed test under alternatives for fixed designs, complementing the asymptotic size control results.
Our contribution lies at the intersection of three strands of literature: 1) debiased inference in high-dimensional linear models, 2) testing many linear restrictions, and 3) heteroskedasticity-robust inference. We start by building on the literature on inference for low-dimensional parameters in high-dimensional linear models \citep[e.g.][]{belloni2012,belloni2014inferenceRES,van-de-geer(2014),javanmard2014confidence,zhang2014confidence,ning2017general}. Our work is positioned within the general framework of ``desparsified'' or ``debiased'' LASSO pioneered by \cite{van-de-geer(2014)}, \cite{javanmard2014confidence}, and \cite{zhang2014confidence}. However, existing results are typically derived under relatively homogeneous designs with relatively light-tailed conditions on regressors or errors. We relax these stringent requirements to allow for heteroskedasticity both in the error terms and independent potentially non-independent covariate rows.
A common approach in the literature is to first establish the validity of a test under a fixed design and then show that an $i.i.d.$ random design satisfies the required fixed-design conditions with high probability. We depart from this approach in two respects. First, we provide high-level conditions for the validity of our test under a fixed design when general LASSO estimators are employed. These conditions are not tied to a particular random-design specification and may therefore be of independent interest. Second, we provide sufficient conditions under which a random design matrix
satisfies these fixed-design requirements. Importantly, our random-design results do not require the observations to be identically distributed. In
particular, we allow the observations to have different marginal distributions, subject to marginal sub-Weibull tail conditions. Our work accommodates the joint testing of multiple linear restrictions, with the number of restrictions allowed to increase with sample size. Finally, we characterize the asymptotic power of the proposed test under alternative hypotheses, complementing our theoretical results on its validity under the null.
Our work is closely related to \cite{caner2018asymptotically}, who also consider inference within the desparsified LASSO framework and allow for heteroskedastic, non-sub-Gaussian covariates and errors under finite-moment conditions. Their asymptotic normality result concerns a single linear contrast, although the number of coefficients entering that contrast is allowed to increase with the sample size. In contrast, we consider a general system of linear restrictions \(R\beta^0=q\), where the number of restrictions itself may diverge. Their formal results are established under i.i.d. covariates, although they note that their arguments could be extended to non-identically distributed covariates. Our random-design theory directly accommodates independent but non-identically distributed covariates under marginal sub-Weibull tail conditions.
Our work also speaks to research on testing many restrictions in high-dimensional linear regressions. \cite{shi2019linear} consider testing for a subvector of the regression coefficients for sparse generalized linear models, including sparse linear models as special cases, where the number of restrictions is smaller than the sample size of order $o(n^{1/3})$. They allow for heteroskedastic error terms but restrict them to be sub-exponential. \cite{zhu2018linear} propose a test statistic for a single contrast of the whole vector of the regression coefficients; \cite{zhang2025novel} extend this test to testing many restrictions where the number of restrictions can be larger than the sample size. Their analysis assumes that the regressors and errors are homoskedastic and the errors are sub-Gaussian, and require the number of restrictions to diverge for their test statistic to be valid. Building on recent work on concentration inequalities under general tail assumptions \citep[e.g.][]{kuchibhotla2022moving}, however, we allow a design matrix with heavier than sub-Gaussian tails and heteroskedastic errors satisfying only polynomial moments. Under favorable sparsity regimes, the restriction dimension can grow nearly at the rate of $n^{1/2}$ up to logarithmic factors. Importantly, our test statistic is valid regardless of whether the number of restrictions is bounded or diverges.
In a different context of moderate dimensional linear regressions, where the number of covariates is smaller than (but can be proportional to) the sample size, inference procedures that do not require a sparsity condition have been proposed under general heteroskedasticity \citep[e.g.][]{cattaneo2018inference,KSS2020,jochmans2022heteroscedasticity,Anatolyev-Solvsten(2023)}. In particular, \cite{Anatolyev-Solvsten(2023)} propose a test statistic for many restrictions that remains valid even when the number of restrictions is proportional to the sample size. We adapt their logic to the high-dimensional case, where the number of covariates can be much larger than the sample size and direct covariance inversion is unavailable, by combining desparsification with regularization of the restriction covariance matrix. Their procedure requires ${\bf X}^\top{\bf X}$ and
\[
\sum_{\ell\neq i,j,k} X_\ell X_\ell^\top
\]
to be invertible for every $i,j,k\in[n]$, where $X_\ell\in\mathbb{R}^p$ denotes the $\ell$-th row of ${\bf X}$. When $p\gg n$, however, these matrices are necessarily rank-deficient and hence non-invertible. Consequently, the OLS estimators underlying their procedure are not available in our setting, requiring us to modify their approach to accommodate high-dimensional covariates. \cite{Anatolyev-Solvsten(2023)} note the possibility of an alternative strategy when ${\bf X}^\top{\bf X}$ is non-invertible. However, as discussed in Footnote~7 of their paper, the resulting heteroskedasticity-robust Wald statistic still requires the inversion of a high-dimensional covariance matrix estimator, which poses substantial challenges for size control.
In the special case of testing many zero restrictions (i.e. $R = I_m$), \cite{dezeure2017high} and \cite{zhang2017simultaneous} propose bootstrap inference procedures based on the high-dimensional Gaussian approximation in \cite{chernozhukov2013}, while \cite{xia2018two} derive an analytical test statistic. These methods allow the number of restrictions to be as large as the number of covariates. In contrast, we set the number of restrictions to be relatively small compared with the sample size, but allow for a more general form of restrictions.
We summarize our contributions in Table \ref{tab:comp}.
The distributional assumptions satisfy the following
hierarchy of implications:
joint Gaussianity $\Rightarrow$ joint sub-Gaussianity
$\Rightarrow$ joint sub-Weibullity
$\Rightarrow$ marginal sub-Weibullity,
where the sub-Weibull order satisfies $0<\nu\leq 2$.
\begin{table}[ht]
\centering
\caption{Comparison with Related Tests}
\label{tab:comp}
\footnotesize
\begin{adjustbox}{max width=\linewidth}
\begin{tabular}{llll}
\toprule
\textbf{Paper}
& \textbf{Design assumptions}
& \textbf{Restrictions}
& \textbf{Error assumptions} \\
\midrule
VBRD14
& $i.i.d.$ joint sub-Gaussian rows
& $R$ selects components of $\beta^0$
& Homoskedastic, sub-Gaussian \\
DBZ17
& $i.i.d.$ joint sub-Gaussian rows
& $R=I_m$; fixed/diverging $r$; $\log m=o(n^{1/7})$
& Heteroskedastic, finite fourth moment \\
ZB18
& $i.i.d.$ joint Gaussian rows
& General $R$; $r=1$; $m=p$
& Homoskedastic, sub-Gaussian \\
SSCL19
& $i.i.d.$ joint sub-Gaussian rows
& General $R$; fixed/diverging $r$; $m=o(n^{1/3})$
& Heteroskedastic, sub-exponential \\
ZYL25
& $i.i.d.$ joint sub-Gaussian rows
& General $R$; $r\rightarrow\infty$
& Homoskedastic, sub-Gaussian \\
\textbf{Our paper}
& Independent, non-identically distributed,
& General $R$; fixed/diverging $r$;
$m=o\!\left(n/(s_0\log p)^2\right)$
& Heteroskedastic, finite fourth moment \\& marginally sub-Weibull rows \\
\bottomrule
\end{tabular}
\end{adjustbox}
\vspace{2mm}
\begin{minipage}{\linewidth}
\footnotesize
\textit{Note: VBRD14: \cite{van-de-geer(2014)}; DBZ17: \cite{dezeure2017high}; ZB18: \cite{zhu2018linear}; SSCL19: \cite{shi2019linear}; ZYL25: \cite{zhang2025novel}}
\end{minipage}
\end{table}
\subsection{Notation}
For a real vector $a$, let $\text{diag} (a)$ denote a diagonal matrix with diagonal entries given by $a$. For a real matrix $A \in \mathbb{R}^{s_1 \times s_2}$, denote its transpose by $A^\top$, its trace by $\text{tr} (A)$, its rank by $\text{rk} (A)$, its Frobenius norm by $\lVert A \rVert_{F}$, its $(i,j)$-th entry by $A_{ij}$. For any $p \in [1,\infty]$, let $\|\cdot\|_p$ denote the usual vector $\ell_p$-norm.
We define the induced (operator) matrix $p$-norm by $\|A\|_{p} \;:= \sup_{\|x\|_{p}=1}\|Ax\|_{p}$. We also define the entrywise max norm by $\|A\|_{\max} \;:=\; \max_{1\le i \le s_1,\; 1\le j \le s_2} |A_{ij}|$. For a square matrix $A$, let $\lambda_{\max} (A)$, $\lambda_i (A)$, $\lambda_{\min} (A)$ denote the largest, $i$-th largest and smallest eigenvalues, respectively. For $a,b\in\mathbb{R}$, we write $a \lesssim b$ if there exists a universal constant $C>0$, independent of $n$, such that $a \le C b$, and $a \asymp b$ if $a \lesssim b$ and $b \lesssim a$. We let $\overline{c}$ and $\underline{c}$ respectively denote constants for upper and lower bound. The generic notation of probability measures is denoted by $\ensuremath{\mathbb{P}}$.
\subsection{Structure of the paper}
The remainder of the paper is organized as follows.
Section~\ref{sec:test_statistic} motivates and constructs our proposed test.
Section~\ref{section:strong_approx} establishes a strong approximation result
for a general class of quadratic-form statistics.
Section~\ref{section:size_control} develops the asymptotic theory underlying
the validity of our test. We first establish, in
Theorem~\ref{thm:size-control_fixed_design}, sufficient conditions for size
control under a fixed design generally stated, and turn to random designs.
Proposition~\ref{prop:sufficient_conditions_for_fixed_design} provides general sufficient conditions under which a random design satisfies the requirements of
Theorem~\ref{thm:size-control_fixed_design}, thereby establishing validity in
the random-design setting. Theorem~\ref{thm:size-control_random_design_extension} extends the results to
the conservative LASSO.
Section~\ref{sec:sim} examines the finite-sample performance of the proposed
test through Monte Carlo simulations. Section~\ref{sec:empirical_application} illustrates its use in an empirical application. Section~\ref{sec:conclusion} suggests future research and concludes.
\section{Test Statistic: Motivation and Construction}
\label{sec:test_statistic}
\subsection{The General Idea} \label{subsec:general-idea}
In this section, we outline the central idea from \cite{van-de-geer(2014)} (see also \cite{javanmard2014confidence} and \cite{zhang2014confidence}), which provides motivation for both the problem and our approach. Their key insight is to invert the Karush-Kuhn-Tucker (KKT) condition characterizing the least absolute shrinkage and selection operator (LASSO) solution of \cite{tibshirani(1996)}. To formalize the problem, recall that the LASSO estimator is defined as
\begin{align}
\widehat{\beta}(\lambda_0) \in \arg \min_{\beta \in \mathbb{R}^p} \frac{1}{n} \lVert Y - {\bf X}\beta \rVert_2^2 + 2\lambda_0 \lVert \beta \rVert_1,
\label{eq:beta-hat}
\end{align}
where $\lambda_0 >0$ is the regularization parameter. For notational simplicity, we write $\widehat{\beta}\equiv\widehat{\beta}(\lambda_0)$ throughout the remainder of the paper, unless the dependence on $\lambda_0$ needs to be made explicit. The {\it desparsified LASSO estimator} in \cite{van-de-geer(2014)} is defined as
\begin{align}
\widehat{b} := \widehat{\beta} + \frac{1}{n} \widehat{\Theta} {\bf X}^\top ( Y - {\bf X}\widehat{\beta}),
\label{eq:b-hat-general}
\end{align}
where $\widehat{\Theta} \in \mathbb{R}^{p \times p}$ is the ``approximate inverse'' of $\widehat{\Sigma} := {\bf X^\top X}/n$, and we have
\begin{align}
\sqrt{n} (\widehat{b} - \beta^0) = W - \Delta, \quad W := \widehat{\Theta} {\bf X}^\top \varepsilon /\sqrt{n},
\label{eq:beta_hat_inference}
\end{align}
where $\Delta := \sqrt{n}(\widehat{\Theta} \widehat{\Sigma} - I_p)(\widehat{\beta} - \beta^0)$ is the remainder term.
To construct $\widehat{\Theta}$, we follow \cite{van-de-geer(2014)} to use the LASSO for the node-wise regression on the design matrix ${\bf X}$. For each $j = 1,...,p$, denote by $X_j\in \mathbb{R}^n$ the $j$-th column vector of ${\bf X}$ and
\begin{align}
\widehat{\gamma}_j(\lambda_j) \in \operatorname*{arg\,min}_{\gamma \in \mathbb{R}^{p-1}} \frac{1}{n} \lVert X_j - {\bf X}_{-j} \gamma \rVert_2^2 + 2\lambda_j \lVert \gamma \rVert_1,
\label{eq:gamma_j}
\end{align}
for some regularization parameter $\lambda_j >0$; for notational simplicity, we write
$\widehat{\gamma}_j \equiv \widehat{\gamma}_j(\lambda_j)$ throughout the remainder of the paper, unless the dependence on $\lambda_j$ needs to be made explicit. We can write the elements of $\widehat{\gamma}_j$ as $\{ \widehat\gamma_{j,k}: k=1,...,p, k \neq j\}$, and define
\begin{align*}
\widehat{C} :=
\begin{bmatrix}
1 & -\widehat{\gamma}_{1,2}& \cdots &-\widehat\gamma_{1,p} \\
-\widehat{\gamma}_{2,1}& 1 & \cdots& -\widehat\gamma_{2,p} \\
\vdots & \vdots & \ddots & \vdots \\
-\widehat{\gamma}_{p,1} & -\widehat{\gamma}_{p,2} & \cdots & 1
\end{bmatrix},
\end{align*}
and
\begin{align*}
\widehat{T}^2 := \text{diag}(\widehat{\tau}_1^2,...,\widehat{\tau}^2_p),
\end{align*}
with $\widehat{\tau}_j^2 := \lVert X_j - {\bf X}_{-j} \widehat{\gamma}_j \rVert_2^2/n + \lambda_j \lVert \widehat{\gamma}_j\rVert_1 $ for each $j = 1,...,p$. We then define
\begin{align*}
\widehat{\Theta} := \widehat{T}^{-2} \widehat{C}.
\end{align*}
In \cite{van-de-geer(2014)}, the properties of $\widehat \beta$ as well as $\widehat \Theta$ are derived when $\bf X$ is sub-Gaussian or bounded under homoskedasticity; they further assume $\varepsilon$ to be i.i.d. sub-exponential when $\bf X$ is sub-Gaussian. Using the concentration inequalities under general tail assumptions (e.g. \cite{kuchibhotla2022moving}), we derive these properties when $\bf X$ is sub-Weibull (see Appendix \ref{sec:appendix_sub_weibull} for definition) and $\varepsilon$ has only polynomial moments under general heteroskedasticity.
\subsection{Our Test Statistic} \label{sec:Test_Statistic}
To test the hypothesis in (\ref{eq:H_0}), we may pre-multiply both sides of \eqref{eq:beta_hat_inference} by the matrix \( R \) to obtain
\begin{align}
\sqrt{n} (R\widehat{b} - R\beta^0) = RW - R \Delta,
\label{eq:beta_hat_inference:R}
\end{align}
where, under the null, we have \( q = R\beta^0 \in \mathbb{R}^r \). Under stated regularity conditions, the remainder term $R \Delta$ is negligible compared to $RW$. If $r$ is fixed, it is then straightforward to derive a Wald test based on the asymptotic normality of $RW$ (e.g. \cite{van-de-geer(2014)}), or a heteroskedasticity-robust Wald test that also takes into account heteroskedasticity in error terms (e.g. \cite{caner2018asymptotically}). However, these test statistics might not be valid when we consider $r$ as diverging. In the context of moderate dimensional linear regressions, a modified Fisher's $F$ test for testing the hypothesis in (\ref{eq:H_0}) is proposed by \cite{Anatolyev-Solvsten(2023)}, which provides size-control under general heteroskedasticity regardless of whether $r$ is bounded or diverging.
This can be attractive in practice, especially when analysts do not want to make a preliminary judgment about whether $r$ is in the low- or high-dimensional regime. Another advantage of the modified $F$ test is that, in contrast with the heteroskedasticity-robust Wald test, it does not involve the inverse of an estimated sandwich matrix whose entries carry the noisy residuals, and is thus expected to be more stable in finite samples as demonstrated by our simulations.
We adapt this $F$ test to the current setting, with two important modifications. First, we replace $\widehat{\Sigma}^{-1}$ in $R \widehat{\Sigma}^{-1} R^\top$, which is not well defined when $p > n$, with $\widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top$, which can be interpreted as an approximation of $\widehat{\Sigma}^{-1}$. Second, we introduce ridge regularization to $R \widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top$, which is motivated by the observation that inverting $R \widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top$ can be unstable or even undefined in finite samples, especially when $r$ is large. We further derive a novel strong approximation result for quadratic forms, which allows us to construct a critical value that leads to valid inference under both fixed and random designs.
We introduce the ridge-regularized projection matrix
\begin{align*}
P_\theta := \frac{1}{n} {\bf X} \widehat{\Theta}^\top R^\top \left( R\widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top + \theta I_r \right)^{-1} R \widehat{\Theta} {\bf X}^\top,
\end{align*}
where the ridge regularization parameter $\theta > 0$ is determined solely by the matrix ${\bf X} \widehat{\Theta}^\top R^\top$. We denote the $(i,j)$-th entry of $P_\theta$ by $P_{\theta,ij}$. The test statistic that we propose in this paper is given by
\begin{align}
\widehat{Q}(q) := \frac{b(q)^\top \left( R\widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top + \theta I_r \right)^{-1} b(q)}{\sum_{i=1}^n P_{\theta,ii} \widehat{\sigma}_i^2}, \quad b(q) := \sqrt{n}(R \widehat{b} - q) \in \mathbb{R}^r,
\label{Q_hat_statistic_definition:1}
\end{align}
where $\widehat{\sigma}_i^2 := (Y_i - X_i \widehat{\beta})^2$ is the estimator for $\sigma_i^2: = \text{Var}(\varepsilon_i|{\bf X})$ and $\widehat{\beta}$ is given in \eqref{eq:beta-hat}. Here we use $\widehat{\beta}$ in $\widehat{\sigma}_i^2$ as we can derive the consistency of $\widehat{\sigma}_i^2$ directly from the properties of $\widehat{\beta}$. We also define the effective degrees of freedom associated with the regularization level $\theta$ as $K_\theta := \sum_{i \in [n]} P_{\theta,ii}$. The regularization parameter $\theta$ is selected to ensure desirable theoretical properties. Specifically, we set
\begin{align}
\theta \equiv \text{the } \lceil \varpi r_n \rceil\text{-th largest eigenvalue of }
R \widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top,
\label{eq:theta}
\end{align}
with $0<\varpi < 1$ that can depend on ${\bf X} \widehat{\Theta}^\top R^\top$, and $r_n := \text{rk}({\bf X} \widehat{\Theta}^\top R^\top)$. Note that $r_n \leq r$, where $r$ is the number of restrictions in our null hypothesis \eqref{eq:H_0}, which we allow to diverge. In our setup, $r_n$ is of the same order as $r$, but we allow for $r_n < r$, which motivates us to use the ridge regularization.
\begin{rem}
When $\theta$ is chosen arbitrarily, it is possible that $K_\theta$ and $r_n$ are not of the same order. This can lead to non-standard limiting distribution and complicate inference (e.g. \cite{KSS2020} and \cite{LWZ2024}). The choice in (\ref{eq:theta}) ensures that the order of $K_\theta$ is the same as $r_n$. In addition, under stated regularity conditions, we show that $\varpi$ can be selected in such a way that $\theta$ is bounded away from zero asymptotically when $\bf X$ is random (see Proposition \ref{prop:ridge_parameter}).
\end{rem}
\begin{rem}
As noted by \cite{caner2018asymptotically}, computing $R \widehat \Theta$ can be much easier than computing the whole $\widehat \Theta$, depending on the structure of $R$. For example, if $R = e_j^\top \in \mathbb{R}^{1 \times p}$ where $e_j$ is the $j$-th unit vector in $\mathbb{R}^p$, then only the $j$-th row of $\widehat \Theta$ needs to be computed. More generally, at most $m$ node-wise LASSO regressions are required for computing $R \widehat \Theta$. This reduces computational demands compared to running the node-wise LASSO regression $p$ times when $m \ll p$.
\end{rem}
\subsection{Construction of Critical Value}\label{sec:CV}
Our analytical test compares the test statistic $\widehat Q(q)$ with a robust critical value $\mathcal{C}_{\alpha}$ that is valid when $r$ is bounded as well as when $r$ diverges, where $\alpha \in (0,1)$ is the significance level. We will reject $H_0: R\beta^0 = q$ at the $\alpha$ significance level whenever
\begin{align}
\widehat Q(q) > \mathcal{C}_{\alpha}.
\label{eq:rejection-criteria}
\end{align}
Inspired by \cite{Anatolyev-Solvsten(2023)}, we define this critical value as
\begin{align}
\mathcal{C}_{\alpha} := 1 + \frac{ \sqrt{ K_\theta \widehat{\Phi}}}{ \sum_{i \in [n]} P_{\theta,ii} \widehat{\sigma}_i^2 } \left( \frac{ q_{1-\alpha} (F_{\widetilde{w}_n})-1}{ \sqrt{2 \sum_{i \in [r_n] } \widetilde{w}_{i,n}^2} } \right),
\label{eq:critical_value}
\end{align}
where $\{\widetilde{w}_{i,n}\}_{i \in [r_n]}$ and $\widehat{\Phi}$ are defined later. Here $q_{1-\alpha}(F_{\widetilde{w}_n})$ is the $(1-\alpha)$-quantile of the weighted chi-square random variable
\[
F_{\widetilde{w}_n} := \sum_{i \in [r_n]} \widetilde{w}_{i,n}\,\chi^2_{1,i}
\]
conditional on the data $(Y, {\bf X})$, where $\{\chi_{1,i}^2\}_{i \in [r_n]}$ are $r_n$ independent chi-squared random variables with one degree of freedom. This quantile can be obtained via simulation given the weights~$\{\widetilde{w}_{i,n}\}_{i \in [r_n]}$.
To see the intuition of the critical value, we need to explain the limiting distribution of our test statistic $\widehat Q(q)$ under the null $R \beta^0= q$. When $r$ is bounded, by using the strong approximation result we derive in this paper, we show that
\begin{align}
\widehat{Q}(q) \overset{d}{\approx} \mathcal{Z}^\top \mathcal{D}_n \mathcal{Z} = \sum_{i \in [r_n]} w_{i,n} \chi^2_{1,i},
\label{eq:heuritsics_fixed_K}
\end{align}
where $\mathcal{Z} \overset{d}{=} \mathcal{N}(0,I_{r_n}) $ and $\mathcal{D}_n:= \text{diag}(w_{1,n},...,w_{r_n,n})$ contains the decreasingly ordered non-zero eigenvalues of
\begin{align*}
\Omega := \frac{ (R \widehat{\Theta} {\bf X}^\top \Lambda {\bf X} \widehat{\Theta}^\top R^\top)^{1/2} (R \widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top + \theta I_r)^{-1} (R \widehat{\Theta} {\bf X}^\top \Lambda {\bf X} \widehat{\Theta}^\top R^\top)^{1/2} }{n \sum_{i \in [n]} P_{\theta,ii} \sigma_i^2},
\end{align*}
with $\Lambda := \text{diag}(\sigma_1^2,\cdots,\sigma_n^2)$. Note that $\text{tr} (\Omega) = 1$ and $\Omega$ is positive semi-definite, so that $w_n := (w_{1,n},\cdots,w_{r_n,n})^\top$ are positive and sum up to one. Using the ``plug-in" rule, we let $\widetilde{w}_n := (\widetilde{w}_{1,n},\cdots,\widetilde{w}_{r_n,n})^\top$ be the first $r_n$ decreasingly ordered eigenvalues of
\begin{align*}
\widehat{\Omega}:= \frac{ (R \widehat{\Theta} {\bf X}^\top \widehat\Lambda {\bf X} \widehat{\Theta}^\top R^\top)^{1/2} (R \widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top + \theta I_r)^{-1} (R \widehat{\Theta} {\bf X}^\top \widehat\Lambda {\bf X} \widehat{\Theta}^\top R^\top)^{1/2} }{n \sum_{i \in [n]} P_{\theta,ii} \widehat{\sigma}_i^2},
\end{align*}
where $\widehat \Lambda := \text{diag}(\widehat{\sigma}_1^2,\cdots ,\widehat{\sigma}_n^2)$. Then, note that the critical value $\mathcal{C}_{\alpha}$ can be rearranged as
\begin{align*}
q_{1-\alpha}(F_{\widetilde{w}_n}) + (q_{1-\alpha}(F_{\widetilde{w}_n})-1) \left( \frac{ {\sqrt{ K_\theta \widehat{\Phi}}} }{\left({\sum_{i \in [n]} P_{\theta,ii} \widehat{\sigma}_i^2}\right)\sqrt{2 \sum_{i \in [r_n]} \widetilde{w}_{i,n}^2 } } -1 \right),
\end{align*}
and we show that the second term is asymptotically negligible. In addition, the first term is close to $q_{1-\alpha}(F_{w_n})$, the $(1-\alpha)$-quantile of the weighted chi-square random variable
\[
F_{w_n} := \sum_{i \in [r_n]} w_{i,n}\,\chi^2_{1,i}
\]
conditional on $\bf X$. This guarantees that our analytical test achieves the correct asymptotic size under bounded $r$.
If, however, $r$ diverges, the weighted chi-square random variable converges in distribution to the standard normal distribution
\begin{align*}
\left( \sum_{i \in [r_n]} 2\tilde w_{i,n}^2 \right)^{-1/2} (F_{\tilde{w}_n}-1) \rightsquigarrow \mathcal{N}(0,1),
\end{align*}
and we show that the standardized $\widehat Q(q)$ also converges in distribution to the standard normal distribution
\begin{align*}
\Phi^{-1/2} \left[\frac{1}{\sqrt{K_\theta}}\sum_{i \in [n]} P_{\theta,ii} \widehat{\sigma}_i^2\right] \left(\widehat Q(q)-1\right) \rightsquigarrow \mathcal{N}(0,1),
\end{align*}
where
\begin{align*}
\Phi := \frac{2}{K_\theta} \sum_{i \in [n]}\sum_{j \neq i}P_{\theta,ij}^2 \sigma_i^2 \sigma_j^2.
\end{align*}
Motivated by the work of \cite{crudu2021} and \cite{dovi2024ridge}, we consider the estimator
\begin{align}
\widehat{\Phi} := \frac{2}{K_\theta} \sum_{i \in [n]} \sum_{j \neq i} P_{\theta,ij}^2 \widehat{\sigma}_i^2 \widehat{\sigma}_j^2 \label{Phihat_def}
\end{align}
for $\Phi$. As we show the consistency of $\widehat \Phi$, the asymptotic size-control of our test is also guaranteed under diverging $r$.
\begin{rem}
In our simulations, we use the following finite-sample adjustment of
$\widehat{\Phi}$, motivated by the diagonal components:
\begin{align*}
\widetilde{\Phi}
:=
\widehat{\Phi}
+
\frac{1}{K_\theta}
\sum_{i\in[n]}
P_{\theta,ii}^2\widehat{\sigma}_i^4.
\end{align*}
The motivation for this adjustment is that, although the diagonal component of
$\varepsilon^\top P_\theta \varepsilon/\sqrt{K_\theta}$ is asymptotically negligible under
our assumptions, it may still contribute non-negligible variability in finite
samples. In particular, the diagonal component is given by
\begin{align*}
\frac{1}{\sqrt{K_\theta}}
\sum_{i\in[n]}P_{\theta,ii}\varepsilon_i^2,
\end{align*}
and the additional term
\begin{align*}
\frac{1}{K_\theta}
\sum_{i\in[n]}P_{\theta,ii}^2\widehat{\sigma}_i^4
\end{align*}
estimates its contribution to the variance. Thus, $\widetilde{\Phi}$ accounts
for a source of finite-sample variability that is omitted from
$\widehat{\Phi}$ but vanishes asymptotically. Since
\begin{align*}
\frac{1}{K_\theta}
\sum_{i\in[n]}P_{\theta,ii}^2\widehat{\sigma}_i^4
= o_p(1)
\end{align*}
under the conditions of our asymptotic analysis, replacing
$\widehat{\Phi}$ by $\widetilde{\Phi}$ does not affect any of our
asymptotic results.
\end{rem}
\section{Strong Approximation of Quadratic Forms} \label{section:strong_approx}
This section focuses on approximating the quadratic form in the numerator of \eqref{Q_hat_statistic_definition:1} by a Gaussian counterpart, without imposing any restriction on the dimension \( r \). As this result may be of independent interest to readers, we present it in a standalone section with greater generality than required for the results of this paper. For this purpose, let $P({\bf X})$ be a generic $n \times n$ symmetric positive semi-definite matrix that depends only on $\bf X$, and satisfies the projection property that $\sum_{j \in [n]} P_{ij}({\bf X})^2 \leq P_{ii}({\bf X}), \forall i \in [n]$. We treat $\bf X$ as fixed throughout the section for simplicity. We also impose the following assumptions on the data generating process (DGP).
\begin{ass}
The following assumptions hold:
\begin{enumerate}
\item The errors $\{ \varepsilon_i \}_{i \in [n]}$ are independent and mean-zero conditional on {\bf X}.
\item There exists a universal constant $\underline{c}>0$ such that $\inf_{i \in [n]}\sigma_i^2 \geq \underline{c}$.
\item There exists a universal constant $\overline{c}<\infty$ such that $\sup_{i \in [n]} \mathbb{E}[|\varepsilon_i|^h|{\bf X}] \leq \overline{c}$ for some $h \geq 4$.
\item The maximum diagonal element of $P({\bf X})$, $\tilde p_n := \max_{i \in [n]} P_{ii}({\bf X})$, satisfies ${\tilde p_n}/{K_n} = O (n^{-1})$ with $K_n := \sum_{i \in [n]} P_{ii}({\bf X})$.
\end{enumerate}
\label{ass:mean_zero}
\end{ass}
\begin{rem}
Assumptions \ref{ass:mean_zero}.1--\ref{ass:mean_zero}.3 are standard regularity conditions and allow for general heteroskedasticity in error terms. Assumption \ref{ass:mean_zero}.4 requires the maximum leverage score $\tilde p_n$ to vanish sufficiently quickly. Consider the special case where $P({\bf X})$ is the ridge-regularized projection matrix $P_{\theta}$ with $\theta = 0$. If the diagonals of $P_\theta$ are well-balanced (e.g. $p_n \leq \overline{c} \, {r_n}/{n}$ for some constant $\overline{c} < \infty$), then noting $K_\theta = \text{tr}(P_\theta) = r_n$, we obtain $p_n/K_\theta = p_n / r_n \lesssim 1/n$.
\end{rem}
\begin{thm}[Strong Approximation for Quadratic Forms]
Suppose Assumption~\ref{ass:mean_zero} holds. Let \( \{ \mathcal{E}_i \}_{i \in [n]} \) be a sequence of independent Gaussian random variables such that \( \mathcal{E}_i \overset{d}{=} \mathcal{N}(0, \sigma_i^2) \), where each \( \mathcal{E}_i \) is independent of \( \{ \varepsilon_j \}_{j \in [n]} \). Suppose in addition that $$\frac{1}{K_n} \sum_{i \in [n]} \sum_{j \neq i} P_{ij}^2({\bf X}) \gtrsim 1.$$ Then for any $\zeta \in (0,1)$
\begin{align}
\sup_{y \in \mathbb{R}}\left| \mathbb{P}_\varepsilon \left( \sum_{i,j \in [n], j \neq i} P_{ij}({\bf X}) \varepsilon_i \varepsilon_j \leq y \right) - \mathbb{P}_\mathcal{E} \left( \sum_{i,j \in [n], j \neq i} P_{ij}({\bf X}) \mathcal{E}_i \mathcal{E}_j \leq y \right) \right| \leq C_\zeta \left( \frac{\tilde p_n}{K_n} \right)^{\frac{1-\zeta}{14-2\zeta}}
\label{eq:strong-approximation:1}
\end{align}
where $C_\zeta$ depends only on $\zeta$. Furthermore, for any $\ell>0$,
\begin{align}
\sup_{y \in \mathbb{R}} \left| \mathbb{P}_\varepsilon \left( \varepsilon^\top P({\bf X}) \varepsilon \leq y \right)
- \mathbb{P}_\mathcal{E} \left( \mathcal{E}^\top P({\bf X}) \mathcal{E} \leq y \right) \right| \leq C_\zeta \left( \frac{\tilde p_n}{K_n} \right)^{\frac{1-\zeta}{14-2\zeta}} + C_\zeta \ell^{\frac{1-\zeta}{2}} + \overline{C}\frac{\tilde p_n}{\ell^2},
\label{eq:strong-approximation:2}
\end{align}
where $\overline{C}$ is a universal constant.
\label{theorem:strong_approx}
\end{thm}
\begin{rem}
The term $C_\zeta \ell^{(1-\zeta)/2}$ in \eqref{eq:strong-approximation:2} arises from an anti-concentration inequality for non-Gaussian errors, which bounds the probability that the off-diagonal sum $\sum_{i \in [n]} \sum_{j \neq i} P_{ij}({\bf X}) \varepsilon_i \varepsilon_j$ falls within a neighborhood of radius $\ell$ centered at $y$, uniformly over $y$; see Lemma \ref{lem:anti-concentration:main} in Appendix \ref{sec:appendix_main_thm}. The term $\tilde p_n / (\ell^2)$ captures the approximation error induced by incorporating the diagonal elements of $P({\bf X})$.
\end{rem}
If we assume that $K_n/n = o(1)$, then $\tilde p_n = O(K_n/n) = o(1)$, and we can set $\ell = (K_n/n)^{1/4} = o(1)$ in the right-hand side of \eqref{eq:strong-approximation:2}, which yields
$\tilde p_n/\ell^2=O(\sqrt{K_n/n})=o(1)$, so that the errors $\varepsilon_i$ can be treated as approximately normal under such quadratic forms. For the ridge-regularized projection matrix $P_{\theta}$, it suffices to require that $r = o(n)$, since then $K_\theta/n = O(r_n/n) = o(1)$ by $r_n \leq r$. The following corollary then follows immediately.
\begin{cor} \label{cor:strong_approx}
Suppose Assumption \ref{ass:mean_zero}.1-\ref{ass:mean_zero}.3 hold and that $p_n = O(K_\theta/n)$ where $p_n := \max_{i \in [n]} P_{\theta,ii}$, and $\sum_{i \in [n]} \sum_{j \neq i} P_{\theta,ij}^2 / K_\theta \gtrsim 1$. If we further assume $r = o(n)$, then
\begin{align*}
\sup_{y \in \mathbb{R}}\left| \mathbb{P}_\varepsilon \left( \sum_{i,j \in [n], j \neq i} P_{\theta,ij} \varepsilon_i \varepsilon_j \leq y \right) - \mathbb{P}_\mathcal{E} \left( \sum_{i,j \in [n], j \neq i} P_{\theta,ij} \mathcal{E}_i \mathcal{E}_j \leq y \right) \right| = o(1),
\end{align*}
and
\begin{align*}
\sup_{y \in \mathbb{R}} \left| \mathbb{P}_\varepsilon \left( \varepsilon^\top P_\theta \varepsilon \leq y \right)
- \mathbb{P}_\mathcal{E} \left( \mathcal{E}^\top P_\theta \mathcal{E} \leq y \right) \right| = o(1).
\end{align*}
\end{cor}
The strong approximation result implies that the distribution of $\widehat{Q}(q)$ can be approximated by a weighted chi-square distribution when $r$ is bounded (see~\eqref{eq:heuritsics_fixed_K}). This approximation is essential for establishing the validity of the critical value defined in~\eqref{eq:critical_value}. Furthermore, when $r$ is diverging, we may approximate the distribution of $\widehat{Q}(q)$ by a standard normal distribution.
\section{Asymptotic Theory} \label{section:size_control}
In this section, we show theoretically that our test has asymptotically exact size control under both fixed and random designs. In particular, our inference procedure is valid regardless of whether the number of restrictions $r$, the number of coefficients involved in the hypothesis being tested $m$, and the number of coefficients $p$ are treated as bounded or diverging, thereby providing uniform inference over a large class of DGPs.
\subsection{Fixed Design}
We first consider the fixed design. The following theorem gives general conditions for our test to have asymptotically exact size control. Random designs can then be deduced from this theorem by verifying the conditions below. Note that these conditions are not restricted to the ordinary LASSO estimator. For example, one can also replace $\widehat \beta$ and $\widehat \gamma_j, j \in [p]$ with alternative LASSO estimators, such as the adaptive LASSO estimator (e.g. \cite{zou2006adaptive}) or the conservative LASSO estimator (e.g. \cite{caner2018asymptotically}), in the definition of $\widehat{Q}(q)$, and the result here can be applied to those inference procedures as well. We provide one such extension below.
\begin{thm}
Suppose that ${\bf X}$ is taken as fixed. Let \( \{ \Lambda_n \}_{n \in \mathbb{N}} \) be a sequence of classes of DGPs under which the following conditions hold uniformly over each \( \Lambda_n \):
\begin{enumerate}
\item Assumption \ref{ass:mean_zero}.1-\ref{ass:mean_zero}.3 hold, $p_n = O(K_\theta/n)$ where $p_n = \max_{i \in [n]} P_{\theta,ii}$, and $r=o(n)$.
\item The matrix \( R \in \mathbb{R}^{r \times p} \) has at most \( m \) non-zero column vectors, and satisfies $\lambda_{\max}(R R^\top) \leq \overline{c} < \infty$.
\item $\left(\sqrt{r} + \lVert {\bf X} (\widehat{\beta} - \beta^0 ) \rVert_{\infty} \right) \times \lVert {\bf X} (\widehat{\beta} - \beta^0 ) \rVert_2 = o_\ensuremath{\mathbb{P}}(\sqrt{n})$ and $\lVert \Delta \rVert_\infty = o_\ensuremath{\mathbb{P}}(1/\sqrt{m})$.
\item There exists $\varpi \in (0,1)$ that is bounded away from zero, such that $\theta \gtrsim 1$ and $$\frac{1}{K_\theta} \sum_{i \in [n]} \sum_{j \neq i} P_{\theta,ij}^2 \gtrsim 1$$
for all $n$ large enough.
\end{enumerate}
Then, under the null,
\begin{align*}
\liminf_{n \rightarrow \infty} \inf_{\lambda_n \in \Lambda_n} \mathbb{P}_{\lambda_n} \left( \widehat{Q}(q) > \mathcal{C}_\alpha \right) = \limsup_{n \rightarrow \infty} \sup_{\lambda_n \in \Lambda_n} \mathbb{P}_{\lambda_n} \left( \widehat{Q}(q) > \mathcal{C}_\alpha \right) = \alpha.
\end{align*}
\label{thm:size-control_fixed_design}
\end{thm}
\begin{rem}
We make three remarks regarding the preceding conditions.
First, the condition on $\lambda_{\max}(RR^\top)$
controls the magnitude of the restriction matrix $R$.
For any full-row-rank restriction matrix
$R\in\mathbb{R}^{r\times p}$, we may define
$\widetilde R=(RR^\top)^{-1/2}R$ and
$\widetilde q=(RR^\top)^{-1/2}q$.
The null hypotheses $R\beta^0=q$ and
$\widetilde R\beta^0=\widetilde q$ are equivalent,
while $\widetilde R\widetilde R^\top=I_r$.
Consequently, any full-row-rank system of linear
restrictions admits an equivalent normalized
representation satisfying
$\lambda_{\max}(\widetilde R\widetilde R^\top)=1$. Throughout the construction of the test statistic, we may therefore adopt this normalized representation of the null hypothesis.
Second, the rate conditions on
$\lVert {\bf X}(\widehat{\beta}-\beta^0)\rVert_\infty$
and
$\lVert {\bf X}(\widehat{\beta}-\beta^0)\rVert_2$
ensure the consistency of the variance estimators
and the estimated weights of the weighted
chi-squared distribution.
The rate condition on $\lVert\Delta\rVert_\infty$
ensures that the remainder term $R\Delta$
is asymptotically negligible in $\widehat Q(q)$.
Third, we impose a condition involving
$\varpi\in(0,1)$ to ensure that
$\left(
R\widehat{\Theta}\widehat{\Sigma}
\widehat{\Theta}^\top R^\top+\theta I_r
\right)^{-1}$
and $\Phi$ are well-defined in the limit.
In the next section, we establish sufficient
conditions under which these requirements hold
for marginally sub-Weibull covariates and
heteroskedastic errors satisfying finite
polynomial moment conditions.
\end{rem}
\begin{rem}
While our primary focus is on establishing the asymptotic validity of the proposed test under the null hypothesis, its behavior under the alternative is also of interest. In Section~\ref{sec:power_analysis} of the Supplemental
Appendix, we establish power properties of the proposed test under a fixed design. To maintain the focus of the main text on asymptotically exact size control, we defer this analysis to the Supplemental Appendix.
\end{rem}
\subsection{Random Design}
We now turn to the random design. Let $\widetilde{X}_i^\top \in \mathbb{R}^p$ for $i \in [n]$ be the row vectors of the design matrix $\bf X$, with $\Sigma_i := \mathbb{E}[\widetilde{X}_i \widetilde{X}_i^\top] \in \mathbb{R}^{p \times p}$. We define the population covariance matrix as $\Sigma := \sum_{i \in [n]} \Sigma_i /n \in \mathbb{R}^{p \times p}$, and $\Theta:= \Sigma^{-1} \in \mathbb{R}^{p \times p}$, which we assume exists. We also define $\eta_j:= X_j - {\bf X}_{-j} \gamma_j \in \mathbb{R}^n$ where $\gamma_j \in \arg \min_{\gamma \in \mathbb{R}^{p-1}} \mathbb{E} \lVert X_j - {\bf X}_{-j} \gamma \rVert_2^2/n$, and we assume that $\max_{i \in [n], j \in [p]} \mathbb E \left| \eta_{j}(i) \right|^l \leq \overline{c} < \infty$ for some $l \geq 4$. Finally, we define $s_0 := \left|\{j: \beta^0_j \neq 0 \}\right| \vee 1$, $s_j := \left| \{ k \neq j:\Theta_{j,k} \neq 0\} \right| $ and $s_0^* := \max_{j \in [p]} s_j \vee 1$.
\begin{ass} \label{ass:size}
\begin{enumerate}
\item Assumption \ref{ass:mean_zero}.1 - \ref{ass:mean_zero}.3 hold and $p_n = O_\ensuremath{\mathbb{P}}(K_\theta/n)$ where $p_n = \max_{i \in [n]} P_{\theta,ii}$.
\item The matrix \( R \in \mathbb{R}^{r \times p} \) has at most \( m \) non-zero column vectors, and satisfies $\lambda_{\max}(R R^\top) \leq \overline{c} < \infty$ and $\text{tr} (R \Theta R^\top)/r \geq \underline{c} > 0$.
\item $(\tilde X_i, \varepsilon_i)^\top \perp\!\!\!\perp (\tilde X_j, \varepsilon_j)^\top$ for $i \neq j$, and conditional on $\bf X$, $\{ \varepsilon_i \}_{i \in [n]}$ are independent and mean-zero. $\tilde X_i$ are marginally sub-Weibull $(\nu, K_{n,p})$ \footnote{A random variable $Z$ is said to be sub-Weibull $(\nu, K_{n,p})$ if $\|Z\|_{\psi_\nu} \le K_{n,p}$ with $\|Z\|_{\psi_\nu} :=\inf\Bigl\{c>0:\ \mathbb{E}\,\psi_\nu\bigl(|Z|/c\bigr)\le 1\Bigr\}$, where $\psi_\nu(t):=\exp(t^\nu)-1$ for $t\ge 0$. See Appendix \ref{sec:appendix_sub_weibull} for further discussion.} for some $0 < \nu \leq 2$ and $K_{n,p} = O(1)$, $1/\lambda_{\min}(\Sigma) = O(1)$ and $\max_{j \in [p]} \Sigma_{jj} = O(1)$.
\item The following rate conditions hold:
\begin{enumerate}[label=(\arabic*)]
\item $\max \{s_0 (\log(np))^{4/\nu}, (s_0^*)^2 \log(np), ms_0^* (\log(np))^{4/\nu}, m s_0^2 \log^2(np), m^2 s_0^* \log(np)\} = o(n)$,
\item $(\log(np))^{2/\nu - 1/2} = o(n^{1/2-1/h})$,
\item $(\log(np))^{2/\nu - 1/2} = o(n^{1/2-1/l})$.
\end{enumerate}
\end{enumerate}
\end{ass}
\begin{rem}
Assumption \ref{ass:size}.1 is the random design version of Assumption \ref{ass:mean_zero} when $\bf X$ is considered as random. The condition $\text{tr} (R \Theta R^\top)/r \geq \underline{c} > 0$ in Assumption \ref{ass:size}.2 ensures that there are sufficiently many large eigenvalues for $R\widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top$, so that $\theta$ can be bounded away from zero; note that this condition is weaker than $\lambda_{\min} (\Theta) \geq \underline{c} > 0$ and $\lambda_{\min} (R R^\top) \geq \underline{c} > 0$. Assumption \ref{ass:size}.3 allows for i.n.i.d. $\tilde X_i^\top$ and quite general tail conditions for the entries of $\tilde X_i^\top$; the conditions on $\Sigma$ are natural extensions of standard conditions on the population covariance matrix of $\bf X$ (e.g. \cite{van-de-geer(2014)}) to i.n.i.d. data, and note that $1/\lambda_{\min}(\Sigma) = O(1)$ is sufficient for $\Sigma$ to satisfy the restricted eigenvalue condition. The rate conditions in Assumption \ref{ass:size}.4 extend the results in \cite{van-de-geer(2014)} for $\nu = 2$. Our results generalize beyond those in \cite{van-de-geer(2014)} by not restricting the joint distribution of $\tilde X_i^\top$ (e.g. joint sub-Gaussian), and we allow $\varepsilon$ to have only polynomial moments. \cite{caner2018asymptotically} consider the even more general case when both $\bf X$ and $\varepsilon$ only have polynomial moments, so that their rate conditions can only allow $p$ to be polynomial of $n$, while we still allow $p$ to be exponential in polynomial of $n$. Finally, we allow $m$ to diverge at a slower rate than $n$. Consider the important special case where $s_0^*=O(1)$, $s_0=O(n^{1/4})$, and $p=an$ for some $a>0$, then we allow $m$ to be at most $o(n^{1/2})$ up to a logarithmic factor.
\end{rem}
\begin{prop}
Suppose Assumption \ref{ass:size} holds and $\lambda_j \asymp \sqrt{\log(np)/n}$ (uniformly in $j \in [p]$) are suitably chosen. There exist $c_{\varpi} \in (0,1]$ and $\underline{c} > 0$ such that, with probability approaching one, $\theta \geq \underline{c}$ and $\Phi \geq \underline{c}$ for all $\varpi \in (0,c_{\varpi})$.
\label{prop:ridge_parameter}
\end{prop}
\begin{rem}
Proposition \ref{prop:ridge_parameter} shows how we can choose $\varpi$ such that $\left( R\widehat{\Theta} \widehat{\Sigma} \widehat{\Theta}^\top R^\top + \theta I_r \right)^{-1}$ and $\Phi$ are well-defined in the limit under random design. We discuss in Appendix \ref{subsec:ridge_parameter} (see Remark \ref{rem:varpi}) how to choose $\varpi$ in a data-driven way. Our simulations suggest that both the size and the power of our test statistic are actually not very sensitive to the choice of $\varpi$.
\end{rem}
Before establishing the validity of our test under general random-design settings, we first provide an intermediate result showing that Assumption~\ref{ass:size}, together with suitable choices of the tuning parameters, is sufficient for the conditions of Theorem~\ref{thm:size-control_fixed_design} to hold.
\begin{prop}[Sufficient conditions for Theorem~\ref{thm:size-control_fixed_design}]
Suppose Assumption~\ref{ass:size} holds,
$\lambda_0 \asymp \sqrt{\log(np)/n}$,
$\lambda_j \asymp \sqrt{\log(np)/n}$ uniformly over $j \in [p]$,
and $\varpi \in (0,1)$ is suitably chosen.
Then all the conditions of
Theorem~\ref{thm:size-control_fixed_design} are satisfied.
\label{prop:sufficient_conditions_for_fixed_design}
\end{prop}
We are now ready to establish the main result for our test under random design, which follows directly from Proposition~\ref{prop:sufficient_conditions_for_fixed_design}.
\begin{thm}
Let \( \{ \Lambda_n \}_{n \in \mathbb{N}} \) be a sequence of classes of DGPs. Suppose that $\lambda_0 \asymp \sqrt{\log(np)/n}$, $\lambda_j \asymp \sqrt{\log(np)/n}$ (uniformly in $j \in [p]$) and $\varpi \in (0,1)$ are suitably chosen. If Assumption \ref{ass:size} holds uniformly over each \( \Lambda_n \), then under the null,
\begin{align*}
\liminf_{n \rightarrow \infty} \inf_{\lambda_n \in \Lambda_n} \mathbb{P}_{\lambda_n} \left( \widehat{Q}(q) > \mathcal{C}_\alpha \right) = \limsup_{n \rightarrow \infty} \sup_{\lambda_n \in \Lambda_n} \mathbb{P}_{\lambda_n} \left( \widehat{Q}(q) > \mathcal{C}_\alpha \right) = \alpha.
\end{align*}
\label{thm:size-control_random_design}
\end{thm}
\subsection{Extension to Conservative LASSO}
\label{subsec:conservative_LASSO}
In practice, it may be desirable to employ alternative LASSO estimators in
the construction of our test statistic to improve its finite-sample performance. In this section, we consider such an extension using the
conservative LASSO estimator proposed by \cite{caner2018asymptotically}. We let $\widehat\beta_j$ denote the $j$th element of $\widehat\beta$ from the original LASSO estimator \eqref{eq:beta-hat}, and consider the conservative LASSO
\begin{align*}
\widehat{\beta}_{\rm C} \in \arg \min_{\beta \in \mathbb{R}^p} \frac{1}{n} \lVert Y - {\bf X}\beta \rVert_2^2 + 2\lambda_0 \sum_{j \in [p]}\hat w_j |\beta_j|,
\end{align*}
where $\hat w_j = \lambda_{{\rm prec},0}/(|\widehat{\beta}_j| \vee \lambda_{{\rm prec},0})$ for $j \in [p]$, and $\lambda_{{\rm prec},0}$ is a prespecified tuning parameter.\footnote{We follow \cite{caner2018asymptotically} and set $\lambda_{\mathrm{prec},j}=\frac{9}{4}\lambda_n$, where $\lambda_n$ denotes the penalty parameter for the ordinary LASSO. In our implementation, $\lambda_n$ is either chosen to satisfy $\lambda_n \asymp \sqrt{\log(np)/n}$ or selected using BIC or GIC; see section \ref{sec:sim:description_of_test} for more details.} Similarly, define
\begin{align*}
\widehat{\gamma}_{{\rm C},j} \in \operatorname*{arg\,min}_{\gamma \in \mathbb{R}^{p-1}} \frac{1}{n} \lVert X_j - {\bf X}_{-j} \gamma \rVert_2^2+ 2\lambda_j \sum_{k\neq j} \hat w_{j,k} |\gamma_k|,
\end{align*}
where $\hat w_{j,k} = \lambda_{{\rm prec},j}/(|\widehat{\gamma}_{j}(k)| \vee \lambda_{{\rm prec},j})$ and we recall that $\widehat{\gamma}_{j}(k)$ is the $k$th element of $\widehat\gamma_j \in \mathbb{R}^{p-1}$ defined in \eqref{eq:gamma_j}.
Following this, we modify $\widehat{\Theta}$, $\widehat{b}$, and $\theta$ in the test statistic as described in Section~\ref{sec:Test_Statistic}, and modify $\{\widehat{\sigma}_i^2\}_{i\in[n]}$, $\{\widetilde{w}_{i,n}\}_{i\in[r_n]}$, and $\widehat{\Phi}$ in the construction of the critical value in Section~\ref{sec:CV} accordingly. We denote the resulting test statistic and critical value by $\widehat{Q}^{\rm C}(q)$ and $\mathcal{C}^{\rm C}_\alpha$, respectively.
\begin{thm}
Let \( \{ \Lambda_n \}_{n \in \mathbb{N}} \) be a sequence of classes of DGPs. Suppose that $\lambda_{{\rm prec},0} \asymp \lambda_0 \asymp \sqrt{\log(np)/n}$, $\lambda_{{\rm prec},j} \asymp \lambda_j \asymp \sqrt{\log(np)/n}$ (uniformly in $j \in [p]$) and $\varpi \in (0,1)$ are suitably chosen. If Assumption \ref{ass:size} holds uniformly over each \( \Lambda_n \) and $s_0^* = O(1)$, then under the null,
\begin{align*}
\liminf_{n \rightarrow \infty} \inf_{\lambda_n \in \Lambda_n} \mathbb{P}_{\lambda_n} \left( \widehat{Q}^{\rm C}(q) > \mathcal{C}^{\rm C}_\alpha \right) = \limsup_{n \rightarrow \infty} \sup_{\lambda_n \in \Lambda_n} \mathbb{P}_{\lambda_n} \left( \widehat{Q}^{\rm C}(q) > \mathcal{C}^{\rm C}_\alpha \right) = \alpha.
\end{align*}
\label{thm:size-control_random_design_extension}
\end{thm}
\section{Monte Carlo Study}
\label{sec:sim}
In this section, we present Monte Carlo simulation results comparing our proposed test (LMTW) with alternative tests in the literature based on the desparsified LASSO framework. In particular, we compare LMTW against the homoskedastic Wald test based on \cite{van-de-geer(2014)} and a heteroskedasticity-robust Wald test based on the conservative LASSO of \cite{caner2018asymptotically}, extended to accommodate $r \geq 1$. For fixed $r$, the theoretical results in these papers can be extended to establish the validity of the corresponding tests under their respective assumptions. It is less clear, however, whether and how their theoretical arguments can be extended to the setting where $r \rightarrow \infty$.
We also consider a variant of LMTW in which the ordinary LASSO estimators are replaced by the conservative LASSO estimators of \cite{caner2018asymptotically}. All reported results are based on $1{,}000$ Monte Carlo replications.
\subsection{Description of Tests}
\label{sec:sim:description_of_test}
We compare the finite-sample performance of our proposed test with two benchmarks, namely homoskedastic and heteroskedasticity-robust Wald tests. Our tests consider the Bayesian Information Crietiron (BIC) and Generalized Information Criterion (GIC).
\begin{enumerate}
\item {\bf Wald-HE}: The heteroskedasticity-robust Wald test of \cite{caner2018asymptotically} is implemented using ordinary LASSO estimation. The tuning parameter selected for the initial regression is also used for the node-wise regression.
\item {\bf Wald-HO}: The homoskedastic Wald test is implemented using the same approximate inverse as Wald-HE and the scaled LASSO estimator of \cite{sun2012scaled} for the initial regression following \cite{van-de-geer(2014)}.
\item {\bf LMTW}: Our proposed test defined in \eqref{eq:rejection-criteria} is implemented using ordinary LASSO estimation for both the initial and node-wise regressions. Tuning parameters are seperately selected using BIC and GIC.
\item {\bf LMTW-CL}: Our proposed test defined in \eqref{eq:rejection-criteria} is implemented using conservative LASSO estimation described in Section \ref{subsec:conservative_LASSO} for both the initial and node-wise regressions. Tuning parameters are seperately selected using BIC and GIC.
\end{enumerate}
Implementation details such as information criteria, tuning parameter selection and procedures, conservative LASSO weights, and critical value computations are found in Appendix \ref{sec:appendix_implement}.
\subsection{Simulation Design}
We consider two dimensional settings: $(n,p)=(200,100)$ and $(n,p)=(200,300)$. In both settings, the number of nonzero elements of $\beta^0$ is fixed at $s_0=5$. Each nonzero coefficient is set equal to $1$, with their locations equally spaced across the coordinates of $\beta^0$.
We consider two restriction matrices. The first imposes restrictions on one nonzero and four zero elements of $\beta^0$, so that $r=m=5$. The second imposes restrictions on all five nonzero and twenty zero elements of $\beta^0$, so that $r=m=25$.
For the power analysis, we consider deviations from the null by testing against
\[
\widetilde q
=
q
+
\delta
\left(
\underbrace{\frac{1}{\sqrt r},\ldots,
\frac{1}{\sqrt r}}_{r}
\right)^\top,
\qquad
\delta\in\{0.2,0.4\}.
\]
We focus on the random-design setting, so that both ${\bf X}$ and $\varepsilon$ are independently resampled in each Monte Carlo replication. For the design matrix ${\bf X}$, we consider Gaussian and exponential designs, both standardized to have zero mean and unit variance, with covariance structure
\[
\Sigma_{ij}=(0.5)^{|i-j|},
\qquad i,j\in[p].
\]
For the errors, we consider Gaussian and Student's $t$ distributions with five degrees of freedom, both standardized to have zero mean and unit variance. We consider both homoskedastic and heteroskedastic specifications. Under heteroskedasticity, the conditional error variance is given by
\[
\sigma_i^2
=
\frac{0.5+|{\bf X}_{i,1}|}
{0.5+\mathbb{E}|{\bf X}_{i,1}|},
\]
so that $\mathbb{E}[\sigma_i^2]=1$, ensuring that the unconditional error variance remains equal to one, as in the homoskedastic specification.
Tables \ref{tab:size_power_n200_p100_r5_homo}-\ref{tab:size_power_n200_p100_r25_hetero} report the simulation results for $(n,p)=(200,100)$. The results for $(n,p)=(200,300)$ are reported in Appendix \ref{sec:additional_simulation}, and the pattern is similar to the one we observe here. The proposed test statistic delivers substantially better size control than the homoskedastic Wald benchmark under heteroskedasticity, especially when the design is non-Gaussian. Under homoskedasticity and small $r$, the homoskedastic Wald test can be competitive, but its size deteriorates more severely as the number of restrictions increases. For the heteroskedasticity-robust Wald test, the size control is competitive with our test statistic when $r$ is small under both homoskedasticity and heteroskedasticity. However, when $r$ is large, it has even more severe size distortion than the homoskedastic Wald test. This suggests that estimating and inverting the full sandwich covariance matrix can be very unstable in this regime. We also note that the use of conservative LASSO in our test statistic can improve size control when $r$ is large. In addition, our test statistic with ordinary LASSO has better power than the homoskedastic Wald test under homoskedasticity and small $r$ where the sizes of all tests are close to each other. Our proposed test with conservative LASSO exhibits higher rejection rates under alternatives than under ordinary LASSO, although these rates remain below those of the heteroskedasticity-robust Wald test under both homoskedasticity and heteroskedasticity.\footnote{Note, however, that the heteroskedasticity-robust Wald test exhibits greater size distortions when $r$ is large. Rejection rates under alternatives are not adjusted for differences in empirical size and therefore should not be interpreted as size-adjusted power comparisons}. Finally, note that the power of our test statistic is slightly better when BIC is used to select the tuning parameters. Given our emphasis on size control, our proposed test obtains smaller size distortions than the benchmark Wald tests in heteroskedastic designs with large $r$, with the conservative LASSO implementation generally reducing distortions further.
\begin{table}[H]
\centering
\begin{threeparttable}
\caption{Empirical Size and Power: $(n,p,r)=(200,100,5)$, Homoskedastic Errors}
\label{tab:size_power_n200_p100_r5_homo}
\scriptsize
\setlength{\tabcolsep}{2.7pt}
\begin{tabular}{lcccccccccccc}
\toprule\toprule
& \multicolumn{3}{c}{D1} & \multicolumn{3}{c}{D2} & \multicolumn{3}{c}{D3} & \multicolumn{3}{c}{D4} \\
\cmidrule(lr){2-4}\cmidrule(lr){5-7}\cmidrule(lr){8-10}\cmidrule(lr){11-13}
& $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ \\
\midrule
Wald-HO (BIC)
& 0.066 & 0.291 & 0.953 & 0.062 & 0.309 & 0.950 & 0.056 & 0.327 & 0.950 & 0.067 & 0.329 & 0.923 \\
Wald-HO (GIC)
& 0.071 & 0.266 & 0.934 & 0.070 & 0.278 & 0.935 & 0.060 & 0.299 & 0.941 & 0.074 & 0.305 & 0.907 \\
\addlinespace
Wald-HE (BIC)
& 0.059 & 0.503 & 0.992 & 0.051 & 0.554 & 0.994 & 0.052 & 0.564 & 0.992 & 0.065 & 0.593 & 0.993 \\
Wald-HE (GIC)
& 0.061 & 0.476 & 0.986 & 0.050 & 0.522 & 0.991 & 0.055 & 0.530 & 0.989 & 0.067 & 0.555 & 0.987 \\
\addlinespace
LMTW (BIC)
& 0.060 & 0.388 & 0.976 & 0.059 & 0.385 & 0.979 & 0.060 & 0.386 & 0.976 & 0.062 & 0.393 & 0.958 \\
LMTW (GIC)
& 0.069 & 0.342 & 0.970 & 0.062 & 0.337 & 0.961 & 0.061 & 0.349 & 0.967 & 0.062 & 0.349 & 0.939 \\
\addlinespace
LMTW-CL (BIC)
& 0.056 & 0.443 & 0.991 & 0.043 & 0.467 & 0.993 & 0.046 & 0.475 & 0.986 & 0.056 & 0.465 & 0.981 \\
LMTW-CL (GIC)
& 0.056 & 0.417 & 0.986 & 0.043 & 0.449 & 0.989 & 0.046 & 0.457 & 0.984 & 0.058 & 0.446 & 0.975 \\
\bottomrule
\end{tabular}
\begin{tablenotes}[flushleft]
\scriptsize
\item \textit{Notes:} The nominal size is $\alpha=0.05$. The column $\delta=0$ reports empirical size, while $\delta=0.2$ and $\delta=0.4$ report empirical power. D1: Gaussian $X$ and Gaussian errors; D2: centered-exponential $X$ and Gaussian errors; D3: Gaussian $X$ and standardized $t_5$ errors; D4: centered-exponential $X$ and standardized $t_5$ errors. The errors are homoskedastic in all four designs. Wald-HO is the homoskedastic Wald test; Wald-HE is the heteroskedasticity-robust Wald test using conservative-LASSO; LMTW and LMTW-CL use ordinary/plain- and conservative-LASSO respectively.
\end{tablenotes}
\end{threeparttable}
\end{table}
\begin{table}[H]
\centering
\begin{threeparttable}
\caption{Empirical Size and Power: $(n,p,r)=(200,100,5)$, Heteroskedastic Errors}
\label{tab:size_power_n200_p100_r5_hetero}
\scriptsize
\setlength{\tabcolsep}{2.7pt}
\begin{tabular}{lcccccccccccc}
\toprule\toprule
& \multicolumn{3}{c}{D1} & \multicolumn{3}{c}{D2} & \multicolumn{3}{c}{D3} & \multicolumn{3}{c}{D4} \\
\cmidrule(lr){2-4}\cmidrule(lr){5-7}\cmidrule(lr){8-10}\cmidrule(lr){11-13}
& $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ \\
\midrule
Wald-HO (BIC)
& 0.118 & 0.373 & 0.936 & 0.231 & 0.470 & 0.951 & 0.100 & 0.385 & 0.945 & 0.246 & 0.504 & 0.926 \\
Wald-HO (GIC)
& 0.124 & 0.345 & 0.915 & 0.247 & 0.451 & 0.936 & 0.106 & 0.355 & 0.924 & 0.260 & 0.482 & 0.907 \\
\addlinespace
Wald-HE (BIC)
& 0.064 & 0.458 & 0.973 & 0.095 & 0.540 & 0.974 & 0.067 & 0.514 & 0.979 & 0.107 & 0.592 & 0.970 \\
Wald-HE (GIC)
& 0.065 & 0.433 & 0.964 & 0.099 & 0.505 & 0.965 & 0.065 & 0.486 & 0.967 & 0.114 & 0.566 & 0.963 \\
\addlinespace
LMTW (BIC)
& 0.062 & 0.307 & 0.934 & 0.083 & 0.287 & 0.835 & 0.056 & 0.333 & 0.934 & 0.091 & 0.313 & 0.821 \\
LMTW (GIC)
& 0.060 & 0.277 & 0.889 & 0.079 & 0.243 & 0.798 & 0.056 & 0.285 & 0.904 & 0.090 & 0.287 & 0.783 \\
\addlinespace
LMTW-CL (BIC)
& 0.051 & 0.355 & 0.952 & 0.087 & 0.314 & 0.872 & 0.045 & 0.392 & 0.955 & 0.088 & 0.331 & 0.838 \\
LMTW-CL (GIC)
& 0.050 & 0.341 & 0.943 & 0.086 & 0.304 & 0.867 & 0.045 & 0.373 & 0.948 & 0.088 & 0.323 & 0.825 \\
\bottomrule
\end{tabular}
\begin{tablenotes}[flushleft]
\scriptsize
\item \textit{Notes:} The nominal size is $\alpha=0.05$. The column $\delta=0$ reports empirical size, while $\delta=0.2$ and $\delta=0.4$ report empirical power. D1: Gaussian $X$ and Gaussian errors; D2: centered-exponential $X$ and Gaussian errors; D3: Gaussian $X$ and standardized $t_5$ errors; D4: centered-exponential $X$ and standardized $t_5$ errors. The errors are heteroskedastic in all four designs. Wald-HO is the homoskedastic Wald test; Wald-HE is the heteroskedasticity-robust Wald test using conservative-LASSO; LMTW and LMTW-CL use ordinary/plain-- and conservative-LASSO respectively.
\end{tablenotes}
\end{threeparttable}
\end{table}
\begin{table}[H]
\centering
\begin{threeparttable}
\caption{Empirical Size and Power: $(n,p,r)=(200,100,25)$, Homoskedastic Errors}
\label{tab:size_power_n200_p100_r25_homo}
\scriptsize
\setlength{\tabcolsep}{2.7pt}
\begin{tabular}{lcccccccccccc}
\toprule\toprule
& \multicolumn{3}{c}{D1} & \multicolumn{3}{c}{D2} & \multicolumn{3}{c}{D3} & \multicolumn{3}{c}{D4} \\
\cmidrule(lr){2-4}\cmidrule(lr){5-7}\cmidrule(lr){8-10}\cmidrule(lr){11-13}
& $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ \\
\midrule
Wald-HO (BIC)
& 0.147 & 0.267 & 0.734 & 0.149 & 0.241 & 0.772 & 0.137 & 0.248 & 0.736 & 0.137 & 0.256 & 0.764 \\
Wald-HO (GIC)
& 0.172 & 0.271 & 0.712 & 0.173 & 0.245 & 0.737 & 0.162 & 0.261 & 0.699 & 0.164 & 0.266 & 0.730 \\
\addlinespace
Wald-HE (BIC)
& 0.199 & 0.521 & 0.976 & 0.192 & 0.551 & 0.975 & 0.185 & 0.581 & 0.976 & 0.177 & 0.590 & 0.990 \\
Wald-HE (GIC)
& 0.214 & 0.510 & 0.963 & 0.202 & 0.526 & 0.970 & 0.195 & 0.571 & 0.973 & 0.192 & 0.567 & 0.988 \\
\addlinespace
LMTW (BIC)
& 0.092 & 0.202 & 0.753 & 0.089 & 0.171 & 0.748 & 0.082 & 0.189 & 0.721 & 0.071 & 0.164 & 0.726 \\
LMTW (GIC)
& 0.096 & 0.177 & 0.695 & 0.087 & 0.151 & 0.675 & 0.086 & 0.164 & 0.649 & 0.069 & 0.137 & 0.642 \\
\addlinespace
LMTW-CL BIC
& 0.071 & 0.209 & 0.840 & 0.069 & 0.193 & 0.841 & 0.056 & 0.184 & 0.805 & 0.050 & 0.171 & 0.807 \\
LMTW-CL (GIC)
& 0.069 & 0.198 & 0.818 & 0.071 & 0.179 & 0.822 & 0.058 & 0.179 & 0.777 & 0.051 & 0.167 & 0.784 \\
\bottomrule
\end{tabular}
\begin{tablenotes}[flushleft]
\scriptsize
\item \textit{Notes:} The nominal size is $\alpha=0.05$. The column $\delta=0$ reports empirical size, while $\delta=0.2$ and $\delta=0.4$ report empirical power. D1: Gaussian $X$ and Gaussian errors; D2: centered-exponential $X$ and Gaussian errors; D3: Gaussian $X$ and standardized $t_5$ errors; D4: centered-exponential $X$ and standardized $t_5$ errors. The errors are homoskedastic in all four designs. Wald-HO is the homoskedastic Wald test; Wald-HE is the heteroskedasticity-robust Wald test using conservative-LASSO; LMTW and LMTW-CL use ordinary/plain-- and conservative-LASSO respectively.
\end{tablenotes}
\end{threeparttable}
\end{table}
\begin{table}[H]
\centering
\begin{threeparttable}
\caption{Empirical Size and Power: $(n,p,r)=(200,100,25)$, Heteroskedastic Errors}
\label{tab:size_power_n200_p100_r25_hetero}
\scriptsize
\setlength{\tabcolsep}{2.7pt}
\begin{tabular}{lcccccccccccc}
\toprule\toprule
& \multicolumn{3}{c}{D1} & \multicolumn{3}{c}{D2} & \multicolumn{3}{c}{D3} & \multicolumn{3}{c}{D4} \\
\cmidrule(lr){2-4}\cmidrule(lr){5-7}\cmidrule(lr){8-10}\cmidrule(lr){11-13}
& $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ & $\delta=0$ & $\delta=0.2$ & $\delta=0.4$ \\
\midrule
Wald-HO (BIC)
& 0.184 & 0.291 & 0.778 & 0.244 & 0.361 & 0.831 & 0.168 & 0.278 & 0.758 & 0.237 & 0.364 & 0.834 \\
Wald-HO (GIC)
& 0.213 & 0.302 & 0.747 & 0.265 & 0.364 & 0.791 & 0.194 & 0.284 & 0.732 & 0.268 & 0.367 & 0.815 \\
\addlinespace
Wald-HE (BIC)
& 0.215 & 0.612 & 0.984 & 0.252 & 0.686 & 0.993 & 0.209 & 0.651 & 0.992 & 0.249 & 0.742 & 0.997 \\
Wald-HE (GIC)
& 0.221 & 0.597 & 0.981 & 0.264 & 0.666 & 0.989 & 0.221 & 0.632 & 0.983 & 0.270 & 0.721 & 0.994 \\
\addlinespace
LMTW (BIC)
& 0.087 & 0.168 & 0.722 & 0.096 & 0.183 & 0.668 & 0.074 & 0.164 & 0.676 & 0.066 & 0.158 & 0.641 \\
LMTW (GIC)
& 0.091 & 0.152 & 0.653 & 0.097 & 0.156 & 0.591 & 0.073 & 0.142 & 0.614 & 0.067 & 0.138 & 0.576 \\
\addlinespace
LMTW-CL (BIC)
& 0.069 & 0.184 & 0.786 & 0.083 & 0.194 & 0.734 & 0.043 & 0.151 & 0.754 & 0.053 & 0.157 & 0.709 \\
LMTW-CL (GIC)
& 0.070 & 0.175 & 0.764 & 0.079 & 0.188 & 0.703 & 0.045 & 0.145 & 0.739 & 0.058 & 0.144 & 0.682 \\
\bottomrule
\end{tabular}
\begin{tablenotes}[flushleft]
\scriptsize
\item \textit{Notes:} The nominal size is $\alpha=0.05$. The column $\delta=0$ reports empirical size, while $\delta=0.2$ and $\delta=0.4$ report empirical power. D1: Gaussian $X$ and Gaussian errors; D2: centered-exponential $X$ and Gaussian errors; D3: Gaussian $X$ and standardized $t_5$ errors; D4: centered-exponential $X$ and standardized $t_5$ errors. The errors are heteroskedastic in all four designs. Wald-HO is the homoskedastic Wald test; Wald-HE is the heteroskedasticity-robust Wald test using conservative-LASSO; LMTW and LMTW-CL use ordinary/plain-- and conservative-LASSO respectively.
\end{tablenotes}
\end{threeparttable}
\end{table}
\section{Empirical Application}
\label{sec:empirical_application}
We illustrate the proposed procedure using a high-dimensional version of the Romer-Romer regression for identifying monetary policy shocks. The empirical problem is to explain changes in the federal funds rate using information available to the Federal Reserve at each Federal Open Market Committee (FOMC) meeting. Omitting information contained in the Fed's own forecasts can confound the regression residuals used to measure policy shocks. The ``narrative identification approach'' of \citet{romer1989does, romer2004new} thus regresses changes in federal funds rate on its prevailing level and Fed's real-time forecasts of macroeconomic conditions, with the residual interpreted as the monetary policy shock.
\cite{aruoba2024identifying} augment these controls with sentiment indicators constructed by scoring economic-concept keywords in the FOMC documents.\footnote{The text consists of Greenbook Parts I and II (merged and renamed Tealbook A in 2010) and the Redbook (renamed Beige Book in 1983). Greenbook Part I summarizes recent economic conditions and presents the Fed staff's economic forecasts. Greenbook part II discusses recent economic and financial developments, whereas the Redbook summarizes economic conditions across the 12 Federal Reserve Districts.} Their data therefore provide a natural high-dimensional setting for illustrating joint inference on the coefficients associated with economically meaningful groups of text-based predictors. We analyze their dataset, which covers $n=210$ FOMC meetings from October 5 1982 to October 29 2008, and contains sentiment indicators for the 296 selected economic-concept keywords.\footnote{The sentiment data is available on the authors' website: \href{https://thomas-drechsel.com/research.html}{https://thomas-drechsel.com/research.html}. The sample terminates in 2008 to avoid the zero-lower-bound period.}
Our benchmark regression is
\begin{align}
\Delta i_t
=
\alpha + \gamma i_{t-1} + \delta^\top N_t + \zeta^\top S_t + \varepsilon_t,
\qquad t=1,\ldots, 210,
\label{RR_reg}
\end{align}
where $i_{t-1}$ denotes the level of the federal funds rate prior to the $t$-th FOMC announcement, $\Delta i_t$ denotes the Fed's policy decision, $N_t$ represents the 132 numerical forecasts in the FOMC documents, $S_t$ is the vector of 296 sentiment indicators, and $\varepsilon_t$ is conceptually the monetary policy shock. Results for $p=858$, constructed by additionally including quadratic terms of $N_t$ and $S_t$, are reported in Appendix \ref{ERR_RC}. Let $\mathbf{X}_t=(1,i_{t-1},N_t^\top, S_t^\top)^\top$, $\mathbf{X}=(\mathbf{X}_1,\ldots,\mathbf{X}_n)^\top$, $\beta^0 = (\alpha, \gamma, \delta^\top, \zeta^\top)^\top$. We assume that $\beta^0$ is sparse, and let $\Delta i = (\Delta i_1, \ldots, \Delta i_n)^\top$, and $\varepsilon = (\varepsilon_1, \ldots, \varepsilon_n)^\top$. Then the regression (\ref{RR_reg}) can be written in the form (\ref{eq:model}):
\begin{align*}
\Delta i = \mathbf{X} \beta^0 + \varepsilon, \qquad \text{with } (n,p) = (210, 430).
\end{align*}
We emphasize that our objective is not to use the high-dimensional linear regression to identify monetary policy shocks, but rather to conduct inference on its regression coefficients under the sparsity assumption. This extends the regression to a high-dimensional setting, from the conventional low-dimensional setting by \cite{romer2004new} and \cite{croushore2018fiscal}, and examine joint hypotheses concerning its coefficients. The error term $\varepsilon_t$ in \eqref{RR_reg}, interpreted as the monetary policy shocks, is conceptually treated as serially uncorrelated (\cite{bauer2023reassessment}). For the purpose of applying our procedure, we assume independence across FOMC meetings, which is stronger than serial uncorrelatedness. Note that substantial heteroskedasticity in $\varepsilon_t$ has been documented by \citet{coibion2012effects}; see also \citet{schlaak2023monetary}. These features provide a natural setting for applying heteroskedasticity-robust high-dimensional inference.
To provide empirical context for the regularity conditions underlying our test, Appendix \ref{cali_ERR} reports tail diagnostics for the standardized covariates and moment diagnostics based on the LASSO residuals. We adopt $(\nu,h)=(1,4)$ as conservative working values, corresponding to marginal sub-exponential covariates and the minimum fourth-moment requirement for the errors. These choices avoid imposing stronger distributional assumptions, such as sub-Gaussian covariates or higher-order error moments. We emphasize that the finite-sample diagnostics do not verify the population assumptions or the asymptotic rate conditions in Assumption \ref{ass:size}. They instead provide descriptive evidence to accompany our empirical application, with the robustness of the findings further examined in Appendix \ref{ERR_RC}.
Let $G$ be one of the six clusters of economic concepts: ``oil'', ``inflation'', ``labor'', ``credit'', ``demand'', and ``housing'', with the constituent keywords of each cluster detailed in Table \ref{tab:lag0-linear-main-concepts}. We apply our test to examine which clusters of sentiment indicators are jointly significant for the Fed's monetary policy decisions. For example, when $G = \text{``oil''}$, the parameter of interest is the subvector $\beta^0_G$ consisting of the regression coefficients associated with the sentiment indicators for the keywords ``crude oil'', ``oil imports'', ``oil prices'', and ``petroleum'', so that $m = |G| = 4$. Recalling the Introduction, our testing problem of interest in subvector form is
\begin{align*}
R_G \beta^0_{G} = 0
\qquad vs. \qquad
R_G \beta^0_{G} \neq 0.
\end{align*}
Equivalently, the null hypothesis can be written as $R\beta^0=0$, where $R$ is obtained by embedding $R_G$ into an $\mathbb{R}^{r \times p}$ matrix and padding the columns corresponding to coefficients outside $G$ with zeros. Our proposed test can then be applied directly. We consider two choices of the restriction matrix:
\begin{align*}
R_G
= I_m \qquad \text{or} \qquad R_G = [ - \mathbb{1}_{m-1}, I_{m-1}]
\end{align*}
where $\mathbb{1}_{m-1}$ is the vector of ones of dimension $m-1$. The first choice, $R_G=I_m$, tests the joint null that all coefficients associated with cluster $G$ are equal to zero. The second choice, $R_G=[-\mathbb{1}_{m-1},I_{m-1}]$, tests whether all coefficients within the cluster are equal, thereby examining whether the sentiment indicators within a given cluster have homogeneous effects on the monetary policy outcome.
\begin{table}[!t]
\centering
\footnotesize
\setlength{\tabcolsep}{5.0pt}
\begin{threeparttable}
\caption{Economic Concept Keywords underlying Sentiment Indicators in Each Cluster $G$}
\label{tab:lag0-linear-main-concepts}
\begin{tabular}{@{}>{\raggedright\arraybackslash}p{1.65cm}>{\raggedright\arraybackslash}p{14.25cm}@{}}
\toprule
Group $G$ & Related concepts \\
\midrule
credit & ``banks'', ``credit'', ``credit standards'', ``financing'', ``lending'', ``liquidity'', ``loan demand'', ``loans'' \\
demand & ``aggregate demand'', ``business activity'', ``consumption'', ``domestic demand'', ``economic growth'', ``final sales'', ``industrial production'', ``investment'', ``manufacturing activity'', ``output gap'', ``retail sales'' \\
housing & ``home prices'', ``home sales'', ``household debt'', ``household net worth'', ``housing'', ``housing activity'', ``housing starts'', ``mortgage'', ``residential investment'' \\
inflation & ``consumer prices'', ``core inflation'', ``headline inflation'', ``inflation'', ``inflation expectations'', ``inflationary pressures'', ``price pressures'', ``producer prices'', ``unit labor cost'', ``wage pressures'' \\
labor & ``employment'', ``jobs'', ``labor force participation'', ``labor market'', ``payrolls'', ``unemployment'', ``wage pressures'', ``wages'', ``workers'' \\
oil & ``crude oil'', ``oil imports'', ``oil prices'', ``petroleum'' \\
\bottomrule
\end{tabular}
\begin{tablenotes}[flushleft]
\scriptsize
\item Notes: The concepts are from \cite{aruoba2024identifying}; they are single-, double-, or triple-word expressions extracted from FOMC documents and scored by a dictionary-based procedure.
\end{tablenotes}
\end{threeparttable}
\end{table}
Table \ref{tab:lag0-linear-main-p429} summarizes our empirical findings. The first $p$-value column shows that the joint null of zero coefficients is rejected at the 5\% level for the credit, demand, and housing clusters across both ordinary/plain and conservative LASSO implementations and under both tuning criteria. These results indicate that sentiment indicators in these clusters are jointly associated with monetary policy decisions conditional on the remaining controls. The evidence is particularly strong for demand, with the 1\% significance level in all four implementations. Viewed in the context of monetary-policy history, this pattern is consistent with the sample period, which spans Volcker's disinflation legacy, Greenspan's risk-management approach, and Bernanke's response to mounting housing and credit stress. Against this backdrop, the joint significance of demand, credit, and housing suggests that qualitative information on real activity and financial conditions contains incremental information beyond that captured by the numerical staff forecasts $N_t$. By contrast, the insignificant results for concept groupings, inflation and oil, need not imply that these factors were unimportant for policy; a possibility is that if they were already central to the Fed's forecasting framework, then less incremental information would remain in the associated text. Labor presents a less robust pattern: BIC selects ``jobs'' but yields no joint significance, whereas GIC selects no individual labor concept yet strongly rejects the joint null. This contrast highlights the distinction between variable selection and joint inference, suggesting that conclusions regarding the incremental information in labor market sentiment are sensitive to the chosen regularization criterion. Importantly, these findings illustrate how joint inference can reveal statistically significant associations within groups of predictors even when individual variables are not selected by the initial penalized estimator.
The coefficient-equality tests in the second $p$-value column complement this analysis by examining whether the coefficients associated with sentiment indicators are equal within each cluster. For demand, the equality null is rejected at the 5\% level in all four implementations. Thus, not only is demand sentiment jointly relevant, but different dimensions of demand appear to carry different policy implications. Credit shows a similar, though less robust, pattern, with equality rejected at the 5\% level under GIC but not under BIC. Housing is more distinctive. The cluster is jointly significant in all four implementations, yet coefficient equality is not rejected. This finding is consistent with the possibility that the housing indicators capture a common underlying housing market signal, although the failure to reject equality in coefficients does not establish that these effects are identical.
More generally, the two tests provide complementary insights into cluster-level relevance and within-cluster heterogeneity. When both the joint-zero and coefficient-equality null hypotheses are rejected, as in the case of demand, the results suggest that the constituent sentiment concepts are jointly relevant but exert heterogeneous effects, highlighting the importance
of their composition for policy interpretation.
In contrast, joint significance accompanied by non-rejection of coefficient equality, as observed for housing, is consistent with a more parsimonious common-factor interpretation, in which the constituent sentiment concepts have similar estimated effects. For inflation and oil, the joint-zero null hypothesis is not rejected at the 5\% significance level, limiting the additional
economic insights provided by the coefficient-equality tests. Labor presents a more nuanced case: under GIC, the results indicate joint significance without detectable coefficient heterogeneity, resembling the findings for housing. However, this conclusion is sensitive to the choice of regularization criterion.
\begin{table}[H]
\centering
\footnotesize
\setlength{\tabcolsep}{3.0pt}
\begin{threeparttable}
\caption{Concept-group tests of zero coefficients and coefficient equality, $(n,p)=(210,430)$}
\label{tab:lag0-linear-main-p429}
\begin{tabular*}{\textwidth}{@{\extracolsep{\fill}}lcllll@{}}
\toprule
& & & & \multicolumn{2}{c}{$p$-value} \\
\cmidrule(l){5-6}
Group $G$ & $m$ & Method & Selected concepts
& $R_G=I_m$
& $R_G=[-\mathbb{1}_{m-1},\,I_{m-1}]$ \\
\midrule
credit & 8 & LMTW (BIC) & ``banks'', ``lending'' & $0.0271^{**}$ & $0.0716^{*}$ \\
& & LMTW (GIC) & ``banks'' & $<0.001^{***}$ & $0.0217^{**}$ \\
& & LMTW-CL (BIC) & ``banks'', ``lending'' & $0.0283^{**}$ & $0.0750^{*}$ \\
& & LMTW-CL (GIC) & ``banks'' & $<0.001^{***}$ & $0.0283^{**}$ \\
\midrule
demand & 11 & LMTW (BIC) & ``business activity'', ``investment'' & $<0.001^{***}$ & $0.0123^{**}$ \\
& & LMTW (GIC) & ``business activity'' & $<0.001^{***}$ & $0.0026^{***}$ \\
& & LMTW-CL (BIC) & ``business activity'', ``investment'' & $<0.001^{***}$ & $0.0123^{**}$ \\
& & LMTW-CL (GIC) & ``business activity'' & $<0.001^{***}$ & $0.0026^{***}$ \\
\midrule
housing & 9 & LMTW (BIC) & ``household net worth'' & $0.0105^{**}$ & $0.4638$ \\
& & LMTW (GIC) & None & $<0.001^{***}$ & $0.6572$ \\
& & LMTW-CL (BIC) & ``household net worth'' & $0.0105^{**}$ & $0.4638$ \\
& & LMTW-CL (GIC) & None & $<0.001^{***}$ & $0.6572$ \\
\midrule
inflation & 10 & LMTW (BIC) & None & $0.5488$ & $0.5927$ \\
& & LMTW (GIC) & None & $0.2133$ & $0.1805$ \\
& & LMTW-CL (BIC) & None & $0.5488$ & $0.5927$ \\
& & LMTW-CL (GIC) & None & $0.2133$ & $0.1805$ \\
\midrule
labor & 9 & LMTW (BIC) & ``jobs'' & $0.1673$ & $0.1720$ \\
& & LMTW (GIC) & None & $0.0025^{***}$ & $0.1475$ \\
& & LMTW-CL (BIC) & ``jobs'' & $0.1673$ & $0.1720$ \\
& & LMTW-CL (GIC) & None & $0.0025^{***}$ & $0.1475$ \\
\midrule
oil & 4 & LMTW (BIC) & None & $0.3089$ & $0.1457$ \\
& & LMTW (GIC) & None & $0.2420$ & $0.0592^{*}$ \\
& & LMTW-CL (BIC) & None & $0.2998$ & $0.1592$ \\
& & LMTW-CL (GIC) & None & $0.2229$ & $0.0623^{*}$ \\
\bottomrule
\end{tabular*}
\begin{tablenotes}[flushleft]
\scriptsize
\item Notes: For coefficient-equality tests $R_G = [-\mathbb{1}_{m-1}, I_{m-1}]$, the reference is the first
concept listed for its group in Table~\ref{tab:lag0-linear-main-concepts}. LMTW and LMTW-CL denote the proposed test implemented with ordinary and conservative LASSO, respectively. BIC and GIC denote the information criteria used to select the penalty parameters in the main and node-wise regressions. ``Selected concepts'' lists the concept selected by the main LASSO within each group. Simulated $p$-values and critical values use $\widetilde\Phi$ and 200,000 weighted chi-square draws.
Stars denote rejection at the 10\% ($^{*}$), 5\% ($^{**}$), and 1\% ($^{***}$) levels.
\end{tablenotes}
\end{threeparttable}
\end{table}
\section{Conclusion}
\label{sec:conclusion}
We propose a novel test for multiple restrictions in high-dimensional linear regression models, allowing for general heteroskedasticity and heavy-tailed covariates and errors. The method combines desparsified estimation with a ridge-regularized quadratic-form statistic and a critical value construction that remains valid under both fixed and diverging restriction dimensions. Our theoretical framework accommodates independent but non-identically distributed covariate rows with marginal sub-Weibull tails and heteroskedastic errors satisfying finite polynomial moment conditions. Simulation results demonstrate substantially improved size control relative to conventional homoskedastic and heteroskedasticity-robust Wald tests, particularly when the number of restrictions is large. Our empirical application illustrates joint inference for groups of text-based predictors in a high-dimensional monetary policy regression.
Several directions for future research remain. One possibility is to extend the proposed framework to generalized linear and nonparametric models. Another is to weaken the conditions relating the effective degrees of freedom of the projection matrix to its rank. Developing alternative critical value approximations, potentially based on resampling methods, may provide a way to relax these conditions, although establishing their validity in high-dimensional settings would require further theoretical investigation.
Overall, our results provide a unified framework for heteroskedasticity-robust joint inference in sparse high-dimensional models, without requiring practitioners to distinguish between fixed and diverging restriction dimensions when implementing the proposed procedure.
\clearpage
\bibliographystyle{chicago}
\singlespacing
\bibliography{reference}
\newpage