Modelo Geometrico - Evacuaciones públicas
Grupo 1
1.Carga de datos
setwd("/cloud/project/")
datos<-read.csv("DerramesEEUU.csv", header = TRUE, sep=";" , dec=".",,na.strings ="-")
str(datos)## 'data.frame': 2760 obs. of 59 variables:
## $ NumeroInforme : int 20100064 20100054 20100092 20100098 20100101 20100102 20100113 20100120 20100039 20100150 ...
## $ NumeroComplementario : int 15072 15114 15120 15127 15130 15132 15146 15162 15197 15205 ...
## $ DiaAccidente : int 8 25 10 28 27 29 11 23 15 11 ...
## $ MesAccidente : int 4 3 5 4 5 5 6 5 3 1 ...
## $ AnioAccidente : int 2010 2010 2010 2010 2010 2010 2010 2010 2010 2010 ...
## $ HoraAccidente : int 6 13 6 24 3 14 7 6 15 2 ...
## $ AmPmAccidente : chr "a. m." "p. m." "a. m." "p. m." ...
## $ IDOperador : int 31684 18779 30829 12105 20160 30003 1248 300 18718 32296 ...
## $ NombreOperador : chr "CONOCOPHILLIPS" "SUNOCO, INC (R&M)" "TEPPCO CRUDE PIPELINE, LLC" "MAGELLAN AMMONIA PIPELINE, L.P." ...
## $ NombreOleoductoInstalacion : chr "GD-03, GOLD LINE" "PHILADELPHIA REFINERY - WEST YARD" "HOBBS TO MIDLAND" "WHITING TO EARLY SEGMENT" ...
## $ UbicacionOleoducto : chr "ONSHORE" "ONSHORE" "ONSHORE" "ONSHORE" ...
## $ TipoOleoducto : chr "ABOVEGROUND" "ABOVEGROUND" "UNDERGROUND" "UNDERGROUND" ...
## $ TipoLiquido : chr "REFINED AND/OR PETROLEUM PRODUCT (NON-HVL), LIQUID" "REFINED AND/OR PETROLEUM PRODUCT (NON-HVL), LIQUID" "CRUDE OIL" "HVL OR OTHER FLAMMABLE OR TOXIC FLUID, GAS" ...
## $ SubtipoLiquido : chr "GASOLINE (NON-ETHANOL)" "OTHER" NA "ANHYDROUS AMMONIA" ...
## $ NombreLiquido : chr NA "VACUUM GAS OIL (VGO)" NA NA ...
## $ CiudadAccidente : chr "GREEN RIDGE" "PHILADELPHIA" "HOBBS" "SCHALLER" ...
## $ CondadoAccidente : chr "PETTIS" "PHILADELPHIA" "LEA" "IDA" ...
## $ EstadoAccidente : chr "MO" "PA" "NM" "IA" ...
## $ LatitudAccidente : chr "38,63064" "39,91934" "32,611" "42,45589" ...
## $ LongitudAccidente : chr "-93,39656" "-75,20447" "-103,0763" "-95,32798" ...
## $ CategoriaCausa : chr "NATURAL FORCE DAMAGE" "MATERIAL/WELD/EQUIP FAILURE" "CORROSION" "MATERIAL/WELD/EQUIP FAILURE" ...
## $ SubcategoriaCausa : chr "TEMPERATURE" "NON-THREADED CONNECTION FAILURE" "EXTERNAL" "CONSTRUCTION, INSTALLATION OR FABRICATION-RELATED" ...
## $ LiberacionInvoluntariaBarriles : chr "0,24" "1700" "2" "0,36" ...
## $ LiberacionIntencionalBarriles : chr "0" "0" NA "0.05" ...
## $ RecuperacionLiquidoBarriles : chr "0,07" "1699" "0,48" "0" ...
## $ PerdidaNetaBarriles : chr "0,17" "1" "1,52" "0,36" ...
## $ IgnicionLiquido : chr "NO" "NO" "NO" "NO" ...
## $ ExplosionLiquido : chr "NO" "NO" "NO" "NO" ...
## $ CierreOleoducto : chr "YES" "YES" "NO" "NO" ...
## $ DiaCierre : int 8 25 NA NA 27 NA NA 23 15 11 ...
## $ MesCierre : int 4 3 NA NA 5 NA NA 5 3 1 ...
## $ AnioCierre : int 2010 2010 NA NA 2010 NA NA 2010 2010 2010 ...
## $ HoraCierre : int 6 18 NA NA 3 NA NA 7 16 2 ...
## $ AmPmCierre : chr "a. m." "p. m." NA NA ...
## $ DiaReinicio : int 9 28 NA NA 27 NA NA 23 15 15 ...
## $ MesReinicio : int 4 3 NA NA 5 NA NA 5 3 1 ...
## $ AnioReinicio : int 2010 2010 NA NA 2010 NA NA 2010 2010 2010 ...
## $ HoraReinicio : int 10 16 NA NA 24 NA NA 9 18 15 ...
## $ AmPmReinicio : chr "a. m." "p. m." NA NA ...
## $ EvacuacionesPublicas : int NA 0 NA NA 0 0 0 0 NA 0 ...
## $ LesionesEmpleadosOperador : int NA NA NA NA NA NA NA NA NA NA ...
## $ LesionesContratistasOperador : int NA NA NA NA NA NA NA NA NA NA ...
## $ LesionesRescatistasEmergencia : int NA NA NA NA NA NA NA NA NA NA ...
## $ OtrasLesiones : int NA NA NA NA NA NA NA NA NA NA ...
## $ LesionesPublico : int NA NA NA NA NA NA NA NA NA NA ...
## $ TodasLesiones : int NA NA NA NA NA NA NA NA NA NA ...
## $ FallecimientosEmpleadosOperador : int NA NA NA NA NA NA NA NA NA NA ...
## $ FallecimientosContratistasOperador : int NA NA NA NA NA NA NA NA NA NA ...
## $ FallecimientosRescatistasEmergencia : int NA NA NA NA NA NA NA NA NA NA ...
## $ OtrosFallecimientos : int NA NA NA NA NA NA NA NA NA NA ...
## $ FallecimientosPublico : int NA NA NA NA NA NA NA NA NA NA ...
## $ TodosFallecimientos : int NA NA NA NA NA NA NA NA NA NA ...
## $ CostosDaniosPropiedad : int 0 0 30000 12000 2720 NA 750 1300 NA 29360 ...
## $ CostosMercanciaPerdidas : int 27 0 100 30 1500 150 300 340 46 136233 ...
## $ CostosDaniosPropiedadesPublicasPrivadas: int 0 0 1000 5000 0 0 0 0 NA NA ...
## $ CostosRespuestaEmergencia : int 0 0 NA 0 1000 NA 400 2445 10999 NA ...
## $ CostosRemediacionAmbiental : int 0 100000 20000 15000 NA NA 6050 3350 452 NA ...
## $ OtrosCostos : int 0 0 NA 0 NA NA 0 2530 NA NA ...
## $ TodosCostos : int 27 100000 51100 32030 5220 150 7500 9965 11497 165593 ...1.1 Extracción de datos
EvacuacionesPublicas <- datos$EvacuacionesPublicas
EvacuacionesPublicas <- na.omit(EvacuacionesPublicas)
# Cantidad de datos de la variable
n <- length(EvacuacionesPublicas)Cantidad total de datos: 2303
1.2 Diagrama de Caja
Dado el gran volumen de datos, se analizo primero el comportamiento general de la variable por ello, se elaboró un diagrama de caja con el fin de identificar el rango en el que se agrupan la mayoría de los valores y obtener así un conjunto representativo de datos.
caja_Evacuaciones<-boxplot(EvacuacionesPublicas , horizontal = TRUE,
col = "tan1",
main = "Gráfica N°1: Distribución de todas las Evacuaciones
públicas de accidentes de oleoductos ocurridos en EE.UU.",
xlab = "Evacuaciones públicas")
En este caso, para lograr un análisis estadístico más representativo, se seleccionó el rango comprendido entre 1 y 90 evacuaciones, que corresponde al intervalo donde se concentra la mayor parte de la distribución.
Evacuaciones_1_90 <- subset(EvacuacionesPublicas, EvacuacionesPublicas >= 1 & EvacuacionesPublicas <= 100)2.Distribución de Frecuencias
Limites e Intervalos
limites <- seq(1, 91, by = 10)
clasificacion <- cut(Evacuaciones_1_90,
breaks = limites,
include.lowest = TRUE,
right = FALSE)Frecuencias Simples
TablaEvacuacion <- as.data.frame(table(clasificacion))
colnames(TablaEvacuacion) <- c("Intervalo", "ni")
TablaEvacuacion$hi <- round((TablaEvacuacion$ni / sum(TablaEvacuacion$ni)) * 100, 2)Frecuencias Acumuladas
# Frecuencia acumulada absoluta
TablaEvacuacion$Niasc <- cumsum(TablaEvacuacion$ni)
TablaEvacuacion$Nidsc <- rev(cumsum(rev(TablaEvacuacion$ni)))
# Frecuencia acumulada relativa (%)
TablaEvacuacion$Hiasc <- cumsum(TablaEvacuacion$hi)
TablaEvacuacion$Hidsc <- rev(cumsum(rev(TablaEvacuacion$hi)))Tabla de distribución de frecuencias
TDFFinalEvacuaciones<- rbind(TablaEvacuacion, data.frame(
Intervalo = "TOTAL",
ni = sum(TablaEvacuacion$ni),
hi = 100,
Niasc = " ",
Hiasc = " ",
Nidsc = " ",
Hidsc = " "
))
library(gt)
tablaEvacuaciones <- TDFFinalEvacuaciones %>%
gt() %>%
cols_label(
Intervalo = md("**Intervalo**"),
ni = md("**ni**"),
hi = md("**hi (%)**"),
Niasc = md("**Ni ↑**"),
Hiasc = md("**Hi ↑ (%)**"),
Nidsc = md("**Ni ↓**"),
Hidsc = md("**Hi ↓ (%)**")
) %>%
tab_header(
title = md("**Tabla N°1**"),
subtitle = md("**Distribución de las Evacuaciones públicas por accidentes ocurridos en EE.UU (2010-2017)**")
) %>%
tab_source_note(
source_note = md("Autor: Grupo 1")
) %>%
tab_options(
table.background.color = "white",
row.striping.background_color = "white",
table.border.top.color = "black",
table.border.bottom.color = "black",
table.border.top.style = "solid",
table.border.bottom.style = "solid",
column_labels.border.top.color = "black",
column_labels.border.bottom.color = "black",
column_labels.border.bottom.width = px(2),
heading.border.bottom.color = "black",
heading.border.bottom.width = px(2),
table_body.hlines.color = "gray",
table_body.border.bottom.color = "black"
) %>%
tab_style(
style = cell_text(weight = "bold"),
locations = cells_body(
rows = Intervalo == "TOTAL"
)
)
tablaEvacuaciones| Tabla N°1 | ||||||
| Distribución de las Evacuaciones públicas por accidentes ocurridos en EE.UU (2010-2017) | ||||||
| Intervalo | ni | hi (%) | Ni ↑ | Ni ↓ | Hi ↑ (%) | Hi ↓ (%) |
|---|---|---|---|---|---|---|
| [1,11) | 23 | 46.94 | 23 | 49 | 46.94 | 99.99 |
| [11,21) | 10 | 20.41 | 33 | 26 | 67.35 | 53.05 |
| [21,31) | 6 | 12.24 | 39 | 16 | 79.59 | 32.64 |
| [31,41) | 4 | 8.16 | 43 | 10 | 87.75 | 20.4 |
| [41,51) | 1 | 2.04 | 44 | 6 | 89.79 | 12.24 |
| [51,61) | 1 | 2.04 | 45 | 5 | 91.83 | 10.2 |
| [61,71) | 2 | 4.08 | 47 | 4 | 95.91 | 8.16 |
| [71,81) | 1 | 2.04 | 48 | 2 | 97.95 | 4.08 |
| [81,91] | 1 | 2.04 | 49 | 1 | 99.99 | 2.04 |
| TOTAL | 49 | 100.00 | ||||
| Autor: Grupo 1 | ||||||
Gráfica
par(mar = c(7, 6, 4, 2))
barplot(
TablaEvacuacion$ni,
main = "Gráfica N°2: Distribución de evacuaciones públicas\npor accidentes en oleoductos de EE.UU.",
ylab = "Cantidad",
names.arg = TablaEvacuacion$Intervalo,
col = "darkseagreen2",
las = 2,
cex.main = 1.1,
cex.lab = 1.1,
cex.axis = 0.9,
cex.names = 0.8,
)
mtext("Evacuaciones públicas", side = 1, line = 4, cex = 1)
3.Conjetura de Modelo
Se conjetura que la variable EvacuacionesPublicas, podría seguir una distribución geometrica, bajo el supuesto de que la probabilidad de ocurrencia es máxima para valores bajos de la variable y disminuye progresivamente a medida que el número de evacuaciones aumenta..
3.1 Definición de Hipótesis
Hipótesis nula(Ho): Las evacuaciones públicas realizadas por accidentes en oleductos siguen una distribución geométrica
Hipótesis alternativa (H1): Las evacuaciones públicas realizadas por accidentes en oleductos NO siguen una distribución geométrica
3.2 Ajuste del modelo geometrico
Se calcula el promedio de evacuaciones registradas. Con este promedio, se estima el parámetro del modelo geométrico, que representa la probabilidad de que ocurra una evacuación en cada intento.
media_obs <- mean(Evacuaciones_1_90)
p <- 1 / (media_obs + 1)Parámetro estimado (p): 0.05
3.2 Extracción de límites de intervalos
Se extraen los valores iniciales y finales de cada rango de evacuaciones para poder analizarlos.
Intervalo_clean <- gsub("\\[|\\]|\\(|\\)", "", TablaEvacuacion$Intervalo)
Inicio <- as.numeric(sub(",.*", "", Intervalo_clean))
Fin <- as.numeric(sub(".*,", "", Intervalo_clean))3.3 Distribución de Frecuencias
3.3.1 Frecuencias Observadas
Fo <- TablaEvacuacion$niProbabilidades geométricas y frecuencias esperadas
k <- length(Inicio)
P_geom <- numeric(k)
Fe <- numeric(k)
for(i in 1:k){
P_geom[i] <- pgeom(Fin[i]-1, prob = p) - pgeom(Inicio[i]-1, prob = p)
Fe[i] <- P_geom[i] * sum(TablaEvacuacion$ni)
}3.4 Gráfica del Modelo
barplot(rbind(Fo, Fe),
main = "Gráfica N°3: Comparación Modelo Geométrico vs Observado",
beside = TRUE,
col = c("darkseagreen2", "#FFC685"),
names.arg = TablaEvacuacion$Intervalo,
ylab = "Cantidad",
las = 2,
cex.names = 0.8,
cex.axis = 1,
cex.main = 1.1)
mtext("Evacuaciones públicas", side = 1, line = 3.5, cex = 1)
legend(x = 27, y = 170,
legend = c("Observado", "Uniforme"),
fill = c("darkseagreen2", "#FFC685"),
bty = "o",
y.intersp = 0.7,
cex = 0.8)
4.Tests
4.1 Test de Pearson
Correlación de frecuencias
Correlacion_geo <- cor(Fo, Fe) * 100La correlación de frecuencias es de = 97.94 %
Gráfica de correlación Fo vs Fe
plot(Fo, Fe,
main = "Gráfica N°4: Correlación de frecuencias en el modelo Geométrico",
xlab = "Frecuencia Observada ",
ylab = "Frecuencia Esperada",
col = "darkseagreen3", pch = 19)
abline(lm(Fe ~ Fo), col = "red", lwd = 2)
4.2 Test de Bondad de ajuste
4.2.1 Cálculo del Estadístico Chi-cuadrado
- Estadístico chi-cuadrado
x2_g <- sum((Fo - Fe)^2 / Fe)El estadistico Chi-cuadrado es: 5.637477
- Cálculo del Umbral de Aceptación
Grados de Libertad
gl_g <- (k - 1) - 1 Definición del nivel de significancia
nivel_significancia <- 0.05Umbral de aceptación
umbral_aceptacion<- qchisq(1 - nivel_significancia, gl_g)El umbral de aceptación es: 14.06714
4.1 Decisión
if (x2_g < umbral_aceptacion) {
cat("Conclusión: No se rechaza H0, las evacuaciones públicas realizadas por accidentes en oleductos podrian seguir una distribución geométrica")
} else {
cat("Conclusión: Se rechaza H0, las evacuaciones públicas realizadas por accidentes en oleductos NO siguen una distribución geométrica")
}Conclusión: No se rechaza H0, las evacuaciones públicas realizadas por accidentes en oleductos podrian seguir una distribución geométrica
4.2 Tabla resumen de test4
Variable <- c("Evacuaciones públicas")
Modelo <- c("Geométrico")
Tabla_resumen <- data.frame(Variable,
Modelo,
Pearson = round(Correlacion_geo,2),
Chi_Cuadrado = round(x2_g,2),
Umbral = round(umbral_aceptacion,2),
TestChi = c("Aprobado"))
colnames(Tabla_resumen) <- c("Variable",
"Modelo",
"Test Pearson (%)",
"Chi-Cuadrado",
"Umbral de aceptación",
"Test de Bondad de ajuste")
library(gt)
Tabla_resumen %>%
gt() %>%
tab_header(
title = md("**Tabla N°2**"),
subtitle = md("**Resumen de los Tests Aplicados al Modelo Geométrico**")
) %>%
tab_source_note(
source_note = md("Autor: Grupo 1")
) %>%
cols_align(
align = "center",
columns = everything()
) %>%
tab_options(
table.border.top.color = "black",
table.border.bottom.color = "black",
table.border.top.style = "solid",
table.border.bottom.style = "solid",
column_labels.font.weight = "bold",
column_labels.border.top.color = "black",
column_labels.border.bottom.color = "black",
column_labels.border.bottom.width = px(2),
heading.border.bottom.color = "black",
heading.border.bottom.width = px(2),
table_body.hlines.color = "grey",
table_body.border.bottom.color = "black"
)| Tabla N°2 | |||||
| Resumen de los Tests Aplicados al Modelo Geométrico | |||||
| Variable | Modelo | Test Pearson (%) | Chi-Cuadrado | Umbral de aceptación | Test de Bondad de ajuste |
|---|---|---|---|---|---|
| Evacuaciones públicas | Geométrico | 97.94 | 5.64 | 14.07 | Aprobado |
| Autor: Grupo 1 | |||||
Los resultados de los tests indican que los datos de Evacuaciones públicas pueden seguir una distribución geométrica.
El estadístico Chi-cuadrado calculado (5.08) es menor que el umbral de aceptación (14.07), lo que indica que no se rechaza la hipótesis nula (H₀) y que el modelo geométrico es adecuado para describir los datos observados.
El Test de Pearson reporta un nivel de ajuste del 97.94%, lo que indica que las frecuencias observadas se aproximan muy bien a las frecuencias esperadas bajo la distribución geométrica, respaldando adicionalmente la validez del modelo.
5.Cálculo de probabilidades
- ¿Qué probabilidad hay de que se registren entre 1 y 11 evacuaciones públicas en un accidente de oleoductos en EE.UU.?
5.1 Estimación de Probabilidad
# Rango de valores
rango_inicio <- 1
rango_fin <- 11
# Probabilidad total para el rango
prob_rango <- pgeom(rango_fin - 1, prob = p) - pgeom(rango_inicio - 1, prob = p)Probabilidad de que se de, de 1 a 11 evacuaciones: 38.06 %
5.2 Gráfica de Probabilidad
TablaEvacuacion$P_geom <- P_geom
colores <- ifelse(Inicio >= rango_inicio & Fin <= rango_fin, "#D75C20", "#FFC685")
# Gráfico
barplot(height = TablaEvacuacion$P_geom,
names.arg = TablaEvacuacion$Intervalo,
main = "Gráfica N°5: Probabilidad Geométrica - Evacuaciones públicas",
ylab = "Probabilidad (%)",
col = colores,
las = 2,
cex.names = 0.8,
ylim = c(0, max(TablaEvacuacion$P_geom)+0.1)
)
mtext("Evacuaciones públicas", side = 1, line = 3.5, cex = 1)
# Leyenda
legend("topright",
legend = c("Evacuaciones fuera del rango", "Rango de Evacuaciones"),
fill = c("#D75C20", "#FFC685"),
border = NA,
cex = 0.8)
6.Teorema del Límite Central
El Teorema del Límite Central establece que, aunque las variables individuales de una población no sigan una distribución normal, la distribución de las medias muestrales tiende a ser normal cuando el tamaño de la muestra es suficientemente grande (n ≥ 30).
Por lo tanto, es posible estimar la media poblacional mediante intervalos de confianza, aun cuando la distribución original de los datos no sea normal.Esto se puede hacer con tres postulados principales:
Donde:
- x es la media aritmética muestral
- e es el margen de error de la media
Media aritmética muestral
x <- mean(Evacuaciones_1_90)La media muestral es de: 19.04082
Desviación estándar muestral
sigma_g<- sd(Evacuaciones_1_90)La desviación estandar muestral es de: 20.90949
Error estándar de la media
n <- length(Evacuaciones_1_90)
e <- sigma_g/ sqrt(n)El error estandar de la media es de: 2.98707
Intervalo de Confianza del 95%
limite_inferior <- x - 2 * e
limite_superior <- x + 2 * eEl limite inferior es: 13.06668
El limite superior es: 25.01496
Tabla
tabla_media_exp <- data.frame(
round(limite_inferior, 2),
round(x, 2),
round(limite_superior, 2),
round(sigma_g, 2)
)
colnames(tabla_media_exp) <- c("Límite inferior", "Media poblacional",
"Límite superior", "Desviación estándar poblacional")
library(gt)
tabla_media_exp%>%
gt() %>%
tab_header(
title = md("**Tabla N°3**"),
subtitle = md("**Parámetros poblacionales de la variable **EvacuacionesPublicas** de los accidentes en oleoductos ocurridos en EE.UU.(2010-2017)**")
) %>%
tab_source_note(
source_note = md("Autor: Grupo 1")
) %>%
cols_align(
align = "center",
columns = everything()
) %>%
tab_options(
table.border.top.color = "black",
table.border.bottom.color = "black",
table.border.top.style = "solid",
table.border.bottom.style = "solid",
column_labels.font.weight = "bold",
column_labels.border.top.color = "black",
column_labels.border.bottom.color = "black",
column_labels.border.bottom.width = px(2),
heading.border.bottom.color = "black",
heading.border.bottom.width = px(2),
table_body.hlines.color = "grey",
table_body.border.bottom.color = "black"
)| Tabla N°3 | |||
| Parámetros poblacionales de la variable EvacuacionesPublicas de los accidentes en oleoductos ocurridos en EE.UU.(2010-2017) | |||
| Límite inferior | Media poblacional | Límite superior | Desviación estándar poblacional |
|---|---|---|---|
| 13.07 | 19.04 | 25.01 | 20.91 |
| Autor: Grupo 1 | |||
7.Conclusión
La variable EvacuacionesPublicas se ajusta a una distribución geometrica, con una media poblacional 19.04 y una desviación estándar poblacional de 20.91.
La variable EvacuacionesPublicas fue considerada inicialmente bajo un modelo geométrico, con una media muestral de 19.04 y una desviación estándar de 20.91.
De acuerdo con este modelo, se estima que la probabilidad de que ocurran entre 1 y 11 evacuaciones públicas es aproximadamente del 38.06%, lo que indica que cerca de 2 de cada 5 accidentes presentarían un número de evacuaciones dentro de este rango específico.
Aplicando el Teorema del Límite Central, se estima que la media
poblacional se encuentra entre 13.07 y 25.01 con un 95% de confianza,
proporcionando un nivel de certeza adecuado sobre el rango probable en
el que se ubican las evacuaciones promedio.