¿Cómo incluir un término lineal y cuadrático cuando también se incluye la interacción con esas variables?

Question

Cuando se añade un predictor numérico con predictores categóricos y sus interacciones, se suele considerar necesario centrar previamente las variables en 0. El razonamiento es que, de lo contrario, los efectos principales son difíciles de interpretar, ya que se evalúan con el predictor numérico en 0.

Mi pregunta ahora es cómo centrar si uno no sólo incluye la variable numérica original (como término lineal) sino también el término cuadrático de esta variable. En este caso, son necesarios dos enfoques diferentes:

1. Centrando ambas variables en su media individual. Esto tiene el desafortunado inconveniente de que el 0 ahora está en una posición diferente para ambas variables considerando la variable original.
2. Centrar ambas variables en la media de la variable original (es decir, restando la media de la variable original para el término lineal y restando el cuadrado de la media de la variable original para el término cuadrático). Con este enfoque, el 0 representaría el mismo valor de la variable original, pero la variable cuadrática no estaría centrada en 0 (es decir, la media de la variable no sería 0).

Creo que el enfoque 2 parece razonable dado el motivo del centrado después de todo. Sin embargo, no encuentro nada al respecto (tampoco en las preguntas relacionadas: a y b ).

¿O es generalmente una mala idea incluir términos lineales y cuadráticos y sus interacciones con otras variables en un modelo?

Answer 1

Cuando se incluyen polinomios e interacciones entre ellos, la multicolinealidad puede ser un gran problema; un enfoque es buscar polinomios ortogonales.

En general, los polinomios ortogonales son una familia de polinomios que son ortogonales con respecto a algún producto interior.

Así, por ejemplo, en el caso de polinomios sobre alguna región con función de peso $w$ El producto interno producto interno es $\int_a^bw(x)p_m(x)p_n(x)dx$ - la ortogonalidad hace que el producto interno $0$ a menos que $m=n$ .

El ejemplo más sencillo de polinomios continuos son los polinomios de Legendre, que tienen función de peso constante sobre un intervalo real finito (comúnmente sobre $[-1,1]$ ).

En nuestro caso, el espacio (las propias observaciones) es discreto, y nuestra función de peso también es constante (normalmente), por lo que los polinomios ortogonales son una especie de equivalente discreto de los polinomios de Legendre. Con la constante incluida en nuestros predictores, el producto interno es simplemente $p_m(x)^Tp_n(x) = \sum_i p_m(x_i)p_n(x_i)$ .

Por ejemplo, considere $x = 1,2,3,4,5$

Comienza con la columna constante, $p_0(x) = x^0 = 1$ . El siguiente polinomio es de la forma $ax-b$ pero no nos preocupamos por la escala en este momento, así que $p_1(x) = x-\bar x = x-3$ . El siguiente polinomio sería de la forma $ax^2+bx+c$ ; resulta que $p_2(x)=(x-3)^2-2 = x^2-6x+7$ es ortogonal a las dos anteriores:

x         p0  p1  p2   
1          1  -2   2   
2          1  -1  -1
3          1   0  -2
4          1   1  -1
5          1   2   2

Con frecuencia, la base también se normaliza (produciendo una familia ortonormal), es decir, la suma de los cuadrados de cada término se fija en alguna constante (por ejemplo, en $n$ o a $n-1$ para que la desviación estándar sea 1, o quizás lo más frecuente, a $1$ ).

Las formas de ortogonalizar un conjunto de predictores polinómicos incluyen la ortogonalización de Gram-Schmidt y la descomposición de Cholesky, aunque hay muchos otros enfoques.

---

Algunas de las ventajas de los polinomios ortogonales:

1) la multicolinealidad no es un problema: estos predictores son todos ortogonales.

2) Los coeficientes de bajo orden no cambian al añadir términos . Si se ajusta a un grado $k$ mediante polinomios ortogonales, se conocen los coeficientes de un ajuste de todos los polinomios de orden inferior sin necesidad de reajustar.

---

Ejemplo en R ( cars datos, distancias de frenado frente a la velocidad):

Aquí consideramos la posibilidad de que un modelo cuadrático sea adecuado:

R utiliza el poly para establecer predictores polinómicos ortogonales:

> p <- model.matrix(dist~poly(speed,2),cars)
> cbind(head(cars),head(p))
  speed dist (Intercept) poly(speed, 2)1 poly(speed, 2)2
1     4    2           1      -0.3079956      0.41625480
2     4   10           1      -0.3079956      0.41625480
3     7    4           1      -0.2269442      0.16583013
4     7   22           1      -0.2269442      0.16583013
5     8   16           1      -0.1999270      0.09974267
6     9   10           1      -0.1729098      0.04234892

Son ortogonales:

> round(crossprod(p),9)
                (Intercept) poly(speed, 2)1 poly(speed, 2)2
(Intercept)              50               0               0
poly(speed, 2)1           0               1               0
poly(speed, 2)2           0               0               1

Aquí hay un gráfico de los polinomios:

Este es el resultado del modelo lineal:

> summary(carsp)

Call:
lm(formula = dist ~ poly(speed, 2), data = cars)

Residuals:
    Min      1Q  Median      3Q     Max 
-28.720  -9.184  -3.188   4.628  45.152 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)       42.980      2.146  20.026  < 2e-16 ***
poly(speed, 2)1  145.552     15.176   9.591 1.21e-12 ***
poly(speed, 2)2   22.996     15.176   1.515    0.136    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 15.18 on 47 degrees of freedom
Multiple R-squared:  0.6673,    Adjusted R-squared:  0.6532 
F-statistic: 47.14 on 2 and 47 DF,  p-value: 5.852e-12

Aquí hay un gráfico del ajuste cuadrático:

Answer 2

No creo que merezca la pena centrar, y el centrado hace más compleja la interpretación de las estimaciones de los parámetros. Si se utiliza un software moderno de álgebra matricial, la colinealidad algebraica no es un problema. Su motivación original de centrar para poder interpretar los efectos principales en presencia de la interacción no es fuerte. Los efectos principales cuando se estiman en cualquier valor elegido automáticamente de un factor continuo que interactúa son algo arbitrarios, y es mejor pensar en esto como un simple problema de estimación comparando los valores predichos. En el R rms paquete contrast.rms por ejemplo, puede obtener cualquier contraste de interés independientemente de la codificación de las variables. He aquí un ejemplo de una variable categórica x1 con niveles "a" "b" "c" y una variable continua x2, ajustada mediante un spline cúbico restringido con 4 nudos por defecto. Se permiten diferentes relaciones entre x2 e y para diferentes x1. Dos de los niveles de x1 se comparan con x2=10.

require(rms)
dd <- datadist(x1, x2); options(datadist='dd')
f <- ols(y ~ x1 * rcs(x2,4))
contrast(f, list(x1='b', x2=10), list(x1='c', x2=10))
# Now get all comparisons with c:
contrast(f, list(x1=c('a','b'), x2=10), list(x1='c', x2=10))
# add type ='joint' to get a 2 d.f. test, or conf.type='simultaneous'
# to get simultaneous individual confidence intervals

Con este enfoque también se pueden estimar fácilmente los contrastes en varios valores del factor o factores que interactúan, por ejemplo

contrast(f, list(x1='b', x2=10:20), list(x1='c', x2=10:20))