Regresión en R: Una comparación

En este artículo, proporciono varios ejemplos de regresión en R. Comparo varios métodos de regresión usando R. Los métodos que presentaré son regresión lineal simple, regresión lineal múltiple, regresión polinomial, regresión de árbol de decisión y regresión de máquina de vectores de soporte.

Paquetes necesarios para la regresión en R

Comienzo cargando los paquetes necesarios (además de la biblioteca base de R):

library(dplyr)         # package for data wrangling
library(ggplot2)       # package for data visualization
library(rpart)         # package for decision tree regression analysis
library(randomForest)  # package for random forest regression analysis
library(e1071)         # package for support vector regression analysis
library(caTools)       # package for e.g. splitting into training and test sets

Lectura de los datos antes de la regresión en R

Ahora leo datos que usaré para implementar mis ejemplos de regresión en R. Los datos contienen datos sobre el producto interno bruto, las inversiones privadas brutas, la esperanza de vida, el tamaño de la población y las inversiones en investigación y desarrollo.

# reading in gross domestic product data
gdp_df <- read.csv("gross domestic product.csv",header=TRUE)

# reading in gross private investment data
investment_df <- read.csv("gross private investments.csv",header=TRUE)

# reading in life expectancy data
lifeexpectancy_df <- read.csv("life expectancy.csv",header=TRUE)

# reading in population data
population_df <- read.csv("population.csv",header=TRUE)

# reading in research and development data
rnd_df <- read.csv("research and development.csv",header=TRUE)

# joining all data into a final input data set
data_df <- gdp_df %>% inner_join(investment_df,by="DATE",fa) %>% 
  inner_join(lifeexpectancy_df,by="DATE") %>%
  inner_join(population_df,by="DATE") %>%
  inner_join(rnd_df,by="DATE")

# rename column headers
colnames(data_df) <-c("date",
                      "gdp",
                      "private_investment",
                      "life_expectancy",
                      "population",
                      "rnd")

# convert date entries from characters into date data types
data_df$date <- as.Date(data_df$date)

Como primer paso quiero visualizar los datos. Para permitir una comparación visual mejorada de las tendencias de crecimiento en los datos, lo normalizo. Para esto defino una función:

normalize_data <- function(x){
  return ((x-min(x))/(max(x)-min(x)))
}

Usando la función, convierto mi data_df en un marco de datos compatible con ggplot que está normalizado. Paso esto a ggplot y visualizo los datos:

dates <- c(data_df$date,
           data_df$date,
           data_df$date,
           data_df$date,
           data_df$date)

values <- c(normalize_data(data_df$gdp),
            normalize_data(data_df$private_investment),
            normalize_data(data_df$life_expectancy),
            normalize_data(data_df$population),
            normalize_data(data_df$rnd))

sources <- c(rep("gdp",times=nrow(data_df)),
             rep("private_investment",times=nrow(data_df)),
             rep("life_expectancy",times=nrow(data_df)),
             rep("population",times=nrow(data_df)),
             rep("rnd",times=nrow(data_df)))

# build ggplot-friendly data frame
ggdata_df <- as.data.frame(matrix(nrow=5*nrow(data_df),ncol=3))
colnames(ggdata_df) <- c("date","value","source")

# populating the ggplot-friendly data frame
ggdata_df$date <- as.Date(dates)
ggdata_df$value <- as.numeric(values)
ggdata_df$source <- sources

# visualize the normalized values (i.e. data set) with ggplot
ggplot(ggdata_df) + 
  geom_point(mapping = aes(x=dates,y=values,color=sources)) +
  ggtitle("Normalized data set values; time-series view") +
  xlab("Time") +
  ylab("Normalized observation values")
resumen de datos para datos de entrenamiento de regresión

Pipeline para regresión en r

Después de haber visualizado el conjunto de datos de interés, pruebo varios modelos de regresión diferentes en los datos.

Para cada enfoque alternativo de regresión en R, el proceso es siempre el mismo:

  1. Divida los datos en conjuntos de entrenamiento y prueba
  2. Entrenar un predictor en el conjunto de entrenamiento
  3. Revisar la bondad de ajuste para el regresor basado en el conjunto de entrenamiento
  4. Probar el regresor en el conjunto de prueba

