Существует множество и множество классических методов обнаружения пиков, любой из которых может работать. Вам нужно будет увидеть, что, в частности, ограничивает качество ваших данных. Вот основные описания:

• Между любыми двумя точками ваших данных (x(0), y(0)) и (x(n), y(n)) добавьте y(i + 1) - y(i) для 0 <= i < n и вызовите это T ( "travel" ) и установите R ( "повышение" ) до y(n) - y(0) + k для подходящего малого k. T/R > 1 указывает пик. Это работает нормально, если большое перемещение из-за шума маловероятно или если шум распределяется симметрично вокруг формы базовой кривой. Для вашего приложения принимайте самый ранний пик со счетом выше заданного порога или анализируйте кривую проезда на значения высоты для более интересных свойств.

• Используйте сопоставленные фильтры для оценки сходства с стандартной формой пика (по существу, используйте нормализованный точечный продукт против некоторой формы, чтобы получить косинус-метрику сходства)

• Деконвертируйте стандартную форму пика и проверьте высокие значения (хотя я часто считаю, что 2 менее чувствителен к шуму для простого выходного сигнала).

• Разгладьте данные и проверьте триплеты одинаково разнесенных точек, где, если x0 < x1 < x2, y1 > 0.5 * (y0 + y2), или проверите евклидовы расстояния следующим образом: D((x0, y0), (x1, y1)) + D((x1, y1), (x2, y2)) > D((x0, y0),(x2, y2)), который опирается на неравенство треугольника. Использование простых коэффициентов снова предоставит вам механизм подсчета очков.

• Настройте очень простую модель из 2-гауссовой смеси для ваших данных (например, Numerical Recipes имеет хороший готовый кусок кода). Возьмите более ранний пик. Это будет иметь дело с перекрывающимися пиками.

• Найдите наилучшее совпадение данных с простой кривой Гаусса, Коши, Пуассона или того, что у вас есть. Оцените эту кривую в широком диапазоне и вычтите ее из копии данных, отметив ее пиковое местоположение. Повторение. Возьмите самый ранний пик, параметры модели которого (возможно, стандартное отклонение, но некоторые приложения могут заботиться о эксцессах или других особенностях) соответствуют некоторому критерию. Следите за отсутствием артефактов, когда пики вычитаются из данных. Лучшее совпадение может быть определено с помощью оценки соответствия, предложенной в № 2 выше.

Я сделал то, что вы делали раньше: находите пики в данных последовательности ДНК, находите пики в производных, оцененные по измеренным кривым, и находите пики на гистограммах.

Я призываю вас внимательно следить за надлежащей базой. Фильтрация Винера или другая фильтрация или простой анализ гистограммы часто являются легким способом базового уровня при наличии шума.

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

Я хотел бы внести в этот поток алгоритм, который Я разработал сам:

Он основан на принципе dispersion: если новый datapoint - это заданное x число стандартных отклонений от некоторого движущегося средние сигналы алгоритма (также называемые z-score). Алгоритм очень надежный, поскольку он создает отдельное движущееся среднее и отклонение, так что сигналы не повреждают порог. Таким образом, будущие сигналы идентифицируются с примерно одинаковой точностью, независимо от количества предыдущих сигналов. Алгоритм принимает 3 входа: lag = the lag of the moving window, threshold = the z-score at which the algorithm signals и influence = the influence (between 0 and 1) of new signals on the mean and standard deviation. Например, a lag из 5 будет использовать последние 5 наблюдений для сглаживания данных. A threshold в 3.5 будет сигнализировать, если точка передачи равна 3,5 стандартным отклонениям от движущегося среднего. А influence 0,5 дает сигналы о половине влияния, которое имеют обычные данные. Точно так же influence of 0 полностью игнорирует сигналы для пересчета нового порога: поэтому значение 0 является наиболее надежным вариантом.

Он работает следующим образом:

псевдокод

