Title:	Joint Estimation of Latent Groups and Group-Specific Coefficients in (Time-Varying) Panel Data Models
Version:	1.1.4
Maintainer:	Paul Haimerl <paul.haimerl@econ.au.dk>
Description:	Latent group structures are a common challenge in panel data analysis. Disregarding group-level heterogeneity can introduce bias. Conversely, estimating individual coefficients for each cross-sectional unit is inefficient and may lead to high uncertainty. This package addresses the issue of unobservable group structures by implementing the pairwise adaptive group fused Lasso (PAGFL) by Mehrabani (2023) <doi:10.1016/j.jeconom.2022.12.002>. PAGFL identifies latent group structures and group-specific coefficients in a single step. On top of that, we extend the PAGFL to time-varying coefficient functions (FUSE-TIME), following Haimerl et al. (2025) <doi:10.48550/arXiv.2503.23165>.
Depends:	R (≥ 3.5.0)
License:	AGPL (≥ 3)
Encoding:	UTF-8
RoxygenNote:	7.3.2
LinkingTo:	Rcpp, RcppArmadillo, RcppParallel, RcppThread
Imports:	Rcpp, lifecycle, ggplot2, RcppParallel
BugReports:	https://github.com/Paul-Haimerl/PAGFL/issues
URL:	https://github.com/Paul-Haimerl/PAGFL
Suggests:	testthat (≥ 3.0.0)
Config/testthat/edition:	3
NeedsCompilation:	yes
Packaged:	2025-11-17 14:42:22 UTC; au772358
Author:	Paul Haimerl [aut, cre], Stephan Smeekes [ctb], Ines Wilms [ctb], Ali Mehrabani [ctb]
Repository:	CRAN
Date/Publication:	2025-11-17 21:51:43 UTC

Fused Unobserved group Spline Estimation of TIME varying coefficients

Description

Estimate a time-varying panel data model subject to a latent group structure using FUSE-TIME–Fused Unobserved group Spline Estimation of TIME varying coefficients–by Haimerl et al. (2025). FUSE-TIME jointly identifies the latent group structure and group-specific time-varying functional coefficients. The time-varying coefficients are approximated as polynomial B-splines. The function supports both static and dynamic panel data models.

Usage

fuse_time(
  formula,
  data,
  index = NULL,
  n_periods = NULL,
  lambda,
  d = 3,
  M = floor(length(y)^(1/7) - log(p)),
  min_group_frac = 0.05,
  const_coef = NULL,
  kappa = 2,
  max_iter = 50000,
  tol_convergence = 1e-10,
  tol_group = 0.001,
  rho = 0.04 * log(N * n_periods)/sqrt(N * n_periods),
  varrho = 1,
  verbose = TRUE,
  parallel = TRUE,
  ...
)

tv_pagfl(
  formula,
  data,
  index = NULL,
  n_periods = NULL,
  lambda,
  d = 3,
  M,
  min_group_frac = 0.05,
  const_coef = NULL,
  kappa = 2,
  max_iter = 50000,
  tol_convergence = 1e-10,
  tol_group = 0.001,
  rho,
  varrho = 1,
  verbose = TRUE,
  parallel = TRUE,
  ...
)

## S3 method for class 'fusetime'
summary(object, ...)

## S3 method for class 'fusetime'
formula(x, ...)

## S3 method for class 'fusetime'
df.residual(object, ...)

## S3 method for class 'fusetime'
print(x, ...)

## S3 method for class 'fusetime'
coef(object, ...)

## S3 method for class 'fusetime'
residuals(object, ...)

## S3 method for class 'fusetime'
fitted(object, ...)

Arguments

Details

Consider the grouped time-varying panel data model subject to a latent group structure

y_{it} = \gamma_i^0 + \bold{\beta}^{0\prime}_{i} (t/T) \bold{x}_{it} + \epsilon_{it}, \quad i = 1, \dots, N, \; t = 1, \dots, T,

where y_{it} is the scalar dependent variable, \gamma_i^0 is an individual fixed effect, \bold{x}_{it} is a p \times 1 vector of explanatory variables, and \epsilon_{it} denotes a zero mean error. The p-dimensional coefficient vector \bold{\beta}_{i}^0 (t/T) contains smooth functions of time and follows the latent group pattern

\bold{\beta}_i^0 \left(\frac{t}{T} \right) = \sum_{k = 1}^K \bold{\alpha}_k^0 \left( \frac{t}{T} \right) \bold{1} \{i \in G_k^0 \},

with \cup_{k = 1}^K G_k^0 = \{1, \dots, N\}, G_k^0 \cap G_j^0 = \emptyset for any k \neq j, k,j = 1, \dots, K.

The time-varying coefficient functions are estimated as polynomial B-splines. To this end, let \bold{b}(v) denote a M + d +1 vector of polynomial basis functions with the polynomial degree d and M interior knots. Then, \bold{\beta}_i^0 (t/T) is approximated by forming linear combinations of these basis functions \bold{\beta}_i^0 (t/T) \approx \bold{\Pi}_i^{0 \prime} \bold{b} (t/T), where \bold{\Pi}_i^{0} is a (M + d + 1) \times p matrix of spline control points.

To estimate \bold{\Pi}_i^{0}, we project the explanatory variables onto the spline basis system, resulting in the (M + d + 1)p \times 1 regressor vector \bold{z}_{it} = \bold{x}_{it} \otimes \bold{b}(v). Subsequently, the DGP can be reformulated as

y_{it} = \gamma_i^0 + \bold{\pi}_{i}^{0 \prime} \bold{z}_{it} + u_{it},

where \bold{\pi}_i^0 = \text{vec}(\bold{\Pi}_i^0), and u_{it} = \epsilon_{it} + \eta_{it} collects the idiosyncratic \epsilon_{it} and the sieve approximation error \eta_{it}.

Following Haimerl et al. (2025, sec. 2), FUSE-TIME jointly estimates the functional coefficients and the group structure by minimizing the criterion

F_{NT} (\bold{\pi}, \lambda) = \frac{1}{NT} \sum^N_{i=1} \sum^{T}_{t=1}(\tilde{y}_{it} - \bold{\pi}_{i}^\prime \tilde{\bold{z}}_{it})^2 + \frac{\lambda}{N} \sum_{i = 1}^{N - 1} \sum_{j = i+1}^N \dot{\omega}_{ij} \| \bold{\pi}_i - \bold{\pi}_j \|_2

with respect to \bold{\pi} = (\bold{\pi}_i^\prime, \dots, \bold{\pi}_N^\prime)^\prime. \tilde{a}_{it} = a_{it} - T^{-1} \sum^{T}_{t=1} a_{it}, a = \{y, \bold{z}\} to concentrate out the individual fixed effects \gamma_i^0 (within-transformation). \lambda is the penalty tuning parameter and \dot{w}_{ij} denotes adaptive penalty weights which are obtained by a preliminary non-penalized estimation. The criterion function is minimized via an iterative alternating direction method of multipliers (ADMM) algorithm (Haimerl et al. 2053, Appendix C).

Two individuals are assigned to the same group if \| \hat{\bold{\pi}}_i - \hat{\bold{\pi}}_j \|_2 \leq \epsilon_{\text{tol}} (and hence \hat{\bold{\xi}}_k = \hat{\bold{\pi}}_i = \hat{\bold{\pi}}_j for some k = 1, \dots, \hat{K}), where \epsilon_{\text{tol}} is determined by tol_group. The time-varying coefficients are then retrieved by taking \hat{\bold{\beta}}_i (t/T) = \hat{\bold{\Pi}}_i^\prime \bold{b}(t/T), where \hat{\bold{\pi}}_i = \text{vec}(\hat{\bold{\Pi}}_i) (analogously \hat{\bold{\alpha}}_k (t/T) = \hat{\bold{\Xi}}_k^\prime \bold{b}(t/T), using \hat{\bold{\xi}}_k = \text{vec}(\hat{\bold{\Xi}}_k)).

Subsequently, the estimated number of groups \hat{K} and group structure follow by examining the number of distinct elements in \hat{\bold{\pi}}. Given an estimated group structure, it is straightforward to obtain post-Lasso estimates \hat{\bold{\alpha}}^p_k (t/T) = \hat{\bold{\Xi}}^{p \prime}_k \bold{b}(t/T) for each k = 1, \dots, \hat{K} using group-wise least squares (see grouped_tv_plm).