Regresión lineal simple en R

Comenzaré mi serie de ejemplos sobre regresión en R con un modelo de regresión lineal simple. Trato de predecir el tamaño de la población a partir de los otros datos de observación proporcionados por el conjunto de datos.

A continuación, implemento la tubería mencionada anteriormente (tenga en cuenta que no tengo que normalizar los datos en este caso).

Dado que el tamaño de la población servirá como variable dependiente, ahora tengo cuatro variables independientes para elegir. Para este ejemplo, elijo la esperanza de vida (al nacer) como variable predictiva independiente. Es decir, para mi ejemplo sobre la regresión en RI, en este caso derivará la relación lineal entre el tamaño de la población y la esperanza de vida.

Mi suposición es que una mayor esperanza de vida va acompañada de un mayor tamaño de la población (para el mismo país y región).

# randomly generate split
set.seed(123)
training_split <- sample.split(data_df$date,SplitRatio = 0.8)

# extract training and test sets
training_df <- subset(data_df, training_split)
test_df <- subset(data_df, !training_split)

# train predictor based on training set
predictor <- lm(formula = population ~ life_expectancy, data=training_df)

# print summary of simple linear regression
summary(predictor)
## 
## Call:
## lm(formula = population ~ life_expectancy, data = training_df)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -15279  -6746   1131   4625  17762 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -840396.8    32055.5  -26.22   <2e-16 ***
## life_expectancy   14629.9      427.4   34.23   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8377 on 44 degrees of freedom
## Multiple R-squared:  0.9638, Adjusted R-squared:  0.963 
## F-statistic:  1172 on 1 and 44 DF,  p-value: < 2.2e-16

Visualizo el resultado de esta regresión en R usando ggplot2. En el siguiente código muestro cómo hacerlo:

ggplot(training_df) + 
  geom_point(mapping=aes(x=life_expectancy,
                         y=population)) +
  geom_line(mapping=aes(x=life_expectancy,
                        y=predict(predictor,training_df)),color="red") +
  ggtitle("Population in dependence of life expectancy; linear regression (training)") +
  xlab("US life expectancy since birth [years]") +
  ylab("US population size [-]")
regresión lineal en R del tamaño de la población de EE. UU. resultante de la esperanza de vida

También quiero ver un histograma de los residuos de predicción:

hist(predictor$residuals,
     main ="Histogram of model residuals",
     xlab="Model residuals [-]")

Además, quiero ver los residuos en función de predecir el tamaño de la observación variable:

plot(x=training_df$life_expectancy,
     y=predictor$residuals, 
     main = "Model residuals in dependence of independent variable",
     xlab="US life expectancy [years]",
     ylab="Model residuals [-]")

Esto parece un sesgo algo sistemático. También puedo intentar mostrar esto en un gráfico qqnorm:

qqnorm(predictor$residuals)

Después de haber entrenado el predictor lineal, pruebo el rendimiento de la predicción en el conjunto de prueba:

# predict the test set values
predictions <- predict(predictor,test_df)

# visualize prediction accuracy
ggplot(test_df) + 
  geom_point(mapping=aes(x=life_expectancy,y=population)) +
  geom_line(mapping=aes(x=life_expectancy,predictions),color="red") +
  ggtitle("Population in dependence of life expectancy; linear regression (test)") +
  xlab("US life expectancy since birth [years]") +
  ylab("US population size [-]")
regresión lineal simple en R

Veamos si podemos obtener predicciones más precisas utilizando un enfoque de regresión lineal multivariable. Pero primero, almacenamos las predicciones de regresión lineal simple para todo el conjunto de datos en un vector:

predictions_slr <- predict(predictor,data_df)

Regresión lineal múltiple

Esta vez empiezo considerando todas las variables en el conjunto de datos al predecir el tamaño de la población de EE. UU.:

# split data in training and test set
set.seed(123)
training_split <- sample.split(data_df$date, SplitRatio = 0.8)
training_df <- subset(data_df, training_split)
test_df <- subset(data_df, !training_split)

# train predictor with mutliple linear regression methodology, on training set
predictor <- lm(formula = population ~ gdp + 
                  private_investment +
                  life_expectancy +
                  rnd, training_df)

