Симулируйте регрессионные модели с ошибками ARMA

Симулируйте модель ошибки AR

В этом примере показано, как симулировать пути расчета из регрессионой модели с ошибками AR, не задавая предварительные нарушения порядка.

Задайте регрессионую модель с AR (2) ошибками:

$\begin{array}{l} y_{t} = 2 + X_{t} [\begin{array}{l} - 2 \\ 1.5 \end{array}] + u_{t} \\ u_{t} = 0.75 u_{t - 1} - 0.5 u_{t - 2} + ε_{t}, \end{array}$

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

Beta = [-2; 1.5];
Intercept = 2;
a1 = 0.75;
a2 = -0.5;
Variance = 1;
Mdl = regARIMA('AR',{a1, a2},'Intercept',Intercept,...
    'Beta',Beta,'Variance',Variance);

Сгенерируйте два ряда предикторов длины T = 50 путем случайного выбора из стандартного Гауссова распределения.

T = 50;
rng(1); % For reproducibility
X = randn(T,2);

Программа рассматривает предикторы как нестохастические ряды.

Сгенерируйте и постройте график одного пути откликов от Mdl.

rng(2);
ySim = simulate(Mdl,T,'X',X);

figure
plot(ySim)
title('{\bf Simulated Response Series}')

$Figure contains an axes. The axes with title {\bf Simulated Response Series} contains an object of type line.$

simulate требует P = 2 предварительных безусловных нарушения порядка ( $u_{t}$ ) для инициализации ряда ошибок. Без них, как в данном случае, simulate устанавливает необходимые предопределительные безусловные нарушения порядка равными 0.

Другой способ - фильтрация ряда случайных инноваций через Mdl использование filter.

rng(2);
e = randn(T,1);
yFilter = filter(Mdl,e,'X',X);

figure
plot(yFilter)
title('{\bf Simulated Response Series Using Filtered Innovations}')

$Figure contains an axes. The axes with title {\bf Simulated Response Series Using Filtered Innovations} contains an object of type line.$

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

Симулируйте 1000 путей отклика от Mdl. Оцените переходные эффекты путем построения графика безусловного нарушения порядка (U) отклонения по моделируемым путям в каждом периоде.

numPaths = 1000;
[Y,~,U] = simulate(Mdl,T,'NumPaths',numPaths,'X',X);

figure
h1 = plot(Y,'Color',[.85,.85,.85]);
title('{\bf 1000 Simulated Response Paths}')
hold on
h2 = plot(1:T,Intercept+X*Beta,'k--','LineWidth',2);
legend([h1(1),h2],'Simulated Path','Mean')
hold off

$Figure contains an axes. The axes with title {\bf 1000 Simulated Response Paths} contains 1001 objects of type line. These objects represent Simulated Path, Mean.$

figure
h1 = plot(var(U,0,2),'r','LineWidth',2);
hold on
theoVarFix = ((1-a2)*Variance)/((1+a2)*((1-a2)^2-a1^2));
h2 = plot([1 T],[theoVarFix theoVarFix],'k--','LineWidth',2);
title('{\bf Unconditional Disturbance Variance}')
legend([h1,h2],'Simulation Variance','Theoretical Variance')
hold off

$Figure contains an axes. The axes with title {\bf Unconditional Disturbance Variance} contains 2 objects of type line. These objects represent Simulation Variance, Theoretical Variance.$

Симулированные отклики следуют их теоретическому среднему значению, $c + X β$ , которая не является постоянной с течением времени (и может выглядеть нестационарно).

Отклонение процесса не постоянно, а выравнивается при теоретическом отклонении к 10-му периоду. Теоретическое отклонение модели ошибки AR (2)

$\frac{(1 - a_{2}) σ_{ε}^{2}}{(1 + a_{2}) [{(1 - a_{2})}^{2} - a_{1}^{2}]} = \frac{(1 + 0.5)}{(1 - 0.5) [{(1 + 0.5)}^{2} - {0.75}^{2}]} = 1.78$

