Profesor: Felipe González

Elaborado por Daniela Pinto Veizaga

library(tidyverse)
library(tidymodels)
library(cowplot)

1 . Definición de la función f_real:

La función real a estimar es una función definida por partes: se comporta como una raíz para \(x<10\) y como una constante para \(x\geq 10\). Misma que para describir en su totalidad necesita solo 2 puntos.

f_real <- function(x){
  ifelse(x < 10, 1000*sqrt(x), 1000*sqrt(10))
}

genera_datos_sin_perturbacion<- function(n = 100){
  x <- runif(n, 0, 25)
  y <- f_real(x)
  tibble(x = x, y = y)
}

Tibble de puntos que se generan aleatoriamente con la función f_real.

Graficando la forma de la función f_real:

  1. Generamos datos perturbados con la función f_real, agregando rnorm con media cero y desviación estándar 500. El resultado: un normal con media f_real(x) y desviación estándar 500.
genera_datos <- function(n = 100){
  x <- runif(n, 0, 25)
  y <- f_real(x) + rnorm(n, 0, 500)
  tibble(x = x, y = y)
}
calcular_grafica <- function(mod, nombre = ""){
  datos_g <- tibble(x = seq(0, 25, 0.01))
  datos_g <- predict(mod, datos_g) %>% 
    bind_cols(datos_g)
  datos_g %>% mutate(nombre = nombre)
}

Graficando la forma de la función f_real con los datos perturbados generados.

  1. Generación de métodos que emplearemos para el ajuste de los datos: regresión lineal, polinomio de grado tres y polinomio de grado ocho.
modelo_lineal <- linear_reg() %>% 
  set_engine("lm")
modelo_svm <- svm_poly() %>% 
  set_engine("kernlab") %>% 
  set_mode("regression")

Simulamos datos de entrenamiento:

set.seed(8181)
datos <- genera_datos(30)
  1. Generamos los objetos de los tres métodos antes mencionados.
# Ajuste
mod_1 <- modelo_lineal %>% fit(y ~ x, datos)
mod_2 <- modelo_svm %>% set_args(cost = 0.1, degree = 3) %>% 
  fit(y ~ x, datos)
mod_3 <- modelo_svm %>% set_args(cost = 100, degree = 8) %>% 
  fit(y ~ x, datos)

Ejemplo de la información guardada en cada uno de los objetos generados anteriormente.

mod_1
parsnip model object

Fit time:  4ms 

Call:
stats::lm(formula = y ~ x, data = data)

Coefficients:
(Intercept)            x  
    2247.55        37.85  
  1. Generamos datos con cada uno de los métodos antes guardados en objetos y graficamos estos datos en una sola gráfica.
datos_1 <- calcular_grafica(mod_1, "regresión lineal")
datos_2 <- calcular_grafica(mod_2, "polinomio de grado 3")
datos_3 <- calcular_grafica(mod_3, "polinomio de grado 8")

# Unimos las observaciones de los tres modelos en un solo dataframe 
datos_g <- bind_rows(datos_1, datos_2, datos_3)

ggplot(datos, aes(x = x)) +
  geom_line(data = datos_g, aes(y = .pred, colour = nombre, group = nombre), size = 1.5) +
    geom_point(aes(y = y)) 

  1. Repetimos con distintas muestras de entrenamiento. Describe cómo se comporta cada ajuste: sesgo y varianza.
datos_f <- tibble(x = seq(0, 25, 0.01)) %>% 
  mutate(.pred = f_real(x)) %>% 
  mutate(nombre = "verdadera f")