# Let y be a vector of timeseries data of at least length lag+2
# Let mean() be a function that calculates the mean
# Let std() be a function that calculates the standard deviaton
# Let absolute() be the absolute value function

# Settings (the ones below are examples: choose what is best for your data)
set lag to 5;          # lag 5 for the smoothing functions
set threshold to 3.5;  # 3.5 standard deviations for signal
set influence to 0.5;  # between 0 and 1, where 1 is normal influence, 0.5 is half

# Initialise variables
set signals to vector 0,...,0 of length of y;   # Initialise signal results
set filteredY to y(1,...,lag)                   # Initialise filtered series
set avgFilter to null;                          # Initialise average filter
set stdFilter to null;                          # Initialise std. filter
set avgFilter(lag) to mean(y(1,...,lag));       # Initialise first value
set stdFilter(lag) to std(y(1,...,lag));        # Initialise first value

for i=lag+1,...,t do
  if absolute(y(i) - avgFilter(i-1)) > threshold*stdFilter(i-1) then
    if y(i) > avgFilter(i-1)
      set signals(i) to +1;                     # Positive signal
    else
      set signals(i) to -1;                     # Negative signal
    end
    # Adjust the filters
    set filteredY(i) to influence*y(i) + (1-influence)*filteredY(i-1);
    set avgFilter(i) to mean(filteredY(i-lag,i),lag);
    set stdFilter(i) to std(filteredY(i-lag,i),lag);
  else
    set signals(i) to 0;                        # No signal
    # Adjust the filters
    set filteredY(i) to y(i);
    set avgFilter(i) to mean(filteredY(i-lag,i),lag);
    set stdFilter(i) to std(filteredY(i-lag,i),lag);
  end
end

Demo

> Оригинальный ответ

Эта проблема была изучена более подробно.

В классах TSpectrum* есть набор очень современных реализаций ROOT ( инструмент анализа физики ядра/частицы). Код работает в одно-трехмерных данных.

Исходный код ROOT доступен, поэтому вы можете захватить эту реализацию, если хотите.

Из документации класса TSpectrum:

Алгоритмы, используемые в этом классе, опубликованы в следующих ссылках:

[1] M.Morhac et al.: предыстория методы исключения для многомерное совпадение гамма-излучения спектры. Ядерные инструменты и Методы исследования физики A 401 (1997) 113- 132.

[2] M.Morhac et al.: эффективное одномерное и двумерное золото деконволюция и ее применение к гамма-спектральное разложение. Ядерные инструменты и методы в Физические исследования A 401 (1997) 385-408.

[3] M.Morhac et al.: идентификация пиков в многомерное совпадение гамма-излучения спектры. Ядерные инструменты и Методы исследования Физика A 443 (2000), 108-125.

Документы связаны с документацией классов для тех из вас, у кого нет подписки на NIM.

Короткий вариант того, что делается, состоит в том, что гистограмма сглаживается для устранения шума, а затем локальные максимумы обнаруживаются грубой силой в сплющенной гистограмме.

Этот метод в основном из книги Дэвида Марра "Видение"

Гауссов размывает ваш сигнал с ожидаемой шириной ваших пиков. это избавляет от шумов, и ваши фазовые данные не повреждены.

Тогда обнаружение края (LOG сделает)

Тогда ваши ребра были краями функций (например, пиками). посмотрите между краями для пиков, сортируйте пики по размеру, и все готово.

Я использовал вариации на этом, и они работают очень хорошо.

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

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

Я не очень разбираюсь в инструментах, поэтому это может быть совершенно непрактичным, но опять же это может быть полезным другим направлением. Если вы знаете, как показания могут потерпеть неудачу, и существует определенный интервал между пиками при таких сбоях, почему бы не сделать градиентный спуск на каждом интервале. Если спуск возвращает вас в область, которую вы искали раньше, вы можете отказаться от нее. В зависимости от формы отбираемой поверхности это также может помочь вам найти пики быстрее, чем поиск.

Вы можете применить стандартное отклонение к своей логике и обратить внимание на пики выше x%.