Сравните симуляцию, более сглаженную со сглаживавшими состояниями
Этот пример показывает, как результаты симуляции модели в пространстве состояний, более сглаженной (simsmooth), выдерживают сравнение со сглаживавшими состояниями (smooth).
Предположим что отношение между изменением в уровне безработицы () и темп роста номинального валового национального продукта (nGNP) () может быть выражен в следующем, форме модели в пространстве состояний.
где:
изменение в уровне безработицы во время t.
фиктивное состояние для MA (1) эффект на .
nGNP темп роста во время t.
фиктивное состояние для MA (1) эффект на .
наблюдаемое изменение в уровне безработицы.
наблюдаемый nGNP темп роста.
и серия Gaussian воздействий состояния, имеющих среднее значение 0 и стандартное отклонение 1.
серия Gaussian инноваций наблюдения, имеющих среднее значение 0 и стандартное отклонение .
серия Gaussian инноваций наблюдения, имеющих среднее значение 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.13996 2.64342 0.05295 0.95778
c(6) | 0.00366 2.45456 0.00149 0.99881
c(7) | 0.00239 2.11321 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.21191 0
x(2) | 0.21778 0.91641 0.23765 0.81216
x(3) | 0.04730 0.00014 329.63031 0
x(4) | 0.03568 0.00015 240.98656 0
EstMdl является моделью ssm, и можно получить доступ к ее свойствам с помощью записи через точку.
Моделируйте пути 1e4 наблюдений от подходящей, модели в пространстве состояний EstMdl с помощью более сглаженной симуляции. Задайте, чтобы моделировать наблюдения в течение каждого периода.
numPaths = 1e4;
SimX = simsmooth(EstMdl,y,'NumPaths',numPaths);SimX является 4 T - 1 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
Моделируемые средние значения состояния практически идентичны сглаживавшим состояниям.