Objetivo

Crear, evaluar y comparar un modelo de regresión lineal simple y un moldeo polinómico para predecir ???

Descripción

Fundamento teórico

La regresión lineal simple apoya a la correlación estimando coeficientes para hacer predicciones.

\[ Y = a + bx \]

\[ Y = \beta_0 + \beta_1\cdot x_i + \epsilon \]

\[ Y = \beta0 + \beta_1\cdot x^2+\beta_2\cdot x^3 + ... \beta_n \cdot x^n \]

Desarrollo

Cargar librerías

library(ggplot2)
# library(plotly) # no se está usando
library(knitr)
library(PerformanceAnalytics) # Para coorelaciones gráficas
library(caret)  # Para particionar
library(Metrics) # Para determinar rmse

Cargar datos

datos <- read.csv("https://raw.githubusercontent.com/rpizarrog/Machine-Learning-con-R/main/datos/Advertising.csv")

Explorar datos

str(datos)
## 'data.frame':    200 obs. of  5 variables:
##  $ X        : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ TV       : num  230.1 44.5 17.2 151.5 180.8 ...
##  $ Radio    : num  37.8 39.3 45.9 41.3 10.8 48.9 32.8 19.6 2.1 2.6 ...
##  $ Newspaper: num  69.2 45.1 69.3 58.5 58.4 75 23.5 11.6 1 21.2 ...
##  $ Sales    : num  22.1 10.4 9.3 18.5 12.9 7.2 11.8 13.2 4.8 10.6 ...
summary(datos)
##        X                TV             Radio          Newspaper     
##  Min.   :  1.00   Min.   :  0.70   Min.   : 0.000   Min.   :  0.30  
##  1st Qu.: 50.75   1st Qu.: 74.38   1st Qu.: 9.975   1st Qu.: 12.75  
##  Median :100.50   Median :149.75   Median :22.900   Median : 25.75  
##  Mean   :100.50   Mean   :147.04   Mean   :23.264   Mean   : 30.55  
##  3rd Qu.:150.25   3rd Qu.:218.82   3rd Qu.:36.525   3rd Qu.: 45.10  
##  Max.   :200.00   Max.   :296.40   Max.   :49.600   Max.   :114.00  
##      Sales      
##  Min.   : 1.60  
##  1st Qu.:10.38  
##  Median :12.90  
##  Mean   :14.02  
##  3rd Qu.:17.40  
##  Max.   :27.00

Visualizar datos

Diagrama de dispersión

source("https://raw.githubusercontent.com/rpizarrog/Machine-Learning-con-R/main/funciones/simular_women")

f_dispersion(data.frame(datos$TV, datos$Sales))

Las variables de interés

  • x la variable independiente o predictora TV
  • y la variable dependiente o resultado es Sales, es decir y depende de x.

Seleccionar las columnas de interes

datos <- select(datos, TV, Sales)

Limpiar datos

En caso necesario

Partir datos

Aleatoriamente se repate las observaciones con el 70% para datos de entrenamiento y el 30% para datos de validación.

Sembrar una semilla con set.seed()

set.seed(2022)
n <- nrow(datos)  # cantidad de observaciones
entrena <- createDataPartition(y = datos$Sales, p = 0.70, list = FALSE, times = 1)

# Datos entrenamiento
datos.entrenamiento <- datos[entrena, ]  # [renglones, columna]

# Datos validación
datos.validacion <- datos[-entrena, ]

Datos de entrenamiento

