Análisis de datos en R

Workbooks del curso de formación en RTVE (2025)

Author

Javier Álvarez Liébana

1 Clase 1

1.1 🐣 Caso práctico I: airquality (manejo de vectores)

En el paquete {datasets} (ya instalado por defecto) tenemos diversos conjuntos de datos y uno de ellos es airquality. Los datos capturan medidas diarias (n = 153 observaciones) de la calidad del aire en Nueva York, de mayo a septiembre de 1973. Se midieron 6 variables: niveles de ozono, radiación solar, viento, temperatura, mes y día.

library(datasets)
airquality
    Ozone Solar.R Wind Temp Month Day
1      41     190  7.4   67     5   1
2      36     118  8.0   72     5   2
3      12     149 12.6   74     5   3
4      18     313 11.5   62     5   4
5      NA      NA 14.3   56     5   5
6      28      NA 14.9   66     5   6
7      23     299  8.6   65     5   7
8      19      99 13.8   59     5   8
9       8      19 20.1   61     5   9
10     NA     194  8.6   69     5  10
11      7      NA  6.9   74     5  11
12     16     256  9.7   69     5  12
13     11     290  9.2   66     5  13
14     14     274 10.9   68     5  14
15     18      65 13.2   58     5  15
16     14     334 11.5   64     5  16
17     34     307 12.0   66     5  17
18      6      78 18.4   57     5  18
19     30     322 11.5   68     5  19
20     11      44  9.7   62     5  20
21      1       8  9.7   59     5  21
22     11     320 16.6   73     5  22
23      4      25  9.7   61     5  23
24     32      92 12.0   61     5  24
25     NA      66 16.6   57     5  25
26     NA     266 14.9   58     5  26
27     NA      NA  8.0   57     5  27
28     23      13 12.0   67     5  28
29     45     252 14.9   81     5  29
30    115     223  5.7   79     5  30
31     37     279  7.4   76     5  31
32     NA     286  8.6   78     6   1
33     NA     287  9.7   74     6   2
34     NA     242 16.1   67     6   3
35     NA     186  9.2   84     6   4
36     NA     220  8.6   85     6   5
37     NA     264 14.3   79     6   6
38     29     127  9.7   82     6   7
39     NA     273  6.9   87     6   8
40     71     291 13.8   90     6   9
41     39     323 11.5   87     6  10
42     NA     259 10.9   93     6  11
43     NA     250  9.2   92     6  12
44     23     148  8.0   82     6  13
45     NA     332 13.8   80     6  14
46     NA     322 11.5   79     6  15
47     21     191 14.9   77     6  16
48     37     284 20.7   72     6  17
49     20      37  9.2   65     6  18
50     12     120 11.5   73     6  19
51     13     137 10.3   76     6  20
52     NA     150  6.3   77     6  21
53     NA      59  1.7   76     6  22
54     NA      91  4.6   76     6  23
55     NA     250  6.3   76     6  24
56     NA     135  8.0   75     6  25
57     NA     127  8.0   78     6  26
58     NA      47 10.3   73     6  27
59     NA      98 11.5   80     6  28
60     NA      31 14.9   77     6  29
61     NA     138  8.0   83     6  30
62    135     269  4.1   84     7   1
63     49     248  9.2   85     7   2
64     32     236  9.2   81     7   3
65     NA     101 10.9   84     7   4
66     64     175  4.6   83     7   5
67     40     314 10.9   83     7   6
68     77     276  5.1   88     7   7
69     97     267  6.3   92     7   8
70     97     272  5.7   92     7   9
71     85     175  7.4   89     7  10
72     NA     139  8.6   82     7  11
73     10     264 14.3   73     7  12
74     27     175 14.9   81     7  13
75     NA     291 14.9   91     7  14
76      7      48 14.3   80     7  15
77     48     260  6.9   81     7  16
78     35     274 10.3   82     7  17
79     61     285  6.3   84     7  18
80     79     187  5.1   87     7  19
81     63     220 11.5   85     7  20
82     16       7  6.9   74     7  21
83     NA     258  9.7   81     7  22
84     NA     295 11.5   82     7  23
85     80     294  8.6   86     7  24
86    108     223  8.0   85     7  25
87     20      81  8.6   82     7  26
88     52      82 12.0   86     7  27
89     82     213  7.4   88     7  28
90     50     275  7.4   86     7  29
91     64     253  7.4   83     7  30
92     59     254  9.2   81     7  31
93     39      83  6.9   81     8   1
94      9      24 13.8   81     8   2
95     16      77  7.4   82     8   3
96     78      NA  6.9   86     8   4
97     35      NA  7.4   85     8   5
98     66      NA  4.6   87     8   6
99    122     255  4.0   89     8   7
100    89     229 10.3   90     8   8
101   110     207  8.0   90     8   9
102    NA     222  8.6   92     8  10
103    NA     137 11.5   86     8  11
104    44     192 11.5   86     8  12
105    28     273 11.5   82     8  13
106    65     157  9.7   80     8  14
107    NA      64 11.5   79     8  15
108    22      71 10.3   77     8  16
109    59      51  6.3   79     8  17
110    23     115  7.4   76     8  18
111    31     244 10.9   78     8  19
112    44     190 10.3   78     8  20
113    21     259 15.5   77     8  21
114     9      36 14.3   72     8  22
115    NA     255 12.6   75     8  23
116    45     212  9.7   79     8  24
117   168     238  3.4   81     8  25
118    73     215  8.0   86     8  26
119    NA     153  5.7   88     8  27
120    76     203  9.7   97     8  28
121   118     225  2.3   94     8  29
122    84     237  6.3   96     8  30
123    85     188  6.3   94     8  31
124    96     167  6.9   91     9   1
125    78     197  5.1   92     9   2
126    73     183  2.8   93     9   3
127    91     189  4.6   93     9   4
128    47      95  7.4   87     9   5
129    32      92 15.5   84     9   6
130    20     252 10.9   80     9   7
131    23     220 10.3   78     9   8
132    21     230 10.9   75     9   9
133    24     259  9.7   73     9  10
134    44     236 14.9   81     9  11
135    21     259 15.5   76     9  12
136    28     238  6.3   77     9  13
137     9      24 10.9   71     9  14
138    13     112 11.5   71     9  15
139    46     237  6.9   78     9  16
140    18     224 13.8   67     9  17
141    13      27 10.3   76     9  18
142    24     238 10.3   68     9  19
143    16     201  8.0   82     9  20
144    13     238 12.6   64     9  21
145    23      14  9.2   71     9  22
146    36     139 10.3   81     9  23
147     7      49 10.3   69     9  24
148    14      20 16.6   63     9  25
149    30     193  6.9   70     9  26
150    NA     145 13.2   77     9  27
151    14     191 14.3   75     9  28
152    18     131  8.0   76     9  29
153    20     223 11.5   68     9  30

1.1.1 Pregunta 1

¿Cómo averiguar qué representan los datos? Piensa algún comando que nos dé información sobre objetos en R, sabiendo que el nombre del dataset es airquality

Code
? airquality

Haciendo uso de ? ... podemos consultar en el panel de ayuda lo que significa el objeto.

1.1.2 Pregunta 2

Accede solo a los 5 primeros registros de temperaturas. Después accede al primero, segundo, quinto y décimo

Code
# secuencia de 1 a 5
airquality$Temp[1:5] 

# otra forma
airquality$Temp[c(1, 2, 3, 4, 5)]

# primero, segundo, quinto y décimo
airquality$Temp[c(1, 2, 5, 10)]

1.1.3 Pregunta 3

Accede solo a los registros de temperaturas de mayo (tienes las temperaturas guardadas, piensa como acceder a ellos pero ahora usando una condición en lugar de índices concretos). Después accede a los elementos de mayo, abril y junio

Code
airquality$Temp[airquality$Month == 5]

# abril, mayo y junio
airquality$Temp[airquality$Month == 4 | airquality$Month == 5 |
                  airquality$Month == 6]

# otra forma más legible: %in% nos comprueba si 
# los valores están dentro de una lista permitida
airquality$Temp[airquality$Month %in% c(4, 5, 6)]

1.1.4 Pregunta 4

¿Cuántos registros tenemos de mayo? ¿Y de abril?

Code
# Una forma para registros de mayo
sum(airquality$Month == 5)

# Otra forma: la longitud de un vector
length(airquality$Temp[airquality$Month == 5])

# ídem en abril
sum(airquality$Month == 4)

1.1.5 Pregunta 5

Construye una nueva variable date con la fecha de cada registro (combinando año, mes y día), sabiendo que todos los datos son de 1973. Pista: para construir una fecha antes debes tener un vector de textos (por ejemplo, “1973-01-01”)

