Что такое анализ мультиразрешения?

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

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

Fs = 1e3;
t = 0:1/Fs:1-1/Fs;
comp1 = cos(2*pi*200*t).*(t>0.7);
comp2 = cos(2*pi*60*t).*(t>=0.1 & t<0.3);
trend = sin(2*pi*1/2*t);
rng default
wgnNoise = 0.4*randn(size(t));
x = comp1+comp2+trend+wgnNoise;
plot(t,x)
xlabel('Seconds')
ylabel('Amplitude')

Сигнал явным образом состоит из трех основных компонентов: локализованное временем колебание с частотой 60 циклов/секунда, локализованное временем колебание с частотой 200 циклов/секунда и термин тренда. Термин тренда здесь является также синусоидальным, но имеет частоту 1/2 цикла в секунду, таким образом, это завершает только 1/2 цикл в одном втором интервале. Эти 60 циклов/секунда или колебание на 60 Гц находятся между 0,1 и 0,3 секундами, в то время как колебание на 200 Гц находится между 0,7 и 1 секундой.

Не все это сразу очевидно из графика необработанных данных, потому что эти компоненты смешаны.

Теперь постройте сигнал с точки зрения частоты.

xdft = fft(x);
N = numel(x);
xdft = xdft(1:numel(xdft)/2+1);
freq = 0:Fs/N:Fs/2;
plot(freq,20*log10(abs(xdft)))
xlabel('Cycles/second'); ylabel('dB')
grid on

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

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

cwt(x,Fs)

Теперь вы видите степени времени компонентов на 200 Гц и на 60 Гц. Однако у нас все еще нет полезной визуализации тренда.

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

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

Разделение компонентов сигнала вовремя

Реальные сигналы являются смесью различных компонентов. Часто вы только интересуетесь подмножеством этих компонентов. Анализ мультиразрешения позволяет вам сужать свой анализ путем разделения сигнала на компоненты в различных разрешениях.

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

Анализируйте и постройте синтетический сигнал с помощью вейвлета MRA. Сигнал анализируется в восьми разрешениях или уровнях.

mra = modwtmra(modwt(x,8));
helperMRAPlot(x,mra,t,'wavelet','Wavelet MRA',[2 3 4 9])

Не объясняя, что обозначения на среднем значении графика, давайте использовать наше знание сигнала и попытаемся изучить то, что этот вейвлет MRA показывает нам. Если при запуске с верхнего графика и продолжаете вниз, пока вы не достигаете графика исходных данных, вы видите, что компонентами прогрессивно становятся более сглаженными. Если вы предпочитаете думать о данных в терминах частоты, частоты, содержавшиеся в компонентах, становятся ниже. Вспомните, что исходный сигнал имел три основных компонента, высокочастотное колебание на уровне 200 Гц, более низкое колебание частоты 60 Гц и термин тренда, все поврежденные аддитивным шумом.

Если вы смотрите на $$\underset{}{\overset{\sim }{D}}2$$ постройте, вы видите, что локализованный временем высокочастотный компонент изолируется там. Вы видите и исследуете эту важную функцию сигнала по существу в изоляции. Следующие два графика содержат более низкое колебание частоты. Это - важный аспект анализа мультиразрешения, а именно, важные компоненты сигнала не могут закончиться изолированные в одном компоненте MRA, но они редко располагаются в больше чем двух. Наконец, мы видим $$S8$$ график содержит термин тренда. Для удобства цвет осей в тех компонентах был изменен, чтобы подсветить их в MRA. Если вы предпочитаете визуализировать этот график или последующие единицы без выделения, не учтите последний числовой вход к helperMRAPlot.

Вейвлет MRA использует стандартные функции, вызвал вейвлеты, чтобы разделить компоненты сигнала. K-ой вейвлет компонент MRA, обозначенный $$\underset{}{\overset{\sim }{D}}k$$ в предыдущем графике, может рассматриваться как фильтрация сигнала в диапазоны частот формы $[\frac{1}{2^{k+1}\Delta t},\frac{1}{2^k\Delta t}]$ где $\Delta t$ период выборки или интервал выборки. Финал сглаживает компонент, обозначенный в графике $$SL$$, где $$L$$ уровень MRA, получает диапазон частот $[0,\frac{1}{2^{L+1}\Delta t}]$. Точность этого приближения зависит от вейвлета, используемого в MRA. См. [4] для подробных описаний вейвлета и пакета вейвлета MRAs.

