INTRODUCCIÓN
Los datos abiertos son una gran fuente de información para aquellos que deseen estudiar las características del entorno en que se encuentran en aras de brindar alternativas que permitan mejorar la situación actual. Dichas alternativas varían desde un sencillo análisis descriptivo para entender los sucesos que se han venido presentando históricamente, hasta la implementación de programas más avanzadas como las herramientas de Aprendizaje de Máquinas que permiten encontrar patrones en los datos que no se pueden visualizar a simple vista por un humano o el entrenamiento de un algoritmo predictivo que ayude a realizar estimaciones de lo que pasará en el futuro, con lo que se pueden apoyar sistemas de toma de decisiones o el desarrollo de estrategias de optimización.
En este trabajo se propone la utilización de herramientas de aprendizaje de máquinas para realizar un agrupamiento y predicción de la accidentalidad en Medellín a partir de los datos recientes de accidentalidad reportados en las bases de datos abiertos de MEData.
Las tablas descargadas, elcódigo y la limpieza de los datos se encuentran alojadas en este repositorio. Se construyó además una aplicación en Shiny para que los usuarios pueda explorar de forma interactiva el estudio de accidentalidad realizado.
Como parte de los lineamientos para la realización de este trabajo, se solicitó considerar fechas especiales dentro de los años de 2014 a 2018. En este caso, se crea una nueva variable en la base de datos que contenga información sobre los días para determinar si es un día laboral o no y adicional para identificar cuáles de los días laborales es un día festivo, a partir de un conjunto de datos encontrado en internet con información de los días festivos en Colombia. Se puede consultar el archivo de limpieza de datos para ver la integración de las fechas especiales al conjunto de datos.
EXPLORACIÓN DE LOS DATOS
Se realiza un análisis exploratorio para descubrir las principales características de la información contenida en la base de datos de accidentalidad a través de gráficos. Con este análisis se busca completar la limpieza de los datos de manera que se adecuen para realizar los procesos posteriores de predicción y agrupamiento.
Tabla de datos y variables
Se muestra a continuación los primeros seis registros de la base de datos “limpia”.
| X | FECHA | DIA | PERIODO | CLASE | GRAVEDAD | DIA_NOMBRE | MES | MES_NOMBRE | LONGITUD | LATITUD | SEMANA | BARRIO | COMUNA | FECHAS_ESP |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 2014-01-01 | 1 | 2014 | CHOQUE | HERIDO | MIERCOLES | 1 | ENERO | -75.60273 | 6.219016 | 1 | LOMA DE LOS BERNAL | BELEN | FESTIVO |
| 2 | 2014-01-01 | 1 | 2014 | ATROPELLO | HERIDO | MIERCOLES | 1 | ENERO | -75.56818 | 6.260009 | 1 | JESUS NAZARENO | LA CANDELARIA | FESTIVO |
| 3 | 2014-01-01 | 1 | 2014 | ATROPELLO | HERIDO | MIERCOLES | 1 | ENERO | -75.54994 | 6.264765 | 1 | MANRIQUE ORIENTAL | MANRIQUE | FESTIVO |
| 4 | 2014-01-01 | 1 | 2014 | ATROPELLO | HERIDO | MIERCOLES | 1 | ENERO | -75.60761 | 6.234327 | 1 | LAS MERCEDES | BELEN | FESTIVO |
| 5 | 2014-01-01 | 1 | 2014 | CAIDA OCUPANTE | HERIDO | MIERCOLES | 1 | ENERO | -75.57969 | 6.299968 | 1 | DOCE DE OCTUBRE NO. 2 | DOCE DE OCTUBRE | FESTIVO |
| 6 | 2014-01-01 | 1 | 2014 | CHOQUE | HERIDO | MIERCOLES | 1 | ENERO | -75.55543 | 6.284499 | 1 | BERLIN | ARANJUEZ | FESTIVO |
Definición de las variables
- FECHA: fecha completa de ocurrencia del accidente.
- DIA: día de ocurrencia del accidene (día del mes) de acuerdo a las fechas.
- PERIODO: año de occurrencia del accidente, de acuerdo a las fechas.
- CLASE: indica el tipo de accidente: choque, atropello, caída ocupante, incendio, volcamiento.
- GRAVEDAD: indica las consecuencias del accidente: herido, muerto, solo daños.
- DIA_NOMBRE: nombre del día de la semana, correspondiente a la fecha de ocurrencia del accidente.
- MES: tiempo de ocurrencia del accidente, mes del año, de acuerdo a las fechas registradas.
- MES_NOMBRE: corresponde al cambio del mes de formato numérico a texto, donde ocurrieron los accidentes.
- LATITUD Y LONGITUD: coordenadas de los accidentes de la base de datos de movilidad de Medellín.
- SEMANA: semana del año de ocurrencia del accidente.
- BARRIO: lugar de la ciudad de Medellín donde ocurrieron los accidentes registrados.
- COMUNA: subdivisión administrativa donde se encuentran contenidos los barrios de una ciudad, y donde ocurrieron los accidentes registrados.
- FECHAS_ESP: es una diferenciación de los días en: día laborar, día no laboral y festivos.
Distribución de los datos por comunas
Con la ayuda de un gráfico de barras se identifica la distribución de los accidentes por comuna durante el periodo de estudio. En este gráfico se muestra el porcentaje de la accidentalidad total que corresponde a cada comuna y corregimiento de Medellín.
La menor ocurrencia de accidentes se presenta en los corregimientos. Este resultado era de esperarse, ya que los corregimientos correspondes a zonas rurales de la ciudad, las cuales no alcanza a llegar al 1% de la accidentalidad total reistrada durante el periodo de estudio, por tanto, se toma la decisión de eliminar los datos de los corregimientos para no afectar los análisis y los modelos predictivos.
Se vuelve a graficar la distribución de la acciedentalidad por comunas para verificar el cambio realizado y los nuevos porcentajes de accidentalidad.
La comuna con mayor porcentaje de accidentes entre los años 2014 y 2018 es La Candelaria, que se compone de los barrios de la zona central de la ciudad y con mayor flujo vehicular, debido a la alta concentración de comercios y entidades administrativas de la ciudad. A esta comuna le sigue la de Laureles - Estadio, en la cual confluyen diferentes vías principales de la ciudad como son la avenida San Juan y las carreras 70 y 80. A estás comunas les sigue la de Castilla, que ocupa cierto tramo de la autopista norte.
Accidentalidad por barrios
Se realiza a continuación un pareto donde se muestran los 10 primeros barrios de acuerdo al número de accidentes ocurridos durante el periodo comprendido entre 2014 y 2018.
Se observa que entre los 265 barrios de la ciudad de Medellín, fueron La Candelaria y Caribe los que presentaron la mayor cantidad de accidentes durante el periodo de estudio. El primer barrio se encuentra en la zona central de la ciudad, con mayor accidentalidad debido a lo explicado antes. El segundo cuenta con una alta afluencia de vehículos ya que por el cruza la carrera 46 C, muy utilizada por los vehículos que se dirigen hacia el norte de la ciudad.
Gráfico clase y tipo de accidente
En las bases de datos la ocurrencia de accidentes se encuentra clasificada de acuerdo a la gravedad y clase de acccidentes. En la siguiente gráfica se muestra cómo están distribuídos los accidentes respecto a estos parámetros.
La mayor cantidad de accidentes presentados son CHOQUES, de los cuales en su mayoría solo ocurren daños, seguido de eventos con heridos y eventos con muertes. La clase OTROS presenta una alta participación en la ocurrencia de accidentes. Habría que realizar un estudio más profundo para identificar qué tipo de eventos se presentan en estos casos.
Los ATROPELLOS y CAIDA DE OCUPANTES presentan una alta cantidad de personas heridas, en ninguno de estos casos se observa solo ocurrencia de daños y ocurren muchas más muertes por los atropellos que por las caídas.
Los VOLCAMIENTOS se presentan en su mayoría heridos, seguido de solo daños y solo 4 muertes ocurridas en el período de estudio. Finalmente, la clase con menos participación es la de INCENDIO, con 15 eventos con daños y 14 con heridos sin ninguna ocurrencia de muertes. Se procede a eliminar esta clase del estudio debido a su baja participación y se muestra la cantidad total de accidentes para las otras clases:
ATROPELLO CAIDA OCUPANTE CHOQUE OTRO VOLCAMIENTO
19919 17791 138229 21243 6351 Accidentalidad por año y mes
El siguiente gráfico muestra la distribución de la accidentalidad por mes cada año.
En general, la menor ocurrencia de accidentes se encuentra en el mes de enero, que es un mes donde no se celebra alguna festividad importante en la ciudad y el flujo de personas es relativamente bajo, debido a los foráneos que regresan a sus lugares de procedencia en el mes de diciembre.
En contraste, en promedio, agosto es el mes con la mayor ocurrencia de accidentes. Se debe tener en cuenta que en este mes se celebra la Feria de las Flores y que la afluencia de personas es alta y la ciudad cuenta con la llegada de personas de todas partes en estas fechas.
En los meses de marzo, mayo, septiembre y octubre también se observan valores altos de accidentalidad en la ciudad, mentras que el resto de los meses la accidentalidad parece mantenerse en un promedio.
Modelamiento
Para la selección del modelo predictivo que se utilizará para construir la aplicación en Shiny a partir de los datos de accidentalidad en Medellín, se procede a comparar diferentes modelos y evaluar su desempeño en entrenamiento y prueba calculando la el error cuadrático medio (MSE). En este caso, variaciones de más del 15% entre los MSE obtenidos en ambos conjuntos, se considerarán sospechosas de sobreentrenamiento del modelo.
Primero se seleccionan las variables que se utilizarán, que para este trabajo serán: DIA, MES, SEMANA, FECHAS_ESP, PERIODO, COMUNA, BARRIO y CLASE. Además de que se agrega la variable DIA_ANNO que muestra el número del día del año (1-365).
Luego se procede a dar un formato adecuado a las variables; se convertirá las variables BARRIO, COMUNAS y CLASE a categóricas (factor). Se tiene un número finito de barrios, sin embargo es bastante alto (271), por lo que se manejará la variable como continua y no como categórica. Las otras variables se convertirán nuevamente a factor luego de convertirlas en numéricas, ya que los modelos pueden trabajar mejor si se ingresan números en lugar de texo. Adicionalmente, las variables MES y FECHAS_ESP se trabajarán como categórica, mientras que se tratarán las variables DIA, SEMANA y PERIODO como continuas.
| FECHA | DIA | MES | SEMANA | FECHAS_ESP | PERIODO | COMUNA | BARRIO | CLASE | DIA_ANNO |
|---|---|---|---|---|---|---|---|---|---|
| 2014-01-01 | 1 | 1 | 1 | 1 | 2014 | 16 | 156 | 3 | 1 |
| 2014-01-01 | 1 | 1 | 1 | 1 | 2014 | 10 | 107 | 1 | 1 |
| 2014-01-01 | 1 | 1 | 1 | 1 | 2014 | 3 | 174 | 1 | 1 |
| 2014-01-01 | 1 | 1 | 1 | 1 | 2014 | 16 | 150 | 1 | 1 |
| 2014-01-01 | 1 | 1 | 1 | 1 | 2014 | 6 | 62 | 2 | 1 |
| 2014-01-01 | 1 | 1 | 1 | 1 | 2014 | 4 | 28 | 3 | 1 |
Se pretende precedir la accidentalidad por mes, semana y día, por lo que se crean tres nuevos dataframes. El primero, realizando un conteo de la accidentalidad por día, en el segundo se cuenta la accidentalidad por semana, y el tercero cuenta la accidentalidad por mes para cada año. Para esto se utilizan las funciones group_by y summarise del paquete dplyr.
Para cada dataframe se procede a entrenar 3 modelos y calcular su desempeño tanto en los datos de entrenamiento como en los datos de validación. Debido a que la cantidad de accidentes es un conteo, se hará uso del modelo lineal generalizado con distribución Poisson. El segundo modelo a entrenar será el de k-Vecinos cercanos, para el cual se planea realizar un proceso previo para seleccionar un k óptimo para cada conjunto de datos. Por último, se ajustará un modelo de Bosque Aleatorio, para un determinado número de árboles.
Se procederá a calcular el MSE de entrenamiento y prueba para cada modelo y calcular la diferencia porcentual entre los mismos, de manera que se pueda establecer si el modelo está sobreentrenado (diferencia mayor al 15%) o no.
Modelamiento de accidentalidad por día
Se realiza el conteo del número de accidentes ocurridos cada día durante el período comprendido entre el 2014 y el 2018.
| PERIODO | DIA_ANNO | FECHAS_ESP | COMUNA | BARRIO | CLASE | CANTIDAD |
|---|---|---|---|---|---|---|
| 2014 | 1 | 1 | 1 | 184 | 1 | 1 |
| 2014 | 1 | 1 | 1 | 206 | 1 | 1 |
| 2014 | 1 | 1 | 1 | 224 | 1 | 1 |
| 2014 | 1 | 1 | 1 | 224 | 3 | 1 |
| 2014 | 1 | 1 | 1 | 236 | 3 | 1 |
| 2014 | 1 | 1 | 3 | 46 | 3 | 1 |
Se observa la distribución de la ocurrencia de accidentes.
En este caso la mayoría de los accidentes son únicos, debido a que el conteo se encuentra discretizado por diferentes variables.
Se dividen los datos en entrenamiento y prueba.
[1] "El conjunto de entrenamiento tiene 129975 filas y 7 columnas"[1] "El conjunto de prueba tiene 30820 filas y 7 columnas"Ahora se procede realizar el entrenamiento y predicción de accidentalidad utilizando los modelos mencionados anteriormente.
Modelo lineal generalizado - Distribución Poisson
Se procede a ajustar el modelo con todas las variables utilizando la función glm, especificando la distribución Poisson como argumento.
Call:
glm(formula = CANTIDAD ~ ., family = "poisson", data = trainSet_dia)
Deviance Residuals:
Min 1Q Median 3Q Max
-0.6043 -0.2795 -0.1560 0.0783 4.1997
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.507e+00 4.458e+00 -0.787 0.43149
PERIODO 1.670e-03 2.212e-03 0.755 0.45023
DIA_ANNO 4.303e-05 2.372e-05 1.814 0.06963 .
FECHAS_ESP2 7.951e-02 1.482e-02 5.366 8.04e-08 ***
FECHAS_ESP3 1.148e-02 1.538e-02 0.747 0.45535
COMUNA2 -1.773e-02 2.923e-02 -0.607 0.54405
COMUNA3 2.394e-02 2.455e-02 0.975 0.32945
COMUNA4 3.341e-02 2.265e-02 1.475 0.14022
COMUNA5 1.838e-01 2.188e-02 8.403 < 2e-16 ***
COMUNA6 2.695e-02 2.501e-02 1.078 0.28125
COMUNA7 4.454e-02 2.251e-02 1.979 0.04782 *
COMUNA8 2.417e-03 2.523e-02 0.096 0.92366
COMUNA9 1.206e-02 2.423e-02 0.498 0.61877
COMUNA10 3.301e-01 2.112e-02 15.632 < 2e-16 ***
COMUNA11 2.075e-01 2.167e-02 9.578 < 2e-16 ***
COMUNA12 3.212e-03 2.441e-02 0.132 0.89534
COMUNA13 -1.373e-02 2.702e-02 -0.508 0.61137
COMUNA14 9.735e-02 2.223e-02 4.380 1.19e-05 ***
COMUNA15 2.727e-01 2.233e-02 12.209 < 2e-16 ***
COMUNA16 6.253e-02 2.256e-02 2.772 0.00557 **
BARRIO 2.733e-05 3.318e-05 0.823 0.41024
CLASE2 -2.017e-02 1.147e-02 -1.758 0.07878 .
CLASE3 2.549e-01 8.474e-03 30.082 < 2e-16 ***
CLASE4 -1.069e-02 1.090e-02 -0.981 0.32679
CLASE5 -4.803e-02 1.591e-02 -3.018 0.00254 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 32478 on 129974 degrees of freedom
Residual deviance: 26697 on 129950 degrees of freedom
AIC: 304929
Number of Fisher Scoring iterations: 4El PERIODO no resulta ser significativo para el modelo, con un valor p bastante alto que no brinda evidencia para rechazar la hipótesis nula, por lo que se concluye que en este modelo ingresar el año no aporta información.
La inclusión de las comunas 5, 10, 11, 14, 15 y 16 es bastante significativa para el modelo. Cabe destacar que las comunas mencionadas son las de mayor ocurrencia de accidentalidad en Medellín durante el periodo de estudio. Entre las clases de accidentes, los CHOQUES (CLASE3) son los de mayor ocurrencia, además de los más significativos para el modelo ajustado, con un valor p que es esencialmente cero, seguido de la clase de accidentes VOLCAMIENTO (CLASE6).
Por último, los días No Laborales no parecen ser significativos para el ajuste, de acuerdo con este modelo. En contraste la inclusión de los días laborales es muy significativa. Además de que la variable clave para la predicción con este modelo (columna DIA_ANNO) muestra un valor p bastante alto, que indica que no puede favorecer la hipótesis alternativa y por tanto dicha variable no presenta relación alguna con la variable de respuesta bajo este escenario. Esto representa un inconveniente, ya que el objetivo es predecir la accidentalidad en Medellín por día.
Se procede a continuación a eliminar variables del modelo de una en una para observar el efecto separado sobre el ajuste del modelo y las otras variables.
Call:
glm(formula = CANTIDAD ~ . - PERIODO, family = "poisson", data = trainSet_dia)
Deviance Residuals:
Min 1Q Median 3Q Max
-0.6013 -0.2795 -0.1559 0.0783 4.1950
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -1.410e-01 2.651e-02 -5.319 1.04e-07 ***
DIA_ANNO 4.310e-05 2.372e-05 1.817 0.06921 .
FECHAS_ESP2 7.951e-02 1.482e-02 5.366 8.03e-08 ***
FECHAS_ESP3 1.146e-02 1.538e-02 0.745 0.45610
COMUNA2 -1.778e-02 2.923e-02 -0.608 0.54308
COMUNA3 2.392e-02 2.455e-02 0.974 0.32986
COMUNA4 3.346e-02 2.265e-02 1.477 0.13964
COMUNA5 1.839e-01 2.187e-02 8.408 < 2e-16 ***
COMUNA6 2.706e-02 2.501e-02 1.082 0.27936
COMUNA7 4.458e-02 2.251e-02 1.981 0.04764 *
COMUNA8 2.436e-03 2.523e-02 0.097 0.92307
COMUNA9 1.212e-02 2.423e-02 0.500 0.61700
COMUNA10 3.302e-01 2.112e-02 15.634 < 2e-16 ***
COMUNA11 2.076e-01 2.167e-02 9.582 < 2e-16 ***
COMUNA12 3.307e-03 2.441e-02 0.135 0.89223
COMUNA13 -1.367e-02 2.702e-02 -0.506 0.61290
COMUNA14 9.749e-02 2.223e-02 4.387 1.15e-05 ***
COMUNA15 2.728e-01 2.233e-02 12.214 < 2e-16 ***
COMUNA16 6.267e-02 2.256e-02 2.778 0.00546 **
BARRIO 2.734e-05 3.318e-05 0.824 0.41005
CLASE2 -2.012e-02 1.147e-02 -1.753 0.07954 .
CLASE3 2.551e-01 8.470e-03 30.118 < 2e-16 ***
CLASE4 -1.048e-02 1.089e-02 -0.962 0.33629
CLASE5 -4.759e-02 1.590e-02 -2.993 0.00276 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 32478 on 129974 degrees of freedom
Residual deviance: 26697 on 129951 degrees of freedom
AIC: 304928
Number of Fisher Scoring iterations: 4La eliminación de la columna PERIODO no altera demasiado la interacción de las demás variables con la variable de respuesta, salvo el INTERCEPTO, que ahora se hace muy significativo para el modelo, con un valor p bastante bajo. La interacción de la columna DIA_ANNO con la variable respuesta no se vio afectada.
Se procede ahora a eliminar la variable CLASE, para observar si la interacción deseada mejora.
Call:
glm(formula = CANTIDAD ~ . - CLASE, family = "poisson", data = trainSet_dia)
Deviance Residuals:
Min 1Q Median 3Q Max
-0.5059 -0.3000 -0.1345 -0.0444 4.3911
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -6.472e+00 4.455e+00 -1.453 0.146309
PERIODO 3.185e-03 2.210e-03 1.441 0.149635
DIA_ANNO 4.637e-05 2.373e-05 1.954 0.050648 .
FECHAS_ESP2 8.676e-02 1.482e-02 5.856 4.75e-09 ***
FECHAS_ESP3 1.946e-02 1.538e-02 1.265 0.205752
COMUNA2 7.475e-03 2.922e-02 0.256 0.798111
COMUNA3 4.660e-02 2.454e-02 1.899 0.057581 .
COMUNA4 8.380e-02 2.261e-02 3.706 0.000211 ***
COMUNA5 2.361e-01 2.181e-02 10.825 < 2e-16 ***
COMUNA6 2.993e-02 2.500e-02 1.197 0.231195
COMUNA7 8.399e-02 2.245e-02 3.741 0.000183 ***
COMUNA8 3.839e-02 2.520e-02 1.523 0.127671
COMUNA9 6.416e-02 2.418e-02 2.654 0.007964 **
COMUNA10 3.963e-01 2.106e-02 18.820 < 2e-16 ***
COMUNA11 2.848e-01 2.158e-02 13.199 < 2e-16 ***
COMUNA12 7.604e-02 2.435e-02 3.123 0.001789 **
COMUNA13 1.041e-02 2.700e-02 0.386 0.699802
COMUNA14 2.004e-01 2.210e-02 9.070 < 2e-16 ***
COMUNA15 3.365e-01 2.226e-02 15.114 < 2e-16 ***
COMUNA16 1.334e-01 2.248e-02 5.933 2.97e-09 ***
BARRIO 5.595e-05 3.319e-05 1.686 0.091879 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 32478 on 129974 degrees of freedom
Residual deviance: 29306 on 129954 degrees of freedom
AIC: 307531
Number of Fisher Scoring iterations: 4Se observa en este caso que el eliminar la variable CLASE tiene un efecto positivo en la interacción de las demás variables con la vaiable de respuesta. En este caso, aumenta el número de comunas que el modelo considera significativas para el ajuste. Además de que el valor p para la columna BARRIO disminuyó considerablemente.
La inclusión de la columna PERIODO provoca que el INTERCEPTO deje de ser significativo para el modelo. Además, la columna DIA_ANNO sigue mostrando un valor p que no permite rechazar la hipótesis nula.
Se ajusta ahora un cuarto modelo eliminando las columnas CLASE y PERIODO.
Call:
glm(formula = CANTIDAD ~ . - PERIODO - CLASE, family = "poisson",
data = trainSet_dia)
Deviance Residuals:
Min 1Q Median 3Q Max
-0.5004 -0.3002 -0.1348 -0.0450 4.3825
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -5.301e-02 2.586e-02 -2.050 0.040386 *
DIA_ANNO 4.649e-05 2.373e-05 1.960 0.050035 .
FECHAS_ESP2 8.677e-02 1.482e-02 5.856 4.73e-09 ***
FECHAS_ESP3 1.941e-02 1.538e-02 1.262 0.206780
COMUNA2 7.424e-03 2.922e-02 0.254 0.799474
COMUNA3 4.657e-02 2.454e-02 1.898 0.057713 .
COMUNA4 8.397e-02 2.261e-02 3.714 0.000204 ***
COMUNA5 2.364e-01 2.181e-02 10.840 < 2e-16 ***
COMUNA6 3.016e-02 2.500e-02 1.206 0.227648
COMUNA7 8.414e-02 2.245e-02 3.748 0.000178 ***
COMUNA8 3.847e-02 2.520e-02 1.527 0.126866
COMUNA9 6.436e-02 2.418e-02 2.662 0.007772 **
COMUNA10 3.965e-01 2.106e-02 18.826 < 2e-16 ***
COMUNA11 2.851e-01 2.158e-02 13.212 < 2e-16 ***
COMUNA12 7.631e-02 2.434e-02 3.135 0.001721 **
COMUNA13 1.057e-02 2.700e-02 0.392 0.695385
COMUNA14 2.008e-01 2.210e-02 9.088 < 2e-16 ***
COMUNA15 3.368e-01 2.226e-02 15.128 < 2e-16 ***
COMUNA16 1.338e-01 2.248e-02 5.950 2.68e-09 ***
BARRIO 5.602e-05 3.319e-05 1.688 0.091436 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 32478 on 129974 degrees of freedom
Residual deviance: 29309 on 129955 degrees of freedom
AIC: 307531
Number of Fisher Scoring iterations: 4En este caso tanto el INTERCEPTO como DIA_ANNO pasan a ser ligeramente significativas para el modelo, con un valor p menor al 5% que no es muy contundente para rechazar la hippótesis nula, pero es mejor que los resultados arrojados por los ajustes anteriores, tomando en cuenta el objetivo del ajuste, que es poder predecir la accidentalidad por día en Medellín con cierto grado de seguridad.
A continuación se calcula el MSE para los datos de entrenamiento y validación con cada uno de los cuatro ajustes anteriores. Además se calcula el porcentaje de variación para descartar sobreentrenamiento.
| Ajuste.Poisson | MSE.Entrenamiento | MSE.Prueba | Porcentaje.Variacion |
|---|---|---|---|
| Ajuste1 | 1.481775 | 1.474705 | 0.4771792 |
| Ajuste2 | 1.481776 | 1.483358 | 0.1067355 |
| Ajuste3 | 1.488708 | 1.488884 | 0.0118282 |
| Ajuste4 | 1.488714 | 1.505348 | 1.1173117 |
Los resultados obtenidos para los cuatro ajustes son buenos, teniendo en cuenta que la variación máxima del MSE caculado para ambos conjuntos debe ser de 15%. El cuarto ajuste presenta una variación de 1.11% y es la más alta. Este ajuste, aunque sacrifica la variable CLASE, permite ganar significancia a la variable DIA_ANNO, y por tanto se toma para realizar la comparación con los modelos que se estudiarán más adelante.
K Vecinos Cercanos
El segundo modelo a ajustar es el K-Vecinos Cercanos para Regresión. Este hace parte de la librería caret.
Antes de realizar el ajuste del modelo, se selecciona un valor de \(K\) óptimo de acuerdo al valor mínimo obtenido para el RMSE. Para hacerlo se tomará una muestra del 5% de los datos de entrenamiento ya que debido a la alta dimensionalidad de los datos conllevaría un alto costo computacional el poder realizar el proceso con los datos completos.
k-Nearest Neighbors
6499 samples
6 predictor
No pre-processing
Resampling: Repeated Train/Test Splits Estimated (15 reps, 75%)
Summary of sample sizes: 4875, 4875, 4875, 4875, 4875, 4875, ...
Resampling results across tuning parameters:
k RMSE Rsquared MAE
1 0.9334002 0.006538715 0.4541598
2 0.8131204 0.008727420 0.4407191
3 0.7737921 0.008535358 0.4381089
4 0.7507789 0.009421495 0.4376415
5 0.7366937 0.009539959 0.4374020
6 0.7258625 0.010976364 0.4357040
7 0.7205006 0.010629708 0.4356655
8 0.7170707 0.010076601 0.4366437
9 0.7130671 0.010679289 0.4365483
10 0.7096793 0.010973084 0.4356434
11 0.7073537 0.011038732 0.4355543
12 0.7053150 0.011019814 0.4352602
13 0.7034602 0.011152480 0.4351307
14 0.7023586 0.010713981 0.4353982
15 0.7011746 0.010602367 0.4356120
16 0.7006561 0.010007158 0.4361337
17 0.6989062 0.010731097 0.4356499
18 0.6977098 0.011028482 0.4357237
19 0.6967296 0.011145812 0.4362060
20 0.6957943 0.011460122 0.4365931
21 0.6951555 0.011520160 0.4369619
22 0.6946478 0.011517317 0.4375259
23 0.6938387 0.011909987 0.4379116
24 0.6936921 0.011621847 0.4385096
25 0.6930351 0.011975907 0.4385341
26 0.6930061 0.011656489 0.4391564
27 0.6924136 0.012069203 0.4391381
28 0.6920373 0.012280088 0.4395216
29 0.6914613 0.012803811 0.4398924
30 0.6913134 0.012765639 0.4401895
31 0.6910443 0.012868980 0.4405082
32 0.6909929 0.012695084 0.4409699
33 0.6908668 0.012614822 0.4412837
34 0.6907494 0.012570141 0.4415893
35 0.6907018 0.012405088 0.4419263
36 0.6907884 0.012047977 0.4424350
37 0.6906026 0.012055169 0.4425381
38 0.6905827 0.011833607 0.4428213
39 0.6904841 0.011709709 0.4430212
40 0.6904212 0.011533032 0.4432252
41 0.6903585 0.011446223 0.4434486
42 0.6902230 0.011383592 0.4435179
43 0.6901213 0.011324710 0.4436829
44 0.6902502 0.010920146 0.4440078
45 0.6902102 0.010837749 0.4441532
46 0.6901407 0.010765209 0.4442043
47 0.6900967 0.010697535 0.4442710
48 0.6900337 0.010674487 0.4443658
49 0.6900974 0.010454928 0.4445716
50 0.6900256 0.010438456 0.4446361
RMSE was used to select the optimal model using the smallest value.
The final value used for the model was k = 50.
De acuerdo con el modelo, el valor óptimo de \(k\) es \(k=50\), sin embargo, se observa que el mse no varía mucho a partir de \(k=15\), por lo que se toma este valor para el entrenamiento posterior.
15-nearest neighbor regression modelSe realiza la predicción para los datos de entrenamiento y se calcula el MSE.
[1] 0.3695411Para el conjunto de prueba:
[1] 0.4236879Se calcula la diferencia porcentual.
[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Vecinos cercanos es de 14.65%"Los valores de MSE obtenidos en este con el modelo de Vecinos Cercanos son más bajos que para el modelo anterior. En este caso se calcula una diferencia en el MSE de 14.9%, que se encuentra en el límite. Se puede decir entonces que el modelo no está sobreentrenado, aunque la diferencia entre ambos valores de MSE es alta.
Bosque Aleatorio
Antes de entrenar el modelo, se debe seleccionar un número óptimo de árboles. Para esto se toma una muestra del 10% de los datos de entrenamiento y se realiza el siguiente proceso, guiado de lo encontrado en este RPubs.
A partir de 100 árboles el error no varía mucho. Se selecciona este valor para entrenar el modelo con los datos de entrenamiento.
Call:
randomForest(formula = CANTIDAD ~ ., data = trainSet_dia, ntree = 100)
Type of random forest: regression
Number of trees: 100
No. of variables tried at each split: 2
Mean of squared residuals: 0.3378411
% Var explained: 22.77El ajuste del modelo con 50 árboles solo explica el 22.94% de la varianza, lo cual es un valor bastamnnte bajo y puede indicar que al modelo le hacen falta variables que aporten más información. A continuación se calcula el mse para los datos de entrenamiento y prueba.
[1] "El MSE para el conjunto de entrenamiento es: 0.308830367506836"[1] "El MSE para el conjunto de prueba es: 0.349864891944701"[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Bosque Aleatorio entrenado con 20 árboles es de 13.29%"La variación del MSE es de 13.42% menor que la variación obtenida para el ajuste anterior y que el límite establecido, por lo que el modelo no se encuentra sobreentrenado. El ajuste del modelo con 100 árboles conlleva un costo computacional bastante alto, puesto que el ajustar el modelo demora al menos 30 minutos.
Se muestra a continuación la importancia de las variables en aras de saber qué tan influyente es la variable DIA_ANNO para la predicción de la accidentalidad.
Se observa que la variable COMUNA es la de mayor influencia, seguida de BARRIO, CLASE y DIA_ANNO. Las variables FECHAS_ESP y PERIODO no son tan inlfuyentes para el modelo. Se realiza un segundo ajuste sin contar estas variables.
Call:
randomForest(formula = CANTIDAD ~ . - FECHAS_ESP - PERIODO, data = trainSet_dia, ntree = 100)
Type of random forest: regression
Number of trees: 100
No. of variables tried at each split: 1
Mean of squared residuals: 0.3595494
% Var explained: 17.8Se observa una disminución considerable en el porcentaje de varianza explicada por el modelo, que para el primer ajuste también era baja. Esto indica que pueden hacer falta variables dentro del ajuste que permitan alcanzar un mayor valor de varianza explicada y permitan mejores predicciones de la accidentalidad. A continuación se calcula el MSE para los conjuntos de entrenamiento y validación.
[1] "El MSE para el conjunto de prueba es: 0.35577046003082"[1] "El MSE para el conjunto de prueba es: 0.373205648304266"[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Bosque Aleatorio entrenado con 20 árboles es de 4.9%"Aunque los valores de MSE en este segundo ajuste son más alto, la diferencia entre el MSE calculado para ambos conjuntos fue más baja que la obtenida con el modelo ajustado con todas las variables. Se toma este segundo ajuste para realizar la comparación entre los tres modelos.
A continuación se presenta un resumen de los resultados obtenidos para los tres modelos.
| Modelo | mse_train | mse_test | Diff_mse |
|---|---|---|---|
| GLM(Poisson) | 1.4887142 | 1.5053478 | 1.12% |
| KNN | 0.3695411 | 0.4236879 | 14.65% |
| Random Forest | 0.3557705 | 0.3732056 | 4.9% |
Se observa que los valores más altos de MSE se obtuvieron para el GLM (Poisson), seguido de KNN y Random Forest. Sin embargo, la diferencia entre dicho valor para ambos conjuntos aumenta, siendo GLM (Poisson) el que arrojó la menor diferencia y KNN el de mayor.
El error cuadrático para los datos de validación es bastante alto con el modelo GLM (Poisson), tomando en cuenta que la cantidad máxima de accidentes de acuerdo a los filtros introducidos sería de 9, y que se tienen muchos datos únicos. Por otro lado, el modelo KNN muestra una diferencia de MSE bastante cercana al límite, por lo que es propenso a estar sobreestimado, además de que el tiempo necesario para el ajuste del modelo es un poco alto. En este caso, el mejor modelo para la predicción de la accidentalidad por días sería el modelo de Bosque Aleatorio, el cual no solo muestra una diferencia en el MSE bastante baja, si no que los valores calculados de MSE son bajos individualmente.
Modelamiento de accidentalidad por semana
Se realiza un proceso similar al hecho para los datos de accidentalidad por día, solo que ahora se excluye la variable DIA del conteo. Ya que la variable FECHAS_ESP hace referencia a un díua puntual, no tiene sentido incluirlo en esta parte del análisis.
| PERIODO | MES | SEMANA | COMUNA | BARRIO | CLASE | CANTIDAD |
|---|---|---|---|---|---|---|
| 2014 | 1 | 1 | 1 | 49 | 4 | 1 |
| 2014 | 1 | 1 | 1 | 66 | 3 | 1 |
| 2014 | 1 | 1 | 1 | 101 | 1 | 1 |
| 2014 | 1 | 1 | 1 | 120 | 2 | 1 |
| 2014 | 1 | 1 | 1 | 184 | 1 | 1 |
| 2014 | 1 | 1 | 1 | 206 | 1 | 1 |
Se observa la distribución de la cantidad de accidentes por semanas.
Nuevamente, se observa un alto número de accidentes únicos, esto debido a que el conteo se encuentra discretizado por diferentes variables. Si se observa la relación entre la accidentalidad con la variable semana.
No se observa un patrón definido al graficar la accidentalidad por semana. Se procede ahora a dividir los datos en entrenamiento (2014 a 2017) y prueba (datos del 2018).
[1] "El conjunto de entrenamiento tiene 79181 filas y 7 columnas"[1] "El conjunto de prueba tiene 18741 filas y 7 columnas"Se ajustan los mismos tres modelos que se entrenaron en la sección anterior, esta vez para predecir la accidentalidad en Medellín por semana.
Modelo lineal generalizado - Distribución Poisson
Se ajusta el modelo con todas las variables para los datos de entrenamiento.
Call:
glm(formula = CANTIDAD ~ ., family = "poisson", data = trainSet_sem)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.4163 -0.6908 -0.0661 0.2533 8.0154
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -1.405e+01 4.454e+00 -3.155 0.001606 **
PERIODO 6.802e-03 2.210e-03 3.078 0.002084 **
MES2 6.021e-02 1.464e-02 4.113 3.90e-05 ***
MES3 7.872e-02 1.933e-02 4.072 4.66e-05 ***
MES4 5.617e-02 2.616e-02 2.147 0.031769 *
MES5 7.764e-02 3.325e-02 2.335 0.019531 *
MES6 2.682e-02 4.041e-02 0.664 0.506871
MES7 3.022e-02 4.828e-02 0.626 0.531389
MES8 6.124e-02 5.585e-02 1.096 0.272907
MES9 7.527e-02 6.347e-02 1.186 0.235665
MES10 5.255e-02 7.165e-02 0.733 0.463294
MES11 1.672e-02 7.917e-02 0.211 0.832742
MES12 1.552e-02 8.687e-02 0.179 0.858194
SEMANA 1.220e-03 1.811e-03 0.674 0.500406
COMUNA2 -6.332e-02 2.923e-02 -2.166 0.030277 *
COMUNA3 1.849e-01 2.454e-02 7.535 4.89e-14 ***
COMUNA4 2.850e-01 2.262e-02 12.599 < 2e-16 ***
COMUNA5 6.267e-01 2.181e-02 28.739 < 2e-16 ***
COMUNA6 1.683e-01 2.500e-02 6.732 1.67e-11 ***
COMUNA7 2.673e-01 2.248e-02 11.891 < 2e-16 ***
COMUNA8 8.600e-02 2.521e-02 3.411 0.000647 ***
COMUNA9 1.375e-01 2.418e-02 5.687 1.29e-08 ***
COMUNA10 1.051e+00 2.106e-02 49.926 < 2e-16 ***
COMUNA11 7.492e-01 2.159e-02 34.706 < 2e-16 ***
COMUNA12 1.630e-01 2.435e-02 6.691 2.21e-11 ***
COMUNA13 -1.670e-02 2.701e-02 -0.618 0.536391
COMUNA14 4.467e-01 2.214e-02 20.173 < 2e-16 ***
COMUNA15 8.652e-01 2.225e-02 38.885 < 2e-16 ***
COMUNA16 3.131e-01 2.248e-02 13.926 < 2e-16 ***
BARRIO 2.003e-04 3.340e-05 5.996 2.02e-09 ***
CLASE2 -6.179e-02 1.145e-02 -5.396 6.83e-08 ***
CLASE3 8.481e-01 8.416e-03 100.773 < 2e-16 ***
CLASE4 -2.307e-02 1.087e-02 -2.123 0.033795 *
CLASE5 -2.192e-01 1.589e-02 -13.792 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 117403 on 79180 degrees of freedom
Residual deviance: 67757 on 79147 degrees of freedom
AIC: 258202
Number of Fisher Scoring iterations: 5En este caso, la variable SEMANA no es significativa para el ajuste del modelo. Se observa, además, que en promedio solo la adición de los meses de febrero y marzo aporta información para el ajuste. En este caso, la mayoría de las comunas son significativas para el ajuste, como muestran sus valores p, además de que para la variable BARRIO dicho valor es bastante bajo, permitiendo rechazar la hipótesis nula.
Solo la clase INCENDIO no es significativa para el modelo. Cabe destacar que la ocurrencia de incendios es muy baja, en comparación con las otras clases de accidentes. La variable PERIODO tiene una influencia sobre la variable de respuesta en este caso, a diferencia del ajuste del modelo para predicción de accidentalidad por días.
Se realiza un segundo ajuste del modelo eliminando la columna MES para observar el comportamiento de la variable SEMANA respecto a la variable de respuesta.
Call:
glm(formula = CANTIDAD ~ . - MES, family = "poisson", data = trainSet_sem)
Deviance Residuals:
Min 1Q Median 3Q Max
-2.3897 -0.6868 -0.0630 0.2570 8.0993
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -1.394e+01 4.454e+00 -3.131 0.001744 **
PERIODO 6.777e-03 2.210e-03 3.067 0.002166 **
SEMANA 8.784e-04 1.653e-04 5.316 1.06e-07 ***
COMUNA2 -6.441e-02 2.923e-02 -2.204 0.027527 *
COMUNA3 1.849e-01 2.454e-02 7.535 4.90e-14 ***
COMUNA4 2.844e-01 2.262e-02 12.572 < 2e-16 ***
COMUNA5 6.263e-01 2.181e-02 28.722 < 2e-16 ***
COMUNA6 1.678e-01 2.500e-02 6.711 1.94e-11 ***
COMUNA7 2.672e-01 2.248e-02 11.888 < 2e-16 ***
COMUNA8 8.560e-02 2.521e-02 3.396 0.000684 ***
COMUNA9 1.373e-01 2.418e-02 5.678 1.37e-08 ***
COMUNA10 1.051e+00 2.106e-02 49.891 < 2e-16 ***
COMUNA11 7.492e-01 2.159e-02 34.707 < 2e-16 ***
COMUNA12 1.634e-01 2.435e-02 6.708 1.97e-11 ***
COMUNA13 -1.703e-02 2.701e-02 -0.630 0.528419
COMUNA14 4.467e-01 2.214e-02 20.172 < 2e-16 ***
COMUNA15 8.651e-01 2.225e-02 38.882 < 2e-16 ***
COMUNA16 3.133e-01 2.248e-02 13.934 < 2e-16 ***
BARRIO 2.006e-04 3.340e-05 6.006 1.90e-09 ***
CLASE2 -6.114e-02 1.145e-02 -5.339 9.34e-08 ***
CLASE3 8.476e-01 8.415e-03 100.720 < 2e-16 ***
CLASE4 -2.291e-02 1.087e-02 -2.108 0.035012 *
CLASE5 -2.186e-01 1.589e-02 -13.755 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 117403 on 79180 degrees of freedom
Residual deviance: 67860 on 79158 degrees of freedom
AIC: 258283
Number of Fisher Scoring iterations: 5Se observa una disminución considerable en el valor p calculado para la variable SEMANA en comparación con el ajuste anterior. Este valor es lo suficientemente contundente para rechazar la hipótesis nula de entrada. Además, el ajuste muestra que las otras variables son bastante significativas para el modelo, por lo que no se realizan más modificaciones.
A continuación se calcula el MSE para los datos de entrenamiento y validación con cada uno de los dos ajustes anteriores. Además se calcula el porcentaje de variación para descartar sobreentrenamiento.
| Ajuste_Pois | MSE_train | MSE_test | variacion_percent |
|---|---|---|---|
| Ajuste1 | 6.326080 | 6.424384 | 1.553946 |
| Ajuste2 | 6.327083 | 6.424770 | 1.543949 |
Se observa que el modelo se desempeña ligeramente mejor en con los datos de entrenamiento para ambos ajustes. La diferencia porcentual entre los valores calculados de MSE es menor para el segundo ajuste del modelo.
Debido a que en el segundo ajuste la variable SEMANA muestra cierto grado de significancia, se toma este como el modelo a comparar con el ajuste de K vecinos cercanos y Random Forest para la predicción de la accidentalidad por semana en Medellín.
K Vecinos Cercanos
Se selecciona un valor de \(K\) óptimo de acuerdo al valor mínimo obtenido para el MSE. Para hacerlo se tomará una muestra del 5% de los datos de entrenamiento ya que debido a la alta dimensionalidad de los datos conllevaría un alto costo computacional el poder realizar el proceso con los datos completos.
k-Nearest Neighbors
3959 samples
6 predictor
No pre-processing
Resampling: Repeated Train/Test Splits Estimated (15 reps, 75%)
Summary of sample sizes: 2970, 2970, 2970, 2970, 2970, 2970, ...
Resampling results across tuning parameters:
k RMSE Rsquared MAE
1 3.001184 0.03769401 1.613398
2 2.635161 0.03047857 1.540418
3 2.487686 0.03576602 1.496311
4 2.414659 0.03753239 1.476593
5 2.385950 0.03411782 1.469324
6 2.363848 0.03308873 1.458918
7 2.346175 0.03196413 1.451412
8 2.334790 0.03108350 1.446250
9 2.325897 0.03034868 1.443707
10 2.319206 0.02985139 1.439485
11 2.312266 0.03023662 1.434445
12 2.305698 0.03072253 1.431016
13 2.299844 0.03111351 1.427637
14 2.295617 0.03134484 1.423914
15 2.291884 0.03123249 1.420012
16 2.291529 0.03049281 1.419608
17 2.288720 0.03074631 1.418588
18 2.287126 0.03054187 1.416569
19 2.281934 0.03218322 1.412781
20 2.279083 0.03264578 1.411426
21 2.276440 0.03336310 1.410526
22 2.274549 0.03376623 1.409412
23 2.273857 0.03367281 1.408302
24 2.272238 0.03406435 1.407398
25 2.270635 0.03480575 1.406761
26 2.270455 0.03436022 1.405604
27 2.270789 0.03376631 1.405356
28 2.270544 0.03351337 1.403815
29 2.272263 0.03222512 1.404756
30 2.273603 0.03100777 1.404880
31 2.275655 0.02953304 1.406156
32 2.275466 0.02933853 1.405573
33 2.275797 0.02877310 1.405611
34 2.275963 0.02826992 1.404891
35 2.275714 0.02813160 1.404392
36 2.275558 0.02805573 1.403971
37 2.275870 0.02756844 1.403130
38 2.276843 0.02672998 1.403217
39 2.276346 0.02687494 1.402184
40 2.275968 0.02689107 1.401569
41 2.275709 0.02699833 1.400305
42 2.275339 0.02705483 1.399660
43 2.275092 0.02693672 1.399269
44 2.274968 0.02683409 1.399151
45 2.275379 0.02625087 1.399127
46 2.275886 0.02582075 1.399230
47 2.276439 0.02529847 1.398847
48 2.277504 0.02450816 1.399565
49 2.277879 0.02413382 1.399504
50 2.278369 0.02372722 1.399600
51 2.278996 0.02325473 1.399779
52 2.279086 0.02311917 1.399344
53 2.279641 0.02264163 1.399624
54 2.279308 0.02270464 1.399034
55 2.279515 0.02251350 1.399162
56 2.280166 0.02194077 1.399214
57 2.279797 0.02205476 1.398599
58 2.279470 0.02221728 1.398362
59 2.279624 0.02207214 1.398176
60 2.279989 0.02170410 1.397964
61 2.279994 0.02163697 1.397752
62 2.280380 0.02127202 1.398006
63 2.281370 0.02053316 1.398300
64 2.281799 0.02017349 1.398601
65 2.282121 0.01993692 1.398967
66 2.282078 0.01984716 1.398619
67 2.282621 0.01951776 1.398872
68 2.282891 0.01924469 1.398714
69 2.282893 0.01914500 1.398668
70 2.283051 0.01898761 1.398860
71 2.283406 0.01872425 1.399191
72 2.283584 0.01843071 1.399114
73 2.283886 0.01824773 1.399090
74 2.283981 0.01817953 1.398908
75 2.284186 0.01794812 1.398969
76 2.284322 0.01780906 1.398784
77 2.284588 0.01747893 1.399108
78 2.284597 0.01737589 1.398915
79 2.284286 0.01753510 1.398804
80 2.284452 0.01736055 1.398969
81 2.284819 0.01705330 1.399265
82 2.284852 0.01700616 1.399495
83 2.284993 0.01686518 1.399682
84 2.284882 0.01693905 1.399845
85 2.284946 0.01683840 1.399849
86 2.285324 0.01653430 1.400438
87 2.285426 0.01642144 1.400683
88 2.285312 0.01642023 1.400874
89 2.285438 0.01628096 1.401063
90 2.285226 0.01637443 1.400968
91 2.285177 0.01631187 1.401018
92 2.285104 0.01628496 1.401298
93 2.285269 0.01608434 1.401547
94 2.285331 0.01601172 1.402109
95 2.284846 0.01633256 1.402115
96 2.284527 0.01655361 1.402001
97 2.284549 0.01650585 1.402134
98 2.284700 0.01638537 1.402132
99 2.284422 0.01653963 1.401916
100 2.284303 0.01662658 1.402127
RMSE was used to select the optimal model using the smallest value.
The final value used for the model was k = 26.
Según el modelo, el valor óptimo de vecinos para este escenario es de 78, sin embargo, la disminución del error no es muy pronunciada a partir de los 10 vecinos, por lo que se toma \(k=20\) para realizar el ajuste del modelo.
20-nearest neighbor regression modelSe calcula el MSE para los datos de entrenamiento:
[1] 3.746555Ahora para los datos de prueba:
[1] 4.306304[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Vecinos cercanos es de 14.94%"Se tienen una diferencia del 23.3% entre el MSE para los datos de entrenamiento y validación. Dicha diferencia está por encima del 15% establecido como límite, por lo que en este caso hay evidencia de sobreentrenamiento del modelo. Se entrena nuevamente el modelo con \(k=40\) vecinos para observar si el desempeño mejora.
40-nearest neighbor regression modelAhora se calcula el MSE para el conjunto de entrenamiento:
[1] 4.112394El MSE para el conjunto de prueba es:
[1] 4.411655[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Vecinos cercanos es de 7.28%"En este caso la diferencia entre el desempeño del modelo en los conjuntos de entrenamiento y prueba es de solo 8.51%, más bajo que el valor obtenido al ajustar el modelo con 20 vecinos y que el valor límite (15%). Los valores calculados de MSE son más altos que los obtenidos con el modelo Poisson en este caso.
Se toma este ajuste para las comparaciones con los otros dos modelos.
Bosque Aleatorio
Se realiza el proceso para seleccionar un número óptimo de árboles para predecir la accidentalidad por semana tomando una muestra del 10% de los datos de entrenamiento.
Alrededor de los 100 árboles no se observa una disminución significativa del error. Se toma este valor y se ajusta el modelo perdictivo.
Call:
randomForest(formula = CANTIDAD ~ ., data = trainSet_sem, ntree = 100)
Type of random forest: regression
Number of trees: 100
No. of variables tried at each split: 2
Mean of squared residuals: 2.029293
% Var explained: 59.68El modelo ajustado con 100 árboles explica un 59.59% de la varianza, puede que hagan falta variables para poder aumentar el porcentaje de varianza explicada por el modelo. Se realiza ahora la predicción para los conjuntos de entrenamiento y prueba y se calcula el MSE.
[1] "El MSE para el conjunto de prueba es: 1.5238620671935"[1] "El MSE para el conjunto de prueba es: 2.08857471076245"[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Bosque Aleatorio entrenado con 150 árboles es de 37.06%"La diferencia entre el MSE para ambos conjuntos está por encima del límite (15%), indicando que el modelo ajustado con todas las variables está sobreentrenado. El modelo se desempeña mejor para el conjunto de entrenamiento. Además, los valores de MSE son más bajos que los obtenidos para el método de Vecinos Cercanos.
Se observa la importancia de las variables en el ajuste.
En este caso las variables MES y PERIODO son las de menor importancia. Se procede a realizar un segundo ajuste sin contar estas variables para observar el desempeño del modelo de Bosque Aleatorio.
Call:
randomForest(formula = CANTIDAD ~ . - PERIODO - MES, data = trainSet_sem, ntree = 100)
Type of random forest: regression
Number of trees: 100
No. of variables tried at each split: 1
Mean of squared residuals: 3.028943
% Var explained: 39.82Este ajuste explica un 41.21% de la varianza, menos que el ajuste anterior. Se procede a calcular su desempeño con los conjuntos de entrenamiento y prueba calculando el MSE.
[1] "El MSE para el conjunto de prueba es: 2.98749733700832"[1] "El MSE para el conjunto de prueba es: 3.17470557598906"[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Bosque Aleatorio sin la variable PERIODO es de 6.27%"En este caso la diferencia entre el MSE obtenido para los conjuntos de entrenamiento y validación es de solo 6.46%, más bajo que el obtenido para el ajuste anterior. Se toma este ajuste para comparar con los mejores ajustes de los otros dos modelos, ya que la eliminación de las dos variables redujo signficativamente la diferencia de MSE.
A continuación se presenta una tabla con los resultados de los ajustes para los tres modelos.
| Modelo | mse_train | mse_test | Diff_mse |
|---|---|---|---|
| GLM(Poisson) | 6.327083 | 6.424770 | 1.54% |
| KNN | 4.112394 | 4.411655 | 7.28% |
| Random Forest | 2.987497 | 3.174706 | 6.27% |
El MSE obtenido para ambos conjuntos con GLM(Poisson) es el más alto. Sin embargo, la diferencia entre ambos valores fue baja (1.54%) en comparación a la obtenida con los otros modelos. Los valores más bajos de MSE se obtuvieron para el modelo de Bosque Aleatorio, el cual se selecciona para realizar la predicción de la accidentalidad por semana en Medellín.
Modelamiento de accidentalidad por mes
Por último, se procede a entrenar un modelo para predecir la accidentalidad por mes. Se crea la tabla al igual que en los casos anteriores y se particionan los datos en los conjuntos de entrenamiento y validación.
| PERIODO | MES | COMUNA | BARRIO | CLASE | CANTIDAD |
|---|---|---|---|---|---|
| 2014 | 1 | 1 | 1 | 3 | 1 |
| 2014 | 1 | 1 | 1 | 4 | 1 |
| 2014 | 1 | 1 | 49 | 2 | 1 |
| 2014 | 1 | 1 | 49 | 4 | 1 |
| 2014 | 1 | 1 | 66 | 3 | 1 |
| 2014 | 1 | 1 | 101 | 1 | 3 |
[1] "El conjunto de entrenamiento tiene 36199 filas y 6 columnas"[1] "El conjunto de prueba tiene 8696 filas y 6 columnas"Se observa cómo está distribuida la accidentalidad en este caso.
La mayoría de los accidentes en la tabla son únicos, debido a que estos se discretizan por diferentes variables, sin embargo, se puede observar que la ocurrencia de accidentes varía entre 1 y 64 en esta ocasión, sebido a que la resolución temporal es mayor.
Modelo lineal generalizado - Distribución Poisson
Se procede a entrenar el modelo con todas las variables inicialmente.
Call:
glm(formula = CANTIDAD ~ ., family = "poisson", data = trainSet_mes)
Deviance Residuals:
Min 1Q Median 3Q Max
-6.5406 -0.9419 -0.2078 0.5757 12.4788
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -1.078e+01 4.450e+00 -2.423 0.015395 *
PERIODO 5.210e-03 2.208e-03 2.360 0.018285 *
MES2 1.016e-01 1.268e-02 8.012 1.13e-15 ***
MES3 1.548e-01 1.246e-02 12.427 < 2e-16 ***
MES4 1.099e-01 1.263e-02 8.704 < 2e-16 ***
MES5 1.707e-01 1.245e-02 13.704 < 2e-16 ***
MES6 8.507e-02 1.270e-02 6.696 2.14e-11 ***
MES7 1.498e-01 1.248e-02 12.000 < 2e-16 ***
MES8 1.858e-01 1.238e-02 15.016 < 2e-16 ***
MES9 1.611e-01 1.242e-02 12.969 < 2e-16 ***
MES10 1.524e-01 1.248e-02 12.218 < 2e-16 ***
MES11 9.417e-02 1.266e-02 7.438 1.02e-13 ***
MES12 1.354e-01 1.256e-02 10.775 < 2e-16 ***
COMUNA2 -1.089e-01 2.922e-02 -3.727 0.000194 ***
COMUNA3 4.968e-01 2.452e-02 20.263 < 2e-16 ***
COMUNA4 7.360e-01 2.261e-02 32.556 < 2e-16 ***
COMUNA5 1.183e+00 2.175e-02 54.383 < 2e-16 ***
COMUNA6 5.133e-01 2.499e-02 20.545 < 2e-16 ***
COMUNA7 6.196e-01 2.244e-02 27.615 < 2e-16 ***
COMUNA8 2.598e-01 2.520e-02 10.307 < 2e-16 ***
COMUNA9 4.425e-01 2.410e-02 18.356 < 2e-16 ***
COMUNA10 1.788e+00 2.103e-02 84.994 < 2e-16 ***
COMUNA11 1.369e+00 2.154e-02 63.546 < 2e-16 ***
COMUNA12 5.188e-01 2.430e-02 21.347 < 2e-16 ***
COMUNA13 -2.993e-02 2.701e-02 -1.108 0.267789
COMUNA14 9.164e-01 2.208e-02 41.500 < 2e-16 ***
COMUNA15 1.511e+00 2.222e-02 68.010 < 2e-16 ***
COMUNA16 7.062e-01 2.242e-02 31.491 < 2e-16 ***
BARRIO 3.366e-04 3.356e-05 10.032 < 2e-16 ***
CLASE2 -1.096e-01 1.144e-02 -9.584 < 2e-16 ***
CLASE3 1.519e+00 8.376e-03 181.308 < 2e-16 ***
CLASE4 1.392e-02 1.085e-02 1.284 0.199301
CLASE5 -5.761e-01 1.588e-02 -36.282 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 231920 on 36198 degrees of freedom
Residual deviance: 88018 on 36166 degrees of freedom
AIC: 190966
Number of Fisher Scoring iterations: 5Para este modelo, salvo la COMUNA13 y la clase de accidente “Otro”, todas las variables son significativas, con valores p lo suficientemente bajos para poder rechazar la hipótesis nula y poder concluir que cada variable aporta información. Sin embargo, el PERIODO parace no ser muy significativo para el modelo. Se ajusta un segundo modelo sin contar dicha variable.
Call:
glm(formula = CANTIDAD ~ . - PERIODO, family = "poisson", data = trainSet_mes)
Deviance Residuals:
Min 1Q Median 3Q Max
-6.5515 -0.9409 -0.2085 0.5752 12.4361
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.812e-01 2.397e-02 -11.734 < 2e-16 ***
MES2 1.017e-01 1.268e-02 8.017 1.08e-15 ***
MES3 1.548e-01 1.246e-02 12.427 < 2e-16 ***
MES4 1.099e-01 1.263e-02 8.704 < 2e-16 ***
MES5 1.707e-01 1.245e-02 13.705 < 2e-16 ***
MES6 8.506e-02 1.270e-02 6.695 2.15e-11 ***
MES7 1.498e-01 1.248e-02 12.003 < 2e-16 ***
MES8 1.859e-01 1.238e-02 15.016 < 2e-16 ***
MES9 1.611e-01 1.242e-02 12.974 < 2e-16 ***
MES10 1.525e-01 1.248e-02 12.224 < 2e-16 ***
MES11 9.426e-02 1.266e-02 7.445 9.72e-14 ***
MES12 1.354e-01 1.256e-02 10.779 < 2e-16 ***
COMUNA2 -1.090e-01 2.922e-02 -3.731 0.00019 ***
COMUNA3 4.967e-01 2.452e-02 20.261 < 2e-16 ***
COMUNA4 7.360e-01 2.261e-02 32.554 < 2e-16 ***
COMUNA5 1.183e+00 2.175e-02 54.387 < 2e-16 ***
COMUNA6 5.135e-01 2.499e-02 20.551 < 2e-16 ***
COMUNA7 6.198e-01 2.244e-02 27.622 < 2e-16 ***
COMUNA8 2.597e-01 2.520e-02 10.306 < 2e-16 ***
COMUNA9 4.425e-01 2.410e-02 18.359 < 2e-16 ***
COMUNA10 1.788e+00 2.103e-02 84.995 < 2e-16 ***
COMUNA11 1.369e+00 2.154e-02 63.550 < 2e-16 ***
COMUNA12 5.188e-01 2.430e-02 21.348 < 2e-16 ***
COMUNA13 -2.985e-02 2.701e-02 -1.105 0.26903
COMUNA14 9.165e-01 2.208e-02 41.503 < 2e-16 ***
COMUNA15 1.511e+00 2.222e-02 68.012 < 2e-16 ***
COMUNA16 7.062e-01 2.242e-02 31.494 < 2e-16 ***
BARRIO 3.368e-04 3.356e-05 10.038 < 2e-16 ***
CLASE2 -1.096e-01 1.144e-02 -9.581 < 2e-16 ***
CLASE3 1.519e+00 8.376e-03 181.343 < 2e-16 ***
CLASE4 1.416e-02 1.085e-02 1.305 0.19177
CLASE5 -5.753e-01 1.588e-02 -36.239 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Dispersion parameter for poisson family taken to be 1)
Null deviance: 231920 on 36198 degrees of freedom
Residual deviance: 88024 on 36167 degrees of freedom
AIC: 190970
Number of Fisher Scoring iterations: 5Se observa una mejora en la interacción de las variables predictoras con la variable de respuesta para este segundo ajuste. Se toma el ajuste y se calcula el MSE para los datos de entrenamiento y prueba.
Se procede a realizar la predicción y realizar del MSE para los datos de entrenamiento:
[1] 62.23592Ahora para los datos de prueba:
[1] 62.56759Se calcula la diferencia porcentual entre ambos resultados.
[1] "La variación entre los mse calculados para los datos entrenamiento y prueba con el GLM con distribución Poisson es de 0.53%"La diferencia obtenida para los MSE en entrenamiento y prueba es mucho menor que el valor límite (15%). El ajuste realizado en este caso no sufre de sobreentrenamiento. Sin embargo, los valores de MSE obtenidos para el ajuste son bastante altos.
K Vecinos Cercanos
Se selecciona un valor de \(K\) óptimo de acuerdo al valor mínimo obtenido para el rmse. Para hacerlo se tomará una muestra del 10% de los datos de entrenamiento ya que debido a la alta dimensionalidad de los datos conllevaría un alto costo computacional el poder realizar el proceso con los datos completos.
k-Nearest Neighbors
3620 samples
5 predictor
No pre-processing
Resampling: Repeated Train/Test Splits Estimated (15 reps, 75%)
Summary of sample sizes: 2717, 2717, 2717, 2717, 2717, 2717, ...
Resampling results across tuning parameters:
k RMSE Rsquared MAE
1 5.881482 0.47689578 2.842152
2 5.650723 0.48197369 2.927454
3 5.901948 0.43194601 3.120463
4 6.175556 0.37780419 3.307284
5 6.413204 0.32965715 3.462802
6 6.600546 0.28965071 3.578220
7 6.700207 0.26766293 3.656508
8 6.764265 0.25423926 3.715264
9 6.840048 0.23726348 3.760408
10 6.906294 0.22230208 3.799990
11 6.948036 0.21247774 3.824835
12 6.983087 0.20414857 3.842669
13 7.011206 0.19766154 3.852327
14 7.040297 0.19102399 3.869068
15 7.064937 0.18512205 3.877648
16 7.071397 0.18411284 3.886201
17 7.089042 0.17988637 3.893887
18 7.119096 0.17321257 3.907959
19 7.150880 0.16537628 3.922706
20 7.164603 0.16219047 3.933974
21 7.182331 0.15813963 3.942950
22 7.199420 0.15414930 3.953537
23 7.217417 0.14979192 3.963063
24 7.238773 0.14449616 3.970477
25 7.248483 0.14220705 3.974324
26 7.264116 0.13864068 3.979141
27 7.281195 0.13435870 3.983443
28 7.289233 0.13241822 3.982845
29 7.301697 0.12943754 3.987380
30 7.313825 0.12645909 3.994797
31 7.321455 0.12438677 3.996652
32 7.335660 0.12077470 4.004292
33 7.348042 0.11762494 4.010268
34 7.354780 0.11589914 4.017496
35 7.363192 0.11377840 4.022122
36 7.372675 0.11150624 4.027065
37 7.389319 0.10732077 4.033379
38 7.396425 0.10554893 4.034875
39 7.399292 0.10480332 4.032829
40 7.403088 0.10393427 4.031510
41 7.407131 0.10298003 4.034567
42 7.410357 0.10230596 4.035781
43 7.416976 0.10066241 4.038557
44 7.418579 0.10038857 4.038284
45 7.422540 0.09948024 4.040323
46 7.425618 0.09887438 4.044156
47 7.429408 0.09798641 4.045104
48 7.435890 0.09643201 4.048512
49 7.440997 0.09510778 4.048771
50 7.448183 0.09322388 4.053025
RMSE was used to select the optimal model using the smallest value.
The final value used for the model was k = 2.
El comportamiento del rmse en este caso es distinto al obtenido para los ajustes realizado para las predicciones de accidentalidad por día y semana. En este caso, el error aumenta a medida que incrementa el número de vecinos y se encuentra que el k óptimo es \(k=2\). Sin embargo, al entrenar el modelo con dicho valor no se obtuvieron resultados concluyentes, por lo que se optó por utilizar un \(k=40\) para entrenar el modelo como sigue, ya que con este valor se logra que el error esté por debajo del 15%.
40-nearest neighbor regression modelSe realiza la predicción para los datos de entrenamiento y se calcula el MSE:
[1] 43.60028El MSE para el conjunto de prueba es:
[1] 48.0726Se calcula la diferencia porcentual.
[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Vecinos cercanos es de 10.26%"El desempeño del modelo en la predicción es mejor para el conjunto de entrenamiento. Los valores obtenidos para el MSE calculado para los conjuntos de entrenamiento y validación son mayores a los obtenidos para el modelo anterior. Sin embargo, la diferencia entre ambos (12.57%) permanece por debajo del límite, por lo que no se sufre de sobreentrenamiento en este caso.
Bosque Aleatorio
Se realiza el proceso para seleccionar un número óptimo de árboles para predecir la accidentalidad por mes tomando una muestra del 10% de los datos de entrenamiento.
Para este modelo se observa que a partir de 100 árboles el error comienza a variar muy poco, por lo que se toma este valor para el ajuste del modelo de predicción de la accidentalidad por mes.
Call:
randomForest(formula = CANTIDAD ~ ., data = trainSet_mes, ntree = 100)
Type of random forest: regression
Number of trees: 100
No. of variables tried at each split: 1
Mean of squared residuals: 28.59936
% Var explained: 50.12El ajuste del modelo con 100 árboles explica el 68.07% de la varianza. A continuación se calcula el MSE para los datos de entrenamiento y prueba para este ajuste.
[1] "El MSE para el conjunto de prueba es: 27.279714049818"[1] "El MSE para el conjunto de prueba es: 29.2152784617111"[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Bosque Aleatorio entrenado con 100 árboles es de 7.1%"La diferencia entre el MSE calculado para los conjuntos de entrenamiento y prueba es más alta que el valor límite, por lo que este ajuste se encuentra sobreentrenado. Se observa a continuación la importancia de las variables para este modelo.
En este caso las variables de menor importancia son MES y PERIODO. Se ajusta a continuación un modelo eliminando la variable PERIODO, para observar su desempeño.
Call:
randomForest(formula = CANTIDAD ~ . - PERIODO, data = trainSet_mes, ntree = 100)
Type of random forest: regression
Number of trees: 100
No. of variables tried at each split: 1
Mean of squared residuals: 25.9904
% Var explained: 54.67El ajuste en este caso explica 55.54% de la varianza, mayor que el ajuste anterior. Se observa el MSE para los conjuntos de entrenamiento y prueba.
[1] "El MSE para el conjunto de prueba es: 25.2200969980806"[1] "El MSE para el conjunto de prueba es: 26.083810971577"[1] "La variación entre los mse de entrenamiento y prueba para el modelo de Bosque Aleatorio sin la variable PERIODO es de 3.42%"En este caso la variación del MSE para los conjuntos de entrenamiento y prueba es de solo 3.08%, además de que los valores calculados de MSE son menores que los obtenidos para el primer ajuste. Se tiene entonces que esta diferencia se encuentra por debajo del límite, por lo que se toma este ajuste como el modelo final de Bosque Aleatorio.
A continuación se presenta una tabla con los resultados de los ajustes para los tres modelos.
| Modelo | mse_train | mse_test | Diff_mse |
|---|---|---|---|
| GLM(Poisson) | 62.23592 | 62.56759 | 0.53% |
| KNN | 43.60028 | 48.07260 | 10.26% |
| Random Forest | 25.22010 | 26.08381 | 3.42% |
Los valores más bajos de MSE se obtuvieron con el modelo Random Forest, que presentó un MSE más bajo que para los dos modelos anteriores. Se tomará como el modelo para predecir la accidentalidad por mes en la aplicación de Shiny.
Clusterización de los Barrios
Se realiza una copia del dataframe original.
| LATITUD | LONGITUD | FECHA | DIA | SEMANA | MES | FECHAS_ESP | PERIODO | DIA_NOMBRE | MES_NOMBRE | COMUNA | BARRIO | CLASE | GRAVEDAD |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 6.219016 | -75.60273 | 2014-01-01 | 1 | 1 | 1 | FESTIVO | 2014 | MIERCOLES | ENERO | BELEN | LOMA DE LOS BERNAL | CHOQUE | HERIDO |
| 6.260009 | -75.56818 | 2014-01-01 | 1 | 1 | 1 | FESTIVO | 2014 | MIERCOLES | ENERO | LA CANDELARIA | JESUS NAZARENO | ATROPELLO | HERIDO |
| 6.264765 | -75.54994 | 2014-01-01 | 1 | 1 | 1 | FESTIVO | 2014 | MIERCOLES | ENERO | MANRIQUE | MANRIQUE ORIENTAL | ATROPELLO | HERIDO |
| 6.234327 | -75.60761 | 2014-01-01 | 1 | 1 | 1 | FESTIVO | 2014 | MIERCOLES | ENERO | BELEN | LAS MERCEDES | ATROPELLO | HERIDO |
| 6.299968 | -75.57969 | 2014-01-01 | 1 | 1 | 1 | FESTIVO | 2014 | MIERCOLES | ENERO | DOCE DE OCTUBRE | DOCE DE OCTUBRE NO. 2 | CAIDA OCUPANTE | HERIDO |
Se realiza el cambio de la ñ por n para la columna gravedad.
Luego, se importan algunas librerias para establecer las mejores alternativas de clusterización.
El siguiente paso consiste en crear 4 tablas que van a ser necesarias para realizar la clusterización. Al final todas se unen con un full joint para crear los indicadores.
La primera tabla que se crea es la de barrio/Gravedad.
| BARRIO | GRAVEDAD | total |
|---|---|---|
| ALDEA PABLO VI | HERIDO | 64 |
| ALDEA PABLO VI | MUERTO | 2 |
| ALDEA PABLO VI | SOLO_DANOS | 17 |
| ALEJANDRIA | HERIDO | 132 |
| ALEJANDRIA | MUERTO | 1 |
La segunda tabla creada contiene la agrupación Barrio/Fechas especiales.
| BARRIO | FECHAS_ESP | total |
|---|---|---|
| ALDEA PABLO VI | FESTIVO | 5 |
| ALDEA PABLO VI | LABORAL | 44 |
| ALDEA PABLO VI | NO LABORAL | 34 |
| ALEJANDRIA | FESTIVO | 18 |
| ALEJANDRIA | LABORAL | 405 |
La tercera tabla contiene la agrupación de barrios y periodos.
| BARRIO | PERIODO | total |
|---|---|---|
| ALDEA PABLO VI | 2014 | 17 |
| ALDEA PABLO VI | 2015 | 15 |
| ALDEA PABLO VI | 2016 | 20 |
| ALDEA PABLO VI | 2017 | 11 |
| ALDEA PABLO VI | 2018 | 20 |
La ultima tabla contiene la agrupación barrio / clase (choque, atropello, caida ocupante, incendio, volcamiento, otro).
| BARRIO | CLASE | total |
|---|---|---|
| ALDEA PABLO VI | ATROPELLO | 27 |
| ALDEA PABLO VI | CAIDA OCUPANTE | 12 |
| ALDEA PABLO VI | CHOQUE | 35 |
| ALDEA PABLO VI | OTRO | 8 |
| ALDEA PABLO VI | VOLCAMIENTO | 1 |
Antes de unir todas las tablas, se transformas las variables de interés de cada una de filas a columnas.
Finalmente se unen todas las tablas para crear el consolidado.
| BARRIO | HERIDO | MUERTO | SOLO_DANOS | FESTIVO | LABORAL | NO LABORAL | 2014 | 2015 | 2016 | 2017 | 2018 | ATROPELLO | CAIDA OCUPANTE | CHOQUE | OTRO | VOLCAMIENTO |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| ALDEA PABLO VI | 64 | 2 | 17 | 5 | 44 | 34 | 17 | 15 | 20 | 11 | 20 | 27 | 12 | 35 | 8 | 1 |
| ALEJANDRIA | 132 | 1 | 366 | 18 | 405 | 76 | 78 | 64 | 105 | 123 | 129 | 14 | 18 | 434 | 25 | 8 |
| ALEJANDRO ECHAVARRIA | 595 | 3 | 243 | 41 | 558 | 242 | 159 | 145 | 187 | 181 | 169 | 98 | 84 | 483 | 125 | 51 |
| ALFONSO LOPEZ | 772 | 7 | 370 | 50 | 779 | 320 | 234 | 244 | 289 | 196 | 186 | 133 | 178 | 607 | 182 | 49 |
| ALTAMIRA | 504 | 2 | 270 | 14 | 598 | 164 | 216 | 162 | 137 | 142 | 119 | 50 | 144 | 466 | 100 | 16 |
Después se realiza el cambio de NA por cero para evitar inconsistencias.
Se crea el primer indicador que establece la proporción de accidentes los fines de semana y festivos vs el total de accidentes: (suma de accidentes en dias festivos y no laborales)/ suma total de accidentes.
El segundo indicador establece la proporción entre el número de accidentes fatales vs el total de accidentes: suma de accidentes faltales / suma total de accidentes.
El tercer indicador permite conocer el número de accidentes promedio por año: suma de los accidentes de todos los años / nro de años.
El ultimo indicador permite establecer que porcentaje de los accidentes que se presentaron fueron atropellos: suma de accidentes con atropellos / suma total de accidentes.
Finalmente se crea un subset para dejar el dataframe unicamente con las columnas necesarias.
| BARRIO | Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello |
|---|---|---|---|---|
| ALDEA PABLO VI | 0.4698795 | 0.0240964 | 16.6 | 0.3253012 |
| ALEJANDRIA | 0.1883768 | 0.0020040 | 99.8 | 0.0280561 |
| ALEJANDRO ECHAVARRIA | 0.3365042 | 0.0035672 | 168.2 | 0.1165279 |
| ALFONSO LOPEZ | 0.3220191 | 0.0060923 | 229.8 | 0.1157528 |
| ALTAMIRA | 0.2293814 | 0.0025773 | 155.2 | 0.0644330 |
Se reemplazan los na del dataframe cluster.
Para poder realizar los calculos se nombran las filas del dataframe con el nombre de los barrios y así poder eliminar esta columna.
| Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello | |
|---|---|---|---|---|
| ALDEA PABLO VI | 0.4698795 | 0.0240964 | 16.6 | 0.3253012 |
| ALEJANDRIA | 0.1883768 | 0.0020040 | 99.8 | 0.0280561 |
| ALEJANDRO ECHAVARRIA | 0.3365042 | 0.0035672 | 168.2 | 0.1165279 |
| ALFONSO LOPEZ | 0.3220191 | 0.0060923 | 229.8 | 0.1157528 |
| ALTAMIRA | 0.2293814 | 0.0025773 | 155.2 | 0.0644330 |
Finalmente se escalan los datos del dataframe para evitar problemas en el algorimo k-means por la dimensionalidad de los datos.
| Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello | |
|---|---|---|---|---|
| ALDEA PABLO VI | 2.5000455 | 2.7311436 | -0.8009166 | 2.0260755 |
| ALEJANDRIA | -1.6626592 | -0.6584602 | -0.3145554 | -1.1405583 |
| ALEJANDRO ECHAVARRIA | 0.5277657 | -0.4186245 | 0.0852897 | -0.1980435 |
| ALFONSO LOPEZ | 0.3135695 | -0.0312061 | 0.4453841 | -0.2063011 |
| ALTAMIRA | -1.0563049 | -0.5704978 | 0.0092957 | -0.7530254 |
Nuevamente se hace un resumen de los indicadores a utilizar la para la clusterización:
- Promedio anual de accidentes: Suma de los accidentes desde el 2014 a 2018 dividido 5 (nro de años)
- Tasa de accidentes fines de semana y festivos: accidentes fines de semana y festivos sobre el total de accidentes
- Tasa de fatalidad: accidentes fatales sobre el total de accidentes
- tasa de atropellos: atropellos sobre el total de accidentes
Una forma sencilla para encontrar el k óptimo de clusters es aplicar el algoritmo de k-means para un rango de valores de K. Luego de esto, se procede a identificar aquel valor a partir del cual la reducción en la suma total de la varianza intra-cluster deja de ser significativa (médoto del códo).
Método del codo: Según este método y la interpretación que se da de la gráfica, el k óptimo es igual a 5.
Metodo Silhouette: el número óptimo de clusters (aquel que maximiza la media del slhouette coeficient de todas las observaciones) se logra con k=2.
Según el método clusgap, el k óptimo estaría en 6.
A continuación se crean diferentes clusters variando K, con el objetivo de ver la representación espacial de diversas formas.
Luego, se genera un dataframe con las agrupaciones con diferentes k.
| Cl1G | Cl2G | Cl3G | Cl4G | Cl5G | Cl6G | Cl7G | Cl8G | Cl9G | C10G | |
|---|---|---|---|---|---|---|---|---|---|---|
| ALDEA PABLO VI | 1 | 1 | 2 | 1 | 4 | 2 | 6 | 6 | 5 | 6 |
| ALEJANDRIA | 2 | 2 | 1 | 4 | 1 | 6 | 3 | 4 | 8 | 4 |
| ALEJANDRO ECHAVARRIA | 2 | 2 | 1 | 4 | 3 | 1 | 7 | 7 | 9 | 7 |
| ALFONSO LOPEZ | 2 | 2 | 1 | 4 | 3 | 1 | 7 | 7 | 9 | 2 |
| ALTAMIRA | 2 | 2 | 1 | 4 | 1 | 6 | 3 | 4 | 8 | 4 |
Después de utlizar varios métodos para encontrar el mejor valor de k, se decide trabajar con k=5.
[1] 111 28 90 16 20Para poder crear el relato de los grupos, se observan las medias de cada una de las variables analizadas.
| cluster5$cluster | Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello |
|---|---|---|---|---|
| 1 | 0.2565408 | 0.0043729 | 135.92432 | 0.0690129 |
| 2 | 0.2320749 | 0.0060653 | 568.87143 | 0.0735121 |
| 3 | 0.3444111 | 0.0057171 | 93.98667 | 0.1709267 |
| 4 | 0.3629875 | 0.0248182 | 34.37500 | 0.2332957 |
| 5 | 0.3968397 | 0.0050750 | 34.09000 | 0.3485666 |
Con el proposito de tener más información, se agrega a la tabla de medias el tamaño de los grupos.
| tamano_grupo | cluster5$cluster | Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello |
|---|---|---|---|---|---|
| 111 | 1 | 0.2565408 | 0.0043729 | 135.92432 | 0.0690129 |
| 28 | 2 | 0.2320749 | 0.0060653 | 568.87143 | 0.0735121 |
| 90 | 3 | 0.3444111 | 0.0057171 | 93.98667 | 0.1709267 |
| 16 | 4 | 0.3629875 | 0.0248182 | 34.37500 | 0.2332957 |
| 20 | 5 | 0.3968397 | 0.0050750 | 34.09000 | 0.3485666 |
Para el mismo análisis por grupo se calculan las desviaciones estándar de cada uno.
| cluster5$cluster | Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello |
|---|---|---|---|---|
| 1 | 0.0419666 | 0.0035937 | 84.69856 | 0.0347746 |
| 2 | 0.0324794 | 0.0026422 | 178.76599 | 0.0362630 |
| 3 | 0.0349268 | 0.0044966 | 62.32880 | 0.0480733 |
| 4 | 0.0581805 | 0.0089293 | 30.75814 | 0.0705911 |
| 5 | 0.0555088 | 0.0058447 | 30.33077 | 0.0743526 |
Al igual que en la tabla de medias, se agrega el tamaño de grupos para interpretar mejor los resultados. Para algunos.
| tamano_grupo | cluster5$cluster | Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello |
|---|---|---|---|---|---|
| 111 | 1 | 0.0419666 | 0.0035937 | 84.69856 | 0.0347746 |
| 28 | 2 | 0.0324794 | 0.0026422 | 178.76599 | 0.0362630 |
| 90 | 3 | 0.0349268 | 0.0044966 | 62.32880 | 0.0480733 |
| 16 | 4 | 0.0581805 | 0.0089293 | 30.75814 | 0.0705911 |
| 20 | 5 | 0.0555088 | 0.0058447 | 30.33077 | 0.0743526 |
Se crea el dataframe cluster_final que contiene las medidas calculadas por barrio y su agrupación.
| BARRIO | Porcentaje_accidentes_finsemyfestivos | Porcentaje_fatalidad | promedio_anual_accidentes | Porcentaje_atropello | Cl5G |
|---|---|---|---|---|---|
| ALDEA PABLO VI | 0.4698795 | 0.0240964 | 16.6 | 0.3253012 | 4 |
| ALEJANDRIA | 0.1883768 | 0.0020040 | 99.8 | 0.0280561 | 1 |
| ALEJANDRO ECHAVARRIA | 0.3365042 | 0.0035672 | 168.2 | 0.1165279 | 3 |
| ALFONSO LOPEZ | 0.3220191 | 0.0060923 | 229.8 | 0.1157528 | 3 |
| ALTAMIRA | 0.2293814 | 0.0025773 | 155.2 | 0.0644330 | 1 |
Relato:
Los cinco grupos están agrupados de mayor a menor promedio anual de accidentalidad.
Los grupos con más accidentes promedio por año tienen la particularidad de tener una menor tasa de accidentalidad los fines de semana y festivo.
Tanto para los grupos de mayor y menor promedio anual de accidentes se presentan variaciones en la tasa de fatalidad de los accidentes:
Grupo 1: Promedio anual de accidentes alto, tasa de accidentes fines de semana y festivos alto, tasa de fatalidad baja y tasa de atropellos baja.
Grupo 2: Promedio anual de accidentes alto, tasa de accidentes fines de semana y festivos bajo, tasa de fatalidad baja y tasa de atropellos baja.
Grupo 3: Promedio anual de accidentes medio, tasa de accidentes fines de semana y festivos media, tasa de fatalidad media y tasa de atropellos media.
Grupo 4: Promedio anual de accidentes bajo, tasa de accidentes fines de semana y festivos alta, tasa de fatalidad alta y tasa de atropellos alta.
Grupo 5: Promedio anual de accidentes bajo, tasa de accidentes fines de semana y festivos alta, tasa de fatalidad baja y tasa de atropellos alta.
Luego se procede a utilizar la libreria leaflet para realizar el análisis espacial. Se une el dataframe que contiene el clustering con el dataframe de poligonos que permite identificar los barrios en el mapa.
Se eliminan las comunas del dataframe med, ya que no se consideraron para este análisis
Se cambian los nombres de algunos barrios en el dataframe med para que coincidan con los del dataframe cluster_final
Se reemplazan los caracteres con tilde por caracteres sin tilde
Se ponen los nombres de los barrios del dataframe med en mayuscula para poder hacer el join con el dataframe cluster_final
Se eliminan filas que contienen barrios repetidos y el centro administrativo la alpujarra, donde no se han presentado accidentes en los años analizados
tibble [265 x 9] (S3: sf/tbl_df/tbl/data.frame)
$ OBJECTID : int [1:265] 14 15 16 17 18 66 67 68 69 70 ...
$ SHAPEAREA : num [1:265] 100044 481662 413759 677361 358612 ...
$ SHAPELEN : num [1:265] 1525 2813 2892 4781 2612 ...
$ COMUNA : chr [1:265] "13" "07" "05" "15" ...
$ BARRIO : chr [1:265] "05" "01" "10" "11" ...
$ CODIGO : chr [1:265] "1305" "0701" "0510" "1511" ...
$ NOMBRE_BAR: chr [1:265] "METROPOLITANO" "UNIVERSIDAD NACIONAL" "TRICENTENARIO" "LA COLINA" ...
$ NOMBRE_COM: chr [1:265] "SAN JAVIER" "ROBLEDO" "CASTILLA" "GUAYABAL" ...
$ geometry :sfc_MULTIPOLYGON of length 265; first list element: List of 1
..$ :List of 1
.. ..$ : num [1:137, 1:2] -75.6 -75.6 -75.6 -75.6 -75.6 ...
..- attr(*, "class")= chr [1:3] "XY" "MULTIPOLYGON" "sfg"
- attr(*, "sf_column")= chr "geometry"
- attr(*, "agr")= Factor w/ 3 levels "constant","aggregate",..: NA NA NA NA NA NA NA NA
..- attr(*, "names")= chr [1:8] "OBJECTID" "SHAPEAREA" "SHAPELEN" "COMUNA" ...Se cambiar el nombre de la columna barrio al dataframe cluster_final
Se fusionan el dataframe med que contiene los polígonos con el dataframe cluster_final que tiene la agrupación de los barrios
tibble [265 x 14] (S3: sf/tbl_df/tbl/data.frame)
$ OBJECTID : int [1:265] 14 15 16 17 18 66 67 68 69 70 ...
$ SHAPEAREA : num [1:265] 100044 481662 413759 677361 358612 ...
$ SHAPELEN : num [1:265] 1525 2813 2892 4781 2612 ...
$ COMUNA : chr [1:265] "13" "07" "05" "15" ...
$ BARRIO : chr [1:265] "05" "01" "10" "11" ...
$ CODIGO : chr [1:265] "1305" "0701" "0510" "1511" ...
$ NOMBRE_BAR : chr [1:265] "METROPOLITANO" "UNIVERSIDAD NACIONAL" "TRICENTENARIO" "LA COLINA" ...
$ NOMBRE_COM : chr [1:265] "SAN JAVIER" "ROBLEDO" "CASTILLA" "GUAYABAL" ...
$ geometry :sfc_MULTIPOLYGON of length 265; first list element: List of 1
..$ :List of 1
.. ..$ : num [1:137, 1:2] -75.6 -75.6 -75.6 -75.6 -75.6 ...
..- attr(*, "class")= chr [1:3] "XY" "MULTIPOLYGON" "sfg"
$ Porcentaje_accidentes_finsemyfestivos: num [1:265] 0.359 0.25 0.284 0.322 0.264 ...
$ Porcentaje_fatalidad : num [1:265] 0 0.0065 0.00514 0.00784 0.00217 ...
$ promedio_anual_accidentes : num [1:265] 7.8 369 194.4 102 276.4 ...
$ Porcentaje_atropello : num [1:265] 0.1538 0.0488 0.0494 0.1 0.0904 ...
$ Cl5G : int [1:265] 3 2 1 3 1 3 1 3 3 1 ...
- attr(*, "sf_column")= chr "geometry"
- attr(*, "agr")= Factor w/ 3 levels "constant","aggregate",..: NA NA NA NA NA NA NA NA NA NA ...
..- attr(*, "names")= chr [1:13] "OBJECTID" "SHAPEAREA" "SHAPELEN" "COMUNA" ...A continuación se importa la librereia leaflet para hacer el gráfico de los barrios
Se carga el mapa base
Se definen los colores de la paleta de colores (Set3), la cual establece un color de grupo distinto para cada uno de los que se encuentra en la columna Cl5G del datraframe.Luego se añaden los poligonos de los barrios, los colores y la leyenda.
Finalmente se definen las propiedades de los polígonos como opacidad, tamaño de las lineas, etc y las propiedades de la leyenda.
Características espaciales:Se encuentra que los barrios con mayor promedio de accidentes por año están ubicados cerca al rio Medellín que es precisamente donde se encuentran las principales vías de la ciudad. El siguiente grupo con mayor promedio de accidentes se ubica un poco más hacia las laderas, sin embargo, este patrón es más visible hacia el sur que hacia el norte. Esto puede deberse a que los barrios del sur por lo general tienen mejores vías y mayor flujo vehicular. Los barrios con menor accidentalidad se encuentran totalmente concentrados en las partes altas de las laderas, con una tendencia más clara hacia el norte. La mayoría de los barrios con menor accidentalidad se encuentran en dos sectores específicos de la ciudad: la comuna 13 y la comuna 1 (Popular).
Retos que el agrupamiento puede resolver y posibles planes de acción:
La estrategia principal se basa en destinar más recursos a los barrios donde más accidentes se presenten. Sin embargo, los demás indicadores servirán para optimizar el uso de los recursos en cada barrio independiente de la cantidad de accidentes que se presenten.
Grupo 1:
Por tener un promedio alto de accidentes, estos barrios deben tener una mayor disponibilidad de agentes de tránsito, ambulancias, cámaras y grúas para atender los incidentes. Los porcentajes de atropellos y fatalidad no son tan altos como los otros grupos, lo que indica que principalmente lo que se presentan son choques. Pueden realizarse mejoras en la señalización, evaluar la necesidad de nuevos policías acostados y todo aquello que ayude a que se presenten menos incidentes menores en horas pico.
Grupo 2: Presenta el promedio mayor de accidentes por año en comparación con el resto del grupo, por lo que son a los que más se les debe prestar atención. Al estar ubicados muy cerca a la avenida regional, las vegas y demás vías céntricas de la ciudad, estos barrios deben de contar con presencia permanente del tránsito en horas pico para mejorar la movilidad. La accidentalidad no es tan alta los fines de semana comparada con los porcentajes de los otros grupos, por lo tanto, la mayor atención debe presentarse en las grandes concentraciones de vehículos que ocurren en los días laborales.
Grupo 3: Este es un grupo con indicadores intermedios, con barrios ubicados principalmente en las laderas occidental, oriental y norte de la ciudad. Aunque no requieren una presencia permanente del tránsito, pueden mejorarse las condiciones de movilidad utilizando herramientas tecnológicas como cámaras.
Grupo 4: Si bien este grupo tiene un promedio anual de accidentes bajo, presenta un porcentaje alto de accidentes los fines de semana en relación con los demás grupos. Sumado a esto tiene tasas de fatalidad y atropello más altas. La estrategia para estos barrios está más enfocada en campañas educativas que permitan disminuir los accidentes fatales, como son aquellas relacionadas con controlar los excesos de velocidad. Al mismo tiempo, pueden realizarse retenes de tránsito con más frecuencia los fines de semana para hacer diferentes tipos de controles, entre ellos el de la alcoholemia, ya que conducir bajo los efectos del alcohol es una de las razones por las que ocurren más accidentes fatales.
Grupo 5: Pueden realizarse acciones muy similares a la del grupo 4, con la diferencia de que este grupo presenta una tasa menor de accidentalidad, razón por la cual no es necesario hacer retenes con tanta frecuencia.