Оцените устойчивость модели в пространстве состояний Используя прокручивающийся анализ окна

Оцените устойчивость модели Используя прокручивающийся анализ окна

В этом примере показано, как использовать прокручивающееся окно, чтобы проверять, независимы ли параметры модели timeseries от времени. Этот пример анализирует два временных рядов:

Timeseries 1: симулированные данные из известной, независимой от времени модели
Timeseries 2: симулированные данные из известной, изменяющейся во времени модели

Полностью задайте эту модель AR (1) для Timeseries 1:

$y_{t} = 0.6 y_{t - 1} + ε_{t}$

где $ε_{t}$ является Гауссовым со средним значением 0 и отклонением 1. Полностью задайте эту изменяющуюся во времени модель для Timeseries 2:

$\begin{array}{cccccccccccccccccccc} y_{t} = 0.2 y_{t - 1} + ε_{t}; t = 1, . . ., 100 \\ y_{t} = 0.75 y_{t - 1} + ε_{t}; t = 101, . . ., 150 \\ y_{t} = - 0.5 y_{t - 1} + ε_{t}; t = 151, . . ., 200, \end{array}$

где $ε_{t}$ является Гауссовым со средним значением 0 и отклонением 1.

Mdl1 = arima('AR',0.6,'Constant',0,'Variance',1);

Mdl2 = cell(3,1); % Preallocate
ARMdl2 = [0.2 0.75 -0.5];
for j = 1:3
    Mdl2{j} = arima('AR',ARMdl2(j),'Constant',0,'Variance',1);
end

Mdl1 arima объекты модели. Можно получить доступ к его свойствам с помощью записи через точку. Mdl2 массив ячеек arima объекты модели. Вы можете вы индексация ячейки и запись через точку, чтобы получить доступ к свойствам моделей в Mdl2. Например, чтобы получить доступ к значению параметров AR третьей модели в Mdl3, введите Mdl2{3}.AR.

Симулируйте T = 200 периодов данных из Mdl1 и Mdl2. Используйте преддемонстрационный ответ 0 для обоих рядов.

rng(1); % For reproducibility
T = 200;
y1 = simulate(Mdl1,T,'Y0',0);
timeMdl2 = [100 50 50]; % Number of observations per model in Mdl2
y2 = 0;
for k = 1:numel(Mdl2)
    y2 = [y2; simulate(Mdl2{k},timeMdl2(k),'Y0',y2(end))];
end

Y = [y1 y2(2:end)];

Задайте пустые модели AR (1) для оценки Mdl1, Mdl2, и Mdl3. Оцените все три модели с помощью соответствующих наборов данных и прокручивающегося размера окна 40 периодов. Кроме того, используйте прокручивающийся шаг окна одного периода. Сохраните авторегрессивные параметры и оцененное инновационное отклонение.

Mdl = arima('ARLags',1,'Constant',0);

m = 100;       % Rolling window size
N = T - m + 1; % Number of rolling windows

EstParams = cell(2,1);     % Preallocate for estimates
EstParamsMat = zeros(N,2);
EstParamsSE = cell(2,1);
EstParamsSEMat = zeros(N,2);

for j = 1:2
    for k = 1:N
        idxRW =  k:(m + k - 1); % In-sample indices
        [EstMdl,EstParamCov] = estimate(Mdl,Y(idxRW,j),'Display','off');
        EstParamsMat(k,:) = [EstMdl.AR{1} EstMdl.Variance];
        EstParamsSEMat(k,:) = sqrt([EstParamCov(2,2) EstParamCov(3,3)]);
    end
    EstParams{j} = EstParamsMat;
    EstParamsSE{j} = EstParamsSEMat;
end

Постройте оценки и их мудрые точкой доверительные интервалы по прокручивающемуся индексу окна.

titleMdls = {'Mdl1','Mdl2'};

for j = 1:2
    figure;
    subplot(2,1,1);
    Estimates = EstParams{j};
    SEs = EstParamsSE{j};
    plot(Estimates(:,1),'LineWidth',2);
    hold on;
    plot(Estimates(:,1) + 2*SEs(:,1),'r:','LineWidth',2);
    plot(Estimates(:,1) - 2*SEs(:,1),'r:','LineWidth',2);
    title(sprintf('%s - AR at Lag 1 Estimate',titleMdls{j}));
    xlabel 'Rolling window index';
    axis tight;
    hold off;
    subplot(2,1,2);
    plot(Estimates(:,2),'LineWidth',2);
    hold on;
    plot(Estimates (:,2) + 2*SEs(:,2),'r:','LineWidth',2);
    plot(Estimates(:,2) - 2*SEs(:,2),'r:','LineWidth',2);
    title(sprintf('%s - Variance Estimate',titleMdls{j}));
    xlabel 'Rolling window index';
    axis tight;
    hold off;
