Ejercicios de Práctica # 3-5
Juan Mendoza
2024-08-30
Actividad Práctica # 3
El análisis se basa en el conjunto de datos “Suicide Rates Overview 1985 to 2016”, disponible en Kaggle. Este conjunto de datos, que abarca 27,820 observaciones, proporciona información socioeconómica y demográfica de diversos países, con el objetivo de identificar señales correlacionadas con el aumento de las tasas de suicidio globales durante el período de 1985 a 2016.
df <- read_csv("C:/Users/jcami/Downloads/archive/master.csv")
df %<>% clean_namesAnálisis inicial
En un primer momento, estandarizamos el nombre de las variables
utilizando la función clean_names y estudiamos las
dimensiones de la base de datos. Tal como se puede observar a
continuación, el dataset cuenta con 27,820 observaciones y 12 variables,
entre las cuales se encuentran: país, year, sexo, edad (o rango de
edad), número de suicidios, población, tasa de suicidio por cada 100,000
habitantes, entre otras.
df %>% glimpse## Rows: 27,820
## Columns: 12
## $ country <chr> "Albania", "Albania", "Albania", "Albania", "Albania…
## $ year <dbl> 1987, 1987, 1987, 1987, 1987, 1987, 1987, 1987, 1987…
## $ sex <chr> "male", "male", "female", "male", "male", "female", …
## $ age <chr> "15-24 years", "35-54 years", "15-24 years", "75+ ye…
## $ suicides_no <dbl> 21, 16, 14, 1, 9, 1, 6, 4, 1, 0, 0, 0, 2, 17, 1, 14,…
## $ population <dbl> 312900, 308000, 289700, 21800, 274300, 35600, 278800…
## $ suicides_100k_pop <dbl> 6.71, 5.19, 4.83, 4.59, 3.28, 2.81, 2.15, 1.56, 0.73…
## $ country_year <chr> "Albania1987", "Albania1987", "Albania1987", "Albani…
## $ hdi_for_year <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ gdp_for_year <dbl> 2156624900, 2156624900, 2156624900, 2156624900, 2156…
## $ gdp_per_capita <dbl> 796, 796, 796, 796, 796, 796, 796, 796, 796, 796, 79…
## $ generation <chr> "Generation X", "Silent", "Generation X", "G.I. Gene…En cuanto al tipo de variable, se observa que la mayoría de las variables tipo factor vienen en forma de cadena de caracteres, por lo cual a continuación se procede a cambiar el tipo de variable para que sea acorde con la descripción del dataset:
selected_columns <- c("country","sex","country_year","generation","age")
df[selected_columns] <- lapply(df[selected_columns], as.factor)Además de las variables categóricas ya mencionadas, el dataset cuenta con variables continuas relacionadas con la población, la tasa de suicidio por cada 100,000 habitantes, el índice de desarrollo humano por year, el producto interno bruto (PIB) por year y el producto interno bruto per cápita. A continuación, se presenta una descripción general de la base de datos:
summary(df)## country year sex age
## Austria : 382 Min. :1985 female:13910 15-24 years:4642
## Iceland : 382 1st Qu.:1995 male :13910 25-34 years:4642
## Mauritius : 382 Median :2002 35-54 years:4642
## Netherlands: 382 Mean :2001 5-14 years :4610
## Argentina : 372 3rd Qu.:2008 55-74 years:4642
## Belgium : 372 Max. :2016 75+ years :4642
## (Other) :25548
## suicides_no population suicides_100k_pop country_year
## Min. : 0.0 Min. : 278 Min. : 0.00 Albania1987: 12
## 1st Qu.: 3.0 1st Qu.: 97498 1st Qu.: 0.92 Albania1988: 12
## Median : 25.0 Median : 430150 Median : 5.99 Albania1989: 12
## Mean : 242.6 Mean : 1844794 Mean : 12.82 Albania1992: 12
## 3rd Qu.: 131.0 3rd Qu.: 1486143 3rd Qu.: 16.62 Albania1993: 12
## Max. :22338.0 Max. :43805214 Max. :224.97 Albania1994: 12
## (Other) :27748
## hdi_for_year gdp_for_year gdp_per_capita generation
## Min. :0.483 Min. :4.692e+07 Min. : 251 Boomers :4990
## 1st Qu.:0.713 1st Qu.:8.985e+09 1st Qu.: 3447 G.I. Generation:2744
## Median :0.779 Median :4.811e+10 Median : 9372 Generation X :6408
## Mean :0.777 Mean :4.456e+11 Mean : 16866 Generation Z :1470
## 3rd Qu.:0.855 3rd Qu.:2.602e+11 3rd Qu.: 24874 Millenials :5844
## Max. :0.944 Max. :1.812e+13 Max. :126352 Silent :6364
## NA's :19456Finalmente, concluimos esta primera sección examinando la presencia
de datos faltantes haciendo uso de la función missmap del
paquete Amelia. Tal como se puede observar a continuación, la variable
asociada al índice de desarrollo humano por year cuenta con una gran
cantidad de datos faltantes:
missmap(df)
Se nota que, aunque en general estos datos faltantes solo representan alrededor del 6% del total de datos, si se analiza solo la proporción de faltantes dentro de la variable, se llega a la conclusión de que estos datos representan alrededor del 70% de las observaciones del índice de desarrollo humano por year:
missing_data <- ifelse(is.na(df$hdi_for_year),1,0)
porcentaje_NA <- round(sum(missing_data)/length(df$hdi_for_year)*100,3)
porcentaje_NA## [1] 69.935Estudio de las varibales
Paises
El primer análisis de los datos no indica que no hay datos faltantes y que tienden a tener un numero consistente de observaciones. A continuación se presentan las primeras 10 observaciones del dataframe con el fin de entender mejor la estructura de esta variable:
head(df,10)## # A tibble: 10 × 12
## country year sex age suicides_no population suicides_100k_pop
## <fct> <dbl> <fct> <fct> <dbl> <dbl> <dbl>
## 1 Albania 1987 male 15-24 years 21 312900 6.71
## 2 Albania 1987 male 35-54 years 16 308000 5.19
## 3 Albania 1987 female 15-24 years 14 289700 4.83
## 4 Albania 1987 male 75+ years 1 21800 4.59
## 5 Albania 1987 male 25-34 years 9 274300 3.28
## 6 Albania 1987 female 75+ years 1 35600 2.81
## 7 Albania 1987 female 35-54 years 6 278800 2.15
## 8 Albania 1987 female 25-34 years 4 257200 1.56
## 9 Albania 1987 male 55-74 years 1 137500 0.73
## 10 Albania 1987 female 5-14 years 0 311000 0
## # ℹ 5 more variables: country_year <fct>, hdi_for_year <dbl>,
## # gdp_for_year <dbl>, gdp_per_capita <dbl>, generation <fct>Ahora se presentan las ultimas 10 observaciones del dataframe.
tail(df,10)## # A tibble: 10 × 12
## country year sex age suicides_no population suicides_100k_pop
## <fct> <dbl> <fct> <fct> <dbl> <dbl> <dbl>
## 1 Uzbekistan 2014 female 15-24 years 347 2992817 11.6
## 2 Uzbekistan 2014 male 55-74 years 144 1271111 11.3
## 3 Uzbekistan 2014 male 15-24 years 347 3126905 11.1
## 4 Uzbekistan 2014 male 75+ years 17 224995 7.56
## 5 Uzbekistan 2014 female 25-34 years 162 2735238 5.92
## 6 Uzbekistan 2014 female 35-54 years 107 3620833 2.96
## 7 Uzbekistan 2014 female 75+ years 9 348465 2.58
## 8 Uzbekistan 2014 male 5-14 years 60 2762158 2.17
## 9 Uzbekistan 2014 female 5-14 years 44 2631600 1.67
## 10 Uzbekistan 2014 female 55-74 years 21 1438935 1.46
## # ℹ 5 more variables: country_year <fct>, hdi_for_year <dbl>,
## # gdp_for_year <dbl>, gdp_per_capita <dbl>, generation <fct>Al examinar más detalladamentes el número de observaciones por país, se evindencia que no se cuenta con los datos de algunos paises en algunos years. Más aún, no encuentran datos faltantes que más alla de lo ya mencionado.
which(is.na(df$country)) # no se encuentran datos faltantes## integer(0)head(table(df$country),20)##
## Albania Antigua and Barbuda Argentina
## 264 324 372
## Armenia Aruba Australia
## 298 168 360
## Austria Azerbaijan Bahamas
## 382 192 276
## Bahrain Barbados Belarus
## 252 300 252
## Belgium Belize Bosnia and Herzegovina
## 372 336 24
## Brazil Bulgaria Cabo Verde
## 372 360 12
## Canada Chile
## 348 372Ahora se construir una tabla de frecuencia con medidas descriptivas:
Tabla_pais <- df %>%
group_by(country) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_pais## # A tibble: 101 × 4
## country n Proporcion Porcentaje
## <fct> <int> <dbl> <chr>
## 1 Mongolia 10 0.000359 0.0359%
## 2 Cabo Verde 12 0.000431 0.0431%
## 3 Dominica 12 0.000431 0.0431%
## 4 Macau 12 0.000431 0.0431%
## 5 Bosnia and Herzegovina 24 0.000863 0.0863%
## 6 Oman 36 0.00129 0.1294%
## 7 Saint Kitts and Nevis 36 0.00129 0.1294%
## 8 San Marino 36 0.00129 0.1294%
## 9 Nicaragua 72 0.00259 0.2588%
## 10 United Arab Emirates 72 0.00259 0.2588%
## # ℹ 91 more rowsFinalmente se culmina que el estudio de esta variable, se presentara un diagrama de barra con los paises con el menor número de observaciones debido a la gran número de categorias.
Tabla_pais <- Tabla_pais %>%
arrange(desc(Proporcion))
ggplot(tail(Tabla_pais,20), aes(x=reorder(
country, -n),y=n))+
geom_bar(stat="identity", position=position_dodge(),fill="darkolivegreen")+
theme(plot.title = element_text(hjust = 0.5, size = 12, vjust = 2),
axis.text.x = element_text(angle = 90, hjust = 0.8, vjust = 0.5, ))+
ggtitle("Observaciones por pais")+
xlab("Pais")+
ylab("Observaciones")
years
En el apartado anterior observamos que cada paises tiene varias observaciones por year y que en principio esta variable no presenta datos faltantes.Con el fin de corroborar esto y de conocer de manera más detallada cuantas observaciones se cuenta por cada year, se presenta un tabla que se resume esta información:
table(df$year)##
## 1985 1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 1997 1998 1999 2000
## 576 576 648 588 624 768 768 780 780 816 936 924 924 948 996 1032
## 2001 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012 2013 2014 2015 2016
## 1056 1032 1032 1008 1008 1020 1032 1020 1068 1056 1032 972 960 936 744 160Dado que no se observan datos faltantes, procedemos a visualizar el número de observaciones por year en un grafico de barras:
Tabla_year <- df %>%
group_by(year) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_year## # A tibble: 32 × 4
## year n Proporcion Porcentaje
## <dbl> <int> <dbl> <chr>
## 1 2016 160 0.00575 0.575%
## 2 1985 576 0.0207 2.070%
## 3 1986 576 0.0207 2.070%
## 4 1988 588 0.0211 2.114%
## 5 1989 624 0.0224 2.243%
## 6 1987 648 0.0233 2.329%
## 7 2015 744 0.0267 2.674%
## 8 1990 768 0.0276 2.761%
## 9 1991 768 0.0276 2.761%
## 10 1992 780 0.0280 2.804%
## # ℹ 22 more rowsggplot(Tabla_year, aes(x = factor(year), y = n)) +
geom_bar(stat = "identity", position = position_dodge(), fill = "cyan4") +
theme(plot.title = element_text(hjust = 0.5, size = 12, vjust = 2),
axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
ggtitle("Observaciones por perido") +
xlab("Periodo") +
ylab("Observaciones") +
scale_x_discrete(breaks = unique(Tabla_year$year))
### Sexo
En lo que concierne al número de observaciones por sexo, el analisis inical muestra que no hay datos faltantes y que las observaciones entre hombres y mujeres es similar. Debido a esto, a continuación se presenta un tabla de con el porcentaje de observaciones por sexo y un diagrama de barras que ayudara a detallar la distribución de esta variables en el dataset:
Tabla_sex <- df %>%
group_by(sex) %>%
count() %>%
ungroup() %>%
mutate(proporcion = `n` / sum(`n`)) %>%
arrange(proporcion) %>%
mutate(Porcentaje = scales::percent(proporcion))
Tabla_sex## # A tibble: 2 × 4
## sex n proporcion Porcentaje
## <fct> <int> <dbl> <chr>
## 1 female 13910 0.5 50%
## 2 male 13910 0.5 50%ggplot(Tabla_sex , aes(x = "", y = proporcion, fill = sex)) +
geom_col(color = "black") +
geom_label(aes(label = Porcentaje),
position = position_stack(vjust = 0.5),
show.legend = FALSE) +
guides(fill = guide_legend(title = "Sexo")) + scale_color_gradient() +
coord_polar(theta = "y") + ggtitle ("Distribucion de sexo")+
theme(plot.title = element_text(hjust = 0.5, size = 12, vjust = 2))
Rango de edad
En lo que respecta al rango de edad, el análisis inicial nos indica
que no hay datos faltantes. Con el fin de examinar de manera más
detallada la variable age a continuación se presenta una
tabla de frecuencias:
table(df$age)##
## 15-24 years 25-34 years 35-54 years 5-14 years 55-74 years 75+ years
## 4642 4642 4642 4610 4642 4642Tal como se noto inicialmente, la tabla indica que nos hay observaciones faltantes y, más aún, todas las categorias de edad tienen el mismo número de observaciones. Por lo cual se procede a construir la tabla de frecuencia con medidas descritivas:
Por lo cual se culmina por el momento en análisis de esta varibale
presentado una diagrama de barras de age:
Tabla_age <- df %>%
group_by(age) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_age## # A tibble: 6 × 4
## age n Proporcion Porcentaje
## <fct> <int> <dbl> <chr>
## 1 5-14 years 4610 0.166 16.57%
## 2 15-24 years 4642 0.167 16.69%
## 3 25-34 years 4642 0.167 16.69%
## 4 35-54 years 4642 0.167 16.69%
## 5 55-74 years 4642 0.167 16.69%
## 6 75+ years 4642 0.167 16.69%ggplot(Tabla_age, aes(x=age,y=n))+
geom_bar(stat="identity", position=position_dodge(),fill="cyan4")+
coord_cartesian(ylim = c(0, 5000))+
theme(plot.title = element_text(hjust = 0.5, size = 12, vjust = 2),
axis.text.x = element_text(angle = 90, hjust = 0.7, vjust = 0.5))+
ggtitle("Observaciones por rango de edad")+
xlab("Rango de edad")+
ylab("Observaciones")
Generación
Tal como en el caso anterior, el análisis inicial indico que no hay
datos faltantes en la variablegeneration. A continuación,
se presente una tabla de frecuencias absolutas:
table(df$generation)##
## Boomers G.I. Generation Generation X Generation Z Millenials
## 4990 2744 6408 1470 5844
## Silent
## 6364Ahora se presenta la tabla con las respectivas medidas descriptivas:
Tabla_gene <- df %>%
group_by(generation) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_gene## # A tibble: 6 × 4
## generation n Proporcion Porcentaje
## <fct> <int> <dbl> <chr>
## 1 Generation Z 1470 0.0528 5.28%
## 2 G.I. Generation 2744 0.0986 9.86%
## 3 Boomers 4990 0.179 17.94%
## 4 Millenials 5844 0.210 21.01%
## 5 Silent 6364 0.229 22.88%
## 6 Generation X 6408 0.230 23.03%La tabla de frecuencias confirma la supocisión en cuanto la ausencia de datos faltantes y que cada una de los generaciones difieren en el número de las observaciones, para ver esto a más detalle de presenta un diagrama de barras:
ggplot(Tabla_gene, aes(x=reorder(
generation, -n),y=n))+
geom_bar(stat="identity", position=position_dodge(),fill="darkolivegreen")+
theme(plot.title = element_text(hjust = 0.5, size = 12, vjust = 2),
axis.text.x = element_text(angle = 90, hjust = 0.8, vjust = 0.5, ))+
ggtitle("Observaciones por generacion")+
xlab("generacion")+
ylab("Observaciones")
Población
Ahora analizaremos las varibale population. Dedido a que
esta es una varibale númerica se comenzara construyendo un histograma y
un diagrama de cajas (Boxplot) para estudiar la distribución de la
varibale.
par(mfrow = c(1, 2))
hist(df$population, xlab = "Poblacion", ylab = "Frecuencia", main = "Poblacion por year")
boxplot <- boxplot(df$population, ylab = "Poblacion por year", main = "Boxplot")
Los graficos muestran que la distribución esta fuertamente sesgada a
la izquierda. Dicho esto, en el contexto de las varibales relacionados
con la economia se suele utilizar la transformación logaritmica para
centrar estas distribuciones. Debido a esto se creara una nueva varibale
population_logy se estudiara su distribución:
df$population_log <- log(df$population)Ahora volveremos a generar el histograma y el diagrama de cajas para estudiar la distribución de la nueva variable:
par(mfrow = c(1, 2))
hist(df$population_log, xlab = "Poblacion por year (Log)", ylab = "Frecuencia", main = "Poblacion por year")
boxplot <- boxplot(df$population_log, ylab = "Poblacion por year", main = "Boxplot")
Tal como se puede apreciar el los graficos la distribución de esta
variable es más simetrica, con lo cual concluimos que la transformación
fue beneficiosa. Sin embargo, el diagrama de cajas aún muestra
observaciones que pueden llegar a ser datos atipicos.
Datos atipicos
Los graficos revelan la existencia de valores atipicos en la variable
population_log . Por lo cual, antes de proseguir
utilizaremos el filtro de Hampel y el método de los percentiles para
estudiar más a profundidad estas observaciones.
# Calcular los límites usando el filtro de Hampel
lower_bound <- median(df$population_log, na.rm = TRUE) - 3 * mad(df$population_log, constant = 1, na.rm = TRUE)
upper_bound <- median(df$population_log, na.rm = TRUE) + 3 * mad(df$population_log, constant = 1, na.rm = TRUE)
# Identificar los índices de los outliers
outlier_ind <- which(df$population_log < lower_bound | df$population_log > upper_bound)
# Crear una columna que marque los outliers
df$Hampel <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)
# Ver el número de outliers
length(outlier_ind)## [1] 1491lower_bound <- quantile(df$population_log, 0.025)
upper_bound <- quantile(df$population_log, 0.975)
outlier_ind <- which(df$population_log < lower_bound | df$population_log > upper_bound)
# Crear una columna que marque los outliers
df$Perc <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)
# Ver el número de outliers
length(outlier_ind)## [1] 1391table(df$Hampel, df$Perc)##
## 0 1
## 0 25723 606
## 1 706 785Vemos que existen discrepancias en entre el filtor de Hampel y el metodo de los percentiles. Por se aplicara la prueba de Rosner para confirmar la existencias de estos datos atipicos:
test <- rosnerTest(df$population_log, k = 1000)
test <- test$all.stats
outlier_ind <- test[test$Outlier == TRUE,]$Obs.Num
df$Rosner <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)Los resultado de la prueba de Rosner indico que ninguno de las
observaciones analisadas son datos atipicos, por lo cual se culmina esta
sección con medidas descriptivas de la variable
population_log:
df %>%
get_summary_stats(population_log, type = "five_number")## # A tibble: 1 × 7
## variable n min max q1 median q3
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 population_log 27820 5.63 17.6 11.5 13.0 14.2Suicidios por year
La variable suicides_no contiene información sobre los
suicidios por year. Debido a su naturaleza procederemos de misma manera
que en la sección anterior:
par(mfrow = c(1, 2))
hist(df$suicides_no, xlab = "Suicidios", ylab = "Frecuencia", main = "Suicidios por year")
boxplot(df$suicides_no, ylab = "Suicidios por year", main = "Boxplot")
La distribución de suicides_no presenta un sesgo
significativo hacia la izquierda. Sin embargo, no es viable aplicar una
transformación debido a la presencia de ceros en los datos. Dado que,
por la naturaleza de la variable, se anticipa una alta frecuencia de
valores bajos en comparación con los valores elevados, no es posible
descartar que los valores que se desvían de la mediana sean errores de
codificación. Por lo tanto, esta sección concluirá presentando algunas
métricas descriptivas.
df%>%
get_summary_stats(suicides_no, type = "five_number")## # A tibble: 1 × 7
## variable n min max q1 median q3
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 suicides_no 27820 0 22338 3 25 131Producto interno bruto por year
En lo que concierne al producto interno producto por year
(gdp_for_year), la exploración inicial de los datos no
revela datos faltanes. Con el fin de estudiar de manera más detallada de
distribucipon de la variables se genera un histograma de frecuencia y un
diagrama de cajas:
par(mfrow = c(1, 2))
hist(df$gdp_for_year, xlab = "GDP por year", ylab = "Frecuencia", main = "Histograma del GDP por year")
boxplot(df$gdp_for_year, ylab = "GDP por year", main = "Boxplot del DGP por year")
Al trabajar con variables economicas es usual utiliar la
transformación logaritmica, por lo cual vamos a calcular el logartimo
del producto interno bruto por year y lo registramos como una nueva
variable en el conjunto de datos (gdp_for_year_log):
df$gdp_for_year_log <- log(df$gdp_for_year)Hecho esto, volvemos volvemos a revisar la distribución de
gdp_for_year_log utlizando un histograma y un diagrama de
cajas:
par(mfrow = c(1, 2))
hist(df$gdp_for_year_log, xlab = "GDP por year (Log)", ylab = "Frecuencia", main = "Histograma del GDP por year")
boxplot <- boxplot(df$gdp_for_year_log, ylab = "GDP por year (Log)", main = "Boxplot del DGP por year")
Del diagrama de cajas, observamos que existen posibles datos atipicos
en la variable gdp_for_year_log. Por lo cual, vamos a
proceder a un estudio detallado de estas observaciones:
Datos atipicos
Los graficos revelan la existencia de valores atipicos en la variable
gdp_for_year_log . Por lo cual, antes de proseguir
utilizaremos el filtro de Hampel y el método de los percentiles para
estudiar más a profundidad estas observaciones.
# Calcular los límites usando el filtro de Hampel
lower_bound <- median(df$gdp_for_year_log, na.rm = TRUE) - 3 * mad(df$gdp_for_year_log, constant = 1, na.rm = TRUE)
upper_bound <- median(df$gdp_for_year_log, na.rm = TRUE) + 3 * mad(df$gdp_for_year_log, constant = 1, na.rm = TRUE)
# Identificar los índices de los outliers
outlier_ind <- which(df$gdp_for_year_log < lower_bound | df$gdp_for_year_log > upper_bound)
# Crear una columna que marque los outliers
df$Hampel <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)
# Ver el número de outliers
length(outlier_ind)## [1] 636lower_bound <- quantile(df$gdp_for_year_log, 0.025)
upper_bound <- quantile(df$gdp_for_year_log, 0.975)
outlier_ind <- which(df$gdp_for_year_log < lower_bound | df$gdp_for_year_log > upper_bound)
# Crear una columna que marque los outliers
df$Perc <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)
# Ver el número de outliers
length(outlier_ind)## [1] 1392table(df$Hampel, df$Perc)##
## 0 1
## 0 26428 756
## 1 0 636Vemos que existen discrepancias en entre el filtor de Hampel y el metodo de los percentiles. Por se aplicara la prueba de Rosner para confirmar la existencias de estos datos atipicos:
test <- rosnerTest(df$gdp_for_year_log, k = 1000)
test <- test$all.stats
outlier_ind <- test[test$Outlier == TRUE,]$Obs.Num
df$Rosner <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)Los resultado de la prueba de Rosner indico que ninguno de las
observaciones analisadas son datos atipicos, por lo cual se culmina esta
sección con medidas descriptivas de la variable
gdp_for_year_log:
nortest::ad.test(df$gdp_for_year_log)##
## Anderson-Darling normality test
##
## data: df$gdp_for_year_log
## A = 38.016, p-value < 2.2e-16nortest::lillie.test(df$gdp_for_year_log)##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: df$gdp_for_year_log
## D = 0.034712, p-value < 2.2e-16Dado que ambas pruebas indican que la distribución de
gdp_for_year_log no se aproxima a una distribución
guassiana, culiminaremos el analisis de la varibale con algunas
estadisticas descriptivas:
df %>%
get_summary_stats(gdp_for_year_log, type = "five_number")## # A tibble: 1 × 7
## variable n min max q1 median q3
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 gdp_for_year_log 27820 17.7 30.5 22.9 24.6 26.3Producto interno bruto per capita
Tal como en el caso anterior, la exploración inicial de los datos no
revela datos faltanes. Con el fin de estudiar de manera más detallada de
distribucipon de la varible gdp_per_capita se genera un
histograma de frecuencia y un diagrama de cajas:
par(mfrow = c(1, 2))
hist(df$gdp_per_capita, xlab = "GDP per capita", ylab = "Frecuencia", main = "Histograma del GDP per capita")
boxplot(df$gdp_per_capita, ylab = "GDP per capita", main = "Boxplot del DGP per capita")
Ahora procedemos a calcular una nueva variable con la función
logaritmica llamada gdp_per_capita_log y utilizamos un
histograma y una diagrama de cajas para estudiar su distribución.
df$gdp_per_capita_log <- log(df$gdp_per_capita)
par(mfrow = c(1, 2))
hist(df$gdp_per_capita_log, xlab = "GDP per capita (Log)", ylab = "Frecuencia", main = "Histograma del GDP per capita")
boxplot(df$gdp_per_capita_log, ylab = "GDP per capita (Log)", main = "Boxplot del DGP per capita")
En este caso, no se evidenicas datos atipicos, por lo cual se procede
a comporbar la distribución de la variable
gdp_per_capita_log.
nortest::ad.test(df$gdp_per_capita_log)##
## Anderson-Darling normality test
##
## data: df$gdp_per_capita_log
## A = 155.57, p-value < 2.2e-16nortest::lillie.test(df$gdp_per_capita_log)##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: df$gdp_per_capita_log
## D = 0.061759, p-value < 2.2e-16Dado que ambas pruebas indican que la distribución de
gdp_per_capita_log no se aproxima a una distribución
guassiana, culiminaremos el analisis de la varibale con algunas
estadisticas descriptivas:
df %>%
get_summary_stats(gdp_per_capita_log, type = "five_number")## # A tibble: 1 × 7
## variable n min max q1 median q3
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 gdp_per_capita_log 27820 5.52 11.7 8.14 9.14 10.1Índice de Desarrollo Humano
Al analizar la variable del Índice de Desarrollo Humano (HDI) por year, se observa una proporción considerable de datos faltantes. Por lo tanto, procederemos a 1) realizar una descripción detallada de la variable y 2) explorar la forma más adecuada de manejar estos valores faltantes.
Descripción de la variable
Primero, se explorará la distribución de la variable hdi_for_year utilizando un histograma y un diagrama de cajas:
par(mfrow = c(1, 2))
hist(df$hdi_for_year, xlab = "HDI por year", ylab = "Frecuencia", main = "Histograma del HDI por year")
boxplot(df$hdi_for_year, ylab = "HDI por year", main = "Boxplot del HDI por year")
Los gráficos sugieren que la distribución muestral del Índice de Desarrollo Humano por year está sesgada a la derecha y presenta al menos un valor atípico. Para complementar esta descripción, se realizaron dos pruebas de normalidad utilizando el paquete nortest para determinar si la variable sigue una distribución gaussiana:
nortest::lillie.test(df$hdi_for_year)##
## Lilliefors (Kolmogorov-Smirnov) normality test
##
## data: df$hdi_for_year
## D = 0.055914, p-value < 2.2e-16nortest::ad.test(df$hdi_for_year)##
## Anderson-Darling normality test
##
## data: df$hdi_for_year
## A = 46.448, p-value < 2.2e-16Dado que ambas pruebas de normalidad indican que la distribución del HDI por year no se aproxima a una distribución normal (probablemente debido a la asimetría y a la presencia de datos atípicos), finalizamos esta parte de la descripción presentando algunas medidas descriptivas robustas:
df%>%
get_summary_stats(hdi_for_year, type = "robust")## # A tibble: 1 × 4
## variable n median iqr
## <fct> <dbl> <dbl> <dbl>
## 1 hdi_for_year 8364 0.779 0.142Imputación de datos
En esta sección, exploraremos dos métodos para tratar los datos faltantes: 1) la imputación simple utilizando la mediana y 2) la imputación utilizando un modelo de regresión lineal. Consideramos inviable la omisión de observaciones con datos faltantes debido a su alta proporción en la base de datos. Además, dado que los datos están agrupados por países, la imputación por mediana podría tener limitaciones debido a posibles influencias de las agrupaciones en el HDI.
A continuación, se estudian las correlaciones entre las variables numéricas de la base de datos como punto de partida para construir el modelo de regresión. Utilizamos el coeficiente de Spearman para explorar las relaciones monotónicas entre variables, dado que el HDI no sigue una distribución normal:
round(cor(na.omit(df[c("suicides_no","population","hdi_for_year","gdp_for_year","gdp_for_year_log","gdp_per_capita","gdp_per_capita_log")]), method = "spearman"),3)## suicides_no population hdi_for_year gdp_for_year
## suicides_no 1.000 0.762 0.185 0.620
## population 0.762 1.000 0.133 0.766
## hdi_for_year 0.185 0.133 1.000 0.612
## gdp_for_year 0.620 0.766 0.612 1.000
## gdp_for_year_log 0.620 0.766 0.612 1.000
## gdp_per_capita 0.084 0.054 0.929 0.591
## gdp_per_capita_log 0.084 0.054 0.929 0.591
## gdp_for_year_log gdp_per_capita gdp_per_capita_log
## suicides_no 0.620 0.084 0.084
## population 0.766 0.054 0.054
## hdi_for_year 0.612 0.929 0.929
## gdp_for_year 1.000 0.591 0.591
## gdp_for_year_log 1.000 0.591 0.591
## gdp_per_capita 0.591 1.000 1.000
## gdp_per_capita_log 0.591 1.000 1.000Como se esperaba, las variables con transformaciones logarítmicas
muestran una relación perfecta con sus respectivas variables originales.
Sin embargo, se observa que el HDI (hdi_for_year) está
moderadamente relacionado con gdp_for_year y fuertemente
relacionado con gdp_per_capita.
# Establece la disposición de 2x2 para los gráficos
par(mfrow = c(2, 2),
mar = c(4, 4, 2, 1)) # Márgenes ajustados para dar más espacio a las etiquetas
# Primer gráfico: GDP por year vs. HDI por year
plot(df$gdp_for_year, df$hdi_for_year,
ylab = "HDI por year", xlab = "GDP por year",
main = "GDP vs. HDI",
pch = 19, cex = 0.8)
# Segundo gráfico: GDP por year (Log) vs. HDI por year
plot(df$gdp_for_year_log, df$hdi_for_year,
ylab = "", xlab = "GDP por year (Log)",
main = "GDP (Log) vs. HDI",
pch = 19,cex = 0.8)
# Tercer gráfico: GDP per capita vs. HDI por year
plot(df$gdp_per_capita, df$hdi_for_year,
ylab = "HDI por year", xlab = "GDP per capita",
main = "GDP per capita vs. HDI",
pch = 19,cex = 0.8)
# Cuarto gráfico: GDP per capita (Log) vs. HDI por year
plot(df$gdp_per_capita_log, df$hdi_for_year,
ylab = "", xlab = "GDP per capita (Log)",
main = "GDP per capita (Log) vs. HDI",
pch = 19,cex = 0.8)
Los gráficos de dispersión a continuación sugieren que las
transformaciones logarítmicas de las variables
gdp_per_capita y gdp_for_year tienen
relaciones lineales considerables con hdi_for_year. Por lo
tanto, procedemos a ajustar un modelo de regresión lineal utilizando
estas dos variables como predictores del HDI:
model <- lm(hdi_for_year~gdp_per_capita_log+gdp_for_year_log, data = df)
summary(model)##
## Call:
## lm(formula = hdi_for_year ~ gdp_per_capita_log + gdp_for_year_log,
## data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.131865 -0.023307 0.004735 0.026218 0.109136
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.0543784 0.0046271 11.75 <2e-16 ***
## gdp_per_capita_log 0.0655451 0.0004095 160.05 <2e-16 ***
## gdp_for_year_log 0.0044041 0.0002238 19.68 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.03724 on 8361 degrees of freedom
## (19456 observations deleted due to missingness)
## Multiple R-squared: 0.8409, Adjusted R-squared: 0.8409
## F-statistic: 2.21e+04 on 2 and 8361 DF, p-value: < 2.2e-16model2 <- lm(hdi_for_year~gdp_per_capita_log, data = df)
summary(model2)##
## Call:
## lm(formula = hdi_for_year ~ gdp_per_capita_log, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.129733 -0.026388 0.004922 0.027833 0.110571
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.1209166 0.0032310 37.42 <2e-16 ***
## gdp_per_capita_log 0.0701741 0.0003429 204.64 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.03809 on 8362 degrees of freedom
## (19456 observations deleted due to missingness)
## Multiple R-squared: 0.8336, Adjusted R-squared: 0.8335
## F-statistic: 4.188e+04 on 1 and 8362 DF, p-value: < 2.2e-16Ambos modelos de regresión explican una parte considerable de la
varianza observada en el HDI por year. La exclusión de
gdp_for_year representa una reducción marginal en la
varianza explicada. Sin embargo, ambos modelos violan los supuestos de
homocedasticidad y normalidad de los residuos.
car::vif(model)## gdp_per_capita_log gdp_for_year_log
## 1.492244 1.492244lmtest::bptest(model)##
## studentized Breusch-Pagan test
##
## data: model
## BP = 770.8, df = 2, p-value < 2.2e-16hist(residuals(model))
nortest::ad.test(residuals(model))##
## Anderson-Darling normality test
##
## data: residuals(model)
## A = 59.082, p-value < 2.2e-16lmtest::bptest(model2)##
## studentized Breusch-Pagan test
##
## data: model2
## BP = 574.2, df = 1, p-value < 2.2e-16hist(residuals(model2))
nortest::ad.test(residuals(model2))##
## Anderson-Darling normality test
##
## data: residuals(model2)
## A = 50.936, p-value < 2.2e-16A partir del segundo modelo, a continuación se presente el grafico de
dispersión entre gdp_per_capita_log y
hdi_for_year, con su respectiva recta de regresión:
plot(df$gdp_per_capita_log,df$hdi_for_year, ylab = "HDI por year", xlab = "GDP per capita",pch = 19,cex = 0.6)
abline(model2, col = "red")
Para evaluar la capacidad predictiva del modelo, se realiza una
validación cruzada con 10 submuestras utlizando la libreria
caret:
train_control <- trainControl(method = "cv", number = 10)
model_cv <- train(hdi_for_year~gdp_per_capita_log, data = na.omit(df), method = "lm", trControl = train_control)
print(model_cv)## Linear Regression
##
## 8364 samples
## 1 predictor
##
## No pre-processing
## Resampling: Cross-Validated (10 fold)
## Summary of sample sizes: 7529, 7527, 7527, 7526, 7528, 7528, ...
## Resampling results:
##
## RMSE Rsquared MAE
## 0.03808609 0.8336506 0.03091247
##
## Tuning parameter 'intercept' was held constant at a value of TRUEDado que los resultados indican que el modelo tiene una capacidad predictiva considerable, y considerando las violaciones a los supuestos de la regresión lineal, se justifica utilizar las estimaciones del modelo para imputar los datos faltantes del HDI:
# Identificar los índices con datos faltantes
missing_indices <- which(is.na(df$hdi_for_year))
# Imputar los datos faltantes
df$hdi_for_year[missing_indices] <- predict(model2, newdata = df[missing_indices,])Finalmente, confirmamos que la base de datos ya no presenta datos faltantes:
missmap(df)## Warning: Unknown or uninitialised column: `arguments`.
## Unknown or uninitialised column: `arguments`.## Warning: Unknown or uninitialised column: `imputations`.
Actividad Practica # 4
Introducción
En esta sección, se filtran los datos para obtener dos conjuntos de
datos específicos: uno para Colombia y otro para los Estados Unidos. Se
crean dos bases de datos llamadas master_col y
master_eu, que contienen solo los datos correspondientes a
cada país.
master_col <- df[df$country == "Colombia",]
master_col$country <- droplevels(master_col$country)
master_eu <- df[df$country == "United States",]
master_eu$country <- droplevels(master_eu$country)
master_col_eu <- rbind(master_col,master_eu)
table(master_col_eu$country)##
## Colombia United States
## 372 372Aquí, se realiza un análisis de la evolución de varios indicadores (suicidios por cada 100,000 habitantes, PIB per cápita e IDH) a lo largo de los años en ambos países. El cálculo se realiza agrupando los datos por año y país, y luego se calcula la tasa de suicidios ajustada por la población.
Suicides_years <- master_col_eu %>%
group_by(year, country) %>%
summarise(total_suicides = sum(suicides_no), population = sum(population), gdp_per_capita = mean(gdp_per_capita), hdi_for_year = mean(hdi_for_year))## `summarise()` has grouped output by 'year'. You can override using the
## `.groups` argument.Suicides_years <- Suicides_years %>%
mutate(suicidios_por_100k = (total_suicides / population) * 100000)
colnames(Suicides_years) <- c("year","Pais","Total_suicidios","Num_habitantes","PIB_per_capita","hdi_por_year","Suicidios_por_1000")
Suicides_years$year <- as.numeric(as.character(Suicides_years$year))Indicadores
Tasa de Suicidios
ggplot(Suicides_years, aes(x = year, y = Suicidios_por_1000,color = Pais, group = Pais)) +
geom_line(size = 1.2) + # Dibuja las líneas
geom_point(size = 2) + # Añade los puntos
labs(title = "Tasa de Suicidio por Cada 100,000 Habitantes (1985-2016)",
x = "year",
y = "Tasa de Suicidio (por 100,000 Habitantes)") +
scale_x_continuous(breaks = seq(min(Suicides_years$year), max(Suicides_years$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
PIB per Capita
ggplot(Suicides_years, aes(x = year, y = PIB_per_capita ,color = Pais, group = Pais)) +
geom_line(size = 1.2) + # Dibuja las líneas
geom_point(size = 2) + # Añade los puntos
labs(title = "PIB per capita de Colombia y Estados Unidos (1985-2016)",
x = "year",
y = "PIB per capita") +
scale_x_continuous(breaks = seq(min(Suicides_years$year), max(Suicides_years$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
IDH
ggplot(Suicides_years, aes(x = year, y = hdi_por_year ,color = Pais, group = Pais)) +
geom_line(size = 1.2) + # Dibuja las líneas
geom_point(size = 2) + # Añade los puntos
labs(title = "HID de Colombia y Estados Unidos (1985-2016)",
x = "year",
y = "HID por year") +
scale_x_continuous(breaks = seq(min(Suicides_years$year), max(Suicides_years$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
Análisis por Género
En esta sección, se examina cómo la evolución de los suicidios, el PIB per cápita y el IDH varía por género en cada país. Los datos se agrupan por año, género y país, y se calcula la tasa de suicidios por cada 100,000 habitantes.
sucides_year_sex <- master_col_eu%>%
group_by(year,sex,country)%>%
summarise(total_suicides = sum(suicides_no), population = sum(population),dp_per_capita = sum(gdp_per_capita), hdi_for_year = sum(hdi_for_year))## `summarise()` has grouped output by 'year', 'sex'. You can override using the
## `.groups` argument.sucides_year_sex <- sucides_year_sex %>%
mutate(suicidios_por_100k = (total_suicides / population) * 100000)
sucides_year_sex$sex <- str_replace(sucides_year_sex$sex, "female", "Mujeres")
sucides_year_sex$sex <- str_replace(sucides_year_sex$sex, "male", "Hombres")
colnames(sucides_year_sex) <- c("year","Sexo","Pais","Total_suicidios","Num_habitantes","PIB_per_capita","hdi_por_year","Suicidios_por_1000")
sucides_year_sex$Pais_Sexo <- paste(sucides_year_sex$Pais, sucides_year_sex$Sexo, sep = "-")Ahora se presentaran los gráficos para cada indicador, diferenciados por género.
Tasa de Suicidios por Género
ggplot(sucides_year_sex, aes(x = year, y = Suicidios_por_1000, color = Pais_Sexo, group = Pais_Sexo)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "Tasa de suicidio por cada 100,000 habitantes por país y sexo",
x = "year",
y = "Tasa de suicidio (por 100,000 habitantes)") +
scale_x_continuous(breaks = seq(min(sucides_year_sex$year), max(sucides_year_sex$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
PIB per Capita por Género
ggplot(sucides_year_sex, aes(x = year, y = PIB_per_capita ,color = Pais_Sexo, group = Pais_Sexo)) +
geom_line(size = 1.2) +
geom_point(size = 2) +
labs(title = "PIB per capita por sexo (1985-2016)",
x = "year",
y = "PIB per capita") +
scale_x_continuous(breaks = seq(min(sucides_year_sex$year), max(sucides_year_sex$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
IDH por Género
ggplot(sucides_year_sex, aes(x = year, y = hdi_por_year ,color = Pais_Sexo, group = Pais_Sexo)) +
geom_line(size = 1.2) + # Dibuja las líneas
geom_point(size = 2) + # Añade los puntos
labs(title = "PIB per capita por sexo (1985-2016)",
x = "year",
y = "PIB per capita") +
scale_x_continuous(breaks = seq(min(sucides_year_sex$year), max(sucides_year_sex$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
Análisis por Grupo de edad
Primero, se agrupan los datos por año, grupo de edad y país. Se calculan los suicidios totales, la población total, el PIB per cápita y el IDH para cada combinación de estos factores. Luego, se calcula la tasa de suicidios por cada 100,000 habitantes.
sucides_year_age <- master_col_eu%>%
group_by(year,age,country)%>%
summarise(total_suicides = sum(suicides_no), population = sum(population),dp_per_capita = sum(gdp_per_capita), hdi_for_year = sum(hdi_for_year))## `summarise()` has grouped output by 'year', 'age'. You can override using the
## `.groups` argument.sucides_year_age <- sucides_year_age %>%
mutate(suicidios_por_100k = (total_suicides / population) * 100000)Se realiza una normalización de los nombres de los rangos de edad para garantizar la consistencia en el análisis.
sucides_year_age$age <- str_replace(sucides_year_age$age , "15-24 years", "15-24 years")
sucides_year_age$age <- str_replace(sucides_year_age$age , "25-34 years", "25-34 years")
sucides_year_age$age <- str_replace(sucides_year_age$age , "35-54 years", "35-54 years")
sucides_year_age$age <- str_replace(sucides_year_age$age , "5-14 years", "5-14 years")
sucides_year_age$age <- str_replace(sucides_year_age$age , "55-74 years", "55-74 years")
sucides_year_age$age <- str_replace(sucides_year_age$age , "75+ years", "75+ years")
colnames(sucides_year_age) <- c("year","Rango_edad","Pais","Total_suicidios","Num_habitantes","PIB_per_capita","hdi_por_year","Suicidios_por_1000")
sucides_year_age$Pais_Age <- paste(sucides_year_age$Pais, sucides_year_age$Rango_edad, sep = "-")Tasa de Suicidios por Edad
Se genera un gráfico que muestra la evolución de la tasa de suicidios por cada 100,000 habitantes, desglosado por rango de edad en Colombia. Se utiliza geom_line para trazar las líneas y geom_point para marcar los puntos en cada año.
ggplot(sucides_year_age[sucides_year_age$Pais == "Colombia",], aes(x = year, y = Suicidios_por_1000, color = Rango_edad, group = Rango_edad)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "Tasa de suicidio por cada 100,000 habitantes por rango de edad en Colombia",
x = "year",
y = "Tasa de suicidio (por 100,000 habitantes)") +
scale_x_continuous(breaks = seq(min(sucides_year_age$year), max(sucides_year_age$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
De manera similar, se genera un gráfico para los Estados Unidos que muestra la tasa de suicidios por grupo de edad.
ggplot(sucides_year_age[sucides_year_age$Pais == "United States",], aes(x = year, y = Suicidios_por_1000, color = Rango_edad, group = Rango_edad)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "Tasa de suicidio por cada 100,000 habitantes por rango de edad en los US",
x = "year",
y = "Tasa de suicidio (por 100,000 habitantes)") +
scale_x_continuous(breaks = seq(min(sucides_year_age$year), max(sucides_year_age$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
PIB per Cápita por Edad
Se crea un gráfico que muestra la evolución del PIB per cápita por grupo de edad en Colombia. Similar al gráfico de suicidios, se usa geom_line y geom_point.
ggplot(sucides_year_age[sucides_year_age$Pais == "Colombia",], aes(x = year, y = PIB_per_capita, color = Rango_edad, group = Rango_edad)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "PIB per capita de Colombia por sexo (1985-2016)",
x = "year",
y = "PIB per capita") +
scale_x_continuous(breaks = seq(min(sucides_year_age$year), max(sucides_year_age$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
El gráfico correspondiente para Estados Unidos se genera de la misma manera.
ggplot(sucides_year_age[sucides_year_age$Pais == "United States",], aes(x = year, y = PIB_per_capita, color = Rango_edad, group = Rango_edad)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "PIB per capita de US por sexo (1985-2016)",
x = "year",
y = "PIB per capita") +
scale_x_continuous(breaks = seq(min(sucides_year_age$year), max(sucides_year_age$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
IDH por Edad
Finalmente, se crea un gráfico para mostrar la evolución del IDH por grupo de edad en Colombia.
ggplot(sucides_year_age[sucides_year_age$Pais == "Colombia",], aes(x = year, y = hdi_por_year, color = Rango_edad, group = Rango_edad)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "HDI por year de Colombia por rango de edad (1985-2016)",
x = "year",
y = "HDI por year") +
scale_x_continuous(breaks = seq(min(sucides_year_age$year), max(sucides_year_age$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
Se genera el gráfico equivalente para Estados Unidos.
ggplot(sucides_year_age[sucides_year_age$Pais == "United States",], aes(x = year, y = hdi_por_year, color = Rango_edad, group = Rango_edad)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "HDI por year de US por rango de edad (1985-2016)",
x = "year",
y = "HDI por year") +
scale_x_continuous(breaks = seq(min(sucides_year_age$year), max(sucides_year_age$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
Actividad práctica # 5
Durante esta actividad, se trabajara con una base titulada
Accidentalidad en Barranquilla la cual fue suministrada por
la Alcaldía
Distrital de Barranquilla y contiene información sobre accidentes de
tránsito ocurridos en el distrito de Barranquilla, obtenida de los
informes policiales de accidentes de tránsito (IPAT). Esta base de datos
incluye registros de accidentes que han sucedido en el área desde el
year 2018 hasta el 2024. Es importante destacar que la información de la
vigencia actual es preliminar y está sujeta a cambios, lo que implica
que los datos pueden ser actualizados o corregidos en el futuro. La base
de datos proporciona un registro detallado de cada accidente,
permitiendo realizar análisis sobre la frecuencia, severidad y ubicación
de los accidentes en la ciudad.
df <- read_csv("C:/Users/jcami/OneDrive/Escritorio/DATAVIZ/Accidentalidad_en_Barranquilla_20240826.csv", col_names = TRUE)
df %<>% clean_namesEl estudio preliminar de la base de datos nos muestra contiene 11 varibales y 25610 observaciones. El estudio superfical de las varibales indica que alguna de ellas necesitaran ser normalizadas, renomabradas, o cambiar su formato para facilitar los analisis posteriores.
dim(df)## [1] 25610 11knitr::kable(descr_df)| Nombre_original | Nombre_recodificado | Tipo_original | Tipo_modificado |
|---|---|---|---|
| fecha_accidente | Fecha | Fecha | Fecha |
| hora_accidente | Hora | Cadena | Hora |
| gravedad_accidente | Gravedad | Cadena | Factor |
| clase_accidente | Clase | Cadena | Factor |
| sitio_exacto_accidente | Sitio | Cadena | Cadena |
| cant_heridos_en_sitio_accidente | Heridos | Númerico | Númerico |
| cant_muertos_en_sitio_accidente | Muertos | Númerico | Númerico |
| cantidad_accidentes | Cantidad | Númerico | Númerico |
| ano_accidente | year | Númerico | Númerico |
| mes_accidente | Mes | Cadena | Factor |
| dia_accidente | Dia | Cadena | Factor |
Con tal de examinar de manera los datos de almacenados en cada una de
las variables y detarminar la presencia de datos faltantes se utiliza la
función Summary de base en R:
summary(df)## fecha_accidente hora_accidente gravedad_accidente
## Min. :2018-01-01 00:00:00.00 Length:25610 Length:25610
## 1st Qu.:2019-02-02 00:00:00.00 Class :character Class :character
## Median :2020-04-23 12:00:00.00 Mode :character Mode :character
## Mean :2020-07-31 19:57:36.05
## 3rd Qu.:2021-12-13 00:00:00.00
## Max. :2024-06-30 00:00:00.00
##
## clase_accidente sitio_exacto_accidente cant_heridos_en_sitio_accidente
## Length:25610 Length:25610 Min. : 1.000
## Class :character Class :character 1st Qu.: 1.000
## Mode :character Mode :character Median : 1.000
## Mean : 1.472
## 3rd Qu.: 2.000
## Max. :42.000
## NA's :15626
## cant_muertos_en_sitio_accidente cantidad_accidentes ano_accidente
## Min. :1.000 Min. :1 Min. :2018
## 1st Qu.:1.000 1st Qu.:1 1st Qu.:2019
## Median :1.000 Median :1 Median :2020
## Mean :1.036 Mean :1 Mean :2020
## 3rd Qu.:1.000 3rd Qu.:1 3rd Qu.:2021
## Max. :2.000 Max. :2 Max. :2024
## NA's :25358
## mes_accidente dia_accidente
## Length:25610 Length:25610
## Class :character Class :character
## Mode :character Mode :character
##
##
##
## Tal como se puede observar, no es posible detarminar la presencia de
datos faltantes o el número de observaciones por categoria en la
varibales tipo carácter por lo que más adelante se procede a cambiar el
tipo de dato y analizarlas más a detalla. Dicho esto, se puede observar
una cantidad considerable de datos faltantes en las varibales
cant_muertos_en_sitio_accidente y
cant_heridos_en_sitio_accidente, por lo cual se planifica
examinar diversos métodos de imputación de datos. Habiendo completado el
estudio superfical de la base de datos, ahora se estudiaran una a una
las varibales implicadas. No sin antes modificar el nombre de las
columnas y cambiar el tipo de datos de algunas varibales:
colnames(df) <- c("Fecha","Hora","Gravedad","Clase","Sitio","Heridos","Muertos","Cantidad","year","Mes","Dia")
selected_columns <- c("Gravedad","Clase","Mes","Dia")
df[selected_columns] <- lapply(df[selected_columns], as.factor)
df$Fecha <- as.Date(df$Fecha)Estudio de las varibales
Fecha
En primera instancia, la variable Fecha no muestra
problemas relacionados con la digitación o los datos faltantes:
knitr::kable(head(df[1:5],5))| Fecha | Hora | Gravedad | Clase | Sitio |
|---|---|---|---|---|
| 2018-01-01 | 01:30:00:am | Con heridos | Atropello | CL 87 9H 24 |
| 2018-01-01 | 02:00:00:pm | Solo daños | Choque | CL 110 CR 46 |
| 2018-01-01 | 04:00:00:am | Solo daños | Choque | AV CIRCUNVALAR CR 9G |
| 2018-01-01 | 04:30:00:am | Solo daños | Choque | CLLE 72 CRA 29 |
| 2018-01-01 | 05:20:00:pm | Solo daños | Choque | VIA 40 CALLE 75 |
which(is.na(as.factor(df$Fecha)))## integer(0)Se concluira visualizando al variable por medio de un histograma:
ggplot(df, aes(x = Fecha)) +
geom_histogram(fill = "cyan4", color = "black") +
theme_minimal() +
labs(
title = "Histograma de los accidentes en Barranquilla por fecha",
x = "Fecha",
y = "Frecuencia de accidentes"
) +
theme(
plot.title = element_text(hjust = 0.5), # Centrar el título
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
Hora
En la vector Hora se almacenan las hora en que
sucedieron los accidentes. Comenzaremos cambiando el tipo de la variable
a hora utilizando el paquete lubridate. Para eso, debemos
normalizar el contenido de la varibale y suprimir caracteres que pueden
dificultar la conversión:
hora_modificada <- gsub(":(am|pm)$", " \\1", df$Hora) # Parsear la hora utilizando el formato correcto
hora_modificada <- parse_date_time(hora_modificada, orders = "I:M:S p") # Extraer solo la parte "HH:MM:SS" en formato de 24 horas
hora_modificada <- format(hora_modificada, "%H:%M:%S")
# Convertir a hms
df$Hora <- hms(hora_modificada)
rm(hora_modificada)Hecho esto podemos asegurarnos que no hayan datos faltantes y calcular estadisticas descriptivas:
which(is.na(df$Hora))## integer(0)summary(df[2])## Hora
## Min. :0S
## 1st Qu.:9H 30M 0S
## Median :13H 30M 0S
## Mean :13H 26M 30.5154236626331S
## 3rd Qu.:17H 30M 0S
## Max. :23H 58M 0SA <- ggplot(df, aes(x = as.numeric(Hora))) +
geom_histogram(binwidth = 3600, fill = "cyan4", color = "black") +
scale_x_continuous(labels = function(x) format(hms::as_hms(x), "%H:%M:%S")) +
labs(title = "Histograma de Horas", x = "Hora del Dia", y = "Frecuencia") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5), # Centrar el título
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15)),
plot.margin = margin(10, 20, 10, 10)
)
B <- ggplot(df, aes(x = "", y = as.numeric(Hora))) +
geom_boxplot(fill = "cyan4", color = "black") +
scale_y_continuous(labels = function(x) format(hms::as_hms(x), "%H:%M:%S")) +
labs(title = "Boxplot de Horas", x = "", y = "Hora del Dia") +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5), # Centrar el título
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15)),
plot.margin = margin(10, 10, 10, 20)
)
grid.arrange(A, B, ncol = 2)
Gravedad
Tal como su nombre lo indica la variable Gravedad indica
la gravedad del accidente. El análisis de esta variable iniciara
construyendo una tabla de frecuencia que nos ayudara a determinar el
número de observaciones en cada condición y la posible existencia de
datos faltantes.
table(df$Gravedad)##
## Con heridos Con muertos Solo daños
## 9901 252 15457Debido a que la tabla de frencuencia muestra que no hay datos faltantes en esta variable, procederemos a realizar una tabla con un conteo y el porcentaje de accidentes según su clase:
Tabla_gravedad <- df %>%
group_by(Gravedad) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_gravedad## # A tibble: 3 × 4
## Gravedad n Proporcion Porcentaje
## <fct> <int> <dbl> <chr>
## 1 Con muertos 252 0.00984 1%
## 2 Con heridos 9901 0.387 39%
## 3 Solo daños 15457 0.604 60%Podemos utilizar esta tabla para crear un diagrama de barras que nos permita visualizar la cantidad de accidentes segun su gravedad:
Tabla_gravedad$Gravedad <- factor(Tabla_gravedad$Gravedad,
levels = Tabla_gravedad$Gravedad[order(Tabla_gravedad$n, decreasing = TRUE)])
ggplot(Tabla_gravedad, aes(x = Gravedad, y = n)) +
geom_bar(width = 0.7, stat = "identity", position = position_dodge(0.7), fill = "cyan4") +
coord_cartesian(ylim = c(0, 20000)) +
ggtitle("Accidentes en Barranquilla segun su gravedad (2018-2024)") +
labs(x = "Gravedad", y = "Frecuencias") +
geom_text(aes(label = Porcentaje),
vjust = -0.5,
color = "black",
hjust = 0.5,
position = position_dodge(0.7),
angle = 0,
size = 4.5) +
scale_y_continuous(labels = scales::label_number(scale = 1e-3, suffix = "K")) +
theme(
axis.text.x = element_text(angle = 0, vjust = 1, hjust = 0.5),
axis.title.x = element_text(margin = margin(t = 20)),
axis.title.y = element_text(margin = margin(r = 20)),
plot.title = element_text(size = 14,vjust = 2,hjust = 0.5)
)+
theme_minimal()
Clase de accidente
La varibale Clase indica la clase de accidente. Tal como
en el caso anterior, construimos una tabla de frecuencias para verificar
el número de observaciones por condición y se cuentan con datos
faltantes:
table(df$Clase)##
## Atropello Caida Ocupante Choque Incendio Otro
## 1344 194 23819 13 123
## Volcamiento
## 117Debido a que la tabla de frencuencia muestra que no hay datos faltantes en esta variable, procederemos a realizar una tabla con un conteo y el porcentaje de accidentes según su tipo:
df$Clase <- str_replace(df$Clase, "Caida Ocupante", "Caida")
Tabla_Clase <- df %>%
group_by(Clase) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_Clase## # A tibble: 6 × 4
## Clase n Proporcion Porcentaje
## <chr> <int> <dbl> <chr>
## 1 Incendio 13 0.000508 0.051%
## 2 Volcamiento 117 0.00457 0.457%
## 3 Otro 123 0.00480 0.480%
## 4 Caida 194 0.00758 0.758%
## 5 Atropello 1344 0.0525 5.248%
## 6 Choque 23819 0.930 93.007%Podemos utilizar esta tabla para crear un diagrama de barras que nos permita visualizar la cantidade de accidentes segun su clase:
Tabla_Clase$Clase <- factor(Tabla_Clase$Clase,
levels = Tabla_Clase$Clase[order(Tabla_Clase$n, decreasing = TRUE)])
ggplot(Tabla_Clase, aes(x = Clase, y = n)) +
geom_bar(width = 0.7, stat = "identity", position = position_dodge(0.7), fill = "cyan4") +
coord_cartesian(ylim = c(0, 30000)) +
ggtitle("Accidentes en Barranquilla segun su clase (2018-2014)") +
labs(x = "Clase", y = "Frecuencias") +
geom_text(aes(label = Porcentaje),
vjust = -0.5,
color = "black",
hjust = 0.5,
position = position_dodge(0.7),
angle = 0,
size = 4.5) +
scale_y_continuous(labels = scales::label_number(scale = 1e-3, suffix = "K")) +
theme(
axis.text.x = element_text(angle = 0, vjust = 1, hjust = 0.5),
axis.title.x = element_text(margin = margin(t = 20)),
axis.title.y = element_text(margin = margin(r = 20)),
plot.title = element_text(size = 14,vjust = 2,hjust = 0.5)
)+
theme_minimal()
Year
En lo que concierne al year del accidente (year)
procederemos a manera similar. Primer construiremos una tabla de
frecuencia para comprobar datos faltantes, calcularemos estadisticos
descriptivos, y finalizaramos vizualizando el número de accidentes a lo
largo de los years.
table(df$year)##
## 2018 2019 2020 2021 2022 2023 2024
## 5898 5645 3281 4700 3683 1662 741Tal que no se encuentran datos faltantes, se procede a construir la tabla de frecuencia con estadisticas descriptivas:
Tabla_year <- df %>%
group_by(year) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_year## # A tibble: 7 × 4
## year n Proporcion Porcentaje
## <dbl> <int> <dbl> <chr>
## 1 2024 741 0.0289 2.89%
## 2 2023 1662 0.0649 6.49%
## 3 2020 3281 0.128 12.81%
## 4 2022 3683 0.144 14.38%
## 5 2021 4700 0.184 18.35%
## 6 2019 5645 0.220 22.04%
## 7 2018 5898 0.230 23.03%Culminamos el estudio de esta varible construyendo su diagrama de barras:
ggplot(Tabla_year, aes(x = as.factor(year), y = n)) +
geom_bar(width = 0.7, stat = "identity", position = position_dodge(1), fill = "cyan4") +
coord_cartesian(ylim = c(0, 10000)) +
ggtitle("Accidentes en Barranquilla (2018-2014)") +
labs(x = "Year", y = "Frecuencias") +
geom_text(aes(label = Porcentaje),
vjust = -0.5,
color = "black",
hjust = 0.5,
position = position_dodge(1),
angle = 0,
size = 4.5) +
scale_y_continuous(labels = scales::label_number(scale = 1e-3, suffix = "K")) +
theme(
axis.text.x = element_text(angle = 0, vjust = 1, hjust = 0.5),
axis.title.x = element_text(margin = margin(t = 20)),
axis.title.y = element_text(margin = margin(r = 20)),
plot.title = element_text(size = 14,vjust = 2,hjust = 0.5)
)+
theme_minimal()
Tambien podemos construir un grafico de lineas que nos permite vizualizar más clara el número de accidentes en la ciudad a lo largo de los years:
ggplot(Tabla_year, aes(x = year, y = n)) +
geom_line(size = 1) +
geom_point(size = 2) +
labs(title = "Numero de accidentes en Barranquilla (2018-2024)",
x = "Year",
y = "Numero de accidentes") +
scale_x_continuous(breaks = seq(min(Tabla_year$year), max(Tabla_year$year), by = 1)) +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1),
axis.title.x = element_text(margin = margin(t = 15)),
axis.title.y = element_text(margin = margin(r = 15))
)
Cantidad de accidentes
Según la información suministrada junto a la base de datos, la
variable Cantidad indica la cantidad de accidentes en un
suceso particular. Debido a que es un conteo procedemos de manera
similar a las secciones anteriores, comenzado por construir una tabla de
frecuencia.
table(df$Cantidad)##
## 1 2
## 25605 5La tabla nos indica que la varibale no presenta datos faltantes y que la la modo es uno. Procedemos a construir una tabla con la proporción y porcentaje de observaciones en cada categoria.
Tabla_cantidad <- df %>%
group_by(Cantidad) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_cantidad## # A tibble: 2 × 4
## Cantidad n Proporcion Porcentaje
## <dbl> <int> <dbl> <chr>
## 1 2 5 0.000195 0%
## 2 1 25605 1.00 100%De manera inmediata que solo cinco observaciones tienen una cantidad de dos, lo cual representa alrededor del 0.01% de las observaciones totales. Es probable que estas observaciones hayan sido el resultado de un error durante la digitación.
knitr::kable(df[which(df$Cantidad == 2),][c(1,2,3,4,6,7,8,9)])| Fecha | Hora | Gravedad | Clase | Heridos | Muertos | Cantidad | year |
|---|---|---|---|---|---|---|---|
| 2018-06-09 | 9H 35M 0S | Solo daños | Choque | NA | NA | 2 | 2018 |
| 2019-12-28 | 12H 28M 0S | Solo daños | Choque | NA | NA | 2 | 2019 |
| 2020-02-06 | 9H 0M 0S | Solo daños | Choque | NA | NA | 2 | 2020 |
| 2020-05-02 | 19H 40M 0S | Solo daños | Choque | NA | NA | 2 | 2020 |
| 2020-11-05 | 17H 25M 0S | Solo daños | Choque | NA | NA | 2 | 2020 |
Dicho esto, se considera que reemplazar dichas observaciones por la
moda no es lo adecuado. Esto es debido a que el resultado seria una
variables sin variación que en principio poco aportaria al ejercicio.
Con esto en mente, se procede a eliminar la variable
Cantidad de la base de datos.
df <- df[-8]Mes
Procedemos de la misma manera para estudiar los meses en donde
ocurrienron accidentes (Mes):
table(df$Mes)##
## April August December February January July June March
## 2010 1918 2189 2477 2349 1932 2103 2446
## May November October September
## 2121 1995 2090 1980Tal que no se identifican datos faltantes se procede a construir la tabla de frecuencias.
Tabla_mes <- df %>%
group_by(Mes) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_mes## # A tibble: 12 × 4
## Mes n Proporcion Porcentaje
## <fct> <int> <dbl> <chr>
## 1 August 1918 0.0749 7.489%
## 2 July 1932 0.0754 7.544%
## 3 September 1980 0.0773 7.731%
## 4 November 1995 0.0779 7.790%
## 5 April 2010 0.0785 7.848%
## 6 October 2090 0.0816 8.161%
## 7 June 2103 0.0821 8.212%
## 8 May 2121 0.0828 8.282%
## 9 December 2189 0.0855 8.547%
## 10 January 2349 0.0917 9.172%
## 11 March 2446 0.0955 9.551%
## 12 February 2477 0.0967 9.672%Dia
En lo que respecta a la varibale Dia utilizamos el mismo
procedimiento anterior. Comenzamos inspeccionando la cantidad de
observaciones por condición:
table(df$Dia)##
## Fri Mon Sat Sun Thu Tue Wed
## 3920 3774 3735 2577 3756 4009 3839Debido a que no se encontraron datos faltantes, se procede a construir una tabla de frecuencia con estadisticos descriptivos:
Tabla_dia <- df %>%
group_by(Dia) %>%
count() %>%
ungroup() %>%
mutate(Proporcion = `n` / sum(`n`)) %>%
arrange(Proporcion) %>%
mutate(Porcentaje = scales::percent(Proporcion))
Tabla_dia## # A tibble: 7 × 4
## Dia n Proporcion Porcentaje
## <fct> <int> <dbl> <chr>
## 1 Sun 2577 0.101 10.062%
## 2 Sat 3735 0.146 14.584%
## 3 Thu 3756 0.147 14.666%
## 4 Mon 3774 0.147 14.736%
## 5 Wed 3839 0.150 14.990%
## 6 Fri 3920 0.153 15.307%
## 7 Tue 4009 0.157 15.654%Sitio
Finalmente la variable sitio contine información de los
lugares en donde ocurrieron accidentes. Sin embargo, tal como se puede
observar a continuación, se suministran las direcciones de los lugares
sin datos de geolocalización que permitan señalar de manera precisisa el
lugar del suceso. En este caso, se considera inviable calcular medidas
descriptivas o estudiar la frecuencia de los datos de manera usual ya
que se encuentran errores de codificación y una gran cantidad de
categoiras.
knitr::kable(tail(df$Sitio,10))| x |
|---|
| CALLE 111 CON CARRERA 37 |
| CALLE 19 CRA 6B |
| CALLE 30 CRA 14 |
| CALLE 30 CON CRA 13C |
| CARRERA 5C CALLE 90 |
| CARRERA 42G CALLE 90 |
| CALLE 116 CARRERA 38 |
| CALLE 56 FRENTE AL # 8E - 48 |
| CARRERA 43 CON CALLE 41 |
| CARRERA 8 CALLE 49 |
Debido a esto, esta varibale se describira apartir de mapas estaticos
que permitiran vizualizar la cantidad de accindentes de manera más
clara. Para esto utilizaremos la libreria de ggmap y las
API de google Could. Baste de decir que debido a gran volumen de datos
obtención de la coordenadas se realizo de manera automatica con las
siguientes funciones:
limpiar_direccion <- function(direccion) {
direccion <- str_to_upper(direccion) # Convertir a mayúsculas
direccion <- str_replace_all(direccion, "\\s+", " ") # Normalizar espacios
direccion <- str_replace_all(direccion, "\\bCR\\b", "CARRERA")
direccion <- str_replace_all(direccion, "\\bCL\\b", "CALLE")
direccion <- str_replace_all(direccion, "\\bAV\\b", "AVENIDA")
direccion <- str_replace_all(direccion, "\\bVIA\\b", "VIA")
direccion <- str_replace_all(direccion, "\\bCRA\\b", "CARRERA")
direccion <- str_replace_all(direccion, "\\s*CON\\s*", " CON ")
direccion <- str_replace_all(direccion, "\\s*#\\s*", " # ")
direccion <- str_trim(direccion)
+
direccion <- paste(direccion, "BARRANQUILLA")
return(direccion)
}
direcciones <- limpiar_direccion(df$Sitio) # Direcciones normalizadas# Inicializar el contador global para iteraciones
contador <- 0
# Función para geocodificar una dirección y manejar errores
geocodificar_direccion <- function(direccion) {
tryCatch({
resultado <- suppressWarnings(suppressMessages(geocode(direccion, source = "google")))
if (is.na(resultado$lat) | is.na(resultado$lon)) {
return(data.frame(lat = NA, lon = NA))
}
return(resultado)
}, error = function(e) {
return(data.frame(lat = NA, lon = NA))
})
}
resultados <- lapply(seq_along(direcciones), function(i) {
contador <<- contador + 1
cat("Procesando iteración:", contador, "de", length(direcciones), "\n")
resultado <- geocodificar_direccion(direcciones[i])
return(resultado)
})
coordenadas_df <- do.call(rbind, resultados)Una vez culminada la geolocalización procedemos a revisar si existen datos faltantes y constuir una base de datos que utilizaremos para crear los mapas:
coordenadas_df %>%
summarize(
Lon_NA = sum(is.na(lon)),
LAt_NA = sum(is.na(lat))
)
# Unir con las direcciones originales
resultados <- cbind(direcciones, coordenadas_df)Apesar de la normalización, aun se encontraron algunoes errores el
geodecodificación. Por lo cual, se procedio revisar las coordenadas más
extremas y determinar la latitud y longitud aproximada de los lugares
utlizando Google Maps:
coordenadas <- data.frame(
indice = c(1823, 4322, 6332, 6935, 10912, 22409, 20209, 13749, 16066, 16255,
17410, 17527, 19494, 18426, 21651, 17980, 18341, 22494, 16299,
13132, 15642, 19094, 15189, 16089, 7393, 15695, 21454, 2836,
7719, 17800, 6122),
lon = c(-74.795060, -74.843081, -74.812223, -74.791578, -74.800250, -74.805110,
-74.821869, -74.802105, -74.836010, -74.770386, -74.765775, -74.844272,
-74.835221, -74.803440, -74.778708, -74.787420, -74.835940, -74.774765,
-74.836026, -74.781519, -74.835707, -74.891751, -74.810851, -74.835752,
-74.761275, -74.761275, -74.761275, -74.773021, -74.769029, -74.904925,
NA),
lat = c(10.918997, 11.015256, 10.922867, 10.989965, 11.022345, 10.923354,
10.925821, 10.960233, 10.981524, 10.985015, 10.969660, 10.997190,
10.978436, 10.922918, 10.987049, 10.963814, 10.980279, 10.972310,
10.980032, 10.976462, 11.018203, 10.956987, 10.951522, 10.982615,
10.950376, 10.950376, 10.950376, 10.955540, 10.952123, 10.951866,
NA)
)
for (i in 1:nrow(coordenadas)) {
resultados$lon[coordenadas$indice[i]] <- coordenadas$lon[i]
resultados$lat[coordenadas$indice[i]] <- coordenadas$lat[i]
}knitr::kable(head(resultados,5))| direcciones | lon | lat |
|---|---|---|
| CALLE 87 9H 24 BARRANQUILLA | -74.82560 | 10.95923 |
| CALLE 110 CARRERA 46 BARRANQUILLA | -74.83095 | 11.03326 |
| AVENIDA CIRCUNVALAR CARRERA 9G BARRANQUILLA | -74.83612 | 10.95671 |
| CLLE 72 CARRERA 29 BARRANQUILLA | -74.81184 | 10.97874 |
| VIA 40 CALLE 75 BARRANQUILLA | -74.79308 | 11.01285 |
Habiendo preparado los datos de geolocalización solo no queda
construir el mapa utilizando a libreria de ggmap y las API
de google Could:
myMap <- get_map(location = "Barranquilla", zoom = 12)
ggmap(myMap) +
stat_density2d(data = resultados,
aes(x = lon, y = lat, fill = ..level.., alpha = ..level..),
size = 0.01, bins = 20, geom = "polygon") +
scale_fill_gradient(low = "green", high = "red",name = "Accidentes") +
scale_alpha(guide = FALSE) +
ggtitle("Mapa de Accidentes en Barranquilla (2018-2024)") +
labs(x = "Longitud", y = "Latitud") +
theme(
plot.title = element_text(hjust = 0.5, size = 12, vjust = 2),
axis.title.x = element_text(margin = margin(t = 10)),
axis.title.y = element_text(margin = margin(r = 10)),
legend.title = element_text(hjust = 1)
)
Heridos
table(df$Heridos)##
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
## 6911 2301 469 155 57 29 14 9 6 7 6 7 2 2 1 1
## 18 19 20 21 22 23 42
## 1 1 1 1 1 1 1Ahora procedemos a examinar por posible datos faltantes:
length(which(is.na(df$Heridos)))## [1] 15626Los resultados nos muestran que hacen faltan 15626 observaciones en
la variable Heridos, alrededor de un 61% del número total
del observaciones. Procedemos a examinar cinco de ellas:
head(df[which(is.na(df$Heridos)),],5)## # A tibble: 5 × 10
## Fecha Hora Gravedad Clase Sitio Heridos Muertos year Mes Dia
## <date> <Period> <fct> <chr> <chr> <dbl> <dbl> <dbl> <fct> <fct>
## 1 2018-01-01 14H 0M 0S Solo daños Choq… CL 1… NA NA 2018 Janu… Mon
## 2 2018-01-01 4H 0M 0S Solo daños Choq… AV C… NA NA 2018 Janu… Mon
## 3 2018-01-01 4H 30M 0S Solo daños Choq… CLLE… NA NA 2018 Janu… Mon
## 4 2018-01-01 17H 20M 0S Solo daños Choq… VIA … NA NA 2018 Janu… Mon
## 5 2018-01-02 14H 30M 0S Solo daños Choq… CARR… NA NA 2018 Janu… TueLa tabla sugiere que los datos faltantes en la variable
Heridos pueden estar relacionadas con la variable
Gravedad, por lo cual se prodera a agrupar los datos según
la gravedad y calcularan algunas estadisticas descriptivas de los datos
faltantes:
df %>%
group_by(Gravedad) %>%
summarise(
total = n(),
Casos = sum(is.na(Heridos)),
Proporcion = round(mean(is.na(Heridos)),2),
Procentaje = round(mean(is.na(Heridos)) * 100,3 ))## # A tibble: 3 × 5
## Gravedad total Casos Proporcion Procentaje
## <fct> <int> <int> <dbl> <dbl>
## 1 Con heridos 9901 0 0 0
## 2 Con muertos 252 169 0.67 67.1
## 3 Solo daños 15457 15457 1 100Resulta evidente que el número de datos faltantes esta relacionado
con la gravedad del accidente. Vemos por ejemplo que no hay datos
faltantes en la categoria de “con heridos”, mientras que no datos sobre
el número de heridos en la categoria de “solo dyears”. Este tipo de
datos faltantes se podría considerar como datos “no al azar”
(MANAR - Missing Not at Random), ya que la ausencia de
datos parece depender directamente de la gravedad del accidente. Esto
implica que el tratamiento de estos datos faltantes debe hacerse con
cuidado, y probablemente no se deba imputar o ignorar sin una estrategia
que considere la estructura de los datos y su relación con la variable
Gravedad.
En este caso comenzaremos por reemplazar los NA con ceros en la
categoría de “Solo dyears”. Esto se debe a que la ausencia de datos en
la variable Heridos podría interpretarse como que no hubo
heridos, lo cual es coherente con la categoría “Solo dyears”.
df <- df %>%
mutate(Heridos = ifelse((Gravedad != "Con heridos" | Gravedad != "Con muertos" )& is.na(Heridos), 0, Heridos))
df %>%
group_by(Gravedad) %>%
summarise(
total = n(),
Casos = sum(is.na(Heridos)),
Proporcion = round(mean(is.na(Heridos)),2),
Procentaje = round(mean(is.na(Heridos)) * 100,3 ))## # A tibble: 3 × 5
## Gravedad total Casos Proporcion Procentaje
## <fct> <int> <int> <dbl> <dbl>
## 1 Con heridos 9901 0 0 0
## 2 Con muertos 252 0 0 0
## 3 Solo daños 15457 0 0 0Dado que, por el momento, no contamos con suficiente información para
tratar con los faltantes en la categoria de “Con muestros” el analisis
de la varibale Heridos se suspendera para analizar el
comportamiento de los datos faltantes en al varibale
Muertos.
Muertos
En la base de datos en cuestión, la varibale Muertos
indica el número de fallecimiento en un determinado accidente.
Procedemos de manera parecida a la variable anterior.
table(df$Muertos)##
## 1 2
## 243 9En un primer momento, observamos que solo en dos casos la cantidad de fallecidos fue de cinco, mientras que los accidentes restantes (de los cuales se tienen datos) solo resultaron en un fallecimiento.Ahora procedemos a revisar el número de datos faltantes:
length(which(is.na(df$Muertos)))## [1] 25358Tal como en la varibale anterior, los resultados nos muestran que
hacen faltan 25358 observaciones en la variable Muertos,
alrededor de un 99% del número total de observaciones. Procedemos a
examinar cinco de ellas:
head(df[which(is.na(df$Muertos)),],5)## # A tibble: 5 × 10
## Fecha Hora Gravedad Clase Sitio Heridos Muertos year Mes Dia
## <date> <Period> <fct> <chr> <chr> <dbl> <dbl> <dbl> <fct> <fct>
## 1 2018-01-01 1H 30M 0S Con herid… Atro… CL 8… 1 NA 2018 Janu… Mon
## 2 2018-01-01 14H 0M 0S Solo daños Choq… CL 1… 0 NA 2018 Janu… Mon
## 3 2018-01-01 4H 0M 0S Solo daños Choq… AV C… 0 NA 2018 Janu… Mon
## 4 2018-01-01 4H 30M 0S Solo daños Choq… CLLE… 0 NA 2018 Janu… Mon
## 5 2018-01-01 17H 20M 0S Solo daños Choq… VIA … 0 NA 2018 Janu… MonEntonces se observa que una gran proporción de los faltantes podrian estar relacionados con la gravedad del accidente. Por lo cual, se construira una tabla de frecuencia agrupada por la gravedad para estudiar la proporción de faltantes en cada condición:
df %>%
group_by(Gravedad) %>%
summarise(
total = n(),
Casos = sum(is.na(Muertos)),
Proporcion = round(mean(is.na(Muertos)),2),
Procentaje = round(mean(is.na(Muertos)) * 100,3 ))## # A tibble: 3 × 5
## Gravedad total Casos Proporcion Procentaje
## <fct> <int> <int> <dbl> <dbl>
## 1 Con heridos 9901 9901 1 100
## 2 Con muertos 252 0 0 0
## 3 Solo daños 15457 15457 1 100Una vez más, los datos faltantes parecen estar fuertemente
relacionados con la gravedad del accidente. En este caso, las categorías
“Solo dyears” y “Con heridos” no tienen datos para la variable Muertos.
Naturalmente, esto indica que los datos faltantes no son
MAR (Missing At Random) y la imputación no es una
estrategia adecuada en este contexto. Además, la estructura de los datos
faltantes sugiere que el registro de heridos y muertos está determinado
por la gravedad del accidente.
Como hipótesis, es posible que en los casos en que la gravedad del accidente fue “Con muertos”, se hayan involucrado personas que, aunque no fallecieron, resultaron heridas, por lo que se registran valores en ambas variables. Sin embargo, si el accidente solo resultó en personas heridas, tiene sentido que solo se reporte este número, y los datos sobre muertos sean faltantes. Asimismo, si el accidente solo involucró dyears materiales, ambas variables Heridos y Muertos deberían estar faltantes. Esto explicaría la estructura encontrada.
df <- df %>%
mutate(Muertos = ifelse(Gravedad != "Con muertos" & is.na(Muertos), 0, Muertos))
df %>%
group_by(Gravedad) %>%
summarise(
total = n(),
Casos = sum(is.na(Muertos)),
Proporcion = round(mean(is.na(Muertos)),2),
Procentaje = round(mean(is.na(Muertos)) * 100,3 ))## # A tibble: 3 × 5
## Gravedad total Casos Proporcion Procentaje
## <fct> <int> <int> <dbl> <dbl>
## 1 Con heridos 9901 0 0 0
## 2 Con muertos 252 0 0 0
## 3 Solo daños 15457 0 0 0summary(df)## Fecha Hora Gravedad
## Min. :2018-01-01 Min. :0S Con heridos: 9901
## 1st Qu.:2019-02-02 1st Qu.:9H 30M 0S Con muertos: 252
## Median :2020-04-23 Median :13H 30M 0S Solo daños :15457
## Mean :2020-07-31 Mean :13H 26M 30.5154236626331S
## 3rd Qu.:2021-12-13 3rd Qu.:17H 30M 0S
## Max. :2024-06-30 Max. :23H 58M 0S
##
## Clase Sitio Heridos Muertos
## Length:25610 Length:25610 Min. : 0.0000 Min. :0.00000
## Class :character Class :character 1st Qu.: 0.0000 1st Qu.:0.00000
## Mode :character Mode :character Median : 0.0000 Median :0.00000
## Mean : 0.5737 Mean :0.01019
## 3rd Qu.: 1.0000 3rd Qu.:0.00000
## Max. :42.0000 Max. :2.00000
##
## year Mes Dia
## Min. :2018 February: 2477 Fri:3920
## 1st Qu.:2019 March : 2446 Mon:3774
## Median :2020 January : 2349 Sat:3735
## Mean :2020 December: 2189 Sun:2577
## 3rd Qu.:2021 May : 2121 Thu:3756
## Max. :2024 June : 2103 Tue:4009
## (Other) :11925 Wed:3839Tratamiento de datos atipicos
En esta sección trataremos con los posibles datos atipicos de que
puden existir en las varibales Heridos y
Muertos. Para esto se seguira el siguiente protocolo 1) se
determinan si existens datos atipicos utilizando graficos, estadisiticos
descriptivos, y pruebas analiticas, 2) se reemplazaran por datos
atipicos en el caso de que sean un baja proporción, y 3) se utilizaran
multiples tecnicas de imputación.
Muertos
df%>%
get_summary_stats(Muertos, type = "five_number")## # A tibble: 1 × 7
## variable n min max q1 median q3
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Muertos 25610 0 2 0 0 0boxplot(df$Muertos)
hist(df$Muertos)
Los digramas de cajas y el histograma muestran que los valores de 1 y 2 pueden ser posibles candidatos atipicos. En particular se estarian considerando un total de 252 observaciones (aproximadamente 1% del total).
length(which(df$Muertos > 0))## [1] 252lower_bound <- median(df$Muertos) - 3 * mad(df$Muertos, constant = 1)
upper_bound <- median(df$Muertos) + 3 * mad(df$Muertos, constant = 1)
outlier_ind <- which(df$Muertos< lower_bound | df$Muertos > upper_bound)
rm(lower_bound)
rm(upper_bound)
length(outlier_ind)## [1] 252rosnerTest(df$Muertos, k = 252)## Warning in rosnerTest(df$Muertos, k = 252): The true Type I error may be larger than assumed.
## Although the help file for 'rosnerTest' has a table with information
## on the estimated Type I error level,
## simulations were not run for k > 10 or k > floor(n/2).##
## Results of Outlier Test
## -------------------------
##
## Test Method: Rosner's Test for Outliers
##
## Hypothesized Distribution: Normal
##
## Data: df$Muertos
##
## Sample Size: 25610
##
## Test Statistics: R.1 = 19.15516
## R.2 = 19.29426
## R.3 = 19.43643
## R.4 = 19.58180
## R.5 = 19.73047
## R.6 = 19.88259
## R.7 = 20.03828
## R.8 = 20.19768
## R.9 = 20.36095
## R.10 = 10.21518
## R.11 = 10.23626
## R.12 = 10.25748
## R.13 = 10.27882
## R.14 = 10.30030
## R.15 = 10.32192
## R.16 = 10.34367
## R.17 = 10.36557
## R.18 = 10.38760
## R.19 = 10.40977
## R.20 = 10.43208
## R.21 = 10.45454
## R.22 = 10.47715
## R.23 = 10.49990
## R.24 = 10.52280
## R.25 = 10.54585
## R.26 = 10.56905
## R.27 = 10.59241
## R.28 = 10.61592
## R.29 = 10.63959
## R.30 = 10.66342
## R.31 = 10.68741
## R.32 = 10.71156
## R.33 = 10.73588
## R.34 = 10.76036
## R.35 = 10.78502
## R.36 = 10.80984
## R.37 = 10.83483
## R.38 = 10.86000
## R.39 = 10.88534
## R.40 = 10.91087
## R.41 = 10.93657
## R.42 = 10.96245
## R.43 = 10.98852
## R.44 = 11.01478
## R.45 = 11.04123
## R.46 = 11.06786
## R.47 = 11.09470
## R.48 = 11.12172
## R.49 = 11.14895
## R.50 = 11.17638
## R.51 = 11.20401
## R.52 = 11.23184
## R.53 = 11.25989
## R.54 = 11.28814
## R.55 = 11.31661
## R.56 = 11.34530
## R.57 = 11.37420
## R.58 = 11.40333
## R.59 = 11.43268
## R.60 = 11.46226
## R.61 = 11.49207
## R.62 = 11.52212
## R.63 = 11.55240
## R.64 = 11.58292
## R.65 = 11.61369
## R.66 = 11.64470
## R.67 = 11.67596
## R.68 = 11.70747
## R.69 = 11.73924
## R.70 = 11.77127
## R.71 = 11.80357
## R.72 = 11.83613
## R.73 = 11.86896
## R.74 = 11.90207
## R.75 = 11.93546
## R.76 = 11.96912
## R.77 = 12.00308
## R.78 = 12.03732
## R.79 = 12.07186
## R.80 = 12.10670
## R.81 = 12.14185
## R.82 = 12.17730
## R.83 = 12.21306
## R.84 = 12.24914
## R.85 = 12.28554
## R.86 = 12.32227
## R.87 = 12.35933
## R.88 = 12.39673
## R.89 = 12.43446
## R.90 = 12.47255
## R.91 = 12.51098
## R.92 = 12.54978
## R.93 = 12.58894
## R.94 = 12.62846
## R.95 = 12.66836
## R.96 = 12.70864
## R.97 = 12.74931
## R.98 = 12.79037
## R.99 = 12.83183
## R.100 = 12.87370
## R.101 = 12.91597
## R.102 = 12.95867
## R.103 = 13.00180
## R.104 = 13.04535
## R.105 = 13.08935
## R.106 = 13.13380
## R.107 = 13.17870
## R.108 = 13.22407
## R.109 = 13.26990
## R.110 = 13.31622
## R.111 = 13.36303
## R.112 = 13.41033
## R.113 = 13.45814
## R.114 = 13.50646
## R.115 = 13.55531
## R.116 = 13.60469
## R.117 = 13.65462
## R.118 = 13.70510
## R.119 = 13.75614
## R.120 = 13.80776
## R.121 = 13.85996
## R.122 = 13.91276
## R.123 = 13.96617
## R.124 = 14.02020
## R.125 = 14.07486
## R.126 = 14.13016
## R.127 = 14.18612
## R.128 = 14.24275
## R.129 = 14.30007
## R.130 = 14.35808
## R.131 = 14.41681
## R.132 = 14.47626
## R.133 = 14.53645
## R.134 = 14.59740
## R.135 = 14.65912
## R.136 = 14.72163
## R.137 = 14.78495
## R.138 = 14.84910
## R.139 = 14.91408
## R.140 = 14.97993
## R.141 = 15.04665
## R.142 = 15.11428
## R.143 = 15.18283
## R.144 = 15.25231
## R.145 = 15.32276
## R.146 = 15.39420
## R.147 = 15.46664
## R.148 = 15.54012
## R.149 = 15.61465
## R.150 = 15.69027
## R.151 = 15.76699
## R.152 = 15.84486
## R.153 = 15.92388
## R.154 = 16.00410
## R.155 = 16.08555
## R.156 = 16.16825
## R.157 = 16.25224
## R.158 = 16.33756
## R.159 = 16.42423
## R.160 = 16.51230
## R.161 = 16.60180
## R.162 = 16.69276
## R.163 = 16.78525
## R.164 = 16.87928
## R.165 = 16.97492
## R.166 = 17.07219
## R.167 = 17.17116
## R.168 = 17.27188
## R.169 = 17.37438
## R.170 = 17.47873
## R.171 = 17.58499
## R.172 = 17.69320
## R.173 = 17.80344
## R.174 = 17.91577
## R.175 = 18.03025
## R.176 = 18.14695
## R.177 = 18.26595
## R.178 = 18.38732
## R.179 = 18.51114
## R.180 = 18.63749
## R.181 = 18.76648
## R.182 = 18.89817
## R.183 = 19.03268
## R.184 = 19.17010
## R.185 = 19.31054
## R.186 = 19.45412
## R.187 = 19.60094
## R.188 = 19.75115
## R.189 = 19.90485
## R.190 = 20.06221
## R.191 = 20.22335
## R.192 = 20.38844
## R.193 = 20.55764
## R.194 = 20.73113
## R.195 = 20.90908
## R.196 = 21.09170
## R.197 = 21.27918
## R.198 = 21.47176
## R.199 = 21.66966
## R.200 = 21.87313
## R.201 = 22.08245
## R.202 = 22.29789
## R.203 = 22.51977
## R.204 = 22.74840
## R.205 = 22.98414
## R.206 = 23.22737
## R.207 = 23.47848
## R.208 = 23.73792
## R.209 = 24.00616
## R.210 = 24.28369
## R.211 = 24.57108
## R.212 = 24.86893
## R.213 = 25.17787
## R.214 = 25.49862
## R.215 = 25.83195
## R.216 = 26.17870
## R.217 = 26.53980
## R.218 = 26.91627
## R.219 = 27.30923
## R.220 = 27.71992
## R.221 = 28.14971
## R.222 = 28.60014
## R.223 = 29.07290
## R.224 = 29.56991
## R.225 = 30.09331
## R.226 = 30.64552
## R.227 = 31.22930
## R.228 = 31.84776
## R.229 = 32.50449
## R.230 = 33.20359
## R.231 = 33.94983
## R.232 = 34.74876
## R.233 = 35.60688
## R.234 = 36.53189
## R.235 = 37.53295
## R.236 = 38.62109
## R.237 = 39.80970
## R.238 = 41.11529
## R.239 = 42.55836
## R.240 = 44.16490
## R.241 = 45.96829
## R.242 = 48.01231
## R.243 = 50.35573
## R.244 = 53.07961
## R.245 = 56.29942
## R.246 = 60.18661
## R.247 = 65.00897
## R.248 = 71.21376
## R.249 = 79.61941
## R.250 = 91.93657
## R.251 = 112.59885
## R.252 = 159.23881
##
## Test Statistic Parameter: k = 252
##
## Alternative Hypothesis: Up to 252 observations are not
## from the same Distribution.
##
## Type I Error: 5%
##
## Number of Outliers Detected: 252
##
## i Mean.i SD.i Value Obs.Num R.i+1 lambda.i+1 Outlier
## 1 0 1.019133e-02 0.103878487 2 1295 19.15516 4.757289 TRUE
## 2 1 1.011363e-02 0.103133613 2 3501 19.29426 4.757281 TRUE
## 3 2 1.003593e-02 0.102383202 2 7088 19.43643 4.757273 TRUE
## 4 3 9.958215e-03 0.101627131 2 11470 19.58180 4.757265 TRUE
## 5 4 9.880497e-03 0.100865273 2 18070 19.73047 4.757257 TRUE
## 6 5 9.802773e-03 0.100097497 2 20074 19.88259 4.757249 TRUE
## 7 6 9.725043e-03 0.099323664 2 20361 20.03828 4.757241 TRUE
## 8 7 9.647307e-03 0.098543632 2 25050 20.19768 4.757233 TRUE
## 9 8 9.569565e-03 0.097757253 2 25342 20.36095 4.757225 TRUE
## 10 9 9.491817e-03 0.096964372 1 69 10.21518 4.757217 TRUE
## 11 10 9.453125e-03 0.096768431 1 288 10.23626 4.757209 TRUE
## 12 11 9.414430e-03 0.096572062 1 293 10.25748 4.757201 TRUE
## 13 12 9.375732e-03 0.096375261 1 349 10.27882 4.757194 TRUE
## 14 13 9.337032e-03 0.096178027 1 522 10.30030 4.757186 TRUE
## 15 14 9.298328e-03 0.095980357 1 624 10.32192 4.757178 TRUE
## 16 15 9.259621e-03 0.095782247 1 813 10.34367 4.757170 TRUE
## 17 16 9.220911e-03 0.095583696 1 989 10.36557 4.757162 TRUE
## 18 17 9.182198e-03 0.095384700 1 1045 10.38760 4.757154 TRUE
## 19 18 9.143482e-03 0.095185257 1 1315 10.40977 4.757146 TRUE
## 20 19 9.104763e-03 0.094985363 1 1520 10.43208 4.757138 TRUE
## 21 20 9.066041e-03 0.094785017 1 1537 10.45454 4.757130 TRUE
## 22 21 9.027316e-03 0.094584214 1 1730 10.47715 4.757122 TRUE
## 23 22 8.988588e-03 0.094382953 1 1760 10.49990 4.757114 TRUE
## 24 23 8.949857e-03 0.094181230 1 1823 10.52280 4.757106 TRUE
## 25 24 8.911123e-03 0.093979042 1 2009 10.54585 4.757098 TRUE
## 26 25 8.872386e-03 0.093776386 1 2104 10.56905 4.757091 TRUE
## 27 26 8.833646e-03 0.093573260 1 2120 10.59241 4.757083 TRUE
## 28 27 8.794903e-03 0.093369660 1 2577 10.61592 4.757075 TRUE
## 29 28 8.756157e-03 0.093165582 1 2767 10.63959 4.757067 TRUE
## 30 29 8.717407e-03 0.092961025 1 2835 10.66342 4.757059 TRUE
## 31 30 8.678655e-03 0.092755984 1 3014 10.68741 4.757051 TRUE
## 32 31 8.639900e-03 0.092550456 1 3238 10.71156 4.757043 TRUE
## 33 32 8.601142e-03 0.092344439 1 3521 10.73588 4.757035 TRUE
## 34 33 8.562380e-03 0.092137928 1 3897 10.76036 4.757027 TRUE
## 35 34 8.523616e-03 0.091930921 1 4079 10.78502 4.757019 TRUE
## 36 35 8.484848e-03 0.091723415 1 4109 10.80984 4.757011 TRUE
## 37 36 8.446078e-03 0.091515405 1 4187 10.83483 4.757003 TRUE
## 38 37 8.407305e-03 0.091306888 1 4294 10.86000 4.756995 TRUE
## 39 38 8.368528e-03 0.091097861 1 4509 10.88534 4.756988 TRUE
## 40 39 8.329749e-03 0.090888321 1 4519 10.91087 4.756980 TRUE
## 41 40 8.290966e-03 0.090678263 1 4633 10.93657 4.756972 TRUE
## 42 41 8.252180e-03 0.090467685 1 4768 10.96245 4.756964 TRUE
## 43 42 8.213392e-03 0.090256582 1 4792 10.98852 4.756956 TRUE
## 44 43 8.174600e-03 0.090044951 1 5268 11.01478 4.756948 TRUE
## 45 44 8.135805e-03 0.089832787 1 5277 11.04123 4.756940 TRUE
## 46 45 8.097008e-03 0.089620089 1 5654 11.06786 4.756932 TRUE
## 47 46 8.058207e-03 0.089406850 1 5906 11.09470 4.756924 TRUE
## 48 47 8.019403e-03 0.089193068 1 6301 11.12172 4.756916 TRUE
## 49 48 7.980596e-03 0.088978739 1 6426 11.14895 4.756908 TRUE
## 50 49 7.941786e-03 0.088763859 1 6541 11.17638 4.756900 TRUE
## 51 50 7.902973e-03 0.088548423 1 6662 11.20401 4.756892 TRUE
## 52 51 7.864157e-03 0.088332427 1 6861 11.23184 4.756884 TRUE
## 53 52 7.825338e-03 0.088115869 1 6980 11.25989 4.756876 TRUE
## 54 53 7.786516e-03 0.087898742 1 7217 11.28814 4.756869 TRUE
## 55 54 7.747691e-03 0.087681044 1 7531 11.31661 4.756861 TRUE
## 56 55 7.708863e-03 0.087462769 1 7535 11.34530 4.756853 TRUE
## 57 56 7.670032e-03 0.087243914 1 8153 11.37420 4.756845 TRUE
## 58 57 7.631198e-03 0.087024474 1 8296 11.40333 4.756837 TRUE
## 59 58 7.592361e-03 0.086804445 1 8598 11.43268 4.756829 TRUE
## 60 59 7.553520e-03 0.086583821 1 8659 11.46226 4.756821 TRUE
## 61 60 7.514677e-03 0.086362600 1 8688 11.49207 4.756813 TRUE
## 62 61 7.475831e-03 0.086140775 1 8698 11.52212 4.756805 TRUE
## 63 62 7.436981e-03 0.085918343 1 8699 11.55240 4.756797 TRUE
## 64 63 7.398129e-03 0.085695298 1 8962 11.58292 4.756789 TRUE
## 65 64 7.359273e-03 0.085471636 1 9347 11.61369 4.756781 TRUE
## 66 65 7.320415e-03 0.085247352 1 9359 11.64470 4.756773 TRUE
## 67 66 7.281553e-03 0.085022440 1 9440 11.67596 4.756765 TRUE
## 68 67 7.242689e-03 0.084796897 1 9645 11.70747 4.756757 TRUE
## 69 68 7.203821e-03 0.084570717 1 9795 11.73924 4.756749 TRUE
## 70 69 7.164950e-03 0.084343894 1 10029 11.77127 4.756742 TRUE
## 71 70 7.126077e-03 0.084116424 1 10088 11.80357 4.756734 TRUE
## 72 71 7.087200e-03 0.083888301 1 10103 11.83613 4.756726 TRUE
## 73 72 7.048320e-03 0.083659521 1 10203 11.86896 4.756718 TRUE
## 74 73 7.009437e-03 0.083430076 1 10220 11.90207 4.756710 TRUE
## 75 74 6.970551e-03 0.083199963 1 10245 11.93546 4.756702 TRUE
## 76 75 6.931662e-03 0.082969175 1 10939 11.96912 4.756694 TRUE
## 77 76 6.892770e-03 0.082737707 1 11174 12.00308 4.756686 TRUE
## 78 77 6.853875e-03 0.082505554 1 11181 12.03732 4.756678 TRUE
## 79 78 6.814977e-03 0.082272708 1 11308 12.07186 4.756670 TRUE
## 80 79 6.776076e-03 0.082039165 1 11310 12.10670 4.756662 TRUE
## 81 80 6.737172e-03 0.081804918 1 11597 12.14185 4.756654 TRUE
## 82 81 6.698265e-03 0.081569961 1 11766 12.17730 4.756646 TRUE
## 83 82 6.659354e-03 0.081334289 1 11783 12.21306 4.756638 TRUE
## 84 83 6.620441e-03 0.081097895 1 11862 12.24914 4.756630 TRUE
## 85 84 6.581525e-03 0.080860772 1 12325 12.28554 4.756622 TRUE
## 86 85 6.542605e-03 0.080622914 1 12495 12.32227 4.756614 TRUE
## 87 86 6.503683e-03 0.080384315 1 12518 12.35933 4.756606 TRUE
## 88 87 6.464757e-03 0.080144968 1 12640 12.39673 4.756599 TRUE
## 89 88 6.425829e-03 0.079904866 1 12643 12.43446 4.756591 TRUE
## 90 89 6.386897e-03 0.079664002 1 12785 12.47255 4.756583 TRUE
## 91 90 6.347962e-03 0.079422370 1 12833 12.51098 4.756575 TRUE
## 92 91 6.309025e-03 0.079179963 1 12836 12.54978 4.756567 TRUE
## 93 92 6.270084e-03 0.078936773 1 13021 12.58894 4.756559 TRUE
## 94 93 6.231140e-03 0.078692793 1 13104 12.62846 4.756551 TRUE
## 95 94 6.192193e-03 0.078448015 1 13170 12.66836 4.756543 TRUE
## 96 95 6.153243e-03 0.078202433 1 13180 12.70864 4.756535 TRUE
## 97 96 6.114290e-03 0.077956038 1 13203 12.74931 4.756527 TRUE
## 98 97 6.075334e-03 0.077708823 1 13237 12.79037 4.756519 TRUE
## 99 98 6.036375e-03 0.077460780 1 13261 12.83183 4.756511 TRUE
## 100 99 5.997413e-03 0.077211901 1 13312 12.87370 4.756503 TRUE
## 101 100 5.958448e-03 0.076962177 1 13355 12.91597 4.756495 TRUE
## 102 101 5.919479e-03 0.076711602 1 13394 12.95867 4.756487 TRUE
## 103 102 5.880508e-03 0.076460165 1 13396 13.00180 4.756479 TRUE
## 104 103 5.841534e-03 0.076207860 1 13452 13.04535 4.756471 TRUE
## 105 104 5.802556e-03 0.075954676 1 13558 13.08935 4.756463 TRUE
## 106 105 5.763576e-03 0.075700605 1 13629 13.13380 4.756455 TRUE
## 107 106 5.724592e-03 0.075445639 1 13993 13.17870 4.756448 TRUE
## 108 107 5.685606e-03 0.075189768 1 14004 13.22407 4.756440 TRUE
## 109 108 5.646616e-03 0.074932982 1 14124 13.26990 4.756432 TRUE
## 110 109 5.607623e-03 0.074675273 1 14181 13.31622 4.756424 TRUE
## 111 110 5.568627e-03 0.074416631 1 14219 13.36303 4.756416 TRUE
## 112 111 5.529629e-03 0.074157046 1 14241 13.41033 4.756408 TRUE
## 113 112 5.490627e-03 0.073896508 1 14396 13.45814 4.756400 TRUE
## 114 113 5.451622e-03 0.073635007 1 14501 13.50646 4.756392 TRUE
## 115 114 5.412614e-03 0.073372532 1 14655 13.55531 4.756384 TRUE
## 116 115 5.373603e-03 0.073109074 1 14670 13.60469 4.756376 TRUE
## 117 116 5.334589e-03 0.072844621 1 14724 13.65462 4.756368 TRUE
## 118 117 5.295571e-03 0.072579163 1 14808 13.70510 4.756360 TRUE
## 119 118 5.256551e-03 0.072312688 1 14828 13.75614 4.756352 TRUE
## 120 119 5.217528e-03 0.072045186 1 14831 13.80776 4.756344 TRUE
## 121 120 5.178501e-03 0.071776644 1 14954 13.85996 4.756336 TRUE
## 122 121 5.139472e-03 0.071507051 1 15631 13.91276 4.756328 TRUE
## 123 122 5.100439e-03 0.071236395 1 15675 13.96617 4.756320 TRUE
## 124 123 5.061404e-03 0.070964665 1 15812 14.02020 4.756312 TRUE
## 125 124 5.022365e-03 0.070691846 1 16117 14.07486 4.756304 TRUE
## 126 125 4.983324e-03 0.070417928 1 16140 14.13016 4.756296 TRUE
## 127 126 4.944279e-03 0.070142897 1 16715 14.18612 4.756288 TRUE
## 128 127 4.905231e-03 0.069866739 1 16797 14.24275 4.756280 TRUE
## 129 128 4.866180e-03 0.069589442 1 16858 14.30007 4.756273 TRUE
## 130 129 4.827126e-03 0.069310991 1 16887 14.35808 4.756265 TRUE
## 131 130 4.788069e-03 0.069031373 1 17141 14.41681 4.756257 TRUE
## 132 131 4.749009e-03 0.068750574 1 17622 14.47626 4.756249 TRUE
## 133 132 4.709946e-03 0.068468579 1 17624 14.53645 4.756241 TRUE
## 134 133 4.670880e-03 0.068185372 1 17935 14.59740 4.756233 TRUE
## 135 134 4.631810e-03 0.067900940 1 18106 14.65912 4.756225 TRUE
## 136 135 4.592738e-03 0.067615266 1 18284 14.72163 4.756217 TRUE
## 137 136 4.553663e-03 0.067328335 1 18370 14.78495 4.756209 TRUE
## 138 137 4.514584e-03 0.067040130 1 18444 14.84910 4.756201 TRUE
## 139 138 4.475503e-03 0.066750635 1 18646 14.91408 4.756193 TRUE
## 140 139 4.436418e-03 0.066459834 1 18689 14.97993 4.756185 TRUE
## 141 140 4.397330e-03 0.066167708 1 18756 15.04665 4.756177 TRUE
## 142 141 4.358239e-03 0.065874240 1 19283 15.11428 4.756169 TRUE
## 143 142 4.319146e-03 0.065579413 1 19362 15.18283 4.756161 TRUE
## 144 143 4.280049e-03 0.065283208 1 19539 15.25231 4.756153 TRUE
## 145 144 4.240949e-03 0.064985605 1 19605 15.32276 4.756145 TRUE
## 146 145 4.201846e-03 0.064686587 1 19684 15.39420 4.756137 TRUE
## 147 146 4.162740e-03 0.064386132 1 19801 15.46664 4.756129 TRUE
## 148 147 4.123630e-03 0.064084221 1 19945 15.54012 4.756121 TRUE
## 149 148 4.084518e-03 0.063780832 1 20009 15.61465 4.756113 TRUE
## 150 149 4.045403e-03 0.063475946 1 20201 15.69027 4.756105 TRUE
## 151 150 4.006284e-03 0.063169540 1 20328 15.76699 4.756097 TRUE
## 152 151 3.967163e-03 0.062861592 1 20344 15.84486 4.756089 TRUE
## 153 152 3.928038e-03 0.062552079 1 20354 15.92388 4.756081 TRUE
## 154 153 3.888911e-03 0.062240977 1 20355 16.00410 4.756073 TRUE
## 155 154 3.849780e-03 0.061928264 1 20393 16.08555 4.756065 TRUE
## 156 155 3.810646e-03 0.061613914 1 20471 16.16825 4.756057 TRUE
## 157 156 3.771509e-03 0.061297901 1 20559 16.25224 4.756050 TRUE
## 158 157 3.732369e-03 0.060980202 1 20575 16.33756 4.756042 TRUE
## 159 158 3.693226e-03 0.060660787 1 20780 16.42423 4.756034 TRUE
## 160 159 3.654080e-03 0.060339632 1 20791 16.51230 4.756026 TRUE
## 161 160 3.614931e-03 0.060016706 1 20804 16.60180 4.756018 TRUE
## 162 161 3.575779e-03 0.059691983 1 20993 16.69276 4.756010 TRUE
## 163 162 3.536624e-03 0.059365432 1 21439 16.78525 4.756002 TRUE
## 164 163 3.497465e-03 0.059037022 1 21599 16.87928 4.755994 TRUE
## 165 164 3.458304e-03 0.058706724 1 21743 16.97492 4.755986 TRUE
## 166 165 3.419139e-03 0.058374504 1 22015 17.07219 4.755978 TRUE
## 167 166 3.379972e-03 0.058040330 1 22020 17.17116 4.755970 TRUE
## 168 167 3.340801e-03 0.057704167 1 22199 17.27188 4.755962 TRUE
## 169 168 3.301627e-03 0.057365982 1 22339 17.37438 4.755954 TRUE
## 170 169 3.262450e-03 0.057025736 1 22380 17.47873 4.755946 TRUE
## 171 170 3.223270e-03 0.056683395 1 22384 17.58499 4.755938 TRUE
## 172 171 3.184087e-03 0.056338919 1 22482 17.69320 4.755930 TRUE
## 173 172 3.144901e-03 0.055992269 1 22497 17.80344 4.755922 TRUE
## 174 173 3.105712e-03 0.055643404 1 22603 17.91577 4.755914 TRUE
## 175 174 3.066520e-03 0.055292283 1 22623 18.03025 4.755906 TRUE
## 176 175 3.027325e-03 0.054938862 1 22632 18.14695 4.755898 TRUE
## 177 176 2.988126e-03 0.054583096 1 22661 18.26595 4.755890 TRUE
## 178 177 2.948925e-03 0.054224940 1 22674 18.38732 4.755882 TRUE
## 179 178 2.909720e-03 0.053864345 1 22675 18.51114 4.755874 TRUE
## 180 179 2.870512e-03 0.053501262 1 22693 18.63749 4.755866 TRUE
## 181 180 2.831302e-03 0.053135641 1 22695 18.76648 4.755858 TRUE
## 182 181 2.792088e-03 0.052767429 1 22730 18.89817 4.755850 TRUE
## 183 182 2.752871e-03 0.052396570 1 22747 19.03268 4.755842 TRUE
## 184 183 2.713651e-03 0.052023008 1 22804 19.17010 4.755834 TRUE
## 185 184 2.674428e-03 0.051646685 1 22850 19.31054 4.755826 TRUE
## 186 185 2.635202e-03 0.051267540 1 22999 19.45412 4.755818 TRUE
## 187 186 2.595972e-03 0.050885510 1 23067 19.60094 4.755810 TRUE
## 188 187 2.556740e-03 0.050500528 1 23089 19.75115 4.755802 TRUE
## 189 188 2.517505e-03 0.050112528 1 23265 19.90485 4.755794 TRUE
## 190 189 2.478266e-03 0.049721439 1 23377 20.06221 4.755786 TRUE
## 191 190 2.439024e-03 0.049327186 1 23395 20.22335 4.755778 TRUE
## 192 191 2.399780e-03 0.048929694 1 23422 20.38844 4.755770 TRUE
## 193 192 2.360532e-03 0.048528883 1 23447 20.55764 4.755762 TRUE
## 194 193 2.321281e-03 0.048124669 1 23473 20.73113 4.755754 TRUE
## 195 194 2.282027e-03 0.047716968 1 23492 20.90908 4.755746 TRUE
## 196 195 2.242770e-03 0.047305687 1 23544 21.09170 4.755739 TRUE
## 197 196 2.203510e-03 0.046890734 1 23587 21.27918 4.755731 TRUE
## 198 197 2.164247e-03 0.046472009 1 23605 21.47176 4.755723 TRUE
## 199 198 2.124980e-03 0.046049411 1 23639 21.66966 4.755715 TRUE
## 200 199 2.085711e-03 0.045622830 1 23667 21.87313 4.755707 TRUE
## 201 200 2.046438e-03 0.045192155 1 23725 22.08245 4.755699 TRUE
## 202 201 2.007163e-03 0.044757267 1 23840 22.29789 4.755691 TRUE
## 203 202 1.967884e-03 0.044318042 1 23879 22.51977 4.755683 TRUE
## 204 203 1.928602e-03 0.043874350 1 23880 22.74840 4.755675 TRUE
## 205 204 1.889317e-03 0.043426054 1 23928 22.98414 4.755667 TRUE
## 206 205 1.850030e-03 0.042973010 1 23950 23.22737 4.755659 TRUE
## 207 206 1.810738e-03 0.042515066 1 24064 23.47848 4.755651 TRUE
## 208 207 1.771444e-03 0.042052062 1 24065 23.73792 4.755643 TRUE
## 209 208 1.732147e-03 0.041583829 1 24217 24.00616 4.755635 TRUE
## 210 209 1.692847e-03 0.041110188 1 24237 24.28369 4.755627 TRUE
## 211 210 1.653543e-03 0.040630950 1 24245 24.57108 4.755619 TRUE
## 212 211 1.614237e-03 0.040145915 1 24263 24.86893 4.755611 TRUE
## 213 212 1.574927e-03 0.039654870 1 24345 25.17787 4.755603 TRUE
## 214 213 1.535614e-03 0.039157588 1 24353 25.49862 4.755595 TRUE
## 215 214 1.496299e-03 0.038653830 1 24362 25.83195 4.755587 TRUE
## 216 215 1.456980e-03 0.038143338 1 24371 26.17870 4.755579 TRUE
## 217 216 1.417658e-03 0.037625838 1 24419 26.53980 4.755571 TRUE
## 218 217 1.378333e-03 0.037101038 1 24472 26.91627 4.755563 TRUE
## 219 218 1.339004e-03 0.036568622 1 24523 27.30923 4.755555 TRUE
## 220 219 1.299673e-03 0.036028254 1 24524 27.71992 4.755547 TRUE
## 221 220 1.260339e-03 0.035479569 1 24531 28.14971 4.755539 TRUE
## 222 221 1.221001e-03 0.034922177 1 24533 28.60014 4.755531 TRUE
## 223 222 1.181661e-03 0.034355652 1 24617 29.07290 4.755523 TRUE
## 224 223 1.142317e-03 0.033779535 1 24630 29.56991 4.755515 TRUE
## 225 224 1.102970e-03 0.033193328 1 24636 30.09331 4.755507 TRUE
## 226 225 1.063620e-03 0.032596485 1 24647 30.64552 4.755499 TRUE
## 227 226 1.024267e-03 0.031988411 1 24650 31.22930 4.755491 TRUE
## 228 227 9.849112e-04 0.031368454 1 24713 31.84776 4.755483 TRUE
## 229 228 9.455520e-04 0.030735893 1 24719 32.50449 4.755475 TRUE
## 230 229 9.061897e-04 0.030089935 1 24729 33.20359 4.755467 TRUE
## 231 230 8.668243e-04 0.029429696 1 24735 33.94983 4.755459 TRUE
## 232 231 8.274558e-04 0.028754194 1 24785 34.74876 4.755451 TRUE
## 233 232 7.880842e-04 0.028062326 1 24810 35.60688 4.755443 TRUE
## 234 233 7.487095e-04 0.027352850 1 24850 36.53189 4.755435 TRUE
## 235 234 7.093317e-04 0.026624358 1 24857 37.53295 4.755427 TRUE
## 236 235 6.699507e-04 0.025875245 1 24864 38.62109 4.755419 TRUE
## 237 236 6.305667e-04 0.025103664 1 24868 39.80970 4.755411 TRUE
## 238 237 5.911796e-04 0.024307476 1 24963 41.11529 4.755403 TRUE
## 239 238 5.517894e-04 0.023484178 1 24990 42.55836 4.755395 TRUE
## 240 239 5.123960e-04 0.022630813 1 25037 44.16490 4.755387 TRUE
## 241 240 4.729996e-04 0.021743839 1 25108 45.96829 4.755379 TRUE
## 242 241 4.336001e-04 0.020818961 1 25250 48.01231 4.755371 TRUE
## 243 242 3.941974e-04 0.019850883 1 25268 50.35573 4.755363 TRUE
## 244 243 3.547917e-04 0.018832944 1 25293 53.07961 4.755355 TRUE
## 245 244 3.153828e-04 0.017756569 1 25401 56.29942 4.755347 TRUE
## 246 245 2.759708e-04 0.016610405 1 25475 60.18661 4.755339 TRUE
## 247 246 2.365557e-04 0.015378853 1 25478 65.00897 4.755331 TRUE
## 248 247 1.971376e-04 0.014039461 1 25485 71.21376 4.755323 TRUE
## 249 248 1.577163e-04 0.012557771 1 25490 79.61941 4.755315 TRUE
## 250 249 1.182919e-04 0.010875778 1 25510 91.93657 4.755307 TRUE
## 251 250 7.886435e-05 0.008880385 1 25527 112.59885 4.755299 TRUE
## 252 251 3.943373e-05 0.006279628 1 25609 159.23881 4.755291 TRUELos graficos, la prueba de Rosner, y el flitro de Hampel han marcado a las 252 observaciones en cuestión como datos atipicos. Sin embargo, por la estrucutra que presentaron los datos, faltantes es necesario revisar su los datos atipicos no estan contenidos en un clase particular.
es_outlier <- ifelse(seq_along(df$Gravedad) %in% outlier_ind, 1, 0)
df_outlier <- as.data.frame(cbind(as.character(df$Gravedad),es_outlier))
colnames(df_outlier) <-c("Gravedad","outlier")
df_outlier$outlier <- as.numeric(df_outlier$outlier)df_outlier %>%
group_by(Gravedad) %>%
summarise(
total = n(),
Casos = sum(outlier),
Proporcion = round(mean(outlier),2),
Procentaje = round(mean(outlier) * 100,3 ))## # A tibble: 3 × 5
## Gravedad total Casos Proporcion Procentaje
## <chr> <int> <dbl> <dbl> <dbl>
## 1 Con heridos 9901 0 0 0
## 2 Con muertos 252 252 1 100
## 3 Solo daños 15457 0 0 0Los resultados muestran que todos los datos atípicos están concentrados en la categoría “con muertos”. Esto sugiere que los valores observados en la variable Muertos probablemente reflejan la gravedad real de los accidentes en lugar de ser errores de digitación. Transformar estos datos en faltantes e imputarlos no es una estrategia adecuada, ya que podría distorsionar el análisis y subestimar la gravedad de ciertos accidentes. Por lo tanto, se ha decidido mantener intactos estos valores, preservando la integridad de los datos y asegurando que el análisis refleje fielmente la realidad de los hechos observados.
Heridos
Procedemos de la misma manera que con la variable
Muertos.
df%>%
get_summary_stats(Heridos, type = "five_number")## # A tibble: 1 × 7
## variable n min max q1 median q3
## <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Heridos 25610 0 42 0 0 1boxplot(df$Heridos)
hist(df$Heridos)
Los graficos revelan la existencia de valores atipicos en la variable
Heridos. Por lo cual, antes de proseguir utilizaremos el
filtro de Hampel y el método de los percentiles para estudiar más a
profundidad estas observaciones.
lower_bound <- quantile(df$Heridos, 0.025)
upper_bound <- quantile(df$Heridos, 0.975)
outlier_ind <- which(df$Heridos < lower_bound | df$Heridos > upper_bound)
# Crear una columna que marque los outliers
df$Perc <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)
# Ver el número de outliers
length(outlier_ind)## [1] 303Hora el flitro de Hampel
# Calcular los límites usando el filtro de Hampel
lower_bound <- median(df$Heridos, na.rm = TRUE) - 3 * mad(df$Heridos, constant = 1, na.rm = TRUE)
upper_bound <- median(df$Heridos, na.rm = TRUE) + 3 * mad(df$Heridos, constant = 1, na.rm = TRUE)
# Identificar los índices de los outliers
outlier_ind <- which(df$Heridos < lower_bound | df$Heridos > upper_bound)
# Crear una columna que marque los outliers
df$Hampel <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)
# Ver el número de outliers
length(outlier_ind)## [1] 9984Revisemos las discrepansias entre los dos criterios:
table(df$Hampel, df$Perc)##
## 0 1
## 0 15626 0
## 1 9681 303Notamos que ambos criterios defieren considerablemente, por lo cual se aplica la prueba de Rosner:
test <- rosnerTest(df$Heridos, k = 1000)## Warning in rosnerTest(df$Heridos, k = 1000): The true Type I error may be larger than assumed.
## Although the help file for 'rosnerTest' has a table with information
## on the estimated Type I error level,
## simulations were not run for k > 10 or k > floor(n/2).test <- test$all.stats
outlier_ind <- test[test$Outlier == TRUE,]$Obs.Num
df$Rosner <- ifelse(1:nrow(df) %in% outlier_ind, 1, 0)df %>%
group_by(Gravedad) %>%
summarise(
total = n(),
Casos_Hampel = sum(Hampel),
Casos_perc = sum(Perc),
Casos_rosenr = sum(Rosner),
Procentaje_Hampel = round((sum(Hampel)/9984) * 100,3 ),
Procentaje_perc = round((sum(Perc)/303) * 100,3 ),
Procentaje_Rosenr = round((sum(Rosner)/148) * 100,3 ))## # A tibble: 3 × 8
## Gravedad total Casos_Hampel Casos_perc Casos_rosenr Procentaje_Hampel
## <fct> <int> <dbl> <dbl> <dbl> <dbl>
## 1 Con heridos 9901 9901 297 145 99.2
## 2 Con muertos 252 83 6 3 0.831
## 3 Solo daños 15457 0 0 0 0
## # ℹ 2 more variables: Procentaje_perc <dbl>, Procentaje_Rosenr <dbl>Notamos que, tal como se reporte en la variable Muertos,
los datos atípicos están concentrados en la categoría “con heridos”.
Esto sugiere que los valores observados en la variable
Heridos probablemente reflejan la gravedad real de los
accidentes en lugar de ser errores de digitación. Por lo tanto, se ha
decidido mantener intactos estos valores.