выбор размера пикселя для карты оценки плотности ядра (density.ppp) в r

#r #spatial #kernel-density #spatstat

Вопрос:

У меня есть шейп-файл полигонов и шейп-файл точек, которые распределены внутри полигонов. Я создал карты оценки плотности ядра (KDE) для каждого полигона на основе содержащихся в нем точек, используя функцию density.ppp из spatstat пакета. Теперь я хочу создать разные kde с разными размерами пикселей, чтобы выбрать тот, который лучше всего подходит для моей модели. Я попытался использовать аргумент eps в as.mask функции, но это изменило только размер пикселя окна, а не самой карты ядра, поэтому результаты не изменились. после того, как я прочитал все руководство по функции density.ppp, все, что я смог найти, что выглядело связанным, — это pixellate.ppp функция из spatstat.geom пакета, но я не уверен, как ее использовать density.ppp .

Есть какие-либо предложения, как изменить размер пикселя kde?

 library(sf)
library(spatstat)

buffer <- st_read("gis/layers/buffers.shp")
pbb<- st_read("gis/layers/points_by_buffer.shp")

for (p in 1:10) {
 if(p %in% pbb$field_id) {
   
   poly123 <- pbb[pbb$field_id == p,]
   
   C <- as.owin(buffer$geometry[p])

   W<- as.mask(C, eps = 100)

   point<- ppp(poly123$X,poly123$Y, window = W)

   sigma <- bw.diggle(point)
   
   d <- density.ppp(point, kernel = "gaussian", sigma=sigma, positive = TRUE, at="pixels", )
   plot(d) 
}
 

Комментарии:

1. Плотность-это не количество точек в пикселе. Плотность-это количество точек на единицу площади — — — количество точек в одной квадратной единице — — — с использованием физической единицы длины, в которой выражены пространственные координаты. Если координаты указаны в метрах, значения плотности равны «количеству точек на квадратный метр». Смотрите более подробный ответ ниже.

Ответ №1:

Три разные вещи:

  • разрешение: расстояние между точками сетки, в которых будет рассчитываться плотность.
  • пропускная способность: величина сглаживания, которая будет применена к входным данным.
  • единица длины: физическая единица длины (например, метр, километр); единицы измерения, в которых заданы пространственные координаты; значения плотности выражаются с использованием одних и тех же единиц (плотность-это количество точек на квадратную единицу).

Если вы измените разрешение (что вы можете сделать с помощью аргументов eps или dimyx в вызове density.ppp ), значения плотности не сильно изменятся; вы просто меняете расстояние между местами, где рассчитывается плотность.

Если вы измените пропускную способность (что вы можете сделать, используя аргументы sigma и adjust в вызове density.ppp ), значения плотности существенно изменятся. При увеличении пропускной способности каждая точка данных эффективно распределяется на большее расстояние, а оценка плотности становится более плавной и плоской.

Если вы измените единицу длины, все изменится с учетом коэффициента масштабирования. Значения плотности-это «количество точек на единицу площади». Если значение плотности равно 3, это означает, что мы ожидаем увидеть в среднем 3 точки на квадратную единицу. Здесь «единица измерения» — это единица длины, в которой выражаются пространственные координаты. Если исходные координаты точек выражены в метрах, то значения плотности равны «количеству точек на квадратный метр». Если вы хотите изменить это, есть два способа:

  • масштабирование входных данных (рекомендуется): например, если исходные координаты данных указаны в метрах, вы можете преобразовать их в километры, разделив координаты на 1000. Мы рекомендуем использовать эту функцию rescale.ppp для этого: Xnew <- rescale(Xold, 1000) . Вы можете объявить имя единицы длины с помощью функции unitname (как в unitname(Xnew) <- "km" ) или с помощью rescale.ppp (как в Xnew <- rescale(Xold, 1000, "km") ). После масштабирования шаблона точек вы вызываете density.ppp (как в случае Dnew <- density(Xnew, ...) с полосой пропускания, выраженной в километрах), и результирующие значения плотности выражаются в точках на квадратный километр.
  • масштабировать результат (не рекомендуется): если ваши координаты были в метрах, и позже вы решите, что хотите узнать плотность, выраженную в точках на квадратный километр, тогда вы можете просто применить коэффициент масштабирования к изображению (введя Dnew <- D * 1000^2 , где D находится исходное изображение плотности), а также масштабировать пространственные координаты изображения (введя Dnew <- rescale(Dnew, 1000, "km") ). Мы не рекомендуем этого делать, потому что это легко испортить.

Дополнительную информацию смотрите в книге spatstat.

Ответ №2:

Изменение размера пикселя изменяет только разрешение результата, но основная математика та же самая. Это просто вопрос о том, насколько точна полученная сетка. Вы можете изменить это, используя аргументы eps или dimyx в своем вызове density.ppp() . Это, конечно, актуально, но в значительной степени ответ прост: более высокое разрешение лучше (или, по крайней мере, обычно не хуже, если вы не впадаете в крайность и не сталкиваетесь с числовой нестабильностью). Более высокое разрешение делает график предполагаемой интенсивности более приятным, но, возможно, это не то, что вы хотите, когда говорите, что хотите «выбрать тот, который лучше всего подходит для моей модели». Обычно гораздо интереснее увидеть эффект выбранной полосы пропускания сглаживания, указанной sigma . Вы только что исправили это, используя правило Диггла, но это никоим образом не является общепринятой «лучшей практикой». На самом деле выбор пропускной способности-это проблема, не имеющая простого универсального решения. Вы можете попробовать различные диапазоны частот и посмотреть, какой из них лучше всего подходит для вашей цели. Разрешение пикселей должно быть достаточно высоким, чтобы график выглядел красиво, и достаточно низким, чтобы расчет не занял слишком много времени.

Комментарии:

1. Спасибо той Эге за ответ. если бы вы могли увидеть комментарий выше и обратиться к нему, я был бы вам очень признателен.

Ответ №3:

Я думаю, вы путаете «разрешение» с «пропускной способностью».

Разрешение пикселей относится к размеру пикселей на выходе — фактически к расстоянию между точками сетки, в которых будет рассчитана оценка интенсивности. Это расстояние не имеет большого значения для большинства целей.

«Пропускная способность» — это степень сглаживания, применяемая к данным. Большая полоса пропускания создаст относительно плоскую поверхность, в то время как небольшая полоса пропускания создаст сильно заостренную поверхность. Пропускная способность действительно оказывает большое влияние, и ее необходимо выбирать соответствующим образом для данных.

Разрешение пикселей регулируется аргументами eps и dimyx может быть передано density.ppp . (По общему признанию , это трудно прочитать из файлов справки , но help(density.ppp) упоминается, что передаются несопоставимые аргументы ... pixellate.ppp , и help(pixellate.ppp) упоминается, что передаются несопоставимые аргументы ... as.mask , и, наконец, в справке as.mask описаны аргументы eps и dimyx . Именно так работает передача аргументов в R.)

Полоса пропускания сглаживания регулируется аргументом sigma . См help(density.ppp) .Информацию о том, как выбрать и управлять полосой пропускания сглаживания. Также прочитайте раздел «См. также».

Комментарии:

1. Спасибо вам за ваш ответ! Я думаю, что понимаю, что в данном случае означает разрешение и пропускная способность. я имел в виду, что хочу изменить способ расчета плотности. Если плотность-это количество точек на единицу площади, я хочу изменить размер единицы площади. возможно ли это вообще?

2. Похоже, что разрешение должно быть сделано путем регулировки «eps», но результаты, которые я получаю (значения kde), остаются прежними.

3. Плотность-это количество точек на единицу площади. Если пространственные координаты точечного рисунка заданы в метрах, то значения плотности-это количество точек на квадратный метр. Чтобы изменить это на количество точек на квадратный километр, вы должны изменить координаты шаблона точек, изменив X значение на Y <- rescale(X, 1000, "km") , а затем применить density(Y) . Пропускная способность также должна быть выражена в километрах. Смотрите справку по rescale.ppp и unitname<-.ppp