Этот пример показывает, как результаты симуляции модели пространства состояний сглаживаются (simsmooth
) сравнить с сглаженными состояниями (smooth
).
Предположим, что связь между изменением уровня безработицы () и номинальные темпы роста валового национального продукта (nGNP) () может быть выражена в следующей форме модели пространства состояний.
где:
- изменение уровня безработицы в момент t.
является фиктивным состоянием для эффекта MA (1) на.
- скорость роста nGNP в момент t.
является фиктивным состоянием для эффекта MA (1) на.
- наблюдаемое изменение уровня безработицы.
- наблюдаемая скорость роста nGNP.
и являются Гауссовым рядом нарушений порядка состояния, имеющих среднее 0 и стандартное отклонение 1.
- гауссов ряд инноваций наблюдений, имеющих среднее 0 и стандартное отклонение .
- гауссов ряд инноваций наблюдений, имеющих среднее 0 и стандартное отклонение .
Загрузите набор данных Нельсона-Плоссера, который содержит, среди прочего, уровень безработицы и серию nGNP.
load Data_NelsonPlosser
Предварительно обработайте данные, взяв естественный логарифм серии nGNP и первое различие каждого. Также удалите стартовую NaN
значения из каждой серии.
isNaN = any(ismissing(DataTable),2); % Flag periods containing NaNs gnpn = DataTable.GNPN(~isNaN); u = DataTable.UR(~isNaN); T = size(gnpn,1); % Sample size y = zeros(T-1,2); % Preallocate y(:,1) = diff(u); y(:,2) = diff(log(gnpn));
Этот пример переходит к использованию ряда без NaN
значения. Однако, используя среду фильтра Калмана, программное обеспечение может включать серии, содержащие отсутствующие значения.
Задайте матрицы коэффициентов.
A = [NaN NaN NaN 0; 0 0 0 0; NaN 0 NaN NaN; 0 0 0 0]; B = [1 0;1 0 ; 0 1; 0 1]; C = [1 0 0 0; 0 0 1 0]; D = [NaN 0; 0 NaN];
Задайте модель пространства состояний используя ssm
. Проверьте, что спецификация модели соответствует модели пространства состояний.
Mdl = ssm(A,B,C,D)
Mdl = State-space model type: ssm State vector length: 4 Observation vector length: 2 State disturbance vector length: 2 Observation innovation vector length: 2 Sample size supported by model: Unlimited Unknown parameters for estimation: 8 State variables: x1, x2,... State disturbances: u1, u2,... Observation series: y1, y2,... Observation innovations: e1, e2,... Unknown parameters: c1, c2,... State equations: x1(t) = (c1)x1(t-1) + (c3)x2(t-1) + (c4)x3(t-1) + u1(t) x2(t) = u1(t) x3(t) = (c2)x1(t-1) + (c5)x3(t-1) + (c6)x4(t-1) + u2(t) x4(t) = u2(t) Observation equations: y1(t) = x1(t) + (c7)e1(t) y2(t) = x3(t) + (c8)e2(t) Initial state distribution: Initial state means are not specified. Initial state covariance matrix is not specified. State types are not specified.
Оцените параметры модели и используйте случайный набор начальных значений параметров для оптимизации. Ограничьте оценку и ко всем положительным, вещественным числам, использующим 'lb'
аргумент пары "имя-значение". Для числовой устойчивости задайте Гессиана, когда программное обеспечение вычисляет ковариационную матрицу параметра, используя 'CovMethod'
аргумент пары "имя-значение".
rng(1); params0 = rand(8,1); [EstMdl,estParams] = estimate(Mdl,y,params0,... 'lb',[-Inf -Inf -Inf -Inf -Inf -Inf 0 0],'CovMethod','hessian');
Method: Maximum likelihood (fmincon) Sample size: 61 Logarithmic likelihood: -199.397 Akaike info criterion: 414.793 Bayesian info criterion: 431.68 | Coeff Std Err t Stat Prob ---------------------------------------------------- c(1) | 0.03387 0.15213 0.22262 0.82383 c(2) | -0.01258 0.05749 -0.21876 0.82684 c(3) | 2.49856 0.22759 10.97827 0 c(4) | 0.77438 2.58648 0.29940 0.76464 c(5) | 0.13994 2.64363 0.05294 0.95778 c(6) | 0.00367 2.45477 0.00149 0.99881 c(7) | 0.00239 2.11325 0.00113 0.99910 c(8) | 0.00014 0.12685 0.00113 0.99910 | | Final State Std Dev t Stat Prob x(1) | 1.40000 0.00239 585.18164 0 x(2) | 0.21778 0.91641 0.23765 0.81216 x(3) | 0.04730 0.00014 329.58394 0 x(4) | 0.03568 0.00015 240.95607 0
EstMdl
является ssm
модель, и вы можете получить доступ к ее свойствам с помощью записи через точку.
Моделируйте 1e4
пути наблюдений от установленной модели пространства состояний EstMdl
использование более плавной симуляции. Задайте, чтобы симулировать наблюдения для каждого периода.
numPaths = 1e4;
SimX = simsmooth(EstMdl,y,'NumPaths',numPaths);
SimX
является T - 1
-by- 4
-by- numPaths
матрица, содержащая моделируемые состояния. Строки SimX
соответствуют периодам, столбцы соответствуют состоянию в модели, а страницы - путям.
Оцените среднее сглаженное состояние, стандартные отклонения и 95% доверительные интервалы.
SmoothBar = mean(SimX,3); SmoothSTD = std(SimX,0,3); SmoothCIL = SmoothBar - 1.96*SmoothSTD; SmoothCIU = SmoothBar + 1.96*SmoothSTD;
Оцените гладкие состояния, используя smooth
.
SmoothX = smooth(EstMdl,y);
Постройте график сглаженных состояний и средств моделируемых состояний и их 95% доверительных интервалов.
figure h = plot(dates(2:T),SmoothBar(:,1),'-r',... dates(2:T),SmoothCIL(:,1),':b',... dates(2:T),SmoothCIU(:,1),':b',... dates(2:T),SmoothX(:,1),':k',... 'LineWidth',3); xlabel 'Period'; ylabel 'Unemployment rate'; legend(h([1,2,4]),{'Simulated, smoothed state mean','95% confidence interval',... 'Smoothed states'},'Location','Best'); title 'Smoothed Unemployment Rate'; axis tight
figure h = plot(dates(2:T),SmoothBar(:,3),'-r',... dates(2:T),SmoothCIL(:,3),':b',... dates(2:T),SmoothCIU(:,3),':b',... dates(2:T),SmoothX(:,3),':k',... 'LineWidth',3); xlabel 'Period'; ylabel 'nGNP'; legend(h([1,2,4]),{'Simulated, smoothed state mean','95% confidence interval',... 'Smoothed states'},'Location','Best'); title 'Smoothed nGNP'; axis tight
Средства моделируемого состояния практически идентичны сглаженным состояниям.