En la siguiente tabla aparece información acerca de ocho automóviles de cuatro cilindros considerados entre los más eficientes en consumo de combustible en 2006. Los tamaños de los motores se dan en volumen total de cilindros, medido en litros (l).
Cargaremos los datos manualmente en dos vectores y luego lo uniremos en la estructura de data frame. También existen estas funciones, por ejemplo:
Para buscar los archivos de la PC y leerlos proveniendo de .txt con “read.delim”
from google.colab import files
uploaded = files.upload()
datos <- read.delim(“datos.txt”, sep = ” “)
Para leer datos de hojas de Spread Sheets
library(googlesheets4)
gs4_deauth() # evita autenticación
url <- “”
datos <- read_sheet(url, sheet = “Hoja1”)
# Vector de datos de variable explicativa: Volumen de los cilindros
volumen<- c(1.8,1.5,2.0,2.5,1.8,2.5,1.6,1.5)
# Vector de datos de variable expicada: Potencia del motor (hp)
potencia<-c(101,91,115,150,126,150,118,106)
# Creamos un data frame a partir de los vectores generados
datos <- data.frame(volumen,potencia)
names(datos) #Nombres de variables
## [1] "volumen" "potencia"
head(datos) #Info de las primeras 5 filas
## volumen potencia
## 1 1.8 101
## 2 1.5 91
## 3 2.0 115
## 4 2.5 150
## 5 1.8 126
## 6 2.5 150
dim(datos) # cantidad de filas y columnas
## [1] 8 2
str(datos) # resumen de la estructura de los datos
## 'data.frame': 8 obs. of 2 variables:
## $ volumen : num 1.8 1.5 2 2.5 1.8 2.5 1.6 1.5
## $ potencia: num 101 91 115 150 126 150 118 106
summary(datos) # resumen estadístico de las variables del conjunto de datos
## volumen potencia
## Min. :1.500 Min. : 91.0
## 1st Qu.:1.575 1st Qu.:104.8
## Median :1.800 Median :116.5
## Mean :1.900 Mean :119.6
## 3rd Qu.:2.125 3rd Qu.:132.0
## Max. :2.500 Max. :150.0
Utilizaremos la librería ggplot2 de R para graficar:
library(ggplot2) #usamos la librería ggplot2 de R para graficar
ggplot(data = datos,
mapping = aes(x = volumen, y = potencia)) +
geom_point(alpha = 0.9, color = "blue") +
labs(x ="Volumen de cilindrada (l)", y = "Potencia del motor (hp)")
Calculamos el coeficiente de correlación, asumiendo que tanto X como Y son variables aleatorias.
cor(volumen,potencia) #Cálculo de coeficiente de correlación lineal
## [1] 0.887626
cor.test(volumen,potencia) #Test para el coeficiente de correlación
##
## Pearson's product-moment correlation
##
## data: volumen and potencia
## t = 4.7208, df = 6, p-value = 0.003255
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.4885078 0.9795833
## sample estimates:
## cor
## 0.887626
Utilizamos la función \(lm\).
modelo <- lm(potencia ~ volumen, datos) ## usamos la función "lm(y~x,datos)"
summary(modelo)
##
## Call:
## lm(formula = potencia ~ volumen, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.918 -9.448 2.134 6.672 12.496
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 30.194 19.321 1.563 0.16914
## volumen 47.069 9.971 4.721 0.00326 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 10.74 on 6 degrees of freedom
## Multiple R-squared: 0.7879, Adjusted R-squared: 0.7525
## F-statistic: 22.29 on 1 and 6 DF, p-value: 0.003255
Se puede hacer de manera automática a partir de la función \(confint\) seteando el nivel de confianza \((1-\alpha)\) del intervalo con el parámetro \(level\). Por defecto, se hacen al 95% de confianza.
confint(modelo,level=0.95) # usamos la función “confint(modelo, level=1-alfa)
## 2.5 % 97.5 %
## (Intercept) -17.08232 77.47025
## volumen 22.67186 71.46607
Se puede hacer “a mano” con la información de la estimación puntual y el error estándar de la tabla de coeficientes. También usamos la función \(qt\) indicando coo argumentos \((1-\alpha/2)\) y los grados de libertad \(df=n-2\)
li_beta<- 47.069-qt(0.975,df=6)*9.971 # IC para Beta, limite inferior
li_beta
## [1] 22.67084
ls_beta<- 47.069+qt(0.975,df=6)*9.971 # IC para Beta, limite superior
ls_beta
## [1] 71.46716
li_alfa<- 30.194-qt(0.975,df=6)*19.321 # IC para Alfa, limite inferior
li_alfa
## [1] -17.08278
ls_alfa<- 30.194+qt(0.975,df=6)*19.321 # IC para Alfa, limite superior
ls_alfa
## [1] 77.47078
Agregaremos una línea de código al diagrama de dispersión con la función “abline”.
library(ggplot2) #usamos la librería ggplot2 de R para graficar
ggplot(data = datos,
mapping = aes(x = volumen, y = potencia)) +
geom_point(alpha = 0.9, color = "blue") +
labs(x ="Volumen de cilindrada (l)", y = "Potencia del motor (hp)") +
geom_smooth(method = "lm",se=FALSE, color = "green")
Comenzaremos calculando los residuos \(e_i=y_i-\hat{y_i}\), los residuos estandarizados \(\frac{e_i}{s}\)y los valores predichos \(\hat{y_i}\).
re<- residuals(modelo) # Residuos
re
## 1 2 3 4 5 6 7
## -13.918103 -9.797414 -9.331897 2.133621 11.081897 2.133621 12.495690
## 8
## 5.202586
zre <- rstandard(modelo) # Residuos estandarizados
zre
## 1 2 3 4 5 6 7
## -1.3924417 -1.0626944 -0.9336129 0.2644092 1.1086923 0.2644092 1.3030748
## 8
## 0.5643080
pre <- predict(modelo) # Predichos
pre
## 1 2 3 4 5 6 7 8
## 114.9181 100.7974 124.3319 147.8664 114.9181 147.8664 105.5043 100.7974
Continuamos haciendo pruebas de normalidad para los mismos (Test de Shapiro-Wilk) y gráficos como qqplot y boxplot.
shapiro.test(re)
##
## Shapiro-Wilk normality test
##
## data: re
## W = 0.91281, p-value = 0.3743
qqnorm(re)
qqline(re, col = "red")
boxplot(re)
Agregamos un gráfico de dispersión de los residuos estandarizados vs. valores predichos de Y.
plot(pre,zre,xlab="Valores predichos", ylab="Residuos estandarizados")
abline(h = 0) # agregamos el eje x con una recta horizontal de ecuación y=0