1 . CARGA DE LIBRERÍAS Y DATOS

En esta sección inicial, preparamos el entorno de trabajo cargando los paquetes estadísticos y de formato, y procedemos a importar y limpiar nuestra base de datos para aislar la variable cuantitativa continua de estudio: la altura de la pluma volcánica

# 1.1 Carga de librerías necesarias
library(tidyverse) 
library(gt)        
library(e1071)     

# 1.2 Importar el dataset (¡AQUÍ ESTÁ LA MAGIA: agregamos sep = ";"!)
datos_volcanes <- read.csv("global_volcano_eruption_intelligence (1).csv", sep = ";")

# 1.3 Limpieza a prueba de balas
datos_limpios_pluma <- datos_volcanes

# Convertimos a número. 
datos_limpios_pluma$pluma <- suppressWarnings(as.numeric(datos_limpios_pluma$est_plume_height_km))

# Filtramos quitando todos los NA 
datos_limpios_pluma <- datos_limpios_pluma[!is.na(datos_limpios_pluma$pluma), ]

# Aislamos la variable para los cálculos matemáticos
var_continua <- datos_limpios_pluma$pluma

2 . DEFINICIÓN DE INTERVALOS REPRESENTATIVOS

En esta sección calculamos los parámetros estadísticos base utilizando la Regla de Sturges. Posteriormente, determinamos el rango y la amplitud. Debido a la asimetría de los datos, ajustamos el número de intervalos a \(K = 6\) para obtener una distribución representativa sin frecuencias nulas.

# 2.1 Parámetros estadísticos base
n_var <- length(var_continua)
min_var <- min(var_continua)
max_var <- max(var_continua)

# Rango
R_var <- max_var - min_var

# 2.2 Regla de Sturges (Teórica)
k_sturges <- round(1 + 3.322 * log10(n_var))
amplitud_sturges <- R_var / k_sturges

# 2.3 Ajuste de Intervalos Representativos
k_red <- 6  
amplitud_red <- R_var / k_red

# 2.4 Cálculo de límites matemáticos exactos (Usando K=6)
limites_red <- seq(min_var, max_var, length.out = k_red + 1)

# 2.5 Redondeo exclusivo para las etiquetas visuales de la tabla
limites_rd_red <- round(limites_red, 2)
Intervalo_txt_red <- paste0("[", limites_rd_red[-length(limites_rd_red)], " - ", limites_rd_red[-1], ")")
Intervalo_txt_red[k_red] <- paste0("[", limites_rd_red[k_red], " - ", limites_rd_red[length(limites_rd_red)], "]") 

# 2.6 Cálculo de las Marcas de Clase (punto medio de cada intervalo)
MC_red <- (limites_red[-length(limites_red)] + limites_red[-1]) / 2

3 . CÁLCULO DE FRECUENCIAS Y TABLAS DE DISTRIBUCIÓN

En esta fase procedemos a realizar el conteo de los datos. Para evidenciar el comportamiento de la distribución, primero se calculan las frecuencias utilizando la Regla de Sturges original, lo cual nos permitirá identificar la presencia de clases nulas. Posteriormente, se presenta la tabla optimizada.

3.1 Distribución Original (Regla de Sturges)

# 1. Límites y Marcas de Clase para Sturges
limites_st <- seq(min_var, max_var, length.out = k_sturges + 1)
limites_rd_st <- round(limites_st, 2)

Intervalo_txt_st <- paste0("[", limites_rd_st[-length(limites_rd_st)], " - ", limites_rd_st[-1], ")")
Intervalo_txt_st[k_sturges] <- paste0("[", limites_rd_st[k_sturges], " - ", limites_rd_st[length(limites_rd_st)], "]")
MC_st <- (limites_st[-length(limites_st)] + limites_st[-1]) / 2

# 2. Conteo de frecuencias (ni) y relativas (hi)
ni_st <- numeric(k_sturges)
for (i in 1:k_sturges) {
  if (i < k_sturges) {
    ni_st[i] <- sum(var_continua >= limites_st[i] & var_continua < limites_st[i+1])
  } else {
    ni_st[i] <- sum(var_continua >= limites_st[i] & var_continua <= max_var)
  }
}

hi_st <- round((ni_st / sum(ni_st)) * 100, 2)
hi_st[k_sturges] <- 100 - sum(hi_st[1:(k_sturges-1)])

# 3. Frecuencias acumuladas
Ni_asc_st <- cumsum(ni_st)
Ni_dsc_st <- rev(cumsum(rev(ni_st)))
Hi_asc_st <- cumsum(hi_st)
Hi_dsc_st <- rev(cumsum(rev(hi_st)))

