0

我是一个相当新的 R 用户，使用嵌入在美国都市区域内的地方的空间数据。我正在为许多变量生成 IDW 地图，但在定义 IDW 的边界时遇到了麻烦。我只想在都市区边界范围内插入数据，而我目前的方法填充了整个图。

我首先ppp从 a 创建一个对象dataframe：

x.points<-SpatialPointsDataFrame(coords =x.data[,c("LON10","LAT10")],
  data=x.data[,c("v1","v2")],
  proj4string=CRS("+proj=longlat +ellps=GRS80 +datum=NAD83 +no_defs"))

x.proj<-spTransform(x.points,CRSobj = CRS("+proj=aea
  +lat_1=29.5 +lat_2=45.5 +lat_0=37.5 +lon_0=-96 +x_0=0 +y_0=0 
  +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0"))    

x.ppp<-as.ppp(x.proj)

使用投影此ppp对象plot(x.ppp)会显示每个变量的图（上面代码中的“v1”和“v2”）。但是，绘图填充了整个矩形框。在绘制对象时，我ppp通过为其分配owin对象来限制此框架，如下所示：

y.spdf<-spTransform(y.spdf,CRSobj = CRS("+proj=aea 
  +lat_1=29.5 +lat_2=45.5 +lat_0=37.5 +lon_0=-96 +x_0=0 +y_0=0 
  +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=0,0,0"))

y.owin<-as.owin(y.spdf)

x.ppp<-x.ppp[y.owin]

现在使用它来绘制这个plot(x.ppp)图，可以让我在他们的都市区的多边形内正确地绘制出我的数据点图。虽然这是一个好的开始，但我想继续将此空间框架传递到 IDW 地图上。查阅 spatstat 手册让我相信我可以将我的owin对象作为as.mask参数传递，但是每当我这样做时，都会出错。示例代码如下：

idw<-idw(x.ppp, power=4, as.mask(y.owin))

这给了我错误“开关错误（在，像素= {：EXPR必须是长度为1的向量”。使用at="points"不起作用，因为我得到的是结果dataframe对象而不是imlist。如果我省略了as.mask参数，我会得到一个正确的 idw，但它插值超出了我之前能够使用该owin对象设置的都市区边界。

我认为这是一项相当常见的任务，但查阅 spatstat 文档让我有点不知所措。这个as.mask论点似乎可行，我不确定我哪里出错了。关于如何解决这个问题的任何意见或建议，或者我可以求助于尝试更好地理解这个问题的来源，都会有很大的帮助。谢谢！

4

1 回答 1

0

我明白这有点令人困惑。基本上idw计算包含数据的整个矩形框（的框架owin）中的结果。如果你这样做，这就是你得到的结果rslt <- idw(x.ppp)。为了将结果限制在多边形窗口，我们需要通过这个窗口对结果图像进行子集化（并设置drop=FALSE为得到结果不规则图像，而不是多边形区域内的像素值）。在代码中解释：

简单案例（一个数字标记值）

library(spatstat)
# Fake data with one numeric mark
X <- humberside
marks(X) <- runif(npoints(X))
# Window is a polygonal region:
W <- Window(X)
W
# IDW in the frame of W
rslt <- idw(X)
# Subset of IDW inside W
rslt <- rslt[W, drop = FALSE]
plot(rslt)

一个变量的 IDW

多变量大小写（两个或多个数字标记值）

library(spatstat)
# Fake data with more numeric marks
X <- humberside
marks(X) <- cbind(x = runif(npoints(X)), y = runif(npoints(X)))
# Window is a polygonal region:
W <- Window(X)
W
# IDW in the frame of W for each mark variable (imlist object)
rslt <- idw(X)
# Subset of IDW inside W for each mark variable (imlist object)
rslt <- solapply(rslt, "[.im", W, drop = FALSE)
plot(rslt)

两个变量的 IDW

于 2016-02-25T12:47:54.237 回答