Обнаружение провалов в 2D-графике - PullRequest
9 голосов
/ 15 июня 2011

Мне нужно автоматически определять провалы на двухмерном графике, например, области, отмеченные красными кружками на рисунке ниже. Меня интересуют только «основные» провалы, это означает, что провалы должны охватывать минимальную длину по оси x. Количество провалов неизвестно, т.е. на разных графиках будет разное количество провалов. Есть идеи?

Dips in a 2D plot

Обновление:

По запросу, вот пример данных вместе с попыткой сгладить их, используя медианную фильтрацию, как предложено лозами. * * 1010

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

y <- c(0.9943,0.9917,0.9879,0.9831,0.9553,0.9316,0.9208,0.9119,0.8857,0.7951,0.7605,0.8074,0.7342,0.6374,0.6035,0.5331,0.4781,0.4825,0.4825,0.4879,0.5374,0.4600,0.3668,0.3456,0.4282,0.3578,0.3630,0.3399,0.3578,0.4116,0.3762,0.3668,0.4420,0.4749,0.4556,0.4458,0.5084,0.5043,0.5043,0.5331,0.4781,0.5623,0.6604,0.5900,0.5084,0.5802,0.5802,0.6174,0.6124,0.6374,0.6827,0.6906,0.7034,0.7418,0.7817,0.8311,0.8001,0.7912,0.7912,0.7540,0.7951,0.7817,0.7644,0.7912,0.8311,0.8311,0.7912,0.7688,0.7418,0.7232,0.7147,0.6906,0.6715,0.6681,0.6374,0.6516,0.6650,0.6604,0.6124,0.6334,0.6374,0.5514,0.5514,0.5412,0.5514,0.5374,0.5473,0.4825,0.5084,0.5126,0.5229,0.5126,0.5043,0.4379,0.4781,0.4600,0.4781,0.3806,0.4078,0.3096,0.3263,0.3399,0.3184,0.2820,0.2167,0.2122,0.2080,0.2558,0.2255,0.1921,0.1766,0.1732,0.1205,0.1732,0.0723,0.0701,0.0405,0.0643,0.0771,0.1018,0.0587,0.0884,0.0884,0.1240,0.1088,0.0554,0.0607,0.0441,0.0387,0.0490,0.0478,0.0231,0.0414,0.0297,0.0701,0.0502,0.0567,0.0405,0.0363,0.0464,0.0701,0.0832,0.0991,0.1322,0.1998,0.3146,0.3146,0.3184,0.3578,0.3311,0.3184,0.4203,0.3578,0.3578,0.3578,0.4282,0.5084,0.5802,0.5667,0.5473,0.5514,0.5331,0.4749,0.4037,0.4116,0.4203,0.3184,0.4037,0.4037,0.4282,0.4513,0.4749,0.4116,0.4825,0.4918,0.4879,0.4918,0.4825,0.4245,0.4333,0.4651,0.4879,0.5412,0.5802,0.5126,0.4458,0.5374,0.4600,0.4600,0.4600,0.4600,0.3992,0.4879,0.4282,0.4333,0.3668,0.3005,0.3096,0.3847,0.3939,0.3630,0.3359,0.2292,0.2292,0.2748,0.3399,0.2963,0.2963,0.2385,0.2531,0.1805,0.2531,0.2786,0.3456,0.3399,0.3491,0.4037,0.3885,0.3806,0.2748,0.2700,0.2657,0.2963,0.2865,0.2167,0.2080,0.1844,0.2041,0.1602,0.1416,0.2041,0.1958,0.1018,0.0744,0.0677,0.0909,0.0789,0.0723,0.0660,0.1322,0.1532,0.1060,0.1018,0.1060,0.1150,0.0789,0.1266,0.0965,0.1732,0.1766,0.1766,0.1805,0.2820,0.3096,0.2602,0.2080,0.2333,0.2385,0.2385,0.2432,0.1602,0.2122,0.2385,0.2333,0.2558,0.2432,0.2292,0.2209,0.2483,0.2531,0.2432,0.2432,0.2432,0.2432,0.3053,0.3630,0.3578,0.3630,0.3668,0.3263,0.3992,0.4037,0.4556,0.4703,0.5173,0.6219,0.6412,0.7275,0.6984,0.6756,0.7079,0.7192,0.7342,0.7458,0.7501,0.7540,0.7605,0.7605,0.7342,0.7912,0.7951,0.8036,0.8074,0.8074,0.8118,0.7951,0.8118,0.8242,0.8488,0.8650,0.8488,0.8311,0.8424,0.7912,0.7951,0.8001,0.8001,0.7458,0.7192,0.6984,0.6412,0.6516,0.5900,0.5802,0.5802,0.5762,0.5623,0.5374,0.4556,0.4556,0.4333,0.3762,0.3456,0.4037,0.3311,0.3263,0.3311,0.3717,0.3762,0.3717,0.3668,0.3491,0.4203,0.4037,0.4149,0.4037,0.3992,0.4078,0.4651,0.4967,0.5229,0.5802,0.5802,0.5846,0.6293,0.6412,0.6374,0.6604,0.7317,0.7034,0.7573,0.7573,0.7573,0.7772,0.7605,0.8036,0.7951,0.7817,0.7869,0.7724,0.7869,0.7869,0.7951,0.7644,0.7912,0.7275,0.7342,0.7275,0.6984,0.7342,0.7605,0.7418,0.7418,0.7275,0.7573,0.7724,0.8118,0.8521,0.8823,0.8984,0.9119,0.9316,0.9512)