# summarize regression outcome 
summary(predictor)
## 
## Call:
## lm(formula = population ~ gdp + private_investment + life_expectancy + 
##     rnd, data = training_df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10472.0  -1969.0    188.3   2421.8   7782.1 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        -3.202e+05  4.488e+04  -7.133 1.07e-08 ***
## gdp                 9.891e+00  2.970e+00   3.330  0.00184 ** 
## private_investment -3.272e+00  5.115e+00  -0.640  0.52598    
## life_expectancy     7.283e+03  6.318e+02  11.527 1.93e-14 ***
## rnd                -1.923e+02  8.065e+01  -2.385  0.02181 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3834 on 41 degrees of freedom
## Multiple R-squared:  0.9929, Adjusted R-squared:  0.9922 
## F-statistic:  1440 on 4 and 41 DF,  p-value: < 2.2e-16

El R-cuadrado ajustado ha mejorado, en comparación con la regresión lineal simple anterior. El volumen de inversión privada no parece ser tan significativo para predecir el tamaño de la población de EE. UU. Por lo tanto, aplico la eliminación de parámetros por pasos en la que elimino el volumen de inversión privada como variable independiente al realizar el análisis de regresión:

# regression without private investment volume as independent variable
predictor <- lm(formula=population~gdp + 
                  life_expectancy +
                  rnd, training_df)

# review model performance
summary(predictor)
## 
## Call:
## lm(formula = population ~ gdp + life_expectancy + rnd, data = training_df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10412.6  -2275.0    247.7   2209.6   7794.9 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -3.225e+05  4.442e+04  -7.259 6.20e-09 ***
## gdp              8.598e+00  2.160e+00   3.980 0.000267 ***
## life_expectancy  7.316e+03  6.252e+02  11.702 8.41e-15 ***
## rnd             -1.673e+02  7.006e+01  -2.388 0.021496 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3807 on 42 degrees of freedom
## Multiple R-squared:  0.9929, Adjusted R-squared:  0.9924 
## F-statistic:  1948 on 3 and 42 DF,  p-value: < 2.2e-16

La eliminación del volumen de inversión privada ha mejorado el R-cuadrado ajustado (ligeramente). Sin embargo, todavía estoy interesado en ver la distribución de los residuos.

Primero, un histograma de residuos:

hist(predictor$residuals)

A continuación, un gráfico qqnorm de residuos:

qqnorm(predictor$residuals)

Ahora evalúo el rendimiento del modelo prediciendo los valores de población del conjunto de prueba:

# predict test set population values
predictions <- predict(predictor,test_df)

# visualize model performance vs test set, along timeline
ggplot(test_df) + 
  geom_point(mapping=aes(x=date,y=population)) +
  geom_line(mapping=aes(x=date,y=predictions),color="red") +
  ggtitle("Multiple linear regression prediction of US population") +
  xlab("Date") +
  ylab("US population size [-]")
regresión lineal múltiple en R que nos predice el tamaño de la población

Antes de continuar con la regresión polinomial, almaceno predicciones para todo el conjunto de datos en un vector:

predictions_mlr <- predict(predictor,data_df)

Regresión polinomial en R

Teniendo en cuenta el resultado de la regresión lineal múltiple después de la eliminación hacia atrás, parece que, después de todo, la esperanza de vida es bastante útil para predecir el tamaño de la población.

Por lo tanto, trato de realizar una regresión polinomial en R, prediciendo el tamaño de la población de EE. UU. con el conocimiento de la esperanza de vida de EE. UU. al nacer.

Mi enfoque es agregar términos polinómicos siempre que esto mejore el R-cuadrado ajustado:

# split in training and test set
set.seed(123)
training_split <- sample.split(data_df$date,SplitRatio = 0.8)
training_df <- subset(data_df, training_split)
test_df <- subset(data_df, !training_split)

# train predictor on first order term only, i.e. simple linear regression
predictor <- lm(formula = population ~ life_expectancy,training_df)

# review model performance
summary(predictor)
## 
## Call:
## lm(formula = population ~ life_expectancy, data = training_df)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -15279  -6746   1131   4625  17762 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -840396.8    32055.5  -26.22   <2e-16 ***
## life_expectancy   14629.9      427.4   34.23   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8377 on 44 degrees of freedom
## Multiple R-squared:  0.9638, Adjusted R-squared:  0.963 
## F-statistic:  1172 on 1 and 44 DF,  p-value: < 2.2e-16

