Generalized Functional Linear Models

• GFLMs are regression models that connect scalar responses to functional predictors using integral operators and known link functions within an RKHS context.
• They unify classical L² models with impact-point and projection methods, supporting applications in spectroscopy, neuroimaging, and spatial analysis.
• Modern estimation techniques include roughness-penalized likelihood, wavelet-based methods, and sparse kernel regressions for enhanced interpretability and accuracy.

Searching arXiv for recent and foundational papers on generalized functional linear models and related variants. Generalized Functional Linear Models (GFLMs) are regression models in which a scalar response is linked to one or more functional covariates through an integral linear predictor and a known link function, typically in exponential-family settings. In the canonical scalar-on-function form, the model is written as

$$g\!\bigl(E[Y\mid X]\bigr)=\alpha+\int_T X(t)\,\beta(t)\,dt,$$

with intercept $\alpha$, slope function $\beta$, and link $g$. The term also appears in a broader RKHS-based formulation that subsumes classical $L^2$ functional linear models, finite-marginal “impact-point” models, and projection models, and it has recently been reused in scientific machine learning for kernelized sparse integral surrogates of neural operators (Shang et al., 2014, Berrendero et al., 2020, Arzani et al., 2023).

1. Foundational formulations and scope

In the contemporary statistical literature, a GFLM usually denotes a model with scalar response $Y$ and square-integrable functional predictor $X(\cdot)$ such that

$$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$

or equivalently

$$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$

This template covers functional linear regression, functional logistic regression, and functional Poisson regression, among others (Chen et al., 13 Nov 2025).

A more general formulation replaces the classical $L^2$ inner product by an RKHS representation built from the covariance kernel $\alpha$0. With centered $\alpha$1 and $\alpha$2, one writes

$\alpha$3

and, by slight abuse of notation,

$\alpha$4

Here $\alpha$5 is Loève’s isometry from the closed linear span of the stochastic process to the RKHS $\alpha$6. In this view, the classical model

$\alpha$7

is only one special case, recovered when the RKHS representer belongs to $\alpha$8 (Berrendero et al., 2020).

This RKHS perspective unifies three families that are often treated separately. First, the standard $\alpha$9 model appears when $\beta$0. Second, finite-marginal or impact-point models arise when $\beta$1, giving

$\beta$2

Third, finite-projection models are obtained from expansions in an orthonormal basis, including the FPCR case when the basis is formed by covariance eigenfunctions. A common misconception is therefore that the functional linear model is intrinsically an $\beta$3 inner-product model; the RKHS definition shows that marginals, projections, and the classical integral form all live inside the same umbrella space $\beta$4 (Berrendero et al., 2020).

Empirical illustrations in this broader framework included Tecator spectra, sugar spectroscopy, and population under-14 trajectories. In those studies, low-dimensional marginal fits with $\beta$5–$\beta$6 were reported to give $\beta$7 as high as, or higher than, $\beta$8-PC fits with $\beta$9–$g$0 components, underscoring the practical relevance of the generalized definition (Berrendero et al., 2020).

2. Penalization, basis reduction, and computational estimators

A central estimation strategy for GFLMs is roughness-penalized likelihood in an RKHS or Sobolev space. One maximizes

$g$1

with $g$2, $g$3. In the Shang–Cheng construction, the covariance operator and roughness penalty are simultaneously diagonalized, a reproducing kernel is obtained explicitly, and asymptotically valid confidence intervals for regression mean, prediction intervals for future response, and several hypothesis-testing procedures follow from a functional Bahadur representation (Shang et al., 2014).

An alternative computational route is wavelet penalized likelihood for generalized functional models with a semiparametric predictor

$g$4

In this formulation, scalar covariates are modeled linearly while the functional effect is represented by an unknown nonparametric function. The estimator uses a maximum penalized likelihood criterion with an $g$5 penalty on wavelet coefficients, combined with backfitting and Fisher scoring. The paper establishes quasi-minimax rates under standard conditions and shows that the LASSO penalty yields adaptive estimation with respect to the regularity of the estimated function (Gannaz, 2011).

Penalized-spline mixed-model representations provide another widely used approximation. Functional predictors are reduced through FPCA, the slope function is expanded in B-splines, and the roughness penalty is rewritten as a variance component. The resulting GLMM form supports restricted maximum likelihood estimation and variance-component tests for nullity, functionality, and linearity. This framework was developed specifically to accommodate dense and sparse functional data with responses from exponential-family distributions (Chen et al., 2019).

