Datoos de Aguacates

Analisis Exploratorio

Combirtiendola en datos espaciales

df_aguacate_geo <- as.geodata(df_aguacate, coords.col = 2:3, data.col = 9) # Variable regionalizada
plot(df_aguacate_geo)

En este grafico podemos observar en azul las temperaturas más bajas de la zona y en rojo las más altas. Amarillo y verde serian las temperaturas más intermedias y son las que predominan en mayor proporción en la zona.

Analizar si hay tendencia o es estacionaria

Construimos el semivariograma para analizar si existe una correlacion espacial significativa, al observar el comportamiento del semivariograma observado con la variable temperatura, vemos que podriamos pensar que si tienen una correlacion significativa dado que sus puntos estan por fuera de la banda, sin esta no pareciera ser muy fuerte dado que en algun punto la linea pareciera estabilizarse dentro de las bandas.

variograma = variog(geodata = df_aguacate_geo, uvec = seq(0,0.0009,0.00005), option = "bin")
## variog: computing omnidirectional variogram
datos_env = variog.mc.env(df_aguacate_geo, obj.variog =variograma)
## variog.env: generating 99 simulations by permutating data values
## variog.env: computing the empirical variogram for the 99 simulations
## variog.env: computing the envelops
plot(variograma, pch=16, envelop = datos_env)

## los puntos negros son el variograma observado
# Las bandas me permiten determinar cuando  el semi-variograma no tiene correlacion si este tuviera su comportamiento dentro de las bandas, entonces podriamos pensar que no hay correlacion espacial

Sin embargo, se procede a analizar que modelo podria ajustarse al comportamiento mostrado por la variable anteriormente.

Ajuste del modelo de semivarianza

Se ajusta el modelo exponencial, gausiano y esferico con la finalizadad de ver cual se ajusta mejor obteniendo la siguiente grafica:

plot(variograma, pch=16, envelop = datos_env)
ini.vals = expand.grid( seq(2,4,l=10), seq(0.0001,0.00025, l = 10))

modelo_mco_exp = variofit(variograma, ini.vals, cov.model = "exponential") # Modleo exponencial
## variofit: covariance model used is exponential 
## variofit: weights used: npairs 
## variofit: minimisation function used: optim 
## variofit: searching for best initial value ... selected values:
##               sigmasq phi   tausq kappa
## initial.value "3.11"  "0"   "0"   "0.5"
## status        "est"   "est" "est" "fix"
## loss value: 3555.92876863577
modelo_mco_gau = variofit(variograma, ini.vals, cov.model = "gaussian") # Modelo gausianico
## variofit: covariance model used is gaussian 
## variofit: weights used: npairs 
## variofit: minimisation function used: optim 
## variofit: searching for best initial value ... selected values:
##               sigmasq phi   tausq kappa
## initial.value "3.11"  "0"   "0"   "0.5"
## status        "est"   "est" "est" "fix"
## loss value: 6633.5155506537
modelo_mco_shp = variofit(variograma, ini.vals, cov.model = "spherical") # Modelo Esferico
## variofit: covariance model used is spherical 
## variofit: weights used: npairs 
## variofit: minimisation function used: optim 
## variofit: searching for best initial value ... selected values:
##               sigmasq phi   tausq kappa
## initial.value "3.11"  "0"   "0"   "0.5"
## status        "est"   "est" "est" "fix"
## loss value: 6005.01775448873
lines(modelo_mco_exp, col = "blue")
lines(modelo_mco_gau, col = "red")
lines(modelo_mco_shp, col = "green")

legend(0.0007, 7, legend=c("Gausiano", "Exponencial", "Esferico"),
       col=c("red", "blue","green"), lty=1, cex=0.8)

Pareciendo que graficamente el modelo que mejor se ajusta es al semivariograma observado, seria el modelo exponencial, sin embargo este dato se corrobora comparando el MSE de de cada uno de los modelos obteniendo lo siguiente:

print(paste("MSE modelo Exponencial =", as.character(modelo_mco_exp$value)))
## [1] "MSE modelo Exponencial = 3314.17633638709"
print(paste("MSE modelo Gausiano =", as.character(modelo_mco_gau$value)))
## [1] "MSE modelo Gausiano = 6339.97796271843"
print(paste("MSE modelo Esferico =", as.character(modelo_mco_shp$value)))
## [1] "MSE modelo Esferico = 5837.58927236167"

Con esto evidenciamos que efectivamente el modelo exponencial es el que mejor describe el comportamiento del semivariograma observado con un MSE de 3.314,18

Es por esto que para realizar la prediccion espacial usaremos este con los siguientes parametros:

## variofit: model parameters estimated by WLS (weighted least squares):
## covariance model is: exponential
## parameter estimates:
##   tausq sigmasq     phi 
##  0.0270  3.1344  0.0001 
## Practical Range with cor=0.05 for asymptotic range: 0.0003183712
## 
## variofit: minimised weighted sum of squares = 3314.176

Con estos parametros se realiza la prediccion espacial, dando lo siguiente:

Prediccion espacial

datos_grid = expand.grid(Latitude = seq(2.392,2.394, l=100),Longitud = seq(-76.7118,-76.71022, l=100))
datos_ko = krige.conv(df_aguacate_geo, loc = datos_grid, krige = krige.control(nugget = 0, trend.d = "cte", trend.l = "cte", cov.pars = c(sigmasq = 3.1344, phi = 0.0001)))
## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhood
image(datos_ko, main = "Kriging Predict - Temperatura", xlab ="Latitude", ylab ="Longitude")

Donde las zonas marcadas en rojos, corresponde a las zonas de mayor temperaturas que observamos en el analisis exploratorio.