We recommend choosing a \lambda tuning parameter by passing a logarithmically spaced grid of candidate values with a lower limit close to 0 and an upper limit that leads to a fully homogeneous panel. A BIC-type information criterion then automatically selects the best fitting \lambda value.

In case of an unbalanced panel data set, the earliest and latest available observations per group define the start and end-points of the interval on which the group-specific time-varying coefficients are defined.

We refer to Haimerl et al. (2025) for more details.

Value

An object of class fusetime holding

An object of class fusetime has print, summary, fitted, residuals, formula, df.residual, and coef S3 methods.

Author(s)

Paul Haimerl

References

Haimerl, P., Smeekes, S., & Wilms, I. (2025). Estimation of latent group structures in time-varying panel data models. arXiv preprint arXiv:2503.23165. doi:10.48550/arXiv.2503.23165.

Examples

# Simulate a time-varying panel with a trend and a group pattern
set.seed(1)
sim <- sim_tv_DGP(N = 10, n_periods = 50, intercept = TRUE, p = 1)
df <- data.frame(y = c(sim$y))

# Run FUSE-TIME
estim <- fuse_time(y ~ ., data = df, n_periods = 50, lambda = 10, parallel = FALSE)
summary(estim)

Grouped Panel Data Model

Description

Estimate a panel data model subject to an observed group structure. Slope parameters are homogeneous within groups but heterogeneous across groups. This function supports both static and dynamic panel data models, with or without endogenous regressors.

Usage

grouped_plm(
  formula,
  data,
  groups,
  index = NULL,
  n_periods = NULL,
  method = "PLS",
  Z = NULL,
  bias_correc = FALSE,
  rho = 0.07 * log(N * n_periods)/sqrt(N * n_periods),
  verbose = TRUE,
  parallel = TRUE,
  ...
)

## S3 method for class 'gplm'
print(x, ...)

## S3 method for class 'gplm'
formula(x, ...)

## S3 method for class 'gplm'
df.residual(object, ...)

## S3 method for class 'gplm'
summary(object, ...)

## S3 method for class 'gplm'
coef(object, ...)

## S3 method for class 'gplm'
residuals(object, ...)

## S3 method for class 'gplm'
fitted(object, ...)

Arguments

Details

Consider the grouped panel data model

y_{it} = \gamma_i^0 + \bold{\beta}^{0 \prime}_{i} \bold{x}_{it} + \epsilon_{it}, \quad i = 1, \dots, N, \; t = 1, \dots, T,

where y_{it} is the scalar dependent variable, \gamma_i^0 is an individual fixed effect, \bold{x}_{it} is a p \times 1 vector of (weakly) exogenous explanatory variables, and \epsilon_{it} denotes a zero mean error. The coefficient vector \bold{\beta}_i^0 follows the observed group pattern

\bold{\beta}_i^0 = \sum_{k = 1}^K \bold{\alpha}_k^0 \bold{1} \{i \in G_k \},

with \cup_{k = 1}^K G_k = \{1, \dots, N\}, G_k \cap G_j = \emptyset and \| \bold{\alpha}_k^0 - \bold{\alpha}_j^0 \| \neq 0 for any k \neq j, k,j = 1, \dots, K. The group structure G_1, \dots, G_K is determined by the argument groups.

Using PLS, the group-specific coefficients of group k, \, k = 1, \dots, K, are obtained by OLS

\hat{\bold{\alpha}}_k = \left( \sum_{i \in G_k} \sum_{t = 1}^T \tilde{\bold{x}}_{it} \tilde{\bold{x}}_{it}^\prime \right)^{-1} \sum_{i \in G_k} \sum_{t = 1}^T \tilde{\bold{x}}_{it} \tilde{y}_{it},

where \tilde{a}_{it} = a_{it} - T^{-1} \sum_{t=1}^T a_{it}, a = \{y, \bold{x}\} to concentrate out the individual fixed effects \gamma_i^0 (within-transformation).

In case of PGMM, the slope coefficients are derived as

\hat{\alpha}_k = \left( \left[ \sum_{i \in G_k} T^{-1} \sum_{t = 1}^T \bold{z}_{it} \Delta \bold{x}_{it} \right]^\prime \bold{W}_k \left[ \sum_{i \in G_k} T^{-1} \sum_{t = 1}^T \bold{z}_{it} \Delta \bold{x}_{it} \right] \right)^{-1}

