Установка выходного разрешения объекта kde

Question

Я пытаюсь выяснить, как правильно установить выходное разрешение объекта kde, созданного с помощью функции kde в библиотеке ks в R. В основном, у меня есть несколько объектов SpatialPoints, которые я подаю в kde и преобразовываю выходные данныев растр.Я хочу, чтобы ячейки этого растра имели определенное разрешение.

Вот пример использования набора данных meuse по запросу.

library(ks)
library(raster)
data(meuse)
points = data.frame(meuse$x,meuse$y)
raster(kde(points,Hlscv(points)))

Выходные данные, которые я получаю из этого кода:

class       : RasterLayer 
dimensions  : 151, 151, 22801  (nrow, ncol, ncell)
resolution  : 31.37394, 46.03558  (x, y) 
extent      : 177628.8, 182366.2, 328186.8, 335138.2  (xmin, xmax, ymin, ymax)
projection  : NA 
values      : in memory
min value   : 0 
max value   : 2.925851e-07 

Я хочу найтиспособ установить выходное разрешение (третий ряд выходных данных) на то, что я хочу.

Теперь я знаю, что у kde есть опции 'gridsize' и 'bgridsize', и для их использования вы устанавливаете количество точек / ячеекВы хотите в каждом измерении.Однако, не зная степени вывода, вы не можете рассчитать количество ячеек, чтобы получить определенное разрешение.

Одна мысль, которая у меня возникла, - это использовать значение H в соответствующем измерении для буферизации мин иМаксимальные координаты каждого измерения и предварительно получить экстент для вывода KDE.Тем не менее, я думаю, что это будет работать только с диагональными матрицами для H, и поэтому я не уверен, что это может быть реализовано с полной матрицей значений H 2x2.

Я также знаю, что вы можете пересчитать растр, но яхочу убедиться, что я не исходил из объекта kde с более низким разрешением.

Answer 1

Как указано в моем первом комментарии, я подумал, что: gridsize = диапазон / разрешение, поэтому проиллюстрировано ниже относительно тонкой сетки (выбранной так, чтобы иметь то, что, по моему мнению, было бы сопоставимым разрешением по умолчанию), а затем продемонстрировал более грубую сетку.Как отмечает blindJesse, единицы измерения «разрешения» сложнее, чем я предполагал.На этом этапе я бы посоветовал bindJesse начать составление письма в список рассылки Spatial-R:

 gridsize.x <- diff(range(meuse@coords[,"x"]))/35
 gridsize.y <- diff(range(meuse@coords[,"y"]))/35
 gridsize.x
#[1] 79.57143
 gridsize.y
# [1] 111.3429
 rimage <- raster(kde(coordinates(meuse),Hlscv(coordinates(meuse)), 
                                  gridsize=c(gridsize.x,gridsize.y) ))
 plot(rimage)

 gridsize.y <- diff(range(meuse@coords[,"y"]))/70
 gridsize.x <- diff(range(meuse@coords[,"x"]))/70
 rimage <- raster(kde(coordinates(meuse),Hlscv(coordinates(meuse)), 
         gridsize=c(gridsize.x,gridsize.y) ))
 plot(rimage)

Answer 2

Я ковырялся и заметил, что существует постоянная связь между увеличением экстента выходного растра kde по сравнению с входными точками и параметрами H.Я не знаю, почему это так (если это связано с реализацией или математикой, лежащей в основе), поэтому может не подходить для всех случаев, но работает для меня последовательно, поэтому я передаю это.В частности, я обнаружил, что (это уравнение, а не код):

(range(coordinates(output)) - range(coordinates(input))) / H = ~7.4 for each dimension.

Рассмотрим следующий код:

for (i in seq(1:10)){
    pts = data.frame(x=rnorm(100,300000,2500),y=rnorm(100,4000000,2500))
    rangePtsX <- diff(range(pts$x))
    rangePtsY <- diff(range(pts$y))
    H  <- Hlscv(pts)
    ras <- raster(kde(pts,H))
    rangeRasX <- xmax(ras) - xmin(ras)
    rangeRasY <- ymax(ras) - ymin(ras)
    rangeDiffX <- rangeRasX - rangePtsX
    rangeDiffY <- rangeRasY - rangePtsY
    print(paste("rangeDiffX/hX:",rangeDiffX/sqrt(H[1,1]),"rangeDiffY/hY:",rangeDiffY/sqrt(H[2,2])))
}

