Unit Averaging
for Heterogeneous Panels

Christian Brownlees^† Vladislav Morozov^, *

In this work we introduce a unit averaging procedure to efficiently recover unit-specific parameters in a heterogeneous panel model.

The procedure consists in estimating the parameter of a given unit using a weighted average of all the unit-specific parameter estimators in the panel.
The weights of the average are determined by minimizing an MSE criterion.
We analyze	the properties of the minimum MSE unit averaging estimator     in a local heterogeneity framework inspired by the literature on frequentist model averaging.
The analysis of the estimator covers both the cases in which the cross-sectional dimension of the panel is fixed and large. In both cases we obtain the local asymptotic       distribution of the minimum MSE unit averaging estimators and of the associated weights.
A GDP nowcasting application for a panel of European countries showcases the benefits of the procedure.

Keywords: heterogeneous panels, frequentist model averaging, prediction

JEL: C33, C52, C53

[0]

$^{\dag}$ Department of Economics and Business, Universitat Pompeu Fabra and Barcelona School of Economics;

e-mail: .

$^{\ddag}$ Department of Economics and Business, Universitat Pompeu Fabra and Barcelona School of Economics;

e-mail: .

$^*$ Corresponding author.

We thank Jan Ditzen, Kirill Evdokimov, Geert Mesters, Luca Neri, Katerina Petrova, Barbara Rossi, Wendun Wang, and the participants at 26th International Panel Data Conference, 7th RCEA Time Series Workshop, 9th SIDE WEEE, and the 2021  Econometric Research in Finance Workshop for comments and discussion.
Christian Brownlees acknowledges support from the Spanish Ministry of Science and Technology (Grant MTM2012-37195)
and the Spanish Ministry of Economy and Competitiveness through the Severo Ochoa Programme for Centres of Excellence in R&D (SEV-2011-0075).

#1 1

§ INTRODUCTION

Estimation of unit-specific parameters in panel data models with heterogeneous parameters is a topic of active research in econometrics <cit.>.
Estimation of unit-specific parameters is relevant, for instance, when interest lies in constructing forecasts for the individual units in the panel <cit.>,
which typically arises in the analysis of international panels of macroeconomic time series <cit.>.   Other unit-specific parameters of interest include individual slopes <cit.> and long-run effects of a change in a covariate <cit.>.

There are three natural strategies for estimating unit-specific parameters <cit.>.
The simplest approach consists in estimating each unit-specific parameter from its individual time series.
While this strategy  typically leads to approximately unbiased estimation, the resulting estimators suffer from large estimation variability when the time dimension is small.
In the second approach, an assumption of parameter homogeneity is imposed and a common panel-wide estimator is used for all unit-specific parameters.
This strategy leads to small variability; however, it suffers from large bias in the presence of heterogeneity.
The third strategy is a compromise between the first two. It combines pooled and individual approaches
in the attempt to obtain an estimator with favorable risk properties <cit.>.
This is appealing when the time dimension is moderate and there is a nontrivial bias-variance trade-off between individual-specific and panel-wide estimation.

In this paper we propose a novel unit-specific compromise estimator that we call the unit averaging estimator.
The estimator is fairly general and is designed for possibly nonlinear panel models estimated by M-estimation.
We are concerned with the estimation of a unit-specific “focus” parameter.
Focus parameters considered include the examples mentioned above as well as other smooth transformations of the unit-specific parameter vector.
The unit averaging estimator for the unit-specific focus parameter is then defined as a weighted average of all the unit-specific focus parameter estimators in the panel.
The weights of the average are chosen by minimizing a unit-specific MSE criterion.
Specifically, we introduce two MSE criteria that differ in whether the cross-sectional dimension of the panel is treated as fixed (fixed-$N$) or large (large-$N$).
The minimum MSE weights solve a quadratic optimization problem that is straightforward to compute in both approaches.

We analyze the theoretical properties of the proposed unit averaging methodology.

The analysis is carried out using a  local heterogeneity assumption, in which we assume that individual coefficients are local in the time dimension to a common mean.
This is inspired by and builds upon the notion of local misspecification used in the frequentist model averaging literature <cit.>.
Local heterogeneity may be interpreted as a theoretical device to emulate panels where the time dimension is moderate.
In such a setting each unit carries information about the other units in the panel and thus may improve estimation accuracy.

Using the local heterogeneity framework we derive the local asymptotic MSE of the unit averaging estimator as well as estimators for this quantity.

The minimum MSE weights  minimize the local asymptotic MSE estimators. As we show, these minimum MSE weights  minimize an appropriately defined notion of the asymptotic MSE contaminated by a noise component that we characterize explicitly.
Finally, we obtain the limiting distribution of the minimum MSE unit averaging   estimator similarly to <cit.>.

In a simulation study, we  assess the finite sample properties of the proposed methodology.
We compare our unit averaging estimator against the unit-specific estimator, the mean group estimator as well as alternative unit averaging estimators based on AIC weights, BIC weights and Mallows weights <cit.>.
The simulation study shows that the proposed methodology performs favorably relative to these benchmarks.
The improvement in MSE is strongest for those units whose individual parameter is sufficiently distant from both the mean of the distribution and the endpoints of the support.

An application to GDP nowcasting for a panel of European countries showcases the methodology <cit.>.

GDP prediction is a natural application of the unit averaging methodology since the literature documents both evidence of heterogeneity between countries and the benefits of pooling data <cit.>.

We find that unit averaging using minimum MSE weights improves prediction accuracy and that the magnitude of the improvement is larger for shorter panels.

This paper is related to different strands of the literature.
First, it is related to the literature on frequentist model averaging.
Important contributions in this area include <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, and <cit.>, among others.

<cit.> deal with model averaging estimators specifically tailored for panel models.
The main difference with respect to these  contributions is that we focus on averaging different units with the same model whereas these papers average different models.
Second, this paper is related to the literature on unit-specific estimation using compromise estimators.
Important contributions in this area include <cit.>, <cit.>, <cit.> and <cit.>.
The main difference with respect to these contributions is that we focus on a setting where the time dimension is moderate (as opposed to either large or small) and that we do not require strict exogeneity.
Moreover, the existing literature largely focuses on linear models <cit.> whereas our framework allows for a large class of nonlinear models.

The rest of the paper is structured as follows.
Section <ref> introduces the unit averaging methodology.
Section <ref> studies the theoretical properties of the procedure.
Section <ref> contains the simulation study.
Section <ref> contains the empirical application.
Concluding remarks follow in section <ref>.
All proofs are collected in the proof appendix. Additional results are collected in the online appendix, available from the authors' websites.

§ METHODOLOGY

We introduce our unit averaging methodology within the framework of a fairly general class of panel data models with heterogeneous parameters.
Let $\{ {\bz_{i\,t}} \}$ with $i=1,\ldots, N$ and $t=1,\ldots, T$ denote a panel where
$\bz_{i\,t}$ denotes a random vector of observations taking values in $\mathcal Z \subset \mathbb R^d$.

For each unit in the panel, we define the unit-specific parameter $\btheta_{i} \in \Theta \subset \mathbb R^p$ as

\begin{equation*}
\btheta_{i} = \argmax_{\btheta\in \Theta} \mathbb E \left(\dfrac{1}{T} \sum_{t=1}^T m(\btheta, \bz_{i\,t}) \right)~,
\end{equation*}

where $m : \Theta \times \mathcal Z \rightarrow \mathbb R $ denotes a smooth criterion function.

Our interest lies in estimating the unit-specific “focus” parameter $\mu(\btheta_{i})$ for a fixed unit $i$ with minimal MSE, where $\mu : \Theta \rightarrow \mathbb R$ denotes a smooth function (similarly to the setup in <cit.>).
For example, $\mu(\btheta_i)$ may denote a component of $\btheta_i$, the conditional mean of a response variable given the covariates, or the long-run effect of a covariate.
To simplify exposition and without loss of generality, we focus on the problem of estimating the focus parameter $\mu(\btheta_{1})$ for unit 1.
In this paper we consider the case in which the focus function $\mu$ takes values in $\mathbb R$. It is straightforward to generalize the framework to a focus function taking values in $\mathbb R^q$ for some $q >1$.

To estimate $\mu(\btheta_1)$ we consider the class of   unit averaging estimators given by

\begin{equation}\label{eqn:avgest}
\hat{\mu}(\bw) = \sum_{i=1}^N w_{i} \mu(\hat{\btheta}_i) ~,
\end{equation}

$\bw  = ( w_{i} ) $ is a $N$-vector such that $ w_{i} \geq 0$ for all $i$ and $\sum_{i=1}^N w_{i}=1$,
and $\hat{\btheta}_i$ for $i=1,\ldots,N$ is the unit-specific estimator given by

\begin{equation}\label{equation:unitSpecificEstimator}
\hat{\btheta}_i = \argmax_{\btheta\in \Theta} \dfrac{1}{T} \sum_{t=1}^T m(\btheta, \bz_{i\,t}) ~.
\end{equation}

The class of estimators in (<ref>) is fairly broad and contains a number of important special cases.
It includes the individual estimator of unit 1 \( \hat{\mu}_1 = \mu(\hat{\btheta}_1) \) and the mean group estimator