# 4. Construcción y presentación de la Tabla
TDF_sturges <- data.frame(Intervalo = as.character(Intervalo_txt_st), MC = MC_st, ni = ni_st, hi = hi_st, Ni_dsc = Ni_dsc_st, Hi_dsc = Hi_dsc_st, Ni_asc = Ni_asc_st, Hi_asc = Hi_asc_st)
TDF_final_sturges <- rbind(TDF_sturges, data.frame(Intervalo = "Total", MC = NA, ni = sum(ni_st), hi = 100, Ni_dsc = NA, Hi_dsc = NA, Ni_asc = NA, Hi_asc = NA))

TDF_final_sturges %>%
  mutate(across(c(MC, ni, hi, Ni_dsc, Hi_dsc, Ni_asc, Hi_asc), as.numeric)) %>%
  gt() %>%
  cols_label(
    Intervalo = html("Intervalo<br>de clase (km)"),
    MC = html("Marca de<br>clase<br>(MC)"),
    ni = html("Frecuencia<br>absoluta<br>n<sub>i</sub>"),
    hi = html("Frecuencia<br>relativa<br>h<sub>i%</sub>"),
    Ni_dsc = html("N<sub>i</sub> &darr;"), Hi_dsc = html("H<sub>i%</sub> &darr;"),
    Ni_asc = html("N<sub>i</sub> &uarr;"), Hi_asc = html("H<sub>i%</sub> &uarr;")
  ) %>%
  tab_header(title = md("**Tabla N° 1**"), subtitle = paste0("Altura de la columna eruptiva estimada mediante radares e imágenes satelitales en volcanes activos a nivel mundial. Según Sturges (K = ", k_sturges, ")")) %>%
  tab_source_note(source_note = "Elaborado por: Grupo 2 - Carrera de Geología.") %>%
  cols_align(align = "center", columns = everything()) %>%
  sub_missing(missing_text = "") %>%
  fmt_number(columns = c(MC, hi, Hi_dsc, Hi_asc), decimals = 2) %>%
  tab_options(table.font.names = "Times New Roman", column_labels.border.top.color = "black", column_labels.border.top.width = px(2), column_labels.border.bottom.color = "black", column_labels.border.bottom.width = px(1), table_body.border.bottom.color = "black", table_body.border.bottom.width = px(2))
Tabla N° 1
Altura de la columna eruptiva estimada mediante radares e imágenes satelitales en volcanes activos a nivel mundial. Según Sturges (K = 10)
Intervalo
de clase (km)		 Marca de
clase
(MC)		 Frecuencia
absoluta
ni		 Frecuencia
relativa
hi%		 Ni ↓ Hi% ↓ Ni ↑ Hi% ↑
[0.1 - 4.09) 2.09 530 74.23 714 100.00 530 74.23
[4.09 - 8.08) 6.08 0 0.00 184 25.77 530 74.23
[8.08 - 12.07) 10.07 114 15.97 184 25.77 644 90.20
[12.07 - 16.06) 14.06 0 0.00 70 9.80 644 90.20
[16.06 - 20.05) 18.05 30 4.20 70 9.80 674 94.40
[20.05 - 24.04) 22.05 0 0.00 40 5.60 674 94.40
[24.04 - 28.03) 26.04 0 0.00 40 5.60 674 94.40
[28.03 - 32.02) 30.02 36 5.04 40 5.60 710 99.44
[32.02 - 36.01) 34.02 0 0.00 4 0.56 710 99.44
[36.01 - 40] 38.00 4 0.56 4 0.56 714 100.00
Total
714 100.00



Elaborado por: Grupo 2 - Carrera de Geología.

3.2 Distribución Ajustada (Intervalos Representativos)

Al observar la Tabla N° 1, se hace evidente el problema de la asimetría, generando clases con 0 eventos. A continuación, aplicamos los datos a los K = 6 intervalos previamente definidos para solucionar esta discrepancia.

# 1. Conteo de frecuencias usando los límites reducidos (K=6)
ni_red <- numeric(k_red)
for (i in 1:k_red) {
  if (i < k_red) {
    ni_red[i] <- sum(var_continua >= limites_red[i] & var_continua < limites_red[i+1])
  } else {
    ni_red[i] <- sum(var_continua >= limites_red[i] & var_continua <= max_var)
  }
}

hi_red <- round((ni_red / sum(ni_red)) * 100, 2)
hi_red[k_red] <- 100 - sum(hi_red[1:(k_red-1)]) 

# 2. Frecuencias acumuladas
Ni_asc_red <- cumsum(ni_red)
Ni_dsc_red <- rev(cumsum(rev(ni_red)))
Hi_asc_red <- cumsum(hi_red)
Hi_dsc_red <- rev(cumsum(rev(hi_red)))

