Документация

Дискретный вейвлет и дискретный пакет вейвлета преобразовывают

Следующие данные показывают дерево DWT для 1D входного сигнала.

Рисунок 1: Дерево DWT вниз к уровню 3 для 1D входного сигнала.

$$\underset{}{\overset{\sim }{G}}(f)$$ масштабирование (lowpass) аналитический фильтр и $$\underset{}{\overset{\sim }{H}}(f)$$ представляет вейвлет (highpass) аналитический фильтр. Метки в нижней части показывают раздел оси частоты [0,1/2] в поддиапазоны.

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

Сравните дерево DWT с полным пакетным деревом вейвлета.

Рисунок 2: Полное пакетное дерево вейвлета вниз к уровню 3.

В вейвлете пакет преобразовывает, операции фильтрации также применяются к вейвлету, или детали, коэффициентам. Результат состоит в том, что пакеты вейвлета обеспечивают фильтрацию поддиапазона входного сигнала на прогрессивно более прекрасные интервалы равной ширины. На каждом уровне, $$j$$, ось частоты [0,1/2] разделена на $$2^j$$ поддиапазоны. Поддиапазоны в герц на уровне j приблизительно

$$[\frac{nFs}{2^{j+1}},\frac{(n+1)Fs}{2^{j+1}})n=0,1,\dots 2^j-1$$ где Фс является частотой дискретизации.

То, насколько хороший это полосовое приближение, зависит от того, насколько локализованный частотой фильтры. Для фильтров Фейера-Коровкина как 'fk18' (18 коэффициентов), приближение довольно хорошо. Для фильтра как Хаар ('Хаар') приближение не точно.

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

Анализ частоты времени с пакетами вейвлета

Поскольку пакеты вейвлета делят ось частоты на более прекрасные интервалы, чем DWT, пакеты вейвлета выше при анализе частоты времени. Как пример, рассмотрите две неустойчивых синусоиды с частотами 150 и 200 Гц в аддитивном шуме. Данные выбираются на уровне 1 кГц. Чтобы предотвратить потерю разрешения времени, свойственного от критически выбранного пакета вейвлета, преобразовывают, используют неподкошенный, преобразовывают.

dt = 0.001;
t = 0:dt:1-dt;
x = ...
cos(2*pi*150*t).*(t>=0.2 & t<0.4)+sin(2*pi*200*t).*(t>0.6 & t<0.9);
y = x+0.05*randn(size(t));
[wpt,~,F] = modwpt(x,'TimeAlign',true);
contour(t,F.*(1/dt),abs(wpt).^2)
grid on
xlabel('Time (secs)')
ylabel('Hz')
title('Time-Frequency Analysis -- Undecimated Wavelet Packet Transform')

Обратите внимание на то, что пакет вейвлета преобразовывает, может разделить компоненты на 150 и 200 Гц. Это не верно для DWT, потому что падение на 150 и 200 Гц той же полосы октавы. Полосы октавы для уровня четыре DWT (в Гц)

Продемонстрируйте, что эти два компонента локализуются временем DWT, но не разделяются в частоте.

mra  = modwtmra(modwt(y,'fk18',4),'fk18');
J = 5:-1:1;
freq = 1/2*(1000./2.^J);
contour(t,freq,flipud(abs(mra).^2))
grid on
ylim([0 500])
xlabel('Time (secs)')
ylabel('Hz')
title('Time-Frequency Analysis -- Undecimated Wavelet Transform')

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

Энергетическое сохранение в пакете вейвлета преобразовывает

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

load kobe
plot(kobe)
grid on
xlabel('Seconds')
title('Kobe Earthquake Data')

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

st = dwtmode('status','nodisplay');
dwtmode('per','nodisplay');
wptreeDecimated = wpdec(kobe,3,'fk18');
wptUndecimated = modwpt(kobe,3);

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

wpcfs = zeros(8,381);
for kk = 7:14
    wpcfs(kk-6,:) = wpcoef(wptreeDecimated,kk);
end

Теперь вычислите полную энергию и на подкошенном и на неподкошенном пакетном уровне вейвлета три коэффициента и сравните с энергией в исходном сигнале.

decimatedEnergy = sum(sum(abs(wpcfs).^2,2))

decimatedEnergy = 2.0189e+11

undecimatedEnergy = sum(sum(abs(wptUndecimated).^2,2))

undecimatedEnergy = 2.0189e+11

signalEnergy = norm(kobe,2)^2

signalEnergy = 2.0189e+11

Доли DWT это важное свойство с пакетом вейвлета преобразовывают.

wt = modwt(kobe,'fk18',3);
undecimatedWTEnergy = sum(sum(abs(wt).^2,2))

undecimatedWTEnergy = 2.0189e+11

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

Классификация пакетов вейвлета - землетрясение или взрыв?

Сейсмические записи поднимают действие со многих источников. Важно смочь классифицировать это действие на основе своего источника. Например, землетрясения и взрывы могут представить подобные подписи временного интервала, но очень важно дифференцироваться между этими двумя событиями. Набор данных EarthQuakeExplosionData содержит 16 записей с 2 048 выборками. Первые 8 записей (столбцы) являются сейсмическими записями землетрясений, последние 8 записей (столбцы) являются сейсмическими записями взрывов. Загрузите данные и постройте землетрясение и запись взрыва для сравнения. Данные взяты из пакета R astsa Стоффер (2014), который добавляет Шумвея и Стоффера (2011). Автор любезно позволил его использование в этом примере.

