Functional Analysis of Variance (FANOVA)

• FANOVA is a comprehensive framework that partitions functional data into interpretable components, including main effects, interactions, and higher-order terms.
• It employs rigorous techniques such as Hilbert-space formulations, basis expansions, and wavelet estimation to achieve precise variance decomposition.
• Applications span testing mean curves, forecasting time series, sensitivity analysis, and interpretable machine learning across diverse disciplines.

Functional Analysis of Variance (FANOVA, often written fANOVA) denotes a class of decompositions and inferential procedures for functions rather than scalars. In one widely used sense, it represents a mapping as a sum of effects of increasing complexity—an overall mean, main effects, pairwise interactions, and higher-order terms. In another, it provides models and test statistics for comparing mean functional responses across groups, treatments, times, or repeated conditions. A recurring theme across the literature is that FANOVA is not distribution-free: the decomposition is unique once a unique distribution is assigned to the covariates, and its orthogonality, variance decomposition, and interpretation depend on that choice (Borgonovo et al., 2018, Hu et al., 2023).

1. Core definitions and canonical decompositions

A standard functional-ANOVA expansion for a prediction function is

$$f(x) = f_0 + \sum_{j=1}^p f_j(x_j) + \sum_{1 \le j < k \le p} f_{jk}(x_j, x_k) + \cdots,$$

with hierarchical centering or orthogonality conditions ensuring that each term captures variation not already assigned to lower-order components. In the recursive form used in explainability work,

$$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$

so the component indexed by a subset $S$ depends explicitly on the distribution $P$ used to integrate out the complementary variables (Fumagalli et al., 2024).

In functional-data analysis proper, FANOVA often starts from a model for functional responses. For group $i$ and replicate $j$, one basic formulation is

$y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$

with the global null

$H_0:\ \mu_1(t)=\mu_2(t)=\cdots=\mu_k(t).$

This formulation treats each observed trajectory as a realization of an underlying smooth stochastic process and shifts the ANOVA question from scalar means to mean curves over a continuum (Vsevolozhskaya et al., 2014).

A closely related usage decomposes arrays of curves into deterministic and stochastic functional components. For transformed Lorenz curves indexed by region $s$ and time $t$,

$$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$0

where $$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$1 is the functional grand effect, $$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$2 is the functional row effect, and $$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$3 is the residual function. In that formulation, the grand and row effects are deterministic, while only the residual functions are modeled dynamically (Shang, 6 Apr 2025).

These formulations share the same structural idea: FANOVA partitions functional variation into interpretable components. What changes across subfields is the mathematical object being decomposed—mean response curves, black-box mappings, stochastic processes in Hilbert spaces, or high-dimensional functional panels.

2. Orthogonality, identifiability, and Hilbert-space formulations

Under the classical Hoeffding–Sobol setting, FANOVA requires a product measure

$$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$4

and effect functions satisfy annihilating conditions such as

$$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$5

This yields orthogonality of distinct components and an additive variance decomposition,

$$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$6

When the input distribution is not uniquely specified, however, the decomposition need not be unique. The literature introduces the notion of a FANOVA core: a set of probability measures that generate identical effect functions. Without a prior over candidate measures, multiple cores imply multiple decompositions; with a prior $$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$7, the mixture decomposition is unique, but orthogonality is generally lost because the mixture usually ceases to be a product measure (Borgonovo et al., 2018).

Identifiability constraints in functional-response FANOVA are typically imposed by centering the effects. In grouped longitudinal models, one common condition is

$$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$8

while median-polish formulations use functional medians,

$$f_S(x) = \int F(x)\, dP(x_{-S}) - \sum_{T \subset S} f_T(x),$$9

Such constraints separate the grand effect from group deviations and residual structure (Ji et al., 2022, Shang, 6 Apr 2025).

A more abstract formulation places FANOVA in a real separable Hilbert space $S$0. In the Hilbert-valued fixed-effect model,

$S$1

the response and coefficients are $S$2-valued, and the error covariance is an operator matrix with compact self-adjoint trace-class blocks. Under a common orthonormal eigenvector system, the infinite-dimensional problem reduces to a sequence of finite-dimensional Gaussian problems. The geometry of the associated reproducing kernel Hilbert space defines the natural quadratic loss,

$S$3

from which generalized least squares, functional sums of squares, and functional linear-hypothesis testing are derived (Ruiz-Medina, 2015).

