---
title: "Calidad del Aire UCA — Exploración de Datos"
output:
  html_notebook:
    toc: true
    toc_float: true
    theme: flatly
---

Notebook de exploración para el sitio de monitoreo de calidad de aire
(`calidadaire.uca.edu.ar`). Fuente: mediciones del sensor exportadas a CSV.


**Requisitos**

```r
install.packages(c("dplyr", "tidyr", "lubridate", "plotly"))
```

## 1. Configuración e imports

```{r setup, message=FALSE, warning=FALSE}
library(dplyr)
library(tidyr)
library(lubridate)
library(plotly)

options(width = 120)

# --- Funciones auxiliares ---------------------------------------------------
# pandas no tiene equivalente directo en base R, así que se replican a mano.

# Equivalente a pandas .rolling(window, min_periods).mean() (ventana trailing)
rolling_mean_min <- function(x, window = 10, min_periods = 2) {
  n <- length(x)
  out <- rep(NA_real_, n)
  for (i in seq_len(n)) {
    ini <- max(1, i - window + 1)
    vals <- x[ini:i]
    vals <- vals[!is.na(vals)]
    if (length(vals) >= min_periods) out[i] <- mean(vals)
  }
  out
}

# Interpolación lineal por posición, análoga a pandas .interpolate()
interpolar_lineal <- function(x) {
  idx <- seq_along(x)
  if (sum(!is.na(x)) < 2) return(x)
  approx(idx[!is.na(x)], x[!is.na(x)], xout = idx, rule = 2)$y
}
```

## 2. Carga de datos

Se busca la carpeta `data/` subiendo por el árbol de directorios desde el
directorio actual, para que la notebook funcione sin importar desde dónde
se la ejecute dentro del repo.

```{r cargar-datos}
file_path <- "./datos_calidad_aire_uca_2026-04-01_2026-07-01.csv"

cat("Archivo de datos  :", file_path, "\n")
```

```{r leer-csv}
df <- read.csv(file_path, stringsAsFactors = FALSE)
df$fecha_reporte <- lubridate::as_datetime(df$fecha_reporte)

head(df)
```

## 3. Preparación de variables

```{r feature-engineering}
df$hora <- lubridate::hour(df$fecha_reporte)
df$anio <- lubridate::year(df$fecha_reporte)

dias_semana_es <- c("Lunes", "Martes", "Miércoles", "Jueves", "Viernes", "Sábado", "Domingo")
# lubridate::wday(..., week_start = 1) -> 1 = Lunes ... 7 = Domingo
df$dia_semana <- dias_semana_es[lubridate::wday(df$fecha_reporte, week_start = 1)]


clasificar_franja <- function(hora) {
  dplyr::case_when(
    hora >= 0  & hora < 6  ~ "Madrugada",
    hora >= 6  & hora < 12 ~ "Mañana",
    hora >= 12 & hora < 14 ~ "Mediodía",
    hora >= 14 & hora < 20 ~ "Tarde",
    TRUE                   ~ "Noche"
  )
}

df$franja_horaria <- clasificar_franja(df$hora)

head(df)
```

## 4. Vista general del dataset

```{r shape-info}
cat(sprintf("Filas x columnas: %d x %d\n", nrow(df), ncol(df)))
str(df)
```

```{r faltantes}
faltantes <- colSums(is.na(df))
faltantes <- sort(faltantes[faltantes > 0], decreasing = TRUE)

if (length(faltantes) == 0) {
  cat("No hay valores faltantes.\n")
} else {
  data.frame(nulos = faltantes)
}
```

```{r describe}
cols_numericas <- names(df)[sapply(df, is.numeric)]
cols_numericas <- setdiff(cols_numericas, c("id", "hora", "anio"))

describe_df <- function(df, cols) {
  do.call(rbind, lapply(cols, function(col) {
    x <- df[[col]]
    data.frame(
      variable = col,
      count = sum(!is.na(x)),
      mean = mean(x, na.rm = TRUE),
      sd = sd(x, na.rm = TRUE),
      min = min(x, na.rm = TRUE),
      p25 = quantile(x, 0.25, na.rm = TRUE),
      p50 = quantile(x, 0.5, na.rm = TRUE),
      p75 = quantile(x, 0.75, na.rm = TRUE),
      max = max(x, na.rm = TRUE),
      row.names = NULL
    )
  }))
}

describe_df(df, cols_numericas)
```

```{r rango-temporal}
rango_inicio <- min(df$fecha_reporte)
rango_fin <- max(df$fecha_reporte)

cat(sprintf("Rango temporal: %s → %s\n", rango_inicio, rango_fin))
cat("Duración:", format(rango_fin - rango_inicio), "\n")
```

