Objetos de datos

En R, toda la información se almacena en objetos. Dependiendo de cómo estén organizados los datos, estos objetos pueden adoptar distintas estructuras. Existen cuatro estructuras principales en R:

I. Vectores

II. Matrices

III. Data frames

IV. Listas

Vamos ahora con la última de las cuatro: las listas.

Listas

Una lista en R es una estructura de datos flexible que permite almacenar objetos de distinto tipo y tamaño dentro de un mismo contenedor. A diferencia de vectores o matrices, los elementos de una lista pueden ser números, texto, vectores, data frames e incluso otras listas.

Para crear una lista se utiliza la función list() que recibe una secuencia de elementos (…) que pueden ser nombrados opcionalmente usando nombre = valor.

# Creación de una lista
mi_lista <- list(
  nombre = "María",
  edad = 25,
  notas = c(8, 9, 10)
)

mi_lista

## $nombre
## [1] "María"
## 
## $edad
## [1] 25
## 
## $notas
## [1]  8  9 10

Usamos corchetes simples [] para acceder a sublistas.

mi_lista[1]

## $nombre
## [1] "María"

class(mi_lista[1]) # sigue siendo una lista, de un solo elemento

## [1] "list"

mi_lista[c(2,3)]

## $edad
## [1] 25
## 
## $notas
## [1]  8  9 10

También utilizando el operador $ cuando las sublistas tienen nombre. Esto es preferido cuando se trabaja con listas nombradas.

mi_lista$nombre      # "María"

## [1] "María"

mi_lista$notas

## [1]  8  9 10

Para acceder al contenido de las sublistas (no a la sublista en sí), usamos corchetes dobles [[ ]].

mi_lista["notas"]     # devuelve una lista con un elemento llamado "notas"

## $notas
## [1]  8  9 10

mi_lista[["notas"]]   # devuelve directamente el vector c(8,9,10)

## [1]  8  9 10

mi_lista[["notas"]][2] # devuelve 9

## [1] 9

Entonces, ¿$ y [[ ]] sirven para lo mismo? En listas nombradas, sí:

mi_lista$nombre

## [1] "María"

mi_lista[["nombre"]]

## [1] "María"

La diferencia aparece cuando el nombre del elemento está guardado en una variable:

x <- "nombre"

mi_lista[[x]]   # funciona: [[ ]] evalúa el contenido de x, es decir busca "nombre"

## [1] "María"

mi_lista$x      # devuelve NULL: $ no evalúa x, busca literalmente un elemento llamado "x"

## NULL

Este último caso no da error, solo devuelve NULL en silencio, lo que puede generar bugs difíciles de detectar.

Para escribir código serio o replicable: [[ ]].
Para exploración rápida con nombres fijos: $ .

Las listas son estructuras dinámicas, por lo que pueden modificarse después de su creación.

Es posible cambiar valores existentes, reemplazar elementos completos o incluso agregar nuevos componentes a la lista sin necesidad de redefinirla por completo.

mi_lista

## $nombre
## [1] "María"
## 
## $edad
## [1] 25
## 
## $notas
## [1]  8  9 10

print("---")

## [1] "---"

# Modificación de elementos
mi_lista$edad <- 26              # cambiar valor
mi_lista$notas <- c(9, 10, 10)   # modificar vector
mi_lista$ciudad <- "CDMX"        # agregar nuevo elemento

mi_lista

## $nombre
## [1] "María"
## 
## $edad
## [1] 26
## 
## $notas
## [1]  9 10 10
## 
## $ciudad
## [1] "CDMX"

Las listas también pueden contener otras listas en su interior, lo que permite construir estructuras jerárquicas más complejas. Este tipo de organización es útil cuando se necesita representar información agrupada, como conjuntos de estudiantes o registros anidados.

clase <- list(e1 = mi_lista,
              e2 = list(nombre = "Luis",
                        edad = 20,
                        notas = c(8, 9, 10))
              )

# Lista clase de dos estudiantes
clase

## $e1
## $e1$nombre
## [1] "María"
## 
## $e1$edad
## [1] 26
## 
## $e1$notas
## [1]  9 10 10
## 
## $e1$ciudad
## [1] "CDMX"
## 
## 
## $e2
## $e2$nombre
## [1] "Luis"
## 
## $e2$edad
## [1] 20
## 
## $e2$notas
## [1]  8  9 10

En estructuras anidadas, el acceso a los elementos requiere encadenar el operador $ tantas veces como niveles tenga la lista. Esto permite llegar a elementos específicos dentro de sublistas.

# Acceso a elementos
clase$e2$nombre             # "Luis"

## [1] "Luis"

clase$e2$notas              # notas de Luis

## [1]  8  9 10

De la misma forma, también es posible modificar elementos dentro de listas anidadas accediendo directamente a su nivel correspondiente, sin necesidad de reconstruir toda la estructura.

# Añadimos ciudad para estudiante 2
clase$e2$ciudad <- "Xalapa"
clase

## $e1
## $e1$nombre
## [1] "María"
## 
## $e1$edad
## [1] 26
## 
## $e1$notas
## [1]  9 10 10
## 
## $e1$ciudad
## [1] "CDMX"
## 
## 
## $e2
## $e2$nombre
## [1] "Luis"
## 
## $e2$edad
## [1] 20
## 
## $e2$notas
## [1]  8  9 10
## 
## $e2$ciudad
## [1] "Xalapa"

Cuando todos los elementos de una lista requieren la misma operación, la función lapply() permite realizarla de forma automática.

En lugar de modificar cada elemento por separado, lapply() aplica una misma función a todos los elementos de la lista y devuelve una nueva lista con los resultados.