\quad \quad \left[ \sum_{i \in G_k} T^{-1} \sum_{t = 1}^T \bold{z}_{it} \Delta \bold{x}_{it} \right]^\prime \bold{W}_k \left[ \sum_{i \in G_k} T^{-1} \sum_{t = 1}^T \bold{z}_{it} \Delta y_{it} \right],

where \bold{W}_k is a q \times q p.d. symmetric weight matrix and \Delta denotes the first difference operator \Delta \bold{x}_{it} = \bold{x}_{it} - \bold{x}_{it-1} (first-difference transformation).

Value

An object of class gplm holding

A gplm object has print, summary, fitted, residuals, formula, df.residual, and coef S3 methods.

Author(s)

Paul Haimerl

References

Dhaene, G., & Jochmans, K. (2015). Split-panel jackknife estimation of fixed-effect models. The Review of Economic Studies, 82(3), 991-1030. doi:10.1093/restud/rdv007.

Mehrabani, A. (2023). Estimation and identification of latent group structures in panel data. Journal of Econometrics, 235(2), 1464-1482. doi:10.1016/j.jeconom.2022.12.002.

Examples

# Simulate a panel with a group structure
set.seed(1)
sim <- sim_DGP(N = 20, n_periods = 80, p = 2, n_groups = 3)
y <- sim$y
X <- sim$X
groups <- sim$groups
df <- cbind(y = c(y), X)

# Estimate the grouped panel data model
estim <- grouped_plm(y ~ ., data = df, groups = groups, n_periods = 80, method = "PLS")
summary(estim)

# Lets pass a panel data set with explicit cross-sectional and time indicators
i_index <- rep(1:20, each = 80)
t_index <- rep(1:80, 20)
df <- data.frame(y = c(y), X, i_index = i_index, t_index = t_index)
estim <- grouped_plm(
  y ~ .,
  data = df, index = c("i_index", "t_index"), groups = groups, method = "PLS"
)
summary(estim)

Grouped Time-varying Panel Data Model

Description

Estimate a time-varying panel data model subject to an observed group structure. Coefficient functions are homogeneous within groups but heterogeneous across groups. Time-varying coefficient functions are approximated as polynomial B-splines. The function supports both static and dynamic panel data models.

Usage

grouped_tv_plm(
  formula,
  data,
  groups,
  index = NULL,
  n_periods = NULL,
  d = 3,
  M = floor(length(y)^(1/7) - log(p)),
  const_coef = NULL,
  rho = 0.04 * log(N * n_periods)/sqrt(N * n_periods),
  verbose = TRUE,
  parallel = TRUE,
  ...
)

## S3 method for class 'tv_gplm'
summary(object, ...)

## S3 method for class 'tv_gplm'
formula(x, ...)

## S3 method for class 'tv_gplm'
df.residual(object, ...)

## S3 method for class 'tv_gplm'
print(x, ...)

## S3 method for class 'tv_gplm'
coef(object, ...)

## S3 method for class 'tv_gplm'
residuals(object, ...)

## S3 method for class 'tv_gplm'
fitted(object, ...)

Arguments

Details

Consider the time-varying panel data model

y_{it} = \gamma_i^0 + \bold{\beta}^{0\prime}_{i} (t/T) \bold{x}_{it} + \epsilon_{it}, \quad i = 1, \dots, N, \; t = 1, \dots, T,

where y_{it} is the scalar dependent variable, \gamma_i^0 is an individual fixed effect, \bold{x}_{it} is a p \times 1 vector of explanatory variables, and \epsilon_{it} is a zero mean error. The p-dimensional coefficient vector \bold{\beta}_{i}^0 (t/T) contains smooth functions of time and follows the observed group pattern

\bold{\beta}_i^0 \left(\frac{t}{T} \right) = \sum_{k = 1}^K \bold{\alpha}_k^0 \left( \frac{t}{T} \right) \bold{1} \{i \in G_k \},

with \cup_{k = 1}^K G_k = \{1, \dots, N\}, G_k \cap G_j = \emptyset for any k \neq j, k,j = 1, \dots, K. The group structure G_1, \dots, G_K is determined by the argument groups.

