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

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

Предположим, что модель AR (2) и MA (1) содержат латентный процесс. Существует 50 периодов, и процесс MA (1) выпадает из модели на последние 25 периодов. Уравнение состояния для первых 25 периодов

$$\begin{array}{l}
{x_{1,t}} = 0.7{x_{1,t - 1}} - 0.2{x_{1,t - 2}} + {u_{1,t}}\\
{x_{2,t}} = {u_{2,t}} + 0.6{u_{2,t - 1}},
\end{array}$$

Для последних 25 периодов уравнение состояния является

$${x_{1,t}} = 0.7{x_{1,t - 1}} - 0.2{x_{1,t - 2}} + {u_{1,t}},$$

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

Сгенерируйте случайную серию из 50 наблюдений от$x_{1,t}$ и, $x_{2,t}$принимая, что серия начинается с 1,5 и 1, соответственно.

T = 50;
ARMdl = arima('AR',{0.7,-0.2},'Constant',0,'Variance',1);
MAMdl = arima('MA',0.6,'Constant',0,'Variance',1);
x0 = [1.5 1; 1.5 1];
rng(1);
x = [simulate(ARMdl,T,'Y0',x0(:,1)),...
    [simulate(MAMdl,T/2,'Y0',x0(:,2));nan(T/2,1)]];

Последние 25 значений для моделируемых данных MA (1) отсутствуют.

Предположим далее, что латентные процессы измеряются с помощью

$${y_t} = 2\left( {{x_{1,t}} + {x_{2,t}}} \right) + {\varepsilon _t},$$

за первые 25 периодов, и

$${y_t} = 2{x_{1,t}} + {\varepsilon _t}$$

за последние 25 периодов.$\varepsilon_t$ является Гауссовым со средним 0 и стандартным отклонением 1.

Сгенерируйте наблюдения, используя процесс случайного скрытого состояния (x) и уравнение наблюдения.

y = 2*sum(x','omitnan')'+randn(T,1);

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

$$\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{{x_{1,t}}}\\
{{x_{2,t}}}\\
{{x_{3,t}}}\\
{{x_{4,t}}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{\phi _1}}&{{\phi _2}}&0&0\\
1&0&0&0\\
0&0&0&{{\theta _1}}\\
0&0&0&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_{1,t - 1}}}\\
{{x_{2,t - 1}}}\\
{{x_{3,t - 1}}}\\
{{x_{4,t - 1}}}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
1&0\\
0&0\\
0&1\\
0&1
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{u_{1,t}}}\\
{{u_{2,t}}}
\end{array}} \right]\\
{y_t} = a({x_{1,t}} + {x_{3,t}}) + {\varepsilon _t}
\end{array}{\rm for\;}t = 1,...,25,$$

$$\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{{x_{1,t}}}\\
{{x_{2,t}}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{\phi _1}}&{{\phi _2}}&0&0\\
1&0&0&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_{1,t - 1}}}\\
{{x_{2,t - 1}}}\\
{{x_{3,t - 1}}}\\
{{x_{4,t - 1}}}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
1\\
0
\end{array}} \right]{u_{1,t}}\\
{y_t} = b{x_{1,t}} + {\varepsilon _t}
\end{array}{\rm for\;}t = 26,$$

$$\begin{array}{l}
\left[ {\begin{array}{*{20}{c}}
{{x_{1,t}}}\\
{{x_{2,t}}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{\phi _1}}&{{\phi _2}}\\
1&0
\end{array}} \right]\left[ {\begin{array}{*{20}{c}}
{{x_{1,t - 1}}}\\
{{x_{2,t - 1}}}
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
1\\
0
\end{array}} \right]{u_{1,t}}\\
{y_t} = b{x_{1,t}} + {\varepsilon _t}
\end{array}{\rm for\;}t = 27,...,50$$

Написание функции, которая задает, как параметры в params сопоставить с матрицами модели пространства состояний, начальными значениями состояний и типом состояния.


% Copyright 2015 The MathWorks, Inc.