# Obtener la mínima nota de cada estudiante
lapply(clase, function(x) min(x$notas))

## $e1
## [1] 9
## 
## $e2
## [1] 8

# Obtener la edad de cada estudiante
lapply(clase, function(x) x$edad)

## $e1
## [1] 26
## 
## $e2
## [1] 20

Para cada estudiante (x) dentro de clase:

toma sus notas: x$notas
encuentra el mínimo: min(...)
devuelve ese valor

Nota: Existe también la función mapply(), que permite aplicar una función utilizando varios objetos al mismo tiempo, por ejemplo, dos o más listas o vectores.

rep(x = "a", times = 3)

## [1] "a" "a" "a"

mapply(rep, x = 4:1, times = 1:4)

## [[1]]
## [1] 4
## 
## [[2]]
## [1] 3 3
## 
## [[3]]
## [1] 2 2 2
## 
## [[4]]
## [1] 1 1 1 1

Las listas se presentan al final porque muchas de las funciones utilizadas en análisis de datos y aprendizaje automático devuelven sus resultados en este formato.

Por ejemplo, al entrenar un modelo estadístico, el objeto resultante suele ser una lista que contiene distintos componentes, como coeficientes, residuos, valores ajustados y métricas de desempeño. Comprender cómo acceder y manipular estos elementos facilitará la interpretación de los modelos que se estudiarán en las siguientes sesiones.

Introducción a Machine Learning

Hasta este momento del curso hemos utilizado R para importar, explorar y describir datos, además de familiarizarnos con sus principales estructuras, como vectores, matrices, data frames, factores y listas.

A partir de esta unidad cambiaremos el enfoque. Nuestro objetivo ya no será únicamente analizar los datos que ya tenemos, sino construir modelos capaces de aprender relaciones presentes en los datos para realizar predicciones sobre nuevas observaciones.

Machine Learning es una rama de la Inteligencia Artificial dedicada al desarrollo de algoritmos capaces de aprender a partir de datos.

En este curso seguiremos el enfoque del aprendizaje estadístico propuesto por el libro Introduction to Statistical Learning with R (ISLR), donde el objetivo es construir modelos que aproximen las relaciones presentes en los datos para realizar predicciones y, cuando sea posible, interpretar esas relaciones.

La idea central puede resumirse mediante la expresión:

\[ Y = f(X) + \varepsilon \] donde:

Y representa la variable que queremos predecir.
X representa la información disponible para hacerlo.
f representa la relación desconocida entre ambas.
ε representa variabilidad no explicada (ruido).

Este último término nos dice que, incluso con la misma X, observaremos distintos valores de Y.

Supongamos que disponemos de información sobre viviendas: tamaño, habitaciones, antigüedad, ubicación y precio. Observando suficientes casos, notaríamos patrones en los datos. La pregunta que nos hacemos es:

¿Podemos construir un algoritmo que descubra automáticamente esas relaciones y las utilice para predecir nuevos casos?

Ese será precisamente el objetivo de los modelos que estudiaremos durante esta unidad. Hoy comenzaremos con el más sencillo: regresión lineal.

Regresión lineal

Exploración del problema

El dataset Advertising contiene datos sobre inversión en publicidad en distintos medios y ventas obtenidas.

Nuestro objetivo será responder:

¿Existe una relación entre la inversión en publicidad en televisión y las ventas?

advertising <- read.csv("https://www.statlearning.com/s/Advertising.csv")

advertising$X <- NULL

head(advertising, 10)

##       TV radio newspaper sales
## 1  230.1  37.8      69.2  22.1
## 2   44.5  39.3      45.1  10.4
## 3   17.2  45.9      69.3   9.3
## 4  151.5  41.3      58.5  18.5
## 5  180.8  10.8      58.4  12.9
## 6    8.7  48.9      75.0   7.2
## 7   57.5  32.8      23.5  11.8
## 8  120.2  19.6      11.6  13.2
## 9    8.6   2.1       1.0   4.8
## 10 199.8   2.6      21.2  10.6

Antes de modelar, necesitamos entender los datos.

¿Cuántas observaciones tenemos?
¿Cuántas variables?
¿Cuál parece ser la variable respuesta?

Esto lo lograremos con la herramienta str(), que resume de forma compacta el tipo de objeto, el número de observaciones y variables si aplica, y un vistazo a los primeros elementos de cada componente.

str(advertising)

## 'data.frame':    200 obs. of  4 variables:
##  $ 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 ...

TV, radio, newspaper: inversión publicitaria en miles de dólares, e.g., TV = 230.1 significa 230.100 dólares invertidos en publicidad televisiva para ese mercado.

sales: ventas del producto en miles de unidades (thousands of units), e.g., sales = 22.1 significa 22.100 unidades vendidas.

Sabemos que en datos reales no podemos esperar relaciones perfectas, porque siempre existe variabilidad no explicada entre observaciones.

Por esta razón, antes de construir un modelo, utilizamos visualizaciones que nos permitan identificar tendencias generales y patrones que no son evidentes al observar solo la tabla.

Cada observación del dataset puede interpretarse como un par (X, Y), lo que en el plano genera una colección de puntos, por lo que usaremos un diagrama de dispersión o “nube de puntos” para visualizar.

plot(advertising$TV,
     advertising$sales,
     pch = 19,
     col = "steelblue",
     xlab = "Publicidad en TV (miles de dólares)",
     ylab = "Ventas (miles de unidades)",
     main = "Ventas vs. inversión en TV")

La nube de puntos sugiere una relación entre variables, pero aún no permite cuantificarla de forma precisa ni utilizarla para predicción.