for (i in 0:6){
  datos_entrena <- genera_datos(100)
  mod_1 <- modelo_lineal %>% fit(y ~ x, datos_entrena)
  mod_2 <- modelo_svm %>% set_args(cost = 0.1, degree = 3) %>% 
    fit(y ~ x, datos_entrena)
  mod_3 <- modelo_svm %>% set_args(cost = 100, degree = 8) %>% 
    fit(y ~ x, datos_entrena)
  datos_1 <- calcular_grafica(mod_1, "regresión lineal")
  datos_2 <- calcular_grafica(mod_2, "polinomio de grado 3")
  datos_3 <- calcular_grafica(mod_3, "polinomio de grado ocho")
  datos_g <- bind_rows(datos_1, datos_2, datos_3, datos_f)
  Object <- ggplot(datos_entrena, aes(x = x)) +
    geom_line(data = datos_g, aes(y = .pred, colour = nombre, group = nombre), size = 1.5)  +
      geom_point(aes(y = y), size = 1.5) +
    ylim(c(-500, 4500))
  assign(paste0("plot", i), Object)
}
plot_grid(plot1,
          plot2,
          plot3,
          plot4,
          plot5,
          plot6,
          labels = 'AUTO',
          label_fontfamily = "serif",
          label_fontface = "plain",
          label_colour = "blue",
          ncol = 2)

Anotaciones: Observaciones en relación al sesgo y la varianza

Intuitivamente, el modelo polinomial de grado 3 es menos variante ante distintas muestras. Esto se debe a la forma de la \(f\) verdadera (como se menciona en el punto 1). Para describir a \(f\) se necesitan 3 puntos para conocerla totalmente. Y el polinomio de grado 3 necesita 4 puntos y el polinomio de grado 8 necesita 9 puntos. Por ello, ante cambios grandes en la muestra el polinomio de grado 8 cambia drásticamente. El polinomio de grado 3 mantiene menos variación. El modelo lineal es sensible a outliers por lo que tampoco es adecuado.

Por los argumentos anteriores, dentro de los modelos presentados, el que puede presentar el mejor redimiento es en polinimio cúbico porque los otros dos modelos se sobreajustan (overfitting) o subajustan (underfitting) a los datos generados.

