Блог

У меня есть куча точек в 2D-пространстве, и я вычислил для них выпуклую оболочку. Теперь я хотел бы «подтянуть» оболочку, чтобы она больше не обязательно охватывала все точки. В типичной аналогии с гвоздями, вставленными в доску с резиновой лентой, я хотел бы добиться того, чтобы иметь возможность настраивать эластичность резиновой ленты и позволять гвоздям сгибаться при давлении выше некоторого предела. Это просто аналогия, здесь нет реальной физики. Это было бы отчасти связано с уменьшением площади корпуса, если бы данная точка была удалена, но не совсем, потому что могут быть две точки, которые очень близки друг к другу. Это не обязательно связано с обнаружением выбросов, потому что вы можете представить себе шаблон, в котором большая часть гвоздей будет гнуться, если они находятся на узкой линии (представьте, например, форму молотка). Все это должно выполняться достаточно быстро для тысяч точек. Любые подсказки, где я должен искать с точки зрения алгоритмов? Реализация в R была бы идеальной, но не нужна.

РЕДАКТИРОВАТЬ ПОСЛЕ КОММЕНТАРИЯ: три точки, которые я обозначил, — это те, которые имеют наибольший потенциал для уменьшения площади корпуса, если их исключить. На участке нет другого набора из трех точек, который привел бы к большему уменьшению площади. Наивная реализация того, что я ищу, возможно, заключалась бы в случайной выборке некоторой доли точек, вычислении площади корпуса, итеративном удалении каждой точки на корпусе, пересчете площади, повторении много раз и удалении точек, которые, как правило, приводят к значительному уменьшению площади. Может быть, это может быть реализовано в каком-то случайном варианте леса? Это не совсем правильно, потому что я бы хотел, чтобы точки удалялись одна за другой, чтобы вы получили следующий результат. Если бы вы посмотрели на все точки сразу, возможно, было бы лучше обрезать от краев «головки молотка».

Предположим, у меня есть набор точек, подобных этому:

 set.seed(69)

x <- runif(20)
y <- runif(20)

plot(x, y)
  

Тогда легко найти точки подмножества, которые находятся на выпуклой оболочке, выполнив:

 ss <- chull(x, y)
  

Это означает, что мы можем построить выпуклую оболочку, выполнив:

 lines(x[c(ss, ss[1])], y[c(ss, ss[1])], col = "red")
  

Теперь мы можем случайным образом удалить одну из точек, расположенных на выпуклой оболочке (т.Е. «Согнуть гвоздь»), выполнив:

 bend <- ss[sample(ss, 1)]

x <- x[-bend]
y <- y[-bend]
  

И затем мы можем повторить процесс нахождения выпуклой оболочки этого нового набора точек:

 ss <- chull(x, y)
lines(x[c(ss, ss[1])], y[c(ss, ss[1])], col = "blue", lty = 2)
  

Чтобы получить точку, которая при удалении приведет к наибольшему уменьшению площади, одним из вариантов будет следующая функция:

 library(sp)

shrink <- function(coords)
{
  ss <- chull(coords[, 1], coords[, 2])
  outlier <- ss[which.min(sapply(seq_along(ss), 
                function(i) Polygon(coords[ss[-i], ], hole = FALSE)@area))]
  coords[-outlier, ]
}
  

Итак, вы могли бы сделать что-то вроде:

 coords <- cbind(x, y)

new_coords <- shrink(coords)

new_chull <- new_coords[chull(new_coords[, 1], new_coords[, 2]),]
new_chull <- rbind(new_chull, new_chull[1,])

plot(x, y)

lines(new_chull[,1], new_chull[, 2], col = "red")
  

Конечно, вы могли бы сделать это в цикле, чтобы оно new_coords возвращалось shrink несколько раз.

Вычислите надежный центр и дисперсию, используя mcd.cov in MASS и расстояние махаланобиса каждой точки от него (используя mahalanobis in psych). Затем мы показываем график квантилей расстояний Махаланобиса с использованием PlotMD from modi, а также показываем связанные выбросы красным цветом на втором графике. (В modi есть и другие функции, которые также могут представлять интерес.)

 library(MASS)
library(modi)
library(psych)

set.seed(69)   
x <- runif(20)
y <- runif(20)
m <- cbind(x, y)

mcd <- cov.mcd(m)
md <- mahalanobis(m, mcd$center, mcd$cov)
stats <- PlotMD(md, 2, alpha = 0.90)
  

предоставление:

(продолжение после скриншота)

и мы показываем выпуклую оболочку с помощью линий, а выбросы выделены красным:

 plot(m)    
ix <- chull(m)
lines(m[c(ix, ix[1]), ])

wx <- which(md > stats$halpha)
points(m[wx, ], col = "red", pch = 20)
  

Спасибо вам обоим! Я пробовал различные методы обнаружения выбросов, но это не совсем то, что я искал. Они плохо работали из-за странных форм моих кластеров. Я знаю, что говорил о площади выпуклой оболочки, но я думаю, что фильтрация по длинам сегментов дает лучшие результаты и ближе к тому, что я действительно хотел. Тогда это выглядело бы примерно так:

 shrink <- function(xy, max_length = 30){
      to_keep <- 1:(dim(xy)[1])
      centroid <- c(mean(xy[,1]), mean(xy[,2]))

      while (TRUE){
        ss <- chull(xy[,1], xy[,2])
        ss <- c(ss, ss[1])
        lengths <- sapply(1:(length(ss)-1), function(i) sum((xy[ss[i 1],] - xy[ss[i],])^2))

        # This gets the point with the longest convex hull segment. chull returns points
        # in clockwise order, so the point to remove is either this one or the one
        # after it. Remove the one furthest from the centroid.
        max_point <- which.max(lengths)
        if (lengths[max_point] < max_length) return(to_keep)

        if (sum((xy[ss[max_point],] - centroid)^2) > sum((xy[ss[max_point   1],] - centroid)^2)){
          xy <- xy[-ss[max_point],]
          to_keep <- to_keep[-ss[max_point]]
        }else{
          xy <- xy[-ss[max_point   1],]
          to_keep <- to_keep[-ss[max_point   1]]
        }
      }
    }
  

Это не оптимально, потому что оно учитывает расстояние до центра тяжести, чего я бы хотел избежать, и есть max_length параметр, который должен вычисляться из данных, а не быть жестко запрограммированным.

Нет фильтра:

Это выглядит так, потому что здесь 500 000 ячеек, и многие из них оказываются «неправильными» при проецировании с ~ 20 000 измерений на 2.

Фильтр:

Обратите внимание, что он отфильтровывает точки на концах некоторых кластеров. Это не оптимально, но нормально. Перекрытие между некоторыми кластерами является истинным и должно быть.