Obtener una regresion lineal tal que :
\[ Y = \beta_0+\sum_{i=1}^{n} \beta_i X_i \] Donde:
\(Y\) : Variable Dependiente.
\(X_i\) : Predictores.
\(\beta_i\) : Coeficientes.
El dueño de criaderos de pescado desea monitorear de mejor manera su producción, para ello está considerando introducir cámaras en los estanques para tomar mediciones de los peces a través de las imágenes y con eso estimar los pesos de cada pez. Nos han contratado como empresa consultora para evaluar la viabilidad de estimar el peso a través de las imágenes.
Dentro de nuestro contexto la variable dependiente es Weight, y los posibles predictores son Length1, Length2, Length3, Height y Width.
Antes de comenzar el análisis de la base de datos, me percató de que había una observación de un pez que tenía como peso 0. Dado que éste fue un error en la transcripción de los datos y no en el cálculo de ellos, eliminó de la base esta observación.
Dentro de la base de datos se cuenta con una variable categórica que toma los siguientes valores Bream, Roach, Whitefish, Parkki, Perch, Pike y Smelt, que corresponden a las especies de los peces. Éstas tienen distintas observaciones y las podemos observar en el siguiente gráfico.
Estudió la distribución de los datos numéricos para tener una noción de qué variables tienen un comportamiento normal, así como también conocer qué varibales poseen valores atípicos.
En la siguiente tabla mostramos las variables numéricas que posee la base de datos.
| Columna | Nombre |
|---|---|
| 1 | Weight |
| 2 | Length1 |
| 3 | Length2 |
| 4 | Length3 |
| 5 | Height |
| 6 | Width |
Ahora veó algunos gráficos de estas variables:
Lo que se observa de estas gráficas es que las tres variables correspondientes a longitud (length) podrían tener la misma distribución, además de que la variable peso tiene una distribución muy distinta a las demás.
Los datos outliers que se muestran son de tres observaciones de la especie Pike que tienen datos más grandes que los demás, pero dentro de su misma especie no son datos anormales, por lo que se decide conservarlos, además de que no afectaron al modelo ni a los supuestos.
Antes de comenzar el desarrollo del modelo, es útil tener una idea visual de las relaciones dentro de los datos.
El gráfico presentado a continuación, proporciona una comparación por pares de todas las observaciones que utilizaremos como posibles predictores.
Se puede ver que las variables referentes a la longitud presentan una correlación muy alta entre sí, por lo que podrían presentar problemas de multicolinealidad. Esto era de esperarse, ya que en los diagramas anteriores notamos que estas observaciones presentaban una distribución idéntica; y este supuesto se refuerza con el siguiente diagrama.
Como se puede observar, las distribuciones de las variables mencionadas son casi idénticas visualmente, por lo que posiblemente más adelante lleguemos a la conclusión de que sólo debemos considerar una de ellas.
Los criterios Best Subset Selection y AIC suguieren usar 5 o 4 predictores, pues éstos resultaron ser los mejores, a pesar de todo, estos modelos sugeridos no lograron cubrir todos los supuestos requeridos para la regresión lineal. Por lo que se decide construir modelos para cada especie; este análisis se resume en el anexo, con el fin de dentificar los predictores comunes para cada especie que tambien cumplen los supuestos necesarios y usarlos para construir un modelo en general.
El modelo es el siguiente y cuenta con las siguientes características, este modelo no es el definitivo pues aún se ajustará alguna transformación para mejorar los supuestos requeridos en la regresión lineal.
\[ Y = -499.73573 + 20.75766 \cdot X_{Height } + 27.15694 \cdot X_{Length1 } \]
##
## Call:
## lm(formula = Weight ~ Length1 + Height, data = FishB)
##
## Residuals:
## Min 1Q Median 3Q Max
## -285.33 -75.27 -24.58 71.12 379.67
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -499.736 28.966 -17.253 < 2e-16 ***
## Length1 27.157 1.291 21.028 < 2e-16 ***
## Height 20.758 3.010 6.896 1.28e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 126.5 on 155 degrees of freedom
## Multiple R-squared: 0.8764, Adjusted R-squared: 0.8749
## F-statistic: 549.8 on 2 and 155 DF, p-value: < 2.2e-16
El modelo es capaz de explicar el 87.64% de la variabilidad observada en el peso (R-squared: 0.8764). El valor de \(R^{2}\)-ajustado es cercano al \(R^{2}\) (Adjusted R-squared: 0.8749 ), lo que indica que el modelo contiene predictores útiles. El p-value del modelo es significativo (< 2.2e-16) por lo que se puede decir que el modelo no es por azar, al menos uno de los coeficientes parciales de regresión es distinto de 0. Las variables que consideramos son significativas por lo que nuestro modelo tiene una buena estructura.
La transformación que ayudará a mejorar el supuesto de normalidad en nuestro modelo es obtenido a partir del diagrama siguiente, donde según el valor de \(\lambda\) se toma la transformación más adecuada.
Transformaciones comunes
| Valor \(\lambda\) | Transformed Data |
|---|---|
| -3 | \(Y^{-3}\) |
| -2 | \(Y^{-2}\) |
| -1 | \(Y^{-1}\) |
| -0.5 | \(Y^{-0.5}\) |
| 0 | \(log(Y)\) |
| 0.5 | \(Y^{0.5}\) |
| 1 | \(Y\) |
| 2 | \(Y^{2}\) |
| 3 | \(Y^{3}\) |
Observamos que \(\lambda\) es cercana a 0.5, por lo cual la transformación que aplicaremos será \(Y^{0.5}\)
##
## Call:
## lm(formula = (Weight)^0.5 ~ Height + Length1, data = Fish_)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.3162 -0.7856 -0.0886 0.7935 5.1670
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.27658 0.37702 -16.65 <2e-16 ***
## Height 0.77971 0.03918 19.90 <2e-16 ***
## Length1 0.65030 0.01681 38.69 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.647 on 155 degrees of freedom
## Multiple R-squared: 0.9679, Adjusted R-squared: 0.9675
## F-statistic: 2339 on 2 and 155 DF, p-value: < 2.2e-16
Se puede observar que bajo esta transformación se mejoró la variabilidad observada en el peso pues ésta mejora el \(R^{2}\) ya que de 87.64% que inicial mente poseía el modelo ha pasado a 96.79%, ambas variables predictoras que hemos contemplado resultan ser significativas incluyendo el intercepto. El p-value del modelo es significativo (< 2.2e-16) por lo que se puede decir que el modelo no es por azar, al menos uno de los coeficientes parciales de regresión es distinto de 0.
\[ Y^{0.5} = -6.2765791 + 0.7797133 \cdot X_{Height } + 0.6502980 \cdot X_{Length1 } \] Para los coeficientes de los predictores tenemos los siguientes intervalos de confianza.
## 2.5 % 97.5 %
## (Intercept) -7.0213387 -5.5318194
## Height 0.7023174 0.8571092
## Length1 0.6170929 0.6835031
Cada una de las pendientes de un modelo de regresión lineal múltiple (coeficientes parciales de regresión de los predictores) se define del siguiente modo: Si el resto de variables se mantienen constantes, por cada unidad que aumenta el predictor en cuestión, la variable (Y) varía en promedio tantas unidades como indica la pendiente.
Para este ejemplo, por cada unidad que aumenta el predictor Length1, \((Weight)^{0.5}\)aumenta en promedio 0.6502980 unidades, manteniéndose constantes el resto de predictores.
Haciendo esta modificación se puede comenzar a probar supuestos:
Se puede notar en el siguiente diagrama que los residuos estandarizados en el modelo con la transformación descrita tienen la forma de una distribución normal, por lo que podemos sospechar que el supuesto de normalidad se está cumpliendo.
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -3.254434 -0.479652 -0.054017 -0.000766 0.484233 3.151633
En este diagrama vemos que la mayoria de los residuales estandarizados se concentran alrededor del 0, lo cual podria indicarnos esperanza igual a 0. Sin embargo, observamos también un patrón en los residuales estandarizados, y que algunos de ellos se salen de las bandas de -2 y 2, lo cual puede traer problemas para los supuestos de no correlación y homocedasticidad, respectivamente.
Se procede a realizar las pruebas de manera formal.
Normalidad Para verificar el supuesto de normalidad, se aplica la prueba de Anderson-Darling de la siguiente forma:
\[Anderson-Darling\] \[ H_0 :\{ X_i \}^{n}_{i=1} \backsim N(\cdot) \] \[ H_a :\{ X_i \}^{n}_{i=1} \not\sim N(\cdot) \]
##
## Anderson-Darling GoF Test
##
## data: pnorm(r_estandarizados, mean = 0, sd = 1)
## AD = 1.996, p-value = 0.09238
## alternative hypothesis: NA
La prueba indica que estadísticamente que el modelo propuesto cumple el supuesto de normalidad, lo cual es muy bueno, pues parte de toda la teoría desarollada está muy ligada con esto.
Homocedasticidad Para verificar este supuesto, se realiza la prueba de Breusch-Pagan:
\[ Breusch-Pagan \] \[ H_0 : \sigma^{2}_j = \sigma^{2} \] \[vs\] \[ H_a : \sigma^{2}_j \neq \sigma^{2} \]
##
## studentized Breusch-Pagan test
##
## data: modelo.trans
## BP = 29.301, df = 2, p-value = 4.339e-07
Por lo que se pude ver que el supuesto de varianza constante no se cumple, tal como lo esperábamos por el análisis gráfico.
No correlación Aplicamos la prueba de Durbin-Watson:
\[ Durbin-Watson \] \[ H_0 : \rho_{j,i}= 0 \] \[ vs \] \[ H_0 : \rho_{j,i} > 0 \]
##
## Durbin-Watson test
##
## data: modelo.trans
## DW = 0.76285, p-value = 7.516e-16
## alternative hypothesis: true autocorrelation is greater than 0
Entonces se tiene que el supuesto de no correlación tampoco se cumple, como se esperaba.
Relación Lineal
Relación lineal entre los predictores numéricos y la variable respuesta:
Esta condición se puede validar bien mediante diagramas de dispersión entre la variable dependiente y cada uno de los predictores o con diagramas de dispersión entre cada uno de los predictores y los residuos del modelo. Si la relación es lineal, los residuos deben de distribuirse aleatoriamente en torno a \(0\) con una variabilidad constante a lo largo del eje X.
Análisis de Inflación de Varianza (VIF)
## [1] 4.324936
## Height Length1
## 1.639038 1.639038
No hay predictores que muestren una correlación lineal muy alta ni inflación de varianza. Es decir, no se tiene problemas de multicolinealidad.
Análisis de Valores Influyentes
Un análisis de valores influyetes nos muestra qué dato tiene influencia considerable en el modelo, la observación 14 tiene un residuo estandarizado >3 (más de 3 veces la desviación estándar de los residuos) por lo que se considera un dato atípico. El siguiente paso es determinar si es influyente
## No Studentized residuals with Bonferroni p < 0.05
## Largest |rstudent|:
## rstudent unadjusted p-value Bonferroni p
## 14 -3.360773 0.00098016 0.15487
## StudRes Hat CookD
## 14 -3.36077281 0.016342622 0.058655327
## 112 3.24721776 0.009175599 0.030661109
## 138 -2.59388283 0.044002504 0.099549626
## 139 -2.57800594 0.043034784 0.096123931
## 143 1.52021825 0.092939706 0.078270435
## 145 0.06316483 0.101538715 0.000151273
El análisis muestra varias observaciones influyentes, aunque ninguna excede los límites de preocupación para los valores de Leverages hat o Distancia Cook (>1).
Por lo que el modelo aunque presenta valores influyentes ninguno supera los criterios mensionados.
El modelo encontrado puede ser vizualizado ya que contamos con 3 variables y podemos trazar las observaciones con respecto al plano descrito por la regresion lineal multiple.
\[ Y^{0.5} = -6.2765791 + 0.7797133 \cdot X_{Height } + 0.6502980 \cdot X_{Length1 } \]
De la base de datos con las que se cuenta, evaluaremos que tan acertivo es el modelo propuesto, por lo que se toma una muestra aleatoria de 5 observaciones y se compara con los valores arrojados por el modelo.
## # A tibble: 5 x 4
## Weight Length1 Height Prediccion
## <dbl> <dbl> <dbl> <dbl>
## 1 8.83 17.5 5.58 9.49
## 2 26.8 32 16.4 27.3
## 3 26.5 31.9 16.2 27.2
## 4 31.9 37 12.4 27.3
## 5 8.83 16.8 5.20 8.74
Una mejor verificacion de esto es el siguiente grafico que compara la observacion del peso y lo contrasta con el valor que arroja el modelo de regresion lineal, esto con la finalidad de observar si el ajuste fue bueno, y lo que se puede observar es que el modelo es bueno al predecir el peso solo en funcion de \(Length1 , Height\)
El modelo propuesto para abordar el problema de estimar el peso de los peces en función de sus medidas observadas se redujo después de un análisis exhaustivo a 2 predictores, Length y Heigth, que eran medidas comunes en los modelos considerados por especies, pues es un buen indicador de que esos predictores son los adeacuados.
El análisis de supuestos nos ayudó a comprobar que se cumple que los residuos tienen distribución normal, además de que la Relación Lineal entre los predictores se cumplió, no se presentó una Inflación de la varianza alta. Existen valores influyentes en nuestras observaciones consideradas para los predictores, pero ninguna excede el criterio \(CookDistance\) por lo que no se considera volver a construir el modelo eliminando estos datos.
No existe una condición establecida para el número mínimo de observaciones pero, para prevenir que una variable resulte muy influyente cuando realmente no lo es, se recomienda que la cantidad de observaciones sea entre 10 y 20 veces el número de predictores. En este caso debería haber como mínimo 40 observaciones y se dispone de 158 por lo que este problema no se presenta en el modelo.
El modelo puede exlicar el 86.79% de la variabilidad observada en el peso de las especies y las pruebas realizadas para contrastar lo observado en la variable respuesta con los valores arrojados por el modelo nos indican que se logra una buena aproximación, lo que nos deja un modelo certero que cumple los supuestos más importantes en una regresión lineal múltiple.
*An Introduction to Statistical Learning: with Applications in R (Springer Texts in Statistics)
*Linear Models with R, Julian J.Faraway
*Points of Significance: Simple linear regression by Naomi Altman & Martin Krzywinski