Можно уменьшить переходные эффекты путем разбиения моделируемых данных на фрагмент записи и фрагмент для анализа. Не используйте фрагмент burn-in для анализа. Включите достаточно периодов в фрагмент горения, чтобы преодолеть переходные эффекты.

burnIn = 1:10;
notBurnIn = burnIn(end)+1:T;
Y = Y(notBurnIn,:);
X = X(notBurnIn,:);
U = U(notBurnIn,:);
figure
h1 = plot(notBurnIn,Y,'Color',[.85,.85,.85]);
hold on
h2 = plot(notBurnIn,Intercept+X*Beta,'k--','LineWidth',2);
title('{\bf 1000 Simulated Response Paths for Analysis}')
legend([h1(1),h2],'Simulated Path','Mean')
hold off

$Figure contains an axes. The axes with title {\bf 1000 Simulated Response Paths for Analysis} contains 1001 objects of type line. These objects represent Simulated Path, Mean.$

figure
h1 = plot(notBurnIn,var(U,0,2),'r','LineWidth',2);
hold on
h2 = plot([notBurnIn(1) notBurnIn(end)],...
    [theoVarFix theoVarFix],'k--','LineWidth',2);
title('{\bf Converged Unconditional Disturbance Variance}')
legend([h1,h2],'Simulation Variance','Theoretical Variance')
hold off

$Figure contains an axes. The axes with title {\bf Converged Unconditional Disturbance Variance} contains 2 objects of type line. These objects represent Simulation Variance, Theoretical Variance.$

Симуляция безусловного нарушения порядка отклонений колебаться вокруг теоретического отклонения из-за ошибки дискретизации Монте-Карло. Имейте в виду, что исключение образца ожога из анализа уменьшает эффективный размер выборки.

Симулируйте модель ошибки MA

Открыть Live Script

В этом примере показано, как моделировать ответы из регрессионой модели с ошибками MA, не задавая предварительную выборку.

Задайте регрессионую модель с ошибками MA (8 ):

$\begin{array}{l} y_{t} = 2 + X_{t} [\begin{array}{l} - 2 \\ 1.5 \end{array}] + u_{t} \\ u_{t} = ε_{t} + 0.4 ε_{t - 1} - 0.3 ε_{t - 4} + 0.2 ε_{t - 8}, \end{array}$

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

Beta = [-2; 1.5];
Intercept = 2;
b1 = 0.4;
b4 = -0.3;
b8 = 0.2;
Variance = 0.5;
Mdl = regARIMA('MA',{b1, b4, b8},'MALags',[1 4 8],...
    'Intercept',Intercept,'Beta',Beta,'Variance',Variance);

Сгенерируйте два ряда предикторов длины T = 100 путем случайного выбора из стандартного Гауссова распределения.

T = 100;
rng(4); % For reproducibility
X = randn(T,2);

Программа рассматривает предикторы как нестохастические ряды.

Сгенерируйте и постройте график одного пути откликов от Mdl.

rng(5);
ySim = simulate(Mdl,T,'X',X);

figure
plot(ySim)
title('{\bf Simulated Response Series}')

$Figure contains an axes. The axes with title {\bf Simulated Response Series} contains an object of type line.$

simulate требует Q = 8 предварительный образец инноваций ( $ε_{t}$ ) для инициализации ряда ошибок. Без них, как в данном случае, simulate устанавливает необходимые предварительные нововведения равными 0.

Кроме того, используйте filter для фильтрации ряда случайных инноваций через Mdl.

rng(5);
e = randn(T,1);
yFilter = filter(Mdl,e,'X',X);

figure
plot(yFilter)
title('{\bf Simulated Response Series Using Filtered Innovations}')

$Figure contains an axes. The axes with title {\bf Simulated Response Series Using Filtered Innovations} contains an object of type line.$

numPaths = 1000;
[Y,~,U] = simulate(Mdl,T,'NumPaths',numPaths,'X',X);

figure
h1 = plot(Y,'Color',[.85,.85,.85]);
title('{\bf 1000 Simulated Response Paths}')
hold on
h2 = plot(1:T,Intercept+X*Beta,'k--','LineWidth',2);
legend([h1(1),h2],'Simulated Path','Mean')
hold off

