EconBase
← Back to paper

Estimator Averaging of Local Projection and VAR Impulse Responses

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

80,436 characters · 19 sections · 35 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Estimator Averaging of Local Projection and VAR Impulse Responses

abstractLocal projections (LP) and vector autoregressions (VAR) are the two standard tools for impulse response analysis, but they often display a finite-sample trade-off: LP is typically less biased but more volatile, while VAR is more precise but can be biased under misspecification. We propose an easy-to-implement estimator-averaging approach that combines LP and VAR at each horizon by minimizing the mean squared error of the impulse response itself, rather than in-sample fit. We derive closed-form oracle weights for this finite-sample risk problem, develop feasible AR-sieve-bootstrap procedures, and compare them against an \(R^2\)-based model-averaging benchmark. For a benchmark class of short-memory linear data generating processes in which LP and VAR are both consistent, we establish the consistency and limiting distribution of the feasible averaged estimator. Monte Carlo results show meaningful risk reductions relative to LP and VAR alone. In an empirical application revisiting BauerSwanson2023, estimator averaging delivers stable and economically intuitive responses for yields, activity, prices, and credit spreads. JEL Classification: C32, C52, C53, E52 \\ Keywords: Local projections; Vector autoregressions; Impulse response functions; Estimator averaging; Model averaging; Monetary policy shocks

\setcounter{page}{1}\parskip0.5em \baselineskip18pt \numberwithin{equation}{section} \doublespacing

Introduction

Local projections (LP) and vector autoregressions (VAR) are the two workhorse approaches to estimating impulse response functions (IRF). In practice, however, the two methods often deliver noticeably different estimates, especially at intermediate and long horizons. This creates a natural problem for empirical researchers: if LP and VAR provide conflicting answers, how should one combine the information in the two estimators?

A useful starting point is that the distinction is not primarily about the population object being estimated. PlagborgMollerWolf2021 show that, with sufficiently rich lag structure, LPs and VARs estimate the same population IRFs. The key difference is therefore one of finite-sample behavior. As emphasized in recent work, LP and VAR exhibit a familiar finite-sample bias--variance trade-off (see, e.g., Lietal2024; MontielOleaetal2025). LPs tend to have lower bias but higher variance, especially at intermediate and long horizons, because they estimate separate horizon-specific regressions. VARs, by imposing parametric dynamic structure, typically achieve lower variance by borrowing strength across horizons, but they can incur higher bias when the DGP deviates from a finite-order VAR (Lietal2024).

This makes combining LP and VAR naturally appealing from a mean-squared-error perspective: an averaged estimator can potentially exploit LP's low-bias properties and VAR's low-variance properties, provided the weights adapt to horizon-specific performance. MontielOleaEtAl2026 sharpen this trade-off in a local-misspecification framework. They show that LP remains first-order robust for inference, whereas VAR can suffer first-order bias that is relevant for coverage. This highlights a practical tension: LP is attractive when robustness to misspecification is the priority, while VAR is attractive when finite-sample precision is the priority.

Existing work on combining IRF estimators has followed several related routes. One strand takes a fit-oriented model-averaging approach. For example, HounyoJung2025 propose a two-stage scheme that averages within LPs and within VARs, and then blends the two classes. Such procedures often yield smooth, interpretable IRFs and stable weights when many horizons are pooled. However, because the objective is in-sample fit, they do not directly target the estimation risk of the structural IRF. A related but distinct approach is developed by NemtyrevBoldea2026, who propose Targeted Local Projections (TLP), a shrinkage estimator that pulls LP impulse responses toward their SVAR counterparts in order to reduce variance at the cost of some bias. Their framework is developed explicitly under a local-misspecification asymptotic setup and is complemented by bootstrap-based inference designed to improve coverage in that setting.

Motivated by these distinctions, we take an estimator-averaging approach that selects weights to minimize the expected error of the IRF itself. Specifically, for each horizon $h$, we choose $w_h$ to minimize the population or estimated risk---variance or MSE---of the convex combination \[ \hat{\theta}_h(w_h) = w_h\,\hat{\theta}_{LP,h}+(1-w_h)\,\hat{\theta}_{VAR,h}. \] Our approach is symmetric between LP and VAR and directly tied to the object of interest. Rather than selecting the best-fitting model (HounyoJung2025) or shrinking one estimator toward the other as a baseline (NemtyrevBoldea2026), we ask how to combine the two estimators so as to minimize IRF estimation risk. This prioritizes precision in the IRF itself, rather than in-sample fit, by choosing weights that directly reflect the LP--VAR bias--variance trade-off and by exploiting the covariance between LP and VAR to reduce sampling noise.

Relative to the existing literature, our contribution is threefold. First, in population, we derive closed-form finite-sample (infeasible) oracle weights that minimize the MSE of the combined estimator and make transparent how the optimal LP share varies with the horizon and with the underlying LP--VAR trade--off. Second, we develop the asymptotic theory for the AR-sieve plug-in averaged estimator under a benchmark short-memory linear DGP in which both LP and VAR are consistent for the same population impulse response. In that benchmark, the limiting risk is variance-based, and the bootstrap bias terms in our Algorithm 1 are a finite-sample refinement that is asymptotically negligible but empirically useful. This places our framework in deliberate complement to NemtyrevBoldea2026, who study a closely-related linear combination of LP and VAR estimators under a local-to-VAR drifting DGP and derive its asymptotic bias--variance trade-off in that regime. We use Monte Carlo simulations to assess its finite-sample performance, compare it with a simple $R^2$-based model-averaging benchmark, and show that our finite-sample MSE criterion and our AR-sieve-bootstrap implementation remain meaningful when the misspecification, if any, is fixed rather than drifting. Third, in an empirical application revisiting the high-frequency monetary policy shocks of BauerSwanson2023, we show that estimator averaging systematically reconciles the often volatile IV-LP and very smooth IV-VAR responses: the estimated weights put more mass on LP at short horizons and on VAR at longer horizons, and the resulting IRFs for the two-year yield, industrial production, consumer prices, and the excess bond premium are reasonably smooth, lie between LP and VAR, and are arguably more economically intuitive than either estimator on its own.

The rest of the paper is organized as follows. Section (ref) presents the estimator-averaging framework, derives the oracle weights, and outlines feasible implementations alongside the $R^2$-based model-averaging benchmark. Section (ref) presents the required assumptions and derives the large-sample properties of the estimator. Section (ref) reports Monte Carlo evidence on small-sample performance and the horizon-specific comparison between estimator and model averaging. Section (ref) revisits BauerSwanson2023 using our methods and documents the empirical pattern of weights and IRFs. Section (ref) concludes. All mathematical proofs are collected in the appendix.

Estimators

In this section, we review the definitions of LP- and VAR-based IRF estimators and introduce our LP--VAR averaged estimator. We derive the optimal weights by minimizing the MSE of the combined estimator---an approach known as estimator averaging (see, e.g., mittelhammer2005combining). We also compare our approach with model averaging, which selects weights by fitting an averaged model to maximize the $R^2$, following HounyoJung2025. The estimator-averaging framework presented here provides the core idea and can be extended to broader settings---for example, to averaging across multiple LP estimators or across multiple VAR estimators.

LP-Based IRF

The local projection method directly regresses the future outcome of a target variable, $y_t$, on a current impulse variable and a set of controls. Let $Y_t$ denote the $n \times 1$ vector of observed macroeconomic variables, which includes $y_t$. The horizon-$h$ LP regression is given by

equation[equation omitted — 124 chars of source]