Code
# variable de tipo date
library(lubridate)
dates <- as_date(glue("{1973}-{airquality$Month}-{airquality$Day}"))

1.1.6 Pregunta 6

Crea una nueva variable temp_celsius con la temperatura en ºC (sabiendo que se calcula como \(celsius = (fahr - 32) * (5/9)\)). Tras ello calcula cuántos días de junio superaron los 30 grados ºC.

Code
# Temperatura en celsius
temp_celsius <- (airquality$Temp - 32) * (5/9)
temp_celsius 

# una forma
sum(temp_celsius[airquality$Month == 6] > 30)

# otra forma
length(temp_celsius[airquality$Month == 6 & temp_celsius > 30])

1.1.7 Pregunta 7

Selecciona aquellos datos que no sean ni de julio ni de agosto.

Code
airquality_tb[airquality_tb$Month != 7 & airquality_tb$Month != 8, ]

# otra opción
airquality_tb[!(airquality_tb$Month %in% c(7, 8)), ]

1.1.8 Pregunta 8

¿Cuál fue la media de temperatura del mes de agosto?

Code
# media en agosto
mean(airquality$Temp[airquality$Month == 8], na.rm = TRUE)
mean(temp_celsius[airquality$Month == 8], na.rm = TRUE)

1.2 🐣 Caso práctico II: covid (tidy data)

En el fichero messy_covid_data.csv tienes un archivo con la cantidad de casos reportados durante la pandemia covid pero en formato messy: el nombre de las columnas codifica el sexo (H hombre, M mujer, NC no consta) y el grupo etario (0-9, 10-19, 20-29, 30-39, 40-49, 50-59, 60-69, 70-79, ≥80 años y NC no consta).

Las soluciones serán asumiendo que no se conoce aún las opciones de tidyverse para modificar y filtrar tablas.

library(readr)
datos <- read_csv(file = "./datos/messy_covid_data.csv")
Rows: 9486 Columns: 31
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr   (1): provincia_iso
dbl  (29): 0-9_H, 10-19_H, 20-29_H, 30-39_H, 40-49_H, 50-59_H, 60-69_H, 70-7...
dttm  (1): fecha

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
datos
# A tibble: 9,486 × 31
   provincia_iso fecha               `0-9_H` `10-19_H` `20-29_H` `30-39_H`
   <chr>         <dttm>                <dbl>     <dbl>     <dbl>     <dbl>
 1 A             2020-03-01 00:00:00       0         0         0         0
 2 AL            2020-03-01 00:00:00       0         0         0         0
 3 B             2020-03-01 00:00:00       0         0         0         0
 4 BA            2020-03-01 00:00:00       0         1         0         0
 5 BI            2020-03-01 00:00:00       0         0         0         0
 6 C             2020-03-01 00:00:00       0         0         0         0
 7 CA            2020-03-01 00:00:00       0         0         0         0
 8 CO            2020-03-01 00:00:00       0         0         0         0
 9 GI            2020-03-01 00:00:00       0         0         0         0
10 GR            2020-03-01 00:00:00       0         0         0         0
# ℹ 9,476 more rows
# ℹ 25 more variables: `40-49_H` <dbl>, `50-59_H` <dbl>, `60-69_H` <dbl>,
#   `70-79_H` <dbl>, `80-Inf_H` <dbl>, `NC-NC_H` <dbl>, `0-9_M` <dbl>,
#   `10-19_M` <dbl>, `20-29_M` <dbl>, `30-39_M` <dbl>, `40-49_M` <dbl>,
#   `50-59_M` <dbl>, `60-69_M` <dbl>, `70-79_M` <dbl>, `80-Inf_M` <dbl>,
#   `NC-NC_M` <dbl>, `0-9_NC` <dbl>, `10-19_NC` <dbl>, `20-29_NC` <dbl>,
#   `30-39_NC` <dbl>, `40-49_NC` <dbl>, `50-59_NC` <dbl>, `60-69_NC` <dbl>, …

1.2.1 Pregunta 1

Razona por qué no es tidydata y conviérte a tidy data en cuatro columnas: provincia_iso, fecha, grupo, casos.

Code
# pivot_longer(): le indicamos que columnas no queremos pivotar,
# la variable a donde irán los nombres y cómo se llamará la variable
# de casos
# le incluimos `values_drop_na = TRUE` para que elimine los ausentes
# si los hubiese
library(tidyr)
tidy_covid <-
  datos |> 
  pivot_longer(cols = -c("provincia_iso", "fecha"),
               names_to = "grupo", values_to = "casos",
               values_drop_na = TRUE)

1.2.2 Pregunta 2

Elimina los registros sin casos.

Code
# pivot_longer(): le indicamos que columnas no queremos pivotar,
# la variable a donde irán los nombres y cómo se llamará la variable
# de casos
# le incluimos `values_drop_na = TRUE` para que elimine los ausentes
# si los hubiese
tidy_covid <- tidy_covid[tidy_covid$casos > 0, ]

1.2.3 Pregunta 3

Una de las columnas (grupo) la tenemos codificada a veces como cosa-cosa_cosa, otras como 80-Inf_cosa (sin cota superior en edad), otras como NC-NC_cosa (NC codifica los “no consta”). Intenta separar dicha columna para generar tres columnas nuevas edad_inf, edad_sup y sexo de manera adecuada. Recuerda que NC es un ausente. Pista: trata antes la variable grupo y reemplaza reemplazar los "Inf" y "NC" por "NA", para que luego al convertirlos pueda pasarlos a numéricos sin problema

Code
# Vamos a reemplazar los `"Inf"` y `"NC"` por `"NA"`, para que luego al
# convertirlos puedas pasarlos a `NA` numéricos
library(stringr)
tidy_covid$grupo <- str_replace_all(tidy_covid$grupo, "Inf|NC", "NA")
tidy_covid2 <-
  tidy_covid |>
  separate(col = "grupo", into = c("edad_inf", "edad_sup", "sexo"), convert = TRUE)

1.2.4 Pregunta 4

Incorpora una nueva variable mes_year a la tabla que codifique conjuntamente el mes y el año (por ejemplo, cualquier día de enero de 2020 será algo similar a “1-2020” y cualquier día de febrero del 2021 será “2-2021”).

Code
library(glue)
library(lubridate)

Attaching package: 'lubridate'
The following objects are masked from 'package:base':

    date, intersect, setdiff, union
Code
tidy_covid2$mes_year <-
  glue("{month(tidy_covid2$fecha)}-{year(tidy_covid2$fecha)}")

1.2.5 Pregunta 5

Haciendo uso de esa variable de grupo mes_year y del vector de provincias permitidas que aparece debajo (codificadas por su código ISO) obtén un resumen que nos devuelva en un tibble, para cada provincia permitida y cada mes-año, la media de casos (sin importar edad ni sexo)

Solo necesitamos recorrer cada provincia permitida y cada valor mes_year distinto, así que creamos ambos vectores. También creamos un dataset vacío resumen_mes_provincia donde, en cada iteración, pegaremos una fila con su resumen. Una vez creados esos objetos, el bucle simplemente deberá recorrer cada valor de provincia y de mes-año: para el mes-año i y la provincia j, filtramos solo los casos para esos valores concretos (de todo el vector de casos nos quedamos con aquellos que cumplan dicha condición, eliminando también provincias ausentes). Tras filtrar guardamos la media ese filtro de casos

Code
# madrid, barcelona, sevilla
prov_permitidas <- c("M", "B", "SE")

filtro_covid <-
  tidy_covid2[!is.na(tidy_covid2$provincia_iso) &
                tidy_covid2$provincia_iso %in% prov_permitidas, ]
resumen_mes_provincia <-
  tibble("M" = mean(filtro_covid$casos[filtro_covid$provincia_iso == "M"],
                    na.rm = TRUE),
         "B" = mean(filtro_covid$casos[filtro_covid$provincia_iso == "B"],
                    na.rm = TRUE),
         "SE" = mean(filtro_covid$casos[filtro_covid$provincia_iso == "SE"],
                    na.rm = TRUE))

1.3 🐣 Caso práctico III: starwars (dplyr)

Vamos a usar al dataset starwars del paquete {dplyr}

library(dplyr)

Attaching package: 'dplyr'
The following objects are masked from 'package:stats':

    filter, lag