Residual-based Alternative Partial Least Squares (RAPLS) extends the alternative PLS algorithm iteratively to accommodate additional scalar covariates and non-continuous outcomes. It combines IRLS working responses, residualization against scalar covariates, and empirical component construction from the residual covariance operator. The paper establishes the convergence rate of the RAPLS estimator for the unknown slope function and, after an additional calibration step, proves asymptotic normality and efficiency of the calibrated estimator for scalar parameters. In an ADNI application predicting Alzheimer’s disease progression from neuroimaging, leave-one-out cross-validation accuracy was reported as $g$6 for RAPLS, compared with $g$7 for FPCR, $g$8 for plsRglm, $g$9 for LDA on top-50 PCs, and $L^2$0 for Random Forest with 1,000 voxels subsampled (Wang et al., 1 May 2025).

3. Inference, hypothesis testing, and goodness-of-fit

The inferential theory for GFLMs has developed beyond point estimation. In the roughness-regularized RKHS framework, asymptotically valid confidence intervals for the regression mean and prediction intervals for future response are available, as are functional-contrast tests and penalized likelihood-ratio tests. One notable result is a new type of Wilks phenomenon for testing functional linear models: after suitable normalization, the null-limit chi-square law is free of nuisance parameters. The same work also constructs an adaptive global test whose null limit is Gumbel and whose separation rate is minimax up to a $L^2$1 factor (Shang et al., 2014).

Model-selection-oriented testing has also been formulated through mixed-effects spline representations. In that setting, nullity corresponds to $L^2$2, functionality to $L^2$3, and linearity to $L^2$4. Each null can be recast as a variance-component hypothesis $L^2$5 against $L^2$6, and approximate RLRT procedures can be used for Gaussian, Bernoulli, Poisson, and Binomial responses. Simulation results reported that the approximate RLRT maintains type I error close to the nominal $L^2$7 even under moderate sparsity and achieves substantially higher power than competing approaches that simultaneously test fixed and random effects (Chen et al., 2019).

Formal goodness-of-fit testing for GFLMs has remained comparatively limited, and a recent proposal addresses this gap through projection averaging. The test targets

$L^2$8

and builds a U-statistic from a Cramér–von–Mises metric integrated over all one-dimensional projections of the functional predictor. The projection-averaging device is intended to mitigate the curse of dimensionality, and the resulting statistic has asymptotic normality under the null and consistency under alternatives. Because the null variance is cumbersome, the authors provide bootstrap procedures for both continuous and discrete responses. Their simulations report empirical size close to nominal at $L^2$9, $Y$0, and $Y$1, together with uniformly higher power than PCvM, $Y$2, and $Y$3 under the alternatives considered (Chen et al., 13 Nov 2025).

4. Sparse observation, irregular design, and measurement error

A major difficulty for GFLMs is that functional covariates are often latent, observed sparsely, or contaminated by noise. MISFIT addresses generalized functional linear regression with sparsely and irregularly sampled data by treating the unobserved trajectories as missing data and using multiple imputation. The model assumes latent smooth curves $Y$4, noisy measurements $Y$5, and a generalized functional linear predictor

$Y$6

Under the Gaussian-process latent-curve model, imputations are drawn from the conditional distribution of $Y$7 or its FPC scores given $Y$8 and the sparse observations, and the completed-data estimates are combined by Rubin’s rules. The paper proves consistency in the linear case without requiring the number of observations per curve to increase with sample size, and emphasizes that this avoids the dense-design assumption. In the motivating macrocephaly study, many EHR curves had $Y$9 as small as $X(\cdot)$0–$X(\cdot)$1, and the method associated pathology with both the overall head circumference and the velocity of head growth (Petrovich et al., 2018).

Measurement error in functional and scalar covariates has motivated two additional strands. One develops SIMEX and regression calibration for generalized functional linear regression models with a mixture of function-valued and scalar-valued covariates prone to classical measurement error. After spline-basis reduction, SIMEX inflates the estimated error covariance over a grid of $X(\cdot)$2 values and extrapolates estimates to $X(\cdot)$3, while regression calibration replaces latent covariates by best linear predictors. In simulations, the oracle had the smallest bias and MSE, SIMEX nearly achieved oracle bias, RC outperformed the average and naïve estimators, and SIMEX was reported to be robust across squared-exponential, AR(1), compound symmetric, and unstructured covariance settings (Luan et al., 2023).

