[R-sig-Geo] Spatial correlation between two 'sf' kriging objects

[Bede-Fazekas Ákos]

Dear Nicola,

Instead of raster::focal(), you can apply the cor() function in 
combination with st_distance() (or st_buffer() and st_within()). E.g. 
let's say that 'grid' is a POINT type sf object containing columns 
'temperature' and 'yield':
distance_threshold <- 40
distance_matrix <- st_distance(grid)
grid$correlation <- vapply(X = 1:nrow(grid), FUN.VALUE = numeric(1), FUN 
= function (point_number) {cor(x = st_set_geometry(grid, 
NULL)[distance_matrix[point_number, distance_matrix[point_number, , drop 
= TRUE] < distance_threshold , drop = TRUE], "temperature"], y = 
st_set_geometry(grid, NULL)[distance_matrix[point_number, 
distance_matrix[point_number, , drop = TRUE] < distance_threshold , drop 
= TRUE], "yield"])})

Or something like this... I have not tested this code, and am sure that 
it is not the most efficient solution.

For large, square grids, raster might be faster than sf. You can convert 
your grid to RasterLayer with function rasterFromXYZ() combined with 
st_coordinates(), st_set_geometry(, NULL) and cbind.data.frame(). There 
might be more straightforward solutions for the conversion...

HTH,
Ákos Bede-Fazekas
Hungarian Academy of Sciences


2020.08.13. 20:11 keltezéssel, Nicola Gambaro írta:
> I have created two ‘sf’ kriging objects (point vectors), one for temperature and another for agricultural yields. To make the grid and carry out the point interpolation, I have remained within the ‘sf’ package.
>
> I would now like to create a spatial local correlation ‘raster’ between these two variables, as shown on this webpage https://statnmap.com/2018-01-27-spatial-correlation-between-rasters/ <https://statnmap.com/2018-01-27-spatial-correlation-between-rasters/>. However, in that example, they use the ‘raster’ package and the ‘focal’ function. I was wondering if there was a way of doing this within ‘sf’, i.e. without having to change classes? If not, what is the best way to convert those objects into raster classes?
>
> Here is an excerpt of my kriging code for reference:
> library(sf)
> sf_data <- st_as_sf(x = data, coords = c("longitude", "latitude"), crs = 4326)
> library(gstat)
> vgm_utci <- variogram(UTCI~1, sf_data)
> utci_fit <- fit.variogram(vgm_utci, vgm("Gau"), fit.kappa = TRUE)
> plot(vgm_utci, utci_fit)
> istria <- read_sf(“./Istria_Boundary.shp")
> istria <- istria$geometry
> istria.grid <- istria %>%
>    st_make_grid(cellsize = 0.05, what = "centers") %>%
>    st_intersection(istria)
> library(ggplot2)
> ggplot() + geom_sf(data = istria) + geom_sf(data = istria.grid)
> library(stars)
> utci_krig <- krige(formula = sf_data$UTCI ~ 1, locations = sf_data,
>                     newdata = istria.grid, model = utci_fit)
>
>
> Thank you very much in advance,
>
> Nicola
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo using r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo