Proyecto estrictamente académico. Datos: EarthChem PetDB (CC BY 4.0), descargados el 13/10/2026. Web del proyecto: https://crisrumiguano-art.github.io/geoquimica-earthchem-ml/

1. Carga de datos y librerías

Se lee la base limpia (159 760 análisis de roca total × 75 variables). La primera vez se descarga del repositorio del proyecto (26 MB).

knitr::opts_chunk$set(
  echo = TRUE,          # muestra el código en el documento
  message = FALSE,      # sin mensajes de carga de paquetes
  warning = FALSE,      # sin advertencias
  fig.align = "center" # gráficos centrados
)
library(e1071)         # asimetría y curtosis (install.packages("e1071"))
if (!file.exists("data/earthchem_limpio.csv.gz")) {
  dir.create("data", showWarnings = FALSE)
  download.file("https://github.com/crisrumiguano-art/geoquimica-earthchem-ml/raw/main/data/earthchem_limpio.csv.gz",
                "data/earthchem_limpio.csv.gz", mode = "wb")
}
base <- read.csv(gzfile("data/earthchem_limpio.csv.gz", encoding = "UTF-8"), check.names = FALSE)
dim(base)              # filas y columnas
## [1] 159760     75

2. Datos

Se toma la columna de la variable y se quitan las filas vacías (vacío = no medido, nunca cero).

x <- base[["Publication Year"]]
x <- x[!is.na(x)]                                      # quita vacíos
n <- length(x)
tipico <- x                                            # aquí se agrupan todos los datos

3. Intervalos (regla de Sturges)

El número de intervalos sale de la regla de Sturges, k = 1 + 3,3·log10(n); la amplitud se redondea a un número cómodo para que los límites sean fáciles de leer.

k <- round(1 + 3.3 * log10(n))                         # número de intervalos
a <- (max(tipico) - min(tipico)) / k                   # amplitud calculada
a <- 10                                               # se redondea hacia arriba a un número cómodo
desde <- floor(min(tipico) / a) * a                    # primer límite: múltiplo de a
k <- floor((max(tipico) - desde) / a) + 1              # intervalos que hacen falta para llegar al máximo
Li <- desde + (0:(k - 1)) * a                          # límites inferiores
Ls <- Li + a                                           # límites superiores
MC <- (Li + Ls) / 2                                    # marca de clase

4. Frecuencias

Tabla de distribución de frecuencias: límites, marca de clase, frecuencia absoluta y relativa, y acumuladas ascendentes y descendentes.

dentro <- x[x >= Li[1] & x <= Ls[k]]
ni <- as.vector(table(cut(dentro, breaks = c(Li, Ls[k]), right = FALSE, include.lowest = TRUE)))  # cuenta por intervalo [Li, Ls[
abajo <- sum(x < Li[1])                                # clase abierta "menor que" (outliers bajos)
arriba <- sum(x > Ls[k])                               # clase abierta "mayor o igual" (outliers altos)
hi <- round(ni / n * 100, 2)                           # frecuencia relativa (%)
todos <- c(abajo, ni, arriba)                          # con las clases abiertas, para acumular
cerradas <- 2:(k + 1)                                  # posiciones de los intervalos normales
Ni_asc <- cumsum(todos)[cerradas]                      # acumulada ascendente
Ni_dsc <- rev(cumsum(rev(todos)))[cerradas]            # acumulada descendente
Hi_asc <- round(cumsum(todos) / n * 100, 2)[cerradas]
Hi_dsc <- round(rev(cumsum(rev(todos))) / n * 100, 2)[cerradas]
TDF <- data.frame(Li, Ls, MC, ni, hi, Ni_asc, Ni_dsc, Hi_asc, Hi_dsc)
knitr::kable(TDF, caption = "Tabla N° 1: distribución de frecuencias")
Tabla N° 1: distribución de frecuencias
Li Ls MC ni hi Ni_asc Ni_dsc Hi_asc Hi_dsc
1930 1940 1935 12 0.01 12 159760 0.01 100.00
1940 1950 1945 9 0.01 21 159748 0.01 99.99
1950 1960 1955 73 0.05 94 159739 0.06 99.99
1960 1970 1965 1707 1.07 1801 159666 1.13 99.94
1970 1980 1975 13472 8.43 15273 157959 9.56 98.87
1980 1990 1985 24880 15.57 40153 144487 25.13 90.44
1990 2000 1995 27049 16.93 67202 119607 42.06 74.87
2000 2010 2005 35834 22.43 103036 92558 64.49 57.94
2010 2020 2015 38786 24.28 141822 56724 88.77 35.51
2020 2030 2025 17938 11.23 159760 17938 100.00 11.23 
c(abajo = abajo, arriba = arriba)                      # datos en las clases abiertas (no se borra ninguno)
##  abajo arriba 
##      0      0

