Interpolación Espacial
Karenth Yuliana Menjura Coy
2025-06-26
#Tabla de contenido ##1. Introducción ##2. Configuración ##3. Lectura de los datos de entrada ##4. Muestreo del mundo ##5. Interpolación ##5.1 Creación del objeto gstat ##5.2 Interpolación IDW ##5.3 Interpolación OK ##6. Evaluación de los resultados ##6.1 Evaluación cualitativa ##6.2 Validación cruzada ##7. Análisis de resultados ##8. Referencias
#1. Introducción La interpolación espacial es una herramienta fundamental en el análisis de datos geográficos, resulta útil para estimar valores en ubicaciones no muestreadas a partir de datos puntuales. En este informe se ilustran dos técnicas ampliamente utilizadas en estudios de suelos y recursos naturales: la ponderación inversa de la distancia (IDW) y el kriging ordinario (OK). La técnica IDW es un método determinista que asigna valores a partir de un promedio ponderado por la inversa de la distancia a los puntos conocidos, bajo la premisa de que los puntos más cercanos tienen mayor influencia que los distantes (Isaaks & Srivastava, 1989). Por otro lado, el kriging ordinario es un método geoestadístico que, además de estimar valores, proporciona una medida de la incertidumbre asociada, utilizando un modelo de variograma que representa la dependencia espacial entre las observaciones (Goovaerts, 1997).
Ambas técnicas se aplican en este estudio para generar una superficie continua del carbono orgánico del suelo (COS) a una profundidad de 15–30 cm, a partir de muestras extraídas de la base de datos SoilGrids con una resolución de 250 metros. La interpolación del COS es relevante debido a su importancia como indicador de la calidad del suelo y su rol clave en el ciclo global del carbono (Minasny et al., 2013). Una representación espacial precisa del COS permite apoyar la planificación agrícola, la conservación del suelo y la evaluación de prácticas de manejo sostenible.
El método IDW, aunque simple y fácil de implementar, no incorpora directamente la estructura espacial del fenómeno, lo que puede limitar su precisión en contextos con variabilidad espacial compleja. En contraste, el kriging incorpora dicha estructura a través del variograma, lo cual suele mejorar las estimaciones, aunque requiere un mayor conocimiento estadístico y una mayor carga computacional, es decir mayor potencia en el equipo(Li & Heap, 2008). Estudios comparativos han demostrado que el desempeño de cada técnica depende en gran medida de las características del terreno, la densidad de muestreo y la variabilidad de los datos (AbdelRahman et al., 2021). Este cuaderno, por tanto, busca no solo aplicar IDW y kriging ordinario al mapeo de COS, sino también comparar sus resultados para evaluar su utilidad en contextos similares. La comparación se enfoca tanto en la calidad de las predicciones como en la interpretación de los mapas generados, con el fin de identificar ventajas y limitaciones de cada enfoque en el análisis espacial de propiedades edáficas.
#2. Configuración Limpieza del espacio de trabajo: La funcion usada limpia completamente la memoria de R al inicio del análisis, para evitar que objetos viejos interfieran con los nuevos resultados.
rm(list=ls())Para leer y escribir datos geoespaciales, se usa terra (para ráster) y sf (para vectores). Para la interpolación espacial, se usa las bibliotecas gstat y automap. Para la visualización (gráficos y mapas), usaremos ggplot, leaflet y tmap.
Primero, se deben instalar las bibliotecas necesarias para cumplir con cada parte de analisis a realizar. (Se usa “{r message= FALSE}” para ocultar los mensajes del sistema que genera el código dentro del bloque R correspondiente. )
library(terra)
library(sf)
library(sp)
library(stars)
library(gstat)
library(automap)
library(RColorBrewer)
library(leaflet)
library(leafem)
library(curl)#3. Lectura de los datos de entrada Se require leer un conjunto de datos que simule datos reales. Por lo tanto, se lee la capa SOC que fue descargada de ISRIC usando la biblioteca Terra:
archivo <- ("C:/Users/menju/Downloads/soil_soc.tif")
(soc <- rast(archivo))## class : SpatRaster
## size : 837, 885, 1 (nrow, ncol, nlyr)
## resolution : 0.002259887, 0.002389486 (x, y)
## extent : -76.2578, -74.2578, 3.0295, 5.0295 (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326)
## source : soil_soc.tif
## name : soil_socSe convierten los datos de SOC a porcentaje. En SoilGrids, los valores de carbono orgánico del suelo (SOC) se almacenan como enteros, y para convertirlos a las unidades convencionales en g kg⁻¹ (‰) se utiliza un factor de escala de 100. Es decir:
Valor almacenado: 100 → Valor real: 100 ÷ 100 = 1 g kg⁻¹
Valor almacenado: 250 → Valor real: 2.5 g kg⁻¹
Este factor se detalla en la documentación de SoilGrids bajo los “Mapped units” y “Conversion factor” de la tabla oficial, donde muestra que el SOC está en centigramos por kilogramo (cg/kg), y para obtener g/kg es necesario dividir por 100,(ISRIC, 2023).
soc.perc <- soc/10La transformación de un Sistema de Referencia de Coordenadas (CRS) es esencial para asegurar que los datos espaciales estén correctamente alineados. El CRS más comúnmente utilizado en el mundo real es el WGS84, identificado también con el código EPSG:4326. Este sistema geodésico se basa en un modelo centrado en la masa de la Tierra (geocéntrico) y define un elipsoide, un datum y un sistema de coordenadas geográficas en grados decimales (Lat/Lon) (Wikipedia, 2024; MapScaping, 2023).
geog ="+proj=longlat +datum=WGS84"
(geog.soc = project(soc.perc, geog))## class : SpatRaster
## size : 861, 861, 1 (nrow, ncol, nlyr)
## resolution : 0.002321978, 0.002321978 (x, y)
## extent : -76.2578, -74.25858, 3.030277, 5.0295 (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +no_defs
## source(s) : memory
## name : soil_soc
## min value : 0.00000
## max value : 17.85256La conversión de un objeto de tipo SpatRaster (del paquete terra) a un objeto de clase stars se realiza para aprovechar las funcionalidades específicas del paquete stars, como el manejo eficiente de datos espaciales multidimensionales y su compatibilidad con funciones de visualización o análisis del paquete sf.
En R, esto se puede hacer mediante la función st_as_stars(),esto convierte la capa raster geog.soc a un formato compatible con la clase stars, facilitando su integración con funciones de análisis espacial o visualización geométrica (Pebesma, 2018). Conversion de la capa SpatRaster en un objeto de estrellas:
stars.soc = st_as_stars(geog.soc)stars.soc## stars object with 2 dimensions and 1 attribute
## attribute(s):
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## soil_soc 0 6.522442 8.990762 8.701489 10.47903 17.85256
## dimension(s):
## from to offset delta refsys x/y
## x 1 861 -76.26 0.002322 GEOGCRS["unknown",\n DA... [x]
## y 1 861 5.03 -0.002322 GEOGCRS["unknown",\n DA... [y]m <- leaflet() %>%
addTiles() %>%
leafem:::addGeoRaster(
stars.soc,
opacity = 0.8,
colorOptions = colorOptions(palette = c("orange", "yellow", "cyan", "green"),
domain = 0:20)
) #4. Muestreo del mundo Se obtiene una muestra de aproximadamente 500 sitios de datos reales utilizando una muestra aleatoria:Esta técnica de muestreo es fundamental en estudios espaciales para reducir el tamaño de los datos manteniendo la representatividad, y es ampliamente utilizada en análisis de suelos, vegetación o clima, el propósito es extraer ubicaciones espaciales representativas de forma aleatoria para análisis posterior. Además, se fija la semilla con set.seed(123456) para garantizar que el muestreo sea reproducible: cada vez que se ejecute el código, se obtendrán los mismos puntos. (Hijmans, 2023).
set.seed(123456)
(samples <- spatSample(geog.soc, 500, "random", as.points=TRUE))## class : SpatVector
## geometry : points
## dimensions : 500, 1 (geometries, attributes)
## extent : -76.25664, -74.27135, 3.03376, 5.026017 (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +no_defs
## names : soil_soc
## type : <num>
## values : 10.8
## 14.93
## 10.16Características principales del objeto samples: Es un data frame con 4 filas y 2 columnas (aunque solo se visualizan tres valores en la columna soil_soc).
La columna soil_soc contiene valores numéricos que representan la cantidad de carbono orgánico del suelo (SOC) en los puntos muestreados.
Es un objeto de clase SpatVector, aunque al imprimirlo como data.frame se muestra como una tabla para facilitar su interpretación.
Este tipo de objeto es útil para análisis posteriores como modelos espaciales, interpolaciones o validaciones de predicciones.
Ahora, se convierte el objeto spatVector en un objeto de características simple:
(muestras <- sf::st_as_sf(samples))## Simple feature collection with 500 features and 1 field
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -76.25664 ymin: 3.03376 xmax: -74.27135 ymax: 5.026017
## Geodetic CRS: GEOGCRS["unknown",
## DATUM["World Geodetic System 1984",
## ELLIPSOID["WGS 84",6378137,298.257223563,
## LENGTHUNIT["metre",1]],
## ID["EPSG",6326]],
## PRIMEM["Greenwich",0,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8901]],
## CS[ellipsoidal,2],
## AXIS["longitude",east,
## ORDER[1],
## ANGLEUNIT["degree",0.0174532925199433,
## ID["EPSG",9122]]],
## AXIS["latitude",north,
## ORDER[2],
## ANGLEUNIT["degree",0.0174532925199433,
## ID["EPSG",9122]]]]
## First 10 features:
## soil_soc geometry
## 1 10.801709 POINT (-74.33404 4.891342)
## 2 14.927185 POINT (-75.03063 4.933138)
## 3 10.160336 POINT (-75.24426 4.738092)
## 4 10.431653 POINT (-75.14906 4.480352)
## 5 6.982718 POINT (-74.63822 4.50125)
## 6 11.838931 POINT (-75.65757 4.858835)
## 7 10.609202 POINT (-75.58791 4.139021)
## 8 12.169822 POINT (-75.37197 4.217969)
## 9 6.077174 POINT (-75.12351 3.491189)
## 10 6.927933 POINT (-74.33868 3.268279)nmuestras <- na.omit(muestras)Visualiczacion de las muestras:
longit <- st_coordinates(muestras)[,1]
latit <- st_coordinates(muestras)[,2]
soc <- muestras$soil_socid <- seq(1,500,1) (sitios <- data.frame(id, longit, latit, soc))## id longit latit soc
## 1 1 -74.33404 4.891342 10.801709
## 2 2 -75.03063 4.933138 14.927185
## 3 3 -75.24426 4.738092 10.160336
## 4 4 -75.14906 4.480352 10.431653
## 5 5 -74.63822 4.501250 6.982718
## 6 6 -75.65757 4.858835 11.838931
## 7 7 -75.58791 4.139021 10.609202
## 8 8 -75.37197 4.217969 12.169822
## 9 9 -75.12351 3.491189 6.077174
## 10 10 -74.33868 3.268279 6.927933
## 11 11 -75.15138 4.645213 12.452711
## 12 12 -75.00277 5.026017 10.486536
## 13 13 -74.97491 3.050013 10.587042
## 14 14 -74.48729 4.603417 7.228143
## 15 15 -74.77289 3.546917 7.711100
## 16 16 -75.41608 3.544595 5.645353
## 17 17 -74.28760 3.551561 9.211179
## 18 18 -74.94937 4.700940 5.652468
## 19 19 -74.82862 3.166112 9.631239
## 20 20 -75.59023 3.256669 8.639798
## 21 21 -75.58327 4.234222 9.897621
## 22 22 -75.23961 3.084843 4.588855
## 23 23 -75.94317 4.027566 11.281161
## 24 24 -75.83404 4.626637 7.988769
## 25 25 -75.78296 4.022922 10.615541
## 26 26 -75.76670 3.574780 11.597811
## 27 27 -75.97104 4.069362 9.152866
## 28 28 -75.45091 4.749702 10.264052
## 29 29 -74.43621 3.730353 9.042409
## 30 30 -75.57398 3.187010 12.712980
## 31 31 -75.21175 3.948619 5.616924
## 32 32 -75.01206 4.575553 7.114019
## 33 33 -75.91763 4.380507 8.799308
## 34 34 -74.28760 3.974161 8.590240
## 35 35 -74.38048 4.552333 12.853364
## 36 36 -76.22878 4.844903 11.007408
## 37 37 -75.45556 3.811622 11.493668
## 38 38 -74.80540 3.930043 6.167068
## 39 39 -74.71949 4.817039 6.653548
## 40 40 -74.90293 4.222612 5.258292
## 41 41 -74.82862 4.208681 5.612529
## 42 42 -75.71562 4.928494 9.303795
## 43 43 -75.58791 4.310848 10.238413
## 44 44 -75.16299 5.005119 12.994623
## 45 45 -75.97104 4.127411 10.580025
## 46 46 -75.47413 3.781437 7.941092
## 47 47 -75.91066 4.642891 6.657237
## 48 48 -76.10803 4.224934 7.025329
## 49 49 -74.39441 4.543046 11.665659
## 50 50 -74.63358 4.714872 5.801264
## 51 51 -75.59720 4.863479 11.333754
## 52 52 -74.57088 3.279889 8.387492
## 53 53 -75.23032 3.514409 6.011607
## 54 54 -75.76902 3.319363 10.561338
## 55 55 -75.92692 3.395988 10.465277
## 56 56 -74.39673 3.649084 8.513721
## 57 57 -75.21639 3.818588 5.491396
## 58 58 -75.30927 3.263635 6.141519
## 59 59 -75.58094 4.004347 12.384007
## 60 60 -75.90370 4.817039 7.792376
## 61 61 -75.13280 3.523697 5.975424
## 62 62 -75.13512 4.847225 12.180776
## 63 63 -75.22104 4.666110 10.222572
## 64 64 -75.23032 4.612705 9.317332
## 65 65 -75.78760 4.127411 10.207938
## 66 66 -74.64983 3.272923 9.442617
## 67 67 -76.02212 3.050013 10.194357
## 68 68 -76.19859 3.405276 8.280443
## 69 69 -75.68311 4.378185 11.730165
## 70 70 -75.74116 4.963324 11.114596
## 71 71 -74.90757 4.640569 6.254510
## 72 72 -75.19085 3.249704 5.142672
## 73 73 -75.79224 3.862706 9.464290
## 74 74 -75.72026 3.869672 10.308852
## 75 75 -75.87351 4.199393 9.423700
## 76 76 -75.50896 4.120445 13.858252
## 77 77 -74.50587 4.238866 5.733439
## 78 78 -75.23497 3.683913 6.640060
## 79 79 -76.13357 3.447072 8.782784
## 80 80 -75.28373 4.598773 9.485847
## 81 81 -75.47181 4.552333 10.863164
## 82 82 -75.50896 4.280662 10.803702
## 83 83 -74.73342 3.335617 10.147791
## 84 84 -74.81237 3.212552 10.010519
## 85 85 -75.83404 3.635152 9.100148
## 86 86 -75.40215 3.084843 6.457154
## 87 87 -74.94704 4.183139 5.324570
## 88 88 -75.72026 3.382056 10.432694
## 89 89 -74.41531 4.145987 10.175943
## 90 90 -76.00587 5.016729 12.320331
## 91 91 -75.96871 3.604966 9.687784
## 92 92 -75.37893 4.510538 10.188061
## 93 93 -75.09797 4.139021 5.882924
## 94 94 -76.05695 3.744285 10.364605
## 95 95 -75.38358 4.371219 12.568676
## 96 96 -76.12429 3.516731 10.923639
## 97 97 -74.65680 4.322458 6.324320
## 98 98 -75.70865 4.350321 11.539417
## 99 99 -75.20478 4.889020 11.290785
## 100 100 -76.11732 3.623542 9.171146
## 101 101 -75.25122 3.646762 6.183482
## 102 102 -75.15370 3.802334 5.806458
## 103 103 -74.79147 3.196298 9.554635
## 104 104 -74.28528 4.208681 9.486403
## 105 105 -74.71949 4.180817 5.537827
## 106 106 -75.10958 4.331745 6.272832
## 107 107 -75.31624 3.121995 5.502897
## 108 108 -75.81779 3.651406 10.804747
## 109 109 -74.45943 4.378185 10.012559
## 110 110 -75.64828 4.721838 13.724036
## 111 111 -75.22568 4.961002 11.733837
## 112 112 -74.51748 4.120445 9.955345
## 113 113 -75.91763 4.782209 7.699980
## 114 114 -76.16376 4.387473 7.218832
## 115 115 -74.33868 3.087165 5.368257
## 116 116 -75.77367 3.611932 13.180048
## 117 117 -75.27444 3.909145 5.345487
## 118 118 -76.23110 3.075555 8.438118
## 119 119 -74.70324 3.354193 8.573370
## 120 120 -75.04921 3.971839 5.778007
## 121 121 -75.97104 3.725709 9.795695
## 122 122 -75.20478 3.296143 6.148403
## 123 123 -74.48961 4.612705 6.932626
## 124 124 -75.68079 4.401405 12.521132
## 125 125 -74.72646 3.725709 10.665610
## 126 126 -75.28373 4.656822 9.622918
## 127 127 -75.10726 4.141343 6.359622
## 128 128 -75.82243 4.350321 7.101722
## 129 129 -74.89364 3.700167 8.076306
## 130 130 -76.16840 3.286855 9.079433
## 131 131 -75.77135 3.693201 10.625616
## 132 132 -75.90138 4.248154 12.064674
## 133 133 -74.33868 4.587163 9.852722
## 134 134 -76.00122 3.584068 9.710753
## 135 135 -74.57553 4.598773 6.439574
## 136 136 -76.08946 4.582519 7.714114
## 137 137 -74.88899 3.614254 6.990876
## 138 138 -75.64596 3.207908 7.617199
## 139 139 -75.53683 3.298465 10.343254
## 140 140 -75.46949 3.774471 6.428236
## 141 141 -74.77986 4.494284 6.198705
## 142 142 -74.93311 4.083294 5.381612
## 143 143 -75.14209 3.372768 5.435059
## 144 144 -75.97800 4.512860 6.244603
## 145 145 -75.67847 3.600322 11.993864
## 146 146 -75.57166 3.070911 12.726831
## 147 147 -75.45788 3.045369 5.991716
## 148 148 -75.53450 3.892892 11.216869
## 149 149 -75.68543 4.628959 10.402876
## 150 150 -75.89441 4.895986 7.947830
## 151 151 -74.94240 4.143665 5.625063
## 152 152 -75.44395 3.505121 6.249770
## 153 153 -75.19549 3.581746 6.250525
## 154 154 -75.37197 3.614254 6.803142
## 155 155 -74.52444 4.621993 5.847041
## 156 156 -75.47181 3.523697 6.726218
## 157 157 -75.04921 4.310848 6.414926
## 158 158 -75.96175 3.990415 10.971891
## 159 159 -74.98652 3.135927 10.995726
## 160 160 -74.83559 5.009763 6.835690
## 161 161 -75.42305 3.217196 6.794930
## 162 162 -74.66144 4.812395 7.247848
## 163 163 -75.91531 3.458682 13.190475
## 164 164 -74.49426 3.414564 10.187453
## 165 165 -74.30618 4.842581 9.327036
## 166 166 -75.79224 3.658372 10.623851
## 167 167 -74.76361 3.755895 6.534280
## 168 168 -75.54147 4.487318 11.332709
## 169 169 -75.09101 4.211003 6.302176
## 170 170 -74.74967 4.025244 6.833659
## 171 171 -76.04998 3.672304 8.815156
## 172 172 -74.89364 3.416886 9.121772
## 173 173 -74.37119 3.124317 5.462683
## 174 174 -74.56624 3.163790 6.466397
## 175 175 -74.67073 4.357287 6.019896
## 176 176 -75.75974 4.080972 11.895828
## 177 177 -75.92227 4.060074 11.077200
## 178 178 -75.59720 4.849547 12.339217
## 179 179 -74.90525 4.905274 5.848634
## 180 180 -74.36423 4.183139 10.724874
## 181 181 -76.25200 3.435462 6.444319
## 182 182 -75.00045 4.285306 6.115094
## 183 183 -74.64054 4.951714 9.979046
## 184 184 -74.49194 4.187783 7.424860
## 185 185 -75.19782 4.048464 5.487835
## 186 186 -75.23497 4.949392 13.024869
## 187 187 -74.63125 5.014407 10.255013
## 188 188 -75.70633 3.869672 12.666335
## 189 189 -76.17537 3.193976 8.939064
## 190 190 -76.19162 3.182366 10.974299
## 191 191 -74.47800 4.847225 9.137359
## 192 192 -75.79457 3.632830 9.705297
## 193 193 -74.49890 4.139021 10.791475
## 194 194 -74.72878 4.002025 5.754549
## 195 195 -74.51051 3.152180 8.081536
## 196 196 -76.07552 4.134377 8.576612
## 197 197 -75.72723 4.905274 12.383555
## 198 198 -75.53450 4.617349 12.816098
## 199 199 -74.98419 4.635925 6.139585
## 200 200 -74.85184 3.119673 9.555013
## 201 201 -75.73187 3.191654 9.675797
## 202 202 -75.24658 4.380507 11.578958
## 203 203 -75.11423 3.679270 5.828213
## 204 204 -75.81082 3.878960 10.100000
## 205 205 -75.70169 3.498155 7.846942
## 206 206 -76.15447 3.558527 7.946811
## 207 207 -75.03528 4.222612 5.897088
## 208 208 -74.47568 3.772149 9.708948
## 209 209 -76.08713 3.050013 9.447316
## 210 210 -74.77522 3.855740 5.935575
## 211 211 -74.94008 4.893664 5.394388
## 212 212 -75.46717 4.675398 10.738787
## 213 213 -75.88048 3.526019 10.015093
## 214 214 -74.35726 3.765183 9.157049
## 215 215 -75.33946 4.429269 12.373422
## 216 216 -74.68002 3.741963 10.100106
## 217 217 -74.60803 3.138249 6.556734
## 218 218 -76.25664 3.279889 7.436344
## 219 219 -75.23497 3.611932 6.342598
## 220 220 -75.55540 5.021373 11.095710
## 221 221 -75.05153 3.567815 5.977952
## 222 222 -74.85649 3.667660 6.476151
## 223 223 -75.18156 4.872766 12.454919
## 224 224 -75.31856 3.228806 6.230792
## 225 225 -75.72955 4.252798 11.377036
## 226 226 -76.17537 4.118123 6.811135
## 227 227 -75.73187 4.847225 10.710344
## 228 228 -75.93621 4.742736 5.821217
## 229 229 -74.54070 3.635152 11.495430
## 230 230 -74.36190 4.076328 9.588242
## 231 231 -75.72258 3.746607 10.531799
## 232 232 -75.07243 4.510538 9.563382
## 233 233 -74.99580 5.002797 11.168987
## 234 234 -75.26283 4.366575 8.309416
## 235 235 -74.47336 4.231900 6.940283
## 236 236 -74.54070 3.256669 10.576662
## 237 237 -75.14673 4.563943 11.976322
## 238 238 -74.56392 4.591807 6.515208
## 239 239 -75.67150 4.399083 12.982393
## 240 240 -75.18388 4.324780 6.743471
## 241 241 -74.82862 4.410693 5.251880
## 242 242 -75.94549 3.688557 10.239090
## 243 243 -75.00741 3.043047 8.074848
## 244 244 -75.85261 3.688557 9.561659
## 245 245 -75.19317 4.689330 11.816648
## 246 246 -75.26748 3.946297 5.066669
## 247 247 -76.24967 3.832520 8.381060
## 248 248 -75.46020 3.927721 11.578142
## 249 249 -74.69163 3.447072 9.162738
## 250 250 -75.20478 4.847225 13.415951
## 251 251 -74.37584 4.475708 11.118632
## 252 252 -75.83172 3.386700 11.403005
## 253 253 -74.89132 3.349549 8.978333
## 254 254 -74.67073 3.941653 7.246947
## 255 255 -74.92382 3.653728 6.220556
## 256 256 -74.43156 4.417659 10.616484
## 257 257 -74.54302 3.523697 9.145444
## 258 258 -74.97026 4.498928 4.611331
## 259 259 -74.35262 3.718743 8.993698
## 260 260 -74.85649 4.879732 5.663822
## 261 261 -75.00045 3.182366 6.984000
## 262 262 -75.65525 3.291499 10.058765
## 263 263 -74.49194 4.807751 9.673142
## 264 264 -74.69163 3.746607 10.636909
## 265 265 -75.07475 3.676948 5.854684
## 266 266 -74.76593 4.155275 5.959586
## 267 267 -75.04689 4.935460 13.719927
## 268 268 -75.75974 3.556205 11.299211
## 269 269 -75.30231 4.608061 10.024953
## 270 270 -74.58249 3.307753 9.898138
## 271 271 -75.64131 4.178495 11.210106
## 272 272 -74.47568 3.247382 6.993076
## 273 273 -76.15679 3.772149 10.319891
## 274 274 -74.70556 3.502799 9.069007
## 275 275 -74.78218 4.176173 6.465563
## 276 276 -74.30850 3.888248 9.619830
## 277 277 -74.28760 4.626637 9.988910
## 278 278 -75.95943 4.385151 6.958313
## 279 279 -76.18001 4.452488 11.335790
## 280 280 -74.80540 4.480352 5.841406
## 281 281 -75.30927 3.180044 5.761398
## 282 282 -75.66453 4.002025 9.895706
## 283 283 -74.77754 4.717194 6.581511
## 284 284 -75.28837 5.002797 11.282930
## 285 285 -74.91918 4.986543 5.843883
## 286 286 -76.11964 3.103419 9.877734
## 287 287 -74.53605 4.972612 11.209393
## 288 288 -74.34565 4.682364 10.306379
## 289 289 -75.42305 4.248154 11.593365
## 290 290 -75.60649 4.503572 12.846023
## 291 291 -74.40602 4.136699 10.600320
## 292 292 -75.52522 4.670754 11.062603
## 293 293 -74.72181 5.016729 5.945249
## 294 294 -74.78218 4.508216 6.264574
## 295 295 -75.55076 3.974161 10.680253
## 296 296 -75.74348 3.509765 10.317769
## 297 297 -76.18234 4.299238 6.379967
## 298 298 -75.07475 4.970290 10.429973
## 299 299 -74.40138 3.997381 9.352732
## 300 300 -76.07552 3.818588 10.356620
## 301 301 -76.16840 4.385151 7.297194
## 302 302 -76.23110 3.809300 11.155326
## 303 303 -75.83636 3.772149 10.252287
## 304 304 -74.88667 3.825554 5.463158
## 305 305 -74.54302 3.423852 8.426931
## 306 306 -75.60649 3.361158 6.450453
## 307 307 -76.04302 3.702489 11.070751
## 308 308 -75.16995 4.524470 11.663638
## 309 309 -74.37351 4.480352 11.347538
## 310 310 -74.63358 3.983449 7.773649
## 311 311 -76.10107 3.465648 10.972401
## 312 312 -74.34797 3.723387 8.933864
## 313 313 -75.82243 3.279889 11.497157
## 314 314 -75.46252 4.357287 11.067483
## 315 315 -75.01206 3.198620 6.275005
## 316 316 -74.42924 3.077877 5.632380
## 317 317 -74.76361 4.900630 5.812014
## 318 318 -74.73342 4.895986 7.277349
## 319 319 -74.32243 3.711777 9.594090
## 320 320 -75.06314 4.698618 12.846326
## 321 321 -76.13590 3.402954 10.061088
## 322 322 -75.02599 3.414564 7.800445
## 323 323 -74.80772 4.350321 5.833250
## 324 324 -74.99813 3.082521 8.132365
## 325 325 -74.58946 3.149858 8.396154
## 326 326 -75.79689 4.289950 11.463914
## 327 327 -75.89673 4.612705 10.654577
## 328 328 -75.11887 4.758990 12.796000
## 329 329 -74.43621 4.269052 9.016408
## 330 330 -75.08172 4.810073 13.628903
## 331 331 -76.10803 4.220291 7.098499
## 332 332 -75.13048 4.385151 6.399830
## 333 333 -76.16144 3.918433 11.088652
## 334 334 -75.59488 3.454038 6.657052
## 335 335 -76.15447 4.847225 12.469318
## 336 336 -75.69472 3.830198 11.906689
## 337 337 -76.15447 4.345677 4.791567
## 338 338 -75.10262 4.113479 5.552363
## 339 339 -74.36423 4.299238 10.398923
## 340 340 -75.86887 4.884376 7.100000
## 341 341 -74.48033 4.422303 9.600532
## 342 342 -75.01438 3.663016 5.071311
## 343 343 -74.50819 3.579424 9.005939
## 344 344 -74.62661 4.937782 6.429270
## 345 345 -75.59488 3.296143 9.450234
## 346 346 -75.47878 4.150631 11.461810
## 347 347 -75.22800 3.047691 5.053403
## 348 348 -74.40370 3.147536 5.520484
## 349 349 -76.22413 4.354965 8.914068
## 350 350 -75.96871 4.605739 5.880092
## 351 351 -76.00819 3.999703 10.678988
## 352 352 -74.67305 3.214874 9.357347
## 353 353 -75.92924 4.076328 10.473546
## 354 354 -74.62661 4.074006 6.313595
## 355 355 -74.87971 3.995059 2.854747
## 356 356 -75.51593 3.751251 6.493718
## 357 357 -74.52212 4.252798 6.600754
## 358 358 -74.95401 4.789175 7.423533
## 359 359 -75.02367 4.013634 6.281922
## 360 360 -74.27831 5.002797 9.683897
## 361 361 -75.63899 4.183139 10.602435
## 362 362 -75.36500 5.005119 12.483944
## 363 363 -76.10571 3.600322 9.947707
## 364 364 -74.78450 3.904501 6.035358
## 365 365 -75.76438 4.043820 10.615919
## 366 366 -75.11423 3.507443 4.031125
## 367 367 -74.73110 3.574780 8.572993
## 368 368 -74.83094 3.776793 5.430199
## 369 369 -74.63358 4.956358 9.289946
## 370 370 -74.96098 3.495833 7.348737
## 371 371 -74.44317 3.289177 9.829047
## 372 372 -74.69859 3.683913 8.982966
## 373 373 -75.36964 4.457132 12.474405
## 374 374 -74.95401 4.494284 5.184835
## 375 375 -75.41840 4.020600 11.446390
## 376 376 -76.16376 4.062396 7.563533
## 377 377 -75.07475 4.429269 6.200038
## 378 378 -74.34333 4.436235 11.660369
## 379 379 -75.99890 4.229578 7.041994
## 380 380 -74.83094 4.347999 5.832593
## 381 381 -74.43156 3.056979 5.317156
## 382 382 -74.80772 4.466420 4.933749
## 383 383 -74.49890 3.909145 8.929915
## 384 384 -74.97491 3.263635 6.218020
## 385 385 -76.03373 4.426947 5.537522
## 386 386 -75.33017 4.756668 9.199310
## 387 387 -75.86655 3.505121 10.149097
## 388 388 -74.70091 3.586390 9.611584
## 389 389 -74.98419 4.436235 4.891315
## 390 390 -74.59875 4.292272 6.863955
## 391 391 -76.11732 3.839486 11.038663
## 392 392 -75.31624 5.014407 10.757126
## 393 393 -75.51593 3.560849 6.430840
## 394 394 -75.63435 3.800012 12.060020
## 395 395 -74.41763 3.660694 8.560393
## 396 396 -75.50200 3.477257 5.979239
## 397 397 -76.14518 3.512087 13.538759
## 398 398 -75.48342 4.984221 12.870726
## 399 399 -74.34565 3.683913 10.215023
## 400 400 -76.19395 3.679270 7.801891
## 401 401 -75.30695 4.501250 13.719677
## 402 402 -75.56933 3.154502 11.570788
## 403 403 -75.45556 3.286855 7.271574
## 404 404 -75.60184 4.118123 9.451626
## 405 405 -75.06082 4.680042 12.772438
## 406 406 -75.51361 4.317814 10.993743
## 407 407 -75.54379 3.844130 12.293010
## 408 408 -74.36887 4.884376 10.198087
## 409 409 -75.88744 4.238866 11.258756
## 410 410 -74.96098 4.020600 5.525968
## 411 411 -75.56005 3.256669 10.516779
## 412 412 -75.15834 4.352643 6.394445
## 413 413 -74.28528 4.213325 9.754179
## 414 414 -76.18930 3.319363 9.015474
## 415 415 -75.62970 3.036081 10.633072
## 416 416 -74.58481 3.474935 9.765560
## 417 417 -75.02135 4.826327 8.752589
## 418 418 -75.89905 3.574780 9.939452
## 419 419 -74.36190 4.241188 10.266441
## 420 420 -75.64596 4.914562 11.981302
## 421 421 -76.15447 4.921528 11.508897
## 422 422 -76.19395 4.440878 12.007215
## 423 423 -74.45711 3.498155 9.471640
## 424 424 -75.79689 4.027566 10.092760
## 425 425 -74.38048 3.070911 5.039159
## 426 426 -74.53141 3.767505 10.377966
## 427 427 -75.59952 4.329423 10.300307
## 428 428 -75.95943 3.776793 9.591829
## 429 429 -75.36732 4.610383 10.340413
## 430 430 -76.03605 3.486545 9.275988
## 431 431 -75.87583 4.501250 7.261450
## 432 432 -74.30153 3.465648 10.823158
## 433 433 -74.88435 4.278340 5.653080
## 434 434 -74.34797 4.533758 10.177104
## 435 435 -75.35803 4.689330 8.880352
## 436 436 -75.98729 4.396761 6.541796
## 437 437 -75.34410 3.052335 5.554998
## 438 438 -75.58327 3.800012 8.128913
## 439 439 -75.45556 3.377412 8.206290
## 440 440 -76.11500 3.145214 11.489865
## 441 441 -75.37661 4.882054 9.840794
## 442 442 -75.32553 4.252798 12.019866
## 443 443 -75.78992 3.444750 11.360051
## 444 444 -75.33017 4.108836 12.692925
## 445 445 -75.63667 3.033760 11.665031
## 446 446 -75.24426 4.415337 9.942396
## 447 447 -75.57166 3.156824 11.836326
## 448 448 -75.78528 4.457132 5.073398
## 449 449 -74.68002 3.370446 10.713280
## 450 450 -74.91221 4.345677 5.889400
## 451 451 -74.28296 3.751251 8.314962
## 452 452 -75.01206 3.070911 8.286225
## 453 453 -74.29225 4.731126 10.057660
## 454 454 -75.45091 4.062396 12.983023
## 455 455 -74.74039 4.231900 4.951744
## 456 456 -75.01206 3.205586 6.555354
## 457 457 -74.51516 4.220291 6.833316
## 458 458 -74.87738 4.770599 6.805016
## 459 459 -74.63358 3.604966 12.100000
## 460 460 -74.27135 4.765955 9.631100
## 461 461 -75.46252 3.611932 5.603315
## 462 462 -76.01051 4.132055 7.930970
## 463 463 -74.75200 3.363480 8.990134
## 464 464 -76.01748 4.596451 6.248755
## 465 465 -75.32553 4.824005 9.749845
## 466 466 -75.37429 4.979577 11.825982
## 467 467 -75.57862 3.967195 11.335144
## 468 468 -75.39751 3.337939 6.525484
## 469 469 -75.32320 4.285306 13.472350
## 470 470 -75.11655 3.971839 5.619627
## 471 471 -76.09874 3.354193 9.878464
## 472 472 -74.99116 3.702489 4.359419
## 473 473 -75.91995 3.665338 10.501062
## 474 474 -76.03141 3.333295 9.991726
## 475 475 -74.60107 4.543046 7.213731
## 476 476 -74.41763 4.700940 10.752914
## 477 477 -74.55231 3.635152 8.435520
## 478 478 -74.91918 3.776793 5.315644
## 479 479 -75.51593 3.365802 10.859879
## 480 480 -74.72181 3.906823 6.419482
## 481 481 -75.97336 3.635152 9.986810
## 482 482 -74.70788 3.964873 6.958068
## 483 483 -74.28064 4.916884 13.770671
## 484 484 -76.18466 4.454810 12.171570
## 485 485 -74.63822 4.373541 6.279824
## 486 486 -75.99193 3.416886 10.671240
## 487 487 -76.11964 3.307753 9.647901
## 488 488 -75.65525 4.531436 13.233150
## 489 489 -75.59720 3.960229 10.222795
## 490 490 -74.72413 4.480352 5.566048
## 491 491 -76.24967 3.289177 7.161597
## 492 492 -75.92692 4.680042 5.848316
## 493 493 -76.05463 3.402954 9.804360
## 494 494 -75.06314 3.904501 5.777020
## 495 495 -74.71485 4.710228 7.997286
## 496 496 -75.33714 4.981899 9.407422
## 497 497 -74.60339 4.217969 6.648452
## 498 498 -75.16763 4.183139 6.625496
## 499 499 -75.83172 4.134377 9.051435
## 500 500 -76.07320 5.005119 12.846138Eliminar los valores NA:La función na.omit() elimina todas las filas que contienen valores faltantes (NA), mientras que drop_na() permite mayor control cuando se trabaja con múltiples columnas.
sitios <- na.omit(sitios)head(sitios)## id longit latit soc
## 1 1 -74.33404 4.891342 10.801709
## 2 2 -75.03063 4.933138 14.927185
## 3 3 -75.24426 4.738092 10.160336
## 4 4 -75.14906 4.480352 10.431653
## 5 5 -74.63822 4.501250 6.982718
## 6 6 -75.65757 4.858835 11.838931Visualización de muestras()
m <- leaflet() %>%
addTiles() %>%
leafem:::addGeoRaster(
stars.soc,
opacity = 0.7,
colorOptions = colorOptions(palette = c("orange", "yellow", "cyan", "green"),
domain = 0:20)
) %>%
addMarkers(lng=sitios$longit,lat=sitios$latit, popup=sitios$soc, clusterOptions = markerClusterOptions())
m # Print the map#Interpretacion Distribución espacial del contenido de carbono orgánico del suelo (SOC) sobre una región del centro de Colombia (centrado al departamento del Tolima).
El fondo está compuesto por un ráster coloreado (stars.soc) que indica el valor de SOC, con una paleta de colores personalizada:
Naranja → Amarillo → Cian → Verde representa valores crecientes de SOC en el dominio de 0 a 20.
Por ejemplo, las zonas más verdes tienen un contenido más alto de SOC. Los puntos agrupados (clusters) sonlos números en los círculos indican cuántas observaciones se encuentran agrupadas en esa zona. Esto facilita la visualización cuando hay muchos puntos cercanos. La mayor concentración de muestras está alrededor del centro del mapa (Quindío, Tolima, Huila), con valores de SOC que varían espacialmente. Algunas regiones muestran valores más bajos (en amarillo), mientras que otras muestran valores más altos (en verde), lo cual puede estar relacionado con factores edáficos, uso del suelo o cobertura vegetal.
Ahora esta todo listo para realizar las tareas de interpolación.
#5. Interpolación Para realizar una interpolación espacial en R, el primer paso consiste en crear un objeto de la clase gstat mediante la función del mismo nombre: gstat(). Este objeto encapsula toda la información necesaria para llevar a cabo el proceso de interpolación.
Un objeto gstat incluye: La definición del modelo de interpolación, los datos de calibración o entrenamiento que se utilizarán, un modelo de variograma y/o covariables, La forma en que se configura el objeto determina el tipo de interpolación que se realizará. Por ejemplo: -Si no se especifica un modelo de variograma, la función ejecuta una interpolación por IDW (Inverse Distance Weighting). -Si se incluye un modelo de variograma pero sin covariables, se aplica Kriging ordinario.
En general, para crear un objeto gstat, se deben definir tres parámetros principales: -formula: especifica la relación entre la variable a predecir (dependiente) y las variables explicativas (independientes o covariables). -data: corresponde al conjunto de datos de calibración (o entrenamiento), generalmente en formato espacial (sf, sp o data.frame). -model: hace referencia al modelo de variograma que describe la estructura espacial de los datos. Esta estructura permite aplicar distintos métodos de interpolación de forma flexible según la disponibilidad de información espacial y la complejidad del análisis deseado.
#5.2 Interpolación IDW Para interpolar el SOC mediante el método IDW, se crea el siguiente objeto gstat, especificando únicamente la fórmula y los datos:
g1 = gstat(formula = soil_soc~ 1, data = nmuestras)EL modelo de interpolación g1 está definido, ahora es posible usar la función de predicción para interpolar, es decir, para estimar los valores de precipitación.
La función de predicción acepta: -Un objeto ráster (stars), como dem -Un modelo (gstat), como g1 El ráster tiene dos funciones: 1. Especificar las ubicaciones donde queremos realizar predicciones (en todos los métodos) 2. Especificar valores de covariables (solo en Kriging Universal)
Se crea un objeto ráster con valores de celda iguales a 4:
# a simple copy
rrr = aggregate(geog.soc, 4)que es rrr
rrr## class : SpatRaster
## size : 216, 216, 1 (nrow, ncol, nlyr)
## resolution : 0.009287914, 0.009287914 (x, y)
## extent : -76.2578, -74.25161, 3.023311, 5.0295 (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +no_defs
## source(s) : memory
## name : soil_soc
## min value : 0.00000
## max value : 16.46579Definir valores
values(rrr) <-3Definir nombres
names(rrr) <- "Nuevo Valor"Ques es rrr ahora
rrr## class : SpatRaster
## size : 216, 216, 1 (nrow, ncol, nlyr)
## resolution : 0.009287914, 0.009287914 (x, y)
## extent : -76.2578, -74.25161, 3.023311, 5.0295 (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +no_defs
## source(s) : memory
## name : Nuevo Valor
## min value : 3
## max value : 3stars.rrr = st_as_stars(rrr)La siguiente expresión interpola valores SOC según el modelo definido en g1 y la plantilla ráster definida en stars.rrr:
z1 = predict(g1, stars.rrr)## [inverse distance weighted interpolation]Que es z1
z1## stars object with 2 dimensions and 2 attributes
## attribute(s):
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## var1.pred 2.878722 7.599568 8.997482 8.848649 9.995632 14.74128 0
## var1.var NA NA NA NaN NA NA 46656
## dimension(s):
## from to offset delta refsys x/y
## x 1 216 -76.26 0.009288 +proj=longlat +datum=WGS8... [x]
## y 1 216 5.03 -0.009288 +proj=longlat +datum=WGS8... [y]Los atributos incluidos en el objeto z1:“x” y “y”, que corresponden a las filas, y de los nombres de las columnas: “from”, “to”, “offset”, “delta” y “refsys”. Estos nombres permitirán identificar y manipular partes específicas del objeto.
Podemos subconjuntar solo el primer atributo y renombrarlo como “soc”:
z1 = z1["var1.pred",,]
names(z1) = "soc"se utiliza una paleta de color
paleta <- colorNumeric(palette = c("orange", "yellow", "cyan", "green"), domain = 10:100, na.color = "transparent")El ráster SOC interpolado, utilizando IDW, se muestra en la siguiente figura:
m <- leaflet() %>%
addTiles() %>%
leafem:::addGeoRaster(
z1,
opacity = 0.7,
colorOptions = colorOptions(palette = c("orange", "yellow", "cyan", "green"),
domain = 0:20)
) %>%
addMarkers(lng=sitios$longit,lat=sitios$latit, popup=sitios$soc, clusterOptions = markerClusterOptions()) %>%
addLegend("bottomright", pal=paleta, values= z1$soc,
title = "IDW SOC interpolation [%]"
)## Warning in pal(c(r[1], cuts, r[2])): Some values were outside the color scale
## and will be treated as NAm Interpretación: El mapa muestra la interpolación del carbono orgánico del suelo (SOC) utilizando el método de ponderación por distancia inversa (IDW). A partir de los valores muestreados en diferentes sitios, se generó un ráster que representa la distribución espacial del SOC, con una escala de colores que va del naranja (valores bajos) al verde (valores altos). Los puntos numerados sobre el mapa indican las ubicaciones de los sitios de muestreo, y la leyenda en la esquina inferior derecha permite interpretar los valores interpolados en porcentaje. La visualización abarca una región del centro-occidente de Colombia, (centrandose en el departamento del Tolima), destacando variaciones en el contenido de SOC a lo largo del territorio.
#5.3 Interpolación OK Los métodos de kriging requieren un modelo de variograma, el cual permite cuantificar de manera objetiva el patrón de autocorrelación espacial en los datos y asignar ponderaciones al realizar predicciones. Como primer paso, se puede calcular y examinar el variograma empírico utilizando la función variogram(). Esta función requiere dos argumentos principales: la fórmula, que especifica la variable dependiente y las covariables (de forma similar a la función gstat), y los datos, que corresponden a la capa de puntos que contiene tanto la variable dependiente como las covariables como atributos. Por ejemplo, la siguiente expresión calcula el variograma empírico de las muestras sin incluir covariables: variogram(SOC ~ 1, data = puntos).
v_emp_ok = variogram(soil_soc ~ 1, data=nmuestras)Representacion del variograma:
plot(v_emp_ok)
El gráfico muestra el variograma empírico, el cual representa la
semivarianza en función de la distancia entre pares de puntos
muestreados. Se observa un aumento progresivo de la semivarianza a
medida que crece la distancia, lo que indica la existencia de
autocorrelación espacial: los puntos cercanos tienden a ser más
similares entre sí. A partir de cierta distancia (alrededor de 80–90
unidades), la semivarianza se estabiliza, lo que sugiere que se ha
alcanzado el rango, es decir, la distancia a partir de la cual las
observaciones dejan de estar correlacionadas.
Hay varias maneras de ajustar un modelo de variograma a un variograma empírico. Se utilizará la más sencilla: el ajuste automático mediante la función autofitVariogram del paquete automap.
v_mod_ok = autofitVariogram(soil_soc ~ 1, as(nmuestras, "Spatial"))La función selecciona el tipo de modelo que mejor se ajusta y ajusta sus parámetros.
El modelo ajustado se puede representar gráficamente con plot:
plot(v_mod_ok)
En el gráfico se muestra el variograma empírico junto con el modelo
teórico ajustado, seleccionado automáticamente por la función. El modelo
utilizado corresponde al tipo Ste (esférico o similar), con los
siguientes parámetros: un nugget de 0.14 que representa la variabilidad
a escala muy pequeña o el error de medición, un sill de 9.4 que indica
la semivarianza máxima alcanzada, y un range de 93, lo que significa que
las observaciones separadas por distancias mayores a 93 unidades ya no
presentan autocorrelación espacial significativa. El ajuste del modelo
permite utilizar el kriging para realizar interpolaciones precisas
basadas en la estructura espacial de los datos.
El objeto resultante es en realidad una lista con varios componentes, incluyendo el variograma empírico y el modelo de variograma ajustado. El componente $var_model del objeto resultante contiene el modelo real:
v_mod_ok$var_model## model psill range kappa
## 1 Nug 0.1367562 0.00000 0.0
## 2 Ste 9.2934731 93.04843 0.4Ahora, el modelo de variograma se puede pasar a la función gstat y asi continuar con la interpolación de Kriging Ordinario:
## [using ordinary kriging]
g2 = gstat(formula = soil_soc ~ 1, model = v_mod_ok$var_model, data = nmuestras)
z2= predict(g2, stars.rrr)## [using ordinary kriging]Nuevamente, subconjunto del atributo de valores predichos y cambio de nombre:
z2 = z2["var1.pred",,]
names(z2) = "soc"Las predicciones del Kriging Ordinario se muestran en la siguiente figura:
m <- leaflet() %>%
addTiles() %>%
leafem:::addGeoRaster(
z2,
opacity = 0.7,
colorOptions = colorOptions(palette = c("orange", "yellow", "cyan", "green"),
domain = 0:20)
) %>%
addMarkers(lng=sitios$longit,lat=sitios$latit, popup=sitios$soc, clusterOptions = markerClusterOptions()) %>%
addLegend("bottomright", pal = paleta, values= z2$soc,
title = "OK SOC interpolation [%]"
)## Warning in pal(c(r[1], cuts, r[2])): Some values were outside the color scale
## and will be treated as NAm El mapa generado mediante Kriging Ordinario (OK) presenta una interpolación más precisa y suavizada del contenido de carbono orgánico del suelo (SOC), expresado en porcentajes. La superficie interpolada muestra una transición continua de colores, que van desde el naranja (valores bajos), pasando por el amarillo y el cyan, hasta llegar al verde (valores altos), de acuerdo con la paleta de colores definida en el código. Esta visualización se logra mediante la función addGeoRaster() de leafem, con una opacidad del 70% que permite observar tanto la interpolación como el mapa base.
En comparación con el mapa anterior generado por interpolación IDW, el mapa por Kriging muestra una superficie mucho más suave y realista, ya que no solo considera la distancia entre puntos, sino también la estructura espacial de los datos, representada a través del variograma. Mientras que el IDW tiende a producir patrones más abruptos y con límites definidos entre zonas, el Kriging genera una representación más coherente con la variación espacial continua del SOC. En conclusión, el método de Kriging Ordinario mejora significativamente la visualización espacial del carbono orgánico del suelo, al incorporar información adicional sobre la autocorrelación espacial, lo que se refleja en una interpolación más refinada y útil para análisis ambientales y agronómicos.
#6. Evaluación de resultados #6.1 Evaluación cualitativa Otra perspectiva de los tres resultados de interpolación:
colores <- colorOptions(palette = c("orange", "yellow", "cyan", "green"), domain = 0:20, na.color = "transparent")m <- leaflet() %>%
addTiles() %>%
addGeoRaster(stars.soc, opacity = 0.8, colorOptions = colores, group="RealWorld") %>%
addGeoRaster(z1, colorOptions = colores, opacity = 0.8, group= "IDW") %>%
addGeoRaster(z2, colorOptions = colores, opacity = 0.8, group= "OK") %>%
# Add layers controls
addLayersControl(
overlayGroups = c("RealWorld", "IDW", "OK"),
options = layersControlOptions(collapsed = FALSE)
) %>%
addLegend("bottomright", pal = paleta, values= z1$soc,
title = "Soil organic carbon [%]"
)## Warning in pal(c(r[1], cuts, r[2])): Some values were outside the color scale
## and will be treated as NAm En el mapa presentado se observa la distribución espacial del carbono orgánico del suelo (SOC) interpolado mediante tres métodos distintos: valores reales del mundo (RealWorld), interpolación por distancia inversa (IDW) y kriging ordinario (OK). Estas capas han sido sobrepuestas en una sola visualización, lo cual permite una comparación directa, aunque con cierta superposición visual por la transparencia compartida.
Se aprecia que los valores de SOC varían entre aproximadamente 4% y 14%, con una gradación cromática que va del naranja (valores bajos) al verde y cian (valores altos). En general, se distingue una distribución heterogénea del SOC en el área de estudio, con zonas de mayor contenido orgánico en suelos ubicadas principalmente hacia el centro y algunas regiones del norte del mapa.
La combinación de capas permite ver diferencias sutiles entre los métodos de interpolación. Por ejemplo, IDW tiende a suavizar los valores entre puntos muestreados, mostrando patrones más uniformes, mientras que el kriging ordinario (OK) ofrece una mayor variabilidad local, al incorporar la estructura espacial de los datos mediante el variograma. La capa RealWorld representa los datos originales o una interpolación directa a partir de ellos, sirviendo como base de comparación.En conjunto, el mapa resalta la importancia del método de interpolación elegido en el análisis espacial del SOC, ya que influye en la representación final del recurso en el territorio. Además, permite identificar zonas potenciales para conservación o manejo diferencial de suelos según su contenido de materia orgánica.
#6.2 Validación cruzada Se han estimado el carbono orgánico del suelo (SOC) utilizando dos métodos diferentes: IDW y Kriging Ordinario. Si bien es útil examinar y comparar los resultados gráficamente, también necesitamos una forma objetiva de evaluar la precisión de la interpolación. De esta manera, podemos elegir objetivamente el método más preciso entre los métodos de interpolación disponibles.
En la Validación Cruzada de Exclusión de Puntos, se procede de la siguiente manera: 1.Extraemos un punto de los datos de calibración. 2.Hacemos una predicción para ese punto. 3.Repetimos el proceso para todos los puntos. 4.Finalmente, obtenemos una tabla con un valor observado y un valor predicho para todos los puntos.
Podemos ejecutar la Validación Cruzada de Exclusión de Puntos utilizando la función gstat.cv, que acepta un objeto gstat.(Esta parte del codigo fue corrida en consola, para evitar errores o perdida de tiempo al cargalos mas de una vez, entonces, todas lan funciones se han puesto como texto para evitar erroresa la momento de la publicacion, recordando que las imagenes posteriosres a esto fueron llamadas desde el ordenador donde se guardaron previamente segun los reasultados obtenidos)
### run the following code from the console -- line-by-line
#cv1 = gstat.cv(g1)
#cv2 = gstat.cv(g2)El resultado de gstat.cv tiene los siguientes atributos:
var1.pred: Valor predicho var1.var: Varianza (solo para Kriging) observed: Valor observado residual: Observado-predicho zscore: Puntuación Z (solo para Kriging) fold: ID de validación cruzada
#cv1 = na.omit(cv1)#cv1Convirtamos el objeto cv1:
#cv1 = st_as_sf(cv1)Ahora, grafiquemos los residuos:
#sp::bubble(as(cv1[, "residual"], "Spatial"))
Ahora, calculamos índices de exactitud de predicción, como el error cuadrático medio (RMSE): Metrica de exactitud: cercania de los valores estimados respecto a los valores verdaderos cuando sus valores sean mayores al 10% del rango de datos o dominio, son muy flojas, por debajo de 10% son aceptables, si es menor al 5% es bueno
# This is the RMSE value for the IDW interpolation
#sqrt(sum((cv1$var1.pred - cv1$observed)^2) / nrow(cv1))Se repite el proceso con los resultados correctos:
Conversión:
#cv2 = na.omit(cv2)
#cv2#cv2 = st_as_sf(cv2)Calcular RSME para obtener resultados OK:
#sp::bubble(as(cv1[, "residual"], "Spatial"))
# This is the RMSE value for the OK interpolation
#sqrt(sum((cv2$var1.pred - cv2$observed)^2) / nrow(cv2))#7. Análisis de resultados Se calcularon los índices de exactitud de predicción utilizando el error cuadrático medio (RMSE) para ambas técnicas de interpolación. El RMSE obtenido para la interpolación mediante IDW fue de 1.4408, mientras que para el Kriging Ordinario (OK) fue de 1.1932. Estas métricas reflejan la cercanía entre los valores estimados y los valores observados, y permiten evaluar la precisión de cada método. Dado que ambos errores se encuentran por debajo del 10% del rango de los datos, se consideran aceptables; sin embargo, el RMSE del Kriging es menor, lo cual indica una mayor precisión en sus predicciones. Por tanto, con base en los resultados de la validación cruzada, se concluye que el Kriging Ordinario es el método de interpolación más preciso entre los dos evaluados.
#7. Referencias Lizarazo, I., 2023. Interpolación espacial del carbono orgánico del suelo. Disponible en https://rpubs.com/ials2un/soc_interp.
AbdelRahman, M. A. E., Saleh, A. M., & Gad, A. (2021). Comparative Study Between Ordinary Kriging and Inverse Distance Weighting for Mapping Soil Properties. Sustainability, 13(1), 194. https://doi.org/10.3390/su13010194
Goovaerts, P. (1997). Geostatistics for natural resources evaluation. Oxford University Press.
Isaaks, E. H., & Srivastava, R. M. (1989). An introduction to applied geostatistics. Oxford University Press.
Li, J., & Heap, A. D. (2008). A review of spatial interpolation methods for environmental scientists. Geoscience Australia, Record, 2008(23), 137.
Minasny, B., Malone, B. P., McBratney, A. B., Angers, D. A., Arrouays, D., Chambers, A., … & Winowiecki, L. (2013). Soil carbon 4 per mille. Geoderma, 292, 59–86. https://doi.org/10.1016/j.geoderma.2017.01.002
ISRIC – World Soil Information. (2023). Frequently asked questions about SoilGrids. Recuperado de https://www.isric.org/explore/soilgrids/faq-soilgrids
MapScaping. (2023, septiembre 11). What is the Difference Between WGS84 and EPSG:4326? Recuperado de https://mapscaping.com/difference-between-wgs84-and-epsg4326/
R-Spatial. (2025). Transform or convert coordinates of simple feature. Recuperado de la documentación del paquete sf.
Wikipedia. (2024). World Geodetic System. Recuperado de https://en.wikipedia.org/wiki/World_Geodetic_System
Pebesma, E. (2018). Simple Features for R: Standardized Support for Spatial Vector Data. The R Journal, 10(1), 439–446. https://doi.org/10.32614/RJ-2018-009
Hijmans, R. J. (2023). terra: Spatial Data Analysis. R package version 1.7-71. https://rspatial.org/terra/