This Hilbert-space perspective makes explicit that FANOVA is not merely a heuristic decomposition. It can be formulated as a rigorously defined inferential framework with operator-valued covariance, exact distributional calculations, and transformed total, regression, and residual sums of squares.

3. Testing mean curves and repeated-measures effects

Inferential FANOVA distinguishes between global tests, which ask whether any group differs anywhere on the domain, and local or interval-specific procedures, which identify where differences occur. For pairwise comparison of treatment levels in functional responses, one approach partitions the domain into prespecified intervals $S$4 and uses the functional $S$5-type numerator statistic

$S$6

A two-level procedure first tests interval-wise global differences and then, within significant intervals, tests pairwise hypotheses $S$7 using permutation $S$8-values and the closure principle to control the family-wise error rate (Vsevolozhskaya et al., 2014).

Repeated-measures FANOVA adapts the classical repeated-measures ANOVA $S$9-statistic pointwise and then aggregates it. With $P$0 repeated visits or conditions, the pointwise statistic is

$P$1

and two global functionals are constructed:

$P$2

Permutation and bootstrap procedures approximate the null distribution, and pairwise post hoc comparisons are handled with Bonferroni-adjusted $P$3-values (Kuryło et al., 2023).

A different route reduces repeated-measures FANOVA to multivariate inference through basis expansion. If each curve is approximated by $P$4 basis coefficients, the two-way functional model becomes a repeated-measures MANOVA on the coefficient vectors. Two frameworks are used: the Doubly Multivariate Model (DMM), which requires $P$5, and the Mixed Multivariate Model (MMM), which requires a multivariate sphericity condition but only $P$6 (Acal et al., 2024).

Setting	Representative formulation	Notable feature
Interval-wise one-way FANOVA	$P$7	Closure-based pairwise follow-up
Functional repeated measures	$P$8, $P$9 built from $i$0	Permutation and bootstrap inference
Basis-expansion repeated measures	Curves $i$1 coefficient MANOVA	DMM or MMM depending on covariance assumptions

For multivariate functional responses, dimension reduction via multivariate FPCA is also used. In the air-pollution study, the multivariate FANOVA problem for independent measures is reduced to testing homogeneity on vectors of principal-component scores, and four principal components explained more than $i$2 of total variability in all situations (Acal et al., 2024).

A common misconception is that FANOVA inference must be either fully global or fully pointwise. The cited methods show a broader spectrum: interval-wise global tests, repeated-measures functionals, score-based multivariate reductions, and post hoc procedures that preserve multiplicity control.

4. Dependence, robustness, phase variation, and spatial geometry

Many functional datasets violate the i.i.d. assumptions of classical FANOVA. One line of work addresses dependent errors by combining wavelet estimation with a Cochrane–Orcutt-type iteration. In the dependent-error model

$i$3

the algorithm alternates between prewhitening, wavelet estimation of the transformed signal, reconstruction of $i$4, and updating $i$5. The paper reports that the iterative correction for dependence improves the $i$6 error norm at each iteration (Kist et al., 2015).

A separate wavelet-domain FANOVA approach uses a Bayesian hierarchical model with spike-and-slab priors and a Markov grove over wavelet-tree indicators. The hidden-state recursion yields an exact posterior through a pyramid algorithm with computational complexity linear in both the number of observations and the number of locations. Posterior marginal alternative probabilities and posterior joint alternative probabilities provide local and global evidence for factor contributions without MCMC (Ma et al., 2016).

Classical FANOVA also conflates amplitude and phase variation when curves differ in timing. Warped FANOVA resolves this by writing the observed trajectory as

$i$7

where $i$8 is a monotone increasing warping function and $i$9 follows a functional random-effects ANOVA decomposition. The model thereby separates amplitude variation from phase variation, accommodates irregularly sampled longitudinal data, and estimates both components jointly by maximum likelihood via an EM algorithm (Gervini et al., 2013).

Robustness is another major extension. RoFANOVA replaces least-squares-type means by functional $j$0-estimators and uses a FuNMAD-based scale estimate inside robust test statistics for main effects and interaction. Inference proceeds by permutation, and the simulation study reports that, under contamination, the redescending RoFANOVA variants outperform competitors in both size control and power. The method is implemented in the R package rofanova (Centofanti et al., 2021).

In spatially structured functional models, domain geometry itself affects FANOVA behavior. For Hilbert-valued fixed-effect models with ARH(1) errors and Dirichlet boundary conditions, rectangles, disks, and circular sectors induce different eigenstructures for the covariance kernels. The rate of decay of the corresponding eigenvalues affects dependence range, GLS stability, truncation, and the derived significance test for functional fixed effects (Álvarez-Liébana et al., 2017).