Вывод этого:

[1] "rangeDiffX/hX: 7.37104456534156 rangeDiffY/hY: 7.37763153209006"
[1] "rangeDiffX/hX: 7.39015683159216 rangeDiffY/hY: 7.39151926375274"
[1] "rangeDiffX/hX: 7.39492769120192 rangeDiffY/hY: 7.39414077521017"
[1] "rangeDiffX/hX: 7.39909462708713 rangeDiffY/hY: 7.39917867776494"
[1] "rangeDiffX/hX: 7.39801448966617 rangeDiffY/hY: 7.39779576937998"
[1] "rangeDiffX/hX: 7.39679742756067 rangeDiffY/hY: 7.39745249174806"
[1] "rangeDiffX/hX: 7.39405975028797 rangeDiffY/hY: 7.39368126656615"
[1] "rangeDiffX/hX: 7.3913950522465 rangeDiffY/hY: 7.38980236385133"
[1] "rangeDiffX/hX: 7.39988585440102 rangeDiffY/hY: 7.39988850314936"
[1] "rangeDiffX/hX: 7.39529001635855 rangeDiffY/hY: 7.39475036015628"

Оттуда было просто написать функцию для вычисления оптимального количества ячеек:

gridSize <- function(pts,H,res){
    sizeX <- (diff(range(pts[,1])) + (7.4 * sqrt(H[1,1]))) / res
    sizeY <- (diff(range(pts[,2])) + (7.4 * sqrt(H[2,2]))) / res
    c(sizeX,sizeY)
}

Это также работает для различных селекторов:

for (hMethod in c("Hlscv","Hlscv.diag","Hscv","Hpi","Hpi.diag")){
    for (i in seq(1:5)){
        pts <- data.frame(x=rnorm(100,300000,2500),y=rnorm(100,4000000,2500))
        ras <- raster(kde(pts,eval(call(hMethod,x=pts)),gridsize=gridSize(pts,H,30)))
        print(paste("xres:",xres(ras),"yres:",yres(ras)))
    }
}

производит:

    [1] "xres: 29.8498137761045 yres: 29.8456700392426"
    [1] "xres: 29.9573491524671 yres: 29.9874090657282"
    [1] "xres: 29.968525344047 yres: 29.9580897162408"
    [1] "xres: 29.9498803408057 yres: 29.964382664777"
    [1] "xres: 29.9711108728299 yres: 29.9773401860409"
    [1] "xres: 29.9950831714231 yres: 29.9658642153949"
    [1] "xres: 29.9905564968586 yres: 29.982738272666"
    [1] "xres: 29.9527729769381 yres: 29.9855591986985"
    [1] "xres: 29.9786535322154 yres: 29.9594421322198"
    [1] "xres: 29.9461263582666 yres: 29.9587155045891"
    [1] "xres: 32.1494143184879 yres: 31.9656619568261"
    [1] "xres: 31.3425525696929 yres: 31.7046601198584"
    [1] "xres: 31.9515186102478 yres: 31.9042586036464"
    [1] "xres: 31.5043829006183 yres: 30.0677630099221"
    [1] "xres: 31.1816325999623 yres: 30.782974425409"
    [1] "xres: 29.4406922338173 yres: 28.70952461795"
    [1] "xres: 31.0945577583419 yres: 31.2995894646556"
    [1] "xres: 29.3060327529272 yres: 29.4708662690499"
    [1] "xres: 29.502732054192 yres: 29.3034911939017"
    [1] "xres: 30.0529058397693 yres: 29.2893540247008"
    [1] "xres: 27.7898933275596 yres: 29.1928490584976"
    [1] "xres: 27.7628745096943 yres: 27.5794864810828"
    [1] "xres: 28.9513504972817 yres: 29.7665791290592"
    [1] "xres: 28.9569389967698 yres: 29.2103688932787"
    [1] "xres: 28.6005579365073 yres: 29.3701912564357"

, который дает результаты, очень близкие к значению 30, которое я хотел (есть ошибка, связанная с выбором целого числа ячеек сетки в kde).Если мне нужно более точное разрешение, я могу отсчитать / дезагрегировать растр отсюда.