Variable Cuantitativa Continua: longitud del Accidente

Variable Original: Accident.Longitude

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.Latitude

2 Justificación de la variable

La variable de estudio Longitud del Accidente se define como una variable cuantitativa continua para el análisis de ubicación geoespacial. Se aplica un Modelo de Mezcla Normal (Bimodal) de manera segmentada sobre los datos, asumiendo que la incidencia de accidentes no es uniforme, sino que presenta una doble concentración en zonas geográficas específicas (Oeste y Este) separadas por una brecha espacial. Este modelo permite ajustar curvas de densidad de probabilidad independientes para cada agrupación, facilitando la identificación de zonas críticas de alta densidad y permitiendo realizar pruebas de bondad de ajuste para validar la precisión de la distribución gaussiana frente a las coordenadas reales de los siniestros observados.

3 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

4 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

5 Cantidad Global de Accidentes (ni)

En esta sección se analiza la dispersión total de los accidentes a lo largo de todas las longitudes registradas.

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)

#Análisis del Rango Principal (Cantidad Local) Se ha filtrado la información para enfocarse en el “Rango Principal”, definido aquí como la zona donde se concentran la mayoría de los datos (entre el percentil 10% y el 90%), eliminando los valores atípicos extremos geográficos para observar mejor la distribución central.

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 del rango principal

Estas secciones transforman las frecuencias absolutas en porcentajes para entender el peso de cada franja longitudinal dentro de la zona de mayor accidentalidad.

# 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 del rango principal

Visualiza la probabilidad de que un accidente caiga dentro de un rango específico de longitud.

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

La ojiva combinada representa la acumulación de accidentes a medida que nos desplazamos de Oeste a Este (ascendente) o viceversa.

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

El diagrama de cajas visualiza la dispersión central de la longitud de los accidentes dentro del rango principal seleccionado.

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 local

Esta tabla proporciona los valores exactos (media, mediana, etc.) para el rango geográfico principal seleccionado.

library(e1071)
library(knitr)

variable_analisis <- variable_box # Local
n <- length(variable_analisis)

# Recálculo básico para la media agrupada local
k_local <- 6
min_val <- min(variable_analisis)
max_val <- max(variable_analisis)
A_local <- (max_val - min_val) / k_local
Li <- seq(min_val, max_val - A_local, length.out = k_local)
Ls <- Li + A_local
MC <- (Li + Ls) / 2

ni <- numeric(k_local)
for(i in 1:k_local){
  if(i == k_local) ni[i] <- sum(variable_analisis >= Li[i] & variable_analisis <= max_val)
  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), nsmall=2), " ; ", format(round(ic_superior, 2), nsmall=2), "]")
ri <- min(variable_analisis)
rs <- max(variable_analisis)
mediana <- median(variable_analisis)

# Moda simple
t <- table(round(variable_analisis, 1)) # Redondeamos para encontrar moda útil
Mo <- as.numeric(names(t)[which.max(t)])

cv <- (desviacion_estandar / abs(media_agrupada)) * 100
As <- skewness(variable_analisis)
K <- kurtosis(variable_analisis)

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

colnames(Tabla_local) <- c("Variable","Min","Max","Media (IC 95%)","Mediana","Moda","Desv. S","CV","As","K")
kable(Tabla_local, format = "markdown", caption = "Tabla 1: Indicadores Locales (Rango Principal).")

Tabla 1: Indicadores Locales (Rango Principal).

Variable	Min	Max	Media (IC 95%)	Mediana	Moda	Desv. S	CV	As	K
Longitud (Local)	-104.9132 °	-83.45348 °	[-95.89 ; -95.48]	-95.489 °	-95.2 °	4.927 °	5.15 %	0.26	-0.27

11 Tabla estadísticos global

Sirve como línea base para comparar la ubicación promedio general contra el rango principal.

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

#Valores Atípicos global Detectamos ubicaciones extremadamente alejadas del centro de los datos (siniestros muy al Oeste o muy al Este).

# --- 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 Valores atípicos del rango principal

En forma de rango principal se diferencio a medida con el anterior modelo.

# --- ANÁLISIS LOCAL ---
variable_local <- variable_box
stats_outliers_local <- boxplot.stats(variable_local)$out
num_outliers_local <- length(stats_outliers_local)

cat("\n--- Análisis de Outliers (LOCAL - Rango Principal) ---\n")

## 
## --- Análisis de Outliers (LOCAL - Rango Principal) ---

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

## Número de valores atípicos internos: 82

13 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.