49 votos

Intuición detrás de las interacciones del producto tensorial en los GAM (paquete MGCV en R)

Los modelos aditivos generalizados son aquellos en los que $$ y = \alpha + f_1(x_1) + f_2(x_2) + e_i $$ Las funciones son suaves y deben ser estimadas. Normalmente mediante splines penalizados. MGCV es un paquete en R que lo hace, y el autor (Simon Wood) escribe un libro sobre su paquete con ejemplos en R. Ruppert, et al. (2003) escriben un libro mucho más accesible sobre versiones más simples de lo mismo.

Mi pregunta es sobre las interacciones dentro de este tipo de modelos. ¿Qué pasa si quiero hacer algo como lo siguiente: $$ y = \alpha + f_1(x_1) + f_2(x_2) + f_3(x_1\times x_2) + e_i $$ si estuviéramos en la tierra de OLS (donde el $f$ es sólo una beta), no tendría ningún problema en interpretar $\hat{f}_3$ . Si estimamos a través de splines penalizados, tampoco tengo problemas de interpretación en el contexto aditivo.

Pero el paquete MGCV en el GAM tiene estas cosas llamadas "suavizaciones de productos tensoriales". Busco en Google "producto tensorial" y se me nublan los ojos inmediatamente al intentar leer las explicaciones que encuentro. O no soy lo suficientemente inteligente o las matemáticas no se explican muy bien, o ambas cosas.

En lugar de codificar

normal = gam(y~s(x1)+s(x2)+s(x1*x2))

un producto tensorial haría lo mismo (?) al

what = gam(y~te(x1,x2))

cuando lo haga

plot(what)

o

vis.gam(what)

Consigo unos resultados realmente geniales. Pero no tengo ni idea de lo que está pasando dentro de la caja negra que es te() ni cómo interpretar la mencionada salida fría. Justo la otra noche tuve una pesadilla en la que daba un seminario. Les mostraba a todos un gráfico muy chulo, me preguntaban qué significaba y yo no lo sabía.

¿Podría alguien ayudarme tanto a mí, como a la posteridad, dando un poco de mecánica e intuición sobre lo que está pasando aquí debajo del capó? Lo ideal sería decir un poco sobre la diferencia entre el caso de interacción aditiva normal y el caso tensorial.

56voto

Tori Puntos 1

Voy a (intentar) responder a esto en tres pasos: primero, vamos a identificar exactamente lo que queremos decir con un suave univariante. A continuación, describiremos un suavizado multivariante (concretamente, un suavizado de dos variables). Por último, haré mi mejor intento de describir un producto tensorial suave.

1) Univariante suave

Digamos que tenemos algunos datos de respuesta $y$ que conjeturamos es una función desconocida $f$ de una variable predictora $x$ más algún error $ε$ . El modelo sería:

$$y=f(x)+ε$$

Ahora, para ajustar este modelo, tenemos que identificar la forma funcional de $f$ . La forma de hacerlo es identificando las funciones base, que se superponen para representar la función $f$ en su totalidad. Un ejemplo muy sencillo es una regresión lineal, en la que las funciones base son simplemente $β_2x$ y $β_1$ La intercepción. Aplicando la expansión de bases, tenemos

$$y=β_1+β_2x+ε$$

En forma de matriz, tendríamos:

$$Y=Xβ+ε$$

Donde $Y$ es un vector de columnas n por 1, $X$ es una matriz modelo de n por 2, $β$ es un vector de columnas de 2 por 1 de los coeficientes del modelo, y $ε$ es un vector de columnas de errores n por 1. $X$ tiene dos columnas porque hay dos términos en nuestra expansión de base: el término lineal y el intercepto.

El mismo principio se aplica a la expansión de bases en MGCV, aunque las funciones de base son mucho más sofisticadas. En concreto, no es necesario definir las funciones de base individuales en todo el dominio de la variable independiente $x$ . Esto es lo que ocurre a menudo cuando se utilizan bases basadas en nudos (véase "ejemplo basado en nudos" ). El modelo se representa entonces como la suma de las funciones base, cada una de las cuales se evalúa en cada valor de la variable independiente. Sin embargo, como he mencionado, algunas de estas funciones de base toman un valor de cero fuera de un intervalo determinado y, por tanto, no contribuyen a la expansión de las bases fuera de ese intervalo. Como ejemplo, consideremos una base spline cúbica en la que cada función base es simétrica respecto a un valor diferente (nudo) de la variable independiente -- en otras palabras, cada función base tiene el mismo aspecto pero es sólo desplazado a lo largo del eje de la variable independiente (esto es una simplificación excesiva, ya que cualquier base práctica incluirá también un intercepto y un término lineal, pero esperamos que se entienda la idea).

Para ser explícitos, una expansión de bases de dimensión $i-2$ podría parecer:

$$y=β_1+β_2x+β_3f_1(x)+β_4f_2(x)+...+β_if_{i-2} (x)+ε$$

donde cada función $f$ es, quizás, una función cúbica de la variable independiente $x$ .