Однако существуют другие методы MRA, чтобы рассмотреть.

Эмпирическое разложение моды (EMD) является адаптивным данными методом мультиразрешения. EMD рекурсивно извлекает различные разрешения из данных без использования стандартных функций или фильтров. EMD рассматривает сигнал как состоящий из быстрого колебания, наложенного на более медленное. После того, как быстрое колебание извлечено, процесс обрабатывает остающийся более медленный компонент как новый сигнал и снова рассматривает его как быстрое колебание, наложенное на более медленное. Процесс продолжается, пока некоторый останавливающийся критерий не достигнут. В то время как EMD не использует стандартные функции как вейвлеты, чтобы извлечь информацию, подход EMD концептуально очень похож на метод вейвлета разделения сигнала в детали и приближения и затем разделение приближения снова в детали и приближение. Компоненты MRA в EMD упоминаются как внутренние функции режима (IMF). См. [3] для подробной обработки EMD.

Постройте анализ EMD того же сигнала.

[imf_emd,resid_emd] = emd(x);
helperMRAPlot(x,imf_emd,t,'emd','Empirical Mode Decomposition',[1 2 3 6])

В то время как количество компонентов MRA отличается, EMD и вейвлет, MRAs производят подобное изображение сигнала. Это не случайно. См. [2] для описания общих черт между вейвлетом, преобразовывают и EMD.

В разложении EMD высокочастотное колебание локализуется к первой внутренней функции режима (МВФ 1). Более низкое колебание частоты локализуется в основном к МВФ 2, но вы видите некоторый эффект также в МВФ 3. Компонент тренда в МВФ 6 очень похож на компонент тренда, извлеченный методом вейвлета.

Еще один метод для адаптивного анализа мультиразрешения является вариационным разложением режима (VMD). Как EMD, VMD пытается извлечь внутренние функции режима или режимы колебания от сигнала, не используя стандартные функции для анализа. Но EMD и VMD определяют режимы совсем другими способами. EMD работает рекурсивно над сигналом области времени извлечь прогрессивно более низкую частоту IMFs. VMD запускается путем идентификации peaks сигнала в частотном диапазоне и извлекает все режимы одновременно. См. [1] для обработки VMD.

[imf_vmd,resid_vmd] = vmd(x);
helperMRAPlot(x,imf_vmd,t,'vmd','Variational Mode Decomposition',[2 4 5])

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

Для реального примера полезного разделения компонента рассмотрите сейсмограф (вертикальное ускорение, $$\text{nm}/{\text{s}}^2$$) из землетрясения Кобе, зарегистрированного в Университете Тасмании, Хобарт, Австралия 16 января 1995, начинаясь в 20:56:51 (GMT) и продолжаясь в течение 51 минуты в 1 втором интервале.

load KobeTimeTable
T = KobeTimeTable.t;
kobe = KobeTimeTable.kobe;
figure
plot(T,kobe)
title('Kobe Earthquake Seismograph')
axis tight
grid on

Получите и постройте вейвлет MRA данных.

mraKobe = modwtmra(modwt(kobe,8));
figure
helperMRAPlot(kobe,mraKobe,T,'Wavelet','Wavelet MRA Kobe Earthquake',[4 5])

График показывает разделение первичных и задержал вторичные волновые компоненты в компонентах MRA $\underset{}{\overset{\sim }{D}}4$ и $\underset{}{\overset{\sim }{D}}5$. Компоненты в сейсмическом перемещении волны при различных скоростях с первичными волнами (сжатия), перемещающимися быстрее, чем вторичный (сдвиг) волны. Методы MRA могут позволить вам изучить эти компоненты в изоляции на исходном масштабе времени.

Восстановление сигналов от MRA

Точка разделения сигналов на компоненты должна часто удалять определенные компоненты или смягчать их последствие на сигнале. Крайне важный для методов MRA способность восстановить исходный сигнал.

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

