Método de interpolacion espacial Kriging
## Importar datos
library(readxl)
fresmo_1_Da <- read_excel("Fresno_1_Da.xlsx",sheet ="Sheet1")
plot (fresmo_1_Da[,3:2])
## Mapa de localizacion de ofertas
require(ggplot2)
## Loading required package: ggplot2

ggplot(fresmo_1_Da,aes(x=Longitude,y=Latitude,color=Temperature))+geom_point()+theme_bw()+
facet_wrap(~Date2)+ scale_color_continuous(type = "viridis")

## Filtramos la base de datos por fecha del 21
fechas=unique(fresmo_1_Da$Date2)
pos=which(fresmo_1_Da$Date2==fechas[4])
datos_m1=fresmo_1_Da[pos,]
datos_m1
## # A tibble: 394 x 9
## id_arbol Latitude Longitude Date2 `Psychro Wet Bu~ `Station Pressu~
## <chr> <dbl> <dbl> <chr> <dbl> <dbl>
## 1 1 2.38 -76.6 1/10~ 17.4 799.
## 2 2 2.38 -76.6 1/10~ 18.4 800.
## 3 3 2.38 -76.6 1/10~ 17.4 800
## 4 4 2.38 -76.6 1/10~ 16.9 800.
## 5 5 2.38 -76.6 1/10~ 17.3 800.
## 6 6 2.38 -76.6 1/10~ 18.4 799.
## 7 7 2.38 -76.6 1/10~ 16.6 799.
## 8 8 2.38 -76.6 1/10~ 17.1 799.
## 9 9 2.38 -76.6 1/10~ 17 800.
## 10 10 2.38 -76.6 1/10~ 16.4 800.
## # ... with 384 more rows, and 3 more variables: `Relative Humidity` <dbl>,
## # Crosswind <dbl>, Temperature <dbl>
require(geoR)
## Loading required package: geoR
## --------------------------------------------------------------
## Analysis of Geostatistical Data
## For an Introduction to geoR go to http://www.leg.ufpr.br/geoR
## geoR version 1.8-1 (built on 2020-02-08) is now loaded
## --------------------------------------------------------------
geo_dato=as.geodata(datos_m1,coords.col=3:2,data.col=9)
plot(geo_dato)

##Predicción Espacial Kriging
min(geo_dato$coords[,1])
## [1] -76.61415
max(geo_dato$coords[,1])
## [1] -76.61254
min(geo_dato$coords[,2])
## [1] 2.380799
max(geo_dato$coords[,2])
## [1] 2.382694
geodatos_grid=expand.grid( lat=seq(-76.61415,-76.61254,l=100),long=seq(2.380799,2.382694,l=100))
plot(geodatos_grid)
points(geo_dato$coords,col="red",pch=16)

geodatos_ko=krige.conv(geo_dato, loc=geodatos_grid,
krige= krige.control(nugget=0,trend.d="cte",
trend.l="cte",cov.pars=c(sigmasq=12.6072, phi=0.0006)))
## krige.conv: model with constant mean
## krige.conv: Kriging performed using global neighbourhood
par(mfrow=c(1,1))
image(geodatos_ko, main="kriging Predict", xlab="Lat", ylab="Long")

image(geodatos_ko, main="kriging StDv Predicted",val=sqrt(geodatos_ko$krige.var), xlab="East", ylab="North")

## Cargamos el ratser
require(raster)
## Loading required package: raster
## Loading required package: sp
mat_imagen=data.frame(geodatos_grid,geodatos_ko$predict)
temp_m1_raster=rasterFromXYZ(mat_imagen)
plot(temp_m1_raster)

finca=shapefile("Fresno_area.shp")
temp_m1_raster_finca=mask(temp_m1_raster,finca)
plot(temp_m1_raster_finca)
plot(finca,add=TRUE)

crs(temp_m1_raster_finca)=crs(finca)
crs(temp_m1_raster_finca)=CRS("+init=epsg:4326")
require(rasterVis)
## Loading required package: rasterVis
## Loading required package: lattice
## Loading required package: latticeExtra
##
## Attaching package: 'latticeExtra'
## The following object is masked from 'package:ggplot2':
##
## layer
levelplot(temp_m1_raster_finca,par.settings=BuRdTheme)

require(leaflet)
## Loading required package: leaflet
leaflet() %>% addTiles() %>% addRasterImage(temp_m1_raster_finca,opacity = 0.7, colors = "Spectral")
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## ellps WGS 84 in CRS definition: +proj=merc +a=6378137 +b=6378137 +lat_ts=0
## +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null +wktext +no_defs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
## datum WGS_1984 in CRS definition