Laboratorio II: Preprocesado

Preprocesamiento

A partir del dataset MPI_subnational.csv (Multidimensional Poverty Measures), se solicita trabajar sobre las siguientes consignas:

Fuente de datos

1. Integración de datos

Para comenzar cargamos las librerías necesarias para esta práctica:

library(pacman)
pacman::p_load(tidyverse, modeest, WVPlots, DT, plotly, GGally, gplots, infotheo)
options(warn=-1)

Luego que nada definimos las funciones ue luego usaremos para graficar:

plot_heatmap <- function(data) {
  p <- plot_ly(
    z = data, 
    y = colnames(data),
    x = rownames(data),
    colors = colorRamp(c("white", "red")),
    type = "heatmap"
  )
  p <- ggplotly(p)
  p <- htmltools::div(p, align="center" )
  p
}

Analice e integre los datasets MPI_subnational.csv y MPI_national.csv1. Tenga en cuenta las cuestiones trabajadas en clase como el método de integración, los nombres de las variables, granularidad, representación, etc.

Luego cargamos el dataset a analizar:

read_csv <- function(name) {
  read.csv(
    paste('https://raw.githubusercontent.com/adrianmarino/dm/master/datasets/', name, '.csv', sep=''), 
    header = TRUE, 
    sep = ','
  )
}

sub_national = read_csv('MPI_subnational')
national = read_csv('MPI_national')

Columnas:

names(sub_national)

## [1] "ISO.country.code"                  "Country"                          
## [3] "Sub.national.region"               "World.region"                     
## [5] "MPI.National"                      "MPI.Regional"                     
## [7] "Headcount.Ratio.Regional"          "Intensity.of.deprivation.Regional"

names(national)

## [1] "ISO"                            "Country"                       
## [3] "MPI.Urban"                      "Headcount.Ratio.Urban"         
## [5] "Intensity.of.Deprivation.Urban" "MPI.Rural"                     
## [7] "Headcount.Ratio.Rural"          "Intensity.of.Deprivation.Rural"

Renombramos las columnas del dataset a nombres mas cortos y reemplazamos valores nulos:

sub_national = sub_national %>%
  rename(
    world_region        = World.region,
    country_code        = ISO.country.code,
    country             = Country,
    country_mpi         = MPI.National,
    country_region      = Sub.national.region,
    country_region_mpi  = MPI.Regional,
    country_region_hc   = Headcount.Ratio.Regional,
    country_region_iod  = Intensity.of.deprivation.Regional
  ) %>% 
  replace_na(list(
    country_mpi = 0,
    country_region_mpi = 0,
    country_region_hc  = 0,
    country_region_iod = 0
  ))

sub_national$world_region   <- trimws(sub_national$world_region,   which = c("both"))
sub_national$country        <- trimws(sub_national$country,        which = c("both"))
sub_national$country_region <- trimws(sub_national$country_region, which = c("both"))

national = national %>%
  rename(
    country_code      = ISO,
    country           = Country,
    country_urban_mpi = MPI.Urban,
    country_urban_hc  = Headcount.Ratio.Urban,
    country_urban_iod = Intensity.of.Deprivation.Urban,
    country_rural_mpi = MPI.Rural,
    country_rural_hc  = Headcount.Ratio.Rural,
    country_rural_iod = Intensity.of.Deprivation.Rural
  ) %>% 
  replace_na(list(
    country_urban_mpi = 0,
    country_urban_hc  = 0,
    country_urban_iod = 0,
    
    country_rural_mpi = 0,
    country_rural_hc  = 0,
    country_rural_iod = 0
  ))

plot_heatmap <- function(data) {
  plot_ly(
    z = data, 
    y = colnames(data),
    x = rownames(data),
    colors = colorRamp(c("white", "red")),
    type = "heatmap"
  )
}

plot_smooth_compare <- function(original, smooths, name, smoth_names, colors) {
  colors <- c(c("black"), colors)
  
  # Plot original curve...
  plot(
    original, 
    type = "l", 
    col  = colors[1], 
    ylab = name,
    xlab = "Observaciones", 
    main = "Original vs suavizado"
  )

  # Plot smooth versions...
  legends <- list("Original")
  i <- 2
  for (smooth in smooths) {
    lines(smooth, type = "l", col=colors[i])
    legends[i] <- smoth_names[i-1]
    i <- i+1
  }
  legend("topleft", legend=unlist(legends), col=colors, lty=1)
}

En el siguiente diagrama se describen los datos y sus relaciones:

Modelo Relacional

Tenemos mas países en el dataset national:

sub_national_countries_count <- sub_national %>% 
  summarise(count = n_distinct(country))

print(paste('Paises en el dataset sub_national:', sub_national_countries_count))