end

Для Mdl1, оценка AR не варьируется очень от 0,6, и оценки не существенно отличаются друг от друга (попарно). Подобные результаты происходят для отклонения Mdl1. Оценка AR Mdl2 растет, и затем падает, который указывает на временную зависимость. Кроме того, на основе доверительных интервалов существует доказательство, что некоторые оценки отличаются от других. Хотя отклонение не изменилось в процессе моделирования, кажется, существует heteroscedasticity, возможно вызванный нестабильностью модели.

Оцените устойчивость неявно созданной модели в пространстве состояний

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

Рассмотрите эту модель в пространстве состояний:

$\begin{matrix} x_{t} = ϕ x_{t - 1} + ε_{t} \\ y_{t} - β z_{t} = x_{t} + u_{t} \end{matrix},$

где _εt и _ut являются Гауссовым процессом со средним значением 0 и отклонением 1. Создайте функциональный rwParamMap.m, который задает как параметры в params сопоставьте с матрицами модели в пространстве состояний, значениями начального состояния и типом состояния, и сохраните его в своей рабочей папке.

function [A,B,C,D,Mean0,Cov0,StateType,deflateY] = rwParamMap(params,y,Z)
%rwParamMap Parameter-to-matrix mapping function for rolling window example
%using ssm and specifying an AR(1) state model
%   The state space model specified by rwParamMap contains a stationary
%   AR(1) state, the observation model includes a regression component, and
%   the variances of the innovation and disturbances are 1. The response y
%   is deflated by the regression component specified by the predictor
%   variables x.
A = params(1);
B = 1;
C = 1;
D = 1;
Mean0 = [];
Cov0 = [];
StateType = 0;
deflateY = y - params(2)*Z;
end

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

Mdl2Sim = ssm(NaN,1,1,1,'StateType',1);

Mdl2Sim неявно заданный ssm объект.

Симулируйте путь с 200 периодами случайных стандартных Гауссовых данных. Затем симулируйте ответы от Mdl2Sim, и раздуйте ответы с компонентом регрессии. В данном примере используйте ϕ = 0.6 и β = 2.

rng(1); % For reproducibility
T = 200;
Z = randn(T,1);
phi = 0.6;
beta = 2;
deflateY = simulate(Mdl2Sim,T,'Params',phi);
y = deflateY + Z*beta;

y расширенный, симулированный путь к ответу и Z симулированный ряд предиктора.

Если вы задаете модель в пространстве состояний неявно и ответ и данные о предикторе (т.е. y и Z) существуйте в MATLAB^® Workspace, затем программное обеспечение создает ссылку из функции отображения параметра к матрице те ряды. Если данные не существуют в рабочем пространстве MATLAB, то программное обеспечение создает модель, но необходимо обеспечить данные с помощью соответствующих аргументов пары "имя-значение", когда вы, например, оцениваете модель.

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

m = 100; 
N = T - m + 1;        % Number of rolling windows

EstParams = nan(N,2); % Preallocation
EstParamSE = nan(N,2);

for j = 1:N;
    idxRW = j:(m + j - 1);
    Mdl = ssm(@(c)rwParamMap(c,y(idxRW),Z(idxRW)));
    [~,EstParams(j,:),EstParamCov] = estimate(Mdl,y(idxRW),[0.5 1]',...
        'Display','off');
    EstParamSE(j,:) = sqrt(diag(EstParamCov));
end

Постройте оценки и мудрые точкой доверительные интервалы для параметра AR и коэффициента регрессии.

figure;
subplot(2,1,1);
plot(EstParams(:,1),'LineWidth',2);
hold on;
plot(EstParams(:,1) + 2*EstParamSE(:,1),':r','LineWidth',2);
plot(EstParams(:,1) - 2*EstParamSE(:,1),':r','LineWidth',2);
title 'State AR Estimate at Lag 1';
xlabel 'Rolling window index';
axis tight;
hold off;
subplot(2,1,2);
plot(EstParams(:,2),'LineWidth',2);
hold on;
plot(EstParams(:,2) + 2*EstParamSE(:,2),':r','LineWidth',2);
plot(EstParams(:,2) - 2*EstParamSE(:,2),':r','LineWidth',2);
title 'Regression Coefficient Estimate';
xlabel 'Rolling window index';
axis tight;
hold off;

Графики показывают, что модель устойчива, поскольку оценка AR не отклоняется очень от своего среднего значения, и при этом коэффициент регрессии не оценивает.

Смотрите также

estimate | forecast | ssm

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