Por eso, nuestro objetivo será que R pueda encontrar una línea que resuma de la mejor manera esa nube de puntos.

Modelo de regresión lineal `lm`

En R, la regresión lineal se ajusta mediante la función lm(). La sintaxis general es:

respuesta ~ predictor

que se lee como “explicar la respuesta utilizando el predictor”.

modelo <- lm(sales ~ TV, data = advertising)

modelo

## 
## Call:
## lm(formula = sales ~ TV, data = advertising)
## 
## Coefficients:
## (Intercept)           TV  
##     7.03259      0.04754

str(modelo)

## List of 12
##  $ coefficients : Named num [1:2] 7.0326 0.0475
##   ..- attr(*, "names")= chr [1:2] "(Intercept)" "TV"
##  $ residuals    : Named num [1:200] 4.13 1.25 1.45 4.27 -2.73 ...
##   ..- attr(*, "names")= chr [1:200] "1" "2" "3" "4" ...
##  $ effects      : Named num [1:200] -198.31 57.57 1.08 3.99 -2.98 ...
##   ..- attr(*, "names")= chr [1:200] "(Intercept)" "TV" "" "" ...
##  $ rank         : int 2
##  $ fitted.values: Named num [1:200] 17.97 9.15 7.85 14.23 15.63 ...
##   ..- attr(*, "names")= chr [1:200] "1" "2" "3" "4" ...
##  $ assign       : int [1:2] 0 1
##  $ qr           :List of 5
##   ..$ qr   : num [1:200, 1:2] -14.1421 0.0707 0.0707 0.0707 0.0707 ...
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:200] "1" "2" "3" "4" ...
##   .. .. ..$ : chr [1:2] "(Intercept)" "TV"
##   .. ..- attr(*, "assign")= int [1:2] 0 1
##   ..$ qraux: num [1:2] 1.07 1.09
##   ..$ pivot: int [1:2] 1 2
##   ..$ tol  : num 1e-07
##   ..$ rank : int 2
##   ..- attr(*, "class")= chr "qr"
##  $ df.residual  : int 198
##  $ xlevels      : Named list()
##  $ call         : language lm(formula = sales ~ TV, data = advertising)
##  $ terms        :Classes 'terms', 'formula'  language sales ~ TV
##   .. ..- attr(*, "variables")= language list(sales, TV)
##   .. ..- attr(*, "factors")= int [1:2, 1] 0 1
##   .. .. ..- attr(*, "dimnames")=List of 2
##   .. .. .. ..$ : chr [1:2] "sales" "TV"
##   .. .. .. ..$ : chr "TV"
##   .. ..- attr(*, "term.labels")= chr "TV"
##   .. ..- attr(*, "order")= int 1
##   .. ..- attr(*, "intercept")= int 1
##   .. ..- attr(*, "response")= int 1
##   .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
##   .. ..- attr(*, "predvars")= language list(sales, TV)
##   .. ..- attr(*, "dataClasses")= Named chr [1:2] "numeric" "numeric"
##   .. .. ..- attr(*, "names")= chr [1:2] "sales" "TV"
##  $ model        :'data.frame':   200 obs. of  2 variables:
##   ..$ sales: num [1:200] 22.1 10.4 9.3 18.5 12.9 7.2 11.8 13.2 4.8 10.6 ...
##   ..$ TV   : num [1:200] 230.1 44.5 17.2 151.5 180.8 ...
##   ..- attr(*, "terms")=Classes 'terms', 'formula'  language sales ~ TV
##   .. .. ..- attr(*, "variables")= language list(sales, TV)
##   .. .. ..- attr(*, "factors")= int [1:2, 1] 0 1
##   .. .. .. ..- attr(*, "dimnames")=List of 2
##   .. .. .. .. ..$ : chr [1:2] "sales" "TV"
##   .. .. .. .. ..$ : chr "TV"
##   .. .. ..- attr(*, "term.labels")= chr "TV"
##   .. .. ..- attr(*, "order")= int 1
##   .. .. ..- attr(*, "intercept")= int 1
##   .. .. ..- attr(*, "response")= int 1
##   .. .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
##   .. .. ..- attr(*, "predvars")= language list(sales, TV)
##   .. .. ..- attr(*, "dataClasses")= Named chr [1:2] "numeric" "numeric"
##   .. .. .. ..- attr(*, "names")= chr [1:2] "sales" "TV"
##  - attr(*, "class")= chr "lm"

La función lm() devuelve un objeto completo que contiene toda la información del modelo ajustado.

Se trata de una lista estructurada y podremos acceder a cada uno de sus componentes utilizando $.

modelo$fitted.values[1:10] # predicciones del modelo para cada observación

##         1         2         3         4         5         6         7         8 
## 17.970775  9.147974  7.850224 14.234395 15.627218  7.446162  9.765950 12.746498 
##         9        10 
##  7.441409 16.530414

modelo$residuals[1:10] # diferencia entre el valor real y la predicción

##          1          2          3          4          5          6          7 
##  4.1292255  1.2520260  1.4497762  4.2656054 -2.7272181 -0.2461623  2.0340496 
##          8          9         10 
##  0.4535023 -2.6414087 -5.9304143

Todo lo que aprendimos sobre listas sigue siendo válido, ahora aplicándolo a una lista particular que fue resultado de un modelo estadístico.

Interpretación visual del modelo

Para este gráfico usaremos ggplot2, aprovechando para practicar cómo agregar títulos y etiquetas de ejes con labs().

Un gráfico sin título ni ejes rotulados no debería considerarse un entregable terminado, ni siquiera cuando es solo para exploración propia.

# install.packages("ggplot2")
library(ggplot2)