The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union
starwars
# A tibble: 87 × 14
   name     height  mass hair_color skin_color eye_color birth_year sex   gender
   <chr>     <int> <dbl> <chr>      <chr>      <chr>          <dbl> <chr> <chr> 
 1 Luke Sk…    172    77 blond      fair       blue            19   male  mascu…
 2 C-3PO       167    75 <NA>       gold       yellow         112   none  mascu…
 3 R2-D2        96    32 <NA>       white, bl… red             33   none  mascu…
 4 Darth V…    202   136 none       white      yellow          41.9 male  mascu…
 5 Leia Or…    150    49 brown      light      brown           19   fema… femin…
 6 Owen La…    178   120 brown, gr… light      blue            52   male  mascu…
 7 Beru Wh…    165    75 brown      light      blue            47   fema… femin…
 8 R5-D4        97    32 <NA>       white, red red             NA   none  mascu…
 9 Biggs D…    183    84 black      light      brown           24   male  mascu…
10 Obi-Wan…    182    77 auburn, w… fair       blue-gray       57   male  mascu…
# ℹ 77 more rows
# ℹ 5 more variables: homeworld <chr>, species <chr>, films <list>,
#   vehicles <list>, starships <list>

1.3.1 Pregunta 1

Selecciona del conjunto de starwars (en el paquete {dplyr}) solo los personajes que sean androides o cuyo valor en species sea desconocido

Code
starwars |>
  filter(species == "Droid" | is.na(species))

1.3.2 Pregunta 2

Selecciona del conjunto de starwars solo los personajes cuyo peso esté entre 65 y 90 kg.

Code
starwars |> filter(between(mass, 65, 90))

1.3.3 Pregunta 3

Tras limpiar de ausentes en todas las variables, selecciona del conjunto de starwars solo los personajes que sean humanos y que vengan de Tatooine

Code
starwars |>
  drop_na() |> 
  filter(species == "Human" & homeworld == "Tatooine")

1.3.4 Pregunta 4

Selecciona solo los personajes que sean humanos y de ojos marrones, para después ordernarlos en altura descendente y peso ascendente.

Code
starwars |>
  filter(eye_color == "brown" & species == "Human") |> 
  arrange(height, desc(mass))

1.3.5 Pregunta 5

Extrae aleatoriamente 10 personajes pero de forma que la probabilidad de que salga cada uno sea proporcional a su peso (más pesados, más probable), teniendo en cuenta que no podemos tener probabilidades ausentes.

Code
starwars |>
  drop_na(mass) |> 
  slice_sample(n = 10, weight_by = mass)

1.3.6 Pregunta 6

De los personajes que son humanos y miden más de 160 cm, elimina duplicados en color de ojos, elimina ausentes en peso, selecciona los 3 más altos, y orden de mayor a menor peso. Devuelve la tabla.

Code
starwars |>
  filter(species == "Human" & height > 160) |> 
  distinct(eye_color, .keep_all = TRUE) |> 
  drop_na(mass) |> 
  slice_max(height, n = 3) |> 
  arrange(desc(mass))

1.3.7 Pregunta 7

Filtra el conjunto de personajes y quédate solo con aquellos que en la variable height no tengan un dato ausente. Con los datos obtenidos del filtro anterior, selecciona solo las variables name, height, así como todas aquellas variables que CONTENGAN la palabra color en su nombre.

Code
starwars_2 <-
  starwars |> 
  drop_na(height) |> 
  select(name, height, contains("color"))

1.3.8 Pregunta 8

Con los datos obtenidos del ejercicio anterior, traduce el nombre de las columnas a castellano.

Code
starwars_2 |> 
  rename(nombre = name, altura = height, color_pelo = hair_color,
         color_piel = skin_color, color_ojos = eye_color)

1.3.9 Pregunta 9

Selecciona solo las variables nombre, altura y así como todas aquellas variables relacionadas con el color, a la vez que te quedas solo con aquellos que no tengan ausente en la altura.

Code
starwars |> 
  select(name, height, contains("color")) |> 
  drop_na(height)

1.3.10 Pregunta 10

Del dataset original, selecciona solo las variables numéricas y de tipo texto. Tras ello define una nueva variable llamada under_18 que nos recategorice la variable de edad: TRUE si es menor de edad y FALSE en caso contrario

Code
starwars |> 
  select(where(is.numeric) | where(is.character)) |> 
  mutate(under_18 = birth_year < 18)

1.3.11 Pregunta 11

Del dataset original, incluye una columna que calcule el IMC. Tras ello, crea una nueva variable que valga NA si no es humano, delgadez por debajo de 18, normal entre 18 y 30, sobrepeso por encima de 30.

Code
starwars |> 
  mutate(IMC = mass / ((height/100)^2),
         IMC_recat = case_when(species != "Human" ~ NA,
                               IMC < 18 ~ "delgadez",
                               IMC < 30 ~ "normal",
                               TRUE ~ "sobrepeso"),
         .after = name)

1.3.12 Pregunta 12

Calcula cuántos personajes hay de cada especie, ordenados de más a menor frecuencia.

Code
starwars |> count(species, sort = TRUE)

1.3.13 Pregunta 13

Tras eliminar ausentes en las variables de peso y estatura, añade una nueva variable que nos calcule el IMC de cada personaje, y determina el IMC medio de nuestros personajes desagregada por sexo

Code
starwars |>
  drop_na(mass, height) |> 
  mutate(IMC = mass / ((height/100)^2)) |> 
  summarise(IMC_medio = mean(IMC), .by = sex)

1.3.14 Pregunta 14

Obtén el personaje más joven por cada sexo.

Code
starwars |>
  slice_min(birth_year, by = sex)

1.3.15 Pregunta 15

Obtén la edad del personaje más joven y más viejo de cada sexo.

Code
starwars |>
  drop_na(birth_year) |>
  summarise(min(birth_year), max(birth_year), .by = sex)

1.3.16 Pregunta 16