datos.entrenamiento
##        TV Sales
## 1   230.1  22.1
## 2    44.5  10.4
## 3    17.2   9.3
## 4   151.5  18.5
## 5   180.8  12.9
## 6     8.7   7.2
## 7    57.5  11.8
## 8   120.2  13.2
## 9     8.6   4.8
## 11   66.1   8.6
## 12  214.7  17.4
## 13   23.8   9.2
## 14   97.5   9.7
## 15  204.1  19.0
## 16  195.4  22.4
## 17   67.8  12.5
## 18  281.4  24.4
## 19   69.2  11.3
## 23   13.2   5.6
## 25   62.3   9.7
## 28  240.1  15.9
## 29  248.8  18.9
## 32  112.9  11.9
## 37  266.9  25.4
## 38   74.7  14.7
## 39   43.1  10.1
## 40  228.0  21.5
## 41  202.5  16.6
## 43  293.6  20.7
## 44  206.9  12.9
## 45   25.1   8.5
## 46  175.1  14.9
## 47   89.7  10.6
## 49  227.2  14.8
## 51  199.8  11.4
## 52  100.4  10.7
## 53  216.4  22.6
## 56  198.9  23.7
## 57    7.3   5.5
## 58  136.2  13.2
## 59  210.8  23.8
## 61   53.5   8.1
## 62  261.3  24.2
## 64  102.7  14.0
## 65  131.1  18.0
## 68  139.3  13.4
## 70  216.8  22.3
## 71  199.1  18.3
## 72  109.8  12.4
## 73   26.8   8.8
## 74  129.4  11.0
## 75  213.4  17.0
## 77   27.5   6.9
## 78  120.5  14.2
## 79    5.4   5.3
## 80  116.0  11.0
## 81   76.4  11.8
## 82  239.8  12.3
## 85  213.5  21.7
## 86  193.2  15.2
## 88  110.7  16.0
## 89   88.3  12.9
## 91  134.3  11.2
## 92   28.6   7.3
## 94  250.9  22.2
## 96  163.3  16.9
## 97  197.6  11.7
## 98  184.9  15.5
## 99  289.7  25.4
## 100 135.2  17.2
## 102 296.4  23.8
## 103 280.2  14.8
## 105 238.2  20.7
## 106 137.9  19.2
## 108  90.4   8.7
## 110 255.4  19.8
## 113 175.7  14.1
## 115  78.2  14.6
## 117 139.2  12.2
## 119 125.7  15.9
## 121 141.3  15.5
## 122  18.8   7.0
## 123 224.0  11.6
## 124 123.1  15.2
## 126  87.2  10.6
## 127   7.8   6.6
## 129 220.3  24.7
## 131   0.7   1.6
## 132 265.2  12.7
## 134 219.8  19.6
## 135  36.9  10.8
## 136  48.3  11.6
## 137  25.6   9.5
## 138 273.7  20.8
## 141  73.4  10.9
## 142 193.7  19.2
## 143 220.5  20.1
## 145  96.2  11.4
## 146 140.3  10.3
## 147 240.1  13.2
## 148 243.2  25.4
## 150  44.7  10.1
## 152 121.0  11.6
## 154 171.3  19.0
## 155 187.8  15.6
## 156   4.1   3.2
## 157  93.9  15.3
## 158 149.8  10.1
## 159  11.7   7.3
## 160 131.7  12.9
## 161 172.5  14.4
## 162  85.7  13.3
## 163 188.4  14.9
## 164 163.5  18.0
## 165 117.2  11.9
## 166 234.5  11.9
## 167  17.9   8.0
## 168 206.8  12.2
## 169 215.4  17.1
## 170 284.3  15.0
## 173  19.6   7.6
## 174 168.4  11.7
## 175 222.4  11.5
## 176 276.9  27.0
## 179 276.7  11.8
## 180 165.6  12.6
## 181 156.6  10.5
## 183  56.2   8.7
## 184 287.6  26.2
## 185 253.8  17.6
## 186 205.0  22.6
## 187 139.5  10.3
## 188 191.1  17.3
## 189 286.0  15.9
## 190  18.7   6.7
## 192  75.5   9.9
## 193  17.2   5.9
## 195 149.7  17.3
## 196  38.2   7.6
## 197  94.2   9.7
## 198 177.0  12.8
## 200 232.1  13.4

Datos de validación

datos.validacion
##        TV Sales
## 10  199.8  10.6
## 20  147.3  14.6
## 21  218.4  18.0
## 22  237.4  12.5
## 24  228.3  15.5
## 26  262.9  12.0
## 27  142.9  15.0
## 30   70.6  10.5
## 31  292.9  21.4
## 33   97.2   9.6
## 34  265.6  17.4
## 35   95.7   9.5
## 36  290.7  12.8
## 42  177.0  17.1
## 48  239.9  23.2
## 50   66.9   9.7
## 54  182.6  21.2
## 55  262.7  20.2
## 60  210.7  18.4
## 63  239.3  15.7
## 66   69.0   9.3
## 67   31.5   9.5
## 69  237.4  18.9
## 76   16.9   8.7
## 83   75.3  11.3
## 84   68.4  13.6
## 87   76.3  12.0
## 90  109.8  16.7
## 93  217.7  19.4
## 95  107.4  11.5
## 101 222.4  11.7
## 104 187.9  14.7
## 107  25.0   7.2
## 109  13.1   5.3
## 111 225.8  13.4
## 112 241.7  21.8
## 114 209.6  15.9
## 116  75.1  12.6
## 118  76.4   9.4
## 120  19.4   6.6
## 125 229.5  19.7
## 128  80.2   8.8
## 130  59.6   9.7
## 133   8.4   5.7
## 139  43.0   9.6
## 140 184.9  20.7
## 144 104.6  10.4
## 149  38.0  10.9
## 151 280.7  16.1
## 153 197.6  16.6
## 171  50.0   8.4
## 172 164.5  14.5
## 177 248.4  20.2
## 178 170.2  11.7
## 182 218.5  12.2
## 191  39.5  10.8
## 194 166.8  19.6
## 199 283.6  25.5

