Используйте преобразование Фурье для очистки данных временных рядов в кратчайшем коде Python

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

Преобразуйте волну временной области в частотную:

Изображение ниже хорошо иллюстрирует преобразование Фурье: разложите сложную волну на множество регулярных синусоид.

Вот полная анимация, объясняющая, что происходит при преобразовании волновых данных во временной области в представление частотной области.

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

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

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

Смешайте чистые данные с шумом

Создайте две синусоидальные волны и объедините их в одну синусоидальную волну, а затем целенаправленно загрязните чистую волну данными, сгенерированными из np.random.randn(len(t)).

import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['figure.figsize'] = [16,10]
plt.rcParams.update({'font.size':18})
#Create a simple signal with two frequencies
data_step   = 0.001
t           = np.arange(start=0,stop=1,step=data_step)
f_clean     = np.sin(2*np.pi*50*t) + np.sin(2*np.pi*120*t)
f_noise     = f_clean + 2.5*np.random.randn(len(t))
plt.plot(t,f_noise,color='c',Linewidth=1.5,label='Noisy')
plt.plot(t,f_clean,color='k',Linewidth=2,label='Clean')
plt.legend()

(Объединение двух сигналов для формирования третьего сигнала также называется сверткой или сверткой сигналов. Это еще одна интересная тема, которая в настоящее время интенсивно используется в моделях нейронных сетей)

У нас есть волны с шумом, черный - это волна, которую мы хотим, зеленые линии - это шум.

Если я скрою цвета на диаграмме, мы едва сможем отделить шум от чистых данных. Здесь может помочь преобразование Фурье, все, что нам нужно сделать, это преобразовать данные в другую перспективу, от представления времени (ось x) до представления частоты (ось x будет частотами волн).

Преобразование из временной области в частотную

Вы можете использовать numpy.fft или scipy.fft. Я обнаружил, что scipy.fft довольно удобен и полностью функционален. Здесь я буду использовать scipy.fft в этой статье, но это ваш выбор, если вы хотите использовать другие модули, или вы даже можете создать свой собственный (см. Код ниже) на основе формулы, которую я представил в начале.

from scipy.fft import rfft,rfftfreq
n    = len(t)
yf   = rfft(f_noise)
xf   = rfftfreq(n,data_step)
plt.plot(xf,np.abs(yf))
  • В коде я использую rfft вместо fft. r означает уменьшение (я думаю), чтобы вычислялись только положительные частоты. Все отрицательные зеркальные частоты будут опущены. и это тоже быстрее. см. больше обсуждения здесь.
  • Результат yf функции rfft - это комплексное число, по форме похожее на a+bj. функция np.abs() вычислит √ (a² + b²) для комплексных чисел.

А вот и волшебный вид наших исходных волн в частотной области. ось абсцисс представляет частоты.

То, что выглядит сложным во временной области, теперь преобразуется в очень простые данные частотной области. Два пика представляют частоту наших двух синусоидальных волн. Одна волна 50 Гц, другая 120 Гц. Еще раз взгляните на код, который генерирует синусоидальные волны.

f_clean     = np.sin(2*np.pi*50*t) + np.sin(2*np.pi*120*t)

Остальные частоты - это шумы, которые будут легко устранены на следующем шаге.

Убрать шумовые частоты

С помощью Numpy мы можем легко установить эти частоты как 0, кроме 50 Гц и 120 Гц.

yf_abs      = np.abs(yf) 
indices     = yf_abs>300   # filter out those value under 300
yf_clean    = indices * yf # noise frequency will be set to 0
plt.plot(xf,np.abs(yf_clean))

Теперь все шумы удалены.

Обратный возврат к данным во временной области

Код:

from scipy.fft import irfft
new_f_clean = irfft(yf_clean)
plt.plot(t,new_f_clean)
plt.ylim(-6,8)

Как показывает результат, все шумовые волны удалены.

Как работает трансформация

Вернемся к уравнению преобразования:

Исходный сигнал во временной области представлен строчными буквами x. x [n] означает точку данных во временной области в позиции n (время). Точка данных в частотной области представлена ​​заглавными буквами X. X [k] означает значение на частоте k.

Скажем, у нас есть 10 точек данных.

x = np.random.random(10)

N должно быть 10, поэтому диапазон n составляет от 0 до 9, 10 точек данных. k представляет частоту #, ее диапазон от 0 до 9, почему? в крайнем случае каждая точка данных представляет собой независимую синусоидальную волну.

В традиционном языке программирования потребуются два цикла for, один цикл для k, другой для n. В Python вы можете векторизовать операцию с помощью 0 явных циклов for (Expressive Python).

А встроенная поддержка комплексных чисел в Python просто потрясающая. Построим функцию преобразования Фурье.

import numpy as np
from scipy.fft import fft
def DFT_slow(x):
    x = np.asarray(x, dtype=float)# ensure the data type
    N = x.shape[0]                # get the x array length
    n = np.arange(N)              # 1d array
    k = n.reshape((N, 1))       # 2d array, 10 x 1, aka column array
    M = np.exp(-2j * np.pi * k * n / N)
    return np.dot(M, x)       # [a,b] . [c,d] = ac + bd, it is a sum
x = np.random.random(1024)
np.allclose(DFT_slow(x), fft(x))

Эта функция относительно медленная по сравнению с функцией из numpy или scipy, но достаточно хороша для понимания того, как работает функция БПФ. Для более быстрой реализации прочтите Понимание алгоритма БПФ Джейка.

Дальнейшие мысли - не резюме

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

Вы можете преобразовывать не только звуковые данные, но и изображения, видео, электромагнитные волны или даже данные биржевой торговли (волна Кондратьева).

Преобразование Фурье также можно интерпретировать с помощью кругов генерации волн.

Большой круг - это наша страна, наша эпоха. Мы как личность - это маленький крошечный внутренний круг. Без большого круга, который движет всем, мы действительно не можем сделать слишком много.

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

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

Чем больше я знаю о преобразовании Фурье, тем больше меня поражает Жозеф Фурье, который придумал это невероятное уравнение в 1822 году. Он никогда не мог знать, что его работа теперь используется повсюду в 21 веке.

Приложение - Четыре вида преобразования Фурье

Все преобразования Фурье, упомянутые в этой статье, относятся к дискретному преобразованию Фурье.

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

  1. Частотная область и преобразования Фурье. Учебник Пола У. Каффа из Принстона.

2. Глава 8 книги Цифровая обработка сигналов - Стивен Смит предоставил онлайн-ссылку: DSP Ch8

использованная литература

  1. Интерактивное руководство по преобразованию Фурье: https://betterexplained.com/articles/an-interactive-guide-to-the-fourier-transform/
  2. Преобразование Фурье с помощью scipy.fft: Обработка сигналов Python: https://realpython.com/python-scipy-fft/
  3. Понимание алгоритма БПФ: http://jakevdp.github.io/blog/2013/08/28/understanding-the-fft/
  4. Частотная область и преобразования Фурье: https://www.princeton.edu/~cuff/ele201/kulkarni_text/frequency.pdf
  5. Удаление шумов с помощью БПФ [Python]: https://www.youtube.com/watch?v=s2K1JfNR7Sc&ab_channel=SteveBrunton
  6. Алгоритм БПФ - простой шаг за шагом: https://www.youtube.com/watch?v=htCj9exbGo0&ab_channel=SimonXu

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