Justificación de la variable

La longitud se parametriza como una variable cuantitativa continua fundamental para cartografiar la varianza Este-Oeste de los siniestros. Su análisis permite validar si los accidentes siguen una distribución de mezcla probabilística dictada por la ubicación de la infraestructura industrial.

1 Cargar datos

Importamos el archivo “database (1).csv” desde una ruta local y lo almacena en el objeto datos, usando espacios o tabulaciones como separador.

datos <- read.csv("database-_1_.csv")
zona<-datos$Accident.Longitude

2 Tabla de frecuencia

Extraemos la variable costos totales, omitimos las celdas en blanco o valores iguales a cero y verificamos el tamaño muestral En la tabla de distribución de frecuencias de la variable Longitud del Accidente, el número de clases se determinó mediante la regla de Sturges y el ancho de clase se calculó a partir del rango geoespacial total de los datos, asegurando una cobertura completa desde la coordenada más occidental hasta la más oriental.

df <- read.csv("database-_1_.csv")
Accident_Longitude <- df$Accident.Longitude
Accident_Longitude <- na.omit(Accident_Longitude)
xmin <- min(Accident_Longitude)
xmax <- max(Accident_Longitude)
R <- xmax - xmin
K <- floor(1 + 3.3 * log10(length(Accident_Longitude)))
A <- R / K
Li <- round(seq(from = xmin, by = A, length.out = K), 2)
Ls <- round(seq(from = xmin + A, by = A, length.out = K), 2)
MC <- round((Li + Ls) / 2, 2)
ni <- numeric(K)
for (i in 1:(K-1)) {
  ni[i] <- sum(Accident_Longitude >= Li[i] & Accident_Longitude < Ls[i])
}
ni[K] <- sum(Accident_Longitude >= Li[K] & Accident_Longitude <= xmax)
hi <- ni / sum(ni) * 100
Ni_asc <- cumsum(ni)
Ni_desc <- rev(cumsum(rev(ni)))
Hi_asc <- cumsum(hi)
Hi_desc <- rev(cumsum(rev(hi)))
TDF <- data.frame(
  Li, Ls, MC, ni, 
  hi_porc = round(hi, 2), 
  Ni_asc, Ni_desc, 
  Hi_asc_porc = round(Hi_asc, 2), 
  Hi_desc_porc = round(Hi_desc, 2)
)
kable(TDF, 
      caption = "Tabla No. 1: Distribución de Frecuencias de longitud de accidente", 
      col.names = c("Lím. Inf.", "Lím. Sup.", "Marca Clase", "ni", "hi (%)", "Ni Asc.", "Ni Desc.", "Hi Asc. (%)", "Hi Desc. (%)"), 
      digits = 2)

Tabla No. 1: Distribución de Frecuencias de longitud de accidente
Lím. Inf.	Lím. Sup.	Marca Clase	ni	hi (%)	Ni Asc.	Ni Desc.	Hi Asc. (%)	Hi Desc. (%)
-158.10	-136.24	-147.17	14	0.50	14	2795	0.50	100.00
-136.24	-114.37	-125.31	169	6.05	183	2781	6.55	99.50
-114.37	-92.51	-103.44	1814	64.90	1997	2612	71.45	93.45
-92.51	-70.65	-81.58	792	28.34	2789	798	99.79	28.55
-70.65	-48.78	-59.72	2	0.07	2791	6	99.86	0.21
-48.78	-26.92	-37.85	0	0.00	2791	4	99.86	0.14
-26.92	-5.05	-15.99	1	0.04	2792	4	99.89	0.14
-5.05	16.81	5.88	0	0.00	2792	3	99.89	0.11
16.81	38.67	27.74	0	0.00	2792	3	99.89	0.11
38.67	60.54	49.61	0	0.00	2792	3	99.89	0.11
60.54	82.40	71.47	0	0.00	2792	3	99.89	0.11
82.40	104.26	93.33	3	0.11	2795	3	100.00	0.11

3 Nueva tabla de frecuencia

Se seleccionó el rango central de la variable Longitud del Accidente para el análisis, omitiendo las coordenadas atípicas (outliers), debido a que en esta zona se concentra la mayor densidad de los datos. Esta elección permite construir la tabla de frecuencia y gráficas más claras y legibles, facilitando la interpretación de la distribución espacial de los siniestros y evitando distorsiones visuales provocadas por ubicaciones extremas con baja frecuencia

