ANÁLISIS ESTADÍSTICO

CARGA DE DATOS Y LIBRERÍAS

CARGA DE DATOS

#Limpiar entorno
rm(list = ls())
#Cargar librerías
if (!require("readr")) install.packages("readr")
if (!require("dplyr")) install.packages("dplyr")
if (!require("knitr")) install.packages("knitr")
if (!require("moments")) install.packages("moments")

library(readr)
library(dplyr)
library(knitr)
library(moments)

#Cargar datos
ruta <- "D:/SIO2_with_depth.csv"
datos <- read_csv(ruta)

cat("✓ Datos cargados exitosamente\n")
## ✓ Datos cargados exitosamente
cat("✓ Dimensiones:", nrow(datos), "observaciones\n")
## ✓ Dimensiones: 2500 observaciones

EXTRAER Y LIMPIAR LA VARIABLE

ASIGNACION DE VARIABLES

#Extraer y limpiar variable de Profundidad
PROFUNDIDAD_raw <- datos$TOP_DEPTH_M
PROFUNDIDAD <- as.numeric(PROFUNDIDAD_raw)
PROFUNDIDAD <- PROFUNDIDAD[!is.na(PROFUNDIDAD)]
# Nota: En profundidad, el valor 0 puede ser válido (superficie), 
# pero siguiendo el esqueleto previo de limpieza:
PROFUNDIDAD <- PROFUNDIDAD[PROFUNDIDAD != 0]

n <- length(PROFUNDIDAD)

TABLA DE DISTRIBUCION DE PARAMETROS POR STURGES

# Parámetros de Sturges
min_val <- min(PROFUNDIDAD)
max_val <- max(PROFUNDIDAD)
R <- max_val - min_val

# Regla de Sturges recomendada (redondeo hacia arriba)
k <- ceiling(1 + 3.322 * log10(n))
A <- R / k

parametros_sturges <- data.frame(
  Parámetro = c("Mínimo", 
                "Máximo",
                "Rango (R)", 
                "Número de datos (n)", 
                "Número de intervalos (k)", 
                "Amplitud de clase (A)"),
  Valor = c(round(min_val, 2), round(max_val, 2), round(R, 2), n, k, round(A, 2))
)

kable(parametros_sturges, 
      caption = "Tabla 1: Parámetros de Sturges para TOP_DEPTH_M")
Tabla 1: Parámetros de Sturges para TOP_DEPTH_M
Parámetro Valor
Mínimo 10.47
Máximo 815.94
Rango (R) 805.47
Número de datos (n) 2499.00
Número de intervalos (k) 13.00
Amplitud de clase (A) 61.96 
#Intervalos y frecuencias
# Intervalos y frecuencias (Sturges)
Li <- seq(min(PROFUNDIDAD), max(PROFUNDIDAD) - A, by = A)
Ls <- seq(min(PROFUNDIDAD) + A, max(PROFUNDIDAD), by = A)
MC <- (Li + Ls) / 2

ni <- numeric(k)
for (i in 1:k) {
  if (i == k) {
    ni[i] <- sum(PROFUNDIDAD >= Li[i] & PROFUNDIDAD <= Ls[i])
  } else {
    ni[i] <- sum(PROFUNDIDAD >= Li[i] & PROFUNDIDAD < Ls[i])
  }
}

# 1. Frecuencia relativa exacta en decimales/porcentaje sin redondear
hi_exacto <- (ni / n) * 100

# 2. Frecuencia relativa individual para la columna (a 2 decimales)
hi <- round(hi_exacto, 2)

# 3. Frecuencias acumuladas absolutas
Niasc <- cumsum(ni)
Nidsc <- rev(cumsum(rev(ni)))

# 4. Frecuencias acumuladas relativas basadas en valores exactos (evita error acumulado)
Hiasc <- round(cumsum(hi_exacto), 2)
Hidsc <- round(rev(cumsum(rev(hi_exacto))), 2)