$Figure contains an axes. The axes with title {\bf 1000 Simulated Response Paths} contains 1001 objects of type line. These objects represent Simulated Path, Mean.$

figure
h1 = plot(var(U,0,2),'r','LineWidth',2);
hold on
theoVarFix = (1+b1^2+b4^2+b8^2)*Variance;
h2 = plot([1 T],[theoVarFix theoVarFix],'k--','LineWidth',2);
title('{\bf Unconditional Disturbance Variance}')
legend([h1,h2],'Simulation Variance','Theoretical Variance')
hold off

$Figure contains an axes. The axes with title {\bf Unconditional Disturbance Variance} contains 2 objects of type line. These objects represent Simulation Variance, Theoretical Variance.$

Моделируемые пути следуют своему теоретическому среднему, $c + X β$ , которая не является постоянной с течением времени (и может выглядеть нестационарно).

Отклонение процесса не постоянно, а выравнивается в теоретическом отклонении к 15-му периоду. Теоретическое отклонение модели ошибки MA (8)

$(1 + b_{1}^{2} + b_{4}^{2} + b_{8}^{2}) σ_{ε}^{2} = (1 + {0.4}^{2} + {(- 0.3)}^{2} + {0.2}^{2}) 0.5 = 0.645 .$

burnIn = 1:15;
notBurnIn = burnIn(end)+1:T;
Y = Y(notBurnIn,:);
X = X(notBurnIn,:);
U = U(notBurnIn,:);
figure
h1 = plot(notBurnIn,Y,'Color',[.85,.85,.85]);
hold on
h2 = plot(notBurnIn,Intercept+X*Beta,'k--','LineWidth',2);
title('{\bf 1000 Simulated Response Paths for Analysis}')
legend([h1(1),h2],'Simulated Path','Mean')
axis tight
hold off

$Figure contains an axes. The axes with title {\bf 1000 Simulated Response Paths for Analysis} contains 1001 objects of type line. These objects represent Simulated Path, Mean.$

figure
h1 = plot(notBurnIn,var(U,0,2),'r','LineWidth',2);
hold on
h2 = plot([notBurnIn(1) notBurnIn(end)],...
    [theoVarFix theoVarFix],'k--','LineWidth',2);
title('{\bf Converged Unconditional Disturbance Variance}')
legend([h1,h2],'Simulation Variance','Theoretical Variance')
axis tight
hold off

$Figure contains an axes. The axes with title {\bf Converged Unconditional Disturbance Variance} contains 2 objects of type line. These objects represent Simulation Variance, Theoretical Variance.$

Симулируйте модель ошибки ARMA

Открыть Live Script

В этом примере показано, как моделировать ответы из регрессионой модели с ошибками ARMA, не задавая предварительную выборку.

Задайте регрессионую модель с ошибками ARMA (2,1):

$\begin{array}{l} y_{t} = 2 + X_{t} [\begin{array}{l} - 2 \\ 1.5 \end{array}] + u_{t} \\ u_{t} = 0.9 u_{t - 1} - 0.1 u_{t - 2} + ε_{t} + 0.5 ε_{t - 1}, \end{array}$

где $ε_{t}$ распределяется с 15 степенями свободы и отклонения 1.

Beta = [-2; 1.5];
Intercept = 2;
a1 = 0.9;
a2 = -0.1;
b1 = 0.5;
Variance = 1;
Distribution = struct('Name','t','DoF',15);
Mdl = regARIMA('AR',{a1, a2},'MA',b1,...
    'Distribution',Distribution,'Intercept',Intercept,...
    'Beta',Beta,'Variance',Variance);

T = 50;
rng(6); % For reproducibility
X = randn(T,2);

Программа рассматривает предикторы как нестохастические ряды.

Сгенерируйте и постройте график одного пути откликов от Mdl.

rng(7);
ySim = simulate(Mdl,T,'X',X);

figure
plot(ySim)
title('{\bf Simulated Response Series}')

$Figure contains an axes. The axes with title {\bf Simulated Response Series} contains an object of type line.$