ggplot(advertising, aes(x = TV, y = sales)) + # define el dataset y mapea TV al eje X, sales al eje Y
  geom_point(color = "steelblue") +           # dibuja los puntos observados (cada fila del dataset)
  geom_smooth(method = "lm", se = FALSE,      # ajusta y dibuja la recta de regresión lineal
              color = "red") +                # sin banda de error (se = FALSE), en color rojo

  labs(
    title = "Relación entre publicidad en TV y ventas",     # título principal del gráfico
    subtitle = "Línea de regresión ajustada y residuos",    # subtítulo aclaratorio
    x = "Inversión en TV",                                  # etiqueta del eje X
    y = "Ventas"                                            # etiqueta del eje Y
  )

## `geom_smooth()` using formula = 'y ~ x'

¿Qué representa exactamente la línea roja? Porque no logra unir los puntos…

Esa línea representa la mejor aproximación lineal de la relación entre variables.

El residuo es la diferencia entre el valor observado y el valor predicho por el modelo.

# install.packages("ggplot2")
library(ggplot2)

ggplot(advertising, aes(x = TV, y = sales)) + # define el dataset y mapea TV al eje X, sales al eje Y
  geom_point(color = "steelblue") +           # dibuja los puntos observados (cada fila del dataset)
  geom_smooth(method = "lm", se = FALSE,      # ajusta y dibuja la recta de regresión lineal
              color = "red") +                # sin banda de error (se = FALSE), en color rojo
  geom_segment(aes(xend = TV,                 # traza un segmento por punto: mismo TV como fin en X...
                   yend = fitted(modelo)),   # ...y como fin en Y, el valor predicho por el modelo (residuo visual)
               alpha = 0.3) +                 # transparencia baja para que no tapen el resto del gráfico
  labs(
    title = "Relación entre publicidad en TV y ventas",     # título principal del gráfico
    subtitle = "Línea de regresión ajustada y residuos",    # subtítulo aclaratorio
    x = "Inversión en TV",                                  # etiqueta del eje X
    y = "Ventas"                                            # etiqueta del eje Y
  )

## `geom_smooth()` using formula = 'y ~ x'

Gráficamente, es la distancia vertical (los segmentos grises) entre el punto observado y la recta.

La regresión lineal busca minimizar estas distancias en conjunto.

La limitante es que si sumáramos todas esas distancias, algunas positivas y otras negativas se cancelarían. Por eso la regresión lineal utiliza los cuadrados de los residuos y elige la recta que minimiza la suma de esos errores al cuadrado.

Podríamos decir entonces que el modelo de regresión lineal buscará una recta que, en promedio, se equivoque lo menos posible.

Evaluación del modelo

La ecuación del modelo es

\[ Y = \beta_{0} + \beta_{1} X \]

La función summary() devuelve los parámetros de la recta que acabamos de ajustar.

summary(modelo)

## 
## Call:
## lm(formula = sales ~ TV, data = advertising)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.3860 -1.9545 -0.1913  2.0671  7.2124 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 7.032594   0.457843   15.36   <2e-16 ***
## TV          0.047537   0.002691   17.67   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.259 on 198 degrees of freedom
## Multiple R-squared:  0.6119, Adjusted R-squared:  0.6099 
## F-statistic: 312.1 on 1 and 198 DF,  p-value: < 2.2e-16

Este resultado trae más información de la que vamos a usar por ahora.

Nota: Por ahora vamos a dejar de lado a las columnas Std. Error, t value, Pr(>|t|) y un F-statistic que se ven al final de la salida. Estos elementos están relacionados con qué tan “confiable” es cada coeficiente, a modo de intuición rápida, nos quedamos con que mientras más chico el valor de la última columna (Pr(>|t|)), más confianza tenemos en que esa variable realmente aporta información y no es solamente ruido.

Nos vamos a quedar únicamente con tres piezas: los coeficientes, el Residual standard error y el Multiple R-squared.

donde:

Intercepto (Intercept): es el valor predicho de las ventas cuando la inversión en TV es 0.
Pendiente (TV): indica cuánto cambian, en promedio, las ventas por cada unidad adicional invertida en publicidad en TV.

Ahora que tenemos un modelo ajustado, podemos preguntarnos qué tan bien representa los datos observados.

Residual Standard Error (RSE)

Encontramos: Residual standard error: 3.259

Esta es una medida del error típico que comete el modelo al hacer predicciones.

En promedio, nuestras predicciones se equivocan aproximadamente 3.259 miles de unidades de ventas.

Mientras menor sea este valor, mejor se ajusta la recta a los datos.

¡OJO! No significa que todas las predicciones se equivoquen exactamente 3.259. Es un error promedio.

Multiple R-squared

Encontraremos también: Multiple R-squared: 0.61

El coeficiente de determinación, o R², indica qué proporción de la variabilidad observada en la variable respuesta logra explicar el modelo.

Por ejemplo, R² = 0.61 se interpreta como:

Aproximadamente el 61% de la variación en las ventas puede explicarse mediante la inversión en publicidad en televisión.

El 39% restante se debe a otros factores que este modelo no está considerando o al ruido inherente de los datos.

Con lm() y summary() logramos ajustar un modelo de regresión lineal a los datos, además de identificar que:

Los coeficientes β0 y β1 nos dicen cuál es la recta que ajustó el modelo.
El Residual Standard Error nos indica qué tan grandes son, en promedio, los errores de predicción.
El R² nos dice qué tan bien esa recta logra explicar la variabilidad observada en los datos.

Predicción

Una vez estimado el modelo, podemos utilizarlo para predecir valores nuevos de la variable respuesta.

Si una empresa invierte 150 unidades monetarias en publicidad en TV, ¿qué predice nuestro modelo?