Construir los modelo

El modelo se construye con datos de entrenamiento ### Modelo de Regresión lineal simple

modelo.rs <- lm(data = datos.entrenamiento, 
                 formula = Sales ~ TV)
summary(modelo.rs)
## 
## Call:
## lm(formula = Sales ~ TV, data = datos.entrenamiento)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.8752 -1.7706 -0.2902  2.1879  6.9256 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 6.801923   0.551158   12.34   <2e-16 ***
## TV          0.050138   0.003278   15.30   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.318 on 140 degrees of freedom
## Multiple R-squared:  0.6256, Adjusted R-squared:  0.623 
## F-statistic:   234 on 1 and 140 DF,  p-value: < 2.2e-16

Modelo de Regresión polinomial segunda potencia

modelo.rpoly2 <- lm(data = datos.entrenamiento, 
                 formula = Sales ~ poly(TV, degree = 2))
summary(modelo.rpoly2)
## 
## Call:
## lm(formula = Sales ~ poly(TV, degree = 2), data = datos.entrenamiento)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.3015 -1.7900 -0.2512  2.1656  6.8913 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)            14.0775     0.2779  50.664   <2e-16 ***
## poly(TV, degree = 2)1  50.7569     3.3111  15.330   <2e-16 ***
## poly(TV, degree = 2)2  -4.1996     3.3111  -1.268    0.207    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.311 on 139 degrees of freedom
## Multiple R-squared:  0.6299, Adjusted R-squared:  0.6246 
## F-statistic: 118.3 on 2 and 139 DF,  p-value: < 2.2e-16

Modelo de Regresión polinomial quinta potencia

modelo.rpoly5 <- lm(data = datos.entrenamiento, 
                 formula = Sales ~ poly(TV, degree = 5))
summary(modelo.rpoly5)
## 
## Call:
## lm(formula = Sales ~ poly(TV, degree = 5), data = datos.entrenamiento)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.2657 -1.7632 -0.1326  1.8232  6.9183 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             14.078      0.276  51.001   <2e-16 ***
## poly(TV, degree = 5)1   50.757      3.289  15.431   <2e-16 ***
## poly(TV, degree = 5)2   -4.200      3.289  -1.277    0.204    
## poly(TV, degree = 5)3    4.261      3.289   1.295    0.197    
## poly(TV, degree = 5)4   -3.294      3.289  -1.002    0.318    
## poly(TV, degree = 5)5    4.845      3.289   1.473    0.143    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.289 on 136 degrees of freedom
## Multiple R-squared:  0.6427, Adjusted R-squared:  0.6295 
## F-statistic: 48.92 on 5 and 136 DF,  p-value: < 2.2e-16

Evaluar el modelo antes de predicciones

En cuanto a a métrica R Square que significa la represntatividad estadística de la variable independiete sobre la variable dependiente (TV sobre Sales), el mejor modelo sería el polinómioc de segndo nivel dado que :

  • Polinomial de segundo nivel: Multiple R-squared: 0.6256
  • Lineal simple: Multiple R-squared: 0.6299

Hacer predicciones

Modelo Reg. Lineal Simple

Se usa la función predcit() para hacer predicciones con los datos de validación luego se compara con los datos reales.

predicciones <- predict(modelo.rs, newdata = datos.validacion)
predicciones
##        10        20        21        22        24        26        27        30 
## 16.819528 14.187274 17.752098 18.704723 18.248466 19.983246 13.966667 10.341677 
##        31        33        34        35        36        42        48        50 
## 21.487391 11.675352 20.118619 11.600145 21.377087 15.676378 18.830069 10.156166 
##        54        55        60        63        66        67        69        76 
## 15.957152 19.973219 17.366034 18.799986 10.261456  8.381275 18.704723  7.649258 
##        83        84        87        90        93        95       101       104 
## 10.577326 10.231373 10.627465 12.307093 17.717001 12.186762 17.952651 16.222884 
##       107       109       111       112       114       116       118       120 
##  8.055377  7.458733 18.123120 18.920317 17.310882 10.567299 10.632478  7.774603 
##       125       128       130       133       139       140       144       149 
## 18.308632 10.823003  9.790157  7.223083  8.957864 16.072469 12.046375  8.707173 
##       151       153       171       172       177       178       182       191 
## 20.875706 16.709224  9.308831 15.049651 19.256243 15.335438 17.757112  8.782380 
##       194       199 
## 15.164969 21.021106

