Остаточный анализ с автокорреляцией
Этот пример показов, как использовать автокорреляцию с доверием интервалом для анализа невязок наименьших квадратов, подгоняемых к зашумленным данным. Невязками являются различия между подобранной моделью и данными. В модели сигнала плюс белый шум, если у вас хорошая подгонка для сигнала, невязки должны быть белым шумом.
Создайте зашумленные данные набор, состоящий из полинома 1-го порядка (прямая линия) в аддитивном белом Гауссовом шуме. Аддитивный шум является последовательностью некоррелированных случайных переменных после N (0,1) распределения. Это означает, что все случайные переменные имеют среднее нулевое и единичное отклонение. Установите генератор случайных чисел в настройки по умолчанию для воспроизводимых результатов.
x = -3:0.01:3;
rng default
y = 2*x+randn(size(x));
plot(x,y)
Использование polyfit чтобы найти линию методом наименьших квадратов для зашумленных данных. Постройте график исходных данных вместе с методом наименьших квадратов.
coeffs = polyfit(x,y,1); yfit = coeffs(2)+coeffs(1)*x; plot(x,y) hold on plot(x,yfit,'linewidth',2)
Найдите невязки. Получаем автокорреляционную последовательность невязок до запаздывания 50.
residuals = y - yfit;
[xc,lags] = xcorr(residuals,50,'coeff');При просмотре автокорреляционной последовательности необходимо определить, есть ли доказательства автокорреляции. Другими словами, вы хотите определить, выглядит ли выборочная автокорреляционная последовательность как автокорреляционная последовательность белого шума. Если автокорреляционная последовательность остатков похожа на автокорреляцию процесса белого шума, вы уверены, что ни один из сигналов не избежал вашей подгонки и оказался в невязках. В этом примере используйте 99% -конфицитный интервал. Чтобы создать доверительный интервал, необходимо знать распределение выборочных значений автокорреляции. Вам также нужно найти критические значения в соответствующем распределении, между которыми лежит 0,99 вероятности. Поскольку распределение в этом случае Гауссово, можно использовать дополнительную функцию обратной ошибки erfcinv. Связь между этой функцией и обратной функцией Гауссова совокупного распределения описана на странице с описанием для erfcinv.
Найдите критическое значение для интервала доверия 99%. Используйте критическое значение, чтобы создать нижнюю и верхнюю доверительные границы.
conf99 = sqrt(2)*erfcinv(2*.01/2); lconf = -conf99/sqrt(length(x)); upconf = conf99/sqrt(length(x));
Постройте график автокорреляционной последовательности вместе с 99% -конфидентными интервалами.
figure stem(lags,xc,'filled') ylim([lconf-0.03 1.05]) hold on plot(lags,lconf*ones(size(lags)),'r','linewidth',2) plot(lags,upconf*ones(size(lags)),'r','linewidth',2) title('Sample Autocorrelation with 99% Confidence Intervals')
За исключением нулевой задержки, значения дискретизации автокорреляции находятся в пределах 99% -доверие для автокорреляции последовательности белого шума. Из этого можно сделать вывод, что невязки являются белым шумом. Более конкретно, вы не можете отвергнуть, что невязки являются реализацией процесса белого шума.
Создайте сигнал, состоящий из синусоиды плюс шум. Данные отбираются с частотой дискретизации 1 кГц. Частота синусоиды составляет 100 Гц. Установите генератор случайных чисел в настройки по умолчанию для воспроизводимых результатов.
Fs = 1000;
t = 0:1/Fs:1-1/Fs;
rng default
x = cos(2*pi*100*t)+randn(size(t));Используйте дискретное преобразование Фурье (ДПФ), чтобы получить аппроксимацию методом наименьших квадратов к синусоиде с частотой 100 Гц. Оценка амплитуды методом наименьших квадратов в 2/N умножает коэффициент ДПФ, соответствующий 100 Гц, где N - длина сигнала. Действительная часть является амплитудой косинуса при 100 Гц и мнимая часть является амплитудой синуса при 100 Гц. Аппроксимация методом наименьших квадратов является суммой косинуса и синуса с правильной амплитудой. В этом примере интервал 101 ДПФ соответствует 100 Гц.
xdft = fft(x); ampest = 2/length(x)*xdft(101); xfit = real(ampest)*cos(2*pi*100*t)+imag(ampest)*sin(2*pi*100*t); figure plot(t,x) hold on plot(t,xfit,'linewidth',2) axis([0 0.30 -4 4]) xlabel('Seconds') ylabel('Amplitude')
Найдите невязки и определите последовательность автокорреляции выборки до задержки 50.
residuals = x-xfit;
[xc,lags] = xcorr(residuals,50,'coeff');Постройте график автокорреляционной последовательности с 99% -конфидентными интервалами.
figure stem(lags,xc,'filled') ylim([lconf-0.03 1.05]) hold on plot(lags,lconf*ones(size(lags)),'r','linewidth',2) plot(lags,upconf*ones(size(lags)),'r','linewidth',2) title('Sample Autocorrelation with 99% Confidence Intervals')
Снова, вы видите, что кроме нулевой задержки, значения дискретизации автокорреляции находятся в пределах 99% -доверие для автокорреляции последовательности белого шума. Из этого можно сделать вывод, что невязки являются белым шумом. Более конкретно, вы не можете отвергнуть, что невязки являются реализацией процесса белого шума.
Наконец, добавьте другую синусоиду с частотой 200 Гц и амплитудой 3/4. Подгонка только синусоида при 100 Гц и обнаружение автокорреляции образца остатков.
x = x+3/4*sin(2*pi*200*t);
xdft = fft(x);
ampest = 2/length(x)*xdft(101);
xfit = real(ampest)*cos(2*pi*100*t)+imag(ampest)*sin(2*pi*100*t);
residuals = x-xfit;
[xc,lags] = xcorr(residuals,50,'coeff');Постройте график автокорреляции выборки вместе с 99% -конфидентными интервалами.
figure stem(lags,xc,'filled') ylim([lconf-0.12 1.05]) hold on plot(lags,lconf*ones(size(lags)),'r','linewidth',2) plot(lags,upconf*ones(size(lags)),'r','linewidth',2) title('Sample Autocorrelation with 99% Confidence Intervals')
В этом случае значения автокорреляции явно превышают 99% -ные границы для автокорреляции белого шума при многих лагах. Здесь можно отвергнуть гипотезу о том, что невязки являются последовательностью белого шума. Подразумевается, что модель не учитывала весь сигнал, и, следовательно, невязки состоят из сигнала плюс шум.