predict(modelo, newdata = data.frame(TV = 150))

##        1 
## 14.16309

De donde interpretamos que el modelo predice aproximadamente 14-15 miles de unidades de ventas para ese nivel de inversión en TV.

Ojo con algo importante: esta predicción y todas las métricas que vimos (RSE, R²) las calculamos evaluando el modelo sobre los mismos datos con los que lo entrenamos. Más adelante en esta clase vamos a ver qué pasa cuando evaluamos el modelo con datos que nunca vio.

Regresión lineal múltiple

Para las ventas en advertising, hicimos el supuesto de que dependen únicamente de la inversión en televisión, aunque el conjunto de datos advertising contiene más información…

En problemas reales, la variable respuesta suele depender de múltiples factores simultáneamente.

Observemos nuevamente cuáles son las otras variables disponibles.

names(advertising)

## [1] "TV"        "radio"     "newspaper" "sales"

Además de la inversión en televisión (TV), también disponemos de información sobre publicidad en radio (radio) y periódicos (newspaper).

Relación marginal de cada predictor

Antes de combinarlas, veamos cómo se relaciona cada una por separado con las ventas. Esto nos va a servir después para entender por qué algunas variables “aportan poco” al modelo conjunto.

ggplot(advertising, aes(x = radio, y = sales)) +   # define el dataset y mapea radio al eje X, sales al eje Y
  geom_point(color = "darkgreen") +                # dibuja los puntos observados en verde
  geom_smooth(method = "lm", se = FALSE,           # ajusta y dibuja la recta de regresión lineal
              color = "black") +                   # sin banda de error, en color negro
  labs(title = "Ventas vs. inversión en radio",    # título del gráfico
       x = "Inversión en radio (miles de dólares)",                   # etiqueta del eje X
       y = "Ventas (miles de unidades)")                                # etiqueta del eje Y

## `geom_smooth()` using formula = 'y ~ x'

ggplot(advertising, aes(x = newspaper, y = sales)) + # define el dataset y mapea newspaper al eje X, sales al eje Y
  geom_point(color = "orange") +                     # dibuja los puntos observados en naranja
  geom_smooth(method = "lm", se = FALSE,             # ajusta y dibuja la recta de regresión lineal
              color = "black") +                     # sin banda de error, en color negro
  labs(title = "Ventas vs. inversión en periódico",  # título del gráfico
       x = "Inversión en periódico (miles de dólares)",                 # etiqueta del eje X
       y = "Ventas (miles de unidades)")                                  # etiqueta del eje Y

## `geom_smooth()` using formula = 'y ~ x'

Podemos facilitar el análisis de estos gráficos a través de la biblioteca patchwork.

# install.packages("patchwork")
library(patchwork) # Combinar gráficos de ggplot2

g1 <- ggplot(advertising, aes(x = radio, y = sales)) +
  geom_point(color = "darkgreen") +
  geom_smooth(method = "lm", se = FALSE, color = "black") +
  labs(title = "Ventas vs. inversión en radio", x = "Inversión en radio (miles de dólares)", y = "Ventas (miles de unidades)")

g2 <- ggplot(advertising, aes(x = newspaper, y = sales)) +
  geom_point(color = "orange") +
  geom_smooth(method = "lm", se = FALSE, color = "black") +
  labs(title = "Ventas vs. inversión en periódico", x = "Inversión en periódico (miles de dólares)", y = "Ventas (miles de unidades)")

g1 + g2 # ambos gráficos en una sola fila.

## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'

g1 / g2 # un gráfico debajo del otro.

## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'

Ya a simple vista, la nube de puntos de newspaper se ve más dispersa y con una tendencia menos clara que la de TV o radio. Guardamos esa observación para más adelante…

Construyendo el modelo paso a paso: de 1 a 3 predictores

Vamos a ampliar el modelo agregando variables de a una, separadas por el operador +, y comparar cómo cambia el R² en cada paso.

modelo_1 <- lm(sales ~ TV, data = advertising)           # 1 predictor: TV
modelo_2 <- lm(sales ~ TV + radio, data = advertising)   # 2 predictores: TV + radio
modelo_3 <- lm(sales ~ TV + radio + newspaper, data = advertising) # 3 predictores: TV + radio + newspaper

summary(modelo_1)$r.squared

## [1] 0.6118751

summary(modelo_2)$r.squared

## [1] 0.8971943

summary(modelo_3)$r.squared

## [1] 0.8972106

Vas a notar que el R² sube (o al menos se mantiene) cada vez que agregamos una variable, sin importar cuál sea. Esto no es casualidad ni buena suerte: es una propiedad matemática del R². Agregar cualquier variable, incluso una que no tenga relación real con las ventas, nunca puede hacer que el R² baje.

Para comprobarlo, agreguemos una columna de puros números aleatorios, sin ninguna relación con las ventas, y veamos qué pasa:

set.seed(50)

advertising$ruido <- rnorm(nrow(advertising), 1, 5)

modelo_ruido <- lm(sales ~ TV + radio + newspaper + ruido, data = advertising)

summary(modelo_3)$r.squared

## [1] 0.8972106

summary(modelo_ruido)$r.squared

## [1] 0.8973091

El R² sube un poquito aunque ruido no tenga ninguna relación real con las ventas. Esto nos deja una pregunta importante: si el R² siempre sube, ¿cómo decidimos si vale la pena agregar una variable?

R² ajustado

Para resolver este problema existe el R² ajustado, una versión del R² que penaliza agregar variables que no aportan información real. A diferencia del R² normal, el R² ajustado puede bajar cuando agregamos una variable que no ayuda.