The time-varying coefficient functions in \bold{\alpha}_k (t/T) and, in turn, \bold{\beta}_i (t/T) are estimated as polynomial B-splines. To this end, let \bold{b}(v) denote a M + d +1 vector of polynomial basis functions with the polynomial degree d and M interior knots. \bold{\alpha}_k^0 (t/T) is approximated by forming linear combinations of the basis functions \bold{\alpha}_k^0 (t/T) \approx \bold{\Xi}_k^{0 \prime} \bold{b} (t/T), where \bold{\Xi}_i^{0} is a group-specific (M + d + 1) \times p matrix of spline control points.

To estimate \bold{\Xi}_k^{0}, we project the explanatory variables onto the spline basis system, resulting in the (M + d + 1)p \times 1 regressor vector \bold{z}_{it} = \bold{x}_{it} \otimes \bold{b}(v). Subsequently, the DGP can be reformulated as

y_{it} = \gamma_i^0 + \bold{\pi}_{i}^{0 \prime} \bold{z}_{it} + u_{it},

where \bold{\pi}_i^0 = \text{vec}(\bold{\Pi}_i^0) and \bold{\Xi}_k^0 = \bold{\Pi}_i^0 if i \in G_k. u_{it} = \epsilon_{it} + \eta_{it} collects the idiosyncratic \epsilon_{it} and the sieve approximation error \eta_{it}. Then, we obtain \hat{\bold{\xi}}_k = \text{vec}( \hat{\bold{\Xi}}_k) by

\hat{\xi}_k = \left( \sum_{i \in G_k} \sum_{t = 1}^T \tilde{\bold{z}}_{it} \tilde{\bold{z}}_{it}^\prime \right)^{-1} \sum_{i \in G_k} \sum_{t = 1}^T \tilde{\bold{z}}_{it} \tilde{y}_{it},

with \tilde{a}_{it} = a_{it} - T^{-1} \sum_{t = 1}^T a_{it}, a = \{y, \bold{z}\} to concentrate out the fixed effect \gamma_i^0 (within-transformation). Lastly, \hat{\bold{\alpha}}_k (t/T) = \hat{\bold{\Xi}}_k^{\prime} \bold{b} (t/T). We refer to Haimerl et al. (2025, sec. 2) for more details.

In case of an unbalanced panel data set, the earliest and latest available observations per group define the start and end-points of the interval on which the group-specific time-varying coefficients are defined.

Value

An object of class tv_gplm holding

An object of class tv_gplm has print, summary, fitted, residuals, formula, df.residual, and coef S3 methods.

Author(s)

Paul Haimerl

References

Haimerl, P., Smeekes, S., & Wilms, I. (2025). Estimation of latent group structures in time-varying panel data models. arXiv preprint arXiv:2503.23165. doi:10.48550/arXiv.2503.23165.

Examples

# Simulate a time-varying panel with a trend and a group pattern
set.seed(1)
sim <- sim_tv_DGP(N = 10, n_periods = 50, intercept = TRUE, p = 2)
df <- data.frame(y = c(sim$y), X = sim$X)
groups <- sim$groups

# Estimate the time-varying grouped panel data model
estim <- grouped_tv_plm(y ~ ., data = df, n_periods = 50, groups = groups)
summary(estim)

Pairwise Adaptive Group Fused Lasso

Description

Estimate a panel data model subject to a latent group structure using the pairwise adaptive group fused Lasso (PAGFL) by Mehrabani (2023). The PAGFL jointly identifies the group structure and group-specific slope parameters. The function supports both static and dynamic panels, with or without endogenous regressors.

Usage

pagfl(
  formula,
  data,
  index = NULL,
  n_periods = NULL,
  lambda,
  method = "PLS",
  Z = NULL,
  min_group_frac = 0.05,
  bias_correc = FALSE,
  kappa = 2,
  max_iter = 10000,
  tol_convergence = 1e-08,
  tol_group = 0.001,
  rho = 0.07 * log(N * n_periods)/sqrt(N * n_periods),
  varrho = max(sqrt(5 * N * n_periods * p)/log(N * n_periods * p) - 7, 1),
  verbose = TRUE,
  parallel = TRUE,
  ...
)

