Generalized Additive Models (GAMs)

Linear regression gives us interpretability. Fully nonparametric methods give us flexibility. Generalized Additive Models (GAMs) give us both—by exploiting a deceptively simple yet profoundly powerful mathematical structure: addition.

The insight is elegant: what if we allowed each feature to have its own arbitrary, potentially complex relationship with the response, but combined these relationships through simple summation? This is the additive structure, and it forms the architectural foundation of one of the most useful modeling frameworks in modern statistical learning.

To appreciate additive models, we must first understand the fundamental problem they solve. Consider the general regression problem: given a feature vector $\mathbf{x} = (x_1, x_2, \ldots, x_p)^\top \in \mathbb{R}^p$, we want to estimate the regression function:

$$E[Y | \mathbf{X} = \mathbf{x}] = f(\mathbf{x}) = f(x_1, x_2, \ldots, x_p)$$

In the fully nonparametric approach, we make no assumptions about the form of $f$. This sounds ideal—let the data speak! But this freedom comes at a devastating cost.

Why does this happen?

Nonparametric regression works by finding observations 'near' the query point. In one dimension, 'near' is intuitive—points along a line. But in high dimensions, points scatter across a vast space. The volume of a $p$-dimensional hypercube grows as $L^p$ where $L$ is the side length. To maintain the same density of points, sample size must grow exponentially.

The statistical consequence:

With finite data in high dimensions, local neighborhoods become sparse. Either we expand neighborhood size (introducing bias by averaging dissimilar points) or accept high variance from small samples. Neither leads to accurate estimation.

This explosion in required sample size is why fully nonparametric methods fail for even moderately high-dimensional problems. We need structural assumptions to make learning tractable—but which assumptions preserve flexibility while enabling estimation?

The simplest structural assumption is linearity. Linear regression assumes:

$$f(x_1, \ldots, x_p) = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \cdots + \beta_p x_p$$

This is remarkably tractable—we have only $p+1$ parameters regardless of sample size. But the assumption is severe: each feature $x_j$ contributes to the prediction through a simple linear term $\beta_j x_j$, and features combine additively.

What linearity gives us:

What linearity takes away:

The assumption $f_j(x_j) = \beta_j x_j$ is restrictive beyond measure. Real relationships are rarely linear:

• Saturation effects: Drug efficacy plateaus at high doses
• Threshold effects: No effect below a critical value, strong effect above
• U-shaped relationships: Stress and performance follow an inverted U
• Periodic patterns: Cyclical behavior in time series

When the true relationship is nonlinear, linear regression suffers specification error—systematic bias that cannot be reduced by collecting more data.

We face a dilemma:

1. Fully nonparametric: Flexible but suffers the curse of dimensionality
2. Linear: Tractable but too restrictive for complex relationships

Is there a middle ground?

The additive model strikes a remarkable balance. Instead of forcing linear effects, we allow arbitrary univariate functions for each feature, but retain additive combination:

$$f(x_1, \ldots, x_p) = \alpha + f_1(x_1) + f_2(x_2) + \cdots + f_p(x_p)$$

where each $f_j: \mathbb{R} \to \mathbb{R}$ is a smooth but otherwise unspecified function, and $\alpha$ is an intercept term.

What the additive assumption preserves:

Mathematical formalization:

Let $(X_1, X_2, \ldots, X_p, Y)$ be random variables. An additive model assumes:

$$Y = \alpha + \sum_{j=1}^{p} f_j(X_j) + \epsilon$$

where:

• $\alpha \in \mathbb{R}$ is the overall intercept (often set to $E[Y]$)
• $f_j: \mathbb{R} \to \mathbb{R}$ are smooth functions with $E[f_j(X_j)] = 0$ for identifiability
• $\epsilon$ is random error with $E[\epsilon | \mathbf{X}] = 0$

The constraint $E[f_j(X_j)] = 0$ prevents ambiguity: without it, we could add a constant to one $f_j$ and subtract it from another, obtaining infinitely many equivalent representations.

The additive structure has a beautiful geometric interpretation. Consider the $p = 2$ case:

$$f(x_1, x_2) = \alpha + f_1(x_1) + f_2(x_2)$$