# 5. Ajuste fino de bordes para asegurar 100.00% exacto
Hiasc[length(Hiasc)] <- 100.00
Hidsc[1] <- 100.00

# Construcción del Data Frame
TDF <- data.frame(
  Clase = 1:k,
  Li = round(Li, 2),
  Ls = round(Ls, 2),
  MC = round(MC, 2),
  ni = ni,
  hi = hi,
  Niasc = Niasc,
  Nidsc = Nidsc,
  Hiasc = Hiasc,
  Hidsc = Hidsc
)

kable(
  TDF,
  caption = "Tabla 2: Distribución de frecuencias de Profundidad en metros"
)
Tabla 2: Distribución de frecuencias de Profundidad en metros
Clase Li Ls MC ni hi Niasc Nidsc Hiasc Hidsc
1 10.47 72.43 41.45 72 2.88 72 2499 2.88 100.00
2 72.43 134.39 103.41 193 7.72 265 2427 10.60 97.12
3 134.39 196.35 165.37 210 8.40 475 2234 19.01 89.40
4 196.35 258.31 227.33 303 12.12 778 2024 31.13 80.99
5 258.31 320.27 289.29 362 14.49 1140 1721 45.62 68.87
6 320.27 382.23 351.25 401 16.05 1541 1359 61.66 54.38
7 382.23 444.18 413.21 361 14.45 1902 958 76.11 38.34
8 444.18 506.14 475.16 285 11.40 2187 597 87.52 23.89
9 506.14 568.10 537.12 169 6.76 2356 312 94.28 12.48
10 568.10 630.06 599.08 90 3.60 2446 143 97.88 5.72
11 630.06 692.02 661.04 42 1.68 2488 53 99.56 2.12
12 692.02 753.98 723.00 8 0.32 2496 11 99.88 0.44
13 753.98 815.94 784.96 3 0.12 2499 3 100.00 0.12

Debido a que los límites de clase obtenidos mediante la regla de Sturges presentan números decimales difíciles de interpretar visualmente se simplifico la tabla.

# 1. Definir cortes exactos de Sturges (13 clases)
cortes_sturges <- seq(min(PROFUNDIDAD), min(PROFUNDIDAD) + k * A, length.out = k + 1)
h <- hist(PROFUNDIDAD, breaks = cortes_sturges, plot = FALSE)

# 2. Extraer límites y conteos
Li <- head(h$breaks, -1)
Ls <- tail(h$breaks, -1)
MC <- h$mids
ni <- h$counts
n <- sum(ni)

# 3. Cálculo de Frecuencias Relativas y Acumuladas
hi_exacto <- (ni / n) * 100
hi <- round(hi_exacto, 2)

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

# Acumular sobre exactos para llegar a 100%
Hiasc <- round(cumsum(hi_exacto), 2)
Hidsc <- round(rev(cumsum(rev(hi_exacto))), 2)

# Asegurar el 100.00% exacto en los extremos acumulados
Hiasc[length(Hiasc)] <- 100.00
Hidsc[1] <- 100.00

# 4. Construcción de la Tabla Resumen
TDF_resumen <- data.frame(
  Clase = seq_along(ni),
  Li = round(Li, 0),
  Ls = round(Ls, 0),
  MC = round(MC, 0),
  ni = ni,
  hi = hi,
  Niasc = Niasc,
  Nidsc = Nidsc,
  Hiasc = Hiasc,
  Hidsc = Hidsc
)