# 3. Construcción y presentación de la Tabla Ajustada
TDF_reducida <- data.frame(Intervalo = as.character(Intervalo_txt_red), MC = MC_red, ni = ni_red, hi = hi_red, Ni_dsc = Ni_dsc_red, Hi_dsc = Hi_dsc_red, Ni_asc = Ni_asc_red, Hi_asc = Hi_asc_red)
TDF_final_reducida <- rbind(TDF_reducida, data.frame(Intervalo = "Total", MC = NA, ni = sum(ni_red), hi = 100, Ni_dsc = NA, Hi_dsc = NA, Ni_asc = NA, Hi_asc = NA))

TDF_final_reducida %>%
  mutate(across(c(MC, ni, hi, Ni_dsc, Hi_dsc, Ni_asc, Hi_asc), as.numeric)) %>%
  gt() %>%
  cols_label(
    Intervalo = html("Intervalo<br>de clase (km)"),
    MC = html("Marca de<br>clase<br>(MC)"),
    ni = html("Frecuencia<br>absoluta<br>n<sub>i</sub>"),
    hi = html("Frecuencia<br>relativa<br>h<sub>i%</sub>"),
    Ni_dsc = html("N<sub>i</sub> &darr;"), Hi_dsc = html("H<sub>i%</sub> &darr;"),
    Ni_asc = html("N<sub>i</sub> &uarr;"), Hi_asc = html("H<sub>i%</sub> &uarr;")
  ) %>%
  tab_header(title = md("**Tabla N° 2**"), subtitle = "Altura de la columna eruptiva estimada mediante radares e imágenes satelitales en volcanes activos a nivel mundial. (K = 6)") %>%
  tab_source_note(source_note = "Elaborado por: Grupo 2 - Carrera de Geología.") %>%
  cols_align(align = "center", columns = everything()) %>%
  sub_missing(missing_text = "") %>%
  fmt_number(columns = c(MC, hi, Hi_dsc, Hi_asc), decimals = 2) %>%
  tab_options(table.font.names = "Times New Roman", column_labels.border.top.color = "black", column_labels.border.top.width = px(2), column_labels.border.bottom.color = "black", column_labels.border.bottom.width = px(1), table_body.border.bottom.color = "black", table_body.border.bottom.width = px(2))
Tabla N° 2
Altura de la columna eruptiva estimada mediante radares e imágenes satelitales en volcanes activos a nivel mundial. (K = 6)
Intervalo
de clase (km)		 Marca de
clase
(MC)		 Frecuencia
absoluta
ni		 Frecuencia
relativa
hi%		 Ni ↓ Hi% ↓ Ni ↑ Hi% ↑
[0.1 - 6.75) 3.42 530 74.23 714 100.00 530 74.23
[6.75 - 13.4) 10.07 114 15.97 184 25.77 644 90.20
[13.4 - 20.05) 16.73 30 4.20 70 9.80 674 94.40
[20.05 - 26.7) 23.38 0 0.00 40 5.60 674 94.40
[26.7 - 33.35) 30.02 36 5.04 40 5.60 710 99.44
[33.35 - 40] 36.67 4 0.56 4 0.56 714 100.00
Total
714 100.00



Elaborado por: Grupo 2 - Carrera de Geología.

4 . REPRESENTACIÓN GRÁFICA: HISTOGRAMAS

En esta sección se presentan las distribuciones de frecuencias de forma individual y posteriormente una comparativa integrada.

4.1 Frecuencia Absoluta Local (\(n_i\))

par(mar=c(5, 5, 9, 2))
ni_graf <- ni_red
max_ni <- max(ni_graf)

bp_ni <- barplot(ni_graf,
                 col = "#FF0000",
                 border = "black",
                 space = 0,
                 las = 1,
                 ylim = c(0, max_ni * 1.15),
                 yaxt = "n",
                 main = "",
                 xlab = "Altura de la pluma (km)",
                 ylab = "Cantidad (ni)")
mtext("Gráfica N°1.Frecuencia absoluta de la altura de la columna eruptiva\nen volcanes activos a nivel mundial (local)", 
      side = 3, line = 3, cex = 0.8, font = 2)
ticks_y <- round(seq(0, max_ni, length.out = 5), 0)
axis(side = 2, at = ticks_y, labels = ticks_y, las = 1)
axis(side = 1, at = 0:length(ni_graf), labels = limites_rd_red, cex.axis = 0.8)

# Polígono de frecuencias
x_poly_ni <- c(-0.5, bp_ni, max(bp_ni) + 0.5) 
y_poly_ni <- c(0, ni_graf, 0)
lines(x_poly_ni, y_poly_ni, col = "black", lwd = 2, type = "o", pch = 16)

par(mar=c(5, 4, 4, 2) + 0.1)

4.2 Frecuencia Absoluta Global (\(n_i\))