La ecuación matricial $Y=Xβ+ε$ puede seguir utilizándose para representar nuestro modelo. La única diferencia es que $X$ es ahora una matriz de n por i; es decir, tiene una columna para cada término de la expansión de bases (incluyendo el intercepto y el término lineal). Como el proceso de expansión de bases nos ha permitido representar el modelo en forma de ecuación matricial, podemos utilizar los mínimos cuadrados lineales para ajustar el modelo y encontrar los coeficientes $β$ .

Este es un ejemplo de regresión no penalizada, y uno de los principales puntos fuertes del MGCV es su estimación de suavidad mediante una matriz de penalización y un parámetro de suavización. En otras palabras, en lugar de:

$$β=(X^TX)^{-1}X^TY$$

que tenemos:

$$β=(X^TX+λS)^{-1}X^TY$$

donde $S$ es un cuadrático $i$ -por- $i$ matriz de penalización y $λ$ es un parámetro de suavización escalar. No voy a entrar en la especificación de la matriz de penalización aquí, pero debería ser suficiente decir que para cualquier expansión de base dada de alguna variable independiente y la definición de una penalización cuadrática de "ondulación" (por ejemplo, una penalización de segunda derivada), se puede calcular la matriz de penalización $S$ .

El MGCV puede utilizar varios medios para estimar el parámetro de suavizado óptimo $λ$ . No voy a entrar en ese tema ya que mi objetivo aquí era dar una amplia visión de cómo se construye un suave univariado, lo que creo que he hecho.

2) Multivariante suave

La explicación anterior puede generalizarse a múltiples dimensiones. Volvamos a nuestro modelo que da la respuesta $y$ como una función $f$ de los predictores $x$ y $z$ . La restricción a dos variables independientes evitará que la explicación se llene de notación arcaica. El modelo es entonces:

$$y=f(x,z)+ε$$

Ahora, debería ser intuitivamente obvio que vamos a representar $f(x,z)$ con una expansión de bases (es decir, una superposición de funciones de base) al igual que hicimos en el caso univariante de $f(x)$ arriba. También debería ser obvio que al menos una, y casi seguramente muchas más, de estas funciones base deben ser funciones de ambos $x$ y $z$ (si no fuera el caso, entonces implícitamente $f$ sería separable de forma que $f(x,z)=f_x(x)+f_z(z)$ ). Se puede encontrar una ilustración visual de una base spline multidimensional ici . Una expansión de base bidimensional completa de dimensión $i-3$ podría ser algo así:

$$y=β_1+β_2x+β_3z+β_4f_1(x,z)+...+β_if_{i-3} (x,z)+ε$$

Creo que está bastante claro que todavía podemos representar esto en forma de matriz con:

$$Y=Xβ+ε$$

simplemente evaluando cada función base en cada combinación única de $x$ y $z$ . La solución sigue siendo:

$$β=(X^TX)^{-1}X^TY$$

El cálculo de la matriz de penalización de la segunda derivada es muy parecido al del caso univariante, salvo que en lugar de integrar la segunda derivada de cada función base con respecto a una sola variable, integramos la suma de todas las segundas derivadas (incluidas las parciales) con respecto a todas las variables independientes. Los detalles de lo anterior no son especialmente importantes: la cuestión es que podemos seguir construyendo la matriz de penalización $S$ y utilizar el mismo método para obtener el valor óptimo del parámetro de suavizado $λ$ y dado ese parámetro de suavizado, el vector de coeficientes sigue siendo

$$β=(X^TX+λS)^{-1}X^TY$$

Ahora, este liso bidimensional tiene un isotrópico pena: esto significa que un solo valor de $λ$ se aplica en ambas direcciones. Esto funciona bien cuando ambos $x$ y $z$ tienen aproximadamente la misma escala, como una aplicación espacial. Pero, ¿y si sustituimos la variable espacial $z$ con la variable temporal $t$ ? Las unidades de $t$ puede ser mucho mayor o menor que las unidades de $x$ y esto puede desviar la integración de nuestras segundas derivadas porque algunas de esas derivadas contribuirán de forma desproporcionada a la integración global (por ejemplo, si medimos $t$ en nanosegundos y $x$ en años luz, la integral de la segunda derivada con respecto a $t$ puede ser ampliamente mayor que la integral de la segunda derivada con respecto a $x$ y, por lo tanto, el "meneo" a lo largo de la $x$ dirección puede quedar en gran medida impune). La diapositiva 15 de la "caja de herramientas suave" que he enlazado tiene más detalles sobre este tema.

Cabe destacar que no hemos descompuesto las funciones de base en bases marginales de $x$ y $z$ . La implicación aquí es que los alisados multivariantes deben construirse a partir de bases que soporten múltiples variables. Los alisados de producto tensorial permiten construir bases multivariantes a partir de bases marginales univariantes, como explico a continuación.

3) El producto tensorial suaviza

Los suavizados del producto tensorial abordan el problema de modelar las respuestas a las interacciones de múltiples entradas con diferentes unidades. Supongamos que tenemos una respuesta $y$ que es una función $f$ de la variable espacial $x$ y la variable temporal $t$ . Nuestro modelo es entonces:

