在 R 中的 2D 直方图中查找高密度区域答案

【问题标题】：Finding high density zones in a 2D Histogram in R在 R 中的 2D 直方图中查找高密度区域
【发布时间】：2018-09-24 16:44:14
【问题描述】：

我使用 MASS 的 kde2d 函数按以下方式创建了一系列 2D 直方图：

    # Loading libraries
    library(MASS)
    library(RcolorBrewer)
    # Loading data
    data <- as.matrix(read.table('data.dat'))
    # Create the 2dhist object      
    hist_2d <- kde2d(data[,1],data[,2],n = 60, lims=c(-180,180,-180,180))
    # Define the color palette
    rf <- colorRampPalette(rev(brewer.pal(11,'Spectral')))
    r <- rf(60)
    # Defining the axis
    at_x = seq(-180,180,by=30)
    at_y = seq(-180,180,by=30)
    # Plot the 2DHistogram
    image(hist_2d,col=r,cex.main=3,main='Q68L',axes=F)
    axis(1,lwd.ticks=2,at=at_x,labels=T,cex.axis=2)
    axis(2,lwd.ticks=2,at=at_y,labels=T,cex.axis=2)

直方图像这样生成looks。如何识别所有高密度区域（我在white squares 内标记）？这个问题的理想解决方案是一个函数，它为每个高密度区域抛出一个 (x,y) 范围，以便它可以应用于多个数据集。

提前致谢，如果您需要更多信息，请告诉我

【问题讨论】：

标签： r mass

【解决方案1】：

通过正确的数据表示，这可以通过聚类分析。由于您不提供数据，我将说明 kde2d 帮助页面上使用的数据 - 间歇泉数据。该数据给出了“高密度”的非常清晰的分离区域（如您的示例图片），所以我将使用一个简单的 k-means 聚类。

library(MASS)
attach(geyser)
f2 <- kde2d(duration, waiting, n = 50, lims = c(0.5, 6, 40, 100),
            h = c(width.SJ(duration), width.SJ(waiting)) )
image(f2, zlim = c(0, 0.05))

我们需要找到“热点”。为了了解有关什么值应该被认为是“高”，我们可以看一个箱线图。

boxplot(as.vector(f2$z))

基于此，我会有些武断地使用点 z 值大于 0.012。你需要调整这个你的特殊问题。

Hot = which(f2$z > 0.012, arr.ind = TRUE)
HotPoints = data.frame(x=f2$x[Hot[,1]], y=f2$y[Hot[,2]])
plot(HotPoints, pch=20, xlim = c(0.5,6), ylim = c(40,100))

现在我们需要对点进行聚类并找到 x 和 y 范围对于集群。首先，我做的很简单，并表明结果是合理的。

KM3 = kmeans(scale(HotPoints), 3)
plot(HotPoints, pch=20, xlim = c(0.5,6), ylim = c(40,100))
for(i in 1:3) {
    Rx = range(HotPoints[KM3$cluster == i,1])
    Ry = range(HotPoints[KM3$cluster == i,2])
    polygon(c(Rx, rev(Rx)), rep(Ry, each=2))
}

我不确定您希望如何将结果呈现给您，但是将它们全部集中在一个地方的一种方法是：

XRanges = sapply(unique(KM3$cluster), 
    function(i) range(HotPoints[KM3$cluster == i,1]))
XRanges
         [,1]     [,2]     [,3]
[1,] 3.979592 3.867347 1.734694
[2,] 4.877551 4.316327 2.071429
YRanges = sapply(unique(KM3$cluster), 
    function(i) range(HotPoints[KM3$cluster == i,2]))
YRanges
         [,1]     [,2]     [,3]
[1,] 47.34694 70.61224 73.06122
[2,] 62.04082 87.75510 95.10204

这给出了三个集群中每个集群的 x 和 y 的最小值和最大值。

但是，我在这里做了一些选择，我想指出我还给你留了一些工作。你还需要做什么：
1.你需要选择一个分界点来决定密度有多高需要得到一个集群。
2. 鉴于你的分界线以上的点，你需要说您要生成多少个集群。

其余的机器都在那里。

【讨论】：