yy <- runmed(y, 41)
plot(y, type="l", ylim=c(0,1), ylab="", xlab="", lwd=0.5)
points(yy, col="blue", type="l", lwd=2)

Median filtering

Ответы [ 4 ]

6 голосов
/ 15 июня 2011

РЕДАКТИРОВАНИЕ: функция очищает регионы, чтобы содержать только низшую часть, если это необходимо.

На самом деле, использование среднего значения проще, чем использование медианы.Это позволяет вам находить регионы, где реальные значения постоянно ниже среднего.Медиана не является достаточно гладкой для простого применения.

Один из примеров функций для этого:

FindLowRegion <- function(x,n=length(x)/4,tol=length(x)/20,p=0.5){
    nx <- length(x)
    n <-  2*(n %/% 2) + 1
    # smooth out based on means
    sx <- rowMeans(embed(c(rep(NA,n/2),x,rep(NA,n/2)),n),na.rm=T)
    # find which series are far from the mean
    rlesx <- rle((sx-x)>0)
    # construct start and end of regions
    int <- embed(cumsum(c(1,rlesx$lengths)),2)
    # which regions fulfill requirements
    id <- rlesx$value & rlesx$length > tol
    # Cut regions to be in general smaller than median
    regions <-
    apply(int[id,],1,function(i){
        i <- min(i):max(i)
        tmp <- x[i]
        id <- which(tmp < quantile(tmp,p))
        id <- min(id):max(id)
        i[id]            
    })
    # return
    unlist(regions)
}

, где

  • n определяет, какмного значений используется для вычисления среднего значения,
  • tol определяет, сколько последовательных значений должно быть ниже, чем среднее значение, чтобы говорить о низкой области, а
  • p определяетотсечка используется (как квантиль) для зачистки областей до их нижней части.Когда p = 1, отображается полная нижняя область.

Функция настроена для работы с данными, которые вы представили, но цифры могут потребоваться немного скорректировать для работы с другими данными.

Эта функция возвращает набор индексов, который позволяет находить нижние регионы.Проиллюстрировано вашим вектором y:

Lows <- FindLowRegion(y)

newx <- seq_along(y)
newy <- ifelse(newx %in% Lows,y,NA)
plot(y, col="blue", type="l", lwd=2)
lines(newx,newy,col="red",lwd="3")

Дает:

enter image description here

3 голосов
/ 15 июня 2011

Вы должны каким-то образом сгладить график.Median filtration весьма полезен для этой цели (см. http://en.wikipedia.org/wiki/Median_filter).. После сглаживания вам просто придется искать минимумы, как обычно (т. Е. Искать точки, где первая производная переключается с отрицательной на положительную).

1 голос
/ 15 июня 2011

Более простой ответ (который также не требует сглаживания) может быть обеспечен путем адаптации функции maxdrawdown() из tseries . просадка обычно определяется как отступление от самого последнего максимума; здесь мы хотим обратного. Такая функция может затем использоваться в скользящем окне над данными или над сегментированными данными.

maxdrawdown <- function(x) {
    if(NCOL(x) > 1)
        stop("x is not a vector or univariate time series")
    if(any(is.na(x)))
        stop("NAs in x")
    cmaxx <- cummax(x)-x
    mdd <- max(cmaxx)
    to <- which(mdd == cmaxx)
    from <- double(NROW(to))
    for (i in 1:NROW(to))
        from[i] <- max(which(cmaxx[1:to[i]] == 0))
    return(list(maxdrawdown = mdd, from = from, to = to))
}

Таким образом, вместо использования cummax(), придется переключиться на cummin() и т. Д.

0 голосов
/ 15 июня 2011

Моя первая мысль была чем-то более грубым, чем фильтрация. Почему бы не искать большие капли, за которыми следуют достаточно длительные стабильные периоды?

span.b <- 20
threshold.b <- 0.2
dy.b <- c(rep(NA, span.b), diff(y, lag = span.b))
span.f <- 10
threshold.f <- 0.05
dy.f <- c(diff(y, lag = span.f), rep(NA, span.f))
down <- which(dy.b < -1 * threshold.b & abs(dy.f) < threshold.f)
abline(v = down)

График показывает, что он не идеален, но не отбрасывает выбросы (я полагаю, это зависит от вашего восприятия данных).

...