database <- read.csv("database-_1_.csv", check.names = FALSE)
Accident_Longitude <- na.omit(database$`Accident Longitude`)
umbral_90 <- quantile(Accident_Longitude, 0.90)
datos_zoom <- Accident_Longitude[Accident_Longitude <= umbral_90]
n_z <- length(datos_zoom)
xmin_z <- min(datos_zoom)
xmax_z <- max(datos_zoom)
K_z <- floor(1 + 3.322 * log10(n_z))
R_z <- xmax_z - xmin_z
A_z <- R_z / K_z
cortes_z <- seq(xmin_z, xmin_z + (K_z * A_z), length.out = K_z + 1)
Li_z <- cortes_z[1:K_z]
Ls_z <- cortes_z[2:(K_z + 1)]
MC_z <- (Li_z + Ls_z) / 2
ni_z <- as.vector(table(cut(datos_zoom, breaks = cortes_z, include.lowest = TRUE)))
hi_z <- (ni_z / n_z) * 100
Ni_asc_z  <- cumsum(ni_z)
Ni_desc_z <- rev(cumsum(rev(ni_z)))
Hi_asc_z  <- cumsum(hi_z)
Hi_desc_z <- rev(cumsum(rev(hi_z)))
TDF_final_zoom <- data.frame(
  Li      = round(Li_z, 2),
  Ls      = round(Ls_z, 2),
  MC      = round(MC_z, 2),
  ni      = ni_z,
  hi_porc = round(hi_z, 2),
  Ni_asc  = Ni_asc_z,
  Ni_desc = Ni_desc_z,
  Hi_asc  = round(Hi_asc_z, 2),
  Hi_desc = round(Hi_desc_z, 2)
)
TDF_final_zoom <- TDF_final_zoom[TDF_final_zoom$ni > 0, ]
kable(TDF_final_zoom, 
      caption = "Tabla No. 2: Distribución de Frecuencias de Accidente Longitud",
      align = 'c',
      row.names = FALSE,
      col.names = c("Lím. Inf.", "Lím. Sup.", "Marca Clase", "ni", "hi (%)", 
                    "Ni Asc.", "Ni Desc.", "Hi Asc. (%)", "Hi Desc. (%)"))

Tabla No. 2: Distribución de Frecuencias de Accidente Longitud
Lím. Inf.	Lím. Sup.	Marca Clase	ni	hi (%)	Ni Asc.	Ni Desc.	Hi Asc. (%)	Hi Desc. (%)
-158.10	-151.88	-154.99	3	0.12	3	2515	0.12	100.00
-151.88	-145.66	-148.77	11	0.44	14	2512	0.56	99.88
-127.00	-120.78	-123.89	37	1.47	51	2501	2.03	99.44
-120.78	-114.56	-117.67	132	5.25	183	2464	7.28	97.97
-114.56	-108.34	-111.45	41	1.63	224	2332	8.91	92.72
-108.34	-102.12	-105.23	339	13.48	563	2291	22.39	91.09
-102.12	-95.89	-99.00	803	31.93	1366	1952	54.31	77.61
-95.89	-89.67	-92.78	847	33.68	2213	1149	87.99	45.69
-89.67	-83.45	-86.56	302	12.01	2515	302	100.00	12.01

4 Cantidad general de longitud del accidente

Esta visualización macroscópica modela la función de masa espacial en su dominio absoluto. El histograma global permite identificar a simple vista la dispersión total de los eventos y la presencia de asimetrías severas generadas por ubicaciones extremas.

library(ggplot2)
library(dplyr)
library(scales)
datos <- read.csv("database-_1_.csv")
variable_interes <- na.omit(datos$Accident.Longitude)
volumen <- variable_interes
k <- 1 + (3.322 * log10(length(volumen)))
# Usamos un ancho de bin basado en Sturges para el plot
bins_calc <- floor(k)

p_global <- ggplot(data.frame(val=volumen), aes(x = val)) +
  geom_histogram(bins = 30, fill = "steelblue", color = "white", alpha = 0.8) +
  scale_x_continuous(labels = number_format(suffix = "°")) +
  scale_y_continuous(expand = expansion(mult = c(0, 0.05))) +
  
  labs(
    title = "Gráfica 1: Distribución Global de Longitud de Accidentes", 
    x = "Longitud Geográfica (°)",
    y = "Cantidad de Accidentes"
  ) +
  
  theme_classic() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    axis.text.x = element_text(angle = 45, hjust = 1)
  )

