Как сделать тесселяцию Вороного из 100 точек на единицу площади несколько раз

Question

Вкратце, я собираюсь выполнить тесселяцию Вороного на 100 баллов, создать 1000 наборов по 100 баллов и выполнить тесселяцию для каждого набора баллов.

Я создал точки как,

x=matrix(runif(100*1000),nrow=100,ncol=1000)
y=matrix(runif(100*1000),nrow=100,ncol=1000)

Базовый код для выполнения тесселяции Вороного с использованием пакета spatstat равен

dir = ppp(x=x, y=y, window = square(c(0,1)))
tess = dirichlet(dir)
plot(tess, main = "")
points(dir, pch=19, cex = 0.5)

Однако мне нужно сделать тесселяцию Вороногодля 1000 образцов и пытается создать петлю.Я хочу выбрать каждый столбец x и y и в итоге получу 1000 dir.Затем выполните тесселяцию «tess» для 1000 «dir».Также мне нужно вычислить площади ячеек вороной, используя функцию area=tile.areas(tess)

Цикл, который я создал, был

for (i in 1000) {
  dir[i] = ppp(x=x[,i], y=y[,i], window = square(c(0,1)))
}

, но все, что я получаю, это ошибки и предупреждения.У вас есть идеи, как это сделать?

Answer 1

Вам необходимо сохранить выходные данные в объекте, в этом случае давайте поместим его в список с именем dirList.Кроме того, вы должны указать последовательность для повторения.for (i in 100) ничего не делает, это должно быть for (i in 1:100)

library(deldir)
library(spatstat)

x <- matrix(runif(100 * 1000), nrow = 100, ncol = 1000)
y <- matrix(runif(100 * 1000), nrow = 100, ncol = 1000)

dirList <- list()

for (i in 1:1000){
  dirList[[i]] <- ppp(x = x[ , i], y = y[ , i], window = square(c(0, 1)))
}

Затем, чтобы построить график, вы получаете доступ к объекту с помощью [[]]

tess = dirichlet(dirList[[1]])
plot(tess, main = "")

Ко второй части вашего вопроса, используя ваш объект tess для одного набора баллов:

Кредит @ DJack для общего процесса

library(deldir)
library(spatstat)
library(maptools)
library(rgeos)

x <- matrix(runif(100 * 1000), nrow = 100, ncol = 1000)
y <- matrix(runif(100 * 1000), nrow = 100, ncol = 1000)

x <- x[ , 1]
y <- y[ , 1]

dir = ppp(x=x, y=y, window = square(c(0,1)))
tess = dirichlet(dir)

t1 <- as(tess, "SpatialPolygons")
t2 <- unionSpatialPolygons(t1, ID=rep(1,length(t1)))
vor_border <- which(gRelate(t1,t2, byid=TRUE) %in% c("2FF11F212"))

par(mfcol=c(1,3))
plot(t1)
plot(t1[vor_border,], col="red")
# you could also invert the mask to get just the inner polygons
plot(t1[-vor_border,], col="lightblue")