function [A,B,C,D,Mean0,Cov0,StateType] = AR2MAParamMap(params,T)
%AR2MAParamMap Time-variant state-space model parameter mapping function
%
% This function maps the vector params to the state-space matrices (A, B,
% C, and D), the initial state value and the initial state variance (Mean0
% and Cov0), and the type of state (StateType). From periods 1 to T/2, the
% state model is an AR(2) and an MA(1) model, and the observation model is
% the sum of the two states. From periods T/2 + 1 to T, the state model is
% just the AR(2) model.
    A1 = {[params(1) params(2) 0 0; 1 0 0 0; 0 0 0 params(3); 0 0 0 0]};
    B1 = {[1 0; 0 0; 0 1; 0 1]}; 
    C1 = {params(4)*[1 0 1 0]};
    Mean0 = ones(4,1);
    Cov0 = 10*eye(4);
    StateType = [0 0 0 0];
    A2 = {[params(1) params(2) 0 0; 1 0 0 0]};
    B2 = {[1; 0]};
    A3 = {[params(1) params(2); 1 0]};
    B3 = {[1; 0]}; 
    C3 = {params(5)*[1 0]};
    A = [repmat(A1,T/2,1);A2;repmat(A3,(T-2)/2,1)];
    B = [repmat(B1,T/2,1);B2;repmat(B3,(T-2)/2,1)];
    C = [repmat(C1,T/2,1);repmat(C3,T/2,1)];
    D = 1;
end

Сохраните этот код в файле с именем AR2MAParamMap и поместите его в путь MATLAB ®.

Создайте модель пространства состояний путем передачи функции AR2MAParamMap как указатель на функцию, чтобы ssm.

Mdl = ssm(@(params)AR2MAParamMap(params,T));

ssm неявно задает модель пространства состояний. Обычно вы не можете проверить неявно определенные модели пространства состояний.

Передайте наблюдаемые отклики (y) к estimate для оценки параметров. Задайте положительные начальные значения для неизвестных параметров.

params0 = 0.1*ones(5,1);
EstMdl = estimate(Mdl,y,params0)
Method: Maximum likelihood (fminunc)
Sample size: 50
Logarithmic  likelihood:     -114.957
Akaike   info criterion:      239.913
Bayesian info criterion:      249.473
      |     Coeff       Std Err   t Stat     Prob  
---------------------------------------------------
 c(1) |  0.47870       0.26634    1.79733  0.07229 
 c(2) |  0.00809       0.27179    0.02975  0.97626 
 c(3) |  0.55735       0.80958    0.68844  0.49118 
 c(4) |  1.62679       0.41622    3.90848  0.00009 
 c(5) |  1.90021       0.49563    3.83391  0.00013 
      |                                            
      |   Final State   Std Dev    t Stat    Prob  
 x(1) | -0.81229       0.46815   -1.73511  0.08272 
 x(2) | -0.31449       0.45918   -0.68490  0.49341 

EstMdl = 


State-space model type: <a href="matlab: doc ssm">ssm</a>

State vector length: Time-varying
Observation vector length: 1
State disturbance vector length: Time-varying
Observation innovation vector length: 1
Sample size supported by model: 50

State variables: x1, x2,...
State disturbances: u1, u2,...
Observation series: y1, y2,...
Observation innovations: e1, e2,...

State equations of period 1, 2, 3,..., 25:
x1(t) = (0.48)x1(t-1) + (8.09e-03)x2(t-1) + u1(t)
x2(t) = x1(t-1)
x3(t) = (0.56)x4(t-1) + u2(t)
x4(t) = u2(t)

State equations of period 26:
x1(t) = (0.48)x1(t-1) + (8.09e-03)x2(t-1) + u1(t)
x2(t) = x1(t-1)

State equations of period 27, 28, 29,..., 50:
x1(t) = (0.48)x1(t-1) + (8.09e-03)x2(t-1) + u1(t)
x2(t) = x1(t-1)


Observation equation of period 1, 2, 3,..., 25:
y1(t) = (1.63)x1(t) + (1.63)x3(t) + e1(t)

Observation equation of period 26, 27, 28,..., 50:
y1(t) = (1.90)x1(t) + e1(t)


Initial state distribution:

Initial state means
 x1  x2  x3  x4 
  1   1   1   1 

Initial state covariance matrix
     x1  x2  x3  x4 
 x1  10  0   0   0  
 x2  0   10  0   0  
 x3  0   0   10  0  
 x4  0   0   0   10 

State types
     x1          x2          x3          x4     
 Stationary  Stationary  Stationary  Stationary 

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

См. также

| | | |

Похожие примеры

Подробнее о