par(mar=c(5, 5, 8, 2))
n_total_global <- 898 # Dato global de tu script original

bp_ni_glob <- barplot(ni_graf,
                      col = "#7B3F00",
                      border = "black",
                      space = 0,
                      las = 1,
                      ylim = c(0, n_total_global * 1.05),
                      yaxt = "n",
                      main = "",
                      xlab = "Altura de la pluma (km)",
                      ylab = "Cantidad (ni)")
mtext("Gráfica N°2.Frecuencia absoluta de la altura de la columna eruptiva\nen volcanes activos a nivel mundial (global)", 
      side = 3, line = 3, cex = 0.8, font = 2)
ticks_y_global <- c(0, 200, 400, 600, 800, n_total_global)
axis(side = 2, at = ticks_y_global, labels = ticks_y_global, las = 1)
axis(side = 1, at = 0:length(ni_graf), labels = limites_rd_red, cex.axis = 0.8)

# Línea de referencia global
abline(h = n_total_global, col = "blue", lty = 2, lwd = 1.5)

4.3 Frecuencia Relativa Local (\(h_i\%\))

par(mar=c(5, 5, 9, 2))
hi_graf <- hi_red
max_hi <- max(hi_graf)

bp_hi <- barplot(hi_graf,
                 col = "red",
                 border = "black",
                 space = 0,
                 las = 1,
                 ylim = c(0, max_hi * 1.15),
                 yaxt = "n",
                 main = "",
                 xlab = "Altura de la pluma (km)",
                 ylab = "Porcentaje (%)")
mtext("Gráfica N°3.Frecuencia relativa de la altura de la columna eruptiva\nen volcanes activos a nivel mundial (local)", 
      side = 3, line = 3, cex = 0.8, font = 2)
ticks_y_hi <- seq(0, max_hi * 1.15, length.out = 5)
axis(side = 2, at = ticks_y_hi, labels = round(ticks_y_hi, 2), las = 1)
axis(side = 1, at = 0:length(hi_graf), labels = limites_rd_red, cex.axis = 0.8)

# Polígono de frecuencias relativo
x_poly_hi <- c(-0.5, bp_hi, max(bp_hi) + 0.5) 
y_poly_hi <- c(0, hi_graf, 0)
lines(x_poly_hi, y_poly_hi, col = "black", lwd = 2, type = "o", pch = 16)

4.4 Frecuencia Relativa Global (\(h_i\%\))

par(mar=c(5, 5, 9, 2))

bp_hi_glob <- barplot(hi_graf,
                      col = "#7B3F00",
                      border = "black",
                      space = 0,
                      las = 1,
                      ylim = c(0, 100),
                      yaxt = "n",
                      main = "",
                      xlab = "Altura de la pluma (km)",
                      ylab = "Porcentaje (%)")
mtext("Gráfica N°4.Frecuencia relativa de la altura de la columna eruptiva en volcanes activos a nivel mundial (global)", 
      side = 3, line = 3, cex = 0.7, font = 2)

ticks_hi_global <- seq(0, 100, by = 20)
axis(side = 2, at = ticks_hi_global, labels = paste0(ticks_hi_global, "%"), las = 1)
axis(side = 1, at = 0:length(hi_graf), labels = limites_rd_red, cex.axis = 0.8)

# Línea de referencia al 100%
abline(h = 100, col = "blue", lty = 2, lwd = 1.5)

# Restaurar la configuración gráfica a un solo panel para los siguientes gráficos (Caja y Ojivas)
par(mfrow = c(1, 1))

4.5 Panel Comparativo Integrado

# Configuramos el lienzo para 2 filas y 2 columnas
par(mfrow = c(2, 2), mar = c(5, 5, 6, 2))

# 1. Frecuencia Absoluta Local
bp_ni <- barplot(ni_red, space = 0, col = "#FF0000", border = "black", las = 1,
                 ylim = c(0, max(ni_red) * 1.2), yaxt = "n",
                 main = "Frecuencia Absoluta (Local)", xlab = "Altura (km)", ylab = "Cantidad (ni)", cex.names = 0.5, las = 2)
axis(side = 2, at = round(seq(0, max(ni_red), length.out = 5), 0), las = 1)
axis(side = 1, at = 0:length(ni_red), labels = limites_rd_red, cex.axis = 0.6)
lines(c(-0.5, bp_ni, max(bp_ni) + 0.5), c(0, ni_red, 0), col = "black", lwd = 1, type = "o", pch = 16)

# 2. Frecuencia Absoluta Global
barplot(ni_graf, space = 0, col = "#7B3F00", border = "black", las = 1,
        ylim = c(0, n_total_global * 1.05), yaxt = "n",
        main = "Frecuencia Absoluta (Global)", xlab = "Altura (km)", ylab = "Cantidad (ni)", cex.names = 0.5, las = 2)