sigrec_wavelet = sum(mra);
sigrec_emd = sum(imf_emd,2)+resid_emd;
sigrec_vmd = sum(imf_vmd,2)+resid_vmd;
figure
subplot(3,1,1)
plot(t,sigrec_wavelet); title('Wavelet reconstruction');
set(gca,'XTickLabel',[]);
ylabel('Amplitude');
subplot(3,1,2);
plot(t,sigrec_emd); title('EMD reconstruction');
set(gca,'XTickLabel',[]);
ylabel('Amplitude');
subplot(3,1,3)
plot(t,sigrec_vmd); title('VMD reconstruction');
ylabel('Amplitude');
xlabel('Time');

Максимальная ошибка реконструкции на основе выборки выборкой для каждого из методов находится на порядке $10^{-12}$ или меньший, указывая они - совершенные методы реконструкции. Поскольку сумма компонентов MRA восстанавливает исходный сигнал, она выдерживает обосновать, что включая или, исключая подмножество компонентов мог произвести полезное приближение.

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

sigWOtrend = sum(mra(1:end-1,:));
figure
plot(t,sigWOtrend)
title('Trend Term Removed')

Чтобы удалить другие компоненты, можно создать логический вектор с false значения в компонентах вы не хотите включать. Здесь мы удаляем тренд и самую высокую частотную составляющую наряду с первым компонентом MRA (который надеется быть в основном шумом). Постройте фактический второй компонент сигнала (60 Гц) наряду с реконструкцией для сравнения.

include = true(size(mra,1),1);
include([1 2 9]) = false;
ts = sum(mra(include,:));
plot(t,comp2,'b')
hold on
plot(t,ts,'r')
title('Trend Term and Highest Frequency Component Removed')
xlabel('Seconds');
legend('Component 2','Partial Reconstruction')
xlim([0.0 0.4])

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

