Proyecto estrictamente académico (Estadística y Probabilidades, FIGEMPA – UCE). 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

Se usa la base completa (159 760 análisis). Para cada par se toman solo las filas que tienen los dos datos: no se inventa ningún valor que falte.

knitr::opts_chunk$set(echo = TRUE, message = FALSE, warning = FALSE, fig.align = "center")
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)
## [1] 159760     75

2. Asignación de variables

Los 8 óxidos mayores de las rocas ígneas con los 8 completos. Todos a la vez: ya no un par, sino la estructura completa de la química de la roca.

oxidos <- c("SiO2 (wt%)", "Al2O3 (wt%)", "MgO (wt%)", "CaO (wt%)", "Na2O (wt%)", "K2O (wt%)", "P2O5 (wt%)", "MnO (wt%)")
ignea <- grepl("^igneous", base[["Sample Type"]])
datos <- na.omit(base[ignea, c(oxidos, "Sample Type")])
X <- datos[, oxidos]
colnames(X) <- sub(" [(]wt%[)]", "", colnames(X))      # nombres cortos: SiO2, Al2O3, ...
tipo <- sub("igneous:", "", datos[["Sample Type"]])   # volcanic:mafic, plutonic:ultramafic, ...
nrow(X)
## [1] 36420

3. Matriz de correlación

R <- cor(X)
knitr::kable(round(R, 2), caption = "Tabla N° 1: matriz de correlación de los óxidos mayores")
Tabla N° 1: matriz de correlación de los óxidos mayores
SiO2 Al2O3 MgO CaO Na2O K2O P2O5 MnO
SiO2 1.00 0.24 -0.49 -0.54 0.48 0.57 -0.13 -0.05
Al2O3 0.24 1.00 -0.67 0.32 0.43 0.15 0.01 -0.05
MgO -0.49 -0.67 1.00 -0.14 -0.58 -0.37 -0.01 -0.04
CaO -0.54 0.32 -0.14 1.00 -0.25 -0.56 0.04 -0.02
Na2O 0.48 0.43 -0.58 -0.25 1.00 0.58 0.03 -0.01
K2O 0.57 0.15 -0.37 -0.56 0.58 1.00 0.05 -0.02
P2O5 -0.13 0.01 -0.01 0.04 0.03 0.05 1.00 0.05
MnO -0.05 -0.05 -0.04 -0.02 -0.01 -0.02 0.05 1.00 
colores <- colorRampPalette(c("#d95926", "white", "#3987e5"))(21)
image(1:8, 1:8, R[, 8:1], col = colores, zlim = c(-1, 1), axes = FALSE, xlab = "", ylab = "",
      main = "Gráfica N° 1: Matriz de correlación")
axis(1, 1:8, colnames(R), las = 2); axis(2, 1:8, rev(colnames(R)), las = 1)
text(rep(1:8, 8), rep(8:1, each = 8), sprintf("%.2f", R), cex = 0.75)

El SiO2 se opone al MgO (r = -0.49), y los álcalis (Na2O, K2O) acompañan a la sílice: es la diferenciación magmática de la Unidad 3, ahora vista con todas las variables juntas.

4. Análisis de componentes principales (PCA)

El PCA reemplaza las 8 variables, muy correlacionadas, por unas pocas combinaciones nuevas (componentes) que no se correlacionan entre sí y resumen la mayor parte de la variación. Se estandariza cada óxido para que todos pesen igual.

pca <- prcomp(X, scale. = TRUE)
varianza <- pca$sdev^2 / sum(pca$sdev^2) * 100        # % de la variación que explica cada componente
knitr::kable(data.frame(componente = paste0("CP", 1:8), varianza = round(varianza, 1), acumulada = round(cumsum(varianza), 1)),
             caption = "Tabla N° 2: variación explicada por cada componente")