## 5. Réplica de los gráficos del sitio

### 5.1 Serie temporal de AQI (últimas 24 hs)

```{r grafico-serie-temporal}
grafico_aqi_serie_temporal <- function(df) {
  g <- df %>% dplyr::arrange(fecha_reporte) %>% utils::tail(24 * 12)  # últimas 24hs = 288 lecturas de 5min
  g$aqi_smooth <- rolling_mean_min(g$aqi, window = 10, min_periods = 2)

  meses <- c("ene", "feb", "mar", "abr", "may", "jun",
             "jul", "ago", "sep", "oct", "nov", "dic")

  g$fecha_str <- sprintf(
    "%d %s %02d:%02d",
    lubridate::day(g$fecha_reporte),
    meses[lubridate::month(g$fecha_reporte)],
    lubridate::hour(g$fecha_reporte),
    lubridate::minute(g$fecha_reporte)
  )
  g$fecha_str <- factor(g$fecha_str, levels = g$fecha_str)

  ticks_idx <- seq(1, nrow(g), by = 24)

  plot_ly(g, x = ~fecha_str, y = ~aqi_smooth, type = "scatter", mode = "lines",
          line = list(color = "#378ADD", width = 2), name = "AQI") %>%
    layout(
      xaxis = list(
        title = "Horario",
        showgrid = FALSE,
        tickmode = "array",
        tickvals = as.character(g$fecha_str[ticks_idx]),
        ticktext = as.character(g$fecha_str[ticks_idx])
      ),
      yaxis = list(title = "AQI", showgrid = FALSE),
      plot_bgcolor = "white",
      paper_bgcolor = "white"
    )
}

grafico_aqi_serie_temporal(df)
```

### 5.2 Heatmap de AQI por día y franja horaria

```{r grafico-heatmap}
grafico_heatmap_aqi <- function(df) {
  df$hora <- lubridate::hour(df$fecha_reporte)

  dias_semana_es <- c("Lunes", "Martes", "Miércoles", "Jueves", "Viernes", "Sábado", "Domingo")
  df$dia_semana <- dias_semana_es[lubridate::wday(df$fecha_reporte, week_start = 1)]

  clasificar_franja <- function(hora) {
    dplyr::case_when(
      hora >= 0  & hora < 6  ~ "Madrugada",
      hora >= 6  & hora < 12 ~ "Mañana",
      hora >= 12 & hora < 14 ~ "Mediodía",
      hora >= 14 & hora < 20 ~ "Tarde",
      TRUE                   ~ "Noche"
    )
  }
  df$franja_horaria <- clasificar_franja(df$hora)

  orden_dias <- c("Lunes", "Martes", "Miércoles", "Jueves", "Viernes", "Sábado", "Domingo")
  orden_franjas <- c("Madrugada", "Mañana", "Mediodía", "Tarde", "Noche")

  tabla <- df %>%
    dplyr::group_by(franja_horaria, dia_semana) %>%
    dplyr::summarise(aqi = mean(aqi, na.rm = TRUE), .groups = "drop") %>%
    tidyr::pivot_wider(names_from = dia_semana, values_from = aqi)

  m <- as.matrix(tabla[, orden_dias])
  rownames(m) <- tabla$franja_horaria
  m <- m[orden_franjas, orden_dias]

  m <- apply(m, 2, interpolar_lineal)        # interpola dentro de cada columna (día)
  m <- t(apply(m, 1, interpolar_lineal))     # interpola dentro de cada fila (franja)
  dimnames(m) <- list(orden_franjas, orden_dias)

  plot_ly(
    x = orden_dias, y = orden_franjas, z = m,
    type = "heatmap",
    colorscale = list(
      list(0, "#d4fcd4"),
      list(0.33, "#66cc66"),
      list(0.66, "#006600"),
      list(1, "#ffff00")
    ),
    text = round(m, 1),
    texttemplate = "%{text}"
  ) %>%
    layout(
      xaxis = list(title = "Día"),
      yaxis = list(title = "Franja horaria"),
      plot_bgcolor = "white",
      paper_bgcolor = "white"
    )
}

grafico_heatmap_aqi(df)
```

### 5.5 Promedio por hora (función auxiliar)

```{r df-promedio-hora}
df_promedio_hora <- function(df) {
  df_hora <- df %>%
    dplyr::group_by(hora) %>%
    dplyr::summarise(dplyr::across(dplyr::where(is.numeric), ~ mean(.x, na.rm = TRUE)), .groups = "drop")

  horas <- data.frame(hora = 0:23)

  horas %>%
    dplyr::left_join(df_hora, by = "hora") %>%
    dplyr::arrange(hora) %>%
    dplyr::mutate(dplyr::across(-hora, interpolar_lineal))
}
```