These developments show that FANOVA is not tied to a single stochastic regime. It has been reformulated for autocorrelation, localized wavelet structure, time warping, outlier contamination, and geometry-driven covariance decay.

5. Forecasting and effect decomposition for complex functional outputs

FANOVA has become a forecasting device as well as an inferential one. In forecasting a time series of Lorenz curves for Italian regions, one-way FANOVA is combined with functional median polish. Lorenz curves are first transformed by the logit,

$j$1

then decomposed as

$j$2

with the residual functions $j$3 modeled by FPCA and ARIMA score forecasts. Point forecasts are assembled as deterministic FANOVA effects plus forecast residuals, prediction intervals are constructed by bootstrap, and isotonic regression is applied when monotonicity must be enforced (Shang, 6 Apr 2025).

For high-dimensional functional time series of Japanese subnational mortality, a two-way FANOVA first separates grand mean, prefecture effect, sex effect, prefecture-by-sex interaction, and residual process,

$j$4

A sex-specific one-way FANOVA is then applied to the interaction-plus-residual term, and the remaining residual process is modeled by a functional factor model. The resulting forecast combines deterministic TWA and OWA effects with stochastic factor-model forecasts. In the reported application, TWA + OWA + FFM achieved the best average RMSFE for both females and males (Shang et al., 30 Mar 2026).

Functional-output computer experiments introduce a different challenge: the output itself is a function indexed by $j$5. The functional-output orthogonal additive Gaussian process (FOAGP) models

$j$6

and enforces conditional orthogonality through

$j$7

This yields exact functional-output effect decomposition, local Sobol’ indices

$j$8

and global expected conditional variance indices

$j$9

without requiring predefined basis functions or strong distributional assumptions (Tan et al., 15 Jun 2025).

A plausible implication is that forecasting-oriented FANOVA increasingly treats interpretable deterministic structure and stochastic residual structure as distinct modeling targets. That pattern appears in Lorenz-curve forecasting, mortality forecasting, and functional-output surrogate modeling alike.

6. Sensitivity analysis, explainability, and interpretable machine learning

In sensitivity analysis and explainable AI, FANOVA functions as a decomposition layer beneath attribution methods. A recent framework defines three decompositions according to how missing features are integrated out: baseline fANOVA, marginal fANOVA, and conditional fANOVA. Their value functions are

$y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$0

The paper’s central point is that many disagreements among feature-based explanations arise from this distributional choice rather than from purely algorithmic differences. Within the same framework, pure, partial, and full interaction summaries correspond to different game-theoretic allocations of higher-order terms (Fumagalli et al., 2024).

Low-order FANOVA models are also central to inherently interpretable machine learning. In the GAMI formulation,

$y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$1

methods such as EBM, GAMI-Net, and GAMI-Lin-T differ mainly in how they estimate main effects and pairwise interactions, how they filter candidate interactions, and whether they explicitly purify interactions to enforce hierarchical orthogonality. In GAMI-Lin-T, purification subtracts an additive fit $y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$2 from the raw interaction and adds those lower-order pieces back to the main effects (Hu et al., 2023).

For FANOVA Gaussian processes, exact Shapley attributions can be computed in quadratic time because the model is already decomposed into orthogonal component functions,

$y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$3

The local value function

$y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$4

and the global variance-based value function

$y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$5

lead to exact local stochastic Shapley values and exact global variance-based Shapley values via closed-form Möbius representations and recursive elementary-symmetric-polynomial calculations (Mohammadi et al., 20 Aug 2025).

A persistent misconception is that FANOVA-based explanations are intrinsic properties of the fitted model alone. The multiple-distribution results in sensitivity analysis show otherwise: effect functions, variance decompositions, and even uniqueness depend on the assumed input law. Under a prior over plausible measures, the variance acquires an additional variability-of-the-mean term,

$y_{ij}(t)=\mu_i(t)+\varepsilon_{ij}(t),$6

so distributional ambiguity contributes directly to the output uncertainty budget (Borgonovo et al., 2018).

Taken together, these strands indicate that FANOVA now operates across several research programs: classical testing of mean curves, Hilbert-space fixed-effect theory, repeated-measures inference, robust and warped functional modeling, forecasting of functional panels, variance-based sensitivity analysis, and interpretable machine learning. The common element is not a single estimator or test statistic, but a disciplined partition of functional variation into components whose inferential, predictive, or explanatory roles can be studied separately.