class: bottom, left <!--- Para correr en ATOM - open terminal, abrir R (simplemente, R y }enter) - rmarkdown::render('12_transform.Rmd', 'xaringan::moon_reader') About macros.js: permite escalar las imágenes como [scale 50%](path to image), hay si que grabar ese archivo js en el directorio. ---> .right[] <br> <br> <br> <br> <br> <br> <br> # Estadística multivariada, 1 sem. 2019 ## Juan Carlos Castillo & Alejandro Plaza ## **Sesión 12**: Transformación de variables --- class: inverse https://cran.r-project.org/web/packages/bestNormalize/vignettes/bestNormalize.html https://cran.r-project.org/web/packages/dlookr/vignettes/transformation.html http://www.sthda.com/english/articles/40-regression-analysis/162-nonlinear-regression-essentials-in-r-polynomial-and-spline-regression-models/ # Contenidos ## 1. Repaso supuestos ## 2. Transformación - Polinomios - Logaritmos - Centrado --- class: inverse, middle, center # 1. Repaso supuestos --- class: roja, middle, center ## ¿Qué tan eficiente es mi modelo para dar cuenta de las asociaciones entre las variables? --- # El cuarteto de Anscombe .center[ <!-- --> ] --- # El cuarteto de Anscombe .small[ <table style="text-align:center"><tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"></td><td colspan="4"><em>Dependent variable:</em></td></tr> <tr><td></td><td colspan="4" style="border-bottom: 1px solid black"></td></tr> <tr><td style="text-align:left"></td><td>y1</td><td>y2</td><td>y3</td><td>y4</td></tr> <tr><td style="text-align:left"></td><td>(1)</td><td>(2)</td><td>(3)</td><td>(4)</td></tr> <tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">x1</td><td>0.500^*** (0.118)</td><td></td><td></td><td></td></tr> <tr><td style="text-align:left">x2</td><td></td><td>0.500^*** (0.118)</td><td></td><td></td></tr> <tr><td style="text-align:left">x3</td><td></td><td></td><td>0.500^*** (0.118)</td><td></td></tr> <tr><td style="text-align:left">x4</td><td></td><td></td><td></td><td>0.500^*** (0.118)</td></tr> <tr><td style="text-align:left">Constant</td><td>3.000^** (1.125)</td><td>3.001^** (1.125)</td><td>3.002^** (1.124)</td><td>3.002^** (1.124)</td></tr> <tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">Observations</td><td>11</td><td>11</td><td>11</td><td>11</td></tr> <tr><td style="text-align:left">R²</td><td>0.667</td><td>0.666</td><td>0.666</td><td>0.667</td></tr> <tr><td style="text-align:left">Adjusted R²</td><td>0.629</td><td>0.629</td><td>0.629</td><td>0.630</td></tr> <tr><td style="text-align:left">Residual Std. Error (df = 9)</td><td>1.237</td><td>1.237</td><td>1.236</td><td>1.236</td></tr> <tr><td style="text-align:left">F Statistic (df = 1; 9)</td><td>17.990^***</td><td>17.966^***</td><td>17.972^***</td><td>18.003^***</td></tr> <tr><td colspan="5" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"><em>Note:</em></td><td colspan="4" style="text-align:right">^*p<0.1; ^**p<0.05; ^***p<0.01</td></tr> <tr><td style="text-align:left"></td><td colspan="4" style="text-align:right">Errores Estándares en Parentesis</td></tr> </table> ] --- class: roja ## Por lo tanto: ### No basta con que un beta de regresión sea significativo para poder contrastar adecuadamente las hipótesis -- ### Es necesario revisar la distribución de las asociaciones entre variables y el cumplimiento de supuestos que se encuentran a la base de una estimación de regresión. --- ## ... Recordando: Medias condicionales .center[  ] --- ## ... Recordando: Mínimos cuadrados .center[  ] --- # Posibles problemas ### - subpredicción ### - sobrepredicción ### - sub o sobre predicción diferenciada según valores de X --- class: roja # Detectando problemas ### Revisar VISUALMENTE la distribución de los residuos alrededor de la recta de regresión. ### La forma más eficiente es analizar la relación entre valor predicho (recta) y residuos. --- ## Gráficos de Diagnóstico .center[  ] --- ## Chequeo de supuestos -- - **Linealidad**: distribución aleatoria de residos alrededor de linea residual=0. . -- - **Homocedasticidad**: además de aleatoria, los residuos deberían estar distribuidos de manera constante en la banda alrededor de la linea residual=0 para los distintos valores estimados. Formas de abanico, o de embudo indicaría que la varianza de los errores no es constante. -- - **Normalidad de residuos** -- - **Observaciones influyentes**: un outlier es una observación que se separa mucho del resto y que puede afectar la estimación -- - **Multicolinealidad**: asociación estrecha entre predictores -- (ver cuadro resumen en Guía 8) --- class: inverse, middle, center # 2. Transformaciones --- # Modelos lineales de relaciones no-lineales ### - Uno de los supuestos principales del modelo de regresión es la relación lineal (recta) entre las variables. -- ### - Esta linealidad se expresa en una distribución aleatoria de los residuos en relación a los valores estimados -- ## - ¿Qué sucede en el caso de relaciones no lineales? -- ### - Para poder capturar eficientemente relaciones no lineales entre las variables y evitar la violación de supuestos una alternativa es la TRANSFORMACIÓN de las variables --- # Ejemplo Relaciones no lineales .center[  ] ## ¿Es posible modelar más apropiadamente esta relación? --- # Tipos usuales de transformación ## - Polinomial ## - Logarítmica --- ## Transformación polinomial - Usa regresores que son transformaciones de potencias sucesivas, tales como `\(X\)`, `\(X^2\)`, `\(X^3\)` .center[  ] --- ## Polinomial `$$Y=\beta_0+\beta_{1}X+\beta_{2}X^2$$` - Polinomio cuadrático o de segundo orden, permite una curva en la relación entre X e Y -- `$$Y=\beta_0+\beta_{1}X+\beta_{2}X^2+\beta_{2}X^3$$` - Polinomio cúbico o de tercer orden, permite dos curvas en la relación entre X e Y -- ## Implica que el aumento de una unidad de X tendrá un efecto diferente en Y dependiendo de dónde se encuentra el valor de X --- # Ejemplo .large[ - Base de datos `Boston`: incluye variables. - `mdev`: valor promedio de las viviendas por barrio - `lstat` : porcentaje de personas de bajo estatus ] - Cargar datos: para realizar los ejercicios de esta sesión, cargar los datos así: ```r load(url(https://juancarloscastillo.github.io/ metsoc-facsouchile/documents/data/boston.rda)) ``` --- ## Scatter medv y lstat ```r ggplot(boston, aes(lstat, medv)) + geom_point() + stat_smooth() ``` <!-- --> --- ## Modelo .pull-left[ ```r # modelo model1 <- lm(medv ~ lstat, data = boston) ``` ] .pull-right[ <table style="text-align:center"><tr><td colspan="2" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"></td><td>medv</td></tr> <tr><td></td><td colspan="1" style="border-bottom: 1px solid black"></td></tr> <tr><td colspan="2" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">lstat</td><td>-0.956^***</td></tr> <tr><td style="text-align:left"></td><td>(0.043)</td></tr> <tr><td style="text-align:left"></td><td></td></tr> <tr><td style="text-align:left">Constant</td><td>34.653^***</td></tr> <tr><td style="text-align:left"></td><td>(0.623)</td></tr> <tr><td style="text-align:left"></td><td></td></tr> <tr><td colspan="2" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">Observations</td><td>407</td></tr> <tr><td style="text-align:left">Adjusted R²</td><td>0.551</td></tr> <tr><td style="text-align:left">F Statistic</td><td>499.157^*** (df = 1; 405)</td></tr> <tr><td colspan="2" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"><em>Note:</em></td><td style="text-align:right">^*p<0.1; ^**p<0.05; ^***p<0.01</td></tr> </table> ] --- ## Modelo ```r ggplot(train.data, aes(lstat, medv) ) + geom_point() + stat_smooth(method = lm, formula = y ~ x) ``` <!-- --> --- ## Análisis linealidad - agregar variables para análisis de supuestos a la base de datos ```r boston_s <- broom::augment_columns(model1, boston) names(boston_s) ``` ``` ## [1] ".rownames" "medv" "lstat" ".fitted" ".se.fit" ## [6] ".resid" ".hat" ".sigma" ".cooksd" ".std.resid" ``` --- ## Análisis linealidad ```r ggplot(boston_s, aes(x=.fitted, y=.resid)) + geom_hline(yintercept=0) + geom_point() + geom_smooth(method='loess', se=TRUE) ``` <!-- --> --- ## Inclusión de término cuadrático ( `\(+\beta_{2}lstat^2\)` ) .pull-left[ ```r model2 <- lm(medv ~ lstat + I(lstat^2), data = boston) ``` ] .pull-right[.medium[ <table style="text-align:center"><tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"></td><td colspan="2">medv</td></tr> <tr><td></td><td colspan="2" style="border-bottom: 1px solid black"></td></tr> <tr><td style="text-align:left"></td><td>(1)</td><td>(2)</td></tr> <tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">lstat</td><td>-0.956^***</td><td>-2.267^***</td></tr> <tr><td style="text-align:left"></td><td>(0.043)</td><td>(0.136)</td></tr> <tr><td style="text-align:left"></td><td></td><td></td></tr> <tr><td style="text-align:left">I(lstat2)</td><td></td><td>0.041^***</td></tr> <tr><td style="text-align:left"></td><td></td><td>(0.004)</td></tr> <tr><td style="text-align:left"></td><td></td><td></td></tr> <tr><td style="text-align:left">Constant</td><td>34.653^***</td><td>42.574^***</td></tr> <tr><td style="text-align:left"></td><td>(0.623)</td><td>(0.965)</td></tr> <tr><td style="text-align:left"></td><td></td><td></td></tr> <tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">Observations</td><td>407</td><td>407</td></tr> <tr><td style="text-align:left">Adjusted R²</td><td>0.551</td><td>0.640</td></tr> <tr><td style="text-align:left">F Statistic</td><td>499.157^*** (df = 1; 405)</td><td>361.940^*** (df = 2; 404)</td></tr> <tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"><em>Note:</em></td><td colspan="2" style="text-align:right">^*p<0.1; ^**p<0.05; ^***p<0.01</td></tr> </table> ] ] --- ## Comparación ajuste modelos .large[ ```r anova(model1,model2) ``` ``` ## Analysis of Variance Table ## ## Model 1: medv ~ lstat ## Model 2: medv ~ lstat + I(lstat^2) ## Res.Df RSS Df Sum of Sq F Pr(>F) ## 1 405 15286 ## 2 404 12224 1 3062.3 101.21 < 2.2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ``` - por lo tanto, se rechaza la hipótesis nula de igualdad de ajuste de los modelos. El modelo con el término cuadrático es significativamente superior al anterior. ] --- ## Análisis residuos vs predichos model2 ```r ggplot(boston_s2, aes(x=.fitted, y=.resid)) + geom_hline(yintercept=0) + geom_point() + geom_smooth(method='loess', se=TRUE) ``` <!-- --> --- ## Inclusión de término cúbico ( `\(+\beta_{3}lstat^3\)` ) .medium[ <table style="text-align:center"><tr><td colspan="4" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"></td><td colspan="3">medv</td></tr> <tr><td></td><td colspan="3" style="border-bottom: 1px solid black"></td></tr> <tr><td style="text-align:left"></td><td>(1)</td><td>(2)</td><td>(3)</td></tr> <tr><td colspan="4" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">lstat</td><td>-0.956^*** (0.043)</td><td>-2.267^*** (0.136)</td><td>-3.841^*** (0.370)</td></tr> <tr><td style="text-align:left">I(lstat2)</td><td></td><td>0.041^*** (0.004)</td><td>0.148^*** (0.024)</td></tr> <tr><td style="text-align:left">I(lstat3)</td><td></td><td></td><td>-0.002^*** (0.0004)</td></tr> <tr><td style="text-align:left">Constant</td><td>34.653^*** (0.623)</td><td>42.574^*** (0.965)</td><td>48.574^*** (1.619)</td></tr> <tr><td colspan="4" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">Observations</td><td>407</td><td>407</td><td>407</td></tr> <tr><td style="text-align:left">Adjusted R²</td><td>0.551</td><td>0.640</td><td>0.657</td></tr> <tr><td style="text-align:left">F Statistic</td><td>499.157^*** (df = 1; 405)</td><td>361.940^*** (df = 2; 404)</td><td>260.030^*** (df = 3; 403)</td></tr> <tr><td colspan="4" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"><em>Note:</em></td><td colspan="3" style="text-align:right">^*p<0.1; ^**p<0.05; ^***p<0.01</td></tr> </table> ] --- ## Comparación ajuste modelos ```r anova(model1,model2,model3) ``` ``` ## Analysis of Variance Table ## ## Model 1: medv ~ lstat ## Model 2: medv ~ lstat + I(lstat^2) ## Model 3: medv ~ lstat + I(lstat^2) + I(lstat^3) ## Res.Df RSS Df Sum of Sq F Pr(>F) ## 1 405 15286 ## 2 404 12224 1 3062.33 106.167 < 2.2e-16 *** ## 3 403 11624 1 599.27 20.776 6.855e-06 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ``` - por lo tanto, se rechaza la hipótesis nula de igualdad de ajuste de los modelos. El modelo con el término cúbico es significativamente superior al anterior. --- ## Análisis residuos vs valores predichos model3 <!-- --> --- ## Graficando predicción polinomial cúbica ```r ggplot(boston, aes(lstat, medv) ) + geom_point() + stat_smooth(method = lm, formula = y ~ poly(x, 3, raw = TRUE)) ``` <!-- --> --- ## Explorando polinomiales ```r lm(medv ~ poly(lstat, 6, raw = TRUE), data = boston) %>% summary() ``` ``` ## ## Call: ## lm(formula = medv ~ poly(lstat, 6, raw = TRUE), data = boston) ## ## Residuals: ## Min 1Q Median 3Q Max ## -13.1962 -3.1527 -0.7655 2.0404 26.7661 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 7.788e+01 6.844e+00 11.379 < 2e-16 *** ## poly(lstat, 6, raw = TRUE)1 -1.767e+01 3.569e+00 -4.952 1.08e-06 *** ## poly(lstat, 6, raw = TRUE)2 2.417e+00 6.779e-01 3.566 0.000407 *** ## poly(lstat, 6, raw = TRUE)3 -1.761e-01 6.105e-02 -2.885 0.004121 ** ## poly(lstat, 6, raw = TRUE)4 6.845e-03 2.799e-03 2.446 0.014883 * ## poly(lstat, 6, raw = TRUE)5 -1.343e-04 6.290e-05 -2.136 0.033323 * ## poly(lstat, 6, raw = TRUE)6 1.047e-06 5.481e-07 1.910 0.056910 . ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 5.188 on 400 degrees of freedom ## Multiple R-squared: 0.6845, Adjusted R-squared: 0.6798 ## F-statistic: 144.6 on 6 and 400 DF, p-value: < 2.2e-16 ``` --- class: inverse # Resumen transformación polinomial .large[ - alternativa en casos de no linealidad - agregar polinomios secuencialmente y realizar test de ajuste - interpretación sustantiva de los efectos ] --- # Transformación logarítmica .large[ - Posibilidad adicional de transformación para asociaciones no lineales - Ej: efecto de ingreso, ¿es lo mismo una diferencia de 100.000 para alguien que gana 400.000 que para alguien que gana 5.000.000? - La transformación logarítmica permite incorporar mayor "peso" a diferencias en el extremo inferior de la escala. - Por ejemplo, los logaritmos de 10, 100, y 1,000 son, respectivamente, 1, 2, y3, porque `\(10^1=10\)`, `\(10^2=100\)`, `\(10^3=1000\)`. ] --- # Transformación logarítmica .center[   ] --- ## Modelo con predictor logarítmico ```r model_log <- lm(medv ~ log(lstat), data = boston) ``` <table style="text-align:center"><tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"></td><td colspan="2">medv</td></tr> <tr><td></td><td colspan="2" style="border-bottom: 1px solid black"></td></tr> <tr><td style="text-align:left"></td><td>(1)</td><td>(2)</td></tr> <tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">lstat</td><td>-0.956^*** (0.043)</td><td></td></tr> <tr><td style="text-align:left">log(lstat)</td><td></td><td>-12.567^*** (0.441)</td></tr> <tr><td style="text-align:left">Constant</td><td>34.653^*** (0.623)</td><td>52.380^*** (1.081)</td></tr> <tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left">Observations</td><td>407</td><td>407</td></tr> <tr><td style="text-align:left">Adjusted R²</td><td>0.551</td><td>0.666</td></tr> <tr><td style="text-align:left">F Statistic (df = 1; 405)</td><td>499.157^***</td><td>812.064^***</td></tr> <tr><td colspan="3" style="border-bottom: 1px solid black"></td></tr><tr><td style="text-align:left"><em>Note:</em></td><td colspan="2" style="text-align:right">^*p<0.1; ^**p<0.05; ^***p<0.01</td></tr> </table> --- ## Residuos vs valores predichos .pull-left[ ### Original <!-- --> ] .pull-right[ ### Log <!-- --> ] --- ## Ajuste regresión log ```r ggplot(boston, aes(lstat, medv) ) + geom_point() + stat_smooth(method = lm, formula = y ~ log(x)) ``` <!-- --> --- class: inverse, middle, center # Notas sobre centrado de variables --- ## Centrado .large[ - El centrado de las variables en regresión simple posee consecuencias en el valor del intercepto - Intercepto: estimación del valor de la variable dependiente cuando los predictores son 0 - La estimación sin centrar en algunos casos puede llevar a algunos problemas en la interpretación de los parámetros del modelo - El centrado se implementa restando una constante del predictor (usualmente el promedio ) ] --- ## Ejemplo centrado regresión simple ### Datos: puntaje en CI y en una escala de felicidad 1-7 .large[ ``` > mydata felicidad ci 1 1 90 2 3 95 3 4 98 4 6 105 lm(formula = felicidad ~ ci) (Intercept) ci -28.5593 0.3305 ``` ] -- ### **¿Qué significa el valor del intercepto?** --- ## Ejemplo centrado regresión simple <br> <br>  --- ## Ejemplo centrado regresión simple .large[ ``` > mydata felicidad ci mean_ci cicent_mean cicent_100 1 1 90 97 -7 -10 2 3 95 97 -2 -5 3 4 98 97 1 -2 4 6 105 97 8 5 lm(formula = felicidad ~ cicent_mean) (Intercept): 3.5000 ; cicent_mean:0.3305 lm(formula = felicidad ~ cicent_100) (Intercept): 4.4915 ; mydata$cicent_100: 0.3305 ``` ] --- ## Ejemplo centrado regresión simple <br> <br>  --- ## Resumen centrado .large[ - `\(\beta\)` constante en todos los modelos, varı́a solo intercepto - Modelo sin centrar: intercepto -28,5, felicidad para alguien con inteligencia 0 (no tiene mucho sentido ...) - Centrado al promedio: 3.5, estimación de la felicidad para alguien con inteligencia promedio - Centrado a 100: 4.49, estimación de la felicidad para alguien con inteligencia 100 ] --- class: bottom, left .right[] <br> <br> <br> <br> <br> <br> <br> # Estadística multivariada, 1 sem. 2019 ## Juan Carlos Castillo & Alejandro Plaza ## **Sesión 12**: Transformación de variables