Нахождение зон высокой плотности в 2D гистограмме в R

Question

У меня есть серия 2D гистограмм, которые я создал с помощью функции MASS для kde2d следующим образом:

    # 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)

Сгенерированная гистограмма выглядит следующим образом .Как я могу определить все зоны с высокой плотностью (которые я отметил внутри белых квадратов )?Идеальным решением этой проблемы была бы функция, которая выбрасывает диапазон (x, y) для каждой зоны высокой плотности, чтобы ее можно было применять в нескольких наборах данных.

Заранее спасибо и дайте мне знать, если вам нужна дополнительная информация

Answer 1

При правильном представлении данных это можно сделать с помощью кластерного анализа.Поскольку вы не предоставляете данные, я проиллюстрирую их данными на странице справки kde2d - данными гейзера.Эти данные дают довольно четкое разделение областей «высокой плотности» (как на ваших примерах рисунков), поэтому я просто буду использовать простую кластеризацию k-средних.

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Принимая во внимание точки над вашей границей, вам нужно будет указать, сколько кластеров вы хотите сгенерировать.

Остальные машины там.