kable(
  TDF_resumen,
  caption = "Tabla 3: Distribución de frecuencias resumen con acumulados exactos"
)
Tabla 3: Distribución de frecuencias resumen con acumulados exactos
Clase Li Ls MC ni hi Niasc Nidsc Hiasc Hidsc
1 10 72 41 72 2.88 72 2499 2.88 100.00
2 72 134 103 193 7.72 265 2427 10.60 97.12
3 134 196 165 210 8.40 475 2234 19.01 89.40
4 196 258 227 303 12.12 778 2024 31.13 80.99
5 258 320 289 362 14.49 1140 1721 45.62 68.87
6 320 382 351 401 16.05 1541 1359 61.66 54.38
7 382 444 413 361 14.45 1902 958 76.11 38.34
8 444 506 475 285 11.40 2187 597 87.52 23.89
9 506 568 537 169 6.76 2356 312 94.28 12.48
10 568 630 599 90 3.60 2446 143 97.88 5.72
11 630 692 661 42 1.68 2488 53 99.56 2.12
12 692 754 723 8 0.32 2496 11 99.88 0.44
13 754 816 785 3 0.12 2499 3 100.00 0.12

GRAFICAS DE DISTRIBUCION DE CANTIDAD

##  GRAFICA ORIGINAL SEGUN STURGES

# Histograma
hist(
  PROFUNDIDAD,
  breaks = h$breaks,
  col = "gray",
  border = "black",
  main = "Gráfica 1: Distribución de cantidad de la profundidad
  en depósitos minerales de Estados Unidos",
  xlab = "Profundidad (m)",
  ylab = "Frecuencia"
)

hist(PROFUNDIDAD,
     breaks = k,
     col = "gray",
     main = "Gráfica 2: Distribución global de la profundidad
     en depósitos minerales de Estados Unidos",
     xlab = "Profundidad (m)",
     ylab = "Frecuencia")

# Gráfica 3: Distribución porcentual
h <- hist(PROFUNDIDAD,
          breaks = k,
          plot = FALSE)

porcentaje <- h$counts / sum(h$counts) * 100

par(mar = c(7,4,4,2))

barplot(
  porcentaje,
  names.arg = paste(round(head(h$breaks,-1),0),
                    round(tail(h$breaks,-1),0),
                    sep = "-"),
  col = "darkorange",
  ylim = c(0, max(porcentaje)*1.15),
  ylab = "Frecuencia Relativa (%)",
  main = "Gráfica 3: Distribución porcentual de la profundidad",
  las = 2,
  space = 0,
  border = "white"
)

title(
  xlab = "Profundidad (m)",
  line = 5.5
)

media_rel_global <- mean(PROFUNDIDAD)
mediana_rel_global <- median(PROFUNDIDAD)

h <- hist(PROFUNDIDAD,
          breaks = k,
          plot = FALSE)

hi <- h$counts / sum(h$counts) * 100

intervalos_graf <- paste(
  round(head(h$breaks,-1),0),
  round(tail(h$breaks,-1),0),
  sep = " - "
)

par(mar = c(7,4,4,2))

barplot(
  hi,
  names.arg = intervalos_graf,
  col = "mediumpurple",
  ylim = c(0,100),
  cex.names = 0.7,
  space = 0,
  ylab = "Frecuencia Relativa (%)",
  xlab = "",
  main = "Gráfica 4: Distribución de cantidad en porcentaje de la
  profundidad en depósitos minerales de Estados Unidos",
  las = 2,
  border = "white"
)

title(
  xlab = "Intervalos de Profundidad (m)",
  line = 6
)

breaks_prof <- seq(
  0,
  max(PROFUNDIDAD) + A,
  by = A
)

histograma <- hist(
  PROFUNDIDAD,
  breaks = breaks_prof,
  col = rgb(0.7,0.7,0.7,0.7),
  border = "black",
  main = "Gráfica Nº5: Histograma con polígono de frecuencias
  de la profundidad de depósitos minerales de Estados Unidos",
  xlab = "Profundidad (m)",
  ylab = "Cantidad"
)

lines(
  histograma$mids,
  histograma$counts,
  type = "o",
  pch = 16,
  lwd = 2
)

