Search code examples
gisrasterspatialr-rasterterra

ELSA package (spatial association): My categorical raster data do not behave as intended

I've started working with the elsa package (https://cran.r-project.org/web/packages/elsa/vignettes/elsa.html) to measure entropy-based spatial association for categorical data, but it revolves a lot around raster data, of which I have little previous experience.
My data comes from a .tif file whose pixel values I have reclassified to 5 categories in QGIS.

library(raster)

r <- new("RasterLayer", file = new(".RasterFile", name = "C:\\Users\\bluew\\Desktop\\stefano\\ELSA\\are_1.81_3003.tif", 
                              datanotation = "FLT4S", byteorder = "little", nodatavalue = -Inf, 
                              NAchanged = FALSE, nbands = 1L, bandorder = "BIL", offset = 0L, 
                              toptobottom = TRUE, blockrows = c(rows = 1L), blockcols = c(cols = 1062L), 
                              driver = "gdal", open = FALSE), data = new(".SingleLayerData", 
                                                                         values = logical(0), offset = 0, gain = 1, inmemory = FALSE, 
                                                                         fromdisk = TRUE, isfactor = FALSE, attributes = list(), haveminmax = TRUE, 
                                                                         min = 1, max = 5, band = 1L, unit = "", names = "are_1.81_3003"), 
    legend = new(".RasterLegend", type = character(0), values = logical(0), 
                 color = logical(0), names = logical(0), colortable = logical(0)), 
    title = character(0), extent = new("Extent", xmin = 1478031.1565, 
                                       xmax = 1563459.1172, ymin = 4933676.6029, ymax = 5009134.9516), 
    rotated = FALSE, rotation = new(".Rotation", geotrans = numeric(0), 
                                    transfun = function () 
                                      NULL), ncols = 1062L, nrows = 662L, crs = new("CRS", 
                                                                                    projargs = "+proj=tmerc +lat_0=0 +lon_0=9 +k=0.9996 +x_0=1500000 +y_0=0 +ellps=intl +units=m +no_defs"), 
    srs = "+proj=tmerc +lat_0=0 +lon_0=9 +k=0.9996 +x_0=1500000 +y_0=0 +ellps=intl +units=m +no_defs", 
    history = list(), z = list())

r <- setMinMax(r,force=T)

I submit this to the elsa function:

library(elsa)

e <- elsa(r,d=500,categorical = T) #the d argument defines the neighbourhood of each pixel

This should result in a RasterLayer object which contains the ELSA values for each pixel, but my result is somewhat strange:

library(tmap)

cl <- colorRampPalette(c('darkblue','yellow','red','black'))(100)

tm_shape(e) + tm_raster(palette= cl, title = "ELSA",style="cont") +
tm_layout(legend.outside = T)

As you can see below, it looks as though the values are only being calculated at the edges of each category, which doesn't make sense.

Example

The examples provided in the vignette I've linked above look nothing like this, but there's no further documentation on how to prepare the raster data, and as I've said I have very little experience on how to handle them.
It might very well have something to do with the way I have reclassified the raster in QGIS, but I really can't tell; I've tried to compare my raster with the ones provided in the package, but I can't spot significant differences.
Could anybody point me in the direction of what might be causing this result? Thank you to anyone who gives advice.

Solution

  • As it turns out, the "border" effect was simply a consequence of the raster resolution being too high, which was not ideal for the type of data which I am working on, because the different categories tend to cover large areas; I changed the resolution of the raster files in QGIS to a more reasonable amount (appr. 2500 meters) and the result is much more in line with what I was expecting. enter image description here