where $c_h$ collects deterministic terms (e.g., a constant or trends), $x_t$ is the impulse variable of interest, and $z_t$ is a vector of control variables. Typically, $z_t$ consists of lagged values of the system variables, $z_t = (Y_{t-1}', Y_{t-2}', \dots, Y_{t-p}')'$. The LP estimator of the structural impulse response is the OLS (or 2SLS) estimate of the scalar coefficient on the impulse variable, denoted $\widehat\theta_{LP,h}$. A key advantage of this approach is that it requires no dynamic parametric assumptions beyond the linear projection itself.

remarkThe definition of $x_t$ depends on the identification strategy. If the structural shock is observed, $x_t$ is the shock itself. If unobserved, $x_t$ is typically an endogenous policy indicator. Because such an indicator may react contemporaneously to other shocks, merely including lagged controls in $z_t$ is insufficient. To isolate the structural shock, researchers generally rely on external proxies (IV-LP) or, when economically credible, a recursive identification scheme that adds the contemporaneous values of slow-moving variables to $z_t$.

VAR-Based IRF

Alternatively, the vector autoregression (VAR) approach (see, e.g., Sims1980) extrapolates the impulse response by iterating on a reduced-form linear dynamic system. A standard VAR($p$) model for the full vector of observables $Y_t$ is given by

equation[equation omitted — 84 chars of source]

where $c$ contains the deterministic terms, $A_l$ are the $n \times n$ matrices of autoregressive coefficients, and $v_t$ is the vector of reduced-form innovations. This system is typically estimated equation-by-equation via OLS. To recover the structural impulse responses, an identification scheme (e.g., recursive ordering via Cholesky decomposition, or external instruments/proxy SVARs) is specified to map the reduced-form innovations to the underlying structural shocks, $v_t = M_0 \epsilon_t$.

From the fitted dynamics and the identified structural impact matrix $M_0$, the implied vector moving average (VMA) representation is computed. The VAR-based estimator, $\widehat\theta_{VAR,h}$, is then defined as the relevant scalar element of the $h$-step-ahead structural VMA matrix $ \widehat\theta_{VAR,h} \;\equiv\; g_h\bigl(\widehat A_1, \dots, \widehat A_p, \widehat M_0\bigr)$ where $g_h(\cdot)$ is the non-linear function mapping the reduced-form VAR parameters and the identification matrix to the horizon-$h$ impulse response of $y_t$.

Combined Estimators

As is well documented, LPs and VARs exhibit a finite-sample bias–variance trade-off, motivating averaging to improve performance. There are two distinct routes: estimator averaging, which chooses weights to minimize the sampling risk (variance/MSE) of the IRF estimator, and model averaging, which chooses weights to optimize either in-sample or predictive fit. To represent both, we implement (i) a risk-based estimator-averaging scheme that selects horizon-specific weights by minimizing the MSE of the combined IRF, and (ii) a fit-based model averaging benchmark (e.g., HounyoJung2025). Our implementation of this risk-based estimator-averaging scheme encompasses two distinct methods: one that approximates the asymptotically optimal weight and a second that directly minimizes the MSE, both using a bootstrap procedure.

Estimator Averaging with MSE-Minimizing Weights

At each horizon $h$, define the averaged IRF

equation[equation omitted — 143 chars of source]

where $\widehat\theta_{LP,h}$ is the LP-based IRF estimate, $\widehat\theta_{VAR,h}$ is the VAR-based IRF estimate, and $w_h$ is the weight on LP at horizon $h$.

We decompose the (population) MSE of $\widehat\theta_h(w_h)$ in (ref) into variance plus squared bias:

align[align omitted — 405 chars of source]

where

align*[align* omitted — 291 chars of source]

Let $a_h=V_{L,h}+b_{L,h}^2$, $d_h=V_{V,h}+b_{V,h}^2$, and $f_h=C_h+b_{L,h}b_{V,h}$. Differentiating (ref) with respect to $w_h$ and solving the first-order condition, assuming an interior solution, yields the oracle weight

equation[equation omitted — 117 chars of source]

Because $w_h^{\star}$ depends on unknown population quantities, we estimate it using a semiparametric AR-sieve-bootstrap, following the time-series bootstrap literature (see, e.g., Buhlmann1997; KreissPaparoditis2003; GoncalvesKilian2004, Goncalves2007). Specifically, we fit an AR($p$) model with $p$ selected by BIC, resample the estimated residuals nonparametrically, simulate pseudo time series, and re-estimate the relevant objects to obtain $(\widehat a_h,\widehat d_h,\widehat f_h)$. This approach avoids committing to a fixed finite-order structural DGP and instead allows the bootstrap to approximate a broader class of short-memory linear processes through a growing-order autoregressive sieve, while still delivering a feasible estimator of $w_h^{\star}$. We refer to the procedure as the AR-sieve-bootstrap plug-in estimator, or simply the plug-in estimator. Algorithm (ref) summarizes the procedure.\footnote{For ease of exposition, Algorithm (ref) is written for a univariate model. It can easily be extended to a multivariate $VAR(p_T)$ sieve.}

The bootstrap bias terms in Algorithm (ref) are used for finite-sample implementation of the MSE criterion. They are not required for the benchmark first-order asymptotic theory in Section (ref). In the benchmark setting, both LP and VAR are asymptotically centered at \(\theta_h\), so the limiting risk is variance-covariance based and the bias contributions are first-order negligible. We therefore treat the bootstrap bias terms as a finite-sample refinement of a variance-based limit rather than as objects that the first-order theory is required to validate: they capture finite-sample features that disappear in the first-order limit but matter empirically. A higher-order justification of these terms would require additional smoothness and moment conditions on the AR-sieve approximation, which we do not impose.\footnote{A formal asymptotic justification for the bias terms used in Algorithm 1 lies outside the variance-based limit considered here; the contemporaneous and complementary work of NemtyrevBoldea2026 develops such a justification in a local-to-VAR drifting framework.}

The weight choice is a finite-sample risk problem. For fixed $T$, the local projection and VAR estimators may exhibit nonzero bias, so the MSE depends on both variance and bias components. In Section (ref), we derive our asymptotic results as $T\to\infty$. Working under the assumption of a short-memory linear process where the LP and VAR estimators are both consistent, these bias terms become asymptotically negligible, and the first-order limit distribution is centered at $\theta_h$. Within this same section, we subsequently establish the consistency of the bootstrap estimator, $\widehat{w}$, and the limiting theory for $\widehat{\theta}(\widehat{w})$.

algorithm[algorithm omitted — 2,384 chars of source]

Our MSE decomposition in (ref) is closely related to the risk criterion in NemtyrevBoldea2026. Both papers minimize the MSE of a linear combination of LP and VAR/SVAR impulse-response estimators, and both deliver weights that are in the form of a ratio in which the numerator captures a variance/covariance gap and the denominator adds a bias-squared term. Beyond this shared algebraic structure, however, the two papers are different both theoretically and methodologically, and we view the contributions as complementary rather than competing. NemtyrevBoldea2026 adopt the local-to-VAR DGP of MontielOleaEtAl2026 in which the misspecification of the VAR shrinks at rate $T^{-\zeta}$ for some $\zeta>1/4$. Under this drifting sequence, the VAR carries an $O(T^{-\zeta})$ asymptotic bias that is first-order relevant, and the weight is determined by a bias--variance trade-off in the asymptotic limit. We instead work under a fixed short-memory Wold DGP in which both LP and VAR are root-$T$ consistent for the same population impulse response. The bias--variance trade-off in our framework is therefore a finite-sample phenomenon, not an asymptotic one, and the limiting oracle weight is a pure variance/covariance ratio. Given their set up, NemtyrevBoldea2026 interpret the estimators combination as shrinkage of LP toward VAR, with the VAR playing the role of a regularization target; their estimator reduces to LP when shrinkage is off and to VAR in the opposite limit. We adopt a symmetric estimator-averaging perspective, treating LP and VAR as two estimators of the same scalar IRF and choosing weights that minimize the MSE of the combined estimator that can be easily extended to $K$ estimators (see Remark 2).

remarkThe estimator-averaging framework extends directly to more than two IRF estimators. Suppose $\widehat{\boldsymbol\theta}_h=(\hat\theta_h^{(1)},\ldots,\hat\theta_h^{(K)})^\top$ collects $K$ estimators of the same scalar IRF $\theta_h$, and define \[ \boldsymbol e_h=\widehat{\boldsymbol\theta}_h-\theta_h\mathbf{1}_K, \qquad \mathbf{b}_h=\mathbb{E}[\boldsymbol e_h], \qquad \boldsymbol\Sigma_h=\mathrm{Var}(\boldsymbol e_h). \] For weights $\mathbf{w}_h\in\mathbb{R}^K$ satisfying $\mathbf{1}_K^\top\mathbf{w}_h=1$, consider the averaged estimator \[ \tilde\theta_h(\mathbf{w}_h)=\mathbf{w}_h^\top \widehat{\boldsymbol\theta}_h. \] Its MSE is \[ \mathrm{MSE}[\tilde\theta_h] = \mathbf{w}_h^\top \big(\boldsymbol\Sigma_h+\mathbf{b}_h\mathbf{b}_h^\top\big) \mathbf{w}_h. \] Hence, minimizing MSE subject to $\mathbf{1}_K^\top\mathbf{w}_h=1$ yields the oracle weight vector \begin{equation} \mathbf{w}_h^\star = \frac{\mathbf{G}_h^{-1}\mathbf{1}_K} {\mathbf{1}_K^\top \mathbf{G}_h^{-1}\mathbf{1}_K}, \quad \mathbf{G}_h=\boldsymbol\Sigma_h+\mathbf{b}_h\mathbf{b}_h^\top. \end{equation} Thus, the same logic applies to averaging across multiple LP- and VAR-based IRF estimators. As before, $\mathbf{w}_h^\star$ depends on unknown population quantities and must be estimated in practice, for example by bootstrap methods.
remarkThe estimator-averaging framework is not tied to OLS estimation. The objects \(\widehat\theta_{LP,h}\) and \(\widehat\theta_{VAR,h}\) can be any two estimators of the same structural impulse response. Thus, under the usual relevance and exogeneity conditions for a valid external instrument, the framework can also combine IV-based estimators. In particular, one may combine a horizon-specific IV-LP estimator with a proxy-SVAR or IV-VAR estimator identified using the same external instrument: $\widehat\theta_h(w_h)=w_h\widehat\theta_{LP,h}^{IV}+(1-w_h)\widehat\theta_{VAR,h}^{IV}$. The MSE formula in (ref) then applies with the bias, variance, and covariance terms interpreted for the corresponding IV-based estimators. The feasible bootstrap implementation proceeds analogously: in each bootstrap sample, both the IV-LP and proxy-SVAR/IV-VAR estimators are re-estimated under the same identification scheme, and the resulting pair of IRF estimates is used to compute the bootstrap risk components.\footnote{The formal first-order theory in Section (ref) is stated for a generic pair of root-\(T\) consistent LP and VAR impulse-response estimators. Extending the primitive conditions to weak instruments, many instruments, or invalid external instruments is beyond the scope of the paper. In the empirical application, we treat the Bauer--Swanson high-frequency surprise as the external instrument and implement the same averaging procedure using the resulting IV-LP and proxy-SVAR estimates.}
remarkFollowing HounyoJung2025, one can adopt a two–stage estimator-averaging scheme: (i) first, average within the LP–based estimators and within the VAR–based estimators separately; (ii) second, average the two class-specific aggregates (LP–avg and VAR–avg) using the same MSE–minimizing weighting rule. The second-stage weight leverages the complementary strengths of LP (lower bias at short horizons) and VAR (lower variance at longer horizons) to reduce the combined estimator’s risk.

As a complementary implementation, we also consider empirically optimal weights that are chosen by directly minimizing the bootstrap MSE of the combined estimator, rather than by first estimating the bias, variance, and covariance components entering (ref). Concretely, for each bootstrap resample we re-estimate the LP and VAR impulse responses and numerically search over weights in $[0,1]$ for the value that minimizes the bootstrap MSE of the combined estimator.

We also study a more flexible specification in which the weight depends on the discrepancy between the two estimators:

equation[equation omitted — 192 chars of source]

where $a$ and $b$ are non-negative parameters chosen numerically. The normalization by the sum ensures scale invariance, and the case $b=0$ nests the constant-weight specification. Intuitively, larger discrepancies between LP and VAR lead this scheme to put less weight on LP, reflecting its typically higher variance. We refer to this extension as flexible estimator averaging. In both cases, implementation relies on the same semiparametric AR-sieve-bootstrap used for the plug-in method. Since these procedures are primarily computational alternatives and do not form the basis of our asymptotic theory, we treat them as complementary rather than central to the paper's main contribution.

Model Averaging with $R^2$-Based Weights

In this subsection, we briefly review model averaging, focusing on achieving the best in--sample fit. In the spirit of HounyoJung2025, we choose the weight between LP- and VAR-based IRF estimates using their respective in-sample coefficients of determination, $R^2$.

First, we calculate the \(R^2\) for the LP regression at each horizon \(h\), denoted as \(R^2_{LP,h}\). Because the LP method estimates a separate regression for each horizon, this measure of fit is horizon-specific. Second, we calculate the \(R^2\) for the VAR model. In the univariate case, this is the in-sample \(R^2\) from the fitted VAR equation. In the multivariate case, we use the in-sample \(R^2\) from the reduced-form VAR equation corresponding to the response variable whose impulse response is being averaged. Since the VAR is estimated once on the fixed sample, this measure is horizon-invariant and is denoted by \(R^2_{VAR}\).

Using these measures of in-sample fit, we assign the weight to the LP estimator for a given horizon \(h\) proportionally to its relative explanatory power. The model averaging weight is constructed as:

equation[equation omitted — 126 chars of source]

We then obtain the model-averaged IRF at horizon \(h\) as \[ \widehat \theta_h^{M} = \widehat w_h^{M} \, \widehat \theta_{LP,h} + \bigl(1 - \widehat w_h^{M}\bigr) \widehat \theta_{VAR,h}. \]

Asymptotic Results

Since the averaged estimator relies on a bootstrap estimate of \(w_h^\star\), this section establishes the asymptotic properties of both \(\hat w_h\) and \(\hat\theta_h(\hat w_h)\). We develop the asymptotic theory for $\widehat w_h$ and $\widehat\theta_h(\widehat w_h)$ by focusing on the baseline environment where both the LP and VAR estimators are consistent. By analyzing the AR-sieve-bootstrap min-MSE procedure under a short-memory linear DGP in this benchmark setting, we provide the rigorous theoretical justification that anchors our methodology. We impose standard regularity conditions below.

Although the estimators in Section (ref) are written as finite-lag LP and VAR regressions, the benchmark asymptotic theory does not assume that the true DGP is a fixed finite-order VAR. Instead, the true process is modeled as a short-memory Wold process. Equivalently, under standard invertibility and summability conditions, this process admits an infinite-order autoregressive representation that can be approximated by a growing-order AR/VAR sieve. The finite-lag LP and VAR specifications used in estimation should therefore be viewed as feasible approximations to this short-memory benchmark, with the sieve order increasing sufficiently slowly relative to the sample size.

assumptionThe vector process \(\{Y_t\}\) is covariance-stationary and purely nondeterministic, and admits the Wold representation \[ Y_t=\sum_{j=0}^{\infty}\Psi_j\varepsilon_{t-j}, \qquad E[\varepsilon_t\mid \mathcal F_{t-1}]=0, \qquad E[\varepsilon_t\varepsilon_t']=\Sigma_\varepsilon\succ 0. \] Assume short memory, \(\sum_{j=0}^{\infty} j\|\Psi_j\|<\infty\), and \(E\|\varepsilon_t\|^{4+\delta}<\infty\) for some \(\delta>0\).
assumptionFix $h$. There exists a (possibly singular) $2\times 2$ matrix $\Omega_h^{(2)}$ such that \begin{align} \sqrt{T}\Big((\widehat\theta_{LP,h},\widehat\theta_{VAR,h})'-(\theta_h,\theta_h)'\Big) \ \Rightarrow\ \mathcal N(0,\Omega_h^{(2)}). \end{align} In addition, for some $\delta>0$, $\sup_T \mathbb{E}\|Z_{T,h}\|^{2+\delta}<\infty$, where $Z_{T,h}$ is defined in Assumption (ref).
assumptionFor any fixed $h$, define \begin{align} Z_{T,h}=\sqrt{T}\Big((\widehat\theta_{LP,h},\widehat\theta_{VAR,h})'-(\theta_h,\theta_h)'\Big), \\ Z_{T,h}^\ast=\sqrt{T}\Big((\widehat\theta_{LP,h}^\ast,\widehat\theta_{VAR,h}^\ast)'- (\widehat\theta_{LP,h},\widehat\theta_{VAR,h})'\Big), \end{align} where $(\widehat\theta_{LP,h}^\ast,\widehat\theta_{VAR,h}^\ast)'$ is generated by the AR-sieve-bootstrap in algorithm (ref). Let $\mathbb{E}^\ast[\cdot]$ denote expectation conditional on the data, and let $BL_1$ be the class of real-valued functions $\varphi$ on $\mathbb{R}^2$ that satisfy $\sup_z|\varphi(z)|\le 1$ and have Lipschitz constant at most $1$. (i) There exists a tight $\mathbb{R}^2$-valued random vector $Z_h$ such that $Z_{T,h}\Rightarrow Z_h$ and \begin{align} \sup_{\varphi\in BL_1}\Big|\mathbb{E}^\ast[\varphi(Z_{T,h}^\ast)]-\mathbb{E}[\varphi(Z_{T,h})]\Big| \xrightarrow{p}0. \end{align} (ii) For some $\delta>0$, $\mathbb{E}^\ast[\|Z_{T,h}^\ast\|^{2+\delta}]=O_p(1)$ and $\sup_T \mathbb{E}[\|Z_{T,h}\|^{2+\delta}]<\infty$. (iii) $B=B(T)\to\infty$ as $T\to\infty$.

To connect the finite-sample MSE weight in (2.6) to the root-\(T\) asymptotic theory, it is useful to rescale the risk components. Define $a_{T,h}=V_{L,T,h}+b_{L,T,h}^2,\quad d_{T,h}=V_{V,T,h}+b_{V,T,h}^2,\quad f_{T,h}=C_{T,h}+b_{L,T,h}b_{V,T,h}$, and let $A_{T,h}=T a_{T,h},\quad D_{T,h}=T d_{T,h},\quad F_{T,h}=T f_{T,h}$. Since multiplying both the numerator and denominator of (2.6) by \(T\) does not change the weight, we can write \[ w_{T,h}^{\star} = \frac{d_{T,h}-f_{T,h}} {a_{T,h}+d_{T,h}-2f_{T,h}} = \frac{D_{T,h}-F_{T,h}} {A_{T,h}+D_{T,h}-2F_{T,h}}. \] In the benchmark case in which both the LP and VAR estimators are root-\(T\) consistent and asymptotically centered at \(\theta_h\), these scaled risk components converge to the corresponding entries of \(\Omega_h^{(2)}\).

assumptionFor any fixed $h$, the oracle min-MSE weight $w_h^\star$ in (ref) is unique after clipping to $[0,1]$, and the asymptotic variance for the averaged estimator evaluated at the oracle weight is finite and well-defined.
assumptionFor any fixed \(h\), the scaled oracle risk components satisfy $(A_{T,h},D_{T,h},F_{T,h}) \to (A_h,D_h,F_h)$, where $A_h+D_h-2F_h \ge c>0$ for some constant \(c>0\). Moreover, the limiting oracle weight $w_h^\star = \frac{D_h-F_h}{A_h+D_h-2F_h}$ lies in the interior of \([0,1]\): \(w_h^\star\in[\underline w,1-\underline w]\) for some \(\underline w\in(0,1/2)\).
assumptionFor any fixed $h$, in addition to Assumption (ref)(ii), assume that conditional on the data $\mathbb{E}^\ast[\|Z_{T,h}^\ast\|^{4+\delta}]=O_p(1)$ for the same $\delta>0$ as in Assumption (ref).

Assumption (ref) assumes a short-memory Wold DGP. Under Assumptions (ref)--(ref) and routine regularity conditions---including nonsingularity of the relevant projection matrices and, for sieve-based VAR/LP approximations, a lag order \(p_T\) that grows sufficiently slowly with \(T\)---the LP and VAR IRF estimators are consistent and jointly asymptotically normal for each fixed horizon \(h\). Moreover, AR-sieve-bootstrap approximations are valid for broad Wold-type processes for statistics whose limits depend on second-order structure. Finally, LP and VAR target the same impulse responses asymptotically when the lag length increases (PlagborgMollerWolf2021).

Assumption (ref) postulates joint root-$T$ asymptotic normality for the LP and VAR IRF estimators at a fixed horizon $h$. This condition is standard for OLS-based LP and VAR estimators under Assumption (ref) and routine rank and weak-dependence assumptions; for VAR IRFs it follows by applying the delta method to the smooth map from VAR OLS coefficients to the horizon-$h$ IRF. Writing $\Omega_h^{(2)}=

pmatrix[pmatrix omitted — 70 chars of source]

$ as in \eqref{eq:jointCLT}, the (infeasible) oracle weight $w_h^\star=\arg\min_{w} e(w)'\Omega_h^{(2)}e(w)$, with $e(w)=(w,1-w)'$, is

align[align omitted — 124 chars of source]

Assumption (ref) is a high-level statement that the AR-sieve-bootstrap reproduces the first-order limit law of the centered/scaled LP--VAR pair at a fixed horizon $h$ (expressed via the bounded-Lipschitz metric). The moment condition in (ii) is included to justify convergence of quadratic functionals used by the MSE weight. The AR-sieve-bootstrap is designed for short-memory linear processes and may provide a poor approximation when the DGP features strong nonlinearity, structural breaks, unit roots or near-unit roots, long memory, heavy tails, or strong conditional heteroskedasticity. In the latter case, a wild or heteroskedasticity-robust bootstrap variant may be more appropriate. Assumption (ref) ensures the oracle problem is well posed (a unique clipped minimizer exists) and that the asymptotic variance of the averaged estimator is finite, so plug-in variance estimation is meaningful. When LP and VAR are both root-$T$ consistent for the same $\theta_h$, the joint asymptotic covariance $\Omega_h^{(2)}$ of $(\widehat\theta_{LP,h},\widehat\theta_{VAR,h})$ may be singular. This does not affect consistency of $\widehat w_h$ or of the plug-in variance for $\widehat\theta_h(\widehat w_h)$; it only implies we should not require $\Omega_h^{(2)}$ to be positive definite.

Assumption (ref) is a nondegeneracy condition for the scaled asymptotic risk problem. The unscaled finite-sample denominator in (2.6) is typically of order \(T^{-1}\) under root-\(T\) asymptotics, so it is the scaled denominator \(A_{T,h}+D_{T,h}-2F_{T,h}\), rather than \(a_{T,h}+d_{T,h}-2f_{T,h}\), that must be bounded away from zero. The lower bound on \(A_{T,h}+D_{T,h}-2F_{T,h}\) makes the mapping from the scaled risk components to the oracle weight locally smooth. The additional interiority condition \(w_h^\star\in[\underline w,1-\underline w]\) ensures that clipping is asymptotically inactive. While it is possible for the denominator $A_h+D_h-2F_h$ to approach zero, for example in the lag-augmented limit of PlagborgMollerWolf2021 when the LP and VAR estimators become asymptotically equivalent, Algorithm (ref) already includes a numerical safeguard that returns a default weight when the estimated denominator falls below a small threshold, so the procedure is operationally well defined for any sample. In addition, when LP and VAR are close to asymptotically equivalent, the weight is poorly identified but the averaged estimator itself is well-behaved, because any convex combination of two nearly-coincident estimators delivers nearly the same value. The Monte Carlo evidence in Section (ref) confirms this pattern.

For the benchmark first-order theory, $\widehat w_h$ denotes the plug-in weight computed from the scaled variance-covariance components $(\widehat A_{T,h},\widehat D_{T,h},\widehat F_{T,h})$; the finite-sample implementation in Algorithm (ref) additionally includes bootstrap bias terms.

theorem[Consistency and rate of AR-sieve-bootstrap weights] For any fixed \(h\), under Assumptions (ref) and (ref), $\widehat w_h \xrightarrow{p} w_h^\star$. If, in addition, \begin{align} \left\| (\widehat A_{T,h},\widehat D_{T,h},\widehat F_{T,h}) - (A_h,D_h,F_h) \right\| = O_p(r_T)+O_p(B^{-1/2}), \end{align} for some deterministic sequence \(r_T\to0\), then $|\widehat w_h-w_h^\star|=O_p(r_T)+O_p(B^{-1/2})$.

Let $\Omega$ be a $2\times 2$ covariance matrix. Define

align[align omitted — 102 chars of source]
theorem[Consistency and asymptotic normality of $\widehat\theta_h(\widehat w_h)$] For any fixed $h$, under Assumptions (ref), (ref), (ref), and (ref), $\widehat\theta_h(\widehat w_h)\xrightarrow{p}\theta_h$, and \begin{align} \sqrt{T}\big(\widehat\theta_h(\widehat w_h)-\theta_h\big) \Rightarrow \mathcal N\!\Big(0,\, V_h(w_h^\star,\Omega_h^{(2)})\Big), \end{align} where $V_h(w_h^\star,\Omega_h^{(2)})$ is defined by equation ((ref)) with $w=w_h^\star$ and $\Omega=\Omega_h^{(2)}$.

Theorem (ref) shows that the AR-sieve-bootstrap plug-in weight $\widehat w_h$ is consistent for the limiting oracle weight $w_h^\star$ at each fixed horizon $h$, provided that the scaled bootstrap risk components consistently estimate their limiting counterparts. The scaling is important because, under the benchmark root-$T$ asymptotics, the unscaled variance and covariance terms entering the finite-sample MSE are of order $T^{-1}$. Thus, the relevant nondegeneracy condition is imposed on the scaled denominator $A_h+D_h-2F_h$, rather than on the raw denominator $a_{T,h}+d_{T,h}-2f_{T,h}$. If a rate statement is desired, it depends on the convergence rate of the scaled risk-component estimator, denoted by $r_T$, together with the Monte Carlo simulation error $B^{-1/2}$. We leave $r_T$ as a high-level rate because its exact value depends on the AR-sieve approximation, the lag-order sequence, and the statistic being bootstrapped. Importantly, Theorem (ref) only requires $\widehat w_h \overset{p}{\longrightarrow} w_h^\star$, since the plug-in weight error is then first-order negligible. Theorem (ref) shows that $\widehat\theta_h(\widehat w_h)$ is consistent for $\theta_h$ and asymptotically normal at the root-$T$ rate, with asymptotic variance $V_h(w_h^\star,\Omega_h^{(2)})$. Although Theorem (ref) applies to correctly specified models, our simulations show that our estimator performs well also under misspecification.

theorem[Consistency and rate of plug-in asymptotic variance] For any fixed $h$, let $\widehat V_h = V_h(\widehat w_h,\widehat\Omega_h^{(2)})$. Suppose that $\widehat\Omega_h^{(2)}\xrightarrow{p}\Omega_h^{(2)}$. Under the conditions of Theorem (ref), $\widehat V_h \xrightarrow{p} V_h(w_h^\star,\Omega_h^{(2)})$. If, in addition, $|\widehat w_h-w_h^\star| = O_p(r_T)+O_p(B^{-1/2})$ and $\|\widehat\Omega_h^{(2)}-\Omega_h^{(2)}\| = O_p(s_T)+O_p(B^{-1/2})$ for some $s_T\to0$, then $|\widehat V_h-V_h(w_h^\star,\Omega_h^{(2)})| = O_p(r_T+s_T)+O_p(B^{-1/2})$.

Theorem (ref) justifies plug-in inference based on $\widehat V_h=V_h(\widehat w_h,\widehat\Omega_h^{(2)})$. In practice, one may estimate the $2\times2$ matrix $\Omega_h^{(2)}$ by forming the stacked vector of estimating equations for $(\widehat\theta_{LP,h},\widehat\theta_{VAR,h})$ and applying a standard HAC estimator (e.g., Newey--West) to the resulting score/influence-function series. Equivalently, one can estimate $\Omega_h^{(2)}$ via the same AR-sieve bootstrap used to construct $\widehat w_h$: compute $Z_{T,h}^{\ast(b)}$ as in (ref) and set $\widehat\Omega_h^{(2)}=B^{-1}\sum_{b=1}^B Z_{T,h}^{\ast(b)}Z_{T,h}^{\ast(b)\prime}$ (with recentering if desired).

As an alternative to plug-in Wald inference, one can use a AR-sieve-bootstrap to obtain confidence bands that are robust to conditional heteroskedasticity. Fit a stable VAR($p$) sieve to $\{y_t\}$, form bootstrap innovations $u_t^{\ast(b)}=\eta_t^{(b)}\widehat u_t$ with i.i.d.\ Rademacher multipliers $\eta_t^{(b)}\in\{-1,+1\}$, and generate $y_t^{\ast(b)}=\sum_{j=1}^p \widehat A_j y_{t-j}^{\ast(b)}+u_t^{\ast(b)}$ recursively. Re-estimate LP and VAR on each pseudo-sample and recompute the averaged estimator to obtain $\widehat\theta_h^{\ast(b)}=\widehat\theta_h^{\ast(b)}(\widehat w_h^{\ast(b)})$. A convenient centered $1-\alpha$ interval is $\big[\widehat\theta_h(\widehat w_h)-\delta_{\alpha,h},\ \widehat\theta_h(\widehat w_h)+\delta_{\alpha,h}\big]$, where $\delta_{\alpha,h}$ is the empirical $(1-\alpha)$ quantile of $\big|\widehat\theta_h^{\ast(b)}-\widehat\theta_h(\widehat w_h)\big|$ across $b=1,\ldots,B$.

Monte Carlo Evidence

To assess the finite-sample performance of our methods, we consider both a simple univariate design and a richer multivariate design. The univariate exercise is useful for transparently illustrating the LP--VAR bias--variance trade-off and for evaluating how well the feasible procedures approximate the oracle weights. The multivariate exercise then assesses the same methods in a more realistic macroeconomic environment. In both the univariate and the multivariate exercises we use a fixed lag-selection rule; thus our estimators are not asymptotically equivalent in the Plagborg-Møller–Wolf sense. Although our theoretical results assume that both methods estimate the true IRF consistently, in the simulations we also explore what happens when the model is misspecified, showing that our proposed estimators still perform relatively well.

An Univariate ARMA Design

We begin with a simple univariate DGP, also considered by MontielOleaetal2025. This design allows us to (i) verify that the infeasible oracle estimator averaging behaves as predicted by the LP--VAR bias--variance trade-off, (ii) evaluate how well the sieve-bootstrap and flexible implementations approximate the oracle weights in finite samples, and (iii) compare the risk properties of estimator averaging with those of $R^2$-based model averaging.

Model and Estimators

We conduct 1,000 Monte Carlo replications of a univariate ARMA(1,1) process

equation[equation omitted — 164 chars of source]

Each replication uses a burn-in period of 200 observations to reach stationarity before retaining a sample of length $T$. The true impulse response to a one-unit innovation at time $t$ is $\theta_0 = 1,$ $\theta_1 = \rho + \alpha$, $\theta_h = \rho\,\theta_{h-1} \;\; \text{for } h\ge 2,$ that is, $\theta_h=\rho^h + \alpha\rho^{h-1}$ for $h\ge1$. We evaluate horizons $h \in \{1,\ldots,H_{\max}\}$ and set $H_{\max}=10$. We consider $\rho\in\{0.5,0.9\}$ and $\alpha\in\{0.5,0.9\}$ for various sample sizes $T$. At each horizon, we report results for the following estimators:(1) the local projection estimator (LP), using the specification in MontielOleaetal2025, (2) the VAR estimator, implemented as an AR(1), (3) the infeasible oracle averaged estimator based on MSE minimization, (4) the sieve-bootstrap plug-in estimator averaging procedure, (5) the sieve-bootstrap flexible estimator averaging procedure, and (6) the model-averaging estimator based on maximizing $R^2$. To approximate the oracle weight and implement the feasible averaging procedures, we use 500 bootstrap draws.

Horizon-Specific Risk and Weights

We first fix $T=240$, an empirically relevant sample size used in MontielOleaetal2025, and study horizons $h=1,\ldots,H_{\max}$. This allows us to trace how the LP--VAR bias--variance trade-off evolves across horizons and to assess whether the oracle and feasible averaging rules shift weight from LP at short horizons toward VAR at longer horizons, as predicted by the MSE decomposition in (ref).

Figures (ref)--(ref) report RMSE and weights across $(\rho,\alpha)$. The patterns line up closely with the bias--variance logic in (ref). Under the ARMA(1,1) DGP in (ref), $\theta_1=\rho+\alpha$ and $\theta_h=\rho^h+\alpha\rho^{h-1}$ for $h\ge2$, whereas the AR(1)-based VAR estimator implies $\widehat\theta_{VAR,h}\approx\widehat\beta^h$. Hence, when $\alpha>0$, the VAR is downward biased at short horizons, with the bias decaying roughly at rate $\alpha\rho^{h-1}$, while LP is approximately unbiased at short horizons but becomes increasingly variable as $h$ rises.

figure[figure omitted — 995 chars of source]

The figures confirm this trade-off. LP performs better at short horizons. When $\alpha$ is sizable, the difference in RMSE between the LP and VAR estimators is noticeably larger. Roughly after 7 periods, the VAR takes over because LP variance rises while VAR bias decays. A higher value of $\rho$ slows the decay of the VAR bias, thus keeping its RMSE elevated longer.

The combined estimator using oracle weights yields an RMSE that is never larger---and is in fact somewhat smaller---than the minimum of the individual LP and VAR estimators at each horizon. The feasible MSE-based estimator averaging approaches closely follow the performance of the oracle, though they naturally cannot match it perfectly. In contrast, the model averaging estimator sometimes exhibits substantially worse performance in terms of RMSE compared to the MSE-based methods, a discrepancy that is especially pronounced at short horizons.

figure[figure omitted — 1,013 chars of source]

The reason for this difference in performance becomes obvious from the weights. The oracle weights in (ref) display the benchmark pattern: they place more mass on LP at short horizons and then gradually shift toward VAR as $h$ increases. The feasible MSE-based estimator averaging broadly reproduces this pattern, and its RMSE remains close to the oracle benchmark.

The under-performance of the model averaging estimator is due to the inflexibility of its weighting scheme based on $R^2$. Because this approach compares the in-sample fit of the VAR with that of the LP, it relies on an $R^2_{VAR}$ that remains constant across all horizons and an $R^2_{LP,h}$ that exhibits very little variation as $h$ increases. Consequently, these $R^2$-based weights are not flexible enough to capture the complex, shifting balance between bias and variance. While the scheme correctly recognizes that the relative advantage of the VAR estimator increases with the horizon, it adjusts too rigidly: the weight assigned to the LP estimator starts near $0.5$ at $h=1$ and declines only gradually as the horizon extends. This excessively smooth trajectory prevents the model averaging approach from adapting to the sharper, horizon-specific changes in the true risk profile.

Overall, the feasible MSE-based procedures successfully capture the bias--variance trade-off predicted by (ref)--(ref): they assign more weight to the LP where bias reduction matters, and shift toward the VAR where variance dominates. Across all horizons, these estimator-averaging approaches closely track the best attainable oracle benchmark, delivering modest to substantial risk reductions, depending on the horizon and the design, relative to relying on a single method. By contrast, the $R^2$-based model averaging is too rigid to adapt to this shifting risk profile, resulting in suboptimal performance throughout the projection horizon.

Finite-Sample Convergence

To complement the horizon-profile evidence above and connect the simulations to our large-sample theory, we vary the sample size $T\in\{200,400,800\}$ and focus on representative horizons $h\in\{1,3,6\}$. For each sample size, we perform 1,000 Monte Carlo replications.

Using the same ARMA(1,1) DGP in (ref) and the same set of estimators, we report the RMSE of the estimated weights in Table (ref) and the RMSE of the corresponding IRF estimators in Table (ref). We again consider the four $(\alpha,\rho)$ designs.

table[table omitted — 1,768 chars of source]
sidewaystable[!htbp] \caption{RMSE of IRF estimates relative to true IRF} \resizebox{\textwidth}{!}{ \begin{tabular}{lcccccccccccccccccc} \toprule & \multicolumn{6}{c}{$h=1$} & \multicolumn{6}{c}{$h=3$} & \multicolumn{6}{c}{$h=6$} \\ \cmidrule(lr){2-7}\cmidrule(lr){8-13}\cmidrule(lr){14-19} $T$ & $\widehat{\theta}_{LP}$ & $\widehat{\theta}_{VAR}$ & $\widehat{\theta}_O$ & $\widehat{\theta}_P$ & $\widehat{\theta}_F$ & $\widehat{\theta}_M$ & $\widehat{\theta}_{LP}$ & $\widehat{\theta}_{VAR}$ & $\widehat{\theta}_O$ & $\widehat{\theta}_P$ & $\widehat{\theta}_F$ & $\widehat{\theta}_M$ & $\widehat{\theta}_{LP}$ & $\widehat{\theta}_{VAR}$ & $\widehat{\theta}_O$ & $\widehat{\theta}_P$ & $\widehat{\theta}_F$ & $\widehat{\theta}_M$ \\ \midrule \multicolumn{19}{c}{$\alpha=0.5,\ \rho=0.5$}\\ \hline $T=200$ & 0.0958 & 0.2972 & 0.0958 & 0.0990 & 0.0979 & 0.1149 & 0.1136 & 0.1204 & 0.0911 & 0.1084 & 0.1116 & 0.1177 & 0.1125 & 0.1070 & 0.0760 & 0.0886 & 0.0973 & 0.1044 \\ $T=400$ & 0.0827 & 0.2929 & 0.0827 & 0.0832 & 0.0830 & 0.0924 & 0.0822 & 0.1155 & 0.0676 & 0.0785 & 0.0801 & 0.0896 & 0.0741 & 0.1029 & 0.0616 & 0.0760 & 0.0778 & 0.0808 \\ $T=800$ & 0.0705 & 0.2881 & 0.0705 & 0.0706 & 0.0706 & 0.0754 & 0.0580 & 0.1166 & 0.0507 & 0.0554 & 0.0564 & 0.0617 & 0.0539 & 0.1035 & 0.0469 & 0.0604 & 0.0600 & 0.0620 \\ \addlinespace \multicolumn{19}{c}{$\alpha=0.5,\ \rho=0.9$}\\ \hline $T=200$ & 0.3446 & 0.6597 & 0.3446 & 0.3446 & 0.3446 & 0.3581 & 0.1525 & 0.0842 & 0.0754 & 0.1202 & 0.1377 & 0.1296 & 0.1238 & 0.1369 & 0.0874 & 0.1185 & 0.1212 & 0.1223 \\ $T=400$ & 0.3425 & 0.6566 & 0.3425 & 0.3425 & 0.3425 & 0.3491 & 0.1266 & 0.0762 & 0.0563 & 0.0923 & 0.1025 & 0.1120 & 0.0820 & 0.1335 & 0.0686 & 0.0891 & 0.0865 & 0.0908 \\ $T=800$ & 0.3346 & 0.6527 & 0.3346 & 0.3346 & 0.3346 & 0.3378 & 0.1034 & 0.0748 & 0.0420 & 0.0659 & 0.0703 & 0.0926 & 0.0605 & 0.1346 & 0.0512 & 0.0663 & 0.0629 & 0.0669 \\ \addlinespace \multicolumn{19}{c}{$\alpha=0.9,\ \rho=0.5$}\\ \hline $T=200$ & 0.1186 & 0.4679 & 0.1186 & 0.1197 & 0.1240 & 0.1304 & 0.1941 & 0.3261 & 0.1941 & 0.2119 & 0.2025 & 0.2367 & 0.2449 & 0.1846 & 0.1835 & 0.1909 & 0.1939 & 0.2063 \\ $T=400$ & 0.1086 & 0.4613 & 0.1086 & 0.1087 & 0.1095 & 0.1145 & 0.1503 & 0.3080 & 0.1503 & 0.1632 & 0.1597 & 0.1882 & 0.1715 & 0.1520 & 0.1454 & 0.1505 & 0.1504 & 0.1569 \\ $T=800$ & 0.0983 & 0.4586 & 0.0983 & 0.0983 & 0.0983 & 0.1012 & 0.1142 & 0.3002 & 0.1142 & 0.1159 & 0.1159 & 0.1365 & 0.1257 & 0.1354 & 0.1193 & 0.1239 & 0.1212 & 0.1287 \\ \addlinespace \multicolumn{19}{c}{$\alpha=0.9,\ \rho=0.9$}\\ \hline $T=200$ & 0.4025 & 0.8608 & 0.4025 & 0.4025 & 0.4029 & 0.4112 & 0.3834 & 0.6300 & 0.3834 & 0.3951 & 0.3882 & 0.4477 & 0.3676 & 0.3802 & 0.3573 & 0.3498 & 0.3424 & 0.3672 \\ $T=400$ & 0.3995 & 0.8550 & 0.3995 & 0.3995 & 0.3995 & 0.4038 & 0.3582 & 0.6143 & 0.3582 & 0.3587 & 0.3585 & 0.4000 & 0.3006 & 0.3530 & 0.3006 & 0.2902 & 0.2832 & 0.3223 \\ $T=800$ & 0.3922 & 0.8526 & 0.3922 & 0.3922 & 0.3922 & 0.3943 & 0.3313 & 0.6078 & 0.3313 & 0.3313 & 0.3313 & 0.3511 & 0.2677 & 0.3410 & 0.2677 & 0.2606 & 0.2582 & 0.2997 \\ \bottomrule \end{tabular} } \begin{minipage}[c]{130 mm} Note: $\widehat{\theta}_{LP}$ and $\widehat{\theta}_{VAR}$ are the LP and VAR IRF estimators. $\widehat{\theta}_O$ uses the infeasible oracle min-MSE weight $w_h^\star$ in (ref). $\widehat{\theta}_P$ uses the sieve-bootstrap plug-in weight $\widehat w_P$ (Algorithm (ref)). $\widehat{\theta}_F$ uses the flexible sieve-bootstrap weight $w_F$ in (ref). $\widehat{\theta}_M$ uses the $R^2$-based model-averaging weight $\widehat w_{M}$ in (ref). \end{minipage}

Two main messages emerge. First, the estimated LP--VAR weights become more accurate as $T$ increases. In Table (ref), the RMSE of both the plug-in weight $\widehat w_P$ and the flexible weight $\widehat w_F$ (measured relative to the oracle weights) generally declines with $T$. It is important to note that while Theorem (ref) establishes the formal consistency of these weights, and a high-level rate under an additional scaled-risk rate condition, for the baseline case where both the LP and VAR estimators are consistent, this specific simulation design features misspecified estimators. Thus, rather than merely illustrating the theorem, these numerical results demonstrate an important broader finding: the feasible weighting procedures maintain strong finite-sample performance---evidenced by a substantially shrinking RMSE---even in the presence of underlying model misspecification. This improvement is especially clear at short horizons, where the LP--VAR trade-off is most informative for risk estimation.

Second, the averaged IRF estimators based on these estimated weights inherit the same convergence and deliver strong risk performance in finite samples. In Table (ref), the estimator-averaging procedures $\widehat{\theta}_P$ and $\widehat{\theta}_F$ move toward the oracle benchmark $\widehat{\theta}_O$ as $T$ increases, in line with Theorem (ref). Across designs, estimator averaging implements the intended bias--variance trade-off: it places more weight on LP where VAR misspecification bias is most relevant, especially at short horizons, and shifts toward VAR as the horizon increases and LP variance becomes more important.

Overall, these additional results reinforce the main findings from the univariate design: bootstrap-based estimator averaging converges toward oracle averaging and yields low RMSE relative to LP, VAR, and $R^2$-based model averaging when the objective is to minimize IRF estimation risk.

A Multivariate SVARMA Design

We next extend the analysis to a multivariate setting to evaluate the performance of the estimators in a more realistic macroeconomic environment. We consider two DGPs: an SVAR(4) model and an SVARMA(4,1) model. In the latter case, the finite-order VAR and LP models used for estimating the impulse response functions are misspecified. In both cases, we assume that the structural shocks are observed so as to abstract from identification issues and focus directly on the estimation of impulse responses.

Model and Estimators

The data are generated from a three-variable system ($n=3$) following a general SVARMA($p,q$) process:

equation[equation omitted — 152 chars of source]

where $M_0$ is the structural impact matrix. We consider two specifications:

enumerate• SVAR(4): $p=4$, $q=0$. In this case, a VAR with sufficient lag length can capture the true dynamics. • SVARMA(4,1): $p=4$, $q=1$. This introduces a moving-average component; consequently, finite-order VAR and LP approximations remain formally misspecified, even when lag lengths are selected via information criteria.

The specific coefficient matrices $(A_j,M_k)$ are reported in Appendix (ref). We simulate 1,000 replications and, for the main finite-sample analysis, set $T=200$.

We evaluate the following estimators:

enumerate• Local Projection (LP): the lag length is set equal to the lag order chosen for the VAR. • VAR: the lag order is selected by AIC, with a maximum lag length of 8. • Oracle: the infeasible averaged estimator based on MSE minimization, computed using 500 bootstrap draws. • Estimator Averaging: the feasible MSE-based averaged estimator, where the weights are estimated via a sieve bootstrap with the DGP approximated by a VAR selected by AIC. We use 500 bootstrap replications. To save space, we report only the plug-in estimator since the flexible estimator yields similar results. • Model Averaging: the averaged estimator with weights chosen to maximize in-sample $R^2$ at each horizon.

Horizon-Specific Risk and Weights

Figure (ref) reports RMSE and estimated weights for the two DGPs. In the case of SVAR(4), the bias--variance trade-off is clearly visible; see Panels (a) and (c). LP delivers lower RMSE for horizons $h=2$ to $h=6$, after which VAR becomes superior because its structured dynamics produces greater efficiency. The oracle estimator closely tracks the lower envelope of the two individual estimators.

figure[figure omitted — 1,098 chars of source]

In terms of RMSE, the feasible estimator-averaging procedure is nearly indistinguishable from the oracle, indicating that the risk-based weighting rule is well estimated in finite samples. By contrast, model averaging performs noticeably worse, with higher RMSE at most horizons. The weight plot in Panel (c) makes the reason clear: estimator averaging tracks the shape of the oracle weights, assigning more weight to LP at short horizons and then shifting toward VAR, whereas model averaging tends to underweight LP early and overweight it later because it is driven by in-sample fit rather than IRF estimation risk.

In the misspecified case (the SVARMA(4,1) DGP), LP retains its advantage at short horizons ($h<4$), while VAR performs better thereafter; see Panels (b) and (d). The oracle again tracks the minimum RMSE, and the estimator averaging remains close to that oracle benchmark. Model averaging, however, performs less well, especially at horizons where the LP--VAR gap is large. In particular, the $R^2$-based rule fails to reproduce the relatively rapid decline in oracle LP weight, leading to suboptimal averaging.

Finite-sample convergence

To examine the large-sample behavior in the multivariate setting, we vary the sample size $T\in\{200,800,2000\}$ and report results for horizons $h\in\{1,6,12\}$.

Table (ref) presents the RMSE of the impulse response estimates based on 1,000 Monte Carlo replications and 500 bootstrap iterations. Consistent with the univariate evidence, the RMSE of estimator averaging declines with sample size and remains below that of model averaging across almost all specifications. Even at $T=2000$, where LP and VAR have both become quite accurate, estimator averaging typically yields a slight improvement by optimally combining their remaining finite-sample differences.

Table (ref) reports the RMSE of the estimated weights relative to the oracle weights. Unlike the IRF estimates, the weights do not converge monotonically in all cases, particularly under SVARMA(4,1). This may be related to a flat or weakly identified weighting problem, which presumably occurs when the VAR and LP estimates converge toward each other. In such cases, estimator averaging can remain close to risk-optimal even when the estimated weights display discernible finite-sample variation.

sidewaystable[!htbp] \caption{RMSE of multivariate IRF estimators} \resizebox{\textwidth}{!}{ \begin{tabular}{lccccccccccccccc} \toprule & \multicolumn{5}{c}{$h=1$} & \multicolumn{5}{c}{$h=6$} & \multicolumn{5}{c}{$h=12$} \\ \cmidrule(lr){2-6}\cmidrule(lr){7-11}\cmidrule(lr){12-16} $T$ & $\widehat{\theta}_{VAR}$ & $\widehat{\theta}_{LP}$ & $\widehat{\theta}_O$ & $\widehat{\theta}_P$ & $\widehat{\theta}_M$ & $\widehat{\theta}_{VAR}$ & $\widehat{\theta}_{LP}$ & $\widehat{\theta}_O$ & $\widehat{\theta}_P$ & $\widehat{\theta}_M$ & $\widehat{\theta}_{VAR}$ & $\widehat{\theta}_{LP}$ & $\widehat{\theta}_O$ & $\widehat{\theta}_P$ & $\widehat{\theta}_M$ \\ \midrule \multicolumn{16}{c}{SVAR(4)}\\ \hline 200 & 0.1299 & 0.1071 & 0.0847 & 0.0906 & 0.1215 & 1.0344 & 1.0445 & 1.0168 & 1.0247 & 1.0320 & 0.8433 & 1.1363 & 0.8433 & 0.8445 & 0.9390 \\ 800 & 0.0330 & 0.0562 & 0.0291 & 0.0272 & 0.0311 & 0.4253 & 0.4439 & 0.4231 & 0.4216 & 0.4240 & 0.4653 & 0.5627 & 0.4653 & 0.4662 & 0.4888 \\ 2000 & 0.0144 & 0.0339 & 0.0135 & 0.0126 & 0.0137 & 0.2649 & 0.2803 & 0.2647 & 0.2639 & 0.2649 & 0.3141 & 0.3428 & 0.3141 & 0.3140 & 0.3167 \\ \addlinespace \multicolumn{16}{c}{SVARMA(4,1)}\\ \hline 200 & 0.0802 & 0.0302 & 0.0289 & 0.0427 & 0.0671 & 0.1955 & 0.2205 & 0.1955 & 0.1964 & 0.1969 & 0.1182 & 0.2449 & 0.1182 & 0.1183 & 0.1656 \\ 800 & 0.0218 & 0.0148 & 0.0124 & 0.0143 & 0.0186 & 0.0985 & 0.1065 & 0.0985 & 0.0985 & 0.0985 & 0.0586 & 0.1174 & 0.0586 & 0.0587 & 0.0642 \\ 2000 & 0.0094 & 0.0091 & 0.0067 & 0.0071 & 0.0081 & 0.0612 & 0.0657 & 0.0612 & 0.0613 & 0.0612 & 0.0383 & 0.0699 & 0.0383 & 0.0383 & 0.0385 \\ \bottomrule \end{tabular} } \begin{minipage}{\textwidth} Note: $\widehat{\theta}_{VAR}$ and $\widehat{\theta}_{LP}$ are the VAR and LP IRF estimators. $\widehat{\theta}_O$ uses the infeasible oracle min-MSE weight $w_h^\star$ in (ref). $\widehat{\theta}_P$ uses the sieve-bootstrap plug-in weight $\widehat w_P$ (Algorithm (ref)). $\widehat{\theta}_M$ uses the $R^2$-based model-averaging weight $\widehat w_{M}$ in (ref). Results are based on 1,000 Monte Carlo replications and 500 bootstrap iterations. \end{minipage}
table[table omitted — 543 chars of source]

Empirical Application

Section V of BauerSwanson2023 reassess the dynamic effects of monetary policy by addressing two key challenges in high-frequency identification: instrument relevance and exogeneity. To improve relevance, they expand the standard set of monetary policy events beyond FOMC announcements to include press conferences, speeches, and testimony by the Federal Reserve Chair, substantially increasing the variation of the surprise series. To ensure exogeneity, they argue that conventional high-frequency surprises suffer from endogeneity because they are systematically correlated with publicly available macroeconomic and financial data predating the announcements---rather than being driven by central bank “information effects.” To address this, they orthogonalize the high-frequency surprises with respect to these pre-announcement variables and use the resulting residual as an external instrument for the monetary policy shock. They then estimate IV-LP and IV-SVAR impulse responses of yields, activity, prices, and financial conditions, finding that this correction produces stronger and more plausible macroeconomic estimates.

We replicate their baseline setup using the same monthly data set, orthogonalized high-frequency surprise, and IV specification. We study the responses of the two-year Treasury yield (GBY), industrial production (IP), consumer prices (CPI), and the excess bond premium (EBP) to a 25-basis-point contractionary monetary policy shock. For each variable, we compute IV-LP and IV-SVAR impulse responses based on their specifications, as well as model-averaging combinations of the two. Finally, we also compute the external-instrument version of the estimator-averaging combination described in Remark (ref). In the empirical application, we present only the plug-in version, as the flexible estimator averaging yields very similar estimates.

For inference, we employ a nested VAR-sieve wild bootstrap procedure. We first approximate the data-generating process by estimating a reduced-form VAR model, where the lag order is selected via the Bayesian information criterion. To handle potential heteroskedasticity and accommodate the missing observations in the instrument, we apply a wild bootstrap to the data. Specifically, to preserve the identifying contemporaneous correlation between the reduced-form VAR residuals and the external instrument, the same sequence of random wild multipliers---drawn from a Rademacher distribution---is applied simultaneously to both the VAR residuals and the instrument during resampling. Our procedure relies on a double bootstrap architecture in which both the outer and inner loops consist of 500 iterations. The inner bootstrap loop is utilized to approximate the finite-sample variances and biases of the individual IV-SVAR and IV-LP estimators, as well as their covariance, which provides the moments necessary to calculate the optimal weights for the combined estimator. Subsequently, the outer bootstrap loop generates the empirical distribution of all the estimators. The confidence bands are computed by extracting the middle 68 percent of the bootstrap estimates centered around the point estimates from the original sample. Figure (ref) displays the resulting impulse responses; Figure (ref) shows the corresponding weights on LP. Appendix (ref) also presents the impulse responses with confidence bands for all estimators.

Our baseline IV-VAR and IV-LP estimates successfully replicate the results documented in BauerSwanson2023 (see their Figures 3 and A2, respectively). The IV-VAR responses (red lines) are smooth and tightly shaped: GBY and EBP jump on impact and then gradually return toward zero, IP exhibits a modest hump-shaped decline, and CPI shows a prolonged disinflation. By contrast, the IV-LP responses (blue lines) are much more volatile. For GBY and EBP, LP IRFs oscillate and change sign several times; for IP the decline is steeper and very persistent; for CPI the response eventually turns positive and drifts upward, implying an implausible long-run rise in price level after a contractionary shock. This stark contrast between low-variance but potentially biased IV-VAR estimates and low-bias but noisy IV-LP estimates mirrors the bias–variance trade-off highlighted in our simulations and motivates the use of averaging estimators.

figure[figure omitted — 312 chars of source]
figure[figure omitted — 265 chars of source]

The averaging procedures stabilize the IRFs while preserving their main qualitative features. In many cases, they tilt toward the more economically plausible trajectory. Estimator averaging (green dashed lines with triangles) pulls the paths toward the smoother IV-VAR trajectories wherever the LP estimates are extremely noisy, but still allows for some deviations where the LP and VAR disagree. For GBY, the estimator-averaged IRF shows a sharp increase in the two-year yield that decays within a year, avoiding the large negative swings of the LP while not over-smoothing the near-term reaction. Furthermore, unlike the VAR response, it never sinks deeply into negative territory, a profile that is economically plausible following a tightening shock. For IP, estimator averaging yields a moderate, hump-shaped decline that lies between the highly persistent LP response and the more muted VAR response. For CPI, the combined estimator produces a smaller and less protracted disinflation than the VAR, while successfully avoiding the “price puzzle” anomaly exhibited by the LP. For EBP, estimator averaging delivers a sharp, short-run rise in credit spreads that quickly mean-reverts, eliminating the extreme volatility of the LP estimator while preserving the intuitive tightening in financial conditions after a monetary contraction.

Model averaging based on $R^2$ (orange lines with squares) behaves differently. As Figure (ref) shows, the $R^2$-based weights on LP are relatively flat across horizons---around one-half for GBY, IP, and CPI, and somewhat lower for EBP. This reflects the fact that LP and VAR achieve similar in-sample fit even when the LP is estimated for longer horizons. Consequently, the model-averaged IRFs remain more heavily influenced by LP at longer horizons. For GBY and EBP, this yields more volatile responses than estimator averaging. The trajectory of the former drops substantially into negative territory after three years; for CPI, the effect disappears entirely after two years; and for IP, the decline is deeper and more persistent. Model averaging therefore improves on raw IV-LP by shrinking its most extreme movements, but it does not fully correct the long-horizon instability generated by the LP estimator and does not eliminate the puzzles.

Figure (ref) also makes clear why estimator averaging tends to deliver the most intuitive IRFs. The estimator-averaging weights on LP are high only on impact and in the very first few months, then quickly decay toward zero as the horizon increases, especially for IP and EBP. This pattern reflects the empirical fact that identification is strongest and LP variance is smallest at very short horizons, whereas the VAR's parametric structure provides more reliable long-run dynamics. The resulting estimator-averaged IRFs are therefore economically appealing: a front-loaded, temporary rise in GBY; a transitory fall in IP; a small and delayed disinflation in CPI; and a pronounced but short-lived increase in EBP. These responses sit between IV-LP and IV-VAR where the two disagree, dampen LP's long-horizon noise, and avoid the overly smooth extremes of VAR, illustrating the practical usefulness of estimator averaging in applied monetary policy analysis.

Conclusion

LP and VAR are the two workhorse methods for impulse response analysis, and their relative appeal is fundamentally a finite-sample question. LP is often attractive because of its robustness and comparatively low bias, whereas VAR is often attractive because of its greater precision. This paper studies how to combine these two estimators through horizon-specific estimator averaging, with weights chosen to minimize the mean squared error of the structural impulse response itself rather than the in-sample fit of the underlying regression.

We derive closed-form oracle weights that make transparent how the optimal LP share depends on the relative bias, variance, and covariance of LP and VAR, and we develop feasible AR-sieve-bootstrap implementations to estimate these weights in practice. The Monte Carlo results show that estimator averaging can deliver meaningful MSE gains relative to LP and VAR alone, especially because its inherent flexibility allows it to precisely track the horizon-specific dynamics of the bias--variance trade-off. In contrast, the fit-based model-averaging approach is less flexible and performs worse in our design.

In an empirical application revisiting the high-frequency IV monetary policy shocks of BauerSwanson2023, estimator averaging yields IRFs for yields, activity, prices, and credit spreads that are stable, economically intuitive, and lie between the often volatile IV-LP estimates and the very smooth IV-VAR estimates. Overall, our results suggest that estimator averaging provides a practical and easy-to-implement complement to existing LP and VAR practice, especially for empirical researchers who want to discipline the finite-sample bias--variance trade-off directly at the level of the impulse response of interest.

\singlespacing

\doublespacing