Estadística Práctica para Científicos de Datos, 2ª Edición

Capítulo 4. Regresión y predicción

Este trabajo se ha traducido utilizando IA. Agradecemos tus opiniones y comentarios: [e-mail]

Quizá el objetivo más común en estadística sea responder a la pregunta "¿Es la variable X (o más probablemente, $X_1,...,X_p$) asociada a una variable Y y, en caso afirmativo, ¿cuál es la relación y podemos utilizarla para predecir Y?".

En ninguna parte es más fuerte el nexo entre la estadística y la ciencia de datos que en el ámbito de la predicción: concretamente, la predicción de una variable de resultado (objetivo) basada en los valores de otras variables "predictoras". Este proceso de entrenamiento de un modelo sobre datos en los que se conoce el resultado, para su posterior aplicación a datos en los que no se conoce el resultado, se denomina aprendizaje supervisado. Otra conexión importante entre la ciencia de los datos y la estadística se da en el ámbito de la detección de anomalías, en el que los diagnósticos de regresión destinados originalmente al análisis de datos y la mejora del modelo de regresión pueden utilizarse para detectar registros inusuales.

Regresión lineal simple

La regresión lineal simple proporciona un modelo de la relación entre la magnitud de una variable y la de una segunda; por ejemplo, a medida que X aumenta, Y también aumenta.O a medida que X aumenta, Y disminuye.¹ La correlación es otra forma de medir cómo se relacionan dos variables -ver el apartado "Correlación"-. La diferencia es que mientras la correlación mide la fuerza de una asociación entre dos variables, la regresión cuantifica la naturaleza de la relación.

Términos clave de la regresión lineal simple

Respuesta

La variable que intentamos predecir.

Sinónimos

variable dependiente, variable Y, objetivo, resultado

Variable independiente

La variable utilizada para predecir la respuesta.

Sinónimos

X variable, característica, atributo, predictor

Grabar

El vector de valores de predictor y resultado para un individuo o caso concreto.

Sinónimos

fila, caso, instancia, ejemplo

Intercepta

El intercepto de la recta de regresión, es decir, el valor previsto cuando $X=0$.

Sinónimos

$b_0$, ${\beta }_0$

Coeficiente de regresión

La pendiente de la recta de regresión.

Sinónimos

pendiente, $b_1$, ${\beta }_1$estimaciones de los parámetros, ponderaciones

Valores ajustados

Las estimaciones ${\widehat{Y}}_i$ obtenidas de la recta de regresión.

Sinónimo

valores previstos

Residuos

La diferencia entre los valores observados y los valores ajustados.

Sinónimo

errores

Mínimos cuadrados

El método de ajustar una regresión minimizando la suma de los residuos al cuadrado.

Sinónimos

mínimos cuadrados ordinarios, MCO

La ecuación de regresión

La regresión lineal simple estima cuánto cambiará Y cuando X cambie en una cantidad determinada.Con el coeficiente de correlación, las variables X e Y son intercambiables. Con la regresión, intentamos predecir la variable Y a partir de X mediante una relación lineal (es decir, una línea):

$$Y=b_0+b_{1}X$$

Leemos esto como "Y es igual a _b1 veces X, más una constante _b0". El símbolo $b_0$ se conoce como intercepto (o constante), y el símbolo $b_1$ como la pendiente de X. Ambos aparecen en la salida de R como coeficientes, aunque en el uso general el término coeficiente suele reservarse para $b_1$. La variable Y se conoce como respuesta o variable dependiente, ya que depende de X.La variable X se conoce como predictor o variable independiente. La comunidad del aprendizaje automático tiende a utilizar otros términos, llamando a Y objetivo y a X vector de características. A lo largo de este libro, utilizaremos los términos predictor y característica indistintamente.

Considera el diagrama de dispersión de la Figura 4-1, que muestra el número de años que un trabajador estuvo expuesto al polvo de algodón (Exposure) frente a una medida de la capacidad pulmonar (PEFR o "tasa de flujo espiratorio máximo"). ¿Qué relación hay entre PEFR y Exposure? Es difícil saberlo basándose sólo en la imagen.

Figura 4-1. Exposición al algodón frente a capacidad pulmonar

La regresión lineal simple intenta encontrar la "mejor" línea para predecir la respuesta PEFR en función de la variable predictora Exposure:

$$\text{PEFR}=b_0+b_1\text{Exposición}$$

La función lm de R puede utilizarse para ajustar una regresión lineal:

model <- lm(PEFR ~ Exposure, data=lung)

lm significa modelo lineal, y el símbolo ~ denota que PEFR es predicho por Exposure. Con esta definición del modelo, el intercepto se incluye y ajusta automáticamente. Si quieres excluir el intercepto del modelo, tienes que escribir la definición del modelo de la siguiente manera:

PEFR ~ Exposure - 1

Al imprimir el objeto model se obtiene el siguiente resultado:

Call:
lm(formula = PEFR ~ Exposure, data = lung)

Coefficients:
(Intercept)     Exposure
    424.583       -4.185

La intercepción, o $b_0$es 424,583 y puede interpretarse como la predicción PEFR para un trabajador con cero años de exposición. El coeficiente de regresión , o $b_1$, puede interpretarse del siguiente modo: por cada año adicional que un trabajador está expuesto al polvo de algodón, la medida PEFR del trabajador se reduce en -4,185.

En Python, podemos utilizar LinearRegression del paquete scikit-learn. (el paquete statsmodels tiene una implementación de regresión lineal más parecida a la de R (sm.OLS); la utilizaremos más adelante en este capítulo):

predictors = ['Exposure']
outcome = 'PEFR'

model = LinearRegression()
model.fit(lung[predictors], lung[outcome])

print(f'Intercept: {model.intercept_:.3f}')
print(f'Coefficient Exposure: {model.coef_[0]:.3f}')

La línea de regresión de este modelo se muestra en la Figura 4-2.

Figura 4-2. Pendiente e intersección de la regresión ajustada a los datos pulmonares

Valores ajustados y residuos

Conceptos importantes en el análisis de regresión son los valores ajustados (las predicciones) y los residuos (errores de predicción).En general, los datos no caen exactamente sobre una recta, por lo que la ecuación de regresión debe incluir un término de error explícito $e_i$:

$$Y_i=b_0+b_1{X}_i+e_i$$