## [1] "Paises en el dataset sub_national: 78"

national_countries_count <- national %>% 
  summarise(count = n_distinct(country))

print(paste('Paises en el dataset national:', national_countries_count))

## [1] "Paises en el dataset national: 102"

Si hacemos un join por país nos quedarían regiones sin países.

Hay dos alternativas:

Podemos integrar ambos datasets por country, y calcular la media de las columnas: mpi_regional, hc_regional y iod_regional por país. De esta manera, nos quedarían todas las métrica a nivel país y podríamos segregar por pases y world_region. Por otro lado, deberíamos llenar los campos faltantes para los 24 países del dataset national que no tienen registros asociados en el dataset sub_national. Finalmente tengamos en cuenta que vamos a tener muchas métricas iguales (mpi, hc y iod), las cuales seguramente van a estar altamente correlacionadas.
Otra alternativa es quedarnos con el dataset nacional, no tenemos datos faltantes y seguimos teniendo las mismas métricas pero a nivel país.

Vamos a optar por la segunda opción:

dataset <- sub_national %>%
  select(world_region, country_code) %>%
  distinct_all() %>%
  left_join(national, by ='country_code')  %>%
  select(
    world_region,
    country,
    country_urban_mpi,
    country_urban_hc,
    country_urban_iod,
    country_rural_mpi,
    country_rural_hc,
    country_rural_iod
  )

Columnas:

names(dataset)

## [1] "world_region"      "country"           "country_urban_mpi"
## [4] "country_urban_hc"  "country_urban_iod" "country_rural_mpi"
## [7] "country_rural_hc"  "country_rural_iod"

head(dataset)

Finalmente nos quedaremos con esta parte del diagrama anterior:

Modelo Relacional Final

Veamos los países por región y que cantidad hay en cada una:

cities.count.by.region <- dataset %>%
  group_by(world_region) %>%
  summarise(count = n_distinct(country)) %>%
  arrange(world_region)

cities.count.by.region.country <- dataset %>%
  group_by(world_region, country) %>%
  summarise(count = n_distinct(country), .groups = "keep") %>%
  arrange(world_region, country)

parent_regions  <- cities.count.by.region$world_region
regions         <- cities.count.by.region.country$world_region
countries       <- cities.count.by.region.country$country

cities.count         <- sum(cities.count.by.region.country$count)
country.cities.count <- cities.count.by.region.country$count
region.cities.count  <- cities.count.by.region$count

labels  <- c(c(parent_regions), c(countries))
parents <- c(rep(" ", length(parent_regions)), regions)
values  <- c(c(region.cities.count), c(country.cities.count))

# print(paste(length(labels), length(parents), length(values), sep=" - "))
p <- plot_ly(labels = labels, parents = parents, values = values, type = 'sunburst', maxdepth=2)
p <- htmltools::div(p, align="center" )
p

2. Atributos redundantes

Verifique si existen atributos (categóricos o numéricos)redundantes en el dataset y actúe en consecuencia de acuerdo a las técnicas abordadas en clase.

En el paso anterior agregamos por world_region y country y nos quedamos con los pares únicos para no tener redundancia. Por cada par world_region-country hay un solo set de métricas quedando una tabla agregada por word_region y country.

Por otro lado veamos que nivel de correlación tiene la variables numéricas.

Matriz de correlación de pearson

cor_matrix = cor(dataset[3:8])
cor_matrix[upper.tri(cor_matrix)] <- NA
cor_matrix

##                   country_urban_mpi country_urban_hc country_urban_iod
## country_urban_mpi         1.0000000               NA                NA
## country_urban_hc          0.9957482        1.0000000                NA
## country_urban_iod         0.9017424        0.9045761         1.0000000
## country_rural_mpi         0.9048242        0.9228277         0.9136211
## country_rural_hc          0.8705625        0.8970252         0.9019570
## country_rural_iod         0.8652554        0.8753630         0.9148407
##                   country_rural_mpi country_rural_hc country_rural_iod
## country_urban_mpi                NA               NA                NA
## country_urban_hc                 NA               NA                NA
## country_urban_iod                NA               NA                NA
## country_rural_mpi         1.0000000               NA                NA
## country_rural_hc          0.9854434        1.0000000                NA
## country_rural_iod         0.9612909        0.9221695                 1

plot_heatmap(cor_matrix)

Observaciones

Se puede apreciar que todas la métricas en general tienen una correlación muy alta, como mínimo es 0.8.

3. Manejo de Ruido

3.a. Verifique en primer lugar la distribución de los datos, utilice algún método gráfico para esto. A su criterio, ¿Cuál es la variable más “ruidosa”?

Veamos las distribuciones de cada variable:

dataset.n <- scale(dataset[,3:8])
dataset.n <- as.data.frame(dataset.n)
summary(dataset.n)