print(p_global)

5 Cantidad local de longitud del accidente

Al aislar probabilísticamente las colas de la distribución, este apartado evalúa exclusivamente el “core” operativo. Esta amplificación visual revela el comportamiento real de la densidad de incidentes en la franja Este-Oeste de mayor criticidad.

library(ggplot2)
library(dplyr)
library(scales)

# Definimos los límites del ZOOM automáticamente basado en tus datos
# (Para simular el "0-2 barriles" pero adaptado a coordenadas)
limit_min_zoom <- quantile(variable_interes, 0.10) # Cortamos el 10% inferior
limit_max_zoom <- quantile(variable_interes, 0.90) # Cortamos el 10% superior

datos_zoom <- datos %>%
  filter(!is.na(Accident.Longitude)) %>%
  filter(Accident.Longitude >= limit_min_zoom & Accident.Longitude <= limit_max_zoom)

p_zoom_5barras <- ggplot(datos_zoom, aes(x = Accident.Longitude)) +
  geom_histogram(bins = 10, fill = "steelblue", color = "white", alpha = 0.8) +
  scale_x_continuous(labels = number_format(accuracy = 0.1, suffix = "°")) +
  
  scale_y_continuous(expand = expansion(mult = c(0, 0.05))) +
  
  labs(
    title = paste("Gráfica No2: Distribución detallada en rango principal\n(", round(limit_min_zoom,1), "° a", round(limit_max_zoom,1), "°)"), 
    x = "Longitud Geográfica (°)",
    y = "Cantidad"
  ) +
  
  theme_classic() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    axis.text.x = element_text(hjust = 0.5)
  )

print(p_zoom_5barras)

6 Cantidad relativa local de longitud del accidente

La transformación estocástica de frecuencias absolutas a proporciones relativas internas permite evaluar la probabilidad de ocurrencia de un evento dentro de los sub-intervalos del núcleo principal, estandarizando la escala de análisis.

# Usamos el mismo dataset 'datos_zoom' filtrado arriba
p_zoom_5barras_pct <- ggplot(datos_zoom, aes(x = Accident.Longitude)) +
  geom_histogram(
    aes(y = after_stat(count / sum(count))), 
    bins = 10, 
    fill = "steelblue", 
    color = "white", 
    alpha = 0.8
  ) +
  scale_x_continuous(labels = number_format(accuracy = 0.1, suffix = "°")) +
  scale_y_continuous(
    labels = scales::percent_format(accuracy = 1), 
    expand = expansion(mult = c(0, 0.05))
  ) +
  
  labs(
    title = "Gráfica No3: Distribución porcentual local del rango principal", 
    x = "Longitud Geográfica (°)",
    y = "Porcentaje (%)"
  ) +
  
  theme_classic() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    axis.text.x = element_text(hjust = 0.5)
  )

print(p_zoom_5barras_pct)

7 Cantidad relativa global de longitud del accidente

Para dimensionar correctamente el peso de la zona de aglomeración, este apartado contrasta las frecuencias del rango principal contra el universo total de los datos. Mide la probabilidad empírica de que un accidente aleatorio pertenezca a este corredor longitudinal específico.

p_zoom_5barras_100 <- ggplot(datos_zoom, aes(x = Accident.Longitude)) +
  geom_histogram(
    aes(y = after_stat(count / sum(count))), 
    bins = 10, 
    fill = "steelblue", 
    color = "white", 
    alpha = 0.8
  ) +
  scale_x_continuous(labels = number_format(accuracy = 0.1, suffix = "°")) +
  scale_y_continuous(
    labels = scales::percent_format(accuracy = 1),
    limits = c(0, 1), # Escala de 0 a 100%
    expand = expansion(mult = c(0, 0.05))
  ) +
  
  labs(
    title = "Gráfica No4: Distribución porcentual global del rango principal", 
    x = "Longitud Geográfica (°)",
    y = "Porcentaje (%)"
  ) +
  
  theme_classic() +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    axis.text.x = element_text(hjust = 0.5)
  )

print(p_zoom_5barras_100)

8 Ojivas combinadas de la cantidad absoluta

El modelado de la función de distribución acumulada (FDA) bidireccional cartografía la progresión de los incidentes espaciales. La intersección de las trayectorias Oeste-Este (ascendente) y Este-Oeste (descendente) señala gráficamente la mediana espacial de la siniestralidad.