# Gráfica 6: Ojiva Absoluta
plot(Ls, Niasc, 
     type = "o", 
     pch = 16, 
     col = "blue",
     main = "Gráfica 6: Ojiva de Frecuencias Absolutas Acumuladas \nde la profundidad de la muestra",
     xlab = "Profundidad (m)",
     ylab = "Frecuencia Acumulada Absoluta (Ni)")

lines(Ls, Nidsc, type = "o", pch = 16, col = "red")

legend("topleft",
       c("Ni Ascendente", "Ni Descendente"),
       col = c("blue", "red"),
       pch = 16)

# Gráfica 7: Ojiva Relativa
plot(Ls, Hiasc, 
     type = "o", 
     pch = 16, 
     col = "blue",
     main = "Gráfica 7: Ojiva de Frecuencias Relativas Acumuladas \nde la profundidad de la muestra",
     xlab = "Profundidad (m)",
     ylab = "Frecuencia Acumulada Relativa (%)")

lines(Ls, Hidsc, type = "o", pch = 16, col = "red")

legend("bottomright", 
       c("Hi Ascendente", "Hi Descendente"),
       col = c("blue", "red"), 
       pch = 16)

# Calcular todas las variables necesarias para el boxplot
media_box <- mean(PROFUNDIDAD)
mediana_box <- median(PROFUNDIDAD)
Q1_box <- quantile(PROFUNDIDAD, 0.25)
Q3_box <- quantile(PROFUNDIDAD, 0.75)
IQR_box <- IQR(PROFUNDIDAD)

# Calcular outliers 
lim_inf_outlier_box <- Q1_box - 1.5 * IQR_box
lim_sup_outlier_box <- Q3_box + 1.5 * IQR_box
outliers_box <- PROFUNDIDAD[PROFUNDIDAD < lim_inf_outlier_box | PROFUNDIDAD > lim_sup_outlier_box]
n_outliers_box <- length(outliers_box)

boxplot(PROFUNDIDAD,
        horizontal = TRUE,
        col = "lightblue",
        main = "Gráfica 8: Distribución de la Profundidad de la muestra \n con detección de valores atípicos",
        xlab = "Profundidad (m)",
        ylab = "",
        border = "darkblue")

# Añadir estadísticos importantes
points(media_box, 1, pch = 23, bg = "red", cex = 1.2)

legend("topright",
       legend = c(paste("Media:", round(media_box, 2)),
                  paste("Mediana:", round(mediana_box, 2)),
                  paste("Q1:", round(Q1_box, 2)),
                  paste("Q3:", round(Q3_box, 2)),
                  paste("Outliers:", n_outliers_box)),
       bty = "n",
       cex = 0.8)

summary(PROFUNDIDAD)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   10.47  231.90  335.97  335.38  438.40  815.94
hist(
  PROFUNDIDAD,
  breaks = h$breaks,
  col = "lightgray",
  border = "white",
  main = "Gráfica Nº9: Histograma con diagrama de caja superpuesto
  de la profundidad de depósitos minerales de Estados Unidos",
  xlab = "Profundidad (m)",
  ylab = "Cantidad"
)

boxplot(
  PROFUNDIDAD,
  horizontal = TRUE,
  at = par("usr")[4] * 0.70,
  add = TRUE,
  xaxt = "n",
  yaxt = "n",
  boxwex = par("usr")[4] * 0.30,
  col = rgb(0.2,0.2,1,0.30),
  border = "blue"
)

INDICADORES ESTADISTICOS Y OUTLIERS

# Calcular indicadores de tendencia central
minimo <- min(PROFUNDIDAD)
maximo <- max(PROFUNDIDAD)
rango <- maximo - minimo
media <- mean(PROFUNDIDAD)
mediana <- median(PROFUNDIDAD)
moda <- as.numeric(names(sort(table(PROFUNDIDAD), decreasing = TRUE)[1]))