simulate требует:

P = 2 предварительная выборка безусловных нарушений порядка для инициализации авторегрессивного компонента ряда ошибок.
Q = 1 предварительный образец инноваций для инициализации компонента скользящего среднего значения ряда ошибок.

Без них, как в данном случае, simulate устанавливает необходимые ошибки предварительного образца равными 0.

Кроме того, используйте filter для фильтрации ряда случайных инноваций через Mdl.

rng(7);
e = randn(T,1);
yFilter = filter(Mdl,e,'X',X);

figure
plot(yFilter)
title('{\bf Simulated Response Series Using Filtered Innovations}')

$Figure contains an axes. The axes with title {\bf Simulated Response Series Using Filtered Innovations} contains an object of type line.$

numPaths = 1000;
[Y,~,U] = simulate(Mdl,T,'NumPaths',numPaths,'X',X);

figure
h1 = plot(Y,'Color',[.85,.85,.85]);
title('{\bf 1000 Simulated Response Paths}')
hold on
h2 = plot(1:T,Intercept+X*Beta,'k--','LineWidth',2);
legend([h1(1),h2],'Simulated Path','Mean')
hold off

$Figure contains an axes. The axes with title {\bf 1000 Simulated Response Paths} contains 1001 objects of type line. These objects represent Simulated Path, Mean.$

figure
h1 = plot(var(U,0,2),'r','LineWidth',2);
hold on
theoVarFix = Variance*(a1*b1*(1+a2)+(1-a2)*(1+a1*b1+b1^2))/...
    ((1+a2)*((1-a2)^2-a1^2));
h2 = plot([1 T],[theoVarFix theoVarFix],'k--','LineWidth',2);
title('{\bf Unconditional Disturbance Variance}')
legend([h1,h2],'Simulation Variance','Theoretical Variance',...
    'Location','Best')
hold off

$Figure contains an axes. The axes with title {\bf Unconditional Disturbance Variance} contains 2 objects of type line. These objects represent Simulation Variance, Theoretical Variance.$

Отклонение процесса не постоянно, а выравнивается при теоретическом отклонении к 10-му периоду. Теоретическое отклонение модели ошибки ARMA (2,1) является :

$\frac{σ_{ε}^{2} [a_{1} b_{1} (1 + a_{2}) + (1 - a_{2}) (1 + a_{1} b_{1} + b_{1}^{2})]}{{(1 + a_{2})}^{2} [{(1 - a_{2})}^{2} - a_{1}^{2}]} = \frac{[0.9 (0.5) (1 - 0.1) + (1 + 0.1) (1 + 0.9 (0.5) + {0.5}^{2})]}{{(1 - 0.1)}^{2} [{(1 + 0.1)}^{2} - {0.9}^{2}]} = 6 . 32 .$

burnIn = 1:10;
notBurnIn = burnIn(end)+1:T;
Y = Y(notBurnIn,:);
X = X(notBurnIn,:);
U = U(notBurnIn,:);
figure
h1 = plot(notBurnIn,Y,'Color',[.85,.85,.85]);
hold on
h2 = plot(notBurnIn,Intercept+X*Beta,'k--','LineWidth',2);
title('{\bf 1000 Simulated Response Paths for Analysis}')
legend([h1(1),h2],'Simulated Path','Mean')
axis tight
hold off

$Figure contains an axes. The axes with title {\bf 1000 Simulated Response Paths for Analysis} contains 1001 objects of type line. These objects represent Simulated Path, Mean.$

figure
h1 = plot(notBurnIn,var(U,0,2),'r','LineWidth',2);
hold on
h2 = plot([notBurnIn(1) notBurnIn(end)],...
    [theoVarFix theoVarFix],'k--','LineWidth',2);
title('{\bf Converged Unconditional Disturbance Variance}')
legend([h1,h2],'Simulation Variance','Theoretical Variance')
axis tight
hold off

$Figure contains an axes. The axes with title {\bf Converged Unconditional Disturbance Variance} contains 2 objects of type line. These objects represent Simulation Variance, Theoretical Variance.$

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