library(ggplot2)
library(dplyr)
library(scales)

# Datos locales (del rango principal definido antes)
datos_local <- variable_interes[variable_interes >= limit_min_zoom & variable_interes <= limit_max_zoom]
k_local <- 10 # Número de cortes para la ojiva
min_val <- min(datos_local)
max_val <- max(datos_local)
R_local <- max_val - min_val
A_local <- R_local / k_local 

Li_num <- seq(min_val, max_val - A_local, length.out = k_local)
Ls_num <- Li_num + A_local
ni_local <- numeric(k_local)

for(i in 1:k_local){
  if(i == k_local){
    ni_local[i] <- sum(datos_local >= Li_num[i] & datos_local <= max_val)
  } else {
    ni_local[i] <- sum(datos_local >= Li_num[i] & datos_local < Ls_num[i])
  }
}

Niasc <- cumsum(ni_local)
Nidsc <- rev(cumsum(rev(ni_local)))

datos_asc <- data.frame(
  x = c(min_val, Ls_num),  
  y = c(0, Niasc),          
  Tipo = "Ascendente"
)

datos_dsc <- data.frame(
  x = c(Li_num, max_val),  
  y = c(Nidsc, 0),          
  Tipo = "Descendente"
)

datos_ojivas_plot <- rbind(datos_asc, datos_dsc)

p_ojiva_cruzada_solida <- ggplot(datos_ojivas_plot, aes(x = x, y = y, color = Tipo, linetype = Tipo)) +
  geom_line(linewidth = 0.8) +
  geom_point(size = 2) +
  
  scale_x_continuous(
    labels = scales::number_format(accuracy = 0.1, suffix = "°"),
    breaks = scales::pretty_breaks(n = 6)
  ) +
  scale_color_manual(values = c("Ascendente" = "black", "Descendente" = "blue")) +
  scale_linetype_manual(values = c("Ascendente" = "solid", "Descendente" = "solid")) +
  
  labs(
    title = "Gráfica 5: Ojivas Geoespaciales del rango principal",
    x = "Longitud Geográfica (°)",
    y = "Cantidad Acumulada",
    color = NULL,
    linetype = NULL
  ) +
  
  theme_bw() + 
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 14),
    legend.position = c(0.85, 0.5),
    legend.background = element_rect(color = "black", fill = "white"),
    axis.text = element_text(color = "black")
  )

print(p_ojiva_cruzada_solida)

9 Diagrama de cajas

Mediante un enfoque no paramétrico, este diagrama de bigotes evalúa la dispersión central y el rango intercuartílico (RIC) de las coordenadas. Es la herramienta visual definitiva para comprobar el apuntamiento de los datos y el sesgo direccional de la mediana frente a la media.

variable_box <- datos_local # Usamos el filtro local
par(mar = c(5, 2, 4, 2)) 
b <- boxplot(variable_box, 
        horizontal = TRUE, 
        col = "skyblue", 
        border = "gray30",          
        medcol = "red",             
        boxwex = 0.6,               
        outline = FALSE,            
        main = "Gráfica 6: Distribución de Longitud (Rango Principal)",
        xlab = "Longitud (°)",
        xaxt = "n",  
        yaxt = "n",  
        frame = FALSE 
)

limite_visible_inf <- min(variable_box)
limite_visible_sup <- max(variable_box)
puntos_eje <- pretty(c(limite_visible_inf, limite_visible_sup))              
axis(1, at = puntos_eje, labels = format(puntos_eje, nsmall = 2), col = "gray30", col.axis = "gray30")
grid(nx = NULL, ny = NA, col = "lightgray", lty = "dotted", lwd = 1)

10 Tabla estadísticos

Como marco de referencia absoluto, se calculan los estimadores univariantes para la totalidad del vector espacial. La comparación de esta tabla con los estadísticos locales permite medir matemáticamente el grado de distorsión que introducen las coordenadas extremas.

library(e1071)
library(knitr)

variable_analisis <- variable_interes # Global
n <- length(variable_analisis)
k_global <- floor(1 + 3.322 * log10(n))
R <- max(variable_analisis) - min(variable_analisis)
A_global <- R / k_global
Li <- seq(min(variable_analisis), max(variable_analisis) - A_global, length.out = k_global)