$$y=f(x,t)+ε$$

Lo que queremos hacer es construir una base bidimensional para las variables $x$ y $t$ . Esto será mucho más fácil si podemos representar $f$ como:

$$f(x,t)=f_x(x)f_t(t)$$

En un sentido algebraico/analítico, esto no es necesariamente posible. Pero recordemos que estamos discretizando los dominios de $x$ y $t$ (imagina un "entramado" bidimensional definido por las ubicaciones de los nudos en el $x$ y $t$ ejes) tal que la función "verdadera" $f$ está representada por la superposición de funciones de base. Al igual que asumimos que una función univariante muy compleja puede ser aproximada por una función cúbica simple en un intervalo específico de su dominio, podemos asumir que la función no separable $f(x,t)$ puede ser aproximado por el producto de funciones más simples $f_x(x)$ y $f_t(t)$ en un intervalo, siempre que nuestra elección de las dimensiones de la base haga que esos intervalos sean suficientemente pequeños.

Nuestra expansión de base, dada una $i$ -base dimensional en $x$ y $j$ -base dimensional en $t$ se vería así:

\begin{align} y = &β_{1} + β_{2}x + β_{3}f_{x1}(x)+β_{4}f_{x2}(x)+...+ \\ &β_{i}f_{x(i-3)}(x)+ β_{i+1}t + β_{i+2}tx + β_{i+3}tf_{x1}(x)+β_{i+4}tf_{x2}(x)+...+ \\ &β_{2i}tf_{x(i-3)}(x)+ β_{2i+1}f_{t1}(t) + β_{2i+2}f_{t1}(t)x + β_{2i+3}f_{t1}(t)f_{x1}(x)+β_{i+4}f_{t1}(t)f_{x2}(x){\small +...+} \\ &β_{2i}f_{t1}(t)f_{x(i-3)}(x)+\ldots+ \\ &β_{ij}f_{t(j-3)}(t)f_{x(i-3)}(x) + ε \end{align}

Que puede interpretarse como un producto tensorial. Imaginemos que evaluamos cada función base en $x$ y $t$ construyendo así las matrices modelo n-por-i y n-por-j $X$ y $T$ respectivamente. A continuación, podemos calcular el $n^2$ -por- $ij$ producto tensorial $X \otimes T$ de estas dos matrices modelo y reorganizarlas en columnas, de manera que cada columna representara una combinación única $ij$ . Recordemos que las matrices marginales del modelo tenían $i$ y $j$ columnas, respectivamente. Estos valores corresponden a las dimensiones de sus bases respectivas. Nuestra nueva base de dos variables debería entonces tener la dimensión $ij$ y, por tanto, el mismo número de columnas en su matriz modelo.

NOTA: Me gustaría señalar que, dado que construimos explícitamente las funciones base del producto tensorial tomando productos de funciones base marginales, las bases del producto tensorial pueden construirse a partir de bases marginales de cualquier tipo. No es necesario que admitan más de una variable, a diferencia de la suavidad multivariante comentada anteriormente.

En realidad, este proceso da lugar a una expansión global de la base de la dimensión $ij-i-j+1$ porque la multiplicación completa incluye la multiplicación de cada $t$ función base por la intersección x $β_{x1}$ (por lo que restamos $j$ ), además de multiplicar cada $x$ función base por la intersección t $β_{t1}$ (por lo que restamos $i$ ), pero debemos sumar el intercepto por sí mismo (por lo que sumamos 1). Esto se conoce como aplicar una restricción de identificabilidad.

Así que podemos representar esto como:

$$y=β_1+β_2x+β_3t+β_4f_1(x,t)+β_5f_2(x,t)+...+β_{ij-i-j+1}f_{ij-i-j-2}(x,t)+ε$$

Donde cada una de las funciones de base multivariante $f$ es el producto de un par de marginales $x$ y $t$ funciones de base. De nuevo, está bastante claro que habiendo construido esta base podemos representar esto con la ecuación matricial:

$$Y=Xβ+ε$$

Que (todavía) tiene la solución:

$$β=(X^TX)^{-1}X^TY$$

Donde la matriz del modelo $X$ tiene $ij-i-j+1$ columnas. En cuanto a las matrices de penalización $J_x$ y $J_t$ se construyen por separado para cada variable independiente de la siguiente manera:

$$J_x=β^T I_j \otimes S_x β$$

y,

$$J_t=β^T S_t \otimes I_i β$$

Esto permite una penalización global anisotrópica (diferente en cada dirección) (Nota: las penalizaciones en la segunda derivada de $x$ se suman en cada nudo de la $t$ y viceversa). Los parámetros de suavizado $λ_x$ y $λ_t$ puede ahora estimarse de forma muy parecida a como se hizo con el parámetro único de suavizado para los suavizados univariantes y multivariantes. El resultado es que la forma general de un producto tensorial suave es invariable al reescalado de sus variables independientes.

Recomiendo la lectura de todas las viñetas en el sitio web de MGCV, así como de " Modelos aditivos generalizados: e introducción con R ." Larga vida a Simon Wood.

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X