Устранение отклонений с помощью идентификатора Hampel

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

Генерировать случайный сигнал, x, содержащий 24 образца. Сбросьте генератор случайных чисел для воспроизводимых результатов.

rng default

lx = 24;
x = randn(1,lx);

Создание окна наблюдения вокруг каждого элемента x. Бери k = 2 соседей по обе стороны образца. Результирующее подвижное окно имеет длину $$2\mathrm{\times }2\mathrm{+}$$ 1 = 5 выборок.

k = 2;

iLo = (1:lx)-k;
iHi = (1:lx)+k;

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

iLo(iLo<1) = 1;
iHi(iHi>lx) = lx;

Запишите медиану каждого окружающего окна. Найдите медиану абсолютного отклонения каждого элемента относительно медианы окна.

for j = 1:lx
    w = x(iLo(j):iHi(j));
    medj = median(w);
    mmed(j) = medj;
    mmad(j) = median(abs(w-medj));
end

Масштабировать среднее абсолютное отклонение с помощью

для получения оценки стандартного отклонения нормального распределения.

sd = mmad/(erfinv(1/2)*sqrt(2));

Найдите образцы, которые отличаются от медианы более чем на nd = 2 стандартных отклонения. Замените каждое из этих отклонений значением медианы окружающего окна. В этом суть алгоритма Хэмпела.

nd = 2;
ki = abs(x-mmed) > nd*sd;

yu = x;
yu(ki) = mmed(ki);

Используйте hampel для вычисления отфильтрованного сигнала и аннотирования отклонений. Наложение отфильтрованных значений, вычисленных в этом примере.

hampel(x,k,nd)

hold on
plot(yu,'o','HandleVisibility','off')
hold off

См. также

hampel