Предположим, что скрытый процесс является случайным обходом. Впоследствии, уравнение состояния

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

Сгенерируйте случайную последовательность 100 наблюдений от $$x_t$$, предположение, что ряд запускается в 1,5.

T = 100;
rng(1); % For reproducibility 
u = randn(T,1);
x = cumsum([1.5;u]);
x = x(2:end);

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

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

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

y = x + randn(T,1);

Вместе, скрытые уравнения процесса и наблюдения составляют модель в пространстве состояний. Примите, что состояние является стационарным AR (1) процесс. Затем модель в пространстве состояний, чтобы оценить

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

A = NaN;
B = NaN;
C = 1;
D = NaN;

Задайте модель в пространстве состояний с помощью содействующих матриц. Укажите, что распределение начального состояния является стационарным использованием аргумента пары "имя-значение" StateType.

StateType = 0;
Mdl = ssm(A,B,C,D,'StateType',StateType);

Mdl является моделью ssm. Программное обеспечение устанавливает значения для среднего значения начального состояния и отклонения. Проверьте, что модель задана правильно с помощью отображения в Командном окне.

Передайте наблюдения estimate, чтобы оценить параметры. Для параметров params0, которые вряд ли будут соответствовать их истинным значениям. Кроме того, задайте ограничения нижней границы 0 для стандартных отклонений.

params0 = [-1e7 1e-6 2000];
EstMdl = estimate(Mdl,y,params0,'lb',[-Inf,0,0]);

Warning: Covariance matrix of estimators cannot be computed precisely due to inversion difficulty. Check parameter identifiability. Also try different starting values and other options to compute the covariance matrix.

Method: Maximum likelihood (fmincon)
Sample size: 100
Logarithmic  likelihood:     -2464.23
Akaike   info criterion:      4934.46
Bayesian info criterion:      4942.27
      |      Coeff           Std Err       t Stat     Prob  
------------------------------------------------------------
 c(1) |   -9.99977e+06      9.99977e+05  -10.00000   0      
 c(2) |    1.23086e+05      1.91567e+13    0.00000  1.00000 
 c(3) | 2006.86501          3.12341e+11    0.00000  1.00000 
      |                                                     
      |     Final State      Std Dev        t Stat    Prob  
 x(1) |   -3.37649       1999.42392       -0.00169  0.99865 

estimate не удалось сходиться, и таким образом, результатами является нежелательный.

Совершенствуйте params0 с помощью refine.

Output = refine(Mdl,y,params0);
logL = cell2mat({Output.LogLikelihood})';
[~,maxLogLIndx] = max(logL)

maxLogLIndx = 2

refinedParams0 = Output(maxLogLIndx).Parameters

refinedParams0 = 1×3

    0.9705   -0.8934    0.9330

Description = Output(maxLogLIndx).Description

Description = 
'Nelder-Mead simplex'

Алгоритмом, который приводит к самому высокому loglikelihood значению, является Loose bound interior point, который является третьим struct в массиве структур Output.

Оцените Mdl с помощью refinedParams0, который является вектором усовершенствованных начальных значений параметров.

EstMdl = estimate(Mdl,y,refinedParams0,'lb',[-Inf,0,0]);

Method: Maximum likelihood (fmincon)
Sample size: 100
Logarithmic  likelihood:     -181.379
Akaike   info criterion:      368.758
Bayesian info criterion:      376.574
      |     Coeff       Std Err   t Stat     Prob  
---------------------------------------------------
 c(1) |  0.97050       0.02863   33.90367   0      
 c(2) |  0.89343       0.18521    4.82401  0.00000 
 c(3) |  0.93303       0.15176    6.14806   0      
      |                                            
      |   Final State   Std Dev    t Stat    Prob  
 x(1) | -3.93007       0.72066   -5.45343   0      

estimate сходился, делая оценки параметра намного более желательными. Коэффициент модели AR в двух стандартных погрешностях 1, который предполагает, что процессы состояния являются случайным обходом.

Предположим, что отношение между уровнем безработицы и номинальным валовым национальным продуктом (nGNP) линейно. Предположим далее, что уровень безработицы является серией AR (1). Символически, и в форме пространства состояний, модель

где:

load Data_NelsonPlosser

isNaN = any(ismissing(DataTable),2);  % Flag periods containing NaNs
gnpn = DataTable.GNPN(~isNaN);
y = DataTable.UR(~isNaN);
y = diff(y);
T = size(y,1);
Z = [ones(T,1) price2ret(gnpn)];

A = NaN;
B = NaN;
C = 1;

Mdl = ssm(A,B,C);

params0 = [150 1000]; % Initial values chosen arbitrarily
Beta0 = [1 -100];
Output = refine(Mdl,y,params0,'Predictors',Z,'Beta0',Beta0);

logL = cell2mat({Output.LogLikelihood})';
[~,maxLogLIndx] = max(logL)
refinedParams0 = Output(maxLogLIndx).Parameters
Description = Output(maxLogLIndx).Description

maxLogLIndx =

     2


refinedParams0 =

    0.0000   -1.3441    1.3477  -24.4336


Description =

    'Nelder-Mead simplex'

EstMdl = estimate(Mdl,y,refinedParams0(1:(end - 2)),'Predictors',Z,...
    'Beta0',refinedParams0((end - 1):end));

Method: Maximum likelihood (fminunc)
Sample size: 61
Logarithmic  likelihood:     -103.321
Akaike   info criterion:      214.642
Bayesian info criterion:      223.085
           |      Coeff       Std Err    t Stat     Prob  
----------------------------------------------------------
 c(1)      |   0.20499       0.12217     1.67793  0.09336 
 c(2)      |  -1.31586       0.08283   -15.88649   0      
 y <- z(1) |   1.38082       0.23315     5.92241   0      
 y <- z(2) | -24.87986       1.76909   -14.06365   0      
           |                                              
           |    Final State   Std Dev     t Stat    Prob  
 x(1)      |   1.19607        0           Inf      0