## S3 method for class 'pagfl'
print(x, ...)

## S3 method for class 'pagfl'
formula(x, ...)

## S3 method for class 'pagfl'
df.residual(object, ...)

## S3 method for class 'pagfl'
summary(object, ...)

## S3 method for class 'pagfl'
coef(object, ...)

## S3 method for class 'pagfl'
residuals(object, ...)

## S3 method for class 'pagfl'
fitted(object, ...)

Arguments

Details

Consider the panel data model

y_{it} = \gamma_i^0 + \bold{\beta}^{0 \prime}_{i} \bold{x}_{it} + \epsilon_{it}, \quad i = 1, \dots, N, \; t = 1, \dots, T,

where y_{it} is the scalar dependent variable, \gamma_i^0 is an individual fixed effect, \bold{x}_{it} is a p \times 1 vector of weakly exogenous explanatory variables, and \epsilon_{it} is a zero mean error. The coefficient vector \bold{\beta}_i^0 follows the latent group pattern

\bold{\beta}_i^0 = \sum_{k = 1}^K \bold{\alpha}_k^0 \bold{1} \{i \in G_k^0 \},

with \cup_{k = 1}^K G_k^0 = \{1, \dots, N\}, G_k^0 \cap G_j^0 = \emptyset and \| \bold{\alpha}_k^0 - \bold{\alpha}_j^0 \| \neq 0 for any k \neq j, k,j = 1, \dots, K.

The PLS method jointly estimates the latent group structure and group-specific coefficients by minimizing the criterion

Q_{NT} (\bold{\beta}, \lambda) = \frac{1}{T} \sum^N_{i=1} \sum^{T}_{t=1}(\tilde{y}_{it} - \bold{\beta}^\prime_i \tilde{\bold{x}}_{it})^2 + \frac{\lambda}{N} \sum_{i = 1}^{N - 1} \sum_{j=i}^N \dot{\omega}_{ij} \| \bold{\beta}_i - \bold{\beta}_j \|

with respect to \bold{\beta} = (\bold{\beta}_1^\prime, \dots, \bold{\beta}_N^\prime)^\prime. \tilde{a}_{it} = a_{it} - T^{-1} \sum_{t = 1}^T a_{it}, a = \{y, \bold{x}\} to concentrate out the individual fixed effects \gamma_i^0 (within-transformation). \lambda is the penalty tuning parameter and \dot{\omega}_{ij} reflects adaptive penalty weights (see Mehrabani, 2023, eq. 2.6). \| \cdot \| denotes the Frobenius norm. The adaptive weights \dot{w}_{ij} are obtained by a preliminary individual least squares estimation. The criterion function is minimized via an iterative alternating direction method of multipliers (ADMM) algorithm (see Mehrabani, 2023, sec. 5.1).

PGMM employs a set of instruments Z to control for endogenous regressors. Using PGMM, \bold{\beta} is estimated by minimizing

Q_{NT}(\bold{\beta}, \lambda) = \sum^N_{i = 1} \left[ \frac{1}{N} \sum_{t=1}^T \bold{z}_{it} (\Delta y_{it} - \bold{\beta}^\prime_i \Delta \bold{x}_{it}) \right]^\prime \bold{W}_i \left[\frac{1}{T} \sum_{t=1}^T \bold{z}_{it}(\Delta y_{it} - \bold{\beta}^\prime_i \Delta \bold{x}_{it}) \right]

\quad + \frac{\lambda}{N} \sum_{i = 1}^{N - 1} \sum_{j=i+1}^N \ddot{\omega}_{ij} \| \bold{\beta}_i - \bold{\beta}_j \|.

\ddot{\omega}_{ij} are obtained by an initial GMM estimation. \Delta gives the first differences operator \Delta y_{it} = y_{it} - y_{i t-1}. \bold{W}_i represents a data-driven q \times q weight matrix. I refer to Mehrabani (2023, eq. 2.10) for more details. Again, the criterion function is minimized using an efficient ADMM algorithm (Mehrabani, 2023, sec. 5.2).