A related multi-level formulation treats repeated longitudinal functional-data replicates $X(\cdot)$4 with heteroscedastic Gaussian-process errors in a two-stage scalable regression-calibration scheme. Stage 1 uses pointwise generalized linear mixed models to estimate latent curves through either UP_MEM or MP_MEM; stage 2 fits the usual scalar-on-function exponential-family model after basis expansion. In extensive simulations, the reported ordering was MP_MEM $X(\cdot)$5 UP_MEM $X(\cdot)$6 PACE $X(\cdot)$7 Average $X(\cdot)$8 Naive in both bias and MSE, with MP_MEM uniformly closest to Oracle. This framework was applied to assess the relationship between physical activity and type 2 diabetes in community-dwelling adults in the United States who participated in NHANES (Luan et al., 2023).

5. Spatial dependence and structured response models

Standard GFLMs assume independent responses, but scalar responses may also exhibit spatial dependence. A spatial generalized functional linear model extends the independence-based formulation by placing a functional covariate $X(\cdot)$9 at each spatial location $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$0 and defining the marginal mean through

$$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$1

while adding a centered neighbor interaction controlled by a scalar spatial dependence parameter $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$2. When $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$3, the model reduces to an ordinary GFLM (Kim et al., 2024).

Dimension reduction is performed by basis expansion and truncation,

$$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$4

followed by a composite likelihood estimating equation to handle the spatial dependence. Under a repeating-lattice asymptotic context, the model yields consistency, quadratic-form normality, a confidence interval for $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$5, and a uniform confidence band for $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$6. In practice, the Godambe information is estimated empirically from composite score and Hessian analogues (Kim et al., 2024).

The binary-conditionals example in the paper used a $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$7 torus with $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$8 sites, $$E\{Y\mid X\}=F\!\Bigl(\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt\Bigr),$$9 independent lattices per Monte Carlo case, and $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$0 cases. Reported findings included nearly unbiased $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$1, empirical CI coverage of about $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$2–$$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$3 as $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$4 grows, and relative FMSE values of $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$5, $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$6, $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$7, and $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$8 of the non-spatial GFLM’s FMSE as $$g\bigl(E[Y\mid X]\bigr)=\alpha_0+\int_0^1 X(t)\beta_0(t)\,dt.$$9. In the Midwestern U.S. corn-yield application, the estimated spatial dependence was $L^2$0 with $L^2$1 CI $L^2$2–$L^2$3, and the estimated $L^2$4 indicated a significantly positive effect of temperature in June and a significantly negative effect of high temperatures from mid-July through mid-August (Kim et al., 2024).

6. Kernelized GFLMs as interpretable surrogates in scientific machine learning

A distinct recent usage of the term appears in physics-based machine learning, where generalized functional linear models are proposed as interpretable surrogates for trained deep learning models. Motivated by functional data analysis, the construction augments the classical linear integral operator with three ingredients: a nonlinear link $L^2$5, a lifting map $L^2$6 on the inputs, and a prespecified library of candidate kernels $L^2$7. For image-to-scalar prediction, the model is

$L^2$8

The library includes Gaussian kernels, exponential kernels, univariate exponential decays, polynomial-weight kernels, and indicator-truncated kernels, sampled over a range of bandwidths; sparse regression is then used to discover a small number of active analytic terms (Arzani et al., 2023).

The estimation problem is linear in the coefficients once the library is fixed. After assembling the design matrix from candidate terms

$L^2$9

the coefficients are selected through a Lasso-type criterion,

$\alpha$00

implemented in practice by a sequentially-thresholded least-squares algorithm described as a SINDy-style sweep. The resulting surrogate is an explicit sparse sum of integral equations and can be used in two modes: post-hoc interpretation of a trained neural network, or direct interpretable operator learning from raw data without a neural network (Arzani et al., 2023).

Across six test cases in solid mechanics, porous-media flow, and fluid-mechanics super-resolution, the paper reports training-error parity with neural networks in most cases, with mean absolute errors within $\alpha$01–$\alpha$02 of the neural network in distribution. It also reports that, under out-of-distribution shifts, deep-network error typically doubles or triples whereas the GFLM error may increase by only $\alpha$03–$\alpha$04. Local-surrogate fitting around a narrow parameter slice of a convolutional autoencoder was further used to recover an analytic model that reproduces, and in the reported case slightly refines, local network behavior. In this literature, the phrase “generalized functional linear model” therefore denotes not only a statistical regression model but also a sparse, kernelized, closed-form surrogate architecture for interpretable operator learning (Arzani et al., 2023).