Поперечная корреляция вейвлета для анализа задержки вывода

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

Поперечная корреляция Вейвлета является просто локализованной в масштабе версией обычной перекрестной корреляции между двумя сигналами. В перекрестной корреляции вы определяете подобие между двумя последовательностями, сдвигая одну относительно другой, умножая элемент сдвинутых последовательностей на элемент и суммируя результат. Для детерминированных последовательностей можно записать это как обычное скалярное произведение: <xn,yn-k>n=nxnyn-k где xn и yn являются последовательностями (сигналами), и шина обозначает комплексное сопряжение. Переменная, k, является переменной задержки и представляет сдвиг, примененный к yn. Если оба xn и yn являются реальными, комплексный сопряженный не является необходимым. Предположим, что yn является такой же последовательностью, как и xn но задерживается на L > 0 выборок, где L является целым числом. Для конкретности примитеyn=xn-10. Если вы выражаете yn с точки зрения xn выше, вы получаете <xn,xn-10-k>n=nxnxn-10-k. По неравенству Коши-Шварца, вышеуказанное максимизируется, когда k=-10. Это означает, что вы оставили сдвиг ( усовершенствование) yn на 10 выборок вы получаете максимальное значение перекрестной корреляционной последовательности. Если xn является L-отсроченной версией yn, xn=yn-L, затем последовательность перекрестной корреляции максимизируется в k=L. Показать это можно при помощи xcorr.

Создайте треугольный сигнал, состоящий из 20 выборок. Создайте шумную сдвинутую версию этого сигнала. Сдвиг пика треугольника составляет 3 выборки. Постройте график x и y последовательности.

n = 20;
x0 = 1:n/2;
x1 = (2*x0-1)/n;
x = [x1 fliplr(x1)]';
rng default;
y = [zeros(3,1);x]+0.3*randn(length(x)+3,1);

subplot(2,1,1)
stem(x,'filled')
axis([0 22 -1 2])
title('Input Sequence')

subplot(2,1,2)
stem(y,'filled')
axis([0 22 -1 2])
title('Output Sequence')

Figure contains 2 axes. Axes 1 with title Input Sequence contains an object of type stem. Axes 2 with title Output Sequence contains an object of type stem.

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

[xc,lags] = xcorr(x,y);
[~,I] = max(abs(xc));
figure
stem(lags,xc,'filled')
legend(sprintf('Maximum at lag %d',lags(I)))
title('Sample Cross-Correlation Sequence')

Figure contains an axes. The axes with title Sample Cross-Correlation Sequence contains an object of type stem. This object represents Maximum at lag -3.

Максимум встречается в лаге -3. Сигнал y является вторым входом для xcorr и это задержанная версия x. Вы должны сдвинуть y 3 выборки влево (отрицательный сдвиг) для максимизации перекрестной корреляции. Если вы сторнируете роли x и y в качестве входов в xcorr, максимальная задержка сейчас происходит в положительной задержке.

[xc,lags] = xcorr(y,x);
[~,I] = max(abs(xc));
figure
stem(lags,xc,'filled')
legend(sprintf('Maximum at lag %d',lags(I)))
title('Sample Cross-Correlation Sequence')

Figure contains an axes. The axes with title Sample Cross-Correlation Sequence contains an object of type stem. This object represents Maximum at lag 3.

x является расширенной версией y и вы задерживаете x тремя выборками для максимизации перекрестной корреляции.

modwtxcorr - масштабируемая версия xcorr. Как использовать modwtxcorr, вы сначала получаете неразрешенные вейвлет.

Примените вейвлет перекрестную корреляцию к двум сигналам, которые являются сдвинутыми версиями друг друга. Создайте два экспоненциально-демпфированных 200-Hz синусоиды с аддитивным шумом. The x сигнал имеет свой временной центр в t=0.2 секунд в то время как y центрировано в t=0.5 секунд.

t = 0:1/2000:1-1/2000;
x = sin(2*pi*200*t).*exp(-50*pi*(t-0.2).^2)+0.1*randn(size(t));
y = sin(2*pi*200*t).*exp(-50*pi*(t-0.5).^2)+0.1*randn(size(t));
figure
plot(t,x)
hold on
plot(t,y)
xlabel('Seconds')
ylabel('Amplitude')
grid on
legend('x','y')

Figure contains an axes. The axes contains 2 objects of type line. These objects represent x, y.

Видишь, что x и y очень похожи, за исключением того, что y задерживается на 0,3 секунды. Получите неразрешенное вейвлет x и y вплоть до уровня 5 с помощью вейвлета Фейера-Коровкина (14). Вейвлет на уровне 3 с частотой дискретизации 2 кГц являются приблизительными [2000/24,2000/23) полосно-пропускающая фильтрация входов. Частотная локализация фильтров Фейера-Коровкина гарантирует, что это полосное приближение довольно хорошо.

wx = modwt(x,'fk14',5);
wy = modwt(y,'fk14',5);

Получите вейвлет последовательности для вейвлет x и y. Постройте график вейвлет корреляционной последовательности уровня 3 для 2000 лагов с нулевой задержкой. Умножьте лаги на период дискретизации, чтобы получить значимую ось времени.