summary(modelo_3)$adj.r.squared

## [1] 0.8956373

summary(modelo_ruido)$adj.r.squared

## [1] 0.8952026

Con la variable ruido agregada, el R² ajustado se mantiene prácticamente igual o incluso baja, mientras que el R² normal seguía subiendo. Esa diferencia es justamente la señal de que ruido no está aportando nada real al modelo.

La regla práctica que nos llevamos es: para comparar modelos con distinta cantidad de variables, conviene mirar el R² ajustado en lugar del R² normal.

names(advertising)

## [1] "TV"        "radio"     "newspaper" "sales"     "ruido"

advertising$ruido <- NULL # la eliminamos, ya cumplió su propósito de demostración

names(advertising)

## [1] "TV"        "radio"     "newspaper" "sales"

Comparando combinaciones de predictores

Ahora que tenemos ambas métricas, comparemos distintas combinaciones de predictores para ver si newspaper realmente aporta algo al modelo, tal como sospechábamos en el gráfico de dispersión.

# modelo_1 - 1 predictor: TV
# modelo_2 - 2 predictores: TV + radio
# modelo_3 - 3 predictores: TV + radio + newspaper

data.frame(
  modelo = c("TV", "TV + radio", "TV + radio + newspaper"),
  R2 = c(
    summary(modelo_1)$r.squared,
    summary(modelo_2)$r.squared, 
    summary(modelo_3)$r.squared),
  R2_ajustado = c(
    summary(modelo_1)$adj.r.squared,
    summary(modelo_2)$adj.r.squared, 
    summary(modelo_3)$adj.r.squared)
)

##                   modelo        R2 R2_ajustado
## 1                     TV 0.6118751   0.6099148
## 2             TV + radio 0.8971943   0.8961505
## 3 TV + radio + newspaper 0.8972106   0.8956373

Observamos que el R² apenas sube al agregar newspaper, mientras que el R² ajustado prácticamente no cambia (o incluso empeora levemente). Esto es consistente con lo que vimos en el gráfico de dispersión: newspaper aporta poca información adicional una vez que ya tenemos TV y radio en el modelo.

Leyendo los coeficientes del modelo final

Antes de pasar al diagnóstico visual, volvamos a summary() pero ahora sobre el modelo completo de 3 predictores, para ver la tabla de coeficientes con su error estándar, t value y Pr(>|t|):

summary(modelo_3)

## 
## Call:
## lm(formula = sales ~ TV + radio + newspaper, data = advertising)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.8277 -0.8908  0.2418  1.1893  2.8292 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  2.938889   0.311908   9.422   <2e-16 ***
## TV           0.045765   0.001395  32.809   <2e-16 ***
## radio        0.188530   0.008611  21.893   <2e-16 ***
## newspaper   -0.001037   0.005871  -0.177     0.86    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.686 on 196 degrees of freedom
## Multiple R-squared:  0.8972, Adjusted R-squared:  0.8956 
## F-statistic: 570.3 on 3 and 196 DF,  p-value: < 2.2e-16

Cuando hay varias variables en el modelo, el coeficiente de cada una ya no se interpreta de forma tan simple y aislada como en el modelo con un solo predictor, porque las variables pueden estar correlacionadas entre sí (esto se conoce como multicolinealidad).

Por ejemplo, si radio y newspaper subieran y bajaran siempre juntos, el modelo tiene dificultad para distinguir cuánto del efecto sobre las ventas corresponde a cada uno por separado, y eso hace que los coeficientes individuales se vuelvan menos confiables o más inestables.

Observamos que la columna Pr(>|t|) de newspaper es notablemente más alta que la de TV y radio, lo cual confirma numéricamente lo que ya sospechábamos por el R² ajustado y el gráfico de dispersión, newspaper no está aportando información estadísticamente significativa una vez que TV y radio ya están en el modelo.

Supuestos para aplicar regresión lineal

La regresión lineal no es válida para cualquier conjunto de datos: hay una serie de supuestos que idealmente deberían cumplirse. Los supuestos clásicos suelen verificarse con tests estadísticos formales, pero en este curso, vamos a acotarnos a verificarlos de forma visual.

R nos da, de manera gratuita, cuatro gráficos de diagnóstico con solo llamar a plot() sobre el modelo:

plot(modelo_3)

Como no es cómodo explorar los gráficos de esta forma, se proponen dos alternativas:

Usar la función par() para visualizar los cuatro gráficos de diagnóstico en una misma ventana.

par(mfrow = c(2, 2)) # para ver los 4 gráficos juntos
plot(modelo_3)

Utilizar la biblioteca ggfortify, que genera automáticamente los mismos gráficos de diagnóstico utilizando ggplot2.

Esto ofrece varias ventajas, como una mejor estética, la posibilidad de personalizar los gráficos con la sintaxis de ggplot2, modificar temas, etiquetas y colores.

# install.packages("ggfortify")
library(ggfortify)

## Warning: package 'ggfortify' was built under R version 4.5.3

autoplot(modelo_3)

## Warning: `fortify(<lm>)` was deprecated in ggplot2 4.0.0.
## ℹ Please use `broom::augment(<lm>)` instead.
## ℹ The deprecated feature was likely used in the ggfortify package.
##   Please report the issue at <https://github.com/sinhrks/ggfortify/issues>.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

