Используйте преобразование Фурье для очистки данных временных рядов в кратчайшем коде 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 - преобразование, которое обсуждается в этой статье. Если вы решили окунуться в болото математики, вот две книги, которые рекомендуются к прочтению:
- Частотная область и преобразования Фурье. Учебник Пола У. Каффа из Принстона.
2. Глава 8 книги Цифровая обработка сигналов - Стивен Смит предоставил онлайн-ссылку: DSP Ch8
использованная литература
- Интерактивное руководство по преобразованию Фурье: https://betterexplained.com/articles/an-interactive-guide-to-the-fourier-transform/
- Преобразование Фурье с помощью scipy.fft: Обработка сигналов Python: https://realpython.com/python-scipy-fft/
- Понимание алгоритма БПФ: http://jakevdp.github.io/blog/2013/08/28/understanding-the-fft/
- Частотная область и преобразования Фурье: https://www.princeton.edu/~cuff/ele201/kulkarni_text/frequency.pdf
- Удаление шумов с помощью БПФ [Python]: https://www.youtube.com/watch?v=s2K1JfNR7Sc&ab_channel=SteveBrunton
- Алгоритм БПФ - простой шаг за шагом: https://www.youtube.com/watch?v=htCj9exbGo0&ab_channel=SimonXu
Если у вас есть какие-либо вопросы, оставьте комментарий, и я постараюсь ответить. Если вы заметили ошибку или ошибку, не стесняйтесь отмечать их. Ваше чтение и поддержка - это большой круг, который побуждает меня продолжать писать больше. Спасибо.