$\hat \mu_{MG} = N^{-1} \sum_{i=1}^N \mu(\hat{\btheta}_i)$.[An alternative mean group estimation approach consists in setting $\hat{\btheta}_{MG} = N^{-1}\sum_{i=1}^N \hat{\btheta}_i$ and defining $\hat{\mu}_{MG} = \mu(\hat{\btheta}^{MG})$. As follows from lemma <ref> and theorem <ref>, the two approaches have identical asymptotic properties in our setup. The two definitions are also numerically identical if $\mu$ is affine.]
Other important special cases include estimators based on smooth AIC/BIC weights <cit.> and Mallows weights <cit.>, as well as a Stein-type estimator <cit.> given by
${\mu}_\mathrm{Stein} = \lambda \mu( \hat \btheta_1 ) + (1-\lambda)  N^{-1}\sum_{i=1}^N \mu(\hat{\btheta}_i)$
where $\lambda \in [0,1]$.

The class of estimators in (<ref>) may be motivated by the following representation for the individual parameters $\btheta_i$.
Assume that $ \btheta_i = \btheta_0 + \bEta_i $ holds for each $i=1,\ldots,N$, where $\btheta_0$ is a common mean component and $\bEta_i$ is a zero-mean random idiosyncratic component.
In such a setup all units in the panel carry information on $\btheta_0$, and so all units may be useful for estimating $\btheta_1=\btheta_0+ \bEta_1$.

The vector of weights $\bw$ controls the balance between bias and variance of   estimator   (<ref>).
Assigning a large weight to unit 1 leads to low bias but may also lead to excessive variability.
Alternatively, assigning large weights to units other than unit 1 induces bias but may substantially reduce variability.

Such a trade-off is particularly relevant in a moderate-$T$ setting, defined  as the range of values of $T$ for which the variability of the individual estimators $\hat{\btheta}_i$ is of the same order of magnitude as $\bEta_i$.

In this work we introduce a weighting scheme called minimum-MSE unit averaging weights.
These weights seek to strike a balance between the bias and variance of the unit averaging estimator.
We introduce two unit averaging schemes which differ in whether the cross-sectional dimension $N$ is treated as fixed or large.

In both regimes the weights are chosen by minimizing an estimator of the  local asymptotic approximation to the MSE (LA-MSE) of the unit averaging estimator.
The LA-MSE provides an approximation to the moderate-$T$ MSE of the unit averaging estimator and is justified in detail in the next section.

In the fixed-$N$ regime we average over a fixed finite  collection of $\bar{N}$ units.
In this regime we allow all units to have a non-negligible weight. Let $\bw^{\bar{N}}=(w_i^{\bar{N}})$ be a $\bar{N}$-vector such that $w_i^{\bar{N}}\geq 0$ for all $i$ and $\sum_{i=1}^{\bar{N}} w_i^{\bar{N}}=1$.
The fixed-$N$ LA-MSE estimator associated with  $\bw^{\bar{N}}$ is given by

\begin{equation}\label{equation:fixedNLAMSEestimatorMethodology}
\widehat{LA\mhyphen MSE}_{\bar{N}}(\bw^{\bar{N}}) = \sum_{i=1}^{\bar{N}} \sum_{j=1}^{\bar{N}} w_{i}^{\bar{N}}[\hat{\bPsi}_{\bar{N}}]_{i\,j} w_{j}^{\bar{N}} ~,
\end{equation}