Determina la cantidad de personajes en cada década (echa un vistazo a round(), primero sin desagregar y luego desagregado por sexo.

Code
starwars |>
  count(birth_decade = round(birth_year, -1))

:::

:::

1.4 🐣 Caso práctico IV: aviones de NYC (join)

Para los ejercicios usaremos las tablas disponibles en el paquete {nycflights13} (echa un vistazo antes)

library(nycflights13)
  • airlines: nombre de aerolíneas (con su abreviatura).
  • airports: datos de aeropuertos (nombres, longitud, latitud, altitud, etc).
  • flights: datos de vuelos.
  • planes: datos de los aviones.
  • weather: datos meteorológicos horarios de las estaciones LGA, JFK y EWR.

1.4.1 Pregunta 1

Del paquete {nycflights13} cruza la tabla flights con airlines. Queremos mantener todos los registros de vuelos, añadiendo la información de las aerolíneas a dicha tabla.

Code
flights_airlines <-
  flights |> 
  left_join(airlines, by = "carrier")

1.4.2 Pregunta 2

A la tabla obtenida del cruce del apartado anterior, cruza después con los datos de los aviones en planes, pero incluyendo solo aquellos vuelos de los que tengamos información de sus aviones (y viceversa).

Code
flights_airlines_planes <- 
  flights_airlines |> 
  inner_join(planes, by = "tailnum")

1.4.3 Pregunta 3

Repite el ejercicio anterior pero conservando ambas variables year (en una es el año del vuelo, en la otra es el año de construcción del avión), y distinguiéndolas entre sí

Code
flights_airlines_planes <- 
  flights_airlines |> 
  inner_join(planes, by = "tailnum",
             suffix = c("_flight", "_build_aircraft"))

1.4.4 Pregunta 4

Al cruce obtenido del ejercicio anterior incluye la longitud y latitud de los aeropuertos en airports, distinguiendo entre la latitud/longitud del aeropuerto en destino y en origen.

Code
flights_airlines_planes %>%
  left_join(airports %>% select(faa, lat, lon),
            by = c("origin" = "faa")) |> 
  rename(lat_origin = lat, lon_origin = lon) |> 
  left_join(airports %>% select(faa, lat, lon),
            by = c("dest" = "faa")) |> 
  rename(lat_dest = lat, lon_dest = lon)
# A tibble: 284,170 × 32
   year_flight month   day dep_time sched_dep_time dep_delay arr_time
         <int> <int> <int>    <int>          <int>     <dbl>    <int>
 1        2013     1     1      517            515         2      830
 2        2013     1     1      533            529         4      850
 3        2013     1     1      542            540         2      923
 4        2013     1     1      544            545        -1     1004
 5        2013     1     1      554            600        -6      812
 6        2013     1     1      554            558        -4      740
 7        2013     1     1      555            600        -5      913
 8        2013     1     1      557            600        -3      709
 9        2013     1     1      557            600        -3      838
10        2013     1     1      558            600        -2      849
# ℹ 284,160 more rows
# ℹ 25 more variables: sched_arr_time <int>, arr_delay <dbl>, carrier <chr>,
#   flight <int>, tailnum <chr>, origin <chr>, dest <chr>, air_time <dbl>,
#   distance <dbl>, hour <dbl>, minute <dbl>, time_hour <dttm>, name <chr>,
#   year_build_aircraft <int>, type <chr>, manufacturer <chr>, model <chr>,
#   engines <int>, seats <int>, speed <int>, engine <chr>, lat_origin <dbl>,
#   lon_origin <dbl>, lat_dest <dbl>, lon_dest <dbl>

1.4.5 Pregunta 5

Filtra de airports solo aquellos aeropuertos de los que salgan vuelos. Repite el proceso filtrado solo aquellos a los que lleguen vuelos

Code
airports |> 
  semi_join(flights, by = c("faa" = "origin"))
airports |> 
  semi_join(flights, by = c("faa" = "dest"))

1.4.6 Pregunta 6

¿De cuántos vuelos no disponemos información del avión? Elimina antes los vuelos que no tengan identificar (diferente a NA) del avión

Code
flights |> 
  drop_na(tailnum) |>
  anti_join(planes, by = "tailnum") |>
  count(tailnum, sort = TRUE) # de mayor a menor ya de paso

2 Clase 2

2.1 🐣 Caso práctico I: renta de municipios

2.1.1 Pregunta 1

En el archivo municipios.csv tenemos guardada la información de los municipios de España a fecha de 2019. La variable LAU_code representa el código como unidad administrativa local según la estandarización de la UE (ver más). La variable codigo_ine está construida uniendo el código de la provincia y el de la comunidad autónoma.

Por otro lado el archivo renta_mun.csv contiene datos de la renta per capita de cada unidad administrativa (municipios, distritos, provincias, comunidades autonónomas, etc) para diferentes años.

Carga ambos de manera adecuada

Code
library(readr)
mun_data <- read_csv(file = "./datos/municipios.csv")
Rows: 8212 Columns: 8
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (8): codauto, ine.ccaa.name, cpro, ine.prov.name, cmun, name, LAU_CODE, ...

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Code
renta_mun <- read_csv2(file = "./datos/renta_mun.csv")
ℹ Using "','" as decimal and "'.'" as grouping mark. Use `read_delim()` for more control.
Rows: 55273 Columns: 7── Column specification ────────────────────────────────────────────────────────
Delimiter: ";"
chr (2): Unidad, codigo_ine
dbl (5): 2019, 2018, 2017, 2016, 2015
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

2.1.2 Pregunta 2

Investiga el paquete {janitor} y haz uso de la función clean_names() para normalizar nombres de variables

Code
mun_data <-
  mun_data |> 
  janitor::clean_names()
renta_mun <-
  renta_mun |> 
  janitor::clean_names()

2.1.3 Pregunta 3

¿Son tidydata? En caso de que no convierte el que toque de manera adecuada (sin ausentes y cada dato del tipo correcto)

Code
renta_mun_tidy <-
  renta_mun |> 
  pivot_longer(cols = contains("x"), names_to = "year",
               values_to = "renta", names_prefix = "x",
               names_transform = list(year = as.numeric),
               values_drop_na = TRUE)

2.1.4 Pregunta 4

Si te fijas en la tabla anterior, tenemos datos de diferentes unidades administrativas que no siempre son municipios. Sabiendo que todos los municipios tienen un código de 5 caracteres (que representan todos ellos números), filtra sólo aquellos registros que correspondan a unidades municipales.

Code
renta_mun_tidy <-
  renta_mun_tidy |>
  filter(str_detect(codigo_ine, pattern = "[0-9]{5}") & 
           str_length(codigo_ine) == 5)

2.1.5 Pregunta 5

A continuación separa adecuadamente la variable de unidad administrativa en dos columnas: una con el código (que ya tiene, por lo que debe eliminar uno de los dos) y el nombre. Elimina los espacios sobrantes (echa un vistazo a las opciones del paquete {stringr}).

Code
renta_mun_tidy <-
  renta_mun_tidy |>
  separate(col = "unidad", into = c("cod_rm", "name"), sep = 5) |> 
  select(-cod_rm) |> 
  mutate(name = str_trim(name)) 

2.1.6 Pregunta 6

¿En qué año fue mayor la renta media? ¿Y más baja? ¿Cuál fue la renta mediana de los municipios de España en 2019?

Code
summary_renta <-
  renta_mun_tidy |> 
  summarise("mean_renta" = mean(renta, na.rm = TRUE),
            .by = year)
summary_renta |>
  slice_min(mean_renta, n = 1)


summary_renta |>
  slice_max(mean_renta, n = 1)

renta_mun_tidy |> 
  filter(year == 2019) |> 
  summarise("median_renta" = median(renta, na.rm = TRUE))

2.1.7 Pregunta 7

Haz lo que consideres para obtener el NOMBRE de la provincia con la renta media más alta en 2019 y la más baja

Code
summary_by_prov <-
  renta_mun_tidy |> 
  filter(year == 2019) |> 
  left_join(mun_data, by = "codigo_ine", suffix = c("", "_rm")) |> 
  select(-contains("rm")) |> 
  summarise("mean_by_prov" = mean(renta, na.rm = TRUE),
            .by = c("cpro", "ine_prov_name"))

summary_by_prov |> 
  slice_max(mean_by_prov, n = 1)
# A tibble: 1 × 3
  cpro  ine_prov_name mean_by_prov
  <chr> <chr>                <dbl>
1 20    Gipuzkoa            15890.

2.1.8 Pregunta 8

Obten de cada ccaa el nombre del municipio con mayor renta en 2019.

Code
renta_mun_tidy |> 
  filter(year == 2019) |> 
  left_join(mun_data, by = "codigo_ine", suffix = c("", "_rm")) |> 
  select(-contains("rm")) |> 
  slice_max(renta, by = "codauto")
# A tibble: 19 × 10
   name   codigo_ine  year renta codauto ine_ccaa_name cpro  ine_prov_name cmun 
   <chr>  <chr>      <dbl> <dbl> <chr>   <chr>         <chr> <chr>         <chr>
 1 Navas… 40904       2019 20242 07      Castilla y L… 40    Segovia       904  
 2 Almar… 10019       2019 13940 11      Extremadura   10    Cáceres       019  
 3 Oleir… 15058       2019 16447 12      Galicia       15    Coruña, A     058  
 4 Laukiz 48053       2019 21672 16      País Vasco    48    Bizkaia       053  
 5 Lumbr… 26091       2019 24007 17      Rioja, La     26    Rioja, La     091  
 6 Murcia 30030       2019 11631 14      Murcia, Regi… 30    Murcia        030  
 7 Rajad… 08178       2019 23401 09      Cataluña      08    Barcelona     178  
 8 Alocén 19023       2019 19700 08      Castilla - L… 19    Guadalajara   023  
 9 Arguis 22037       2019 19219 02      Aragón        22    Huesca        037  
10 Alcud… 04009       2019 14196 01      Andalucía     04    Almería       009  
11 Arang… 31023       2019 15517 15      Navarra, Com… 31    Navarra       023  
12 Santa… 35021       2019 15982 05      Canarias      35    Palmas, Las   021  
13 Rocaf… 46216       2019 17872 10      Comunitat Va… 46    Valencia/Val… 216  
14 Pozue… 28115       2019 26367 13      Madrid, Comu… 28    Madrid        115  
15 Valld… 07063       2019 19565 04      Balears, Ill… 07    Balears, Ill… 063  
16 Valde… 39093       2019 15574 06      Cantabria     39    Cantabria     093  
17 Tever… 33072       2019 15409 03      Asturias, Pr… 33    Asturias      072  
18 Ceuta  51001       2019 12038 18      Ceuta         51    Ceuta         001  
19 Melil… 52001       2019 11463 19      Melilla       52    Melilla       001  
# ℹ 1 more variable: lau_code <chr>

2.2 🐣 Caso práctico II: simulación

2.2.1 Pregunta 1

Define una variable llamada precio que empiece en 100. Diseña un bucle de 20 iteraciones donde, en cada iteración, el precio se reduzca a la mitad de su valor. Piensa que tipo de estructura de bucle deberías usar. El valor final de precio debería ser 0.00009536743 (aprox).

Code
# Usamos un for ya que sabemos el número de iteraciones 
# de manera prefijada (y no depende de nada)

# definimos inicialmente precio en 100
precio <- 100 

# para el bucle usamos por ejemplo i como índice, que va de 1 a 20
for (i in 1:20) {
  
  # el código fíjate que es el mismo y no depende de i
  precio <- precio/2
  
}
precio

2.2.2 Pregunta 2

Diseña una estructura de bucle de manera que encuentres la iteración en la que precio es menor que 0.001 por primera vez. Una vez encontrado guarda el número de iteración y para el bucle.

Code
# dos formas de hacerlo: for y while

# con for
precio <- 100 

# ya sabemos que en 20 es menor que 0.001 así que podemos poner
# dicha cantidad como tope sabiendo que no llegará
for (i in 1:20) {
  
  # si todavía no es menor, seguimos dividiendo
  if (precio >= 0.001) { precio <- precio/2 }
  else {
    
    # si ya es menor, guardamos la iteración (piensa por qué i - 1)
    iter <- i - 1 
    
    # y paramos
    break
  }
}

# con while
precio <- 100 

iter <- 0 # debemos inicializar las iteraciones

# no sabemos cuantas iteraciones, solo que debe parar cuando
# precio esté por debajo de dicha cantidad
while (precio >= 0.001) {
  
  precio <- precio/2
    
  # estructura clásica de while: si se corre iteración
  # actualizamos un valor (en este caso que cuente una iteración)
  iter <- iter + 1
}

2.2.3 Pregunta 3

En R la función sample(x = ..., size = ...) va sernos muy útil: de una colección de elementos x, selecciona un número size al azar de ellos. Por ejemplo, si queremos simular 3 veces el lanzamiento de un dado tenemos 6 elementos posibles (x = 1:6) y lo seleccionamos 3 veces (size = 3)

sample(x = 1:6, size = 3)
[1] 1 6 4

Al ser aleatorio, cada vez que lo ejecutas saldrá algo distinto

sample(x = 1:6, size = 3)
[1] 1 5 4

¿Y si queremos tirar 10 veces? Al tener solo 6 elementos posibles y elegir 10, no puede, así que le tenemos que indicar que queremos un sample (muestreo) con reemplazamiento (como sucede en el dado, cada cara puede repetirse al volver a tirarlo)

sample(x = 1:6, size = 10, replace = TRUE)
 [1] 3 6 3 2 1 6 4 6 6 3

Con lo anterior, imagina que estás en un concurso de televisión donde te dan a elegir 3 puertas: en una hay un premio millonario y en las otras 2 una galleta oreo. Diseña el estudio de simulación con bucles for para aproximar la probabilidad de que te toque el premio (obviamente tiene que darte aprox 0.3333333). Realiza el experimento para 1000 intentos. ¿Qué observas?

Code
# Definimos las puertas posibles
puertas <- c(1, 2, 3)

# Definimos los intentos
intentos <- 1000

# Definimos las veces que hemos ganado (al inicio empieza en 0 claro)
exitos <- 0

# Simulación de una eleccion de puerta y un premio
for (i in 1:intentos) {
  
  # premio: de 3 puertas, solo está en una
  premio <- sample(x = puertas, size = 1)
    
  # puerta que seleccionas como concursante: de 3 puertas, te quedas con una
  eleccion <- sample(x = puertas, size = 1)
    
  # si la puerta seleccionada coincide con la que tiene el premio
  # sumas un éxito, sino te quedas como estás
  exitos<- dplyr::if_else(eleccion == premio, exitos + 1, exitos)
    
}
# Tras jugar, lo dividimos entre el número de veces que has jugado
exitos <- exitos / intentos
exitos

2.2.4 Pregunta 4

Realiza el experimento para 10, 50 intentos, 100 intentos, 500 intentos, 1000 intentos, 10 000 intentos y 25 000 intentos (pista: necesitas un bucle dentro de otro). ¿Qué observas?

Code
# Definimos las puertas posibles
puertas <- c(1, 2, 3)

# Definimos los intentos
intentos <- c(10, 50, 100, 500, 1000, 10000, 25000)

# Definimos las veces que hemos ganado (al inicio empieza en 0 claro)
exitos <- rep(0, length(intentos))

# Simulación de una eleccion de puerta y un premio
for (j in 1:length(intentos)) {
  for (i in 1:intentos[j]) {
    
    # premio: de 3 puertas, solo está en una
    premio <- sample(x = puertas, size = 1)
      
    # puerta que seleccionas como concursante: de 3 puertas, te quedas con una
    eleccion <- sample(x = puertas, size = 1)
      
    # si la puerta seleccionada coincide con la que tiene el premio
    # sumas un éxito, sino te quedas como estás
    exitos[j] <- dplyr::if_else(eleccion == premio, exitos[j] + 1, exitos[j])
      
  }
  # Tras jugar, lo dividimos entre el número de veces que has jugado
  exitos[j] <- exitos[j] / intentos[j]
}
exitos

2.2.5 Pregunta 5

¿Y si en cada ronda el presentador (tras tu elegir una puerta) te abriese una de las puertas no premiadas y que no has elegido: cambiarías de puerta o te mantendrías? Simula ambos casos y descubre cuál es la estrategia correcta (este problema se conoce como problema de Monty Hall y aparece incluso en películas como 21 Black Jack)

Code
# Definimos las puertas posibles
puertas <- c(1, 2, 3)

# Definimos los intentos
intentos <- c(10, 50, 100, 500, 1000, 10000, 25000)

# Definimos las veces que hemos ganado (al inicio empieza en 0 claro)
# ahora dos estrategias: mantener o cambiar
exitos_mantengo <- exitos_cambio <- rep(0, length(intentos))

# Simulación de una eleccion de puerta y un premio
for (j in 1:length(intentos)) {
  for (i in 1:intentos[j]) {
    
    # premio: de 3 puertas, solo está en una
    premio <- sample(x = puertas, size = 1)
      
    # puerta que seleccionas inicialmente
    eleccion_inicial <- sample(x = puertas, size = 1)
    
    # puerta que te podría abrir el presentador
    puertas_posibles <- puertas[!(puertas %in% c(eleccion_inicial, premio))]
    
    # si solo te podría abrir una (porque tu elección y el premio son distintas)
    # te abre esa directamente. si te puede abrir dos (porque tu elección
    # y el premio son la misma) decide una de las dos aleatorias
    if (length(puertas_posibles) == 1) {
      puerta_abierta <- puertas_posibles
      
    } else {
      puerta_abierta <- sample(x = puertas_posibles, size = 1)
    }
    
      
    # si mantienes: solo comprobar que eleccion_inicial == premio
    exitos_mantengo[j] <- 
      dplyr::if_else(eleccion_inicial == premio, exitos_mantengo[j] + 1, exitos_mantengo[j])
    
    # si cambias: tu elegiste una, el presentador te abrio otra, solo 
    # puedes cambiar a la que no conoces
    eleccion_cambio <- puertas[!(puertas %in% c(eleccion_inicial, puerta_abierta))]
    exitos_cambio[j] <-
      dplyr::if_else(eleccion_cambio == premio, exitos_cambio[j] + 1, exitos_cambio[j])
      
  }
  # Tras jugar, lo dividimos entre el número de veces que has jugado
  exitos_cambio[j] <- exitos_cambio[j] / intentos[j]
  exitos_mantengo[j] <- exitos_mantengo[j] / intentos[j]
}
exitos_cambio
exitos_mantengo

2.2.6 Pregunta 6

Diseña una función propia con 3 argumentos: puertas, puerta_elegida y cambio. La primera podrá tomar los valores c(1, 2, 3, ..., n) (simulando un concurso de n puertas, n al menos 3); la segunda un solo número de puerta elegida entre las posibles; la tercera TRUE/FALSE. La función debe replicar una simulación de la pregunta 5 y devolver TRUE (el concursante ganó el premio) o FALSE. Comprueba dentro de la función que los argumentos son correctos; en caso contrario devuelve un mensaje de error. Tras diseñar la función aplícala y revisa que funciona cuando toca y que devuelve error cuando toca.

Code
# defino la función con unos argumentos por defecto (los que tomará
# si al llamarle no se le asigna nada)
concurso <- 
  function(puertas = 1:3, puerta_elegida = 1, cambio = TRUE) {
  
    # primer check: ¿puertas_posibles es un vector numérico de al menos 3 puertas?
    if (length(puertas) < 3 | !is.numeric(puertas)) {
      stop("puertas debe ser un vector numérico de al menos longitud 3")
    }
    
    # segundo check: si puerta elegida está dentro de las posibles
    if (!(puerta_elegida %in% puertas)) {
      stop("puerta_elegida debe ser una opción dentro de puertas")
    }
    
    # tercer check: si cambio es una variable lógica sin ausentes
    if (is.na(cambio) | !is.logical(cambio)) {
      stop("cambio debe ser una variable binaria sin ausentes")
    }
    
    # premio: de 3 puertas, solo está en una
    premio <- sample(x = puertas, size = 1)
      
    # puerta que seleccionas inicialmente
    eleccion_inicial <- puerta_elegida
    
    # puerta que te podría abrir el presentador
    puertas_posibles <- puertas[!(puertas %in% c(eleccion_inicial, premio))]
    
    if (length(puertas_posibles) == 1) {
      puerta_abierta <- puertas_posibles
    } else {
      puerta_abierta <- sample(x = puertas_posibles, size = 1)
    }
    
    eleccion_cambio <- puertas[!(puertas %in% c(eleccion_inicial, puerta_abierta))]
    
    if (length(eleccion_cambio) > 1) {
      eleccion_cambio <- sample(x = eleccion_cambio, size = 1)
    }
    exito <- dplyr::if_else(cambio, eleccion_cambio == premio, eleccion_inicial == premio)
  
    return(exito)
  }

# ejemplos que deberían funcionar
concurso()
concurso(cambio = FALSE)
concurso(puertas = 1:7)
concurso(puertas = 1:7, puerta_elegida = 4, cambio = FALSE)

# ejemplos que no
concurso(puertas = 1:2)
concurso(puertas = c("a", "b", "c"))
concurso(puertas = 1:5, puerta_elegida = 7)
concurso(puertas = 1:5, puerta_elegida = 3, cambio = "si")

3 Clase 3

3.1 🐣 Caso práctico I: brexit

Veamos un estudio de caso real basado en los artículos «The Brexit Vote: A Divided Nation, a Divided Continent» (S. B. Hobolt, 2016) y «Who Voted for Brexit? A Comprehensive District-Level Analysis» (S. O. Becker et al., 2017). Los datos se extraerán del repositorio de Github de Elena Llaudet

3.1.1 Pregunta 1

Carga del repositorio el .csv

Code
library(readr)
brexit_data <- read_csv(file = "https://raw.githubusercontent.com/ellaudet/DSS/refs/heads/master/BES.csv")
Rows: 30895 Columns: 4
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (1): vote
dbl (3): leave, education, age

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

3.1.2 Pregunta 2

Escribe el código que consideres para responder a las siguientes preguntas

  • ¿Cuál es el tamaño muestral de la encuesta?
  • ¿Cuántas variables se han recopilado?
  • ¿Cuántos ausentes tiene cada variable?
Code
# tamaño muestral
brexit_data |> 
  nrow()

# variables
brexit_data |> 
  ncol()

# ausentes
brexit_data |> 
  summarise(across(everything(), function(x) { sum(is.na(x)) }))

3.1.3 Pregunta 3

Si te fijas hay 2851 ausentes en leave (variable binaria). Comprueba que todos ellos se deben a registros en los que la variable vote (variable cualitativa) toma valores don't known o won't vote. Hazlo tanto en modo tidyverse (count) como en R base

Code
brexit_data |>
  count(leave, vote)

table(brexit_data$leave, brexit_data$vote)

3.1.4 Pregunta 4

Debajo tienes el código de como se construiría una tabla de frecuencias bidimensional en R base

table(brexit_data$vote, brexit_data$education)
            
                1    2    3    4    5
  don't know  157  512  453  696  144
  leave      1356 3388 2685 3783  631
  stay        498 1763 3014 6081 1898
  won't vote   34  118  120  116   23

Diseña lo mismo pero haciendo uso de tidyverse (pista: necesitas pivotar)

Code
brexit_data |> 
  count(vote, education) |> 
  pivot_wider(names_from = "education", values_from = "n")

Moraleja: aunque casi siempre será más legible y eficiente (y coherente) el código en tidyverse, a veces nos tocará echar mano de R base. Por eso es importante distinguir ambas formas de trabajar.

Te dejo por cierto debajo cómo se calcularían las frecuencias relativas (margin = 1 normaliza por filas y margin = 2 por columnas)

freq_abs <- table(brexit_data$vote, brexit_data$education)
prop.table(freq_abs)
            
                       1           2           3           4           5
  don't know 0.005715326 0.018638515 0.016490717 0.025336731 0.005242082
  leave      0.049362941 0.123334547 0.097742992 0.137713870 0.022970513
  stay       0.018128868 0.064179104 0.109719694 0.221368766 0.069093557
  won't vote 0.001237714 0.004295595 0.004368402 0.004222788 0.000837277
prop.table(freq_abs, margin = 1)
            
                      1          2          3          4          5
  don't know 0.08002039 0.26095821 0.23088685 0.35474006 0.07339450
  leave      0.11449802 0.28607616 0.22671620 0.31942920 0.05328042
  stay       0.03757356 0.13301645 0.22740305 0.45880489 0.14320205
  won't vote 0.08272506 0.28710462 0.29197080 0.28223844 0.05596107
prop.table(freq_abs, margin = 2)
            
                       1           2           3           4           5
  don't know 0.076772616 0.088565992 0.072225765 0.065192956 0.053412463
  leave      0.663080685 0.586057775 0.428093112 0.354346197 0.234050445
  stay       0.243520782 0.304964539 0.480548469 0.569595354 0.704005935
  won't vote 0.016625917 0.020411693 0.019132653 0.010865493 0.008531157

3.1.5 Pregunta 5

Usando las tablas relativas anteriores responder a

  • ¿Cuántas personas de los que votaron leave tiene un alto nivel educativo?

  • De la gente con un bajo nivel educativo, ¿qué % voto leave?

  • ¿Cuál fue el nivel educativo es el menos propenso a votar?

3.1.6 Pregunta 6

¿Qué porcentaje estimó la encuesta que votó por permanecer en la Unión Europea (el resultado real después del referéndum fue del 51.89%)?

Code
prop.table(table(brexit_data$vote))
# 46.45 % sobre el total de encuestas

prop.table(table(brexit_data$vote[brexit_data$vote != "won't vote"]))
# 47.27 si obviamos los que decían no ir a votar

3.1.7 Pregunta 7

Calcula la media de edad para cada opción de la variable vote.

Code
brexit_data |> 
  summarise("avg_age" = mean(age), .by = vote)

3.2 🐣 Caso práctico II: encuesta de satisfacción

Vamos a seguir poniendo en práctica lo aprendido el dataset SatisfaccionPacientes.csv

library(readr)
datos <-
  read_csv(file = "./datos/SatisfaccionPacientes.csv") |> 
  janitor::clean_names() |>
  mutate(estado_civil = factor(estado_civil),
         genero = factor(genero))
Rows: 100 Columns: 8
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (3): Genero, EstadoCivil, EstadoSalud
dbl (5): ID, Edad, TiempoEspera, GradoSatisfaccion, NumeroVisitas

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
datos
# A tibble: 100 × 8
      id  edad genero    estado_civil tiempo_espera grado_satisfaccion
   <dbl> <dbl> <fct>     <fct>                <dbl>              <dbl>
 1     1    60 Masculino Casado                  28                  8
 2     2    44 Femenino  Soltero                 22                  8
 3     3    43 Masculino Soltero                  8                  9
 4     4    32 Masculino Soltero                 21                  8
 5     5    66 Masculino Divorciado               7                 10
 6     6    43 Masculino Divorciado              20                  8
 7     7    54 Masculino Casado                  18                  6
 8     8    55 Masculino Soltero                 29                  6
 9     9    56 Masculino Viudo                   17                  9
10    10    34 Femenino  Casado                  34                  8
# ℹ 90 more rows
# ℹ 2 more variables: numero_visitas <dbl>, estado_salud <chr>

3.2.1 Pregunta 1

Convierte de manera adecuada la variable estado_salud a cualitativa ORDINAL

Code
datos <-
  datos |>
  mutate(estado_salud =
           factor(estado_salud, levels = c("Malo", "Regular", "Bueno", "Excelente"),
                  ordered = TRUE))

3.2.2 Pregunta 2

Haz uso de table() para calcular la tabla de frecuencias de genero y estado_civil

Code
table(datos$genero)
table(datos$estado_civil)

3.2.3 Pregunta 3

Calcula la tabla de frecuencias de las ORDINALES y piensa si ahora puedes añadir algo más a la tabla de frecuencias). Tras ello usa el código más sencillo para responder a: ¿cuántas personas tienen un estado de salud regular (o peor)?