## Warning: `aes_string()` was deprecated in ggplot2 3.0.0.
## ℹ Please use tidy evaluation idioms with `aes()`.
## ℹ See also `vignette("ggplot2-in-packages")` for more information.
## ℹ The deprecated feature was likely used in the ggfortify package.
##   Please report the issue at <https://github.com/sinhrks/ggfortify/issues>.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## ℹ The deprecated feature was likely used in the ggfortify package.
##   Please report the issue at <https://github.com/sinhrks/ggfortify/issues>.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Residuals vs Fitted: si la nube de puntos no tiene una forma clara (una especie de nube dispersa sin patrón), es buena señal. Si en cambio ven una curva marcada o un “embudo” que se abre o se cierra, es señal de que el modelo lineal podría no ser adecuado, o de que la variabilidad de los errores no es constante.
Q-Q plot (Normal Q-Q): si los puntos siguen aproximadamente la línea diagonal, los residuos son razonablemente normales. Si se despegan mucho en los extremos, hay que sospechar de este supuesto.
Scale-Location: es una variante del primer gráfico, pero pensada específicamente para detectar si la dispersión de los residuos es constante a lo largo de todos los valores predichos. Se lee igual que Residuals vs Fitted: buscamos una nube pareja de arriba a abajo; si la línea sube o baja de forma marcada, es la misma señal de variabilidad no constante, solo que más fácil de detectar visualmente.
Residuals vs Leverage: este identifica observaciones individuales que podrían estar influyendo demasiado en el modelo (puntos “raros” que, si los sacáramos, cambiarían bastante la recta ajustada).

¿Y en nuestro caso, se cumplen los supuestos? Corriendo plot(modelo_3) sobre Advertising, no del todo:

En Residuals vs Fitted aparece una curva notoria, no una nube pareja sin forma. Esto es una señal real de no-linealidad: el modelo aditivo (TV + radio + newspaper) se está perdiendo algo. En este dataset en particular, lo que falta es un efecto de interacción entre TV y radio (invertir en ambos a la vez genera más ventas que la suma de invertir en cada uno por separado). Un modelo puramente aditivo no puede capturar esta relación entre predictores.
En el Q-Q plot los puntos se despegan de la diagonal en la cola inferior, sugiriendo que los residuos no son perfectamente normales (hay algunas observaciones con ventas bastante más bajas de lo que el modelo predice).

La conclusión honesta para transmitirles: el modelo lineal simple no es perfecto para estos datos aunque funciona razonablemente bien y es un buen punto de partida, pero el diagnóstico visual ya nos avisa que hay estructura en los datos (la interacción) que este modelo no está capturando.

Este ejercicio sirve de ejemplo en donde “ajustar un modelo” y “que el modelo cumpla sus supuestos” son cosas distintas, y que el diagnóstico visual nos advierte acerca de qué tanto podemos confiar.

Dividiendo los datos: entrenamiento y prueba

Hasta ahora evaluamos todos los modelos con los mismos datos que usamos para entrenarlos. Esto es un problema: un modelo puede parecer muy bueno simplemente porque “memorizó” las particularidades de esos datos puntuales, y fallar al recibir datos nuevos. A este fenómeno se le llama sobreajuste (overfitting), y es una de las ideas centrales en Machine Learning.

La solución más simple es dividir el dataset en dos partes:

Conjunto de entrenamiento (train): con el que ajustamos el modelo.
Conjunto de prueba (test): que el modelo nunca ve durante el ajuste, y que usamos solo para evaluar qué tan bien predice datos nuevos.

¿En qué proporción dividimos?

No existe una única proporción “correcta”, las convenciones prácticas son 70/30 y 80/20 (entrenamiento/prueba). La lógica general es: cuantos más datos totales tengamos, menos necesitamos reservar para prueba (con datasets muy grandes, incluso se llega a usar 90/10).

Nuestro dataset Advertising tiene 200 filas, por lo que una división 80/20 es razonable: deja 160 filas para entrenar y 40 para probar.

Dividiendo los datos en R

set.seed(123) # fija la semilla aleatoria para que la partición sea reproducible

n <- nrow(advertising)
indices_entrenamiento <- sample(1:n, size = 0.8 * n)

train <- advertising[indices_entrenamiento, ]
test <- advertising[-indices_entrenamiento, ]

nrow(train)

## [1] 160

nrow(test)

## [1] 40

¿Por qué set.seed()? Sin fijar una semilla, cada vez que se corra el código obtendríamos una partición distinta, y por lo tanto métricas distintas porque estarían comparando particiones diferentes. set.seed() asegura que el “azar” sea siempre el mismo, para que los resultados sean reproducibles en el contexto de la clase.

Entrenando y evaluando con datos nunca vistos

Ahora sí, ajustamos el modelo únicamente con train…

modelo_train <- lm(sales ~ TV + radio, data = train)

summary(modelo_train)$r.squared # R² sobre los datos de entrenamiento

## [1] 0.8896639

…y lo evaluamos prediciendo sobre test, datos que el modelo nunca vio durante el ajuste:

predicciones_test <- predict(modelo_train, newdata = test)

# Calculamos el error de estas predicciones "a mano"
errores <- test$sales - predicciones_test
rmse_test <- sqrt(mean(errores^2))

rmse_test

## [1] 1.405747

Este número RMSE (Root Mean Squared Error, o raíz del error cuadrático medio) es una de las métricas más utilizadas para evaluar modelos de regresión porque resume, en una sola cantidad, el tamaño típico del error de predicción.

En nuestro ejemplo, el RMSE es cercano a 1.4, lo que significa que, en promedio, las predicciones del modelo difieren de las ventas reales en aproximadamente 1.4 unidades (recordemos que sales está medida en miles de unidades).

El RMSE se calcula como la raíz cuadrada del promedio de los errores al cuadrado:

\[ RMSE=\sqrt{\frac{1}{n}\sum_{i=1}^{n}(y_i-\hat y_i)^2} \]