LS0tCnRpdGxlOiAiQXByZW5kaXphamUgZGUgTcOhcXVpbmEuIFRhcmVhIDE6IFNlc2dvIHkgdmFyaWFuemEiCm91dHB1dDogaHRtbF9ub3RlYm9vawotLS0KKlByb2Zlc29yOiBGZWxpcGUgR29uesOhbGV6KgoKKkVsYWJvcmFkbyBwb3IgRGFuaWVsYSBQaW50byBWZWl6YWdhKgpgYGB7ciBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQpsaWJyYXJ5KHRpZHl2ZXJzZSkKbGlicmFyeSh0aWR5bW9kZWxzKQpsaWJyYXJ5KGNvd3Bsb3QpCmBgYAoKPiAxIC4gRGVmaW5pY2nDs24gZGUgbGEgZnVuY2nDs24gYGZfcmVhbGA6CgoKTGEgZnVuY2nDs24gcmVhbCBhIGVzdGltYXIgZXMgdW5hIGZ1bmNpw7NuIGRlZmluaWRhIHBvciBwYXJ0ZXM6IHNlIGNvbXBvcnRhIGNvbW8gdW5hIHJhw616IHBhcmEgJHg8MTAkIHkgY29tbyB1bmEgY29uc3RhbnRlIHBhcmEgJHhcZ2VxIDEwJC4gTWlzbWEgcXVlIHBhcmEgZGVzY3JpYmlyIGVuIHN1IHRvdGFsaWRhZCBuZWNlc2l0YSBzb2xvIDIgcHVudG9zLgoKYGBge3J9CmZfcmVhbCA8LSBmdW5jdGlvbih4KXsKICBpZmVsc2UoeCA8IDEwLCAxMDAwKnNxcnQoeCksIDEwMDAqc3FydCgxMCkpCn0KCmBgYAoKCmBgYHtyfQoKZ2VuZXJhX2RhdG9zX3Npbl9wZXJ0dXJiYWNpb248LSBmdW5jdGlvbihuID0gMTAwKXsKICB4IDwtIHJ1bmlmKG4sIDAsIDI1KQogIHkgPC0gZl9yZWFsKHgpCiAgdGliYmxlKHggPSB4LCB5ID0geSkKfQoKCmBgYAoKClRpYmJsZSBkZSBwdW50b3MgcXVlIHNlIGdlbmVyYW4gYWxlYXRvcmlhbWVudGUgY29uIGxhIGZ1bmNpw7NuIGBmX3JlYWxgLgpgYGB7ciBlY2hvPUZBTFNFfQpnZW5lcmFfZGF0b3Nfc2luX3BlcnR1cmJhY2lvbigKICAKKQpgYGAKCkdyYWZpY2FuZG8gbGEgZm9ybWEgZGUgbGEgZnVuY2nDs24gKipmX3JlYWwqKjoKCmBgYHtyIGVjaG89RkFMU0UsIHdhcm5pbmc9RkFMU0UsIGZpZy5hbGlnbiA9ICdjZW50ZXInfQpnZ3Bsb3QoZ2VuZXJhX2RhdG9zX3Npbl9wZXJ0dXJiYWNpb24oKSwgYWVzKHggPSB4LCB5ID0geSkpICsgZ2VvbV9wb2ludCh3aWR0aCA9IDAuMywgY29sb3VyID0gIiM3Q0FFMDAiKSAgKyBnZ3RpdGxlKCJGIHJlYWwiKQpgYGAKCgoKPiAyLiBHZW5lcmFtb3MgZGF0b3MgcGVydHVyYmFkb3MgY29uIGxhIGZ1bmNpw7NuICoqZl9yZWFsKiosIGFncmVnYW5kbyBybm9ybSBjb24gbWVkaWEgY2VybyB5IGRlc3ZpYWNpw7NuIGVzdMOhbmRhciA1MDAuIEVsIHJlc3VsdGFkbzogdW4gbm9ybWFsIGNvbiBtZWRpYSBmX3JlYWwoeCkgeSBkZXN2aWFjacOzbiBlc3TDoW5kYXIgNTAwLgoKCmBgYHtyfQpnZW5lcmFfZGF0b3MgPC0gZnVuY3Rpb24obiA9IDEwMCl7CiAgeCA8LSBydW5pZihuLCAwLCAyNSkKICB5IDwtIGZfcmVhbCh4KSArIHJub3JtKG4sIDAsIDUwMCkKICB0aWJibGUoeCA9IHgsIHkgPSB5KQp9CmNhbGN1bGFyX2dyYWZpY2EgPC0gZnVuY3Rpb24obW9kLCBub21icmUgPSAiIil7CiAgZGF0b3NfZyA8LSB0aWJibGUoeCA9IHNlcSgwLCAyNSwgMC4wMSkpCiAgZGF0b3NfZyA8LSBwcmVkaWN0KG1vZCwgZGF0b3NfZykgJT4lIAogICAgYmluZF9jb2xzKGRhdG9zX2cpCiAgZGF0b3NfZyAlPiUgbXV0YXRlKG5vbWJyZSA9IG5vbWJyZSkKfQpgYGAKCgpHcmFmaWNhbmRvIGxhIGZvcm1hIGRlIGxhIGZ1bmNpw7NuICoqZl9yZWFsKiogY29uIGxvcyBkYXRvcyBwZXJ0dXJiYWRvcyBnZW5lcmFkb3MuCgpgYGB7ciBlY2hvPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQpnZ3Bsb3QoZ2VuZXJhX2RhdG9zKCksIGFlcyh4ID0geCwgeSA9IHkpKSArIGdlb21fcG9pbnQod2lkdGggPSAwLjMsIGNvbG91ciA9ICIjN0NBRTAwIikgICsgZ2d0aXRsZSgiRiByZWFsIGNvbiBkYXRvcyBwZXJ0dXJiYWRvcyIpCmBgYAoKCj4gMy4gR2VuZXJhY2nDs24gZGUgbcOpdG9kb3MgcXVlIGVtcGxlYXJlbW9zIHBhcmEgZWwgYWp1c3RlIGRlIGxvcyBkYXRvczogcmVncmVzacOzbiBsaW5lYWwsIHBvbGlub21pbyBkZSBncmFkbyB0cmVzIHkgcG9saW5vbWlvIGRlIGdyYWRvIG9jaG8uCgpgYGB7ciwgbWVzc2FnZSA9IEZBTFNFLCB3YXJuaW5nPUZBTFNFfQptb2RlbG9fbGluZWFsIDwtIGxpbmVhcl9yZWcoKSAlPiUgCiAgc2V0X2VuZ2luZSgibG0iKQptb2RlbG9fc3ZtIDwtIHN2bV9wb2x5KCkgJT4lIAogIHNldF9lbmdpbmUoImtlcm5sYWIiKSAlPiUgCiAgc2V0X21vZGUoInJlZ3Jlc3Npb24iKQpgYGAKClNpbXVsYW1vcyBkYXRvcyBkZSBlbnRyZW5hbWllbnRvOgoKYGBge3J9CnNldC5zZWVkKDgxODEpCmRhdG9zIDwtIGdlbmVyYV9kYXRvcygzMCkKYGBgCgoKYGBge3IgZWNobz1GQUxTRX0KZGF0b3MKYGBgCgo+IDQuIEdlbmVyYW1vcyBsb3Mgb2JqZXRvcyBkZSBsb3MgdHJlcyBtw6l0b2RvcyBhbnRlcyBtZW5jaW9uYWRvcy4KCmBgYHtyLCBtZXNzYWdlID0gRkFMU0UsIHdhcm5pbmc9RkFMU0V9CiMgQWp1c3RlCm1vZF8xIDwtIG1vZGVsb19saW5lYWwgJT4lIGZpdCh5IH4geCwgZGF0b3MpCm1vZF8yIDwtIG1vZGVsb19zdm0gJT4lIHNldF9hcmdzKGNvc3QgPSAwLjEsIGRlZ3JlZSA9IDMpICU+JSAKICBmaXQoeSB+IHgsIGRhdG9zKQptb2RfMyA8LSBtb2RlbG9fc3ZtICU+JSBzZXRfYXJncyhjb3N0ID0gMTAwLCBkZWdyZWUgPSA4KSAlPiUgCiAgZml0KHkgfiB4LCBkYXRvcykKYGBgCgpFamVtcGxvIGRlIGxhIGluZm9ybWFjacOzbiBndWFyZGFkYSBlbiBjYWRhIHVubyBkZSBsb3Mgb2JqZXRvcyBnZW5lcmFkb3MgYW50ZXJpb3JtZW50ZS4KCmBgYHtyfQptb2RfMQpgYGAKPiA1LiBHZW5lcmFtb3MgZGF0b3MgY29uIGNhZGEgdW5vIGRlIGxvcyBtw6l0b2RvcyBhbnRlcyBndWFyZGFkb3MgZW4gb2JqZXRvcyB5IGdyYWZpY2Ftb3MgZXN0b3MgZGF0b3MgZW4gdW5hIHNvbGEgZ3LDoWZpY2EuCgpgYGB7ciwgZmlnLndpZHRoPTcsIGZpZy5hc3A9MC43fQpkYXRvc18xIDwtIGNhbGN1bGFyX2dyYWZpY2EobW9kXzEsICJyZWdyZXNpw7NuIGxpbmVhbCIpCmRhdG9zXzIgPC0gY2FsY3VsYXJfZ3JhZmljYShtb2RfMiwgInBvbGlub21pbyBkZSBncmFkbyAzIikKZGF0b3NfMyA8LSBjYWxjdWxhcl9ncmFmaWNhKG1vZF8zLCAicG9saW5vbWlvIGRlIGdyYWRvIDgiKQoKIyBVbmltb3MgbGFzIG9ic2VydmFjaW9uZXMgZGUgbG9zIHRyZXMgbW9kZWxvcyBlbiB1biBzb2xvIGRhdGFmcmFtZSAKZGF0b3NfZyA8LSBiaW5kX3Jvd3MoZGF0b3NfMSwgZGF0b3NfMiwgZGF0b3NfMykKCmdncGxvdChkYXRvcywgYWVzKHggPSB4KSkgKwogIGdlb21fbGluZShkYXRhID0gZGF0b3NfZywgYWVzKHkgPSAucHJlZCwgY29sb3VyID0gbm9tYnJlLCBncm91cCA9IG5vbWJyZSksIHNpemUgPSAxLjUpICsKICAgIGdlb21fcG9pbnQoYWVzKHkgPSB5KSkgCmBgYAoKPiA2LiBSZXBldGltb3MgY29uIGRpc3RpbnRhcyBtdWVzdHJhcyBkZSBlbnRyZW5hbWllbnRvLiBEZXNjcmliZSBjw7NtbyBzZSBjb21wb3J0YSBjYWRhIGFqdXN0ZTogc2VzZ28geSB2YXJpYW56YS4KCmBgYHtyIGVjaG89VFJVRSwgZmlnLmFsaWduPSdjZW50ZXInLCBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQpkYXRvc19mIDwtIHRpYmJsZSh4ID0gc2VxKDAsIDI1LCAwLjAxKSkgJT4lIAogIG11dGF0ZSgucHJlZCA9IGZfcmVhbCh4KSkgJT4lIAogIG11dGF0ZShub21icmUgPSAidmVyZGFkZXJhIGYiKQpmb3IgKGkgaW4gMDo2KXsKICBkYXRvc19lbnRyZW5hIDwtIGdlbmVyYV9kYXRvcygxMDApCiAgbW9kXzEgPC0gbW9kZWxvX2xpbmVhbCAlPiUgZml0KHkgfiB4LCBkYXRvc19lbnRyZW5hKQogIG1vZF8yIDwtIG1vZGVsb19zdm0gJT4lIHNldF9hcmdzKGNvc3QgPSAwLjEsIGRlZ3JlZSA9IDMpICU+JSAKICAgIGZpdCh5IH4geCwgZGF0b3NfZW50cmVuYSkKICBtb2RfMyA8LSBtb2RlbG9fc3ZtICU+JSBzZXRfYXJncyhjb3N0ID0gMTAwLCBkZWdyZWUgPSA4KSAlPiUgCiAgICBmaXQoeSB+IHgsIGRhdG9zX2VudHJlbmEpCiAgZGF0b3NfMSA8LSBjYWxjdWxhcl9ncmFmaWNhKG1vZF8xLCAicmVncmVzacOzbiBsaW5lYWwiKQogIGRhdG9zXzIgPC0gY2FsY3VsYXJfZ3JhZmljYShtb2RfMiwgInBvbGlub21pbyBkZSBncmFkbyAzIikKICBkYXRvc18zIDwtIGNhbGN1bGFyX2dyYWZpY2EobW9kXzMsICJwb2xpbm9taW8gZGUgZ3JhZG8gb2NobyIpCiAgZGF0b3NfZyA8LSBiaW5kX3Jvd3MoZGF0b3NfMSwgZGF0b3NfMiwgZGF0b3NfMywgZGF0b3NfZikKICBPYmplY3QgPC0gZ2dwbG90KGRhdG9zX2VudHJlbmEsIGFlcyh4ID0geCkpICsKICAgIGdlb21fbGluZShkYXRhID0gZGF0b3NfZywgYWVzKHkgPSAucHJlZCwgY29sb3VyID0gbm9tYnJlLCBncm91cCA9IG5vbWJyZSksIHNpemUgPSAxLjUpICArCiAgICAgIGdlb21fcG9pbnQoYWVzKHkgPSB5KSwgc2l6ZSA9IDEuNSkgKwogICAgeWxpbShjKC01MDAsIDQ1MDApKQogIGFzc2lnbihwYXN0ZTAoInBsb3QiLCBpKSwgT2JqZWN0KQp9CgpgYGAKCgpgYGB7ciBlY2hvPVRSVUUsIGZpZy5oZWlnaHQ9MTAsIGZpZy53aWR0aD0xNSwgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KcGxvdF9ncmlkKHBsb3QxLAogICAgICAgICAgcGxvdDIsCiAgICAgICAgICBwbG90MywKICAgICAgICAgIHBsb3Q0LAogICAgICAgICAgcGxvdDUsCiAgICAgICAgICBwbG90NiwKICAgICAgICAgIGxhYmVscyA9ICdBVVRPJywKICAgICAgICAgIGxhYmVsX2ZvbnRmYW1pbHkgPSAic2VyaWYiLAogICAgICAgICAgbGFiZWxfZm9udGZhY2UgPSAicGxhaW4iLAogICAgICAgICAgbGFiZWxfY29sb3VyID0gImJsdWUiLAogICAgICAgICAgbmNvbCA9IDIpCmBgYAoKKipBbm90YWNpb25lczogT2JzZXJ2YWNpb25lcyBlbiByZWxhY2nDs24gYWwgc2VzZ28geSBsYSB2YXJpYW56YSoqCgpJbnR1aXRpdmFtZW50ZSwgZWwgbW9kZWxvIHBvbGlub21pYWwgZGUgZ3JhZG8gMyBlcyBtZW5vcyB2YXJpYW50ZSBhbnRlIGRpc3RpbnRhcyBtdWVzdHJhcy4gRXN0byBzZSBkZWJlIGEgbGEgZm9ybWEgZGUgbGEgJGYkIHZlcmRhZGVyYSAoY29tbyBzZSBtZW5jaW9uYSBlbiBlbCBwdW50byAxKS4gUGFyYSBkZXNjcmliaXIgYSAkZiQgc2UgbmVjZXNpdGFuIDMgcHVudG9zIHBhcmEgY29ub2NlcmxhIHRvdGFsbWVudGUuIFkgZWwgcG9saW5vbWlvIGRlIGdyYWRvIDMgbmVjZXNpdGEgNCBwdW50b3MgeSBlbCBwb2xpbm9taW8gZGUgZ3JhZG8gOCBuZWNlc2l0YSA5IHB1bnRvcy4gUG9yIGVsbG8sIGFudGUgY2FtYmlvcyBncmFuZGVzIGVuIGxhIG11ZXN0cmEgZWwgcG9saW5vbWlvIGRlIGdyYWRvIDggY2FtYmlhIGRyw6FzdGljYW1lbnRlLiBFbCBwb2xpbm9taW8gZGUgZ3JhZG8gMyBtYW50aWVuZSBtZW5vcyB2YXJpYWNpw7NuLiBFbCBtb2RlbG8gbGluZWFsIGVzIHNlbnNpYmxlIGEgX291dGxpZXJzXyBwb3IgbG8gcXVlIHRhbXBvY28gZXMgYWRlY3VhZG8uCgpQb3IgbG9zIGFyZ3VtZW50b3MgYW50ZXJpb3JlcywgZGVudHJvIGRlIGxvcyBtb2RlbG9zIHByZXNlbnRhZG9zLCBlbCBxdWUgcHVlZGUgcHJlc2VudGFyIGVsIG1lam9yIHJlZGltaWVudG8gZXMgZW4gcG9saW5pbWlvIGPDumJpY28gcG9ycXVlIGxvcyBvdHJvcyBkb3MgbW9kZWxvcyBzZSBzb2JyZWFqdXN0YW4gKGBvdmVyZml0dGluZ2ApIG8gIHN1YmFqdXN0YW4gKGB1bmRlcmZpdHRpbmdgKSAgYSBsb3MgZGF0b3MgZ2VuZXJhZG9zLgoK