Code
freq_estado_salud <-
  datos |> 
  count(estado_salud) |> 
  rename(frecuencia_abs = n) |> 
  mutate(frecuencia_rel = frecuencia_abs/sum(frecuencia_abs),
         frecuencia_acum_abs = cumsum(frecuencia_abs),
         frecuencia_acum_rel = cumsum(frecuencia_rel))
# Se ve dentro de la tabla. Hay 44+15 = 59 personas con un estado de salud malo o regular. 
# Con código
datos |>
  count(estado_salud <= "Regular")

3.2.4 Pregunta 4

Si te fijas una de las modalidades es totalmente anecdótica (solo 1 Excelente). Sería conveniente recategorizar la categoría Excelente: siempre que detecte Excelente, lo debe recategorizar a Bueno (criterio general: las categorías deben contener al menos un 5% de los individuos de toda la muestra).

Code
datos <- 
  datos |> 
  mutate(estado_salud  = if_else(estado_salud  == "Excelente", "Bueno", estado_salud),
         # ojo: hay que redefinir los niveles de la cualitativa
         # ya que ha dejado de ser factor (veremos un día el paquete forcats para esto)
         estado_salud =
           factor(estado_salud, levels = c("Malo", "Regular", "Bueno"),
                  ordered = TRUE))

