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.
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 fusetime
summary(object, ...)
# S3 method for fusetime
formula(x, ...)
# S3 method for fusetime
df.residual(object, ...)
# S3 method for fusetime
print(x, ...)
# S3 method for fusetime
coef(object, ...)
# S3 method for fusetime
residuals(object, ...)
# S3 method for fusetime
fitted(object, ...)
An object of class fusetime holding
modela data.frame containing the dependent and explanatory variables as well as cross-sectional and time indices,
coefficientslet \(p^{(1)}\) denote the number of time-varying coefficients and \(p^{(2)}\) the number of time constant parameters. A list holding (i) a \(T \times p^{(1)} \times \hat{K}\) array of the post-Lasso group-specific functional coefficients and (ii) a \(K \times p^{(2)}\) matrix of time-constant post-Lasso estimates.
groupsa list containing (i) the total number of groups \(\hat{K}\) and (ii) a vector of estimated group memberships \((\hat{g}_1, \dots, \hat{g}_N)\), where \(\hat{g}_i = k\) if \(i\) is assigned to group \(k\),
residualsa vector of residuals of the demeaned model,
fitteda vector of fitted values of the demeaned model,
argsa list of additional arguments,
ICa list containing (i) the value of the IC, (ii) the employed tuning parameter \(\lambda\), and (iii) the MSE,
convergencea list containing (i) a logical variable if convergence was achieved and (ii) the number of executed ADMM algorithm iterations,
callthe function call.
An object of class fusetime has print, summary, fitted, residuals, formula, df.residual, and coef S3 methods.
a formula object describing the model to be estimated.
a data.frame or matrix holding a panel data set. If no index variables are provided, the panel must be balanced and ordered in the long format \(\bold{Y}=(Y_1^\prime, \dots, Y_N^\prime)^\prime\), \(Y_i = (Y_{i1}, \dots, Y_{iT})^\prime\) with \(Y_{it} = (y_{it}, \bold{x}_{it}^\prime)^\prime\). Conversely, if data is not ordered or not balanced, data must include two index variables that declare the cross-sectional unit \(i\) and the time period \(t\) of each observation.
a character vector holding two strings. The first string denotes the name of the index variable identifying the cross-sectional unit \(i\) and the second string represents the name of the variable declaring the time period \(t\). The data is automatically sorted according to the variables in index, which may produce errors when the time index is a character variable. In case of a balanced panel data set that is ordered in the long format, index can be left empty if the number of time periods n_periods is supplied.
the number of observed time periods \(T\). If an index character vector is passed, this argument can be left empty. Default is NULL.
the tuning parameter determining the strength of the penalty term. Either a single \(\lambda\) or a vector of candidate values can be passed. If a vector is supplied, a BIC-type IC automatically selects the best fitting \(\lambda\) value.
the polynomial degree of the B-splines. Default is 3.
the number of interior knots of the B-splines. If left unspecified, the default heuristic \(M = \text{floor}((NT)^{\frac{1}{7}} - \log(p))\) following Haimerl et al. (2025) is used.
the minimum group cardinality as a fraction of the total number of individuals \(N\). In case a group falls short of this threshold, each of its members is allocated to one of the remaining groups according to the MSE. Default is 0.05.
a character vector containing the variable names of explanatory variables that enter with time-constant coefficients.
the a non-negative weight used to obtain the adaptive penalty weights. Default is 2.
the maximum number of iterations for the ADMM estimation algorithm. Default is \(5*10^4\).
the tolerance limit for the stopping criterion of the iterative ADMM estimation algorithm. Default is \(1*10^{-10}\).
the tolerance limit for within-group differences. Two individuals are assigned to the same group if the Frobenius norm of their coefficient vector difference is below this threshold. Default is \(1*10^{-3}\).
the tuning parameter balancing the fitness and penalty terms in the IC that determines the penalty parameter \(\lambda\). If left unspecified, the heuristic \(\rho = 0.07 \frac{\log(NT)}{\sqrt{NT}}\) of Haimerl et al. (2025) is used. We recommend the default.
the non-negative Lagrangian ADMM penalty parameter. For the employed penalized sieve estimation PSE, the \(\varrho\) value is not very influential. We recommend the default 1.
logical. If TRUE, helpful warning messages are shown. Default is TRUE.
logical. If TRUE, certain operations are parallelized across multiple cores. Default is TRUE.
ellipsis
of class fusetime.
of class fusetime.
Paul Haimerl
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.
Haimerl, P., Smeekes, S., & Wilms, I. (2025). Estimation of latent group structures in time-varying panel data models. arXiv preprint arXiv:2503.23165. tools:::Rd_expr_doi("10.48550/arXiv.2503.23165").
# 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)
Run the code above in your browser using DataLab