Construir un data frame para comparar

comparaciones <- data.frame(datos.validacion, predicciones)
comparaciones
##        TV Sales predicciones
## 10  199.8  10.6    16.819528
## 20  147.3  14.6    14.187274
## 21  218.4  18.0    17.752098
## 22  237.4  12.5    18.704723
## 24  228.3  15.5    18.248466
## 26  262.9  12.0    19.983246
## 27  142.9  15.0    13.966667
## 30   70.6  10.5    10.341677
## 31  292.9  21.4    21.487391
## 33   97.2   9.6    11.675352
## 34  265.6  17.4    20.118619
## 35   95.7   9.5    11.600145
## 36  290.7  12.8    21.377087
## 42  177.0  17.1    15.676378
## 48  239.9  23.2    18.830069
## 50   66.9   9.7    10.156166
## 54  182.6  21.2    15.957152
## 55  262.7  20.2    19.973219
## 60  210.7  18.4    17.366034
## 63  239.3  15.7    18.799986
## 66   69.0   9.3    10.261456
## 67   31.5   9.5     8.381275
## 69  237.4  18.9    18.704723
## 76   16.9   8.7     7.649258
## 83   75.3  11.3    10.577326
## 84   68.4  13.6    10.231373
## 87   76.3  12.0    10.627465
## 90  109.8  16.7    12.307093
## 93  217.7  19.4    17.717001
## 95  107.4  11.5    12.186762
## 101 222.4  11.7    17.952651
## 104 187.9  14.7    16.222884
## 107  25.0   7.2     8.055377
## 109  13.1   5.3     7.458733
## 111 225.8  13.4    18.123120
## 112 241.7  21.8    18.920317
## 114 209.6  15.9    17.310882
## 116  75.1  12.6    10.567299
## 118  76.4   9.4    10.632478
## 120  19.4   6.6     7.774603
## 125 229.5  19.7    18.308632
## 128  80.2   8.8    10.823003
## 130  59.6   9.7     9.790157
## 133   8.4   5.7     7.223083
## 139  43.0   9.6     8.957864
## 140 184.9  20.7    16.072469
## 144 104.6  10.4    12.046375
## 149  38.0  10.9     8.707173
## 151 280.7  16.1    20.875706
## 153 197.6  16.6    16.709224
## 171  50.0   8.4     9.308831
## 172 164.5  14.5    15.049651
## 177 248.4  20.2    19.256243
## 178 170.2  11.7    15.335438
## 182 218.5  12.2    17.757112
## 191  39.5  10.8     8.782380
## 194 166.8  19.6    15.164969
## 199 283.6  25.5    21.021106

Evaluar predicciones

rmse Root Mean Stándar Error, este valor normalmente se compara contra otro modelo y el que esté mas cerca de cero es mejor.

RMSE es una forma útil de ver qué tan bien un modelo de regresión puede ajustarse a un conjunto de datos.

Cuanto mayor sea el rmse, mayor será la diferencia entre los valores predichos y reales, lo que significa que peor se ajusta un modelo de regresión a los datos. Por el contrario, cuanto más pequeño sea el rmse, mejor podrá un modelo ajustar los datos.

\[ rmse = \sqrt{\frac{\sum(predicho_i - real_i)^{2}}{n}} \]

rmse <- rmse(actual = comparaciones$TV, predicted = comparaciones$predicciones)
rmse
## [1] 160.3604

Modelo Reg. Polinomial de segundo nivel

Se usa la función predcit() para hacer predicciones con los datos de validación luego se compara con los datos reales.

