Classificació atès al motiu de desplaçament dels conductors ferits en accidents de trànsit a Barcelona (2023)
Joan Manel Ramírez Jávega
27 Agost del 2024
1 Introducció i presentació de les dades
Donant continuïtat al treball ja realitzat en el marc de l’assignatura Mineria de dades, cursada durant el segon semestre del curs 2023-2024 del Màster de Ciència de Dades (Universitat Oberta de Catalunya) i podent-se consultar tant la primera part com la segona part en el repositori RPubs, ens proposem millorar-ne la modelització ja realitzada i escollir el model més addient per realitzar la imputació en la variable depenent “Es_mon_treball”. Finalment, també prepararem les dades agregant-hi les seves respectives geocoordenades perque així puguin ser utilitzades per programaris d’anàlisi de dades geoespacials. Atès a que part de les tasques ja s’han presentat i comentat en les dues parts anteriors, evitarem repetir aquelles seccions amb codi que ja s’han detallat en els dos articles anteriorment referits.
En aquest article realitzarem les tasques d’estandarització, neteja i codificació pel conjunt de dades corresponent a aquells conductors dels qui constava en l’expedient obert per la Guàrdia Urbana de Barcelona que el seu motiu de desplaçament com Es desconeix. Atès a que la naturalesa repetitiva d’aquest fragment amb el ja detallat en la primera part, s’evitarà repetir-lo de nou.
maskDistDescon <- dades$Nom_districte == "Desconegut"
dadesDistDescon <- dades[maskDistDescon, ]
dades_ <- dades[!maskDistDescon, ]
maskDistDescon <- dadesDesconegut$Nom_districte == "Desconegut"
dadesDistDescon <- dadesDesconegut[maskDistDescon, ]
dadesDesconegut_ <- dadesDesconegut[!maskDistDescon, ]
nMotiuConegut <- nrow(dades_)
nMotiuDesconegut <- nrow(dadesDesconegut_)
nConductors <- nMotiuConegut + nMotiuDesconegut
propMotiuConegut <- round((nMotiuConegut /nConductors) * 100, 3)
propMotiuDesconegut <- round((nMotiuDesconegut /nConductors) * 100, 3)
print(glue("Es prenen en consideració {nConductors} conductors en total, dels qui {nMotiuConegut} ({propMotiuConegut} %) sí
es coneix el motiu del seu desplaçament i d'altres {nMotiuDesconegut} ({propMotiuDesconegut} %) se'n
desconeix el motiu del seu desplaçament."))## Es prenen en consideració 5019 conductors en total, dels qui 2823 (56.246 %) sí
## es coneix el motiu del seu desplaçament i d'altres 2196 (43.754 %) se'n
## desconeix el motiu del seu desplaçament.En canvi sí resulta rellevant observar que, a diferència de la tasca de codificació realitzada en la segona part de la variable “Es_ocupacional”, que identifica si el moment que va tenir lloc l’accident era durant una franja horària amb un elevada proporció de desplaçament relacionats amb el món del treball o no en dies laborables, ara actualitzarem amb els resultats de l’Enquesta de mobilitat en dia feiner 2023 corresponents a l’any 2023 que fan referència a les franjes horàries i el motiu del desplaçament (INSTITUT METRÒPOLI, 2023: 49), tot i que a partir d’aquesta evidència en contrast a les conclusions de l’estudi esmentat, trobem rellevant el nombre de desplaçaments relacionats amb el món del treball durant la franja des de les 5 fins les 16 hores - ambdues franges incloses-:
dades_["Es_ocupacional"] <- "No"
horaOcupacional <- c(5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16)
maskHorari <- dades_$Hora_dia %in% horaOcupacional
maskLaborable <- dades_$Es_laborable == "Si"
dades_$Es_ocupacional[maskHorari & maskLaborable] <- "Si"
aggEs_ocupacional <- aggregate(dades_$Numero_expedient,
by=dades_["Es_ocupacional"], FUN=length)
ggplot(aggEs_ocupacional, aes(x="", y=x, fill=Es_ocupacional)) +
geom_bar(width = 1, stat = "identity", color = "black") +
ggtitle(glue("Fig. 1.1. Diagrama de barres apilades del recompte dels
conductors en funció de si l'accident va tenir lloc durant
horari ocupacional o no i també es coneix el motiu del
desplaçament.")) +
xlab(glue("Es_ocupacional")) +
ylab("Nombre")
nAccidents <- length(unique(dades_$Numero_expedient))
nNo_Es_ocupacional <- aggEs_ocupacional[1, 2]
propNo_Es_ocupacional <- round((nNo_Es_ocupacional / nMotiuConegut) * 100, 3)
nSi_Es_ocupacional <- aggEs_ocupacional[2, 2]
propSi_Es_ocupacional <- round((nSi_Es_ocupacional / nMotiuConegut) * 100, 3)
print(glue("Del conjunt de conductors dels qui sí es coneix el motiu del seu desplaçament,
{nSi_Es_ocupacional} ({propSi_Es_ocupacional} %) es desplaçaven durant horari ocupacional i altres {nNo_Es_ocupacional}
({propNo_Es_ocupacional} %), comptant-se {nAccidents} accidents en total."))## Del conjunt de conductors dels qui sí es coneix el motiu del seu desplaçament,
## 1272 (45.058 %) es desplaçaven durant horari ocupacional i altres 1551
## (54.942 %), comptant-se 2561 accidents en total.dadesDesconegut_["Es_ocupacional"] <- "No"
horaOcupacional <- c(5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16)
maskHorari <- dadesDesconegut_$Hora_dia %in% horaOcupacional
maskLaborable <- dadesDesconegut_$Es_laborable == "Si"
dadesDesconegut_$Es_ocupacional[maskHorari & maskLaborable] <- "Si"
aggEs_ocupacional <- aggregate(dadesDesconegut_$Numero_expedient,
by=dadesDesconegut_["Es_ocupacional"], FUN=length)
ggplot(aggEs_ocupacional, aes(x="", y=x, fill=Es_ocupacional)) +
geom_bar(width = 1, stat = "identity", color = "black") +
ggtitle(glue("Fig. 1.2. Diagrama de barres apilades del recompte dels
conductors en funció de si l'accident va tenir lloc durant
horari ocupacional o no i no es coneix el motiu del
desplaçament.")) +
xlab(glue("Es_ocupacional")) +
ylab("Nombre")
nAccidents <- length(unique(dadesDesconegut_$Numero_expedient))
nNo_Es_ocupacional <- aggEs_ocupacional[1, 2]
propNo_Es_ocupacional <- round((nNo_Es_ocupacional / nMotiuDesconegut) * 100, 3)
nSi_Es_ocupacional <- aggEs_ocupacional[2, 2]
propSi_Es_ocupacional <- round((nSi_Es_ocupacional / nMotiuDesconegut) * 100, 3)
print(glue("Del conjunt de conductors dels qui sí es coneix el motiu del seu desplaçament,
{nSi_Es_ocupacional} ({propSi_Es_ocupacional} %) es desplaçaven durant horari ocupacional i altres {nNo_Es_ocupacional}
({propNo_Es_ocupacional} %) no, comptant-se {nAccidents} accidents en total."))## Del conjunt de conductors dels qui sí es coneix el motiu del seu desplaçament,
## 1150 (52.368 %) es desplaçaven durant horari ocupacional i altres 1046
## (47.632 %) no, comptant-se 2029 accidents en total.Atès al resultat observable en les figures 1.1 i 1.2, es constata que la proporció de conductors ferits en accidents que van tenir lloc durant una franja horària ocupacional o no s’inverteix en el cas que es coneix o no també el motiu del seu desplaçament.
També comprovem si té lloc algun canvi rellevant en el tamany de les classes en la variable de la descripció de victimització com mostra, atès que en les Figures 15.1 i 15.2 de la primera part ja en fem la comprovació que la distribució territorial no es veia fectada per la neteja de les dades:
dades_["Victimitzacio_est_"] <- ""
dades_$Victimitzacio_est_ <- sapply(strsplit(dades_$Descripcio_victimitzacio,
": "), "[", 2)
maskMorts <- grepl("Mort", dades_$Descripcio_victimitzacio)
dades_$Victimitzacio_est_[maskMorts] <- "Mort"
dades_$Victimitzacio_est_ <- str_to_sentence(dades_$Victimitzacio_est_)
aggVictmitzacio <- aggregate(dades_$Numero_expedient,
by=dades_["Victimitzacio_est_"],
FUN=length)
ggplot(aggVictmitzacio, aes(x="", y=x, fill=Victimitzacio_est_)) +
geom_bar(width = 1, stat = "identity", color = "black") +
ggtitle(glue("Fig. 2.1. Diagrama de barres apilades del recompte dels
conductors en funció del tipus de victimització i que sí es
coneixia el motiu del desplaçament.")) +
xlab(glue("Victimitzacio_est_")) +
ylab("Nombre")
dadesDesconegut_["Victimitzacio_est_"] <- ""
dadesDesconegut_$Victimitzacio_est_ <- sapply(strsplit(dadesDesconegut_$Descripcio_victimitzacio,
": "), "[", 2)
maskMorts <- grepl("Mort", dadesDesconegut_$Descripcio_victimitzacio)
dadesDesconegut_$Victimitzacio_est_[maskMorts] <- "Mort"
dadesDesconegut_$Victimitzacio_est_ <- str_to_sentence(dadesDesconegut_$Victimitzacio_est_)
aggVictmitzacio <- aggregate(dadesDesconegut_$Numero_expedient,
by=dadesDesconegut_["Victimitzacio_est_"],
FUN=length)
ggplot(aggVictmitzacio, aes(x="", y=x, fill=Victimitzacio_est_)) +
geom_bar(width = 1, stat = "identity", color = "black") +
ggtitle(glue("Fig. 2.2. Diagrama de barres apilades del recompte dels
conductors en funció del tipus de victimització i que sí es
coneixia el motiu del desplaçament.")) +
xlab(glue("Victimitzacio_est_")) +
ylab("Nombre")
Tot apreciant-se al comparar les Figures 2.1 i 2.2 que la proporció dels grups és pràcticament idèntica tant el cas dels conductors que sí es coneixia el motiu del desplaçament com en el cas dels que es desconeixia.
També procedim a codificar una nova variable a partir del mes en que es produí l’accident, tot classificant-los pels trimestres corresponents a les estacions de l’any: Hivern, Primavera, Estiu i Tardor.
dades_['Estacio_any'] <- "Estiu"
dades_$Estacio_any[dades_$Mes_any <= 3] <- "Hivern"
dades_$Estacio_any[dades_$Mes_any > 3 & dades_$Mes_any <= 6] <- "Primavera"
dades_$Estacio_any[dades_$Mes_any > 9] <- 'Tardor'
aggEstacions <- aggregate(dades_$Numero_expedient,
by=dades_["Estacio_any"],
FUN=length)
ggplot(aggEstacions, aes(x="", y=x, fill=Estacio_any)) +
geom_bar(width = 1, stat = "identity", color = "black") +
ggtitle(glue("Fig. 3.1. Diagrama de barres apilades del recompte dels
conductors en funció de l'estació de l'any en la que es produí l'accident
i que sí es coneixia el motiu del desplaçament.")) +
xlab(glue("Estacio_any")) +
ylab("Nombre")
dades_$Tipus_vehicle_estandaritzat <- gsub("d'", "",
dades_$Tipus_vehicle_estandaritzat)
dades_$Victimitzacio_est_ <- gsub("d'", "",
dades_$Victimitzacio_est_)
dades_$Descripcio_causa_mediata <- gsub("d'", "",
dades_$Descripcio_causa_mediata)
dades_$Nom_mes <- gsub("ç", "s", dades_$Nom_mes)
dades_$Tipus_vehicle_estandaritzat <-
stri_trans_general(str = dades_$Tipus_vehicle_estandaritzat,
id = "Latin-ASCII")
dades_$Victimitzacio_est_ <-
stri_trans_general(str = dades_$Victimitzacio_est_,
id = "Latin-ASCII")
dades_$Descripcio_causa_mediata <-
stri_trans_general(str = dades_$Descripcio_causa_mediata,
id = "Latin-ASCII")
dades_$Nom_districte <-
stri_trans_general(str = dades_$Nom_districte,
id = "Latin-ASCII")
varsSup <- c("Descripcio_sexe", "Tipus_vehicle_estandaritzat",
"Victimitzacio_est_", "Nom_mes", "Nom_districte",
"Descripcio_causa_mediata", "Estacio_any",
"Victimes_CODIF", "Numero_morts",
"Numero_lesionats_lleus", "Numero_lesionats_greus",
"Numero_victimes", "Numero_vehicles_implicats",
"Vehicles_CODIF", "VM2R_CODIF", "VM4R_CODIF",
"V_no_permis_CODIF", "VUP_CODIF", "Interve_conductor_novell",
"Edat_CODIF", "Edat", "Es_laborable", "Es_ocupacional",
"Lleus_CODIF", "Greus_CODIF", "Morts_CODIF", "Es_atropellament",
"Es_mon_treball")
dadesSup <- dades_[varsSup]dadesDesconegut_['Estacio_any'] <- "Estiu"
dadesDesconegut_$Estacio_any[dadesDesconegut_$Mes_any <= 3] <- "Hivern"
dadesDesconegut_$Estacio_any[dadesDesconegut_$Mes_any > 3 & dadesDesconegut_$Mes_any <= 6] <- "Primavera"
dadesDesconegut_$Estacio_any[dadesDesconegut_$Mes_any > 9] <- 'Tardor'
aggEstacions <- aggregate(dadesDesconegut_$Numero_expedient,
by=dadesDesconegut_["Estacio_any"],
FUN=length)
ggplot(aggEstacions, aes(x="", y=x, fill=Estacio_any)) +
geom_bar(width = 1, stat = "identity", color = "black") +
ggtitle(glue("Fig. 3.2. Diagrama de barres apilades del recompte dels
conductors en funció de l'estació de l'any en la que es produí l'accident
i que no es coneixia el motiu del desplaçament.")) +
xlab(glue("Estacio_any")) +
ylab("Nombre")
dadesDesconegut_$Tipus_vehicle_estandaritzat <- gsub("d'", "",
dadesDesconegut_$Tipus_vehicle_estandaritzat)
dadesDesconegut_$Victimitzacio_est_ <- gsub("d'", "",
dadesDesconegut_$Victimitzacio_est_)
dadesDesconegut_$Descripcio_causa_mediata <- gsub("d'", "",
dadesDesconegut_$Descripcio_causa_mediata)
dadesDesconegut_$Nom_mes <- gsub("ç", "s", dadesDesconegut_$Nom_mes)
dadesDesconegut_$Tipus_vehicle_estandaritzat <-
stri_trans_general(str = dadesDesconegut_$Tipus_vehicle_estandaritzat,
id = "Latin-ASCII")
dadesDesconegut_$Victimitzacio_est_ <-
stri_trans_general(str = dadesDesconegut_$Victimitzacio_est_,
id = "Latin-ASCII")
dadesDesconegut_$Descripcio_causa_mediata <-
stri_trans_general(str = dadesDesconegut_$Descripcio_causa_mediata,
id = "Latin-ASCII")
dadesDesconegut_$Nom_districte <-
stri_trans_general(str = dadesDesconegut_$Nom_districte,
id = "Latin-ASCII")
varsSup <- c("Descripcio_sexe", "Tipus_vehicle_estandaritzat",
"Victimitzacio_est_", "Nom_mes", "Nom_districte",
"Descripcio_causa_mediata", "Estacio_any",
"Victimes_CODIF", "Numero_morts",
"Numero_lesionats_lleus", "Numero_lesionats_greus",
"Numero_victimes", "Numero_vehicles_implicats",
"Vehicles_CODIF", "VM2R_CODIF", "VM4R_CODIF",
"V_no_permis_CODIF", "VUP_CODIF", "Interve_conductor_novell",
"Edat_CODIF", "Edat", "Es_laborable", "Es_ocupacional",
"Lleus_CODIF", "Greus_CODIF", "Morts_CODIF", "Es_atropellament",
"Es_mon_treball")
dadesDesconegutSup <- dadesDesconegut_[varsSup]Una vegada més, s’observa que entre ambdós subconjunts de dades no hi ha una diferència substancial en la distribució dels grups i sent en ambdós casos durant la primavera quan es comptabilitzen més conductors ferits, tot i que en el cas dels conductors dels que no es coneixia el motiu del seu desplaçament s’observa un lleu increment de la proporció per l’estació hivernal.
2 Modelització
En les conclusions de la segona part d’aquesta sèrie vam constatar tot un seguit de mancances en les modelitzacions i la seva acpacitat d’encert, sent possiblement un dels motius el fet que hi havia un nombre molt elevat de variables en joc. Per aquest motiu en aquest apartat cercarem, en primer lloc, bastir un model de classificació amb l’algoritme de la Regressió Logística que ens permeti identificar quines de les variables del total de 27 variables utilitzables són estadísticament rellevants per identificar la variable dependent “Es_mon_treball”. Posteriorment, a partir del resultat obtingut d’aquest primer algoritme, bastirem diversos models amb diversos tipus algoritmes com l’arbre de decisió sense poda o amb poda i amb diverses parametritzacions a més dels models del bosc aleatori i i el XGBoost.
2.1 Regressió logística
De de la descomposició dels elements principlas realitzada en l’apartat 5 de la primera part ja coneixem que les variables numèriques, amb l’excepció de “Numero_victimes” i “Numero_lesionats_greus”, presenten una escassa correlació entre sí. En tot cas, cercarem introduir-les en el model tant en la seva versió contínua com la codificada. I a continuació procedim a realitzar l’anàlisi de les correlacions entre variables mitjançant l’algorisme de la Regressió logística:
dadesSup["EMT_recodificat"] <- ifelse(dadesSup$Es_mon_treball=="Si", 1, 0)
cols_ <- colnames(dadesSup)
df_dtype_ <- as.data.frame(sapply(dadesSup, class))
for (col in cols_){
data_type <- paste(df_dtype_[col, 1])
n_domini <- length(unique(dadesSup[, col]))
if (data_type == "character"){
dadesSup[, col] <- as.factor(dadesSup[, col])
}
else if (n_domini <= 2) {
dadesSup[, col] <- as.factor(dadesSup[, col])
}
}
dadesSup$Numero_morts <- as.integer(dadesSup$Numero_morts)cols_ <- colnames(dadesDesconegutSup)
df_dtype_ <- as.data.frame(sapply(dadesDesconegutSup, class))
for (col in cols_){
data_type <- paste(df_dtype_[col, 1])
n_domini <- length(unique(dadesDesconegutSup[, col]))
if (data_type == "character"){
dadesDesconegutSup[, col] <- as.factor(dadesDesconegutSup[, col])
}
else if (n_domini <= 2) {
dadesDesconegutSup[, col] <- as.factor(dadesDesconegutSup[, col])
}
}
dadesDesconegutSup$Numero_morts <- as.integer(dadesDesconegutSup$Numero_morts)model <- glm(EMT_recodificat ~ Descripcio_sexe +
Edat_CODIF + Edat + Tipus_vehicle_estandaritzat +
Victimitzacio_est_ +
Descripcio_causa_mediata + Estacio_any +
Numero_morts + Morts_CODIF +
Numero_lesionats_lleus + Lleus_CODIF +
Numero_lesionats_greus + Greus_CODIF +
Numero_vehicles_implicats + Vehicles_CODIF +
Numero_victimes + Victimes_CODIF +
Nom_districte + VM2R_CODIF + VM4R_CODIF + V_no_permis_CODIF +
VUP_CODIF + Interve_conductor_novell + Es_atropellament +
Es_laborable + Es_ocupacional,
data=dadesSup, family=binomial(link=logit), na.action = NULL)
summary(model)##
## Call:
## glm(formula = EMT_recodificat ~ Descripcio_sexe + Edat_CODIF +
## Edat + Tipus_vehicle_estandaritzat + Victimitzacio_est_ +
## Descripcio_causa_mediata + Estacio_any + Numero_morts + Morts_CODIF +
## Numero_lesionats_lleus + Lleus_CODIF + Numero_lesionats_greus +
## Greus_CODIF + Numero_vehicles_implicats + Vehicles_CODIF +
## Numero_victimes + Victimes_CODIF + Nom_districte + VM2R_CODIF +
## VM4R_CODIF + V_no_permis_CODIF + VUP_CODIF + Interve_conductor_novell +
## Es_atropellament + Es_laborable + Es_ocupacional, family = binomial(link = logit),
## data = dadesSup, na.action = NULL)
##
## Coefficients: (2 not defined because of singularities)
## Estimate
## (Intercept) -3.141604
## Descripcio_sexeHome 0.057154
## Edat_CODIFConductors entre 30 i 37 anys edat -0.099422
## Edat_CODIFConductors fins 29 anys edat -0.595047
## Edat_CODIFConductors majors de 49 anys edat 0.201110
## Edat -0.018628
## Tipus_vehicle_estandaritzatVehicles motoritzats de quatre rodes 0.167033
## Tipus_vehicle_estandaritzatVehicles sense permis de conduccio -0.190535
## Tipus_vehicle_estandaritzatVehicles Us Professional 4.530266
## Victimitzacio_est_Hospitalitzacio fins a 24h 0.159140
## Victimitzacio_est_Hospitalitzacio superior a 24h 0.611889
## Victimitzacio_est_Mort -0.896369
## Victimitzacio_est_Rebutja assistencia sanitaria 0.284140
## Descripcio_causa_mediataCanvi de carril sense precaucio -0.227889
## Descripcio_causa_mediataDesobeir altres senyals -0.402497
## Descripcio_causa_mediataDesobeir semafor -0.126338
## Descripcio_causa_mediataEnvair calcada contraria -1.181084
## Descripcio_causa_mediataFallada mecanica o avaria -10.738271
## Descripcio_causa_mediataGir indegut o sense precaucio -0.339264
## Descripcio_causa_mediataManca atencio a la conduccio -0.277985
## Descripcio_causa_mediataManca precaucio efectuar marxa enrera -0.285622
## Descripcio_causa_mediataManca precaucio incorporacio circulacio -0.220823
## Descripcio_causa_mediataNo cedir la dreta -0.181772
## Descripcio_causa_mediataNo respectar distancies -0.497367
## Descripcio_causa_mediataNo respectat pas de vianants -0.022101
## Estacio_anyHivern 0.007894
## Estacio_anyPrimavera 0.233228
## Estacio_anyTardor 0.038313
## Numero_morts 1.657122
## Morts_CODIFSi NA
## Numero_lesionats_lleus -0.020486
## Lleus_CODIFSi 0.552134
## Numero_lesionats_greus 1.202418
## Greus_CODIFSi -1.876713
## Numero_vehicles_implicats -0.033703
## Vehicles_CODIFSi 0.799640
## Numero_victimes NA
## Victimes_CODIFSi -0.383877
## Nom_districteEixample 0.070432
## Nom_districteGracia 0.020947
## Nom_districteHorta-Guinardo -0.344914
## Nom_districteLes Corts 0.450337
## Nom_districteNou Barris -0.086124
## Nom_districteSant Andreu -0.328615
## Nom_districteSant Marti 0.113460
## Nom_districteSants-Montjuic 0.345712
## Nom_districteSarria-Sant Gervasi -0.011230
## VM2R_CODIFSi 0.103620
## VM4R_CODIFSi 0.128787
## V_no_permis_CODIFSi 0.058489
## VUP_CODIFSi 0.347619
## Interve_conductor_novellSi -0.031369
## Es_atropellamentSi 0.471575
## Es_laborableSi 0.621579
## Es_ocupacionalSi 0.872729
## Std. Error
## (Intercept) 1.549508
## Descripcio_sexeHome 0.092090
## Edat_CODIFConductors entre 30 i 37 anys edat 0.152061
## Edat_CODIFConductors fins 29 anys edat 0.210195
## Edat_CODIFConductors majors de 49 anys edat 0.177586
## Edat 0.009312
## Tipus_vehicle_estandaritzatVehicles motoritzats de quatre rodes 0.124465
## Tipus_vehicle_estandaritzatVehicles sense permis de conduccio 0.135211
## Tipus_vehicle_estandaritzatVehicles Us Professional 0.736989
## Victimitzacio_est_Hospitalitzacio fins a 24h 0.101444
## Victimitzacio_est_Hospitalitzacio superior a 24h 0.743511
## Victimitzacio_est_Mort 1.733621
## Victimitzacio_est_Rebutja assistencia sanitaria 0.169696
## Descripcio_causa_mediataCanvi de carril sense precaucio 0.221669
## Descripcio_causa_mediataDesobeir altres senyals 0.246703
## Descripcio_causa_mediataDesobeir semafor 0.227335
## Descripcio_causa_mediataEnvair calcada contraria 0.978785
## Descripcio_causa_mediataFallada mecanica o avaria 196.968101
## Descripcio_causa_mediataGir indegut o sense precaucio 0.220003
## Descripcio_causa_mediataManca atencio a la conduccio 0.206952
## Descripcio_causa_mediataManca precaucio efectuar marxa enrera 0.436302
## Descripcio_causa_mediataManca precaucio incorporacio circulacio 0.274093
## Descripcio_causa_mediataNo cedir la dreta 0.371908
## Descripcio_causa_mediataNo respectar distancies 0.213964
## Descripcio_causa_mediataNo respectat pas de vianants 0.607251
## Estacio_anyHivern 0.121827
## Estacio_anyPrimavera 0.115624
## Estacio_anyTardor 0.119853
## Numero_morts 1.320503
## Morts_CODIFSi NA
## Numero_lesionats_lleus 0.073968
## Lleus_CODIFSi 0.601838
## Numero_lesionats_greus 0.578543
## Greus_CODIFSi 0.780793
## Numero_vehicles_implicats 0.090138
## Vehicles_CODIFSi 0.283757
## Numero_victimes NA
## Victimes_CODIFSi 0.140993
## Nom_districteEixample 0.220617
## Nom_districteGracia 0.284376
## Nom_districteHorta-Guinardo 0.262054
## Nom_districteLes Corts 0.257007
## Nom_districteNou Barris 0.278184
## Nom_districteSant Andreu 0.261826
## Nom_districteSant Marti 0.234503
## Nom_districteSants-Montjuic 0.242471
## Nom_districteSarria-Sant Gervasi 0.241102
## VM2R_CODIFSi 0.128622
## VM4R_CODIFSi 0.108049
## V_no_permis_CODIFSi 0.102270
## VUP_CODIFSi 0.227924
## Interve_conductor_novellSi 0.096416
## Es_atropellamentSi 0.415103
## Es_laborableSi 0.118900
## Es_ocupacionalSi 0.095891
## z value
## (Intercept) -2.027
## Descripcio_sexeHome 0.621
## Edat_CODIFConductors entre 30 i 37 anys edat -0.654
## Edat_CODIFConductors fins 29 anys edat -2.831
## Edat_CODIFConductors majors de 49 anys edat 1.132
## Edat -2.000
## Tipus_vehicle_estandaritzatVehicles motoritzats de quatre rodes 1.342
## Tipus_vehicle_estandaritzatVehicles sense permis de conduccio -1.409
## Tipus_vehicle_estandaritzatVehicles Us Professional 6.147
## Victimitzacio_est_Hospitalitzacio fins a 24h 1.569
## Victimitzacio_est_Hospitalitzacio superior a 24h 0.823
## Victimitzacio_est_Mort -0.517
## Victimitzacio_est_Rebutja assistencia sanitaria 1.674
## Descripcio_causa_mediataCanvi de carril sense precaucio -1.028
## Descripcio_causa_mediataDesobeir altres senyals -1.632
## Descripcio_causa_mediataDesobeir semafor -0.556
## Descripcio_causa_mediataEnvair calcada contraria -1.207
## Descripcio_causa_mediataFallada mecanica o avaria -0.055
## Descripcio_causa_mediataGir indegut o sense precaucio -1.542
## Descripcio_causa_mediataManca atencio a la conduccio -1.343
## Descripcio_causa_mediataManca precaucio efectuar marxa enrera -0.655
## Descripcio_causa_mediataManca precaucio incorporacio circulacio -0.806
## Descripcio_causa_mediataNo cedir la dreta -0.489
## Descripcio_causa_mediataNo respectar distancies -2.325
## Descripcio_causa_mediataNo respectat pas de vianants -0.036
## Estacio_anyHivern 0.065
## Estacio_anyPrimavera 2.017
## Estacio_anyTardor 0.320
## Numero_morts 1.255
## Morts_CODIFSi NA
## Numero_lesionats_lleus -0.277
## Lleus_CODIFSi 0.917
## Numero_lesionats_greus 2.078
## Greus_CODIFSi -2.404
## Numero_vehicles_implicats -0.374
## Vehicles_CODIFSi 2.818
## Numero_victimes NA
## Victimes_CODIFSi -2.723
## Nom_districteEixample 0.319
## Nom_districteGracia 0.074
## Nom_districteHorta-Guinardo -1.316
## Nom_districteLes Corts 1.752
## Nom_districteNou Barris -0.310
## Nom_districteSant Andreu -1.255
## Nom_districteSant Marti 0.484
## Nom_districteSants-Montjuic 1.426
## Nom_districteSarria-Sant Gervasi -0.047
## VM2R_CODIFSi 0.806
## VM4R_CODIFSi 1.192
## V_no_permis_CODIFSi 0.572
## VUP_CODIFSi 1.525
## Interve_conductor_novellSi -0.325
## Es_atropellamentSi 1.136
## Es_laborableSi 5.228
## Es_ocupacionalSi 9.101
## Pr(>|z|)
## (Intercept) 0.04261 *
## Descripcio_sexeHome 0.53484
## Edat_CODIFConductors entre 30 i 37 anys edat 0.51322
## Edat_CODIFConductors fins 29 anys edat 0.00464 **
## Edat_CODIFConductors majors de 49 anys edat 0.25744
## Edat 0.04546 *
## Tipus_vehicle_estandaritzatVehicles motoritzats de quatre rodes 0.17959
## Tipus_vehicle_estandaritzatVehicles sense permis de conduccio 0.15878
## Tipus_vehicle_estandaritzatVehicles Us Professional 7.90e-10 ***
## Victimitzacio_est_Hospitalitzacio fins a 24h 0.11671
## Victimitzacio_est_Hospitalitzacio superior a 24h 0.41052
## Victimitzacio_est_Mort 0.60512
## Victimitzacio_est_Rebutja assistencia sanitaria 0.09405 .
## Descripcio_causa_mediataCanvi de carril sense precaucio 0.30392
## Descripcio_causa_mediataDesobeir altres senyals 0.10278
## Descripcio_causa_mediataDesobeir semafor 0.57839
## Descripcio_causa_mediataEnvair calcada contraria 0.22755
## Descripcio_causa_mediataFallada mecanica o avaria 0.95652
## Descripcio_causa_mediataGir indegut o sense precaucio 0.12305
## Descripcio_causa_mediataManca atencio a la conduccio 0.17919
## Descripcio_causa_mediataManca precaucio efectuar marxa enrera 0.51270
## Descripcio_causa_mediataManca precaucio incorporacio circulacio 0.42044
## Descripcio_causa_mediataNo cedir la dreta 0.62501
## Descripcio_causa_mediataNo respectar distancies 0.02010 *
## Descripcio_causa_mediataNo respectat pas de vianants 0.97097
## Estacio_anyHivern 0.94834
## Estacio_anyPrimavera 0.04368 *
## Estacio_anyTardor 0.74922
## Numero_morts 0.20951
## Morts_CODIFSi NA
## Numero_lesionats_lleus 0.78182
## Lleus_CODIFSi 0.35893
## Numero_lesionats_greus 0.03768 *
## Greus_CODIFSi 0.01623 *
## Numero_vehicles_implicats 0.70847
## Vehicles_CODIFSi 0.00483 **
## Numero_victimes NA
## Victimes_CODIFSi 0.00648 **
## Nom_districteEixample 0.74953
## Nom_districteGracia 0.94128
## Nom_districteHorta-Guinardo 0.18811
## Nom_districteLes Corts 0.07973 .
## Nom_districteNou Barris 0.75687
## Nom_districteSant Andreu 0.20945
## Nom_districteSant Marti 0.62850
## Nom_districteSants-Montjuic 0.15393
## Nom_districteSarria-Sant Gervasi 0.96285
## VM2R_CODIFSi 0.42046
## VM4R_CODIFSi 0.23329
## V_no_permis_CODIFSi 0.56739
## VUP_CODIFSi 0.12722
## Interve_conductor_novellSi 0.74491
## Es_atropellamentSi 0.25594
## Es_laborableSi 1.72e-07 ***
## Es_ocupacionalSi < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 3900.6 on 2822 degrees of freedom
## Residual deviance: 3432.0 on 2770 degrees of freedom
## AIC: 3538
##
## Number of Fisher Scoring iterations: 10Si establim que el llindar per acceptar la variable és 0.05 > valor-p (Pr(>|z|)) per considerar-les que tenen un valor explicatiu estadísticament rellevant, atès als resultats hauríem d’acceptar només 11 de les variables tingudes en consideració en l’algorisme de Regressió logística, comptant-se entre aquestes 2 variables numèriques contínues i 9 variables categòriques que, tot seguit, detallem en aquesta taula:
| Nom de la variable | Tipus de dada | Descripció de la variable |
|---|---|---|
| “Tipus_ve hicle_esta ndaritzat” | Categòrica | Estandarització en només quatre categories de la variable “Tipus_vehicle”: “Vehicles d’ús professional”, “Vehicles motoritzats de quatre rodes”, “Vehicles motoritzats de dues rodes” i “Vehicles sense permís de conducció” |
| “Desc ripcio_caus a_mediata” | Categòrica | Detalla la causa mediata identificada per la corresponent patrulla de la Guàrdia Urbana que fer l’atestat de l’accident, fent referència al tipus de maniobra o circumstància immediata que va causar l’accident. |
| “N umero_lesio nats_greus” | Numèrica contínua | Nombre de persones ferides en l’accident que van requerir d’un temps d’hospitalització superior a les 24 hores. |
| “G reus_CODIF” | Categòrica | Estableix si va haver-hi un o més ferits greus en l’accident. |
| “Nom_mes” | Categòrica | Mes de l’any en que va tenir lloc l’accident. |
| “Vict imes_CODIF” | Categòrica | Estableix si va haver-hi dos o més víctimes en l’accident. |
| “Vehi cles_CODIF” | Categòrica | Estableix si va haver-hi dos o més vehicles implicats en l’accident. |
| “Edat” | Numèrica contínua | Anys d’edat complets del conductor en la data de l’accident. |
| ” Edat_CODIF” | Categòrica | Codificació de grups d’edat a partir dels valors interquartils: menors de 29 anys, més dels 29 fins els 37 anys, més dels 37 fins els 49 anys i majors de 49 anys. |
| “Es_ laborable” | Categòrica | Estableix si la data de l’accident era laborable durant l’any 2023. |
| “Es_o cupacional” | Categòrica | Estableix si l’hora de l’accident va ser dins d’una franja horària ocupacional, des de les 5 fins les 16h, en un dia laboral. |
També comprovem quin grau d’encert o precisió té aquest model al comprovar amb les mateixes dades amb les que s’ha entrenat, mesurant també la seva sensitivitat - o grau per encertar que sí es desplaçava per un motiu relacionat amb el món del treball- i especificitat - o grau per encertar que no es desplaçava per un motiu relacionat amb el món del treball:
pred <- predict(model, newdata=dadesSup, type='response')
dadesSup$prediction <- ifelse(pred < 0.5 ,0, 1)
prediccio <- as.factor(dadesSup$prediction)
valor_esperat <- dadesSup$EMT_recodificat
tb <- table(valor_esperat, prediccio); tb## prediccio
## valor_esperat 0 1
## 0 1061 446
## 1 508 808# TN: tb[1]
# FP: tb[2]
# FN: tb[3]
# TP: tb[4]
# accuracy
precisio <- (tb[4] + tb[1]) / (tb[1] + tb[2] + tb[3] + tb[4])
print(glue("\n\nLa precisió de la predicció és {round(precisio, 4)}."))##
## La precisió de la predicció és 0.6621.# sensitivity
sensitivitat <- tb[4] / (tb[4] + tb[3])
print(glue("\n\nLa sensitivitat de la predicció és {round(sensitivitat, 4)}."))##
## La sensitivitat de la predicció és 0.6443.# specifity
especificitat <- tb[1] / (tb[1] + tb[2])
print(glue("\n\nLa especificitat de la predicció és {round(especificitat, 4)}."))##
## La especificitat de la predicció és 0.6762.Observem que té una moderada capacitat predictiva tot i que l’especificitat és superior a la seva sensitivitat, tret que acostuma ser un indicador de sobreajust del model. Per aquest motiu, com també pel fet que l’objectiu d’emprar aquest algoritme de regressió logística era poder escollir, a partir de la seva significància estadística, quines variables tenen potencialment capacitat predictora i, per tant, no s’ha testejat amb una part de les dades que no hagi format part del model, el descartarem.
També pel cas d’aquelles variables que resulten en NA i que denotaria que el seu valors és zero o, atès a altres interpretacions, suggeriria multicolinealitat amb altres variables o classes, les interpretarem que tampoc no resulten explicatives.
Per altra banda, també observem que quatre de les variables, com és el cas tant de “Numero_lesionats_greus” i “Greus_CODIF” com també el de “Edat” i “Edat_CODIF”, deriven l’una de l’altra. Aquestes quatre variables es tindran en consideració a continuació per la modelització amb l’agorisme XGBoost però, pel cas dels arbres de decisió, només es podran tenir en consideració les dues variables categòriques “Greus_CODIF” i “Edat_CODIF”.
2.2 XGBoost (eXtreme Gradient Boosting)
A continuació posem en pràctica l’algorisme de XGBoost tot seguint l’exemple descrit en l’obra de P. Bruce, A. Bruce i P. Gedeck per treballar amb algorismes de boosting (BRUCE; BRUCE i GEDECK, 2022: 263-265), tal i com proposàvem en les conclusions de la segona part, aplicant-hi el fluxe de treball per s’hi executin 500 models diferents amb diverses parametritzacions aleatòries, tot realitzent el fluxe automatitzat el control del resultat per sel·leccionar el model més óptim a més de realitzar també aquestes comprovacions amb diverses porcions d’entrenament del conjunt de dades, tot seguint el l’exemple detallat per T. Pham en RPubs. Posteriorment, farem la classificació corresponent amb el conjunt de dades dels conductors accidentats dels qui desconeixem el motiu del seu desplaçament.
També bastirem el model amb les 11 variables que, atès a l’anàlisi realitzat en l’apartat anterior, sí resultarien explicatives per conèixer el resultat de la variable dependent “Es_mon_treball”; successivament també anirem comprovant quin joc de les variables relacionades resulten en un millor resultat en la classificació sobre els subconjunts de dades de testeig creat durant el fluxe.
En primer lloc, tindrem en consideració totes 11 variables:
data.xgb <- dadesSup[, c("Tipus_vehicle_estandaritzat",
"Descripcio_causa_mediata",
"Nom_mes",
"Numero_lesionats_greus",
"Greus_CODIF",
"Victimes_CODIF",
"Vehicles_CODIF",
"Edat_CODIF", "Edat", "Es_laborable", "Es_ocupacional",
"Es_mon_treball")]
set.seed(123)
cust_split <- data.xgb %>%
initial_split(prop = 0.8, strata = Es_mon_treball)
train <- training(cust_split)
test <- testing(cust_split)
# Cross validation folds from training dataset
set.seed(234)
folds <- vfold_cv(train, strata = Es_mon_treball)
cust_rec <- recipe(Es_mon_treball ~., data = train) %>%
# update_role(customerID, new_role = "ID") %>%
# step_corr(all_numeric()) %>%
step_corr(all_numeric(), threshold = 0.7, method = "spearman") %>%
step_zv(all_numeric()) %>% # filter zero variance
step_normalize(all_numeric()) %>%
step_dummy(all_nominal_predictors())# Setup a model specification
xgb_spec <-boost_tree(
trees = 500,
tree_depth = tune(),
min_n = tune(),
loss_reduction = tune(), ## first three: model complexity
sample_size = tune(), mtry = tune(), ## randomness
learn_rate = tune() ## step size
) %>%
set_engine("xgboost") %>%
set_mode("classification")
xgb_wf <- workflow() %>%
add_formula(Es_mon_treball ~.) %>%
add_model(xgb_spec)
xgb_grid <- grid_latin_hypercube(
tree_depth(),
min_n(),
loss_reduction(),
sample_size = sample_prop(),
finalize(mtry(), train),
learn_rate(),
size = 20
)
doParallel::registerDoParallel()
set.seed(234)
xgb_res <-tune_grid(
xgb_wf,
resamples = folds,
grid = xgb_grid,
control = control_grid(save_pred = TRUE)
)
xgb_res %>%
collect_metrics() %>%
filter(.metric == "roc_auc") %>%
select(mean, mtry:sample_size) %>%
pivot_longer(mtry:sample_size,
names_to = "parameter",
values_to = "value") %>%
ggplot(aes(value, mean, color = parameter)) +
geom_point(show.legend = FALSE)+
facet_wrap(~parameter, scales = "free_x")
show_best(xgb_res)## # A tibble: 5 × 12
## mtry min_n tree_depth learn_rate loss_reduction sample_size .metric
## <int> <int> <int> <dbl> <dbl> <dbl> <chr>
## 1 9 4 7 0.00000000628 0.0985 0.690 roc_auc
## 2 7 3 5 0.0000366 0.000384 0.386 roc_auc
## 3 8 14 2 0.0500 0.0000364 0.492 roc_auc
## 4 10 7 4 0.0224 2.45 0.136 roc_auc
## 5 11 18 6 0.000269 0.00000428 0.902 roc_auc
## # ℹ 5 more variables: .estimator <chr>, mean <dbl>, n <int>, std_err <dbl>,
## # .config <chr>best_auc <- select_best(xgb_res)
final_xgb <- finalize_workflow(xgb_wf, best_auc)
final_rs1 <- last_fit(final_xgb, cust_split,
metrics = metric_set(accuracy, roc_auc, sens,spec))
final_rs1 %>%
collect_metrics()## # A tibble: 4 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 accuracy binary 0.671 Preprocessor1_Model1
## 2 sens binary 0.725 Preprocessor1_Model1
## 3 spec binary 0.610 Preprocessor1_Model1
## 4 roc_auc binary 0.708 Preprocessor1_Model1En segon lloc, bastirem el model amb 9 variables, prescindint de les variables “Greus_CODIF” i d’“Edat_CODIF”:
best_auc <- select_best(xgb_res)
final_xgb <- finalize_workflow(xgb_wf, best_auc)
final_rs2 <- last_fit(final_xgb, cust_split,
metrics = metric_set(accuracy, roc_auc, sens,spec))
final_rs2 %>%
collect_metrics()## # A tibble: 4 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 accuracy binary 0.673 Preprocessor1_Model1
## 2 sens binary 0.699 Preprocessor1_Model1
## 3 spec binary 0.644 Preprocessor1_Model1
## 4 roc_auc binary 0.710 Preprocessor1_Model1En tercer lloc, bastirem el model amb 10 variables, prescindint de la variable “Numero_lesionats_greus”:
best_auc <- select_best(xgb_res)
final_xgb <- finalize_workflow(xgb_wf, best_auc)
final_rs3 <- last_fit(final_xgb, cust_split,
metrics = metric_set(accuracy, roc_auc, sens,spec))
final_rs3 %>%
collect_metrics()## # A tibble: 4 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 accuracy binary 0.673 Preprocessor1_Model1
## 2 sens binary 0.699 Preprocessor1_Model1
## 3 spec binary 0.644 Preprocessor1_Model1
## 4 roc_auc binary 0.710 Preprocessor1_Model1En quart lloc, bastirem el model amb 9 variables, prescindint de les variables “Numero_lesionats_greus” i “Edat”:
best_auc <- select_best(xgb_res)
final_xgb <- finalize_workflow(xgb_wf, best_auc)
final_rs4 <- last_fit(final_xgb, cust_split,
metrics = metric_set(accuracy, roc_auc, sens,spec))
final_rs4 %>%
collect_metrics()## # A tibble: 4 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 accuracy binary 0.666 Preprocessor1_Model1
## 2 sens binary 0.679 Preprocessor1_Model1
## 3 spec binary 0.652 Preprocessor1_Model1
## 4 roc_auc binary 0.702 Preprocessor1_Model1En cinqué lloc, bastirem el model amb 10 variables, prescindint de la variable “Greus_CODIF”:
best_auc <- select_best(xgb_res)
final_xgb <- finalize_workflow(xgb_wf, best_auc)
final_rs5 <- last_fit(final_xgb, cust_split,
metrics = metric_set(accuracy, roc_auc, sens,spec))
final_rs5 %>%
collect_metrics()## # A tibble: 4 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 accuracy binary 0.668 Preprocessor1_Model1
## 2 sens binary 0.712 Preprocessor1_Model1
## 3 spec binary 0.617 Preprocessor1_Model1
## 4 roc_auc binary 0.705 Preprocessor1_Model1best_auc <- select_best(xgb_res)
final_xgb <- finalize_workflow(xgb_wf, best_auc)
final_rs6 <- last_fit(final_xgb, cust_split,
metrics = metric_set(accuracy, roc_auc, sens,spec))
final_rs6 %>%
collect_metrics()## # A tibble: 4 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 accuracy binary 0.666 Preprocessor1_Model1
## 2 sens binary 0.679 Preprocessor1_Model1
## 3 spec binary 0.652 Preprocessor1_Model1
## 4 roc_auc binary 0.699 Preprocessor1_Model1I, en darrer lloc, bastirem el model amb 10 variables, prescindint de la variable “Edat_CODIF”:
best_auc <- select_best(xgb_res)
final_xgb <- finalize_workflow(xgb_wf, best_auc)
final_rs7 <- last_fit(final_xgb, cust_split,
metrics = metric_set(accuracy, roc_auc, sens,spec))
final_rs7 %>%
collect_metrics()## # A tibble: 4 × 4
## .metric .estimator .estimate .config
## <chr> <chr> <dbl> <chr>
## 1 accuracy binary 0.670 Preprocessor1_Model1
## 2 sens binary 0.699 Preprocessor1_Model1
## 3 spec binary 0.636 Preprocessor1_Model1
## 4 roc_auc binary 0.712 Preprocessor1_Model1I, a continuació, sintetitzem els resultats obtinguts en aquesta taula i que inclouran a més de la precisió, la sensitivitat i l’especificitat, també hi indicarem l’anomenada Àrea sota la característica operativa del receptor - Area under the receiver operating characteristics o (AUROC)- que ens descriu la proporció de l’àrea del conjunt de dades que estaria explicada per la funció de probabilitat del model:
| Model | Variables excloses | Precisió | Sensitivitat | Especificitat | AUROC |
|---|---|---|---|---|---|
| Model #1 | Cap | 0.6714 | 0.7252 | 0.6098 | 0.7077 |
| Model #2 | “Greus_CODIF”, “Edat_CODIF” |
0.6731 | 0.6987 | 0.6439 | 0.7099 |
| Model #3 | “Numero_lesionats_greus” | 0.6731 | 0.6987 | 0.6439 | 0.7104 |
| Model #4 | “Numero_lesionats_greus”, “Edat” |
0.6661 | 0.6788 | 0.6515 | 0.7016 |
| Model #5 | “Greus_CODIF” | 0.6678 | 0.7119 | 0.6174 | 0.7049 |
| Model #6 | “Edat” | 0.6661 | 0.6788 | 0.6515 | 0.6994 |
| Model #7 | “Edat_CODIF” | 0.6696 | 0.6987 | 0.6364 | 0.7120 |
En total s’han bastit fins set models de classificació diferents que es diferencien en l’ús o no de les dues variables numèriques contínues “Edat” i “Numerio_lesionats_greus” i les corresponents variables codificades “Edat_CODIF” i “Greus_CODIF”. En general, tots set models presenten una precisió similar i també en l’AUROC, observant-se en canvi les diferències més rellevants en la sensitivitat i l’especificitat corresponents, excepte en el cas dels Models #2 i #3 en els quals només s’observa la diferència en la seva precisió o en l’AUROC. tambés és constatable en tots els casos que la sensitivitat és superior a l’especificitat. Entre aquests, en el cas dels Models #3 s’observa que només es prescindeix en cadascuna d’una de les variables numèriques, mentres que en el cas del Model #2 es prescindeix de les dues variables categòriques codificades.
Finalment, escollirem el Model #2 atès que és el model de classificació que combina una elevada de propoció d’encerts positius i negatius i, per altra banda, perque pot resultar més fàcilment explicable atès que no utilitza dues vegades la mateixa variable, tant en la seva expressió numèrica com també la categorica.
final_wf <- final_rs2 %>%
extract_workflow()
desconegutsPred <- dadesDesconegutSup[, c("Tipus_vehicle_estandaritzat",
"Descripcio_causa_mediata",
"Nom_mes",
"Numero_lesionats_greus",
"Victimes_CODIF",
"Vehicles_CODIF",
"Edat", "Es_laborable", "Es_ocupacional")]
prediction.xgb <- predict(final_wf, desconegutsPred)
table(prediction.xgb)## .pred_class
## No Si
## 1113 1083En els resultats de la classificació del conjunt de dades de conductors del qui es desconeix el motiu per que es desplaçaven hi observem qu s’ha assignat pràcticament de forma equitativa ambdues classes.
2.3 Arbre de decisió (Decision tree)
A continuació, bastirem un model de classificació mitjançant l’algorisme de l’Arbre de decisió del módul C50 de R sense realitzar-ne cap poda ni cap canvi en els seu paràmetres per defecte:
En primer lloc, ho comprovarem incloent-hi també la variable “Nom_mes”:
dummy[] <- lapply(dummy, factor)
colsEscollides <- c("Tipus_vehicle_estandaritzat",
"Descripcio_causa_mediata",
"Nom_mes",
"Greus_CODIF",
"Victimes_CODIF",
"Vehicles_CODIF",
"Edat_CODIF", "Es_laborable", "Es_ocupacional")
y <- dummy$Es_mon_treball
colsSup <- colnames(dummy) %in% colsEscollides
x <- dummy[colsSup]
split_prop <- 5
indexes = sample(1:nrow(dummy),
size=floor(((split_prop-1)/split_prop)*nrow(dummy)))
train_x <- x[indexes, ]
train_y <- y[indexes]
test_x <- x[-indexes, ]
test_y <- y[-indexes]model1 <- C50::C5.0(train_x, train_y, rules=TRUE)
summary(model1)##
## Call:
## C5.0.default(x = train_x, y = train_y, rules = TRUE)
##
##
## C5.0 [Release 2.07 GPL Edition] Tue Aug 27 15:26:46 2024
## -------------------------------
##
## Class specified by attribute `outcome'
##
## Read 2258 cases (10 attributes) from undefined.data
##
## Rules:
##
## Rule 1: (1217/386, lift 1.3)
## Tipus_vehicle_estandaritzat in {Vehicles motoritzats de 2 rodes,
## Vehicles motoritzats de quatre rodes,
## Vehicles sense permis de conduccio}
## Es_ocupacional = No
## -> class No [0.683]
##
## Rule 2: (61/1, lift 2.1)
## Tipus_vehicle_estandaritzat = Vehicles Us Professional
## -> class Si [0.968]
##
## Rule 3: (999/391, lift 1.3)
## Es_ocupacional = Si
## -> class Si [0.608]
##
## Default class: No
##
##
## Evaluation on training data (2258 cases):
##
## Rules
## ----------------
## No Errors
##
## 3 778(34.5%) <<
##
##
## (a) (b) <-classified as
## ---- ----
## 831 392 (a): class No
## 386 649 (b): class Si
##
##
## Attribute usage:
##
## 98.14% Es_ocupacional
## 56.60% Tipus_vehicle_estandaritzat
##
##
## Time: 0.0 secsmodel1 <- C50::C5.0(train_x, train_y)
plot(model1, type="s", title="Fig. 4. Arbre de decisió sense podar.")
predicted_model1 <- predict(model1, test_x, type="class",
threshold=0.7)
print(sprintf("La precisió de l'arbre sense podar és del %.4f %%.",
100*sum(predicted_model1 == test_y) / length(predicted_model1)))## [1] "La precisió de l'arbre sense podar és del 65.1327 %."mat_conf <- table(test_y, Predicted=predicted_model1); mat_conf## Predicted
## test_y No Si
## No 183 101
## Si 96 185# TN: mat_conf[1]
# FP: mat_conf[2]
# FN: mat_conf[3]
# TP: mat_conf[4]
# sensitivity
sensitivitat <- mat_conf[4] / (mat_conf[4] + mat_conf[3])
print(glue("\n\nLa sensitivitat de la predicció és {round(sensitivitat, 4)}."))##
## La sensitivitat de la predicció és 0.6469.# specifity
especificitat <- mat_conf[1] / (mat_conf[1] + mat_conf[2])
print(glue("\n\nL'especificitat de la predicció és {round(especificitat, 4)}."))##
## L'especificitat de la predicció és 0.6559.En aquest cas, amb l’arbre de decisió sense podar s’observa una notable millora en la seva intel·legibilitat en contrast als models d’arbre de decisió bastits en la segona part, constatant-se que el model primordialment fa servir les variables “Tipus_vehicle_estandarditzat” i “Es_ocupacional” per realitzar la classificació. Tot i així, constatem també que la tendència continua sent que la sensitivitat sigui inferior a l’especificitat.
És necessari recordar que l’arbre de decisió aquí implementat presenta el defecte de la poca robustesa dels seus resultats. Com ja vam observar en l’anterior segona part, aquest conjunt de dades presenta grups de classes molt petits que, durant la creació del conjunt de dades d’entrenament mitjançant una sel·lecció aleatòria, poden resultar especialment infrarrepresentats o inclús deixats de banda en la seva totalitat. Aquest fenomen pot comportar que els resultats obtinguts en una iteració resultin diferents als de la següent iteració i, per tant, la seva fiabilitat es torna potencialment qüestionable. I, finalment, constatem que, en el cas de l’arbre de decisió, atès a que l’algorisme prioritza l’ús òptim de les variables amb més capacitat explicativa com són les del tipus de vehicle i si l’accident va tenir lluc o no durant un horari ocupacional, no hi tenen pràcticament cap efecte les altres set variables categòriques.
desconegutsPred <- dadesDesconegutSup[, c("Tipus_vehicle_estandaritzat",
"Descripcio_causa_mediata",
"Nom_mes",
"Greus_CODIF",
"Victimes_CODIF",
"Vehicles_CODIF",
"Edat_CODIF", "Es_laborable", "Es_ocupacional")]
predictions <- predict(model1, desconegutsPred, type="class", threshold=0.7)
table(predictions)## predictions
## No Si
## 1043 1153Per altra banda, crida l’atenció el fet que la proporció de classificacions negatives i positives és pràcticament idèntica a la proporció de classes positiva i negativa de la variable “Es_ocupacional” dels conductors dels qui es desconeix el motiu del seu desplaçament (v. Fig. 1.2 supra). Aquest fet ens suggereix que el model està emprant, en la pràctica, només la variable “Es_ocupacional” per fer aquesta classificació i, possiblement, sense aconseguir classificar als conductors que es desplacessin en dies no laborables o durant els horaris que no eran ocupacionals.
2.4 Arbre de decisió amb parametrització diferent
Realitzarem a continuació el mateix arbre de decisió però forçant que realitzi fins 99 iteracions a fi que escolleixi les prediccions més òptimes en cadascuna:
nTrials <-99
model_nTrials <- C50::C5.0(train_x, train_y, trials = nTrials)
predicted_modelnTrials <- predict(model_nTrials, test_x, type="class")
print(sprintf("La precisió de l'arbre amb 99 iteracions és del %.4f %%.",
100*sum(predicted_modelnTrials == test_y) / length(predicted_modelnTrials)))## [1] "La precisió de l'arbre amb 99 iteracions és del 64.6018 %."mat_conf <- table(test_y, Predicted=predicted_modelnTrials); mat_conf## Predicted
## test_y No Si
## No 199 85
## Si 115 166# TN: mat_conf[1]
# FP: mat_conf[2]
# FN: mat_conf[3]
# TP: mat_conf[4]
# sensitivity
sensitivitat <- mat_conf[4] / (mat_conf[4] + mat_conf[3])
print(glue("\n\nLa sensitivitat de la predicció és {round(sensitivitat, 4)}."))##
## La sensitivitat de la predicció és 0.6614.# specifity
especificitat <- mat_conf[1] / (mat_conf[1] + mat_conf[2])
print(glue("\n\nL'especificitat de la predicció és {round(especificitat, 4)}."))##
## L'especificitat de la predicció és 0.6338.desconegutsPred <- dadesDesconegutSup[, c("Tipus_vehicle_estandaritzat",
"Descripcio_causa_mediata",
"Nom_mes",
"Greus_CODIF",
"Victimes_CODIF",
"Vehicles_CODIF",
"Edat_CODIF", "Es_laborable", "Es_ocupacional")]
predictions <- predict(model_nTrials, desconegutsPred, type="class", threshold=0.7)
table(predictions)## predictions
## No Si
## 1190 1006En aquest model s’observa una millora en la sensitivitat, estant també aquest indicador per sobre de l’especificitat. Per altra banda, al fer la classificació de les dades dels conductors dels qui es desconeix el motiu del desplaçament, observem que la gran majoria de resultats són negatius.
2.5 Bosc aleatori (Random forest)
Tot seguit, també bastim un model de classificació del Bosc aleatoriq ue realitzarà 500 iteracions sobre el conjunt de dades d’entrenament:
library(rpart)##
## Adjuntando el paquete: 'rpart'## The following object is masked from 'package:dials':
##
## prunedummy[] <- lapply(dummy, factor)
colsEscollides <- c("Tipus_vehicle_estandaritzat",
"Descripcio_causa_mediata",
"Nom_mes",
"Greus_CODIF",
"Victimes_CODIF",
"Vehicles_CODIF",
"Edat_CODIF", "Es_laborable", "Es_ocupacional")
y <- dummy$Es_mon_treball
colsSup <- colnames(dummy) %in% colsEscollides
x <- dummy[colsSup]
y <- dummy$Es_mon_treball
rf <- randomForest(y ~., data = x, threshold=0.7); rf##
## Call:
## randomForest(formula = y ~ ., data = x, threshold = 0.7)
## Type of random forest: classification
## Number of trees: 500
## No. of variables tried at each split: 3
##
## OOB estimate of error rate: 37.27%
## Confusion matrix:
## No Si class.error
## No 962 545 0.3616457
## Si 507 809 0.3852584mat_conf <- rf$confusion
# sensitivity
sensitivitat <- mat_conf[4] / (mat_conf[4] + mat_conf[3])
print(glue("\n\nLa sensitivitat de la predicció és {round(sensitivitat, 4)}."))##
## La sensitivitat de la predicció és 0.5975.especificitat <- mat_conf[1] / (mat_conf[1] + mat_conf[2])
print(glue("\n\nLa especificitat de la predicció és {round(especificitat, 4)}."))##
## La especificitat de la predicció és 0.6549.predictions.rf <- predict( rf, desconegutsPred, type="class")
table(predictions.rf)## predictions.rf
## No Si
## 1064 1132Es constata que el resultat del model és molts similar a l’anterior on havíem modificat la parametrització amb 99 iteracions, tot i que el resultat de les classificacions resulta manifestament molt diferent.
2.6 Arbre de decisió podat (Prune tree)
Finalment, també provarem amb el mètode de l’arbre podat que ja vam possar en pràctica en la segona part. en aquest cas, ens limitarem a bastir el model mitjançant la funció tree() del módul Rpart de tal manera que seleccioni les tres branques que dinin un resultat més òptim.
dummy[] <- lapply(dummy, factor)
colsEscollides <- c("Tipus_vehicle_estandaritzat",
"Descripcio_causa_mediata",
"Nom_mes",
"Greus_CODIF",
"Victimes_CODIF",
"Vehicles_CODIF",
"Edat_CODIF", "Es_laborable", "Es_ocupacional",
"Es_mon_treball")
y <- dummy$Es_mon_treball
colsSup <- colnames(dummy) %in% colsEscollides
dummy_ <- dummy[colsSup]
split <- createDataPartition(y=dummy_$Es_mon_treball, p=4/5, list=FALSE)
train <- dummy_[split,]
test <- dummy_[-split,]
trees <- tree(Es_mon_treball~., train)prune.trees <- prune.tree(trees, best=3)tree.pred <- predict(prune.trees, test, type='class', threshold=0.7)
confusionMatrix(tree.pred, test$Es_mon_treball, positive = "Si")## Confusion Matrix and Statistics
##
## Reference
## Prediction No Si
## No 197 100
## Si 104 163
##
## Accuracy : 0.6383
## 95% CI : (0.5971, 0.678)
## No Information Rate : 0.5337
## P-Value [Acc > NIR] : 3.21e-07
##
## Kappa : 0.274
##
## Mcnemar's Test P-Value : 0.8336
##
## Sensitivity : 0.6198
## Specificity : 0.6545
## Pos Pred Value : 0.6105
## Neg Pred Value : 0.6633
## Prevalence : 0.4663
## Detection Rate : 0.2890
## Detection Prevalence : 0.4734
## Balanced Accuracy : 0.6371
##
## 'Positive' Class : Si
## prediction.prune.tree <- predict(prune.trees, desconegutsPred, type='class',
threshold=0.7)
table(prediction.prune.tree)## prediction.prune.tree
## No Si
## 1043 1153S’observa que en el cas del model de l’arbre podat els resultats són similars als del model bastit amb l’algorisme de l’arbre de decisió sense podar, sent idèntica també la classe positiva de la variable dependent “Es_mon_treball” al classificar el subconjunt de dades dels conductors dels qui es desconeix el motiu del seu desplaçament.
3 Elecció del model
En primer lloc presentarem la Taula amb els resultats més rellevants d’aquesta tercera part per, a continuació, justificar l’elecció del resultat del model XGBoost.
3.1 Taula de resultats
Només es detallaran els indicadors de precisió dels models de regressió logística, XGBoost i de l’arbre podat atès que en tots tres casos el model resultant és prou robust, es a dir, no mostra una variabilitat apreciable en els seus resultats.
| Índex | Descripció del resultat | Tipus de dades | Resultats |
|---|---|---|---|
| #1 | Total de conductors** | Enter | 5019 |
| #2 | Subtotal de conductors dels qui sí es coneix el motiu del seu d e splaçament, proporció sobre el total | Enter, percentil | 2823, 56.246 % |
| #3 | Nombre d’accidents amb conductors dels qui sí es coneix el motiu del seu d esplaçament | Enter | 2561 |
| #4 | Subtotal de conductors dels qui no es coneix el motiu del seu d e splaçament, proporció sobre el total | Enter, percentil | 2196, 43.754 % |
| #5 | Nombre d’accidents amb conductors dels qui no es coneix el motiu del seu d esplaçament | Enter | 2029 |
| #6 | Conductors accidenats en horari o cupacional, proporció sobre el total de conductors del que sí es coneixia el motiu del d esplaçament | Enter, percentil | 1272, 45.058 % |
| #7 | Conductors accidenats en horari no o cupacional, proporció sobre el total de conductors del que sí es coneixia el motiu del d esplaçament | Enter, percentil | 1551, 54.942 % |
| #8 | Conductors accidenats en horari o cupacional, proporció sobre el total de conductors del que no es coneixia el motiu del d esplaçament | Enter, percentil | 1150, 52.368 % |
| #9 | Conductors accidenats en horari no o cupacional, proporció sobre el total de conductors del que no es coneixia el motiu del d esplaçament | Enter, percentil | 1046, 47.632 % |
| #10 | Nombre de variables independents valorades per incloure en el model de classificació | Enter | 29 |
| #11 | Nombre de variables independents acceptades per incloure en els posteriors models de classificació | Enter | 11 |
| #12 | Llindar d’acceptació del resultat de la classe positiva | Decimal | 0.7 |
| #13 | En el model de Regressió Logística, la Sensitivitat era superior a l’ Especificitat | Booleà | Fals |
| #14 | Precisió, Sensitivitat i Especificitat del model Regressió Logística | Decimal, decimal, decimal | 0.6645, 0.6468, 0.6788 |
| #14 | En el model de XGBoost, la Sensitivitat era superior a l ’ Especificitat | Booleà | Cert |
| #15 | Precisió, Sensitivitat i Especificitat del model XGBoost | Decimal, decimal, decimal | 0.6731, 0.6987, 0.6439 |
| #16 | En el model d’Arbre de decisió sense pa r ametritzar, la Sensitivitat era superior a l ’ Especificitat | Booleà | Fals |
| #17 | En el model d’Arbre de decisió p a rametritzat per executar 99 cicles, la Sensitivitat era superior a l ’ Especificitat | Booleà | Fals |
| #18 | En el model del Bosc aleatori, la Sensitivitat era superior a l ’ Especificitat | Booleà | Fals |
| #19 | En el model de l’Arbre podat, la Sensitivitat era superior a l ’ Especificitat | Booleà | Fals |
| #20 | Precisió, Sensitivitat i Especificitat del model de l’Arbre podat | Decimal, decimal, decimal | 0.6525, 0.6084, 0.6910 |
3.2 Justificació
En aquesta tercera part d’aquesta tasca de classificació dels conductors ferits en accidents de trànsit durant l’any 2023, en certa mesura ens hem centrat a realitzar en la pràctica una imputació per aquells 2196 conductors dels qui es desconeixia el motiu del seu desplaçament. La rellevància de realitzar aquesta imputació la trobem en el fet que aquests representaven fins el 43% del total de conductors feris i dels qui disposàvem d’un joc complet de dades vàlides (v. fila #3 en Taula de resultats). Una proporció tant elevada tornava com poc fiable la variable del motiu dels desplaçament pel seu anàlisi però, amb la tasca realitzada fins ara valorem que sí podria ser ja utilitzable com una variable més per enriquir l’estudi les tendències en els accidents de trànsit amb ferits en els que hi intervingueren agents de la Guàrdia Urbana de Barcelona.
També cal valorar com rellevants les variables acceptades per la seva significància estadística per explicar el resultat que, en aquest cas, és indicar “Si”/“No” per cada conductor en relació a si el motiu del seu desplaçament en el moment de ser partícep en el corresponent accident documentat tenia relació amb el món del treball - v. detall del significat d’aquesta variable en la primera part d’aquest estudi-, no només per que l’agorisme corresponent el valori estadísticament òptim, si no també per que resulten coherents amb allò que coneixem del món del treball. Per això, resulta congruent que el tipus de vehicle - especialment els vehicles d’ús professional- o el fet que l’accident tingui lloc durant un dia laborable i un horari ocupacional siguin rellevants per considerar que el desplaçament era per un motiu laboral. Per altra banda, ja de forma més intuïtiva també resulta coherent que hi siguin rellevants l’edat o en quin mes o estació de l’any va tenir lloc, a més de tenir ja una clara rellevància analítica. I, en aquesta línia i, per tant, també per la seva posterior interpretació resulta d’interès que el nombre de ferits greus, si el nombre de víctimes o vehicles implicats fos superior o no a dos, a més de la causa mediata de l’accident també presentin relacions estadísticament relllevants amb un motiu de desplçament relacionat, tot i que no s’ha d’entendre necessàriament també com la seva causa. També resulta d’interès que ni la variable referent al sexe biològic ni tampoc si algun dels conductors implicats tenia menys de 5 anys d’experiència no hagin resultat significatives des del punt de vista estadístic a l’hora de classificar els resultats amb l’algorisme de Regressió logística. També s’ha descartat el cas del Districte on va tenir lloc, tot suggerint la conclusió que les divisions administratives del territori no resulten significatives per explicar la classificació per si el conductor es desplaçava per un motiu relacionat amb el món del treball o no.
En relació al model de XGBoost de classificació finalment escollit, la classificació implementada resulta robusta pel fet que no deriva d’una sola partició de dades per obtenir el conjunt de dades d’entrenament o una aprametrització escollida a l’atzar, si no que és l’escollida per que és la que explica una major proporció de resultats de la variable dependent “Es_mon_treball” atès a l’indicador conegut com Area under a curve (AUC) després de realitzar-ne 500 proves diferents, tot cercant optimitzar successivament el darrer resultat en relació a l’anteriorment obtingut, tot aplicant la ponderació als grups de classes més reduït, tot reduïnt el biaix ja detectat en el conjunt de dades en la segona part. Un tret que hi suma fiabilitat al resultat el trobem en el fet que l’únic dels models de classificació que, al comprovar la classificació proposada, demostra també una major proporció d’encerts pel cas del resultat positiu - Sí que el motiu del desplaçament està relacionat amb el món del treball- en comparació amb el cas del resultat negatiu.
4 Bibliografia
- BRUCE, P.; BRUCE, A.; i GEDECK, P. 2022. Estadística práctica para ciencia de datos con R y Python, Barcelona: Marcombo.
- INSTITUT METRÒPOLI. 2023. Enquesta de mobilitat en dia feiner 2023. (EMEF 2023), Barcelona: Autoritat del Transport Metropolità. URL (darrera visita el 13-VIII-2024).