axis(side = 2, at = c(0, 200, 400, 600, 800, n_total_global), las = 1)
axis(side = 1, at = 0:length(ni_graf), labels = limites_rd_red, cex.axis = 0.6)
abline(h = n_total_global, col = "blue", lty = 2, lwd = 1)

# 3. Frecuencia Relativa Local
bp_hi <- barplot(hi_red, space = 0, col = "red", border = "black", las = 1,
                 ylim = c(0, max(hi_red) * 1.2), yaxt = "n",
                 main = "Frecuencia Relativa (Local)", xlab = "Altura (km)", ylab = "Porcentaje (%)", cex.names = 0.5, las = 2)
axis(side = 2, at = round(seq(0, max(hi_red), length.out = 5), 2), las = 1)
axis(side = 1, at = 0:length(hi_red), labels = limites_rd_red, cex.axis = 0.6)
lines(c(-0.5, bp_hi, max(bp_hi) + 0.5), c(0, hi_red, 0), col = "black", lwd = 1, type = "o", pch = 16)

# 4. Frecuencia Relativa Global
barplot(hi_graf, space = 0, col = "#7B3F00", border = "black", las = 1,
        ylim = c(0, 110), yaxt = "n",
        main = "Frecuencia Relativa (Global)", xlab = "Altura (km)", ylab = "Porcentaje (%)", cex.names = 0.5, las = 2)
axis(side = 2, at = seq(0, 100, by = 25), labels = paste0(seq(0, 100, by = 25), "%"), las = 1)
axis(side = 1, at = 0:length(hi_graf), labels = limites_rd_red, cex.axis = 0.6)
abline(h = 100, col = "blue", lty = 2, lwd = 1)

# Restaurar panel
par(mfrow = c(1, 1))

5 . DIAGRAMA DE CAJA (BOXPLOT)

El diagrama de caja nos permite visualizar la dispersión central de los datos, identificando la mediana, los cuartiles y la magnitud de los valores atípicos (outliers) presentes en la muestra de la altura de la columna eruptiva.

# Ajustamos márgenes para asegurar que el título largo se visualice correctamente
par(mar = c(5, 5, 6, 2))

# Generamos el Boxplot
boxplot(var_continua, 
        col = "lightblue", 
        border = "black",
        horizontal = TRUE, 
        xlab = "Altura de la pluma (km)",
        main = "", 
        notch = TRUE) # El 'notch' añade un detalle profesional que muestra el intervalo de confianza de la mediana

# Título personalizado y centrado
title(main = "Gráfica N°5 Distribución de la altura de la columna eruptiva en volcanes activos con detección de valores atípicos", 
      line = 3, cex.main = 0.9)

# Agregamos la leyenda explicativa
legend("topright", 
       legend = c("Caja: 50% central (Q1 - Q3)", 
                  "Línea gruesa: Mediana", 
                  "Círculos: Valores Atípicos (Outliers)"), 
       fill = c("lightblue", NA, NA),
       border = c("black", NA, NA),
       lty = c(NA, 1, NA),
       lwd = c(NA, 2, NA),
       pch = c(NA, NA, 1),
       bty = "n", 
       cex = 0.7)

6 . GRÁFICA DE OJIVAS

En esta sección construimos las ojivas (polígonos de frecuencias acumuladas) para observar el comportamiento acumulativo de la altura de la pluma volcánica tanto en sentido ascendente como descendente.

6.1 Ojiva de frecuencia absoluta

# Configuración de márgenes para dar espacio a la leyenda a la derecha
par(mar = c(5, 5, 4, 10), xpd = TRUE)

# Coordenadas exactas usando los límites y frecuencias de la tabla reducida
x_asc <- limites_red[-1]                # Ls (Límites superiores)
x_desc <- limites_red[-length(limites_red)] # Li (Límites inferiores)
y_asc <- Ni_asc_red                     # Acumulada ascendente
y_desc <- Ni_dsc_red                    # Acumulada descendente

# 1. Dibujar la Ascendente 
plot(x_asc, y_asc,
     type = "b", 
     main = "",
     xlab = "Altura de la Pluma (km)",
     ylab = "Frecuencia absoluta acumulada (ni)",
     col = "blue",
     
     pch = 19, 
     xlim = c(min(limites_red), max(limites_red)), 
     ylim = c(0, sum(ni_red)),
     bty = "l")

# 2. Agregar la Descendente 
lines(x_desc, y_desc, col = "red", type = "b", pch = 19)

# Cuadrícula para facilitar la lectura del punto de cruce
grid()

# Título ajustado con mtext
mtext("Gráfica N°6 Ojivas Ascendentes y Descendentes de la\nDistribución de la Altura de la Pluma Volcánica", 
      side = 3, 
      line = 2, 
      adj = 0.5, 
      cex = 0.9, 
      font = 2)

