Создайте частично заданную модель для оценки

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

P = NaN(2);
mc = dtmc(P,'StateNames',["Expansion" "Recession"]);
mdl = arima(0,0,0);
Mdl = msVAR(mc,[mdl; mdl]);

Mdl частично заданный msVAR объект. NaN- ценные элементы Switch и SubModels свойства указывают на допускающие оценку параметры.

Создайте полностью заданную модель, содержащую начальную стоимость

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

P0 = 0.5*ones(2);
mc0 = dtmc(P0,'StateNames',Mdl.StateNames);
mdl01 = arima('Constant',1,'Variance',1);
mdl02 = arima('Constant',-1,'Variance',1);
Mdl0 = msVAR(mc0,[mdl01; mdl02]);

Mdl0 полностью заданный msVAR объект.

Загрузите и предварительно обработайте данные

Загрузите набор данных GDP США.

load Data_GDP

Data содержит ежеквартальные измерения США действительный GDP в период 1947:Q1–2005:Q2. Период оценки в [1] является 1947:Q2–2004:Q2. Для получения дополнительной информации о наборе данных введите Description в командной строке.

Преобразуйте данные к ряду годового показателя:

qrate = diff(Data(2:230))./Data(2:229); % Quarterly rate
arate = 100*((1 + qrate).^4 - 1);       % Annualized rate

Оценочная модель

Соответствуйте модели Mdl к серии arate годового показателя. Задайте Mdl0 как модель, содержащая начальные допускающие оценку значения параметров.

EstMdl = estimate(Mdl,Mdl0,arate);