3.2.5 Pregunta 5

Calcula la media, mediana, rango intercuartílico y desviación típica de grado de satisfacción desagregado por sexo.

Code
resumen <-
  datos |>
  summarise(media_grado_satisfaccion = mean(grado_satisfaccion),
           sd_grado_satisfaccion = sd(grado_satisfaccion),
           mediana_grado_satisfaccion = median(grado_satisfaccion),
           IQR_grado_satisfaccion =
             quantile(grado_satisfaccion, probs = 0.75) -
             quantile(grado_satisfaccion, probs = 0.25),
           .by = genero)

3.2.6 Pregunta 6

Exporta los resultados anteriores (resumen) en un archivo resumen.csv. En lugar de un read_csv() vamos a usar write_csv(tabla, file = "ruta")

Code
# importado como csv
write_csv(resumen, file = "./resumen.csv")

3.2.7 Pregunta 7

Crear un diagrama de barras para la variable género. ¿Cómo podríamos decirle que cada barra (es decir, para cada modalidad de género) sea de un color (de relleno)?

Code
library(ggplot2)
ggplot(datos) +
  # dentro de aes() para que dependa de la tabla
  geom_bar(aes(x = genero, fill = genero))

3.2.8 Pregunta 8

Crear desde cero un diagrama de barras, con ajustes personalizados para la variable estado de salud

