Geoestadística
Brayan Gómez
2025-12-06
Predicción espacial de variables climáticas mediante Kriging en Cordoba, Colombia
El presente ejercicio pretende realizar una predicción espacial mediante el método geoestadístico Kriging de las variables de temperatura media diaria, húmedad relativa media diaria y precipitación acumulada diaria para el periodo 2011-2020.
En primer lugar, se procede a importar las librerias necesarias para el procedimiento.
## Warning: package 'dplyr' was built under R version 4.5.2## Warning: package 'readr' was built under R version 4.5.2## Warning: package 'lubridate' was built under R version 4.5.2## Warning: package 'geoR' was built under R version 4.5.2## Warning: package 'akima' was built under R version 4.5.2## Warning: package 'scatterplot3d' was built under R version 4.5.2## Warning: package 'tidyverse' was built under R version 4.5.2## Warning: package 'ggplot2' was built under R version 4.5.2## Warning: package 'tidyr' was built under R version 4.5.2## Warning: package 'purrr' was built under R version 4.5.2## Warning: package 'forcats' was built under R version 4.5.2## Warning: package 'DT' was built under R version 4.5.2## Warning: package 'sf' was built under R version 4.5.2## Warning: package 'fields' was built under R version 4.5.2## Warning: package 'spam' was built under R version 4.5.2## Warning: package 'viridisLite' was built under R version 4.5.2## Warning: package 'RColorBrewer' was built under R version 4.5.2## Warning: package 'metR' was built under R version 4.5.2## Warning: package 'terra' was built under R version 4.5.2En archivos individuales de formato ‘Variable’@‘Número_de_estación’.data se encuentran las observaciones disponibles por estación y variable en la región. De este modo, se procede a leer los archivos por variable y agregarlos a un data.frame en el que se organiza en la primera columna la fecha y en las demás los respectivos valores registrados por estación. Además, se crea un archivo .csv en donde es compilan los datos para que una vez creado no se tenga que volver a correr este bloque de código que puede ser demorado.
ruta = "C:/Users/gomez/Desktop/Estadistica/datos_est/"
for (variable in c("HR_CAL_MEDIA_D", "TSSM_MEDIA_D", "PTPM_CON")) {
archivos = list.files(ruta, pattern = paste0("^",variable,".*\\.data$"), full.names = TRUE)
datos = map_df(archivos, function(archivo) {
est = str_extract(basename(archivo), "(?<=@)\\d+")
df = read.table(archivo, header = TRUE, sep = "|", stringsAsFactors = FALSE)
df$Fecha = as.Date(df$Fecha)
df = df %>% rename(!!est := Valor)
return(df)
})
if (variable=="HR_CAL_MEDIA_D") {
df_hr = datos %>%
group_by(Fecha) %>%
summarise(across(everything(), ~ first(na.omit(.x)), .names = "{.col}")) %>%
ungroup() %>%
arrange(Fecha)
}
if (variable=="TSSM_MEDIA_D") {
df_temp = datos %>%
group_by(Fecha) %>%
summarise(across(everything(), ~ first(na.omit(.x)), .names = "{.col}")) %>%
ungroup() %>%
arrange(Fecha)
}
if (variable=="PTPM_CON") {
df_prec = datos %>%
group_by(Fecha) %>%
summarise(across(everything(), ~ first(na.omit(.x)), .names = "{.col}")) %>%
ungroup() %>%
arrange(Fecha)
}
}
write.csv(df_hr,file=paste0(ruta,"datos_hr.csv"), row.names=F)
write.csv(df_temp,file=paste0(ruta,"datos_temp.csv"), row.names=F)
write.csv(df_prec,file=paste0(ruta,"datos_prec.csv"), row.names=F)
rm(variable,archivos)Una vez creado el .csv se hace la lectura de los datos y se filtra para los registros entre 2011 y 2020 que posean al menos el 80% de los datos para el periodo en cuestión. A partir de allí se crea un data.frame por cada variable.
Por otro lado, mediante la base de datos CNE del IDEAM se agrega las coordenadas geográficas a cada estación y posteriormente se convierten a coordenadas planas.
ruta="C:/Users/gomez/Desktop/Estadistica/datos_est/"
df_prec <- read.csv(paste0(ruta,"datos_prec.csv"), check.names = F)
df_temp = read.csv(paste0(ruta,"datos_temp.csv"), check.names = F)
df_hr = read.csv(paste0(ruta,"datos_hr.csv"), check.names = F)
#--------------Se escogen las estaciones con más del 80% de datos disponibles entre 2011 y 2020--------
fecha_inicio = as.Date("2011-01-01")
fecha_fin = as.Date("2020-12-31")
##--------Precipitacion-------------
df_prec$Fecha = as.Date(df_prec$Fecha)
df_prec_1120 = df_prec[, c("Fecha", names(colMeans(is.na(df_prec[df_prec$Fecha >= fecha_inicio & df_prec$Fecha <= fecha_fin, ][, -1])) * 100)[colMeans(is.na(df_prec[df_prec$Fecha >= fecha_inicio & df_prec$Fecha <= fecha_fin, ][, -1])) * 100 <= 20])]
##--------Temperatura-------------
df_temp$Fecha = as.Date(df_temp$Fecha)
df_temp_1120 = df_temp[, c("Fecha", names(colMeans(is.na(df_temp[df_temp$Fecha >= fecha_inicio & df_temp$Fecha <= fecha_fin, ][, -1])) * 100)[colMeans(is.na(df_temp[df_temp$Fecha >= fecha_inicio & df_temp$Fecha <= fecha_fin, ][, -1])) * 100 <= 100])]
##--------Humedad-------------
df_hr$Fecha = as.Date(df_hr$Fecha)
df_hr_1120 = df_hr[, c("Fecha", names(colMeans(is.na(df_hr[df_hr$Fecha >= fecha_inicio & df_hr$Fecha <= fecha_fin, ][, -1])) * 100)[colMeans(is.na(df_hr[df_hr$Fecha >= fecha_inicio & df_hr$Fecha <= fecha_fin, ][, -1])) * 100 <= 100])]
#------------------------Metadata Y TRANSFORMACIÓN DE COORDENADAS-------------------------
cne = read_excel(paste0(ruta,"/Catálogo_Nacional_de_Estaciones_del_IDEAM_20251125.xlsx"))
cne_sf = st_as_sf(cne,
coords = c("LONGITUD", "LATITUD"),
crs = 4326)
cne_plano = st_transform(cne_sf, crs = 3116)
cne_plano = cne_plano %>%
mutate(
COORD_X = st_coordinates(.)[,1],
COORD_Y = st_coordinates(.)[,2]
)
cne_final = cne_plano %>%
st_drop_geometry() %>%
dplyr::select(Codigo, Altitud, COORD_X, COORD_Y)
#-----------------Promedios y union con metadata------------
df_precprom = data.frame(estacion = names(colMeans(df_prec_1120[,-1], na.rm = TRUE)),promedio = as.numeric(colMeans(df_prec_1120[,-1], na.rm = TRUE)))
df_precprom = df_precprom %>% left_join(cne_final, by = c("estacion"="Codigo"))
df_tempprom = data.frame(estacion = names(colMeans(df_temp_1120[,-1], na.rm = TRUE)),promedio = as.numeric(colMeans(df_temp_1120[,-1], na.rm = TRUE)))
df_tempprom = df_tempprom %>% left_join(cne_final, by = c("estacion"="Codigo"))
df_hrprom = data.frame(estacion = names(colMeans(df_hr_1120[,-1], na.rm = TRUE)),promedio = as.numeric(colMeans(df_hr_1120[,-1], na.rm = TRUE)))
df_hrprom = df_hrprom %>% left_join(cne_final, by = c("estacion"="Codigo"))
rm(df_hr_1120,df_prec_1120,df_temp_1120,cne,cne_final,cne_plano,cne_sf,fecha_fin,fecha_inicio,df_hr,df_temp,df_prec)Ahora se procede con la visualización de la ubicación de las estaciones y el valor medio diario de cada variable en la región.
departamentos = st_read(paste0(ruta,"Deps3116.shp"))## Reading layer `Deps3116' from data source
## `C:\Users\gomez\Desktop\Estadistica\datos_est\Deps3116.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 33 features and 7 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 165290.2 ymin: 23255.83 xmax: 1805760 ymax: 1984871
## Projected CRS: MAGNA-SIRGAS / Colombia Bogota zonecordoba = st_read(paste0(ruta,"Cordoba.shp"))## Reading layer `Cordoba' from data source
## `C:\Users\gomez\Desktop\Estadistica\datos_est\Cordoba.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1 feature and 7 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: 730992.2 ymin: 1304596 xmax: 922502.3 ymax: 1537076
## Projected CRS: MAGNA-SIRGAS / Colombia Bogota zonecontornos = st_read(paste0(ruta,"Contornos.shp"))## Reading layer `Contornos' from data source
## `C:\Users\gomez\Desktop\Estadistica\datos_est\Contornos.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 77 features and 3 fields
## Geometry type: LINESTRING
## Dimension: XY
## Bounding box: xmin: 4553270 ymin: 2289029 xmax: 4862458 ymax: 2665367
## Projected CRS: MAGNA_CTM12precg=as.geodata(df_precprom,coords.col = c(4,5),data.col=2)
tempg=as.geodata(df_tempprom,coords.col = c(4,5),data.col=2)
hrg=as.geodata(df_hrprom,coords.col = c(4,5),data.col=2)
prec_sf = st_as_sf(df_precprom, coords = c("COORD_X", "COORD_Y"), crs = st_crs(departamentos))
temp_sf = st_as_sf(df_tempprom, coords = c("COORD_X", "COORD_Y"), crs = st_crs(departamentos))
hr_sf = st_as_sf(df_hrprom, coords = c("COORD_X", "COORD_Y"), crs = st_crs(departamentos))
ggplot() +
geom_sf(data = departamentos, fill = "gray90", color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
geom_sf(data = prec_sf, aes(color = promedio), size = 3) +
scale_color_viridis_c(option = "C") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Precipitación promedio por estación",
color = "mm/día") +
theme(plot.title = element_text(hjust = 0.5))
ggplot() +
geom_sf(data = departamentos, fill = "gray90", color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
geom_sf(data = temp_sf, aes(color = promedio), size = 3) +
scale_color_viridis_c(option = "C") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Temperatura promedio por estación",
color = "°C") +
theme(plot.title = element_text(hjust = 0.5))
ggplot() +
geom_sf(data = departamentos, fill = "gray90", color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
geom_sf(data = hr_sf, aes(color = promedio), size = 3) +
scale_color_viridis_c(option = "C") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Humedad promedio por estación",
color = "%") +
theme(plot.title = element_text(hjust = 0.5))
Como se observa en los mapas, existen algunas estaciones que se encuentran sobre grandes elevaciones o en la vertiente contraria de la serranía. Es por ello que se procede a descartar dichas estaciones para el procedimiento de Kriging.
df_precprom = df_precprom[df_precprom$COORD_X>730000 & df_precprom$COORD_Y>1350000,]
df_tempprom = df_tempprom[df_tempprom$COORD_X>730000,]
df_hrprom = df_hrprom[df_hrprom$COORD_X>730000,]
departamentos = st_read(paste0(ruta,"Deps3116.shp"))## Reading layer `Deps3116' from data source
## `C:\Users\gomez\Desktop\Estadistica\datos_est\Deps3116.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 33 features and 7 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 165290.2 ymin: 23255.83 xmax: 1805760 ymax: 1984871
## Projected CRS: MAGNA-SIRGAS / Colombia Bogota zonecordoba = st_read(paste0(ruta,"Cordoba.shp"))## Reading layer `Cordoba' from data source
## `C:\Users\gomez\Desktop\Estadistica\datos_est\Cordoba.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1 feature and 7 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: 730992.2 ymin: 1304596 xmax: 922502.3 ymax: 1537076
## Projected CRS: MAGNA-SIRGAS / Colombia Bogota zonecontornos = st_read(paste0(ruta,"Contornos.shp"))## Reading layer `Contornos' from data source
## `C:\Users\gomez\Desktop\Estadistica\datos_est\Contornos.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 77 features and 3 fields
## Geometry type: LINESTRING
## Dimension: XY
## Bounding box: xmin: 4553270 ymin: 2289029 xmax: 4862458 ymax: 2665367
## Projected CRS: MAGNA_CTM12precg=as.geodata(df_precprom,coords.col = c(4,5),data.col=2)
tempg=as.geodata(df_tempprom,coords.col = c(4,5),data.col=2)
hrg=as.geodata(df_hrprom,coords.col = c(4,5),data.col=2)
prec_sf = st_as_sf(df_precprom, coords = c("COORD_X", "COORD_Y"), crs = st_crs(departamentos))
temp_sf = st_as_sf(df_tempprom, coords = c("COORD_X", "COORD_Y"), crs = st_crs(departamentos))
hr_sf = st_as_sf(df_hrprom, coords = c("COORD_X", "COORD_Y"), crs = st_crs(departamentos))
ggplot() +
geom_sf(data = departamentos, fill = "gray90", color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
geom_sf(data = prec_sf, aes(color = promedio), size = 3) +
scale_color_viridis_c(option = "C") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Precipitación promedio por estación",
color = "mm/día") +
theme(plot.title = element_text(hjust = 0.5))
ggplot() +
geom_sf(data = departamentos, fill = "gray90", color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
geom_sf(data = temp_sf, aes(color = promedio), size = 3) +
scale_color_viridis_c(option = "C") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Temperatura promedio por estación",
color = "°C") +
theme(plot.title = element_text(hjust = 0.5))
ggplot() +
geom_sf(data = departamentos, fill = "gray90", color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
geom_sf(data = hr_sf, aes(color = promedio), size = 3) +
scale_color_viridis_c(option = "C") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Humedad promedio por estación",
color = "%") +
theme(plot.title = element_text(hjust = 0.5))
Teniendo definidas las estaciones a utilizar, se elabora una matriz de diagramas de dispersión para identificar posibles tendencias entre los valores y las coordenadas
plot(df_precprom %>% dplyr::select(promedio, Altitud, COORD_Y, COORD_X), main="Matriz Precipitacion")
plot(df_tempprom %>% dplyr::select(promedio, Altitud, COORD_Y, COORD_X), main="Matriz Temperatura")
plot(df_hrprom %>% dplyr::select(promedio, Altitud, COORD_Y, COORD_X), main="Matriz Humedad")
Una vez visualizada las posibles tendencias se construyen los modelos de regresión lineal apropiados que mejor representen el comportamiento espacial de los valores. No obstante, como se utilizará kriging ordinario no se utilizarán los residuales, interceptos o coeficientes de los modelos.
##
## Call:
## lm(formula = promedio ~ COORD_X + COORD_Y, data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1160 -0.8393 -0.3056 0.4610 6.4729
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.534e+01 4.524e+00 5.601 2.19e-07 ***
## COORD_X 2.457e-05 2.863e-06 8.584 2.17e-13 ***
## COORD_Y -2.808e-05 2.675e-06 -10.493 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.417 on 92 degrees of freedom
## Multiple R-squared: 0.6581, Adjusted R-squared: 0.6506
## F-statistic: 88.53 on 2 and 92 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ COORD_X * COORD_Y, data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.0354 -0.6512 -0.2486 0.4810 5.3487
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.706e+02 6.298e+01 -5.884 6.59e-08 ***
## COORD_X 4.860e-04 7.331e-05 6.630 2.34e-09 ***
## COORD_Y 2.447e-04 4.337e-05 5.642 1.88e-07 ***
## COORD_X:COORD_Y -3.179e-10 5.048e-11 -6.298 1.05e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.189 on 91 degrees of freedom
## Multiple R-squared: 0.7619, Adjusted R-squared: 0.754
## F-statistic: 97.05 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) + COORD_Y, data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.0680 -0.8101 -0.2842 0.3940 6.3908
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.537e+01 4.012e+00 8.816 7.07e-14 ***
## I(COORD_X^2) 1.473e-11 1.671e-12 8.816 7.08e-14 ***
## COORD_Y -2.795e-05 2.643e-06 -10.576 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.4 on 92 degrees of freedom
## Multiple R-squared: 0.6662, Adjusted R-squared: 0.6589
## F-statistic: 91.8 on 2 and 92 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ COORD_X + I(COORD_Y^2), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1039 -0.8052 -0.3100 0.4484 6.4683
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.776e+00 3.073e+00 1.554 0.124
## COORD_X 2.469e-05 2.875e-06 8.591 2.1e-13 ***
## I(COORD_Y^2) -9.615e-12 9.231e-13 -10.417 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.423 on 92 degrees of freedom
## Multiple R-squared: 0.6553, Adjusted R-squared: 0.6478
## F-statistic: 87.46 on 2 and 92 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) * I(COORD_Y^2), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.9478 -0.6545 -0.1866 0.4406 5.1770
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -8.848e+01 1.561e+01 -5.669 1.67e-07 ***
## I(COORD_X^2) 1.550e-10 2.104e-11 7.363 7.79e-11 ***
## I(COORD_Y^2) 3.945e-11 7.382e-12 5.343 6.70e-07 ***
## I(COORD_X^2):I(COORD_Y^2) -6.643e-23 9.953e-24 -6.675 1.90e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.158 on 91 degrees of freedom
## Multiple R-squared: 0.774, Adjusted R-squared: 0.7666
## F-statistic: 103.9 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * COORD_Y, data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8992 -0.6211 -0.1438 0.4334 5.0353
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.035e+02 2.029e+01 -5.101 1.83e-06 ***
## I(COORD_X^3) 2.349e-16 3.150e-17 7.457 5.02e-11 ***
## COORD_Y 7.002e-05 1.396e-05 5.015 2.61e-06 ***
## I(COORD_X^3):COORD_Y -1.537e-22 2.168e-23 -7.089 2.81e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.117 on 91 degrees of freedom
## Multiple R-squared: 0.7898, Adjusted R-squared: 0.7829
## F-statistic: 114 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ COORD_X * I(COORD_Y^3), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.9896 -0.6777 -0.2454 0.4588 5.3273
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.354e+02 2.140e+01 -6.328 9.18e-09 ***
## COORD_X 1.802e-04 2.492e-05 7.230 1.46e-10 ***
## I(COORD_Y^3) 3.907e-17 6.947e-18 5.624 2.03e-07 ***
## COORD_X:I(COORD_Y^3) -5.066e-23 8.088e-24 -6.264 1.22e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.201 on 91 degrees of freedom
## Multiple R-squared: 0.7571, Adjusted R-squared: 0.7491
## F-statistic: 94.55 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * I(COORD_Y^2), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8776 -0.6349 -0.1382 0.4306 5.0215
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -5.317e+01 1.024e+01 -5.193 1.25e-06 ***
## I(COORD_X^3) 1.240e-16 1.590e-17 7.798 9.99e-12 ***
## I(COORD_Y^2) 2.432e-11 4.839e-12 5.025 2.50e-06 ***
## I(COORD_X^3):I(COORD_Y^2) -5.318e-29 7.518e-30 -7.074 3.01e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.123 on 91 degrees of freedom
## Multiple R-squared: 0.7877, Adjusted R-squared: 0.7807
## F-statistic: 112.5 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) * I(COORD_Y^3), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.9247 -0.6769 -0.1824 0.4404 5.1667
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.122e+01 1.052e+01 -5.822 8.65e-08 ***
## I(COORD_X^2) 1.088e-10 1.419e-11 7.672 1.82e-11 ***
## I(COORD_Y^3) 1.822e-17 3.410e-18 5.342 6.74e-07 ***
## I(COORD_X^2):I(COORD_Y^3) -3.061e-29 4.599e-30 -6.656 2.07e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.164 on 91 degrees of freedom
## Multiple R-squared: 0.7716, Adjusted R-squared: 0.7641
## F-statistic: 102.5 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * I(COORD_Y^3), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8548 -0.6489 -0.1318 0.4274 5.0102
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.642e+01 6.899e+00 -5.279 8.78e-07 ***
## I(COORD_X^3) 8.712e-17 1.073e-17 8.122 2.14e-12 ***
## I(COORD_Y^3) 1.125e-17 2.235e-18 5.032 2.43e-06 ***
## I(COORD_X^3):I(COORD_Y^3) -2.451e-35 3.475e-36 -7.056 3.28e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.129 on 91 degrees of freedom
## Multiple R-squared: 0.7854, Adjusted R-squared: 0.7783
## F-statistic: 111 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(1/COORD_X) * COORD_Y, data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1493 -0.7101 -0.3013 0.4729 5.6741
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.287e+02 6.548e+01 6.548 3.40e-09 ***
## I(1/COORD_X) -3.269e+08 5.595e+07 -5.842 7.90e-08 ***
## COORD_Y -2.781e-04 4.512e-05 -6.164 1.91e-08 ***
## I(1/COORD_X):COORD_Y 2.137e+02 3.855e+01 5.544 2.86e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.263 on 91 degrees of freedom
## Multiple R-squared: 0.7312, Adjusted R-squared: 0.7223
## F-statistic: 82.5 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ COORD_X * I(1/COORD_Y), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.0759 -0.6479 -0.2531 0.4540 5.3803
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.358e+02 6.219e+01 5.399 5.29e-07 ***
## COORD_X -4.326e-04 7.237e-05 -5.977 4.37e-08 ***
## I(1/COORD_Y) -5.089e+08 9.011e+07 -5.647 1.84e-07 ***
## COORD_X:I(1/COORD_Y) 6.624e+02 1.049e+02 6.317 9.65e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.18 on 91 degrees of freedom
## Multiple R-squared: 0.7657, Adjusted R-squared: 0.7579
## F-statistic: 99.1 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ I(1/COORD_X) * I(1/COORD_Y), data = df_precprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1913 -0.6947 -0.2867 0.4587 5.7015
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.758e+02 6.469e+01 -5.810 9.11e-08 ***
## I(1/COORD_X) 2.910e+08 5.529e+07 5.264 9.34e-07 ***
## I(1/COORD_Y) 5.809e+08 9.368e+07 6.201 1.62e-08 ***
## I(1/COORD_X):I(1/COORD_Y) -4.458e+14 8.007e+13 -5.568 2.59e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.254 on 91 degrees of freedom
## Multiple R-squared: 0.7353, Adjusted R-squared: 0.7266
## F-statistic: 84.27 on 3 and 91 DF, p-value: < 2.2e-16##
## Call:
## lm(formula = promedio ~ COORD_X + COORD_Y, data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.44697 -0.23080 -0.01564 0.19158 0.48245
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.423e+01 3.232e+00 7.497 7.26e-06 ***
## COORD_X 3.201e-06 2.449e-06 1.307 0.216
## COORD_Y 5.991e-07 1.570e-06 0.382 0.709
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3074 on 12 degrees of freedom
## Multiple R-squared: 0.1277, Adjusted R-squared: -0.01773
## F-statistic: 0.878 on 2 and 12 DF, p-value: 0.4407##
## Call:
## lm(formula = promedio ~ COORD_X * COORD_Y, data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43483 -0.25464 0.07437 0.19639 0.47751
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -5.439e+00 5.182e+01 -0.105 0.918
## COORD_X 3.942e-05 6.318e-05 0.624 0.545
## COORD_Y 2.110e-05 3.577e-05 0.590 0.567
## COORD_X:COORD_Y -2.503e-11 4.363e-11 -0.574 0.578
##
## Residual standard error: 0.3164 on 11 degrees of freedom
## Multiple R-squared: 0.153, Adjusted R-squared: -0.078
## F-statistic: 0.6623 on 3 and 11 DF, p-value: 0.5923##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) + COORD_Y, data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4472 -0.2289 -0.0155 0.1931 0.4839
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.550e+01 2.613e+00 9.758 4.66e-07 ***
## I(COORD_X^2) 1.984e-12 1.490e-12 1.332 0.208
## COORD_Y 6.113e-07 1.567e-06 0.390 0.703
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3067 on 12 degrees of freedom
## Multiple R-squared: 0.1317, Adjusted R-squared: -0.01296
## F-statistic: 0.9104 on 2 and 12 DF, p-value: 0.4284##
## Call:
## lm(formula = promedio ~ COORD_X + I(COORD_Y^2), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.44639 -0.23046 -0.01714 0.19153 0.48240
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.467e+01 2.441e+00 10.109 3.19e-07 ***
## COORD_X 3.200e-06 2.450e-06 1.306 0.216
## I(COORD_Y^2) 2.025e-13 5.404e-13 0.375 0.714
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3075 on 12 degrees of freedom
## Multiple R-squared: 0.1273, Adjusted R-squared: -0.01817
## F-statistic: 0.8751 on 2 and 12 DF, p-value: 0.4418##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) * I(COORD_Y^2), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43500 -0.25319 0.07736 0.19823 0.47936
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.844e+01 1.302e+01 1.416 0.185
## I(COORD_X^2) 1.315e-11 1.926e-11 0.683 0.509
## I(COORD_Y^2) 3.785e-12 6.180e-12 0.612 0.553
## I(COORD_X^2):I(COORD_Y^2) -5.323e-24 9.155e-24 -0.581 0.573
##
## Residual standard error: 0.3156 on 11 degrees of freedom
## Multiple R-squared: 0.1573, Adjusted R-squared: -0.07259
## F-statistic: 0.6842 on 3 and 11 DF, p-value: 0.5801##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * COORD_Y, data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43518 -0.25152 0.07747 0.20029 0.48185
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.603e+01 1.743e+01 0.920 0.377
## I(COORD_X^3) 1.950e-17 3.121e-17 0.625 0.545
## COORD_Y 7.446e-06 1.202e-05 0.619 0.548
## I(COORD_X^3):COORD_Y -1.234e-23 2.155e-23 -0.573 0.578
##
## Residual standard error: 0.3149 on 11 degrees of freedom
## Multiple R-squared: 0.1608, Adjusted R-squared: -0.06809
## F-statistic: 0.7025 on 3 and 11 DF, p-value: 0.5701##
## Call:
## lm(formula = promedio ~ COORD_X * I(COORD_Y^3), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43493 -0.25492 0.07605 0.19627 0.47697
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.467e+01 1.737e+01 0.844 0.416
## COORD_X 1.558e-05 2.116e-05 0.736 0.477
## I(COORD_Y^3) 3.425e-18 5.661e-18 0.605 0.557
## COORD_X:I(COORD_Y^3) -4.069e-24 6.902e-24 -0.589 0.567
##
## Residual standard error: 0.3163 on 11 degrees of freedom
## Multiple R-squared: 0.1537, Adjusted R-squared: -0.07716
## F-statistic: 0.6657 on 3 and 11 DF, p-value: 0.5904##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * I(COORD_Y^2), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43506 -0.25166 0.07666 0.20025 0.48162
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.137e+01 8.738e+00 2.445 0.0325 *
## I(COORD_X^3) 1.068e-17 1.562e-17 0.684 0.5082
## I(COORD_Y^2) 2.595e-12 4.141e-12 0.627 0.5437
## I(COORD_X^3):I(COORD_Y^2) -4.312e-30 7.421e-30 -0.581 0.5730
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3149 on 11 degrees of freedom
## Multiple R-squared: 0.1611, Adjusted R-squared: -0.06769
## F-statistic: 0.7041 on 3 and 11 DF, p-value: 0.5692##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) * I(COORD_Y^3), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43545 -0.25331 0.07749 0.19817 0.47910
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.103e+01 8.722e+00 2.411 0.0345 *
## I(COORD_X^2) 9.516e-12 1.288e-11 0.739 0.4756
## I(COORD_Y^3) 1.757e-18 2.837e-18 0.619 0.5483
## I(COORD_X^2):I(COORD_Y^3) -2.475e-30 4.201e-30 -0.589 0.5678
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3155 on 11 degrees of freedom
## Multiple R-squared: 0.1575, Adjusted R-squared: -0.07221
## F-statistic: 0.6857 on 3 and 11 DF, p-value: 0.5792##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * I(COORD_Y^3), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43576 -0.25178 0.07573 0.20019 0.48138
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.314e+01 5.853e+00 3.954 0.00226 **
## I(COORD_X^3) 7.736e-18 1.044e-17 0.741 0.47426
## I(COORD_Y^3) 1.204e-18 1.901e-18 0.633 0.53954
## I(COORD_X^3):I(COORD_Y^3) -2.004e-36 3.405e-36 -0.589 0.56802
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3148 on 11 degrees of freedom
## Multiple R-squared: 0.1614, Adjusted R-squared: -0.06733
## F-statistic: 0.7056 on 3 and 11 DF, p-value: 0.5684##
## Call:
## lm(formula = promedio ~ I(1/COORD_X) * COORD_Y, data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43378 -0.25804 0.07012 0.19302 0.47376
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.922e+01 5.188e+01 1.141 0.278
## I(1/COORD_X) -2.643e+07 4.237e+07 -0.624 0.546
## COORD_Y -2.003e-05 3.585e-05 -0.559 0.587
## I(1/COORD_X):COORD_Y 1.684e+01 2.927e+01 0.575 0.577
##
## Residual standard error: 0.3179 on 11 degrees of freedom
## Multiple R-squared: 0.1451, Adjusted R-squared: -0.08811
## F-statistic: 0.6221 on 3 and 11 DF, p-value: 0.6153##
## Call:
## lm(formula = promedio ~ COORD_X * I(1/COORD_Y), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43514 -0.25430 0.07268 0.19646 0.47800
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.484e+01 5.202e+01 1.054 0.314
## COORD_X -3.203e-05 6.350e-05 -0.505 0.624
## I(1/COORD_Y) -4.298e+07 7.511e+07 -0.572 0.579
## COORD_X:I(1/COORD_Y) 5.090e+01 9.165e+01 0.555 0.590
##
## Residual standard error: 0.3165 on 11 degrees of freedom
## Multiple R-squared: 0.1522, Adjusted R-squared: -0.079
## F-statistic: 0.6583 on 3 and 11 DF, p-value: 0.5945##
## Call:
## lm(formula = promedio ~ I(1/COORD_X) * I(1/COORD_Y), data = df_tempprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.43399 -0.25769 0.06852 0.19313 0.47434
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.104e+00 5.215e+01 0.040 0.969
## I(1/COORD_X) 2.163e+07 4.255e+07 0.508 0.621
## I(1/COORD_Y) 4.064e+07 7.523e+07 0.540 0.600
## I(1/COORD_X):I(1/COORD_Y) -3.422e+13 6.140e+13 -0.557 0.588
##
## Residual standard error: 0.318 on 11 degrees of freedom
## Multiple R-squared: 0.1442, Adjusted R-squared: -0.08917
## F-statistic: 0.6179 on 3 and 11 DF, p-value: 0.6177##
## Call:
## lm(formula = promedio ~ COORD_X + COORD_Y, data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2135 -0.9350 -0.0266 2.0452 3.3961
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.073e+01 2.776e+01 2.908 0.0131 *
## COORD_X -9.830e-06 2.103e-05 -0.467 0.6486
## COORD_Y 6.579e-06 1.349e-05 0.488 0.6345
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.64 on 12 degrees of freedom
## Multiple R-squared: 0.04175, Adjusted R-squared: -0.118
## F-statistic: 0.2614 on 2 and 12 DF, p-value: 0.7742##
## Call:
## lm(formula = promedio ~ COORD_X * COORD_Y, data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0179 -1.0171 -0.0841 1.4930 3.2557
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -7.352e+02 3.785e+02 -1.943 0.0781 .
## COORD_X 9.863e-04 4.614e-04 2.137 0.0559 .
## COORD_Y 5.704e-04 2.613e-04 2.183 0.0516 .
## COORD_X:COORD_Y -6.885e-10 3.187e-10 -2.160 0.0537 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.311 on 11 degrees of freedom
## Multiple R-squared: 0.3272, Adjusted R-squared: 0.1437
## F-statistic: 1.783 on 3 and 11 DF, p-value: 0.2084##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) + COORD_Y, data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2088 -0.9303 -0.0300 2.0472 3.3970
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.670e+01 2.250e+01 3.409 0.00518 **
## I(COORD_X^2) -5.947e-12 1.283e-11 -0.464 0.65128
## COORD_Y 6.561e-06 1.349e-05 0.486 0.63557
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.641 on 12 degrees of freedom
## Multiple R-squared: 0.04146, Adjusted R-squared: -0.1183
## F-statistic: 0.2595 on 2 and 12 DF, p-value: 0.7756##
## Call:
## lm(formula = promedio ~ COORD_X + I(COORD_Y^2), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2299 -0.9297 -0.0174 2.0447 3.3802
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.537e+01 2.095e+01 4.075 0.00154 **
## COORD_X -9.795e-06 2.103e-05 -0.466 0.64966
## I(COORD_Y^2) 2.310e-12 4.638e-12 0.498 0.62741
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.639 on 12 degrees of freedom
## Multiple R-squared: 0.04254, Adjusted R-squared: -0.117
## F-statistic: 0.2666 on 2 and 12 DF, p-value: 0.7704##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) * I(COORD_Y^2), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0057 -0.9776 -0.0712 1.5069 3.2394
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.270e+02 9.453e+01 -1.343 0.2062
## I(COORD_X^2) 3.037e-10 1.398e-10 2.173 0.0525 .
## I(COORD_Y^2) 1.016e-10 4.485e-11 2.265 0.0447 *
## I(COORD_X^2):I(COORD_Y^2) -1.476e-22 6.644e-23 -2.222 0.0482 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.29 on 11 degrees of freedom
## Multiple R-squared: 0.339, Adjusted R-squared: 0.1587
## F-statistic: 1.88 on 3 and 11 DF, p-value: 0.1914##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * COORD_Y, data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.9998 -0.9547 -0.0941 1.4952 3.2351
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.057e+02 1.264e+02 -1.627 0.1320
## I(COORD_X^3) 5.034e-16 2.263e-16 2.224 0.0480 *
## COORD_Y 2.006e-04 8.717e-05 2.302 0.0419 *
## I(COORD_X^3):COORD_Y -3.512e-22 1.563e-22 -2.247 0.0461 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.284 on 11 degrees of freedom
## Multiple R-squared: 0.3428, Adjusted R-squared: 0.1635
## F-statistic: 1.912 on 3 and 11 DF, p-value: 0.1861##
## Call:
## lm(formula = promedio ~ COORD_X * I(COORD_Y^3), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0116 -0.9996 -0.0483 1.5179 3.2428
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.879e+02 1.262e+02 -1.489 0.1646
## COORD_X 3.255e-04 1.537e-04 2.118 0.0578 .
## I(COORD_Y^3) 9.132e-17 4.111e-17 2.221 0.0483 *
## COORD_X:I(COORD_Y^3) -1.101e-22 5.013e-23 -2.197 0.0503 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.297 on 11 degrees of freedom
## Multiple R-squared: 0.3351, Adjusted R-squared: 0.1538
## F-statistic: 1.848 on 3 and 11 DF, p-value: 0.1969##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * I(COORD_Y^2), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.9959 -0.9454 -0.0754 1.5076 3.2286
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.116e+01 6.319e+01 -0.968 0.3539
## I(COORD_X^3) 2.502e-16 1.129e-16 2.216 0.0487 *
## I(COORD_Y^2) 6.953e-11 2.995e-11 2.322 0.0405 *
## I(COORD_X^3):I(COORD_Y^2) -1.216e-28 5.367e-29 -2.265 0.0447 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.277 on 11 degrees of freedom
## Multiple R-squared: 0.3467, Adjusted R-squared: 0.1685
## F-statistic: 1.946 on 3 and 11 DF, p-value: 0.1807##
## Call:
## lm(formula = promedio ~ I(COORD_X^2) * I(COORD_Y^3), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0022 -0.9682 -0.0524 1.5185 3.2331
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -5.676e+01 6.312e+01 -0.899 0.3878
## I(COORD_X^2) 2.015e-10 9.322e-11 2.161 0.0536 .
## I(COORD_Y^3) 4.690e-17 2.053e-17 2.284 0.0432 *
## I(COORD_X^2):I(COORD_Y^3) -6.813e-29 3.041e-29 -2.241 0.0466 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.284 on 11 degrees of freedom
## Multiple R-squared: 0.3429, Adjusted R-squared: 0.1638
## F-statistic: 1.914 on 3 and 11 DF, p-value: 0.1859##
## Call:
## lm(formula = promedio ~ I(COORD_X^3) * I(COORD_Y^3), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.9920 -0.9357 -0.0561 1.5189 3.2223
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.308e+01 4.221e+01 -0.310 0.7625
## I(COORD_X^3) 1.660e-16 7.529e-17 2.205 0.0497 *
## I(COORD_Y^3) 3.210e-17 1.371e-17 2.342 0.0390 *
## I(COORD_X^3):I(COORD_Y^3) -5.607e-35 2.455e-35 -2.284 0.0432 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.27 on 11 degrees of freedom
## Multiple R-squared: 0.3507, Adjusted R-squared: 0.1736
## F-statistic: 1.98 on 3 and 11 DF, p-value: 0.1755##
## Call:
## lm(formula = promedio ~ I(1/COORD_X) * COORD_Y, data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0330 -1.0740 -0.0721 1.4900 3.2717
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.542e+02 3.816e+02 2.238 0.0468 *
## I(1/COORD_X) -6.389e+08 3.117e+08 -2.050 0.0650 .
## COORD_Y -5.394e-04 2.637e-04 -2.046 0.0655 .
## I(1/COORD_X):COORD_Y 4.463e+02 2.153e+02 2.073 0.0625 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.338 on 11 degrees of freedom
## Multiple R-squared: 0.3112, Adjusted R-squared: 0.1234
## F-statistic: 1.657 on 3 and 11 DF, p-value: 0.2332##
## Call:
## lm(formula = promedio ~ COORD_X * I(1/COORD_Y), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0243 -1.0329 -0.1173 1.4644 3.2691
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.097e+02 3.819e+02 2.382 0.0364 *
## COORD_X -9.997e-04 4.662e-04 -2.144 0.0552 .
## I(1/COORD_Y) -1.184e+09 5.515e+08 -2.147 0.0550 .
## COORD_X:I(1/COORD_Y) 1.430e+03 6.729e+02 2.125 0.0571 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.324 on 11 degrees of freedom
## Multiple R-squared: 0.3195, Adjusted R-squared: 0.1339
## F-statistic: 1.721 on 3 and 11 DF, p-value: 0.2202##
## Call:
## lm(formula = promedio ~ I(1/COORD_X) * I(1/COORD_Y), data = df_hrprom)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.038 -1.088 -0.103 1.460 3.285
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -7.010e+02 3.855e+02 -1.818 0.0963 .
## I(1/COORD_X) 6.473e+08 3.146e+08 2.058 0.0641 .
## I(1/COORD_Y) 1.119e+09 5.561e+08 2.012 0.0694 .
## I(1/COORD_X):I(1/COORD_Y) -9.252e+14 4.539e+14 -2.038 0.0663 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.351 on 11 degrees of freedom
## Multiple R-squared: 0.3036, Adjusted R-squared: 0.1137
## F-statistic: 1.599 on 3 and 11 DF, p-value: 0.2458En este punto, se procede a convertir los datos de cada data.frame en formato geodata para la manipulación con la libreria geoR.
precg=as.geodata(df_precprom,coords.col=4:5,data.col=2)
tempg=as.geodata(df_tempprom,coords.col=4:5,data.col=2)
hrg=as.geodata(df_hrprom,coords.col=4:5,data.col=2)
plot(precg)
mtext("Precipitación", side=2, line=1, cex=1.2)
plot(tempg)
mtext("Temperatura", side=2, line=1, cex=1.2)
plot(hrg)
mtext("Humedad", side=2, line=1, cex=1.2)
Posteriormente, se procede a graficar los variogramas por variable y en las 4 direcciones posibles para observar la anisometría. En el caso de la precipitación y habiendo visualizado los diagramas de dispersión, se le asigna una tendencia de segundo orden y se utiliza el estimador Cressie-Modulus para disminuir el efecto de los datos atípicos.
Se puede observar que existe anisotropía tanto geométrica como zonal. Para lidiar con esta condición es necesario extraer el angulo con respecto al eje Y de mayor rango presentado y la razón de anisotropía que es el cociente entre el máximo rango y el mínimo (hallado a 90° del eje de mayor rango).
vario_0_prec = variog(precg, trend="2nd", estimator.type="modulus", direction = 0, tolerance = pi/8, pairs.min = 20, max.dist=180000)## variog: computing variogram for direction = 0 degrees (0 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_45_prec = variog(precg, trend="2nd", estimator.type="modulus", direction = pi/4, tolerance = pi/8, pairs.min = 20, max.dist=180000)## variog: computing variogram for direction = 45 degrees (0.785 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_90_prec = variog(precg, trend="2nd", estimator.type="modulus", direction = pi/2, tolerance = pi/8, pairs.min = 20, max.dist=180000)## variog: computing variogram for direction = 90 degrees (1.571 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_135_prec = variog(precg, trend="2nd", estimator.type="modulus", direction = 3*pi/4, tolerance = pi/8, pairs.min = 20, max.dist=180000)## variog: computing variogram for direction = 135 degrees (2.356 radians)
## tolerance angle = 22.5 degrees (0.393 radians)par(mfrow=c(2,2), mar=c(3,3,2,2))
plot(vario_0_prec, main="Direccion 0")
plot(vario_45_prec, main="Direccion 45")
plot(vario_90_prec, main="Direccion 90")
plot(vario_135_prec, main="Direccion 135")
angle_aniso_prec = 3*pi/4
ratio_aniso_prec = 46759/32728
vario_prec=variog(precg, trend="2nd", estimator.type="modulus", pairs.min = 20, max.dist=180000)## variog: computing omnidirectional variogram#---------Mejores ajustes----------
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 wave 1.42 46759.07 0.26 <NA> <NA> 139877.720796716
#2 cubic 0.79 275878.49 0.17 <NA> <NA> 275878.49
#3 gaussian 0.82 126249.48 0.17 <NA> <NA> 218514.92077798
fitvar_prec1 = variofit(vario_prec, cov.model="wave",ini.cov.pars = c(1.42,46759.07),nugget=0.17,wei="npairs")## variofit: covariance model used is wave
## variofit: weights used: npairs
## variofit: minimisation function used: optimfitvar_prec2 = variofit(vario_prec, cov.model="cubic",ini.cov.pars = c(0.79,275878.49),nugget=0.17,wei="npairs")## variofit: covariance model used is cubic
## variofit: weights used: npairs
## variofit: minimisation function used: optim## Warning in variofit(vario_prec, cov.model = "cubic", ini.cov.pars = c(0.79, :
## unreasonable initial value for phi (too high)fitvar_prec3 = variofit(vario_prec, cov.model="gaussian",ini.cov.pars = c(0.82,126249.48),nugget=0.17,wei="npairs")## variofit: covariance model used is gaussian
## variofit: weights used: npairs
## variofit: minimisation function used: optimpar(mfrow=c(1,1))
plot(vario_prec,
xlab = "h",
ylab = "semivarianza",
cex.lab = 1.3,
cex.axis = 1.2,
main = "Estimación teórica del modelo de semivariograma de precipitación",
col.main = 4, cex.main =1.3)
lines(fitvar_prec1, col = 1)
lines(fitvar_prec2, col = 2)
lines(fitvar_prec3, col = 3)
legend(1, 1,
c("Wave","Cubic","Gaussian"),
lwd = 2,
lty = 1:3,
col = 1:3,
box.col = 9,
text.col = 1:3)
Una vez definidos los posibles parámetros de las funciones de semivariograma y los coeficientes de anisometría se continua con el Kriging evaluando el mejor modelo posible mediante boxplot de su varianza. Se evita utilizar el modelo Wave dado que genera patrones oscilatorios poco coherentes.
grilla=st_sample(cordoba,size=c(1000,1000),type="regular")
coords = st_coordinates(grilla)[,1:2]
kc_prec_wave=krige.conv(precg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = fitvar_prec1$cov.pars,nugget=fitvar_prec1$nugget,cov.model = "wave", aniso.pars=c(angle_aniso_prec, ratio_aniso_prec)))## krige.conv: model with constant mean
## krige.conv: anisotropy correction performed
## krige.conv: Kriging performed using global neighbourhoodkc_prec_cubic=krige.conv(precg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = fitvar_prec2$cov.pars,nugget=fitvar_prec2$nugget,cov.model = "cubic", aniso.pars=c(angle_aniso_prec, ratio_aniso_prec)))## krige.conv: model with constant mean
## krige.conv: anisotropy correction performed
## krige.conv: Kriging performed using global neighbourhoodkc_prec_gaus=krige.conv(precg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = fitvar_prec3$cov.pars,nugget=fitvar_prec3$nugget,cov.model = "gaussian", aniso.pars=c(angle_aniso_prec, ratio_aniso_prec)))## krige.conv: model with constant mean
## krige.conv: anisotropy correction performed
## krige.conv: Kriging performed using global neighbourhoodboxplot(kc_prec_wave$krige.var,
kc_prec_cubic$krige.var,
kc_prec_gaus$krige.var,
names=c("Wave","Cubic","Gaussian"),
main="Comparación de modelos",
ylab="Varianza de kriging")
datos_grafica = data.frame(coords, pred = kc_prec_gaus$predict, var=kc_prec_gaus$krige.var)
grilla_sf = st_as_sf(
datos_grafica,
coords = c("X", "Y"),
crs = st_crs(cordoba)
)
r_prec = rast(datos_grafica[, c("X", "Y", "pred", "var")],type = "xyz",crs = st_crs(cordoba)$wkt)
r_prec = crop(r_prec, vect(cordoba))
r_prec = mask(r_prec, vect(cordoba))
par(mfrow = c(2, 1),
mar = c(1,1,1,1))
r_df = as.data.frame(r_prec, xy = TRUE)
ggplot() +
geom_raster(data = r_df,
aes(x = x, y = y, fill = pred)) +
geom_sf(data = departamentos, fill = NA, color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
scale_fill_viridis_c(name = "Precipitación") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Precipitación en Cordoba mediante kriging")
ggplot() +
geom_raster(data = r_df,
aes(x = x, y = y, fill = var)) +
geom_sf(data = departamentos, fill = NA, color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
scale_fill_viridis_c(name = "Varianza") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Varianza del kriging")
Se realizan los procedimientos analogos para temperatura teniendo en cuenta que no mostró tendencia.
vario_0_temp = variog(tempg, trend="cte", estimator.type="modulus", direction = 0, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 0 degrees (0 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_45_temp = variog(tempg, trend="cte", estimator.type="modulus", direction = pi/4, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 45 degrees (0.785 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_90_temp = variog(tempg, trend="cte", estimator.type="modulus", direction = pi/2, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 90 degrees (1.571 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_135_temp = variog(tempg, trend="cte", estimator.type="modulus", direction = 3*pi/4, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 135 degrees (2.356 radians)
## tolerance angle = 22.5 degrees (0.393 radians)par(mfrow=c(2,2), mar=c(3,3,2,2))
plot(vario_0_temp, main="Direccion 0")
plot(vario_45_temp, main="Direccion 45")
plot(vario_90_temp, main="Direccion 90")
plot(vario_135_temp, main="Direccion 135")
vario_temp=variog(tempg, trend="cte", estimator.type="modulus", pairs.min = 3, max.dist=130000)## variog: computing omnidirectional variogram#---------Mejores ajustes----------
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 exponential 0.04 36921.62 0.05 <NA> <NA> 110607.288623314
#2 gaussian 0.04 67130.23 0.05 <NA> <NA> 116190.236110056
#3 cubic 0.04 134260.45 0.05 <NA> <NA> 134260.45
fitvar_temp1 = variofit(vario_temp, cov.model="exponential",ini.cov.pars = c(0.04,36921.32),nugget=0.05,wei="npairs")## variofit: covariance model used is exponential
## variofit: weights used: npairs
## variofit: minimisation function used: optimfitvar_temp2 = variofit(vario_temp, cov.model="gaussian",ini.cov.pars = c(0.04,67130.23),nugget=0.05,wei="npairs")## variofit: covariance model used is gaussian
## variofit: weights used: npairs
## variofit: minimisation function used: optimfitvar_temp3 = variofit(vario_temp, cov.model="cubic",ini.cov.pars = c(0.04,134260.45),nugget=0.05,wei="npairs")## variofit: covariance model used is cubic
## variofit: weights used: npairs
## variofit: minimisation function used: optimpar(mfrow=c(1,1))
plot(vario_temp,
xlab = "h",
ylab = "semivarianza",
cex.lab = 1.3,
cex.axis = 1.2,
main = "Estimación teórica del modelo de semivariograma de temperatura",
col.main = 4, cex.main =1.3)
lines(fitvar_temp1, col = 1)
lines(fitvar_temp2, col = 2)
lines(fitvar_temp3, col = 3)
legend(1, 1,
c("Expo","Gaussian","Cubic"),
lwd = 2,
lty = 1:3,
col = 1:3,
box.col = 9,
text.col = 1:3)
En este caso, seguramente por los pocos números de pares disponibles, la función variofit parece no dar buenos resultados, así que se utilizarán los modelos a sentimiento seleccionados. Además, por la falta de información en los semivariogramas direccionados, se opta por asumir isotropía.
kc_temp_expo=krige.conv(tempg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = c(0.04,36921.32),nugget=0.05,cov.model = "exponential"))## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhoodkc_temp_gaus=krige.conv(tempg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = c(0.04,67130.23),nugget=0.05,cov.model = "gaussian"))## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhoodkc_temp_cubic=krige.conv(tempg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = c(0.04,134260.45),nugget=0.05,cov.model = "cubic"))## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhoodboxplot(kc_temp_expo$krige.var,
kc_temp_gaus$krige.var,
kc_temp_cubic$krige.var,
names=c("Expo","Gaussian","Cubic"),
main="Comparación de modelos",
ylab="Varianza de kriging")
datos_grafica_temp = data.frame(coords, pred = kc_temp_gaus$predict, var=kc_temp_gaus$krige.var)
grilla_sf = st_as_sf(
datos_grafica_temp,
coords = c("X", "Y"),
crs = st_crs(cordoba)
)
r_temp = rast(datos_grafica_temp[, c("X", "Y", "pred", "var")],type = "xyz",crs = st_crs(cordoba)$wkt)
r_temp = crop(r_temp, vect(cordoba))
r_temp = mask(r_temp, vect(cordoba))
par(mfrow = c(2, 1),
mar = c(1,1,1,1))
r_df = as.data.frame(r_temp, xy = TRUE)
ggplot() +
geom_raster(data = r_df,
aes(x = x, y = y, fill = pred)) +
geom_sf(data = departamentos, fill = NA, color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
scale_fill_viridis_c(name = "Temperatura") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Temperatura en Cordoba mediante kriging")
ggplot() +
geom_raster(data = r_df,
aes(x = x, y = y, fill = var)) +
geom_sf(data = departamentos, fill = NA, color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
scale_fill_viridis_c(name = "Varianza") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Varianza del kriging")
Se realiza el procedimiento analogo para la humedad. En este caso, se asigna una tendencia de primer grado.
vario_0_hr = variog(hrg, trend="1st", estimator.type="modulus", direction = 0, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 0 degrees (0 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_45_hr = variog(hrg, trend="1st", estimator.type="modulus", direction = pi/4, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 45 degrees (0.785 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_90_hr = variog(hrg, trend="1st", estimator.type="modulus", direction = pi/2, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 90 degrees (1.571 radians)
## tolerance angle = 22.5 degrees (0.393 radians)vario_135_hr = variog(hrg, trend="1st", estimator.type="modulus", direction = 3*pi/4, tolerance = pi/8, max.dist=130000)## variog: computing variogram for direction = 135 degrees (2.356 radians)
## tolerance angle = 22.5 degrees (0.393 radians)par(mfrow=c(2,2), mar=c(3,3,2,2))
plot(vario_0_temp, main="Direccion 0")
plot(vario_45_temp, main="Direccion 45")
plot(vario_90_temp, main="Direccion 90")
plot(vario_135_temp, main="Direccion 135")
vario_hr=variog(hrg, trend="1st", estimator.type="modulus", pairs.min = 5, max.dist=130000)## variog: computing omnidirectional variogram#---------Mejores ajustes----------
# cov.model sigmasq phi tausq kappa kappa2 practicalRange
#1 exponential 3.59 38085.22 3 <NA> <NA> 114093.12269674
#2 gaussian 3.4 58671.82 3 <NA> <NA> 101550.264587064
#3 circular 3.4 111167.66 3 <NA> <NA> 111167.66
fitvar_hr1 = variofit(vario_hr, cov.model="exponential",ini.cov.pars = c(3.59,38085.22),nugget=3,wei="npairs")## variofit: covariance model used is exponential
## variofit: weights used: npairs
## variofit: minimisation function used: optimfitvar_hr2 = variofit(vario_hr, cov.model="gaussian",ini.cov.pars = c(3.4,58671.82),nugget=3,wei="npairs")## variofit: covariance model used is gaussian
## variofit: weights used: npairs
## variofit: minimisation function used: optimfitvar_hr3 = variofit(vario_hr, cov.model="circular",ini.cov.pars = c(3.4,111167.66),nugget=3,wei="npairs")## variofit: covariance model used is circular
## variofit: weights used: npairs
## variofit: minimisation function used: optimpar(mfrow=c(1,1))
plot(vario_hr,
xlab = "h",
ylab = "semivarianza",
cex.lab = 1.3,
cex.axis = 1.2,
main = "Estimación teórica del modelo de semivariograma de humedad",
col.main = 4, cex.main =1.3)
lines(fitvar_hr1, col = 1)
lines(fitvar_hr2, col = 2)
lines(fitvar_hr3, col = 3)
legend(1, 1,
c("Expo","Gaussian","Circular"),
lwd = 2,
lty = 1:3,
col = 1:3,
box.col = 9,
text.col = 1:3)
Por la falta de información en los semivariogramas direccionados, se opta por asumir isotropía.
kc_hr_expo=krige.conv(hrg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = fitvar_hr1$cov.pars,nugget=fitvar_hr1$nugget,cov.model = "exponential"))## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhoodkc_hr_gaus=krige.conv(hrg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = fitvar_hr2$cov.pars,nugget=fitvar_hr2$nugget,cov.model = "gaussian"))## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhoodkc_hr_circular=krige.conv(hrg,
locations = coords,
krige = krige.control(type.krige="OK", cov.pars = fitvar_hr3$cov.pars,nugget=fitvar_hr3$nugget,cov.model = "circular"))## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhoodboxplot(kc_hr_expo$krige.var,
kc_hr_gaus$krige.var,
kc_hr_circular$krige.var,
names=c("Expo","Gaussian","Circular"),
main="Comparación de modelos",
ylab="Varianza de kriging")
datos_grafica_hr = data.frame(coords, pred = kc_hr_expo$predict, var=kc_hr_expo$krige.var)
grilla_sf = st_as_sf(
datos_grafica_hr,
coords = c("X", "Y"),
crs = st_crs(cordoba)
)
r_hr = rast(datos_grafica_hr[, c("X", "Y", "pred", "var")],type = "xyz",crs = st_crs(cordoba)$wkt)
r_hr = crop(r_hr, vect(cordoba))
r_hr = mask(r_hr, vect(cordoba))
par(mfrow = c(2, 1),
mar = c(1,1,1,1))
r_df = as.data.frame(r_hr, xy = TRUE)
ggplot() +
geom_raster(data = r_df,
aes(x = x, y = y, fill = pred)) +
geom_sf(data = departamentos, fill = NA, color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
scale_fill_viridis_c(name = "Humedad relativa") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Humedad relativa en Cordoba mediante kriging")
ggplot() +
geom_raster(data = r_df,
aes(x = x, y = y, fill = var)) +
geom_sf(data = departamentos, fill = NA, color = "black") +
geom_sf(data = contornos, color = "brown", lwd = 0.2, alpha=0.6) +
scale_fill_viridis_c(name = "Varianza") +
coord_sf(xlim = c(700000, 950000), ylim = c(1280000, 1580000)) +
theme_minimal() +
labs(title = "Varianza del kriging")