Ahora me pregunto: ¿Puedo mejorar el R-cuadrado ajustado agregando un término polinomial de segundo grado? Déjame probar esto:

# polynomial two term regression
training_df$LE2 <- training_df$life_expectancy^2
predictor <- lm(formula = population ~ life_expectancy + LE2, training_df)

# review model performance
summary(predictor)
## 
## Call:
## lm(formula = population ~ life_expectancy + LE2, data = training_df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10785.9  -4081.8   -715.1   4186.9   8913.4 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     4014088.8   604959.5   6.635 4.35e-08 ***
## life_expectancy -116186.0    16295.0  -7.130 8.34e-09 ***
## LE2                 879.9      109.6   8.029 4.31e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 5360 on 43 degrees of freedom
## Multiple R-squared:  0.9855, Adjusted R-squared:  0.9848 
## F-statistic:  1463 on 2 and 43 DF,  p-value: < 2.2e-16

Sí. Se mejoró la R cuadrada ajustada. ¿Qué sucede si agrego un tercer término?

# polynomial three term regression
training_df$LE3 <- training_df$life_expectancy^3
predictor <- lm(formula = population ~ life_expectancy + LE2 + LE3, training_df)

# review model performance
summary(predictor)
## 
## Call:
## lm(formula = population ~ life_expectancy + LE2 + LE3, data = training_df)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -10908  -3739   -862   4079   9125 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)
## (Intercept)      1.196e+07  2.240e+07   0.534    0.596
## life_expectancy -4.367e+05  9.037e+05  -0.483    0.631
## LE2              5.187e+03  1.214e+04   0.427    0.671
## LE3             -1.927e+01  5.432e+01  -0.355    0.725
## 
## Residual standard error: 5415 on 42 degrees of freedom
## Multiple R-squared:  0.9856, Adjusted R-squared:  0.9845 
## F-statistic: 955.5 on 3 and 42 DF,  p-value: < 2.2e-16

Agregar un tercer término no mejoró el R-cuadrado ajustado. ¿Qué sucederá si agrego un cuarto término a la regresión en R?

# polynomial four term regression
training_df$LE4 <- training_df$life_expectancy^4
predictor <- lm(formula = population ~ life_expectancy + LE2 + LE3 + LE4, training_df)

# review model performance
summary(predictor)
## 
## Call:
## lm(formula = population ~ life_expectancy + LE2 + LE3 + LE4, 
##     data = training_df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -11422.7  -2990.7     -4.8   2420.7  11496.0 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -2.345e+09  6.041e+08  -3.882 0.000369 ***
## life_expectancy  1.264e+08  3.250e+07   3.889 0.000362 ***
## LE2             -2.553e+06  6.555e+05  -3.895 0.000355 ***
## LE3              2.290e+04  5.873e+03   3.900 0.000350 ***
## LE4             -7.697e+01  1.972e+01  -3.903 0.000346 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4680 on 41 degrees of freedom
## Multiple R-squared:  0.9895, Adjusted R-squared:  0.9884 
## F-statistic: 963.4 on 4 and 41 DF,  p-value: < 2.2e-16

Agregar un cuarto término mejoró la regresión en R, ya que aumentó R-cuadrado ajustado. Procedo sumando términos, agregando un quinto y un sexto término de regresión.

# polynomial six term regression
training_df$LE5 <- training_df$life_expectancy^5
training_df$LE6 <- training_df$life_expectancy^6
predictor <- lm(formula = population ~ life_expectancy + LE2 + LE3 + 
                  LE4 + LE5 + LE6, training_df)

