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.
305,703 characters · 81 sections · 246 citation commands
Optimal Estimation Methodologies for Panel Data Regression Models
Optimal estimation methodologies (e.g., see hilborn1969optimal and dreze1976bayesian) in the case of multivariate data series has been a research area of interest both in the econometrics as well as in the statistics literature the past 7 decades (see haavelmo1943statistical, marschak1944random, anderson1949estimation, koopmans1950identification, james1954normal, zellner1974time, sargan1976econometric, espasa1977spectral and forchini2003conditional). Moving into an era of ultra-high dimensional data structures where tools such as machine learning and deep learning techniques are employed for statistical learning purposes; it is of paramount importance to have a deep understanding of the classical optimal estimation methodologies in various econometric environments (see, white1996estimation, white2014asymptotic). In previous lecture series, katsouris2023high, katsouris2023limit, katsouris2023quantile discusses recent developments and several open problems in the time series and network econometrics literature are mentioned, with a special interest in nonstationary regression models and quantile regressions. The purpose of this lecture series is to present relevant issues on optimal estimation methodologies for panel data regression models (see also chamberlain1982multivariate, chamberlain1984chapter and baltagi2008econometric).
We shall begin by discussing the so-called identification problem in linear and nonlinear econometric models (e.g., see rothenberg1971identification and dreze1976bayesian) which are commonly presented in the literature based on suitable distributional conditions such as the Gaussianity assumption of structural disturbances (e.g., see phillips1976iterated). However, in order to narrow down the relevant literature we only focus on estimation and inference methodologies for panel data regression models which is considered as a statistical learning mechanism for various applications found in economics, finance, biostatistics and climate sciences among many other fields. Specifically, in the structural econometrics literature, weak identification is commonly discussed as a model specification issue. Some relevant questions of interest:
Overall, system identification plays an important role in revealing the unknown mechanisms of underlying complex phenomena. System identification includes detection of the model structure and estimation of the associated parameters. Moreover, a system identification problem can be thought of as an optimization problem where the optimal model is searched from a large predefined candidate model, given a criterion. The criterion is used to evaluate the performance of each model by measuring the discrepancy between the observed data and the model predictions (see, guo2016ultra). Good criteria result to not only better parameter estimation but also a good search path along which the search process converges quickly to the optimal solution. Different criteria have been used in system identification such as the $L^2$ norm in least squares regression and the $L^1$ norm in least absolute value regression (guo2016ultra).
Among these criteria, the least squares criterion is the most used because of its excellent properties, for example, least squares estimation can be configured to give estimates which are unbiased and efficient when the noise satisfies some basic assumptions. Thus, the least squares problem has analytic solutions and can be easily solved using the QR decomposition technique. However, the least squares technique does not allow to capture some special characteristics in system identification such as the interconnectedness in data points; so-called network dependence (see, katsouris2023limit).
A more accurate estimation methodology, which particularly overcomes the over-parametrization problem commonly found when using the least squares algorithm, is to employ an alternative criterion so-called ultra-least squares. The ULS criterion enables to characterise the model fitness more accurately. To put things into perspective, the ULS criterion considers model fitting in a smaller space, more specifically, the Sobolev space $H^m ( [0,T] )$. In other words, the ULS criterion takes into consideration not only the residuals but also the associate weak derivaties to measure the model fitness. Notice that the $L^2$ norm only emphasizes the similarity of two functions as a whole but disregards the closeness or shape (e.g., shape restrictions). Therefore, system identification can be interpreted as discovering unknown rules form a set of observations.
A useful norm to develop metric properties in the Sobolev space is defined as below:
where $D^j$ represents the $j-$th differentation operator. Based on the $\left\lVert \cdot \right\rVert_{ H^m }$ norm, a new criterion can then be defined as below:
Due to the fact differentiation is a linear operator, the above criterion can be written as below:
Thus, the $\mathcal{J}_H$ criterion consists of two parts: the first term corresponds to the standard least squares criterion which is based on evaluating the statistical distance (agreements) over the sample points; while the second term corresponds to the agreement of the weak derivatives which essentially emphases the agreement in shape (shape restrictions). Lastly, an interesting aspect worth investigating further is an ULS criterion obtained by adapting the $\mathcal{J}_H$ criterion to the nonlinear system identification problem\footnote{Nonlinear system identification involves both the estimation of the parameters and more importantly the problem of how to detect the structure of the unknown model. Model structure detection for linear systems is relatively easy and usually involves determining the order and time delay in a linear model. On the other hand, model detection can be complicated when the system is nonlinear due to the presence of many potential model terms and complex dynamics.}.
In order to evaluate the contribution of the unknown weak derivatives in the $\mathcal{J}_H$ criterion, we consider distributional assumptions. Consider the signal process $y(t)$, then the associated distribition $T_y$ is defined as a functional $T_y: C_0^{\infty} ( [0,T] ) \to \mathbb{R}$ such that
for all $\varphi \in C_0^{\infty} ( [0,T] )$. Therefore, the distribution $T_y$ has weak derivatives which are defined:
Similarly , the distributions that correspond to $x_i$ are defined as:
Thus, the regression problem is now solved in terms of conditional distribution projections. In other words, the system identification problem involves fitting the distribution $T_y$ by the combination of a set of distributions $T_{x_i}$. Therefore, the ultra-least squares problem becomes:
Notice that above we have used interchangeable the terms ultra-least squares and conditional density projections to distinguish to the case where conditional density projections refer explicitly to the least squares problem. For the remainder of this study we focus on $L^2$ estimation techniques and leave the extensions to the Sobolev space for future discussions (see, mahan2021nonclosedness, abdeljawad2022approximations and duzgun2023clustering).
When estimating the parameters of a correctly specified model, identification of the true parameters is a necessary condition for the consistent estimation and robust inference. However, identification is not a sufficient condition for consistency because the estimator may be constructed in such a way as not to be consistent for the true parameters, but for some other parameters which are nevertheless identifiable. Therefore, a model may fail to be identified, but estimation of a misspecified version of the model may yield identifiable parameters (see, bates1985unified).
Suppose that a correctly specified model has the following functional form
where $\mathbb{E} \left( X_t^{\prime} \epsilon_t \right) \neq 0$, and no instrumental variables are available for $X_t$. In general, we consider related conditions for existence and unique characterization of optimal estimation methodologies in panel data regression models (e.g., see andrews2001consistent).
\paragraph{Structure} The remainder of this study is organized as follows. In Section (ref) we discuss commonly used GMM estimation techniques and their asymptotic properties. In Section (ref) we present key aspects related to estimation and inference for panel data regression models under time series stationarity while in Section (ref) we consider the case of nonstationary panel data regressions. Section (ref) presents recent developments in relation to panel data model estimation from the network econometrics literature. Section (ref) and Section (ref) present some further applications. Network connectivity is a driving force for the risk transmission of various economic, financial and social phenomena. Related statistical problems include the modelling of financial contagion in stock markets, cointegration dynamics and market exuberance as well as the spread of epidemic diseases and the monitoring of climate change and biodiversity loss (see, gove2023coral). We present relevant econometric methodologies which can be employed when considering the empirical and theoretical implications of such topics.
We follow the framework proposed by Martinez2020asymptotic who consider the econometric estimation using the GMM methodology, as briefly described below. Let $\theta \in \Theta$ denote a $p-$dimensional vector of parameters partitioned into $\theta = \left( \vartheta^{\prime}, \psi^{\prime} \right)^{\prime}$ of dimensions of $p_{\vartheta}$ and $p_{\psi}$, respectively.
Denote with
to represent the sample moments, where $f_t( \theta )$ is a $k-$dimensional vector-valued function of data and parameters with $k \geq p$ and $\mathbb{E} \big[ f_t( \theta ) \big] = 0$ as the true value of $\theta$. Moreover, we let with $r( \theta )$ to be a known function of the parameters such that $r : \Theta \to \mathbb{R}^q, q \leq p_{\psi}$.
Suppose that $f_t ( \vartheta, . )$ and $r ( \vartheta, . )$ are continuously differentiable with respect to $\psi$, and let
Moreover, we denote with $\hat{V}_f( \theta )$ the $(k \times k)$ matrix that is positive definite almost surely, and define the GMM objective function as below
Furthermore, suppose that the constrained GMM estimator of $\psi$ given $\vartheta$ exists and is given by the following expression
We also simplify the notation as below
In addition we consider $\hat{C} ( \vartheta )$ be an almost surely full-rank $k \times ( k - p_{\psi} )$ matrix that spans the null-space of $\tilde{V}_f ( \vartheta )^{-1/2} \hat{J}_T ( \vartheta )$ such that
In particular, Martinez2020asymptotic develops an asymptotic theory framework based on fixed-smoothing asymptotics for the test statistics in order to account for the estimation uncertainty in the underlying LRV estimators. Consider the following long-run variance estimator
Therefore, a non-parametric estimator of the LRV takes the quadratic form below
such that $\omega(.,.)$ is a weighting function, and $h$ is the smoothing parameter indicating the amount of nonparametric smoothing. For example, we can estimate the kernel density the following way
for some kernel function $k(.)$, leading to the usual kernel LRV estimator. Thus, by substituting the smoothing estimator of the particular kernel function, we obtain the following test statistic
Then, the K statistic is based on the first-order derivative of $Q_T (\theta)$. Define as below the gradients
Taking the first-order and second-order derivatives of $\hat{V}_{ff}(\theta)$ with respect to $\theta_j$, we obtain
Then, it follows that
Denote with $D_T (\theta) = \big[ D_{T,1}(\theta),...., D_{T,d}(\theta) \big] \in \mathbb{R}^{ m \times d}$, such that
Then, the $K$ statistic for testing the null hypothesis $H_0: \theta = \theta_0$ against the alternative hypothesis given by $H_1 : \theta \neq \theta_0$ is given by
where for any concave function $\phi(\theta)$, $\partial \phi(\theta_0) / \partial \theta$ is defined to be
Thus, to consider fixed-smoothing asymptotics, we employ the orthonormal series LRV estimator
where $G$ is a smoothing parameter for this estimator and $\Phi_{\ell} (.)$ is a set of a basis functions on $L^2 [0,1]$. The weighting function is expressed with respect to a set of basis functions on the space of $L^2 [0,1]$. Therefore, the LRV estimator takes the following form
Then, an updated estimator for $\mathcal{J}_T (\theta_0)$ needs to be obtained from the sample such that
Therefore, the CLT to hold the following asymptotic distribution to hold
where $\psi_f \in \mathbb{R}^{ m \times 1 }$ and $\psi_{\mathsf{g}} \in \mathbb{R}^{ m d \times 1 }$. Therefore, it holds that we have a sequence of matrices
In other words, we consider the convergence rate of the $m-$system equations. By multiplying with $\textcolor{red}{T^ {\kappa} }$ we ensure that we take into account the different convergence rate depending on the type of functional form specification for the model under investigation. Specifically, when $\kappa = 0$, then the $m-$moment conditions doesn't include the correct rate of convergence which implies that the matrix $\Pi$ has a full column rank and therefore the parameter under the null hypothesis $\theta_0$, can be estimated at the usual parametric $\sqrt{T}$-rate. In other words, the case which corresponds to the weak identification of the model specification occurs when $\kappa = 1/2$, since it asymptotically converge into the null matrix, such that, $\Pi = 0$ and therefore, $\theta_0$ cannot be consistently estimated (see, Martinez2020asymptotic). Therefore, the estimation procedure for the case of fixed autocorrelation is given as below
Now, to estimate the above moment conditions the important component of the estimation procedure is to obtain unbiased estimators for the LRV covariance matrices which are computed based on a set of basis functions. Let $\ell \in \left\{ 1,..., G \right\}$, then the basis functions shall satisfy: (i). $\Phi_{\ell} (.)$ are piesewise monotonic, continuously differentiable, and (ii). $\Phi_{\ell} (.)$ are orthonormal in the space of $L^2 [ 0,1]$ functions and satisfy $\displaystyle \int_0^1 \Phi_{\ell} (x) dx = 0$. Therefore, the corresponding estimators are obtained as below
Following gospodinov2012local, consider that a univariate process is strictly stationary and geometically ergodic and denote the conditional moment restrictions imposed by economic theory
where $u : \mathbb{R}^{p+1} \times \Theta$ is a known function up to a vector of unknown parameters $\theta_0 \in \Theta$.
Consider the AR(1) model with martingale difference errors
In this case, the moment function is specified as $u( r_{t+1}, r_t, \theta_0 ) = r_{t+1} - \gamma_0 - \gamma_1 r_t$, with $\theta_0 = \left( \gamma_0, \gamma_1 \right)$. Let $x_t = ( r_t,..., r_{t-p+1} )^{\prime}$ and $y_{t+1} = ( r_{t+1}, x_t^{\prime} )^{\prime}$. Then, the conditional moment restriction model (ref) is estimated by the GMM estimator based on the unconditional moment restrictions
with a matrix of instruments $\mathcal{A} ( x_t, \theta_0 )$, which is implied by the original model. The continuously updated GMM estimator is defined as
is an optimal weight matrix to estimate the parameters from the unconditional moment restrictions $\mathbb{E} \big[ \mathsf{g} ( y_{t+1}, \theta_0 \big] = 0$. In particular, gospodinov2012local pursue an alternative approach and use a localized version of the GMM estimator that operates directly on the conditional moment restriction.
Therefore, the kernel estimator of the conditional moment condition $\mathbb{E} \big[ u ( y_{t+1}, \theta ) | x_t \big]$ is defined as
Then, the LGMM estimator minimizes its quadratic form as below
Moreover, the objective function of the LGMM estimator and its population counterpart are written respectively as below
where $V( x, \theta )^{-1}$ exists for each $x \in \mathbb{R}^p$ and $\theta \in \Theta$.
Define with
Then, it holds that
Then, by a change of variables $a = \frac{ x_j - x }{h}$ and an expansion around $a = 0$, whe get that
Combining these results we obtain that
\paragraph{Proof of Part (b).} By expanding the first-order condition $\partial \mathcal{Q}_n / \partial \theta = 0$ around $\theta_0$ we obtain
Moreover, consider the score function such that
for $\ell = 1,...,k$.
Another relevant methodology to GMM estimation especially for irregular data structures with nonlinear dynamics is the estimation approach of Efficient Method of Moments (EMM), proposed by newey1987hypothesis (see, also ortelli2005robust). Specifically, the simulation-based EMM technique provides a systematic way for generating moment conditions for simulated method of moments (SMM) estimation. This approach is also related to the indirect inference estimator. Both EMM and II use a first-stage auxiliary statistical model to generate moment conditions. However, the EMM mimics the first-order conditions for estimation of the auxiliary model, which is computationally much more tractable than mimicking the optimization problem itself, as is done in the II approach (see, li2009simulation, li2010indirect). Furthermore, EMM is useful in situations in which analytical characterization and evaluation of the likelihood function is infeasible. Thus, EMM selects moments based on the score function of an auxiliary model, called the score generator, to define a criterion function for SMM estimation (see, chung2001testing). In the case of panel data estimation, asymptotically EMM estimators for dynamic panel data regressions provide bias corrections when the number of time periods is fixed or tends to infinity with the number of panel units (see, breitung2022bias).
In particular, the weak instrument problem of the system GMM estimator in dynamic panel data models is studied by bhargava1991identification and bun2010weak. Moreover estimation and inference in panel data with cross-sectional dependence is discussed in bai2004estimating.
Consider again the following dynamic panel data regression model
where $\alpha_i$ represents the time-invariant unobserved heterogeneity (relevant references include among others the studies of huang2020identifying and bonhomme2015grouped.
Consider the model of interest as below:
where it is assumed that $\mathbb{E} \big( z_{it} \varepsilon_{it} \big) = 0$ for a given vector of instruments $z_{it}$ of dimension $R \geq K$, where $K$ is the number of elements of $\beta$ (that is, the number of regressors in the model). Then, the set of population moment conditions can be written as below:
In other words, these $R$ conditions can help to estimate the $K$ unknown parameters in $\beta$. Thus, the identification assumption is satisfied only for the true parameter values and is nonzero otherwise. On the other hand, since these expectations are unobservable in practice we rely on sample moments for statistical inference purposes which is given by
Moreover, the GMM estimator for $\beta$ is obtained by minimizing a quadratic form in the sample averages
where $W_{NT}$ is an $( R \times R )$ positive definite weighting matrix, which might depend upon the observed sample and thus needs to be estimated. When implementing the GMM estimation approach, we usually adjust the weighting matrix to obtain an asymptotically more efficient estimator. In particular, if the error term is heteroscedastic, but there is no correlation between different error terms, an empirical weighting matrix is given by the following expression
where $\hat{\varepsilon}$ is the residual given by $y_{it} - x_{it}^{\prime} \hat{\beta}_1$ such that $\hat{\beta}_1$ denotes an initial consistent estimator for $\beta$. Therefore, this makes the optimal GMM estimator a two-step estimator. During the first step, a consistent estimator for $\beta$ is obtained, which is used to calculate residuals and construct the estimated optimal weighting matrix. During the second step, an asymptotically efficient estimator is obtained. Thus, the optimal estimator can be obtained as below:
In this section, we discuss the framework proposed by gonccalves2015bootstrap that corresponds to modeling panel data with cross-sectional dependence (see, also phillips2003dynamic, bond2002projection and pesaran2021general). A relevant issue for identification and estimation is examined by a large stream of literature which develops econometric methodologies for capturing cross sectional dependence and heterogeneity via the use of dynamic panel models, based on the seminal contributions of pesaran2006estimation. Moreover, kapetanios2014nonlinear present a framework for nonlinear panel models with cross-sectional dependence. Recently, in the spatial econometrics literature various methodologies have been proposed to model both spatial dependence and cross-sectional effects such as li2020spatial. Moreover, Olmo2023 propose a network regression model with an estimated interaction matrix which incorporates both the cross-sectional as well as the network dependence in the form of a metric distance between the set of regressors (see, also kapar2022dynamic).
We shall denote with $\mathsf{cum} ( w_0 ) = \mathbb{E} ( w_0 )$ and $\mathsf{cum} ( w_0, w_{t_1} ) = \mathsf{Cov} ( w_0, w_{t_1} )$. Notice that for a given time series $\left\{ w_t \right\}$ and for $j \in \mathbb{N}$, we let $\mathsf{cum} \left( w_0, w_{t_1},..., w_{t_j - 1} \right)$ to denote the $j-$th order joint cumulant of $\left( w_0, w_{t_1},..., w_{t_j - 1} \right)$, where $t_1,...,t_{j-1}$ are integers. In particular, gonccalves2015bootstrap impose the assumption of a martingale difference sequence restriction on $\left\{ \epsilon_{it}, t = 1,2,... \right\}$ for each $i \in \left\{ 1,..., n \right\}$. Therefore, the m.d.s assumption implies that the model for the conditional mean of $y_{it}$ given $\mathcal{F}_{i,t-1}$ is correctly specified.
Denote with $Z_{nt}^{*}$ be a sequence of bootstrap statistics. These convergence modes hold
More precisely, the particular assumption of a correctly specified model for the conditional mean, allow us to obtain results for the recursive-design bootstrap based on the wild bootstrap. We use the following notation for the bootstrap asymptotics.
Consider the stationary linear dynamic panel model with fixed effects
where $| \theta_0 | < 1$ and $\alpha_i$ are individual specific fixed effects that capture the unobserved individual heterogeneity. The standard fixed effects OLS estimator of $\theta_0$ is given by (see, gonccalves2015bootstrap)
where
Therefore, the main goal of this section is to provide a set of assumptions under which we can prove the bootstrap results that will follow and at the same present the asymptotic theory of the fixed effects estimator under these assumptions.
Thus, we consider the joint asymptotic theory of $\hat{\theta}$ as $N, T \to \infty$. Then, the fixed effects OLS estimator can be represented as below
and $\mu_i = \mathbb{E} \left( y_{it-1} \right) = \alpha_i / (1 - \theta_0)$. Therefore, we obtain that
since it can be shown that $A_{NT} \overset{ d }{ \to } A$.
Furthermore, the following decomposition holds for the normalized score,
Notice that we can investigate the stochastic behaviour of the two terms above separately. The above result has two implications for the validity of the proposed bootstrap procedure. First, the bootstrap needs to mimic the asymptotic variance of $\hat{\theta}$ by $C = A^{-1} B A^{-1}$. More precisely, the variance has the usual sandwich form under conditional heteroscedasticity. In particular, it depends on the long run variance of the score process which is defined as below
In other words, the bootstrap validity depends on replicating the properties of the cross sectional average of the fourth order cumulants of $\varepsilon_{it}$. Second, the bootstrap needs to capture the asymptotic bias term $D$ created by the estimation of the fixed effects. More specifically, as the decomposition above shows, this noncentrality parameter results from the correlation between the averaged error terms $\bar{\epsilon}_i$ and the demeaned regressors $\left( y_{it-1} - \mu_i \right)$ and is non zero when $\rho = \mathsf{lim} \frac{N}{T} \neq 0$.
\paragraph{Recursive-design wild bootstrap}
[gonccalves2015bootstrap] The recursive-design bootstrap can generate a panel of pseudo observations $\left\{ y_{it}^{*}, i = 1,..., n ; t = 1,..., T \right\}$ recursively from the panel AR(1) model with estimated parameters,
where $\hat{\alpha}_i = \frac{1}{T} \sum_{t=1}^T \left( y_{it} - \hat{\theta} y_{it-1} \right)$ and $\hat{\theta}$ is a fixed effects OLS consistent estimator. Moreover, the initial condition is given by $y_{i0}^{*} = \frac{ \hat{\alpha_i } }{ 1 - \hat{\theta} }$, which is equivalent to setting $y_{it-1}^{*}$ to the stationary mean in the bootstrap world. In particular, the bootstrap residuals are obtained with the wild bootstrap $\epsilon_{it}^{*} = \hat{\epsilon}_{it} \eta_{it}$, where $\eta_{it} \overset{ \textit{i.i.d} }{ \sim } (0,1)$ over $(i,t)$ such that $\hat{\epsilon}_{it} = y_{it} - \hat{\alpha}_i - \hat{\theta} y_{it-1}$ are the estimated residuals.
Based on the above definitions, gonccalves2015bootstrap consider the bootstrap analogue of $\hat{\theta}$ for the recursive-design wild bootstrap OLS estimator, denoted by $\hat{\theta}^{*}_{rd}$ as below
where $\bar{y}^{*}_{i(t-1)}$ and $\bar{y}^{*}_{i (t)}$ are defined analogously to $\bar{y}_{i(t-1)}$ and $\bar{y}_{i (t)}$. Then, the following theorem, provides a result related to the asymptotic bootstrap validity of the recursive design estimator.
Notice that the proof for the above theorem needs to account for the incidental parameter bias generated by the estimation of the fixed effects due to the particular panel data structure.
\paragraph{Pairs Bootstrap}[gonccalves2015bootstrap]
An alternative bootstrap approach which is found to be robust to conditional heteroscedasticity of unknown form in the error term of a pure time series autoregressive model is the pairs bootstrap, where one resamples with replacement the vector that collects the dependent variable and its lagged values.
\paragraph{Proof of Lemma B1}[gonccalves2015bootstrap]
We can write the following
since $\varepsilon_{it}^{*2} = \hat{\varepsilon}^2_{it} \otimes \eta_{it}^2$. Moreover, the residual term can be expressed as below
We can also write $\big( \alpha_i - \hat{\alpha}_i \big) = \big( \hat{\varepsilon}_{it} - \varepsilon_{it} \big) + \big( \theta_0 - \hat{\theta} \big) y_{it-1}$. Therefore to show that $F_2 = o_p(1)$ which implies that $\frac{1}{NT} \sum_{i=1}^N \sum_{t=1}^T \hat{\varepsilon}_{it}^2 \overset{ p }{ \to } \sigma^2$, and thus we need to show that $\underset{ 1 \leq i \leq n }{ \mathsf{sup} } \left| \hat{\alpha}_i - \alpha_i \right| = o_p(1)$ under the assumptions above (see, gonccalves2015bootstrap for further details). Moreover, it holds that
which implies that $\sum_{t=1}^T \varepsilon_{it} = \mathcal{O}_p \left( \sqrt{T} \right)$, uniformly in $i$, and thus $\frac{1}{ \sqrt{T} } \sum_{t=1}^T \varepsilon_{it} = \mathcal{O}_p \left( 1 \right)$ uniformly in $i$. Furthermore, given that $\frac{1}{T} \sum_{t=1}^T y_{it-1} = \mathcal{O}_p(1)$, uniformly in $i$ and $\left( \hat{\theta} - \theta_0 \right) = o_p(1)$, we have that $\underset{ 1 \leq i \leq n }{ \mathsf{sup} } \left| \hat{\alpha}_i - \alpha_i \right| = o_p(1)$, that is, the fixed effect estimator is bounded in probability almost surely.
The common correlated effects estimation approach proposed by pesaran2006estimation, provides a sufficiently general setting for panel data models with cross-sectional dependence and thus renders a variety of panel model specifications as special cases. In the panel data literature with $T$ small and $n$ being large, the primary parameters of interest are the means of the individual specific slope coefficients, $\boldsymbol{\beta}_i$. Let $\boldsymbol{M}_{ \widehat{\boldsymbol{F}}_x } = \boldsymbol{I}_T - \bar{\boldsymbol{X}} \left( \bar{\boldsymbol{X}}^{\prime} \bar{\boldsymbol{X}} \right) \bar{\boldsymbol{X}}^{\prime}$ denote a projection matrix. Then, the modified estimator is given by
where $N$ are the number of cross-sectional units and $T$ are the number of time-series observations. The quantity of interest here is the asymptotic distribution of the above estimator. It can be proved that the asymptotic variance of $\widehat{\boldsymbol{\beta}}_{x}$ is identical to that of $\widehat{\boldsymbol{\beta}}$, so no asymptotic efficiency is lost by omitting $\bar{\boldsymbol{y}}$, although the bias term of the particular expression still remains computational intractable.
Due to the reasons explained above in order to correct the bias term that appears in the asymptotic expression, we need to employ a bootstrap approximation which can ensure a coordinate-wise convergence in probability to the true value of the population parameter. The asymptotic validity of this expression is obtained with the use of uniform coordinatewise convergence as shown below
Therefore, the above expression establishes the consistency of the bootstrap for the distribution of the estimator $\widehat{\boldsymbol{\beta}}_x$ for general $m \leq k$, and hence validates the construction of bootstrap confidence intervals. Therefore, to establish the asymptotic validity of the bootstrap $t-$intervals, define with $\boldsymbol{\Theta} = \boldsymbol{\Sigma}^{-1} \boldsymbol{\Psi} \boldsymbol{\Sigma}^{-1} $, and let $\widehat{\boldsymbol{\Theta}}^{*}$ be the bootstrap world equivalent of the corresponding variance estimator as below
Next we concentrate on the following sample variance estimator
Therefore, under the assumption of homogeneous slopes $\boldsymbol{\beta}_i = \boldsymbol{\beta}$, we establish its asymptotic distribution as $(N, T) \to \infty$ such that $T / N \to \tau < \infty$ in the case of general $m \leq (k+1)$ (see, harding2020common).
Statistical inference techniques to panel data with or without cross-sectional dependence include slope homogeneity testing. Several studies have extended these methods to nonstationary panel data models as in kapetanios2011panels and huang2021nonstationary. However, no statistical testing methodology exists for slope homogeneity that covers these cases, which is currently a topic worth investigating further.
Consider the formulation of the CCE estimator using vector notation such that
where we have that $\boldsymbol{Y}_i = ( y_{i1}, ..., y_{iT} )^{\prime}$, $\boldsymbol{X}_i = ( x_{i1},..., x_{iT} )$, $\boldsymbol{U}_i = ( u_{i1},..., u_{in} )$ and $\boldsymbol{F} = ( \boldsymbol{f}_1,..., \boldsymbol{f}_T )^{\prime}$. Define the projection matrices $\boldsymbol{P}_{A} = \boldsymbol{A} ( \boldsymbol{A}^{\prime} \boldsymbol{A} )^{-1} \boldsymbol{A}^{\prime}$ and $\boldsymbol{M}_A = ( \boldsymbol{I} - \boldsymbol{P}_A )$. Then, the transformed equation can be written as $\boldsymbol{M} \boldsymbol{Y}_i = \boldsymbol{M} \boldsymbol{X}_i \boldsymbol{\beta}_i + \boldsymbol{M} \boldsymbol{U}_i$. Therefore, the CCE pool estimator is defined as below
Under the alternative hypothesis, we have that the CCE estimator deviates from the true parameters at least for a non-zero fraction of individual units. Therefore, the particular model parametrization can be employed to construct tests statistics for slope homogeneity in panel data models with multifactor error structure. Define the weighted average CCE estimator as below
Then, the proposed test statistic is constructed as below
The modified PY test is interpreted as the weighted average distance between $\hat{\boldsymbol{\beta}}_{i, cce}$ and $\tilde{\boldsymbol{\beta}}_{i, wcce}$.
Combining the two equations we obtain
where
Therefore, based on the above reparametrizations we have that
where $\bar{\boldsymbol{P} } = \frac{1}{n} \sum_{ i = 1 }^n \boldsymbol{P}_i \bar{\boldsymbol{M}}$ is the residual maker of $\boldsymbol{G} \bar{\boldsymbol{P}}$.
Consider the following specification
where $\boldsymbol{G} = ( \boldsymbol{D}, \boldsymbol{F} )$ is the $T \times m + n$ matrix of integrated factors and $\boldsymbol{V}_i$ is a stationary error matrix. Moreover, we denote the OLS residuals of the multiple regression as $\hat{\boldsymbol{V}}_i = \boldsymbol{X}_i - \boldsymbol{G} \hat{\boldsymbol{\Pi}}_i$, where $\hat{\boldsymbol{\Pi}}_i = \left( \boldsymbol{G}^{\prime} \boldsymbol{G} \right)^{-1}\boldsymbol{G}^{\prime} \boldsymbol{X}_i$. Moreover, observe that $\hat{\boldsymbol{V} }_i = \boldsymbol{M}_{ \mathsf{g} } \boldsymbol{X}_i$. Then, we can write
since it holds that $\boldsymbol{M}_{ \mathsf{g} } \boldsymbol{G} = \boldsymbol{0}$.
Consider the following autoregressive distributed lag, ARDL(1,0), panel data model with homogeneous slopes and a multifactor error structure (see, norkute2021instrumental) such that
where the multifactor error structure is captured with the following equations
where $| \rho | < 1$ and $\boldsymbol{\beta} = ( \beta_1, \beta_2,..., \beta_k )^{ \prime }$ such that at least one of $\left\{ \beta_{\ell} \right\}_{ \ell = 1 }^k$ is non-zero and $\boldsymbol{x}_{it} = ( x_{1it},..., x_{kit} )^{\prime}$ is a $( k \times 1 )$ vector of regressors and $\boldsymbol{f}_{y,t}^0 = ( f_{x,1t}^0, f_{x,2t}^0,..., f_{x, m_x t}^0 )$ denoters an $( m_x \times 1)$ vector of true factors, and $\boldsymbol{v}_{it} = ( v_{1it}, v_{2it},..., v_{kit} )^{\prime}$
Although, under the presence of fixed effects, the number of fixed-effects $(\eta_i)$ approaches infinity at the same rate as $N \to \infty$. In other words, for each case we add the MLE estimation adds a parameter to be estimated. Therefore, we cannot rely on asymptotics as $N \to \infty$ since the application of maximum likelihood leads to inconsistent estimates and thus an alternative estimation or transformation approach is required for robust statistical estimation and inference purposes. The key is to consider the orthogonal reparametrization\footnote{The orthogonal reparametrization proposed by lancaster2002orthogonal it changes the meaning of the parameters representing the individual effects but not the meaning of the other parameters. This approach is particularly useful when due to the functional form of the model orthogonality cannot be achieved but information orthogonality can.} (OPM) approach such that we are not actually interested in estimates of the $( \eta_i )$ (as these are incidental parameters). In particular, we are interested in estimates of the common parameters such as $\left\{ \beta_1, \beta_2, \alpha_1, \sigma^2 \right\}$. Then the OPM approach implies a reparametrization of the incidental parameters so that the incidental and common parameters are information orthogonal. The lancaster2002orthogonal reparametrization approach allows us to write the likelihood in which the incidental parameters are informationally orthogonal from the other parameters.
Some important terms which we will need to obtain relevant results for their asymptotic behaviour can be obtained as below
Next, we consider expanding the following sample moments
Next, we can consider the kernel density estimates of the distribution of the Mahalanobis distance given by the following expression
compared to the theoretical $\chi^2$ densities.
Moreover notice that the above process is asymptotically stationary, and the specification $\boldsymbol{A} = \boldsymbol{M}$ means that it will also be asymptotically unidentified. Therefore, such a lack of identification is well known to manifest itself in $\sum_{t=1}^T \boldsymbol{W}_t^{\prime} \boldsymbol{\Sigma}_{\varepsilon}^{-1} \boldsymbol{W}_t$ having less than full rank and $Q_T$ having fewer degree of freedom that might be anticipated on the basis of conventional asymptotic theory.
Given a quantile $\tau \in (0,1)$, consider the following QR model defined by galvao2013estimation such that
where $\boldsymbol{x}_{it}$ is a $( p \times 1 )$ vector of regressors, $\boldsymbol{\beta}_0 (\tau)$ is a $(p \times 1)$ vector of parameters and $\alpha_{i0}(\tau)$ is a scalar individual effect for each $i$, and $u_{it}$ is the innovation term whose $\tau-$th conditional quantile is zero. Notice that the quantile-specific individual effect, $\alpha_{i0}(\tau)$, is intended to capture individual specific sources of variability, or unobserved heterogeneity that was not adequately controlled by other covariates. In general, each $\alpha_{i0}(\tau)$ and $\boldsymbol{\beta}_0(\tau)$ can depend on $\tau$, but we assume $\tau$ to be fixed throughout the framework here. Moreover, the model is semiparametric in the sence that the functional form of the conditional distribution of $y_{it}^{*}$ given $\left( \boldsymbol{x}_{it}, \alpha_{i0} \right)$ is left unspecified and no parametric assumption is made on the relation between $\boldsymbol{x}_{it}$ and $\alpha_{i0}$. Thus, the QR model can be written as below
The main problem of the above estimator is caused by its low frequency of convergence. Furthermore, additional regressors, large proportions of censored observations, and large samples only worsen the problem. Due to censored effects we consider the equivalent minimizer (see, galvao2013estimation)
Therefore, we denote with $\delta_{it} = \boldsymbol{1} \left( y_{it}^{*} > C_{it} \right)$ to indicate uncensored observations.
We define with
whose $\tau-$th conditional quantile given $( \boldsymbol{x}_{it}, \alpha_i, C_{it} )$ equals zero. Furthermore, it holds that
In other words, the restriction set selects those observations $(i,t)$ where the conditional quantile line is above the censoring point $C_{it}$. Then, the objective function is equivalent to the following
We investigate the asymptotic properties of the proposed two-step estimator. A particular issue we impose is that the individual fixed effects parameter $\boldsymbol{\alpha}$ whose dimension tends to infinity. However, it has been noted in the literature that leaving the individual heterogeneity unrestricted in a nonlinear or dynamic panel model generally results in inconsistent estimators of the common parameters due to the incidental parameters problem. In other words, noise in the estimation of the individual specific effects leads to inconsistent estimates of the common parameters due to the nonlinearity of the problem. Therefore, to overcome this problem it has become standard in the panel QR literature to employ a large $N$ and $T$ asymptotics (as joint limits). Denote with $\left\lVert \pi - \pi_0 \right\rVert_{\infty} = \underset{ w }{ \mathsf{sup} } \left| \pi(\boldsymbol{w}) - \pi_0(\boldsymbol{w}) \right|$ for a given $\pi(.)$ and a generic vector $\boldsymbol{w}$ (see, galvao2013estimation).
In this section, we focus on the asymptotic validity of statistical procedures for cluster-robust bootstrap inference and cluster-robust confidence intervals in quantile regression models (see, galvao2011quantile, hagemann2017cluster, galvao2020unbiased, galvao2023bootstrap and galvao2023hac among others). We consider the recentered population objective function given by the following expression
Notice that the map $\beta \mapsto M_n ( \beta, \tau )$ is differentiable with derivative given by $M^{\prime}_n ( \beta, \tau ) := \partial M_n ( \beta, \tau ) \big/ \partial \beta^{\top}$. Specifically, the first-order condition of the QR objective function can be written as
where $\psi_{\tau} ( \mathsf{z} ) = \big( \tau - \boldsymbol{1} \left\{ \mathsf{z} < 0 \right\} \big)$. Then, the sample analogue of this condition is
can be thought of as nearly solved by the QR estimate $\beta = \hat{\beta}_n ( \tau )$. Notice that to ensure that the bootstrap counterparts of the above quantities, that correspond to the QR estimate, accurately reflect the within-cluster dependence, the resampling scheme perturbs the gradient condition at the cluster level. In particular, the bootstrap resampling is approximated using the bootstrap gradient process $\mathbb{W}_n ( \tau ) := \mathbb{W}_n \big( \hat{\beta}_n( \tau ) , \tau \big)$ evaluated at the original QR estimate to construct the new objective function
and define the process $\tau \hat{\beta}^{*}_n ( \tau )$ as any solution to $\mathsf{min}_{ \beta \in B } \mathbb{M}^{*}_n ( \beta, \tau )$. Then, $\hat{\beta}^{*}_n ( \tau )$ can be interpreted as the $\beta$ that nearly solves the corresponding first-order solution based on the proposed bootstrap resampling procedure. Then, the distributional convergence occurs both in the standard sense and with probability approaching one, conditional on the sample data $D_n := \left\{ ( Y_{ik}, X_{ik}^{\top} )^{\top}: 1 \leq k \leq c_i, 1 \leq i \leq n \right\}$.
A block bootstrap algorithm in a longitudinal model is proposed by ju2015moving. In particular, assume that the data generating process is based on a longitudinal data model specification.
Thus, we can add up to $n_0$ individuals and plug this into the model and the results is a pseudo sample series $y_{11}^{*},..., y_{n m}^{*}$. Then, from the model $y_{ij}^{*} = \hat{\beta}_0 + \hat{\beta}_1 x_{ij} + \hat{e}^{*}_{ij}$, we fit the regression model and produce the new parameters $\hat{\beta}_0^{*}$ and $\hat{\beta}_1^{*}$. As a result, the asymptotic validity and justification of Moving Block Bootstrap in Longitudinal Data can be established by carefully considering analytical expressions using the robust regression M-estimator which solves the following optimization problem
in relation to the mixing properties of innovation and the bootstrapping scheme. Non-asymptotic theory and related probability bound results such as the Hoeffding's inequality can be found to be useful for these derivations (see, praestgaard1993exchangeably and bentkus2004hoeffding).
In this section we consider relevant aspects to specification testing in panel data regression models. Relevant studies include metcalf1996specification, su2013nonparametric and su2015specification among others. Thus, in order to correctly define the estimation and inference procedure we first need related regularity conditions regarding the dependence structure across the panel data. Specifically, based on existing results in the literature we can assume that we may have independence across cross-sectional units and strong mixing over time. In other words, we may assume that the innovation sequences in the given setting have bounded higher order moments using results such that Berneisten's inequalities for strong mixing processes. This allow us to study the asymptotic properties of related test statistics and estimators without worrying about the existence of cross-sectional dependence since we decompose that effect into conditional independence within a small neighborhood of values, similar to the meaning of near-epoch dependence in related econometric models.
We consider the example below which represents a panel data model where cross-sectional dependence is captured by the presence of common factor loadings. In other words, using individual fixed effects in the panel, facilitates the presence of heterogeneity of shocks across the cross-sectional units. As a matter of fact, this is a more realistic assumption since shocks such as technology shocks, oil price shocks and financial crises are more likely to have unequal effect across the cross-section. A small economy for example, tends to be more vulnerable to such shocks than a large economy.
Therefore our main objective is to construct a test for linearity for the proposed specification form. In other words, we are interested in testing the null hypothesis
Under the alternative hypothesis we have that
Furthermore, to facilitate the local power analysis, we define a sequence of Pitman local alternatives
where the function $\Delta (.) \equiv \Delta_{NT} (.)$ is a measurable nonlinear function, $\gamma_{Nt} \to 0$, as $( N,T ) \to \infty$. To do this we use the following notation. Define with $e_{it} \equiv Y_{it} - X_{it}^{\prime} \beta^0 - F_t^{0 \prime} \lambda_i^0$. Define the probability density function of the covariates with $f_{it}(.)$ which satisfies related regularity conditions that ensure its validity. Moreover, since we have that $e_{it} = \varepsilon_{it}$ and it holds that $\mathbb{E} \left( e_{it} | X_{it} \right) = 0$ under $H_0$, such that
Under the null hypothesis we have that the following relation holds: $e_{it} = \varepsilon_{it} + m(X_{it}) - X_{it}^{\prime} \beta^0$.
implying that $\mathbb{E} \big[ e_{it} \mathbb{E} \left( e_{it} | X_{it} \right) f_i (X_{it}) \big] > 0$, under $H_1$ (see, su2015specification).
Based on the above notation we can proceed with the introduction of the consistent test for the correct specification of the linear panel data model based on this observation. Specifically, in order to construct the test statistic, we need to estimate the model under the null hypothesis and obtain the restricted residuals $\hat{\epsilon}_i = \left( \hat{\epsilon}_{i1},..., \hat{\epsilon}_{iT} \right)^{\prime}$ for $i \in \left\{ 1,..., N \right\}$. Then, we can obtain the sample analog of $J$ such that
where $K_h (x) = \prod_{\ell=1}^p h_{\ell}^{-1} k \left( \frac{x_{\ell}}{ h_{\ell} } \right)$ is a univariate kernel function with the vector $h = \left( h_1,..., h_p \right)$ is a bandwidth parameter, and $\mathcal{K}_{ij}$ is an $\left( T \times T \right)$ matrix whose $(t,s)-$th element is given by definition $\mathcal{K}_{ij} = K_h \left( X_{it} - X_{js} \right)$.
Notice that Assumption 1 rules out conditional heteroscedasticity that depends on the past information at time $t-1$. However, it does allow for unconditional heteroscedasticity that depends on cross-sectional units and the scaled time index $\tau_t$ and is therefore, less restrictive than conditional homoscedasticity.
Based on the above regularity conditions, su2015specification presents the exact estimation procedure to construct a consistent specification testing procedure.
The framework proposed by wu2023testing considers testing for trend specification in panel data regression models. In particular, the asymptotic distributions of the proposed test statistic are established under the assumption of cross-sectional dependence, although by restricting to the case that the error components to follow a martingale difference sequence (MDS) and thus, rule out serial dependence. Therefore, the panel data trend model and the hypotheses of interest are presented below. Suppose that we observe the panel data of $\left\{ y_{it}, i = 1,...,N, t = 1,..., T \right\}$, where $y_{it}$ is a scalar dependent variable of interest, $N$ the number of panel individuals and $T$ the number of periods.
Thus, the model becomes as below:
where $\alpha_i$ represents the unobserved individual-specific effect that satisfies $\sum_{i=1}^N \alpha_i = 0$ and $u_{it}$ is the error component. More flexible error structure can be also allowed using a suitable specification for heteroscedasticity, cross-sectional and serial dependence in $u_{it}$ (see, wu2023testing).
An example of an application, is when $y_{it}$ represents the total rainfall or temperature across the United Kingdom, $\alpha_i$ is the unobserved region-specific effect and $\mathsf{g}_t$ represents the common climate change trend, and $u_{it}$ is the region specific error. Notice that the classical panel models often assume i.i.d disturbances. This assumption is likely to be violated as the dynamic effect of exogenous shocks to the dependent variable is often distributed over several time periods. Additionally, we assume that spillover effects, competition and global shocks can all induce disturbances that display cross-sectional dependence. Therefore, in order to allow for $y_{it}$ to be general enough to accommodate both cross-sectional and serial dependence, we assume $u_{it}$ to follow an AR$(p)$ process such that
where $A(L) = \left( 1 - \sum_{j=1}^p \rho_j L^j \right)$ with $p \geq 1$ a fixed integer, such that the polynomial operator has all roots strictly outside the unit circle. Furthermore, we assume that the dynamic structure of $u_{it}$ is homogeneous across units. In particular, the homogenous panel autoregressive models are widely used to capture the dynamics of macroeconomic and financial variables. However, most of the studies in the literature consider the case in which innovations are i.i.d over time and and across individuals. Here, we assume that that the innovation $\varepsilon_{it}$ is assumed to follow an MDS such that $\mathbb{E} \left( \varepsilon_{it} | \mathcal{F}_{t-1} \right) = 0$ almost surely for each $i$ and allow for cross-sectional dependence and heteroscedasticity, where $\mathcal{F}_{t-1}$ is the information set available at time $t-1$.
Several studies consider specification testing in panel data regressions (e.g., see lee2012hahn). Let $\chi_t = ( y_t, x_t ) \in \mathbb{R}^2$ be a strictly stationary $\beta-$mixing process and define with $\mathsf{g}(x) = \mathbb{E} \big[ y_t | x_t \big]$.
Consider testing the hypothesis that $\mathsf{g} (x) = \beta_0 + \beta_1 x$ against the alternative that $\mathsf{g}(x)$ is non-linear function of $x$. Let $\theta = \left( \theta_{l}, \theta_{nl} \right)$ where $\theta_l$ is the average partial effect under the linear specification and $\theta_{nl} = \mathbb{E} \left[ \frac{ \partial \mathsf{g} (x_t) }{ \partial x } \right]$ is the average partial effect under the non-linear specification. An estimator for $\theta$ based on a $Z-$estimator using a plug in non-parametric estimate $\hat{\mathsf{g}}_k = \hat{\mathsf{g}}_k(x)$. For this purpose we define the moment function below
and let $m_n ( \theta ) = \frac{1}{n} \sum_{t=1}^n \hat{m} \left( \chi_t, \theta, \hat{\mathsf{g}}_k \right)$.
The limiting distribution of the test statistic is analyzed for the following data-generating mechanism under local alternatives $\mathsf{g}_h (x)$,
where $u_t = y_t - \mathbb{E} \left[ y_t | x_t \right]$ is such that $\mathbb{E} [ u_t | x_t ] = 0$. Let $\theta_0 = ( \psi_1, \theta_{nl} )^{\prime}$ be the value of $\theta$ for the true data generating process under local alternatives.
Under regularity conditions it follows (from newey1994asymptotic), that for $h$ fixed
The correction term $\gamma ( \chi_t )$ accounts for non-parametric estimation of the nuisance parameter $\mathsf{g}_h$ and can be derived using the methods developed in Newey (1994). It is given by
such that $\zeta_x(x)$ is the marginal density of $x_t$. Define the empirical process
Notice that the stochastic equicontinuity properties of the empirical process given above can be used to verify regularity conditions. Furthermore, the functional central limit theorem delivers a stochastic process representation of the limiting distribution of $\hat{\theta}_{\kappa}$ over the class of local alternatives. To obtain the limiting distribution of the above empirical process we consider the following auxiliary vector
and the corresponding long-run covariance matrix given by
In summary, by expanding the following components separately, the convergence in probability of the estimator $\hat{\theta}_{\kappa}$ from its true parameter value can be expressed as
The asymptotic orthogonality between the stationary and the integrated regressors that is, $x_{it}$ and $y_{i,t-1}$ can be examined in the context of nonstationary panel data model specification. In particular, related studies to unit roots and cointegration in panels can be found in breitung2005parametric, breitung2008unit and han2010gmm. Further applications include aspects of estimation and inference for panel VAR models (see, juodis2018first, hayakawa2016improved and camehl2023penalized).
Consider the panel AR(1) model as below (see, juodis2021backward)
where the data observed over $i = 1,...,N$ cross-sectional units in $t = 1,...,T$ time periods.
As is well-known, the conventional Fixed Effects (FE) estimator suffers from a sizeable finite sample bias for small values of $T$. However, the bias is general more noticeable in case of persistent data which is a common pattern for most applications involving macroeconomic panels. Furthermore, we assume that idiosyncratic errors are $\varepsilon_{i,t}$ are independent over $i$, while the initial conditions $y_{i,0}$ are assumed to be observed. Therefore, an alternative estimator to mitigate the finite sample bias we consider the LS estimator of $\rho$ from the following augmented regression
with the new composite error term given by
where $\bar{y}_{i \bullet } = \frac{1}{T} \sum_{t=1}^T y_{i,t-1}$. The inclusion of $\bar{y}_{i \bullet }$ in the regression model, while at the same time ignoring the presence of $\eta_i$, ensures that the LS estimator, which is numerical equivalent to the FE estimator, is consistent as $T \to \infty$. On the other hand, using the full sample mean such that $\bar{y}_{i \bullet }$ creates other problems as it is correlated with all $\left\{ \varepsilon_{i,t} \right\}_{t=1}^{T-1}$. In particular, the sequence of combined error terms $\left\{ \tilde{\varepsilon}_{i,t}, \tilde{\varepsilon}_{i,t-1},... \right\}$ is not a MDS even when $\eta_i = 0$. However, this issue can be fixed easily by considering in the estimation of $\bar{y}_{i \bullet }$ is replaced by the backward (recursive) mean of $y_{i,t-1}$ such that $\bar{y}_{i,t-1}^b = \frac{1}{t} \sum_{k=0}^{t - 1 } y_{i,k}$. Therefore, unlike the full sample mean, the backward mean by construction is not correlated with the current and future values of $\varepsilon_{i,t}$. However, similar to the standard FE estimator, the LS estimator of this type is not consistent for any fixed $T$ but is consistent for $T$ large if the data is stationary (see, juodis2018first and juodis2021backward). In particular, showed that in a model with the autoregressive parameter equal to unity both estimators have a substantially smaller asymptotic variance than the FE estimator. These results are complementary to those provided by some other authors who study asymptotic and finite sample results under stationarity.
In this section, we examine estimation and inference in panel data predictive regression systems with Cross-Sectional Dependence (CSD) and heterogeneous degree of persistence. Specifically, we focus on the IVX instrumentation proposed by kostakis2015Robust, but modifying the framework to account for panel data structure with Cross-Sectional Dependence. In terms of unobserved common factors in the panel structure, we assume that the proposed econometric specification identifies a set of common factors which impose the cross-sectional dependence structure in the panel data predictive regression system. We consider the identification and estimation of bias-corrected pooled IVX estimators, for which we investigate both the finite and asymptotic distributional properties in relation to the correct implementation of the instrumentation methodology using available information in the cross-sectional set of predictors which are observable for each cross-sectional observation.
Panel data models are employed to account for unobserved heterogeneity and dependence. The standard fixed effects model captures only time-invariant heterogeneity, however in many time series applications heterogeneity and latent dynamics are time-varying effects which affect the parameter stability and robust econometric inference. Moreover, predictive regression systems are employed by econometricians to examine the joint predictability of a set of time series while accounting for the existence of local to unity dynamics in predictors, which allows to model the degree of persistence appeared in their time series. In its basic form, the degree of persistence is an unobserved random variable that describes the stochastic behaviour of the time series captured by the common localizing parameter $(c)$.
When one is interested to model the dynamic persistence of a cross-section of observations the homogeneous persistence dynamics may not adequately reflect the corresponding stochastic processes. For example, the framework proposed by katsouris2023statistical can be extended to provide a unified framework for examining identification and estimation aspects related to panel predictive regression systems with network structure. The proposed modelling approach provides a robust representation of the predictive-generating mechanism which allows to examine aspects of financial connectedness via the use of predictive regression systems in panel data structures.
Our interest is the construction of a suitable IVX estimator (see, kostakis2015Robust) to accommodate the cross-sectional time series structure of the data (see, hjalmarsson2006predictive). We discuss the effectiveness of a pooled panel IVX estimator which can smooth the persistence levels across all cross-sectional observations $i$. Persistence homogeneity implies that the localizing coefficient of persistence remains fixed across panels although in a multivariate setting persistence levels is permitted to be different but of the same class. In the literature, various ways are presented regarding the analysis of cross section dependence. Specifically, we are interested for the case where both $N$ and $n$ are large and of the same order of magnitude, that is, $(N,n) \to \infty$ jointly which requires the notion of sequential asymptotic theory.
We follow the econometric specification proposed by hjalmarsson2006predictive, which consists of a cross-sectional representation of panel predictive regression systems with common factors. This provides a natural way of examining the implementation of the IVX instrumentation when the econometric identification induces cross-sectional dependence structure. The IVX instrumentation can provide a suitable methodology of smoothing out persistence effects across panel data model which accommodate cross-sectional dependence, a reasonable assumption especially when considering common factors (such as macroeconomic conditions, volatility spillovers, industrial factors etc.) affecting the degree of persistence of the cross-sectional observations. The proposed framework has various applications in the examination of long-run economic relations driven by stochastic processes in which effects of shocks are propagated across the cross-sectional observations.
Suppose the information set up to time $t$, that is, $\mathcal{F}_t$ includes a large number of predictors $x_{it}$ for $i=1,...,N$ and $t=1,...,n$. We are interested to examine the predictive ability of a panel data predictive regression system using predictors of heterogeneous persistence. Specifically, consider a panel data structure with the pair of random variables $( y_{i,t}, x_{i,t} )$ where $x_{i,t}$ is an $m \times 1$ dimensional vector. Thus, we are constructing a predictive regression system which accounts for the predictability of a set of variables $y_{i,t}$, accounting this way for the idiosyncratic persistence of the other firms in the network.
The particular specification is given by the following system of equations.
where $\boldsymbol{R}_i = \displaystyle \left( \boldsymbol{I} - \frac{ \boldsymbol{C}_i }{n^{\alpha_i}} \right)$ is the $m \times m$ matrix of degree of persistence, $f_t$ is a $k \times 1$ vector capturing common factors in the error terms of the predictive regression system equations (see, kostakis2015Robust).
Firstly, the specification builds on the current frameworks of multivariate predictive regression systems, by considering instead of modelling simultaneously the joint predictability of stock returns with a panel data structure. Secondly, allowing for cross-sectional dependence which in terms of the econometric specification can be represented via the use of common factors affecting the panel of firms, can capture spillovers and network effects. Moreover, cross-sectional dependence which has an economic interpretation as well in terms of macroeconomic conditions and shocks provides a suitable framework for examining heterogeneous persistence in panel structures accounting for such cross-sectional effects across firms. Such common factors allow the inclusion of a baseline persistence across the set of predictands. Thirdly, panel estimators are derived using sequential limits (hjalmarsson2006predictive), which usually implies first keeping the cross-sectional dimensions, $N$, fixed and letting the time-series dimension, $n$, go to infinity, and then letting $n$ go to infinity. We denote such sequential convergence denoted as $(N,n \to \infty)_{\text{seq}}$. We denote with $BM( \boldsymbol{\Omega} )$ the Brownian montion with covariance matrix $\boldsymbol{\Omega}$.
In particular, moon2000estimation, examine the case of inference in panel data autoregressive models with near to unity roots. In particular, the proposed framework allows for sequential asymptotic theory in order to examine the asymptotic properties of such panel data estimators which can accommodate for near to unity cases. We build on the theory of near to unity for panel data autoregressive estimators by focusing on the panel data predictive regression econometric specification. Assumption 1 below provides the necessary conditions for the identification of the proposed econometric specification. The innovation processes of the system indicates the stochastic behaviour of a system which describes a panel data structure with cross-sectional dependence.
Our aim is to investigate whether the IVX instrumentation across the panel smooths out the abstract degree of persistence across the cross-sectional observation $i$ and under the sequential asymptotic framework. Specifically, assuming the existence of common effects for the cross-sectional observations $i$, then using the IVX estimator instead of a pooled estimator can achieve the mixed normality assumption even under abstract degree of persistence. Furthermore, we impose the following assumption which provides a necessary condition for the degree of cross-sectional dependence among the observations of the cross-section. In other words, we impose an assumption which ensures a maximum bound for the cross-sectional dependence which also ensures that there is a limited amount of cross-sectional interactions and financial interconnectedness. Such assumption provides also conditions for the asymptotic efficiency of the IVX estimator in the case of panel data predictive regression specification. Take for example, the case of $m-$dependence in time series. Then, imposing such a condition for the panel data models, induces an asymptotic independence condition across the cross-sectional observations.
Even though this is a strong assumption, we consider that the cross-sectional constructed instruments are strictly exogenous which ensures consistency and asymptotic normality and mixed normality results. Moreover, the assumption of weak exogeneity, provides conditions under which we can perform efficient inference on the conditional model without imposing any regulatory assumptions on the functional form of the system hatanaka1996time. Therefore from the above econometric setting we can see that the main challenge with the robust estimation and identification is that usually the literature on large linear models focuses on the $\textit{i.i.d}$ case, while in our setting we have some type of cross-sectional dependence which could distort the asymptotic theory, especially with persistence and endogenous regressors.
\paragraph{Cross-Sectional Independence}
\paragraph{Cross-Sectional Dependence}
This is the main section of the paper. We aim to investigate whether the assumption of cross-sectional dependence as seen via the existence of common factors, along with the assumption of heterogeneous degree of persistence, induces a consistent IVX estimator. We argue that even though we include a set of common factors in the model, there is still the possibility of existence of heterogeneous degree of persistence across the predictors of the cross-sectional observations.
In this section, we examine the use of the IVX instrumentation methodology as in kostakis2015Robust, to the framework of the panel data predictive regressions with cross-sectional dependence. We consider the first order difference of the corresponding state equation, for the cross-sectional observation $i$ is
where $\alpha_i$ the cross-sectional rate of convergence of the first difference equation. Note that, the particular first difference is not an innovation process unless the regressor belongs to the persistence class of integrated processes. However, it behaves asymptotically as an innovation after linear filtering by a matrix consisting of near-stationary roots. The mildly integrated instrument for each cross-sectional observation $i$ is given by (see, kostakis2015Robust)
The artificial matrix $\widetilde{\mathbf{R}}_{i}$ has the following form
The particular artificial matrix facilitates a way of smoothing the degree of persistence in the cross-section by constructing instruments derived from the existing information in the regressor of the data and therefore the induced instrumentation method has milder degree of persistence (see, kostakis2015Robust). For the proposed panel data specification, the aim is to investigate whether using instruments corresponding to each of the cross-sectional observations $i$ induces a mixed Gaussian asymptotic distribution for the IVX estimator.
Consider the model
Moreover, consider the differencing estimators proposed by camponovo2015differencing for autoregressive models to the predictive regression models. Then, we consider the differenced observations
Moreover, we consider also differenced response variables such that
Thus to determine the stationary instruments, $w_t$, $t = \ell + 1,..., n$, which are strongly correlated with stationary differenced predictors $\Delta x_{t-\ell-1}$, those should satisfy the following moment conditions
Using $w_t := \Delta x_{t-\ell-1}$, we get that $\mathbb{E} \big[ \left( \Delta y_{t-\ell} - \beta \Delta x_{t-\ell-1} \right)\Delta x_{t-\ell-1} \big] = \mathbb{E} \big[ \left( u_t - u_{t-\ell} \right) \Delta x_{t-\ell-1} \big]$, so
Unless $\rho = 0$ or $\sigma_{uv} = 0$. Therefore, $\Delta x_{t-\ell-1}, t = \ell + 1,...,n$, are not valid instruments. However, with slight modifications and based on the following moment equalities
Then, we can prove that the instruments
satisfy the moment conditions below
Therefore, based on the above moment conditions, for a fixed value of $\ell \geq 2$, we define the new class of estimators $\beta_{n,\rho}^{(\rho)}$ of the parameter $\beta$ such that
The null hypothesis of exogeneity has also been reported in the cointegration literature before so this is an important aspect of consideration especially when considering cointegrating panel data regression models. Relevant studies on panel cointegration with respect to cross-sectional dependence and factor dynamics include quintos1998analysis, kao1999international, bai2004estimating, bai2009panel, kapetanios2011panels and westerlund2022factor among many others.
To study the distributional properties of such tests, we will describe the DGP in terms of the partitioned vector $z_{it}^{\prime} \equiv \big( y_{it}, X_{it}^{\prime} \big)$ such that the true process $z_{it}$ is generated as
Consider the following data generating process
Stacking the error process defines $\eta_t = \left[ u_t , \ v_t^{\prime} \right]^{\prime}$. Furthermore, it is assumed that $\eta_t$ is a vector of $I(0)$ processes in which case $x_t$ is a non-cointegrating vector of $I(1)$ processes and there exists a cointegrating relationship among $[ y_t, x_t^{\prime} ]^{\prime}$ with cointegrating vector $[ 1, - \beta^{\prime} ]^{\prime}$. To review existing theory and to obtain the key theoretical results in the paper, assumptions about $\eta_t$ are required. It is sufficient to assume that $\eta_t$ satisfies a functional central limit theorem (FCLT) of the form given below
Define the partial sum process such that $\widehat{S}_t = \sum_{j=1}^t \widehat{\eta}_t$. We start by establishing an invariance principle
Using the definition of $\widehat{\eta} = \left[ \widehat{u}_t, v_t^{\prime} \right]^{\prime}$ and stacking now leads to the following asymptotic theory result
Under the stated assumption it holds that
We sketch the result for the Bartlett kernel only. The proposition is established by showing
To begin with, showing that $\mathsf{plim}_{ b \to 0 } Q_b \left( B_v, B_{u.v} \right) = 0$ is trivial. Furthermore, it is well-known in the fixed$-b$ literature that as $b \to 0$, fixed$-b$ limiting random variables converge to the long run variance being estimated. Specifically, the long-run covariance between $B_v$ and $B_{u.v}$ are independent and so it follows that $\mathsf{plim}_{ b \to 0 } Q_b \left( B_v, B_{u.v} \right) = 0$.
To correctly define the limiting distributions we consider a sequence of rolling statistics. In particular, we first estimate the subsample statistics $t_{zx} ( \tau, \tau + \Delta \tau )$ for $t = \left\{ \floor{ \tau T} + 1,..., \floor{ \tau T} + \floor{ T \Delta \tau} \right\}$, where the window width is $\floor{ T \Delta \tau}$. Thus, the framework proposed by Magdalinos (2020) provides the first instance of standard Gaussian and chi-squared asymptotics applying respectively to the OLS estimator and the Wald statistic in a vector autoregression or predictive regression model with conditionally heteroscedastic innovations.Notice that because of the divergence of these partial sums we need to use an appropriate normalization factor which depends on the exponent rate of persistence in the LUR specification. Check also the paper: bruggemann2016inference.
where the model includes a $k-$dimensional vector of non-stationary regressors and the regression error $e_{it}$ is stationary and i.i.d across $i$. Then, we can show that the pooled OLS estimator of $\beta$ is defined as
The limiting distribution is shifted away from zero due to an asymptotic bias induced by the long-run correlation between $e_{it}$ and $\varepsilon_{it}$. The exception is when $x_{it}$ is strictly exogenous, in which case the estimator is $\sqrt{n} T$ consistent. The asymptotic bias can be estimated and a panel FM estimator can be developed along the lines of phillips1990statistical to achieve $\sqrt{n} T$ consistency and asymptotic normality. Furthermore, the cross-section independence assumption is restrictive and difficult to justify when the data under investigation are economic time series. In view of co-movements of economic variables and shocks, we model the cross-section dependence by imposing a factor structure on $e_{it}$,
where $F_t$ is an $r \times 1$ vector of latent common factors, $\lambda_i$ is an $r \times 1$ vector of factor loadings and $u_{it}$ is the idiosyncratic error. If $F_t$ and $u_{it}$ are both stationary, then $e_{it}$ is also stationary. In that case, a consistent estimator of the regression coefficients can still be obtained even when the cross-section dependence is ignored. In the first step, pooled OLS is used to obtain a consistent estimate of $\beta$. The residuals are then used to construct a FM estimator. In other words, nuisance parameters induced by cross-section correlation are dealt similar to the case of serial correlation by suitable estimation of the long-run covariance matrices. An alternative estimator can be developed by rewriting the equation as
Moving $F_t$ from the error term to the regression function (i.e., treated as parameters) is desirable for the following reason. If some components of $x_{it}$ are actually $I(0)$, treating $F_t$ as part of the error process will yield an inconsistent estimate for $\beta$ when $F_t$ and $x_{it}$ are correlated. Under the presence of global stochastic trends, that is, $F_t$, which are shared by each cross-sectional unit, a new methodology needs to be developed. Denote with $(n,T) \to \infty$ as the joint limit and with $( n, T )_{ \mathsf{sq} } \to \infty$ as the sequential limit which implies that $T \to \infty$ first and $n \to \infty$ later. Moreover, we denote with $\mathcal{MN} (0, V)$ the mixed normal distribution with variance $V$.
Relevant applications of panel cointegrating regressions is the framework proposed by wagner2020fully who develop a SUR system with cointegrating dynamics. Although, the SUR cointegration literature differs in several respects from the SUR literature, where the former case implies that regressors are assumed to be strictly exogenous and stationary and the stationary errors serially uncorrelated.
The econometric estimation of the SUR representation proposed by wagner2020fully is based on the assumption of serially correlated errors which involves the estimation of long-run variance matrices rather than estimates of contemporaneous variance matrices. Thus, the presence of regressor endogeneity in a cointegration setting necessitates the usage of modified least squares estimators to allow for asymptotically normal or chi-squared inference. The proposed framework is then employed for the analysis of the environmental Kuznets curve (EKC), specifically for carbon dioxide (CO$_2$) emissions\footnote{Specifically, the EKC hypothesis postulates an inverted U-shaped relationship between the level of economic development and pollution of emissions. Furthermore, there is significant use of unit root and cointegration techniques both in (single) time series and panel data settings.}. Specifically, wagner2020fully consider the case of panel data with small cross-sectional dimension and develop estimation and inference techniques to combine cointegrating polynomial regressions to a system which allows to test general hypotheses concerning group-wise pooling\footnote{Notice that pooled estimation, when appropriate, leads to considerable efficiency gains, but can lead to misleading results when it is implemented incorrectly. In particular, wagner2020fully conduct a simulation study to assess the finite sample performance of our estimators and tests based upon them. The FM-SUR estimator outperforms the FM-SOLS estimator in terms of bias and root mean squared error (RMSE). However, test statistics based on FM-SOLS exhibit in many configurations lower size distortions than tests based on FM-SUR. The latter have, however, higher size-corrected power than the former. Therefore, the evidence concerning hypothesis testing is mixed. Pooling, illustrated in the simulations by pooling the coefficients for the integrated regressors and its square over all cross-section members, leads to major performance improvements, in particular for the performance of tests. Moreover, the simulations also indicate the limitations of an unrestricted SUR approach in case of large $N$ and small $T$. Thus, it turns out that the estimation of large unrestricted long-run and half long-run covariance matrices is the key reason for the relatively poor performance in these constellations.}.
Consider the system of equations with observations available for $i = 1,..., N$ and $t = 1,...,T$ such that
Then the $N$ equations can be written in matrix form as a system of seemingly unrelated cointegrating regressions defined as $y_t := Z_t^{\prime} \theta + u_t$.
To establish the limiting distributions of model estimators we consider that the following functional central limit theorem (FCLT) to hold for $\left\{ \xi_t \right\}_{ t \in \mathbb{Z} } := \left\{ \big[ u_t^{\prime}, v_t^{\prime} \big]^{\prime} \right\}_{ t \in \mathbb{Z} }$
with the long-run covariance matrix $\boldsymbol{\Omega} := \sum_{ j = - \infty}^{ \infty } \mathbb{E} \left( \xi_t \xi_{t-j}^{\prime} \right)$ and $W(r)$ is a $2N$ standard BM. Regarding the estimation methodology, wagner2020fully examine the statistical properties of the OLS estimator given by $\widehat{\boldsymbol{\theta}}_{OLS} = \left( Z^{\prime} Z \right)^{-1} Z^{\prime} y$ which is consistent but its limiting distribution is contaminated by second order bias terms. Moreover, wagner2020fully consider a feasible GLS-type SUR estimator in which rather than the estimated error covariance matrix $\widehat{\boldsymbol{\Sigma}}_{uu}$, the estimated long-run covariance matrix $\widehat{\boldsymbol{\Omega}}_{uu}$ is used as a weighting matrix. The modified SUR estimator is defined as below
Based on the above formulations, wagner2020fully show that the MSUR estimator is consistent with a nuisance parameter dependent limiting distribution. Consequently, to construct estimators with a zero mean Gaussian mixture limiting distribution that allow for asymptotic chi-square inference, we propose fully modified type corrections to these two estimators. In particular, we consider the FM-OLS estimator of phillips1990statistical which is based on a two-part transformation. The first transformation changes the dependent variable such that $y_t^{+} := y_t - \widehat{\boldsymbol{\Omega}}_{uv} \widehat{\boldsymbol{\Omega}}_{vv}^{-1} \Delta x_t$. In particular, the second transformation consists of subtracting an appropriately constructed correction term to remove bias terms otherwise present in the limiting distributions. This transformation depends upon the estimator considered as starting point OLS or MSUR and the specification of the equation system.
\paragraph{Proof of Proposition 4}
We follow the framework proposed by madsen2005estimating. Consider the variables $Y_{it}$, $X_{1it}$ and $X_{2it}$ where $i = 1,...,N, t = 1,...,$ such that these variables are of dimensions $k_0 \times 1$, $k_1 \times 1$ and $k_2 \times 1$ respectively. For every cross-section unit we assume that the variables are generated by the following model
where $\gamma_1$ and $\gamma_2$ are $k_1 \times k_0$ and $k_2 \times k_0$ matrices of parameters, respectively , and where the time-series processes $\eta_{0it}$, $\eta_{1it}$ and $\eta_{2it}$ are weakly stationary for every cross-section unit for $i = 1,...,N$.
In terms of the identification and model specification, we employ the dynamic model for the variables at the individual level. Suppose that a cross section obtained at some point in time is available. The cross-section is sampled consisting of observations of $N$ cross-section units at time $t \in \mathbb{N}$. In particular, we consider the cross-section demeaned quantities
For notational convenience the stacked $( k_1 + k_2 )-$dimensional stochastic variable $X_{it}^{*}$ is defined such that $X_{it}^{*} = \left( X_{1it}^{*\prime}, X_{2it}^{*\prime} \right)^{\prime}$ and the corresponding $( k_1 + k_2 ) \times k_0$ parameter matrix as $\gamma = \left( \gamma_1^{\prime}, \gamma_2^{\prime} \right)^{\prime}$. We can now specify the regression equation describing the cointegrating relations can be expressed in terms of the demeaned variables below
where $\eta_{0it}^{*} = \eta_{0it} - \frac{1}{N} \sum_{i=1}^N \eta_{0it}$. The corresponding cross-section OLS estimator is defined as below
Moreover, by relevant assumptions the regressor $X_{it}^{*}$ is independent of the regression error $\eta_{0it}^{*}$ since the aggregate shocks have been removed from the variables. This immediately implies that $\hat{ \gamma }_{N,t}$ is an unbiased estimator of $\gamma$, that is, $\mathbb{E} \left( \hat{ \gamma }_{N,t} \right) = \gamma$. Then, the asymptotic behaviour as $N \to \infty$ of the cross-section estimator $\hat{ \gamma }_{N,t}$ is given by the following proposition.
Notice that since the regressor $X_{2it}$ is nonstationary when viewed as a time series, the asymptotic variance of $\sqrt{N} \left( \hat{ \gamma }_{N,t} - \gamma \right)$ depends on the point in time where the cross-section is obtained. We consider the cross-section estimator of $\gamma_2$ defined as the submatrix of $\hat{ \gamma }_{N,t}$ corresponding to the regressor $X_{2it}$. To be more specific let this estimator denoted by $\hat{ \gamma }_{2N,t}$ be the last $k_2$ rows in $\hat{ \gamma }_{N,t}$ and let $\Sigma^{22t}$ be the lower $k_2 \times k_2$ diagonal block matrix of $ \Sigma_t^{-1}$, that is,
where $\Sigma_{t}$ is decomposed according to $X_{1it}$ and $X_{2it}$ as below
Then according to Proposition 1 the limiting distribution of $\hat{ \gamma }_{2,N,t}$ is given by
Consider the fixed effects functional coefficient panel data model (see, phillips2022functional)
where $x_{it}$ is a $p-$vector of regressors, $z_{it}$ is a $q-$vector of covariates that determine the (random) coefficients $\beta( z_{it} ) = \big( \beta_1 ( z_{it} ),..., \beta_p ( z_{it} ) \big)$, the $\alpha_i$ are individual fixed effects, and the error $u_{it}$ has zero mean and finite variance $\sigma_u^2$. Thus, we focus on the case where both $x_{it}$ and $z_{it}$ are exogenous.
Moreover, phillips2022functional developed asymptotic theory for the estimator $\hat{\beta}$ in various settings depending on whether $N$ is fixed or $N \to \infty$ and whether $\eta =0$ or $\eta \neq 0$. In all cases, it is presumed that $T \to \infty$. Results with $N \to \infty$ include both sequential limit $( N, T )_{ \mathsf{seq} } \to \infty$ theory, where $T \to \infty$ followed by $N \to \infty$, and joint limit $( N, T ) \to \infty$ theory, where $T, N$ pass to infinity together. Joint limit theory is obtained by following the double indexed limit theory but our results provide an important extension that covers cases of multiple convergence rates and possibly degenerate limit distributions. In consequence, the divergence rate of the cross-section sample size $N$ needs to be controlled in order to control the random bias contributed by estimation of $\beta_0 ( z_t )$ (see, phillips2022functional).
Statistical inference for the econometric specification of phillips2022functional corresponds to testing specific parametric forms of functional coefficients. Thus, the relevant hypothesis concerning the functional coefficient $\beta(z)$ is whether this vector of coefficient functions can be treated as a constant vector. In particular, tests of such hypotheses can be constructed by examining the discrepancy between the nonparametric estimate of $\beta_0$ and parametric estimate of $\beta_0$. To distinguish the alternative from the null we further require some conditions, so that the function $g(z)$ is not a constant function. This kind of local alternative is commonly used in the study of nonparametric and semiparametric inference involving stationary and nonstationary data. Moreover, to establish joint asymptotics for $\beta$ as $( N, T ) \to \infty$ when $\eta \neq 0$, phillips2022functional, take into account of the singularity that arises in the limiting signal matrix in the passage to joint asymptotics.
As expected, under general weak dependence assumptions on $u_t$, the simple reduced rank regression models are susceptible to the effects of potential misspecification in the transient dynamics. These effects bear on the stationary components in the system. In particular, due to the centering term, both the OLS estimator and the shrinkage estimator are asymptotically biased (see, phillips2022functional).
Using an empirical illustration we can demonstrate the need to implementing a different cluster-based structure when regressors are near-unit root; in the case when there are persistent data and there is no correction due to endogeneity and high persistence. Assuming that there is a known cluster structure a common approach in the literature is to use a clustering algorithm. Current methodologies presented in the literature can be compared with new panel tests that make use of clustering individual time series into common groups reveal the additional discriminatory power obtained by grouping. Specifically, the clustered panel $t-$test introduced by liu2022panel help to diagnose mildly explosive price behaviour in US city housing markets where individual time series tests reveal no evidence of such behaviour, confirming this way the discriminatory power gains that arise from cross section aggregation.
In order to capture explosive and mildly explosive behaviour in panels we used the following data generating process based on the time series model of PM (2007).
Notice that the exponent rate $\gamma \in (0,1)$ and the scale coefficients $c_{ \mathsf{g}_i }$, both influence the extend of departure of the autoregressive coefficients $\rho_{ \mathsf{g}_i }$ from unity, and $\mathsf{g}_i$ denotes the group membership of individual $i$, for which the group structure is defined late. For nonstationary data of each cluster, it holds that the innovations $u_{it}$, follow a stationary linear process for each $i$ and the various variance estimates
Latent group membership of the $\rho_{ \mathsf{g}_i }$ arises through the localizing scale parameters $\left\{ c_{ \mathsf{g}_i } \right\}_{ i = 1}^n$. The framework we adopt lies between a homogeneous panel; in which case $c_{ \mathsf{g}_i } \equiv c$ for all $i$, and a fully heterogeneous panel; in which case $c_{ \mathsf{g}_i } \neq c_{ \mathsf{g}_{\ell} }$ for all $i \neq \ell$. In the paper of PCB et al, the authors assume a group structure involving a fixed number $G < n$ of unknown separate groups that are classified according to the scale parameter $c_{ \mathsf{g}_i }$. The group membership variables are given by the $\left\{ \mathsf{g}_i \right\}_{i=1}^n$ which maps individual units such that $i \in \left\{ 1,..., n \right\}$ into specific groups for which $j \in \left\{ 1,..., G \right\}$ with $G < n$ and allows for several possible midly explosive and mildly integrated groups together with a unit root group.
Consider the following VAR-type model with explosive roots (see, chen2023seemingly)
where $\boldsymbol{X}_t$ is a $d-$dimensional vector with $\boldsymbol{X}_t = \left[ x_{1,t},..., x_{d,t} \right]^{\top}$, while the initial value is set to $x_{i,0} = 0$ for $i \in \left\{ 1,..., k \right\}$ for simplicity. The residual sequence $\boldsymbol{u}_t = \left[ u_{1,t},..., u_{d,t} \right]^{\top}$ is assumed to be a martingale difference sequence with resepect to $\mathcal{F}_t = \sigma \big( \boldsymbol{u}_t, \boldsymbol{u}_{t-1},... \big)$ satisfying
with $\mathsf{Cov} ( u_{i,t}, u_{j,t} ) = \sigma_{i,j}$ for $i,j \in \left\{ 1,..., d \right\}$. Then, the autoregressive coefficient matrix is defined as $\boldsymbol{R} = \mathsf{diag} ( \rho_1,..., \rho_d )$. Moreover, we consider two cases:
Notice that the OLS estimator of the above model with a common explosive root is inconsistent. In particular, the standardized sample variance matrix $\sum_{t=1}^n X_t X_t^{\top}$ is asymptotically singular.
Let the $i-$th regression model be
Furthermore, define with $X_i= \big[ x_{i,1},..., x_{i,n} \big]^{\top}$ which is an $( n \times 1 )$ vector. Let $A = [ \rho_1,..., \rho_d ]^{\top}$ and $U = [ U_1,..., U_d ]^{\top}$ such that
The dependence structure of the model is given by $\mathsf{Var} ( \boldsymbol{U} ) = \boldsymbol{\Sigma}_u \otimes \boldsymbol{I}_n$ is a $\left( nd \times nd \right)$ matrix.
Therefore, in order to facilitate the development of the asymptotic theory, we begin by considering the relevant assumptions for the regressors and the error terms.
The aim of this section is to present an asymptotic theory analysis of the performance of fixed $T$ consistent estimation techniques for PVARX$(1)$ model-based on observations in first differences which can be found in the framework proposed by juodis2018first.
Consider the following PVAR$(1)$ model specification defined as below:
where $\boldsymbol{y}_{i,t}$ is an $( m \times 1 )$ vector and $\boldsymbol{\Phi}_{ m \times m}$ matrix of parameters to be estimated, where $\boldsymbol{\eta}_i$ is an $( m \times 1 )$ vector of fixed effects and $\boldsymbol{\epsilon}_{i,t}$ is an $(m \times 1)$ vector of innovations independent across $i$, with zero mean and constant covariance matrix $\boldsymbol{\Sigma}$. For various empirical applications the PVAR$(1)$ model specification might be two restrictive and incomplete. Therefore, in that case the original model specification can be extended by including strictly exogenous variables, the so-called PVARX$(1)$ model
where $\boldsymbol{x}_{i,t}$ is a $( k \times 1 )$ vector of strictly exogenous regressors and $\boldsymbol{B}$ is an $( m \times k )$ unknown parameter matrix. Furthermore, econometric models with group specific spatial dependence can be formulated as a reduced form PVARX$(1)$ (see, kripfganz2014unconditional and verdier2016estimation).
The FD transformation can be employed in order to remove any individual effects that are present in the original model in levels. Therefore, the econometric specification is as below:
Define the following variables as below:
and denote with
Therefore, after pooling observations for all $t$ and $i$, we define the pooled panel FD estimator (FD-OLS)
In the econometric model without exogenous regressors the FD-OLS estimator is given as below:
Moreover, assuming that $\boldsymbol{y}_{i,0}$ is covariance stationary and as a consequence it holds that
A large stream of literature has proposed econometric methodologies for capturing cross sectional dependence and heterogeneity via the use of dynamic panel models. Firstly, the particular literature has been significantly developed with the seminal paper of pesaran2006estimation. Moreover, kapetanios2014nonlinear present a framework for nonlinear panel models with cross-sectional dependence. Recently, in the spatial econometrics literature various methodologies have been proposed to model both spatial dependence and cross-sectional effects. Olmo2023 propose a network regression model with an estimated interaction matrix which incorporates both the cross-sectional dependence as well as the network dependence in the form of a metric distance between the set of regressors.
We assume a sample of $T$ observations for $N$ agents. Then, we specify the following model
for $t = 2,...,T$ and $i = 1,...,N$, where
This specification implies that $x_{i,t}$ is influenced by the cross-sectional average of a selection of $x_{j,t-1}$ and that in particular that the relevant $x_{j,t-1}$ are those that lie closest to $x_{i,t-1}$. The model involves a $K$ nearest neighbour mechanism, however all neighbours $x_{j,t-1}$ within a given threshold $r$ contribute equally. The formulation aims to capture the intuition that people are affected by those with whom they share common views or behaviour, reflecting the fact that similar agents are affected by similar effects (kapetanios2014nonlinear, moon2017dynamic).
In this section, we discuss the estimation methodology of the aformentioned nonlinear model proposed by kapetanios2014nonlinear. We consider the standard estimation procedure for a threshold model, whereby a grid of values for $r$ is constructed. Then, for all values on that grid the model is estimated by least squares to obtain estimates of the autoregression parameter, $\rho$. Specifically, denoting with
with $\widetilde{x}_{i} = \left( \widetilde{x}_{i,1}, ...., \widetilde{x}_{i,T - 1} \right)^{ \prime }$. The value of $r$ that minimizes the sum of of squared residuals is
A panel data model with intercept is given by
where $\nu_i \sim \textit{i.i.d} \left( 0, \sigma_v^2 \right)$. A more general version is given by
for $r \times 1$ vectors of observable variables, $\zeta_t$, and coefficients $\kappa_i = \left( \kappa_{i,1}, ...., \kappa_{i,r} \right)$, where $\kappa_{i,j} \sim \textit{i.i.d} \left( 0 , \sigma^2_{ \kappa_j } \right)$, for $j = 1,...,r$. It is worth noticing, given our interest in the persistence properties of our class of models, that it is known from the literature on linear dynamic panel data models that high persistence can be generated by moderate values of $\rho$ in combination with a large variance for the individual specific effects. We now examine the properties of the least squares estimator above. The presence of $\nu_i$ induces endogeneity in standard panel AR models, leading to biased estimation of the autoregressive parameter for finite $T$, when standard panel least squares estimators, such as the within group estimator, are used. Endogeneity arises because consistency for least squares estimators requires that
Moreover, it is straightforward to allow for higher order, $p$, lags such that
Similarly, we define with
are the cross-section averages associated with the group of neighbours and non-neighbours respectively. The particular model is more relevant in the case where we are interested to model heterogenous interactions. Furthermore, another important issue is how best to modify the basic model to decompose the slope parameter, $\rho$, into an own effect and a neighbour effect. This case, can be captured with
Similarly a time-space recursive model is formulated as below
where the weights are given by
The estimation of the last model can be conducted consists of a two step estimation procedure. First, the consistent estimate of $r$ is obtained from (2), then we construct the weights and estimate the model by least squares. Notice that the above modelling approaches involved threshold mechanisms for constructing the unit-specific cross-sectional averages. But as discussed in Section 1, the class of models we wish to propose is much more general. In particular, we consider for example models of the form
where $\gamma$ is a finite-dimensional vector of parameters and $w \left( x ; \gamma \right)$ is a positive twice differentiable integrable function such as the exponential function $\text{exp} \left( - \gamma x^2 \right)$.
Therefore, it can be also shown that
which implies that the within estimator is valid for estimating (ref) when fixed effects are incorporating in (ref). For example, the model relates to a univariate process $y_t$, $t = 1,...., T$ and the associated finite-dimensional covariates, denoted by $p_t$. Moreover, if the data are ordered, as in the case of time series, then their model is given by
The above model has the property that it can incorporate the similarity/distance $w$. Another set of extensions to the above models, arises by introducing other variables or lags to the model, either linearly
where $z_{i,t} = \left( z_{1,i,t},..., z_{k,i,t} \right)^{\prime}$ is a set of exogenous stationary variables, that is, $\beta = \left( \beta_1,...., \beta_k \right)^{\prime}$, or nonlinearly as below
\paragraph{ Proof of Theorem (ref)}
Consider the model
for $t= 2,...,T$, $i = 1,...,N$ where
and $\left\{ \epsilon_{i,t} \right\}_{ t = 1}^T$ is an error process.
\paragraph{Proof of Theorem (ref)}
We prove that the NLS estimator of $\left( \rho^0, \gamma^0 \right)$ denoted by $\left( \rho^0, \gamma^0 \right)$ is consistent and asymptotically normal. Notice that for consistency we need to establish the following two conditions:
Condition C2:
where $B( \alpha, \beta )$ is an open ball of radius $b$ centered around $a$, is satisfied. These three conditions together imply the uniform convergence of the objective function given by
to the limit objective function which is the key to establishing consistency.
In this section, we discuss in details the framework proposed by Olmo2023. Consider the cross-sectional regression model with $N$ the number of units
where $Y$ is the demeaned outcome variable, and $\mathbb{X} = \left[ X_1,..., X_L \right]$ is a vector of exogenous demeaned covariates, and $\gamma$ the vector of slope coefficients associated to $\mathbb{X}$. More precisely, Olmo2023 develop a novel econometric framework for network dependence that allows for exogenous spillover effects between the cross-sectional units. To do this, Olmo2023 define a distance measure $d: \mathbb{Z} \times \mathbb{Z} \to \mathbb{R}^{+}$, with $\mathbb{Z}$ a set of elements that reflect the network features of the model. For example, let $\mathbf{z}_i = \left( x_{1i},..., x_{Li} \right)$ and $\mathbf{z}_j = \left( x_{1j},..., x_{Lj} \right)$ be elements of this set. Then, the distance between these two elements satisfy that $d \left( \mathbf{z}_i, \mathbf{z}_j \right) \geq 0$ and $d \left( \mathbf{z}_i, \mathbf{z}_j \right) = 0$ if and only if $ \mathbf{z}_i = \mathbf{z}_j$. For example, a suitable metric to capture this is given by the Euclidean distance,
Another choice is the distance characterized by the $l_1$ norm: $\displaystyle d \left( \mathbf{z}_i, \mathbf{z}_j \right) = \sum_{l=1}^L \left| x_{l,i} - x_{l,j} \right|$.
Consider $\delta_K = \left( \gamma, \Gamma_K \right)^{\prime}$ be the OLS regression of $Y$ on $\mathbb{U}_K = \left( \mathbb{X}, \mathbb{X}_K \right)$, which is as
Using the partitioned inverse one can show that
with $\widehat{\mathbb{X}}_u = \mathbb{X} - \widehat{\mathbb{X}}$, where $\widehat{\mathbb{X}} = \mathbb{X}_K \left( \mathbb{X}_K ^{\prime} \mathbb{X}_K \right)^{-1} \mathbb{X}_K^{\prime} \mathbb{X}$. Similarly, we have that $\widehat{Y} = \mathbb{X}_K \left( \mathbb{X}_K^{\prime} \mathbb{X}_K \right)^{-1} \mathbb{X}_K^{\prime} Y$ is the projection of Y on $\mathbb{X}_K$. Therefore, the network parameters are estimated from the partitioned regression proposed by Olmo2023
Another important quantity to make statistical inference about the network parameters is the variance of $\widehat{\delta}_K$. Let $\mathbb{X}_u = \mathbb{X} - \mathbf{E} \left[ \mathbb{X} | \mathbb{X}_K \right]$, $\Phi = \mathbf{E} \left[ \mathbb{X}^{\prime} \mathbb{X} \right]$, and $\Psi = \mathbf{E} \left[ \left( \mathbb{X}_u \epsilon \right) \left( \mathbb{X}_u \epsilon \right)^{\prime} \right]$, with $\epsilon = Y - \mathbb{X} \gamma - \mathbb{X}_K \Gamma_K$.
The asymptotic variance of the standardized estimator of $\gamma$ is $\Phi^{-1} \Psi \Phi^{-1}$, which can be estimated as
with $\widehat{\Phi} = \displaystyle \frac{1}{N} \sum_{i=1}^N \widehat{ \mathbb{X} }_{i,u}^{\prime} \widehat{ \mathbb{X} }_{i,u}$ and $\widehat{\Psi} = \displaystyle \frac{1}{N} \sum_{i=1}^N e_i^2 \widehat{ \mathbb{X} }_{i,u}^{\prime} \widehat{ \mathbb{X} }_{i,u}$, where $e = Y - \mathbb{X} \widehat{\gamma} - \mathbb{X}_K \widehat{\Gamma}_K$. Under homoscedasticity of the error term, the asymptotic variance of the standardized estimator is $V \left( \widehat{\gamma} \right) = \Phi^{-1} \sigma_{\epsilon}^2 = \mathbf{E} \left[ \epsilon^2 \right]$, and the corresponding estimator is $\widehat{V} \left( \widehat{\gamma} \right) = \widehat{\Phi}^{-1} \widehat{\sigma}_{\epsilon}^2 / N$, with $ \widehat{\sigma}_{\epsilon}^2 = \frac{1}{N} \sum_{i=1}^N \displaystyle \epsilon_i^2$. Olmo2023 derive the asymptotic variance of $\widehat{\Gamma}_{K}$ is formally derived in the proof of Proposition 2 below.
and let $Q = \mathbf{E} \left[ \mathbb{X}_K^{\prime} \mathbb{X}_K \right]$. Then, it follows that,
since $\mathbf{E} \left[ \mathbb{X}_K^{\prime} \mathbb{X} \left( \gamma - \widehat{\gamma} \right) \epsilon^{\prime} \mathbb{X}_K^{\prime} \right] = 0$.
Under homoscedasticity of the error term $\epsilon$ the variance of the standardized estimator satisfies
Therefore, a suitable estimator of the variance $\widehat{\Gamma}_K$ is
where $\widehat{Q}_K = \displaystyle \frac{1}{N} \sum_{ i=1 }^N \mathbb{X}_{iK}^{\prime} \mathbb{X}_{iK}$ and $\mathbb{X}_{iK}$ are rows of the matrix $\mathbb{X}_{K}$. The random quantities of interest for estimating the presence of network effects are the parameters, $\beta_l (d)$, for $\ell = 1,...,L$ which are interpreted as realizations of the continuous and differentiable functional coefficients $\beta_{\ell} (d)$, with $d \in (0, C] \subset \mathbb{R}^{+}$. Then, the estimator of $\beta_{\ell} ( d)$ is defined as
with $v( d ) = \big[ v_1( d )^{\prime},..., v_K( d )^{\prime} \big]^{\prime}$, where
\paragraph{Proof of Theorem 1}[Olmo2023] The proof of this result follows from noting
with $\widehat{X}_u = M_{ \mathbb{X} } X$, $M_{ \mathbb{X} } = \left[ I_N - \mathbb{X}_K \left( \mathbb{X}_K^{\prime} \mathbb{X}_K \right)^{-1} \mathbb{X}_K \right]$, and $I_N$ is the identity matrix, $Y - \widehat{Y}= M_{ \mathbb{X} } Y$. Then, given that $M_{ \mathbb{X} }$ and $\mathbb{X}_K$ are orthogonal, we obtain
Now, under the assumption that $\mathbb{E} \left[ \epsilon | X \right] = 0$, applying the law of large numbers, it follows that $\frac{ X^{\prime} M_{\mathbb{X} } \epsilon }{N} \overset{ p }{ \to } \Phi$, with $\Phi$ a positive definite matrix. Therefore, $\widehat{\gamma} - \gamma = o_p(1)$, as $N \to \infty$. Then, the asymptotic convergence of the standardized estimator immediately follows by slightly adapting the proof of nonlinear additive partial models for series estimators.
\paragraph{Proof of Proposition 1}[Olmo2023] In order to prove the result in Proposition 1, it is sufficient to show that $\mathbf{E} \left[ \left\lVert \widehat{Q}_K - I_{\widetilde{K} } \right\rVert^2 \right] = \mathcal{O} \left( \frac{ \xi_0(h) }{ Nh} \right)$, as $N \to \infty$, with $I_{\widetilde{K} }$ the $\widetilde{K} \times \widetilde{K}$ identity matrix for $\widehat{K} = K \left( q + 1 \right)$. Thus, for a square symmetric square matrix $Q^{-1/2}$ of $Q^{-1}$, with $Q = I_{ \widehat{K} }$, the vector $\mathbb{X}_K (x) Q^{-1 / 2}$ is a nonsingular transformation of $\mathbb{X}_K (x)$, and thus it can be shown that
with $\bar{C}$ some positive constant. Next, we assume that $\mathbb{X}_K$ is a standardized version of our regression matrix. Let $\mathbb{X}_{ij, K}$ denote the element $( i, j)$ of the matrix $\mathbb{X}_{K}$, and $\delta_{ij}$ denote the element $( j, \ell)$ of the matrix $I_{ \widetilde{K} }$. Then, the assumption $Q = I_{ \widetilde{K} }$ implies that $\mathbf{E} \left[ \mathbb{X}_{ij, K}^{\prime} \mathbb{X}_{ij, K} \right] = \delta_{ij}$, and we have that
We have that, under Assumption A.5, $\underset{ x \in \mathcal{X}_X }{ \text{sup} } \ \leq \xi_0 \left( h \right)$, with $h \in 0$. Furthermore, $\mathbb{E} \left[ \text{trace} \big( \mathbb{X}_{K}^{\prime} \mathbb{X}_{K} \big) \right] = \text{trace} \left( I_{ \widetilde{K} } \right)$, with $\widetilde{K} = (q+1)K$. Then, $\mathbf{E} \big[ \left\lVert \widehat{Q}_K - I_{ \widetilde{K}} \right\rVert \big] \leq \xi_0 \left( h \right)^2 C ( q + 1)/ ( 2hN )$. Therefore, for $q$ and $C$ fixed, $\left\lVert \widehat{Q}_K - I_{ \widetilde{K}} \right\rVert = \mathcal{O} \left( \xi_0 \left( h \right) \left( N h \right)^{-1/2} \right)$.
Furthermore, since the smallest eigenvalue of $\widehat{Q}_K - I_{ \widetilde{K}}$ is bounded by $\left\lVert \widehat{Q}_K - I_{ \widetilde{K}} \right\rVert$, this implies that the smallest eigenvalue of $\widehat{Q}_K$ converges to one in probability. Letting $1_N$ be the indicator function for the smallest eigenvalue of $\widehat{Q}_K$ being greater than $1 / 2$, then $\mathbb{P} \left( 1_N = 1 \right) = 1$.
\paragraph{Proof of Lemma 1}[Olmo2023]
Let $X = \mathbf{E} \left[ X | \mathbb{X}_K \right] + X_u$, with $X_u$ the error term of the projection of $X$ on $\mathbb{X}_K$, and let $\widehat{X} = \mathbb{X}_K \big( \mathbb{X}_K^{\prime} \mathbb{X}_K \big) \mathbb{X}_K^{\prime} X$ be the linear projection of $X$ on $\mathbb{X}_K$. Then, $\widehat{X} = \mathbb{X}_K^{\prime} \widehat{\beta}_X$, with $\widehat{\beta}_X = \big( \mathbb{X}_K^{\prime} \mathbb{X}_K \big)^{-1} \mathbb{X}_K^{\prime} X$ the slope coefficient of the approximating regression given by $\widetilde{K}$ regressors. Then, $\widehat{X}_u = X - \widehat{X}_u$,
Simple algebra shows that
with $\beta_X$ the vector of coefficients associated to the linear prediction model $\mathbb{X}_K \beta_X$. Therefore, it follows
where $\mathbb{X}_{K,i}$ represents the $i$th row of matrix $\mathbb{X}_{K}$. The convergence in probability of the matrix $\widehat{\Phi}$ implies $\frac{1}{N} \sum_{ i = 1}^N X_{i,u}^{\prime} X_{i,u} \overset{ p }{ \to } \Phi$, where $\Phi = \mathbb{E} \left[ X_{u}^{\prime} X_{u} \right]$ by applying the law of large numbers. Next we show that
From Assumption A.6, it follows that $\left\lVert \mathbb{E} \left[ X | \mathbb{X}_K \right] - \mathbb{X}_{K} \beta_X \right\rVert^2 = \mathcal{O} \left( K^{-2 \kappa } \right)$. Now, it holds that
by Proposition 1. We now prove that $\left\lVert \beta_X - \widehat{\beta}_X \right\rVert^2 = \mathcal{O} \left( K^{- \kappa } \right)$. Notice that,
\paragraph{Proof of Theorem 1}[Olmo2023]
We prove first the consistency of the estimator $\gamma$. The estimator is defined as
where $\widehat{Y} = \mathbb{X}_K \big( \mathbb{X}_K^{\prime} \mathbb{X}_K \big)^{-1} \mathbb{X}_K^{\prime}$, the linear projection of Y on $\mathbb{X}_K$. The infeasible estimator is
Recently, the literature is considering cluster-robust methods are widely used to account for cross-sectional dependence. The standard model of cluster dependence partitions the set of observations into may independent clusters. Usually researchers use HAC variance estimators, which account for spatial or temporal dependence. According to leung2023network such simulation evidence show that for spatially or temporally dependent data, tests using HAC estimators can exhibit size distortion in smaller samples, unlike cluster-robust inference methods. Specifically, leung2023network develop a novel framework with simulated and theoretical evidence for applying cluster-robust methods to network-dependent data. Within the proposed setting of leung2023network a main econometric challenge to ensure robust estimation and inference is to obtain obtain HAC robust estimation techniques which have some special characteristics especially in the choice of the bandwidth for data is network, rather than spatially, dependent. Moreover, the simulation results of leung2023network show that there are advantages to using cluster-robust methods for network data, in comparison to to classical HAC estimators. In particular, leung2023network find that the randomization test, a leading method for cluster-robust inference with a small number of clusters, better controls size in smaller samples, provided clusters have low conductance. However, when no such clusters exist, the test, when naively applied to the output of spectral clustering, can exhibit substantial size distortion even in large samples, unlike the HAC estimator. This is because clusters in this case cannot generally satisfy the requirement of asymptotic independence, so we expect all existing cluster-robust methods to exhibit similar size distortion.
Following the framework proposed by leung2023network, observe a set of units $\mathcal{N}_n = \left\{ 1,..., n \right\}$, data $W_i \in \mathbb{R}^{d_w}$ associated with each unit $i \in \mathcal{N}_n$, and an undirected network or graph $\boldsymbol{A}$ on $\mathcal{N}_n$. We represent $\boldsymbol{A}$ as a binary, symmetric adjacency matrix with $ij-$th entry $A_{ij}$, where $A_{ij} = 1$ signifies a link between $i$ and $j$. There are no self-links, meaning $A_{ii} = 0$ for all $i$. Moreover, the settings of leung2023network treats $\boldsymbol{A}$ as fixed (conditional upon), whereas $\left\{ W_i \right\}$ is random and not necessarily indetically distributed. Let $\theta_0 \in \mathbb{R}^{d_{\theta}}$ be the estimand of interest and $g : \mathbb{R}^{d_{\theta}} \to \mathbb{R}^{d_{\theta}}$ a moment function such that
Denote with $\displaystyle G(\theta) = \frac{1}{n} \sum_{i=1}^n g( W_i, \theta )$, and $\Psi_n$ be a weighted matrix. Define the generalized method of moments (GMM) estimator
Various studies in the literature develop cluster-robust methods for GMM when $\left\{ W_i \right\}$ satisfies weak temporal or spatial dependence. However, we instead employ a notion of weak network dependence, which is analogous to mixing conditions used in time series and spatial econometrics.
Cluster-robust methods take as input a partition of $\mathcal{N}_n$ into $L$ clusters, which we denote by $\left\{ \mathcal{C}_{\ell} \right\}_{\ell = 1}^L$. Being a partition, the clusters satisfy $\cup_{\ell = 1}^L \mathcal{C}_{\ell} = \mathcal{N}_n$ and $\mathcal{C}_{\ell} \cap \mathcal{C}_m = \varnothing$ for all $\ell \neq m$. Notice that $\mathcal{C}_{\ell}$ depends on $\boldsymbol{A}$ since different networks may be partitioned differently. Furthermore, the number of "quality" clusters in a network is small, so we develop inference procedures robust to a small number of clusters. We refer to these procedures as “small-L cluster robust methods”. In particular such methods have been shown to exhibit substantially improved size control relative to conventional cluster-robust procedures. Moreover, under weak network dependence, observations in different components are independent, so components may therefore be treated as separate clusters. This implies that, if a network consists of many components, standard many-cluster asymptotics are applicable, and one can simply use standard errors on the components (see, leung2023network).
Let $\mathcal{N}_k = \left( \mathcal{C}_K, \boldsymbol{X}_K \right)$ be a $\mathcal{C}-$stationary network and $\mathcal{C}_K^* \in \mathbb{G}_k$ be its realization. Then a natural estimator of $\mathsf{Var} \left( \boldsymbol{X}_K | \mathcal{C}_K = \mathcal{C}_K^{*} \right)$ is given by
Consider formal conditions under which a given set of clusters can be used for asymptotically valid cluster-robust inference. In particular, we consider a sequence of networks with associated clusters indexed by the network size $n$, taking $n$ to infinity while keeping the number of clusters $L$ fixed. Usually the correct way to consider asymptotics is to employ a sequence of networks. Under weak network dependence and standard regularity conditions, we can show that
where $\rho_{\ell} = \mathsf{lim}_{n \to \infty} n_{\ell} / n$ and
The above asymptotic result can ensure that the vector of GMM estimates $\left\{ \sqrt{n} \left( \hat{\theta}_{\ell} - \theta_0 \right) \right\}_{\ell = 1}^L$ is asympotically normal. Conventional cluster-robust methods require independent clusters. Then small-L cluster robust methods exploit the weaker requirement of asymptotic independence, that $\boldsymbol{\Sigma}_{\ell m} = \boldsymbol{0}$ for all $\ell \neq m$. Furthermore, the symmetry of the limit distribution which, under the group of transformations corresponds to having off-diagonal blocks equal to zero. In other words, we interpret the zero off-diagonal blocks $\boldsymbol{\Sigma}_{\ell m}$ as the key requirement for the validity of cluster-robust methods.
In other words, the main assumption imposed by leung2023network regarding dependence between observations is that they are weakly dependent, meaning that random variables approach independence as the distance between their locations grows. Moreover, our methods involve partitioning the data into groups defined by the researcher. We define $G_N$ to be the total number of groups and index them by $g = 1,..., G_N$. The OLS estimator can be written as below
using group-level notation. Therefore, the most common approach to inference with weakly dependent data is to use a plug-in estimator, call it $\widetilde{V}_N$, of the variance matrix of $x_{s_i} \varepsilon_{s_i}$, along with the usual large-sample approximation for the distribution of $\widehat{\beta}_N$ such that
where the long-run covariance matrix is defined as below
where $Q$ is the limit of the second moment matrix for $x$. The typical method uses the sample average of $x_{s_i} x_{s_i}^{\prime}$ to estimate $Q$ and plugs-in a consistent estimators, $\widetilde{V}_N$, of $V$ to arrive at the approximation
We consider a sequence of networks and associated clusters, both implicitly indexed by the network size $n$. Recall that $n_{\ell} = \left| \mathcal{C}_{\ell} \right|$, the size of cluster $\ell$.
Furthermore, based on the notions proposed by kojevnikov2021bootstrap, leung2023network consider a formal notion of weak network dependence called $\psi-$dependence. The particular notion of dependence is analogous to familiar notions of temporal or spatial weak dependence, except distance between observations is measured using path distance. In other words, weak dependence simply means that the correlation between two sets of observations decays as the network distance between the sets grows.
Therefore, for any $H, H^{\prime} \subset \mathcal{N}_n$, define with
the distance between the two sets. Let $\mathcal{L}_d$ be the set of bounded $\mathbb{R}-$valued Lipschitz functions on $\mathbb{R}^d$, on $\left\lVertf\right\rVert_{\infty} = \mathsf{sup}_x | f(x) |$, $\mathsf{Lip} (f)$, the Lipschitz constant of $f \in \mathcal{L}_d$, and
the set of pairs of sets $H, H^{\prime}$ with respective sizes $h, h^{\prime}$ that are at least distances $s$ apart in the network.
Moreover, leung2023network define the $i'$s $s-$neighbourhood boundary $\mathcal{N}^{\partial}_{\boldsymbol{A}} (i,s) = \left\{ j \in \mathcal{N}_n : \ell_{\boldsymbol{A}} (i,j) = s \right\}$, and its $k-$th moment is given by
\paragraph{Open Problems} The last few decades the time series econometrics literature has considered applications of moderate deviations from unit root in univariate autoregressive models and multivariate predictive regression models as well as in settings such as panel data . In particular, predictive regression models with regressors generated as stable or unstable autoregressive processes has recently seen a growing attention in the econometrics and statistics literature. A related open problem include the development of a formal econometric framework for cluster-based predictive regression models. In this direction, our objective is to develop an identification and estimation method that allows us to develop the asymptotic theory for a suitable system estimator under the presence of both nonstationarity and network dependence. To establish a robust estimation and inference procedure, developing central limit theory that explicitly accounts for the dependence between the cross-sectional and time series data will be essential, such as the notation of stable dependence (see, anatolyev2021limit). On the other hand, in our framework we assume that the network induced dependence corresponds to the nodes of a network and the corresponding time series observations. We remain agnostic regarding the structure of the network as well as the form of the distance between nodes which will require to impose further metric space assumptions. As a result, the proposed approach allows to further generalize to a data structure that allows for both a cross-sectional network type dependence in the time series dimension.
Consider a random sample $\mathcal{X}_n = \left\{ \boldsymbol{X}_1,..., \boldsymbol{X}_n \right\}$ from an unknown $p-$dimensional distribution depending on a scalar parameter $\theta$. In particular, the aim is to construct a $( 1 - \alpha )-$level upper confidence bound for $\theta$, based on some appropriate point estimator $\hat{\theta}_n$ for $\theta$. Denote by $\mathcal{X}^{*}_n = \left\{ \boldsymbol{X}^{*}_1,..., \boldsymbol{X}^{*}_n \right\}$ a random sample of size $n$ taken with replacement from $\mathcal{X}_n$ and let $\hat{\theta}_n^{*}$ be the equivalent function of $\mathcal{X}^{*}_n$ as $\hat{\theta}_n$ is to $\mathcal{X}_n$. Let $\hat{\sigma}_n^{*}$ be the bootstrap version of $\hat{\sigma}_n$. Then, the standard percentile $( 1 - \alpha )-$level bootstrap confidence bounds for $\theta$ can be written as below:
where $\xi_{n, \alpha}$ is the $\alpha-$quantile of the bootstrap distribution of the standardized $\hat{\theta}_n$ such that
where $\mathbb{P}^{*}$ refers to the conditional probability law of $\mathcal{X}_n^{*}$ given $\mathcal{X}_n$. Generally, it holds that
Following the framework proposed by yan2022factor, the forecast error is given by $\left( \hat{y}_{T+h|T} - y_{T+h|T} \right)$. Notice that when forecasting, one would be more interested in the distribution of the forecast error. Since $Y_{T+h} = y_{T+h|T} + \varepsilon_{T+h}$, it follows that the forecasting error is given by
Thus, if $\hat{\varepsilon}_{T+h}$ is asymptotically normal with
with $\mathsf{var} \left( \hat{y}_{T+h|T} \right) = B_T^2$.
Therefore, a confidence interval can be obtained by replacing $\sigma^2$ by its consistent estimate, $\frac{1}{T-h} \sum_{t=1}^{T-h} \hat{\varepsilon}^2_{T+h}$. Therefore, we can get the $95\%$ confidence interval for the conditional mean $y_{T+h|T}$ such that
Therefore, the $95\%$ confidence interval for the forecasting variable $y_{T+h}$ is given by
Moreover, the framework of yan2022factor allows to test for the presence of threshold effects of our model. The question of interest is whether the nonlinear term
enters the regression model, that is, whether $\delta_T = 0$. If $\gamma_0$ were known, the traditional Lagrange multiplier statistic and Wald statistic would be good choices to solve this issue. However, $\gamma_0$ is usually unknown and not identified under the null hypothesis. Therefore, we propose a sup-Wald statistic which extends the seminal work of Hansen to test the linearity of the model, as it does not require prior knowledge of $\gamma_0$. Therefore, to facilitate the establishment of distributional theory, we consider a local-to-null reparametrisation: $\delta_T = \frac{c}{\sqrt{T}}$. Thus, based on this model specification the null hypothesis is $H_0: c = 0$ with alternative $H_1: c \neq 0$. Therefore, for each $\gamma \in \Gamma = \left[ \underline{\gamma}, \bar{\gamma} \right]$, we obtain the estimator $\hat{\beta}(\gamma)$ and $\hat{\delta}(\gamma)$. Then, we build a sup-Wald statistic to test the presence of threshold effects such that
where $R = \big[ 0, I_q \big]$ and $q$ denotes the dimension of $z_t$. Thus, we obtain that
Moreover, the endogenous threshold regression model (ETR) has attracted much attention in recent econometric practice. This is due to the fact that economic relationships may shift over time.
Suppose the first-stage regression is given by
where the instruments $\boldsymbol{z}$ contain both exogenous regressors such as 1 and $q$, and excluded exogenous regressors, $\mathbb{E} \left( \boldsymbol{v} | \boldsymbol{z} \right) = 0$ and $\mathbb{E} \left( u | \boldsymbol{z} \right) = 0$. Then, by taking the conditional expectation we obtain the following expression
where $\theta = \left( \theta^{\prime}, \Pi^{\prime} \right)^{\prime}$. Then, the estimator proposed by caner2004instrumental for the endogenous threshold variable $\gamma$ minimizes the sample analogue of the following unconditional condition
Consider the following GMM estimators which use the moment conditions given below
to identify $\gamma$. Although GMM estimators are essential in handling endogeneity, compared with $M-$estimators, they suffer from at least three drawbacks. First GMM changes the nature of $\gamma$ from a threshold point (which is nonregular) to a quantile of $q$ (which is regular), which implies that the convergence rate of $\widehat{\gamma}$ is $n^{1/2}$, much slower than the convergence rate $n$ of $M-$estimators. Second, $\gamma$ is not always identiable by GMM, for example, when $q$ is independent of the rest of the system such as the time index in structural change models, $\gamma$ cannot be identified by GMM. Third, GMM requires more instruments than our CF estimators for identification, which implies that GMM may have less applicability since good instruments are hard to find in practice. Specifically, the model becomes nonlinear threshold regression, and $\gamma$ can be estimated by minimizing the objective function below
where $\widehat{\mathsf{g}}_{xi} = \widehat{\Pi}_x^{\prime} \boldsymbol{z}_i$ and $\widehat{\mathsf{g}}_{qi} = \widehat{\pi}^{\prime} \boldsymbol{z}_i$, such that the estimators $\widehat{\Pi}_x$ and $\widehat{\pi}$ are obtained from a first-stage regression. The estimation procedure of the parameter $\gamma$ requires to regress $y_i$ on $\mathsf{g}_{xi} \boldsymbol{1} \left\{ q_i \leq \gamma \right\}$ and $\mathsf{g}_{xi} \boldsymbol{1} \left\{ q_i > \gamma \right\}$ where
to obtain $\widehat{\beta}_1 (\gamma), \widehat{\beta}_2 (\gamma)$ and $\widehat{\kappa} (\gamma)$. Then, $\gamma$ can be estimated by the extremum problem as below
Given $\widehat{\gamma}$, the parameter $\beta$ can be estimated by 2SLS/GMMs as in CH. In particular, in the small-threshold-effect framework of hansen2000sample, KST show that $\widehat{\gamma}$ is $n^{1 - 2 \alpha}-$consistent and its asymptotic distribution is based on a functional of two-sided Brownian motion, under the assumption that both $\delta_{\beta}$ and $\kappa$ are $\mathcal{O}\left( n^{ 1 - 2 \alpha} \right)$ with $\alpha \in (0, 1/2)$. Now assume that the unknown parameters $( \gamma, \delta, \kappa )^{\prime}$ lie a compact set with their true value in the interior. Moreover, we define centered versions of $S_n$ and $S$ as
Assume that $y_{t+h}$ follows a factor-augmented regression model (see, bai2006confidence) given by
where $W_t$ is a vector of observed regressors (including for instance lags of $y_t$), which jointly with $F_t$ help forecast $y_{t+h}$. Then, the $k-$dimensional vector $F_t$, describes the common latent factors in the panel factor model, $X_{it} = \lambda_i^{\prime} F_t + e_{it}, \ \ \ i = 1,...,N, \ t = 1,...,T$, where the $r \times 1$ vector $\lambda_i$ contains the factor loadings and $e_{it}$ is an idiosyncratic error term. Thus, we can forecast $y_{t+h}$ or its conditional mean $y_{T+h|T} = \alpha^{\prime} F_T + \beta^{\prime} W_T$ using the pair $\big\{ ( y_t, X_t, W_t ) : t = 1,...,T \big\}$, the available data at time $T$. Since factors are not observed, the diffusion index forecast approach typically involves a two-step procedure (see, gonccalves2017bootstrap).
Notice that the above result can be employed to construct point-wise confidence intervals of $f(x)$ at a fixed $x$. To assess shapes of density functions so that one can perform goodness-of-fit tests, however, one needs to construct uniform or simultaneous confidence bands (SCB). Therefore, to do this we need to deal with the maximum absolute deviation over some interval $[ \ell, u ]$:
In practice we first need to study the asymptotic uniform distributional theory for the NW estimator $\mu_n(x)$> Specifically, one needs to find the asymptotic distribution for $\underset{ x \in T }{ \mathsf{sup} } \left| \mu_n(x) - \mu(x) \right|, \ \ \text{where} \ T = \left[ \ell, u \right] $. Then, building on the aforementioned result one can construct an asymptotic $( 1 - \alpha )$ SCB, where $0 < \alpha < 1$, by finding two functions $\mu_n^{ \mathsf{lower} } (x)$ and $\mu_n^{ \mathsf{upper} } (x)$, such that the following holds
In the standard case the SCB can be used for model validation: one can test whether $\mu(.)$ is of certain parametric functional form by checking whether the fitted parametric form lies in the SCB.
Let $m = \floor{ n^{\tau} }$, where $\delta_1 / \gamma < \tau < 1 - \delta_1$, and
Consequently, we can show that
Set the following random quantity
Constructing robust confidence intervals and simultaneous confidence bands has important applications. In particular, it is related to the statistical problem of understanding the trends of extremes and variability of climate variables. By interpreting climate extremes as upper and lower quantiles and climate variability as interpercentile ranges, the nonparametric quantile estimation provides a simple and effective means to address the latter problem (e.g., see li2022simultaneous). The simultaneous confidence band (SCB) is a classical tool for nonparametric inference.
To construct a 100 $(1 - \beta) \%$ SCB for $Q_{\alpha}(.)$, one finds two functions $\ell$ and $u$ depending on $\left( X_i \right)_{i=1}^n$, such that the following condition holds
In this section we briefly discuss the framework of lee2021factor who proposed a novel two-regime regression model where regime switching is driven by a vector of possibly unobservable factors. When the factors are latent, these are estimated by the principle component analysis of a panel data set. Then, lee2021factor show that the optimization problem can be reformulated as mixed integer optimization, and we present two alternative computational algorithms. Moreover, they derive the asymptotic distribution of the resulting estimator under the scheme that the threshold effect shrinks to zero.
Suppose that $y_t$ is generated from
where $x_t$ and $f_t$ are adapted to the filtration $\mathcal{F}_{t-1}$, $( \beta_0, \delta_0, \gamma_0 )$ is a vector of unknown parameters and the unobserved random variable $\varepsilon_t$ satisfies the conditional mean restriction. Thus the above model specification is related to the literature on thereshold models with unknown change points.
In particular, the regression function can be written as:
When the factor $f_t$ is latent, we estimate it using PCA from a potentially much larger dataset, whose dimension is $N$. We illustrate that the asymptotic distribution for the estimator $\alpha_0 \equiv \left( \beta_0^{\prime}, \delta_0^{\prime} \right)^{\prime}$ is identical to that when $\gamma_0$ were known regardless of whether factors are directly observable or not; therefore, the estimator of $\alpha_0$ enjoys an oracle property (see, lee2021factor).
Furthermore, the asymptotic properties of the proposed estimator are established by adopting a diminishing threshold effect. In other words, we assume that $\delta_0 = T^{ - \varphi } d_0$ for some unknown $\phi \in (0, 1/2)$ and unknown nondiminishing vector $d_0$. Specifically, the unknown parameter $\varphi$ reflects the difficulty of estimating $\gamma_0$ and affects the identification and estimation of the change-point $\gamma_0$. Both the rate of convergence and the asymptotic distribution depends on $\varphi$. In particular, when factors are directly observable, then the distribution of the estimator of $\gamma_0$ is given by
Consider the linear dynamic panel data model given by
as a special case of the econometric specification studied by liu2020forecasting. Then, a suitable estimator for $\hat{\rho}$ is given by the truncated instrumental variable (IV) estimator such that
where $M_N$ is a sequence that slowly diverges to infinity. Define the residuals $\hat{U}_{it} = Y_{it} - \hat{\lambda} \hat{\rho}_{IV} - \hat{\rho}_{IV} \hat{Y}_{it-1}$
\paragraph{Pooled-OLS Predictor}
Ignoring the heterogeneity in the $\lambda_i$'s and imposing that $\lambda_i = \lambda$ for all $i$, we can define that
A key point here is to ensure that the empirical model is able to accurately predict bank revenues and balance sheet characteristics under observed macroeconomic conditions.
\paragraph{Compound Risk}
\paragraph{Ratio Optimality}
The convergence is uniform with respect to the correlated random effects distributions $\pi$ in some set $\Pi$. The uniformity holds in the neighborhood of point masses for which the prior and posterior variances of $\lambda_i$ are zero. Thus, the convergence statement covers the case of $\lambda_i$ being homogenous across $i$. Then, the autoregressive coefficient in the basic dynamic panel model can be $\sqrt{N}-$consistently estimated, which suggests that $\sum_{i=1}^N \left( \hat{\rho} - \rho \right)^2 Y_{iT}^2 = \mathcal{O}_p(1)$. Then, the discrepancy between the predictor $\widehat{Y}_{iT+1}$ and $\lambda_i + \rho Y_{it}$ can be decomposed into three terms is given by
Additional discussion on some aspects to identification, estimation and inference for high dimensional models is presented in katsouris2023high. Further related studies from the perspective of shrinkage-based estimation using GMM methods include among others cheng2015select. In this section we focus on some specific applications to panel data regression models such as the framework proposed by lu2016shrinkage. The literature on discovering latent structures in panel data models include su2016identifying as well as the study of bing2022inference.
The particular study considers the problem of determining the number of factors and selecting the proper regressors in linear dynamic panel data models with interactive fixed effects (see, lu2016shrinkage). The authors propose a methodology for simultaneous selection of regressors and factors and estimation through the method of adaptive group Lasso. In particular, we show that with probability approaching one, the proposed method correctly select all relevant regressors and factors and shrink the coefficients of irrelevant regressors and redundant factors to zero.
Therefore, the model in matrix form can be written as below:
Without loss of generality we assume that only the first $K_0$ elements of $X_{it}$ have nonzero slope coefficients, and express the vector of regressors $X_{it} = \big( X_{(1)it}^{\prime}, X_{(2)it}^{\prime} \big)^{\prime}$, where $X_{(1)it}^{\prime}$ and $X_{(2)it}^{\prime}$ are $( K_0 \times 1 )$ and $( K - K_0) \times 1$ vectors respectively, and the true coefficients of $X_{(1)it}^{\prime}$ are nonzero while those of $X_{(2)it}^{\prime}$ are assumed to be all zero. Acoordingly, we decompose the parameter vector $\beta^0$ as $\beta^0 = \left( \beta^{0 \prime}_{ (1) }, \beta^{0 \prime}_{ (2) } \right)^{\prime} = \left( \beta^{0 \prime}_{ (1) }, 0 \right)$. Consider the Gaussian QMLE $\left( \tilde{\beta}, \tilde{\lambda}, \tilde{F} \right)$ of $\left( \beta^0, \lambda^0, F^0 \right)$ which is given by
where
where $\beta \equiv \left( \beta_1,..., \beta_K \right)^{\prime}$ is a $(K \times 1)$ vector and $F \equiv \left( F_1,..., F_T \right)$ is a $(T \times R)$ and $\lambda \equiv \left( \lambda_1,..., \lambda_N \right)^{\prime}$ is an $(N \times R)$ matrix. Therefore, further aspects of the statistical optimization methodology will need to tackle the estimation of the unknown factors in the high dimensional regression model.
In particular, a possible solution is the use of group Lasso for selection of the number of factors. Thus, considering the simultaneous variable and factor selection in a dynamic panel data model along with relevant asymptotic theory analysis is the aim of this section.
Consider the following expression
Notice that we define with $\left\lVertA\right\rVert^2_{\mathsf{sp}} \equiv \lambda_1 \left( A^{\prime} A \right)$, which implies that $\left\lVertA\right\rVert_{\mathsf{sp}} \equiv \sqrt{ \lambda_1 \left( A^{\prime} A \right) }$, where $\lambda_1$ denotes the largest eigenvalue of a real symmetric matrix $A$. Moreover, we use $\lambda_{ \mathsf{min} } (A)$ and $\lambda_{ \mathsf{max} } (A)$ to denote the smallest and largest eigenvalues of symmetric matrix $A$, respectively.
Denote with
We also partition the variance matrix such that $V_{NT} \equiv \mathsf{diag} \big( V_{(11), NT}, V_{(22), NT} \big)$.
In this section we study estimation and variable selection methodologies for high dimensional linear regression models under the presence of endogeneity. In particular, this setting requires to apply an IV estimation approach to tackle the aspect of endogeneity. Several studies in the literature consider the construction of confidence intervals for estimators of interest, under the assumption that all the IVs are valid after controlling for the said covariates. Moreover, in invalid IV settings different statistical frameworks were developed to provide robust inferential methods. In other words, a propose statistical procedure under the presence of endogeneity especially in a high dimensional estimation setting should be valid even under the presence of invalid instruments.
In terms of the asymptotic theory intially we need to examine the stability conditions of the system such that $\mathsf{max}_{ i = 1,..., v } \left| \lambda_i \left( \boldsymbol{C} \right) \right| < 1$, then for any $\boldsymbol{d}$ and all $t = 1,2,3,...$ the norm
Therefore, the above theorem can be applied such that $\boldsymbol{V}_{\theta}$ and $\boldsymbol{V}_{\gamma}$ can be consistently estimated by replacing the unknown parameters of the estimator $\boldsymbol{\beta}$ by their Gaussian maximum likelihood values and substituting a consistent estimate for $\boldsymbol{\Sigma}_{\varepsilon}$, while dropping the expectation from $\boldsymbol{V}_{\theta}$. In practice we are interested to obtain simple estimates of the covariance matrix $\boldsymbol{\Sigma}_{\varepsilon}$.