Interpolación de datos del departamento de Caquetá
Sandra Katteryne Rodríguez Hurtado
Junio de 2020
Introducción
En el presente trabajo se abordará el tema de interpolación de datos, en esta ocasión se realizarán los métodos y análisis correspondientes para un conjunto de datos de precipitación que cubren el departamento de Caquetá. Se hará la recopilación de varios cuadernos realizados durante el primer semestre de 2020, cuyas gráficas serán un apoyo para describir la zona de estudio, que es el departamento de Caquetá. Se hizo uso de la herramienta QGIS, la cual se trabajó en algunas de las clases del curso de geomática básica, además de consultar datos en el DANE o CHIRPS.
La interpolación de datos es muy importante para distintas carreras, entre ellas la Ingeniería Agronómica, donde el conocimiento de herramientas y métodos permite realizar importantes estudios de distintas variables como la degradación de suelos en una zona determinada, índice de vegetación, humedad del suelo, precipitación, entre otras. Estos datos permiten no sólo entender en contexto del momento, también predecir qué puede ocurrir en un futuro y, de ser posible, tomar medidas al respecto. Este trabajo fue elaborado para la materia de geomática básica, cuyo objetivo es entender distintas funcionalidades geoespaciales, generadas por el ambiente, en el programa de RStudio. Se hablará sobre qué es la interpolación, algunos conceptos básicos de la misma, explicación de tres técnicas diferentes para interpolación de datos de precipitación y su respectivo análisis, y algunas conclusiones.
Descripción de la zona
El departamento de Caquetá fue fundado el 02 de Mayo de 1845, se encuentra localizado al Sur de Colombia, tiene una extensión de 88.965 km2 aproximadamente y juega un papel importante en la biodiversidad debido a que un porcentaje importante de su territorio pertenece al sistema Amazónico, piedemonte Amazónico, se encuentran asentamientos humanos, y por ende la dinámica sociopolítica y económica (Gobernación de Caquetá, s.f).
Ubicación
Se encuentra ubicado hacia el sur del territorio colombiano, Su capital es Florencia y limita por el norte con el departamento del Huila y Meta, por el oriente con Vaupés y Guaviare, por el sur con Amazonas y Putumayo, y por el Occidente con Cauca y Huila (IGAC, 2008). Está dividido en diesiseis municipios, donde se encuentran muchos caseríos y sitios poblados.
Figura 1. Ubicación del departamento de Caquetá en Colombia. Mapa realizado en QGIS con datos de diva-GIS. Autoría propia. Se puede apreciar la ubicación del departamento de Caquetá, en color amarilo, al sur del país. 
Figura 2. Ubicación del departamento de Caquetá en Colombia. Mapa realizado en QGIS con datos de diva-GIS. Autoría propia. Se pueden apreciar los diesiseis municipios que conforman el departamento de Caquetá: Albania, Belén de los Andaquíes, Cartagena de Chaira, Curillo, El Doncello, El Paujil, Florencia, Milán, Montañita, Morelia, Puerto Rico, San José del Fragua, San Vicente del Caguán, Solano, Solita y Valparaiso.
Clima y población
La mayoría del departamento presenta un clima cálido que oscila entre los 27°C y los 29°C, alta pluviosidad, donde el verano generalmente se presenta en los meses de Diciembre a Febrero y los meses restantes son de invierno (López, 2006).
Necesidades Básicas Insatisfechas
Las Necesidades Básicas Insatisfechas (NBI) que presentan los habitantes de cada municipio y producción en el año 2018 de algunos cultivos en Caquetá.
Tabla 1. Datos de NBI que presentan los municipios de Caqueta. Los municipios se encuentran en orden descendente acorde con su respectivo porcentaje de NBI, donde el municipio de Milán es el que mayor porcentaje de NBI presenta, donde aproximadamente el 47% de la población tiene Necesidades Basicas Instatisfechas, mientras que el municipio con menor porcentaje de NBI es Florencia, donde cerca del 14.5% de su población presenta Necesidades Basicas Instafisfechas. Fuente: Autoria propia con datos recopilados por el DANE.
Reading layer `MGN_MPIO_POLITICO' from data source `D:\Documentos\Geomática\ADMINISTRATIVO1\MGN_MPIO_POLITICO.shp' using driver `ESRI Shapefile'
Simple feature collection with 16 features and 9 fields
geometry type: POLYGON
dimension: XY
bbox: xmin: -76.30622 ymin: -0.70584 xmax: -71.25385 ymax: 2.964148
geographic CRS: WGS 84Simple feature collection with 6 features and 3 fields
geometry type: POLYGON
dimension: XY
bbox: xmin: -76.1027 ymin: 0.9764735 xmax: -74.89527 ymax: 2.326755
geographic CRS: WGS 84
MUNICIPIO MPIO_CCDGO NBI geometry
1 FLORENCIA 18001 14.44902 POLYGON ((-75.42074 2.19413...
2 ALBANIA 18029 19.80289 POLYGON ((-75.89506 1.36569...
3 BELÉN DE LOS ANDAQUIES 18094 28.28503 POLYGON ((-75.78705 1.74982...
4 EL DONCELLO 18247 20.37901 POLYGON ((-75.36167 2.32142...
5 EL PAUJIL 18256 23.33489 POLYGON ((-75.36691 2.21234...
6 LA MONTAÑITA 18410 32.51475 POLYGON ((-75.40404 1.76944...Simple feature collection with 6 features and 21 fields
geometry type: POLYGON
dimension: XY
bbox: xmin: 774588.2 ymin: 599751.7 xmax: 908977.8 ymax: 749113.8
projected CRS: MAGNA-SIRGAS / Colombia Bogota zone
DPTO_CCDGO MPIO_CCDGO MPIO_CNMBR
1 18 18001 FLORENCIA
2 18 18029 ALBANIA
3 18 18094 BELÉN DE LOS ANDAQUÍES
4 18 18247 EL DONCELLO
5 18 18256 EL PAUJIL
6 18 18410 LA MONTAÑITA
MPIO_CRSLC MPIO_NAREA MPIO_NANO
1 Decreto 642 de Junio 17 de 1912 2547.638 2017
2 Ordenanza 3 de Noviembre 12 de 1985 414.122 2017
3 Decreto 963 de Marzo 14 de 1950 1191.619 2017
4 Decreto 1678 de Septiembre 7 de 1967 1105.803 2017
5 Decreto 1678 de Septiembret 7 de 1967 1234.743 2017
6 Decreto 83 de Julio 6 de 1955 1701.052 2017
DPTO_CNMBR Shape_Leng Shape_Area COD_DEPTO DEPTO COD_MUN
1 CAQUETÁ 2.942508 0.20692777 18 CAQUETÁ 001
2 CAQUETÁ 1.112829 0.03361758 18 CAQUETÁ 029
3 CAQUETÁ 2.234657 0.09674460 18 CAQUETÁ 094
4 CAQUETÁ 3.154370 0.08986744 18 CAQUETÁ 247
5 CAQUETÁ 3.529316 0.10030928 18 CAQUETÁ 256
6 CAQUETÁ 3.402939 0.13817351 18 CAQUETÁ 410
MUNICIPIO NBI MISERIA VIVIENDA
1 FLORENCIA 14.44902 2.774899 5.224497
2 ALBANIA 19.80289 4.400642 9.832684
3 BELÉN DE LOS ANDAQUIES 28.28503 7.091305 12.143430
4 EL DONCELLO 20.37901 3.898247 10.300478
5 EL PAUJIL 23.33489 5.170800 13.718609
6 LA MONTAÑITA 32.51475 8.972579 17.086081
SERVICIOS HACINAMIENTO INASISTENCIA ECONOMIA
1 1.488778 4.434199 2.077731 4.529297
2 3.208801 4.606922 1.260601 6.119642
3 3.906519 5.602016 4.926108 10.654141
4 1.473936 3.704758 2.025950 7.608696
5 2.823273 5.186398 2.737482 5.357978
6 7.150295 5.128428 2.724748 11.315515
geometry Pobreza
1 POLYGON ((850567.7 734450.4... Intermedia
2 POLYGON ((797693.6 642855.1... Intermedia
3 POLYGON ((809754.1 685340, ... Alta
4 POLYGON ((857152.4 748524, ... Intermedia
5 POLYGON ((856558 736458.9, ... Alta
6 POLYGON ((852387.3 687475.8... Altaopar <- par(mar = c(0,0,1.2,0))
par(bg="white")
plot(st_geometry(nbi_munic_2), col = "#e4e9de", border = "darkseagreen4",
bg = "grey85", lwd = 0.5)
choroLayer(
x = nbi_munic_new,
var = "NBI",
method = "geom",
nclass=5,
col = carto.pal(pal1 = "orange.pal", n1 = 5),
border = "Black",
lwd = 0.2,
legend.pos = "topleft",
legend.title.txt = "NBI",
add = TRUE
)
labelLayer(
x = nbi_munic_2,
txt = "MUNICIPIO",
col= "white",
cex = 0.4,
font = 4,
halo = TRUE,
bg = "black",
r = 0.1,
overlap = FALSE,
show.lines = FALSE
)
layoutLayer(
title = "Municipios del departamento de Caqueta",
sources = "Fuente: DANE, 2018",
author = "Sandra Katteryne Rodríguez Hurtado",
frame = TRUE,
north = TRUE,
tabtitle = TRUE,
theme = "taupe.pal"
) 
Figura 3. Mapa del porcentaje de NBI de los municipios de Caquetá. El mapa refleja los datos de la tabla 1, donde cada municipio lleva la etiqueta de su nombre correspondiente, cuando el color se torna más oscuro, mayor será el porcentade de NBI que presenta esa población. Fuente: Autoria propia con datos recopilados por el DANE.
Economía
Dentro de las principales actividades económicas que se encuentran en el departamento sobresale la ganadería y la agricultura, en esta última se encuentran cereales, flores y follajes, forestales, frutales, hortalizas, leguminosas, oleaginosas, otros permanentes, tubérculos y plátanos (DANE, 2017). Se debe tener en cuenta que, pese a que en las tablas no se dan unidades, las de área son hectáreas (Ha), y las unidades de rendimientos son toneladas/Ha.
Parsed with column specification:
cols(
COD_DEP = [32mcol_double()[39m,
DEPARTAMENTO = [31mcol_character()[39m,
COD_MUN = [32mcol_double()[39m,
MUNICIPIO = [31mcol_character()[39m,
GRUPO = [31mcol_character()[39m,
SUBGRUPO = [31mcol_character()[39m,
CULTIVO = [31mcol_character()[39m,
YEAR = [32mcol_double()[39m,
Area_Siembra = [32mcol_double()[39m,
Area_cosechada = [32mcol_double()[39m,
Produccion = [32mcol_double()[39m,
Rendimiento = [32mcol_double()[39m,
ESTADO = [31mcol_character()[39m,
CICLO = [31mcol_character()[39m
)`summarise()` regrouping output by 'GRUPO' (override with `.groups` argument)Tabla 2. Municipios con mayor área cosechada en el año 2018. Se puede apreciar que el municipio San Vicente del Caguán fue el de mayor área sembrada en el grupo de cereales, seguido de Cartagena del Chaira para el grupo de Tubérculos y plátanos, y Florencia, para el grupo de otros permanentes. Fuente: Autoria propia con datos recopilados por el DANE.
Tabla 3. Producción del cultivo de Maiz, perteneciente al grupo de cereales, para el Municipio de San Vicente del Caguán (en toneladas).
g + ggtitle("Evolución de la producción de Maiz en San Vicente del Caguán de 2007 a 2018") + labs(caption= "Basado en los datos de EMA (DANE, 2018)")
Gráfico 1. Evolución de la producción de Maiz en San Vicente del Caguán de 2007 a 2018. Se puede apreciar que para el año 2018 se obtuvo una producción de Maiz, en el municipio de San Vicente del Caguán, de cerca de cinco mil toneladas. La producción máxima de maiz fue para el año de 2013, con cerca de siete mil toneladas, lo anterior se puede corroborar con la información de la tabla 3.
Tabla 4. Producción del cultivo de plátano, perteneciente al grupo de Tubérculos y plátanos, para el Municipio de Cartagena del Chaira (en toneladas).
g + ggtitle("Evolución de la producción de Plátano en Milán de 2007 a 2018") + labs(caption= "Basado en los datos de EMA (DANE, 2018)")
Gráfico 2. Evolución de la producción de Plátano en Milán de 2007 a 2018. Se puede apreciar que la producción de plátano en el municipio de Milán ha ido aumentando en los últimos años, donde en este último año tuvo una producción de cerca de 22 toneladas de plátano.
`summarise()` ungrouping output (override with `.groups` argument)Tabla 5. Rendimientos para cultivos del departamento de Caquetá en 2018. La siguiente tabla muestra los grupos de cultivo con el total de rendimientos (Tonelada/ha) de manera descendente, donde los frutales, tubérculos y plátanos, y hortalizas ocuán los primeros lugares.
g+ ggtitle("Rendimiento por grupos de cultivos en Caquetá para 2018") + theme(plot.title = element_text(hjust = 0.5)) +
labs(caption= "Basado en los datos de EMA (DANE, 2018)")
Gráfico 3. Rendimientos por grupos de cultivos en Caquetá para el año 2018. Se puede apreciar que los grupos de cultivos que presentaron mayores rendimientos para el año 2018 fueron los frutales, tubérculos y plátanos, hortalizas y otros permanentes.
`summarise()` regrouping output by 'CULTIVO' (override with `.groups` argument)Tabla 6. Municipios con mayores rendimientos en Frutales. Se pueden apreciar los municipios con mayores rendimientos (Tonelada/ha) en los cultivos que componen el grupo de Frutales, en el departamento de Caquetá. Se realizó con el grupo de frutales debido a que es el grupo que presentó los mayores rendimientos en el año 2018 en Caquetá. Los municipios de San Vicente del Caguán y Valparaiso tuvieron los mayores rendimientos con los cultivos de Fresa y Arazá respectivamente.
datos %>%
filter(YEAR==2018 & GRUPO=="TUBERCULOS Y PLATANOS") %>%
group_by(CULTIVO, MUNICIPIO) %>%
summarize(Rendimiento = max(Rendimiento, na.rm = TRUE)) %>%
slice(which.max(Rendimiento)) -> Rendimiento`summarise()` regrouping output by 'CULTIVO' (override with `.groups` argument)RendimientoTabla 7. Municipios con mayores rendimientos en Tubérculos y plátanos. Se puede apreciar que los municipios de Cartagena del Chaira, Milán y Solano fueron los que tuvieron mayores rendimientos en los cultivos de yuca, plátano y Ñame respectivamente.
png("D:/Documentos/Geomática/platano_2018p.png", width = 1024, height = 763)
opar <- par(mar = c(2,2,5,5))
plot(st_geometry(rep_platano), col="darkseagreen", border="darkseagreen4",
bg = "white", lwd = 0.6)
propSymbolsChoroLayer(x = rep_platano, var = "Produccion", var2 = "Rendimiento",
col = carto.pal(pal1 = "turquoise.pal", n1 = 3,
pal2 = "pink.pal", n2 = 3),
inches = 0.5, method = "q6",
border = "black", lwd = 1,
legend.title.cex = 1.3,
legend.values.cex = 1.3,
legend.var.pos = "topright",
legend.var2.pos = "topleft",
legend.var2.values.rnd = 2,
legend.var2.title.txt = "Rendimiento\n(en Ton/Ha)",
legend.var.title.txt = "Producción de Plátano, 2018",
legend.var.style = "e")
labelLayer(
x = rep_platano,
txt = "MPIO_CNMBR",
col= "black",
cex = 1,
font = 18,
halo = FALSE,
bg = "white",
r = 0.1,
overlap = FALSE,
show.lines = FALSE
)
layoutLayer(title="Producción y rendimiento de Plátano en Caquetá, 2018",
author = "Sandra Katteryne Rodríguez Hurtado",
sources = "Fuentes: MADR & DANE, 2018",
scale = 50, tabtitle = FALSE, frame = TRUE)
north(pos = "top")
title(main="Producción y rendimiento de Plátano en Caquetá, 2018", cex.main=3,
sub= "Fuentes: MADR & DANE, 2018", cex.sub=2)
graticule = TRUE
par(opar)
dev.off()
Producción y rendimiento de plátano en Caquetá para el año 2018
Figura 4. Producción y rendimiento de plátano en Caquetá para el año 2018 por municipios. Se pueden apreciar que los municipios de mayor producción de plátano son Milán, San Vicente del Caguán y Cartagena del Chaira.
Descripción de datos y métodos
Datos
Carreras tales como la ingeniería Agronómica necesitan recolectar datos en puntos específicos y partiendo de ellos hacer representaciones continuas de una zona determinada, y para ello es necesario hacer interpolación (Lizarazo, 2020). La interpolación o modelación, toma datos recolectados en puntos finitos específicos para volverlos datos continuos, representaciones de mapas útiles y con sentido, estos últimos puntos predicen características similares a las de los puntos muestreados para áreas no muestreadas, predicen valores para celdas de un ráster (ArcGIS, s.f). Las superficies ráster presentan una forma de cuadrícula, donde a menor tamaño de las celdas, mayor será su precisión.
Para poder interpolar datos, estos deben ser cuantitativos, analizando variables como pH, precipitación, humedad, índice de vegetación, no cuantitativos. Los datos pueden ser recolectados en una estación meteorológica, en el caso de Colombia estas estaciones no son muchas y no abarcan todo el territorio nacional, pues algunas se encuentran en los principales aeropuertos y otras pertenecen a gremios, donde el acceso no es público (Lizarazo, 2020).
En este trabajo se usaron los datos del CHIRPS (Climate Hazards Group InfraRed Precipitation with Station Data), que son datos de estimaciones globales de precipitación, basadas en satélites de la NASA y NOAA, y observaciones sobre la tierra. Los datos tienen una resolución de 0,05 y 0,25 grados (5,5 y 28 km) en formatos BIL, TIF o NetCDF (Grupo Kraken, 2019). Se expanden de 50°S a 50°N y a lo largo de todas las latitudes, comprendiendo datos desde 1981 hasta la actualidad, tienen alta resolución y un largo periodo de registro.
Tienen la capacidad de incorporar imágenes satelitales a datos obtenidos de estaciones meteorológicas para crear series de cuadrículas pluviométricas con unidades en milímetros, que permiten realizar estudios como tendencias o sequías (IRI, 2015). CHIRPS incorpora datos de estimaciones de precipitación diarias, pentadales y mensuales; hay seis pentads en un mes calendario, cinco pentads de 5 días y un pentad con los 3 a 5 días restantes del mes (Funk et al., 2015). Los datos fueron descargados del siguiente link http://data.chc.ucsb.edu/products/CHIRPS-2.0/, directorio global_pentad, posteriormente entrar a tifs y seleccionar la fecha más reciente, que fue el 18/05/2020.
Métodos
Existen dos tipos de técnicas para generar superficies de interpolación: modelos deterministas, que usan funciones matemáticas para predecir valores desconocidos, y técnicas estadísticas, que producen límites para determinar la precisión de una predicción. Se explicarán tres métodos diferentes para interpolar datos, siendo IDW, primer orden polinomial y Kriging.
Inverse Distance Weighting (IDW)
El método de IDW (Inverse Distance Weighting) pertenece a la técnica determinista, es matemático, donde los valores más cercanos están más relacionados con la variable que se esté estudiando, es decir que los puntos más cercanos comparten características similares a diferencia de los puntos lejanos (GIS Geography, s.f.). Por ejemplo, al llover en un punto determinado, se puede predecir que cerca de este punto también está lloviendo, a diferencia de un punto lejano, donde no se podría hacer esta afirmación (autocorrelación espacial). Se pueden estimar valores desconocidos de celda raster a partir de valores conocidos tomados en puntos de muestra.
Al tener un área de estudio y determinados puntos de muestra con valores conocidos, si se quisiera medir un punto en específico, es necesario tomar los puntos más cercanos dentro del radio de la zona de interés para poder realizar la interpolación, puede ser un número fijo o variable de puntos (GIS Geography, s.f.). Los puntos desconocidos se estiman partiendo de los valores conocidos de los puntos muestreados, el valor de potencia le permite controlar la significancia de los puntos desconocidos (ArcGIS, s.f.), donde el valor de 2 es el más usado. Su fórmula es: 
Figura 5. Fórmula de IDW. Esta fórmula tiene el número sigma, que es el número de puntos que se interpolarán, z: Variable de estudio, zp: punto desconocido, n: número de puntos que están dentro del radio de búsqueda, zi: valor conocido del punto en el radio de búsqueda, di^p: Distancia entre el punto y el punto desconocido, p: Potencia. Se considera determinístico ya que está basado en los valores medidos y en fórmula matemática, que determinan la suavidad de la superficie generada (ArcGIS, s.f).
Kriging
Es un método de interpolación geoestadístico para crear superficies estimadas en base a unos puntos con valores conocidos de la variable de estudio, se basa en modelos estadísticos, es decir que los puntos muestreados se relacionan con la estadística. Lo anterior permite que resulte una superficie de predicción y medidas de precisión de las estimaciones (ArcGIS, s.f).
Al igual que IDW, el kriging utiliza valores de puntos muestreados dentro de un radio, ajustándolos a una función matemática. El kriging se basa en un variograma que determinan la autocorrelación, acorde con el geógrafo Tobler: “Todo está relacionado con todo lo demás, pero las cosas cercanas están más relacionadas que las distantes”, esto quiere decir la autocorrelación es el grado de similitud entre un objeto y los objetos cercanos al mismo (GIS Geography, 2020), representando gráficamente la varianza de los datos acorde con su distancia. Si los puntos están más cercanos, tendrán una varianza pequeña, mientras que si los puntos están más alejados, probablemente se tendrá una varianza alta. De lo anterior se puede hacer una gráfica que tiene un rango y un umbral, donde el último se alcanza al nivelarse la variación, no habría autocorrelación espacial entre los datos (GIS Geography, 2020).
Figura 6. Fórmula de Kriging. z(si): Valor medido en la ubicación i, Lamda i: ponderación desconocida para el valor medido en la ubicación i, s0: Ubicación de la predicción, N: cantidad de valores medidos. Tomado de: https://desktop.arcgis.com/es/arcmap/10.3/tools/spatial-analyst-toolbox/how-kriging-works.htm
Un semivariograma puede ser superficial, exponencial, circular, gaussiano o lineal; esto después de haber ajustado la superficie como un polinomio. Para usar este método, los datos deben tener una distribución normal, deben ser estacionarios y no tener ninguna tendencia (ArcGIS, s.f). Se dice que los datos tienen una distribución normal si los puntos están sobre una línea de 45°, los datos son estacionarios si la variación local de estos no cambia en diferentes áreas del mapa, es decir que no hay cambios abruptos o líneas discontinuas, se puede verificar mediante un mapa de voronoi, que refleja la desviación estándar. Según las predicciones que la técnica de Kriging haga, generará también una medida de incertidumbre o error, estimando la confianza de la superficie.
Los gráficos de los semivariogramas pueden ser calculados con la ecuación:
Semivariograma(distanciah) = 0.5 * promedio((valori – valorj)^2)dando como resultado un gráfico de los valores de semivariograma promediados ubicados en el eje Y, y la distancia o intervalo en el eje X. Entonces, la ubicación de los puntos proyectados en este gráfico permite conocer si son valores similares o no, si los puntos están uno cerca del otro, significa que comparten características similares, a diferencia de si se encuentran separados uno del otro, significando que son distintos (ArcGIS, s.f).
Polígonos de Thiessen o Diagramas de Voronoi
Son uno de los métodos para interpolar más sencillos ya que hay, están basados en una construcción geométrica que permiten hacer particiones de un plano que representará el área de interés. Su nombre fue puesto en honor al meteorólogo Alfred Thiessen, aunque también fueron estudiados por el matemático Georgy Voronoi, método denominado también Diagramas de Voronoi en honor al matemático (Bravo, s.f).
Este método tiene múltiples funciones, dentro de ellas se encuentra el análisis de datos meteorológicos o determinar áreas de influencia en una zona, por ejemplo iglesias o centros comerciales. Para el caso de precipitación, es necesario saber dónde se encuentran ubicadas las estaciones pluviométricas para identificar el área de influencia de dicha estación. Al unirse los puntos entre sí, se van trazando a su vez mediatrices (líneas perpendiculares a un segmento que pasa por su punto medio) de los segmentos de unión, creando una serie de polígonos en un espacio de dos dimensiones alrededor de un conjunto finito de puntos (Bravo, s.f). Forman triángulos que unen las estaciones, por medio de líneas rectas sin que estas se entrecortan (Gimond, 2019). Una limitante que presenta este método es que no se puede estimar el valor de incertidumbre o error debido a que un solo polígono se obtiene de un solo punto. El tamaño y la forma de cada polígono depende de la distribución de los puntos. Estas estructuras son de vital importancia en la geometría computacional, son capaces de guardar información referente a una proximidad entre puntos.
Presentación de resultados
Una vez descargados los datos y cargadas las librerías, se procedió a leerlos, presentaron una resolución de celda de 0.05 y fue una capa ráster. Son datos globales, tienen un sistema de coordenadas geográficas WGS84.
class : RasterLayer
dimensions : 2000, 7200, 14400000 (nrow, ncol, ncell)
resolution : 0.05, 0.05 (x, y)
extent : -180, 180, -50, 50 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
source : D:/Escritorio/chirps-v2.0.2020.04.6.tif
names : chirps.v2.0.2020.04.6 Se importará un archivo de formato shapefile para el área que se desea estudiar, en este caso el departamento de Caquetá, el cual se compone de 16 municipios, sistema de coordenadas WGS84 y nueve variables.
class : SpatialPolygonsDataFrame
features : 16
extent : -76.30622, -71.25385, -0.70584, 2.964148 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 9
names : DPTO_CCDGO, MPIO_CCDGO, MPIO_CNMBR, MPIO_CRSLC, MPIO_NAREA, MPIO_NANO, DPTO_CNMBR, Shape_Leng, Shape_Area
min values : 18, 18001, ALBANIA, Decreto 1678 de Septiembre 7 de 1967, 403.35827433, 2017, CAQUETÃ, 1.1128288978, 0.0327393194813
max values : 18, 18860, VALPARAÃSO, Ordenanza 3 de Noviembre 12 de 1985, 42312.7997376, 2017, CAQUETÃ, 20.0458589874, 3.43579026269 Para dar valores nulos a las áreas que se encuentran fuera de los límites del departamento se usó la función mask, siendo una capa raster.
class : RasterLayer
dimensions : 19, 26, 494 (nrow, ncol, ncell)
resolution : 0.2, 0.2 (x, y)
extent : -76.4, -71.2, -0.8000008, 2.999999 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
source : memory
names : chirps.v2.0.2020.04.6
values : 18.95572, 73.44761 (min, max)plot(precip.mask, main= "CHIRPS de Precipitación en Caquetá desde 18/05/2020 hasta 23/05/2020 [mm]")
plot(aoi, add=TRUE)
Figura 7. Mapa de vista del precipitación del departamento de Caquetá.
pal <- colorNumeric(c("red", "orange", "yellow", "blue", "darkblue"), values(precip.mask),
na.color = "transparent")
leaflet() %>% addTiles() %>%
addRasterImage(precip.mask, colors = pal, opacity = 0.6) %>%
addLegend(pal = pal, values = values(precip.mask),
title = "CHIRPS de Precipitación en Caquetá\n desde 18/05/2020 hasta 23/05/2020 [mm]")Some values were outside the color scale and will be treated as NAFigura 8. Mapa interactivo de precipitación en Caquetá.
class : SpatialPointsDataFrame
features : 183
extent : -76.1, -71.3, -0.5000008, 2.899999 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 1
names : chirps.v2.0.2020.04.6
min values : 18.9557183384895
max values : 73.4476075172424 class : SpatialPointsDataFrame
features : 183
extent : -76.1, -71.3, -0.5000008, 2.899999 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 1
names : Lluvia
min values : 18.9557183384895
max values : 73.4476075172424 Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
..@ data :'data.frame': 183 obs. of 1 variable:
.. ..$ Lluvia: num [1:183] 34.7 37 44.5 47.8 49.4 ...
..@ coords.nrs : num(0)
..@ coords : num [1:183, 1:2] -74.9 -75.1 -74.9 -74.7 -75.1 ...
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : NULL
.. .. ..$ : chr [1:2] "x" "y"
..@ bbox : num [1:2, 1:2] -76.1 -0.5 -71.3 2.9
.. ..- attr(*, "dimnames")=List of 2
.. .. ..$ : chr [1:2] "x" "y"
.. .. ..$ : chr [1:2] "min" "max"
..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
.. .. ..@ projargs: chr "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"Plotear los puntos para llevar a cabo la interpolación.
plot(precip.mask, main= "CHIRPS Lluvia desde 18/05/2020 hasta 23/05/2020 [mm]")
plot(aoi, add=TRUE)points(precip.points$x, precip.points$y, col = "purple", cex = .1)
Figura 9.Departamento de Caquetá con los puntos para la interpolación.
class : SpatialPointsDataFrame
features : 183
extent : -76.1, -71.3, -0.5000008, 2.899999 (xmin, xmax, ymin, ymax)CRS object has comment, which is lost in outputcrs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 1
names : Lluvia
min values : 18.9557183384895
max values : 73.4476075172424 Es necesario pasar los datos a coordenadas planas porque no es recomendable hacer interpolación con los datos en coordenadas geográficas (datumWGS84).
class : SpatialPolygonsDataFrame
features : 16
extent : -76.30622, -71.25385, -0.70584, 2.964148 (xmin, xmax, ymin, ymax)
crs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
variables : 9
names : DPTO_CCDGO, MPIO_CCDGO, MPIO_CNMBR, MPIO_CRSLC, MPIO_NAREA, MPIO_NANO, DPTO_CNMBR, Shape_Leng, Shape_Area
min values : 18, 18001, ALBANIA, Decreto 1678 de Septiembre 7 de 1967, 403.35827433, 2017, CAQUETÃ, 1.1128288978, 0.0327393194813
max values : 18, 18860, VALPARAÃSO, Ordenanza 3 de Noviembre 12 de 1985, 42312.7997376, 2017, CAQUETÃ, 20.0458589874, 3.43579026269 class : SpatialPolygonsDataFrame
features : 1
extent : -76.30622, -71.25385, -0.70584, 2.964148 (xmin, xmax, ymin, ymax)CRS object has comment, which is lost in outputcrs : +proj=longlat +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +no_defs
variables : 1
names : area
value : 90103.00816646 Polígonos
Simple feature collection with 16 features and 2 fields
geometry type: POLYGON
dimension: XY
bbox: xmin: -76.30622 ymin: -0.70584 xmax: -71.25385 ymax: 2.964148
CRS: +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
First 10 features:
MUNIC CODIGO
1 FLORENCIA 18001
2 ALBANIA 18029
3 BELÉN DE LOS ANDAQUÃ\u008dES 18094
4 EL DONCELLO 18247
5 EL PAUJIL 18256
6 LA MONTAÑITA 18410
7 MILÃ\u0081N 18460
8 MORELIA 18479
9 SAN JOSÉ DEL FRAGUA 18610
10 VALPARAÃ\u008dSO 18860
geometry
1 POLYGON ((-75.42074 2.19413...
2 POLYGON ((-75.89506 1.36569...
3 POLYGON ((-75.78705 1.74982...
4 POLYGON ((-75.36167 2.32142...
5 POLYGON ((-75.36691 2.21234...
6 POLYGON ((-75.40404 1.76944...
7 POLYGON ((-75.39362 1.35738...
8 POLYGON ((-75.77185 1.57991...
9 POLYGON ((-76.16722 1.58752...
10 POLYGON ((-75.73128 1.32740...Una intersección de datos es una operación geométrica entre dos elementos distintos, cuyo resultado busca que por cada elemento que se intercepte y arroje un elemento, donde queden los atributos de ambos objetos en uno solo. En este caso, en un elemento quedaron consignados los elementos municipio (polígonos) y la precipitación (puntos.)
although coordinates are longitude/latitude, st_intersection assumes that they are planar
attribute variables are assumed to be spatially constant throughout all geometriesSimple feature collection with 183 features and 3 fields
geometry type: POINT
dimension: XY
bbox: xmin: -76.1 ymin: -0.5000008 xmax: -71.3 ymax: 2.899999
CRS: +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
First 10 features:
MUNIC CODIGO Lluvia
12 SAN VICENTE DEL CAGUÃ\u0081N 18753 34.66168
12.1 SAN VICENTE DEL CAGUÃ\u0081N 18753 36.99089
12.2 SAN VICENTE DEL CAGUÃ\u0081N 18753 44.50369
12.3 SAN VICENTE DEL CAGUÃ\u0081N 18753 47.77253
16 PUERTO RICO 18592 49.42595
12.4 SAN VICENTE DEL CAGUÃ\u0081N 18753 48.35983
12.5 SAN VICENTE DEL CAGUÃ\u0081N 18753 48.97882
16.1 PUERTO RICO 18592 32.69428
16.2 PUERTO RICO 18592 51.43533
12.6 SAN VICENTE DEL CAGUÃ\u0081N 18753 51.31482
geometry
12 POINT (-74.9 2.899999)
12.1 POINT (-75.1 2.699999)
12.2 POINT (-74.9 2.699999)
12.3 POINT (-74.7 2.699999)
16 POINT (-75.1 2.499999)
12.4 POINT (-74.9 2.499999)
12.5 POINT (-74.7 2.499999)
16.1 POINT (-75.3 2.299999)
16.2 POINT (-75.1 2.299999)
12.6 POINT (-74.9 2.299999)Dos tareas de reprotección: Dos maneras para reproyección, convirtiendo un objeto de característica simple (simplefeature), a un objeto espacial, pues las librerías requieren estos objetos.
class : SpatialPointsDataFrame
features : 183
extent : 774851.6, 1309312, 436453.5, 812468.8 (xmin, xmax, ymin, ymax)CRS object has comment, which is lost in outputcrs : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
variables : 3
names : MUNIC, CODIGO, Lluvia
min values : ALBANIA, 18001, 18.9557183384895
max values : VALPARAÃSO, 18860, 73.4476075172424 Eliminar cifras decimales de la variable Lluvia, que representa la preciítación en Caquetá.
#Esta instrucción es útil para que los valores de lluvia no reflejen tantas cifras decimales en el mapa.
precip2$Lluvia <- round(precip2$Lluvia, 1)class : SpatialPointsDataFrame
features : 183
extent : 774851.6, 1309312, 436453.5, 812468.8 (xmin, xmax, ymin, ymax)CRS object has comment, which is lost in outputcrs : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
variables : 3
names : MUNIC, CODIGO, Lluvia
min values : ALBANIA, 18001, 19
max values : VALPARAÃSO, 18860, 73.4 class : SpatialPolygonsDataFrame
features : 16
extent : 751894.5, 1314457, 413692, 819562.9 (xmin, xmax, ymin, ymax)CRS object has comment, which is lost in outputcrs : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
variables : 2
names : MUNIC, CODIGO
min values : ALBANIA, 18001
max values : VALPARAÃSO, 18860 Trazar los datos con tmap: La libreía tmap permite hacer mapas, en este caso con los datos de precipitación y municipio del departamento de Caquetá.
tm_shape(caqueta2) + tm_polygons() +There were 13 warnings (use warnings() to see them) tm_shape(precip2) +
tm_dots(col="Lluvia", palette = "RdBu", midpoint = 3.0,
title="Precipitación muestreada \n(en mm)", size=0.2) +
tm_text("Lluvia", just="center", xmod=.1, size = 0.5) +
tm_legend(legend.outside=TRUE)CRS object has comment, which is lost in outputCRS object has comment, which is lost in output
Figura 10. Mapa con los puntos y sus respectivos valores de precipitación, en este caso variando en intervalos de a 10 mm. ### Métodos para interpolar datos de precipitación A continuación se representarán algunos de los métodos trabajados para la interpolación de datos de precipitación. #### 1. Polígonos de Thiessen Estos polígonos son (o polígnos de proximidad) se pueden crear utilizando la función dirichlet de spatstat. Arman unos polígonos que contienen todas las zonas que están más cerca de un punto que de otro
crs(th) <- crs(precip2)
crs(caqueta2) <- crs(precip2)CRS arguments:
+proj=tmerc +lat_0=4.59620041666667
+lon_0=-74.0775079166667 +k=1 +x_0=1000000
+y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0
+units=m +no_defs CRS object has comment, which is lost in outputCRS arguments:
+proj=tmerc +lat_0=4.59620041666667
+lon_0=-74.0775079166667 +k=1 +x_0=1000000
+y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0
+units=m +no_defs Mapa
Figura 11. Mapa de polígonos de Thiessen. #### Interpolación ponderada de distancia inversa (IDW). La función IDW está presente en dos paquetes spatstat y gstat. Se debe crear primero una grilla con celdas que aún no tienen valores, esa grilla tiene 100.000 celdas, las cuales van a cubrir un territorio determinado, este último será el ráster de salida.
CRS object has comment, which is lost in outputCRS object has comment, which is lost in output[inverse distance weighted interpolation]class : RasterLayer
dimensions : 224, 311, 69664 (nrow, ncol, ncell)
resolution : 1806.052, 1806.052 (x, y)
extent : 752544.8, 1314227, 414359.5, 818915.1 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
source : memory
names : var1.pred
values : 19.12434, 72.96179 (min, max)class : SpatialPolygonsDataFrame
features : 16
extent : 751894.5, 1314457, 413692, 819562.9 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
variables : 2
names : MUNIC, CODIGO
min values : ALBANIA, 18001
max values : VALPARAÃSO, 18860 class : RasterLayer
dimensions : 224, 311, 69664 (nrow, ncol, ncell)
resolution : 1806.052, 1806.052 (x, y)
extent : 752544.8, 1314227, 414359.5, 818915.1 (xmin, xmax, ymin, ymax)
crs : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
source : memory
names : var1.pred
values : 19.12434, 72.96179 (min, max)Plotear.
tm_shape(r.m) +
tm_raster(n=8,palette = "RdBu", auto.palette.mapping = FALSE,
title="Distancia Inversa Poderada\nPrecipitación prevista\n(en mm)") +
tm_shape(precip2) + tm_dots(size=0.2) +
tm_legend(legend.outside=TRUE)The argument auto.palette.mapping is deprecated. Please use midpoint for numeric data and stretch.palette for categorical data to control the palette mapping.CRS object has comment, which is lost in output
Figura 12. Mapa de IDW sin la vista de los valores.
pal <- colorNumeric(c("red", "orange", "yellow", "blue", "darkblue"), values(precip.mask),
na.color = "transparent")
leaflet() %>% addTiles() %>%
addRasterImage(r.m, colors = pal, opacity = 0.6) %>%
addLegend(pal = pal, values = values(r.m),
title = "Precipitación interpolada IDW en Caquetá\n desde 18/05/2020 hasta 23/05/2020 en [mm]")Figura 13. Mapa de interpolación IDW de precipitación.
class : SpatialPointsDataFrame
features : 183
extent : 751894.5, 1314457, 413692, 819562.9 (xmin, xmax, ymin, ymax)CRS object has comment, which is lost in outputcrs : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
variables : 3
names : MUNIC, CODIGO, Lluvia
min values : ALBANIA, 18001, 19
max values : VALPARAÃSO, 18860, 73.4 Cálcular el error en los datos interpolados.
OP <- par(pty="s", mar=c(4,3,0,0))
plot(IDW.out ~ P$Lluvia, asp=1, xlab="Observada", ylab="Prevista", pch=16,
col=rgb(0,0,0,0.5))
abline(lm(IDW.out ~ P$Lluvia), col="red", lw=2,lty=2) abline(0,1)
Gráfico 3. Semivariograma
Error medio cuadrático.
sqrt( sum((IDW.out - P$Lluvia)^2) / length(P))[1] 6.650464Validación cruzada.
r <- raster(img.sig, layer="v")
r.m <- raster::mask(r, caqueta2)
tm_shape(r.m) + tm_raster(n=8,title="IDW\n95% intervalo de confianza \n(en mm)") +
tm_shape(P) + tm_dots(size=0.2) +
tm_legend(legend.outside=TRUE)CRS object has comment, which is lost in output
Figura 14. Mapa de intervalo de confianza de 95% del modelo de interpolación.
The argument auto.palette.mapping is deprecated. Please use midpoint for numeric data and stretch.palette for categorical data to control the palette mapping.CRS object has comment, which is lost in output
Kriging
# The following plot allows us to assess the fit
plot(var.smpl, dat.fit, xlim=c(0,130000))
Gráfico 4. Semivariograma.
r.mMapa.
# Plot the map
tm_shape(r.m) +
tm_raster(n=10, palette="RdBu", auto.palette.mapping=FALSE,
title="Kriging Universal\nPrecipitación prevista \n(en mm)") +
tm_shape(P) + tm_dots(size=0.2) +
tm_legend(legend.outside=TRUE)The argument auto.palette.mapping is deprecated. Please use midpoint for numeric data and stretch.palette for categorical data to control the palette mapping.CRS object has comment, which is lost in output
Figura 15. Mapa de kriging de precipitación.
library(leaflet)
library(RColorBrewer)
pal <- colorNumeric(c("red", "orange", "yellow", "blue", "darkblue"), values(precip.mask),
na.color = "transparent")
leaflet() %>% addTiles() %>%
addRasterImage(r.m, colors = pal, opacity = 0.6) %>%
addLegend(pal = pal, values = values(r.m),
title = "Kriging Precipitación Interpolada en Caquetá\n desde 18/05/2020 hasta 23/05/2020 [mm]")Some values were outside the color scale and will be treated as NASome values were outside the color scale and will be treated as NAFigura 16. Mapa interactivo Kriging de la precipitación interpolada. Generar los mapas de varianza e intervalo de confianza.
r <- raster(dat.krg, layer="var1.var")
r.m <- raster::mask(r, caqueta2)
tm_shape(r.m) +
tm_raster(n=5, palette ="Reds",
title="Interpolación de Kriging\nMapa de variación \n(en squared mm)") +tm_shape(P) + tm_dots(size=0.2) +
tm_legend(legend.outside=TRUE)CRS object has comment, which is lost in output
Figura 17. Mapa de variación en Kriging.
r <- sqrt(raster(dat.krg, layer="var1.var")) * 1.96There were 23 warnings (use warnings() to see them)r.m <- raster::mask(r, caqueta2)
tm_shape(r.m) +
tm_raster(n=5, palette ="Reds",
title="Interpolación Kriging\n Mapa de intervalo de confianza al 95% \n(en mm)") +tm_shape(P) + tm_dots(size=0.2) +
tm_legend(legend.outside=TRUE)CRS object has comment, which is lost in output
Figura 18. Mapa de intervalo de confianza de kriging.
sessionInfo()Análisis de resultados
En base a los resultados anteriores, se puede ver que los métodos para la interpolación de datos presentan diferencias y unos son mejores que otros dependiendo del tipo de variable que se quiera estudiar. También presentan diferencias en cuanto a la manera de calcularlos, de tomar las medidas, creación de las estructuras y resolución de los mapas.
En las figuras 7 y 8 se puede apreciar una vista del departamento de caquetá, se puede ver una capa ráster compuesta de varias celdas de diferentes colores, donde cada color representa la cantidad de precipitación en milímetros (mm) que tuvo el departamento, donde los municipios de Cartagena del Chaira, San José del Fragua, Valparaíso y Milan fueron los que presentaron mayor precipitación, a comparación de los municipios de El Paujil, Florencia y San Vicente del Caguán, donde sus niveles de precipitación fueron menores. En ambos mapas se pueden ver una serie de celdas de forma cuadrada cuya resolución es de 0.2 * 0.2. Algunas diferencias entre ambos mapas es que el de la figura 7 tiene una vista del departamento son los límites de cada uno de sus municipios, no es interactivo y su color de celdas varía de rojo a verde; mientras que el mapa de la figura 8 es interactivo (OpenStreetMap) y más detallado, permite identificar el departamento desde una vista del mapa de Colombia, en el que se pueden apreciar con cuáles departamentos limita, una vista simple de las zonas en las que hay presencia de vegetación, maneja una paleta de color que va del rojo al azul y no muestra los límites entre sus municipios. La resolución de las celdas que presentan, no es muy precisa, sin embargo se trabajó con esta resolución debido a que más adelante se requerirá, para la vista de los mapas, puntos no tan concentrados. Requirió hacer varias operaciones para volver los datos espaciales. En la figura 9 se pueden apreciar las celdas de la capa raster con un punto color morado en el centro de cada una, estos serán los puntos para hacer la interpolación.
Los datos también se trabajaron por medio del formato Geojson, el cual está diseñado para generar, visualizar y editar datos geográficos (geodatos)GeoJson.io es una herramienta pública y virtual que permite importar y exportar geodatos, y una de sus funciones es reproyectar capas de geodatos (Viloria, s.f). Los puntos fueron pasados a polígonos con el fin de convertir la característica precipitación en un marco de datos espacial. Vale aclarar que los valores de los datos de precipitación presentaron muchas cifras decimales, lo que dificulta la visualización de los puntos en el mapa, para ello fue necesario redondear las cifras para que apareciera solo un decimal.
La figura 9 representa el mapa de Caquetá con intervalos de precipitación de 10 mm. Cada punto representa el nivel de precipitación muestreada en cada uno de los puntos, variando de blanco claro a oscuro, donde a mayor nivel de precipitación, la tonalidad azul de cada punto se va tornando más oscura. Los polígonos de Thiessen se pueden apreciar en la figura 11 polígonos cuadrados, esto debido a que la posición de sus puntos tiene la misma distancia en el eje X y Y. Los polígonos formados en este mapa no cumplen el propósito de este este método ya que los polígonos están distribuidos de forma regular, entonces no proporciona una información clara respecto de la lluvia en el departamento. Algo que se puede deducir del mapa es que este departamento presenta precipitaciones que cubren toda la zona, esto en concordancia con los mapas anteriores.
La interpolación de distancia inversa (IDW) se puede apreciar en la figura 12, donde se ven celdas ráster que están sobrepuestas en los puntos de muestreo, sin embargo no se ve la presencia de los valores (expresados en números) sobre los puntos. Presenta una paleta de colores que varía entre el rojo y el azul, donde los valores más cercanos, están más relacionados con la variable, que es la precipitación, es decir que los puntos más cercanos comparten características similares a diferencia de los puntos lejanos. Lo anterior se puede ver reflejado en el agrupamiento de puntos por colores que varían entre tonalidades azules y rojas, unas más claras claras y otras oscuras, cada agrupamiento significa que los puntos que se encuentran dentro del mismo comparten características similares. Por ejemplo se pueden ver tonalidades azules oscuras en la zona a la que pertenecen los municipios de Cartagena del Chaira, Solano y Milán, presentando mayores precipitaciones en esta zona, el los puntos que se encuentran agrupados en el color azul oscuro, tienen en común que hay precipitaciones fuerte allí, y a medida que se van alejando, los puntos se van agrupando en coloraciones más claras, lo que significa que comparten la característica de que hay precipitaciones, pero más leves. Los agrupamientos de color rojizo, significan que los puntos contenidos ahí, comparten la característica de una baja precipitación, perteneciendo a los municipios de Florencia, Belén de los Andaquíes, una parte de San Vicente del Caguán y el sur del municipio de Solano. La figura 13 presenta una visión más detallada del mapa, con las características de la figura 8, pero esta vez representando el método IDW.
El método IDW permite calcular el valor de error cuadrático medio (RMSE) de los datos interpolados. El gráfico 3, muestra cuál es el valor de precipitación observada en el eje X y cuál fue la prevista en el eje Y. Si los datos coincidieran perfectamente, estarían alineados con la línea de color negro, que es de 45°; sin embargo los resultados muestran que los puntos están sobre la línea roja, la cual puede variar si el exponente fuera distinto a dos (2.2, 2,4), pero no muy lejano a dos. Se trabajó con una potencia IDW de 2, que es la que generalmente se usa.
Al calcular el error medio cuadrático de los datos, arrojó un resultado de 6.650, y proviene de la fórmula:
(sum((IDW.out - P$Lluvia)^2) / length(P))Sum: Sumatoria, IDW.out: Longitud del vector, P$Lluvia: valor de precipitación en el punto length(P): Longitud.
La validación cruzada puede determinar cuál es la incertidumbre de los resultados, donde se creó un mapa de intervalo de confianza de 95% del modelo de interpolación, que corresponde a la figura 14. Esta figura permite apreciar que en la mayoría de los puntos interpolados del departamento de Caquetá, el rango de incertidumbre es muy bajo, tal como lo indica el color. Entonces, se puede decir en base al mapa, que los intervalos de error más alto se presentan en una pequeña área del municipio de Solano.
En el método de Kriging, se calculó un semivariograma lineal porque se usa la fórmula de primer orden polinomial, siendo precip=línea+aX+bY. Para ajustar los valores se realiza el semivariograma, dando como resultado un gráfico donde los valores de semivarianza promediados, los cuales están ubicados en el eje Y, y la distancia o intervalo está ubicada en el eje X, además tiene un nugget igual a cero. Se pueden ver algunos puntos sobre una recta, habiendo una correlación espacial entre dichos puntos, donde están más relacionados entre sí los puntos que se encuentran uno más cerca de otro. Las diferencias no presenta variaciones altas, entonces ese modelo de variograma se ajusta con los puntos de muestreo.
La figura 15 presenta un mapa con mejor resolución, una paleta de color que varía entre el rojo y el azul, reflejando zonas con menor precipitación, y zonas con mayor precipitación respectivamente. Los datos están agrupados por colores según la similitud de precipitación que estos presenten, como se había mencionado en anteriormente, las agrupaciones en tonalidades de color azul son similares en cuanto al nivel de precipitación en el departamento, y las zonas con coloraciones rojizas son aquellas en las que los puntos son parecidos en cuanto a su nivel de precipitación, presentando intervalos menores que los de la zona azul.
La figura 16 muestra un mapa interactivo de la precipitación interpolada por el método de Kriging, presentando representando la información de la figura 15, pero en esta ocasión es el mapa con mejor resolución en comparación de todos los anteriores, esto se dice porque no se ven las celdas tan marcadas, presenta una coloración que varía de una forma adecuada y coherente.
La figura 17 permite observar el intervalo de variación respecto a los datos obtenidos, donde se puede ver que este intervalo presenta valores mayores en los bordes del departamento, permitiendo deducir que sus valores tienen cierto grado de incertidumbre, mientras que los datos cubren el interior del departamento presentan un intervalo de variación entre 10 y 20 mm, pero menor a comparación de los límites.
La figura 18 representa el intervalo de confianza de los datos de precipitación a través de la interpolación por el método de kriging. Se puede apreciar que la mayoría de sus resultados presentan coloraciones correspondientes a intervalos entre 6 y 12 mm, donde el nivel de confianza es favorable.
Conclusiones
La interpolación es una herramienta de vital importancia para muchas carreras debido a que permite realizar estudios de diferentes variales como la precipitación, que fue la que se trabajó en el informe, sin tener que ir directamente a la zona que se quiere estudiar, representa también importancia en cuanto a costos ya que se pueden obtener y representar datos de diferentes variables partiendo de un número limitado de puntos que termiten crear estructuras ráster.
Existen varios métodos para poder realizar interpolación, presentando diferencias como la manera de calcularlos, las estructuras que forman, las variables a las que cada uno se acomoda mejor para estudiar, la manera de interpretar y leer los datos, algunos proporcionan información adicional como la incertidumbre de los datos.
La herramienta RStudio ofrece a sus usuarios una amplia gama de librerías con diferentes funciones, brindando diferentes alternativas para representar, en base a datos descargados de plataformas confiables, tablas, gráficos, realizar mapas que representen gráficamente valores de los gráficos, tales como mapas temáticos, interactivos, entre otros, permitiendo hacer análisis de todo tipo, por ejemplo análisis de las necesidades básicas insatisfechas de una región, nivel de escolaridad, nivel de natalidad, niveles de producción y rendimiento de determinados cultivos, entre otros.
Dentro de las herramientas de R se encuentra la interpolación de datos, permitiendo realizar análisis de precipitación y así poder realizar, en base a literatura y estadísticas, estimaciones sobre impactos de sequías o inundaciones.
Bibliografía
Bravo, L. s.f. Diagrama de Voronoi. Recuperado de: http://matematicas.unex.es/~trinidad/mui/voronoi.pdf;consulta: Junio de 2020.
Funk, C. P, Peterson. M, Landsfeld. D, Pedreros. J, Verdin. S, Shukla. G, Husak.J, Rowland. L, Harrison. A, Hoell y J Michaelsen. 2015. The climate hazards infrared precipitation with stations—a new environmental record for monitoring extremes. Recuperado de: https://www.nature.com/articles/sdata201566; consulta: Junio de 2020.
Gimond, M. 2019. Intro to GIS and Spatial Analysis. Recuperado de: https://mgimond.github.io/Spatial/interpolation-in-r.htm; consulta: Junio de 2020.
GISGeography. s.f.Inverse Distance Weighting (IDW) Interpolation. Recuperado de: https://gisgeography.com/inverse-distance-weighting-idw-interpolation/; consulta: Junio de 2020.
Grupo Kraken IDE. 2019. Global Precipitation Data (CHIRPS). Recuperado de: http://ide.unex.es/conocimiento/index.php?/article/AA-31745/0/Datos-globales-de-precipitacin-CHIRPS.htm; Consulta: Junio de 2020
Ministerio de Agricultura. 2014. Principales Cultivos por Área sembrada en el año 2014. Recuperado de: http://www.agronet.gov.co/Documents/Caquet%C3%A1.pdf; consulta Junio de 2020.
Viloria, J. s.f.GeoJson.io, herramienta para generar geodatos.Recuperado de: https://geoinnova.org/blog-territorio/geojson-io-herramienta-para-generar-geodatos/; consulta: Junio de 2020.