Code
# Estado de salud (ahora el orden importa)
ggplot(datos) +
  geom_bar(aes(x = estado_salud, fill = estado_salud), alpha = 0.75) +
  ggthemes::scale_fill_colorblind() +
  labs(title = "Diagrama de barras de la variable estado salud", 
       x = "Categoría",  y = "Frecuencia absoluta",
       fill = "Categoría") +
  theme_minimal() 

Fíjate que si no tuviésemos la variable como cuali ordinal, las barras van por orden alfabético, no por jerarquía real

ggplot(datos) +
  geom_bar(aes(x = as.character(estado_salud), fill = as.character(estado_salud)),
           alpha = 0.75) +
  ggthemes::scale_fill_colorblind() +
  labs(title = "Diagrama de barras de la variable estado salud", 
       x = "Categoría",  y = "Frecuencia absoluta",
       fill = "Categoría") +
  theme_minimal() 

3.2.9 Pregunta 9

Crea el histograma inferior para las variable edad y tiempo de espera.

Code
ggplot(datos) +
  geom_histogram(aes(x = edad), bins = 30, fill = "darkorange", alpha = 0.75) + 
  labs(title = "Histograma de edad", subtitle = "Bins = 30",
       x = "Valores", y = "Frecuencia absoluta") +
  theme_minimal()

Code
ggplot(datos) +
  # Define el ancho de las barras y colores
  geom_histogram(aes(x = tiempo_espera), bins = 30, fill = "orchid", alpha = 0.75) + 
  labs(title = "Histograma de tiempo de espera", subtitle = "Bins = 30",
       x = "Valores", y = "Frecuencia absoluta") +
  theme_minimal()

3.2.10 Pregunta 10

Crea el gráfico de densidad inferior para las variable edad y tiempo de espera.

Code
ggplot(datos) +
  geom_density(aes(x = edad), color = "darkorange", 
               fill = "darkorange", alpha = 0.75) + 
  labs(title = "Gráfico de densidad de edad",
       x = "Valores", y = "Frecuencia relativa") +
  theme_minimal()

Code
ggplot(datos) +
  geom_density(aes(x = tiempo_espera), color = "orchid", 
               fill = "orchid", alpha = 0.75) + 
  labs(title = "Gráfico de densidad de tiempo de espera",
       x = "Valores", y = "Frecuencia relativa") +
  theme_minimal()

3.2.11 Pregunta 11

Realiza un boxplot para edad y un boxplot para numero de visitas PERO por género (dos variables, piensa cómo)

Code
ggplot(datos) +
  geom_boxplot(aes(y = edad), fill = "lightblue", alpha = 0.75) +  
  labs(title = "Boxplot de edad",  y = "Edad") +
  theme_minimal()

Code
ggplot(datos) +
  geom_boxplot(aes(x = genero, y = tiempo_espera, fill = genero),
               alpha = 0.75) +
  labs(title = "Boxplot de tiempo de espera por género", 
       x = "Género", y = "Tiempo de Espera") +
  theme_minimal()

Haciendo uso del gráfico anterior:

  1. ¿La variable edad tiene outliers? ¿Qué edad tienen esos pacientes?

  2. ¿Quién ha esperado más los hombres o las mujeres?

3.3 🐣 Caso práctico III: bronquitis y tabaco

Vamos a cargar el archivo de datos fumadores.csv donde tenemos datos de 96 pacientes sobre sí o fuman y quienes han desarrollado o no bronquitis.

datos <- read_csv(file = "./datos/fumadores.csv")
Rows: 96 Columns: 3
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (2): fumador, estado
dbl (1): id

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
datos
# A tibble: 96 × 3
      id fumador estado    
   <dbl> <chr>   <chr>     
 1     1 sí      bronquitis
 2     2 sí      bronquitis
 3     3 sí      bronquitis
 4     4 sí      bronquitis
 5     5 sí      bronquitis
 6     6 sí      bronquitis
 7     7 sí      bronquitis
 8     8 sí      bronquitis
 9     9 sí      bronquitis
10    10 sí      bronquitis
# ℹ 86 more rows

3.3.1 Pregunta 1

Realiza la tabla de contigencia de manera absoluta y relativa y responde a las siguientes preguntas

  1. ¿Cuántas personas fumaoras tienen bronquitis?

  2. ¿Qué % de los fumadores está sano?

  3. ¿Qué % del total son a la vez no fumadores y enfermos de bronquitis?

  4. ¿Qué % de los enfermos son fumadores?

Code
table(datos$fumador, datos$estado)
prop.table(table(datos$fumador, datos$estado))
prop.table(table(datos$fumador, datos$estado), margin = 1)
prop.table(table(datos$fumador, datos$estado), margin = 2)
# a) 32 personas
# b) 38.46%
# c) 16%
# d) 61.53%

3.3.2 Pregunta 2

Visualiza ambas variables (fumador y estado) a la vez de manera adecuada que nos permita comparar

Code
ggplot(datos) +
  geom_bar(aes(x = fumador, fill = estado), alpha = 0.6, position = "fill") +
  labs(x = "fumador", y = "Frec relativa", fill = "Estado") +
  theme_minimal()

3.3.3 Pregunta 3

¿Existen evidencias en la muestra de una asociación entre ambas variables?

Code
datos |> 
  summarise("sig_chisq" = chisq.test(datos$fumador, datos$estado)$p.value,
            "sig_fisher" = fisher.test(datos$fumador, datos$estado)$p.value)
# p-valor < alpha --> hay evidencias para rechazar la hip nula
# hay evidencias (no muy fuertes, quizás aumentar tamaño muestral?) de
# que las variables son dependientes y existe una asociación

3.3.4 Pregunta 4

Si hubiera asociación, cuantifica la fuerza de dicha asociación (y el sentido) y calcula el riesgo relativo de los fumadores a contraer bronquitis (respecto a no fumadores)

Code
fisher.test(datos$fumador, datos$estado)
# OR estimado = 0.3611 ==> piensa que tenemos no fumar primero y luego fumar ==>
# 1/0.3611 = 2.769316 > 1 ==> hay una asociación positiva entre fumar y tener 
# bronquitis 
# La bronquitis en pacientes que fuman es 2.77 veces más frecuente
# que en los pacientes que no fuman

# RR ratio
a <- 32  # Expuestos con evento
b <- 16  # Expuestos sin evento
c <- 20  # No expuestos con evento
d <- 28  # No expuestos sin evento

RR <- (a / (a + b)) / (c / (c + d))
# El grupo que fuma tiene un riesgo 1.6 veces mayor de que desarrollar bronquitis en comparación con el grupo que no fuma.

3.4 🐣 Caso práctico IV: salud mental

Esta la base de datos datos_salud_mental.csv tenemos información recopilada de 100 pacientes que acuden a un centro de salud mental. Se quiere realizar un estudio para ver el impacto que tienen distintas características sobre la ansiedad y depresión en estos 100 pacientes. Los datos incluyen una variedad de variables relacionadas con la salud mental, así como características demográficas y de estilo de vida.

datos <-
  read_csv(file = "./datos/datos_salud_mental.csv") |> 
  janitor::clean_names()
Rows: 100 Columns: 10
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (3): Genero, UsoDrogasRecreativas, TipoDrogas
dbl (7): ID, Edad, Ansiedad, Depresion, SesionesTerapia, ActividadFisica, Ho...

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
datos
# A tibble: 100 × 10
      id  edad genero    ansiedad depresion sesiones_terapia actividad_fisica
   <dbl> <dbl> <chr>        <dbl>     <dbl>            <dbl>            <dbl>
 1     1    48 Masculino        5         5               60                6
 2     2    32 Masculino        5         6               30               NA
 3     3    68 Masculino        9        NA               32                7
 4     4    31 NoBinario        9        NA              100                5
 5     5    20 Femenino         9        NA               27                6
 6     6    59 Masculino        7         5               38                7
 7     7    67 Femenino        10        NA               49                7
 8     8    NA Femenino         6         6               55                7
 9     9    54 Femenino         6         6               55                6
10    10    69 Masculino        4         6               32                7
# ℹ 90 more rows
# ℹ 3 more variables: horas_sueno <dbl>, uso_drogas_recreativas <chr>,
#   tipo_drogas <chr>

Las variables son:

  • id: identificador único del paciente.
  • edad: edad del paciente en años.
  • Genero: género del paciente.
  • ansiedad: nivel de ansiedad del paciente en una escala del 1 al 10.
  • depresión: nivel de depresión del paciente en una escala del 1 al 10.
  • sesiones_terapia: número de sesiones de terapia asistidas en el último año.
  • actividad_fisica: número de días a la semana que el paciente realiza actividad física.
  • horas_sueno: número de horas promedio de sueño por noche.
  • uso_drogas_recreativas: indicador de si el paciente ha usado drogas recreativas en el último año.
  • tipo_drogas: tipo de drogas que ha consumido el paciente.

3.4.1 Pregunta 1

¿De qué tipo es cada variable? Convierte las que consideres a cualis nominales y a cualis ordinales, y si hay alguna variable que deba ser lógica

Code
# id: en realidad esto tendría ser un factor (un texto) ya que no cuenta nada
# Cuantitativas: edad, horas_sueno
# Cualitativas nominales: genero, tipo_drogas
# Cuanitativas discretas: sesiones_terapia y actividad_fisica
# Cuantitativas discretas pero que deberíamos tratarlas como cualis ordinales
# ya que son escalas: ansiedad, depresión
# Binarias (cualis ordinales muy concretas): uso_drogas_recreativas
datos <-
  datos |> 
  mutate("id" = as.character(id),
         "genero" = factor(genero), "tipo_drogas" = factor(tipo_drogas),
         "ansiedad" = factor(ansiedad, levels = as.character(1:10), ordered = TRUE),
         "depresion" = factor(depresion, levels = as.character(1:10), ordered = TRUE),
         "uso_drogas_recreativas" = (uso_drogas_recreativas == "Si"))

3.4.2 Pregunta 2

Calcula la tabla de frecuencias absolutas y relativas de género.

Code
tabla_freq_abs <- table(datos$genero)
tabla_freq_rel <- prop.table(tabla_freq_abs)

3.4.3 Pregunta 3

Calcula la media de las 4 variables cuantitativas que tenemos desagregado por género. Exporta dicho resumen en un .csv

Code
resumen <- 
  datos |> 
  drop_na(where(is.numeric)) |> 
  summarise(across(where(is.numeric), mean), .by = genero)
write_csv(resumen, file = "./datos/resumen.csv")

3.4.4 Pregunta 4

Calcula la tabla de contigencia de las variables ansiedad vs genero. Calcula otra para ansiedad vs depresion. Usa useNA = "always" como argumento para incluir los NA

Code
tabla_freq_genero_ansiedad <- table(datos$genero, datos$ansiedad, useNA = "always")
tabla_freq_depresion_ansiedad <- table(datos$depresion, datos$ansiedad, useNA = "always")

3.4.5 Pregunta 5

Realiza un gráfico adecuado para la variable edad. Piensa como adaptarlo para tenerlo desagregado por genero.

Code
# Densidades
ggplot(datos) +
  geom_density(aes(x = edad), fill = "#459191", alpha = 0.4) +
  labs(x = "Edad (años)", y = "Densidad (frec relativa)") +
  theme_minimal()

library(ggridges)
ggplot(datos) +
  geom_density_ridges(aes(x = edad, y = genero, fill = genero), alpha = 0.4) +
  ggthemes::scale_fill_colorblind() +
  labs(x = "Edad (años)", y = "Sexo", fill = "Género") +
  theme_minimal()

# Histograma (mala idea con pocos datos)
ggplot(datos) +
  geom_histogram(aes(x = edad), bins = 15, fill = "#459191", alpha = 0.4) +
  labs(x = "Edad (años)", y = "Frec absoluta") +
  theme_minimal()
ggplot(datos) +
    geom_histogram(aes(x = edad, fill = genero), bins = 15, alpha = 0.25) +
    labs(x = "Edad (años)", y = "Frec absoluta") +
    theme_minimal()

# Boxplot
ggplot(datos) +
  geom_boxplot(aes(y = edad), fill = "#459191", alpha = 0.4,
               outlier.size = 3, outlier.alpha = 0.9,
               outlier.color = "#991293", outlier.shape = 18) +
  labs(y = "Edad") +
  theme_minimal()

ggplot(datos, aes(x = genero, y = edad, fill = genero, color = genero)) +
  geom_boxplot(alpha = 0.4, outlier.size = 3, outlier.alpha = 0.9,
               outlier.color = "#991293", outlier.shape = 18) +
  ggthemes::scale_color_colorblind() +
  ggthemes::scale_fill_colorblind() +
  guides(color = "none") +
  labs(x = "Género", y = "Edad", fill = "Género") +
  theme_minimal()

3.4.6 Pregunta 6

Realiza un gráfico adecuado para visualizar a la vez depresión y ansiedad.

Code
# fíjate que aunque sean números, dado que son variables discretas
# de una escala, no permite una correcta visualización un diagrama
# de dispersión ya que hay muchos puntos iguales que se solapan
ggplot(datos) + 
  geom_point(aes(x = depresion, y = ansiedad)) +
  theme_minimal()

# una opción: se ve un patrón (tipo "recta ascedente")
ggplot(datos |> count(depresion, ansiedad)) + 
  geom_point(aes(x = depresion, y = ansiedad, color = n, size = n)) +
  scale_color_viridis_c() +
  guides(size = "none") +
  theme_minimal()

3.4.7 Pregunta 7

¿Existe asociación entre genero y uso de drogas? ¿Y entre depresión y ansiedad? Cuantifica la respuesta todo lo que puedas.

Code
resumen_p_valores_1 <-
  datos |> 
  summarise("sig_chisq" = chisq.test(datos$genero, datos$uso_drogas_recreativas)$p.value,
            "sig_fisher" = fisher.test(datos$genero, datos$uso_drogas_recreativas)$p.value)
Warning: There was 1 warning in `summarise()`.
ℹ In argument: `sig_chisq = chisq.test(datos$genero,
  datos$uso_drogas_recreativas)$p.value`.
Caused by warning in `chisq.test()`:
! Chi-squared approximation may be incorrect
Code
# p-valores >> alpha --> no evidencias para rechazar la independencia -->
# no hay evidencias para afirmar la dependencia

# con tantas categorías Fisher no funciona
chisq.test(datos$depresion, datos$ansiedad)
Warning in chisq.test(datos$depresion, datos$ansiedad): Chi-squared
approximation may be incorrect

    Pearson's Chi-squared test

data:  datos$depresion and datos$ansiedad
X-squared = 174.39, df = 63, p-value = 2.148e-12
Code
# p-valores << alpha --> sí hay evidencias para rechazar la independencia -->
# sí hay evidencias para afirmar la dependencia