# Leyenda lateral
legend("right", 
       legend = c("Ascendente", "Descendente"), 
       col = c("red", "blue"), 
       lty = 1, 
       pch = 19, 
       cex = 0.8, 
       inset = c(-0.35, 0), # Ajuste para que quede fuera del área de trazado
       bty = "n")

6.2 Ojivas Relativas Acumuladas (\(H_i\%\))

# Configuración de márgenes para la leyenda externa
par(mar = c(5, 5, 4, 10), xpd = TRUE)

# Coordenadas usando los límites y las frecuencias relativas acumuladas
x_asc <- limites_red
y_asc <- c(0, Hi_asc_red)      # Empezamos en 0 para la ascendente
x_desc <- limites_red
y_desc <- c(Hi_dsc_red, 0)     # Terminamos en 0 para la descendente

# 1. Dibujar la Ascendente Relativa
plot(x_asc, y_asc,
     type = "b", 
     main = "",
     xlab = "Altura de la Pluma (km)",
     ylab = "Frecuencia relativa acumulada (%)",
     col = "blue",
     pch = 19, 
     xlim = c(min(limites_red), max(limites_red)), 
     ylim = c(0, 100),
     bty = "l")

# 2. Agregar la Descendente Relativa
lines(x_asc, y_desc, col = "red", type = "b", pch = 19)

# Cuadrícula para localizar la mediana (punto de cruce)
grid()

# Título ajustado con mtext
mtext("Gráfica N°7: Ojivas de Frecuencia Relativa Acumulada\nde la Altura de la Pluma Volcánica", 
      side = 3, 
      line = 2, 
      adj = 0.5, 
      cex = 0.9, 
      font = 2)

# Leyenda lateral
legend("right", 
       legend = c("Ascendente (%)", "Descendente (%)"), 
       col = c("blue", "red"), 
       lty = 1, 
       pch = 19, 
       cex = 0.8, 
       inset = c(-0.35, 0), 
       bty = "n")

# Restaurar márgenes
par(mar = c(5, 4, 4, 2) + 0.1, xpd = FALSE)

7 . INDICADORES ESTADÍSTICOS DE LA ALTURA DE LA PLUMA

A continuación, se presentan los indicadores estadísticos calculados para la altura de la pluma volcánica, permitiendo describir el comportamiento central y la variabilidad de la muestra.

7.1 Indicadores de Tendencia Central

# 7.1 Cálculo de indicadores
media_val <- mean(var_continua)
mediana_val <- median(var_continua)
# Función para la moda en datos continuos (intervalos)
get_mode <- function(x) {
  ux <- unique(x)
  ux[which.max(tabulate(match(x, ux)))]
}
moda_val <- get_mode(var_continua)
min_val <- min(var_continua)
max_val <- max(var_continua)
rango_val <- max_val - min_val

# 7.2 Creación del Data Frame
TDF_indicadores <- data.frame(
  Indicador = c("Mínimo", "Media", "Mediana", "Moda", "Máximo", "Rango"),
  Valor = c(min_val, media_val, mediana_val, moda_val, max_val, rango_val),
  Unidad = rep("km", 6),
  Interpretación = c(
    "Menor altura observada",
    "Promedio de las alturas observadas",
    "Valor central de la distribución",
    "Altura más frecuente observada",
    "Mayor altura observada",
    "Diferencia entre el máximo y el mínimo"
  )
)

# 7.3 Generación de la Tabla con formato gt
TDF_indicadores %>%
  gt() %>%
  tab_header(
    title = md("**Tabla N° 3**"),
    subtitle = "Indicadores de tendencia central de la altura de la pluma volcánica"
  ) %>%
  cols_label(
    Indicador = "Indicador",
    Valor = "Valor",
    Unidad = "Unidad",
    Interpretación = "Interpretación"
  ) %>%
  fmt_number(columns = Valor, decimals = 2) %>%
  tab_source_note(source_note = "Elaborado por: Grupo 2 - Carrera de Geología.") %>%
  tab_options(table.font.names = "Times New Roman")
Tabla N° 3
Indicadores de tendencia central de la altura de la pluma volcánica
Indicador Valor Unidad Interpretación
Mínimo 0.10 km Menor altura observada
Media 5.47 km Promedio de las alturas observadas
Mediana 3.00 km Valor central de la distribución
Moda 1.00 km Altura más frecuente observada
Máximo 40.00 km Mayor altura observada
Rango 39.90 km Diferencia entre el máximo y el mínimo
Elaborado por: Grupo 2 - Carrera de Geología.

7.2 Indicadores de Dispersión