if(max(Li) + A_global < max(variable_analisis)) { 
   Li <- c(Li, tail(Li, 1) + A_global) 
   k_global <- k_global + 1
}
Ls <- Li + A_global
MC <- (Li + Ls) / 2

ni <- numeric(length(MC))
for(i in 1:length(MC)){
  if(i == length(MC)) ni[i] <- sum(variable_analisis >= Li[i] & variable_analisis <= (max(variable_analisis) + 0.001))
  else ni[i] <- sum(variable_analisis >= Li[i] & variable_analisis < Ls[i])
}
media_agrupada <- sum(MC * ni) / sum(ni)
desviacion_estandar <- sd(variable_analisis)
error_estandar <- desviacion_estandar / sqrt(n)
margen_error <- 1.96 * error_estandar
ic_inferior <- media_agrupada - margen_error
ic_superior <- media_agrupada + margen_error
texto_media_intervalo <- paste0("[", format(round(ic_inferior, 2), big.mark=","), " ; ", format(round(ic_superior, 2), big.mark=","), "]")
ri <- min(variable_analisis)
rs <- max(variable_analisis)
mediana <- median(variable_analisis)
t <- table(round(variable_analisis, 1))
Mo <- as.numeric(names(t)[which.max(t)])
cv <- (desviacion_estandar / abs(media_agrupada)) * 100
As <- skewness(variable_analisis)
K <- kurtosis(variable_analisis)

Tabla_global <- data.frame(
  "Longitud (Global)",
  paste(format(ri, nsmall=2), "°"),
  paste(format(rs, nsmall=2), "°"),
  texto_media_intervalo, 
  paste(format(round(mediana, 2), big.mark=","), "°"),
  paste(format(round(Mo, 2), big.mark=","), "°"),
  paste(format(round(desviacion_estandar, 2), big.mark=","), "°"),
  paste(round(cv, 2), "%"),
  round(As, 2),
  round(K, 2)
)

colnames(Tabla_global) <- c("Variable","Min","Max","Media (IC 95%)","Mediana","Moda","Desv. S","CV","As","K")
kable(Tabla_global, format = "markdown", caption = "Tabla 2: Indicadores Globales.")

Tabla 2: Indicadores Globales.
Variable	Min	Max	Media (IC 95%)	Mediana	Moda	Desv. S	CV	As	K
Longitud (Global)	-158.0999 °	104.2634 °	[-98.97 ; -98.06]	-95.49 °	-95.2 °	12.33 °	12.51 %	3.66	64.99

11 Valores Atípicos

Aplicando el criterio analítico de Tukey sobre los cuartiles globales, esta sección identifica las coordenadas espaciales extremas. Estos puntos representan anomalías geográficas severas (siniestros ubicados a distancias marginales hacia el Este o el Oeste).

# --- ANÁLISIS GLOBAL ---
variable_global <- variable_interes
stats_outliers_global <- boxplot.stats(variable_global)$out
num_outliers_global <- length(stats_outliers_global)
minimooutliers_global <- if(num_outliers_global > 0) min(stats_outliers_global) else NA
maximooutliers_global <- if(num_outliers_global > 0) max(stats_outliers_global) else NA

cat("\n--- Análisis de Outliers (GLOBAL) ---\n")

## 
## --- Análisis de Outliers (GLOBAL) ---

cat("Número de valores atípicos:", num_outliers_global, "\n")

## Número de valores atípicos: 355

cat("Ubicación más extrema (Min):", if(!is.na(minimooutliers_global)) paste(format(minimooutliers_global, big.mark=","), "°") else "Ninguno", "\n")

## Ubicación más extrema (Min): -158.0999 °

cat("Ubicación más extrema (Max):", if(!is.na(maximooutliers_global)) paste(format(maximooutliers_global, big.mark=","), "°") else "Ninguno", "\n")

## Ubicación más extrema (Max): 104.2634 °

12 Conclusión

El análisis de la variable Accident.Longitude revela la estructura geoespacial de la siniestralidad. Al comparar la distribución global con el rango principal, podemos identificar si los accidentes se concentran en una franja longitudinal específica (zona central de operaciones) o si están dispersos geográficamente. Los valores de curtosis y asimetría indicarán si hay una tendencia a la aglomeración en ciertas coordenadas (moda) o si la distribución es uniforme a lo largo del territorio.

Variable Cuantitativa Continua: longitud del Accidente

Variable Original: Accident.Longitude

Brandon