Los valores ajustados, también denominados valores predichos en , se suelen denotar como ${\widehat{Y}}_i$ (y vienen dados por:

$${\widehat{Y}}_i={\widehat{b}}_0+{\widehat{b}}_1{X}_i$$

La notación ${\widehat{b}}_0$ y ${\widehat{b}}_1$ indica que los coeficientes son estimados frente a conocidos.

Notación del sombrero: Estimaciones frente a valores conocidos

La notación "sombrero" se utiliza para diferenciar entre estimaciones y valores conocidos.Así que el símbolo $\widehat{b}$ ("sombrero b") es una estimación del parámetro desconocido $b$¿Por qué los estadísticos distinguen entre la estimación y el valor real? La estimación tiene incertidumbre, mientras que el valor real es fijo.²

Calculamos los residuos ${\widehat{e}}_i$ restando los valores predichos de los datos originales:

$${\widehat{e}}_i=Y_i-{\widehat{Y}}_i$$

En R, podemos obtener los valores ajustados y los residuos utilizando las funciones predict y residuals:

fitted <- predict(model)
resid <- residuals(model)

Con el modelo LinearRegression de scikit-learn, utilizamos el método predict en los datos de entrenamiento para obtener los valores de fitted y, posteriormente, los de residuals. Como veremos, éste es un patrón general que siguen todos los modelos de scikit-learn:

fitted = model.predict(lung[predictors])
residuals = lung[outcome] - fitted

La figura 4-3 ilustra los residuos de la recta de regresión ajustada a los datos pulmonares. Los residuos son la longitud de las líneas discontinuas verticales desde los datos hasta la recta.

Figura 4-3. Residuos de una recta de regresión (para acomodar todos los datos, la escala del eje y difiere de la Figura 4-2, de ahí la pendiente aparentemente diferente)

Mínimos cuadrados

¿Cómo se ajusta el modelo a los datos? Cuando existe una relación clara, podrías imaginar el ajuste de la recta a mano. En la práctica, la recta de regresión es la estimación que minimiza la suma de valores residuales al cuadrado, también llamada suma residual de cuadrados o RSS:

$$\begin{array}{ccc}\hfill RSS& =& \sum _{i=1}^n{\left(Y_i-{\widehat{Y}}_i\right)}^2\hfill \\ & =& \sum _{i=1}^n{\left(Y_i-{\widehat{b}}_0-{\widehat{b}}_1{X}_i\right)}^2\hfill \end{array}$$

Las estimaciones ${\widehat{b}}_0$ y ${\widehat{b}}_1$ son los valores que minimizan la RSS.

El método de minimizar la suma de los residuos al cuadrado se denomina regresión por mínimos cuadrados, o regresión por mínimos cuadrados ordinarios (MCO ).A menudo se atribuye a Carl Friedrich Gauss, matemático alemán, pero fue publicado por primera vez por el matemático francés Adrien-Marie Legendre en 1805.La regresión por mínimos cuadrados puede calcularse rápida y fácilmente con cualquier software estadístico estándar.

Históricamente, la comodidad informática es una de las razones del uso generalizado de los mínimos cuadrados en la regresión. Con la llegada de los grandes datos, la velocidad informática sigue siendo un factor importante. Los mínimos cuadrados, al igual que la media (véase "Mediana y estimaciones robustas"), son sensibles a los valores atípicos, aunque esto tiende a ser un problema significativo sólo en conjuntos de datos de tamaño pequeño o moderado. Véase "Valores atípicos" para un análisis de los valores atípicos en la regresión.

Terminología de la regresión

Cuando los analistas e investigadores utilizan el término regresión por sí mismo, suelen referirse a la regresión lineal; la atención suele centrarse en el desarrollo de un modelo lineal para explicar la relación entre las variables predictoras y una variable numérica de resultado. En su sentido estadístico formal, la regresión también incluye los modelos no lineales que producen una relación funcional entre las variables predictoras y las variables de resultado. En la comunidad del aprendizaje automático, el término también se utiliza ocasionalmente en sentido amplio para referirse al uso de cualquier modelo predictivo que produzca un resultado numérico predicho (por oposición a los métodos de clasificación que predicen un resultado binario o categórico).

Predicción frente a explicación (elaboración de perfiles)

Históricamente, un uso primario de la regresión ha sido iluminar una supuesta relación lineal entre las variables predictoras y una variable de resultado. El objetivo ha sido comprender una relación y explicarla utilizando los datos a los que se ajustaba la regresión. En este caso, la atención se centra principalmente en la pendiente estimada de la ecuación de regresión, $\widehat{b}$Los economistas quieren conocer la relación entre el gasto de los consumidores y el crecimiento del PIB. Los responsables de la sanidad pública pueden querer saber si una campaña de información pública es eficaz para promover prácticas sexuales seguras. En estos casos, la atención no se centra en predecir casos individuales, sino en comprender la relación global entre las variables.

Con la llegada de los big data, la regresión se utiliza mucho para formar un modelo que prediga los resultados individuales de los nuevos datos (es decir, un modelo predictivo) en lugar de explicar los datos en cuestión. En este caso, los principales elementos de interés son los valores ajustados $\widehat{Y}$En marketing, la regresión puede utilizarse para predecir el cambio en los ingresos en respuesta al tamaño de una campaña publicitaria. Las universidades utilizan la regresión para predecir el GPA de los estudiantes basándose en sus puntuaciones SAT.

Un modelo de regresión que se ajusta bien a los datos se establece de forma que los cambios en X conducen a cambios en Y.Sin embargo, por sí sola, la ecuación de regresión no prueba la dirección de la causalidad. Las conclusiones sobre la causalidad deben proceder de una comprensión más amplia sobre la relación. Por ejemplo, una ecuación de regresión podría mostrar una relación definida entre el número de clics en un anuncio web y el número de conversiones. Es nuestro conocimiento del proceso de marketing, y no la ecuación de regresión, lo que nos lleva a la conclusión de que los clics en el anuncio conducen a las ventas, y no al revés.

Ideas clave

• La ecuación de regresión modela la relación entre una variable de respuesta Y y una variable de predicción X como una recta.

• Un modelo de regresión produce valores ajustados y residuos: predicciones de la respuesta y los errores de las predicciones.

• Los modelos de regresión suelen ajustarse por el método de los mínimos cuadrados.

• La regresión se utiliza tanto para predecir como para explicar.

Otras lecturas

Para un tratamiento en profundidad de la predicción frente a la explicación, véase el artículo de Galit Shmueli "¿Explicar o predecir?".

Regresión lineal múltiple

Cuando hay múltiples predictores, simplemente se amplía la ecuación para acomodarlos:

$$Y=b_0+b_1{X}_1+b_2{X}_2+...+b_p{X}_p+e$$

En lugar de una línea, ahora tenemos un modelo lineal: la relación entre cada coeficiente y su variable (característica) es lineal.

Términos clave de la regresión lineal múltiple

Error cuadrático medio

La raíz cuadrada del error cuadrático medio de la regresión (es la métrica más utilizada para comparar modelos de regresión).

Sinónimo

RMSE

Error típico residual

Lo mismo que el error cuadrático medio, pero ajustado a los grados de libertad.

Sinónimo

RSE

R-cuadrado

La proporción de varianza explicada por el modelo, de 0 a 1.

Sinónimos

coeficiente de determinación, $R^2$

estadístico t

El coeficiente de un predictor, dividido por el error típico del coeficiente, da una métrica para comparar la importancia de las variables en el modelo. Ver "Pruebas t".

Regresión ponderada

Regresión con los registros con pesos diferentes.

Todos los demás conceptos de la regresión lineal simple, como ajuste por mínimos cuadrados y la definición de valores ajustados y residuos, se extienden al entorno de la regresión lineal múltiple.Por ejemplo, los valores ajustados vienen dados por:

$${\widehat{Y}}_i={\widehat{b}}_0+{\widehat{b}}_1{X}_{1,i}+{\widehat{b}}_2{X}_{2,i}+...+{\widehat{b}}_p{X}_{p,i}$$

Ejemplo: Datos de vivienda del condado de King

Un ejemplo de uso de la regresión lineal múltiple es la estimación del valor de las casas.Los tasadores del condado deben estimar el valor de una casa para calcular los impuestos. Los profesionales inmobiliarios y los compradores de viviendas consultan sitios web populares como Zillow para determinar un precio justo. Aquí tienes algunas filas de datos sobre viviendas del condado de King (Seattle), Washington, del sitio house data.frame:

head(house[, c('AdjSalePrice', 'SqFtTotLiving', 'SqFtLot', 'Bathrooms',
               'Bedrooms', 'BldgGrade')])
Source: local data frame [6 x 6]

  AdjSalePrice SqFtTotLiving SqFtLot Bathrooms Bedrooms BldgGrade
         (dbl)         (int)   (int)     (dbl)    (int)     (int)
1       300805          2400    9373      3.00        6         7
2      1076162          3764   20156      3.75        4        10
3       761805          2060   26036      1.75        4         8
4       442065          3200    8618      3.75        5         7
5       297065          1720    8620      1.75        4         7
6       411781           930    1012      1.50        2         8

El método head del marco de datos pandas enumera las filas superiores:

subset = ['AdjSalePrice', 'SqFtTotLiving', 'SqFtLot', 'Bathrooms',
          'Bedrooms', 'BldgGrade']
house[subset].head()

El objetivo es predecir el precio de venta a partir de las demás variables. La función lm se ocupa del caso de regresión múltiple simplemente incluyendo más términos en el lado derecho de la ecuación; el argumento na.action=na.omit hace que el modelo elimine los registros que tienen valores perdidos:

house_lm <- lm(AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
               Bedrooms + BldgGrade,
               data=house, na.action=na.omit)

scikit-learn's LinearRegression también puede utilizarse para la regresión lineal múltiple:

predictors = ['SqFtTotLiving', 'SqFtLot', 'Bathrooms', 'Bedrooms', 'BldgGrade']
outcome = 'AdjSalePrice'

house_lm = LinearRegression()
house_lm.fit(house[predictors], house[outcome])

La impresión del objeto house_lm produce el siguiente resultado:

house_lm

Call:
lm(formula = AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
    Bedrooms + BldgGrade, data = house, na.action = na.omit)

Coefficients:
  (Intercept)  SqFtTotLiving        SqFtLot      Bathrooms
   -5.219e+05      2.288e+02     -6.047e-02     -1.944e+04
     Bedrooms      BldgGrade
   -4.777e+04      1.061e+05

Para un modelo LinearRegression, el intercepto y los coeficientes son los campos intercept_ y coef_ del modelo ajustado:

print(f'Intercept: {house_lm.intercept_:.3f}')
print('Coefficients:')
for name, coef in zip(predictors, house_lm.coef_):
    print(f' {name}: {coef}')

La interpretación de los coeficientes es como en la regresión lineal simple: el valor predicho $\widehat{Y}$ cambia según el coeficiente $b_j$ por cada cambio unitario en $X_j$ asumiendo todas las demás variables, $X_k$ para $k\ne j$Por ejemplo, añadir un pie cuadrado acabado más a una casa aumenta el valor estimado en unos 229 $; añadir 1.000 pies cuadrados acabados implica que el valor aumentará en 228.800 $.

Evaluación del modelo

La métrica de rendimiento más importante desde la perspectiva de la ciencia de datos es el error cuadrático medio, o RMSE. El RMSE es la raíz cuadrada del error cuadrático medio en la predicción de ${\widehat{y}}_i$ predichos:

$$RMSE=\sqrt{\frac{{\sum }_{i=1}^n{\left(y_i-{\widehat{y}}_i\right)}^2}{n}}$$

Mide la precisión global del modelo y sirve de base para compararlo con otros modelos (incluidos los modelos ajustados mediante técnicas de aprendizaje automático). Similar al RMSE es el error estándar residual, o RSE. En este caso tenemos p predictores, y el RSE viene dado por:

$$RSE=\sqrt{\frac{{\sum }_{i=1}^n{\left(y_i-{\widehat{y}}_i\right)}^2}{\left(n-p-1\right)}}$$

La única diferencia es que el denominador son los grados de libertad, en lugar del número de registros (véase " Grados de libertad"). En la práctica, para la regresión lineal, la diferencia entre RMSE y RSE es muy pequeña, sobre todo para las aplicaciones de big data.

La función summary de R calcula el RSE y otras métricas de un modelo de regresión:

summary(house_lm)

Call:
lm(formula = AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
    Bedrooms + BldgGrade, data = house, na.action = na.omit)

Residuals:
     Min       1Q   Median       3Q      Max
-1199479  -118908   -20977    87435  9473035

Coefficients:
                Estimate Std. Error t value Pr(>|t|)
(Intercept)   -5.219e+05  1.565e+04 -33.342  < 2e-16 ***
SqFtTotLiving  2.288e+02  3.899e+00  58.694  < 2e-16 ***
SqFtLot       -6.047e-02  6.118e-02  -0.988    0.323
Bathrooms     -1.944e+04  3.625e+03  -5.363 8.27e-08 ***
Bedrooms      -4.777e+04  2.490e+03 -19.187  < 2e-16 ***
BldgGrade      1.061e+05  2.396e+03  44.277  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 261300 on 22681 degrees of freedom
Multiple R-squared:  0.5406,	Adjusted R-squared:  0.5405
F-statistic:  5338 on 5 and 22681 DF,  p-value: < 2.2e-16

scikit-learn proporciona una serie de métricas para la regresión y la clasificación. Aquí, utilizamos mean_squared_error para obtener el RMSE yr2_score para el coeficiente de determinación:

fitted = house_lm.predict(house[predictors])
RMSE = np.sqrt(mean_squared_error(house[outcome], fitted))
r2 = r2_score(house[outcome], fitted)
print(f'RMSE: {RMSE:.0f}')
print(f'r2: {r2:.4f}')

Utiliza statsmodels para obtener un análisis más detallado del modelo de regresión en Python:

model = sm.OLS(house[outcome], house[predictors].assign(const=1))
results = model.fit()
results.summary()

El método pandas assign , tal y como se utiliza aquí, añade una columna constante con valor 1 a los predictores. Esto es necesario para modelizar el intercepto.

Otra métrica útil que verás en la salida del software es el coeficiente de determinación, también llamado estadístico R-cuadrado o $R^2$El R-cuadrado oscila entre 0 y 1 y mide la proporción de variación de los datos que se tiene en cuenta en el modelo. Es útil sobre todo en los usos explicativos de la regresión, en los que quieres evaluar lo bien que se ajusta el modelo a los datos. La fórmula de $R^2$ es

$$R^2=1-\frac{{\sum }_{i=1}^n{\left(y_i-{\widehat{y}}_i\right)}^2}{{\sum }_{i=1}^n{\left(y_i-\overline{y}\right)}^2}$$

El denominador es proporcional a la varianza de Y.La salida de R también informa de un R-cuadrado ajustado, que se ajusta a los grados de libertad, penalizando efectivamente la adición de más predictores a un modelo. Rara vez es esto significativamente diferente de R-cuadrado en regresión múltiple con grandes conjuntos de datos.

Junto con los coeficientes estimados, R y statsmodels indican el error estándar de los coeficientes (SE) y un estadístico t:

$$t_b=\frac{\widehat{b}}{\text{SE}\left(\widehat{b}\right)}$$

El estadístico t -y su imagen especular, el valor p- mide hasta qué punto un coeficiente es "estadísticamente significativo", es decir, está fuera del rango de lo que podría producir una disposición aleatoria del predictor y la variable objetivo. Cuanto mayor sea el estadístico t (y menor el valor p), más significativo será el predictor. Dado que la parsimonia es una valiosa característica del modelo, resulta útil disponer de una herramienta como ésta para orientar la elección de las variables que se incluirán como predictores (véase "Selección de modelos y regresión por pasos").

Advertencia

Además del estadístico t, R y otros paquetes informarán a menudo de un valor p(Pr(>|t|) en la salida de R ) y del estadístico F. Los científicos de datos no suelen implicarse demasiado en la interpretación de estos estadísticos, ni en la cuestión de la significación estadística.Los científicos de datos se centran principalmente en el estadístico t como guía útil para decidir si incluir o no un predictor en un modelo.Los estadísticos t altos (que van acompañados de valores p cercanos a 0) indican que un predictor debe mantenerse en un modelo, mientras que los estadísticos t muy bajos indican que un predictor podría descartarse. Para más información, consulta "Valor p".

Validación cruzada

Las métricas clásicas de regresión estadística(^R2, estadísticos F y valores p) son todas métricas "en la muestra": se aplican a los mismos datos que se utilizaron para ajustar el modelo. Intuitivamente, puedes ver que tendría mucho sentido apartar algunos de los datos originales, no utilizarlos para ajustar el modelo y, a continuación, aplicar el modelo a los datos apartados (retenidos) para ver su eficacia. Normalmente, utilizarías la mayoría de los datos para ajustar el modelo y una parte más pequeña para probarlo.

Esta idea de la validación "fuera de la muestra" no es nueva, pero no se impuso realmente hasta que se generalizaron los conjuntos de datos más grandes; con un conjunto de datos pequeño, los analistas suelen querer utilizar todos los datos y ajustar el mejor modelo posible.

Sin embargo, el uso de una muestra "holdout" te deja sujeto a cierta incertidumbre que surge simplemente de la variabilidad de la pequeña muestra "holdout". ¿Cómo de diferente sería la evaluación si seleccionaras una muestra "holdout" diferente?

La validación cruzada amplía la idea de una muestra de retención a múltiples muestras de retención secuenciales. El algoritmo para la validación cruzada k-fold básica es el siguiente:

1. Reserva 1/k de los datos como muestra de reserva.

2. Entrena el modelo con los datos restantes.

3. Aplica (puntúa) el modelo a la retención 1/k, y registra las métricas de evaluación del modelo necesarias.

4. Restaura el primer 1/k de los datos, y aparta el siguiente 1/k (excluyendo los registros que se escogieron la primera vez).

5. Repite los pasos 2 y 3.

6. Repítelo hasta que se haya utilizado cada registro en la parte de retención.

7. Promediar o combinar de otro modo las métricas de evaluación del modelo.

La división de los datos en la muestra de entrenamiento y la muestra retenida también se denomina un pliegue.

Selección de modelos y regresión por pasos

En algunos problemas, podrían utilizarse muchas variables como predictores en una regresión.Por ejemplo, para predecir el valor de una casa, podrían utilizarse variables adicionales como el tamaño del sótano o el año de construcción. En R, es fácil añadirlas a la ecuación de regresión:

house_full <- lm(AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
                 Bedrooms + BldgGrade + PropertyType + NbrLivingUnits +
                 SqFtFinBasement + YrBuilt + YrRenovated +
                 NewConstruction,
               data=house, na.action=na.omit)

En Python, tenemos que convertir las variables categóricas y booleanas en números:

predictors = ['SqFtTotLiving', 'SqFtLot', 'Bathrooms', 'Bedrooms', 'BldgGrade',
              'PropertyType', 'NbrLivingUnits', 'SqFtFinBasement', 'YrBuilt',
              'YrRenovated', 'NewConstruction']

X = pd.get_dummies(house[predictors], drop_first=True)
X['NewConstruction'] = [1 if nc else 0 for nc in X['NewConstruction']]

house_full = sm.OLS(house[outcome], X.assign(const=1))
results = house_full.fit()
results.summary()

Sin embargo, añadir más variables no significa necesariamente que tengamos un modelo mejor.Los estadísticos utilizan el principio de la navaja de Occam para guiar la elección de un modelo: en igualdad de condiciones, debe utilizarse un modelo más sencillo con preferencia a un modelo más complicado.

La inclusión de variables adicionales siempre reduce el RMSE y aumenta $R^2$ para los datos de entrenamiento. Por lo tanto, no son apropiados para ayudar a guiar la elección del modelo. Un enfoque para incluir la complejidad del modelo consiste en utilizar el valor ajustado de $R^2$:

${\displaystyle R_{adj}^2=1-<!--\; -\; -->(1-<!--\; -\; -->R^2)\frac{n-<!--\; -\; -->1}{n-<!--\; -\; -->P-<!--\; -\; -->1}}$

Aquí, n es el número de registros y P es el número de variables del modelo.

En la década de 1970, Hirotugu Akaike, eminente estadístico japonés, desarrolló una métrica denominada AIC (Criterio de Información de Akaike) que penaliza la adición de términos a un modelo.En el caso de la regresión, el AIC tiene la forma:

• AIC = 2P + n log(RSS/n)

donde P es el número de variables y n es el número de registros. El objetivo es encontrar el modelo que minimice el AIC; los modelos con k variables adicionales más se penalizan con 2k.

AIC, BIC y Cp de Mallows

La fórmula del AIC puede parecer un poco misteriosa, pero en realidad se basa en resultados asintóticos de la teoría de la información. Hay varias variantes del AIC:

AICc

Una versión del AIC corregida para tamaños de muestra pequeños.

BIC o criterio de información bayesiano

Similar al AIC, con una penalización mayor por incluir variables adicionales al modelo.

Malvas Cp

Una variante del AIC desarrollada por Colin Mallows.

Normalmente se presentan como métricas dentro de la muestra (es decir, sobre los datos de entrenamiento), y los científicos de datos que utilizan datos retenidos para la evaluación de modelos no necesitan preocuparse por las diferencias entre ellos ni por la teoría subyacente que los sustenta.

¿Cómo hallamos el modelo que minimiza el AIC o maximiza el ajustado $R^2$Una forma es buscar en todos los modelos posibles, un enfoque llamado regresión de todos los subconjuntos.Esto es caro desde el punto de vista informático y no es factible para problemas con grandes datos y muchas variables. Una alternativa atractiva es utilizar la regresión por pasos. Se podría empezar con un modelo completo e ir eliminando sucesivamente las variables que no contribuyen de forma significativa. Esto se denomina eliminación hacia atrás. Otra posibilidad es empezar con un modelo constante e ir añadiendo variables sucesivamente(selección hacia delante). Como tercera opción, también podemos añadir y descartar sucesivamente predictores para encontrar un modelo que reduzca el AIC o el ajustado $R^2$El paquete MASS en R de Venebles y Ripley ofrece una función de regresión por pasos llamada stepAIC:

library(MASS)
step <- stepAIC(house_full, direction="both")
step

Call:
lm(formula = AdjSalePrice ~ SqFtTotLiving + Bathrooms + Bedrooms +
    BldgGrade + PropertyType + SqFtFinBasement + YrBuilt, data = house,
    na.action = na.omit)

Coefficients:
              (Intercept)              SqFtTotLiving
                6.179e+06                  1.993e+02
                Bathrooms                   Bedrooms
                4.240e+04                 -5.195e+04
                BldgGrade  PropertyTypeSingle Family
                1.372e+05                  2.291e+04
    PropertyTypeTownhouse            SqFtFinBasement
                8.448e+04                  7.047e+00
                  YrBuilt
               -3.565e+03

scikit-learn no tiene ninguna implementación para la regresión por pasos. Hemos implementado las funciones stepwise_selection, forward_selection, y backward_elimination en nuestro paquete dmba:

y = house[outcome]

def train_model(variables): 
    if len(variables) == 0:
        return None
    model = LinearRegression()
    model.fit(X[variables], y)
    return model

def score_model(model, variables): 
    if len(variables) == 0:
        return AIC_score(y, [y.mean()] * len(y), model, df=1)
    return AIC_score(y, model.predict(X[variables]), model)

best_model, best_variables = stepwise_selection(X.columns, train_model,
                                                score_model, verbose=True)

print(f'Intercept: {best_model.intercept_:.3f}')
print('Coefficients:')
for name, coef in zip(best_variables, best_model.coef_):
    print(f' {name}: {coef}')

Define una función que devuelva un modelo ajustado para un conjunto dado de variables.

Define una función que devuelva una puntuación para un modelo y un conjunto de variables dados. En este caso, utilizamos el AIC_score implementado en el paquete dmba.

La función eligió un modelo en el que se eliminaron varias variables de house_full:SqFtLot, NbrLivingUnits, YrRenovated, y NewConstruction.

Más sencillas aún son la selección hacia delante y la selección hacia atrás. En la selección hacia delante, empiezas sin predictores y los vas añadiendo uno a uno, añadiendo en cada paso el predictor que tenga la mayor contribución a $R^2$y te detienes cuando la contribución deja de ser estadísticamente significativa. En la selección hacia atrás, o eliminación hacia atrás, empiezas con el modelo completo y quitas los predictores que no son estadísticamente significativos hasta que te quedas con un modelo en el que todos los predictores son estadísticamente significativos.

La regresiónpenalizada es similar en espíritu al AIC.En lugar de buscar explícitamente entre un conjunto discreto de modelos, la ecuación de ajuste del modelo incorpora una restricción que penaliza al modelo por demasiadas variables (parámetros). En lugar de eliminar por completo las variables predictoras -como ocurre con la selección por pasos, hacia delante y hacia atrás-, la regresión penalizada aplica la penalización reduciendo los coeficientes, en algunos casos hasta casi cero. Los métodos habituales de regresión penalizada son la regresión ridge y la regresión lasso.

La regresión escalonada y la regresión de todos los subconjuntos son métodos dentro de la muestra para evaluar y afinar los modelos.Esto significa que la selección del modelo está posiblemente sujeta a sobreajuste (ajuste al ruido de los datos) y puede no funcionar tan bien cuando se aplica a nuevos datos. Un enfoque habitual para evitarlo es utilizar la validación cruzada para validar los modelos. En la regresión lineal, el sobreajuste no suele ser un problema importante, debido a la sencilla estructura global (lineal) impuesta a los datos. Para tipos de modelos más sofisticados, en particular los procedimientos iterativos que responden a la estructura local de los datos, la validación cruzada es una herramienta muy importante; para más detalles, consulta "Validación cruzada".

Regresión ponderada

Los estadísticos utilizan la regresión ponderada para diversos fines; en particular, es importante para el análisis de encuestas complejas. Los científicos de datos pueden encontrar útil la regresión ponderada en dos casos:

• Ponderación de varianza inversa cuando las distintas observaciones se han medido con distinta precisión; las de mayor varianza reciben ponderaciones más bajas.

• Análisis de datos en el que las filas representan casos múltiples; la variable de peso codifica cuántas observaciones originales representa cada fila.

Por ejemplo, con los datos de vivienda, las ventas más antiguas son menos fiables que las más recientes. Utilizando el DocumentDate para determinar el año de la venta, podemos calcular un Weight como el número de años transcurridos desde 2005 (el inicio de los datos):

R

library(lubridate)
house$Year = year(house$DocumentDate)
house$Weight = house$Year - 2005

Python

house['Year'] = [int(date.split('-')[0]) for date in house.DocumentDate]
house['Weight'] = house.Year - 2005

Podemos calcular una regresión ponderada con la función lm utilizando el argumento weight:

house_wt <- lm(AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
                 Bedrooms + BldgGrade,
               data=house, weight=Weight)
round(cbind(house_lm=house_lm$coefficients,
            house_wt=house_wt$coefficients), digits=3)

                 house_lm    house_wt
(Intercept)   -521871.368 -584189.329
SqFtTotLiving     228.831     245.024
SqFtLot            -0.060      -0.292
Bathrooms      -19442.840  -26085.970
Bedrooms       -47769.955  -53608.876
BldgGrade      106106.963  115242.435

Los coeficientes de la regresión ponderada son ligeramente diferentes de los de la regresión original.

La mayoría de los modelos de scikit-learn aceptan pesos como argumento de palabra clave sample_weight en la llamada al método fit:

predictors = ['SqFtTotLiving', 'SqFtLot', 'Bathrooms', 'Bedrooms', 'BldgGrade']
outcome = 'AdjSalePrice'

house_wt = LinearRegression()
house_wt.fit(house[predictors], house[outcome], sample_weight=house.Weight)

Ideas clave

• La regresión lineal múltiple modela la relación entre una variable de respuesta Y y múltiples variables predictoras $X_1,...,X_p$.

• Las métricas más importantes para evaluar un modelo son el error cuadrático medio (RMSE) y R-cuadrado(^R2).

• El error típico de los coeficientes puede utilizarse para medir la fiabilidad de la contribución de una variable a un modelo.

• La regresión por pasos es una forma de determinar automáticamente qué variables deben incluirse en el modelo.

• La regresión ponderada se utiliza para dar más o menos peso a determinados registros a la hora de ajustar la ecuación.

Otras lecturas

Puedes encontrar un excelente tratamiento de la validación cruzada y el remuestreo en An Introduction to Statistical Learning, de Gareth James, Daniela Witten, Trevor Hastie y Robert Tibshirani (Springer, 2013).

Predicción mediante regresión

El objetivo principal de la regresión en la ciencia de datos es la predicción.Es útil tener esto en cuenta, ya que la regresión, al ser un método estadístico antiguo y establecido, viene con un bagaje que es más relevante para su papel tradicional como herramienta de modelización explicativa que para la predicción.

Términos clave para la predicción mediante regresión

Intervalo de predicción

Un intervalo de incertidumbre en torno a un valor individual predicho.

Extrapolación

Extensión de un modelo más allá del alcance de los datos utilizados para ajustarlo.

Los peligros de la extrapolación

Los modelos de regresión no deben utilizarse para extrapolar más allá del rango de los datos (dejando a un lado el uso de la regresión para la previsión de series temporales). El modelo sólo es válido para los valores predictores para los que los datos tienen valores suficientes (incluso en el caso de que se disponga de datos suficientes, podría haber otros problemas-ver "Diagnóstico de la regresión"). Como caso extremo, supongamos que model_lm se utiliza para predecir el valor de un solar vacío de 5.000 metros cuadrados. En tal caso, todos los predictores relacionados con el edificio tendrían un valor de 0, y la ecuación de regresión arrojaría una predicción absurda de -521.900 + 5.000 × -,0605 = -522.202 $. ¿Por qué ocurre esto? Los datos sólo contienen parcelas con edificios; no hay registros correspondientes a terrenos baldíos. En consecuencia, el modelo no dispone de información que le indique cómo predecir el precio de venta de un terreno baldío.

Intervalos de confianza y predicción

Gran parte de la estadística implica comprender y medir la variabilidad (incertidumbre). Los estadísticos t y los valores p que aparecen en los resultados de la regresión se ocupan de esto de manera formal, lo que a veces resulta útil para la selección de variables (consulta "Evaluación del modelo"). Más útiles son los intervalos de confianza, que son intervalos de incertidumbre situados alrededor de los coeficientes de regresión y las predicciones. Una forma fácil de entender esto es mediante el bootstrap (consulta "El bootstrap" para obtener más detalles sobre el procedimiento general de bootstrap). Los intervalos de confianza de regresión más comunes que se encuentran en los resultados del software son los de los parámetros de regresión (coeficientes).Aquí tienes un algoritmo de bootstrap para generar intervalos de confianza para los parámetros de regresión (coeficientes) de un conjunto de datos con P predictores y n registros (filas):

1. Considera cada fila (incluida la variable de resultado) como un único "billete" y coloca todos los n billetes en una caja.

2. Saca un billete al azar, anota los valores y colócalo en la casilla.

3. Repite el paso 2 n veces; ahora tienes una remuestra bootstrap.

4. Ajusta una regresión a la muestra bootstrap y registra los coeficientes estimados.

5. Repite los pasos 2 a 4, digamos, 1.000 veces.

6. Ahora tienes 1.000 valores bootstrap para cada coeficiente; encuentra los percentiles adecuados para cada uno (por ejemplo, el 5º y el 95º para un intervalo de confianza del 90%).

Puedes utilizar la función Boot de R para generar intervalos de confianza bootstrap reales para los coeficientes, o puedes utilizar simplemente los intervalos basados en fórmulas que son una salida rutinaria de R. El significado conceptual y la interpretación son los mismos, y no tienen una importancia central para los científicos de datos, porque se refieren a los coeficientes de regresión. De mayor interés para los científicos de datos son los intervalos en torno a los valores predichos de y (${\widehat{Y}}_i$). La incertidumbre en torno a ${\widehat{Y}}_i$ procede de dos fuentes:

• Incertidumbre sobre cuáles son las variables predictoras relevantes y sus coeficientes (véase el algoritmo bootstrap anterior)

• Error adicional inherente a los puntos de datos individuales

El error de punto de datos individual puede considerarse del siguiente modo: aunque supiéramos con certeza cuál es la ecuación de regresión (p. ej, si tuviéramos un gran número de registros para ajustarla), los valores reales de los resultados para un conjunto determinado de valores predictores variarán. Por ejemplo, varias casas -cada una con 8 habitaciones, un terreno de 6.500 pies cuadrados, 3 baños y un sótano- podrían tener valores diferentes. Podemos modelizar este error individual con los residuos de los valores ajustados. El algoritmo bootstrap para modelizar tanto el error del modelo de regresión como el error de los puntos de datos individuales sería el siguiente:

1. Toma una muestra bootstrap de los datos (explicada antes con más detalle).

2. Ajusta la regresión y predice el nuevo valor.

3. Toma un único residuo al azar del ajuste de regresión original, añádelo al valor predicho y anota el resultado.

4. Repite los pasos del 1 al 3, digamos, 1.000 veces.

5. Halla los percentiles 2,5 y 97,5 de los resultados.

Ideas clave

• La extrapolación más allá del rango de los datos puede conducir a error.

• Los intervalos de confianza cuantifican la incertidumbre en torno a los coeficientes de regresión.

• Los intervalos de predicción cuantifican la incertidumbre en las predicciones individuales.

• La mayoría de los programas informáticos, incluido R, producirán predicciones e intervalos de confianza en la salida predeterminada o especificada, utilizando fórmulas.

• El bootstrap también puede utilizarse para producir intervalos de predicción y de confianza; la interpretación y la idea son las mismas.

¿Intervalo de predicción o intervalo de confianza?

Un intervalo de predicción se refiere a la incertidumbre en torno a un único valor, mientras que un intervalo de confianza se refiere a una media u otra estadística calculada a partir de múltiples valores.Por tanto, un intervalo de predicción será normalmente mucho más amplio que un intervalo de confianza para el mismo valor. Modelizamos este error de valor individual en el modelo bootstrap seleccionando un residuo individual para añadirlo al valor predicho. ¿Cuál deberías utilizar? Eso depende del contexto y del objetivo del análisis, pero, en general, los científicos de datos están interesados en predicciones individuales específicas, por lo que un intervalo de predicción sería más apropiado. Utilizar un intervalo de confianza cuando deberías utilizar un intervalo de predicción subestimará en gran medida la incertidumbre de un valor predicho determinado.

Variables factoriales en la regresión

Las variablesfactoriales, también denominadas variables categóricas, adoptan un número limitado de valores discretos. Por ejemplo, la finalidad de un préstamo puede ser "consolidación de deudas", "boda", "coche", etcétera. La variable binaria (sí/no), también llamada variable indicadora, es un caso especial de variable factorial.La regresión requiere entradas numéricas, por lo que las variables factoriales deben recodificarse para utilizarlas en el modelo.El enfoque más habitual consiste en convertir una variable en un conjunto de variables binarias ficticias.

Términos clave de las variables factoriales

Variables ficticias

Variables binarias 0-1 derivadas de la recodificación de datos factoriales para su uso en modelos de regresión y otros.

Codificación de referencia

Es el tipo de codificación más utilizado por los estadísticos, en el que un nivel de un factor se utiliza como referencia y los demás factores se comparan con ese nivel.

Sinónimo

codificación del tratamiento

Un codificador en caliente

Un tipo común de codificación utilizado en la comunidad del aprendizaje automático en el que se conservan todos los niveles de los factores.Aunque es útil para determinados algoritmos de aprendizaje automático, este enfoque no es apropiado para la regresión lineal múltiple.

Codificación de desviaciones

Tipo de codificación que compara cada nivel con la media global, en contraposición al nivel de referencia.

Sinónimo

suma contrastes

Variables ficticias Representación

En los datos de vivienda del condado de King, hay una variable factorial para el tipo de propiedad; a continuación se muestra un pequeño subconjunto de seis registros:

R:

head(house[, 'PropertyType'])
Source: local data frame [6 x 1]

   PropertyType
         (fctr)
1     Multiplex
2 Single Family
3 Single Family
4 Single Family
5 Single Family
6     Townhouse

Pitón:

house.PropertyType.head()

Hay tres valores posibles: Multiplex, Single Family, y Townhouse. Para utilizar esta variable factorial, tenemos que convertirla en un conjunto de variables binarias. Para ello, creamos una variable binaria para cada valor posible de la variable factorial. Para hacerlo en R, utilizamos la función model.matrix:³

prop_type_dummies <- model.matrix(~PropertyType -1, data=house)
head(prop_type_dummies)
  PropertyTypeMultiplex PropertyTypeSingle Family PropertyTypeTownhouse
1                     1                         0                     0
2                     0                         1                     0
3                     0                         1                     0
4                     0                         1                     0
5                     0                         1                     0
6                     0                         0                     1

La función model.matrix convierte un marco de datos en una matriz adecuada para un modelo lineal. La variable factorial PropertyType, que tiene tres niveles distintos, se representa como una matriz con tres columnas.En la comunidad del aprendizaje automático, esta representación se denomina codificación en caliente (véase "Codificador en caliente").

En Python, podemos convertir variables categóricas en ficticias utilizando el método pandas get_dummies :

pd.get_dummies(house['PropertyType']).head() 
pd.get_dummies(house['PropertyType'], drop_first=True).head() 

Por defecto, devuelve una codificación en caliente de la variable categórica.

El argumento de la palabra clave drop_first devolverá P - 1 columnas. Utilízalo para evitar el problema de la multicolinealidad.

En determinados algoritmos de aprendizaje automático, como los vecinos más próximos y los modelos en árbol, una codificación en caliente es la forma estándar de representar las variables factoriales (por ejemplo, consulta "Modelos en árbol").

En el ámbito de la regresión, una variable factorial con P niveles distintos suele representarse mediante una matriz con sólo P - 1 columnas. Esto se debe a que un modelo de regresión suele incluir un término de intercepción. Con una intercepción, una vez que has definido los valores de P - 1 binarias, el valor de la Pthes conocido y podría considerarse redundante. Si añades la columna Pth, se producirá un error de multicolinealidad (véase "Multicolinealidad").

La representación por defecto en R es utilizar el primer nivel del factor como referencia e interpretar los niveles restantes en relación con ese factor:

lm(AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
       Bedrooms + BldgGrade + PropertyType, data=house)

Call:
lm(formula = AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
    Bedrooms + BldgGrade + PropertyType, data = house)

Coefficients:
              (Intercept)              SqFtTotLiving
               -4.468e+05                  2.234e+02
                  SqFtLot                  Bathrooms
               -7.037e-02                 -1.598e+04
                 Bedrooms                  BldgGrade
               -5.089e+04                  1.094e+05
PropertyTypeSingle Family      PropertyTypeTownhouse
               -8.468e+04                 -1.151e+05

El método get_dummies toma como argumento opcional la palabra clave drop_first para excluir el primer factor como referencia:

predictors = ['SqFtTotLiving', 'SqFtLot', 'Bathrooms', 'Bedrooms',
              'BldgGrade', 'PropertyType']

X = pd.get_dummies(house[predictors], drop_first=True)

house_lm_factor = LinearRegression()
house_lm_factor.fit(X, house[outcome])

print(f'Intercept: {house_lm_factor.intercept_:.3f}')
print('Coefficients:')
for name, coef in zip(X.columns, house_lm_factor.coef_):
    print(f' {name}: {coef}')

La salida de la regresión R muestra dos coeficientes correspondientes a PropertyType: PropertyTypeSingle Family y PropertyTypeTownhouse. No hay coeficiente de Multiplex, ya que está definido implícitamente cuando PropertyTypeSingle Family == 0 y PropertyTypeTownhouse == 0. Los coeficientes se interpretan como relativos a Multiplex, de modo que una vivienda que es Single Family vale casi 85.000 $ menos, y una vivienda que es Townhouse vale más de 150.000 $ menos.⁴

Diferentes codificaciones de factores

Hay varias formas distintas de codificar las variables factoriales, conocidas como sistemas de codificación de contrastes.Por ejemplo, la codificación de desviaciones, también conocida como contrastes de sumas, compara cada nivel con la media general. Otro contraste es la codificación polinómica, que es adecuada para los factores ordenados; véase el apartado "Variables factoriales ordenadas". Con la excepción de los factores ordenados, los científicos de datos no suelen encontrar ningún tipo de codificación aparte de la codificación de referencia o de un codificador en caliente.

Variables factoriales con muchos niveles

Algunas variables factoriales pueden producir un número enorme de variables binarias ficticias -los códigos postales son una variable factorial, y hay 43.000 códigos postales en EE.UU.-. En tales casos, es útil explorar los datos y las relaciones entre las variables predictoras y el resultado, para determinar si las categorías contienen información útil. En caso afirmativo, debes decidir además si es útil conservar todos los factores, o si deben consolidarse los niveles.

En el condado de King hay 80 códigos postales con venta de casas:

table(house$ZipCode)

98001 98002 98003 98004 98005 98006 98007 98008 98010 98011 98014 98019
  358   180   241   293   133   460   112   291    56   163    85   242
98022 98023 98024 98027 98028 98029 98030 98031 98032 98033 98034 98038
  188   455    31   366   252   475   263   308   121   517   575   788
98039 98040 98042 98043 98045 98047 98050 98051 98052 98053 98055 98056
   47   244   641     1   222    48     7    32   614   499   332   402
98057 98058 98059 98065 98068 98070 98072 98074 98075 98077 98092 98102
    4   420   513   430     1    89   245   502   388   204   289   106
98103 98105 98106 98107 98108 98109 98112 98113 98115 98116 98117 98118
  671   313   361   296   155   149   357     1   620   364   619   492
98119 98122 98125 98126 98133 98136 98144 98146 98148 98155 98166 98168
  260   380   409   473   465   310   332   287    40   358   193   332
98177 98178 98188 98198 98199 98224 98288 98354
  216   266   101   225   393     3     4     9

El método value_counts de pandas marcos de datos devuelve la misma información:

pd.DataFrame(house['ZipCode'].value_counts()).transpose()

ZipCode es una variable importante, ya que es una aproximación al efecto de la ubicación sobre el valor de una casa. Incluir todos los niveles requiere 79 coeficientes correspondientes a 79 grados de libertad. El modelo original house_lm sólo tiene 5 grados de libertad; véase "Evaluación del modelo". Además, varios códigos postales sólo tienen una venta. En algunos problemas, puedes consolidar un código postal utilizando los dos o tres primeros dígitos, correspondientes a una región geográfica submetropolitana. En el caso del condado de King, casi todas las ventas se producen en 980xx o 981xx, así que esto no ayuda.

Un enfoque alternativo es agrupar los códigos postales según otra variable, como el precio de venta. Aún mejor es formar grupos de códigos postales utilizando los residuos de un modelo inicial. El siguiente código dplyr en R consolida los 80 códigos postales en cinco grupos basados en la mediana del residuo de la regresión house_lm:

zip_groups <- house %>%
  mutate(resid = residuals(house_lm)) %>%
  group_by(ZipCode) %>%
  summarize(med_resid = median(resid),
            cnt = n()) %>%
  arrange(med_resid) %>%
  mutate(cum_cnt = cumsum(cnt),
         ZipGroup = ntile(cum_cnt, 5))
house <- house %>%
  left_join(select(zip_groups, ZipCode, ZipGroup), by='ZipCode')

Se calcula la mediana del residuo para cada código postal, y se utiliza la función ntile para dividir los códigos postales, ordenados por la mediana, en cinco grupos. Consulta "Variables de confusión" para ver un ejemplo de cómo se utiliza este término en una regresión que mejora el ajuste original.

En Python podemos calcular esta información del siguiente modo:

zip_groups = pd.DataFrame([
    *pd.DataFrame({
        'ZipCode': house['ZipCode'],
        'residual' : house[outcome] - house_lm.predict(house[predictors]),
    })
    .groupby(['ZipCode'])
    .apply(lambda x: {
        'ZipCode': x.iloc[0,0],
        'count': len(x),
        'median_residual': x.residual.median()
    })
]).sort_values('median_residual')
zip_groups['cum_count'] = np.cumsum(zip_groups['count'])
zip_groups['ZipGroup'] = pd.qcut(zip_groups['cum_count'], 5, labels=False,
                                 retbins=False)

to_join = zip_groups[['ZipCode', 'ZipGroup']].set_index('ZipCode')
house = house.join(to_join, on='ZipCode')
house['ZipGroup'] = house['ZipGroup'].astype('category')

El concepto de utilizar los residuos para ayudar a guiar el ajuste de la regresión es un paso fundamental en el proceso de modelización; véase "Diagnóstico de la regresión".

Variables factoriales ordenadas

Algunas variables factoriales reflejan niveles de un factor; éstas se denominan variables factoriales ordenadas o variables categóricas ordenadas. Por ejemplo, la calificación del préstamo podría ser A, B, C, y así sucesivamente: cada calificación conlleva más riesgo que la calificación anterior. A menudo, las variables factoriales ordenadas pueden convertirse en valores numéricos y utilizarse tal cual. Por ejemplo, la variable BldgGrade es una variable factorial ordenada. En la Tabla 4-1 se muestran varios tipos de grados. Aunque los grados tienen un significado específico, el valor numérico está ordenado de menor a mayor, lo que corresponde a las viviendas de mayor grado. Con el modelo de regresión house_lm, ajustado en "Regresión lineal múltiple",BldgGrade se trató como una variable numérica.

Tabla 4-1. Grados de construcción y equivalentes numéricos

Valor	Descripción
1	Cabina
2	Subestándar
5	Feria
10	Muy buena
12	Lujo
13	Mansión

Tratar los factores ordenados como una variable numérica conserva la información contenida en la ordenación, que se perdería si se convirtiera en un factor.

Ideas clave

• Las variables factoriales deben convertirse en variables numéricas para utilizarlas en una regresión.

• El método más habitual para codificar una variable factorial con P valores distintos es representarlos mediante P - 1 variables ficticias.

• Una variable factorial con muchos niveles, incluso en conjuntos de datos muy grandes, puede tener que consolidarse en una variable con menos niveles.

• Algunos factores tienen niveles ordenados y pueden representarse como una única variable numérica.

Interpretar la ecuación de regresión

En la ciencia de datos, el uso más importante de la regresión es predecir alguna variable dependiente (resultado).En algunos casos, sin embargo, puede ser útil obtener información de la propia ecuación para comprender la naturaleza de la relación entre los predictores y el resultado. Esta sección proporciona orientación para examinar la ecuación de regresión e interpretarla.

Términos clave para interpretar la ecuación de regresión

Variables correlacionadas

Variables que tienden a moverse en la misma dirección: cuando una sube también lo hace la otra, y viceversa (con correlación negativa, cuando una sube la otra baja). Cuando las variables predictoras están muy correlacionadas, es difícil interpretar los coeficientes individuales.

Multicolinealidad

Cuando las variables predictoras tienen una correlación perfecta, o casi perfecta, la regresión puede ser inestable o imposible de calcular.

Sinónimo

colinealidad

Variables de confusión

Predictor importante que, cuando se omite, da lugar a relaciones espurias en una ecuación de regresión.

Efectos principales

La relación entre un predictor y la variable de resultado, independiente de otras variables.

Interacciones

Una relación interdependiente entre dos o más predictores y la respuesta.

Predictores correlacionados

En la regresión múltiple, las variables predictoras suelen estar correlacionadas entre sí. A modo de ejemplo, examina los coeficientes de regresión del modelo step_lm, ajustado en "Selección de modelos y regresión por pasos".

R:

step_lm$coefficients
              (Intercept)             SqFtTotLiving                 Bathrooms
             6.178645e+06              1.992776e+02              4.239616e+04
                 Bedrooms                 BldgGrade PropertyTypeSingle Family
            -5.194738e+04              1.371596e+05              2.291206e+04
    PropertyTypeTownhouse           SqFtFinBasement                   YrBuilt
             8.447916e+04              7.046975e+00             -3.565425e+03

Pitón:

print(f'Intercept: {best_model.intercept_:.3f}')
print('Coefficients:')
for name, coef in zip(best_variables, best_model.coef_):
    print(f' {name}: {coef}')

El coeficiente de Bedrooms es negativo, lo que implica que añadir un dormitorio a una casa reducirá su valor. ¿Cómo puede ser? Porque las variables predictoras están correlacionadas: las casas más grandes suelen tener más dormitorios, y es el tamaño lo que determina el valor de la casa, no el número de dormitorios. Considera dos casas exactamente del mismo tamaño: es razonable esperar que una casa con más dormitorios, pero más pequeños, se considere menos deseable.

Tener predictores correlacionados puede dificultar la interpretación del signo y el valor de los coeficientes de regresión (y puede inflar el error típico de las estimaciones). Las variables de dormitorios, tamaño de la casa y número de baños están correlacionadas. Esto se ilustra con el siguiente ejemplo en R, que ajusta otra regresión eliminando de la ecuación las variables SqFtTotLiving, SqFtFinBasement y Bathrooms:

update(step_lm, . ~ . - SqFtTotLiving - SqFtFinBasement - Bathrooms)

Call:
lm(formula = AdjSalePrice ~ Bedrooms + BldgGrade + PropertyType +
    YrBuilt, data = house, na.action = na.omit)

Coefficients:
              (Intercept)                   Bedrooms
                  4913973                      27151
                BldgGrade  PropertyTypeSingle Family
                   248998                     -19898
    PropertyTypeTownhouse                    YrBuilt
                   -47355                      -3212

La función update puede utilizarse para añadir o eliminar variables de un modelo. Ahora el coeficiente de los dormitorios es positivo, en línea con lo que cabría esperar (aunque en realidad está actuando como sustituto del tamaño de la casa, ahora que se han eliminado esas variables).

En Python no existe un equivalente a la función update de R. Tenemos que volver a ajustar el modelo con la lista de predictores modificada:

predictors = ['Bedrooms', 'BldgGrade', 'PropertyType', 'YrBuilt']
outcome = 'AdjSalePrice'

X = pd.get_dummies(house[predictors], drop_first=True)

reduced_lm = LinearRegression()
reduced_lm.fit(X, house[outcome])

Las variables correlacionadas son sólo un problema a la hora de interpretar los coeficientes de regresión. En house_lm, no hay ninguna variable que tenga en cuenta la ubicación de la vivienda, y el modelo está mezclando tipos de regiones muy diferentes.La ubicación puede ser una variable de confusión; para más información, véase "Variables de confusión".

Multicolinealidad

Un caso extremo de variables correlacionadas produce multicolinealidad, una condición en la que hay redundancia entre las variables predictoras. La multicolinealidad perfecta se produce cuando una variable predictora puede expresarse como una combinación lineal de otras. La multicolinealidad se produce cuando:

• Una variable se incluye varias veces por error.

• Las dummiesP, en lugar de las dummies P - 1, se crean a partir de una variable factorial (véase "Variables factoriales en la regresión").

• Dos variables están casi perfectamente correlacionadas entre sí.

La multicolinealidad en la regresión debe abordarse: deben eliminarse variables hasta que desaparezca la multicolinealidad. Una regresión no tiene una solución bien definida en presencia de multicolinealidad perfecta. Muchos paquetes de software, incluidos R y Python, manejan automáticamente ciertos tipos de multicolinealidad. Por ejemplo, si SqFtTotLiving se incluye dos veces en la regresión de los datos de house, los resultados son los mismos que para el modelo house_lm. En caso de multicolinealidad no perfecta, el software puede obtener una solución, pero los resultados pueden ser inestables.

Nota

La multicolinealidad no es tan problemática en los métodos de regresión no lineal como los árboles, la agrupación y los vecinos más próximos, y en esos métodos puede ser aconsejable mantener P variables ficticias (en lugar de P - 1). Dicho esto, incluso en esos métodos, la no redundancia en las variables predictoras sigue siendo una virtud.

Variables de confusión

Con las variables correlacionadas, el problema es de comisión: incluir distintas variables que tienen una relación predictiva similar con la respuesta.Con las variables de confusión, el problema es de omisión: no se incluye una variable importante en la ecuación de regresión. La interpretación ingenua de los coeficientes de la ecuación puede llevar a conclusiones no válidas.

Tomemos, por ejemplo, la ecuación de regresión del condado de King house_lm de "Ejemplo: Datos de vivienda del condado de King". Los coeficientes de regresión de SqFtLot, Bathrooms y Bedrooms son todos negativos. El modelo de regresión original no contiene una variable que represente la ubicación, un factor de predicción muy importante del precio de la vivienda. Para modelizar la ubicación, incluye una variable ZipGroup que clasifique el código postal en uno de cinco grupos, del menos caro (1) al más caro (5):⁵

lm(formula = AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
    Bedrooms + BldgGrade + PropertyType + ZipGroup, data = house,
    na.action = na.omit)

Coefficients:
              (Intercept)              SqFtTotLiving
               -6.666e+05                  2.106e+02
                  SqFtLot                  Bathrooms
                4.550e-01                  5.928e+03
                 Bedrooms                  BldgGrade
               -4.168e+04                  9.854e+04
PropertyTypeSingle Family      PropertyTypeTownhouse
                1.932e+04                 -7.820e+04
                ZipGroup2                  ZipGroup3
                5.332e+04                  1.163e+05
                ZipGroup4                  ZipGroup5
                1.784e+05                  3.384e+05

El mismo modelo en Python:

predictors = ['SqFtTotLiving', 'SqFtLot', 'Bathrooms', 'Bedrooms',
              'BldgGrade', 'PropertyType', 'ZipGroup']
outcome = 'AdjSalePrice'

X = pd.get_dummies(house[predictors], drop_first=True)

confounding_lm = LinearRegression()
confounding_lm.fit(X, house[outcome])

print(f'Intercept: {confounding_lm.intercept_:.3f}')
print('Coefficients:')
for name, coef in zip(X.columns, confounding_lm.coef_):
    print(f' {name}: {coef}')

ZipGroup es claramente una variable importante: se estima que una vivienda del grupo de códigos postales más caros tiene un precio de venta superior en casi 340.000 $. Los coeficientes de SqFtLot y Bathrooms son ahora positivos, y añadir un cuarto de baño aumenta el precio de venta en 5.928 $.

El coeficiente de Bedrooms sigue siendo negativo. Aunque esto no es intuitivo, se trata de un fenómeno bien conocido en el sector inmobiliario. Para viviendas de la misma superficie habitable y número de baños, tener más dormitorios y, por tanto, más pequeños, se asocia a viviendas menos valiosas.

Interacciones y efectos principales

A los estadísticos les gusta distinguir entre efectos principales, o variables independientes, y las interacciones entre los efectos principales. Los efectos principales son lo que a menudo se denominan variables predic toras en la ecuación de regresión. Una suposición implícita cuando sólo se utilizan efectos principales en un modelo es que la relación entre una variable predictora y la respuesta es independiente de las demás variables predictoras. A menudo no es así.

Por ejemplo, el modelo ajustado a los datos de vivienda del condado de King en "Variables de confusión" incluye varias variables como efectos principales, entre ellas ZipCode. La ubicación en el sector inmobiliario lo es todo, y es natural suponer que la relación entre, por ejemplo, el tamaño de la casa y el precio de venta depende de la ubicación. Una casa grande construida en un barrio de renta baja no va a conservar el mismo valor que una casa grande construida en una zona cara. Las interacciones entre variables se incluyen en R utilizando el operador *. Para los datos del condado de King, a continuación se ajusta una interacción entre SqFtTotLiving y ZipGroup:

lm(formula = AdjSalePrice ~ SqFtTotLiving * ZipGroup + SqFtLot +
    Bathrooms + Bedrooms + BldgGrade + PropertyType, data = house,
    na.action = na.omit)

Coefficients:
              (Intercept)              SqFtTotLiving
               -4.853e+05                  1.148e+02
                ZipGroup2                  ZipGroup3
               -1.113e+04                  2.032e+04
                ZipGroup4                  ZipGroup5
                2.050e+04                 -1.499e+05
                  SqFtLot                  Bathrooms
                6.869e-01                 -3.619e+03
                 Bedrooms                  BldgGrade
               -4.180e+04                  1.047e+05
PropertyTypeSingle Family      PropertyTypeTownhouse
                1.357e+04                 -5.884e+04
  SqFtTotLiving:ZipGroup2    SqFtTotLiving:ZipGroup3
                3.260e+01                  4.178e+01
  SqFtTotLiving:ZipGroup4    SqFtTotLiving:ZipGroup5
                6.934e+01                  2.267e+02

El modelo resultante tiene cuatro nuevos términos:SqFtTotLiving:ZipGroup2, SqFtTotLiving:ZipGroup3, etc.

En Python, necesitamos utilizar el paquete statsmodels para entrenar modelos de regresión lineal con interacciones. Este paquete se diseñó de forma similar a R y permite definir modelos utilizando una interfaz de fórmulas:

model = smf.ols(formula='AdjSalePrice ~ SqFtTotLiving*ZipGroup + SqFtLot + ' +
     'Bathrooms + Bedrooms + BldgGrade + PropertyType', data=house)
results = model.fit()
results.summary()

El paquete statsmodels se ocupa de las variables categóricas (por ejemplo, ZipGroup[T.1], PropertyType[T.Single Family]) y de los términos de interacción (por ejemplo, SqFtTotLiving:ZipGroup[T.1]).

La ubicación y el tamaño de la casa parecen tener una fuerte interacción. Para una casa en el ZipGroup más bajo, la pendiente es la misma que la pendiente del efecto principal SqFtTotLiving, que es de 118 $ por pie cuadrado (esto se debe a que R utiliza la codificación de referencia para las variables factoriales; véase "Variables factoriales en la regresión"). En el caso de una vivienda en el código postal más caro ZipGroup, la pendiente es la suma del efecto principal más SqFtTotLiving:ZipGroup5, o sea 115 $ + 227 $ = 342 $ por pie cuadrado. En otras palabras, añadir un pie cuadrado en el grupo de códigos postales más caros aumenta el precio de venta previsto en un factor de casi tres, en comparación con el aumento medio de añadir un pie cuadrado.

Selección del modelo con términos de interacción

En los problemas en los que intervienen muchas variables, puede resultar difícil decidir qué términos de interacción deben incluirse en el modelo. Se suelen adoptar varios enfoques diferentes:

• En algunos problemas, el conocimiento previo y la intuición pueden guiar la elección de qué términos de interacción incluir en el modelo.

• La selección por pasos (véase "Selección de modelos y regresión por pasos") puede utilizarse para cribar los distintos modelos.

• La regresión penalizada puede ajustarse automáticamente a un gran conjunto de posibles términos de interacción.

• Quizá el enfoque más habitual sea utilizar modelos de árbol, así como sus descendientes, los árboles de bosque aleatorio y los árboles potenciados por gradiente. Esta clase de modelos busca automáticamente los términos de interacción óptimos; véase "Modelos de árbol".

Ideas clave

• Debido a la correlación entre los predictores, hay que tener cuidado en la interpretación de los coeficientes en la regresión lineal múltiple.

• La multicolinealidad puede provocar inestabilidad numérica en el ajuste de la ecuación de regresión.

• Una variable de confusión es un predictor importante que se omite en un modelo y puede dar lugar a una ecuación de regresión con relaciones espurias.

• Se necesita un término de interacción entre dos variables si la relación entre las variables y la respuesta es interdependiente.

Diagnóstico de regresión

En la modelización explicativa (es decir, en un contexto de investigación), se siguen varios pasos, además de las métricas mencionadas anteriormente (véase "Evaluación del modelo"), para evaluar lo bien que el modelo se ajusta a los datos; la mayoría se basan en el análisis de los residuos.Estos pasos no abordan directamente la precisión predictiva, pero pueden proporcionar una perspectiva útil en un entorno predictivo.

Términos clave para el diagnóstico de regresión

Residuos normalizados

Residuos divididos por el error típico de los residuos.

Valores atípicos

Registros (o valores de resultado) que se alejan del resto de los datos (o del resultado previsto).

Valor de influencia

Un valor o registro cuya presencia o ausencia marca una gran diferencia en la ecuación de regresión.

Aprovecha

El grado de influencia que tiene un único registro en una ecuación de regresión.

Sinónimo

valor del sombrero

Residuos no normales

Los residuos con distribución no normal pueden invalidar algunos requisitos técnicos de la regresión, pero no suelen ser motivo de preocupación en la ciencia de datos.

Heteroscedasticidad

Cuando algunos rangos del resultado experimentan residuos con mayor varianza (puede indicar que falta un predictor en la ecuación).

Parcelas de residuos parciales

Un gráfico de diagnóstico para iluminar la relación entre la variable de resultado y un único predictor.

Sinónimo

diagrama de variables añadidas

Valores atípicos

En términos generales, un valor extremo, también llamado valor atípico, es aquel que está alejado de la mayoría de las demás observaciones.Del mismo modo que hay que tratar los valores atípicos para las estimaciones de localización y variabilidad (ver "Estimaciones de localización" y "Estimaciones de variabilidad"), los valores atípicos pueden causar problemas con los modelos de regresión. En regresión, un valor atípico es un registro cuyo valor y real está alejado del valor predicho. Puedes detectar los valores atípicos examinando el residuo estandarizado, que es el residuo dividido por el error estándar de los residuos.

No existe una teoría estadística que separe los valores atípicos de los que no lo son. Por ejemplo, con el diagrama de caja, los valores atípicos son los puntos de datos que están demasiado por encima o por debajo de los límites de la caja (véase "Percentiles y diagramas de caja"), donde "demasiado" = "más de 1,5 veces el rango intercuartílico".5 veces el rango intercuartílico". En regresión, el residuo estandarizado es la métrica que se suele utilizar para determinar si un registro se clasifica como atípico. El residuo estandarizado puede interpretarse como "el número de errores estándar que se alejan de la recta de regresión".

Vamos a ajustar una regresión a los datos de ventas de viviendas del condado de King para todas las ventas del código postal 98105 en R:

house_98105 <- house[house$ZipCode == 98105,]
lm_98105 <- lm(AdjSalePrice ~ SqFtTotLiving + SqFtLot + Bathrooms +
                 Bedrooms + BldgGrade, data=house_98105)

En Python:

house_98105 = house.loc[house['ZipCode'] == 98105, ]

predictors = ['SqFtTotLiving', 'SqFtLot', 'Bathrooms', 'Bedrooms', 'BldgGrade']
outcome = 'AdjSalePrice'

house_outlier = sm.OLS(house_98105[outcome],
                       house_98105[predictors].assign(const=1))
result_98105 = house_outlier.fit()

Extraemos los residuos estandarizados en R utilizando la función rstandard y obtenemos el índice del residuo más pequeño utilizando la función order:

sresid <- rstandard(lm_98105)
idx <- order(sresid)
sresid[idx[1]]
    20429
-4.326732

En statsmodels, utiliza OLSInfluence para analizar los residuos:

influence = OLSInfluence(result_98105)
sresiduals = influence.resid_studentized_internal
sresiduals.idxmin(), sresiduals.min()

La mayor sobreestimación del modelo está más de cuatro errores estándar por encima de la línea de regresión, lo que corresponde a una sobreestimación de 757.754 $. El registro de datos original correspondiente a este valor atípico es el siguiente en R:

house_98105[idx[1], c('AdjSalePrice', 'SqFtTotLiving', 'SqFtLot',
              'Bathrooms', 'Bedrooms', 'BldgGrade')]

AdjSalePrice SqFtTotLiving SqFtLot Bathrooms Bedrooms BldgGrade
         (dbl)         (int)   (int)     (dbl)    (int)     (int)
20429   119748          2900    7276         3        6         7

En Python:

outlier = house_98105.loc[sresiduals.idxmin(), :]
print('AdjSalePrice', outlier[outcome])
print(outlier[predictors])

En este caso, parece que hay algo erróneo en el registro: una casa de ese tamaño suele venderse por mucho más de 119.748 $ en ese código postal.La Figura 4-4 muestra un extracto de la escritura legal de esta venta: está claro que la venta sólo implicaba un interés parcial en la propiedad. En este caso, el valor atípico corresponde a una venta que es anómala y no debería incluirse en la regresión. Los valores atípicos también podrían ser el resultado de otros problemas, como la introducción de datos con "dedos gordos" o un desajuste de unidades (por ejemplo, declarar una venta en miles de dólares en lugar de simplemente en dólares).

Figura 4-4. Escritura de garantía legal para el mayor residuo negativo

Para los problemas de big data, los valores atípicos no suelen ser un problema a la hora de ajustar la regresión que se utilizará para predecir nuevos datos. Sin embargo, los valores atípicos son fundamentales para la detección de anomalías, en la que encontrar valores atípicos es el objetivo principal. El valor atípico también podría corresponder a un caso de fraude o a una acción accidental. En cualquier caso, detectar valores atípicos puede ser una necesidad empresarial crítica.

Valores influyentes

Un valor cuya ausencia cambiaría significativamente la ecuación de regresión se denomina observación influyente.En regresión, no es necesario que un valor de este tipo esté asociado a un residuo grande. Como ejemplo, considera las líneas de regresión de la Figura 4-5. La línea continua corresponde a la regresión con todos los datos, mientras que la línea discontinua corresponde a la regresión con el punto de la esquina superior derecha eliminado. La línea continua corresponde a la regresión con todos los datos, mientras que la línea discontinua corresponde a la regresión con el punto de la esquina superior derecha eliminado. Es evidente que ese valor de dato tiene una gran influencia en la regresión, aunque no esté asociado a un valor atípico grande (de la regresión completa). Se considera que ese valor de dato tiene una gran influencia en la regresión.

Además de los residuos normalizados (véase "Valores atípicos"), los estadísticos han desarrollado varias métricas para determinar la influencia de un único registro en una regresión.Una medida común de influencia es el valor sombrero; los valores superiores a $2(P+1)/n$ indican un valor de datos de alto apalancamiento.⁶

Figura 4-5. Ejemplo de punto de datos influyente en la regresión

Otra métrica es la distancia de Cook, que define la influencia como una combinación de apalancamiento y tamaño residual.Una regla empírica es que una observación tiene una influencia alta si la distancia de Cook supera $4/(n-P-1)$.

Un gráfico de influencia o gráfico de burbujas combina los residuos estandarizados, el valor sombrero y la distancia de Cook en un único gráfico.La Figura 4-6 muestra el gráfico de influencia para los datos de viviendas del condado de King y puede crearse mediante el siguiente código R:

std_resid <- rstandard(lm_98105)
cooks_D <- cooks.distance(lm_98105)
hat_values <- hatvalues(lm_98105)
plot(subset(hat_values, cooks_D > 0.08), subset(std_resid, cooks_D > 0.08),
     xlab='hat_values', ylab='std_resid',
     cex=10*sqrt(subset(cooks_D, cooks_D > 0.08)), pch=16, col='lightgrey')
points(hat_values, std_resid, cex=10*sqrt(cooks_D))
abline(h=c(-2.5, 2.5), lty=2)

Aquí tienes el código Python para crear una figura similar:

influence = OLSInfluence(result_98105)
fig, ax = plt.subplots(figsize=(5, 5))
ax.axhline(-2.5, linestyle='--', color='C1')
ax.axhline(2.5, linestyle='--', color='C1')
ax.scatter(influence.hat_matrix_diag, influence.resid_studentized_internal,
           s=1000 * np.sqrt(influence.cooks_distance[0]),
           alpha=0.5)
ax.set_xlabel('hat values')
ax.set_ylabel('studentized residuals')

Aparentemente, hay varios puntos de datos que muestran una gran influencia en la regresión. La distancia de Cook puede calcularse utilizando la función cooks.distance, y puedes utilizar hatvalues para calcular los diagnósticos. Los valores del sombrero se representan en el eje x, los residuos en el eje y, y el tamaño de los puntos está relacionado con el valor de la distancia de Cook.

Figura 4-6. Un gráfico para determinar qué observaciones tienen una gran influencia; los puntos con una distancia de Cook superior a 0,08 están resaltados en gris

La Tabla 4-2 compara la regresión con el conjunto de datos completo y con los puntos de datos muy influyentes eliminados (distancia de Cook > 0,08).

El coeficiente de regresión de para Bathrooms cambia bastante.⁷

Tabla 4-2. Comparación de los coeficientes de regresión con los datos completos y con los datos influyentes eliminados

	Original	Influyente eliminado
(Interceptar)	-772,550	-647,137
SqFtTotVivir	210	230
Parcela	39	33
Baños	2282	-16,132
Dormitorios	-26,320	-22,888
EdificioGrado	130,000	114,871

Para ajustar una regresión que prediga de forma fiable los datos futuros, identificar las observaciones influyentes sólo es útil en conjuntos de datos pequeños. En las regresiones que incluyen muchos registros, es poco probable que una sola observación tenga el peso suficiente para causar una influencia extrema en la ecuación ajustada (aunque la regresión puede seguir teniendo grandes valores atípicos). Sin embargo, para detectar anomalías, identificar las observaciones influyentes puede ser muy útil.

Heteroscedasticidad, no normalidad y errores correlacionados

Los estadísticos prestan mucha atención a la distribución de los residuos.Resulta que los mínimos cuadrados ordinarios (véase "Mínimos cuadrados") son insesgados, y en algunos casos son el estimador "óptimo", bajo una amplia gama de supuestos distribucionales. Esto significa que, en la mayoría de los problemas, los científicos de datos no necesitan preocuparse demasiado por la distribución de los residuos.

La distribución de los residuos es relevante principalmente para la validez de la inferencia estadística formal (pruebas de hipótesis y valores p), que tiene una importancia mínima para los científicos de datos preocupados principalmente por la precisión predictiva. Los errores distribuidos normalmente son una señal de que el modelo está completo; los errores que no se distribuyen normalmente indican que al modelo le falta algo. Para que la inferencia formal sea totalmente válida, se supone que los residuos se distribuyen normalmente, tienen la misma varianza y son independientes. Un área en la que esto puede preocupar a los científicos de datos es el cálculo estándar de los intervalos de confianza para los valores predichos, que se basan en las suposiciones sobre los residuos (véase "Intervalos de confianza y predicción").

La heteroscedasticidad es la falta de varianza residual constante en todo el intervalo de los valores predichos.En otras palabras, los errores son mayores en algunas partes del intervalo que en otras. Visualizar los datos es una forma cómoda de analizar los residuos.

El siguiente código en R traza los residuos absolutos frente a los valores predichos para el ajuste de regresión lm_98105 en "Valores atípicos":

df <- data.frame(resid = residuals(lm_98105), pred = predict(lm_98105))
ggplot(df, aes(pred, abs(resid))) + geom_point() + geom_smooth()

La Figura 4-7 muestra el gráfico resultante. Utilizando geom_smooth, es fácil superponer un suavizado de los residuos absolutos. La función llama al método loess (suavizado de dispersión estimado localmente) para producir una estimación suavizada de la relación entre las variables del eje x y del eje y en un gráfico de dispersión (consulta "Suavizadores de dispersión").

En Python, el paquete seaborn tiene la función regplot para crear una figura similar:

fig, ax = plt.subplots(figsize=(5, 5))
sns.regplot(result_98105.fittedvalues, np.abs(result_98105.resid),
            scatter_kws={'alpha': 0.25}, line_kws={'color': 'C1'},
            lowess=True, ax=ax)
ax.set_xlabel('predicted')
ax.set_ylabel('abs(residual)')

Figura 4-7. Gráfico del valor absoluto de los residuos frente a los valores predichos

Evidentemente, la varianza de los residuos tiende a aumentar para las viviendas de mayor valor, pero también es grande para las de menor valor.Este gráfico indica que lm_98105 tiene errores heteroscedásticos.

¿Por qué un científico de datos se preocupa por la heteroscedasticidad?

La heteroscedasticidad indica que los errores de predicción difieren para distintos rangos del valor predicho, y puede sugerir un modelo incompleto. Por ejemplo, la heteroscedasticidad en lm_98105 puede indicar que la regresión ha dejado algo sin contabilizar en los hogares de rango alto y bajo.

La Figura 4-8 es un histograma de los residuos estandarizados de la regresión lm_98105.La distribución tiene colas decididamente más largas que la distribución normal y muestra una leve asimetría hacia los residuos más grandes.

Figura 4-8. Histograma de los residuos de la regresión de los datos de vivienda

Los estadísticos también pueden comprobar el supuesto de que los errores son independientes. Esto es especialmente cierto para los datos que se recogen a lo largo del tiempo o del espacio. El estadístico de Durbin-Watson puede utilizarse para detectar si existe una autocorrelación significativa en una regresión que incluya datos de series temporales. Si los errores de un modelo de regresión están correlacionados, esta información puede ser útil para hacer previsiones a corto plazo y debe incorporarse al modelo. Consulta Practical Time Series Forecasting with R, 2ª ed., de Galit Shmueli y Kenneth Lichtendahl (Axelrod Schnall, 2018) para saber más sobre cómo incorporar información de autocorrelación en modelos de regresión para datos de series temporales. Si el objetivo son las previsiones a largo plazo o los modelos explicativos, el exceso de datos autocorrelacionados a nivel micro puede distraer. En ese caso, puede ser conveniente suavizar los datos o recopilarlos de forma menos granular.

Aunque una regresión pueda incumplir uno de los supuestos de distribución, ¿debería importarnos? La mayoría de las veces, en la ciencia de datos, el interés se centra principalmente en la precisión predictiva, por lo que puede ser conveniente revisar un poco la heteroscedasticidad. Puede que descubras que hay alguna señal en los datos que tu modelo no ha captado. Sin embargo, satisfacer los supuestos de distribución simplemente para validar la inferencia estadística formal (valores p, estadísticos F, etc.) no es tan importante para el científico de datos.

Suavizadores de dispersión

La regresión consiste en modelizar la relación entre la respuesta y las variables predictoras.Al evaluar un modelo de regresión, es útil utilizar un suavizador de dispersión para resaltar visualmente las relaciones entre dos variables.

Por ejemplo, en la Figura 4-7, un suavizado de la relación entre los residuos absolutos y el valor predicho muestra que la varianza de los residuos depende del valor del residuo. En este caso, se utilizó la función loess; loess funciona ajustando repetidamente una serie de regresiones locales a subconjuntos contiguos para obtener un suavizado. Aunque loess es probablemente el suavizador más utilizado, en R hay disponibles otros suavizadores de gráficos de dispersión, como super smooth (supsmu) y kernel smoothing (ksmooth). En Python, podemos encontrar otros suavizadores en scipy (wiener o sav) y statsmodels (kernel_regression). Para evaluar un modelo de regresión, no suele ser necesario preocuparse por los detalles de estos suavizadores de dispersión.

Gráficos de residuos parciales y no linealidad

Los gráficos de residuos parciales son una forma de visualizar hasta qué punto el ajuste estimado explica la relación entre un predictor y el resultado.La idea básica de un gráfico de residuos parciales es aislar la relación entre una variable predictora y la respuesta, teniendo en cuenta todas las demás variables predictoras.Un residuo parcial podría considerarse como un valor de "resultado sintético", que combina la predicción basada en un único predictor con el residuo real de la ecuación de regresión completa. Un residuo parcial para el predictor $X_i$ es el residuo ordinario más el término de regresión asociado a $X_i$:

$$\text{Parcial}\text{residual}=\text{Residual}+{\widehat{b}}_i{X}_i$$

donde ${\widehat{b}}_i$ es el coeficiente de regresión estimado. La función predict de R tiene una opción para devolver los términos de regresión individuales ${\widehat{b}}_i{X}_i$:

terms <- predict(lm_98105, type='terms')
partial_resid <- resid(lm_98105) + terms

El gráfico de residuos parciales muestra las $X_i$ en el eje de abscisas y los residuos parciales en el eje de ordenadas. Utilizar ggplot2 facilita la superposición de un alisado de los residuos parciales:

df <- data.frame(SqFtTotLiving = house_98105[, 'SqFtTotLiving'],
                 Terms = terms[, 'SqFtTotLiving'],
                 PartialResid = partial_resid[, 'SqFtTotLiving'])
ggplot(df, aes(SqFtTotLiving, PartialResid)) +
  geom_point(shape=1) + scale_shape(solid = FALSE) +
  geom_smooth(linetype=2) +
  geom_line(aes(SqFtTotLiving, Terms))

El paquete statsmodels tiene el método sm.graphics.plot_ccpr que crea un gráfico de residuos parciales similar:

sm.graphics.plot_ccpr(result_98105, 'SqFtTotLiving')

Los gráficos de R y Python se diferencian por un desplazamiento constante. En R, se añade una constante para que la media de los términos sea cero.

El gráfico resultante se muestra en la Figura 4-9. El residuo parcial es una estimación de la contribución que SqFtTotLiving añade al precio de venta. La relación entre SqFtTotLiving y el precio de venta es evidentemente no lineal (línea discontinua). La línea de regresión (línea continua) subestima el precio de venta de las viviendas de menos de 1.000 pies cuadrados y sobreestima el precio de las viviendas de entre 2.000 y 3.000 pies cuadrados. Hay muy pocos puntos de datos por encima de 4.000 pies cuadrados para sacar conclusiones sobre esas viviendas.

Figura 4-9. Diagrama de residuos parciales de la variable SqFtTotLiving

Esta no linealidad de tiene sentido en este caso: añadir 500 pies en una casa pequeña supone una diferencia mucho mayor que añadir 500 pies en una casa grande. Esto sugiere que, en lugar de un simple término lineal para SqFtTotLiving, debería considerarse un término no lineal (ver "Regresión polinómica y Spline").

Ideas clave

• Aunque los valores atípicos pueden causar problemas en conjuntos de datos pequeños, el interés principal de los valores atípicos es identificar problemas en los datos, o localizar anomalías.

• Los registros individuales (incluidos los valores atípicos de la regresión) pueden tener una gran influencia en una ecuación de regresión con datos pequeños, pero este efecto se desvanece en los datos grandes.

• Si el modelo de regresión se utiliza para realizar inferencias formales (valores p y similares), habrá que comprobar ciertos supuestos sobre la distribución de los residuos. En general, sin embargo, la distribución de los residuos no es crítica en la ciencia de datos.

• El gráfico de residuos parciales puede utilizarse para evaluar cualitativamente el ajuste de cada término de regresión, lo que posiblemente conduzca a una especificación alternativa del modelo.

Regresión polinómica y por splines

La relación entre la respuesta y una variable predictora no es necesariamente lineal.La respuesta a la dosis de un fármaco suele ser no lineal: duplicar la dosis no suele conducir a duplicar la respuesta.La demanda de un producto no es una función lineal de los dólares gastados en marketing; en algún momento, es probable que la demanda se sature. Hay muchas formas de ampliar la regresión para captar estos efectos no lineales.

Términos clave de la regresión no lineal

Regresión polinómica

Añade términos polinómicos (cuadrados, cubos, etc.) a una regresión.

Regresión Spline

Ajustar una curva suave con una serie de segmentos polinómicos.

Nudos

Valores que separan los segmentos de la spline.

Modelos aditivos generalizados

Modelos spline con selección automatizada de nudos.

Sinónimo

GAM

Regresión no lineal

Cuando los estadísticos hablan de regresión no lineal, se refieren a los modelos que no pueden ajustarse mediante mínimos cuadrados. ¿Qué tipo de modelos son no lineales?Esencialmente, todos los modelos en los que la respuesta no puede expresarse como una combinación lineal de los predictores o alguna transformación de los predictores. Los modelos de regresión no lineal son más difíciles y computacionalmente más intensivos de ajustar, ya que requieren una optimización numérica. Por este motivo, generalmente se prefiere utilizar un modelo lineal si es posible.

Polinomio

La regresión polinómica consiste en incluir términos polinómicos en una ecuación de regresión.El uso de la regresión polinómica se remonta casi al desarrollo de la propia regresión, con un artículo de Gergonne en 1815. Por ejemplo, una regresión cuadrática entre la respuesta Y y el predictor X tendría la forma:

$$Y=b_0+b_{1}X+b_2{X}^2+e$$

La regresión polinómica puede ajustarse en R mediante la función poly. Por ejemplo, a continuación se ajusta un polinomio cuadrático para SqFtTotLiving con los datos de vivienda del condado de King:

lm(AdjSalePrice ~  poly(SqFtTotLiving, 2) + SqFtLot +
                BldgGrade + Bathrooms + Bedrooms,
                    data=house_98105)

Call:
lm(formula = AdjSalePrice ~ poly(SqFtTotLiving, 2) + SqFtLot +
   BldgGrade + Bathrooms + Bedrooms, data = house_98105)

Coefficients:
           (Intercept)  poly(SqFtTotLiving, 2)1  poly(SqFtTotLiving, 2)2
            -402530.47               3271519.49                776934.02
               SqFtLot                BldgGrade                Bathrooms
                 32.56                135717.06                 -1435.12
              Bedrooms
              -9191.94

En statsmodels, añadimos el término cuadrado a la definición del modelo utilizando I(SqFtTotLiving**2):

model_poly = smf.ols(formula='AdjSalePrice ~  SqFtTotLiving + ' +
                '+ I(SqFtTotLiving**2) + ' +
                'SqFtLot + Bathrooms + Bedrooms + BldgGrade', data=house_98105)
result_poly = model_poly.fit()
result_poly.summary()  

El intercepto y los coeficientes polinómicos son diferentes en comparación con R. Esto se debe a las diferentes implementaciones. Los demás coeficientes y las predicciones son equivalentes.

Ahora hay dos coeficientes asociados a SqFtTotLiving: uno para el término lineal y otro para el término cuadrático.

El gráfico de residuos parciales (ver " Gráficos de residuos parciales y no linealidad") indica cierta curvatura en la ecuación de regresión asociada a SqFtTotLiving. La línea ajustada se ajusta más a la suavidad (ver "Splines") de los residuos parciales, en comparación con un ajuste lineal (ver Figura 4-10).

La implementación de statsmodels sólo funciona para términos lineales. El código fuente adjunto ofrece una implementación que también funcionará para la regresión polinómica.

Figura 4-10. Un ajuste de regresión polinómica para la variable SqFtTotLiving (línea continua) frente a una línea suave (línea discontinua; véase la sección siguiente sobre splines)

Estrías

La regresión polinómica sólo capta una cierta cantidad de curvatura en una relación no lineal.Añadir términos de orden superior, como un polinomio cúbico cuártico, a menudo conduce a una "ondulación" indeseable en la ecuación de regresión. Un enfoque alternativo, y a menudo superior, para modelar relaciones no lineales es utilizar splines.Los splines proporcionan una forma de interpolar suavemente entre puntos fijos. Los splines fueron utilizados originalmente por los dibujantes para dibujar una curva suave, sobre todo en la construcción naval y aeronáutica.

Las estrías se crearon doblando una fina pieza de madera mediante pesos, denominados "patos"; véase la Figura 4-11.

Figura 4-11. Las splines se crearon originalmente utilizando madera flexible y "patos" y se utilizaban como herramienta de delineante para ajustar curvas (foto cortesía de Bob Perry)

La definición técnica de spline es una serie de polinomios continuos a trozos. Fueron desarrollados por primera vez durante la Segunda Guerra Mundial en el campo de pruebas estadounidense de Aberdeen por I. J. Schoenberg, un matemático rumano. Las piezas polinómicas se conectan suavemente en una serie de puntos fijos de una variable predictora, denominados nudos. La formulación de splines es mucho más complicada que la regresión polinómica; el software estadístico suele encargarse de los detalles del ajuste de un spline. El paquete R splines incluye la función bs para crear un término b-spline (spline de base) en un modelo de regresión.Por ejemplo, lo siguiente añade un término b-spline al modelo de regresión de la casa:

library(splines)
knots <- quantile(house_98105$SqFtTotLiving, p=c(.25, .5, .75))
lm_spline <- lm(AdjSalePrice ~ bs(SqFtTotLiving, knots=knots, degree=3) +
  SqFtLot + Bathrooms + Bedrooms + BldgGrade,  data=house_98105)

Hay que especificar dos parámetros: el grado del polinomio y la ubicación de los nodos. En este caso, el predictor SqFtTotLiving se incluye en el modelo mediante un spline cúbico (degree=3). Por defecto, bs coloca los nodos en los límites; además, también se colocaron nodos en el cuartil inferior, el cuartil medio y el cuartil superior.

La interfaz de fórmulas statsmodels admite el uso de splines de forma similar a R. Aquí, especificamos la b-spline utilizando df, los grados de libertad. Esto creará df - degree = 6 - 3 = 3 nudos internos con posiciones calculadas del mismo modo que en el código R anterior:

formula = 'AdjSalePrice ~ bs(SqFtTotLiving, df=6, degree=3) + ' +
          'SqFtLot + Bathrooms + Bedrooms + BldgGrade'
model_spline = smf.ols(formula=formula, data=house_98105)
result_spline = model_spline.fit()

A diferencia de un término lineal, para el que el coeficiente tiene un significado directo, los coeficientes de un término spline no son interpretables. En su lugar, es más útil utilizar la representación visual para revelar la naturaleza del ajuste spline.La Figura 4-12 muestra el gráfico de residuos parciales de la regresión. En contraste con el modelo polinómico, el modelo spline se ajusta más al liso, lo que demuestra la mayor flexibilidad de los splines. En este caso, la línea se ajusta más a los datos. ¿Significa esto que la regresión spline es un modelo mejor? No necesariamente: no tiene sentido económico que las viviendas muy pequeñas (menos de 1.000 pies cuadrados) tengan un valor más alto que las viviendas ligeramente más grandes. Esto es posiblemente un artefacto de una variable de confusión; véase "Variables de confusión".

Figura 4-12. Un ajuste de regresión spline para la variable SqFtTotLiving (línea continua) comparado con un ajuste suave (línea discontinua)

Modelos aditivos generalizados

Supongamos que sospechas que existe una relación no lineal entre la respuesta y una variable predictora, ya sea por conocimiento a priori o examinando los diagnósticos de regresión. Los términos polinómicos pueden no ser lo bastante flexibles para captar la relación, y los términos spline requieren especificar los nudos.Los modelos aditivos generalizados, o GAM, son una técnica de modelización flexible que puede utilizarse para ajustar automáticamente una regresión spline. El paquete mgcv de R puede utilizarse para ajustar un modelo GAM a los datos de la vivienda:

library(mgcv)
lm_gam <- gam(AdjSalePrice ~ s(SqFtTotLiving) + SqFtLot +
                    Bathrooms +  Bedrooms + BldgGrade,
                    data=house_98105)

El término s(SqFtTotLiving) indica a la función gam que busque los "mejores" nudos para un término spline (ver Figura 4-13).

Figura 4-13. Ajuste de una regresión GAM para la variable SqFtTotLiving (línea continua) comparada con una suave (línea discontinua)

En Python, podemos utilizar el paquete pyGAM. Proporciona métodos para la regresión y la clasificación. Aquí utilizamos LinearGAM para crear un modelo de regresión:

predictors = ['SqFtTotLiving', 'SqFtLot', 'Bathrooms',  'Bedrooms', 'BldgGrade']
outcome = 'AdjSalePrice'
X = house_98105[predictors].values
y = house_98105[outcome]

gam = LinearGAM(s(0, n_splines=12) + l(1) + l(2) + l(3) + l(4))  
gam.gridsearch(X, y)

El valor por defecto de n_splines es 20. Esto conduce a un ajuste excesivo para valores mayores de SqFtTotLiving. Un valor de 12 conduce a un ajuste más razonable.

Ideas clave

• Los valores atípicos en una regresión son registros con un residuo grande.

• La multicolinealidad puede provocar inestabilidad numérica en el ajuste de la ecuación de regresión.

• Una variable de confusión es un predictor importante que se omite en un modelo y puede dar lugar a una ecuación de regresión con relaciones espurias.

• Se necesita un término de interacción entre dos variables si el efecto de una variable depende del nivel o la magnitud de la otra.

• La regresión polinómica puede ajustar relaciones no lineales entre los predictores y la variable de resultado.

• Las splines son series de segmentos polinómicos encadenados, que se unen en nudos.

• Podemos automatizar el proceso de especificar los nudos en splines utilizando modelos aditivos generalizados (GAM).

Otras lecturas

• Para más información sobre los modelos spline y los GAM, consulta The Elements of Statistical Learning, 2ª ed., de Trevor Hastie, Robert Tibshirani y Jerome Friedman (2009), y su primo más breve basado en R, An Introduction to Statistical Learning, de Gareth James, Daniela Witten, Trevor Hastie y Robert Tibshirani (2013); ambos son libros de Springer.

• Para saber más sobre el uso de modelos de regresión para la previsión de series temporales, consulta Practical Time Series Forecasting with R, de Galit Shmueli y Kenneth Lichtendahl (Axelrod Schnall, 2018).

Resumen

Tal vez ningún otro método estadístico se haya utilizado más a lo largo de los años que la regresión: el proceso de establecer una relación entre múltiples variables predictoras y una variable de resultado. La forma fundamental es lineal: cada variable predictora tiene un coeficiente que describe una relación lineal entre el predictor y el resultado. Las formas más avanzadas de regresión, como la regresión polinómica y la spline, permiten que la relación no sea lineal. En la estadística clásica, el énfasis se pone en encontrar un buen ajuste a los datos observados para explicar o describir algún fenómeno, y la fuerza de este ajuste es la forma en que se utilizan las métricas tradicionales en la muestra para evaluar el modelo. En cambio, en la ciencia de datos, el objetivo suele ser predecir valores para nuevos datos, por lo que se utilizan métricas basadas en la precisión predictiva para datos fuera de muestra. Los métodos de selección de variables se utilizan para reducir la dimensionalidad y crear modelos más compactos.

¹ Esta sección y las siguientes de este capítulo © 2020 Datastats, LLC, Peter Bruce, Andrew Bruce y Peter Gedeck; utilizadas con permiso.

² En la estadística bayesiana, se supone que el valor verdadero es una variable aleatoria con una distribución determinada. En el contexto bayesiano, en lugar de estimaciones de parámetros desconocidos, existen distribuciones a posteriori y a priori.

³ El argumento -1 en model.matrix produce una representación de codificación en caliente (eliminando el intercepto, de ahí el "-"). De lo contrario, el valor por defecto en R es producir una matriz con P - 1 columnas con el primer nivel de factor como referencia.

⁴ Esto es poco intuitivo, pero puede explicarse por el impacto de la ubicación como variable de confusión; véase "Variables de confusión".

⁵ Hay 80 códigos postales en el condado de King, varios de ellos con sólo un puñado de ventas. Como alternativa al uso directo del código postal como variable factorial, ZipGroup agrupa códigos postales similares en un único grupo. Para más detalles, consulta "Variables factoriales con muchos niveles".

⁶ El término valor sombrero procede de la noción de matriz sombrero en regresión. La regresión lineal múltiple puede expresarse mediante la fórmula $\widehat{Y}=HY$ donde $H$ es la matriz sombrero. Los valores del sombrero corresponden a la diagonal de $H$.

⁷ El coeficiente de Bathrooms se vuelve negativo, lo que es poco intuitivo. No se ha tenido en cuenta la ubicación, y el código postal 98105 contiene zonas de tipos de viviendas dispares. Véase "Variables de confusión " para un análisis de las variables de confusión.