我正在使用kde2d (MASS) 对经纬度数据创建密度图。我想知道原始数据中哪些点在特定轮廓线内。
我使用两种方法创建90%和50%的轮廓线。我想知道哪些点在90%轮廓线内,哪些点在50%轮廓线内。90%轮廓线内的点将包含所有50%轮廓线内的点。最后一步是找到在90%轮廓线内但不在50%轮廓线内的点(我不一定需要这一步的帮助)。
我创建了包含所有纬度和经度点以及90%和50%轮廓的绘图。我只想知道如何提取在90%和50%轮廓内的确切点。
我尝试查找与每行纬度和经度值相关联的z值(从kde2d的密度图中获取的高程),但没有成功。我还考虑过为数据添加一个ID列来标记每一行,然后在使用melt()之后以某种方式将其转移过去。然后,我可以简单地对具有与我想要的每个轮廓匹配的z值的数据进行子集,并根据ID列比较它们与原始BW数据中的哪些纬度和经度。
这是我所说的图片:
我想知道哪些红色点在50%轮廓(蓝色)内,哪些在90%轮廓(红色)内。
注意:很多代码都来自其他问题。向所有做出贡献的人致以最大的赞扬!
我使用两种方法创建90%和50%的轮廓线。我想知道哪些点在90%轮廓线内,哪些点在50%轮廓线内。90%轮廓线内的点将包含所有50%轮廓线内的点。最后一步是找到在90%轮廓线内但不在50%轮廓线内的点(我不一定需要这一步的帮助)。
# bw = data of 2 cols (lat and lon) and 363 rows
# two versions to do this:
# would ideally like to use the second version (with ggplot2)
# version 1 (without ggplot2)
library(MASS)
x <- bw$lon
y <- bw$lat
dens <- kde2d(x, y, n=200)
# the contours to plot
prob <- c(0.9, 0.5)
dx <- diff(dens$x[1:2])
dy <- diff(dens$y[1:2])
sz <- sort(dens$z)
c1 <- cumsum(sz) * dx * dy
levels <- sapply(prob, function(x) {
approx(c1, sz, xout = 1 - x)$y
})
plot(x,y)
contour(dens, levels=levels, labels=prob, add=T)
这里是第二个版本 - 使用 ggplot2。我希望使用这个版本来找到在90%和50%轮廓线内的点。
# version 2 (with ggplot2)
getLevel <- function(x,y,prob) {
kk <- MASS::kde2d(x,y)
dx <- diff(kk$x[1:2])
dy <- diff(kk$y[1:2])
sz <- sort(kk$z)
c1 <- cumsum(sz) * dx * dy
approx(c1, sz, xout = 1 - prob)$y
}
# 90 and 50% contours
L90 <- getLevel(bw$lon, bw$lat, 0.9)
L50 <- getLevel(bw$lon, bw$lat, 0.5)
kk <- MASS::kde2d(bw$lon, bw$lat)
dimnames(kk$z) <- list(kk$x, kk$y)
dc <- melt(kk$z)
p <- ggplot(dc, aes(x=Var1, y=Var2)) + geom_tile(aes(fill=value))
+ geom_contour(aes(z=value), breaks=L90, colour="red")
+ geom_contour(aes(z=value), breaks=L50, color="yellow")
+ ggtitle("90 (red) and 50 (yellow) contours of BW")
我创建了包含所有纬度和经度点以及90%和50%轮廓的绘图。我只想知道如何提取在90%和50%轮廓内的确切点。
我尝试查找与每行纬度和经度值相关联的z值(从kde2d的密度图中获取的高程),但没有成功。我还考虑过为数据添加一个ID列来标记每一行,然后在使用melt()之后以某种方式将其转移过去。然后,我可以简单地对具有与我想要的每个轮廓匹配的z值的数据进行子集,并根据ID列比较它们与原始BW数据中的哪些纬度和经度。
这是我所说的图片:
我想知道哪些红色点在50%轮廓(蓝色)内,哪些在90%轮廓(红色)内。
注意:很多代码都来自其他问题。向所有做出贡献的人致以最大的赞扬!