r - 如何提取等高线图的特定边界？

Question

考虑以下代码提供的计数图:

# input data x,y,z:
x<-rep(c(1:100),times=10)
y<-rep(c(1:10),each=100)
z<-NULL
for(i in 1:10){
  n<-sample(c(10:30),1)
  m<-sample(c(50:70),1)
  z<-c(z,c(1:n,sample(c(50:100),100-m-n),c(m:1)))}

# bilinear interpolation of irregular data:
library('akima')
interpo<-interp(x=x,y=y,z=z,xo=seq(min(x),max(x),length.out=1000),yo=seq(min(y),max(y),length.out=1000))
image(interpo)

# convert data format for ggplot:
x<-rep(interpo$x,times=1000)
y<-rep(interpo$y,each=1000)
r<-NULL
for(i in 1:1000){
  r<-c(r,interpo$z[,i])
}
plo<-data.frame(x,y,r)

# plot the countour plot:
library(ggplot2) 
library(scales)
ggplot(plo, aes(y = y, x = x,  fill = r)) +
  geom_raster()+
  scale_fill_gradient(low="blue",high="red",limits=c(min(r),max(r)))

剧情大概是这样的:

现在我的问题是:

使用 R，是否可以提取和可视化 r 值大于 70 的某些边界(例如 r=70)？

最好能推导出边界的具体位置(z(x,y)的x、y坐标值不小于70)。

Answer 1

不太确定您的预期输出(如果您添加一些会很好)，但以下可能接近:

# visualize the points where r = 70
ggplot(plo, aes(y = y, x = x,  fill = r)) +
  geom_raster()+
  scale_fill_gradient(low="blue",high="red",limits=c(min(r),max(r))) + 
  geom_point(data=plo[plo$r == 70,], col='white')      

或

# visualize the points where r > 70
ggplot(plo, aes(y = y, x = x,  fill = r)) +
  geom_raster()+
  scale_fill_gradient(low="blue",high="red",limits=c(min(r),max(r))) + 
  geom_point(data=plo[plo$r > 70,])

如果我们想(从 ggplot)取回 r > 70 的数据，我们可以尝试以下操作:

p <- ggplot(plo, aes(y = y, x = x,  fill = r)) +
  geom_raster()+
  scale_fill_gradient(low="blue",high="red",limits=c(min(r),max(r))) + 
  geom_point(data=plo[plo$r > 70,])

pg <- ggplot_build(p)
str(pg)
head(pg$data[[2]])
# fill        x y PANEL group shape colour size alpha stroke
#1 #E50056 26.96396 1     1    -1    19  black  1.5    NA    0.5
#2 #E70052 27.06306 1     1    -1    19  black  1.5    NA    0.5
#3 #E70050 27.16216 1     1    -1    19  black  1.5    NA    0.5
#4 #E8004F 27.26126 1     1    -1    19  black  1.5    NA    0.5
#5 #E9004D 27.36036 1     1    -1    19  black  1.5    NA    0.5
#6 #E9004C 27.45946 1     1    -1    19  black  1.5    NA    0.5

dplyr 的另一次尝试:

data <- data.frame(x=x, y=y, z=z)
dim(data)
#[1] 1000000       3
library(dplyr)
data <- plo %>% 
  inner_join(data, by=c('x'='x', 'y'='y')) %>% 
  filter(z >= 70 & r >= 70) # change the filter condition if needed
dim(data)
#[1] 10058     4
ggplot(plo, aes(y = y, x = x,  fill = r)) +
  geom_raster()+
  scale_fill_gradient(low="blue",high="red",limits=c(min(r),max(r))) + 
  geom_point(data=data)