Para comenzar, revisamos si tenemos todas la librerías necesarias con estos comandos en la consola:

library(performance)
library(skimr)
library(dplyr)

En caso de no tenerlos instalados, utilizamos estas líneas de código en la consola (removiendo el numeral # del inicio de cada línea):

# install.packages("performance", dependencies=TRUE)
# install.packages("see", dependencies=TRUE)
# install.packages("skimr", dependencies=TRUE)
# install.packages("dplyr", dependencies=TRUE)

Cargar datos

En esta data, tenemos los siguiente datos de arboles:

Diámetro: pulgadas
Altura en pies.
Volumen en pies cúbicos.

library(dplyr)
data_arboles <- trees
names(data_arboles) <- c("diametro","altura","volumen")
head(data_arboles)

diametro	altura	volumen
8.3	70	10.3
8.6	65	10.3
8.8	63	10.2
10.5	72	16.4
10.7	81	18.8
10.8	83	19.7

Realizamos un análisis exploratorio univariado para cada variable.

library(skimr)
skim(data_arboles)

Data summary
Name	data_arboles
Number of rows	31
Number of columns	3
_______________________
Column type frequency:
numeric	3
________________________
Group variables	None

Variable type: numeric

skim_variable	complete_rate	mean	sd	p0	p25	p50	p75	p100	hist
diametro	1	13.25	3.14	8.3	11.05	12.9	15.25	20.6	▃▇▃▅▁
altura	1	76.00	6.37	63.0	72.00	76.0	80.00	87.0	▃▃▆▇▃
volumen	1	30.17	16.44	10.2	19.40	24.2	37.30	77.0	▇▅▁▂▁

Exploración bivariada.

Observaremos aquí la relación entre dos variables cuantitativas utilizando la correlación.

Correlación entre dos variables: diámetro y volumen.

cor(data_arboles$diametro,data_arboles$volumen)

## [1] 0.9671194

¿Existe una alta correlación entre estas variables? Si.

Ahora, exploramos las relaciones entre pares de variables con el comando cor y pairs aplicado en todo el conjunto de datos.

# Matriz de correlaciones:
cor(data_arboles)

##           diametro    altura   volumen
## diametro 1.0000000 0.5192801 0.9671194
## altura   0.5192801 1.0000000 0.5982497
## volumen  0.9671194 0.5982497 1.0000000

# Gráfico de dispersión para pares de variables.
pairs(data_arboles)

¿Qué tipo de relaciones observamos entre las variables? ¿Entre qué variables encontramos alta correlación?

Observamos relaciones lineales positivas entre las variables.
Alta correlación entre volumen y diámetro. (\(r=0.97\))
Moderada correlación entre volumen y altura. (\(r=0.52\))

Modelos de regresión lineal.

Modelo 1

Modelo de regresión lineal simple para el volumen utilizando el diámetro como variable explicativa.

Obtenemos los coeficientes del modelo con el código:

m1_diametro <- lm(volumen~diametro,data_arboles)
m1_diametro

## 
## Call:
## lm(formula = volumen ~ diametro, data = data_arboles)
## 
## Coefficients:
## (Intercept)     diametro  
##     -36.943        5.066

Modelo ajustado:

\[Volumen = -36.943 + 0.0564*Diametro\]

Interpretación de coeficientes:

El intercepto (-36.943) no tiene interpretación en este caso.
Por cada pulgada de aumento en el diámetro, el volumen aumenta en 5.066 pies cúbicos.

summary(m1_diametro)

## 
## Call:
## lm(formula = volumen ~ diametro, data = data_arboles)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -8.065 -3.107  0.152  3.495  9.587 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -36.9435     3.3651  -10.98 7.62e-12 ***
## diametro      5.0659     0.2474   20.48  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.252 on 29 degrees of freedom
## Multiple R-squared:  0.9353, Adjusted R-squared:  0.9331 
## F-statistic: 419.4 on 1 and 29 DF,  p-value: < 2.2e-16

Revisión de supuestos.

Por lo tanto, para ver si un modelo de regresión lineal ajustado es válido, debemos comprobar que se cumplen estos supuestos sobre los residuos:

Normalidad: Los residuos \(e_i\) tienen distribución normal con media cero para cada valor de X.
Independencia: los residuos deben ser independientes entre sí.
Homocedasticidad (igualdad de varianzas): para cada valor de la variable X, la varianza de los residuos (\(e_i=\hat{Y_i}-Y_i\)) debe ser la misma (es decir, que el ajuste es igual de preciso independientemente de los valores que tome X)

# Normalidad
performance::check_normality(m1_diametro) # Test Shapiro-Wilks

## OK: residuals appear as normally distributed (p = 0.718).

# Independencia
performance::check_autocorrelation(m1_diametro) # Test Durbin-Watson

## OK: Residuals appear to be independent and not autocorrelated (p = 0.088).

# Homocedasticidad
performance::check_heteroscedasticity(m1_diametro) #Test de Breusch-Pagan

## Warning: Heteroscedasticity (non-constant error variance) detected (p = 0.043).

Gráficamente:

performance::check_model(m1_diametro)

Si no se cumplen algunos de sus hipótesis básicas:

Efectuar una transformación de los datos de manera que los datos ya cumplan todas las hipótesis del modelo.
- Transformaciones de la familia Box-Cox son una alternativa.
Probar con otro tipo de modelos para los que si cumplan los supuestos.

Modelo 2

Modelo de regresión lineal simple para el volumen utilizando la altura como variable explicativa.

m2_altura <- lm(volumen~altura,data_arboles)
m2_altura

## 
## Call:
## lm(formula = volumen ~ altura, data = data_arboles)
## 
## Coefficients:
## (Intercept)       altura  
##     -87.124        1.543

Modelo ajustado

\[Volumen = -87.124 + 1.543*Altura\]

Interpretación de coeficientes:

El intercepto (-87.124) no tiene interpretación en este caso.
Por un aumento en 1 pie en la altura, el volumen aumenta en 1.543 pies cúbicos.

summary(m2_altura)

## 
## Call:
## lm(formula = volumen ~ altura, data = data_arboles)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -21.274  -9.894  -2.894  12.068  29.852 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -87.1236    29.2731  -2.976 0.005835 ** 
## altura        1.5433     0.3839   4.021 0.000378 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 13.4 on 29 degrees of freedom
## Multiple R-squared:  0.3579, Adjusted R-squared:  0.3358 
## F-statistic: 16.16 on 1 and 29 DF,  p-value: 0.0003784

Revisión de supuestos.

performance::check_model(m2_altura)

# Normalidad
performance::check_normality(m2_altura) # Test Shapiro-Wilks

## OK: residuals appear as normally distributed (p = 0.168).

# Independencia
performance::check_autocorrelation(m2_altura) # Test Durbin-Watson

## Warning: Autocorrelated residuals detected (p < .001).

# Homocedasticidad
performance::check_heteroscedasticity(m2_altura) #Test de Breusch-Pagan

## Warning: Heteroscedasticity (non-constant error variance) detected (p = 0.006).

¿Cuál de los modelos es mejor?

Veamos sus estadísticos y gráficos:

library(performance)
compare_performance(m1_diametro, m2_altura)

Name	Model	AIC	BIC	R2	R2_adjusted	RMSE	Sigma
m1_diametro	lm	181.6447	185.9467	0.9353199	0.9330895	4.11254	4.251988
m2_altura	lm	252.7986	257.1005	0.3579026	0.3357614	12.95762	13.396982

plot(compare_performance(m1_diametro, m2_altura, rank = TRUE))

Modelo 3

Modelo de regresión lineal múltiple para el volumen utilizando el diámetro y altura como variables explicativas.

m3_diam_altura <- lm(volumen~altura+diametro,data_arboles)
m3_diam_altura

## 
## Call:
## lm(formula = volumen ~ altura + diametro, data = data_arboles)
## 
## Coefficients:
## (Intercept)       altura     diametro  
##    -57.9877       0.3393       4.7082

Observamos sus resultados:

summary(m3_diam_altura)

## 
## Call:
## lm(formula = volumen ~ altura + diametro, data = data_arboles)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -6.4065 -2.6493 -0.2876  2.2003  8.4847 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -57.9877     8.6382  -6.713 2.75e-07 ***
## altura        0.3393     0.1302   2.607   0.0145 *  
## diametro      4.7082     0.2643  17.816  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.882 on 28 degrees of freedom
## Multiple R-squared:  0.948,  Adjusted R-squared:  0.9442 
## F-statistic:   255 on 2 and 28 DF,  p-value: < 2.2e-16

Revisión de supuestos

performance::check_model(m3_diam_altura)

performance::check_autocorrelation(m3_diam_altura)

## Warning: Autocorrelated residuals detected (p = 0.016).

# Revisión de multicolinealidad
performance::check_collinearity(m3_diam_altura)

Term	VIF	SE_factor
altura	1.36921	1.170132
diametro	1.36921	1.170132

Comparación de varios modelos:

Podemos comparar rappidamente sus resultados utilizando compare_performance del paquete performance.

compare_performance(m1_diametro, m2_altura, m3_diam_altura)

Name	Model	AIC	BIC	R2	R2_adjusted	RMSE	Sigma
m1_diametro	lm	181.6447	185.9467	0.9353199	0.9330895	4.112540	4.251988
m2_altura	lm	252.7986	257.1005	0.3579026	0.3357614	12.957617	13.396982
m3_diam_altura	lm	176.9100	182.6459	0.9479500	0.9442322	3.689223	3.881832

plot(compare_performance(m1_diametro, 
                         m2_altura, 
                         m3_diam_altura,
                         rank = TRUE))

Caso: árboles

Timteo Torres Pinchi

13 Mayo, 2021