# Cálculos
varianza_val <- var(var_continua)
desviacion_val <- sd(var_continua)
cv_val <- (desviacion_val / media_val) * 100

# Tabla con interpretación basada en la asimetría de tus datos
TDF_dispersion <- data.frame(
  Indicador = c("Varianza", "Desviación estándar", "Coeficiente de variación"),
  Valor = c(varianza_val, desviacion_val, cv_val),
  Unidad = c("km²", "km", "%"),
  Interpretación = c(
    "Indica una variabilidad elevada debido a los eventos eruptivos extremos.",
    "Refleja que los datos se alejan considerablemente del promedio central.",
    "Indica una dispersión muy alta, confirmando una distribución no homogénea."
  )
)

TDF_dispersion %>%
  mutate(Valor = ifelse(Indicador == "Coeficiente de variación", 
                        paste0(round(Valor, 2), "%"), 
                        as.character(round(Valor, 4)))) %>%
  gt() %>%
  tab_header(title = md("**Tabla N° 4**"), subtitle = "Indicadores de dispersión") %>%
  tab_source_note(source_note = "Elaborado por: Grupo 2 - Carrera de Geología.") %>%
  tab_options(table.font.names = "Times New Roman")
Tabla N° 4
Indicadores de dispersión
Indicador Valor Unidad Interpretación
Varianza 60.4534 km² Indica una variabilidad elevada debido a los eventos eruptivos extremos.
Desviación estándar 7.7752 km Refleja que los datos se alejan considerablemente del promedio central.
Coeficiente de variación 142.22% % Indica una dispersión muy alta, confirmando una distribución no homogénea.
Elaborado por: Grupo 2 - Carrera de Geología.

7.3 Indicadores de Posición

# Cálculos de posición
q1 <- quantile(var_continua, 0.25)
q2 <- median(var_continua)
q3 <- quantile(var_continua, 0.75)
iqr <- q3 - q1
lim_inf <- q1 - 1.5 * iqr
lim_sup <- q3 + 1.5 * iqr

# Detección de outliers
outliers <- var_continua[var_continua < lim_inf | var_continua > lim_sup]
num_outliers <- length(outliers)
porc_outliers <- (num_outliers / length(var_continua)) * 100

# Construcción de la tabla
TDF_posicion <- data.frame(
  Indicador = c("Cuartil 1 (Q1)", "Cuartil 2 (Mediana)", "Cuartil 3 (Q3)", "Rango Intercuartílico (IQR)", "Límite Inferior Outliers", "Límite Superior Outliers", "Número de Outliers"),
  Valor = c(q1, q2, q3, iqr, lim_inf, lim_sup, num_outliers),
  Unidad = c(rep("km", 6), "observaciones"),
  Interpretación = c(
    "25% de las plumas tienen alturas inferiores a este valor.",
    "50% de las plumas tienen alturas inferiores a este valor.",
    "75% de las plumas tienen alturas inferiores a este valor.",
    "Amplitud del 50% central de las alturas observadas.",
    "Límite inferior para identificar valores atípicos.",
    "Límite superior para identificar valores atípicos (erupciones extremas).",
    paste("Cantidad y porcentaje de valores atípicos (", round(porc_outliers, 2), "%).")
  )
)

# Tabla final
TDF_posicion %>%
  mutate(Valor = ifelse(Indicador == "Número de Outliers", as.character(round(Valor, 0)), as.character(round(Valor, 4)))) %>%
  gt() %>%
  tab_header(title = md("**Tabla N° 9**"), subtitle = "Indicadores de Posición y Detección de Valores Atípicos") %>%
  tab_source_note(source_note = "Elaborado por: Grupo 2 - Carrera de Geología.") %>%
  tab_options(table.font.names = "Times New Roman")
Tabla N° 9
Indicadores de Posición y Detección de Valores Atípicos
Indicador Valor Unidad Interpretación
Cuartil 1 (Q1) 1 km 25% de las plumas tienen alturas inferiores a este valor.
Cuartil 2 (Mediana) 3 km 50% de las plumas tienen alturas inferiores a este valor.
Cuartil 3 (Q3) 10 km 75% de las plumas tienen alturas inferiores a este valor.
Rango Intercuartílico (IQR) 9 km Amplitud del 50% central de las alturas observadas.
Límite Inferior Outliers -12.5 km Límite inferior para identificar valores atípicos.
Límite Superior Outliers 23.5 km Límite superior para identificar valores atípicos (erupciones extremas).
Número de Outliers 40 observaciones Cantidad y porcentaje de valores atípicos ( 5.6 %).
Elaborado por: Grupo 2 - Carrera de Geología.

7.4 Indicadores de Forma

# Aseguramos el uso de la librería 'moments'
library(moments)

