Без усиления шума!

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

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

Хорошо известная проблема для зашумленного сигнала заключается в том, что:

  • Дифференциация усиливает шум.
  • Интеграция вносит дрейф.

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

Недавно я наткнулся на Total Variation Regularization (TVR) — возможно, вы уже сталкивались с этим раньше в области шумоподавления изображений, похожее на описание, которое вы найдете в Википедии. Здесь я покажу вам, как это применимо к различению зашумленных сигналов по этой статье и по этой книге (глава 8). Я расскажу вам об основной идее, а затем перейду к некоторым примерам и коду на Python и Mathematica.

Вы можете найти полный код этого проекта на GitHub здесь.

Полная вариационная регуляризация (TVR)

Рассмотрим некоторый вектор данных f, который содержит зашумленные данные f_i для i=1,...,N-1, то есть длиной N. Здесь мы рассматриваем случай, когда данные равномерно распределены с интервалом \Delta x — «мы оставляем обобщение данных с неравномерным интервалом в качестве упражнения для читателя» (вы когда-нибудь читали более страшные слова?). Здесь мы будем делать все в дискретной настройке; конечно, вы также можете сделать это в непрерывном пространстве.

Взятие производной наивно определяется как:

Их всего N-1. Вы можете рассматривать g_i как производную в средней точке, т.е. если f_i соответствует x_i , то g_i является центральной разностью, дающей производную в 0.5*(x_i + x_{i+1}) . Мы также можем сформулировать это как матричную задачу:

Здесь D имеет размер (N-1) x N. Антипроизводный оператор определить немного сложнее:

Важной частью является матрица A, реализующая правило трапеций; остальное для начальной точки. Обратите внимание, что здесь A имеет размер (N-1) x (N-1) .

При регуляризации полной вариации вместо этого мы формулируем производную как задачу оптимизации с целевой функцией S(u) для желаемого вектора производной u длины N-2:

Посмотрим на первый термин:

  • \alpha — это настраиваемый гиперпараметр. Большее \alpha приводит к более гладкой производной u; меньшее \alpha приведет к более точному воспроизведению данных.
  • D — это матрица, которая реализует оператор дифференцирования — решать вам, но давайте просто воспользуемся предыдущей.
  • Мы эффективно регуляризируем вторую производную, что способствует тому, чтобы результирующая производная не имела разрывов.

Во втором сроке:

  • A — это матрица, реализующая оператор антидифференцирования. Мы будем следовать примеру использования правила трапеций (см. код ниже), так как это лучше, чем просто необработанная формула Эйлера выше. Мы можем отбросить первую точку, поскольку она дана.
  • Потери L2 во втором члене (термин достоверности данных) — хороший выбор, если вы считаете, что шум в ваших данных является гауссовым, что также хорошо работает в целом, но это следует изменить, если у вас есть конкретное известное распределение шума (см. "Эта бумага").

Как это обычно бывает в задачах оптимизации, метод оптимизатора имеет решающее значение. Хотя вы можете попытать счастья с обычным градиентным спуском или методом Ньютона, победителем, похоже, станет итерация фиксированной точки диффузии с задержкой. Все три описаны в Главе 8 этой книги, поэтому мы оставим это описание и вместо этого объясним соответствующие уравнения для обновления с u_k до u_{k+1}:

Эта реализация записывает гигантские матрицы дифференцирования и интегрирования D,A . Можно разработать более эффективную версию для больших задач, когда запись гигантских матриц слишком дорога.

Если вам нужна собственная реализация, вы можете использовать эти уравнения. Кроме того, вы можете найти реализации Mathematica и Python в этом репозитории.

Простой пример: функция абсолютного значения

В качестве простого примера мы возьмем довольно удивительный пример из following this paper. Это функция абсолютного значения, которая значительно искажена гауссовским шумом:

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

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

Подробный пример: затухающий осциллятор

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

который снова будет испорчен гауссовым шумом:

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

Теперь производная выглядит более гладкой, но (1) начальные точки выглядят не очень хорошо, и (2) интегральный сигнал отсутствует в пиках. Возможно, вы считаете, что выбранный порог слишком низкий. Давайте попробуем еще раз, установив более высокий порог:

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

На самом деле кажется, что у нас есть худшее из обоих миров — шумная производная идрейф в интеграции.

Мы можем попробовать еще две стратегии:

  • Взяв производную необработанного зашумленного сигнала, а затем фильтрация нижних частот производной:

Производная выглядит лучше, но интегрированный сигнал показывает некоторый дрейф от истинного решения.

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

Все равно не хорошо! Шум в конце производной, неточность в начале и дрейф при интегрировании.

ТВР снова

Давайте посмотрим, как TVR справляется с этой задачей.

Вы можете быть разочарованы — но подождите.

Не суди слишком рано! Производная выглядит немного прерывистой, но интегрированный сигнал выглядит потрясающе! Без дрейфа и гладкости. Это задумано TVR — задача оптимизации направлена ​​на минимизацию ошибки интегрированной производной. Одна из стратегий сделать производную более гладкой состоит в том, чтобы увеличить \alpha при оптимизации, но это принесет в жертву точность интегрированного сигнала. Настоящая магия начинается с раунда фильтрации нижних частот производной TVR:

Это победитель! Плавная и точная производная, а также плавный и точный интегрированный сигнал без дрейфа.

Сравнение вместе

Вы можете сравнить все подходы:

Я включил оценки ошибки в производную и сигнал, которые вычисляются из вектора x как:

Обратите внимание, что мы делим на диапазон, а не на истинный сигнал, так как это усилит ошибку, если истинный сигнал близок к нулю.

В одной категории (производная или ошибка сигнала) TVR не имеет самой низкой ошибки; однако он имеет низкую ошибку по обеим категориям. Кроме того, он не борется с начальными точками, которые очень неверны для многих подходов, и не вносит дрейф в интеграцию.

Последние мысли

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

Вы можете просмотреть код этого проекта на GitHub здесь.