使用 spatstat 在 window 内进行插值会遗漏一个点

Question

我有一个映射森林图的数据框，其中所有树干都有 X、Y 坐标、胸高直径（以厘米为单位）和生存率 (0,1)（在此处找到，名为 "MFP14_surv_forSO.csv": https://www.dropbox.com/sh/t10b53qcobvxlzg/AACZyASgudtFLiZ79QRIjHH_a?dl=0).

我创建了一个空间点模式，然后创建了一个平滑的大茎死亡内核（im 对象）作为树冠 "gappiness"（或光穿透量）的代理。我认为大茎是指直径 >9cm 的茎，“0”存活表示茎死亡。

我想将这个 "gappiness" 的度量值插值到所有树干上。我已经知道模式中有几个点位于 window 之外，因此我使用 inside.owin() 将这些点排除在分析之外以对数据框进行子集化。

surv14 <- read.csv("MFP14_surv_forSO.csv")
win14 <- owin(poly=list(x=c(0,250,250,225,225,0),y=c(0,0,50,50,100,100)))  #specifying window extent
surv14 <- surv14[inside.owin(surv14[,1],surv14[,2], win14)==TRUE,]  #removing points outside of window
death <- surv14[!is.na(surv14$diam90) & surv14$diam90>9,]  #subsetting only large stems
death <- death[death$surv==0,]  #subsetting only the large stems that died
death.pp <- as.ppp(death,win14)  #creating point pattern from large stem death
death.fun <- Smoothfun(death.pp,sigma=10,edge=TRUE)  #smoothed kernel of large stem death
im <- as.im(death.fun)  #converting smoothed kernel into im object
ext <- im[surv14[,1:2]]  #yields pixel values of gappiness for each stem
surv14 <- cbind(surv14,ext)  #adding gappiness measure to data frame

但是当我插值时，它遗漏了一个点 - surv14 中有 4873 个观测值，而插值 ext 中只有 4872 个观测值。因此，当我尝试使用 cbind 将间隙度量绑定到我的数据框时，我收到以下错误消息：

Error in data.frame(..., check.names = FALSE) : 
arguments imply differing number of rows: 4873, 4872

我不知道如何弄清楚它遗漏了哪一点以及为什么。任何指导将不胜感激！

Answer 1

这个失败的原因有点复杂

我将尝试检查您的代码并解释发生了什么：

library(spatstat)
surv14 <- read.csv("MFP14_surv_forSO.csv")
win14 <- owin(poly=list(x=c(0,250,250,225,225,0),y=c(0,0,50,50,100,100)))  #specifying window extent
surv14 <- surv14[inside.owin(surv14[,1],surv14[,2], win14)==TRUE,]  #removing points outside of window
death <- surv14[!is.na(surv14$diam90) & surv14$diam90>9,]  #subsetting only large stems
death <- death[death$surv==0,]  #subsetting only the large stems that died
death.pp <- as.ppp(death,win14)  #creating point pattern from large stem death

此时你有一个带有两个不同标记的标记点图案： diam90 和 surv，为了避免并发症，我只会保留 diam90（因为我们知道 surv 对于所有点都是 0）：

marks(death.pp) <- marks(death.pp)$diam90

继续您的代码，您接下来创建 death.fun 这是一个 R 函数两个参数 x 和 y 并且它 returns 内核平滑树位置 x,y 处的直径必须在多边形 window 内 win14 以上定义。

death.fun <- Smoothfun(death.pp, sigma=10, edge=TRUE)  #smoothed kernel of large stem death

现在您将此函数转换为具有默认分辨率的像素图像封闭矩形上方 128x128 像素。

im <- as.im(death.fun)  #converting smoothed kernel into im object

在一个简单的图上用原来的 window 覆盖了数字近似值看起来不错：

plot(im, main = "")
plot(win14, add = TRUE, border = "green", lwd = 3)

但是，你不能完美地代表原来的window 解析度。 surv14 的点号 1850 在原来的 window 里面但不是数字近似：

inside.owin(surv14$X[1850], surv14$Y[1850], win14)
#> [1] TRUE
inside.owin(surv14$X[1850], surv14$Y[1850], Window(im))
#> [1] FALSE

然后你通过两列 data.frame 对像素图像进行子集化被解释为点列表和每个点的像素值 在图像的 window 内返回 。因此只有4872分返回：

ext <- im[surv14[,1:2]]  #yields pixel values of gappiness for each stem
length(ext)
#> [1] 4872

解决这个问题的简单方法是直接使用 death.fun 并跳过像素近似：

ext <- death.fun(surv14$X, surv14$Y)
#> Warning in ppp(X$x, X$y, window = win, check = check): data contain
#> duplicated points
length(ext)
#> [1] 4873

那么你的最后一行代码应该可以工作：

surv14 <- cbind(surv14,ext)  #adding gappiness measure to data frame

请注意，我还没有考虑使用这里的间隙测量值。我会把它留给你。

Answer 2

AJ 布朗问：

I don't know how to figure out which point it's leaving out, and why.

在命令 ext <- im[surv14[,1:2]] 中，您使用 [.im 在坐标 surv14[,1:2] 的位置提取像素图像值。 [.im 的帮助说它有一个参数 drop 指定位于像素化 window 之外的位置是应该省略（drop=TRUE）还是应该 return NA 值 (drop=FALSE)。默认值为 drop=TRUE，因此如果您在像素化 window 之外有任何位置，则生成的向量 ext 将比预期的短。所以如果你这样做

ext <- im[surv14[,1:2], drop=FALSE]

然后您可以使用 is.na(ext) 找到违规点。

在此特定示例中您不需要像素化，因此最好只使用 Smoothfun，正如 Ege 所建议的。

使用 spatstat 在 window 内进行插值会遗漏一个点

Interpolation within window using spatstat leaves out one point

r

spatstat