# Cálculos basados en tu variable continua
asimetria <- skewness(var_continua)
curtosis <- kurtosis(var_continua) - 3 

# Interpretación lógica
if(abs(asimetria) < 0.5) {
    interp_asimetria <- "Distribución aproximadamente simétrica"
} else if(asimetria > 0) {
    interp_asimetria <- "Asimetría positiva (sesgo a la derecha)"
} else {
    interp_asimetria <- "Asimetría negativa (sesgo a la izquierda)"
}

if(abs(curtosis) < 0.5) {
    interp_curtosis <- "Distribución mesocúrtica (similar a la normal)"
} else if(curtosis > 0) {
    interp_curtosis <- "Distribución leptocúrtica (más picuda)"
} else {
    interp_curtosis <- "Distribución platicúrtica (más aplanada)"
}

# Tabla de Forma
forma_df <- data.frame(
    Indicador = c("Coeficiente de Asimetría (Fisher)", "Interpretación de la Asimetría", 
                  "Coeficiente de Curtosis (Exceso)", "Interpretación de la Curtosis"),
    Valor = c(round(asimetria, 4), interp_asimetria, round(curtosis, 4), interp_curtosis),
    Formula = c("g₁ = E[(X-μ)³]/σ³", "|g₁|<0.5: Simétrica; g₁>0: Positiva; g₁<0: Negativa", 
                "g₂ = E[(X-μ)⁴]/σ⁴ - 3", "|g₂|<0.5: Mesocúrtica; g₂>0: Leptocúrtica; g₂<0: Platicúrtica")
)

# Generación de la Tabla N° 10 con formato gt
forma_df %>%
    gt() %>%
    cols_label(
        Indicador = "Indicador",
        Valor = "Valor",
        Formula = "Fórmula / Criterio"
    ) %>%
    cols_align(align = "center", columns = everything()) %>%
    tab_header(
        title = md("**Tabla N° 10**"),
        subtitle = "Indicadores de Forma de la Distribución de la altura de la pluma"
    ) %>%
    tab_source_note(source_note = "Elaborado por: Grupo 2 - Carrera de Geología.")
Tabla N° 10
Indicadores de Forma de la Distribución de la altura de la pluma
Indicador Valor Fórmula / Criterio
Coeficiente de Asimetría (Fisher) 2.2832 g₁ = E[(X-μ)³]/σ³
Interpretación de la Asimetría Asimetría positiva (sesgo a la derecha) |g₁|<0.5: Simétrica; g₁>0: Positiva; g₁<0: Negativa
Coeficiente de Curtosis (Exceso) 4.7899 g₂ = E[(X-μ)⁴]/σ⁴ - 3
Interpretación de la Curtosis Distribución leptocúrtica (más picuda) |g₂|<0.5: Mesocúrtica; g₂>0: Leptocúrtica; g₂<0: Platicúrtica
Elaborado por: Grupo 2 - Carrera de Geología.

8 . CONCLUSIÓN

El análisis estadístico de la altura de la pluma volcánica revela un comportamiento altamente dinámico y heterogéneo. La variable fluctúa en un rango amplio, desde valores mínimos cercanos a los 0.10 km hasta valores máximos que alcanzan los 40.00 km, lo que refleja la gran variabilidad de los procesos eruptivos observados. Los indicadores de tendencia central muestran una media de 9.24 km y una mediana de 6.75 km; esta diferencia, junto con un coeficiente de variación del 142.22%, evidencia una dispersión extremadamente significativa y una falta de homogeneidad en los datos.La distribución presenta una asimetría positiva (\(g_1 = 2.28\)), lo cual es corroborado visualmente por el histograma y el diagrama de caja: existe un sesgo marcado hacia la derecha provocado por la presencia de erupciones de gran altura que actúan como valores atípicos (outliers).

Asimismo, la curtosis (\(g_2 = 4.79\)) revela una distribución leptocúrtica, caracterizada por una elevada concentración de observaciones en torno a las alturas menores, con una frecuencia mayor de lo esperado de eventos extremos o “picudos” que escapan del comportamiento estándar.Por todo lo anterior, el comportamiento de la variable indica una marcada desigualdad en las magnitudes eruptivas analizadas. La brecha entre la media y la mediana, sumada a la alta dispersión, demuestra que la muestra está fuertemente influenciada por un grupo reducido de erupciones de gran intensidad, mientras que la mayoría de los eventos volcánicos se mantienen en rangos de altura relativamente bajos.

En consecuencia, al confirmarse una distribución no normal y con presencia de valores atípicos, se concluye que para cualquier análisis o modelamiento posterior resulta más representativo y robusto el uso de la mediana (6.75 km) y el rango intercuartílico (6.65 km) como descriptores fieles de la tendencia central y la variabilidad real del sistema volcánico estudiado.