Tatiana
Petukhova,
“INTRODUCCIÓN A R”
9. Regresión
En esta sección vamos a utilizar datos, determinar la línea de regresión de una variable dependiente dadas
variables independientes (variables explicativas), realizar predicciones y estimaciones a partir de la recta de
regresión mediante valores puntuales.
Directorio de trabajo
Estableceremos nosotros mismos un directorio donde queremos guardar los resultados.
>
setwd(
‘/Users/Documents/Curso')
Datos
En esta sección vamos a utilizar datos “mtcars” que están disponibles en R. Los datos fueron extraídos de la revista
Motor Trend US 1974 y comprenden el rendimiento y 8 aspectos del diseño y rendimiento del automóvil para 32
automóviles (modelos 1973-74).
Cargaremos el conjunto de datos que queremos “mtcars”, usando la función data()
>
data(mtcars)
>
?mtcars
>
attach(mtcars)
Formato:
Un
marco
de
datos
con
32
observaciones
sobre
11
variables.
[,
1]
mpg
-‐
el
rendimiento
en
millas
/
galón
(EE.
UU.)
[,
2]
cyl
-‐
Número
de
cilindros
[,
3]
disp.
-‐
Desplazamiento
([Link].)
[,
4]
hp
-‐
Potencia
bruta
[,
5]
drat
-‐
Radio
del
eje
trasero
[,
6]
wt
-‐
Peso
(1000
lbs)
[,
7]
qsec
-‐
Tiempo
de
1/4
milla
[,
8]
vs
-‐
El
tipo
de
motor
:
0=V-‐motor
y
1=S-‐motor
(motor
recto)
[,
9]
am
-‐
Transmisión
(0
=
automático,
1
=
manual)
[,
10]gear
-‐
Engranaje:
número
de
marchas
adelante
[,
11]carb
-‐
Número
de
carburadores
>
str(mtcars)
'[Link]':
32
obs.
of
11
variables:
$
mpg
:
num
21
21
22.8
21.4
18.7
18.1
14.3
24.4
22.8
19.2
...
$
cyl
:
num
6
6
4
6
8
6
8
4
4
6
...
$
disp:
num
160
160
108
258
360
...
$
hp
:
num
110
110
93
110
175
105
245
62
95
123
...
$
drat:
num
3.9
3.9
3.85
3.08
3.15
2.76
3.21
3.69
3.92
3.92
...
$
wt
:
num
2.62
2.88
2.32
3.21
3.44
...
$
qsec:
num
16.5
17
18.6
19.4
17
...
$
vs
:
num
0
0
1
1
0
1
0
1
1
1
...
$
am
:
num
1
1
1
0
0
0
0
0
0
0
...
$
gear:
num
4
4
4
3
3
3
3
4
4
4
...
$
carb:
num
4
4
1
1
2
1
4
2
2
4
...
Convertiremos las variables numéricas en factores:
>
mtcars$cyl=[Link](mtcars$cyl)
Tatiana
Petukhova,
“INTRODUCCIÓN A R”
>
mtcars$vs=[Link](mtcars$vs)
>
mtcars$am=[Link](mtcars$am)
>
mtcars$gear=[Link](mtcars$gear)
>
mtcars$carb=[Link](mtcars$carb)
Diagramas de dispersión
Realizaremos la matriz con los diagramas de dispersión mediante la función pairs
>
attach(mtcars)
>
pairs(mtcars)
La función rcorr(x,
type=c("pearson","spearman)) del paquete Hmisc produce correlaciones y niveles de significancia.
Argumentos:
x
-‐
una
matriz
numérica
con
al
menos
5
filas
y
al
menos
2
columnas
type
-‐
especifica
el
tipo
de
correlaciones
Se comprueban las siguientes hipótesis:
Ho: la correlación verdadera es 0
Ha: la correlación verdadera no es 0
Seleccionaremos las columnas con algunas variables cuantitativas y obtendremos correlaciones entre ellas.
>
[Link]("Hmisc")
>
library(Hmisc)
>
rcorr([Link](coches[,c(1,3:7)]))
mpg
disp
hp
drat
wt
qsec
mpg
1.00
-‐0.85
-‐0.78
0.68
-‐0.87
0.42
Tatiana
Petukhova,
“INTRODUCCIÓN A R”
disp
-‐0.85
1.00
0.79
-‐0.71
0.89
-‐0.43
hp
-‐0.78
0.79
1.00
-‐0.44
0.66
-‐0.71
drat
0.68
-‐0.71
-‐0.44
1.00
-‐0.71
0.09
wt
-‐0.87
0.89
0.66
-‐0.71
1.00
-‐0.17
qsec
0.42
-‐0.43
-‐0.71
0.09
-‐0.17
1.00
n=
32
P
mpg
disp
hp
drat
wt
qsec
mpg
0.0000
0.0000
0.0000
0.0000
0.0171
disp
0.0000
0.0000
0.0000
0.0000
0.0131
hp
0.0000
0.0000
0.0132
0.0000
0.0000
drat
0.0000
0.0000
0.0132
0.0000
0.6196
wt
0.0000
0.0000
0.0000
0.0000
0.3389
qsec
0.0171
0.0131
0.0000
0.6196
0.3389
La mayoría de la las variables están correlacionadas. Entre todas las variables, drat y qsec, wt y qsec no están
correlacionadas, pero drat y wt sí están correlacionadas. Por lo tanto, drat y wt no pueden estar incluidas en el
mismo modelo.
Pruebas de normalidad
Comprobaremos la hipótesis nula de que la muestra ha sido extraída de una población con una distribución de
probabilidad normal. A modo de ejemplo, las variables mpg (el
rendimiento
en
millas
/
galón
(EE.
UU)
y drat (Radio
del
eje
trasero).
Prueba de Kolmogorov-Smirnov
La prueba se realiza con la función [Link](), la cuál contiene varias opciones predeterminadas.
Usaremos algunas de las ellas
[Link](x,
y,
…,
alternative
=
c("[Link]",
"less",
"greater),
exact
=
NULL)
Argumentos:
x
-‐un
vector
numérico
de
valores
de
datos.
y
-‐un
vector
numérico
de
valores
de
datos.
…
-‐
parámetros
de
la
distribución
especificada
(como
una
cadena
de
caracteres)
por
y.
alternative
-‐
indica
la
hipótesis
alternativa
y
debe
ser
uno
de
"[Link]"
(predeterminado),
"less"
(menos)
o
"greater"
(mayor).
exact
-‐
un
argumento
lógico
que
indica
si
se
debe
calcular
un
valor
de
p
exacto.
#
la
variable
mpg
(el
rendimiento
en
millas
/
galón
(EE.
UU)
>
[Link](x=mpg,
"pnorm",
mean
=
mean(mpg),
sd
=
sd(mpg))
One-‐sample
Kolmogorov-‐Smirnov
test
data:
mpg
D
=
0.1263,
p-‐value
=
0.687
alternative
hypothesis:
two-‐sided
#
la
variable
drat
(Radio
del
eje
trasero).
>
[Link](x=drat,"pnorm",mean
=
mean(drat),
sd
=
sd(drat))
One-‐sample
Kolmogorov-‐Smirnov
test
Tatiana
Petukhova,
“INTRODUCCIÓN A R”
data:
drat
D
=
0.15976,
p-‐value
=
0.3876
alternative
hypothesis:
two-‐sided
Prueba de normalidad Shapiro-Wilk
La prueba se realiza para muestras pequeñas (n<30) con la función [Link]()
>
[Link](x=mpg)
Shapiro-‐Wilk
normality
test
data:
mpg
W
=
0.94756,
p-‐value
=
0.1229
El análisis de regresión lineal múltiple
Determinaremos la línea de regresión de una variable dependiente dadas las variables independientes, usando la
función lm(), la cuál se basa en el método de mínimos cuadrados. Usaremos la ayuda en línea de R para consultar
como utilizar la función lm().
>?lm
La función se determina como
lm(formula,
data,...)
Argumentos:
formula
la
expresión
que
define
el
modelo:
la
variable
dependiente
(numérica)
seguida
por
el
símbolo
tilda
“~”
y
la
variable
independiente.
data
los
datos
en
el
formato
de
marco
de
datos.
…
argumentos
adicionales
predeterminados;
los
dejamos
sin
modificar.
El algoritmo de la línea que ajuste mejor los datos y las funciones genéricas aplicadas al objeto es el siguiente
>
attach(mtcars)
>
ajuste
=
lm(mpg~cyl+drat+qsec+vs+am+gear,
data=mtcars)
>
summary(ajuste)
Call:
lm(formula
=
mpg
~
cyl
+
drat
+
qsec
+
vs
+
am
+
gear,
data
=
mtcars)
Residuals:
Min
1Q
Median
3Q
Max
-‐6.3851
-‐1.5661
0.0514
1.5051
5.5417
Coefficients:
Estimate
Std.
Error
t
value
Pr(>|t|)
(Intercept)
11.9928
21.3463
0.562
0.5797
cyl6
-‐4.3795
2.5905
-‐1.691
0.1044
cyl8
-‐7.1791
4.0838
-‐1.758
0.0921
.
drat
1.4201
2.4069
0.590
0.5609
qsec
0.3208
0.8395
0.382
0.7059
Tatiana
Petukhova,
“INTRODUCCIÓN A R”
vs1
1.5841
2.5156
0.630
0.5351
am1
4.6510
2.6219
1.774
0.0893
.
gear4
-‐2.2463
2.8736
-‐0.782
0.4424
gear5
-‐2.4115
3.0254
-‐0.797
0.4336
-‐-‐-‐
Signif.
codes:
0
'***'
0.001
'**'
0.01
'*'
0.05
'.'
0.1
'
'
1
Residual
standard
error:
3.287
on
23
degrees
of
freedom
Multiple
R-‐squared:
0.7793,
Adjusted
R-‐squared:
0.7026
F-‐statistic:
10.15
on
8
and
23
DF,
p-‐value:
5.862e-‐06
#
un
vector
nombrado
de
coeficientes
>
ajuste$coefficients
(Intercept)
cyl
drat
qsec
vs
am
gear
28.5543163
-‐1.9638545
0.9806355
0.1506520
1.1196981
4.2695900
-‐1.2888397
Análisis de residuos
Para ver si el modelo de regresión lineal es adecuado, analizaremos los residuos.
Residuos contra valores predictivos
La gráfica de los residuos contra valores predictivos permite detectar:
1. Heterocedasticidad – la varianza de los residuos no es constante: se deben transformar los datos (la
variable de respuesta 𝑌).
2. Error en el análisis: se ha realizado mal el ajuste y se verifica que los residuos negativos se corresponden
con los valores predictivos pequeños y los residuos positivos se corresponden con los valores predictivos
grandes, o al revés.
3. El modelo es inadecuado: falta de linealidad y se deben de transformar los datos o introducir nuevas
variables que pueden ser polinomios.
4. Existencia de observaciones atípicas o puntos extremos.
Las gráficas de los residuos se obtiene con la función plot(). La sintaxis de esta función es
>
par(mfrow=c(2,2))
>
plot(ajuste)
>
par(mfrow=c(1,1))
Si después de ejecutar el código para hacer la gráfica, aparece un mensaje
>
plot(ajuste)
Error
in
[Link]()
:
figure
margins
too
large
¡tiene que aumentar el tamaño de la ventana de la gráfica!
Veremos como resultado las gráficas de los residuos que aparecen abajo
Tatiana
Petukhova,
“INTRODUCCIÓN A R”
La gráfica “Residuals
vs
Fitted” (Residuos versus Valores Ajustados (pronosticados)) – ¿Hay algún patrón notable?
Los valores atípicos están etiquetados por su número de observación.
La gráfica “Normal
Q-‐Q” (quantile-quantile de los residuos estandarizados) – la evaluación de la normalidad de los
residuos. Los residuos siguen una línea recta. Las desviaciones de la línea recta significarían que los errores no
siguen una distribución normal.
La gráfica “Scale-‐Location” (Escala-Ubicación) - ¿Hay algún patrón notable?
La gráfica “Residuals vs Factor Levels” (Residuos versus niveles del factor) – ¿Hay observaciones influyentes?
Residuos contra una variable explicativa
La gráfica de residuos contra una variable explicativa permite detectar si la existencia de heterocedasticidad o la
falta de linealidad en el modelo son debidas a la variable explicativa.
Primero, tendremos que extraer los residuos del objeto de la clase lm con el comando resid(). El algoritmo para
obtener la gráfica es el siguiente
>
par(mfrow=c(2,3))
>
plot(cyl,
residuos,
xlab="cilindros")
#residuos
vs
número
de
cilindros
>
plot(drat,
residuos)
#residuos
vs
radio
del
eje
trasero
>
plot(qsec,
residuos)
#residuos
vs
tiempo
de
1/4
milla
>
plot(vs,
residuos,
xlab
=
"vs")
#residuos
vs
el
tipo
de
motor
>
plot(am,
residuos,
xlab
=
"am")
#residuos
vs
transmisión
>
plot(gear,
residuos,
xlab
=
"gear")
#residuos
vs
engranaje
>
par(mfrow=c(1,1))
Tatiana
Petukhova,
“INTRODUCCIÓN A R”
Evaluación sobre si los errores están correlacionados entre sí.
Verificaremos si los errores están correlacionados entre sí utilizando la prueba de errores autocorrelacionados -
Prueba de Durbin-Watson con la función lmtest() del paquete lmtest.
Hipótesis nula: la correlación verdadera = 0
Hipótesis alternativa: la correlación verdadera ≠ 0 - una prueba de dos colas
>
[Link]("lmtest")
>
library(lmtest)
>
dwtest(ajuste)
Durbin-‐Watson
test
data:
ajuste
DW
=
1.8821,
p-‐value
=
0.2111
alternative
hypothesis:
true
autocorrelation
is
greater
than
0
Predicciones
Si lo que queremos es ajustar el modelo para poder usarlo para predecir datos, utilizaremos la función predict(), la
cual contiene varias opciones predeterminadas. Para obtener la información sobre las variaciones disponibles de la
función predict(),
usaremos la ayuda en línea de R con el código
>
?[Link]
La función se determina como
predict(object,
newdata,
interval=c("none",
"confidence",
"prediction"),
level
=
0.95)
Argumentos:
object
Objeto
de
clase
de
"lm"
newdata
Marco
de
datos
con
los
valores
nuevos
de
las
variables
explicativas
interval
Tipo
de
cálculo
de
intervalo:
"none"=ninguno,
"confidence"=confianza,
"prediction"=predicción
level
Nivel
de
confianza,
el
valor
predeterminado
es
0.95
Utilizaremos la función predict()
modificando algunas opciones
Tatiana
Petukhova,
“INTRODUCCIÓN A R”
>
[Link]=[Link](cyl
=6,
drat=4,
qsec=20,
vs=0,
am=0,
gear=5)
>
predict(object=ajuste,
newdata=[Link],
interval="prediction")
fit
lwr
upr
1
17.26257
7.470713
27.05443
donde
fit
es
el
valor
predicho
de
la
variable
de
respuesta;
lwr
es
el
límite
inferior
del
intervalo
de
confianza,
upr
es
el
límite
superior
del
intervalo
de
confianza
