Referencias
Las notas se han elaborado a partir de vario material:
- Joaquín Amat Rodrigo, Métodos de regresión no lineal: regresión polinómica, regression splines, smooth splines y GAMs1.
- G. James et al., An Introduction to Statistical Learning, 2023, cap.7 (Disponible en statlearning.com).2
- Varias aportaciones personales.
Base teórica
Definición
La regresión spline es una técnica no paramétrica utilizada para modelar relaciones no lineales entre una variable dependiente (\(Y\)) y una independiente (\(Y\)). A diferencia de la regresión polinómica tradicional que ajusta una única curva a todo el conjunto de datos, la regresión spline divide el rango de \(X\) en segmentos más pequeños y ajusta polinomios de bajo grado (generalmente cúbicos) en cada segmento, unidos suavemente en puntos específicos llamados nudos (knots).
Piecewise polinomial
Modelo de regresión
Se trata de una combinación directa de la regresión polinómica y de step functions. En lugar de ajustar un polinomio de grado alto a todo el rango del predictor X, este se divide en subintervalos y en cada uno de ellos se ajusta un polinomio de menor grado. Un piecewise polinomial de grado 3 emplea la función:
\[ y_i = \beta_0 + \beta_1x_i + \beta_2x^2_i + \beta_3x^3_i + \epsilon_i \]
donde los coeficientes \(β_0, β_1, β_2, β_3\) toman diferente valor en cada una de las regiones en las que se divide el rango de \(X\). A los puntos de corte que delimitan cada región o subintervalo se les denomina knots (\(c_1,..., c_k\)).
Si el modelo solo tiene un único punto de corte (knot) en el valor \(c\), las funciones que describen el modelo son:
\[ y_i = \begin{cases} \beta_{01} + \beta_{11} x_i + \beta_{21} x^2_i + \beta_{31}x^3_i + \epsilon_i & \text{ si } x_i<c \\ \beta_{02} + \beta_{12} x_i + \beta_{22} x^2_i + \beta_{32}x^3_i + \epsilon_i & \text{ si } x_i \geq c \end{cases} \]
Imagen de un Piecewise Polinomial de grado 3 obtenida del libro ISLR
Cubic spline
Para evitar la discontinuidad y exceso de flexibilidad de los piecewise polinomial se pueden imponer restricciones a los polinomios de cada región de forma que el modelo final sea una curva continua. La restricción más básica consiste en forzar a cada polinomio a pasar por los puntos de corte que lo flanquean, de esta forma, el polinomio de la región \(i\) termina en el mismo punto donde empieza el polinomio de la región \(i+1\).
Para conseguir que la transición de una región a otra, además de continua, sea suave, se imponen restriccioness: que las \(d\)-derivadas de los polinomios sean continuas en los puntos de corte (que pasen por ellos), siendo d el grado del polinomio menos 1. En el caso de un piecewise polinomial de grado 3, la máxima continuidad de obtiene aplicando las siguientes restricciones:
-Los polinomios deben ser continuos. Cada polinomio debe pasar por los puntos de corte que lo delimitan.
-La primera derivada de los polinomios debe de ser continua. La primera derivada de cada polinomio debe pasar por los puntos de corte que lo delimitan.
-La segunda derivada de los polinomios debe de ser continua. La segunda derivada de cada polinomio debe pasar por los puntos de corte que lo delimitan3.
A la curva final obtenida al imponer las restricciones de continuidad sobre un piecewise polinomial de grado 3 se le conoce como cubic spline o spline de grado 3. Los cubic spline se emplean con frecuencia ya que a partir de la segunda derivada el ojo humano no aprecia la discontinuidad en los knots. A modo general, un spline de grado \(d\) se define como un piecewise polinomial de grado \(d\) con la condición de continuidad en las derivadas \(\text{d-1}\) en cada punto de corte (knots).
Para poder obtener el ajuste de un cubic spline con \(k\) knots, se emplea una regresión por mínimos cuadrados sobre una ecuación formada por la intersección \(y\) \(3 + k\) predictores:
\[ \beta_0, \ X, \ X^2, \ X^3, \ h(X,\xi_1), \ h(X,\xi_2),..., h(X,\xi_K) \]
donde \(\xi_1,...,\xi_k\) son los puntos de corte (knots).
El modelo resultante implica estimar un total de \(4 + k\) coeficientes de regresión a estimar, por lo que el ajuste de un cubic spline con \(k\) knots tiene \(4 + k\) coeficiente a estimar.
Natural spline
Las regression splines pueden tener mucha varianza en los extremos superior e inferior del predictor, lo que genera intervalos de confianza muy amplios. Esto es así debido a que la primera y última región carecen de restricción de continuidad en uno de sus extremos, por lo que tienen un exceso de flexibilidad. **Los *natural splines evitan este problema incorporando una restricción extra a las regression splines: la función tiene que ser lineal en los extremos. Esto significa que las regiones para las que el predictor \(X\) es menor que el menor de los knots o mayor que el mayor de los knots siguen un ajuste lineal.
Se ajustan dos funciones polinómicas por mínimos cuadrados. Para la primera se emplean las observaciones en las que \(x_i<c\) obteniendo los coeficientes \(β_{01}, β_{11}, β_{21}, β_{31}\) y para la segunda aquellas en las que \(x_i≥c\), obteniendo los coeficientes \(β_{02}, β_{12}, β_{22}, β_{32}\).
Ejemplo Regression Cubic Spline
El set de datos \(\tt Wage\) del paquete \(\tt ISRL\) contiene información sobre 3000 trabajadores. Entre las 12 variables registradas se encuentra el salario (\(\tt wage\)) y la edad (\(\tt age\)). Dada la relación no lineal existente entre estas dos variables, se recurre , se pretende crear un modelo de regression splines que permita predecir el salario de un trabajador en función de la edad que tiene.
Para ajustar regression splines en R se emplea la librería \(\tt splines\), que viene instalada por defecto. La función \(\tt bs()\), cuyo nombre viene de \(\text{B-spline}\), genera la matriz de ecuaciones necesaria para crear las splines acorde a los knots indicados o de los \(df\). Por defecto, \(\tt bs()\) se generan polinomios de grado 3 (cubic splines) pero esto puede cambiarse con el argumento \(\tt degree\).
Regresión decidiendo el grado de libertad \(df\)
Para que los knots estén equidistribuidos, en lugar de indicar las posiciones en el argumento \(\tt knots\), se indican los grados de libertad del spline que, en este caso son los coeficientes a estimar menos el intercept.
Ajuste del modelo
library(ISLR)
data('Wage')
library(splines)
# CUBIC SPLINE CON 3 KNOTS EQUIDISTRIBUIDOS
# -----------------------------------------------------------------------------
modelo_splines <- lm(wage ~ bs(age, df = 6, degree = 3), data = Wage)
modelo_splines##
## Call:
## lm(formula = wage ~ bs(age, df = 6, degree = 3), data = Wage)
##
## Coefficients:
## (Intercept) bs(age, df = 6, degree = 3)1
## 56.31 27.82
## bs(age, df = 6, degree = 3)2 bs(age, df = 6, degree = 3)3
## 54.06 65.83
## bs(age, df = 6, degree = 3)4 bs(age, df = 6, degree = 3)5
## 55.81 72.13
## bs(age, df = 6, degree = 3)6
## 14.75Como se ve, con \(\tt df=6\), el modelo nos da 7 coeficientes
Interpolación de puntos dentro del rango del predictor
Se crea un set de edades continuo, sin repeticiones, sobre el que ajustar el modelo.
Predicción de la variable respuesta y del error estándar
En este caso, además de querer visualizar la función, queremos también calcular el error estándar que tiene ciertas peculiaridades.
Cálculo del intervalo de confianza superior e inferior
En este caso podemos visualizar la curva de la función y los intervalos de confianza
# -----------------------------------------------------------------------------
intervalo_conf <- data.frame(
inferior = predicciones$fit - 1.96*predicciones$se.fit,
superior = predicciones$fit + 1.96*predicciones$se.fit)
attach(Wage)
plot(x = age, y = wage, pch = 20, col = "darkgrey")
title("Cubic Spline: df = 6, 3 knots")
lines(x = nuevos_puntos$age, predicciones$fit, col = "red", lwd = 2)
lines(x = nuevos_puntos$age, intervalo_conf$inferior, col = "blue",
lwd = 2, lty = 2)
lines(x = nuevos_puntos$age, intervalo_conf$superior, col = "blue",
lwd = 2, lty = 2)
Ejemplo Regression Natural Spline
Para ajustar natural splines en R, el proceso es equivalente al seguido en el ajuste de regression splines pero empleando la función ns(). Dado que las natural spline introducen dos restricciones adicionales (linealidad en los extremos), se necesitan dos grados de libertad menos en comparación con las regression splines para introducir el mismo número de knots.
Ajuste del modelo
# AUSTE DEL MODELO
# ------------------------------------------------------------------------------
library(splines)
attach(Wage)
modelo_nsplines <- lm(wage ~ ns(age, df = 4), data = Wage)
modelo_nsplines##
## Call:
## lm(formula = wage ~ ns(age, df = 4), data = Wage)
##
## Coefficients:
## (Intercept) ns(age, df = 4)1 ns(age, df = 4)2 ns(age, df = 4)3
## 58.556 60.462 41.963 97.020
## ns(age, df = 4)4
## 9.773Como vemos los cortes permanecen en los mismos puntos:
## [1] 33.75 42.00 51.00Predicción de la variable respuesta y del error estándar
Cálculo del intervalo de confianza superior e inferior
Representación gráfica y comparación con Cubic Spline
par(mfrow= c(1,2))
plot(x = age, y = wage, pch = 20, col = "darkgrey")
title("Natural Spline con 3 knots (df=4)")
lines(x = nuevos_puntos$age, prediccionesns$fit, col = "red", lwd = 2)
lines(x = nuevos_puntos$age, intervalo_confns$inferior, col = "blue",
lwd = 2, lty = 2)
lines(x = nuevos_puntos$age, intervalo_confns$superior, col = "blue",
lwd = 2, lty = 2)
plot(x = age, y = wage, pch = 20, col = "darkgrey")
title("Cubic Spline: df = 6, 3 knots")
lines(x = nuevos_puntos$age, predicciones$fit, col = "red", lwd = 2)
lines(x = nuevos_puntos$age, intervalo_conf$inferior, col = "blue",
lwd = 2, lty = 2)
lines(x = nuevos_puntos$age, intervalo_conf$superior, col = "blue",
lwd = 2, lty = 2)
La diferencia entre ambos modelos se identifica entre los valores extremos de la variable independiente, con unos niveles mayores del error estándar.
Análiside residuals
par(mfrow= c(1,2))
plot(resid(modelo_nsplines), main = "Residuals Natural Splines", ylab = "Residuals",pch=16,cex=0.75)
abline(h = 0, col = "red",lwd=2)
plot(resid(modelo_splines), main = "Residuals Cubic Splines", ylab = "Residuals",pch=16,cex=0.75)
abline(h = 0, col = "red",lwd=2)
A nivel visual, no se identifican diferencias entre ambos modelos.
Joaquín Amat Rodrigo, Métodos de regresión no lineal: regresión polinómica, regression splines, smooth splines y GAMs. Sobre este tema, ver tambien Michael Brenndoerfer↩︎
De hecho, tanto J. Amat como M. Brenndoerfer, parecen haber utilizado este texto para elaborar sus notas.↩︎
Conceptos Clave de Primera y Segunda Derivada Continua
- Primera Derivada \((f'(x))\): Analiza la pendiente de la recta tangente. Si \(f'(x)>0\), la función crece; si \(f'(x)<0\), decrece.
- Segunda Derivada \((f''(x))\): Analiza la curvatura. Si \((f''(x))>0\), la función es cóncava hacia arriba (convexo); y si \((f''(x))<0\), es cóncava hacia abajo.
- Continuidad de Derivadas: Se dice que una función es de clase \(C_2\) si su primera y segunda derivada son continuas, lo que implica una curva suave sin esquinas ni cambios abruptos en la pendiente.
Relación con Puntos Críticos:
- Máximos/Mínimos: Si \(f'(c)=0\) y si \(f''(x)<0\), hay un máximo local. Si \(f'(c)=0\) \(f''(c)>0\) y, hay un mínimo local.
Puntos de Inflexión: Si \(f''(c)=0\)y la segunda derivada cambia de signo, existe un punto de inflexión donde la curva cambia de concavidad.↩︎