Two individuals are assigned to the same group if \| \hat{\bold{\beta}}_i - \hat{\bold{\beta}}_j \| \leq \epsilon_{\text{tol}} (and hence \hat{\bold{\alpha}}_k = \hat{\bold{\beta}}_i = \hat{\bold{\beta}}_j for some k = 1, \dots, \hat{K}), where \epsilon_{\text{tol}} is determined by tol_group. Subsequently, the estimated number of groups \hat{K} and group structure follows by examining the number of distinct elements in \hat{\bold{\beta}}. Given an estimated group structure, it is straightforward to obtain post-Lasso estimates using group-wise least squares or GMM (see grouped_plm).

We recommend identifying a suitable \lambda parameter by passing a logarithmically spaced grid of candidate values with a lower limit close to 0 and an upper limit that leads to a fully homogeneous panel. A BIC-type information criterion then selects the best fitting \lambda value.

Value

An object of class pagfl holding

A pagfl object has print, summary, fitted, residuals, formula, df.residual, and coef S3 methods.

Author(s)

Paul Haimerl

References

Dhaene, G., & Jochmans, K. (2015). Split-panel jackknife estimation of fixed-effect models. The Review of Economic Studies, 82(3), 991-1030. doi:10.1093/restud/rdv007.

Mehrabani, A. (2023). Estimation and identification of latent group structures in panel data. Journal of Econometrics, 235(2), 1464-1482. doi:10.1016/j.jeconom.2022.12.002.

Examples

# Simulate a panel with a group structure
set.seed(1)
sim <- sim_DGP(N = 20, n_periods = 80, p = 2, n_groups = 3)
y <- sim$y
X <- sim$X
df <- cbind(y = c(y), X)

# Run the PAGFL procedure
estim <- pagfl(y ~ ., data = df, n_periods = 80, lambda = 0.5, method = "PLS")
summary(estim)

# Lets pass a panel data set with explicit cross-sectional and time indicators
i_index <- rep(1:20, each = 80)
t_index <- rep(1:80, 20)
df <- data.frame(y = c(y), X, i_index = i_index, t_index = t_index)
estim <- pagfl(
  y ~ .,
  data = df, index = c("i_index", "t_index"), lambda = 0.5, method = "PLS"
)
summary(estim)

Simulate a Panel With a Group Structure in the Slope Coefficients

Description

Construct a static or dynamic, exogenous or endogenous panel data set subject to a group structure in the slope coefficients with optional AR(1) or GARCH(1,1) innovations.

Usage

sim_DGP(
  N = 50,
  n_periods = 40,
  p = 2,
  n_groups = 3,
  group_proportions = NULL,
  error_spec = "iid",
  dynamic = FALSE,
  dyn_panel = lifecycle::deprecated(),
  q = NULL,
  alpha_0 = NULL
)

Arguments

Details

The scalar dependent variable y_{it} is generated according to the panel data model

y_{it} = \gamma_i + \bold{\beta}_i^\prime \bold{x}_{it} + u_{it}, \quad i = 1, \dots, N, \quad t = 1, \dots, T.

\gamma_i represents individual fixed effects and \bold{x}_{it} a p \times 1 vector of regressors. The individual slope coefficient vectors \bold{\beta}_i follow the group structure

\bold{\beta}_i = \sum_{k = 1}^K \bold{\alpha}_k \bold{1} \{i \in G_k\},

with \cup_{k = 1}^K G_k = \{1, \dots, N\}, G_k \cap G_j = \emptyset and \| \bold{\alpha}_k - \bold{\alpha}_j \| \neq 0 for any k \neq j, k,j = 1, \dots, K. The total number of groups K is determined by n_groups.

If a panel data set with exogenous regressors is generated (set q = NULL), the explanatory variables are simulated according to

x_{it,j} = 0.2 \gamma_i + e_{it,j}, \quad \gamma_i,e_{it,j} \sim i.i.d. N(0, 1), \quad j = 1, \dots, p,

where e_{it,j} denotes a series of innovations. \gamma_i and e_i are independent of each other.

In case alpha_0 = NULL, the group-level slope parameters \bold{\alpha}_{k} are drawn from \sim Unif[-2, 2].

