Análisis de Regresión Lineal: Guía Completa
Análisis de Regresión Lineal: Guía Completa
Contents
1 Introducción 1
7 Ejercicios 25
1 Introducción
El objetivo del análisis de regresión lineal es desarrollar o establecer un modelo lineal que explica el
comportamiento de una variable (la cual llamaremos dependiente o explicada) a través de un grupo de una o
más variables independientes (o explicativas), de esta manera la regresión lineal se refiere a una ecuación
lineal que indica la relación que existe entre la variable dependiente con las independientes.
La regresión lineal se refiere a una función lineal con la siguiente forma:
y = β0 + β1 x 1 + u (1)
y = β0 + β1 x1 + β2 x2 + u (2)
Ambas ecuaciones, son funciones lineales que describen el comportamiento o la relación que existe de x1 sobre
y en la primera ecuación y la relación existente entre x1 y x2 sobre y en la segunda ecuación. y es la variable
dependiente y x1 y x2 las variables independientes.
El problema consiste en tomar los datos observados: y, x1 y x2 para estimar la mejor función lineal que
describa el comportamiento de la variable dependiente. Nota que ambas ecuaciones incluyen un término u,
este es un error estócastico, el cual surge debido a que los modelos de regresión no son modelos deterministas,
sino que más bien es probabilístico y es necesario establecer algunos supuestos sobre el comportamiento de
este término de error.
1
En caso donde existe solo una variable independiente, se llama regresión lineal simple, y si tenemos más una
variable independiente se le llama regresión múltiple, en todos los casos siempre tendremos solo una variable
independiente. Si queremos analizar diferentes variables como función de otras, entonces debemos desarrollar
una regresión lineal para cada una de las variables independientes.
Como parte de la notación o terminología que usaremos en este apartado, cuando hablamos de regresión,
ecuación de regresión, regresión lineal o modelo de regresión, nos referimos al mismo concepto, es decir, la
ecuación que describe el comportamiento lineal de una variable dependiente como función de una o más
variables independientes, como la ecuación 1 y 2.
Las variables observadas en una regresión linean son la variable dependiente y las independientes, mientras
que los parámetros β así como el termino de error no son observables, el problema de consiste en estimar
estos coeficientes los cuales son los que determinan la relación que existe entre las variables.
Considera el siguiente ejemplo. Las ventas de una papelería en una colonia dependen de la población que vive
en esa colonia, entre más personas vivan en esa zona mayor será la demanda de sus productos. La persona
encargada de la tienda le interesa estimar el efecto que tendría en promedio sobre sus ventas si se aumenta (o
disminuye) la población de la colonia y propuso el siguiente modelo de regresión:
ventas = β0 + β1 población + u
Ventas, es una variable que se refiere al monto en pesos de las ventas realizadas y población es una variable
que indica el número de personas que viven en la colonia, β0 y β1 , son los coeficientes que deben ser estimados.
Tomando datos de las ventas de su papelería y la población por año, al estimar la regresión obtuvo el siguiente
resultado:
Los resultados de esta regresión indican que β0 = 100 y β1 = 0.5. El signo de estos coeficientes determina la
relación entre la población y las ventas; ya que estos resultados nos indican que sí la población en la colonia
aumenta en una persona las ventas aumentan en 0.5 pesos, dado que β1 es positivo, la relación entre la
población y las ventas es positiva, si la población aumenta, también las ventas, si la población disminuye, las
ventas disminuyen, ¿en cuánto? en 0.5 pesos por persona.
Nota que β0 es también positivo e igual a 100, en una regresión este coeficiente representa una constante,
ya que no depende de ninguna variable, en ese contexto indica que si no hay población en esta colonia,
la papelería espera tener unas ventas mínimas de 100 pesos, ¿por qué? es muy probable que aunque las
principales ventas sean por la población que viven en la misma colonia, una parte de las ventas se debe
a compras que realizan personas que van de paso, de tal forma que la constante, representa estas ventas.
Dependiendo del problema que se esté analizando, el coeficiente que representa la constante en una regresión
lineal puede o no tener un significado especifico.
El signo de los coeficientes asociados a un variable independiente determinan la relación entre las variables:
Un coeficiente positivo indica una relación positiva (sí x aumenta y también, si x disminuye y también),
mientras que un coeficiente negativo indica una relación negativa o inversa (si x aumenta y disminuye, si x
disminuye y aumenta). La magnitud del impacto lo determina el valor del coeficiente, de hecho, nos interesa
poder estimar un modelo cuyos coeficientes sean estadísticamente significativos, de lo contrario no serían
útiles para explicar el comportamiento de la variable dependiente.
A manera de resumen, el análisis de regresión lo que nos indica es el efecto promedio esperado de una variable
independiente sobre una variable dependiente, lo cual nos permite identificar la dependencia estadística que
existe entre las variables, sin embargo, al ser un método estadístico, se necesitan cumplir ciertos supuestos
sobre el comportamiento de las variables y los residuales para que los resultados sean válidos.
En este capítulo resumiremos los conceptos y herramientas necesarias para realizar un análisis de regresión
lineal básico, considera que detrás de muchas explicaciones se encuentran detalles técnicos, y algunos pueden
llegar a ser complicados, por lo que no estaremos desarrollando detalladamente los aspectos matemáticos
2
detrás del análisis de regresión, a favor de mostrar cómo realizar este análisis usando R, sin embargo, notarás
cómo algunos conceptos como pruebas de hipótesis y distribuciones muestrales se vuelven relevantes para
realizar este ejercicio.
y = β0 + β1 x 1 + u
Es la recta de regresión poblacional, al igual que en los casos cuando se analiza la media o la proporción de
una variable, en muchas ocasiones necesitamos tomar una muestra y utilizar la media muestra o la proporción
muestral como un estimador puntual de los valores poblacionales, en el análisis de regresión sucede algo
similar, dado que no conocemos los coeficientes de la población utilizamos los valores muestral de la variable y
y x para estimar una recta de regresión muestral, es decir, utiliamos los estadísticos βˆ0 y βˆ1 como estimadores
de los coeficientes poblaciones β0 y β1 , y con estos podemos construir la recta de regresión estimada:
ŷ = βˆ0 + βˆ1 x
Donde x es la variable dependiente observada pero ahora ŷ es la variable dependiente estimada, el objetivo
del análisis de regresión lineal es encontrar los mejores estadísticos βˆ0 y βˆ1 de tal manera que la diferencia
entre la variable dependiente estimada y la observada sea la mínima posible, es decir, se busca minimizar la
diferencia y − ŷ.
Esta diferencia es el residual:
û = y − ŷ
û = y − [βˆ0 + βˆ1 x]
Vale la pena añadir un subíndice i a la expresión anterior, ya que esto indica que existe un residual para
cada observación en mi tabla de datos, así como para cada observación corresponde un par y y x. Donde
i = 1, 2, ..., n, y n es el número de observaciones en nuestra muestra:
Sin embargo, necesitamos minimizar en conjunto todos los residuales, para ello simplemente tomamos la
suma de los residuales y añadiendo un operador de suma a cada lado de la igualdad, obtenemos la siguiente
expresión:
Por último, necesitamos elevar al cuadrado ambos lados de la igualdad, ya que en general suma de los
residuales es muy pequeña, un problema similar al que ocurren cuando se calcula la varianza de una variable,
es necesario elevarla al cuadrado. Finalmente obtenemos la siguiente expresión y para poder minimizar la
suma al cuadrado de esta diferencia, se utiliza el método de mínimos cuadrados ordinarios, de ahí el nombre
del método.
3
No entraremos en los detalles técnicos del método, ya que dependeremos de R para que realice todos los
cálculos necesarios y nos concentraremos en las funciones del software para estimar modelos de regresión
simple y múltiple así como poder determinar qué tan bueno o malo es el modelo estimado y como interpretar
los resultados.
2.1 Residuales
En un análisis de regresión se necesita verificar que se cumplen los siguientes supuestos antes de continuar
con la interpretación de los resultados:
• Los residuales tienen media cero y son normales
• Los residuales son homocedásticos, es decir, la varianza de los residuales tiene una varianza constante
para los diferentes valores de x, cuando se rompe este supuesto se dice que los errores son heterocedásticos
• Los valores de los residuales son independientes
Es muy común que para verificar estos supuestos se utilicen métodos gráficos, sin embargo, aunque son muy
didácticos cuando se tienen muestras pequeñas o cuando los datos exhiben un patrón muy obvio que se aleja
de los supuestos, en muestras grandes resulta poco práctico el uso de métodos gráficos y es por ello por lo
que existen un gran número de pruebas de hipótesis para verificar la validez de estos supuestos, aunque para
un ojo entrenado puede sí puede resultar muy útil. En este capítulo estudiaremos ambos enfoques y se deja
al estudiante que utilice el método que mejor se adapte a sus necesidades y experiencias.
Algunos de los métodos gráficos utilizados con frecuencia para analizar estos supuestos son los siguientes:
• Gráfica de residuales contra la variable independiente x: Esta gráfica muestra en el eje vertical los
residuales y en el eje horizontal la variable independiente, sirve para identificar si la varianza de los
residuales es constante.
• Gráfica de residuales contra la variable estimada ŷ: Esta gráfica muestra en el eje vertical los residuales
y en eje horinzontal la variable estimada ŷ, sirve para identificar si la varianza de los residuales es
constante, suele utilizarse esta gráfica de los residuales contra la variable estimada cuando se realiza
una regresión lineal múltiple.
• Gráfica de residuales estandarizados: Es común también que se realicen gráficas con los residuales
estandarizados. Una variable estandarizada se calcular restando su media y dividiendo el resultado entre
la desviación estándar, dado que en la regresión lineal la media de los residuales es cero, los residuales
estandarizados se calcula solo dividiendo entre la desviación estándar. Una gráfica de los residuales
estandarizados en el eje vertical y la variable independiente (o la variable estimada ŷ) permite ver si el
error tiene distribución normal.
• Gráfica de probabilidad normal: Esta gráfica muestra en el eje vertical los residuales estandarizados y
en el eje horizontal puntos normales, es decir, la distribución teórica que debería seguir, esta gráfica
permite verificar si los residuales se distribuyen normalmente.
y = β0 + β1 x + u
Comenzamos cargando los datos de manera habitual, para esta sección será necesario descargar la librería
modelr la cual nos permitirá predecir la variable dependiente estimada (ŷ) y los residuales de la regresión (û).
A continuación mostraremos las funciones para realizar un análisis de regresión lineal simple, describiendo los
resultados obtenidos y su interpretación, las funciones para realizar el análisis de residuales (gráfico y formal).
Posteriormente se mostrará el mismo procedimiento para una regresión lineal múltiple, donde las funciones
son las mismas y solo cambia la cantidad de variables del modelo.
4
Tras cargar los datos podemos observar que solo contiene dos variables con nombres genericos:
setwd("~/Dropbox/Curso Est R/Cap_6/capitulo 6")
rm(list=ls())
library(readxl)
library(dplyr)
library(modelr)
library(tseries)
library(lmtest)
glimpse(base1)
## Rows: 500
## Columns: 2
## $ y <dbl> 68.17749, 53.91755, 76.74119, 67.18822, 59.30210, 80.44753, 100.3023~
## $ x <dbl> 4.652587, 3.105337, 5.634285, 4.505145, 4.095304, 5.820374, 7.976000~
Para estimar el modelo de regresión debemos utilizar la función lm la cual permite estimar modelos lineales,
la sintaxis es la siguiente, debemos guardar este modelo en un objeto, ya que genera información que
necesitaremos más adelante. La sintaxis de la función lm requiere indicar primero la ecuación que se desea
estimar, indicando primero la variable dependiente y después del signo ~ las variables independientes, en este
caso al ser una regresión simple, solo se indica una variable, finalmente se debe indicar en la opción data el
nombre del objeto donde tomará los datos en este caso el objeto que llamamos base1 :
modelo <- lm(y ~ x, data=base1)
Para ver los resultados solo indicamos indicados el nombre del objeto, el cual nos muestra la fórmula utilizada
y los coeficientes, el intercepto se refiere a la constante y x al coeficiente asociado la variable independiente:
modelo
##
## Call:
## lm(formula = y ~ x, data = base1)
##
## Coefficients:
## (Intercept) x
## 20.79 10.07
Con estos datos, podemos decir que la ecuación de regresión estimada es:
y = 20.79 + 10.07x
Otra forma de ver los resultados de manera completa es pedir un summary al objeto de la regresión, de esta
forma nos reporta mucha más información para ver el ajuste del modelo:
summary(modelo)
##
## Call:
## lm(formula = y ~ x, data = base1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.7267 -0.7023 0.0306 0.7625 3.4174
##
5
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 20.7859 0.2410 86.26 <2e-16 ***
## x 10.0708 0.0471 213.81 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.044 on 498 degrees of freedom
## Multiple R-squared: 0.9892, Adjusted R-squared: 0.9892
## F-statistic: 4.572e+04 on 1 and 498 DF, p-value: < 2.2e-16
En la parte superior muestra la regresión estimada. Debajo de esto los residuales, indicando algunas medidas
de dispersión. Seguido una tabla de coeficientes donde nos muestra la estimación puntual, el error estándar,
el estadístico t y su valor-p. En la parte inferior nos muestra el erro estándar de los residuales, una medida
de ajuste que se llama r2 y r2 ajustado, y un estadístico F con sus grados de libertad y su valor-p. Esta
información es necesaria para validar el ajuste del modelo, es decir, tener una medida de que tan bueno es la
regresión estimada.
La primera parte de los resultados muestra algunas medias de dispersión de los residuales, recuerda que se
busca que sean normales con media cero.
Los coeficientes estimados al ser estadísticos tienen por lo tanto una distribución muestral, siguiendo una
distribución t, y podemos realizar pruebas de hipótesis e intervalos de confianza sobre los valores de estos
coeficientes, para ello la función reporta no solo la estimación del coeficiente sino también su error estándar.
De hecho, la prueba de hipótesis “estándar” que sea realiza la función es para probar la hipótesis nula de que
el coeficiente es igual a cero.
H0 : β1 = 0
Ha : β1 ̸= 0
Si no podemos rechazar esta hipótesis nula, indica que el efecto que tiene la variable independiente sobre y es
estadísticamente igual a cero, por lo tanto, no es una variable que explique comportamiento de la dependiente.
El estadístico t y el valor-p que se reporta en la función de regresión estan asociados a esta prueba de hipótesis,
los cuales se pueden verificar con algunos cálculos. Una variable independiente será significativa si es posible
rechazar esta hipótesis nula y no significativa si su coeficiente resulta estadísticamente igual a cero.
La r2 es el coeficiente de determinación y es una medida de bondad de ajuste del modelo, lo que indica es
el porcentaje de la variación en y que es explicada por la regresión o por las variables independientes. La
varianza de la variable dependiente y puede seperarse en dos partes, una parte explicada por el modelo y una
parte no explicada, la variación no explicada se refiere a los residuales.
V ariación Explicada
r2 =
V ariación de y
O también
V ariación no explicada
r2 = 1 −
V ariación de y
En general entre mayor sea r2 mayor poder explicativo tiene el modelo sobre la variable dependiente, y añadir
más variables al modelo aumenta este coeficiente, independientemente si las variables resultan significativas o
no. Es por ello que también se presenta el r2 ajustado, este es el coeficiente de determinación que se ajusta
penalizando el número de variables independientes que se introducen al modelo.
6
Finalmente, la última parte de la regresión muestra un estadístico F, este estadístico es el resultado de una
prueba de hipótesis sobre la significancia total de las variables, y es una prueba complementaria sobre todo
cuando se realiza un análisis de regresión lineal múltiple. Este estadístico analiza la siguiente prueba de
hipótesis:
H0 : β1 = β2 = ... = βk = 0
3.1 Residuales
Antes de continuar con la interpretación de los coeficientes, es necesario verificar que se cumplan los supuestos
de los residuales, una forma sencilla de hacer esto es mediante el análisis gráfico, si utilizamos la función
plot() sobre el objeto que contiene el modelo de regresión obtendremos cuatro gráficos que permiten analizar
los residuales, utilizamos primero la instrucción “par(mfrow=c(2,2))” para juntar los cuatro gráficos en uno:
par(mfrow=c(2,2))
plot(modelo)
420 420
214
3
388 214
388
Residuals
1
−1
−2
−3
50 60 70 80 90 100 −3 −2 −1 0 1 2 3
420
Standardized residuals
110
1.0
−1
0.5
246
Cook's distance
−3
0.0
7
• En la parte superior izquierda: Muestra la gráfica de residuales contra la variable estimada ŷ. Sí los
residuales están dispersos alrededor del cero, sin mostrar ningún patrón o relación entre los residuales y
la variable dependiente estimada entonces nuestro modelo es apropiado. Sí se llega a observar algún
patrón entre ambas variables, indica que la variaza de los errores no es constante para los diferentes
valores de x, es decir son heterocedásticos. La etiqueta “fitted values” se refiere a los valores estimados
de la variable dependiente (ŷ).
• En la parte superior derecha: Muestra la gráfica de probabilidad normal, sirve para comprobar que
los residuales se distribuyen normalmente utilizando los residuales estandarizados, si los residuales
forman algo parecido a una recta de 45°, entonces los residuales son normales, si forman algún patrón,
o cambian de estar arriba de la recta hacia debajo de la recta de 45°, entonces los residuales no son
normales.
• En la parte inferior izquierda: Muestra la gráfica de la raíz cuadrada del valor absoluto de los residuales
estandarizados contra la variable dependiente estimada ŷ, la cuál sirve para determinar si los residuales
son de varianza constante, al igual que la primera gráfica, se espera encontrar que los valores están
dispersos sin mostrar ningún patrón, de lo contrario el modelo no es adecuado.
• En la parte inferior derecha: Muestra la gráfica de distancia de Cook, la cual sirve para identificar si
existen algunos valores grandes que tengan una gran influencia sobre la estimación de los coeficientes,
una gran distancia indica que esas observaciones tienen un gran peso sobre la estimación. La existencia
de valores atípicos son una de las causas de heterocedasticidad.
¿Qué podemos decir sobre el ajuste de nuestro modelo? Las gráficas de los residuales no muestran ningún
patrón claro e indican que son normales, sin embargo, se recomienda el uso de pruebas formales.
Para realizar pruebas formales necesitamos estimar los residuales de la regresión y añadirlo a nuestra tabla de
datos inicial, para estimar o predecir los valores residuales solo debemos utilizar la función add_residuals() la
cual forma parte de la librería modelr, con la siguiente instrucción estamos añadiendo a nuestra tabla de
datos una variable que se llama resid, la cual contiene los residuales estimados de la regresión:
base1<- base1 %>%
add_residuals(modelo)
Para determinar si son media cero y se distribuyen normalmente, basta con pedir la media de la variable, o
hacer un histograma, como se muestra a continuación:
summary(base1$resid)
Histogram of base1$resid
80
Frequency
40
0
−3 −2 −1 0 1 2 3
base1$resid
Sin embargo, una prueba formal para validar si los residuales de una regresión son normales o no sería realizar
8
la prueba de normalidad de Jarque-Bera, la cual prueba la hipótesis nula de que la variable se distribuye
normalmente, contra la alternativa de que la variable no es normal, para utilizar esta prueba es necesario
llamar a la librería tseries, y utilizar la función jarque.bera.test(), simplemente debemos indicar la variable a
la cual realiza la prueba:
jarque.bera.test(base1$resid)
##
## Jarque Bera Test
##
## data: base1$resid
## X-squared = 0.44563, df = 2, p-value = 0.8003
La prueba de Jarque-Bera utiliza un estadístico χ2 , el cual muestra un valor pequeño y por lo tanto un
valor-p grande, de 0.8003, lo cual no permite rechazar la hipótesis nula de que los residuales son normales.
Aunque existen diferentes pruebas de normalidad, la de Jarque-Bera es bastante utilizada.
Aunque los residuales son normales, debemos verificar que sean homocedásticos, para determinar esto es
necesario realizar una prueba de heterocedasticidad, entre las diversas pruebas que existen, la de Breuch-Pagan
es una de más comunes, para realizar esta prueba es necesario cargar la librería lmtest y utilizar la función
bptest() simplemente indicando el objeto que contiene la estimación realizada:
bptest(modelo)
##
## studentized Breusch-Pagan test
##
## data: modelo
## BP = 1.2066, df = 1, p-value = 0.272
La hipótesis nula de la prueba es que los errores son homocedásticos, contra la hipótesis alternativa de que
los errores son heterocedásticos. El estadístico de prueba es un estadístico χ2 , al ser pequeño su valor-p es
alto y por lo tanto no es posible rechazar la hipótesis nula, lo cual es deseable.
3.2 Interpretación
Con todos estos elementos podemos revisar los resultados de la regresión, ya que no hay ningún problema
con los residuales de la regresión pues: son normales, con media cero y homocedásticos.
summary(modelo)
##
## Call:
## lm(formula = y ~ x, data = base1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.7267 -0.7023 0.0306 0.7625 3.4174
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 20.7859 0.2410 86.26 <2e-16 ***
## x 10.0708 0.0471 213.81 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.044 on 498 degrees of freedom
## Multiple R-squared: 0.9892, Adjusted R-squared: 0.9892
9
## F-statistic: 4.572e+04 on 1 and 498 DF, p-value: < 2.2e-16
Los coeficientes son positivos y significativos, esto debido a que podemos rechazar la prueba de hipótesis de
que ambos coeficientes son iguales a cero (nota que el estadístico t de ambos coeficientes es muy grande y
por lo tanto su valor-p es prácticamente cero, mostrando evidencia en contra de la hipótesis nula). La r2 es
de 0.9892 o 98.92%, lo que indica que el modelo de regresión explica el 99% de la variación de y, al ser un
modelo de regresión simple, la r2 ajustada es igual.
La regresión muestra que existe una relación positiva y significativa al 1% entre x y y, de tal manera que si x
aumenta en una unidad, y aumenta en promedio 10 unidades. Cuando x es igual cero, el valor esperado de y
es de 20.78, esto debido al valor de la constante.
Nota que los coeficientes asociados a la variable independiente nos indican el cambio esperado en y por un
cambio en x; es decir, supongamos que x pasa de x1 a x2 , donde x1 = 0 a x2 = 1, segun la ecuación de
regresión:
∆y = β1 ∆x
Por simplicidad, podemos asumir que ∆x = 1, y el coeficiente representaria el efecto de una unidad de x sobre
y, sin embargo, también puede usarse la expresión anterior para analizar impactos de magnitud diferente, por
ejemplo, que pasa si repentinamente x aumenta en 30 unidades, ¿cuál es el efecto esperado sobre y?
Nota que, estamos considerando el efecto o cambio esperado sobre y, no el valor final de y, para predecir el
valor final de y tendriamos que utilizar la ecuación completa considerando el valor de la constante.
10
3.4 Variables categóricas
Dentro del análisis de regresión, es posible añadir variables categóricas como variables independientes. Una
variable categórica es aquella que describe algunas características de los individuos (o unidades) que se están
analizando, por ejemplo, si hablamos de personas, podemos tener una variable categórica que indique el sexo,
la religión, la preferencia sexual, la entidad de residencia, etc.
Para añadir este tipo de variables al análisis se debe hacer en forma de variables binarias (es decir, aquellas
que solo toman el valor de cero y uno) y son de manera comparativa. Una variable categórica que tiene m
categorias, y la cual deseamos introducir en el análisis de regresión, debemos generar m − 1 variables binarias,
la categoría que sea omitida es la base con la cual se comparan las demás variables.
Afortunadamente el software nos ayudará a generar estas variables binarias, solo debemos asegurarnos de que
nuestra variable categórica esta de tipo factor.
y = β0 + β1 x1 + ... + βk xk + u
A continuación estimaremos un modelo de regresión múltiple, utilizando los datos simulados que se encuentran
en el archivo de excel “datos_múltiple.xlsx”, este archivo contiene cuatro variables una dependiente (y) y
tres variables independientes (x1 , x2 , x3 ). Con estos datos podemos estimar la siguiente ecuación:
y = β 0 + β1 x 1 + β2 x 2 + β3 x 3 + u
glimpse(base2)
## Rows: 120
## Columns: 4
## $ y <dbl> 165.19530, 63.91434, 114.58400, 151.51250, 174.84470, 202.94170, 18~
## $ x1 <dbl> 9.305174, 6.210675, 11.268570, 9.010290, 8.190608, 11.640750, 15.95~
## $ x2 <dbl> 0.1464079, 1.9276310, 0.0064795, 2.0569030, 4.4861890, 1.8561320, 2~
## $ x3 <dbl> 3.4417380, -0.4453995, -0.0658304, 0.3598472, 5.7223590, 10.3247100~
Las funciones y librerías son las mismas que en el caso de una regresión simple, solo cambia un poco la
sintaxis al momento de estimar la regresión, debemos añadir las variables independientes añadiendo el signo
de +,
modelo_multiple <- lm(y ~ x1 + x2 + x3, data=base2)
modelo_multiple
##
## Call:
## lm(formula = y ~ x1 + x2 + x3, data = base2)
##
## Coefficients:
## (Intercept) x1 x2 x3
## 36.302 9.596 8.783 6.027
11
summary(modelo_multiple)
##
## Call:
## lm(formula = y ~ x1 + x2 + x3, data = base2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -70.098 -12.628 -3.554 16.089 75.800
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 36.3019 11.6999 3.103 0.00241 **
## x1 9.5959 1.1237 8.540 6.08e-14 ***
## x2 8.7831 1.0896 8.061 7.60e-13 ***
## x3 6.0269 0.5651 10.666 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 23.93 on 116 degrees of freedom
## Multiple R-squared: 0.6794, Adjusted R-squared: 0.6711
## F-statistic: 81.95 on 3 and 116 DF, p-value: < 2.2e-16
Le ecuación de regresión estimada seria entonces:
## Residuales
Recuerda que antes de interpretar los resultados necesitamos verificar si se cumplen los supuestos de los
residuales o no. Utilizando el método gráfico no parece observase ningún problema con el modelo.
par(mfrow=c(2,2))
plot(modelo_multiple)
12
Residuals vs Fitted Normal Q−Q
Standardized residuals
60 60
1 2 3
50
Residuals
−1
−50
117 117
73
−3
73
4
Standardized residuals
60 0.5
Standardized residuals
73
60
1.5
117
2
1.0
0
0.5
−2 87
Cook's distance 73
0.0
summary(base2$resid)
13
Histogram of base2$resid
40
Frequency
20
0
−50 0 50
base2$resid
jarque.bera.test(base2$resid)
##
## Jarque Bera Test
##
## data: base2$resid
## X-squared = 2.7957, df = 2, p-value = 0.2471
bptest(modelo_multiple)
##
## studentized Breusch-Pagan test
##
## data: modelo_multiple
## BP = 2.413, df = 3, p-value = 0.4912
Además, cuando se realiza un análisis de regresión múltiple, se debe verificar que las variables independientes
no están fuertemente correlacionadas entre sí, de lo contrario, tendremos un problema de multicolinelidad,
es de decir, que la relación lineal entre las variables independientes es muy alta y esto genera problemas
en la estimación de la regresión. Una forma sencilla de verificar que existe multicolinealidad es verificar la
correlación entre las variables independientes, si existe una alta correlación entre ellas, y además el modelo
de regresión presenta un alto coeficiente de determinación y pocas variables significativas, son un indicio de
multicolinealidad. Se puede optar por eliminar una de las variables que estén fuertemente correlacionadas.
En nuestro ejemplo, al realizar un análisis de correlación observamos que la correlación entre las variables
independientes es baja:
base2 %>%
select(x1,x2,x3) %>%
cor()
## x1 x2 x3
## x1 1.00000000 -0.05872322 0.08086442
## x2 -0.05872322 1.00000000 -0.08905364
## x3 0.08086442 -0.08905364 1.00000000
4.1 Interpretación
Ahora que verificamos que el modelo de regresión múltiple estimado cumple con los supuestos, podemos
interpretar los resultados:
14
summary(modelo_multiple)
##
## Call:
## lm(formula = y ~ x1 + x2 + x3, data = base2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -70.098 -12.628 -3.554 16.089 75.800
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 36.3019 11.6999 3.103 0.00241 **
## x1 9.5959 1.1237 8.540 6.08e-14 ***
## x2 8.7831 1.0896 8.061 7.60e-13 ***
## x3 6.0269 0.5651 10.666 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 23.93 on 116 degrees of freedom
## Multiple R-squared: 0.6794, Adjusted R-squared: 0.6711
## F-statistic: 81.95 on 3 and 116 DF, p-value: < 2.2e-16
Todas las variables tienen coeficientes positivos, y todos resultaron ser estadísticamente significativos, es
decir, en todos los casos se puede rechazar la hipótesis nula de que coeficiente estimado es igual a cero.
La r2 es de 67.94%, mientras que la r2 ajustada es de 67.11%, recuerda que este último penaliza la cantidad
de variables que se introducen al modelo.
El estadístico F es de 81.95 y su valor-p es prácticamente cero, por lo tanto podemos rechazar la hipótesis
nula de que los coeficientes (β1 , β2 yβ3 ) son iguales cero, concluyendo que el modelo en general es significativo.
La variable con el coeficiente más alto es la variable x1, pues por un incremento en una unidad de este, y
aumenta 9.56 unidades, mientras que por un aumento en una unidad de x2 aumenta 8.78 unidades y por un
aumento en x3 , aumenta 6.02.
## 2.5 % 97.5 %
## (Intercept) 13.128726 59.475158
## x1 7.370299 11.821541
## x2 6.625072 10.941202
## x3 4.907739 7.146062
El uso de los intervalos de confianza para analizar los coeficientes de regresión es útil, pues nos da una idea
de los valores que puede tener el coeficiente de la variable, entre mayor sea el error estándar más amplio será
el intervalo.
15
5 Ejemplo: Violación de los supuestos
Para mostrar un caso en que no se cumplan los supuestos, podemos utilizar los datos simulados que se
encuentran en el archivo de excel “datos_múltiple2.xlsx”, con el cual podemos estimar la siguiente regresión:
y = β0 + β1 x 1 + β2 x 2 + β3 x 3 + β4 x 4 + u
summary(modelo_multiple2)
##
## Call:
## lm(formula = y ~ x1 + x2 + x3 + x4, data = base3)
##
## Residuals:
## Min 1Q Median 3Q Max
## -34253 -17807 -4835 15842 52616
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21478.5 2162.9 9.930 < 2e-16 ***
## x1 -184.2 1367.0 -0.135 0.892932
## x2 2465.5 696.0 3.542 0.000475 ***
## x3 1423.5 416.5 3.418 0.000739 ***
## x4 580.2 342.8 1.692 0.091849 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 21480 on 245 degrees of freedom
## Multiple R-squared: 0.1025, Adjusted R-squared: 0.08789
## F-statistic: 6.999 on 4 and 245 DF, p-value: 2.379e-05
base3<- base3 %>%
add_residuals(modelo_multiple2)
par(mfrow=c(2,2))
plot(modelo_multiple2)
16
Residuals vs Fitted Normal Q−Q
Standardized residuals
249 249
40000 250 244 244
250
2
Residuals
1
0
0
−1
−40000
3
Standardized residuals
249
1.5
Standardized residuals
250 244 249
244
2
214
1.0
1
0
0.5
−2 −1
Cook's distance
0.0
10000 20000 30000 40000 50000 0.00 0.02 0.04 0.06 0.08
summary(base3$resid)
##
## Jarque Bera Test
##
## data: base3$resid
## X-squared = 18.441, df = 2, p-value = 9.898e-05
bptest(modelo_multiple2)
##
## studentized Breusch-Pagan test
##
## data: modelo_multiple2
## BP = 18.767, df = 4, p-value = 0.0008734
Puede observarse tanto gráficamente, como formalmente por las pruebas de hipótesis que en este caso los
residuales no son normales (aunque si de media cero) y no son homocedásticos, por lo tanto, los errores
estándar de los coeficientes no son útiles y no podemos utilizarlos para analizar los coeficientes. Naturalmente,
surgen dos preguntas:
• ¿Qué podemos hacer si los residuales no son normales? Si no son normales pero la muestra es grande,
entonces la violación al supuesto de normalidad no es grave, y podemos continuar con el análisis, si la
muestra es pequeña, el supuesto de normalidad es importante, y se debe verificar que no está omitiendo
17
alguna variable importante, errores de medición en las variables o transformar las variables; por ejemplo,
utilizar una variable en logaritmo natural en vez de nivel.
• ¿Qué podemos hacer si los residuales son heterocedásticos? Si esto ocurren los errores estándar de los
coeficientes no son válidos, una solución consiste en verificar la calidad de los datos, es decir, se debe
verificar que no existan datos atípicos en las variables, errores de medición, o transformar variables.
Otra solución es calcular errores estándar robustos con heterocedasticidad con el método de White.
Incluso algunos autores recomiendan siempre utilizar errores estándar robustos.
Para calcular los errores estándar robustos es necesario instalar la paquetería sandwich para utilizar la función
vcovHC y junto con la función coeftest de la librería lmtest, poder estimar errores estándar robustos con
heterocedasticidad.
library(sandwich)
coeftest(modelo_multiple2, vcov = vcovHC(modelo_multiple2, type="HC1"))
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21478.48 1698.17 12.6480 < 2.2e-16 ***
## x1 -184.18 1253.62 -0.1469 0.8833169
## x2 2465.50 713.19 3.4570 0.0006438 ***
## x3 1423.48 415.41 3.4267 0.0007165 ***
## x4 580.18 328.62 1.7655 0.0787179 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
glimpse(calificaciones)
## Rows: 420
## Columns: 18
## $ observation_number <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, ~
## $ dist_cod <dbl> 75119, 61499, 61549, 61457, 61523, 62042, 68536, 63~
## $ county <chr> "Alameda", "Butte", "Butte", "Butte", "Butte", "Fre~
## $ district <chr> "Sunol Glen Unified", "Manzanita Elementary", "Ther~
## $ gr_span <chr> "KK-08", "KK-08", "KK-08", "KK-08", "KK-08", "KK-08~
## $ enrl_tot <dbl> 195, 240, 1550, 243, 1335, 137, 195, 888, 379, 2247~
## $ teachers <dbl> 10.90, 11.15, 82.90, 14.00, 71.50, 6.40, 10.00, 42.~
18
## $ calw_pct <dbl> 0.5102, 15.4167, 55.0323, 36.4754, 33.1086, 12.3188~
## $ meal_pct <dbl> 2.0408, 47.9167, 76.3226, 77.0492, 78.4270, 86.9565~
## $ computer <dbl> 67, 101, 169, 85, 171, 25, 28, 66, 35, 0, 86, 56, 2~
## $ testscr <dbl> 690.80, 661.20, 643.60, 647.70, 640.85, 605.55, 606~
## $ comp_stu <dbl> 0.3435898, 0.4208333, 0.1090323, 0.3497942, 0.12808~
## $ expn_stu <dbl> 6384.911, 5099.381, 5501.955, 7101.831, 5235.988, 5~
## $ str <dbl> 17.88991, 21.52466, 18.69723, 17.35714, 18.67133, 2~
## $ avginc <dbl> 22.690000, 9.824000, 8.978000, 8.978000, 9.080333, ~
## $ el_pct <dbl> 0.000000, 4.583333, 30.000000, 0.000000, 13.857680,~
## $ read_scr <dbl> 691.6, 660.5, 636.3, 651.9, 641.8, 605.7, 604.5, 60~
## $ math_scr <dbl> 690.0, 661.9, 650.9, 643.5, 639.9, 605.4, 609.0, 61~
La unidad de medida son los distritos escolares de California y las variables representan los promedios de las
escuelas de estos distritos. La descripción de las variables es la siguiente:
• dist_code: Código del distrito escolar
• read_scr: Calificación promedio en la sección de lectura
• math_scr: Calificación promedio en la sección de matemáticas
• county : Condado
• district: Distrito escolar
• gr_span: Grados escolares en el distrito
• enrl_tot: Total de alumnos
• teachers: Número de profesores
• computer: Número de computadoras
• testscr: Calificación promedio (read_scr+math_scr)/2
• comp_stu: Computadotas por estudiante (computer/enrl_tot)
• expn_stu: Gasto por estudiante (dólares)
• str: Alumnos por profesor (enrl_tot/teachers)
• el_pct: Porcentaje de alumnos que estudian inglés, porcentaje
• meal_pct: Porcentaje de alumnos que reciben el programa de almuerzo gratis, porcentaje
• calw_pct: Porcentaje de alumnos que reciben el programa de asistencia calworks, porcentaje
• avginc: Ingreso promedio del distrito (miles de dólares)
El problema consiste en determinar si existe alguna relación entre el desempeño de los alumnos, medido a
través de la variable testscr, y el cociente str (alumnos por profesor), para contestar esta pregunta podemos
proponer la siguiente regresión lineal. Se omite el análisis descriptivo para efectos didácticos y dedicar este
espacio al tema de capítulo, sin embargo, se recomienda realizar un análisis descriptivo o gráfico de las
variables de interés utilizando lo aprendido en capítulos previos, este tipo de análisis generalmente sirve para
darnos una idea de las relaciones esperadas entre las variables.
testscr = β0 + β1 str + u
Nota: Observa como en este caso la constante β0 no tiene un significado real, y representa solo un intercepto
de la ecuación de regresión. Si la variable str es igual cero, la calificación esperada en ese distrito es de β0 , lo
cual no tiene sentido.
Podemos estimar una regresión lineal simple:
regresion <- lm(testscr ~ str,data=calificaciones)
regresion
##
## Call:
## lm(formula = testscr ~ str, data = calificaciones)
##
## Coefficients:
19
## (Intercept) str
## 698.93 -2.28
Sin embargo, contamos con más datos que pueden servir para explicar o determinar mejor si la relación que
existe entre las variables de interes es significativa o no. Por ello, podemos utilizar las siguientes variables:
summary(regresion)
##
## Call:
## lm(formula = testscr ~ str + el_pct + meal_pct + calw_pct + comp_stu +
## avginc, data = calificaciones)
##
## Residuals:
## Min 1Q Median 3Q Max
## -30.6193 -5.1734 -0.0294 5.2517 26.5672
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 671.59606 5.68803 118.072 < 2e-16 ***
## str -0.45978 0.23622 -1.946 0.0523 .
## el_pct -0.19708 0.03331 -5.917 6.88e-09 ***
## meal_pct -0.37051 0.03585 -10.335 < 2e-16 ***
## calw_pct -0.06372 0.05675 -1.123 0.2622
## comp_stu 13.58963 6.84235 1.986 0.0477 *
## avginc 0.66980 0.08329 8.042 9.42e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.411 on 413 degrees of freedom
## Multiple R-squared: 0.8079, Adjusted R-squared: 0.8051
## F-statistic: 289.5 on 6 and 413 DF, p-value: < 2.2e-16
Recuerda que antes de interpretar los resultados, es necesario verificar los supuestos del modelo:
20
ajuste<- calificaciones %>%
add_residuals(regresion)
par(mfrow=c(2,2))
plot(regresion)
4
30
Standardized residuals
367 367
2
10
Residuals
0
−10
−2
−30
6 180
1806
6 180
Standardized residuals
367
1.5
2
1.0
10
−2
0.5
Cook's
1806 distance
−4
0.0
0.5
620 640 660 680 700 0.00 0.02 0.04 0.06 0.08 0.10 0.12 0.14
summary(ajuste$resid)
##
## Jarque Bera Test
##
## data: na.omit(ajuste$resid)
## X-squared = 10.027, df = 2, p-value = 0.006648
bptest(regresion)
##
## studentized Breusch-Pagan test
##
## data: regresion
## BP = 7.1131, df = 6, p-value = 0.3105
Los resultados indican que los residuales no son normales, y no existe hetercodeasticidad.
21
Podemos tratar de mejorar los resultados, si realizamos alguna transformación de la variable de ingreso
(avginc), muchas de las variables monetarias tienen distribuciones no normales, pero al transformarlas en
logaritmo natural, se comportan de manera normal, podemos realizar esta transformación para evitar tener
mucha variabilidad causada por la dispersión del ingreso.
Generamos la variable lningreso, la cual es el logaritmo natural del ingreso promedio del distrito, y con esto
realizamos de nuevo la estimación de la regresión:
calificaciones <- calificaciones %>%
mutate(lningreso=log(avginc))
summary(regresion)
##
## Call:
## lm(formula = testscr ~ str + el_pct + meal_pct + calw_pct + comp_stu +
## lningreso, data = calificaciones)
##
## Residuals:
## Min 1Q Median 3Q Max
## -31.789 -5.238 -0.283 4.972 25.545
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 652.71255 7.98174 81.776 < 2e-16 ***
## str -0.59260 0.23903 -2.479 0.0136 *
## el_pct -0.17642 0.03361 -5.248 2.46e-07 ***
## meal_pct -0.37128 0.03852 -9.637 < 2e-16 ***
## calw_pct -0.05844 0.05791 -1.009 0.3136
## comp_stu 17.21103 6.96801 2.470 0.0139 *
## lningreso 11.68209 1.73174 6.746 5.15e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.585 on 413 degrees of freedom
## Multiple R-squared: 0.7999, Adjusted R-squared: 0.797
## F-statistic: 275.1 on 6 and 413 DF, p-value: < 2.2e-16
Verificamos de nuevo el ajuste del modelo:
ajuste<- calificaciones %>%
add_residuals(regresion)
par(mfrow=c(2,2))
plot(regresion)
22
Residuals vs Fitted Normal Q−Q
30
Standardized residuals
367 367
2
10
Residuals
0
−10
−2
−30
6 180
1806
−4
620 640 660 680 −3 −2 −1 0 1 2 3
180
Standardized residuals
6
367
1.5
2
1.0
0
10
−2
0.5
Cook's
180 6 distance
−4
0.0
0.5
summary(ajuste$resid)
##
## Jarque Bera Test
##
## data: na.omit(ajuste$resid)
## X-squared = 9.8367, df = 2, p-value = 0.007311
bptest(regresion)
##
## studentized Breusch-Pagan test
##
## data: regresion
## BP = 5.7049, df = 6, p-value = 0.457
Finalmente, si seguimos la recomendación de utilizar errores estándar robustos, recuerda que se requiere de la
librería sandwich y lmtest. Aunque los residuales no resultaron normales, la muestra es relativamente grande,
y no debemos preocuparnos por esto, y solo debemos verificar que los errores si sean homocedásticos, caso
que se resuelve considerando los errores estándar robustos.
coeftest(regresion, vcov = vcovHC(regresion, type="HC2"))
##
23
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 652.712550 9.224675 70.7572 < 2.2e-16 ***
## str -0.592597 0.267035 -2.2192 0.02702 *
## el_pct -0.176421 0.038060 -4.6353 4.786e-06 ***
## meal_pct -0.371276 0.046906 -7.9153 2.294e-14 ***
## calw_pct -0.058435 0.063632 -0.9183 0.35898
## comp_stu 17.211025 7.905303 2.1771 0.03003 *
## lningreso 11.682090 1.811451 6.4490 3.150e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Además, podemos también calcular los intervalos de confianza, utilizando los errores estándar robustos:
coefci(regresion ,vcov = vcovHC(regresion, type="HC2"))
## 2.5 % 97.5 %
## (Intercept) 634.5793797 670.84572101
## str -1.1175136 -0.06768101
## el_pct -0.2512357 -0.10160545
## meal_pct -0.4634802 -0.27907213
## calw_pct -0.1835186 0.06664809
## comp_stu 1.6713775 32.75067315
## lningreso 8.1212774 15.24290358
Utilizando los resultados de la regresión, observamos que las variables independientes explican en 79.7% la
variación de las calificaciones de las pruebas estandarizadas, el estadístico F también indica que el modelo es
significativo, pues en conjunto los coeficientes son distintos de cero. La variable de interés (str) es negativa y
significativa, es decir, reducir la relación de estudiantes por profesor, aumenta la calificación de los estudiantes,
la variable es significativa al 5%. Reducir el número de estudiantes por profesor en un alumno, indica que la
calificación promedio aumenta 0.59 puntos en promedio.
El resto de las variables resultaron significativas, con excepción de la variable “calw_pct”. Las variables
“el_pct”, “meal_pct” tiene efectos negativos sobre la calificación, es decir, al igual que cociente de estudiantes
por alumnos, guardan una relación inversa con la calificación, y la variable “comp_stu”, tiene un efecto
positivo (una relación directa con la calificación).
Dicho de otra forma, si las autoridades escolares desean aumentar las calificaciones de los estudiantes en los
distritos escolares, necesitan tomar en cuenta las variables que son significativas, ¿cuál variable tiene el mayor
impacto?
Nota que algunas variables son porcentajes, y otra es un logaritmo natural, la interpretación de los coeficientes
puede cambiar según la forma funcional de la regresión, además que es necesario tener en cuenta las unidades
de la medida de la variable dependiente y las independientes:
• en la variable dependiente, la calificación promedio se mide en puntos del examen
• en el caso de la variable str, la unidad de medida es alumnos por profesor, un aumento de un alumno
por profesor en el distrito disminuye en promedio la calificación en 0.5925 puntos
• en el caso de la variable el_pct, la unidad de medida es un porcentaje, de 0% a 100%, un aumento en
1% en el porcentaje de alumnos que asisten a clases de inglés en el distrito disminuye en promedio la
calificación en 0.1764 puntos
• en el caso de la variable meal_pct, la unidad de medida es un porcentaje, de 0% a 100%, un aumento
en 1% en el porcentaje de alumnos que reciben este beneficio en el distrito disminuye en promedio la
calificación en 0.3712 puntos
• en el caso de la variable comp_stu, la unidad de medida es computadores por estudiante, un aumento
de una computadora por estudiante en el distrito, aumenta en promedio la calificación en 17.21 puntos
24
• en el caso de la variable lningreso, la unidad de medida es una escala logarítmica, un aumento en 1%
en el ingreso en el distrito, aumenta en promedio la calificación en 0.1168 puntos
Estos resultados son mejores que los de una regresión lineal simple, y le dan a los encargados de tomar
decisiones un panorama más amplio sobre las condiciones que afectan el desempeño escolar de los estudiantes.
7 Ejercicios
1. El archivo “ventas.csv” contiene datos sobre las ventas de una compañía y los montos asignados a
publicidad en distintos medios (TV, radio y redes sociales).
Las variables que contiene la base de datos son:
• TV: Presupuesto designado a promoción en (millones de dólares)
• Radio: Presupuesto designado a promoción en redes sociales (millones de dólares)
• Social media: Presupuesto designado a promoción en radio (millones de dólares)
• Influencer: Indica si la promoción incluye la promoción de un influencer según su tamaño: Mega, Macro,
Nano, Micro
• Sales: Ventas de la empresa (millones de dólares)
Utiliza esta información para estimar la siguiente regresión lineal:
25