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

[Nicola Gambaro]

Dear Ákos,

Thank you very much for your help. 

I ran your code but unfortunately it returns the error:

Error in Ops.units(distance_matrix[as_units(point_number), , drop = TRUE],  : 
  		both operands of the expression should be "units" objects 

Any ideas how to fix this?

Cheers!

Nicola Gambaro
BSc Environmental Geoscience, First Class
Durham University 

> Message: 1
> Date: Fri, 14 Aug 2020 14:01:43 +0200
> From: =?UTF-8?Q?Bede-Fazekas_=c3=81kos?= <bfalevlist using gmail.com>
> To: r-sig-geo using r-project.org
> Subject: Re: [R-sig-Geo]  Spatial correlation between two 'sf' kriging
> 	objects
> Message-ID: <043f13cc-20e1-a768-3419-03301f8d8f7b using gmail.com>
> Content-Type: text/plain; charset="utf-8"; Format="flowed"
> 
> 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


	[[alternative HTML version deleted]]