This describes a surface in 3D (the response $f$ over the $(x_1, x_2)$ plane) that can be built by combining two curves:

1. Fix $x_2$ and vary $x_1$: the profile of the surface is just $f_1(x_1)$ shifted by constant $\alpha + f_2(x_2)$
2. Fix $x_1$ and vary $x_2$: the profile is $f_2(x_2)$ shifted by constant $\alpha + f_1(x_1)$

Key geometric property: All horizontal slices (fixing $x_1$) have the same shape, just shifted vertically. All vertical slices (fixing $x_2$) have the same shape, just shifted vertically.

The 'no interaction' implication:

In an additive model, the effect of changing $x_1$ from $a$ to $b$ is:

$$f(b, x_2) - f(a, x_2) = f_1(b) - f_1(a)$$

This difference does not depend on $x_2$. The effect of moving $x_1$ from $a$ to $b$ is the same regardless of where we are in the $x_2$ dimension.

Contrast with an interaction model where $f(x_1, x_2) = x_1 \cdot x_2$. Here: $$f(b, x_2) - f(a, x_2) = (b - a) \cdot x_2$$

The effect of changing $x_1$ depends on $x_2$—this is an interaction, which additive models cannot capture.

Let's make precise how additivity circumvents the curse of dimensionality.

Problem setup:

• We have $n$ observations and $p$ features
• We want to estimate the regression function to some fixed accuracy $\delta$

Fully nonparametric approach:

To estimate a $p$-dimensional function with accuracy $\delta$ using kernel or local polynomial smoothing, the optimal rate of convergence is:

$$\text{MSE} \sim n^{-4/(4+p)}$$

For $p = 10$, achieving $\text{MSE} = 0.01$ requires $n \approx 10^{7}$ observations—often unavailable.

Additive approach:

With the additive structure, we estimate $p$ one-dimensional functions. Each $f_j$ can be estimated at the univariate rate:

$$\text{MSE}_j \sim n^{-4/5}$$

The overall rate depends on how we combine these estimates, but crucially, the rate does not degrade exponentially with $p$. For smooth additive models:

$$\text{MSE} \sim n^{-4/5}$$

This is the one-dimensional optimal rate, achieved regardless of $p$!

Intuition:

Each observation $(x_{i1}, \ldots, x_{ip}, y_i)$ provides information about all $p$ component functions simultaneously. When estimating $f_1$, we use the residual $y_i - \alpha - \sum_{j eq 1} f_j(x_{ij})$ evaluated at $x_{i1}$. All $n$ observations contribute to estimating each $f_j$.

Contrast with fully nonparametric estimation, where only observations with all features near the query point contribute—and these become exponentially rare in high dimensions.

Additive models are a direct generalization of linear regression. To see this clearly, observe that linear regression is a special case where each component function is linear:

$$f_j(x_j) = \beta_j x_j$$

The additive model generalizes by allowing: $$f_j(x_j) = \text{arbitrary smooth function}$$

Hierarchy of models:

The additive model sits between polynomial regression (fixed functional forms) and fully nonparametric regression (no assumptions). It offers adaptive flexibility—the data determines the shape of each $f_j$—while maintaining structure.

Why not always use additive models?

If additive models are more flexible than linear regression and more tractable than nonparametric regression, why ever use anything else?

1. Interpretability vs flexibility tradeoff: Linear models give a single coefficient per variable; additive models give an entire curve. Sometimes simpler is better.

2. Interactions matter: When interactions genuinely exist, additive models will be biased. Testing for interactions can guide model selection.

3. Sample size limitations: Even though additive models have favorable rates, very small samples may still struggle with smooth function estimation.

4. Domain knowledge: In some domains, we know the relationship is linear (or another parametric form). Imposing this knowledge improves efficiency.

Before proceeding to estimation, we must address a subtle but important issue: the additive decomposition is not unique without additional constraints.

The problem:

Suppose $f(x_1, x_2) = \alpha + f_1(x_1) + f_2(x_2)$. Then: $$f(x_1, x_2) = (\alpha + c) + (f_1(x_1) - c) + f_2(x_2)$$

also holds for any constant $c$. We can shift constants between the intercept and component functions arbitrarily.

Worse, if $E[f_1(X_1)] = \mu_1 eq 0$, we can redefine:

• $\tilde{\alpha} = \alpha + \mu_1$
• $\tilde{f}_1(x_1) = f_1(x_1) - \mu_1$

This gives the same function $f$ but different components.

The solution: centering constraints

We impose the constraint: $$E[f_j(X_j)] = \int f_j(x) , dP_{X_j}(x) = 0 \quad \forall j$$

where $P_{X_j}$ is the marginal distribution of $X_j$.

With this constraint:

• The intercept $\alpha = E[Y]$ (the overall mean response)
• Each $f_j$ represents the deviation from the mean attributable to feature $j$
• The decomposition is unique

Formal statement of uniqueness:

Under the centering constraint, if $f = \alpha + \sum_j f_j$ and $f = \tilde{\alpha} + \sum_j \tilde{f}_j$ with both satisfying $E[f_j] = E[\tilde{f}_j] = 0$, then:

• $\alpha = \tilde{\alpha}$
• $f_j = \tilde{f}_j$ for all $j$

The additive decomposition is unique given the centering constraint.

Implication for interpretation:

The constraint $E[f_j(X_j)] = 0$ means $f_j$ represents deviations from average for feature $j$. A positive value $f_j(x_j) > 0$ means this value of $X_j$ is associated with above-average $Y$, and vice versa. This interpretation is crisp and actionable.

Additive structures appear naturally in many domains. Understanding where additivity is reasonable—and where it breaks down—is key to successful modeling.

Additive functions have rich mathematical structure that underlies both estimation methods and theoretical guarantees.

Property 1: Projection decomposition

Let $\mathcal{H}$ be the Hilbert space of square-integrable functions of $(X_1, \ldots, X_p)$. Define: $$\mathcal{A} = \left{ g \in \mathcal{H} : g = \alpha + \sum_{j=1}^p g_j(X_j), \ E[g_j(X_j)] = 0 \right}$$

This is the set of additive functions with centered components. $\mathcal{A}$ is a closed linear subspace of $\mathcal{H}$.

Implication: For any $f \in \mathcal{H}$, there exists a unique additive projection $f^* \in \mathcal{A}$ that minimizes $E[(f - g)^2]$ over all $g \in \mathcal{A}$. This $f^*$ is the best additive approximation to $f$ in mean square sense.

Property 2: Component characterization

The components of the additive projection have a clean form. If $f^* = \alpha^* + \sum_j f_j^$, then: $$\alpha^ = E[f(\mathbf{X})]$$ $$f_j^(x_j) = E[f(\mathbf{X}) - \alpha^ | X_j = x_j] - \sum_{k eq j} E[f_k^*(X_k) | X_j = x_j]$$

For independent features ($X_j$ independent of $X_k$ for $j eq k$), this simplifies dramatically: $$f_j^*(x_j) = E[f(\mathbf{X}) | X_j = x_j] - E[f(\mathbf{X})]$$

The component $f_j^*$ is just the conditional expectation minus the overall mean—the marginal effect of $X_j$.

Property 3: ANOVA decomposition

Additive models connect to the classical ANOVA (Analysis of Variance) decomposition. For any integrable function $f$:

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

where each term is centered (integrates to zero over any of its arguments) and orthogonal to lower-order terms.

An additive model assumes all interaction terms ($f_{jk}, f_{jkl}, \ldots$) are zero. The additive structure is the zeroth-degree interaction assumption in the ANOVA hierarchy.

Property 4: Variance decomposition

Under independence and additivity: $$\text{Var}(f(\mathbf{X})) = \sum_{j=1}^p \text{Var}(f_j(X_j))$$

Total variance decomposes into additive contributions from each feature—interpretable variance attribution.

We have established the theoretical foundation of additive models, understanding both why additivity matters and what it means mathematically.

What's next:

With the additive structure understood, we turn to the component functions themselves. How do we represent and estimate each $f_j$? The next page explores component functions—the building blocks of additive models, including splines, local polynomials, and other smooth function representations that make GAM estimation practical.