figure
plot(t,trend,'LineWidth',2)
hold on
plot(t,[mra(end,:)' imf_vmd(:,end) imf_emd(:,end)])
grid on
legend('Trend','Wavelet','VMD','EMD')
title('Trend in three MRAs')

Обратите внимание на то, что тренд более сглажен и наиболее точно полученный методом вейвлета. EMD находит сглаженный термин тренда, но она смещена относительно истинной амплитуды тренда, в то время как метод VMD кажется по сути более смещенным к нахождению колебаний, чем методы EMD и вейвлет. Последствия этого далее обсуждены в разделе MRA Techniques — Advantages и Disadvantages.

Обнаружение переходных изменений Используя MRA

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

Чтобы проиллюстрировать это, рассмотрите ежеквартальные взвешенные цепью американские действительные данные о валовом внутреннем продукте (ВВП) для 1947Q1 к 2011Q4. Ежеквартальные выборки соответствуют частоте дискретизации 4 выборок/год. Вертикальная черная линия отмечает начало "Большого Модерирования" выражение периода уменьшенной макроэкономической энергозависимости в США, начинающихся в середине 1980-х. Обратите внимание на то, что это затрудняет, чтобы различить рассмотрение необработанных данных.

load GDPData
figure
plot(year,realgdp)
hold on
plot([year(146) year(146)],[-0.06 0.14],'k')
title('GDP Data'); xlabel('Year');
hold off

Получите вейвлет MRA данных о GDP. Постройте самое прекрасное разрешение компонент MRA с периодом большого отмеченного модерирования. Поскольку вейвлет, MRA получен с помощью зафиксированных фильтров, мы можем сопоставить самую прекрасную шкалу компонент MRA с частотами 1 цикла в год к 2 циклам в год. Компонент, который колебался с периодом двух четвертей, имеет частоту 2 циклов в год. В этом случае, самое прекрасное разрешение, компонент MRA получает изменения в GDP, которые находятся между смежными интервалами 2D четверти к изменениям, которые происходят от четверти до четверти.

mra = modwtmra(modwt(realgdp,'db2'));
figure
plot(year,mra(1,:))
hold on
plot([year(146) year(146)],[-0.015 0.015],'k')
title('Wavelet MRA - Quarter-to-Quarter Changes');
hold off

Сокращение изменчивости, или экономической энергозависимости условий, намного с большей готовностью очевидно в самом прекрасном разрешении компонент MRA, чем в необработанных данных. Методы, чтобы обнаружить изменения в отклонении как findchangepts MATLAB (Signal Processing Toolbox) часто работает лучше над компонентами MRA, чем на необработанных данных.

Методы MRA — преимущества и недостатки

В этом примере мы обсудили вейвлет и адаптивные данными методы для анализа мультиразрешения. Что такое некоторые преимущества и недостатки различных методов? Другими словами, поскольку, какие приложения я могу выбрать один по другому? Давайте запустимся с вейвлетов. Методы вейвлета в этом использовании в качестве примера зафиксировали фильтры, чтобы получить MRA. Это означает, что вейвлет, MRA имеет четко определенное математическое объяснение и мы можем предсказать поведение MRA. Мы также можем связать события в MRA к определенным масштабам времени в данных, как был сделан в примере GDP. Недостаток - то, что вейвлет преобразовывает, делит сигнал на полосы октавы (сокращение центральной частоты 1/2 в каждом компоненте) так, чтобы на высоких центральных частотах пропускная способность была намного больше, чем те на более низких частотах. Это означает, что два близко расположенных высокочастотных колебания могут легко закончиться в том же компоненте MRA с методом вейвлета.

Повторите первый синтетический пример, но переместите два колебательных компонента в одной октаве друг друга.

Fs = 1e3;
t = 0:1/Fs:1-1/Fs;
comp1 = cos(2*pi*150*t).*(t>=0.1 & t<0.3);
comp2 = cos(2*pi*200*t).*(t>0.7);
trend = sin(2*pi*1/2*t);
rng default;
wgnNoise = 0.4*randn(size(t));
x = comp1+comp2+trend+wgnNoise;
plot(t,x)
xlabel('Seconds')
ylabel('Amplitude')

Повторите и постройте вейвлет MRA.

mra = modwtmra(modwt(x,8));
helperMRAPlot(x,mra,t,'wavelet','Wavelet MRA',[2 9])

Теперь мы видим это $\underset{}{\overset{\sim }{D}}2$ содержит и компоненты на 200 Гц и на 150 Гц. Если вы повторяете это аналитическое использование EMD, вы видите тот же результат.

Теперь давайте использовать VMD.

[imf_vmd,~,info_vmd] = vmd(x);
helperMRAPlot(x,imf_vmd,t,'vmd','VMD',[1 2 3 5]);

VMD может разделить эти два компонента. Высокочастотное колебание локализовало к МВФ 1, в то время как второй компонент распространен по двум смежным IMFs. Если вы смотрите на предполагаемые центральные частоты режимов VMD, метод локализовал первые два компонента приблизительно 200 и 150 Гц. Третий МВФ имеет центральную частоту близко к 150 Гц, которая является, почему мы видим второй компонент в двух компонентах MRA.

info_vmd.CentralFrequencies*Fs

ans = 5×1

  202.7204
  153.3278
  148.8022
   84.2802
    0.2667

VMD может сделать это, потому что он запускается путем идентификации частот центра кандидата для IMFs путем рассмотрения анализа частотного диапазона данных.

В то время как вейвлет, MRA не может разделить два высокочастотных компонента, существует дополнительный пакет вейвлета MRA, который может.

wpt = modwptdetails(x,3);
helperMRAPlot(x,flip(wpt),t,'wavelet','Wavelet Packet MRA',[5 6 8]);

Теперь вы видите, что эти два колебания разделяются в $\underset{}{\overset{\sim }{D}}5$ и $\underset{}{\overset{\sim }{D}}6$. Из этого примера мы можем извлечь общее правило. Если начальное разложение вейвлета или EMD показывает компоненты с по-видимому различными уровнями колебания в том же компоненте, считайте VMD или пакет вейвлета MRA. Если вы будете подозревать, что ваши данные имеют высокочастотные компоненты близко друг к другу в частоте, VMD или пакетный подход вейвлета будут обычно работать лучше, чем подход EMD или вейвлет.

Вспомните проблему извлечения сглаженного тренда. Повторите и вейвлет MRA и EMD.

mra = modwtmra(modwt(x,8));
helperMRAPlot(x,mra,t,'wavelet','Wavelet MRA',[2 9])

imf_emd = emd(x);
helperMRAPlot(x,imf_emd,t,'EMD','Empirical Mode Decomposition',[1 2 6])

Тренды, извлеченные вейвлетом и методами EMD, ближе к истинному тренду, чем извлеченные VMD и пакетным методом вейвлета. VMD по сути смещается к нахождению узкополосных колебательных компонентов. Это - сила VMD в обнаружении близко расположенных колебаний, но недостатка при извлечении сглаженных трендов в данных. То же самое верно для пакета вейвлета MRAs, когда разложение идет вне нескольких уровней. Это приводит к второй общей рекомендации. Если вы интересуетесь охарактеризованием сглаженного тренда в ваших данных для идентификации или удаления, попробуйте метод EMD или вейвлет.

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

[imf_vmd,~,vmd_info] = vmd(realgdp);
figure
subplot(2,1,1)
plot(year,realgdp)
title('Real GDP');
hold on
plot([year(146) year(146)],[-0.1 0.15],'k')
hold off
subplot(2,1,2)
plot(year,imf_vmd(:,1));
title('First VMD IMF');
hold on
plot([year(146) year(146)],[-0.02 0.02],'k')
hold off

В то время как самая высокая частота, которую компонент VMD также, кажется, показывает некоторому сокращению изменчивости, начинающейся в середине 1980-х, является им не так с готовностью очевидный как в вейвлете MRA. Поскольку метод VMD смещается к нахождению колебаний, существует существенное звенение в первом МВФ VMD, который затеняет изменения энергозависимости.

Повторите это аналитическое использование EMD.

imf_emd = emd(realgdp);
figure
subplot(2,1,1)
plot(year,realgdp)
title('Real GDP');
hold on
plot([year(146) year(146)],[-0.1 0.15],'k')
hold off
subplot(2,1,2)
plot(year,imf_emd(:,1));
title('First EMD IMF');
hold on
plot([year(146) year(146)],[-0.06 0.05],'k')
hold off

Метод EMD менее полезен в нахождении изменения в энергозависимости (отклонение). В этом случае использование стандартных функций в вейвлете MRA было более выгодным, чем данные адаптивные методы.

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

Заключения

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

С Signal Multiresolution Analyzer можно выполнить анализ мультиразрешения сигнала, получить метрики на различных компонентах MRA, эксперименте с частичными реконструкциями, и сгенерировать скрипты MATLAB, чтобы воспроизвести анализ в командной строке.

Ссылки

[1] Dragomiretskiy, Константин и Доминик Зоссо. “Вариационное Разложение Режима”. Транзакции IEEE на Обработке сигналов 62, № 3 (февраль 2014): 531–44. https://doi.org/10.1109/TSP.2013.2288675.

[2] Фландрен, P. G. Вытекание струей, и П. Гонсалвеш. “Эмпирическое разложение моды как Набор фильтров”. Буквы Обработки сигналов IEEE 11, № 2 (февраль 2004): 112–14. https://doi.org/10.1109/LSP.2003.821662.

[3] Хуан, Норден Э., Чжен Шен, Стивен Р. Лонг, Манли К. Ву, Син Х. Ши, Куэнэн Чжен, Най-Чюань Янь, Ши Чао Туньг и Генри Х. Лю. “Эмпирическое разложение моды и Гильбертов Спектр для Нелинейного и Неустановившегося Анализа Временных рядов”. Продолжения Королевского общества Лондона. Серии А: Математические, Физические и Технические науки 454, № 1971 (8 марта 1998): 903–95. https://doi.org/10.1098/rspa.1998.0193.

[4] Персиваль, Дональд Б. и Эндрю Т. Уолден. Методы вейвлета для анализа временных рядов. Кембриджский ряд в статистической и вероятностной математике. Кембридж  ; Нью-Йорк: Издательство Кембриджского университета, 2000.