EstMdl предполагаемое (полностью заданный) переключающая Маркова модель динамической регрессии. EstMdl.Switch предполагаемая модель дискретной цепи Маркова (dtmc объект), и EstMdl.Submodels вектор из предполагаемых одномерных моделей VAR (0) (varm объекты.

Отобразите предполагаемые динамические модели состояния специфичные.

EstMdlExp = EstMdl.Submodels(1)

EstMdlExp = 
  varm with properties:

     Description: "ARIMA(0,0,0) Model (Gaussian Distribution)"
     SeriesNames: "Y1" 
       NumSeries: 1
               P: 0
        Constant: 4.90146
              AR: {}
           Trend: 0
            Beta: [1×0 matrix]
      Covariance: 12.087
 
   1-Dimensional VAR(0) Model

EstMdlRec = EstMdl.Submodels(2)

EstMdlRec = 
  varm with properties:

     Description: "ARIMA(0,0,0) Model (Gaussian Distribution)"
     SeriesNames: "Y1" 
       NumSeries: 1
               P: 0
        Constant: 0.0084884
              AR: {}
           Trend: 0
            Beta: [1×0 matrix]
      Covariance: 12.6876
 
   1-Dimensional VAR(0) Model

Отобразите предполагаемую матрицу переходов.

EstP = EstMdl.Switch.P

EstP = 2×2

    0.9088    0.0912
    0.2303    0.7697

Создайте модель для DGP

Создайте полностью заданную, модель дискретной цепи Маркова с двумя состояниями для переключающегося механизма.

P = [0.7 0.3; 0.1 0.9];
mc = dtmc(P);

Для каждого состояния создайте полностью заданную модель AR (1) для процесса ответа.

% Constants
C1 = 5;
C2 = -2;

% Autoregression coefficients
AR1 = 0.4;
AR2 = 0.2;

% Variances
V1 = 4;
V2 = 2;

% AR Submodels
dgp1 = arima('Constant',C1,'AR',AR1,'Variance',V1);
dgp2 = arima('Constant',C2,'AR',AR2,'Variance',V2);

Создайте полностью заданную переключающую Маркова модель динамической регрессии для DGP.

DGP = msVAR(mc,[dgp1,dgp2]);

Симулируйте пути к ответу от DGP

Сгенерируйте 10 случайных путей к ответу длины 1000 от DGP.

rng(1); % For reproducibility
N = 10;
n = 1000;
Data = simulate(DGP,n,'Numpaths',N);

Data 1000 10 матрица симулированных откликов.

Создайте модель для оценки

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

PEst = NaN(2);
mcEst = dtmc(PEst);
mdl = arima(1,0,0);
Mdl = msVAR(mcEst,[mdl; mdl]);

Создайте модель, содержащую начальную стоимость

Создайте полностью заданную переключающую Маркова модель динамической регрессии, которая имеет ту же структуру как Mdl, но установленный все допускающие оценку параметры на начальные значения.

P0 = 0.5*ones(2);
mc0 = dtmc(P0);
mdl01 = arima('Constant',1,'AR',0.5,'Variance',2);
mdl02 = arima('Constant',-1,'AR',0.5,'Variance',1);
Mdl0 = msVAR(mc0,[mdl01,mdl02]);

Оцените модели

Подбирайте модель к каждому симулированному пути. Для каждого пути постройте логарифмическую правдоподобность в каждой итерации алгоритма EM.

c1 = zeros(N,1);
c2 = zeros(N,1);
v1 = zeros(N,1);
v2 = zeros(N,1);
ar1 = zeros(N,1);
ar2 = zeros(N,1);
PStack = zeros(2,2,N);

figure
hold on

for i = 1:N

    EstModel = estimate(Mdl,Mdl0,Data(:,i),'IterationPlot',true);
    
    c1(i) = EstModel.Submodels(1).Constant;
    c2(i) = EstModel.Submodels(2).Constant;
    v1(i) = EstModel.Submodels(1).Covariance;
    v2(i) = EstModel.Submodels(2).Covariance;
    ar1(i) = EstModel.Submodels(1).AR{1};
    ar2(i) = EstModel.Submodels(2).AR{1};
    PStack(:,:,i) = EstModel.Switch.P;
    
end

hold off

Оцените точность

Вычислите среднее значение Монте-Карло каждого предполагаемого параметра.

c1Mean = mean(c1);
c2Mean = mean(c2);
v1Mean = mean(v1);
v2Mean = mean(v2);
ar1Mean = mean(ar1);
ar2Mean = mean(ar2);
PMean = mean(PStack,3);

Сравните параметры населения с соответствующими оценками Монте-Карло.

DGPvsEstimate = [...
    C1  c1Mean
    C2  c2Mean
    V1  v1Mean
    V2  v2Mean
    AR1 ar1Mean
    AR2 ar2Mean]

DGPvsEstimate = 6×2

    5.0000    5.0260
   -2.0000   -1.9615
    4.0000    3.9710
    2.0000    1.9903
    0.4000    0.4061
    0.2000    0.2017

P

P = 2×2

    0.7000    0.3000
    0.1000    0.9000

PEstimate = PMean

PEstimate = 2×2

    0.7065    0.2935
    0.1023    0.8977

Создайте модель для DGP

Создайте полностью заданную, модель дискретной цепи Маркова с тремя состояниями для переключающегося механизма.

P = [0.8 0.1 0.1; 0.2 0.6 0.2; 0 0.1 0.9];
mc = dtmc(P);

Для каждого состояния создайте полностью заданную модель VARX(1) для процесса ответа. Задайте ту же константу модели и изолируйте 1 матрицу коэффициентов AR для всех подмоделей. Для каждой модели задайте различное 2 1 коэффициент регрессии для одной внешней переменной.

% Constants
C = [1;-1];

% Autoregression coefficients
AR = {[0.6 0.1; 0.4 0.2]};

% Regression coefficients
Beta1 = [0.2;-0.4];
Beta2 = [0.6;-1.0];
Beta3 = [0.9;-1.3];

% VAR Submodels
dgp = varm('Constant',C,'AR',AR,'Covariance',5*eye(2));
dgp1 = dgp;
dgp1.Beta = Beta1;
dgp2 = dgp;
dgp2.Beta = Beta2;
dgp3 = dgp;
dgp3.Beta = Beta3;

Создайте полностью заданную переключающую Маркова модель динамической регрессии для DGP.

DGP = msVAR(mc,[dgp1; dgp2; dgp3]);

Симулируйте данные из DGP

Симулируйте данные для внешнего ряда путем генерации 1 000 наблюдений от Распределения Гаусса со средним значением 0 и отклонением 100.

rng(1); % For reproducibility
X = 10*randn(1000,1);

Сгенерируйте случайный путь к ответу длины 1000 от DGP. Задайте симулированные внешние данные для компонентов регрессии подмодели.

Data = simulate(DGP,1000,'X',X);

Data 1000 1 вектор из симулированных откликов.

Создайте модель для оценки

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

PEst = NaN(3);
mcEst = dtmc(PEst);
mdl = varm(2,1);
mdl.Constant = C;
mdl.AR = AR;
mdl.Beta = NaN(2,1);
Mdl = msVAR(mcEst,[mdl; mdl; mdl]);

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

Создайте модель, содержащую начальную стоимость

Создайте полностью заданную переключающую Маркова модель динамической регрессии, которая имеет ту же структуру как Mdl, но установленный все допускающие оценку параметры на начальные значения и установленные параметры с ограничениями равенства к их значениям заданы в Mdl.

P0 = (1/3)*ones(3);
mc0 = dtmc(P0);
mdl0 = varm('Constant',C,'AR',AR,'Covariance',eye(2));
mdl01 = mdl0;
mdl01.Beta = [0.1;-0.1];
mdl02 = mdl0;
mdl02.Beta = [0.5;-0.5];
mdl03 = mdl0;
mdl03.Beta = [1;-1];
Mdl0 = msVAR(mc0,[mdl01; mdl02; mdl03]);

Оценочная модель

Подбирайте модель к симулированным данным. Задайте внешние данные для компонента регрессии. Постройте логарифмическую правдоподобность в каждой итерации алгоритма EM.

figure
EstMdl = estimate(Mdl,Mdl0,Data,'X',X,'IterationPlot',true);

Оцените точность

Сравните предполагаемые векторы коэффициента регрессии и матрицу перехода к их истинным значениям.

Beta1

Beta1 = 2×1

    0.2000
   -0.4000

Beta1Estimate = EstMdl.Submodels(1).Beta

Beta1Estimate = 2×1

    0.1596
   -0.4040

Beta2

Beta2 = 2×1

    0.6000
   -1.0000

Beta2Estimate = EstMdl.Submodels(2).Beta

Beta2Estimate = 2×1

    0.5888
   -0.9771

Beta3

Beta3 = 2×1

    0.9000
   -1.3000

Beta3Estimate = EstMdl.Submodels(3).Beta

Beta3Estimate = 2×1

    0.8987
   -1.2991

P

P = 3×3

    0.8000    0.1000    0.1000
    0.2000    0.6000    0.2000
         0    0.1000    0.9000

PEstimate = EstMdl.Switch.P

PEstimate = 3×3

    0.7787    0.0856    0.1357
    0.1366    0.6906    0.1727
    0.0086    0.0787    0.9127