Tabla N° 2: variación explicada por cada componente
componente varianza acumulada
CP1 36.8 36.8
CP2 21.9 58.7
CP3 13.5 72.1
CP4 12.0 84.1
CP5 6.0 90.2
CP6 4.1 94.3
CP7 3.9 98.1
CP8 1.9 100.0 
barplot(varianza, names.arg = paste0("CP", 1:8), col = "#86b6ef", border = "white", las = 1,
        main = "Gráfica N° 2: Variación explicada", ylab = "%")

knitr::kable(round(pca$rotation[, 1:2], 3), caption = "Tabla N° 3: peso de cada óxido en CP1 y CP2")
Tabla N° 3: peso de cada óxido en CP1 y CP2
PC1 PC2
SiO2 0.467 -0.207
Al2O3 0.310 0.548
MgO -0.437 -0.395
CaO -0.239 0.643
Na2O 0.480 0.070
K2O 0.455 -0.276
P2O5 -0.012 0.068
MnO -0.013 -0.012 
set.seed(2026)
m <- sample(nrow(X), 6000)                             # 6000 rocas para que se vea
grupo <- ifelse(grepl("ultramafic", tipo), "ultramáfica", ifelse(grepl("mafic", tipo), "máfica",
                ifelse(grepl("intermediate", tipo), "intermedia", ifelse(grepl("felsic", tipo), "félsica", "otra"))))
paleta <- c("ultramáfica" = "#1b5e20", "máfica" = "#3987e5", "intermedia" = "#e3a008", "félsica" = "#d95926", "otra" = "gray80")
plot(pca$x[m, 1], pca$x[m, 2], pch = 16, cex = 0.5, col = paleta[grupo[m]], las = 1,
     xlim = quantile(pca$x[, 1], c(0.002, 0.998)), ylim = quantile(pca$x[, 2], c(0.002, 0.998)),   # sin el 0,4 % extremo
     main = "Gráfica N° 3: Las rocas en el plano CP1–CP2", xlab = paste0("CP1 (", round(varianza[1], 1), " %)"),
     ylab = paste0("CP2 (", round(varianza[2], 1), " %)"))
flechas <- pca$rotation[, 1:2] * 3                     # dirección de cada óxido
arrows(0, 0, flechas[, 1], flechas[, 2], length = 0.08, lwd = 2)
text(flechas[, 1] * 1.15, flechas[, 2] * 1.15, rownames(flechas), font = 2)
legend("topright", names(paleta), col = paleta, pch = 16, bty = "n")

5. Agrupamiento k-means

Sin decirle el tipo de roca, k-means reparte las rocas en k grupos de química parecida. Luego se compara con el tipo real.

Z <- scale(X)
set.seed(2026)
inercia <- sapply(1:8, function(k) kmeans(Z, k, nstart = 5, iter.max = 50)$tot.withinss)
plot(1:8, inercia / 1000, type = "b", pch = 16, col = "#3987e5", lwd = 2, las = 1,
     main = "Gráfica N° 4: Método del codo", xlab = "Número de grupos (k)", ylab = "Variación dentro de los grupos (miles)")

grupos <- kmeans(Z, 4, nstart = 10, iter.max = 50)$cluster
tabla <- table(grupo_kmeans = grupos, tipo_de_roca = grupo)
knitr::kable(tabla, caption = "Tabla N° 4: grupos de k-means frente al tipo de roca")
Tabla N° 4: grupos de k-means frente al tipo de roca
félsica intermedia máfica otra ultramáfica
116 2046 3987 990 71
943 1254 332 3323 53
23 647 18619 586 322
1 8 945 1 2153

6. Conclusión

Los ocho óxidos mayores están muy correlacionados: los dos primeros componentes principales resumen el 58.7 % de la variación. El CP1 opone el magnesio a la sílice y los álcalis: es el eje de la diferenciación magmática y ordena las rocas de ultramáficas a félsicas. Sin conocer el tipo de roca, k-means forma grupos que coinciden en buena parte con los tipos reales, lo que muestra que la química sola separa las rocas; ese es el fundamento de la red neuronal de la Unidad 4.