# Crear tabla
tendencia_central <- data.frame(
  Indicador = c("Mínimo", "Media", "Mediana", "Moda", "Máximo", "Rango"),
  Valor = c(
    round(minimo, 4),
    round(media, 4),
    round(mediana, 4),
    round(moda, 4),
    round(maximo, 4),
    round(rango, 4)
  ),
  Unidad = rep("m", 6),
  Interpretación = c(
    "Valor mínimo observado",
    "Promedio de todos los valores",
    "Valor que divide la muestra en dos partes iguales",
    "Valor más frecuente en la muestra",
    "Valor máximo observado",
    "Diferencia entre el máximo y el mínimo"
  )
)

kable(tendencia_central,
      caption = "Tabla 3: Indicadores de Tendencia Central para la variable TOP_DEPTH_M",
      align = "l")
Tabla 3: Indicadores de Tendencia Central para la variable TOP_DEPTH_M
Indicador Valor Unidad Interpretación
Mínimo 10.4700 m Valor mínimo observado
Media 335.3845 m Promedio de todos los valores
Mediana 335.9700 m Valor que divide la muestra en dos partes iguales
Moda 233.2000 m Valor más frecuente en la muestra
Máximo 815.9400 m Valor máximo observado
Rango 805.4700 m Diferencia entre el máximo y el mínimo 
# Calcular indicadores de dispersión
varianza <- var(PROFUNDIDAD)
desv_est <- sd(PROFUNDIDAD)
CV <- (desv_est / media) * 100

# Interpretación del CV
if(CV < 15) {
  interpretacion_CV <- "BAJA (CV < 15%)"
} else if(CV < 30) {
  interpretacion_CV <- "MODERADA (15% ≤ CV < 30%)"
} else {
  interpretacion_CV <- "ALTA (CV ≥ 30%)"
}

# Crear tabla
dispersion <- data.frame(
  Indicador = c("Varianza", "Desviación Estándar", "Coeficiente de Variación"),
  Valor = c(
    round(varianza, 4),
    round(desv_est, 4),
    paste0(round(CV, 2), "%")
  ),
  Unidad = c("m²", "m", ""),
  Interpretación = c(
    "Medida de dispersión al cuadrado",
    "Dispersión promedio respecto a la media",
    paste("Dispersión relativa:", interpretacion_CV)
  )
)

kable(dispersion,
      caption = "Tabla 4: Indicadores de Dispersión para la variable TOP_DEPTH_M",
      align = "l")
Tabla 4: Indicadores de Dispersión para la variable TOP_DEPTH_M
Indicador Valor Unidad Interpretación
Varianza 21428.8446 m² Medida de dispersión al cuadrado
Desviación Estándar 146.3859 m Dispersión promedio respecto a la media
Coeficiente de Variación 43.65% Dispersión relativa: ALTA (CV ≥ 30%) 
# Calcular indicadores de posición
cuartiles <- quantile(PROFUNDIDAD)
Q1 <- cuartiles[2]
Q2 <- cuartiles[3]
Q3 <- cuartiles[4]
IQR_val <- IQR(PROFUNDIDAD)

# Detección de outliers
lim_inf_outlier <- Q1 - 1.5 * IQR_val
lim_sup_outlier <- Q3 + 1.5 * IQR_val
outliers <- PROFUNDIDAD[PROFUNDIDAD < lim_inf_outlier | PROFUNDIDAD > lim_sup_outlier]
n_outliers <- length(outliers)
porc_outliers <- round((n_outliers / n) * 100, 2)