[xc,~,lags] = modwtxcorr(wx,wy,'fk14');
lev = 3;
zerolag = floor(numel(xc{lev})/2+1);
tlag = lags{lev}(zerolag-999:zerolag+1000).*(1/2000);
figure
plot(tlag,xc{lev}(zerolag-999:zerolag+1000))
title('Wavelet Cross-Correlation Sequence (level 3)')
xlabel('Time')
ylabel('Cross-Correlation Coefficient')

Figure contains an axes. The axes with title Wavelet Cross-Correlation Sequence (level 3) contains an object of type line.

Последовательность перекрестной корреляции достигает пика с задержкой -0,3 секунды. Вейвлет y является вторым входом для modwtxcorr. Потому что второй вход modwtxcorr сдвинут относительно первого, пиковая корреляция происходит с отрицательной задержкой. Для выравнивания временных рядов необходимо сдвинуть (продвинуть) последовательность перекрестной корреляции влево. Если вы сторнируете роли входов на modwtxcorr, вы получаете пиковую корреляцию при положительной задержке.

[xc,~,lags] = modwtxcorr(wy,wx,'fk14');
lev = 3;
zerolag = floor(numel(xc{lev})/2+1);
tlag = lags{lev}(zerolag-999:zerolag+1000).*(1/2000);
figure
plot(tlag,xc{lev}(zerolag-999:zerolag+1000))
title('Wavelet Cross-Correlation Sequence (level 3)')
xlabel('Time')
ylabel('Cross-Correlation Coefficient')

Figure contains an axes. The axes with title Wavelet Cross-Correlation Sequence (level 3) contains an object of type line.

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

lev = 1;
zerolag = floor(numel(xc{lev})/2+1);
tlag = lags{lev}(zerolag-999:zerolag+1000).*(1/2000);
plot(tlag,xc{lev}(zerolag-999:zerolag+1000))
title('Wavelet Cross-Correlation Sequence (level 1)')
xlabel('Time')
ylabel('Cross-Correlation Coefficient')

Figure contains an axes. The axes with title Wavelet Cross-Correlation Sequence (level 1) contains an object of type line.

figure
lev = 5;
zerolag = floor(numel(xc{lev})/2+1);
tlag = lags{lev}(zerolag-999:zerolag+1000).*(1/2000);
plot(tlag,xc{lev}(zerolag-999:zerolag+1000))
title('Wavelet Cross-Correlation Sequence (level 5)')
xlabel('Time')
ylabel('Cross-Correlation Coefficient')

Figure contains an axes. The axes with title Wavelet Cross-Correlation Sequence (level 5) contains an object of type line.

Вейвлет последовательности на уровнях 1 и 5 не показывают никаких доказательств экспоненциально-взвешенных синусоидов из-за характера полосы пропускания вейвлет.

С финансовыми данными часто существует ведущая или отстающая связь между переменными. В этих случаях полезно изучить последовательность перекрестной корреляции, чтобы определить, максимизирует ли отставание одной переменной по отношению к другой их перекрестную корреляцию. Для того чтобы проиллюстрировать это, рассмотрим взаимосвязь между двумя компонентами ВВП -- личными расходами на потребление и валовыми частными внутренними инвестициями. Данные представляют собой ежеквартальные взвешенные по цепочке данные по реальному ВВП США за 1974Q1- 2012Q4 годы. Данные были преобразованы путем сначала взятия естественного логарифма, а затем вычисления различия за год. Посмотрите на корреляцию между двумя компонентами ВВП - личные расходы на потребление, pcи валовые частные внутренние инвестиции, privateinvest.

load GDPcomponents
piwt = modwt(privateinvest,'fk8',5);
pcwt = modwt(pc,'fk8',5);
figure
modwtcorr(piwt,pcwt,'fk8')

Figure contains an axes. The axes with title Correlation by Scale -- Wavelet Coefficients contains 2 objects of type errorbar, line.

Личные расходы и личные инвестиции отрицательно коррелируют в течение 2-4 кварталов. В более длительных шкалах наблюдается сильная положительная корреляция между личными расходами и личными инвестициями. Исследуйте вейвлет последовательность в шкале, представляющей 2-4 квартальных цикла. Постройте график последовательности перекрестной корреляции вместе с 95% доверительными интервалами.

[xcseq,xcseqci,lags] = modwtxcorr(piwt,pcwt,'fk8');
zerolag = floor(numel(xcseq{1})/2)+1;
figure
plot(lags{1}(zerolag:zerolag+20),xcseq{1}(zerolag:zerolag+20))
hold on
plot(lags{1}(zerolag:zerolag+20),xcseqci{1}(zerolag:zerolag+20,:),'r--')
xlabel('Lag (Quarters)')
ylabel('Cross-Correlation')
grid on
title({'Wavelet Cross-Correlation Sequence -- [2Q,4Q)'; ...
    'Personal Consumption and Private Investment'})

Figure contains an axes. The axes with title Wavelet Cross-Correlation Sequence -- [2Q,4Q) Personal Consumption and Private Investment contains 3 objects of type line.

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