## 6. Análisis exploratorio adicional

Gráficos pensados para explorar patrones que los gráficos del
sitio no muestran: perfil horario promedio, correlación entre
contaminantes, distribución del AQI por categoría y continuidad de las
mediciones en el tiempo.

### 6.1 Perfil horario promedio (24 hs)

Usa `df_promedio_hora`, la misma lógica que probablemente alimenta vistas tipo "día típico" del sitio.

```{r perfil-horario}
df_h <- df_promedio_hora(df)

colores <- c(aqi = "#378ADD", pm25 = "#2c6e49", pm10 = "#4c956c")

p <- plot_ly(df_h, x = ~hora)
for (col in c("aqi", "pm25", "pm10")) {
  if (col %in% names(df_h)) {
    p <- p %>% add_trace(
      y = df_h[[col]], name = toupper(col), type = "scatter", mode = "lines+markers",
      line = list(color = colores[[col]], width = 2)
    )
  }
}

p %>% layout(
  title = "Promedio por hora del día",
  xaxis = list(title = "Hora", dtick = 1, showgrid = FALSE),
  yaxis = list(title = "Valor promedio", showgrid = FALSE),
  plot_bgcolor = "white",
  paper_bgcolor = "white"
)
```

### 6.2 Correlación entre variables numéricas

```{r correlacion}
corr_cols <- c("aqi", "pm25", "pm10", "temperatura_c", "humedad", "punto_rocio_c", "indice_calor_c")

corr <- cor(df[, corr_cols], use = "pairwise.complete.obs")

plot_ly(
  x = corr_cols, y = corr_cols, z = corr,
  type = "heatmap",
  colorscale = "RdBu", reversescale = TRUE,
  zmin = -1, zmax = 1,
  text = round(corr, 2),
  texttemplate = "%{text}"
) %>%
  layout(title = "Matriz de correlación")
```

### 6.3 Distribución de AQI por categoría

Se usan las bandas estándar de AQI (EPA) para ver qué proporción del tiempo el aire estuvo en cada categoría.

```{r categorias-aqi}
categorizar_aqi <- function(valor) {
  dplyr::case_when(
    is.na(valor)  ~ NA_character_,
    valor <= 50   ~ "Buena",
    valor <= 100  ~ "Moderada",
    valor <= 150  ~ "Dañina (grupos sensibles)",
    valor <= 200  ~ "Dañina",
    valor <= 300  ~ "Muy dañina",
    TRUE          ~ "Peligrosa"
  )
}

orden_categorias <- c("Buena", "Moderada", "Dañina (grupos sensibles)", "Dañina", "Muy dañina", "Peligrosa")
colores_categorias <- c(
  "Buena" = "#00e400",
  "Moderada" = "#ffff00",
  "Dañina (grupos sensibles)" = "#ff7e00",
  "Dañina" = "#ff0000",
  "Muy dañina" = "#8f3f97",
  "Peligrosa" = "#7e0023"
)

df$categoria_aqi <- categorizar_aqi(df$aqi)

conteo <- df %>%
  dplyr::mutate(categoria_aqi = factor(categoria_aqi, levels = orden_categorias)) %>%
  dplyr::count(categoria_aqi, .drop = FALSE) %>%
  dplyr::mutate(porcentaje = n / sum(n) * 100)

plot_ly(
  conteo, x = ~categoria_aqi, y = ~porcentaje, type = "bar",
  marker = list(color = colores_categorias[as.character(conteo$categoria_aqi)])
) %>%
  layout(
    title = "% del tiempo en cada categoría de AQI",
    xaxis = list(title = "Categoría", categoryorder = "array", categoryarray = orden_categorias, showgrid = FALSE),
    yaxis = list(title = "% del tiempo", showgrid = FALSE),
    showlegend = FALSE
  )
```

### 6.4 Distribución de AQI por día de la semana

```{r boxplot-dia-semana}
orden_dias <- c("Lunes", "Martes", "Miércoles", "Jueves", "Viernes", "Sábado", "Domingo")

plot_ly(
  df, x = ~factor(dia_semana, levels = orden_dias), y = ~aqi, type = "box",
  marker = list(color = "#378ADD"), line = list(color = "#378ADD")
) %>%
  layout(
    title = "Distribución de AQI por día de la semana",
    xaxis = list(title = "Día", showgrid = FALSE),
    yaxis = list(title = "AQI", showgrid = FALSE)
  )
```