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 : num 38.6 39.9 32.6 42.5 30.2 ...
## $ LongitudAccidente : num -93.4 -75.2 -103.1 -95.3 -91.2 ...
## $ 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 : num 0.24 1700 2 0.36 1.31 ...
## $ LiberacionIntencionalBarriles : chr "0" "0" NA "0.05" ...
## $ RecuperacionLiquidoBarriles : num 0.07 1699 0.48 0 0 ...
## $ PerdidaNetaBarriles : num 0.17 1 1.52 0.36 1.31 ...
## $ 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 ...
2.Histograma de la variable
options(scipen = 999)
hist(variable_original,
freq = FALSE,
ylim = c(0,0.08),
col = "lightgreen",
main = "Histograma de Latitud de los accidentes",
xlab = "Latitud del accidente",
ylab = "Frecuencia relativa (Probabilidad)")
# Ajuste del modelo lognormal
log_var1 <- log(variable_original)
ulog1 <- mean(log_var1)
sigmalog1 <- sd(log_var1)
curve(dlnorm(x, meanlog = ulog1, sdlog = sigmalog1),
col = "skyblue4",
lwd = 3,
add = TRUE)
3.Conjetura del Modelo
3.1 Filtrar datos
caja <- boxplot(variable_original, plot = FALSE)
variable_modelo <- variable_original[variable_original < min(caja$out)]
3.2 Parametros Lognormal
log_var <- log(variable_modelo)
ulog <- mean(log_var)
sigmalog <- sd(log_var)
3.3 GRAFICO 2 — ZOOM CON MODELO LOGNORMAL
xmin <- 20
xmax <- 50
HistZoom <- hist(variable_modelo,
freq = FALSE,
breaks = 7,
xlim = c(xmin, xmax),
ylim = c(0,0.08),
col = "darkseagreen3",
main = "Histograma de Latitud con modelo Log-Normal",
xlab = "Latitud (°)",
ylab = "Density")
curve(dlnorm(x, meanlog = ulog, sdlog = sigmalog),
from = xmin,
to = xmax,
col = "red",
lwd = 3,
add = TRUE)
legend("topright",
legend = "Modelo Log-Normal",
col = "red",
lwd = 3,
bty = "n")
3.4 Cálculo de Frecuencias
Frecuencias Observadas y Esperadas
Fo <- HistZoom$counts
h <- length(Fo)
P <- numeric(h)
for(i in 1:h){
P[i] <- plnorm(HistZoom$breaks[i+1], ulog, sigmalog) -
plnorm(HistZoom$breaks[i], ulog, sigmalog)
}
Fe <- P * length(variable_modelo)
4.Test
4.1 Test de Pearson
Correlación de frecuencias
Correlacion_log <- cor(Fo, Fe) * 100
4.2 Test de Bondad de ajuste
Estadístico chi-cuadrado
n <- length(variable_modelo)
Fo_pct <- (Fo/n) * 100
Fe_pct <- P * 100
x2 <- sum((Fe_pct - Fo_pct)^2 / Fe_pct)
gl <- (h - 1) - 2
critico <- qchisq(0.95, gl)
4.3 Tabla resumen de test
Tabla_resumen <- data.frame(Variable = "Latitud",
Modelo = "Log-Normal",
Pearson = round(Correlacion_log,2),
Chi_Cuadrado = round(x2,2),
Umbral = round(critico,2),
Resultado = decision
)
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°1**"),
subtitle = md("**Resumen de los Tests Aplicados al Modelo Normal**")
) %>%
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°1 |
| Resumen de los Tests Aplicados al Modelo Normal |
| Variable |
Modelo |
Test Pearson (%) |
Chi-Cuadrado |
Umbral de aceptación |
Test de Bondad de ajuste |
| Latitud |
Log-Normal |
93.81 |
7.42 |
7.81 |
Aprobado |
| Autor: Grupo 1 |