##  country_urban_mpi country_urban_hc  country_urban_iod  country_rural_mpi
##  Min.   :-1.0081   Min.   :-1.0990   Min.   :-1.90359   Min.   :-1.3474  
##  1st Qu.:-0.8183   1st Qu.:-0.8518   1st Qu.:-0.78204   1st Qu.:-0.9367  
##  Median :-0.1963   Median :-0.1453   Median : 0.02796   Median :-0.1651  
##  Mean   : 0.0000   Mean   : 0.0000   Mean   : 0.00000   Mean   : 0.0000  
##  3rd Qu.: 0.3836   3rd Qu.: 0.4694   3rd Qu.: 0.78085   3rd Qu.: 0.7805  
##  Max.   : 3.8207   Max.   : 3.3375   Max.   : 2.68643   Max.   : 2.1159  
##  country_rural_hc   country_rural_iod
##  Min.   :-1.55311   Min.   :-1.7279  
##  1st Qu.:-0.94045   1st Qu.:-0.8535  
##  Median :-0.01497   Median :-0.2110  
##  Mean   : 0.00000   Mean   : 0.0000  
##  3rd Qu.: 0.90724   3rd Qu.: 0.7386  
##  Max.   : 1.53790   Max.   : 2.5470

Comparemos las distribuciones:

p <- dataset.n %>% select(
      country_urban_mpi, 
      country_rural_mpi, 
      country_urban_hc, 
      country_rural_hc, 
      country_urban_iod, 
      country_rural_iod
    ) %>%
    pivot_longer(
      ., 
      cols = c(
        country_urban_mpi, 
        country_rural_mpi, 
        country_urban_hc, 
        country_rural_hc, 
        country_urban_iod, 
        country_rural_iod
      ), 
      names_to = "Variables", 
      values_to = "Frecuencia"
    ) %>%
    ggplot(aes(x = Variables, y = Frecuencia, fill = Variables)) +
    geom_violin(trim=FALSE, fill = "transparent") + 
    geom_boxplot(width=0.1) + 
    geom_jitter(shape=16, position=position_jitter(0.2))

ggplotly(p)

library(statip)
apply(dataset[,3:8], 2, statip::cv, na_rm=TRUE)

## country_urban_mpi  country_urban_hc country_urban_iod country_rural_mpi 
##         0.9817094         0.8923289         0.1125864         0.7225289 
##  country_rural_hc country_rural_iod 
##         0.6202180         0.1621852

Observaciones

La realidad es que necesitamos tener conocimiento del dominio para entender que variable tiene ruido o no. No podemos asumir que tener outliers o alta variabilidad es un indicador de ruido. La única forma de saber si una variable tiene ruido es conocer sobre el dominio del problema y de esa manera si se puede detectar cuando una medida fue mal tomada o se le suma o resta algún ruido a un valor.

En este caso vamos a elegir una variable cualquiera para hacer el suavizado.

3.b. Realice un suavizado utilizando binning por frecuencias iguales, estime la cantidad de bins óptima y calcule el valor para bin por el cálculo de medias. Grafique las series resultantes y comente los resultados observados.

3.c. Ahora aplique la técnica por anchos iguales del mismo modo que en el punto anterior. ¿Qué observa para esta técnica ante la presencia de outliers?

3.d. Compare los métodos de suavizado de los puntos b. y c.

smooth <- function(
  values,
  criteria = "equalfreq",
  bin_size  = 5,
  fn        = mean
) {
  bin_df   <- discretize(values, criteria, bin_size)

  bin_data <- data.frame(bin = bin_df$X, value = values)
  
  for(bin in 1:bin_size) {
    indexes = bin_data$bin == bin
    bin_data$smooth[indexes] = fn(bin_data$value[indexes])
  }
  bin_data
}

df_smooth_1 <- smooth(dataset$country_urban_iod, bin=3)
df_smooth_2 <- smooth(dataset$country_urban_iod, criteria="equalwidth", bin=3)

plot_smooth_compare(
  original    = dataset$country_urban_iod,
  smooths     = list(df_smooth_1$smooth, df_smooth_2$smooth),
  smoth_names = c('equalfreq', 'equalwidth'),
  colors      = c("yellow", "blue"),
  name        = "country_urban_iod"
)

Conclusiones

Al parecer el método de igual frecuencia parece copiar con mas fidelidad la señal original sin recortarla demasiado los picos.
El método de igual ancho recorta mucho mas los picos pero también en otros casos copia más outliers con mas facilidad.
Pareciera que el método de igual frecuencia suaviza mejor los picos. Se puede apreciar que la señal se recorta de forma mas homogénea.

Laboratorio II: Preprocesado - Parte I

Adrian Marino

4/01/2021

Introducción