where $\hat{\bPsi}_{\bar{N}} \in \mathbb R^{\bar{N} \times \bar{N}}$   with entries
$[\hat{\bPsi}_{\bar{N}}]_{i\,i} = \nabla \mu(\hat{\btheta}_1)' (  T ( \hat{\btheta}_i -\hat{\btheta}_1 )( \hat{\btheta}_i -\hat{\btheta}_1 )' + \hat{\bV}_i ) \nabla \mu(\hat{\btheta}_1)$ and $[\hat{\bPsi}_{N}]_{i\,j} = \nabla \mu(\hat{\btheta}_1)' T ( \hat{\btheta}_i -\hat{\btheta}_1 )( \hat{\btheta}_j -\hat{\btheta}_1 )' \nabla \mu(\hat{\btheta}_1)$ when $i\neq j$,
and $\hat \bV_i$ is an estimator of the asymptotic variance of $\hat \btheta_i$.
We remark that \( \nabla \mu(\hat{\btheta}_1)'T (\hat{\btheta}_i -\hat{\btheta}_1 )( \hat{\btheta}_i -\hat{\btheta}_1 )'  \nabla \mu(\hat{\btheta}_1)\) and \( \nabla \mu(\hat{\btheta}_1)'\hat{\bV}_i \nabla \mu(\hat{\btheta}_1)\) are estimators of, respectively, the squared bias and variance of $\mu (\hat \btheta_i )$ as estimators of $\mu (\btheta_1 )$.
The fixed-$N$ minimum MSE weights are

\begin{equation} \label{equation:fixedNweights}
\hat{\bw}^{\bar{N}} = \argmin_{\bw\in\Delta^{\bar{N}}}   \widehat{LA\mhyphen MSE}_{\bar{N}}(\bw) ~,
\end{equation}

where $\Delta^{\bar{N}} = \curl{\bw\in \R^{\bar{N}}: \sum_{i=1}^{\bar{N}} w_i=1, w_i\geq 0, i=1,\dots, \bar{N} }$.

In the large-$N$ regime we average an arbitrarily large collection of units, and we model $N$   as diverging to infinity.
When $N$ is large, some units are mechanically constrained to have a small weight since weights are non-negative and sum to unity.
Accordingly, in our framework we partition units into two sets, a set of $\bar{N} \geq 0$ unrestricted units that are allowed to have a non-negligible weight, and a set of remaining $N-\bar{N}$ units that are restricted to have a negligible weight in the limit.

In the large-$N$ regime we show that the LA-MSE is determined   by the weights assigned to the $\bar{N}$ unrestricted units.

Let $\bw^{N,\infty}=(w_{i}^{N,\infty})$ be an $N$-vector and assume, without loss of generality, that the weights of the unrestricted units are placed in the first $\bar{N}$ positions.
The vector of weights $\bw^{N,\infty}$ is such that $w_i^{N,\infty}\geq 0$ for all $i$, $\sum_{i=1}^N w_i^{N,\infty}= 1$, and the weights of the restricted units  ($i > \bar N$) are given by \( w_i^{N,\infty} = ( 1 - \sum_{j=1}^{\bar N} w^{N,\infty}_i )/( N - \bar N) \).

We remark that other weighting schemes are allowed for the restricted units. However the negligibility restriction implies that any admissible sequence of weights leads to the same LA-MSE and thus, for simplicity, we opt for equal weights here.

Let $\bw^{\bar{N},\infty}=(w_i^{\bar{N},\infty})$ be a $\bar{N}$-vector such that $w_i^{\bar{N},\infty} = w_i^{N,\infty}$ for $i=1,\ldots,\bar N$.
The large-$N$ LA-MSE estimator associated with $\bw^{N,\infty}$ is controlled by $\bw^{\bar{N},\infty}$ and  given by

\begin{align}
& \widehat{LA\mhyphen MSE}_{\infty}(\bw^{\bar{N}, \infty}) \label{equation:largeNLAMSEestimatorMethodology} \\
& = \sum_{i=1}^{{\bar{N}}} \sum_{j=1}^{{\bar{N}}} w_{i}^{\bar{N}, \infty} [\hat{\bPsi}_{{\bar{N}}}]_{i\,j} w_{j}^{\bar{N}, \infty} + \left[   \left(1-\sum_{i=1}^{\bar{N}}w_{i}^{\bar{N}, \infty} \right) \left( 	\sqrt{T}\nabla \mu(\hat{\btheta}_1)'\left(\hat{\btheta}_1- \frac{1}{N} \sum_{i=1}^N \hat{\btheta}_i  \right)\right)
\right.
\\ & \quad
\left. - 2\sum_{i=1}^{\bar{N}}w_{i}^{\bar{N}, \infty}\nabla \mu(\hat{\btheta}_1) \sqrt{T}\left(\hat{\btheta}_i-\hat{\btheta}_1\right)  \right]
\left( 1-\sum_{i=1}^{\bar{N}}w_{i}^{\bar{N}, \infty} \right)\hspace{-4pt}
\left( 	\sqrt{T}\nabla \mu(\hat{\btheta}_1)'\left(\hat{\btheta}_1- \frac{1}{N} \sum_{i=1}^N \hat{\btheta}_i  \right) \hspace{-4pt}
\right).
\end{align}

Notice that the restricted units may have an impact on the bias but not the variance of the unit averaging estimator.

The large-$N$ minimum MSE weights $\hat{\bw}^{N, \infty} = (\hat{w}^{N, \infty}_i)$ are given by

\begin{equation}\label{equation:largeNweights}
\hat{w}^{N, \infty}_i =
\begin{cases}
\hat{w}^{\bar{N}, \infty}_i & i \leq \bar N \\
{\left(1-\sum_{j=1}^{\bar N} w^{\bar{N}, \infty}_j  \right)(N - \bar N )^{-1}} & i > \bar N
\end{cases}
\end{equation}

\[
\hat{\bw}^{\bar{N}, \infty} = \argmin_{\bw\in\tilde{\Delta}^{\bar{N}}} \widehat{LA\mhyphen MSE}_{\infty}(\bw)
\]

with  $\tilde \Delta^{\bar{N}}= \curl{\bw\in \R^{\bar{N}}: w_i\geq 0, \sum_{i=1}^{N} w_i\leq 1}$.
It is important to emphasize that the optimization problem defining $\hat{\bw}^{\bar{N}, \infty}$ is $\bar N$-dimensional and can be solved by standard quadratic programming methods.

Two remarks are in order before we proceed.
First, the fixed- and large-$N$ regimes cover the cases of practical importance.[We also derive  LA-MSE in case $N$ is arbitrarily large and an infinite number of units have a non-zero weight in the limit. However, this case appears to be of limited interest in practice, and we do not study properties of the data-dependent weights in this case.]     If each unit is potentially important and we do not wish to restrict any weight, then we can apply the fixed-$N$ approximation. The fixed-$N$ regime is agnostic in this sense. Alternatively, if some units can only make an individually negligible contribution to the average, we can apply the large-$N$ regime. Using the large-$N$ approximation requires choosing $\bar{N}$. In principle, $\bar{N}$ can be chosen arbitrarily, though choosing $\bar{N}>N$  effectively results in a fixed-$N$ approximation. An appropriate  choice of $\bar{N}$  might be implied by economic logic. For example, in a macroeconomic application using country-level data,  we might order the countries so that immediate neighbors of the country of interest, its  important trading partners, and the large economies of the world are placed in the first $\bar{N}$ positions. The other countries are judged to contribute relatively little individually, and are placed in positions beyond $\bar{N}$.

Second, the fixed- and large-$N$ LA-MSE estimators have the appealing property of being applicable  both   when the amount of time series information in the panel is moderate or large.
When the amount of time series information is moderate, the LA-MSE approximates the infeasible population problem of minimizing the MSE, along with uncertainty about individual parameters (see the discussion following theorem <ref>).
When the amount of time series information is large, the bias term in the MSE dominates, and the unit averaging estimator based on the minimum MSE weights converges to the individual estimator $\mu(\hat{\btheta}_1)$.

§ THEORY

§.§ Assumptions

Our focus is on a moderate-$T$ setting. To emulate it and the trade-off between unit-specific and panel-wide information, we introduce what we call the local heterogeneity assumption.
This is inspired and is analogous to the local misspecification assumption used in the frequentist model averaging literature <cit.>.

[Local Heterogeneity]
The sequence of unit-specific parameters $\{ \btheta_i \}$ is such that

\begin{equation}
\btheta_{i} = \btheta_0 + \frac{\bEta_i}{\sqrt{T}} ~,
\end{equation}

where $\{\bEta_i\}$ is
a sequence of i.i.d. random vectors  such that $\E_{\bEta}( \bEta_i )=\bm 0$ and $\E_{\bEta}\norm{\bEta_{i} }^{12} < \infty$.[Here and below $\norm{\cdot} $ means the 2-norm, unless specifically labeled otherwise.] All analysis is done conditional on $\sigma(\bEta_1, \bEta_2, \dots)$ and all distributional statements below are conditional on $\sigma(\bEta_1, \bEta_2, \dots)$ unless specifically stated otherwise.

A number of remarks on this assumption are in order.

First, the scaling by $\sqrt{T}$ allows us to approximate a moderate-$T$ setting using asymptotic theory techniques by creating a nontrivial asymptotic bias-variance trade-off.
Intuitively, as $T$ becomes larger, the signal strength becomes proportionally weaker so that the amount of information in each time series remains constant. This is a standard technique to approximate finite-sample properties of model selection and averaging estimators (see, for example, <cit.>).
Second, since the focus of the analysis lies in recovering individual effects, all probability statements are implicitly conditional on $\sigma(\bEta_1, \bEta_2, \dots)$. Such a conditioning is natural, given our focus on estimating $\mu(\btheta_1)$ and typical when individual parameters are of interest <cit.>.
Importantly, all the results we establish are shown to hold with $\bEta$-probability 1, that is, for almost any realization of $\{ \bEta_i \}$.

In this paper we assume that the cross-sectional units are independent.

For each $i, j_1, \dots, j_k, k$ such that $i\neq j_1, \dots, j_k$ $\{ \{ \bz_{i\,t} \}_{t=0}^\infty , \bEta_i \}$ and $\curl*{\{ \{ \bz_{j_1\,t} \}_{t=0}^\infty , \bEta_{j_1} \}, \dots, \{ \{ \bz_{j_k\,t} \}_{t=0}^\infty , \bEta_{j_k} \} }$ are independent.

The unit-specific estimators are assumed to satisfy a number of regularity conditions.

[Individual Objective Function]

* The parameter space $\Theta$ is convex.

* The function $m(\btheta, \bz): \Theta\times \mathcal{Z}\to \R$ is twice continuously differentiable in $\btheta$ for each value of $\bz$. $m(\btheta, \bz)$ is measurable as a function of $\bz$ for every value of $\btheta$.

* There exists a positive finite constant $T_0$ (which does not depend on $i$) such that for all $i$ and $T>T_0$ it holds that
the unit-specific estimator satisfies $\hat{\btheta}_i \in \mathrm{int}( \Theta )$ a.s..

* The gradient of the unit-specific objective function satisfies
\[ \frac{1}{\sqrt{T}}\sum_{t=1}^T \nabla m(\btheta_i, \bz_{i\,t})\Rightarrow N(0, \bSigma_i) ~,\]

where $ \bSigma_i  =\lim_{T\to\infty} T^{-1}\sum_{t=1}^T \E\left[ \left( \sum_{t=1}^T \nabla m(\btheta_i,\bz_{i\,t}) \right)\left( \sum_{t=1}^T \nabla m(\btheta_i,\bz_{i\,t}) \right)'\right]$.

* There exist a positive finite constant $C_{\nabla m}$ (which does not depend on $i$ or $T$) such that, for all $i$ and all $T>T_0$ and for some $\delta \in (0,1)$, it holds that

\[
\E\norm{  \frac{1}{\sqrt{T}}\sum_{t=1}^T \nabla  m({\btheta}_i, \bz_{i\,t})  }^{2(1+\delta)} \leq  C_{\nabla m} ~.
\]

* The Hessian of the unit-specific objective function satisfies

\[
\sup_{ \btheta \in [\btheta_i, \hat{\btheta}_i]} \norm{ \dfrac{1}{T}\sum_{t=1}^T \nabla^2  m({\btheta}, \bz_{i\,t})  - \bH_i }  \xrightarrow{p} 0 ~,
\]

where $ \bH_i = \lim_{T\to\infty}\E(T^{-1} \sum_{t=1}^T \nabla^2 m( \btheta_i, \bz_{i\,t}))$.

* Let $ D_{i\, T} = \sup_{ \btheta \in [\btheta_i, \hat{\btheta}_i]}  \norm{ \left(T^{-1}\sum_{t=1}^T \nabla^2 m( {\btheta}, \bz_{i\, t}) \right)\bH^{-1}_i-\bI   }_{\infty}$.
$D_{i\, T}<1$ a.s. for all $i$ and all $T>T_0$.
There exists a positive constant $C_{\nabla^2 m}$ such that, for all $i$ and all $T>T_0$ and for $\delta$ as in (v), it holds that

\begin{equation}
\E\left[\left( \dfrac{D_{i\, T}}{1-D_{i\, T}} \right)^{\frac{2(2+\delta)(1+\delta)}{\delta} }\right]  \leq C_{\nabla^2 m}.
\end{equation}

* The matrices $\bSigma_i$ and $\bH_i$ satisfy
$\underline{\lambda}_{\bSigma} \leq \lambda_{\min}( \bSigma_i ) \leq \lambda_{\max}(\bSigma_i) \leq \overline \lambda_{\bSigma} $ and
$\underline \lambda_{\bH} \leq \lambda_{\min}( \bH_i ) \leq \lambda_{\max}(\bH_i) \leq \overline \lambda_{\bH} $ where
$\underline \lambda_{\bSigma}$, $\overline \lambda_{\bSigma}$, $\underline \lambda_{\bH}$ and $\overline \lambda_{\bH}$ are positive constants that do not depend on $i$.

* Let $\bV_i = \bH_i^{-1} \bSigma_i \bH_i^{-1}$.
Then, there is a sequence of estimators $\{ \hat{\bV}_i \}$ such that, for all $i$, $\hat{\bV}_i$ is consistent for $\bV_i$, and, for all $T>T_0$, $\lambda_{\min}(\hat{\bV}_i) > 0$ holds almost surely.

Assumption <ref> requires the unit-specific estimators to be consistent, asymptotically normal and to satisfy a number of regularity conditions.

We remark that this assumption allows for a fair amount of dependence and heterogeneity in the unit-specific observations.[A classic reference on M-estimation for dependent and heterogeneous data is the book by <cit.>.]

Assumption <ref>$(iii)$ states that the unit-specific estimator lies in the interior of the parameter space almost surely.

If the problem is linear or  defined by a convex smooth objective function and continuous covariates, the parameter space can be taken to be $\R^p$, and the condition holds automatically.

Assumption <ref>$(iv)$ is standard in the M-estimation literature, it requires the gradient of the objective function evaluated at $\btheta_i$ to satisfy a CLT.

Assumption <ref>$(v)$ is a moment condition on the gradient of the objective function.
In an i.i.d. setting such an assumption translates into a moment condition on the individual gradients.
More generally, this would be implied by appropriate moment and dependence assumption on the individual gradients.

Assumption <ref>$(vi)$ is also standard in the M-estimation literature; it requires the Hessian to satisfy a uniform law of large numbers.

Assumption <ref>$(vii)$ effectively  requires that the sample Hessian is nonsingular in a small enough neighborhood of $\btheta_i$. In a scalar problem,   $(vii)$ restricts the possible range of the second derivative as $\btheta$ ranges over a shrinking interval around $\btheta_i$. In addition, $(vii)$ places an assumption on the moments of   deviation from the population limit Hessian.

In case of linear regression, the sample and population Hessians do  not depend on the slope parameters and $(vii)$ is an assumption on moments of covariates.

Assumption <ref>$(viii)$ restricts the spectrum of the matrices $\bSigma_i$ and $\bH_i$. In particular, it implies a uniform restriction on the asymptotic variance $\bV_i$ of the individual estimators.

Assumption <ref>$(ix)$ states that there exists a sequence of nonsingular estimators $\{ \hat \bV_i \}$ for the asymptotic variance-covariance matrix of the individual estimator.

We remark that Assumptions <ref>$(iii)$ and $(vii)$ state that the sequence of unit-specific estimation problems satisfies approprite uniformity conditions.
Such conditions allow us to distill the key arguments relevant to our averaging theory and, in a sense, should be intrepreted as a simplifying approximation.
In general, $(iii)$ and $(vii)$  would hold with probability approaching one for each unit.
In this case all our results would still hold, though under appropriate rate conditions on $(N, T)$ and trimming to ensure certain well-behavedness of individual estimators.
We further note that assumptions $(iii)$ and $(vii)$ might hold in practice in certain special cases  regardless (such as linear or nonlinear models with a convex and smooth objective function and continuous covariates).

[Unit-specific Bias]
There exists a constant $C_{Bias}$, which does not depend on $i$, such that $\norm{  \E(\hat{\btheta}_i-\btheta_i) }_1 \leq C_{Bias}/T$ for all $T>T_0$.

Assumption <ref>  requires  that the bias of individual estimators for  their own parameters is bounded uniformly in $i$.

The order of the bias is  consistent with the results obtained by <cit.> and <cit.>. The  higher order terms can be subsumed into the $T^{-1}$ term for a sufficiently large $C_{Bias}$.[Explicitly including terms of order $-3/2$ and higher   does not change the analysis, as long as all the constants do not depend on $i$ ]
Assumption <ref>  is satisfied for linear models under assumption <ref>. For nonlinear models it is sufficient that for all $s$ and $i$ it holds that $ T\E(\norm{\nabla^s m(\btheta_i, z_{i\,t})}^2)\leq C_s<\infty$ <cit.>.

[Focus Parameter]
The focus function $\mu:\Theta\to \R$ is twice-differentiable.
There exists a constant $C_{\nabla \mu}$ such that for all $\btheta \in \Theta$ it holds that $\norm{\nabla \mu(\btheta)}<C_{\nabla \mu}$.
There exists a constant $C_{\nabla^2 \mu}$ such that for all $\btheta \in \Theta$ the largest and smallest eigenvalues of the Hessian $\nabla^2 \mu(\btheta)$ are bounded in absolute value by $C_{\nabla^2 \mu}$.
Let $\bd_0 = \nabla \mu(\btheta_0)$ be the gradient of $\mu$ at $\btheta_0$. Then $\bd_0 \neq 0$.

Assumption <ref> lays out mild smoothness assumptions on $\mu$.
For simplicity we assume that $\mu$ is a scalar focus parameter. However, all our results can be extended to the case in which $\mu$ is a vector focus parameter.

§.§ Asymptotic Properties of the Minimum MSE Unit Averaging Estimator

We begin with an auxiliary lemma that establishes the distribution of the unit-specific estimators in the local asymptotic framework of assumption <ref>.

Assume that assumptions <ref>–<ref> are satisfied. Let the unit-specific estimators  $\hat{\btheta}_i$ for $i=1,2,\ldots$ be defined as in eq. (<ref>).

\begin{align}
\sqrt{T}\left(\hat{\btheta}_i-\btheta_1 \right)	 &  \Rightarrow N(\bEta_i-\bEta_1, \bV_i) \eqqcolon \bZ_i ~, \\
\sqrt{T}\left(\mu(\hat{\btheta}_i) -\mu(\btheta_1) \right) & \Rightarrow N(\bd_0'\left(\bEta_i-\bEta_1\right), \bd_0'\bV_i\bd_0) \eqqcolon \Lambda_i
\end{align}

holds as $T\to\infty$ for $i=1,2,\ldots$. Convergence is joint (that is, with respect to the product topology), and all $\bZ_i$ and $\Lambda_i$ are independent across $i$.

Lemma <ref> characterizes the local asymptotic distribution of $\mu(\hat{\btheta}_i)$. In particular, the mean and the  variance of the limit distribution provide a local asymptotic approximation to the exact moderate-$T$ bias and variance  of $\mu(\hat{\btheta}_i)$ as an estimator of $\mu( {\btheta}_1)$. By adding together the square mean and the variance of the limit distribution,  we obtain a local  asymptotic approximation  to the MSE (LA-MSE) of each $\mu(\hat{\btheta}_i)$.

We now introduce a local asymptotic approximation to the MSE (LA-MSE) of the unit averaging estimator (<ref>).
Let $\{ \bw_N \} = \{ \bw_1, \bw_2 , \ldots \}$ be a (non-random) sequence where $\bw_k$ is a $k$-vector of weights. Suppose that $\bw_N$ converges to some $\bw \in \R^{\infty}$ in the sense defined below.  Consider the unit averaging estimator $\hat \mu( \bw_N )$ associated with such a sequence.
The following theorem derives $\lim_{N, T\to\infty} T\times MSE\left(\hat{\mu}(\bw_N) \right) $ where $N, T\to\infty$ jointly. This is a standard notion of risk in a local asymptotic framework (see e.g. <cit.>).
A remark on notation is in order before we proceed. In what follows we treat $\bw_N = (w_{i\,N})$ as   an element both in  $\R^{N}$ and in $\R^{\infty}$ (with coordinates $i>N$ restricted to zero).
This duality will not cause any confusion.

Let assumptions <ref>–<ref> be satisfied. Let $\{ \bw_N \}$ be such that
$(i)$ for each $N$, $\bw_N$ is measurable with respect to $\sigma(\bEta_1, \dots, \bEta_N)$,
$(ii)$ for each $N$,  $w_{i\,N}\geq 0$ for all $i$, $\sum_{i=1}^N w_{i\,N}=1$, $w_{j\,N}=0$ for $j>N$,
$(iii)$ $\sup_{i} \abs*{w_{i\,N}-w_i} = o(N^{-1/2})$ where $\bw=(w_i) \in \mathbb R^{\infty}$ is a vector such that $w_i\geq 0$ and $\sum_{i=1}^{\infty} w_i\leq1$. Let $T_0$ be as in assumption <ref>.

Then $(i)$ $\sum_{i=1}^{\infty} w_i\bd_0'\bEta_i$ and $\sum_{i=1}^{\infty}w_i^2\bd_0'\bV_i\bd_0$ exist;
$(ii)$ for any $N$ and $T>T_0$ the MSE of the averaging estimator is finite; and
$(iii)$  as $N, T\to\infty$ jointly it holds that

\begin{equation}
T\times MSE\left(\hat{\mu}(\bw_N) \right) \to  \left( \sum_{i=1}^{\infty} w_i \bd_0'\bEta_i-\bd_0'\bEta_1 \right)^2 + \sum_{i=1}^{\infty}w_i^2 \bd_0'\bV_i\bd_0. %\eqqcolon LA\text{-}MSE(\bw).
\end{equation}

Two remarks on theorem <ref> are in order before we proceed.
First, the theorem provides a local asymptotic approximation to the MSE of the averaging estimator (LA-MSE).  It highlights the bias-variance trade-off associated with the choice of the weights.
The two extremes of the trade-off correspond to the individual estimator of unit 1 $\mu(\hat{\btheta}_1)$ and the mean group estimator $\hat{\mu}_{MG}$.
The individual estimator of unit 1 is obtained by setting $w_{1\,N}=1$ for all $N$.
It is straightforward to see that $\mu(\hat{\btheta}_1)$ is asymptotically unbiased and that its LA-MSE is equal to $\bd_0'\bV_i\bd_0$, which is its asymptotic variance.

The mean group estimator is obtained by setting $w_{i\,N} = (N)^{-1}\I_{i\leq N}$ for $i=1,\ldots,N$ for all $N$.
$\hat{\mu}_{MG}$ has zero asymptotic variance and  its LA-MSE is equal to $(\bd_0'\bEta_1)^2$.[   This corresponds to the estimator of <cit.>. In a  genuine large-$T$ setting  it is feasible to estimate a bias of such of type and correct for it.  However, in a moderate-$T$ setting this quantity cannot be consistently estimated.]
Second, the theorem requires that the  sequence  $\{ \bw_N \}$ of weight vectors   converge uniformly to some limit $\bw$.[The assumption is needed to ensure convergence of the bias and the variance series and to prevent a “sliding hump” given by weighting structures like $w_{i\,N} = \mathbb{I}_{i=N}$.] In addition,   convergence must occur	 at a rate faster than $N^{-1/2}$.
Notice that this condition is trivially satisfied by the mean group estimator.
We also emphasize that the sum of the limit $\bw$ can be less than one.

In order to study the properties of the two averaging regimes introduced in section <ref>,
we provide a refined version of theorem <ref>.
It imposes  a stronger uniform convergence condition on the weights that reflects the assumptions of the large-$N$ regime.

Assume that assumptions <ref>–<ref> are satisfied.
Let $\{ \bw_N \}$ be such that
$(i)$ for each $N$, $\bw_N$ is measurable w.r.t. $\sigma(\bEta_1, \dots, \bEta_N)$,
$(ii)$ for each $N$, $w_{i\,N}\geq 0$ for all $i$, $\sum_{i=1}^N w_{i\,N}=1$, $w_{j\,N}=0$ for $j>N$,
$(iii)$ for some $\bar{N}\geq 0$ it holds that $\sup_{i> \bar{N}} w_{i\,N} = o(N^{-1/2})$, and
$(iv)$ $\curl{w_{i\,N}}_{i=1}^{\bar{N}}\to \curl{w_i}_{i=1}^{\bar{N}}$.

Then  as $N, T\to\infty$ jointly it holds that

\begin{align}
\sqrt{T}\left(	\hat{\mu}(\bw_N) - \mu(\btheta_1)\right)   & \Rightarrow   %N\left(  \sum_{i=1}^{\bar{N}} w_{i}\bd_0'(\bEta_i-\bEta_1) -\left(1-\sum_{i=1}^{\bar{N}}w_i \right)\bd_0'\bEta_1,   \sum_{i=1}^{\bar{N}} w_{i}^2 \bd_0'\bV_i\bd_0   \right)\\
%& =
N\left(\sum_{i=1}^{\bar{N}} w_{i}\bd_0'\bEta_i - \bd_0'\bEta_1, \sum_{i=1}^{\bar{N}}w_i^2\bd_0'\bV_i\bd_0\right)  .
\end{align}

Theorem <ref>    shows that the unit averaging estimator is asymptotically normal. The mean and variance are given by the weighted sums of biases and variances, respectively.

The theorem can be applied in both the fixed-$N$ and large-$N$ regimes introduced in Section <ref>.
In the fixed-$N$ regime, suppose only the first $\bar{N}<\infty$ units are being averaged (some of them potentially with zero weights). Then we set $w_{i\,N}=0$ for all $N$ and  $i>\bar{N}$, and  condition $(iii)$ holds automatically. The condition that $N\to\infty$ becomes superfluous. Conditions $(ii)$-$(iv)$ reduce to the requirement that  $\bw_N$ converge (pointwise) to a vector of weights $\bw$, where for $i>\bar{N}$ the weights satisfy $w_{i\,N}=w_i=0$ and $\sum_{i=1}^{\bar{N}} w_{i\,N}=\sum_{i=1}^{\bar{N}} w_i=1$. The limit distribution is normal with mean $\sum_{i=1}^{\bar N} w_i\bd_1'(\bEta_i-\bEta_1)$ and variance $\sum_{i=1}^{\bar N} w_{i}^2 \bd_1'\bV_i\bd_1$.
In the large-$N$ regime,  order the units with potentially large weights to be the first $\bar{N}$ units where $\bar{N}$ is chosen in advance.
Theorem <ref> shows that when $N$ is large, the units beyond $\bar{N}$ contribute no variance and approximate bias $-\left(1-\sum_{i=1}^{\bar{N}}w_i \right)\bd_0'\bEta_1$.

Importantly, the limit distribution of the averaging estimator does not depend on the actual shape of the weight vector beyond $\bar{N}$. The limit distribution is controlled  by finitely many weights, $\bar{N}$, and the total weight mass placed beyond $\bar{N}$. The weights in $\bw_N$ beyond $\bar{N}$ can display   strong variations in orders of magnitude, with some weights decaying like $N^{-1/2-\varepsilon}$, and some at a faster rate.
The total weight mass placed beyond  $\bar{N}$ in the limit is also not restricted, and may approach 1, like in the case of the mean group estimator.   As an application of this observation, if condition (iii) holds with $\bar{N}=0$ $\sqrt{T}\left(\hat{\mu}(\bw_N) - \mu(\btheta_1) \right)\xrightarrow{p} -\bd_0'\bEta_1$ for any weight vector sequence satisfying this condition.

Using theorem <ref>, we can obtain expressions for the population LA-MSE of the unit averaging estimator (<ref>) in the fixed-$N$ and large-$N$ regimes.  Consider the fixed-$N$ case in which we average $\bar{N}$  units. Let  $\bw^{\bar{N}}$ be a $\bar{N}$-vector, the limit vector of theorem <ref>.
In the fixed-$N$ regime the population LA-MSE is

\begin{equation}\label{equation:mseFixed}
LA\mhyphen MSE_{\bar{N}}(\bw^{\bar{N}}) = \sum_{i=1}^{\bar{N}} \sum_{j=1}^{\bar{N}} w_{i}^{\bar{N}} [{\bPsi}_{\bar{N}}]_{i\,j} w_{j}^{\bar{N}} ~,
\end{equation}

where $\bPsi_{\bar{N}}$ is an $\bar{N} \times \bar{N}$ matrix with elements
$[\bPsi_{\bar{N}}]_{i\,i} =\bd_0'\left((\bEta_i-\bEta_1)\left(\bEta_i-\bEta_1 \right)' + \bV_i\right)\bd_0 $  and
$[\bPsi_{\bar{N}}]_{i\,j}   = \bd_0'(\bEta_i-\bEta_1)\left(\bEta_j-\bEta_1 \right)'\bd_0$ when $ i\neq j$. Now consider the large-$N$ regime. As theorems <ref> and <ref> show, the LA-MSE is controlled by a $\bar{N}$-vector $\bw^{\bar{N}, \infty}$ such that $w_i^{\bar N,\infty} \geq 0$ for all $i$ and $\sum_{i=1}^{\bar{N}}w_{i}^{\bar{N}, \infty}\leq 1$.
In the large-$N$ regime, the population LA-MSE is

\begin{align}
& LA\mhyphen MSE_{\infty}(\bw^{\bar{N},\infty}) =  \sum_{i=1}^{\bar{N}} \sum_{j=1}^{\bar{N}} w_{i}^{\bar{N}, \infty} [{\bPsi}_{\bar{N}}]_{i\,j} w_{j}^{\bar{N}, \infty} \\
& \quad + \left( \left( 1-\sum_{i=1}^{\bar{N}}w_i^{\bar{N}, \infty}\right)\bd_0'\bEta_1 - 2\sum_{i=1}^{\bar{N}}w_i^{\bar{N}, \infty}\bd_0'(\bEta_i-\bEta_1)  \right) \left( 1-\sum_{i=1}^{\bar{N}}w_i^{\bar{N}, \infty}\right)\bd_0'\bEta_1~.
\end{align}

The quantities $\widehat{LA\mhyphen MSE}_{\bar{N}}$ and $\widehat{LA\mhyphen MSE}_{\infty}$ used to define the minimum MSE weights introduced in Section <ref>
are estimators of the population LA-MSE given above. In the rest of the section we focus on the properties of these estimators as well as the   optimal weights (<ref>) and (<ref>) associated with them.

We begin by noting that in our framework the population LA-MSE  cannot be consistently estimated. Under local heterogeneity the idiosyncratic components $\bEta_i$  cannot be consistently estimated <cit.>.
Instead, following <cit.>, we form $\widehat{LA\mhyphen MSE}_{\bar{N}}$ and $\widehat{LA\mhyphen MSE}_{\infty}$ by plugging in asymptotically unbiased estimators for  $\bEta_i-\bEta_1$ and $\bEta_1$.
Such estimators are provided by $\sqrt{T}(\hat{\btheta}_i-\hat{\btheta})$ and $\sqrt{T}(\hat{\btheta}_1-N^{-1} \sum_{i=1}^N \hat{\btheta}_i)$, respectively, as the following lemma establishes.

Let assumptions <ref>-<ref> hold.
Then as $N, T\to \infty$ jointly, it holds that

\begin{align}
\sqrt{T}\left( \hat{\btheta}_i -\hat{\btheta}_1 \right) &  \Rightarrow N\left(\bEta_i-\bEta_1, \bV_i+\bV_1 \right) = \bZ_i- \bZ_1, \\
\sqrt{T}\left(\hat{\btheta}_1- \frac{1}{N} \sum_{i=1}^N \hat{\btheta}_i  \right)  & \Rightarrow       N(\bEta_1, \bV_1) = \bZ_1  + \bEta_1 .
\end{align}

Convergence is joint for all $i$.

We remark that the matrix  $\hat{\bPsi}_{\bar{N}}$ of equations (<ref>) and (<ref>)   is a biased   estimator of $\bPsi_{\bar{N}}$.
Such a bias   ensures that $\widehat{LA\mhyphen MSE}$ is nonnegative for all admissible weight vectors.  An  asymptotically unbiased estimator  instead would have diagonal elements
$\tilde{\bPsi}_{i\,i} = \hat{\bd}'_1\left( T \left( \hat{\btheta}_i -\hat{\btheta}_1 \right)\left( \hat{\btheta}_i -\hat{\btheta}_1 \right)' - \left(\hat{\bV}_i+\hat{\bV}_1 \right)\right)\hat{\bd}_1$ and   off-diagonal elements $\tilde{\bPsi}_{i\,j} = \hat{\bd}'_1\left( T \left( \hat{\btheta}_i -\hat{\btheta}_1 \right)\left( \hat{\btheta}_j -\hat{\btheta}_1 \right)' - \hat{\bV}_1 \right)	\hat{\bd}_1$ where   $\hat{\bd}_1=\nabla \mu\left(\hat{\btheta}_1 \right)$.
However, this can easily fail to be positive definite,

as it involves a difference of   positive definite matrices. This would lead to the undesirable possibility of obtaining negative estimates of the LA-MSE.[An additional argument in favor of focusing on $\hat{\bPsi}_{\bar{N}}$ is given by <cit.> who examines the behavior of both $\hat{\bPsi}_{\bar{N}}$ and $\tilde{\bPsi}_{\bar{N}}$ in the context of model averaging. <cit.> finds that $\hat{\bPsi}_{\bar{N}}$ leads to superior performance of resulting weights.]

The following  two theorems   establish the properties of our LA-MSE estimators  and the associated minimum MSE weights  (<ref>) and (<ref>).
The theorem also characterizes the asymptotic distribution of the minimum MSE unit averaging estimator $\hat{\mu}(\hat{\bw}^{\bar{N}})$  in this regime.
First, we state a result for the fixed-$N$ regime. Recall that $\Delta^{\bar{N}} = \curl{\bw\in \R^{\bar{N}}: \sum_{i=1}^{\bar{N}} w_i=1, w_i\geq 0, i=1,\dots, \bar{N} }$.

Let assumptions <ref>-<ref> hold.

* For any $\bw^{\bar{N}} \in \Delta^{\bar{N}}$ it holds that
\( \widehat{LA\mhyphen MSE}_{\bar{N}}(\bw^{\bar{N}})\Rightarrow \overline{LA\mhyphen MSE}_{\bar{N}}(\bw^{\bar{N}}) \coloneqq \bw^{\bar{N}'}\overline{\bPsi}_{\bar{N}}\bw^{\bar{N}} \)
as $T\to\infty$,
where ${\overline{\bPsi}}_{\bar{N}}$ is an $\bar{N}\times \bar{N}$ matrix with elements $[{\overline{\bPsi}}_{\bar{N}}]_{i\,j} = \bd_0'((\bZ_i-\bZ_1)(\bZ_i-\bZ_1)' +\bV_i)\bd_0$ when $i=j$ and $\bd_0'((\bZ_i-\bZ_1)(\bZ_j-\bZ_1)' )\bd_0$ when $i\neq j$; and $\bZ_i$ is as in Lemma <ref>.

* As $T\to\infty$, the minimum MSE weights satisfy

\begin{equation}
\hat{\bw}^{\bar{N}} = \argmin_{\bw^{\bar{N}}\in\Delta^{\bar{N}}}  \widehat{LA\mhyphen MSE}_{\bar{N}}(\bw^{\bar{N}}) \Rightarrow \overline{\bw}^{\bar{N}}=  \argmin_{\bw^{\bar{N}}\in\Delta^{\bar{N}}}  \overline{LA\mhyphen MSE}_{\bar{N}}(\bw^{\bar{N}}) .
\end{equation}

* As $T\to\infty$, for $\Lambda_i$ of lemma <ref>, the minimum MSE unit averaging estimator satisfies

\[ \sqrt{T}\left(\hat{\mu}(\hat{\bw}^{\bar{N}}) - \mu(\btheta_1) \right)\Rightarrow  \sum_{i=1}^{\bar N} \overline{w}_{i}^{\bar{N}} \Lambda_i.  \]

A number of remarks on theorem <ref> are in order.
First, $\overline{LA\mhyphen MSE}_{\bar{N}}$ plays the same role to $\widehat{LA\mhyphen MSE}_{\bar{N}}$  as $\bZ_i$ does to $\sqrt{T}(\hat{\btheta}_i-\btheta_1)$ in lemma <ref>. $\overline{LA\mhyphen MSE}_{\bar{N}}$ uses  a local approximation to express   $\widehat{LA\mhyphen MSE}_{\bar{N}}$ in terms of the components of the MSE and the approximate distribution of the individual estimators.

We can see that $\overline{LA\mhyphen MSE}_{\bar{N}}$ is composed of the population LA-MSE, a bias term, and a noise component.
In fact, entries of the matrix $\overline{\bPsi}_{\bar{N}}$   may be expressed as

\begin{align}%\label{key}
[\overline{\bPsi}_{\bar{N}}]_{i\, i} &
=[{\bPsi}_{\bar{N}}]_{i\, i}  + \bd_0'(\bV_1+\bV_i)\bd_0 + \bd_0'\be_{i\,i}\bd_0 ~,\\
[\overline{\bPsi}_{\bar{N}}]_{i\, j}  &  = [ {\bPsi}_{\bar{N}}]_{i\, j}  + \bd_0'\bV_1\bd_0 + \bd_0'\be_{i\,j}\bd_0, \quad i\neq j ~,
\end{align}

where $\be_{i\,j} = (\bZ_i-\bZ_1)(\bZ_j-\bZ_1)' - \E\left((\bZ_i-\bZ_1)(\bZ_j-\bZ_1)' \right)$.

The noise terms $\be_{i\,j}$ may be interpreted as the  result of the fact that in a moderate-$T$ setting there is limited  information about the idiosyncratic components $\bEta_i$. These terms are mean zero and independent conditional on unit 1.
The bias terms guarantee that $\overline{\bPsi}_{\bar{N}}$ is positive definite  and arise as a consequence of using the biased positive definite estimator $\hat{\bPsi}_{\bar{N}}$.  The bias  can be split into two components.
The $\bd'_0\bV_1\bd_0$ is common for all elements of $\overline{\bPsi}_{\bar{N}}$ and does not affect the solution of the MSE minimization problem.
The second component  $\bd_0'\bV_i\bd_0$ only affects the diagonal of $\overline{\bPsi}_{\bar{N}}$ by inflating the variance associated to each unit by the  variance of the corresponding individual estimator.
This component does not modify the ordering of the estimators in terms of their variances. If all the estimators were unbiased, optimal weights would be inversely proportional to individual variances, and this resulting ordering of the weights would be  preserved.

Second, the fixed-$N$ minimum MSE unit averaging estimator has a nonstandard asymptotic distribution in the local heterogeneity framework.
Assertion 3 of theorem <ref>  shows that the limit distribution is a randomly weighted sum of independent normal random variables.
In the online appendix, we show how to construct confidence intervals based on this result.

Third, minimizing  $\widehat{LA\mhyphen MSE}_N$ is natural even in a non-local setting where we drop assumption
<ref> and allow the amount of information in each time series to grow as $T\to\infty$.  For all $i$ such that $\btheta_i\neq \btheta_1$, the bias estimators $\sqrt{T}(\hat{\btheta}_i-\hat{\btheta}_1)$  will diverge, while all variance terms remain bounded. Then the procedure will place  zero weight asymptotically on all such units $i$ and on the tail term.
The minimum MSE weights will only use those units that share the same $\btheta_1$.
In particular, if  the distribution of $\bEta$ is continuous, asymptotically the procedure will use only information from unit 1, and the difference between the averaging estimator with   minimum MSE weights and the individual estimator will converge to zero in probability. In this case the averaging estimator is asymptotically normal with mean zero and variance $\bd_0' \bV_1 \bd_0$.
Such a result has a parallel in fixed parameter asymptotics for model averaging. As <cit.>  show in a recent paper, the weight assigned to the true (or least wrong) model tends to one as sample size increases.

The following result establishes an analogous result for the minimum MSE weights (<ref>) in the large-$N$ regime.  Recall that  $\tilde \Delta^{\bar{N}}= \curl{\bw\in \R^{\bar{N}}: w_i\geq 0, \sum_{i=1}^{N} w_i\leq 1}$.

Let assumptions <ref>-<ref> hold.

* For any $\bw^{\bar{N}, \infty}\in  \tilde{\Delta}^{\bar{N}}$
it holds that
$\widehat{LA\mhyphen MSE}_{\infty}(\bw^{\bar{N}, \infty}) 	\Rightarrow \overline{LA\mhyphen MSE}_{\infty}(\bw^{\bar{N}, \infty}) $
as $N, T\to \infty$ jointly

\begin{align}
\overline{LA\mhyphen MSE}_{\infty}(\bw^{\bar{N}, \infty}) & = \bw^{\bar{N}, \infty'} \overline{\bPsi}_{\bar{N}}\bw^{{\bar{N}, \infty}} +  \Bigg[   \left(1-\sum_{i=1}^{\bar{N}}w^{\bar{N}, \infty}_i \right) \bd_0'\left(\bEta_1+ \bZ_1\right) \\
& \quad 	- 2\sum_{i=1}^{\bar{N}}w^{\bar{N}, \infty}_i\bd_0'  \left(\bZ_i-\bZ_1\right)  \Bigg]
\left( 1-\sum_{i=1}^{\bar{N}}w_i^{\bar{N}, \infty} \right)
\bd_0'\left(\bEta_1 + \bZ_1
\right).
\end{align}

* As $N, T\to\infty$,  the   minimum MSE weights satisfy
\begin{equation}
\hat{\bw}^{\bar{N}, \infty} = \argmin_{\bw\in\tilde{\Delta}^{\bar{N}}}  \widehat{LA\mhyphen MSE}_{\infty}(\bw)\Rightarrow \overline{\bw}^{\bar{N},\infty} =  \argmin_{\bw\in\tilde{\Delta}^{\bar{N}}} \overline{LA\mhyphen MSE}_\infty(\bw) .
\end{equation}

* Let $\bv_{N-\bar N}=(v_{\bar{N}\, N}, \dots, v_{N\, N})$ be a $(N-\bar{N})$-vector such that $\sup_{i} v_{i\, N-\bar N}=o(N^{-1/2})$, $v_{i\,N-\bar{N}}\geq 0$, for each $N$ it holds that $\sum_{i=N-\bar{N}}^{N} v_{i\,N - \bar N}=1$.
Then  as $N, T\to\infty$ jointly

\begin{align}
& \sqrt{T}\left(\sum_{i=1}^{\bar{N}}\hat{w}_{i}^{\bar{N}, \infty} \mu\left(\hat{\btheta}_i \right) + \left(1-\sum_{i=1}^{\bar{N}}\hat{w}_{i}^{\bar{N}, \infty} \right)\sum_{j=N-\bar{N}}^{N} v_{j\,N-\bar N}\mu(\hat{\btheta}_j)  - \mu(\btheta_1) \right)\\ &\Rightarrow
\sum_{i=1}^{\bar{N}} \overline{w}_{i}^{\bar{N},\infty} \Lambda_i - \left(1-\sum_{i=1}^{\bar{N}}\overline{w}_{i}^{\bar{N},\infty} \right)\bd_0'\bEta_1  . \label{equation:optimalWeightsInfinite}
\end{align}

Observe that the estimator in equation (<ref>) is a valid averaging estimator, as the weights  sum to unity by construction.

The exact way $\bv_N$ is picked does not matter, as long as the decay condition holds. All admissible choices lead to the same limit. In particular,  we may pick equal weights $v_{i\,N} = 1/(N-\bar{N})$, as we do in eq. (<ref>).

§ SIMULATION STUDY

We study the finite sample performance of our minimum MSE unit averaging methodology via a simulation exercise.
We consider a model similar to the one we use in our empirical application – a linear dynamic heterogeneous panel model defined as

\begin{align}
y_{i\,t} & = \beta_i x_{i\,t} + \lambda_{i} y_{i\,t-1} + u_{i\,t}~, \hspace{1em} \varepsilon_{i\,t} \stackrel{i.i.d.}{\sim} N( 0 , \sigma^2_i )~,
\end{align}

where we assume $\E(\varepsilon_{i\,t}|y_{i\,t-1}, x_{i\,t})=0$. The prediction error $\varepsilon_{i\,t}$ is cross-sectionally heteroskedastic, with variance $\sigma^2_i$ drawn independently from  a standard exponential distribution.

The exogenous variable $x_{i\,t}$ is independently drawn from from a standard normal distribution.

The initial conditions for $y_{i\,0}$ are independently drawn from a normal distribution with mean zero and variance $(1+\sigma^2_i)/(1-\lambda_i^2)$ to ensure that $\{ y_{it} \}_t$ is covariance stationary.

The heterogeneous parameter $\btheta_i=(\beta_i,\lambda_i)'$ is specified as $\beta_i = 1 + \eta_{i\, \beta}/\sqrt{T}$ and $\lambda_i = \eta_{i\, \lambda}/\sqrt{T}$.
The sequences of idiosyncratic components $\{ \eta_{i,\, \beta} \}$ and $\{ \eta_{i\, \lambda} \}$ are independently drawn from, respectively, a standard normal distribution and a uniform distribution with support $[-4,4]$.

We simulate panels with a cross-sectional dimension $N$ of $10, 25$, and $50$ units and a time dimension $T$ of 60 periods. We remark that $T=60$ matches the time dimension of the estimation sample in the empirical application. As $T=60$, the support of $\lambda_i$ is approximately $[-0.5, 0.5]$.

We consider two focus parameters of interest:  $\lambda_{1}$ and the forecast of $y_{1\, T+1}$  given by the conditional mean $\E(y_{1\,T+1}|y_{1\,T}, x_{1\,T}=1)= \lambda_1 y_{1\,T}+ \beta_1$.[The value of $y_{1\,T}$ is determined in each sample by the dynamic process governing $y$.][In the online appendix, we also report results for estimating $\beta_1$ and the long-run effect $\beta_1/(1-\lambda_1)$.]

The minimum MSE unit averaging estimators corresponding to these focus parameters are used for estimation.

We consider weights based on both the fixed-$N$ and large-$N$ regime. In the large-$N$ regime we consider weights based on the choices $\bar{N}=10$ and $\bar{N}=20$.[When the panel cross-sectional dimension is $N=10$ all three approximations are identical.]

The performance of the minimum MSE unit averaging estimator is benchmarked against a number of alternative approaches.
We consider the individual estimators of unit 1, the mean group estimator, as well as the unit averaging estimator based on AIC/BIC weights <cit.> and MMA weights <cit.>.  Note that AIC and BIC generate the same weights, since each unit has the same number of coefficients.[More formally, our weights  effectively correspond to conditional AIC weights <cit.> where we treat the random slopes as fixed parameters of interest. See also <cit.>] Similarly, the MMA weights reduce to minimal BIC model selection.

MSE performance of the unit averaging estimator relative to the individual estimator. The focus parameters are the autoregressive parameter $\lambda_1$ and conditional mean (forecast). $T=60$, $N=10, 25, 50$

Figure <ref> summarizes the main results of the simulation study based on 2500 replications.
We report performance of different averaging strategies in terms of MSE relative to the individual estimator for a range of values of $\lambda_1$. As the distribution of $\lambda_1$ is symmetric, we report performance for half of the support of $\lambda_1$.

In each plot, we plot the MSE of different averaging estimators relative to the MSE of the individual estimator of unit 1. A value larger than 1 indicates that the individual estimator is more efficient.

The minimum MSE unit averaging estimator performs favorably throughout most  of the parameter space. Large-$N$ approximations work better for $\lambda_1$ closer to $\E(\lambda_1)=0$. However, for larger values of $\lambda_1$  the flexibility of the fixed-$N$ regime is an advantage, as it can freely choose units most similar to unit 1.

None of the averaging estimators dominate the other.
As expected, the mean group estimator performs very well close to the mean of $\lambda_{1}$, but bias starts to dominate as $\lambda_{1}$ is increased. AIC weights offers a risk profile somewhat similar to that of the MG estimator for values closer mean of $\lambda_{1}$. As $\lambda_1$ becomes larger, AIC moves towards the individual estimator.
As noted above, MMA is effectively minimal BIC model selection in the context of unit averaging. MMA almost always selects unit 1, and correspondingly MMA does not offer improvements over the individual estimator.

§ EMPIRICAL APPLICATION

We illustrate our methodology with a pseudo-real time nowcasting exercise for   quarterly GDP for a panel of European countries.
GDP prediction is a natural application of our unit averaging methodology.
There is evidence of considerable heterogeneity between countries, yet at the same time pooling the data at least partially improves prediction accuracy <cit.>.

The design of our application follows standard practices in the nowcasting literature <cit.>.

The literature on nowcasting is vast and we do not to cover it here. We refer to <cit.> for a survey.

We use quarterly GDP data  from 1995Q1 to 2019Q4 for 12 European countries: the 11 founding Eurozone economies and the UK.
We enrich our dataset with a set of 162 monthly GDP predictors for each country. The set of predictors include both real, price, and survey data.
Table OA.1 in the online appendix contains the complete list of variables and descriptions.[Not all variables are available for all countries at a given time. This only impacts the precision in estimating country-specific factors.]
All non-survey data is available from Eurostat whereas the survey data is available from the DG ECFIN.
We use final data releases incorporating all revisions, making our study a pseudo-real time one.

Our empirical design takes into account both  the
delays in publication of monthly data (“ragged-edge problem”) and  the impact of timing on the information set available (“vintages” of data).

First, the predictor  variables are typically released with different  delays  after the end of the corresponding month, which is known as the “ragged-edge” problem <cit.>.[For example, industrial production data is released 6 weeks after the end of the month, while  survey data is released at the end of the month without delay.] We adopt a stylized release calendar of bimonthly releases to account for this (table OA.1 in the online appendix lists the release delay for all variables).

Second,  as the quarter goes by, more data becomes available.[For example, nowcasting  Q4 GDP can be done at any moment between October 1 when no data on Q4 is available yet
up to the middle of the following February, when GDP data for Q4 is released. The amount of data available increases monotonically between these two dates.

Each possible position in time determines a data “vintage".  We assume that a month has 4 weeks; in accordance with our release calendar, we nowcast every two weeks starting from the first day of the quarter at $-12$ weeks (relative to the quarter end)
until $+4$ weeks after the end of the quarter (GDP is released at +6 weeks).
Formally,    let $t$ index months. Then $v = -3, -5/2,  -2, \ldots, 1/2, +1$ is a fractional value that  describes the monthly position (or vintage) relative to the end of the quarter, in increments of two weeks.[For example, if $v=0$, nowcasting uses all information that is available at the  end of quarter $3t$. If $v=+1$, nowcasting uses all the data available  +4 weeks after  the end of quarter $3t$, the last weekly position we consider. Each step of $-1/2$ corresponds to stepping back 2 weeks until $v=-3$ corresponds to the position of $-12$ weeks.]

We nowcast GDP in quarter $3t$  using all information available at time $3t+v$, separately for each value of $v \in \curl{-3, -5/2, \dots, +1}$.
As we have a large number of predictors available at monthly frequency, we opt for  factor unrestricted MIDAS (U-MIDAS) <cit.>.

Given $v$, for each country we estimate monthly factors $f_{i\,t}$ with  $\hat{f}_{i\,t|v}$ for all $t= 1, \dots, \floor{T+v}$ using the full dataset available at $T+v$.[For example, suppose we wish to nowcast Q4 GDP. If $v=1$, we estimate factors up to January of the following year using information available at the end of January. If $v=1/2$, we  estimate  factors up to  December  using all the information available in  the middle of January.]
The GDP is modeled as

\begin{equation*}
y_{i\,3t} = \alpha_{i|v} + \sum_{k=0}^{11} \beta_{i\,k|v} \hat{f}_{i\,\floor{3t+v-k}|v}+ \lambda_{i|v} y_{i\,3(t-1)} +  \varepsilon_{i\, 3t|v},
\end{equation*}

where $y_{i\,3t}$ is GDP of country $i$ in quarter $3t$  and $\varepsilon_{i \,3t|v}$ is the prediction error.

The country factors estimates $\hat{f}_{i\,t|v}$ are extracted from the large set of predictor variables  using the  EM-PCA method <cit.>. If GDP of quarter $3(t-1)$ is not available at moment $v$, we use $y_{i\, 3(t-2)}$ instead.[Quarterly GDP is released six weeks after the end of the relevant quarter, which corresponds to $v=-3/2$. For $v=-3, -5/2, -2$, we use $y_{i\, 3(t-2)}$ in place of $y_{i\,3(t-1)}$.]

We use only one factor for prediction following <cit.> and we include the lag of GDP following <cit.>.

We nowcast GDP for each country using the conditional mean of GDP implied by the U-MIDAS specification.
Parameter estimation is carried out using a rolling-windows of sizes 44, 60 and 76 quarters.[Forecast evaluation begins in 2006Q1 for window size 44, 2010Q1 for T=60 and 2014Q1 for T=76.]
Factors are also re-estimated every two weeks using the all the data available at each point in time.

We estimate the conditional mean using the fixed-$N$ minimum MSE unit averaging estimator, since cross-sectional dimension is not large and each unit is potentially relevant.

The performance of our minimum MSE unit averaging estimator is benchmarked against the individual, mean group, and averaging estimators using AIC and Mallows weights.

In table <ref> we provide  a summary of forecasting performance results for GDP nowcasting.
The table reports the MSE of the individual estimator as well as the MSE relative to the individual estimator for all other strategies.
The table reports results for the five largest economies in our sample along with the GDP-weighted mean.[Weighing by GDP as in <cit.> emulates forecasting the Eurozone GDP using individual forecasts.]
We select the vintages that correspond to $-6, 0, +4$ weeks relative to the end of the quarter (corresponding to $v=-3/2, 0, +1$).
Full results for all vintages and countries  are provided in table OA.7 in the online appendix,
and they are   similar to the ones reported here.

3c|$-6$ weeks     3c|$0$ weeks     3c$+4$ weeks

Averaging       44q      60q      76q      44q      60q      76q      44q      60q      76q

Mean     Individual     1.113     0.986     1.167     0.973     1.010     1.196     0.933     0.914     1.124

minMSE     0.916     0.936     0.907     0.889     0.936     0.910     0.881     0.928     0.901

AIC     0.933     0.962     0.980     0.908     0.960     0.974     0.878     0.949     0.955

MMA     1.119     1.178     1.099     1.011     1.036     1.114     1.000     0.921     0.916

Mean group     1.417     1.570     1.524     1.635     1.696     1.505     1.879     1.704     1.489

DE     Individual     0.661     0.546     0.537     0.509     0.421     0.434     0.565     0.449     0.456

minMSE     0.793${}^{*}$     0.822     0.815     0.787${}^{*}$     0.818     0.775     0.821${}^{*}$     0.809     0.794

AIC     0.963${}^{*}$     0.977${}^{*}$     0.989     0.974     0.982${}^{*}$     0.973     0.978     0.989${}^{*}$     0.978

MMA     1.002     1.257     1.098     0.860     0.970     0.742     0.741     0.830     0.834

Mean group     0.987     0.937     0.773     1.069     1.157     0.742     0.957     1.153     0.849

FR     Individual     0.194     0.154     0.129     0.143     0.100     0.086     0.155     0.121     0.098

minMSE     0.988     1.067     1.037     0.971     1.059     1.159     0.916     0.978     1.069

AIC     0.883${}^{*}$     0.934     0.975     0.833${}^{*}$     0.978     1.049     0.828${}^{*}$     0.935${}^{*}$     0.999

MMA     1.223     1.343${}^{*}$     1.139     1.082     1.372${}^{*}$     1.659${}^{*}$     1.164     1.048     1.337

Mean group     2.125${}^{*}$     2.068${}^{*}$     1.348     2.736${}^{*}$     2.942${}^{*}$     2.169${}^{*}$     2.473${}^{*}$     2.652${}^{*}$     2.156${}^{*}$

IT     Individual     0.591     0.253     0.156     0.279     0.178     0.116     0.232     0.131     0.082

minMSE     0.893${}^{*}$     0.908${}^{*}$     0.852${}^{*}$     0.973     0.974     0.858${}^{*}$     1.046     1.025     0.857${}^{*}$

AIC     0.945     0.972${}^{*}$     0.980     0.955     0.951${}^{*}$     0.976${}^{*}$     0.947     0.969${}^{*}$     0.975${}^{*}$

MMA     0.907     0.710${}^{*}$     0.729${}^{*}$     1.323     0.650${}^{*}$     0.704     1.351     0.917     0.688${}^{*}$

Mean group     0.895     0.901     0.822     1.289     1.042     0.719     1.491${}^{*}$     1.595${}^{*}$     1.239

ES     Individual     0.288     0.198     0.147     0.233     0.121     0.106     0.253     0.114     0.102

minMSE     0.919     0.909${}^{*}$     0.856${}^{*}$     0.955     0.951     0.927     0.957     0.919     0.889

AIC     0.958     0.961${}^{*}$     0.974${}^{*}$     0.860     0.940     0.940${}^{*}$     0.813${}^{*}$     0.934     0.933${}^{*}$

MMA     0.928     0.900     1.000     0.922     0.906     1.002     1.312     0.821     1.000

Mean group     1.237${}^{*}$     1.427${}^{*}$     1.225     1.011     1.561${}^{*}$     1.352     0.946     1.886${}^{*}$     1.248

UK     Individual      0.281     0.116     0.044     0.254     0.142     0.047     0.244     0.142     0.047

minMSE     0.928     0.988     0.953     0.840     0.984     0.868     0.743     1.034     0.913

AIC     0.871${}^{*}$     0.933     0.953     0.876     0.958     0.941     0.726${}^{*}$     0.917${}^{*}$     0.898${}^{*}$

MMA     1.457     1.375${}^{*}$     1.350     1.240     1.318${}^{*}$     1.580     1.523     1.069     0.854

Mean group     1.714     2.444${}^{*}$     3.688${}^{*}$     2.530${}^{*}$     2.445${}^{*}$     3.121${}^{*}$     4.330${}^{*}$     2.214${}^{*}$     2.650${}^{*}$

Nowcasting  MSE. For individual estimator: absolute value. For averaging estimators: MSE relative to individual estimator. For different estimation window sizes (44, 60, 76); selected weekly horizons relative to quarter end ($-6, 0, +4$ weeks).  * – forecasting performance difference significant at 10% in Diebold-Mariano test <cit.>

Distribution of relative MSEs across countries. Split by different averaging strategies and estimation window size. Same weekly positions as reported in table <ref>

Our key finding is that averaging with  smooth data-dependent criteria –  minimum MSE or AIC weights – generally leads to  improved forecasting performance, though the degree of improvement varies with the country in question.
This is clear from   table <ref>, as the vast majority of entries corresponding to those weights display relative MSE smaller than one, with improvements reaching up to 20%. The average gain in performance is on the scale of about 9% for minimum MSE weights and 5% for AIC weights. We also observe that minimum MSE weights and AIC weights do not uniformly dominate each other.

Figure <ref>   provides a box plot for relative MSEs for nowcasting GDP for all the countries in the panel for the vintages considered in table <ref>.
The figure illustrates that the favorable  performance is robust across countries and not limited to the five biggest economies reported in table <ref>.
Both minimum MSE and  AIC weights generally lead to an improvement in performance, as both rarely have relative MSE above one.  There is some evidence that the  minimum MSE weights have a   greater upside, at the price of potentially some more variability in the results, while AIC leads to smaller, but more tightly concentrated improvements.

Further, we find that   averaging is more attractive for the smallest sample size of $T=44$, with relative MSE generally approaching one as $T$ increases.
This can be clearly seen in figure <ref>, as the improvement range for AIC and minimum MSE weights becomes more concentrated and to closer to one. As previously remarked, as $T$ increases, the minimum MSE  estimator converges to the individual estimator; a similar point applies to AIC weights if the log likelihood is not divided by samples size and allowed to diverge as sample size grows.

Averaging methods that are data-independent (mean group) or  not smooth (MMA)  lead to generally poor results. While capable of offering improvements, the two approaches often lead to significantly worse performance than the individual estimator, regardless of $T$. This is clear from figure  <ref>, where the distribution of the  MSE of both averaging strategies has substantial mass above one.

We remark that the online appendix contains a number of additional robustness checks.
First, we consider nowcasting using bridge equations instead of U-MIDAS.
Second, we consider one- and two-quarter-ahead GDP forecasting.
The evidence emerging from these additional robustness checks on the performance of minimum MSE unit averaging estimator is consistent with the overall evidence reported here.

§ CONCLUSIONS

In this work we introduce a unit averaging estimator to recover unit-specific parameters in a general class of panel data models with heterogeneous parameters.
The procedure consists in estimating the parameter of a given unit using a weighted average of all the unit-specific parameter estimators in the panel. The weights of the average are determined by minimizing an MSE criterion.
The paper studies the properties of the procedures using a local heterogeneity framework that builds upon the literature on frequentist model averaging <cit.>.
An application to GDP nowcasting for a panel of European countries shows that the procedure performs favorably for prediction relative to a number of alternative procedures.