Сгенерируйте 100 выборок синусоидального сигнала. Замените шестые и двадцатые выборки на скачки.

x = sin(2*pi*(0:99)/100);
x(6) = 2;
x(20) = -2;

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

[y,i,xmedian,xsigma] = hampel(x);

Постройте отфильтрованный сигнал и аннотируйте выбросы.

n = 1:length(x);
plot(n,x)
hold on
plot(n,xmedian-3*xsigma,n,xmedian+3*xsigma)
plot(find(i),x(i),'sk')
hold off
legend('Original signal','Lower limit','Upper limit','Outliers')

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

hampel(x,1)

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

rng('default')

n = 59;
x = sin(pi./[15 10]'*(1:n)+pi/3)';

spk = randi(2*n,9,1);
x(spk) = x(spk)*2;
x(randi(2*n,6,1)) = NaN;

plot(x)

Отфильтруйте сигнал с помощью hampel с настройками по умолчанию.

y = hampel(x);
plot(y)

Увеличьте длину движущегося окна и уменьшите порог, чтобы обработать выборку как выброс.

y = hampel(x,4,2);
plot(y)

Выведите рабочую медиану для каждого канала. Наложите медианы на графике сигнала.

[y,j,xmd,xsd] = hampel(x,4,2);
plot(x)
hold on
plot(xmd,'--')

Сгенерируйте многоканальный сигнал, который состоит из двух синусоид различных частот, встроенных в белый Гауссов шум модульного отклонения.

rng('default')

t = 0:60;
x = sin(pi./[10;2]*t)'+randn(numel(t),2);

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

k = 4;
nsig = 2;

[y,h] = hampel(x,k,nsig);

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

for k = 1:2
    hk = h(:,k);
    ax = subplot(2,1,k);
    plot(t,x(:,k))
    hold on
    plot(t,y(:,k))
    plot(t(hk),x(hk,k),'*')
    hold off
    ax.XTick = t(hk);
end

Сгенерируйте 100 выборок синусоидального сигнала. Замените шестые и двадцатые выборки на скачки.

n = 1:100;
x = sin(2*pi*n/100);
x(6) = 2;
x(20) = -2;

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

Постройте результат.

[y,i,xmedian,xsigma] = hampel(x);

plot(n,x)
hold on
plot(n,[1;1]*xmedian+3*[-1;1]*xsigma)
plot(find(i),x(i),'sk')
hold off
legend('Signal','Lower','Upper','Outliers')

Повторите вычисление с помощью размера окна $$2\times 10+1=21$$ и два стандартных отклонения как критерии идентификации выбросов.

sds = 2;
adj = 10;
[y,i,xmedian,xsigma] = hampel(x,adj,sds);

plot(n,x)
hold on
plot(n,[1;1]*xmedian+sds*[-1;1]*xsigma)
plot(find(i),x(i),'sk')
hold off
legend('Signal','Lower','Upper','Outliers')