# review model performance
summary(predictor)
## 
## Call:
## lm(formula = population ~ life_expectancy + LE2 + LE3 + LE4 + 
##     LE5 + LE6, data = training_df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -12579.2  -2669.2   -306.2   1781.2  12719.2 
## 
## Coefficients: (1 not defined because of singularities)
##                   Estimate Std. Error t value Pr(>|t|)  
## (Intercept)     -3.842e+10  1.696e+10  -2.265   0.0290 *
## life_expectancy  2.459e+09  1.097e+09   2.243   0.0305 *
## LE2             -6.145e+07  2.768e+07  -2.220   0.0322 *
## LE3              7.276e+05  3.312e+05   2.197   0.0339 *
## LE4             -3.632e+03  1.671e+03  -2.174   0.0357 *
## LE5                     NA         NA      NA       NA  
## LE6              4.285e-02  2.013e-02   2.128   0.0395 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4491 on 40 degrees of freedom
## Multiple R-squared:  0.9905, Adjusted R-squared:  0.9894 
## F-statistic: 837.9 on 5 and 40 DF,  p-value: < 2.2e-16

Agregar un quinto y sexto término de regresión mejoró el R-cuadrado ajustado. Decido continuar y agregar un séptimo y un octavo término:

# polynomial seven term regression
training_df$LE7 <- training_df$life_expectancy^7
training_df$LE8 <- training_df$life_expectancy^8
predictor <- lm(formula = population ~ life_expectancy + LE2 + LE3 + LE4 + 
                  LE5 + LE6 + LE7 + LE8, training_df)

# review model performance
summary(predictor)
## 
## Call:
## lm(formula = population ~ life_expectancy + LE2 + LE3 + LE4 + 
##     LE5 + LE6 + LE7 + LE8, data = training_df)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -12579.2  -2669.2   -306.2   1781.2  12719.2 
## 
## Coefficients: (3 not defined because of singularities)
##                   Estimate Std. Error t value Pr(>|t|)  
## (Intercept)     -3.842e+10  1.696e+10  -2.265   0.0290 *
## life_expectancy  2.459e+09  1.097e+09   2.243   0.0305 *
## LE2             -6.145e+07  2.768e+07  -2.220   0.0322 *
## LE3              7.276e+05  3.312e+05   2.197   0.0339 *
## LE4             -3.632e+03  1.671e+03  -2.174   0.0357 *
## LE5                     NA         NA      NA       NA  
## LE6              4.285e-02  2.013e-02   2.128   0.0395 *
## LE7                     NA         NA      NA       NA  
## LE8                     NA         NA      NA       NA  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4491 on 40 degrees of freedom
## Multiple R-squared:  0.9905, Adjusted R-squared:  0.9894 
## F-statistic: 837.9 on 5 and 40 DF,  p-value: < 2.2e-16

El R-cuadrado ajustado no mejoró. Decido detener la regresión polinomial aquí y visualizar el resultado de la regresión:

ggplot(training_df) + 
  geom_point(mapping=aes(x=life_expectancy,y=population)) +
  geom_line(mapping=aes(x=life_expectancy,y=predict(predictor,training_df)),
            color="red") +
  ggtitle("Trained polynomial regression model; degree 6") + 
  xlab("US life expectancy at birth [years]") + 
  ylab("US population size [-]")
Regresión polinomial de sexto grado en R

Ahora pruebo el modelo polinomial contra el conjunto de prueba:

# add terms to test data frame
test_df$LE2 <- test_df$life_expectancy^2
test_df$LE3 <- test_df$life_expectancy^3
test_df$LE4 <- test_df$life_expectancy^4
test_df$LE5 <- test_df$life_expectancy^5
test_df$LE6 <- test_df$life_expectancy^6
test_df$LE7 <- test_df$life_expectancy^7
test_df$LE8 <- test_df$life_expectancy^8

# visualize predictor performance on test set
predictions <- predict(predictor,test_df)
ggplot(test_df) + 
  geom_point(mapping=aes(x=life_expectancy,y=population)) +
  geom_line(mapping=aes(x=life_expectancy,y=predictions), color = "red") +
  ggtitle("Polynomial regression on test data; degree 6") + 
  xlab("US life expectancy at birth [years]") + 
  ylab("US population size [-]")
Regresión polinomial en R

Antes de pasar a la regresión del árbol de decisiones, almaceno predicciones para todo el conjunto de datos en un vector:

data_df$LE2 <- data_df$life_expectancy^2
data_df$LE3 <- data_df$life_expectancy^3
data_df$LE4 <- data_df$life_expectancy^4
data_df$LE5 <- data_df$life_expectancy^5
data_df$LE6 <- data_df$life_expectancy^6
data_df$LE7 <- data_df$life_expectancy^7
data_df$LE8 <- data_df$life_expectancy^8
predictions_poly <- predict(predictor,data_df)

Regresión del árbol de decisión en R

Otra metodología de regresión es la regresión del árbol de decisión. En este tipo de regresión, el espacio de soluciones se divide en subconjuntos, para los cuales el valor de observación promedio medido en el conjunto de entrenamiento constituye una predicción.

En la siguiente sección, implemento la regresión del árbol de decisiones en R:

# split in training and test set
set.seed(123)
training_split <- sample.split(data_df$date,SplitRatio = 0.8)
training_df <- subset(data_df, training_split)
test_df <- subset(data_df, !training_split)

# conduct decision tree regression
predictor <- rpart(formula = population ~ gdp +
                     private_investment +
                     life_expectancy + 
                     population + 
                     rnd,
                   data = training_df)

# generate predictions from test set
predictions <- predict(predictor,test_df)
residuals <- test_df$population - predictions

#  visualize prediction accuracy along the timeline
ggplot(test_df) + 
  geom_point(mapping = aes(x=date,y=population)) + 
  geom_point(mapping = aes(x=date,y=predictions), color = "red") + 
  ggtitle("Decision tree regression performance on test set (red: predictions)") + 
  xlab("Date") +
  ylab("US population size [-]")
Regresión del árbol de decisión en R

Antes de pasar a admitir la regresión vectorial, almaceno mis predicciones para todo el conjunto de datos en otro vector:

predictions_dtr <- predict(predictor,data_df)

Regresión de vectores de soporte

Eventualmente, realizo una regresión de vectores de soporte. Solo introduzco las variables que la regresión lineal múltiple encontró que son estadísticamente significativas en la predicción del tamaño de la población de EE. UU.:

# split in training and test set
set.seed(123)
training_split <- sample.split(data_df$date,SplitRatio = 0.8)
training_df <- subset(data_df, training_split)
test_df <- subset(data_df, !training_split)

# construct predictor
predictor <- svm(formula = population ~ gdp + 
                   life_expectancy + 
                   rnd,
                 data = training_df,
                 type = "eps-regression")

# visualize prediction accuracy on test set along time line
ggplot(test_df) + 
  geom_point(mapping = aes(x=date, y = population)) + 
  geom_line(mapping = aes(x=date, y = predict(predictor,test_df)), color = "red") + 
  ggtitle("Support vector regression (eps methode) on test set") + 
  xlab("Date") + 
  ylab("US population size [-]")
Soporte de regresión vectorial en R

Al igual que en las otras partes de la regresión, almaceno la predicción sobre todo el conjunto de datos como un espacio de solución en un vector separado:

predictions_svr <- predict(predictor,data_df)

Comparación de las predicciones del modelo

Usando los vectores con valores de predicción almacenados, podemos comparar los diferentes modelos de predicción a lo largo de la línea de tiempo:

# construct ggplot friendly data frame
ggdata_df <- as.data.frame(matrix(nrow=nrow(data_df)*6,ncol=3))
colnames(ggdata_df) <- c("date","value","source")
ggdata_df$date <- rep(data_df$date,times=6)
ggdata_df$value <- c(predictions_slr,
                     predictions_mlr,
                     predictions_poly,
                     predictions_dtr,
                     predictions_svr,
                     data_df$population)
ggdata_df$source <- c(rep("simple linear", times=nrow(data_df)),
                      rep("multiple linear", times=nrow(data_df)),
                      rep("polynomial", times=nrow(data_df)),
                      rep("decision tree", times=nrow(data_df)),
                      rep("support vector", times=nrow(data_df)),
                      rep("original observations", times=nrow(data_df)))

# visualize
ggplot(ggdata_df) + 
  geom_point(mapping = aes(x=date, y=value, color=source)) +
  ggtitle("Comparison of model predictions against real observations values") + 
  xlab("Date") +
  ylab("US population size [-]")
Comparando varios métodos de regresión en R

Leave a Reply

Leave a Reply

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *

Este sitio usa Akismet para reducir el spam. Aprende cómo se procesan los datos de tus comentarios.