5. Gráficos: histograma + polígono, ojiva y caja

Histograma con polígono de frecuencias, ojivas ascendente y descendente (se cruzan en la mediana) y diagrama de caja.

etiq <- paste0("[", Li, "-", Ls, "[")                  # texto de cada intervalo
nota <- c(if (abajo > 0) paste0("< ", Li[1], ": ", abajo, " (", round(abajo / n * 100, 2), " %)"),
          if (arriba > 0) paste0("≥ ", Ls[k], ": ", arriba, " (", round(arriba / n * 100, 2), " %)"))
layout(matrix(c(1, 1, 2, 3), 2, 2, byrow = TRUE), heights = c(1.15, 1))   # arriba 1 gráfico, abajo 2
par(mar = c(8, 5, 4, 1))
h <- barplot(hi, names.arg = etiq, space = 0, las = 2, col = "#3987e5", border = "white",   # histograma
             ylab = "Porcentaje (%)", cex.names = .75, ylim = c(0, max(hi) * 1.15),
             main = "Gráfica No. 1: Histograma y polígono de frecuencias de Año de publicación (año)")
lines(c(h[1] - 1, h, h[length(h)] + 1), c(0, hi, 0), col = "#d95926", lwd = 2, type = "b", pch = 16, cex = .6)  # polígono
mtext(if (length(nota) > 0) paste0("Clases abiertas (fuera del gráfico): ", paste(nota, collapse = " · ")) else "",
      side = 1, line = 6.5, cex = .7, col = "gray30")  # nota bajo el histograma
par(mar = c(5, 5, 4, 1))
plot(Ls, Hi_asc, type = "b", pch = 16, col = "#3987e5", ylim = c(0, 100),            # ojiva ascendente
     xlab = "Año de publicación (año)", ylab = "Porcentaje acumulado (%)",
     main = "Gráfica No. 2: Ojiva ascendente y descendente")
lines(Li, Hi_dsc, type = "b", pch = 17, col = "#d95926")                              # ojiva descendente
abline(h = 50, lty = 3, col = "gray50")                # línea del 50 %: donde se cruzan está la mediana
legend("right", c("ascendente (≤)", "descendente (≥)"), col = c("#3987e5", "#d95926"), pch = c(16, 17), bty = "n", cex = .85)
boxplot(x, horizontal = TRUE, col = "#86b6ef", xlab = "Año de publicación (año)",       # diagrama de caja
        main = "Gráfica No. 3: Diagrama de caja (todos los datos)")

6. Indicadores (con los datos sin agrupar)

Medidas de tendencia central, dispersión y forma con los datos sin agrupar (asimetría y curtosis con el paquete e1071).

ind <- c(n = n, media = mean(x), mediana = median(x), desviacion = sd(x), varianza = var(x),
  CV = sd(x) / mean(x) * 100, asimetria = skewness(x), curtosis = kurtosis(x))
as.numeric(names(which.max(table(x))))                # moda: el valor que más se repite
## [1] 2015
knitr::kable(data.frame(indicador = names(ind), valor = round(ind, 4)), row.names = FALSE,
             caption = "Tabla N° 2: indicadores estadísticos")
Tabla N° 2: indicadores estadísticos
indicador valor
n 159760.0000
media 2001.2396
mediana 2004.0000
desviacion 15.0637
varianza 226.9158
CV 0.7527
asimetria -0.2863
curtosis -0.9556

7. Conclusión

Se analizaron 159 760 casos con dato. La mediana es 2004 años; la asimetría es -0,29 (cola hacia los valores bajos) y el coeficiente de variación es 0,75 %. La moda es 2015. La mayor parte de los datos es posterior a 1990: la producción creció con las técnicas modernas de laboratorio (ICP-MS). El año es discreto, pero se agrupa en intervalos porque toma muchos valores distintos.