# Crear tabla
posicion <- data.frame(
  Indicador = c("Cuartil 1 (Q1)", "Cuartil 2 (Q2 - Mediana)", "Cuartil 3 (Q3)", 
                "Rango Intercuartílico (IQR)", "Límite Inferior Outliers", 
                "Límite Superior Outliers", "Número de Outliers"),
  Valor = c(
    round(Q1, 4),
    round(Q2, 4),
    round(Q3, 4),
    round(IQR_val, 4),
    round(lim_inf_outlier, 4),
    round(lim_sup_outlier, 4),
    paste0(n_outliers, " (", porc_outliers, "%)")
  ),
  Unidad = c(rep("m", 6), "observaciones"),
  Interpretación = c(
    "25% de datos por debajo de este valor",
    "50% de datos por debajo de este valor (coincide con mediana)",
    "75% de datos por debajo de este valor",
    "Rango del 50% central de datos (Q3 - Q1)",
    "Límite inferior para detección de valores atípicos",
    "Límite superior para detección de valores atípicos",
    "Cantidad y porcentaje de valores atípicos detectados"
  )
)

kable(posicion,
      caption = "Tabla 5: Indicadores de Posición y detección de outliers en TOP_DEPTH_M",
      align = "l")
Tabla 5: Indicadores de Posición y detección de outliers en TOP_DEPTH_M
Indicador Valor Unidad Interpretación
Cuartil 1 (Q1) 231.895 m 25% de datos por debajo de este valor
Cuartil 2 (Q2 - Mediana) 335.97 m 50% de datos por debajo de este valor (coincide con mediana)
Cuartil 3 (Q3) 438.4 m 75% de datos por debajo de este valor
Rango Intercuartílico (IQR) 206.505 m Rango del 50% central de datos (Q3 - Q1)
Límite Inferior Outliers -77.8625 m Límite inferior para detección de valores atípicos
Límite Superior Outliers 748.1575 m Límite superior para detección de valores atípicos
Número de Outliers 4 (0.16%) observaciones Cantidad y porcentaje de valores atípicos detectados 
# Calcular coeficiente de asimetría de Fisher
asimetria <- moments::skewness(PROFUNDIDAD)

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

# Curtosis
curtosis <- moments::kurtosis(PROFUNDIDAD) - 3

if(abs(curtosis) < 0.5) {
  interpretacion_curtosis <- "Distribución mesocúrtica (normal)"
} else if(curtosis > 0) {
  interpretacion_curtosis <- "Distribución leptocúrtica (picuda)"
} else {
  interpretacion_curtosis <- "Distribución platicúrtica (aplanada)"
}

# Crear tabla
forma <- data.frame(
  Indicador = c("Coeficiente de Asimetría (Fisher)", "Interpretación Asimetría",
                "Coeficiente de Curtosis (Exceso)", "Interpretación Curtosis"),
  Valor = c(
    round(asimetria, 4),
    interpretacion_asimetria,
    round(curtosis, 4),
    interpretacion_curtosis
  ),
  Fórmula = 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"
  )
)

kable(forma,
      caption = "Tabla 6: Indicadores de forma de la distribución de Profundidad",
      align = "l")
Tabla 6: Indicadores de forma de la distribución de Profundidad
Indicador Valor Fórmula
Coeficiente de Asimetría (Fisher) 0.0728 g₁ = E[(X-μ)³]/σ³
Interpretación Asimetría Distribución simétrica |g₁| < 0.5: Simétrica; g₁ > 0: Positiva; g₁ < 0: Negativa
Coeficiente de Curtosis (Exceso) -0.5337 g₂ = E[(X-μ)⁴]/σ⁴ - 3
Interpretación Curtosis Distribución platicúrtica (aplanada) |g₂| < 0.5: Mesocúrtica; g₂ > 0: Leptocúrtica; g₂ < 0: Platicúrtica

Conclusión

La variable profundidad de la muestra presenta valores que fluctúan entre 10.47 m y 815.94 m, con valores en torno a la mediana de 335.97 m. La desviación estándar de 146.39 m y el coeficiente de variación de 43.65% indican un conjunto heterogéneo. La distribución es simétrica y presenta únicamente 4 valores atípicos, por lo que el comportamiento de la variable resulta beneficioso para el análisis minero, ya que existe una adecuada variabilidad en las profundidades evaluadas.