If a dynamic panel is specified (dynamic = TRUE), the AR coefficients \beta^{\text{AR}}_i are drawn from a uniform distribution with support (-1, 1) and x_{it,j} = e_{it,j}. Moreover, the individual fixed effects enter the dependent variable via (1 - \beta^{\text{AR}}_i) \gamma_i to account for the autoregressive dependency. We refer to Mehrabani (2023, sec 6) for details.

When specifying an endogenous panel (set q to q \geq p), the e_{it,j} correlate with the cross-sectional innovations u_{it} by a magnitude of 0.5 to produce endogenous regressors (\text{E}(u | \bold{X} ) \neq 0). However, the endogenous regressors can be accounted for by exploiting the q instruments in \bold{Z}, for which \text{E}(u| \bold{Z}) = 0 holds. The instruments and the first stage coefficients are generated in the same fashion as \bold{X} and \bold{\alpha} when q = NULL.

The function nests, among other, the DGPs employed in the simulation study of Mehrabani (2023, sec. 6).

Value

A list holding

Author(s)

Paul Haimerl

References

Mehrabani, A. (2023). Estimation and identification of latent group structures in panel data. Journal of Econometrics, 235(2), 1464-1482. doi:10.1016/j.jeconom.2022.12.002.

Examples

# Simulate DGP 1 from Mehrabani (2023, sec. 6)
set.seed(1)
alpha_0_DGP1 <- matrix(c(0.4, 1, 1.6, 1.6, 1, 0.4), ncol = 2)
DGP1 <- sim_DGP(
  N = 50, n_periods = 20, p = 2, n_groups = 3,
  group_proportions = c(.4, .3, .3), alpha_0 = alpha_0_DGP1
)

Simulate a Time-varying Panel With a Group Structure in the Slope Coefficients

Description

Construct a time-varying panel data set subject to a group structure in the slope coefficients with optional AR(1) innovations.

Usage

sim_tv_DGP(
  N = 50,
  n_periods = 40,
  intercept = TRUE,
  p = 1,
  n_groups = 3,
  d = 3,
  dynamic = FALSE,
  group_proportions = NULL,
  error_spec = "iid",
  locations = NULL,
  scales = NULL,
  polynomial_coef = NULL,
  sd_error = 1
)

Arguments

Details

The scalar dependent variable y_{it} is generated according to the time-varying panel data model

y_{it} = \gamma_i + \bold{\beta}^\prime_{i} (t/T) \bold{x}_{it} + u_{it}, \quad i = 1, \dots, N, \; t = 1, \dots, T,

where \gamma_i is an individual fixed effect and \bold{x}_{it} is a p \times 1 vector of explanatory variables. The p \times 1 coefficient vector \bold{\beta}_i (t/T) follows the group pattern

\bold{\beta}_i \left( \frac{t}{T} \right) = \sum_{k = 1}^K \bold{\alpha}_k \left( \frac{t}{T} \right) \bold{1} \{i \in G_k \},

with \cup_{k = 1}^K G_k = \{1, \dots, N\} and G_k \cap G_j = \emptyset for any k \neq j, k,j = 1, \dots, K. The total number of groups K is determined by n_groups.

The predictors are simulated as:

x_{it,j} = 0.2 \gamma_i + e_{it,j}, \quad \gamma_i,e_{it,j} \sim i.i.d. N(0, 1), \quad j = \{1, \dots, p\},

where e_{it,j} denotes a series of innovations. \gamma_i and e_i are independent of each other.

The errors u_{it} feature a iid standard normal distribution.

In case locations = NULL, the location parameters are drawn from \sim Unif[0.3, 0.9]. In case scales = NULL, the scale parameters are drawn from \sim Unif[0.01, 0.09]. In case polynomial_coef = NULL, the polynomial coefficients are drawn from \sim Unif[-20, 20] and normalized so that all coefficients of one polynomial sum up to 1. The final coefficient function follows as \bold{\alpha}_k (t/T) = 3 * F(t/T, location, scale) + \sum_{j=1}^d a_j (t/T)^j, where F(\cdot, location, scale) denotes a cumulative logistic distribution function and a_j reflects a polynomial coefficient.

Value

A list holding

Author(s)

Paul Haimerl

Examples

# Simulate a time-varying panel subject to a time trend and a group structure
set.seed(1)
sim <- sim_tv_DGP(N = 20, n_periods = 50, p = 1)
y <- sim$y
X <- sim$X