Es inevitable notar una aprobación polarizada junto a la diferencia de aprobación entre municipios y no preguntarse si hay una correlación entre las diferencias sociodemográficas de los distintos municipios y la percepción ciudadana sobre las corporaciones policiales.
Aún y con los avances en materia de derechos humanos, es una realidad que en México se siguen vulnerando los derechos de algunos ciudadanos por motivos de raza, color de piel, género, orientación sexual, estrato socioeconómico, apariencia física, etc. “Según la Encuesta Nacional sobre Discriminación del 2017, el 20.2% de la población de 18 años y más ha sido víctima de algún acto discriminatorio […] 30% de las personas ha Sido marginado por su forma de vestir, arreglo personal o por usar tatuajes. Además, el 29.1% ha sido discriminado por su peso o estatura”. (UNAM Global, 2023)
La discriminación y el perfilamiento racial en encuentros civiles con policías es un tema con vastos estudios en países como Estados Unidos y Brasil, sin embargo poca literatura ni documentación oficial se encuentra sobre este mismo tema en México, aún así en años recientes se ha iniciado más la discusión sobre temas de discriminación en encuentros policiales. Ángel Placenta, rescata un estudio del 2017 por Amnistía Internacional (AI) en cual se concluye con que la discriminación policiál afecta principalmente a hombres jóvenes que viven en pobreza.
Según datos de la encuesta Cómo Vamos Nuevo León 2023, casi la mitad de la población confía en la policía de su colonia y considera que la toma de decisiones de la policía es justa, con 45% y 42.7%, respectivamente.
Así mismo, la encuesta señala que existe una gran diferencia de aprobación entre zonas: Pues San Pedro cuenta con una aprobación de más del 80% mientras que en municipios de la periferia esta está alrededor del 50%
Inicialmente, se quería investigar si estas diferencias en la aprobación podían deberse a algún tipo de discriminación o perfilamiento que experimentaban las personas al momento de tener algún encuentro con un cuerpo policial de su colonia. Placenta (2017) rescata un estudio de amnistía internacional que señala que los hombres jóvenes racializados en condiciones de pobreza son más propensos a vivir discriminación policial.
Sin embargo, la encuesta cómo vamos no cuenta con ninguna pregunta que haga referencia a la raza de una persona, auto adscripción indígena, o tonalidad de la piel. Debido a esto se buscó encontrar otros factores que fueran relevantes para responder ¿Cómo evaluará un individuo a su policía municipal?
Esta pregunta es relevante ya que la percepción de los individuos sobre sus cuerpos policiacos tiene diversas consecuencias en materias de seguridad. En primer lugar, Azfar y Gurgur (2008) señalan que una mala percepción de la policía se traduce en una baja denuncia ciudadana de delitos. Shelley (2001), indica que esta relación negativa puede resultar en la penetración de grupos relacionados con el crimen organizado dentro de la estructura estatal.
Otras consecuencias de la erosión en la relación entre el individuo y el cuerpo policial son un estado de insatisfacción con la democracia (Morales, 2009), y la reducción del capital social (Prats, 2008).
Debido a esto, es importante identificar los elementos que impactan en la calificación que le da un individuo a su policía municipal, ya que estos pueden representar un primer avance a sanar estas relaciones, mejorar las percepciones sobre democracia, e incluso aumentar la seguridad.
¿Cómo evaluará un individuo a su policía municipal según sus vivencias y percepciones de los cuerpos de seguridad?
Se ingresó a la base de datos de la encuesta Cómo Vamos Nuevo León 2023. Tras inspeccionar dicha base de datos se extrayeron los datos de preguntas que se creían relevantes para responder la pregunta de investigación. Las preguntas extraídas de la base de datos son las siguientes:
NOM_MUN_MV <- Municipio (queremos sólo AMM y periferia)
Sex_Persona_Entrevistada_MV <- Sexo de la persona entrevistada
CP4_1 <- años cumplidos
P3 <- Durante la semana pasada, ¿cuál fue su principal actividad laboral?
P74 <- ¿Considera que tiene sobrepeso u obesidad?
P93 <- ¿Qué tan seguro se siente en su municipio?
P94 <- ¿En qué lugar se siente más inseguro(a)?
P95 <- En el último año, usted o algún miembro del hogar ¿ha sido víctima de algún delito?
P99 <- (Para todos los encuestados), ¿Qué tan probable considera que los delitos denunciados sean investigados y castigados por las autoridades?
P100_1 <- ¿La presencia de policías en su colonia es suficiente?
P100_2 <- ¿La policía de su colonia lo hace sentir más seguro?
P100_3 <- ¿Confía en la policía de su colonia?
P100_4 <-¿El trato del policía al ciudadano es de forma respetuosa?
P100_5 <- ¿Considera que las decisiones y actuar de la policía benefician a las personas en su colonia?
P100_6 <- ¿La policía es justa al tomar decisiones?
P101_1 <- ¿Cómo evalúa el desempeño de los servicios de seguridad? En una escala del 1 al 10, en donde 1 es “muy malo” y 10 “muy bueno” Policía de su Municipio y Barrio
P101_2 <- Idem. Fuerza Civil
P101_3 <- Idem. Guardia Nacional
P101_4 <- Idem. Ejército y Marina
P101_5 <- Idem. Tránsito
P102 <- Durante el último año, ¿ha tenido contacto con la policía o ha requerido de ella?
P103 <- Si respondió que sí a la pregunta anterior, ¿con qué institución tuvo el contacto?
P104_1 <- La policía cumplió con su labor.
P104_2 <- La policía solicitó una mordida o moche.
P104_3 <- Si requirió de la policía, ¿ésta acudió al llamado?
P144 <- FINALMENTE: ¿Cuál es el ingreso mensual total del hogar?
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
setwd("C:/Users/Fercs/OneDrive/Tec de Monterrey/Uni/S4/Ciencia de datos/Victor")
eav<-read.csv("Base-de-datos_Encuesta-Asi-Vamos-2023 - csv.csv",encoding="latin1") #El encoding es para que lea los acentos
attach(eav)
df_eav<-data.frame(NOM_MUN_MV,Sex_Persona_Entrevistada_MV,CP4_1,P3,P74,P93,P94,P95,P96,P97,P98,P99,P100_1,P100_2,P100_3,P100_4,P100_5,P100_6,P101_1,P101_2,P101_3,P101_4,P101_5,P102,P103,P104_1,P104_2,P104_3,P144,Factor_CVNL,stringsAsFactors = TRUE) #se crea un nuevo dataframe que contenga las columnas de nuestro interés. Datos generales de la población, y variables que podrían estar relacionadas al tema de investigación
head(df_eav) #observar que tipo de variables tenemos
## NOM_MUN_MV Sex_Persona_Entrevistada_MV CP4_1 P3 P74 P93 P94 P95 P96 P97 P98
## 1 Apodaca Mujer 45 5 1 3 1 0 NA NA NA
## 2 Apodaca Hombre 68 7 0 1 7 0 NA NA NA
## 3 Apodaca Mujer 35 5 1 2 2 0 NA NA NA
## 4 Apodaca Hombre 35 1 0 4 1 1 0 NA NA
## 5 Apodaca Mujer 57 5 0 4 3 0 NA NA NA
## 6 Apodaca Hombre 25 8 0 1 2 1 0 NA NA
## P99 P100_1 P100_2 P100_3 P100_4 P100_5 P100_6 P101_1 P101_2 P101_3 P101_4
## 1 2 0 0 0 0 0 1 8 7 6 6
## 2 2 1 0 1 0 0 0 7 6 3 4
## 3 2 1 0 0 0 0 0 9 7 7 4
## 4 1 0 0 0 0 0 0 8 9 9 9
## 5 3 1 0 0 0 0 0 9 9 8 7
## 6 1 0 0 0 1 0 0 8 6 6 4
## P101_5 P102 P103 P104_1 P104_2 P104_3 P144 Factor_CVNL
## 1 6 0 NA NA NA NA 4 1359.894
## 2 4 0 NA NA NA NA 4 1428.555
## 3 7 0 NA NA NA NA 4 1359.894
## 4 7 0 NA NA NA NA 4 1428.555
## 5 7 0 NA NA NA NA 5 1359.894
## 6 6 1 3 1 0 0 5 1428.555
base<-df_eav%>% #aquí se emmpieza a hacer la limpieza de la base de datos
mutate(Factor_CVNL==as.numeric(Factor_CVNL))%>%
filter(NOM_MUN_MV%in%c("Apodaca","Cadereyta Jiménez","General Escobedo","García","Guadalupe","Juárez","Monterrey","San Nicolás de los Garza","San Pedro Garza García","Santa Catarina","Santiago","El Carmen","Ciénega de Flores","General Zuazua","Pesquería","Salinas Victoria"))%>%
mutate(municipio=as.factor(NOM_MUN_MV))%>% #se seleccionaron los muinicipios del AMM y periferia, y se cambió el nombre de la variable a "municipio" en forma de factor para que sea más fácil trabajar con ella
mutate(zona=as.factor(ifelse(municipio=="Apodaca","AMM",
ifelse(municipio=="Cadereyta Jiménez","AMM",
ifelse(municipio=="General Escobedo","AMM",
ifelse(municipio=="García","AMM",
ifelse(municipio=="Guadalupe","AMM",
ifelse(municipio=="Juárez","AMM",
ifelse(municipio=="Monterrey","AMM",
ifelse(municipio=="San Nicolás de los Garza","AMM",
ifelse(municipio=="San Pedro Garza García","AMM",
ifelse(municipio=="Santa Catarina","AMM",
ifelse(municipio=="Santiago","AMM",
ifelse(municipio=="El Carmen","Periferia",
ifelse(municipio=="Ciénega de Flores","Periferia",
ifelse(municipio=="General Zuazua","Periferia",
ifelse(municipio=="Pesquería","Periferia",
ifelse(municipio=="Salinas Victoria","Periferia","else"))))))))))))))))))%>%
filter(zona!="else")%>% #Tras haber filtrado por municipios, a cada municipio se le asigna si pertenece al AMM o a la periferia. Esto para evitar que existan muchas variables por cada municipio, eliminar el sesgo de la cantidad de observaciones, y establecer una comparación más simple entre el AMM y la periferia as wholes
mutate(sexo=as.factor(Sex_Persona_Entrevistada_MV))%>% # se renombra la variable a sexo y como factor
mutate(edad=as.numeric(CP4_1))%>% #Se crea la variable edad como dato numérico que sale de la pregunta CP4_1
mutate(empleo=as.factor(ifelse(P3=="1","empleado",
ifelse(P3=="2","Buscando empleo",
ifelse(P3=="3","Estudiante",
ifelse(P3=="4","Independiente",
ifelse(P3=="5","Trabajo doméstico no renumerado",
ifelse(P3=="6","Trabajo doméstico renumerado",
ifelse(P3=="7","Jubilado",
ifelse(P3=="8","NiNi",
ifelse(P3=="9999","NA","NA")))))))))))%>%
filter(empleo!="NA")%>% #De la pregunta P3 se obtiene la variable empleo, donde las respuestas registradas como numéricas se convierten a factores y se les asigna la etiqueta que les corresponde
mutate(sobrepeso=as.numeric(P74))%>%
filter(sobrepeso%in%c(0,1))%>% #La pregunta 74 sobre el sobrepeso, se conserva como dato numérico donde 0 es no y 1 es sí
mutate(sobrepeso=as.factor(ifelse(sobrepeso=="0","No","Sí")))%>%
mutate(nivel.inseguridad=as.numeric(P93))%>%
filter(nivel.inseguridad%in%c(1,2,3,4))%>% #El nivel de inseguridad se conserva como dato numérico y se conservan los registros 1,2,3,4 que corresponden a la percepción de inseguridad
mutate(lugar.inseguridad=as.factor(ifelse(P94=="1","Transporte Público",
ifelse(P94=="2","Vía Pública",
ifelse(P94=="3","Casa",
ifelse(P94=="4","Trabajo",
ifelse(P94=="6","Centros de entretenimiento",
ifelse(P94=="7","No me siento insegurx",
ifelse(P94=="8","En todos lados","NA")))))))))%>%
filter(lugar.inseguridad!="NA")%>% #como en el caso del empleo, se convierte en factor los lugares en donde las personas reportan sentirse más inseguras
mutate(victim=as.numeric(P95))%>%
filter(victim%in%c(0,1))%>% #Se registra si la persona ha sido víctima de un delito como dato numérico, donde 0 es no y 1 es sí
mutate(victim=as.factor(ifelse(victim=="0","No","Sí")))%>%
mutate(percepcion.delito.castigado=as.numeric(P99))%>%
filter(percepcion.delito.castigado%in%c(1,2,3,4))%>% #La pregunta 99 se mantiene como numérica y se conservan los valores 1,2,3,4 que corresponden a que tanto cree una persona que un delito cometido sea castigado
mutate(policia.suficiente=as.numeric(P100_1))%>%
filter(policia.suficiente%in%c(0,1))%>% # Como variable numérica, si la persona considera que la policía es suficiente con 0 es no y 1 es sí
mutate(policia.suficiente=as.factor(ifelse(policia.suficiente=="0","No","Sí")))%>%
mutate(policia.seguro=as.numeric(P100_2))%>%
filter(policia.seguro%in%c(0,1))%>% # Como variable numérica, si la persona considera que la policía le hace sentir seguro con 0 es no y 1 es sí
mutate(policia.seguro=as.factor(ifelse(policia.seguro=="0","No","Sí")))%>%
mutate(policia.confia=as.numeric(P100_3))%>%
filter(policia.confia%in%c(0,1))%>% # Como variable numérica, si la persona confía en la policía con 0 es no y 1 es sí
mutate(policia.confia=as.factor(ifelse(policia.confia=="0","No","Sí")))%>%
mutate(policia.respetuoso=as.numeric(P100_4))%>%
filter(policia.respetuoso%in%c(0,1))%>% # Como variable numérica, si la persona considera que la policía es respetuosa con 0 es no y 1 es sí
mutate(policia.respetuoso=as.factor(ifelse(policia.respetuoso=="0","No","Sí")))%>%
mutate(policia.beneficio=as.numeric(P100_5))%>%
filter(policia.beneficio%in%c(0,1))%>% # Como variable numérica, si la persona considera que la policía es un beneficio para la comunidad con 0 es no y 1 es sí
mutate(policia.beneficio=as.factor(ifelse(policia.beneficio=="0","No","Sí")))%>%
mutate(policia.justo=as.numeric(P100_6))%>%
filter(policia.justo%in%c(0,1))%>% # Como variable numérica, si la persona considera que la policía es justa con 0 es no y 1 es sí
mutate(policia.justo=as.factor(ifelse(policia.justo=="0","No","Sí")))%>%
mutate(policia.evaluacion=as.numeric(P101_1))%>%
filter(policia.evaluacion%in%c(1,2,3,4,5,6,7,8,9,10))%>% #La evaluación de la policía municipal se da como variable numérica que va del 1 al 10
mutate(fc.evaluacion=as.numeric(P101_2))%>%
filter(fc.evaluacion%in%c(1,2,3,4,5,6,7,8,9,10))%>% #La evaluación de la fuerza civil se da como variable numérica que va del 1 al 10
mutate(gn.evaluacion=as.numeric(P101_3))%>%
filter(gn.evaluacion%in%c(1,2,3,4,5,6,7,8,9,10))%>% #La evaluación de la guardia nacional se da como variable numérica que va del 1 al 10
mutate(ejercito.evaluacion=as.numeric(P101_4))%>%
filter(ejercito.evaluacion%in%c(1,2,3,4,5,6,7,8,9,10))%>% #La evaluación del ejército/marina se da como variable numérica que va del 1 al 10
mutate(transito.evaluacion=as.numeric(P101_5))%>%
filter(transito.evaluacion%in%c(1,2,3,4,5,6,7,8,9,10))%>% #La evaluación de tránsito se da como variable numérica que va del 1 al 10
mutate(contacto.con.policia=as.numeric(P102))%>%
filter(contacto.con.policia%in%c(0,1))%>% #Variable numérica que define si la persona ha tenido contacto con la policía, donde 0 es no y 1 es sí
mutate(contacto.con.policia=as.factor(ifelse(contacto.con.policia=="0","No","Sí")))%>%
mutate(maximo.ingreso.miles=as.numeric(ifelse(P144=="1",0,
ifelse(P144=="2",6.223,
ifelse(P144=="3",12.446,
ifelse(P144=="4",18.670,
ifelse(P144=="5",24.893,
ifelse(P144=="6",31.116,
ifelse(P144=="7",37.339,
ifelse(P144=="8",43.562,
ifelse(P144=="9",49.786,
ifelse(P144=="11",62.232,"NA"))))))))))))%>%
filter(maximo.ingreso.miles!="NA")%>% #Se define el ingreso máximo que podría tener el individuo de una observación según los rangos dados por el tabulado de Cómo vamos NL
select(municipio,zona,sexo,edad,empleo,sobrepeso,nivel.inseguridad,lugar.inseguridad,victim,
percepcion.delito.castigado,policia.suficiente,policia.seguro,policia.confia,policia.respetuoso,
policia.beneficio,policia.justo,policia.evaluacion,fc.evaluacion,gn.evaluacion,ejercito.evaluacion,
transito.evaluacion,contacto.con.policia,maximo.ingreso.miles,Factor_CVNL)
## Warning: There was 1 warning in `mutate()`.
## ℹ In argument: `maximo.ingreso.miles = as.numeric(...)`.
## Caused by warning:
## ! NAs introducidos por coerción
Se buscaron variables que tuvieran relación con la percepción de seguridad, la percepción de un individuo sobre las corporaciones policiales, y características demográficas que podrían brindar cruces importantes sobre dichas percepciones.
El municipio, sexo, edad, empleo, apariencia física, y nivel socioeconómico de una persona podrían causar que la experiencia del encuentro con corporaciones policiales difiera entre individuos.
Las percepciones de seguridad, el haber sido víctima de un delito, y el proceso de denuncia de dicho delito, también influyen en la relación de una persona con esta rama del poder ejecutivo. De igual forma, para poder determinar si una persona tendrá una percepción positiva del cuerpo policial es importante saber si se confía en este, se le considera respetuoso y justo, y se percibe que tiene un buen desmpeño a lo largo de todas las ramas.
Finalmente, es importante observar si al momento de tener contacto con un cuerpo policial este cumplió con su labor en tiempo y forma, y cómo esto puede afectar la percepción de un individuo sobre la corporación.
lm.fit<-lm(policia.evaluacion~.,base,weight=base$Factor_CVNL)
summary(lm.fit)
##
## Call:
## lm(formula = policia.evaluacion ~ ., data = base, weights = base$Factor_CVNL)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -170.63 -13.94 0.58 15.93 136.47
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value
## (Intercept) 1.9329030 1.3657970 1.415
## municipioCadereyta Jiménez 0.4867221 1.1009798 0.442
## municipioCiénega de Flores -0.1610397 0.7421961 -0.217
## municipioEl Carmen 0.1156508 0.7181380 0.161
## municipioGarcía 0.2245646 0.6088262 0.369
## municipioGeneral Escobedo 0.1379785 0.3908645 0.353
## municipioGeneral Zuazua 0.2058757 0.7235054 0.285
## municipioGuadalupe -0.0189805 0.1038253 -0.183
## municipioJuárez 0.3763354 0.4499132 0.836
## municipioMonterrey -0.9219890 1.1392501 -0.809
## municipioPesquería -0.0520044 0.7169792 -0.073
## municipioSalinas Victoria 0.3808056 0.7232934 0.526
## municipioSan Nicolás de los Garza 0.4281774 0.4148114 1.032
## municipioSan Pedro Garza García 0.7060981 1.0618305 0.665
## municipioSanta Catarina 0.3582859 0.7182617 0.499
## municipioSantiago 0.6766350 1.2577533 0.538
## zonaPeriferia NA NA NA
## sexoMujer -0.1074712 0.0552447 -1.945
## edad -0.0024557 0.0017675 -1.389
## empleoempleado -0.0909420 0.1445664 -0.629
## empleoEstudiante -0.0275391 0.2078568 -0.132
## empleoIndependiente 0.0679972 0.1547833 0.439
## empleoJubilado 0.1113733 0.1604613 0.694
## empleoNiNi 0.0263177 0.1542096 0.171
## empleoTrabajo doméstico no renumerado 0.1507369 0.1522999 0.990
## empleoTrabajo doméstico renumerado 0.2406390 0.2802848 0.859
## sobrepesoSí -0.1240557 0.0600692 -2.065
## nivel.inseguridad -0.0695736 0.0351056 -1.982
## lugar.inseguridadCentros de entretenimiento -0.7533222 0.1356079 -5.555
## lugar.inseguridadEn todos lados -0.8353163 0.9943868 -0.840
## lugar.inseguridadNo me siento insegurx -0.6330492 0.1076874 -5.879
## lugar.inseguridadTrabajo -0.6878224 0.2969350 -2.316
## lugar.inseguridadTransporte Público -0.3541045 0.1211387 -2.923
## lugar.inseguridadVía Pública -0.5335114 0.0959323 -5.561
## victimSí 0.0224361 0.0770265 0.291
## percepcion.delito.castigado -0.0953603 0.0245151 -3.890
## policia.suficienteSí 0.1999867 0.0658837 3.035
## policia.seguroSí -0.0249672 0.0948688 -0.263
## policia.confiaSí 0.3350774 0.0926628 3.616
## policia.respetuosoSí 0.1423872 0.0827082 1.722
## policia.beneficioSí 0.1346959 0.0855979 1.574
## policia.justoSí 0.2377044 0.0822923 2.889
## fc.evaluacion 0.5059118 0.0205504 24.618
## gn.evaluacion 0.1011377 0.0276756 3.654
## ejercito.evaluacion -0.1865313 0.0236258 -7.895
## transito.evaluacion 0.2640346 0.0116766 22.612
## contacto.con.policiaSí -0.1340136 0.0951668 -1.408
## maximo.ingreso.miles 0.0176458 0.0033287 5.301
## Factor_CVNL 0.0006213 0.0009608 0.647
## Pr(>|t|)
## (Intercept) 0.157130
## municipioCadereyta Jiménez 0.658469
## municipioCiénega de Flores 0.828244
## municipioEl Carmen 0.872073
## municipioGarcía 0.712272
## municipioGeneral Escobedo 0.724112
## municipioGeneral Zuazua 0.776010
## municipioGuadalupe 0.854961
## municipioJuárez 0.402975
## municipioMonterrey 0.418423
## municipioPesquería 0.942184
## municipioSalinas Victoria 0.598596
## municipioSan Nicolás de los Garza 0.302069
## municipioSan Pedro Garza García 0.506124
## municipioSanta Catarina 0.617948
## municipioSantiago 0.590645
## zonaPeriferia NA
## sexoMujer 0.051844 .
## edad 0.164840
## empleoempleado 0.529363
## empleoEstudiante 0.894607
## empleoIndependiente 0.660478
## empleoJubilado 0.487696
## empleoNiNi 0.864504
## empleoTrabajo doméstico no renumerado 0.322399
## empleoTrabajo doméstico renumerado 0.390671
## sobrepesoSí 0.039006 *
## nivel.inseguridad 0.047607 *
## lugar.inseguridadCentros de entretenimiento 3.07e-08 ***
## lugar.inseguridadEn todos lados 0.400972
## lugar.inseguridadNo me siento insegurx 4.69e-09 ***
## lugar.inseguridadTrabajo 0.020617 *
## lugar.inseguridadTransporte Público 0.003497 **
## lugar.inseguridadVía Pública 2.96e-08 ***
## victimSí 0.770863
## percepcion.delito.castigado 0.000103 ***
## policia.suficienteSí 0.002427 **
## policia.seguroSí 0.792437
## policia.confiaSí 0.000305 ***
## policia.respetuosoSí 0.085274 .
## policia.beneficioSí 0.115710
## policia.justoSí 0.003904 **
## fc.evaluacion < 2e-16 ***
## gn.evaluacion 0.000263 ***
## ejercito.evaluacion 4.31e-15 ***
## transito.evaluacion < 2e-16 ***
## contacto.con.policiaSí 0.159198
## maximo.ingreso.miles 1.25e-07 ***
## Factor_CVNL 0.517901
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 30.05 on 2479 degrees of freedom
## Multiple R-squared: 0.7235, Adjusted R-squared: 0.7183
## F-statistic: 138 on 47 and 2479 DF, p-value: < 2.2e-16
lm.fit2<-lm(policia.evaluacion~.-municipio,base,weight=base$Factor_CVNL) #se elimina la variable municipio para no ser redundante con zona, lo cual permite calcular pvalues y betas para las dos zonas
summary(lm.fit2)
##
## Call:
## lm(formula = policia.evaluacion ~ . - municipio, data = base,
## weights = base$Factor_CVNL)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -167.582 -14.158 0.037 15.042 145.341
##
## Coefficients:
## Estimate Std. Error t value
## (Intercept) 2.853e+00 2.405e-01 11.862
## zonaPeriferia -3.187e-01 7.423e-02 -4.293
## sexoMujer -1.032e-01 5.451e-02 -1.893
## edad -1.979e-03 1.749e-03 -1.132
## empleoempleado -9.206e-02 1.442e-01 -0.639
## empleoEstudiante -1.047e-02 2.069e-01 -0.051
## empleoIndependiente 6.397e-02 1.547e-01 0.414
## empleoJubilado 1.204e-01 1.604e-01 0.751
## empleoNiNi 5.453e-02 1.535e-01 0.355
## empleoTrabajo doméstico no renumerado 1.521e-01 1.514e-01 1.005
## empleoTrabajo doméstico renumerado 1.884e-01 2.803e-01 0.672
## sobrepesoSí -1.160e-01 6.005e-02 -1.932
## nivel.inseguridad -7.344e-02 3.462e-02 -2.121
## lugar.inseguridadCentros de entretenimiento -7.877e-01 1.351e-01 -5.830
## lugar.inseguridadEn todos lados -6.973e-01 9.950e-01 -0.701
## lugar.inseguridadNo me siento insegurx -6.636e-01 1.062e-01 -6.248
## lugar.inseguridadTrabajo -6.913e-01 2.966e-01 -2.330
## lugar.inseguridadTransporte Público -3.718e-01 1.192e-01 -3.121
## lugar.inseguridadVía Pública -5.700e-01 9.460e-02 -6.025
## victimSí 2.108e-02 7.713e-02 0.273
## percepcion.delito.castigado -1.070e-01 2.411e-02 -4.438
## policia.suficienteSí 2.273e-01 6.561e-02 3.464
## policia.seguroSí -2.157e-02 9.475e-02 -0.228
## policia.confiaSí 3.439e-01 9.268e-02 3.711
## policia.respetuosoSí 1.259e-01 8.245e-02 1.527
## policia.beneficioSí 1.191e-01 8.541e-02 1.394
## policia.justoSí 2.368e-01 8.238e-02 2.874
## fc.evaluacion 5.045e-01 2.055e-02 24.548
## gn.evaluacion 9.832e-02 2.768e-02 3.552
## ejercito.evaluacion -1.875e-01 2.357e-02 -7.955
## transito.evaluacion 2.665e-01 1.161e-02 22.961
## contacto.con.policiaSí -1.285e-01 9.519e-02 -1.350
## maximo.ingreso.miles 1.968e-02 3.230e-03 6.092
## Factor_CVNL -5.182e-05 3.486e-05 -1.486
## Pr(>|t|)
## (Intercept) < 2e-16 ***
## zonaPeriferia 1.83e-05 ***
## sexoMujer 0.058504 .
## edad 0.257841
## empleoempleado 0.523134
## empleoEstudiante 0.959644
## empleoIndependiente 0.679158
## empleoJubilado 0.452869
## empleoNiNi 0.722444
## empleoTrabajo doméstico no renumerado 0.315233
## empleoTrabajo doméstico renumerado 0.501491
## sobrepesoSí 0.053457 .
## nivel.inseguridad 0.033986 *
## lugar.inseguridadCentros de entretenimiento 6.27e-09 ***
## lugar.inseguridadEn todos lados 0.483497
## lugar.inseguridadNo me siento insegurx 4.89e-10 ***
## lugar.inseguridadTrabajo 0.019868 *
## lugar.inseguridadTransporte Público 0.001825 **
## lugar.inseguridadVía Pública 1.94e-09 ***
## victimSí 0.784621
## percepcion.delito.castigado 9.47e-06 ***
## policia.suficienteSí 0.000540 ***
## policia.seguroSí 0.819927
## policia.confiaSí 0.000211 ***
## policia.respetuosoSí 0.127006
## policia.beneficioSí 0.163358
## policia.justoSí 0.004084 **
## fc.evaluacion < 2e-16 ***
## gn.evaluacion 0.000390 ***
## ejercito.evaluacion 2.69e-15 ***
## transito.evaluacion < 2e-16 ***
## contacto.con.policiaSí 0.177207
## maximo.ingreso.miles 1.29e-09 ***
## Factor_CVNL 0.137308
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 30.16 on 2493 degrees of freedom
## Multiple R-squared: 0.72, Adjusted R-squared: 0.7163
## F-statistic: 194.3 on 33 and 2493 DF, p-value: < 2.2e-16
Finalmente, el modelo que decidimos usar fue tomando las zonas geográficas como zona metropolitana y periferia en vez de hacerlo por cada municipio. A continuación las variables que resultaron significativas y sus coeficientes:
El Intercepto (2.762) representa que de inicio, una persona evaluaría la policía municipal con un 2.7 de 10 puntos posibles. Sin embargo si la persona vive en un municipio en periferia la evaluación a la policía municipal disminuye -2.695e-01, si la persona es mujer también representa una relación negativa de -1.156e-01. El sobrepeso en una persona también significa que es probable que la evaluación al policía disminuya(-1.230e-01). Hablando sobre percepción de inseguridad, las siguientes variables representan un decremento en la evaluación a la policía municipal: nivel.inseguridad -1.124e-01 lugar.inseguridadCentros de entretenimiento -9.766e-01 lugar.inseguridadNo me siento insegurx -7.779e-01 lugar.inseguridadTransporte Público -4.302e-01 lugar.inseguridadVía Pública -6.450e-01 percepcion.delito.castigado -1.040e-01 Entre mayor sea la percepción de la inseguridad de una persona en los ámbitos previos, más disminuirá la evaluación que le da a la policía municipal en el grado del coeficiente de cada uno.
Las últimas dos variables que muestran una relación negativa con la evaluación al policía municipal son contacto.con.policia (-1.984e-01) y ejercito.evaluacion (-1.741e-01), es decir, si han tenido algún contacto con la policía y si evaluaron bien al ejército.
Por el otro lado, aquellas varaibels que representarían un aumento en la evaluación de un ciudadano hacia la policía municipal son: que haya suficientes policías (2.144e-01), que tenga confianza a la policía (4.216e-01), que perciba a la policía como justa en sus decisiones (2.659e-01), que evalúe bien al policía de tránsito (2.660e-01), entre mayor sea el ingreso máximo del hogar irá aumentando la evaluación en 2.243e-05, por útlimo fc.evaluacion (5.043e-01) y gn.evaluacion (1.034e-01)
library (boot) #Se llama a la librería boot
glm.fit<-glm(policia.evaluacion~.-municipio-empleo-lugar.inseguridad,data=base,weight=base$Factor_CVNL) #Se crea un modelo con glm, pero lineal, para poder usar la función cv.glm
base_sin<-subset(base,select=-c(lugar.inseguridad))#Se crea una base sin lugar de inseguridad para corregir el error de nuevo nivel en el factor
glm.fit2<-glm(policia.evaluacion~.-municipio-empleo,data=base_sin) #Se crea nueva función glm con la nueva base de datos que corrige el error de los niveles
cv.err<-cv.glm(base_sin, glm.fit2) #Se corre la función cv.glm para obtener la lista de varios componentes
cv.err$delta #Se obtienen los números del vector delta para obtener los MSE
## [1] 1.252437 1.252432
En este caso, para nuestra Validación Cruzada LOOCV los valores de los errores cuadráticos medios obtenidos fueron de 1.252268 y 1.252263. Estos representan una estimación del error de predicción promedio del modelo cuando se deja una observación fuera en cada iteración. En estos casos, un error más alto indicaría una mayor discrepancia entre los valores predichos por el modelo y los valores reales, por lo que en nuestra validación los valores sugieren que el modelo tiene un error de predicción promedio de alrededor de 1.25.
Esto significa que las predicciones del modelo difieren en aproximadamente 1.25 unidades de los valores reales de la variable de respuesta (policia.evaluacion), lo cual nos indica una buena precisión y confiabilidad, pues estas unidades no son muy altas y por ende son muy cercanas a los valores reales.
cv.error.10<-rep(0,10)
for (i in 1:10) {
glm.fitK<-glm(policia.evaluacion~.-municipio-empleo,data=base_sin)
cv.error.10[i]<-cv.glm(base_sin,glm.fitK,K=50)$delta[1]
}
cv.error.10
## [1] 1.252014 1.250735 1.252750 1.254858 1.252564 1.253506 1.251784 1.251960
## [9] 1.252567 1.254497
# Se hace un método de K-fold para valores del 1 al 10, usando la misma regresión logística que se usó en el LOOCV. Se utilizó una K de 50
En nuestra Validación Cruzada k-fold los valores de los errores cuadráticos medios obtenidos con K=50 tienen un promedio de 1.2530263. Aquí podemos destacar que estos son sumamente similares a los obtenidos con LOOCV, lo cual sugiere una estimación consistente del error de predicción del modelo.
Al igual que en LOOCV, estos valores indican que el modelo tiene un error de predicción promedio de alrededor de 1.25 unidades en la escala de la variable de respuesta por lo que podemos afirmar que la consistencia entre los diferentes enfoques de validación cruzada aumenta la confianza en la estimación del error de predicción del modelo, sobre todo con los valores bajos que tenemos, los cuales nos indican mayor precisión.
library(boot)
boot.fn<-function(data,index){
return(coef(lm(policia.evaluacion~.-municipio-empleo-lugar.inseguridad,base,weights=base$Factor_CVNL,subset=index)))
} # Se eliminaron las variables que estaban en forma de factores para que al hacer el bootstrap los items a reemplazar fueran múltiplos de la longitud del reemplazo
boot.fn(base,1:2527)
## (Intercept) zonaPeriferia
## 2.1691085136 -0.3230965321
## sexoMujer edad
## 0.0158225199 0.0008230182
## sobrepesoSí nivel.inseguridad
## -0.1031373116 -0.0292313481
## victimSí percepcion.delito.castigado
## 0.0287180062 -0.1354244000
## policia.suficienteSí policia.seguroSí
## 0.2213406603 -0.0520961334
## policia.confiaSí policia.respetuosoSí
## 0.3249655213 0.0797599523
## policia.beneficioSí policia.justoSí
## 0.1403039661 0.2225236741
## fc.evaluacion gn.evaluacion
## 0.5096322861 0.0930886013
## ejercito.evaluacion transito.evaluacion
## -0.2067522665 0.2827762539
## contacto.con.policiaSí maximo.ingreso.miles
## -0.1075892014 0.0225467893
## Factor_CVNL
## -0.0000500272
boot(base,boot.fn,1000)
##
## ORDINARY NONPARAMETRIC BOOTSTRAP
##
##
## Call:
## boot(data = base, statistic = boot.fn, R = 1000)
##
##
## Bootstrap Statistics :
## original bias std. error
## t1* 2.1691085136 -1.875325e-02 2.355788e-01
## t2* -0.3230965321 6.111232e-03 8.517414e-02
## t3* 0.0158225199 7.219316e-04 4.778517e-02
## t4* 0.0008230182 4.263344e-05 1.450504e-03
## t5* -0.1031373116 -3.810777e-03 7.058079e-02
## t6* -0.0292313481 9.644863e-05 4.223690e-02
## t7* 0.0287180062 -6.791095e-03 1.048571e-01
## t8* -0.1354244000 2.384917e-04 2.763776e-02
## t9* 0.2213406603 2.463531e-03 7.820508e-02
## t10* -0.0520961334 -9.749909e-03 1.297922e-01
## t11* 0.3249655213 5.270341e-03 1.137934e-01
## t12* 0.0797599523 3.045622e-03 9.791678e-02
## t13* 0.1403039661 -5.891288e-03 1.087351e-01
## t14* 0.2225236741 1.818576e-03 1.086738e-01
## t15* 0.5096322861 6.895063e-04 3.138758e-02
## t16* 0.0930886013 1.919365e-03 3.699168e-02
## t17* -0.2067522665 -1.042293e-03 2.969408e-02
## t18* 0.2827762539 -3.195061e-04 1.587322e-02
## t19* -0.1075892014 9.545082e-04 1.332165e-01
## t20* 0.0225467893 2.843013e-04 3.756342e-03
## t21* -0.0000500272 1.429607e-06 3.766646e-05
Finalmente, nuestro bootstrap nos proporciona información sobre la incertidumbre y la precisión de las estimaciones de los diferentes coeficientes del modelo. Por ejemplo, para el coeficiente de la variable “zona” (t2), el valor original es -2.681423e-01, con un sesgo estimado de 5.846232e-04 y un error estándar de 8.492677e-02. Esto significa que, aunque el valor original del coeficiente es -0.2681423, el procedimiento de bootstrap sugiere que podría haber un pequeño sesgo de 0.0005846232 en esta estimación y que el verdadero valor del coeficiente podría estar dentro de un intervalo de aproximadamente -0.2681423 ± 1.96 * 0.08492677 con un 95% de confianza.
Podemos interpretar los sesgos y errores estándar de manera similar para los demás coeficientes y construir intervalos de confianza para evaluar la precisión de las estimaciones y cómo podemos observar, en su gran mayoría los sesgos de las variables son sumamente pequeños, indicando así una buena precisión en la mayoría de variables.
library(leaps)
#Mejor selección de subconjuntos con nivel máximo de 19 debido a las 19 variables predictoras
regfit.full <- regsubsets (policia.evaluacion~.-municipio-empleo-lugar.inseguridad,base ,weights=base$Factor_CVNL,nvmax =19)
reg.summary <- summary (regfit.full)
reg.summary
## Subset selection object
## Call: regsubsets.formula(policia.evaluacion ~ . - municipio - empleo -
## lugar.inseguridad, base, weights = base$Factor_CVNL, nvmax = 19)
## 20 Variables (and intercept)
## Forced in Forced out
## zonaPeriferia FALSE FALSE
## sexoMujer FALSE FALSE
## edad FALSE FALSE
## sobrepesoSí FALSE FALSE
## nivel.inseguridad FALSE FALSE
## victimSí FALSE FALSE
## percepcion.delito.castigado FALSE FALSE
## policia.suficienteSí FALSE FALSE
## policia.seguroSí FALSE FALSE
## policia.confiaSí FALSE FALSE
## policia.respetuosoSí FALSE FALSE
## policia.beneficioSí FALSE FALSE
## policia.justoSí FALSE FALSE
## fc.evaluacion FALSE FALSE
## gn.evaluacion FALSE FALSE
## ejercito.evaluacion FALSE FALSE
## transito.evaluacion FALSE FALSE
## contacto.con.policiaSí FALSE FALSE
## maximo.ingreso.miles FALSE FALSE
## Factor_CVNL FALSE FALSE
## 1 subsets of each size up to 19
## Selection Algorithm: exhaustive
## zonaPeriferia sexoMujer edad sobrepesoSí nivel.inseguridad victimSí
## 1 ( 1 ) " " " " " " " " " " " "
## 2 ( 1 ) " " " " " " " " " " " "
## 3 ( 1 ) " " " " " " " " " " " "
## 4 ( 1 ) " " " " " " " " " " " "
## 5 ( 1 ) " " " " " " " " " " " "
## 6 ( 1 ) " " " " " " " " " " " "
## 7 ( 1 ) " " " " " " " " " " " "
## 8 ( 1 ) "*" " " " " " " " " " "
## 9 ( 1 ) "*" " " " " " " " " " "
## 10 ( 1 ) "*" " " " " " " " " " "
## 11 ( 1 ) "*" " " " " " " " " " "
## 12 ( 1 ) "*" " " " " "*" " " " "
## 13 ( 1 ) "*" " " " " "*" " " " "
## 14 ( 1 ) "*" " " " " "*" " " " "
## 15 ( 1 ) "*" " " " " "*" " " " "
## 16 ( 1 ) "*" " " " " "*" "*" " "
## 17 ( 1 ) "*" " " " " "*" "*" " "
## 18 ( 1 ) "*" " " "*" "*" "*" " "
## 19 ( 1 ) "*" " " "*" "*" "*" "*"
## percepcion.delito.castigado policia.suficienteSí policia.seguroSí
## 1 ( 1 ) " " " " " "
## 2 ( 1 ) " " " " " "
## 3 ( 1 ) " " " " " "
## 4 ( 1 ) " " " " " "
## 5 ( 1 ) " " " " " "
## 6 ( 1 ) "*" " " " "
## 7 ( 1 ) "*" " " " "
## 8 ( 1 ) "*" " " " "
## 9 ( 1 ) "*" "*" " "
## 10 ( 1 ) "*" "*" " "
## 11 ( 1 ) "*" "*" " "
## 12 ( 1 ) "*" "*" " "
## 13 ( 1 ) "*" "*" " "
## 14 ( 1 ) "*" "*" " "
## 15 ( 1 ) "*" "*" " "
## 16 ( 1 ) "*" "*" " "
## 17 ( 1 ) "*" "*" "*"
## 18 ( 1 ) "*" "*" "*"
## 19 ( 1 ) "*" "*" "*"
## policia.confiaSí policia.respetuosoSí policia.beneficioSí
## 1 ( 1 ) " " " " " "
## 2 ( 1 ) " " " " " "
## 3 ( 1 ) "*" " " " "
## 4 ( 1 ) "*" " " " "
## 5 ( 1 ) "*" " " " "
## 6 ( 1 ) "*" " " " "
## 7 ( 1 ) "*" " " " "
## 8 ( 1 ) "*" " " " "
## 9 ( 1 ) "*" " " " "
## 10 ( 1 ) "*" " " " "
## 11 ( 1 ) "*" " " "*"
## 12 ( 1 ) "*" " " "*"
## 13 ( 1 ) "*" " " "*"
## 14 ( 1 ) "*" " " "*"
## 15 ( 1 ) "*" "*" "*"
## 16 ( 1 ) "*" "*" "*"
## 17 ( 1 ) "*" "*" "*"
## 18 ( 1 ) "*" "*" "*"
## 19 ( 1 ) "*" "*" "*"
## policia.justoSí fc.evaluacion gn.evaluacion ejercito.evaluacion
## 1 ( 1 ) " " "*" " " " "
## 2 ( 1 ) " " "*" " " " "
## 3 ( 1 ) " " "*" " " " "
## 4 ( 1 ) " " "*" " " "*"
## 5 ( 1 ) " " "*" " " "*"
## 6 ( 1 ) " " "*" " " "*"
## 7 ( 1 ) "*" "*" " " "*"
## 8 ( 1 ) "*" "*" " " "*"
## 9 ( 1 ) "*" "*" " " "*"
## 10 ( 1 ) "*" "*" "*" "*"
## 11 ( 1 ) "*" "*" "*" "*"
## 12 ( 1 ) "*" "*" "*" "*"
## 13 ( 1 ) "*" "*" "*" "*"
## 14 ( 1 ) "*" "*" "*" "*"
## 15 ( 1 ) "*" "*" "*" "*"
## 16 ( 1 ) "*" "*" "*" "*"
## 17 ( 1 ) "*" "*" "*" "*"
## 18 ( 1 ) "*" "*" "*" "*"
## 19 ( 1 ) "*" "*" "*" "*"
## transito.evaluacion contacto.con.policiaSí maximo.ingreso.miles
## 1 ( 1 ) " " " " " "
## 2 ( 1 ) "*" " " " "
## 3 ( 1 ) "*" " " " "
## 4 ( 1 ) "*" " " " "
## 5 ( 1 ) "*" " " "*"
## 6 ( 1 ) "*" " " "*"
## 7 ( 1 ) "*" " " "*"
## 8 ( 1 ) "*" " " "*"
## 9 ( 1 ) "*" " " "*"
## 10 ( 1 ) "*" " " "*"
## 11 ( 1 ) "*" " " "*"
## 12 ( 1 ) "*" " " "*"
## 13 ( 1 ) "*" " " "*"
## 14 ( 1 ) "*" "*" "*"
## 15 ( 1 ) "*" "*" "*"
## 16 ( 1 ) "*" "*" "*"
## 17 ( 1 ) "*" "*" "*"
## 18 ( 1 ) "*" "*" "*"
## 19 ( 1 ) "*" "*" "*"
## Factor_CVNL
## 1 ( 1 ) " "
## 2 ( 1 ) " "
## 3 ( 1 ) " "
## 4 ( 1 ) " "
## 5 ( 1 ) " "
## 6 ( 1 ) " "
## 7 ( 1 ) " "
## 8 ( 1 ) " "
## 9 ( 1 ) " "
## 10 ( 1 ) " "
## 11 ( 1 ) " "
## 12 ( 1 ) " "
## 13 ( 1 ) "*"
## 14 ( 1 ) "*"
## 15 ( 1 ) "*"
## 16 ( 1 ) "*"
## 17 ( 1 ) "*"
## 18 ( 1 ) "*"
## 19 ( 1 ) "*"
names(reg.summary )
## [1] "which" "rsq" "rss" "adjr2" "cp" "bic" "outmat" "obj"
reg.summary$rsq
## [1] 0.4965480 0.6445734 0.6732337 0.6879734 0.6968848 0.7021893 0.7056240
## [8] 0.7077304 0.7094659 0.7108760 0.7112655 0.7115955 0.7118253 0.7119565
## [15] 0.7120628 0.7121488 0.7121865 0.7122222 0.7122384
# gráfica de escala de grises para escoger las mejores variables para la construcción del modelo
plot(regfit.full ,scale ="r2")
plot(regfit.full ,scale ="adjr2")
plot(regfit.full ,scale ="Cp")
plot(regfit.full ,scale ="bic")
A través del método de mejor selección de subconjuntos se analizaron las mejoras en la R2, R2 ajustada, el CP y el BIC que resultan de la inclusión de una variable más en el análisis.
A partir de lo visto en la gráfica de escala de grises, se toma la decisión de incluir 11 variables en el modelo. Los coeficientes de estas 11 variables se obtienen a continuación:
coef(regfit.full,11) #se obtiene los coeficientes de las 11 variables más relevante
## (Intercept) zonaPeriferia
## 2.04349245 -0.30157116
## percepcion.delito.castigado policia.suficienteSí
## -0.13520413 0.21954175
## policia.confiaSí policia.beneficioSí
## 0.33898888 0.14600239
## policia.justoSí fc.evaluacion
## 0.23463549 0.50858032
## gn.evaluacion ejercito.evaluacion
## 0.09317186 -0.20507643
## transito.evaluacion maximo.ingreso.miles
## 0.28490519 0.02257519
Así, se tiene que el modelo obtenido a partir de la mejor selección de subconjuntos es el siguiente:
\(policia.evaluacion=2.04-0.3(zona)-0.14(delito.castigado)+0.22(policia.suficiente)+0.34(policia.confia)+0.15(policia.beneficio)+0.24(policia.justo)+0.51(fc.evaluacion)+0.1(gn.evaluacion)-0.21(ejrcto.evaluacion)+0.29(trnsto.evaluacion)+0.02(max.ingreso.miles)\)
# Se construyen las matrices x & y que servirán para la refresión ridge y lasso
x <- model.matrix(policia.evaluacion~zona + policia.confia + policia.justo + policia.beneficio + policia.suficiente +percepcion.delito.castigado+fc.evaluacion+ejercito.evaluacion+transito.evaluacion+gn.evaluacion+maximo.ingreso.miles,base,weights=base$Factor_CVNL)
x <- model.matrix(policia.evaluacion~zona + policia.confia + policia.justo + policia.beneficio + policia.suficiente +percepcion.delito.castigado+fc.evaluacion+ejercito.evaluacion+transito.evaluacion+gn.evaluacion+maximo.ingreso.miles,base,weights=base$Factor_CVNL)[ ,-1]
y <- base$policia.evaluacion
library(glmnet) # insertar librería glmnet
## Loading required package: Matrix
## Loaded glmnet 4.1-8
grid <- 10^seq(10,-2,length=100) #se define la grid que servirá para probar diferentes valores de lambda
ridge.mod <- glmnet(x,y,alpha=0,lambda=grid) # se realiza la regresión ridge
set.seed(777) #semilla de matrícula y se construyen los segmentos de entrenamiento y de prueba a usarse en las regresiones
train <- sample(1:nrow(x), nrow(x)/2)
test <- (-train)
y.test <- y[test]
set.seed(777) # semilla en matrícula
cv.out <- cv.glmnet(x[train ,],y[train],alpha=0) #cross validation de librería glmnet con argumento family=binomial para que funcione en la regresión logística
plot(cv.out)
bestlam <- cv.out$lambda.min #se calcula el mejor valor de lambda para la regresión ridge
bestlam
## [1] 0.1538951
# a partir del best lambda se van a obtener los valores de los coeficientes para las 11 variables e intercepto
out <- glmnet(x,y,alpha=0)
predict(out,type="coefficients",s=bestlam)[1:12,]
## (Intercept) zonaPeriferia
## 2.02168586 -0.26653312
## policia.confiaSí policia.justoSí
## 0.37494910 0.26664615
## policia.beneficioSí policia.suficienteSí
## 0.13319181 0.21120384
## percepcion.delito.castigado fc.evaluacion
## -0.14105125 0.43349223
## ejercito.evaluacion transito.evaluacion
## -0.12837354 0.27273788
## gn.evaluacion maximo.ingreso.miles
## 0.08759236 0.02542701
Tras realizar una regresión Ridge se estimo una mejor lambda con valor de 0.15. Se observa que sí cambiaron los coeficientes con respecto al modelo estimado a partir de la mejor selección de subconjuntos. Sin embargo, a ninguna variable se le asignó un coeficiente que tendiera a o estuviera muy cerca de 0.
Debido a esto, se puede teorizar que el poco cambio en los coeficientes de las variables en la mejor lambda, con justifican la penalización que esta realiza sobre los predictores. A pesar de esto, más adelante se realizará una comparación de todos los modelos estimados para identificar cuál es el mejor modelo para predecir la evaluación de un policía municipal.
#Lasso
lasso.mod <- glmnet(x[train ,],y[train],alpha=1,lambda=grid) #se realiza la regresión lasso (alpha=1) especificanco familia binomial ya que se cuenta con una regresión logística
plot(lasso.mod)
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
set.seed(777) # semilla en matrícula
cv.out <- cv.glmnet(x[train ,],y[train],alpha=1)
plot(cv.out)
bestlamL <- cv.out$lambda.min # se encuentra la mejor lambda
bestlamL
## [1] 0.002508107
lasso.pred <- predict(lasso.mod,s=bestlamL ,newx=x[test,])
out <- glmnet(x,y,alpha=1,lambda=grid)
lasso.coef <- predict(out,type="coefficients",s=bestlamL)[1:12,] # se obtienen los coeficientes de las 11 variables e intercepto de la regresión lasso a partir del mejor lambda
lasso.coef
## (Intercept) zonaPeriferia
## 1.88729236 -0.25388210
## policia.confiaSí policia.justoSí
## 0.37468925 0.24665702
## policia.beneficioSí policia.suficienteSí
## 0.09056762 0.18797715
## percepcion.delito.castigado fc.evaluacion
## -0.13423755 0.50392644
## ejercito.evaluacion transito.evaluacion
## -0.15145146 0.28347259
## gn.evaluacion maximo.ingreso.miles
## 0.05943157 0.02415346
Similar a la regresión ridge, la regresión lasso también alteró los coeficientes de las variables predictoras del modelo. Para lasso, la lambda óptima toma un valor de 0.0025.
Sin embargo, lasso no ha logrado transformar ninguno de los coeficientes de las variables en 0. Debido a esto, no es posible afirmar que la regresión lasso adapte de mejor forma el modelo.
# Regresión polinómica
model1 <- lm(policia.evaluacion ~ zona + policia.confia + policia.justo + policia.beneficio + policia.suficiente +percepcion.delito.castigado+fc.evaluacion+ejercito.evaluacion+transito.evaluacion+gn.evaluacion+maximo.ingreso.miles,data=base,weights=base$Factor_CVNL)
model2 <- lm(policia.evaluacion ~ zona + policia.confia + policia.justo + policia.beneficio + policia.suficiente + poly(percepcion.delito.castigado,2)+poly(fc.evaluacion,2)+poly(ejercito.evaluacion,2)+poly(transito.evaluacion,2)+poly(gn.evaluacion,2)+poly(maximo.ingreso.miles,2),data=base,weights=base$Factor_CVNL)
transito2<-(base$transito.evaluacion)^2
ingreso2<-(base$maximo.ingreso.miles)^2
modelp<-lm(policia.evaluacion ~ zona + policia.confia + policia.justo + policia.beneficio + policia.suficiente +percepcion.delito.castigado+fc.evaluacion+ejercito.evaluacion+transito2+gn.evaluacion+ingreso2,data=base,weights=base$Factor_CVNL)
summary(model1)
##
## Call:
## lm(formula = policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito.evaluacion +
## gn.evaluacion + maximo.ingreso.miles, data = base, weights = base$Factor_CVNL)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -166.52 -13.99 0.22 15.36 148.44
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.043492 0.137135 14.901 < 2e-16 ***
## zonaPeriferia -0.301571 0.071091 -4.242 2.30e-05 ***
## policia.confiaSí 0.338989 0.078297 4.330 1.55e-05 ***
## policia.justoSí 0.234635 0.081841 2.867 0.004179 **
## policia.beneficioSí 0.146002 0.079270 1.842 0.065616 .
## policia.suficienteSí 0.219542 0.061648 3.561 0.000376 ***
## percepcion.delito.castigado -0.135204 0.023865 -5.665 1.63e-08 ***
## fc.evaluacion 0.508580 0.020499 24.810 < 2e-16 ***
## ejercito.evaluacion -0.205076 0.023264 -8.815 < 2e-16 ***
## transito.evaluacion 0.284905 0.011102 25.662 < 2e-16 ***
## gn.evaluacion 0.093172 0.027794 3.352 0.000813 ***
## maximo.ingreso.miles 0.022575 0.003165 7.133 1.28e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 30.49 on 2515 degrees of freedom
## Multiple R-squared: 0.7113, Adjusted R-squared: 0.71
## F-statistic: 563.2 on 11 and 2515 DF, p-value: < 2.2e-16
summary(model2)
##
## Call:
## lm(formula = policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + poly(percepcion.delito.castigado,
## 2) + poly(fc.evaluacion, 2) + poly(ejercito.evaluacion, 2) +
## poly(transito.evaluacion, 2) + poly(gn.evaluacion, 2) + poly(maximo.ingreso.miles,
## 2), data = base, weights = base$Factor_CVNL)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -156.403 -14.608 0.199 15.838 158.541
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.54847 0.03734 175.381 < 2e-16 ***
## zonaPeriferia -0.28890 0.07061 -4.092 4.42e-05 ***
## policia.confiaSí 0.38298 0.07778 4.924 9.04e-07 ***
## policia.justoSí 0.29513 0.08136 3.627 0.000292 ***
## policia.beneficioSí 0.12852 0.07892 1.628 0.103565
## policia.suficienteSí 0.21209 0.06099 3.477 0.000515 ***
## poly(percepcion.delito.castigado, 2)1 -5.86868 1.13841 -5.155 2.73e-07 ***
## poly(percepcion.delito.castigado, 2)2 1.01700 1.09534 0.928 0.353248
## poly(fc.evaluacion, 2)1 46.59406 1.90722 24.430 < 2e-16 ***
## poly(fc.evaluacion, 2)2 -0.29232 1.51632 -0.193 0.847143
## poly(ejercito.evaluacion, 2)1 -18.97324 2.31256 -8.204 3.66e-16 ***
## poly(ejercito.evaluacion, 2)2 -0.44755 1.74675 -0.256 0.797803
## poly(transito.evaluacion, 2)1 33.97428 1.38665 24.501 < 2e-16 ***
## poly(transito.evaluacion, 2)2 -6.77249 1.18599 -5.710 1.26e-08 ***
## poly(gn.evaluacion, 2)1 9.38805 2.59866 3.613 0.000309 ***
## poly(gn.evaluacion, 2)2 -3.10451 1.89784 -1.636 0.102005
## poly(maximo.ingreso.miles, 2)1 7.64044 1.10414 6.920 5.72e-12 ***
## poly(maximo.ingreso.miles, 2)2 -2.39769 1.08002 -2.220 0.026505 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 30.13 on 2509 degrees of freedom
## Multiple R-squared: 0.7186, Adjusted R-squared: 0.7167
## F-statistic: 376.9 on 17 and 2509 DF, p-value: < 2.2e-16
summary(modelp)
##
## Call:
## lm(formula = policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito2 + gn.evaluacion +
## ingreso2, data = base, weights = base$Factor_CVNL)
##
## Weighted Residuals:
## Min 1Q Median 3Q Max
## -182.798 -14.451 0.114 16.252 138.016
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.012e+00 1.300e-01 23.165 < 2e-16 ***
## zonaPeriferia -3.068e-01 7.320e-02 -4.192 2.87e-05 ***
## policia.confiaSí 3.350e-01 8.086e-02 4.143 3.54e-05 ***
## policia.justoSí 1.922e-01 8.450e-02 2.275 0.023007 *
## policia.beneficioSí 2.094e-01 8.164e-02 2.565 0.010390 *
## policia.suficienteSí 2.548e-01 6.357e-02 4.008 6.31e-05 ***
## percepcion.delito.castigado -1.518e-01 2.464e-02 -6.161 8.38e-10 ***
## fc.evaluacion 5.369e-01 2.105e-02 25.500 < 2e-16 ***
## ejercito.evaluacion -2.370e-01 2.388e-02 -9.925 < 2e-16 ***
## transito2 2.273e-02 1.038e-03 21.893 < 2e-16 ***
## gn.evaluacion 1.030e-01 2.867e-02 3.591 0.000336 ***
## ingreso2 4.813e-04 8.039e-05 5.987 2.44e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 31.47 on 2515 degrees of freedom
## Multiple R-squared: 0.6924, Adjusted R-squared: 0.6911
## F-statistic: 514.7 on 11 and 2515 DF, p-value: < 2.2e-16
Se realizó un modelo completamente lineal, seguido de un modelo que incorporaba polinomios cuadráticos a todas aquellas que no fueran variables dummy. A partir de la comparación de ambos modelos, se decidió incorporar la calificación de tránsito y el máximo ingreso en miles como variables cuadráticas ya que fueron las únicas que resultaron significativas al ser elevadas al cuadrado en el segundo modelo.
anova(model1,model2,modelp)
## Analysis of Variance Table
##
## Model 1: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito.evaluacion +
## gn.evaluacion + maximo.ingreso.miles
## Model 2: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + poly(percepcion.delito.castigado,
## 2) + poly(fc.evaluacion, 2) + poly(ejercito.evaluacion, 2) +
## poly(transito.evaluacion, 2) + poly(gn.evaluacion, 2) + poly(maximo.ingreso.miles,
## 2)
## Model 3: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito2 + gn.evaluacion +
## ingreso2
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 2515 2337839
## 2 2509 2278247 6 59593 10.938 4.659e-12 ***
## 3 2515 2490474 -6 -212228 38.954 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
La prueba Anova: Tras realizar una prueba de ANOVA (análisis de varianzas) que determina si hay diferencia significativa en los modelos.Se pued eobservar que el menor RSS es el del modelo 2 por lo tanto este es el mejor.
library(splines)
model.sp <- lm(policia.evaluacion ~ zona + policia.confia + policia.justo + policia.beneficio+ policia.suficiente+
bs(percepcion.delito.castigado, df = 3) +
bs(fc.evaluacion, df = 5) +
bs(ejercito.evaluacion, df = 5) +
bs(transito.evaluacion, df = 5)+
bs(gn.evaluacion,df=5)+
bs(maximo.ingreso.miles,df=5),
data = base,weights = base$Factor_CVNL)
# Grids de valores para las variables numéricas
percepcion.delito.castigado.grid <- seq(min(base$percepcion.delito.castigado), max(base$percepcion.delito.castigado))
fc.evaluacion.grid <- seq(min(base$fc.evaluacion), max(base$fc.evaluacion))
ejercito.evaluacion.grid <- seq(min(base$ejercito.evaluacion), max(base$ejercito.evaluacion))
transito.evaluacion.grid <- seq(min(base$transito.evaluacion), max(base$transito.evaluacion))
gn.evaluacion.grid<-seq(min(base$gn.evaluacion),max(base$gn.evaluacion))
maximo.ingreso.miles.grid<-seq(min(base$maximo.ingreso.miles),max(base$maximo.ingreso.miles))
newdata <- expand.grid(zona = unique(base$zona), # pq variables dummy
policia.confia = unique(base$policia.confia),
policia.justo = unique(base$policia.justo),
percepcion.delito.castigado = percepcion.delito.castigado.grid,
fc.evaluacion = fc.evaluacion.grid,
ejercito.evaluacion = ejercito.evaluacion.grid,
transito.evaluacion = transito.evaluacion.grid,
gn.evaluacion=gn.evaluacion.grid,
maximo.ingreso.miles=maximo.ingreso.miles.grid,
policia.beneficio=unique(base$policia.beneficio),
policia.suficiente=unique(base$policia.suficiente))
#predicciones <- predict(model.sp, newdata = newdata, se.fit = TRUE)
#Predicciones <- predicciones$fit
#SE <- predicciones$se.fit
# predicciones como comentario porque debido al tamaño no deja convetir al archivo en HTML
anova(model1,model2,modelp,model.sp)
## Analysis of Variance Table
##
## Model 1: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito.evaluacion +
## gn.evaluacion + maximo.ingreso.miles
## Model 2: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + poly(percepcion.delito.castigado,
## 2) + poly(fc.evaluacion, 2) + poly(ejercito.evaluacion, 2) +
## poly(transito.evaluacion, 2) + poly(gn.evaluacion, 2) + poly(maximo.ingreso.miles,
## 2)
## Model 3: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito2 + gn.evaluacion +
## ingreso2
## Model 4: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + bs(percepcion.delito.castigado,
## df = 3) + bs(fc.evaluacion, df = 5) + bs(ejercito.evaluacion,
## df = 5) + bs(transito.evaluacion, df = 5) + bs(gn.evaluacion,
## df = 5) + bs(maximo.ingreso.miles, df = 5)
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 2515 2337839
## 2 2509 2278247 6 59593 11.283 1.808e-12 ***
## 3 2515 2490474 -6 -212228 40.182 < 2.2e-16 ***
## 4 2493 2194551 22 295923 15.280 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
En esta segunda prueba ANOVA se incluyó el modelo splines y según los resultados de la prueba se ve que sí hjay diferencia significativa entre los cuatro modelos. Ahora, comparando los RSS demuestra que el mejor modelo es el cuarto, es decir el modelo spline.
library(gam)
## Loading required package: foreach
## Loaded gam 1.22-3
gam<-gam(policia.evaluacion~policia.confia + policia.justo+ ns(percepcion.delito.castigado,2)+ns(fc.evaluacion,3)+ns(ejercito.evaluacion,3)+ns(transito.evaluacion,3)+ns(gn.evaluacion,3)+ns(maximo.ingreso.miles,6)+policia.beneficio+policia.suficiente, data=base,weights=base$Factor_CVNL)
plot.Gam(gam, se=TRUE, col ="purple")
summary(gam)
##
## Call: gam(formula = policia.evaluacion ~ policia.confia + policia.justo +
## ns(percepcion.delito.castigado, 2) + ns(fc.evaluacion, 3) +
## ns(ejercito.evaluacion, 3) + ns(transito.evaluacion, 3) +
## ns(gn.evaluacion, 3) + ns(maximo.ingreso.miles, 6) + policia.beneficio +
## policia.suficiente, data = base, weights = base$Factor_CVNL)
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -155.85475 -13.75417 0.09271 15.30672 146.65357
##
## (Dispersion Parameter for gaussian family taken to be 898.9461)
##
## Null Deviance: 8096846 on 2526 degrees of freedom
## Residual Deviance: 2249163 on 2502 degrees of freedom
## AIC: 8075.234
##
## Number of Local Scoring Iterations: 2
##
## Anova for Parametric Effects
## Df Sum Sq Mean Sq F value Pr(>F)
## policia.confia 1 2728494 2728494 3035.2140 < 2.2e-16 ***
## policia.justo 1 190599 190599 212.0251 < 2.2e-16 ***
## ns(percepcion.delito.castigado, 2) 2 143553 71776 79.8450 < 2.2e-16 ***
## ns(fc.evaluacion, 3) 3 1727657 575886 640.6231 < 2.2e-16 ***
## ns(ejercito.evaluacion, 3) 3 262482 87494 97.3295 < 2.2e-16 ***
## ns(transito.evaluacion, 3) 3 692351 230784 256.7269 < 2.2e-16 ***
## ns(gn.evaluacion, 3) 3 13918 4639 5.1609 0.0014781 **
## ns(maximo.ingreso.miles, 6) 6 71806 11968 13.3131 6.652e-15 ***
## policia.beneficio 1 5635 5635 6.2684 0.0123542 *
## policia.suficiente 1 11188 11188 12.4461 0.0004264 ***
## Residuals 2502 2249163 899
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
gam2<-gam(policia.evaluacion~policia.confia + policia.justo+ percepcion.delito.castigado+fc.evaluacion+ejercito.evaluacion+transito.evaluacion+gn.evaluacion+ns(maximo.ingreso.miles,3)+policia.beneficio+policia.suficiente, data=base,weights=base$Factor_CVNL)
plot.Gam(gam2, se=TRUE, col ="purple")
summary(gam2)
##
## Call: gam(formula = policia.evaluacion ~ policia.confia + policia.justo +
## percepcion.delito.castigado + fc.evaluacion + ejercito.evaluacion +
## transito.evaluacion + gn.evaluacion + ns(maximo.ingreso.miles,
## 3) + policia.beneficio + policia.suficiente, data = base,
## weights = base$Factor_CVNL)
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -168.5217 -14.1577 0.5771 15.1015 141.4341
##
## (Dispersion Parameter for gaussian family taken to be 930.5123)
##
## Null Deviance: 8096846 on 2526 degrees of freedom
## Residual Deviance: 2339308 on 2514 degrees of freedom
## AIC: 8150.537
##
## Number of Local Scoring Iterations: 2
##
## Anova for Parametric Effects
## Df Sum Sq Mean Sq F value Pr(>F)
## policia.confia 1 2728494 2728494 2932.2488 < 2.2e-16 ***
## policia.justo 1 190599 190599 204.8324 < 2.2e-16 ***
## percepcion.delito.castigado 1 128897 128897 138.5222 < 2.2e-16 ***
## fc.evaluacion 1 1719033 1719033 1847.4046 < 2.2e-16 ***
## ejercito.evaluacion 1 176169 176169 189.3251 < 2.2e-16 ***
## transito.evaluacion 1 710775 710775 763.8529 < 2.2e-16 ***
## gn.evaluacion 1 9477 9477 10.1851 0.0014333 **
## ns(maximo.ingreso.miles, 3) 3 77628 25876 27.8084 < 2.2e-16 ***
## policia.beneficio 1 4893 4893 5.2579 0.0219290 *
## policia.suficiente 1 11574 11574 12.4384 0.0004281 ***
## Residuals 2514 2339308 931
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#Hacer comparación GAM de las diferentes interacciones polinómicas de las variables.
Los resultados de las regresiones realizadas con los modelo aditivos generalizados señalan que el modelo se beneficia de relaciones no lineales para todas la variables excepto sobrepeso. Al tener un valor p de 0.05, la variable sobrepeso está correctamente ajustada con una relación lineal.
library(tree)
tree1<-tree(policia.evaluacion~zona + policia.beneficio+policia.suficiente+policia.confia + policia.justo + percepcion.delito.castigado + fc.evaluacion + gn.evaluacion+ejercito.evaluacion + transito.evaluacion+maximo.ingreso.miles, data=base,weights=base$Factor_CVNL) #árbol ajustado a conjunto de entrenamiento
summary(tree1)
##
## Regression tree:
## tree(formula = policia.evaluacion ~ zona + policia.beneficio +
## policia.suficiente + policia.confia + policia.justo + percepcion.delito.castigado +
## fc.evaluacion + gn.evaluacion + ejercito.evaluacion + transito.evaluacion +
## maximo.ingreso.miles, data = base, weights = base$Factor_CVNL)
## Variables actually used in tree construction:
## [1] "fc.evaluacion" "maximo.ingreso.miles"
## [3] "percepcion.delito.castigado" "transito.evaluacion"
## [5] "gn.evaluacion" "ejercito.evaluacion"
## Number of terminal nodes: 64
## Residual mean deviance: 1.424 = 2972000 / 2088000
## Distribution of residuals:
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -7.55000 -0.54410 0.07982 0.00649 0.75260 4.82200
#trazo del árbol
plot(tree1)
text(tree1, pretty=0,cex=0.3)
#uso de validación cruzada
cv.base<-cv.tree(tree1)
plot(cv.base$size, cv.base$dev, type="b")
#predicción
tree1
## node), split, n, deviance, yval
## * denotes terminal node
##
## 1) root 2087596.48 8097000 6.970
## 2) fc.evaluacion < 8.5 1509856.96 5227000 6.341
## 4) maximo.ingreso.miles < 21.7815 1222888.1 4273000 6.119
## 8) percepcion.delito.castigado < 3.5 976631.01 3134000 6.336
## 16) transito.evaluacion < 8.5 951700.06 3048000 6.298
## 32) transito.evaluacion < 7.5 813817.8 2628000 6.089
## 64) gn.evaluacion < 9.5 767742.18 2484000 6.046
## 128) fc.evaluacion < 5.5 159026.02 575700 4.229
## 256) fc.evaluacion < 4.5 79045.18 371500 3.585
## 512) ejercito.evaluacion < 3.5 34950.3 134300 2.881
## 1024) fc.evaluacion < 2.5 16729.29 56500 2.178 *
## 1025) fc.evaluacion > 2.5 18221.01 61980 3.526 *
## 513) ejercito.evaluacion > 3.5 44094.88 206100 4.142
## 1026) percepcion.delito.castigado < 2.5 24590.63 124500 4.751
## 2052) gn.evaluacion < 4.5 15351.61 66510 4.537 *
## 2053) gn.evaluacion > 4.5 9239.02 56150 5.107 *
## 1027) percepcion.delito.castigado > 2.5 19504.25 61010 3.375 *
## 257) fc.evaluacion > 4.5 79980.84 139000 4.865
## 514) ejercito.evaluacion < 7.5 50909.04 96520 5.016
## 1028) ejercito.evaluacion < 5.5 36450.29 65550 5.053 *
## 1029) ejercito.evaluacion > 5.5 14458.75 30790 4.920 *
## 515) ejercito.evaluacion > 7.5 29071.8 39290 4.601 *
## 129) fc.evaluacion > 5.5 608716.17 1246000 6.521
## 258) fc.evaluacion < 7.5 397269.16 737800 6.228
## 516) transito.evaluacion < 4.5 152055.08 225900 5.325
## 1032) gn.evaluacion < 7.5 74625.7 170600 5.459
## 2064) transito.evaluacion < 3.5 49574.15 91590 4.919
## 4128) transito.evaluacion < 2.5 37619.86 65500 4.758 *
## 4129) transito.evaluacion > 2.5 11954.28 22050 5.425 *
## 2065) transito.evaluacion > 3.5 25051.55 35930 6.529 *
## 1033) gn.evaluacion > 7.5 77429.38 52650 5.195 *
## 517) transito.evaluacion > 4.5 245214.08 311100 6.787
## 1034) maximo.ingreso.miles < 15.558 125814.96 162800 6.562
## 2068) gn.evaluacion < 7.5 87109.98 110500 6.561
## 4136) transito.evaluacion < 6.5 40982 51030 6.221 *
## 4137) transito.evaluacion > 6.5 46127.98 50460 6.864 *
## 2069) gn.evaluacion > 7.5 38704.98 52300 6.564 *
## 1035) maximo.ingreso.miles > 15.558 119399.12 135200 7.024
## 2070) percepcion.delito.castigado < 2.5 79325.96 97630 7.144
## 4140) gn.evaluacion < 6.5 34400.62 35890 7.047 *
## 4141) gn.evaluacion > 6.5 44925.35 61170 7.218 *
## 2071) percepcion.delito.castigado > 2.5 40073.16 34210 6.788 *
## 259) fc.evaluacion > 7.5 211447.01 409600 7.073
## 518) maximo.ingreso.miles < 15.558 127071.9 200000 6.994
## 1036) transito.evaluacion < 5.5 65617.97 137700 6.494
## 2072) percepcion.delito.castigado < 2.5 39196.28 65710 6.737 *
## 2073) percepcion.delito.castigado > 2.5 26421.69 66220 6.133 *
## 1037) transito.evaluacion > 5.5 61453.93 28360 7.528 *
## 519) maximo.ingreso.miles > 15.558 84375.1 207600 7.192
## 1038) transito.evaluacion < 5.5 44663.77 158700 6.693
## 2076) transito.evaluacion < 4.5 23585.33 125500 6.122
## 4152) transito.evaluacion < 2.5 9609.86 54270 4.435 *
## 4153) transito.evaluacion > 2.5 13975.48 25130 7.282 *
## 2077) transito.evaluacion > 4.5 21078.43 16830 7.332 *
## 1039) transito.evaluacion > 5.5 39711.34 25320 7.753 *
## 65) gn.evaluacion > 9.5 46075.61 118600 6.799
## 130) transito.evaluacion < 5.5 25945.76 65480 6.317 *
## 131) transito.evaluacion > 5.5 20129.85 39260 7.421 *
## 33) transito.evaluacion > 7.5 137882.26 174400 7.533
## 66) percepcion.delito.castigado < 2.5 83921.78 105600 7.657
## 132) maximo.ingreso.miles < 15.558 49004.35 56170 7.544 *
## 133) maximo.ingreso.miles > 15.558 34917.43 47960 7.816 *
## 67) percepcion.delito.castigado > 2.5 53960.49 65410 7.339 *
## 17) transito.evaluacion > 8.5 24930.95 32470 7.788 *
## 9) percepcion.delito.castigado > 3.5 246257.1 910300 5.258
## 18) transito.evaluacion < 5.5 163503.91 482400 4.438
## 36) gn.evaluacion < 7.5 87087.41 276600 4.098
## 72) ejercito.evaluacion < 6.5 37085.41 140900 4.018
## 144) transito.evaluacion < 4.5 20469.55 91750 3.234
## 288) gn.evaluacion < 4.5 12034.32 48770 2.832 *
## 289) gn.evaluacion > 4.5 8435.23 38250 3.809 *
## 145) transito.evaluacion > 4.5 16615.86 21050 4.984 *
## 73) ejercito.evaluacion > 6.5 50002 135300 4.157
## 146) ejercito.evaluacion < 8.5 30364.39 89600 4.457
## 292) fc.evaluacion < 6.5 15379.76 44350 3.479 *
## 293) fc.evaluacion > 6.5 14984.64 15460 5.460 *
## 147) ejercito.evaluacion > 8.5 19637.61 38800 3.695 *
## 37) gn.evaluacion > 7.5 76416.5 184300 4.825
## 74) fc.evaluacion < 7.5 44611.46 78550 4.247 *
## 75) fc.evaluacion > 7.5 31805.03 69990 5.635 *
## 19) transito.evaluacion > 5.5 82753.19 100300 6.879
## 38) gn.evaluacion < 7.5 41763.48 41420 6.548 *
## 39) gn.evaluacion > 7.5 40989.71 49670 7.217 *
## 5) maximo.ingreso.miles > 21.7815 286968.85 636100 7.288
## 10) gn.evaluacion < 9.5 280311.22 620400 7.302
## 20) fc.evaluacion < 7.5 180220.83 451200 7.038
## 40) ejercito.evaluacion < 7.5 125951.39 232500 7.274
## 80) transito.evaluacion < 6.5 80254.02 208200 7.243
## 160) fc.evaluacion < 6.5 51273.81 144000 7.012
## 320) percepcion.delito.castigado < 2.5 37654.91 112700 7.223
## 640) transito.evaluacion < 5.5 25957.13 87860 6.805
## 1280) gn.evaluacion < 5.5 15694.74 57160 6.362 *
## 1281) gn.evaluacion > 5.5 10262.4 22920 7.482 *
## 641) transito.evaluacion > 5.5 11697.78 10210 8.151 *
## 321) percepcion.delito.castigado > 2.5 13618.9 25040 6.429 *
## 161) fc.evaluacion > 6.5 28980.21 56620 7.652 *
## 81) transito.evaluacion > 6.5 45697.37 24110 7.329 *
## 41) ejercito.evaluacion > 7.5 54269.44 195400 6.490
## 82) gn.evaluacion < 8.5 41039.04 117900 6.723
## 164) percepcion.delito.castigado < 2.5 21749.54 42250 7.282 *
## 165) percepcion.delito.castigado > 2.5 19289.5 61200 6.092 *
## 83) gn.evaluacion > 8.5 13230.4 68330 5.767 *
## 21) fc.evaluacion > 7.5 100090.38 134000 7.776
## 42) percepcion.delito.castigado < 2.5 61145.42 65020 7.838 *
## 43) percepcion.delito.castigado > 2.5 38944.96 68390 7.679 *
## 11) gn.evaluacion > 9.5 6657.64 13490 6.719 *
## 3) fc.evaluacion > 8.5 577739.53 711400 8.614
## 6) transito.evaluacion < 6.5 153799.34 293300 7.859
## 12) transito.evaluacion < 5.5 80161.26 194000 7.631
## 24) gn.evaluacion < 9.5 44756.31 87330 7.704
## 48) transito.evaluacion < 4.5 18946.84 47390 7.266 *
## 49) transito.evaluacion > 4.5 25809.47 33630 8.026 *
## 25) gn.evaluacion > 9.5 35404.95 106100 7.539
## 50) percepcion.delito.castigado < 2.5 23745.77 28140 7.935 *
## 51) percepcion.delito.castigado > 2.5 11659.18 66730 6.734 *
## 13) transito.evaluacion > 5.5 73638.08 90590 8.107
## 26) gn.evaluacion < 9.5 50966.72 57810 7.934 *
## 27) gn.evaluacion > 9.5 22671.36 27830 8.496 *
## 7) transito.evaluacion > 6.5 423940.19 298600 8.888
## 14) maximo.ingreso.miles < 15.558 159110.08 128200 8.774
## 28) gn.evaluacion < 9.5 83501.1 56390 8.550 *
## 29) gn.evaluacion > 9.5 75608.98 62970 9.021 *
## 15) maximo.ingreso.miles > 15.558 264830.11 167100 8.957
## 30) ejercito.evaluacion < 9.5 155350.15 85440 8.786
## 60) transito.evaluacion < 8.5 67343.77 41970 8.450 *
## 61) transito.evaluacion > 8.5 88006.38 30060 9.043 *
## 31) ejercito.evaluacion > 9.5 109479.96 70680 9.200 *
train <- sample(1:nrow(base), nrow(base)/2) #conjunto de entrenamiento
test<-base[-train, ]
base.pred<-predict(tree1, newdata=test)
base.test1<-base[-train,"policia.evaluacion"]
sqrt(mean((base.pred-base.test1)^2)) #error
## [1] 1.315022
Como se puede observar, el árbol de rgresión contiene 64 nodos terminales que emanan de 6 variables predictoras. Sin embargo, como se observa en la gráfica de validación cruzada, el modelo minimiza el error con tan solo 13 nodos terminales. Debido a esto, más adelante se realizará un proceso de podado del árbol para conservar únicamente los nodos necesarios y no tener una sobresaturación en el árbol.
library(randomForest)
## randomForest 4.7-1.1
## Type rfNews() to see new features/changes/bug fixes.
##
## Attaching package: 'randomForest'
## The following object is masked from 'package:dplyr':
##
## combine
set.seed(777)
attach(base)
## The following object is masked from eav:
##
## Factor_CVNL
train <- sample(1:nrow(base), nrow(base)/2) #conjunto de entrenamiento
test<-base[-train, ]
ev.test<-base[-train,"policia.evaluacion"]
#proceso bagging
bag.model<-randomForest(policia.evaluacion~ zona + policia.beneficio+policia.suficiente+policia.confia + policia.justo + percepcion.delito.castigado + fc.evaluacion + gn.evaluacion+ejercito.evaluacion + transito.evaluacion+maximo.ingreso.miles, data=base,subset=train,mtry=11,importance=TRUE,ntree=80)
bag.model
##
## Call:
## randomForest(formula = policia.evaluacion ~ zona + policia.beneficio + policia.suficiente + policia.confia + policia.justo + percepcion.delito.castigado + fc.evaluacion + gn.evaluacion + ejercito.evaluacion + transito.evaluacion + maximo.ingreso.miles, data = base, mtry = 11, importance = TRUE, ntree = 80, subset = train)
## Type of random forest: regression
## Number of trees: 80
## No. of variables tried at each split: 11
##
## Mean of squared residuals: 1.296351
## % Var explained: 69.82
plot(bag.model)
bag.pred <- predict(bag.model, newdata=test)
mean((bag.pred-ev.test)^2) #error del modelo
## [1] 1.324765
set.seed(777)
bag.model.rf<-randomForest(policia.evaluacion~ zona + policia.beneficio+policia.suficiente+policia.confia + policia.justo + percepcion.delito.castigado + fc.evaluacion + gn.evaluacion+ejercito.evaluacion + transito.evaluacion+maximo.ingreso.miles, data=base,subset=train,mtry=11,ntree=25)
bag.pred.rf <- predict(bag.model.rf, newdata=test)
mean((bag.pred.rf-ev.test)^2)
## [1] 1.322576
set.seed(777)
rf.ev<-randomForest(policia.evaluacion~ zona + policia.beneficio+policia.suficiente+policia.confia + policia.justo + percepcion.delito.castigado + fc.evaluacion + gn.evaluacion+ejercito.evaluacion + transito.evaluacion+maximo.ingreso.miles,data=base,subset=train,mtry=11,importance =TRUE)
ev.pred.rf <- predict(rf.ev ,newdata=test)
mean((ev.pred.rf-ev.test)^2)
## [1] 1.300508
importance(rf.ev)
## %IncMSE IncNodePurity
## zona 2.366636 55.53813
## policia.beneficio 20.625978 72.39982
## policia.suficiente 20.003454 58.08238
## policia.confia 28.781068 390.02983
## policia.justo 23.751242 201.60180
## percepcion.delito.castigado 24.927564 184.88492
## fc.evaluacion 68.756753 1525.00545
## gn.evaluacion 27.712954 224.07692
## ejercito.evaluacion 31.323072 206.80114
## transito.evaluacion 74.579878 2018.88590
## maximo.ingreso.miles 21.556725 230.26339
varImpPlot(rf.ev)
Se realizó un proceso de bagging y random forests donde se encontró que a partir de los 40 árboles cultivados se minimiza el error medio cuadrático del modelo. Sin embargo, de los diferentes modelos realizados con esta metodología ninguno de ellos alcanzó un error menor al obtenido por el árbol de regresión original.
Adicionalmente, se llevó a cabo un análisis de importancia en la varianza del modelo de random forest. Este arrojó que las variables más relevantes para explicar la varianza en la calificación otorgada a la policía municipal son la evaluación al equipo de tránsito, la evaluación a la fuerza civil, y, en menor medida, la confianza en la policía. Estos resultados coinciden con las variables de las que derivaban los nodos en el árbol de regresión.
library(gbm)
## Loaded gbm 2.1.9
## This version of gbm is no longer under development. Consider transitioning to gbm3, https://github.com/gbm-developers/gbm3
set.seed(777)
boost.base<-gbm(policia.evaluacion~ zona + policia.beneficio+policia.suficiente+policia.confia + policia.justo + percepcion.delito.castigado + fc.evaluacion + gn.evaluacion+ejercito.evaluacion + transito.evaluacion+maximo.ingreso.miles,data=base[train,],distribution="gaussian",n.trees=5000,interaction.depth=8)
summary(boost.base)
## var rel.inf
## transito.evaluacion transito.evaluacion 24.559484
## fc.evaluacion fc.evaluacion 20.938639
## gn.evaluacion gn.evaluacion 11.310356
## ejercito.evaluacion ejercito.evaluacion 10.136452
## percepcion.delito.castigado percepcion.delito.castigado 8.755268
## maximo.ingreso.miles maximo.ingreso.miles 8.183007
## policia.confia policia.confia 5.527859
## zona zona 2.960420
## policia.suficiente policia.suficiente 2.883544
## policia.beneficio policia.beneficio 2.481071
## policia.justo policia.justo 2.263900
par(mfrow=c(1,2))
plot(boost.base, i="fc.evaluacion")
plot(boost.base, i="transito.evaluacion")
Otra manera para observar las variables con mayor efecto en la predicción es usar una función de boosting. Esta, al igual que random forest, señala que las dos variables más importantes en la calificación del policía municipal son las evaluaciones al tránsito y a la fuerza civil.
Adicionalmente, permite observar el efecto marginal que los diferentes cambios en estas dos variables tienen sobre el resultado final. Es natural que, a mayor calificación que un individuo le de al oficial de tránsito y/o el de fuerza civil, tendrá una mejor percepción de su oficial municipal.
prune.tree<-prune.tree(tree1, best=13)
plot(prune.tree)
text(prune.tree, pretty=0,cex=0.5)
tree.prune.pred<-predict(prune.tree, test)
table(tree.prune.pred, base.test1)
## base.test1
## tree.prune.pred 1 2 3 4 5 6 7 8 9 10
## 3.58453731539177 1 1 3 5 9 5 15 9 10 6
## 4.43780685751249 3 0 3 5 20 11 13 16 16 2
## 4.86496615964418 2 0 0 4 10 4 11 10 7 2
## 5.32473988267226 0 0 0 8 9 7 24 18 16 4
## 6.78737277368505 4 6 2 3 22 10 32 32 23 11
## 6.79948472090269 0 0 3 0 1 1 5 6 3 4
## 6.87943820557236 3 2 2 3 8 10 12 15 3 3
## 7.07304385158353 6 2 1 3 10 11 36 28 14 5
## 7.28806886213783 3 5 2 6 20 14 20 38 29 6
## 7.53265224433698 2 4 2 6 8 9 10 26 8 2
## 7.78839619841333 1 0 0 0 5 2 5 4 2 0
## 7.85899277453261 8 3 3 2 15 7 12 24 20 8
## 8.88823039868771 9 7 4 11 30 24 41 86 50 27
Como se puede ver, después de haber llevado a cabo un proceso de podado, únicamente se cuentan con 13 nodos terminales. Este número se escogió porque como se vió anteriormente, es la cantidad que minimiza el error. El mejor árbol utiliza los predictores de evaluación a las fuerzas civiles, de tránsito y guardia nacional, la percepción del delito castigado, y el máximo ingreso del individuo.
#PRUEBA ANOVA
anova(model1,model2,modelp)
## Analysis of Variance Table
##
## Model 1: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito.evaluacion +
## gn.evaluacion + maximo.ingreso.miles
## Model 2: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + poly(percepcion.delito.castigado,
## 2) + poly(fc.evaluacion, 2) + poly(ejercito.evaluacion, 2) +
## poly(transito.evaluacion, 2) + poly(gn.evaluacion, 2) + poly(maximo.ingreso.miles,
## 2)
## Model 3: policia.evaluacion ~ zona + policia.confia + policia.justo +
## policia.beneficio + policia.suficiente + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito2 + gn.evaluacion +
## ingreso2
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 2515 2337839
## 2 2509 2278247 6 59593 10.938 4.659e-12 ***
## 3 2515 2490474 -6 -212228 38.954 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(gam,gam2)
## Analysis of Deviance Table
##
## Model 1: policia.evaluacion ~ policia.confia + policia.justo + ns(percepcion.delito.castigado,
## 2) + ns(fc.evaluacion, 3) + ns(ejercito.evaluacion, 3) +
## ns(transito.evaluacion, 3) + ns(gn.evaluacion, 3) + ns(maximo.ingreso.miles,
## 6) + policia.beneficio + policia.suficiente
## Model 2: policia.evaluacion ~ policia.confia + policia.justo + percepcion.delito.castigado +
## fc.evaluacion + ejercito.evaluacion + transito.evaluacion +
## gn.evaluacion + ns(maximo.ingreso.miles, 3) + policia.beneficio +
## policia.suficiente
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 2502 2249163
## 2 2514 2339308 -12 -90145 4.911e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Como se observa en los procesos de ANOVA realizados, hay evidencia estadística para decir que el model2 es diferente del model1; y que a su vez, modelp es diferente de model2.
Adicionalmente, gam2 también es estadísticamente diferente de gam.
Sin embargo, como se observa en los summary de cada uno de los modelos, el agregar interacciones polinómicas no mejora en una cantidad significativa su error ni ajuste a los datos. Adicionalmente, los splines naturales de los modelos gam son complejos pero no lograr dar un modelo muchísimo mejor.
Debido a esto, se ha optado por usar el modelo lineal model1 para predecir la calificación que un individuo dará a su policía local. de forma que el modelo que se obtiene es: \(policia.evaluacion=2.04-0.3*zona-0.1315*delito.castigado+0.339*policia.confia +0.235*policia.justa+0.146*policia.beneficio+0.51*fc.evaluacion-0.21*ejercito.evaluacion+0.29*transito.evaluacion+0.09*gn.evaluacion+0.02*ing.maximo.miles\)
Para llegar al modelo final que explica la percepción ciudadana sobre los cuerpos policiales en Nuevo León, seguimos un proceso meticuloso que incluyó varias etapas. Primero, seleccionamos las variables iniciales relevantes basadas en estudios previos y su potencial influencia en la percepción de la policía y estas variables incluyeron el municipio, sexo, edad, actividad laboral, percepción de seguridad, y experiencias con delitos y cuerpos policiales. Luego, realizamos una limpieza exhaustiva de los datos, lo que implicó filtrar datos irrelevantes, manejar valores faltantes, y crear un nuevo data frame con las variables seleccionadas.
Después, aplicamos métodos de análisis estadístico, específicamente regresión lineal y análisis de correlación, para identificar la relación entre las variables independientes (como la percepción de seguridad y experiencias con delitos) y la variable dependiente (evaluación de la policía), estos métodos fueron seleccionados porque permiten cuantificar la influencia de cada variable en la percepción de la policía.
Durante la validación del modelo, empleamos pruebas ANOVA para comparar diferentes modelos y seleccionar el más adecuado. La prueba ANOVA nos permitió evaluar la significancia estadística de las diferencias entre los modelos y determinar cuál explicaba mejor la variabilidad en la percepción de la policía. Nuestro modelo final mostró una significancia estadística superior y un mayor poder explicativo en comparación con los otros modelos probados, confirmando su idoneidad.
Al realizar la comparación entre los modelos, se llegó a un modelo lineal final. Este modelo utiliza las variables de zona, percepciones de la policía, calificaciones a otros cuerpos policiacos, y el ingreso máximo del individuo, para predecir la calificación que éste dará a su policía local.
Se decidió conservar el modelo lineal debido a que es fácil de explicar, al mismo tiempo que el resto de los modelos no llevaba a mejoras tan significativas que justificaran su complicación en el análisis.
El modelo final fue elegido basado en su capacidad para explicar la variabilidad en la percepción de la policía y la significancia estadística de los coeficientes, validamos este modelo utilizando técnicas de validación cruzada para asegurar que los resultados no fueran producto del overfitting. Los hallazgos clave del análisis indicaron que la percepción de la policía varía significativamente entre municipios, siendo más alta en San Pedro y más baja en la periferia. Además, factores sociodemográficos como sexo, edad, y actividad laboral tienen un impacto notable en la evaluación de la policía, al igual que las experiencias personales relacionadas con la seguridad y los encuentros con cuerpos policiales.
Basándonos en las problemáticas que enfrentamos a la hora de usar la encuesta y en los resultados obtenidos, hemos formulado varias recomendaciones. Primero, consideramos que es esencial fortalecer la encuesta incluyendo preguntas sobre raza y auto adscripción indígena, pues esta inclusión permitirá identificar discriminaciones específicas y mejorar la precisión del análisis, proporcionando una visión más completa de la percepción ciudadana. Proponemos añadir preguntas directas sobre estos temas en futuras ediciones de la encuesta.
Además, pensamos que pudiera ser muy valioso detallar más las experiencias personales con la policía en la encuesta, pues actualmente, la falta de detalles puede limitar la riqueza y precisión de la información obtenida. Sugerimos incluir preguntas adicionales sobre el tipo de interacciones con la policía, la frecuencia de estas interacciones, y la percepción de profesionalismo durante las mismas.
Identificamos algunas problemáticas durante el análisis que deben ser abordadas. La primera es la posible presencia de sesgos en las respuestas debido a la percepción socialmente deseable y para combatir esto, proponemos implementar técnicas de validación cruzada y métodos de detección de sesgos que aseguren la autenticidad de las respuestas. La segunda problemática es la ambigüedad en algunas preguntas, lo cual puede llevar a respuestas inconsistentes, por lo mismo recomendamos revisar y mejorar la redacción de las preguntas para asegurar que sean claras y específicas.