Nos recuerda a otra métrica muy similar: el RSE (Residual Standard Error, o error estándar residual), que obtuvimos al analizar el resumen del modelo. Ambas miden, aproximadamente, el tamaño promedio de los errores y utilizan las mismas unidades que la variable respuesta.

La diferencia está en que:

RSE: se calcula sobre los datos utilizados para entrenar el modelo e incorpora una corrección por los grados de libertad del modelo.
RMSE: se calcula sobre cualquier conjunto de datos (normalmente el conjunto de prueba) y mide directamente el error de predicción sobre observaciones que el modelo nunca vio durante el entrenamiento.

Por esta razón decimos que el RMSE es un “primo” del RSE: ambos cuantifican el tamaño de los errores, pero el primero evalúa la capacidad predictiva del modelo mientras que el segundo describe qué tan bien se ajustó a los datos con los que fue entrenado.

Podemos recuperar ambos valores y compararlos fácilmente:

# Error sobre el conjunto de entrenamiento
rse_train <- summary(modelo_train)$sigma

# Error sobre el conjunto de prueba
rmse_test <- sqrt(mean((test$sales - predicciones_test)^2))

rse_train; rmse_test

## [1] 1.745841

## [1] 1.405747

Comparar el valor del RMSE contra el RSE calculado sobre train es la forma más honesta de saber si nuestro modelo realmente generaliza, o si solo se ajustó bien a los datos que ya conocía.

Si ambos valores son similares, el modelo parece generalizar bien: el error al predecir datos nuevos es parecido al error obtenido durante el entrenamiento.

Si, por el contrario, el RMSE es mucho mayor que el RSE, es una señal de que el modelo podría estar sobreajustado (overfitting), es decir, aprendió demasiado bien las particularidades de los datos de entrenamiento y perdió capacidad para predecir observaciones nuevas.

Con la semilla 123, este código da:

RSE en train: ≈ 1.75
RMSE en test: ≈ 1.41

¿Es bueno este resultado? Hay dos formas de leerlo:

Comparando train vs. test: el error en test es igual o incluso levemente menor que el de train. Esa es justo la señal que buscábamos: si el modelo hubiera memorizado el train, veríamos un error mucho más chico ahí y uno mucho más grande en test. Que el error en datos nunca vistos no sea peor nos dice que el modelo generaliza bien, por lo que no hay señales de sobreajuste con este modelo tan simple.

Nota: que el error de test salga igual o menor que el de train no es lo típico ni lo esperable en general, normalmente el error de test es un poco mayor. En este ejercicio pasa porque el modelo es simple y el split de 40 filas para test, al ser chico, puede tener casos “fáciles” de predecir. No hay que generalizar esta observación puntual como una regla.

Comparando el error contra la escala de la variable:

summary(advertising$sales)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1.60   10.38   12.90   14.02   17.40   27.00

Las ventas (sales) en este dataset van de 1.6 a 27, con una media de aproximadamente 14. Un error promedio de ~1.4 unidades representa, entonces, un error de más o menos el 10% respecto a la media, que un margen razonable para un modelo con solo dos predictores.

En resumen: para este split en particular, el modelo generaliza bien (train y test dan errores similares), aunque como vimos en la sección de supuestos, no es un modelo perfecto, faltaría capturar la interacción entre TV y radio. Ambas cosas pueden ser ciertas a la vez: un modelo puede generalizar razonablemente bien y, al mismo tiempo, tener margen de mejora en cómo está especificado.

Hasta ahora hemos trabajado con predicción de variables numéricas.

Sin embargo, muchos problemas de interés práctico no buscan un número continuo, sino una categoría, como:

clasificar un correo como spam o no
predecir si un cliente comprará o no
determinar si un paciente tiene una enfermedad

Esto nos lleva a la necesidad de modelos de clasificación.

En la siguiente clase estudiaremos la regresión logística, que adapta estas ideas al caso donde la variable respuesta es categórica.

Data Analysis with R Studio

Sesión 6

Act. Casandra Sofía Pejcovich

24 septiembre, 2026

Objetos de datos

I. Vectores

II. Matrices

III. Data frames

IV. Listas

Listas

Introducción a Machine Learning

Regresión lineal

Exploración del problema

Modelo de regresión lineal `lm`

Interpretación visual del modelo

Evaluación del modelo

Residual Standard Error (RSE)

Multiple R-squared

Predicción

Regresión lineal múltiple

Relación marginal de cada predictor

Construyendo el modelo paso a paso: de 1 a 3 predictores

R² ajustado

Comparando combinaciones de predictores

Leyendo los coeficientes del modelo final

Supuestos para aplicar regresión lineal

Dividiendo los datos: entrenamiento y prueba

¿En qué proporción dividimos?

Dividiendo los datos en R

Entrenando y evaluando con datos nunca vistos

Data Analysis with R Studio

Sesión 6

Act. Casandra Sofía Pejcovich

24 septiembre, 2026

Objetos de datos

I. Vectores

II. Matrices

III. Data frames

IV. Listas

Listas

Introducción a Machine Learning

Regresión lineal

Exploración del problema

Modelo de regresión lineal lm

Interpretación visual del modelo

Evaluación del modelo

Residual Standard Error (RSE)

Multiple R-squared

Predicción

Regresión lineal múltiple

Relación marginal de cada predictor

Construyendo el modelo paso a paso: de 1 a 3 predictores

R² ajustado

Comparando combinaciones de predictores

Leyendo los coeficientes del modelo final

Supuestos para aplicar regresión lineal

Dividiendo los datos: entrenamiento y prueba

¿En qué proporción dividimos?

Dividiendo los datos en R

Entrenando y evaluando con datos nunca vistos

Modelo de regresión lineal `lm`