predicciones <- predict(modelo.rpoly2, newdata = datos.validacion)
predicciones
##        10        20        21        22        24        26        27        30 
## 17.028317 14.569565 17.829207 18.609391 18.240507 19.596234 14.350206 10.451326 
##        31        33        34        35        36        42        48        50 
## 20.668829 11.950308 19.696681 11.867778 20.593419 15.996468 18.709193 10.236869 
##        54        55        60        63        66        67        69        76 
## 16.255018 19.588763 17.502112 18.685301 10.358766  8.111550 18.609391  7.196249 
##        83        84        87        90        93        95       101       104 
## 10.721648 10.323986 10.778861 12.634129 17.799731 12.505177 17.996642 16.496651 
##       107       109       111       112       114       116       118       120 
##  7.706848  6.954308 18.137625 18.780640 17.454870 10.710193 10.784577  7.354584 
##       125       128       130       133       139       140       144       149 
## 18.289654 11.000977  9.809488  6.652944  8.816571 16.360244 12.353960  8.511765 
##       151       153       171       172       177       178       182       191 
## 20.244170 16.931159  9.238838 15.407334 19.043556 15.678037 17.833413  8.603486 
##       194       199 
## 15.516980 20.346545

Construir un data frame para comparar

comparaciones <- data.frame(datos.validacion, predicciones)
comparaciones
##        TV Sales predicciones
## 10  199.8  10.6    17.028317
## 20  147.3  14.6    14.569565
## 21  218.4  18.0    17.829207
## 22  237.4  12.5    18.609391
## 24  228.3  15.5    18.240507
## 26  262.9  12.0    19.596234
## 27  142.9  15.0    14.350206
## 30   70.6  10.5    10.451326
## 31  292.9  21.4    20.668829
## 33   97.2   9.6    11.950308
## 34  265.6  17.4    19.696681
## 35   95.7   9.5    11.867778
## 36  290.7  12.8    20.593419
## 42  177.0  17.1    15.996468
## 48  239.9  23.2    18.709193
## 50   66.9   9.7    10.236869
## 54  182.6  21.2    16.255018
## 55  262.7  20.2    19.588763
## 60  210.7  18.4    17.502112
## 63  239.3  15.7    18.685301
## 66   69.0   9.3    10.358766
## 67   31.5   9.5     8.111550
## 69  237.4  18.9    18.609391
## 76   16.9   8.7     7.196249
## 83   75.3  11.3    10.721648
## 84   68.4  13.6    10.323986
## 87   76.3  12.0    10.778861
## 90  109.8  16.7    12.634129
## 93  217.7  19.4    17.799731
## 95  107.4  11.5    12.505177
## 101 222.4  11.7    17.996642
## 104 187.9  14.7    16.496651
## 107  25.0   7.2     7.706848
## 109  13.1   5.3     6.954308
## 111 225.8  13.4    18.137625
## 112 241.7  21.8    18.780640
## 114 209.6  15.9    17.454870
## 116  75.1  12.6    10.710193
## 118  76.4   9.4    10.784577
## 120  19.4   6.6     7.354584
## 125 229.5  19.7    18.289654
## 128  80.2   8.8    11.000977
## 130  59.6   9.7     9.809488
## 133   8.4   5.7     6.652944
## 139  43.0   9.6     8.816571
## 140 184.9  20.7    16.360244
## 144 104.6  10.4    12.353960
## 149  38.0  10.9     8.511765
## 151 280.7  16.1    20.244170
## 153 197.6  16.6    16.931159
## 171  50.0   8.4     9.238838
## 172 164.5  14.5    15.407334
## 177 248.4  20.2    19.043556
## 178 170.2  11.7    15.678037
## 182 218.5  12.2    17.833413
## 191  39.5  10.8     8.603486
## 194 166.8  19.6    15.516980
## 199 283.6  25.5    20.346545

Evaluar predicciones

rmse Root Mean Stándar Error, este valor normalmente se compara contra otro modelo y el que esté mas cerca de cero es mejor.

RMSE es una forma útil de ver qué tan bien un modelo de regresión puede ajustarse a un conjunto de datos.

Cuanto mayor sea el rmse, mayor será la diferencia entre los valores predichos y reales, lo que significa que peor se ajusta un modelo de regresión a los datos. Por el contrario, cuanto más pequeño sea el rmse, mejor podrá un modelo ajustar los datos.

rmse <- rmse(actual = comparaciones$TV, predicted = comparaciones$predicciones)
rmse
## [1] 160.4068

Dispersión con lineas de tendecias

Unicamente se reflejan las tendencias lineal simple y polinómico de segundo nivel.

ggplot(data = datos.entrenamiento) +
  geom_point(aes(x = TV, y = Sales), col='blue') +
  geom_line(aes(x = TV, y = modelo.rs$fitted.values), col='red') +
  geom_line(aes(x = TV, y = modelo.rpoly2$fitted.values), col='green')

Predicciones con datos nuevos

Interpretración

Bibliografía