Постройте временные ряды от каждой группы для сравнения.

load EarthQuakeExplosionData
subplot(2,1,1)
plot(EarthQuakeExplosionData(:,3))
xlabel('Time')
title('Earthquake Recording')
grid on
subplot(2,1,2)
plot(EarthQuakeExplosionData(:,9))
xlabel('Time')
title('Explosion Recording')
grid on

Сформируйте пакетный характеристический вектор вейвлета путем разложения каждых временных рядов вниз, чтобы выровняться, три использования 'fk6' вейвлета с неподкошенным пакетом вейвлета преобразовывают. Это приводит к 8 поддиапазонам с аппроксимированной шириной 1/16 циклов/выборки. Используйте относительную энергию в каждом поддиапазоне, чтобы создать характеристический вектор. Как пример, получите пакетный энергетический характеристический вектор родственника вейвлета для первой записи землетрясения.

[wpt,~,F,E,RE] = modwpt(EarthQuakeExplosionData(:,1),3,'fk6');

RE содержит относительную энергию в каждом из этих 8 поддиапазонов. Если вы суммируете все элементы в RE, это равно 1. Обратите внимание на то, что это - значительное сокращение данных. Временные ряды длины 2048 уменьшаются до вектора с 8 элементами, каждым элементом, представляющим относительную энергию в пакетных узлах вейвлета на уровне 3. Функция помощника helperEarthQuakeExplosionClassifier получает относительные энергии для каждой из этих 16 записей на уровне 3 с помощью 'fk6' вейвлета. Это приводит к 16 8 матрица и использует эти характеристические векторы в качестве входных параметров к kmeans классификатору. Классификатор использует критерий статистической величины Разрыва, чтобы определить оптимальное количество кластеров для характеристических векторов между 1 и 6 и классифицирует каждый вектор. Кроме того, классификатор выполняет ту же самую классификацию на неподкошенном вейвлете, преобразовывают коэффициенты на уровне 3, полученном с 'fk6' вейвлетом и спектрами мощности для каждых из временных рядов. Неподкошенный вейвлет преобразовывает результаты в 16 4 матрица (3 поддиапазона вейвлета и 1 масштабирующийся поддиапазон). Спектры мощности приводят к 16 1025 матрица. У вас должны быть Statistics and Machine Learning Toolbox™ и Signal Processing Toolbox™, чтобы запустить классификатор.

Level = 3;
[WavPacketCluster,WtCluster,PspecCluster] = ...
helperEarthQuakeExplosionClassifier(EarthQuakeExplosionData,Level)

WavPacketCluster = 
  GapEvaluation with properties:

    NumObservations: 16
         InspectedK: [1 2 3 4 5 6]
    CriterionValues: [-0.0039 0.0779 0.1314 0.0862 0.0296 -0.0096]
           OptimalK: 2

WtCluster = 
  GapEvaluation with properties:

    NumObservations: 16
         InspectedK: [1 2 3 4 5 6]
    CriterionValues: [-0.0095 0.0474 0.0706 0.0260 -0.0280 -0.0346]
           OptimalK: 1

PspecCluster = 
  GapEvaluation with properties:

    NumObservations: 16
         InspectedK: [1 2 3 4 5 6]
    CriterionValues: [0.3804 0.3574 0.3466 0.2879 0.3256 0.3326]
           OptimalK: 1

dwtmode(st,'nodisplay')

WavPacketCluster является кластеризацией вывод для неподкошенных пакетных характеристических векторов вейвлета. WtCluster является кластеризацией вывод с помощью неподкошенных характеристических векторов DWT, и PspecCluster является вывод для кластеризации на основе спектров мощности. Классификация пакетов вейвлета идентифицировала два кластера (OptimalK: 2) как оптимальный номер. Исследуйте результаты пакетной кластеризации вейвлета.

res = [ WavPacketCluster.OptimalY(1:8)' ; WavPacketCluster.OptimalY(9:end)']

res = 2×8

     1     2     2     2     2     2     2     2
     1     1     1     1     1     1     1     2

Вы видите, что были совершены только две ошибки. Из восьми записей землетрясения, семь классифицируются вместе (группа 2), и каждый по ошибке классифицируется как принадлежащий группе 1. То же самое верно для 8 записей взрыва - 7, классифицируются вместе, и 1 неправильно классифицируется. Классификация на основе уровня три неподкошенных DWT и спектры мощности возвращает один кластер как оптимальное решение.

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

Заключения

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

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

Ссылки

Wickerhauser, Т-х "Адаптировал анализ вейвлета от теории до программного обеспечения". Нажатие IEEE, 1994.

Персиваль, Д.Б. и Уолден, A.T. "Методы вейвлета для анализа временных рядов". Издательство Кембриджского университета, 2000.

Shumway, Р.Х. и Стоффер, D.H. "Анализ временных рядов и его приложения с примерами R". Спрингер, 2011.

Стоффер, D.H. "astsa: Прикладной Статистический Анализ Временных рядов". R версия пакета 1.3. https://CRAN.R-project.org/package=astsa, 2014.

Приложение

Следующая функция помощника используется в этом примере: